1 /* *****************************************************************
2     MESQUITE -- The Mesh Quality Improvement Toolkit
3 
4     Copyright 2004 Sandia Corporation and Argonne National
5     Laboratory.  Under the terms of Contract DE-AC04-94AL85000
6     with Sandia Corporation, the U.S. Government retains certain
7     rights in this software.
8 
9     This library is free software; you can redistribute it and/or
10     modify it under the terms of the GNU Lesser General Public
11     License as published by the Free Software Foundation; either
12     version 2.1 of the License, or (at your option) any later version.
13 
14     This library is distributed in the hope that it will be useful,
15     but WITHOUT ANY WARRANTY; without even the implied warranty of
16     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17     Lesser General Public License for more details.
18 
19     You should have received a copy of the GNU Lesser General Public License
20     (lgpl.txt) along with this library; if not, write to the Free Software
21     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
22 
23     diachin2@llnl.gov, djmelan@sandia.gov, mbrewer@sandia.gov,
24     pknupp@sandia.gov, tleurent@mcs.anl.gov, tmunson@mcs.anl.gov
25 
26   ***************************************************************** */
27 // -*- Mode : c++; tab-width: 3; c-tab-always-indent: t; indent-tabs-mode: nil; c-basic-offset: 3 -*-
28 //
29 //   SUMMARY:
30 //     USAGE:
31 //
32 //    AUTHOR: Michael Brewer
33 //       ORG: Sandia National Labs
34 //    E-MAIL: mbrewer@sandia.gov
35 //
36 // ORIG-DATE: Jan. 29, 2003
37 //  LAST-MOD: 25-Feb-04 at 10:49:04 by Thomas Leurent
38 //
39 // DESCRIPTION:
40 // ============
41 /*! \file SphericalGeometryTest.cpp
42 
43 Regression testing using the spherical geometry capabilities in
44 SimplifiedGeometryEngine.
45  */
46 // DESCRIP-END.
47 //
48 
49 #include "meshfiles.h"
50 
51 #include "PatchDataInstances.hpp"
52 #include "cppunit/extensions/HelperMacros.h"
53 #include <math.h>
54 
55 #include "Mesquite.hpp"
56 #include "MsqError.hpp"
57 #include "Vector3D.hpp"
58 #include "InstructionQueue.hpp"
59 #include "PatchData.hpp"
60 //#include "StoppingCriterion.hpp"
61 #include "QualityAssessor.hpp"
62 
63 #include "IdealWeightInverseMeanRatio.hpp"
64 #include "ConditionNumberQualityMetric.hpp"
65 #include "LPtoPTemplate.hpp"
66 #include "EdgeLengthQualityMetric.hpp"
67 #include "LaplacianSmoother.hpp"
68 #include "SmartLaplacianSmoother.hpp"
69 #include "LInfTemplate.hpp"
70 #include "SteepestDescent.hpp"
71 #include "ConjugateGradient.hpp"
72 #include "AspectRatioGammaQualityMetric.hpp"
73 #include "UntangleBetaQualityMetric.hpp"
74 #include "MultiplyQualityMetric.hpp"
75 #include "SphericalDomain.hpp"
76 
77 #include "MeshImpl.hpp"
78 
79 #include <iostream>
80 using std::cout;
81 using std::endl;
82 using namespace MBMesquite;
83 
84 class SphericalGeometryTest : public CppUnit::TestFixture
85 {
86 private:
87   CPPUNIT_TEST_SUITE(SphericalGeometryTest);
88     //run cg on the quad sphere mesh with condition number L2
89   CPPUNIT_TEST (test_cg_mesh_cond_sphere);
90     //run cg on the quad sphere mesh with condition number L2
91   CPPUNIT_TEST (test_smart_lapl_sphere);
92     //run laplacian smoothing on the geo tri mesh
93   CPPUNIT_TEST (test_lapl_geo_sphere);
94 
95   CPPUNIT_TEST_SUITE_END();
96 
97 private:
98   double qualTol;//double used for double comparisons
99   int pF;//PRINT_FLAG
100 public:
setUp()101   void setUp()
102   {
103       //pF=1;//PRINT_FLAG IS ON
104       pF=0;//PRINT_FLAG IS OFF
105         //tolerance double
106       qualTol=MSQ_MIN;
107   }
108 
tearDown()109   void tearDown()
110   {
111   }
112 
113 public:
SphericalGeometryTest()114   SphericalGeometryTest()
115     {}
116 
test_cg_mesh_cond_sphere()117    void test_cg_mesh_cond_sphere()
118    {
119      MBMesquite::MeshImpl mesh;
120      MBMesquite::MsqPrintError err(cout);
121      mesh.read_vtk(MESH_FILES_DIR "2D/vtk/quads/untangled/quads_on_sphere_529.vtk", err);
122      CPPUNIT_ASSERT(!err);
123 
124        //create geometry: sphere, center (2,2,0), radius 3
125      Vector3D center(2,2,0);
126      SphericalDomain msq_geom(center, 3.0);
127 
128        // creates an intruction queue
129      InstructionQueue queue1;
130 
131        // creates a mean ratio quality metric ...
132      ConditionNumberQualityMetric shape;
133      UntangleBetaQualityMetric untan;
134 
135        // ... and builds an objective function with it
136      LPtoPTemplate obj_func(&shape, 2, err);
137        //Make sure no errors
138      CPPUNIT_ASSERT(!err);
139        // creates the steepest descent optimization procedures
140      ConjugateGradient pass1( &obj_func, err );
141        //SteepestDescent* pass2 = new SteepestDescent( obj_func );
142      pass1.use_global_patch();
143        //Make sure no errors
144      CPPUNIT_ASSERT(!err);
145      QualityAssessor qa=QualityAssessor( &shape );
146 
147        //**********Set stopping criterion  5 iterates ****************
148        //StoppingCriterion sc5(StoppingCriterion::NUMBER_OF_PASSES,5);
149        //pass1->set_stopping_criterion(&sc5);
150      TerminationCriterion sc5;
151      sc5.add_iteration_limit( 5 );
152      pass1.set_inner_termination_criterion(&sc5);
153        //CG's debugging print, increase integer to get more print info
154      pass1.set_debugging_level(0);
155 
156        //  queue1.add_preconditioner(pass2, err); CPPUNIT_ASSERT(!err);
157      queue1.set_master_quality_improver(&pass1, err); CPPUNIT_ASSERT(!err);
158        //Make sure no errors
159      CPPUNIT_ASSERT(!err);
160        // launches optimization on mesh_set1
161      MeshDomainAssoc mesh_and_domain = MeshDomainAssoc(&mesh, &msq_geom);
162      double orig_qa_val=qa.loop_over_mesh(&mesh_and_domain, 0, err);
163        //Make sure no errors
164      CPPUNIT_ASSERT(!err);
165      queue1.run_instructions(&mesh_and_domain, err); CPPUNIT_ASSERT(!err);
166        //Make sure no errors
167      CPPUNIT_ASSERT(!err);
168      double fin_qa_val=qa.loop_over_mesh(&mesh_and_domain, 0, err);
169        //Make sure no errors
170      CPPUNIT_ASSERT(!err);
171        //make sure 'quality' improved
172      CPPUNIT_ASSERT( fin_qa_val <= orig_qa_val );
173    }
test_smart_lapl_sphere()174    void test_smart_lapl_sphere()
175      {
176        MBMesquite::MeshImpl mesh;
177        MBMesquite::MsqPrintError err(cout);
178        mesh.read_vtk(MESH_FILES_DIR "2D/vtk/quads/untangled/quads_on_sphere_529.vtk", err);
179 
180          //create geometry sphere:  ratius 1, centered at (0,0,0)
181        Vector3D center(2,2,0);
182        SphericalDomain msq_geom(center, 3.0);
183 
184          // creates an intruction queue
185        InstructionQueue queue1;
186 
187          // creates an edge length metric ...
188        IdealWeightInverseMeanRatio shape_metric(err);
189        LInfTemplate shape_func(&shape_metric);
190 
191          //create the smart laplacian smoother
192        SmartLaplacianSmoother s_lapl(&shape_func);
193          //Make sure no errors
194        CPPUNIT_ASSERT(!err);
195 
196          //*******Set stopping criterion 5 iterates  ***********
197        TerminationCriterion sc5;
198        sc5.add_iteration_limit( 5 );
199        s_lapl.set_outer_termination_criterion(&sc5);
200          //qa, qi, qa
201        queue1.set_master_quality_improver(&s_lapl, err); CPPUNIT_ASSERT(!err);
202          //Make sure no errors
203        CPPUNIT_ASSERT(!err);
204          // launches optimization on mesh_set1
205        QualityAssessor qa=QualityAssessor( &shape_metric );
206        MeshDomainAssoc mesh_and_domain = MeshDomainAssoc(&mesh, &msq_geom);
207        double orig_val=qa.loop_over_mesh(&mesh_and_domain, 0, err);
208 
209          //Make sure no errors
210        CPPUNIT_ASSERT(!err);
211        queue1.run_instructions(&mesh_and_domain, err); CPPUNIT_ASSERT(!err);
212          //Make sure no errors
213        CPPUNIT_ASSERT(!err);
214 
215        double final_val= qa.loop_over_mesh(&mesh_and_domain, 0, err);
216 
217          //Make sure no errors
218        CPPUNIT_ASSERT(!err);
219          //make sure 'quality' improved
220        CPPUNIT_ASSERT( final_val < orig_val );
221      }
222 
test_lapl_geo_sphere()223   void test_lapl_geo_sphere()
224      {
225        MBMesquite::MeshImpl mesh;
226        MBMesquite::MsqPrintError err(cout);
227 
228        mesh.read_vtk(MESH_FILES_DIR "2D/vtk/tris/untangled/Mesquite_geo_10242.vtk", err);
229 
230          //create geometry sphere:  ratius 1, centered at (0,0,0)
231        Vector3D center(0,0,0);
232        MBMesquite::SphericalDomain msq_geom(center, 1.0);
233 
234          // creates an intruction queue
235        InstructionQueue queue1;
236 
237          // creates an edge length metric ...
238        EdgeLengthQualityMetric edg_len;
239 
240          //create the laplacian smoother
241        LaplacianSmoother lapl;
242          //Make sure no errors
243        CPPUNIT_ASSERT(!err);
244 
245          //create a quality assessor
246        QualityAssessor qa=QualityAssessor( &edg_len );
247 
248          //*******Set stopping criterion 10 iterates  ***********
249          //StoppingCriterion sc10(StoppingCriterion::NUMBER_OF_PASSES,10);
250          //lapl->set_stopping_criterion(&sc10);
251        TerminationCriterion sc10;
252        sc10.add_iteration_limit( 10 );
253        lapl.set_outer_termination_criterion(&sc10);
254          //qa, qi, qa
255        queue1.set_master_quality_improver(&lapl, err);
256          //Make sure no errors
257        CPPUNIT_ASSERT(!err);
258          // launches optimization on mesh_set1
259        MeshDomainAssoc mesh_and_domain = MeshDomainAssoc(&mesh, &msq_geom);
260        double orig_qa_val=qa.loop_over_mesh(&mesh_and_domain, 0, err);
261          //Make sure no errors
262        CPPUNIT_ASSERT(!err);
263        queue1.run_instructions(&mesh_and_domain, err); CPPUNIT_ASSERT(!err);
264          //Make sure no errors
265        CPPUNIT_ASSERT(!err);
266        double fin_qa_val=qa.loop_over_mesh(&mesh_and_domain, 0, err);
267          //Make sure no errors
268        CPPUNIT_ASSERT(!err);
269          //make sure 'quality' improved
270        CPPUNIT_ASSERT( fin_qa_val <= orig_qa_val );
271      }
272 
273 
274 };
275 
276 
277 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION(SphericalGeometryTest, "SphericalGeometryTest");
278 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION(SphericalGeometryTest, "Regression");
279