Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
#include <iostream>
#include <set>
#include <limits>
#include <ctime>
#include <vector>
#include <cstdlib>
#include <cstdio>
#include <cassert>
#include <unistd.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <fcntl.h>
#include "moab/Interface.hpp"
#include "MBTagConventions.hpp"
#include "moab/Core.hpp"
#include "moab/Range.hpp"
#include "moab/Skinner.hpp"
#include "moab/AdaptiveKDTree.hpp"
#include "moab/CN.hpp"
#include <sys/resource.h>
Go to the source code of this file.
Defines | |
#define | CHKERROR(A) |
Functions | |
static void | get_time_mem (double &tot_time, double &tot_mem) |
static ErrorCode | merge_duplicate_vertices (Interface &, double epsilon) |
static ErrorCode | min_edge_length (Interface &, double &result) |
static void | usage (const char *argv0, bool help=false) |
int | main (int argc, char *argv[]) |
Variables | |
const char | DEFAULT_FIXED_TAG [] = "fixed" |
const int | MIN_EDGE_LEN_DENOM = 4 |
#define CHKERROR | ( | A | ) |
do \ { \ if( MB_SUCCESS != ( A ) ) \ { \ std::cerr << "Internal error at line " << __LINE__ << std::endl; \ return 3; \ } \ } while( false )
Definition at line 41 of file skin.cpp.
Referenced by main().
void get_time_mem | ( | double & | tot_time, |
double & | tot_mem | ||
) | [static] |
Definition at line 465 of file skin.cpp.
Referenced by main().
{ struct rusage r_usage; getrusage( RUSAGE_SELF, &r_usage ); double utime = (double)r_usage.ru_utime.tv_sec + ( (double)r_usage.ru_utime.tv_usec / 1.e6 ); double stime = (double)r_usage.ru_stime.tv_sec + ( (double)r_usage.ru_stime.tv_usec / 1.e6 ); tot_time = utime + stime; tot_mem = 0; if( 0 != r_usage.ru_maxrss ) { tot_mem = (double)r_usage.ru_maxrss; } else { // this machine doesn't return rss - try going to /proc // print the file name to open char file_str[4096], dum_str[4096]; int file_ptr = open( "/proc/self/stat", O_RDONLY ); int file_len = read( file_ptr, file_str, sizeof( file_str ) - 1 ); if( file_len <= 0 ) { close( file_ptr ); return; } close( file_ptr ); file_str[file_len] = '\0'; // read the preceding fields and the ones we really want... int dum_int; unsigned int dum_uint, vm_size, rss; int num_fields = sscanf( file_str, "%d " // pid "%s " // comm "%c " // state "%d %d %d %d %d " // ppid, pgrp, session, tty, tpgid "%u %u %u %u %u " // flags, minflt, cminflt, majflt, cmajflt "%d %d %d %d %d %d " // utime, stime, cutime, cstime, counter, priority "%u %u " // timeout, itrealvalue "%d " // starttime "%u %u", // vsize, rss &dum_int, dum_str, dum_str, &dum_int, &dum_int, &dum_int, &dum_int, &dum_int, &dum_uint, &dum_uint, &dum_uint, &dum_uint, &dum_uint, &dum_int, &dum_int, &dum_int, &dum_int, &dum_int, &dum_int, &dum_uint, &dum_uint, &dum_int, &vm_size, &rss ); if( num_fields == 24 ) tot_mem = ( (double)vm_size ); } }
int main | ( | int | argc, |
char * | argv[] | ||
) |
Definition at line 82 of file skin.cpp.
References moab::Interface::add_entities(), moab::Range::begin(), CHKERROR, moab::Interface::create_meshset(), DEFAULT_FIXED_TAG, dim, moab::dum, moab::Range::empty(), entities, ErrorCode, moab::Skinner::find_skin(), moab::Interface::get_adjacencies(), moab::Interface::get_entities_by_dimension(), moab::Interface::get_entities_by_type_and_tag(), get_time_mem(), iface, input_file, moab::Interface::list_entities(), moab::Interface::load_mesh(), MATERIAL_SET_TAG_NAME, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_DENSE, MB_TAG_SPARSE, MB_TYPE_INTEGER, MBENTITYSET, moab::Range::merge(), merge_duplicate_vertices(), MESHSET_SET, MIN_EDGE_LEN_DENOM, min_edge_length(), NEUMANN_SET_TAG_NAME, moab::Range::num_of_dimension(), moab::Range::size(), moab::Interface::tag_get_handle(), moab::Interface::tag_set_data(), moab::Interface::UNION, usage, and moab::Interface::write_mesh().
{ int i = 1; std::vector< int > matsets; int neuset_num = -1; bool write_tag = false, write_whole_mesh = false; bool print_perf = false; bool use_vert_elem_adjs = false; bool merge_vertices = false; double merge_epsilon = -1; bool list_skin = false; bool use_scd = false; const char* fixed_tag = DEFAULT_FIXED_TAG; const char *input_file = 0, *output_file = 0; bool no_more_flags = false; char* endptr = 0; long block = 0; // initialize to eliminate compiler warning while( i < argc ) { if( !no_more_flags && argv[i][0] == '-' ) { const int f = i++; for( int j = 1; argv[f][j]; ++j ) { switch( argv[f][j] ) { case 'a': use_vert_elem_adjs = true; break; case 'p': print_perf = true; break; case 't': write_tag = true; break; case 'w': write_whole_mesh = true; break; case 'm': merge_vertices = true; break; case '-': no_more_flags = true; break; case 'h': usage( argv[0], true ); break; case 'l': list_skin = true; break; case 'S': use_scd = true; break; case 'b': if( i == argc || 0 >= ( block = strtol( argv[i], &endptr, 0 ) ) || *endptr ) { std::cerr << "Expected positive integer following '-b' flag" << std::endl; usage( argv[0] ); } matsets.push_back( (int)block ); ++i; break; case 's': if( i == argc || 0 >= ( neuset_num = strtol( argv[i], &endptr, 0 ) ) || *endptr ) { std::cerr << "Expected positive integer following '-s' flag" << std::endl; usage( argv[0] ); } ++i; break; case 'T': if( i == argc || argv[i][0] == '-' ) { std::cerr << "Expected tag name following '-T' flag" << std::endl; usage( argv[0] ); } fixed_tag = argv[i++]; break; case 'M': if( i == argc || 0.0 > ( merge_epsilon = strtod( argv[i], &endptr ) ) || *endptr ) { std::cerr << "Expected positive numeric value following '-M' flag" << std::endl; usage( argv[0] ); } merge_vertices = true; ++i; break; default: std::cerr << "Unrecognized flag: '" << argv[f][j] << "'" << std::endl; usage( argv[0] ); break; } } } else if( input_file && output_file ) { std::cerr << "Extra argument: " << argv[i] << std::endl; usage( argv[0] ); } else if( input_file ) { output_file = argv[i++]; } else { input_file = argv[i++]; } } if( !input_file ) { std::cerr << "No input file specified" << std::endl; usage( argv[0] ); } ErrorCode result; Core mbimpl; Interface* iface = &mbimpl; if( print_perf ) { double tmp_time1, tmp_mem1; get_time_mem( tmp_time1, tmp_mem1 ); std::cout << "Before reading: cpu time = " << tmp_time1 << ", memory = " << tmp_mem1 / 1.0e6 << "MB." << std::endl; } // read input file result = iface->load_mesh( input_file ); if( MB_SUCCESS != result ) { std::cerr << "Failed to load \"" << input_file << "\"." << std::endl; return 2; } std::cerr << "Read \"" << input_file << "\"" << std::endl; if( print_perf ) { double tmp_time2, tmp_mem2; get_time_mem( tmp_time2, tmp_mem2 ); std::cout << "After reading: cpu time = " << tmp_time2 << ", memory = " << tmp_mem2 / 1.0e6 << "MB." << std::endl; } if( merge_vertices ) { if( merge_epsilon < 0.0 ) { if( MB_SUCCESS != min_edge_length( *iface, merge_epsilon ) ) { std::cerr << "Error determining minimum edge length" << std::endl; return 1; } merge_epsilon /= MIN_EDGE_LEN_DENOM; } if( MB_SUCCESS != merge_duplicate_vertices( *iface, merge_epsilon ) ) { std::cerr << "Error merging duplicate vertices" << std::endl; return 1; } } // get entities of largest dimension int dim = 4; Range entities; while( entities.empty() && dim > 1 ) { dim--; result = iface->get_entities_by_dimension( 0, dim, entities );CHKERROR( result ); } Range skin_ents; Tag matset_tag = 0, neuset_tag = 0; result = iface->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, matset_tag ); if( MB_SUCCESS != result ) return 1; result = iface->tag_get_handle( NEUMANN_SET_TAG_NAME, 1, MB_TYPE_INTEGER, neuset_tag ); if( MB_SUCCESS != result ) return 1; if( matsets.empty() ) skin_ents = entities; else { // get all entities in the specified blocks if( 0 == matset_tag ) { std::cerr << "Couldn't find any material sets in this mesh." << std::endl; return 1; } for( std::vector< int >::iterator vit = matsets.begin(); vit != matsets.end(); ++vit ) { int this_matset = *vit; const void* this_matset_ptr = &this_matset; Range this_range, ent_range; result = iface->get_entities_by_type_and_tag( 0, MBENTITYSET, &matset_tag, &this_matset_ptr, 1, this_range ); if( MB_SUCCESS != result ) { std::cerr << "Trouble getting material set #" << *vit << std::endl; return 1; } else if( this_range.empty() ) { std::cerr << "Warning: couldn't find material set " << *vit << std::endl; continue; } result = iface->get_entities_by_dimension( *this_range.begin(), dim, ent_range, true ); if( MB_SUCCESS != result ) continue; skin_ents.merge( ent_range ); } } if( skin_ents.empty() ) { std::cerr << "No entities for which to compute skin; exiting." << std::endl; return 1; } if( use_vert_elem_adjs ) { // make a call which we know will generate vert-elem adjs Range dum_range; result = iface->get_adjacencies( &( *skin_ents.begin() ), 1, 1, false, dum_range ); if( MB_SUCCESS != result ) return 1; } double tmp_time = 0.0, tmp_mem = 0.0; if( print_perf ) { get_time_mem( tmp_time, tmp_mem ); std::cout << "Before skinning: cpu time = " << tmp_time << ", memory = " << tmp_mem / 1.0e6 << "MB." << std::endl; } // skin the mesh Range forward_lower, reverse_lower; Skinner tool( iface ); if( use_scd ) result = tool.find_skin( 0, skin_ents, false, forward_lower, NULL, false, true, true ); else result = tool.find_skin( 0, skin_ents, false, forward_lower, &reverse_lower ); Range boundary; boundary.merge( forward_lower ); boundary.merge( reverse_lower ); if( MB_SUCCESS != result || boundary.empty() ) { std::cerr << "Mesh skinning failed." << std::endl; return 3; } if( list_skin ) { Range skin_verts; result = iface->get_adjacencies( boundary, 0, true, skin_verts, Interface::UNION ); std::cout << "Skin has "; if( skin_ents.num_of_dimension( 3 ) ) std::cout << boundary.num_of_dimension( 2 ) << " faces and "; else if( skin_ents.num_of_dimension( 2 ) ) std::cout << boundary.num_of_dimension( 1 ) << " edges and "; std::cout << skin_verts.size() << " vertices." << std::endl; } if( write_tag ) { // get tag handle Tag tag; int zero = 0; result = iface->tag_get_handle( fixed_tag, 1, MB_TYPE_INTEGER, tag, MB_TAG_DENSE | MB_TAG_CREAT, &zero );CHKERROR( result ); // Set tags std::vector< int > ones; Range bverts; result = iface->get_adjacencies( boundary, 0, false, bverts, Interface::UNION ); if( MB_SUCCESS != result ) { std::cerr << "Trouble getting vertices on boundary." << std::endl; return 1; } ones.resize( bverts.size(), 1 ); result = iface->tag_set_data( tag, bverts, &ones[0] );CHKERROR( result ); } if( -1 != neuset_num ) { // create a neumann set with these entities if( 0 == neuset_tag ) { result = iface->tag_get_handle( "NEUMANN_SET_TAG_NAME", 1, MB_TYPE_INTEGER, neuset_tag, MB_TAG_SPARSE | MB_TAG_CREAT ); if( MB_SUCCESS != result || 0 == neuset_tag ) return 1; } // always create a forward neumann set, assuming we have something in the set EntityHandle forward_neuset = 0; result = iface->create_meshset( MESHSET_SET, forward_neuset ); if( MB_SUCCESS != result || 0 == forward_neuset ) return 1; result = iface->tag_set_data( neuset_tag, &forward_neuset, 1, &neuset_num ); if( MB_SUCCESS != result ) return 1; if( !forward_lower.empty() ) { result = iface->add_entities( forward_neuset, forward_lower ); if( MB_SUCCESS != result ) return 1; } if( !reverse_lower.empty() ) { EntityHandle reverse_neuset = 1; result = iface->create_meshset( MESHSET_SET, reverse_neuset ); if( MB_SUCCESS != result || 0 == forward_neuset ) return 1; result = iface->add_entities( reverse_neuset, reverse_lower ); if( MB_SUCCESS != result ) return 1; Tag sense_tag; int dum_sense = 0; result = iface->tag_get_handle( "SENSE", 1, MB_TYPE_INTEGER, sense_tag, MB_TAG_SPARSE | MB_TAG_CREAT, &dum_sense ); if( result != MB_SUCCESS ) return 1; int sense_val = -1; result = iface->tag_set_data( neuset_tag, &reverse_neuset, 1, &sense_val ); if( MB_SUCCESS != result ) return 0; result = iface->add_entities( forward_neuset, &reverse_neuset, 1 ); if( MB_SUCCESS != result ) return 0; } } if( NULL != output_file && write_whole_mesh ) { // write output file result = iface->write_mesh( output_file ); if( MB_SUCCESS != result ) { std::cerr << "Failed to write \"" << output_file << "\"." << std::endl; return 2; } std::cerr << "Wrote \"" << output_file << "\"" << std::endl; } else if( NULL != output_file ) { // write only skin; write them as one set EntityHandle skin_set; result = iface->create_meshset( MESHSET_SET, skin_set ); if( MB_SUCCESS != result ) return 1; result = iface->add_entities( skin_set, forward_lower ); if( MB_SUCCESS != result ) return 1; result = iface->add_entities( skin_set, reverse_lower ); if( MB_SUCCESS != result ) return 1; int dum = 10000; result = iface->tag_set_data( matset_tag, &skin_set, 1, &dum ); if( MB_SUCCESS != result ) return 1; result = iface->write_mesh( output_file, &skin_set, 1 ); if( MB_SUCCESS != result ) { std::cerr << "Failed to write \"" << output_file << "\"." << std::endl; return 2; } std::cerr << "Wrote \"" << output_file << "\"" << std::endl; } if( print_perf ) { double tot_time, tot_mem; get_time_mem( tot_time, tot_mem ); std::cout << "Total cpu time = " << tot_time << " seconds." << std::endl; std::cout << "Total skin cpu time = " << tot_time - tmp_time << " seconds." << std::endl; std::cout << "Total memory = " << tot_mem / 1024 << " MB." << std::endl; std::cout << "Total skin memory = " << ( tot_mem - tmp_mem ) / 1024 << " MB." << std::endl; std::cout << "Entities: " << std::endl; iface->list_entities( 0, 0 ); } return 0; }
ErrorCode merge_duplicate_vertices | ( | Interface & | moab, |
double | epsilon | ||
) | [static] |
Definition at line 558 of file skin.cpp.
References moab::Range::begin(), moab::AdaptiveKDTree::build_tree(), moab::AdaptiveKDTree::distance_search(), moab::Range::end(), moab::Range::erase(), ErrorCode, moab::Range::find(), moab::Interface::get_coords(), moab::Interface::get_entities_by_type(), MB_SUCCESS, MBVERTEX, moab::Range::merge(), and moab::Interface::merge_entities().
Referenced by main().
{ ErrorCode rval; Range verts; rval = moab.get_entities_by_type( 0, MBVERTEX, verts ); if( MB_SUCCESS != rval ) return rval; AdaptiveKDTree tree( &moab ); EntityHandle root; rval = tree.build_tree( verts, &root ); if( MB_SUCCESS != rval ) { fprintf( stderr, "Failed to build kD-tree.\n" ); return rval; } std::set< EntityHandle > dead_verts; std::vector< EntityHandle > leaves; for( Range::iterator i = verts.begin(); i != verts.end(); ++i ) { double coords[3]; rval = moab.get_coords( &*i, 1, coords ); if( MB_SUCCESS != rval ) return rval; leaves.clear(); ; rval = tree.distance_search( coords, epsilon, leaves, epsilon, epsilon ); if( MB_SUCCESS != rval ) return rval; Range near; for( std::vector< EntityHandle >::iterator j = leaves.begin(); j != leaves.end(); ++j ) { Range tmp; rval = moab.get_entities_by_type( *j, MBVERTEX, tmp ); if( MB_SUCCESS != rval ) return rval; near.merge( tmp.begin(), tmp.end() ); } Range::iterator v = near.find( *i ); assert( v != near.end() ); near.erase( v ); EntityHandle merge = 0; for( Range::iterator j = near.begin(); j != near.end(); ++j ) { if( *j < *i && dead_verts.find( *j ) != dead_verts.end() ) continue; double coords2[3]; rval = moab.get_coords( &*j, 1, coords2 ); if( MB_SUCCESS != rval ) return rval; coords2[0] -= coords[0]; coords2[1] -= coords[1]; coords2[2] -= coords[2]; double dsqr = coords2[0] * coords2[0] + coords2[1] * coords2[1] + coords2[2] * coords2[2]; if( dsqr <= epsilon * epsilon ) { merge = *j; break; } } if( merge ) { dead_verts.insert( *i ); rval = moab.merge_entities( merge, *i, false, true ); if( MB_SUCCESS != rval ) return rval; } } if( dead_verts.empty() ) std::cout << "No duplicate/coincident vertices." << std::endl; else std::cout << "Merged and deleted " << dead_verts.size() << " vertices " << "coincident within " << epsilon << std::endl; return MB_SUCCESS; }
ErrorCode min_edge_length | ( | Interface & | moab, |
double & | result | ||
) | [static] |
Definition at line 513 of file skin.cpp.
References moab::Range::begin(), moab::Range::end(), entities, moab::Range::erase(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_handle(), moab::Range::lower_bound(), MB_SUCCESS, MBENTITYSET, MBVERTEX, moab::CN::NumSubEntities(), moab::CN::SubEntityVertexIndices(), t, moab::Interface::type_from_handle(), and moab::Range::upper_bound().
Referenced by main().
{ double sqr_result = std::numeric_limits< double >::max(); ErrorCode rval; Range entities; rval = moab.get_entities_by_handle( 0, entities ); if( MB_SUCCESS != rval ) return rval; Range::iterator i = entities.upper_bound( MBVERTEX ); entities.erase( entities.begin(), i ); i = entities.lower_bound( MBENTITYSET ); entities.erase( i, entities.end() ); std::vector< EntityHandle > storage; for( i = entities.begin(); i != entities.end(); ++i ) { EntityType t = moab.type_from_handle( *i ); const EntityHandle* conn; int conn_len, indices[2]; rval = moab.get_connectivity( *i, conn, conn_len, true, &storage ); if( MB_SUCCESS != rval ) return rval; int num_edges = CN::NumSubEntities( t, 1 ); for( int j = 0; j < num_edges; ++j ) { CN::SubEntityVertexIndices( t, 1, j, indices ); EntityHandle v[2] = { conn[indices[0]], conn[indices[1]] }; if( v[0] == v[1] ) continue; double c[6]; rval = moab.get_coords( v, 2, c ); if( MB_SUCCESS != rval ) return rval; c[0] -= c[3]; c[1] -= c[4]; c[2] -= c[5]; double len_sqr = c[0] * c[0] + c[1] * c[1] + c[2] * c[2]; if( len_sqr < sqr_result ) sqr_result = len_sqr; } } result = sqrt( sqr_result ); return MB_SUCCESS; }
static void usage | ( | const char * | argv0, |
bool | help = false |
||
) | [static] |
Definition at line 54 of file skin.cpp.
References DEFAULT_FIXED_TAG, help(), and MIN_EDGE_LEN_DENOM.
{ std::ostream& str = help ? std::cout : std::cerr; str << "Usage: " << argv0 << " [-b <block_num> [-b ...] ] [-l] [-m] [-M <n>] [-p] [-s <sideset_num>] [-S] [-t|-T " "<name>] [-w] [-v|-V <n>]" << " <input_file> [<output_file>]" << std::endl; str << "Help : " << argv0 << " -h" << std::endl; if( !help ) exit( 1 ); str << "Options: " << std::endl; str << "-a : Compute skin using vert-elem adjacencies (more memory, less time)." << std::endl; str << "-b <block_num> : Compute skin only for material set/block <block_num>." << std::endl; str << "-p : Print cpu & memory performance." << std::endl; str << "-s <sideset_num> : Put skin in neumann set/sideset <sideset_num>." << std::endl; str << "-S : Look for and use structured mesh information to speed up skinning." << std::endl; str << "-t : Set '" << DEFAULT_FIXED_TAG << "' tag on skin vertices." << std::endl; str << "-T <name> : Create tag with specified name and set to 1 on skin vertices." << std::endl; str << "-w : Write out whole mesh (otherwise just writes skin)." << std::endl; str << "-m : consolidate duplicate vertices" << std::endl; str << "-M <n> : consolidate duplicate vertices with specified tolerance. " "(Default: min_edge_length/" << MIN_EDGE_LEN_DENOM << ")" << std::endl; str << "-l : List total numbers of entities and vertices in skin." << std::endl; exit( 0 ); }
const char DEFAULT_FIXED_TAG[] = "fixed" |
const int MIN_EDGE_LEN_DENOM = 4 |