Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
#include <iostream>
#include <vector>
#include <string>
#include <iomanip>
#include <fstream>
#include <ctime>
#include <cmath>
#include <cassert>
#include <cfloat>
#include "moab/Core.hpp"
#include "moab/Interface.hpp"
#include "moab/verdict/VerdictWrapper.hpp"
#include "moab/NestedRefine.hpp"
Go to the source code of this file.
Defines | |
#define | SUCCESS 0 |
#define | USAGE_ERROR 1 |
#define | NOT_IMPLEMENTED 2 |
Functions | |
static void | print_usage (const char *name, std::ostream &stream) |
static void | usage_error (const char *name) |
bool | parse_id_list (const char *string, int dim, int nval, std::vector< int > &results) |
bool | make_opts_string (std::vector< std::string > options, std::string &opts) |
ErrorCode | get_degree_seq (Core &mb, EntityHandle fileset, int dim, double desired_vol, int &num_levels, std::vector< int > &level_degs) |
ErrorCode | get_max_volume (Core &mb, EntityHandle fileset, int dim, double &vmax) |
int | main (int argc, char *argv[]) |
#define NOT_IMPLEMENTED 2 |
#define USAGE_ERROR 1 |
Definition at line 30 of file umr.cpp.
Referenced by main(), print_usage(), and usage_error().
ErrorCode get_degree_seq | ( | Core & | mb, |
EntityHandle | fileset, | ||
int | dim, | ||
double | desired_vol, | ||
int & | num_levels, | ||
std::vector< int > & | level_degs | ||
) |
Definition at line 411 of file umr.cpp.
References moab::error(), ErrorCode, get_max_volume(), MB_CHK_ERR, and MB_SUCCESS.
Referenced by main().
{ // Find max volume double vmax_global; ErrorCode error = get_max_volume( mb, fileset, dim, vmax_global );MB_CHK_ERR( error ); int init_nl = num_levels; num_levels = 0; level_degs.clear(); // Find a sequence that reduces the maximum volume to desired. // desired_vol should be a fraction or absolute value ? // double remV = vmax_global*desired_vol/dim; double Vdesired = desired_vol; double try_x; double remV = vmax_global; int degs[3][3] = { { 5, 3, 2 }, { 25, 9, 4 }, { 0, 27, 8 } }; if( dim == 1 || dim == 2 ) { while( remV - Vdesired >= 0 ) { try_x = degs[dim - 1][0]; if( ( remV / try_x - Vdesired ) >= 0 ) { level_degs.push_back( 5 ); num_levels += 1; remV /= try_x; } else { try_x = degs[dim - 1][1]; if( ( remV / try_x - Vdesired ) >= 0 ) { level_degs.push_back( 3 ); num_levels += 1; remV /= try_x; } else { try_x = degs[dim - 1][2]; if( ( remV / try_x - Vdesired ) >= 0 ) { level_degs.push_back( 2 ); num_levels += 1; remV /= try_x; } else break; } } } } else { while( remV - Vdesired >= 0 ) { try_x = degs[dim - 1][1]; if( ( remV / try_x - Vdesired ) >= 0 ) { level_degs.push_back( 3 ); num_levels += 1; remV /= try_x; } else { try_x = degs[dim - 1][2]; if( ( remV / try_x - Vdesired ) >= 0 ) { level_degs.push_back( 2 ); num_levels += 1; remV /= try_x; } else break; } } } if( init_nl != 0 && init_nl < num_levels ) { for( int i = level_degs.size(); i >= init_nl; i-- ) level_degs.pop_back(); num_levels = init_nl; } return MB_SUCCESS; }
ErrorCode get_max_volume | ( | Core & | mb, |
EntityHandle | fileset, | ||
int | dim, | ||
double & | vmax | ||
) |
Definition at line 506 of file umr.cpp.
References moab::Range::begin(), moab::Range::clear(), moab::Range::end(), moab::error(), ErrorCode, moab::ParallelComm::filter_pstatus(), moab::Core::get_entities_by_handle(), moab::ParallelComm::get_pcomm(), moab::MB_AREA, MB_CHK_ERR, moab::MB_LENGTH, MB_SUCCESS, moab::MB_VOLUME, PSTATUS_NOT, PSTATUS_NOT_OWNED, moab::VerdictWrapper::quality_measure(), size, and moab::Range::subset_by_dimension().
Referenced by get_degree_seq(), and main().
{ ErrorCode error; VerdictWrapper vw( &mb ); QualityType q; switch( dim ) { case 1: q = MB_LENGTH; break; case 2: q = MB_AREA; break; case 3: q = MB_VOLUME; break; default: return MB_FAILURE; break; } // Get all entities of the highest dimension which is passed as a command line argument. Range allents, owned; error = mb.get_entities_by_handle( fileset, allents );MB_CHK_ERR( error ); owned = allents.subset_by_dimension( dim );MB_CHK_ERR( error ); // Get all owned entities #ifdef MOAB_HAVE_MPI int size = 1; MPI_Comm_size( MPI_COMM_WORLD, &size ); int mpi_err; if( size > 1 ) { // filter the entities not owned, so we do not process them more than once ParallelComm* pcomm = moab::ParallelComm::get_pcomm( &mb, 0 ); Range current = owned; owned.clear(); error = pcomm->filter_pstatus( current, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1, &owned ); if( error != MB_SUCCESS ) { MPI_Finalize(); return MB_FAILURE; } } #endif double vmax_local = 0; // Find the maximum volume of an entity in the owned mesh for( Range::iterator it = owned.begin(); it != owned.end(); it++ ) { double volume; error = vw.quality_measure( *it, q, volume );MB_CHK_ERR( error ); if( volume > vmax_local ) vmax_local = volume; } // Get the global maximum double vmax_global = vmax_local; #ifdef MOAB_HAVE_MPI mpi_err = MPI_Reduce( &vmax_local, &vmax_global, 1, MPI_DOUBLE, MPI_MAX, 0, MPI_COMM_WORLD ); if( mpi_err ) { MPI_Finalize(); return MB_FAILURE; } #endif vmax = vmax_global; return MB_SUCCESS; }
int main | ( | int | argc, |
char * | argv[] | ||
) |
Definition at line 89 of file umr.cpp.
References moab::Range::clear(), moab::Interface::create_meshset(), dim, moab::error(), ErrorCode, moab::NestedRefine::generate_mesh_hierarchy(), get_degree_seq(), moab::Interface::get_entities_by_handle(), get_max_volume(), moab::Interface::load_file(), make_opts_string(), mb, MB_CHK_ERR, MESHSET_SET, output, parse_id_list(), print_usage(), size, moab::Range::size(), moab::Range::subset_by_dimension(), SUCCESS, moab::NestedRefine::timeall, moab::NestedRefine::codeperf::tm_refine, moab::NestedRefine::codeperf::tm_resolve, moab::NestedRefine::codeperf::tm_total, USAGE_ERROR, usage_error(), and moab::Interface::write_file().
{ int proc_id = 0, size = 1; #ifdef MOAB_HAVE_MPI MPI_Init( &argc, &argv ); MPI_Comm_rank( MPI_COMM_WORLD, &proc_id ); MPI_Comm_size( MPI_COMM_WORLD, &size ); #endif int num_levels = 0, dim = 0; std::vector< int > level_degrees; bool optimize = false; bool do_flag = true; bool print_times = false, print_size = false, output = false; bool parallel = false, resolve_shared = false, exchange_ghosts = false; bool printusage = false, parselevels = true; bool qc_vol = false, only_quality = false; double cvol = 0; std::string infile; int i = 1; while( i < argc ) { if( !argv[i][0] && 0 == proc_id ) { usage_error( argv[0] ); #ifdef MOAB_HAVE_MPI MPI_Finalize(); #endif } if( do_flag && argv[i][0] == '-' ) { switch( argv[i][1] ) { // do flag arguments: case '-': do_flag = false; break; case 't': print_times = true; break; case 's': print_size = true; break; case 'h': case 'H': print_usage( argv[0], std::cerr ); printusage = true; break; case 'd': dim = atol( argv[i + 1] ); ++i; break; case 'n': num_levels = atol( argv[i + 1] ); ++i; break; case 'L': if( dim != 0 && num_levels != 0 ) { parselevels = parse_id_list( argv[i + 1], dim, num_levels, level_degrees ); ++i; } else { print_usage( argv[0], std::cerr ); printusage = true; } break; case 'V': qc_vol = true; cvol = strtod( argv[i + 1], NULL ); ++i; break; case 'q': only_quality = true; break; case 'o': output = true; break; case 'O': optimize = true; break; #ifdef MOAB_HAVE_MPI case 'p': parallel = true; resolve_shared = true; if( argv[i][1] == '1' ) exchange_ghosts = true; break; #endif default: ++i; switch( argv[i - 1][1] ) { // case 'O': read_opts.push_back(argv[i]); break; default: std::cerr << "Invalid option: " << argv[i] << std::endl; } } ++i; } // do file names else { infile = argv[i]; ++i; } } if( infile.empty() && !printusage ) print_usage( argv[0], std::cerr ); if( !parselevels ) exit( USAGE_ERROR ); #ifdef MOAB_HAVE_MPI parallel = true; resolve_shared = true; #endif ErrorCode error; Core* moab = new Core; Interface* mb = (Interface*)moab; EntityHandle fileset; // Create a fileset error = mb->create_meshset( MESHSET_SET, fileset );MB_CHK_ERR( error ); // Set the read options for parallel file loading std::vector< std::string > read_opts, write_opts; std::string read_options, write_options; if( parallel && size > 1 ) { read_opts.push_back( "PARALLEL=READ_PART" ); read_opts.push_back( "PARTITION=PARALLEL_PARTITION" ); if( resolve_shared ) read_opts.push_back( "PARALLEL_RESOLVE_SHARED_ENTS" ); if( exchange_ghosts ) read_opts.push_back( "PARALLEL_GHOSTS=3.0.1" ); /* if (output) { write_opts.push_back(";;PARALLEL=WRITE_PART"); std::cout<<"Here"<<std::endl; }*/ } if( !make_opts_string( read_opts, read_options ) || !make_opts_string( write_opts, write_options ) ) { #ifdef MOAB_HAVE_MPI MPI_Finalize(); #endif return USAGE_ERROR; } // Load file std::cout << "READ OPTS=" << (char*)read_options.c_str() << std::endl; error = mb->load_file( infile.c_str(), &fileset, read_options.c_str() );MB_CHK_ERR( error ); // Create the nestedrefine instance #ifdef MOAB_HAVE_MPI ParallelComm* pc = new ParallelComm( mb, MPI_COMM_WORLD ); NestedRefine* uref = new NestedRefine( moab, pc, fileset ); #else NestedRefine* uref = new NestedRefine( moab ); #endif std::vector< EntityHandle > lsets; // std::cout<<"Read input file"<<std::endl; if( only_quality ) { double vmax; error = get_max_volume( *moab, fileset, dim, vmax );MB_CHK_ERR( error ); #ifdef MOAB_HAVE_MPI int rank = 0; MPI_Comm_rank( MPI_COMM_WORLD, &rank ); if( rank == 0 ) std::cout << "Maximum global volume = " << vmax << std::endl; #else std::cout << "Maximum volume = " << vmax << std::endl; #endif exit( SUCCESS ); } // If a volume constraint is given, find an optimal degree sequence to reach the desired volume // constraint. if( qc_vol ) { error = get_degree_seq( *moab, fileset, dim, cvol, num_levels, level_degrees );MB_CHK_ERR( error ); if( dim == 0 ) print_usage( argv[0], std::cerr ); } if( num_levels == 0 ) num_levels = 1; int* ldeg = new int[num_levels]; // std::cout<<"level_degrees.size = "<<level_degrees.size()<<std::endl; if( level_degrees.empty() ) { for( int l = 0; l < num_levels; l++ ) ldeg[l] = 2; } else { for( int l = 0; l < num_levels; l++ ) ldeg[l] = level_degrees[l]; } std::cout << "Starting hierarchy generation" << std::endl; std::cout << "opt = " << optimize << std::endl; error = uref->generate_mesh_hierarchy( num_levels, ldeg, lsets, optimize );MB_CHK_ERR( error ); if( print_times ) { std::cout << "Finished hierarchy generation in " << uref->timeall.tm_total << " secs" << std::endl; if( parallel ) { std::cout << "Time taken for refinement " << uref->timeall.tm_refine << " secs" << std::endl; std::cout << "Time taken for resolving shared interface " << uref->timeall.tm_resolve << " secs" << std::endl; } } else std::cout << "Finished hierarchy generation " << std::endl; if( print_size ) { Range all_ents, ents[4]; error = mb->get_entities_by_handle( fileset, all_ents );MB_CHK_ERR( error ); for( int k = 0; k < 4; k++ ) ents[k] = all_ents.subset_by_dimension( k ); std::cout << std::endl; if( qc_vol ) { double volume; error = get_max_volume( *moab, fileset, dim, volume );MB_CHK_ERR( error ); std::cout << "Mesh size for level 0" << " :: nverts = " << ents[0].size() << ", nedges = " << ents[1].size() << ", nfaces = " << ents[2].size() << ", ncells = " << ents[3].size() << " :: Vmax = " << volume << std::endl; } else std::cout << "Mesh size for level 0" << " :: nverts = " << ents[0].size() << ", nedges = " << ents[1].size() << ", nfaces = " << ents[2].size() << ", ncells = " << ents[3].size() << std::endl; for( int l = 0; l < num_levels; l++ ) { all_ents.clear(); ents[0].clear(); ents[1].clear(); ents[2].clear(); ents[3].clear(); error = mb->get_entities_by_handle( lsets[l + 1], all_ents );MB_CHK_ERR( error ); for( int k = 0; k < 4; k++ ) ents[k] = all_ents.subset_by_dimension( k ); // std::cout<<std::endl; if( qc_vol ) { double volume; error = get_max_volume( *moab, lsets[l + 1], dim, volume );MB_CHK_ERR( error ); std::cout << "Mesh size for level " << l + 1 << " :: nverts = " << ents[0].size() << ", nedges = " << ents[1].size() << ", nfaces = " << ents[2].size() << ", ncells = " << ents[3].size() << " :: Vmax = " << volume << std::endl; } else std::cout << "Mesh size for level " << l + 1 << " :: nverts = " << ents[0].size() << ", nedges = " << ents[1].size() << ", nfaces = " << ents[2].size() << ", ncells = " << ents[3].size() << std::endl; } } if( output ) { for( int l = 0; l < num_levels; l++ ) { std::string::size_type idx1 = infile.find_last_of( "\\/" ); std::string::size_type idx2 = infile.find_last_of( "." ); std::string file = infile.substr( idx1 + 1, idx2 - idx1 - 1 ); std::stringstream out; if( parallel ) out << "_ML_" << l + 1 << ".h5m"; else out << "_ML_" << l + 1 << ".vtk"; file = file + out.str(); const char* output_file = file.c_str(); #ifdef MOAB_HAVE_MPI error = mb->write_file( output_file, 0, ";;PARALLEL=WRITE_PART", &lsets[l + 1], 1 );MB_CHK_ERR( error ); #else error = mb->write_file( output_file, 0, NULL, &lsets[l + 1], 1 );MB_CHK_ERR( error ); #endif // const char* output_file = file.c_str(); // error = mb->write_file(output_file, 0, write_options.c_str(), &lsets[l+1], // 1);MB_CHK_ERR(error); // mb->list_entity(lsets[l+1]); // mb->write_file(output_file, 0, "PARALLEL=WRITE_PART", &lsets[l+1], 1); } } delete uref; #ifdef MOAB_HAVE_MPI delete pc; #endif delete[] ldeg; delete moab; #ifdef MOAB_HAVE_MPI MPI_Finalize(); #endif exit( SUCCESS ); }
bool make_opts_string | ( | std::vector< std::string > | options, |
std::string & | opts | ||
) |
bool parse_id_list | ( | const char * | string, |
int | dim, | ||
int | nval, | ||
std::vector< int > & | results | ||
) |
Definition at line 578 of file umr.cpp.
{ bool okay = true; char* mystr = strdup( string ); for( const char* ptr = strtok( mystr, "," ); ptr; ptr = strtok( 0, "," ) ) { char* endptr; int val = strtol( ptr, &endptr, 0 ); if( dim == 1 || dim == 2 ) { if( val != 2 && val != 3 && val != 5 ) { std::cerr << "Not a valid degree for the passed dimension" << std::endl; okay = false; break; } } else if( dim == 3 ) { if( val != 2 && val != 3 ) { std::cerr << "Not a valid degree for the passed dimension" << std::endl; okay = false; break; } } if( endptr == ptr || val <= 0 ) { std::cerr << "Not a valid id: " << ptr << std::endl; okay = false; break; } results.push_back( val ); } if( (int)results.size() < nval ) { for( int i = results.size(); i <= nval - 1; i++ ) results.push_back( results[0] ); } else if( (int)results.size() > nval ) { for( int i = results.size(); i > nval; i-- ) results.pop_back(); } free( mystr ); return okay; }
static void print_usage | ( | const char * | name, |
std::ostream & | stream | ||
) | [static] |
Definition at line 35 of file umr.cpp.
References USAGE_ERROR.
{ stream << "Usage: " << name << " <options> <input_file> [<input_file2> ...]" << std::endl << "Options: " << std::endl << "\t-h - Print this help text and exit." << std::endl << "\t-d - Dimension of the mesh." << std::endl << "\t-n - Exact or a maximum number of levels for the hierarchy. Default 1." << std::endl << "\t-L - Degree of refinement for each level. Pass an array or a number. It " "is mandatory to pass dimension and num_levels before to use this option. If this flag " "is not used, a default of deg 2 refinement is used. " << std::endl << "\t-V - Pass a desired volume (absolute) . This will generate a hierarchy " "such that the maximum volume is reduced to the given amount approximately. The length " "of the hierarchy can be constrained further if a maximum number of levels is passed. " "It is mandatory to pass the dimension for this option. " << std::endl << "\t-q - Prints out the maximum volume of the mesh and exits the program. " "This option can be used as a guide to volume constrainted mesh hierarchy later. " << std::endl << "\t-t - Print out the time taken to generate hierarchy." << std::endl << "\t-s - Print out the mesh sizes of each level of the generated hierarchy." << std::endl << "\t-o - Specify true for output files for the mesh levels of the hierarchy." << std::endl //<< "\t-O option - Specify read option." << std::endl #ifdef MOAB_HAVE_MPI << "\t-p[0|1|2] - Read in parallel[0], optionally also doing resolve_shared_ents (1) " "and exchange_ghosts (2)" << std::endl #endif ; exit( USAGE_ERROR ); }
static void usage_error | ( | const char * | name | ) | [static] |
Definition at line 67 of file umr.cpp.
References print_usage(), and USAGE_ERROR.
Referenced by main().
{ print_usage( name, std::cerr ); #ifdef MOAB_HAVE_MPI MPI_Finalize(); #endif exit( USAGE_ERROR ); }