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