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 // ORIG-DATE: 19-Feb-02 at 10:57:52
33 // LAST-MOD: 10-Feb-04 at 22:44:58 by Thomas Leurent
34 //
35 //
36 // DESCRIPTION:
37 // ============
38 /*! \file main.cpp
39
40 describe main.cpp here
41
42 */
43 // DESCRIP-END.
44 //
45
46 #include "Mesquite.hpp"
47 #include "MsqIMesh.hpp"
48 #include "MeshImpl.hpp"
49 #include "MsqError.hpp"
50 #include "InstructionQueue.hpp"
51 #include "TerminationCriterion.hpp"
52 #include "QualityAssessor.hpp"
53 #include "PlanarDomain.hpp"
54 #include "MeshWriter.hpp"
55 #include "TestUtil.hpp"
56
57 // algorithms
58 #include "IdealWeightInverseMeanRatio.hpp"
59 #include "EdgeLengthQualityMetric.hpp"
60 #include "LPtoPTemplate.hpp"
61 #include "FeasibleNewton.hpp"
62 #include "ConjugateGradient.hpp"
63 #include "SmartLaplacianSmoother.hpp"
64
65 #include <iostream>
66 using std::cerr;
67 using std::cout;
68 using std::endl;
69
70 #include "iBase.h"
71
72 using namespace MBMesquite;
73
74 std::string default_file_name = TestDir + "/3D/vtk/large_box_hex_1000.vtk";
75
usage()76 void usage()
77 {
78 cout << "main [-N] [filename]" << endl;
79 cout << " -N : Use native representation instead of TSTT implementation\n";
80 cout << " If no file name is specified, will use \""
81 << default_file_name << '"' << endl;
82 exit (1);
83 }
84
85 // Construct a MeshTSTT from the file
86 Mesh* get_imesh_mesh( const char* file_name );
87
88 // Construct a MeshImpl from the file
89 Mesh* get_native_mesh( const char* file_name );
90
91 // Run FeasibleNewton solver
92 int run_global_smoother( Mesh* mesh, MsqError& err );
93
94 // Run SmoothLaplacian solver
95 int run_local_smoother( Mesh* mesh, MsqError& err );
96
main(int argc,char * argv[])97 int main(int argc, char* argv[])
98 {
99 MBMesquite::MsqPrintError err(cout);
100
101 // command line arguments
102 const char* file_name = 0;
103 bool use_native = false, opts_done = false;
104 for (int arg = 1; arg < argc; ++arg)
105 {
106 if (!opts_done && argv[arg][0] == '-')
107 {
108 if (!strcmp( argv[arg], "-N"))
109 use_native = true;
110 else if(!strcmp( argv[arg], "--"))
111 opts_done = true;
112 else
113 usage();
114 }
115 else if (!file_name)
116 file_name = argv[arg];
117 else
118 usage();
119 }
120 if (!file_name)
121 {
122 file_name = default_file_name.c_str();
123 cout << "No file specified: using default: " << default_file_name << endl;
124 }
125
126 // Try running a global smoother on the mesh
127 Mesh* mesh = use_native ?
128 get_native_mesh(file_name) :
129 get_imesh_mesh(file_name);
130 if (!mesh) {
131 std::cerr << "Failed to load input file. Aborting." << std::endl;
132 return 1;
133 }
134
135 MeshWriter::write_vtk(mesh, "original.vtk", err);
136 if (err) return 1;
137 cout << "Wrote \"original.vtk\"" << endl;
138 run_global_smoother( mesh, err );
139 if (err) return 1;
140
141 // Try running a local smoother on the mesh
142 mesh = use_native ?
143 get_native_mesh(file_name) :
144 get_imesh_mesh(file_name);
145 if (!mesh) {
146 std::cerr << "Failed to load input file. Aborting." << std::endl;
147 return 1;
148 }
149
150 run_local_smoother( mesh, err );
151 if (err) return 1;
152
153 return 0;
154 }
155
156
run_global_smoother(Mesh * mesh,MsqError & err)157 int run_global_smoother( Mesh* mesh, MsqError& err )
158 {
159 double OF_value = 0.0001;
160
161 // creates an intruction queue
162 InstructionQueue queue1;
163
164 // creates a mean ratio quality metric ...
165 IdealWeightInverseMeanRatio* mean_ratio = new IdealWeightInverseMeanRatio(err);
166 if (err) return 1;
167 mean_ratio->set_averaging_method(QualityMetric::SUM, err);
168 if (err) return 1;
169
170 // ... and builds an objective function with it
171 LPtoPTemplate* obj_func = new LPtoPTemplate(mean_ratio, 1, err);
172 if (err) return 1;
173
174 // creates the feas newt optimization procedures
175 FeasibleNewton* pass1 = new FeasibleNewton( obj_func );
176 pass1->use_global_patch();
177 if (err) return 1;
178
179 QualityAssessor stop_qa( mean_ratio );
180
181 // **************Set stopping criterion****************
182 TerminationCriterion tc_inner;
183 tc_inner.add_absolute_vertex_movement( OF_value );
184 if (err) return 1;
185 TerminationCriterion tc_outer;
186 tc_outer.add_iteration_limit( 1 );
187 pass1->set_inner_termination_criterion(&tc_inner);
188 pass1->set_outer_termination_criterion(&tc_outer);
189
190 queue1.add_quality_assessor(&stop_qa, err);
191 if (err) return 1;
192
193 // adds 1 pass of pass1 to mesh_set1
194 queue1.set_master_quality_improver(pass1, err);
195 if (err) return 1;
196
197 queue1.add_quality_assessor(&stop_qa, err);
198 if (err) return 1;
199
200 // launches optimization on mesh_set
201 queue1.run_instructions(mesh, err);
202 if (err) return 1;
203
204 MeshWriter::write_vtk(mesh, "feasible-newton-result.vtk", err);
205 if (err) return 1;
206 cout << "Wrote \"feasible-newton-result.vtk\"" << endl;
207
208 //print_timing_diagnostics( cout );
209 return 0;
210 }
211
run_local_smoother(Mesh * mesh,MsqError & err)212 int run_local_smoother( Mesh* mesh, MsqError& err )
213 {
214 double OF_value = 0.0001;
215
216 // creates an intruction queue
217 InstructionQueue queue1;
218
219 // creates a mean ratio quality metric ...
220 IdealWeightInverseMeanRatio* mean_ratio = new IdealWeightInverseMeanRatio(err);
221 if (err) return 1;
222 mean_ratio->set_averaging_method(QualityMetric::SUM, err);
223 if (err) return 1;
224
225 // ... and builds an objective function with it
226 LPtoPTemplate* obj_func = new LPtoPTemplate(mean_ratio, 1, err);
227 if (err) return 1;
228
229 // creates the smart laplacian optimization procedures
230 SmartLaplacianSmoother* pass1 = new SmartLaplacianSmoother( obj_func );
231
232 QualityAssessor stop_qa( mean_ratio );
233
234 // **************Set stopping criterion****************
235 TerminationCriterion tc_inner;
236 tc_inner.add_absolute_vertex_movement( OF_value );
237 TerminationCriterion tc_outer;
238 tc_outer.add_iteration_limit( 1 );
239 pass1->set_inner_termination_criterion(&tc_inner);
240 pass1->set_outer_termination_criterion(&tc_outer);
241
242 queue1.add_quality_assessor(&stop_qa, err);
243 if (err) return 1;
244
245 // adds 1 pass of pass1 to mesh_set
246 queue1.set_master_quality_improver(pass1, err);
247 if (err) return 1;
248
249 queue1.add_quality_assessor(&stop_qa, err);
250 if (err) return 1;
251
252 // launches optimization on mesh_set
253 queue1.run_instructions(mesh, err);
254 if (err) return 1;
255
256 MeshWriter::write_vtk(mesh, "smart-laplacian-result.vtk", err);
257 if (err) return 1;
258 cout << "Wrote \"smart-laplacian-result.vtk\"" << endl;
259
260 //print_timing_diagnostics( cout );
261 return 0;
262 }
263
264
265
get_imesh_mesh(const char * file_name)266 Mesh* get_imesh_mesh( const char* file_name )
267 {
268 int ierr;
269 iMesh_Instance imesh_mesh = 0;
270 iMesh_newMesh( NULL, &imesh_mesh, &ierr, 0 );
271 if (iBase_SUCCESS != ierr) {
272 return 0;
273 }
274
275 iBase_EntitySetHandle root_set;
276 iMesh_getRootSet( imesh_mesh, &root_set, &ierr );
277 if (iBase_SUCCESS != ierr) {
278 iMesh_dtor( imesh_mesh, &ierr );
279 return 0;
280 }
281
282 iMesh_load( imesh_mesh, root_set, file_name, 0, &ierr, strlen(file_name), 0 );
283 if (iBase_SUCCESS != ierr) {
284 std::cerr << file_name << ": failed to load file." << std::endl;
285 iMesh_dtor( imesh_mesh, &ierr );
286 return 0;
287 }
288
289 iBase_TagHandle fixed_tag;
290 iMesh_getTagHandle( imesh_mesh, "fixed", &fixed_tag, &ierr, strlen("fixed") );
291 if (iBase_SUCCESS != ierr) {
292 iMesh_dtor( imesh_mesh, &ierr );
293 return 0;
294 }
295
296 MsqError err;
297 Mesh* result = new MBMesquite::MsqIMesh( imesh_mesh, root_set, iBase_REGION, err, &fixed_tag );
298 if (MSQ_CHKERR(err)) {
299 delete result;
300 cerr << err << endl;
301 return 0;
302 }
303
304 return result;
305 }
306
307
308
get_native_mesh(const char * file_name)309 Mesh* get_native_mesh( const char* file_name )
310 {
311 MsqError err;
312 MeshImpl* mesh = new MeshImpl;
313 mesh->read_vtk( file_name, err );
314 if (err)
315 {
316 cerr << err << endl;
317 exit(3);
318 }
319
320 return mesh;
321 }
322
323
324