MOAB: Mesh Oriented datABase
(version 5.2.1)
|
00001 /* ***************************************************************** 00002 MESQUITE -- The Mesh Quality Improvement Toolkit 00003 00004 Copyright 2004 Sandia Corporation and Argonne National 00005 Laboratory. Under the terms of Contract DE-AC04-94AL85000 00006 with Sandia Corporation, the U.S. Government retains certain 00007 rights in this software. 00008 00009 This library is free software; you can redistribute it and/or 00010 modify it under the terms of the GNU Lesser General Public 00011 License as published by the Free Software Foundation; either 00012 version 2.1 of the License, or (at your option) any later version. 00013 00014 This library is distributed in the hope that it will be useful, 00015 but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00017 Lesser General Public License for more details. 00018 00019 You should have received a copy of the GNU Lesser General Public License 00020 (lgpl.txt) along with this library; if not, write to the Free Software 00021 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00022 00023 diachin2@llnl.gov, djmelan@sandia.gov, mbrewer@sandia.gov, 00024 pknupp@sandia.gov, tleurent@mcs.anl.gov, tmunson@mcs.anl.gov 00025 00026 ***************************************************************** */ 00027 // -*- Mode : c++; tab-width: 3; c-tab-always-indent: t; indent-tabs-mode: nil; c-basic-offset: 3 00028 // -*- 00029 // 00030 // SUMMARY: 00031 // USAGE: 00032 // 00033 // ORIG-DATE: 19-Feb-02 at 10:57:52 00034 // LAST-MOD: 10-Feb-04 at 22:44:58 by Thomas Leurent 00035 // 00036 // 00037 // DESCRIPTION: 00038 // ============ 00039 /*! \file main.cpp 00040 00041 describe main.cpp here 00042 00043 */ 00044 // DESCRIP-END. 00045 // 00046 00047 // #define USE_IMESH 00048 00049 #include "Mesquite.hpp" 00050 #ifdef USE_IMESH 00051 #include "MsqIMesh.hpp" 00052 #else 00053 #include "MsqMOAB.hpp" 00054 #endif 00055 #include "MeshImpl.hpp" 00056 #include "MsqError.hpp" 00057 #include "InstructionQueue.hpp" 00058 #include "TerminationCriterion.hpp" 00059 #include "QualityAssessor.hpp" 00060 #include "PlanarDomain.hpp" 00061 #include "MeshWriter.hpp" 00062 #include "TestUtil.hpp" 00063 00064 // algorithms 00065 #include "IdealWeightInverseMeanRatio.hpp" 00066 #include "EdgeLengthQualityMetric.hpp" 00067 #include "LPtoPTemplate.hpp" 00068 #include "FeasibleNewton.hpp" 00069 #include "ConjugateGradient.hpp" 00070 #include "SmartLaplacianSmoother.hpp" 00071 00072 #include <iostream> 00073 using std::cerr; 00074 using std::cout; 00075 using std::endl; 00076 00077 #include "iBase.h" 00078 00079 using namespace MBMesquite; 00080 00081 std::string default_file_name = TestDir + "/3D/vtk/large_box_hex_1000.vtk"; 00082 00083 void usage() 00084 { 00085 cout << "main [-N] [filename]" << endl; 00086 cout << " -N : Use native representation instead of TSTT implementation\n"; 00087 cout << " If no file name is specified, will use \"" << default_file_name << '"' << endl; 00088 exit( 1 ); 00089 } 00090 00091 // Construct a MeshTSTT from the file 00092 Mesh* get_imesh_mesh( const char* file_name ); 00093 00094 // Construct a MeshImpl from the file 00095 Mesh* get_native_mesh( const char* file_name ); 00096 00097 // Run FeasibleNewton solver 00098 int run_global_smoother( Mesh* mesh, MsqError& err ); 00099 00100 // Run SmoothLaplacian solver 00101 int run_local_smoother( Mesh* mesh, MsqError& err ); 00102 00103 int main( int argc, char* argv[] ) 00104 { 00105 MBMesquite::MsqPrintError err( cout ); 00106 00107 // command line arguments 00108 const char* file_name = 0; 00109 bool use_native = false, opts_done = false; 00110 for( int arg = 1; arg < argc; ++arg ) 00111 { 00112 if( !opts_done && argv[arg][0] == '-' ) 00113 { 00114 if( !strcmp( argv[arg], "-N" ) ) 00115 use_native = true; 00116 else if( !strcmp( argv[arg], "--" ) ) 00117 opts_done = true; 00118 else 00119 usage(); 00120 } 00121 else if( !file_name ) 00122 file_name = argv[arg]; 00123 else 00124 usage(); 00125 } 00126 if( !file_name ) 00127 { 00128 file_name = default_file_name.c_str(); 00129 cout << "No file specified: using default: " << default_file_name << endl; 00130 } 00131 00132 // Try running a global smoother on the mesh 00133 Mesh* mesh = use_native ? get_native_mesh( file_name ) : get_imesh_mesh( file_name ); 00134 if( !mesh ) 00135 { 00136 std::cerr << "Failed to load input file. Aborting." << std::endl; 00137 return 1; 00138 } 00139 00140 MeshWriter::write_vtk( mesh, "original.vtk", err ); 00141 if( err ) return 1; 00142 cout << "Wrote \"original.vtk\"" << endl; 00143 run_global_smoother( mesh, err ); 00144 if( err ) return 1; 00145 00146 // Try running a local smoother on the mesh 00147 mesh = use_native ? get_native_mesh( file_name ) : get_imesh_mesh( file_name ); 00148 if( !mesh ) 00149 { 00150 std::cerr << "Failed to load input file. Aborting." << std::endl; 00151 return 1; 00152 } 00153 00154 run_local_smoother( mesh, err ); 00155 if( err ) return 1; 00156 00157 return 0; 00158 } 00159 00160 int run_global_smoother( Mesh* mesh, MsqError& err ) 00161 { 00162 double OF_value = 0.0001; 00163 00164 // creates an intruction queue 00165 InstructionQueue queue1; 00166 00167 // creates a mean ratio quality metric ... 00168 IdealWeightInverseMeanRatio* mean_ratio = new IdealWeightInverseMeanRatio( err ); 00169 if( err ) return 1; 00170 mean_ratio->set_averaging_method( QualityMetric::SUM, err ); 00171 if( err ) return 1; 00172 00173 // ... and builds an objective function with it 00174 LPtoPTemplate* obj_func = new LPtoPTemplate( mean_ratio, 1, err ); 00175 if( err ) return 1; 00176 00177 // creates the feas newt optimization procedures 00178 FeasibleNewton* pass1 = new FeasibleNewton( obj_func ); 00179 pass1->use_global_patch(); 00180 if( err ) return 1; 00181 00182 QualityAssessor stop_qa( mean_ratio ); 00183 00184 // **************Set stopping criterion**************** 00185 TerminationCriterion tc_inner; 00186 tc_inner.add_absolute_vertex_movement( OF_value ); 00187 if( err ) return 1; 00188 TerminationCriterion tc_outer; 00189 tc_outer.add_iteration_limit( 1 ); 00190 pass1->set_inner_termination_criterion( &tc_inner ); 00191 pass1->set_outer_termination_criterion( &tc_outer ); 00192 00193 queue1.add_quality_assessor( &stop_qa, err ); 00194 if( err ) return 1; 00195 00196 // adds 1 pass of pass1 to mesh_set1 00197 queue1.set_master_quality_improver( pass1, err ); 00198 if( err ) return 1; 00199 00200 queue1.add_quality_assessor( &stop_qa, err ); 00201 if( err ) return 1; 00202 00203 // launches optimization on mesh_set 00204 queue1.run_instructions( mesh, err ); 00205 if( err ) return 1; 00206 00207 MeshWriter::write_vtk( mesh, "feasible-newton-result.vtk", err ); 00208 if( err ) return 1; 00209 cout << "Wrote \"feasible-newton-result.vtk\"" << endl; 00210 00211 // print_timing_diagnostics( cout ); 00212 return 0; 00213 } 00214 00215 int run_local_smoother( Mesh* mesh, MsqError& err ) 00216 { 00217 double OF_value = 0.0001; 00218 00219 // creates an intruction queue 00220 InstructionQueue queue1; 00221 00222 // creates a mean ratio quality metric ... 00223 IdealWeightInverseMeanRatio* mean_ratio = new IdealWeightInverseMeanRatio( err ); 00224 if( err ) return 1; 00225 mean_ratio->set_averaging_method( QualityMetric::SUM, err ); 00226 if( err ) return 1; 00227 00228 // ... and builds an objective function with it 00229 LPtoPTemplate* obj_func = new LPtoPTemplate( mean_ratio, 1, err ); 00230 if( err ) return 1; 00231 00232 // creates the smart laplacian optimization procedures 00233 SmartLaplacianSmoother* pass1 = new SmartLaplacianSmoother( obj_func ); 00234 00235 QualityAssessor stop_qa( mean_ratio ); 00236 00237 // **************Set stopping criterion**************** 00238 TerminationCriterion tc_inner; 00239 tc_inner.add_absolute_vertex_movement( OF_value ); 00240 TerminationCriterion tc_outer; 00241 tc_outer.add_iteration_limit( 1 ); 00242 pass1->set_inner_termination_criterion( &tc_inner ); 00243 pass1->set_outer_termination_criterion( &tc_outer ); 00244 00245 queue1.add_quality_assessor( &stop_qa, err ); 00246 if( err ) return 1; 00247 00248 // adds 1 pass of pass1 to mesh_set 00249 queue1.set_master_quality_improver( pass1, err ); 00250 if( err ) return 1; 00251 00252 queue1.add_quality_assessor( &stop_qa, err ); 00253 if( err ) return 1; 00254 00255 // launches optimization on mesh_set 00256 queue1.run_instructions( mesh, err ); 00257 if( err ) return 1; 00258 00259 MeshWriter::write_vtk( mesh, "smart-laplacian-result.vtk", err ); 00260 if( err ) return 1; 00261 cout << "Wrote \"smart-laplacian-result.vtk\"" << endl; 00262 00263 // print_timing_diagnostics( cout ); 00264 return 0; 00265 } 00266 00267 Mesh* get_imesh_mesh( const char* file_name ) 00268 { 00269 #ifdef USE_IMESH 00270 int ierr; 00271 iMesh_Instance imesh_mesh = 0; 00272 iMesh_newMesh( NULL, &imesh_mesh, &ierr, 0 ); 00273 if( iBase_SUCCESS != ierr ) { return 0; } 00274 00275 iBase_EntitySetHandle root_set; 00276 iMesh_getRootSet( imesh_mesh, &root_set, &ierr ); 00277 if( iBase_SUCCESS != ierr ) 00278 { 00279 iMesh_dtor( imesh_mesh, &ierr ); 00280 return 0; 00281 } 00282 00283 iMesh_load( imesh_mesh, root_set, file_name, 0, &ierr, strlen( file_name ), 0 ); 00284 if( iBase_SUCCESS != ierr ) 00285 { 00286 std::cerr << file_name << ": failed to load file." << std::endl; 00287 iMesh_dtor( imesh_mesh, &ierr ); 00288 return 0; 00289 } 00290 00291 iBase_TagHandle fixed_tag; 00292 iMesh_getTagHandle( imesh_mesh, "fixed", &fixed_tag, &ierr, strlen( "fixed" ) ); 00293 if( iBase_SUCCESS != ierr ) 00294 { 00295 iMesh_dtor( imesh_mesh, &ierr ); 00296 return 0; 00297 } 00298 00299 MsqError err; 00300 Mesh* result = new MBMesquite::MsqIMesh( imesh_mesh, root_set, iBase_REGION, err, &fixed_tag ); 00301 if( MSQ_CHKERR( err ) ) 00302 { 00303 delete result; 00304 cerr << err << endl; 00305 return 0; 00306 } 00307 00308 #else 00309 00310 moab::Core* mb = new( std::nothrow ) moab::Core; 00311 if( NULL == mb ) return 0; 00312 00313 moab::ErrorCode rval; 00314 // This file is in the mesh files directory 00315 rval = mb->load_file( file_name );MB_CHK_SET_ERR_RET_VAL( rval, "Failed to read", 0 ); 00316 00317 moab::Tag fixed_tag; 00318 rval = mb->tag_get_handle( "fixed", fixed_tag );MB_CHK_SET_ERR_RET_VAL( rval, "Failed to create fixed tag", 0 ); 00319 00320 moab::EntityHandle root_set = 0; 00321 MsqError err; 00322 Mesh* result = new MBMesquite::MsqMOAB( mb, root_set, moab::MBHEX, err, &fixed_tag ); 00323 if( MSQ_CHKERR( err ) ) 00324 { 00325 delete result; 00326 cerr << err << endl; 00327 return 0; 00328 } 00329 00330 #endif 00331 00332 return result; 00333 } 00334 00335 Mesh* get_native_mesh( const char* file_name ) 00336 { 00337 MsqError err; 00338 MeshImpl* mesh = new MeshImpl; 00339 if( !mesh ) 00340 { 00341 cerr << "Failed during MeshImpl construction.\n"; 00342 exit( 2 ); 00343 } 00344 mesh->read_vtk( file_name, err ); 00345 if( err ) 00346 { 00347 cerr << err << endl; 00348 exit( 3 ); 00349 } 00350 00351 return mesh; 00352 }