MOAB: Mesh Oriented datABase  (version 5.2.1)
test_mesh_refiner.cpp
Go to the documentation of this file.
00001 #include "moab/Core.hpp"
00002 #include "EdgeSizeSimpleImplicit.hpp"
00003 #include "SimplexTemplateRefiner.hpp"
00004 #include "MeshRefiner.hpp"
00005 #include "moab/Interface.hpp"
00006 #include "MBParallelConventions.h"
00007 
00008 #ifdef MOAB_HAVE_MPI
00009 #include "moab/ParallelComm.hpp"
00010 #include "moab_mpi.h"
00011 #endif  // MOAB_HAVE_MPI
00012 
00013 #include "TestUtil.hpp"
00014 
00015 #include <iostream>
00016 #include <sstream>
00017 #include <map>
00018 
00019 #include <ctime>
00020 
00021 using namespace moab;
00022 
00023 int TestMeshRefiner( int argc, char* argv[] )
00024 {
00025     int nprocs, rank;
00026 #ifdef MOAB_HAVE_MPI
00027     MPI_Init( &argc, &argv );
00028     MPI_Comm_size( MPI_COMM_WORLD, &nprocs );
00029     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
00030 #else   // MOAB_HAVE_MPI
00031     nprocs = 1;
00032     rank   = 0;
00033 #endif  // MOAB_HAVE_MPI
00034     // sleep(20);
00035 
00036     // Create the input mesh and, if -new-mesh is specified, an output mesh
00037     std::string ifname   = argc > 1 ? argv[1] : TestDir + "fourVolsBare.cub";
00038     bool input_is_output = false, do_output = false;
00039     std::string output_filename;
00040     if( argc > 2 )
00041     {
00042         if( !strcmp( argv[2], "-new-mesh" ) )
00043             input_is_output = true;
00044         else
00045         {
00046             do_output       = true;
00047             output_filename = std::string( argv[2] );
00048         }
00049     }
00050 
00051     Interface* imesh = new Core;                            // ( rank, nprocs );
00052     Interface* omesh = input_is_output ? imesh : new Core;  // ( rank, nprocs );
00053 
00054 #ifdef MOAB_HAVE_MPI
00055     // Use an ParallelComm object to help set up the input mesh
00056     ParallelComm* ipcomm = new ParallelComm( imesh, MPI_COMM_WORLD );
00057     // ReadParallel* readpar = new ReadParallel( imesh, ipcomm );
00058 #endif  // MOAB_HAVE_MPI
00059 
00060     EntityHandle set_handle;
00061     std::ostringstream parallel_options;
00062 #ifdef MOAB_HAVE_MPI
00063     if( nprocs > 1 )
00064     {
00065         parallel_options << "PARALLEL=READ_DELETE"
00066                          << ";"  // NB: You can use BCAST_DELETE or READ_DELETE here.
00067                          //<< "PARALLEL=BCAST_DELETE" << ";" // NB: You can use BCAST_DELETE or READ_DELETE here.
00068                          << "PARTITION=MATERIAL_SET"
00069                          << ";"
00070                          << "PARTITION_DISTRIBUTE"
00071                          << ";"
00072                          << "PARALLEL_RESOLVE_SHARED_ENTS"
00073                          << ";"
00074                          << "CPUTIME";
00075     }
00076 #endif
00077     ErrorCode rval = imesh->create_meshset( MESHSET_SET, set_handle );
00078     if( MB_SUCCESS != rval )
00079     {
00080         std::cout << "Trouble creating set, exiting." << std::endl;
00081         return 1;
00082     }
00083 
00084     rval = imesh->load_file( ifname.c_str(), &set_handle, parallel_options.str().c_str() );
00085     if( MB_SUCCESS != rval )
00086     {
00087         std::cout << "Trouble reading mesh file " << ifname << ", exiting." << std::endl;
00088         return 1;
00089     }
00090 
00091     // Print out what we have so far, one process at a time
00092     for( int i = 0; i < nprocs; ++i )
00093     {
00094         MPI_Barrier( MPI_COMM_WORLD );
00095         if( i == rank )
00096         {
00097             std::cout << "\n************** Rank: " << ( rank ) << " of: " << nprocs << "\n";
00098             imesh->list_entities( 0, 0 );
00099             std::cout << "**************\n\n";
00100         }
00101         MPI_Barrier( MPI_COMM_WORLD );
00102     }
00103 
00104     // The refiner will need an implicit function to be used as an indicator function for
00105     // subdivision:
00106     EdgeSizeSimpleImplicit* eval = new EdgeSizeSimpleImplicit();
00107     eval->set_ratio( 2. );
00108     // Refine the mesh
00109     MeshRefiner* mref            = new MeshRefiner( imesh, omesh );
00110     SimplexTemplateRefiner* eref = new SimplexTemplateRefiner;
00111     mref->set_entity_refiner( eref );
00112     // mref->add_vertex_tag( tag_floatular );
00113     // mref->add_vertex_tag( tag_intular );
00114     // (We don't add tag_gid to the refiner's tag manager because it is special)
00115     eref->set_edge_size_evaluator( eval );
00116     Range ents_to_refine;
00117     imesh->get_entities_by_type( set_handle, MBTET, ents_to_refine );  // refine just the tets
00118     // ents_to_refine.insert( set_handle ); // refine everything multiple times (because subsets are
00119     // not disjoint)
00120 #ifdef _WIN32
00121     struct timeval tic, toc;
00122     gettimeofday( &tic, 0 );
00123     mref->refine( ents_to_refine );
00124     gettimeofday( &toc, 0 );
00125     std::cout << "\nTime: " << ( ( toc.tv_sec - tic.tv_sec ) * 1000 + ( toc.tv_usec - tic.tv_usec ) / 1000. )
00126               << " ms\n\n";
00127 #endif
00128 
00129     if( do_output )
00130     {
00131         parallel_options.clear();
00132         if( nprocs > 1 ) parallel_options << "PARALLEL=WRITE_PART";
00133         omesh->write_file( output_filename.c_str(), NULL, parallel_options.str().c_str() );
00134     }
00135 
00136     // Print out the results, one process at a time
00137 #ifdef MOAB_HAVE_MPI
00138     for( int i = 0; i < nprocs; ++i )
00139     {
00140         MPI_Barrier( MPI_COMM_WORLD );
00141         if( i == rank )
00142         {
00143             std::cout << "\n************** Rank: " << ( rank ) << " of: " << nprocs << "\n";
00144             omesh->list_entities( 0, 0 );
00145             std::cout << "**************\n\n";
00146         }
00147         MPI_Barrier( MPI_COMM_WORLD );
00148     }
00149 #else   // MOAB_HAVE_MPI
00150     omesh->list_entities( 0, 1 );
00151 #endif  // MOAB_HAVE_MPI
00152 
00153     // Clean up
00154 #ifdef MOAB_HAVE_MPI
00155     // delete readpar;
00156     delete ipcomm;
00157 #endif  // MOAB_HAVE_MPI
00158     if( omesh != imesh ) delete omesh;
00159     delete imesh;
00160     delete mref;  // mref will delete eref
00161 
00162 #ifdef MOAB_HAVE_MPI
00163     MPI_Barrier( MPI_COMM_WORLD );
00164     MPI_Finalize();
00165 #endif  // MOAB_HAVE_MPI
00166 
00167     return 0;
00168 }
00169 
00170 int main( int argc, char* argv[] )
00171 {
00172     return TestMeshRefiner( argc, argv );
00173 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines