Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
#include "moab/GeomQueryTool.hpp"
#include "moab/GeomTopoTool.hpp"
#include "InitCGMA.hpp"
#include "CGMApp.hpp"
#include "moab/CN.hpp"
#include "moab/Core.hpp"
#include "moab/CartVect.hpp"
#include "moab/FileOptions.hpp"
#include "moab/Skinner.hpp"
#include "quads_to_tris.hpp"
#include <limits>
#include <cstdlib>
#include <sstream>
#include <ctime>
Go to the source code of this file.
Defines | |
#define | GF_CUBIT_FILE_TYPE "CUBIT" |
#define | GF_STEP_FILE_TYPE "STEP" |
#define | GF_IGES_FILE_TYPE "IGES" |
#define | GF_ACIS_TXT_FILE_TYPE "ACIS_SAT" |
#define | GF_ACIS_BIN_FILE_TYPE "ACIS_SAB" |
#define | GF_OCC_BREP_FILE_TYPE "OCC" |
Functions | |
void | tokenize (const std::string &str, std::vector< std::string > &tokens, const char *delimiters) |
ErrorCode | get_group_names (Interface *MBI, const EntityHandle group_set, const Tag nameTag, std::vector< std::string > &grp_names) |
ErrorCode | summarize_cell_volume_change (Interface *MBI, const EntityHandle cgm_file_set, const Tag categoryTag, const Tag dimTag, const Tag sizeTag, const Tag nameTag, const Tag idTag, const bool conserve_mass, const bool debug) |
ErrorCode | remove_empty_cgm_surfs_and_vols (Interface *MBI, const EntityHandle cgm_file_set, const Tag idTag, const Tag dimTag, const bool) |
ErrorCode | build_new_surface (Interface *MBI, EntityHandle &new_surf, const EntityHandle forward_parent_vol, const EntityHandle reverse_parent_vol, const int new_surf_id, const Tag dimTag, const Tag idTag, const Tag categoryTag, const Tag senseTag) |
ErrorCode | orient_faces_outward (Interface *MBI, const Range &faces, const bool) |
void | plot_histogram (const std::string &title, const std::string &x_axis_label, const std::string &y_axis_label, const int n_bins, const double data[], const int n_data) |
void | generate_plots (const double orig[], const double defo[], const int n_elems, const std::string &time_step) |
static double | tet_volume (const CartVect &v0, const CartVect &v1, const CartVect &v2, const CartVect &v3) |
double | measure (Interface *MBI, const EntityHandle element) |
ErrorCode | get_signed_volume (Interface *MBI, const EntityHandle surf_set, const CartVect &offset, double &signed_volume) |
ErrorCode | fix_surface_senses (Interface *MBI, const EntityHandle cgm_file_set, const EntityHandle cub_file_set, const Tag idTag, const Tag dimTag, const Tag senseTag, const bool debug) |
ErrorCode | replace_faceted_cgm_surfs (Interface *MBI, const EntityHandle cgm_file_set, const EntityHandle cub_file_set, const Tag idTag, const Tag dimTag, const bool debug) |
ErrorCode | add_dead_elems_to_impl_compl (Interface *MBI, const EntityHandle cgm_file_set, const EntityHandle cub_file_set, const Tag idTag, const Tag dimTag, const Tag categoryTag, const Tag senseTag, const bool debug) |
const char * | get_geom_file_type (const char *filename) |
const char * | get_geom_fptr_type (FILE *file) |
int | is_cubit_file (FILE *file) |
int | is_step_file (FILE *file) |
int | is_iges_file (FILE *file) |
int | is_acis_txt_file (FILE *file) |
int | is_acis_bin_file (FILE *file) |
int | is_occ_brep_file (FILE *file) |
int | main (int argc, char *argv[]) |
Variables | |
double | DEFAULT_DISTANCE = 0.001 |
double | DEFAULT_LEN = 0.0 |
int | DEFAULT_NORM = 5 |
#define GF_ACIS_BIN_FILE_TYPE "ACIS_SAB" |
Definition at line 20 of file cub2h5m.cpp.
Referenced by get_geom_fptr_type().
#define GF_ACIS_TXT_FILE_TYPE "ACIS_SAT" |
Definition at line 19 of file cub2h5m.cpp.
Referenced by get_geom_fptr_type().
#define GF_CUBIT_FILE_TYPE "CUBIT" |
Definition at line 16 of file cub2h5m.cpp.
Referenced by get_geom_fptr_type().
#define GF_IGES_FILE_TYPE "IGES" |
Definition at line 18 of file cub2h5m.cpp.
Referenced by get_geom_fptr_type().
#define GF_OCC_BREP_FILE_TYPE "OCC" |
Definition at line 21 of file cub2h5m.cpp.
Referenced by get_geom_fptr_type().
#define GF_STEP_FILE_TYPE "STEP" |
Definition at line 17 of file cub2h5m.cpp.
Referenced by get_geom_fptr_type().
ErrorCode add_dead_elems_to_impl_compl | ( | Interface * | MBI, |
const EntityHandle | cgm_file_set, | ||
const EntityHandle | cub_file_set, | ||
const Tag | idTag, | ||
const Tag | dimTag, | ||
const Tag | categoryTag, | ||
const Tag | senseTag, | ||
const bool | debug | ||
) |
Definition at line 860 of file cub2h5m.cpp.
References moab::Interface::add_entities(), moab::Range::begin(), build_new_surface(), moab::Range::empty(), moab::Range::end(), ErrorCode, moab::Skinner::find_skin(), moab::Range::front(), moab::Interface::get_child_meshsets(), moab::Interface::get_entities_by_type(), moab::Interface::get_entities_by_type_and_tag(), idTag, moab::intersect(), MB_INVALID_SIZE, MB_SUCCESS, MBENTITYSET, MBHEX, MBQUAD, orient_faces_outward(), moab::Interface::remove_entities(), moab::Range::size(), moab::subtract(), and moab::Interface::tag_get_data().
Referenced by main().
{ // Get the cgm surfaces ErrorCode result; const int two = 2; const void* const two_val[] = { &two }; Range cgm_surfs; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, &dimTag, two_val, 1, cgm_surfs ); if( MB_SUCCESS != result ) return result; // Get the maximum surface id. This is so that new surfaces do not have // duplicate ids. int max_surf_id = -1; for( Range::const_iterator i = cgm_surfs.begin(); i != cgm_surfs.end(); ++i ) { int surf_id; result = MBI->tag_get_data( idTag, &( *i ), 1, &surf_id ); if( MB_SUCCESS != result ) return result; if( max_surf_id < surf_id ) max_surf_id = surf_id; } std::cout << " Maximum surface id=" << max_surf_id << std::endl; // For each cgm volume, does a cub volume with the same id exist? const int three = 3; const void* const three_val[] = { &three }; Range cgm_vols; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, &dimTag, three_val, 1, cgm_vols ); if( MB_SUCCESS != result ) return result; // get the corresponding cub volume for( Range::iterator i = cgm_vols.begin(); i != cgm_vols.end(); ++i ) { int vol_id; result = MBI->tag_get_data( idTag, &( *i ), 1, &vol_id ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; std::cout << " Volume " << vol_id; // Find the meshed vol set with the same id Range cub_vol; const Tag tags[] = { idTag, dimTag }; const void* const tag_vals[] = { &vol_id, &three }; result = MBI->get_entities_by_type_and_tag( cub_file_set, MBENTITYSET, tags, tag_vals, 2, cub_vol ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; if( 1 != cub_vol.size() ) { std::cout << ": no meshed representation found" << std::endl; continue; } else { std::cout << std::endl; } // get the mesh elements of the volume. Range elems; result = MBI->get_entities_by_type( cub_vol.front(), MBHEX, elems ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; if( debug ) std::cout << " found " << elems.size() << " hex elems" << std::endl; // skin the volumes Skinner tool( MBI ); Range skin_faces; result = tool.find_skin( 0, elems, 2, skin_faces, true ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // Reconcile the difference between faces of the cub file surfaces and skin // faces of the 3D mesh. The difference exists because dead elements have been // removed. Faces are divided into: // cub_faces (in) - the faces in the cub file surface // skin_faces (in) - the faces of the 3D mesh elements // common_faces (out)- the faces common to both the cub file surface and 3D mesh // they are still adjacent to this volume // old_faces (out)- the faces of the cub surface not on the 3D mesh skin // they are no longer adjacent to this vol // new_faces (out)- the faces of the 3D mesh skin not on the cub surface // they are now adjacent to this volume // get cub child surfaces. Range cub_surfs; result = MBI->get_child_meshsets( cub_vol.front(), cub_surfs ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; for( Range::iterator j = cub_surfs.begin(); j != cub_surfs.end(); ++j ) { int surf_id; result = MBI->tag_get_data( idTag, &( *j ), 1, &surf_id ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // get the quads on each surface Range cub_faces; result = MBI->get_entities_by_type( *j, MBQUAD, cub_faces ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // Meshed volumes must have meshed surfaces if( cub_faces.empty() ) { std::cout << " Surface " << surf_id << ": contains no meshed faces" << std::endl; // return MB_ENTITY_NOT_FOUND; } // get the faces common to both the skin and this surface Range common_faces = intersect( cub_faces, skin_faces ); // find the surface faces not on the skin - these are old and need removed Range old_faces = subtract( cub_faces, common_faces ); result = MBI->remove_entities( *j, old_faces ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // remove the common faces from the skin faces skin_faces = subtract( skin_faces, common_faces ); // If no old faces exist we are done if( old_faces.empty() ) continue; std::cout << " Surface " << surf_id << ": " << old_faces.size() << " old faces removed" << std::endl; // Place the old faces in a new surface, because they may still be adjacent // to 3D mesh in another volume. Get the parent vols of the surface. Range cgm_surf; const Tag tags2[] = { idTag, dimTag }; const void* const tag_vals2[] = { &surf_id, &two }; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, tags2, tag_vals2, 2, cgm_surf ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; if( 1 != cgm_surf.size() ) { std::cout << "invalid size" << std::endl; return MB_INVALID_SIZE; } EntityHandle cgm_vols2[2], cub_vols[2]; result = MBI->tag_get_data( senseTag, cgm_surf, &cgm_vols2 ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; result = MBI->tag_get_data( senseTag, &( *j ), 1, &cub_vols ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // for the new surf, replace the current volume with the impl compl vol. // This is because the faces that no longer exist will become adjacent to // the impl compl if( *i == cgm_vols2[0] ) { cgm_vols2[0] = 0; cub_vols[0] = 0; } if( *i == cgm_vols2[1] ) { cgm_vols2[1] = 0; cub_vols[1] = 0; } // If both sides of the surface are the impl comp, do not create the surface. if( 0 == cgm_vols2[0] && 0 == cgm_vols2[1] ) { std::cout << " New surface was not created for old faces because both parents " "are impl_compl volume " << std::endl; continue; } // build the new surface. EntityHandle new_cgm_surf, new_cub_surf; ++max_surf_id; result = build_new_surface( MBI, new_cgm_surf, cgm_vols2[0], cgm_vols2[1], max_surf_id, dimTag, idTag, categoryTag, senseTag ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; result = build_new_surface( MBI, new_cub_surf, cub_vols[0], cub_vols[1], max_surf_id, dimTag, idTag, categoryTag, senseTag ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // add the new surface to the file set and populate it with the old faces result = MBI->add_entities( cgm_file_set, &new_cgm_surf, 1 ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; assert( MB_SUCCESS == result ); result = MBI->add_entities( cub_file_set, &new_cub_surf, 1 ); if( MB_SUCCESS != result ) return result; assert( MB_SUCCESS == result ); result = MBI->add_entities( new_cub_surf, old_faces ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; std::cout << " Surface " << max_surf_id << ": created for " << old_faces.size() << " old faces" << std::endl; } // the remaining skin faces are newly exposed faces Range new_faces = skin_faces; // new skin faces must be assigned to a surface if( new_faces.empty() ) continue; std::cout << " Surface " << max_surf_id + 1 << ": created for " << new_faces.size() << " new faces" << std::endl; // Ensure that faces are oriented outwards result = orient_faces_outward( MBI, new_faces, debug ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // Create the new surface. EntityHandle new_cgm_surf, new_cub_surf; ++max_surf_id; result = build_new_surface( MBI, new_cgm_surf, *i, 0, max_surf_id, dimTag, idTag, categoryTag, senseTag ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; result = build_new_surface( MBI, new_cub_surf, cub_vol.front(), 0, max_surf_id, dimTag, idTag, categoryTag, senseTag ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; // Insert the new surf into file sets and populate it with faces. result = MBI->add_entities( cgm_file_set, &new_cgm_surf, 1 ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; result = MBI->add_entities( cub_file_set, &new_cub_surf, 1 ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; result = MBI->add_entities( new_cub_surf, new_faces ); assert( MB_SUCCESS == result ); if( MB_SUCCESS != result ) return result; } return MB_SUCCESS; }
ErrorCode build_new_surface | ( | Interface * | MBI, |
EntityHandle & | new_surf, | ||
const EntityHandle | forward_parent_vol, | ||
const EntityHandle | reverse_parent_vol, | ||
const int | new_surf_id, | ||
const Tag | dimTag, | ||
const Tag | idTag, | ||
const Tag | categoryTag, | ||
const Tag | senseTag | ||
) |
Definition at line 300 of file cub2h5m.cpp.
References moab::Interface::add_parent_child(), CATEGORY_TAG_SIZE, moab::Interface::create_meshset(), ErrorCode, moab::geom_category, MB_SUCCESS, and moab::Interface::tag_set_data().
Referenced by add_dead_elems_to_impl_compl().
{ ErrorCode result; result = MBI->create_meshset( 0, new_surf ); if( MB_SUCCESS != result ) return result; if( 0 != forward_parent_vol ) { result = MBI->add_parent_child( forward_parent_vol, new_surf ); if( MB_SUCCESS != result ) return result; } if( 0 != reverse_parent_vol ) { result = MBI->add_parent_child( reverse_parent_vol, new_surf ); if( MB_SUCCESS != result ) return result; } const int two = 2; result = MBI->tag_set_data( dimTag, &new_surf, 1, &two ); if( MB_SUCCESS != result ) return result; result = MBI->tag_set_data( idTag, &new_surf, 1, &new_surf_id ); if( MB_SUCCESS != result ) return result; const char geom_category[CATEGORY_TAG_SIZE] = { "Surface\0" }; result = MBI->tag_set_data( categoryTag, &new_surf, 1, &geom_category ); if( MB_SUCCESS != result ) return result; EntityHandle vols[2] = { forward_parent_vol, reverse_parent_vol }; result = MBI->tag_set_data( senseTag, &new_surf, 1, vols ); if( MB_SUCCESS != result ) return result; return MB_SUCCESS; }
ErrorCode fix_surface_senses | ( | Interface * | MBI, |
const EntityHandle | cgm_file_set, | ||
const EntityHandle | cub_file_set, | ||
const Tag | idTag, | ||
const Tag | dimTag, | ||
const Tag | senseTag, | ||
const bool | debug | ||
) |
Definition at line 694 of file cub2h5m.cpp.
References moab::Interface::add_entities(), moab::Range::begin(), moab::Range::end(), ErrorCode, moab::Range::front(), moab::Interface::get_entities_by_type(), moab::Interface::get_entities_by_type_and_tag(), get_signed_volume(), idTag, make_tris_from_quads(), MB_SUCCESS, MBENTITYSET, MBQUAD, moab::Range::size(), moab::Interface::tag_get_data(), and moab::Interface::tag_set_data().
Referenced by main().
{ ErrorCode result; const int two = 2; const void* const two_val[] = { &two }; Range cgm_surfs; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, &dimTag, two_val, 1, cgm_surfs ); if( MB_SUCCESS != result ) return result; for( Range::iterator i = cgm_surfs.begin(); i != cgm_surfs.end(); ++i ) { int surf_id; result = MBI->tag_get_data( idTag, &( *i ), 1, &surf_id ); if( MB_SUCCESS != result ) return result; // Find the meshed surface set with the same id Range cub_surf; const Tag tags[] = { idTag, dimTag }; const void* const tag_vals[] = { &surf_id, &two }; result = MBI->get_entities_by_type_and_tag( cub_file_set, MBENTITYSET, tags, tag_vals, 2, cub_surf ); if( MB_SUCCESS != result ) return result; if( 1 != cub_surf.size() ) { std::cout << " Surface " << surf_id << ": no meshed representation found, using CAD representation instead" << std::endl; continue; } // Get tris that represent the quads of the cub surf Range quads; result = MBI->get_entities_by_type( cub_surf.front(), MBQUAD, quads ); if( MB_SUCCESS != result ) return result; Range cub_tris; result = make_tris_from_quads( MBI, quads, cub_tris ); if( MB_SUCCESS != result ) return result; // Add the tris to the same surface meshset as the quads are inside. result = MBI->add_entities( cub_surf.front(), cub_tris ); if( MB_SUCCESS != result ) return result; // get the signed volume for each surface representation. Keep trying until // The signed volumes are much greater than numerical precision. Planar // surfaces will have a signed volume of zero if the plane goes through the // origin, unless we apply an offset. const int n_attempts = 100; const int max_random = 10; const double min_signed_vol = 0.1; double cgm_signed_vol, cub_signed_vol; for( int j = 0; j < n_attempts; ++j ) { cgm_signed_vol = 0; cub_signed_vol = 0; CartVect offset( std::rand() % max_random, std::rand() % max_random, std::rand() % max_random ); result = get_signed_volume( MBI, *i, offset, cgm_signed_vol ); if( MB_SUCCESS != result ) return result; result = get_signed_volume( MBI, cub_surf.front(), offset, cub_signed_vol ); if( MB_SUCCESS != result ) return result; if( debug ) std::cout << " surf_id=" << surf_id << " cgm_signed_vol=" << cgm_signed_vol << " cub_signed_vol=" << cub_signed_vol << std::endl; if( fabs( cgm_signed_vol ) > min_signed_vol && fabs( cub_signed_vol ) > min_signed_vol ) break; if( n_attempts == j + 1 ) { std::cout << "error: signed volume could not be calculated unambiguously" << std::endl; return MB_FAILURE; } } // If the sign is different, reverse the cgm senses so that both // representations have the same signed volume. if( ( cgm_signed_vol < 0 && cub_signed_vol > 0 ) || ( cgm_signed_vol > 0 && cub_signed_vol < 0 ) ) { EntityHandle cgm_surf_volumes[2], reversed_cgm_surf_volumes[2]; result = MBI->tag_get_data( senseTag, &( *i ), 1, cgm_surf_volumes ); if( MB_SUCCESS != result ) return result; if( MB_SUCCESS != result ) return result; reversed_cgm_surf_volumes[0] = cgm_surf_volumes[1]; reversed_cgm_surf_volumes[1] = cgm_surf_volumes[0]; result = MBI->tag_set_data( senseTag, &( *i ), 1, reversed_cgm_surf_volumes ); if( MB_SUCCESS != result ) return result; } } return MB_SUCCESS; }
void generate_plots | ( | const double | orig[], |
const double | defo[], | ||
const int | n_elems, | ||
const std::string & | time_step | ||
) |
Definition at line 570 of file cub2h5m.cpp.
References plot_histogram().
Referenced by main().
{ // find volume ratio then max and min double* ratio = new double[n_elems]; for( int i = 0; i < n_elems; ++i ) ratio[i] = ( defo[i] - orig[i] ) / orig[i]; plot_histogram( "Predeformed Element Volume", "Num_Elems", "Volume [cc]", 10, orig, n_elems ); std::string title = "Element Volume Change Ratio at Time Step " + time_step; plot_histogram( title, "Num_Elems", "Volume Ratio", 10, ratio, n_elems ); delete[] ratio; }
const char * get_geom_file_type | ( | const char * | filename | ) |
Definition at line 1440 of file cub2h5m.cpp.
References get_geom_fptr_type().
Referenced by main().
{ FILE* file; const char* result = 0; file = fopen( name, "r" ); if( file ) { result = get_geom_fptr_type( file ); fclose( file ); } return result; }
const char * get_geom_fptr_type | ( | FILE * | file | ) |
Definition at line 1455 of file cub2h5m.cpp.
References GF_ACIS_BIN_FILE_TYPE, GF_ACIS_TXT_FILE_TYPE, GF_CUBIT_FILE_TYPE, GF_IGES_FILE_TYPE, GF_OCC_BREP_FILE_TYPE, GF_STEP_FILE_TYPE, is_acis_bin_file(), is_acis_txt_file(), is_cubit_file(), is_iges_file(), is_occ_brep_file(), and is_step_file().
Referenced by get_geom_file_type().
{ static const char* CUBIT_NAME = GF_CUBIT_FILE_TYPE; static const char* STEP_NAME = GF_STEP_FILE_TYPE; static const char* IGES_NAME = GF_IGES_FILE_TYPE; static const char* SAT_NAME = GF_ACIS_TXT_FILE_TYPE; static const char* SAB_NAME = GF_ACIS_BIN_FILE_TYPE; static const char* BREP_NAME = GF_OCC_BREP_FILE_TYPE; if( is_cubit_file( file ) ) return CUBIT_NAME; else if( is_step_file( file ) ) return STEP_NAME; else if( is_iges_file( file ) ) return IGES_NAME; else if( is_acis_bin_file( file ) ) return SAB_NAME; else if( is_acis_txt_file( file ) ) return SAT_NAME; else if( is_occ_brep_file( file ) ) return BREP_NAME; else return 0; }
ErrorCode get_group_names | ( | Interface * | MBI, |
const EntityHandle | group_set, | ||
const Tag | nameTag, | ||
std::vector< std::string > & | grp_names | ||
) |
Definition at line 41 of file cub2h5m.cpp.
References ErrorCode, MB_SUCCESS, MB_TAG_NOT_FOUND, NAME_TAG_SIZE, and moab::Interface::tag_get_data().
Referenced by summarize_cell_volume_change().
{ // get names char name0[NAME_TAG_SIZE]; std::fill( name0, name0 + NAME_TAG_SIZE, '\0' ); ErrorCode result = MBI->tag_get_data( nameTag, &group_set, 1, &name0 ); if( MB_SUCCESS != result && MB_TAG_NOT_FOUND != result ) return MB_FAILURE; if( MB_TAG_NOT_FOUND != result ) grp_names.push_back( std::string( name0 ) ); return MB_SUCCESS; }
ErrorCode get_signed_volume | ( | Interface * | MBI, |
const EntityHandle | surf_set, | ||
const CartVect & | offset, | ||
double & | signed_volume | ||
) |
Definition at line 653 of file cub2h5m.cpp.
References moab::Range::begin(), moab::Range::end(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_type(), MB_INVALID_SIZE, MB_SUCCESS, and MBTRI.
Referenced by fix_surface_senses().
{ ErrorCode rval; Range tris; rval = MBI->get_entities_by_type( surf_set, MBTRI, tris ); if( MB_SUCCESS != rval ) return rval; signed_volume = 0.0; const EntityHandle* conn; int len; CartVect coords[3]; for( Range::iterator j = tris.begin(); j != tris.end(); ++j ) { rval = MBI->get_connectivity( *j, conn, len, true ); if( MB_SUCCESS != rval ) return rval; if( 3 != len ) return MB_INVALID_SIZE; rval = MBI->get_coords( conn, 3, coords[0].array() ); if( MB_SUCCESS != rval ) return rval; // Apply offset to avoid calculating 0 for cases when the origin is in the // plane of the surface. for( unsigned int k = 0; k < 3; ++k ) { coords[k][0] += offset[0]; coords[k][1] += offset[1]; coords[k][2] += offset[2]; } coords[1] -= coords[0]; coords[2] -= coords[0]; signed_volume += ( coords[0] % ( coords[1] * coords[2] ) ); } return MB_SUCCESS; }
int is_acis_bin_file | ( | FILE * | file | ) |
Definition at line 1498 of file cub2h5m.cpp.
References buffer.
Referenced by get_geom_fptr_type().
{ char buffer[15]; return !fseek( file, 0, SEEK_SET ) && fread( buffer, 15, 1, file ) && !memcmp( buffer, "ACIS BinaryFile", 9 ); }
int is_acis_txt_file | ( | FILE * | file | ) |
Definition at line 1504 of file cub2h5m.cpp.
References buffer, and length().
Referenced by get_geom_fptr_type().
{ char buffer[5]; int version, length; if( fseek( file, 0, SEEK_SET ) || 2 != fscanf( file, "%d %*d %*d %*d %d ", &version, &length ) ) return 0; if( version < 1 || version > 0xFFFF ) return 0; // Skip application name if( fseek( file, length, SEEK_CUR ) ) return 0; // Read length of version string followed by first 5 characters if( 2 != fscanf( file, "%d %4s", &length, buffer ) ) return 0; return !strcmp( buffer, "ACIS" ); }
int is_cubit_file | ( | FILE * | file | ) |
Definition at line 1480 of file cub2h5m.cpp.
References buffer.
Referenced by get_geom_fptr_type().
{ unsigned char buffer[4]; return !fseek( file, 0, SEEK_SET ) && fread( buffer, 4, 1, file ) && !memcmp( buffer, "CUBE", 4 ); }
int is_iges_file | ( | FILE * | file | ) |
Definition at line 1492 of file cub2h5m.cpp.
References buffer.
Referenced by get_geom_fptr_type().
{ unsigned char buffer[10]; return !fseek( file, 72, SEEK_SET ) && fread( buffer, 10, 1, file ) && !memcmp( buffer, "S 1\r\n", 10 ); }
int is_occ_brep_file | ( | FILE * | file | ) |
Definition at line 1522 of file cub2h5m.cpp.
References buffer.
Referenced by get_geom_fptr_type().
{ unsigned char buffer[6]; return !fseek( file, 0, SEEK_SET ) && fread( buffer, 6, 1, file ) && !memcmp( buffer, "DBRep_", 6 ); }
int is_step_file | ( | FILE * | file | ) |
Definition at line 1486 of file cub2h5m.cpp.
References buffer.
Referenced by get_geom_fptr_type().
{ unsigned char buffer[9]; return !fseek( file, 0, SEEK_SET ) && fread( buffer, 9, 1, file ) && !memcmp( buffer, "ISO-10303", 9 ); }
int main | ( | int | argc, |
char * | argv[] | ||
) |
Definition at line 1129 of file cub2h5m.cpp.
References add_dead_elems_to_impl_compl(), moab::Range::begin(), CATEGORY_TAG_NAME, CATEGORY_TAG_SIZE, moab::Core::create_meshset(), moab::debug, moab::Range::end(), ErrorCode, fix_surface_senses(), generate_plots(), GEOM_DIMENSION_TAG_NAME, moab::Core::get_entities_by_dimension(), moab::Core::get_entities_by_type_and_tag(), get_geom_file_type(), moab::Core::get_last_error(), moab::Core::globalId_tag(), idTag, moab::Core::load_file(), MB_SUCCESS, MB_TAG_CREAT, MB_TAG_DENSE, MB_TAG_SPARSE, MB_TYPE_DOUBLE, MB_TYPE_HANDLE, MB_TYPE_INTEGER, MB_TYPE_OPAQUE, MBENTITYSET, moab::measure(), moab::GeomQueryTool::measure_volume(), MESHSET_TRACK_OWNER, NAME_TAG_NAME, NAME_TAG_SIZE, nameTag, remove_empty_cgm_surfs_and_vols(), replace_faceted_cgm_surfs(), moab::Core::set_meshset_options(), size, moab::Range::size(), start_time, summarize_cell_volume_change(), moab::Core::tag_get_handle(), moab::Core::tag_set_data(), and moab::Core::write_mesh().
{ clock_t start_time = clock(); const bool debug = false; const char* file_type = NULL; const char* cub_name = 0; const char* exo_name = 0; const char* out_name = 0; const char* time_step = 0; const char* sat_name = 0; double dist_tol = 0.001, len_tol = 0.0; int norm_tol = 5; if( 6 != argc && 9 != argc ) { std::cerr << "To read meshed geometry for MOAB:" << std::endl; std::cerr << "$> <cub_file.cub> <acis_file.sat> <facet_tol> <output_file.h5m> conserve_mass<bool>" << std::endl; std::cerr << "To read meshed geometry for MOAB and update node coordinates:" << std::endl; std::cerr << "$> <cub_file.cub> <acis_file.sat> <facet_tol> <output_file.h5m> " "<deformed_exo_file.e> time_step<int> check_vol_change<bool> conserve_mass<bool>" << std::endl; exit( 4 ); } // check filenames for proper suffix std::string temp; cub_name = argv[1]; temp.assign( cub_name ); if( std::string::npos == temp.find( ".cub" ) ) { std::cerr << "cub_file does not have correct suffix" << std::endl; return 1; } sat_name = argv[2]; // needed because the cub file's embedded sat file does not have groups temp.assign( sat_name ); if( std::string::npos == temp.find( ".sat" ) ) { std::cerr << "sat_file does not have correct suffix" << std::endl; return 1; } out_name = argv[4]; temp.assign( out_name ); if( std::string::npos == temp.find( ".h5m" ) ) { std::cerr << "out_file does not have correct suffix" << std::endl; return 1; } // get the facet tolerance dist_tol = atof( argv[3] ); if( 0 > dist_tol || 1 < dist_tol ) { std::cout << "error: facet_tolerance=" << dist_tol << std::endl; return 1; } // Should the nodes be updated? bool update_coords = false; if( 9 == argc ) { exo_name = argv[5]; temp.assign( exo_name ); if( std::string::npos == temp.find( ".e" ) ) { std::cerr << "e_file does not have correct suffix" << std::endl; return 1; } time_step = argv[6]; update_coords = true; } // Should the volume change be determined? bool determine_volume_change = false; if( 9 == argc ) { temp.assign( argv[7] ); if( std::string::npos != temp.find( "true" ) ) determine_volume_change = true; } // Should densities be changed to conserve mass? bool conserve_mass = false; if( 9 == argc ) { temp.assign( argv[8] ); if( std::string::npos != temp.find( "true" ) ) conserve_mass = true; } else if( 6 == argc ) { temp.assign( argv[5] ); if( std::string::npos != temp.find( "true" ) ) conserve_mass = true; } // Get CGM file type if( !file_type ) { file_type = get_geom_file_type( cub_name ); if( !file_type ) { std::cerr << cub_name << " : unknown file type, try '-t'" << std::endl; exit( 1 ); } } // Read the mesh from the cub file with Tqcdfr Core* MBI = new Core(); ErrorCode result; EntityHandle cub_file_set; result = MBI->create_meshset( 0, cub_file_set ); if( MB_SUCCESS != result ) return result; // Do not ignore the Cubit file version. In testing, a cub file from Cubit12 // did not work. // char cub_options[256] = "120"; // char cub_options[256] = "IGNORE_VERSION"; // result = MBI->load_file(cub_name, &cub_file_set, cub_options, NULL, 0, 0); result = MBI->load_file( cub_name, &cub_file_set, 0, NULL, 0, 0 ); if( MB_SUCCESS != result ) { std::cout << "error: problem reading cub file" << std::endl; return result; } std::cout << "Mesh file read." << std::endl; // Read the ACIS file with ReadCGM char cgm_options[256]; std::cout << " facet tolerance=" << dist_tol << std::endl; sprintf( cgm_options, "CGM_ATTRIBS=yes;FACET_DISTANCE_TOLERANCE=%g;FACET_NORMAL_TOLERANCE=%d;MAX_FACET_EDGE_" "LENGTH=%g;", dist_tol, norm_tol, len_tol ); EntityHandle cgm_file_set; result = MBI->create_meshset( 0, cgm_file_set ); if( MB_SUCCESS != result ) return result; result = MBI->load_file( sat_name, &cgm_file_set, cgm_options, NULL, 0, 0 ); if( MB_SUCCESS != result ) { std::cout << "error: problem reading sat file" << std::endl; return result; } std::cout << "CAD file read." << std::endl; // Create tags Tag dimTag, idTag, categoryTag, senseTag, sizeTag, nameTag; result = MBI->tag_get_handle( GEOM_DIMENSION_TAG_NAME, 1, MB_TYPE_INTEGER, dimTag, MB_TAG_SPARSE | MB_TAG_CREAT ); if( MB_SUCCESS != result ) return result; idTag = MBI->globalId_tag(); result = MBI->tag_get_handle( CATEGORY_TAG_NAME, CATEGORY_TAG_SIZE, MB_TYPE_OPAQUE, categoryTag, MB_TAG_SPARSE | MB_TAG_CREAT ); if( MB_SUCCESS != result ) return result; result = MBI->tag_get_handle( "GEOM_SENSE_2", 2, MB_TYPE_HANDLE, senseTag, MB_TAG_SPARSE | MB_TAG_CREAT ); if( MB_SUCCESS != result ) return result; result = MBI->tag_get_handle( "GEOM_SIZE", 1, MB_TYPE_DOUBLE, sizeTag, MB_TAG_DENSE | MB_TAG_CREAT ); if( MB_SUCCESS != result ) return result; result = MBI->tag_get_handle( NAME_TAG_NAME, NAME_TAG_SIZE, MB_TYPE_OPAQUE, nameTag, MB_TAG_SPARSE | MB_TAG_CREAT ); if( MB_SUCCESS != result ) return result; // Create triangles from the quads of the cub surface sets and add them to the // cub surface sets. Get the signed volume of each surface for both cgm and // cub representations. Change the sense of the cgm representation to match // the cub representation. result = fix_surface_senses( MBI, cgm_file_set, cub_file_set, idTag, dimTag, senseTag, debug ); if( MB_SUCCESS != result ) return result; std::cout << "Fixed CAD surface senses to match meshed surface senses." << std::endl; // Get the 3D elements in the cub file and measure their volume. Range orig_elems; std::vector< double > orig_size; if( determine_volume_change ) { result = MBI->get_entities_by_dimension( 0, 3, orig_elems ); if( MB_SUCCESS != result ) return result; orig_size.resize( orig_elems.size() ); for( unsigned int i = 0; i < orig_elems.size(); ++i ) { orig_size[i] = measure( MBI, orig_elems[i] ); } } // Before updating the nodes and removing dead elements, force the cub volume // sets to track ownership so that dead elements will be deleted from the sets. const int three = 3; const void* const three_val[] = { &three }; Range cub_vols; result = MBI->get_entities_by_type_and_tag( cub_file_set, MBENTITYSET, &dimTag, three_val, 1, cub_vols ); if( MB_SUCCESS != result ) return result; for( Range::const_iterator i = cub_vols.begin(); i != cub_vols.end(); ++i ) { result = MBI->set_meshset_options( *i, MESHSET_TRACK_OWNER ); if( MB_SUCCESS != result ) return result; } // Tag volume sets with the undeformed size of each volume. Range vol_sets; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, &dimTag, three_val, 1, vol_sets ); if( MB_SUCCESS != result ) return result; moab::GeomTopoTool gtt = moab::GeomTopoTool( MBI, false ); moab::GeomQueryTool gqt = moab::GeomQueryTool( >t ); for( Range::const_iterator i = vol_sets.begin(); i != vol_sets.end(); ++i ) { double size; result = gqt.measure_volume( *i, size ); if( MB_SUCCESS != result ) return result; result = MBI->tag_set_data( sizeTag, &( *i ), 1, &size ); if( MB_SUCCESS != result ) return result; } // Update the coordinates if needed. Do not do this before checking surface // sense, because the coordinate update could deform the surfaces too much // to make an accurate comparison. // The cub node ids are unique because cgm vertex ids are tagged on the vertex // meshset and not the vertex itself. // result = MBI->delete_entities( &cub_file_set, 1 ); // if(MB_SUCCESS != result) return result; // Assume dead elements exist until I think of something better. bool dead_elements_exist = true; if( update_coords ) { // ReadNCDF my_ex_reader(MBI); char exo_options[120] = "tdata=coord,"; strcat( exo_options, time_step ); strcat( exo_options, ",set" ); // FileOptions exo_opts(exo_options) ; // opts = "tdata=coord, 100, sum, temp.exo"; // result = my_ex_reader.load_file(exo_name, cgm_file_set, exo_opts, NULL, 0 , 0); // result = my_ex_reader.load_file(exo_name, cub_file_set, exo_opts, NULL, 0 , 0); // result = my_ex_reader.load_file(exo_name, &cub_file_set, exo_opts, NULL, 0 , 0); MBI->load_file( exo_name, &cub_file_set, exo_options ); if( MB_SUCCESS != result ) { std::string last_error; MBI->get_last_error( last_error ); std::cout << "coordinate update failed, " << last_error << std::endl; return result; } std::cout << "Updated mesh nodes with deformed coordinates from exodus file." << std::endl; } if( determine_volume_change ) { // Dead elements have been removed by the deformation. Get the elements that // still exist. Range defo_elems; result = MBI->get_entities_by_dimension( 0, 3, defo_elems ); if( MB_SUCCESS != result ) return result; // Determine the volume of the elements now that a deformation has been // applied. Condense the original array by removing dead elements. double* orig_size_condensed = new double[defo_elems.size()]; double* defo_size_condensed = new double[defo_elems.size()]; int j = 0; for( unsigned int i = 0; i < orig_elems.size(); ++i ) { if( orig_elems[i] == defo_elems[j] ) { orig_size_condensed[j] = orig_size[i]; defo_size_condensed[j] = measure( MBI, defo_elems[j] ); ++j; } } generate_plots( orig_size_condensed, defo_size_condensed, defo_elems.size(), std::string( time_step ) ); delete[] orig_size_condensed; // can't use the same delete[] for both delete[] defo_size_condensed; } // Deal with dead elements. For now, add them to the impl_compl volume. // Extra surfaces are created to do this. if( update_coords && dead_elements_exist ) { result = add_dead_elems_to_impl_compl( MBI, cgm_file_set, cub_file_set, idTag, dimTag, categoryTag, senseTag, debug ); if( MB_SUCCESS != result ) return result; std::cout << "Placed dead elements in implicit complement volume and added required surfaces." << std::endl; } // The quads in the cub_file_set have been updated for dead elements. For each // cgm_surf, if there exists a cub_surf with the same id, replace the cgm tris // with cub_tris (created from the quads). Note the a surface that is not // meshed (in cub file) will not be effected. result = replace_faceted_cgm_surfs( MBI, cgm_file_set, cub_file_set, idTag, dimTag, debug ); if( MB_SUCCESS != result ) return result; std::cout << "Replaced faceted CAD surfaces with meshed surfaces of triangles." << std::endl; result = remove_empty_cgm_surfs_and_vols( MBI, cgm_file_set, idTag, dimTag, debug ); if( MB_SUCCESS != result ) return result; std::cout << "Removed surfaces and volumes that no longer have any triangles." << std::endl; // For each material, sum the volume. If the coordinates were updated for // deformation, summarize the volume change. result = summarize_cell_volume_change( MBI, cgm_file_set, categoryTag, dimTag, sizeTag, nameTag, idTag, conserve_mass, debug ); if( MB_SUCCESS != result ) return result; std::cout << "Summarized the volume change of each material, with respect to the solid model." << std::endl; result = MBI->write_mesh( out_name, &cgm_file_set, 1 ); if( MB_SUCCESS != result ) { std::cout << "write mesh failed" << std::endl; return result; } std::cout << "Saved output file for mesh-based analysis." << std::endl; clock_t end_time = clock(); std::cout << " " << (double)( end_time - start_time ) / CLOCKS_PER_SEC << " seconds" << std::endl; std::cout << std::endl; return 0; }
double measure | ( | Interface * | MBI, |
const EntityHandle | element | ||
) |
Definition at line 591 of file cub2h5m.cpp.
References ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), length(), MB_SUCCESS, MBEDGE, MBHEX, MBPOLYGON, MBPRISM, MBPYRAMID, MBQUAD, MBTET, MBTRI, moab::sum(), tet_volume(), and moab::Interface::type_from_handle().
{ EntityType type = MBI->type_from_handle( element ); const EntityHandle* conn; int num_vertices; ErrorCode result = MBI->get_connectivity( element, conn, num_vertices ); if( MB_SUCCESS != result ) return result; std::vector< CartVect > coords( num_vertices ); result = MBI->get_coords( conn, num_vertices, coords[0].array() ); if( MB_SUCCESS != result ) return result; switch( type ) { case MBEDGE: return ( coords[0] - coords[1] ).length(); case MBTRI: return 0.5 * ( ( coords[1] - coords[0] ) * ( coords[2] - coords[0] ) ).length(); case MBQUAD: num_vertices = 4; case MBPOLYGON: { CartVect mid( 0, 0, 0 ); for( int i = 0; i < num_vertices; ++i ) mid += coords[i]; mid /= num_vertices; double sum = 0.0; for( int i = 0; i < num_vertices; ++i ) { int j = ( i + 1 ) % num_vertices; sum += ( ( mid - coords[i] ) * ( mid - coords[j] ) ).length(); } return 0.5 * sum; } case MBTET: return tet_volume( coords[0], coords[1], coords[2], coords[3] ); case MBPYRAMID: return tet_volume( coords[0], coords[1], coords[2], coords[4] ) + tet_volume( coords[0], coords[2], coords[3], coords[4] ); case MBPRISM: return tet_volume( coords[0], coords[1], coords[2], coords[5] ) + tet_volume( coords[3], coords[5], coords[4], coords[0] ) + tet_volume( coords[1], coords[4], coords[5], coords[0] ); case MBHEX: return tet_volume( coords[0], coords[1], coords[3], coords[4] ) + tet_volume( coords[7], coords[3], coords[6], coords[4] ) + tet_volume( coords[4], coords[5], coords[1], coords[6] ) + tet_volume( coords[1], coords[6], coords[3], coords[4] ) + tet_volume( coords[2], coords[6], coords[3], coords[1] ); default: return 0.0; } }
ErrorCode orient_faces_outward | ( | Interface * | MBI, |
const Range & | faces, | ||
const bool | |||
) |
Definition at line 341 of file cub2h5m.cpp.
References moab::Range::begin(), moab::CN::Dimension(), moab::Range::end(), ErrorCode, moab::Range::front(), moab::Interface::get_adjacencies(), moab::Interface::get_connectivity(), MB_INVALID_SIZE, MB_SUCCESS, moab::Interface::set_connectivity(), moab::CN::SideNumber(), moab::Range::size(), and moab::Interface::type_from_handle().
Referenced by add_dead_elems_to_impl_compl().
{ ErrorCode result; for( Range::const_iterator i = faces.begin(); i != faces.end(); ++i ) { Range adj_elem; result = MBI->get_adjacencies( &( *i ), 1, 3, false, adj_elem ); if( MB_SUCCESS != result ) return result; if( 1 != adj_elem.size() ) return MB_INVALID_SIZE; // get connectivity for element and face const EntityHandle* elem_conn; int elem_n_nodes; result = MBI->get_connectivity( adj_elem.front(), elem_conn, elem_n_nodes ); if( MB_SUCCESS != result ) return result; const EntityHandle* face_conn; int face_n_nodes; result = MBI->get_connectivity( *i, face_conn, face_n_nodes ); if( MB_SUCCESS != result ) return result; // Get the sense of the face wrt the element EntityType elem_type = MBI->type_from_handle( adj_elem.front() ); EntityType face_type = MBI->type_from_handle( *i ); int face_num, offset; int sense = 0; const int face_dim = CN::Dimension( face_type ); int rval = CN::SideNumber( elem_type, elem_conn, face_conn, face_n_nodes, face_dim, face_num, sense, offset ); if( 0 != rval ) return MB_FAILURE; // If the face is not oriented outward wrt the element, reverse it if( -1 == sense ) { EntityHandle new_face_conn[4] = { face_conn[3], face_conn[2], face_conn[1], face_conn[0] }; result = MBI->set_connectivity( *i, new_face_conn, 4 ); if( MB_SUCCESS != result ) return result; } } return MB_SUCCESS; }
void plot_histogram | ( | const std::string & | title, |
const std::string & | x_axis_label, | ||
const std::string & | y_axis_label, | ||
const int | n_bins, | ||
const double | data[], | ||
const int | n_data | ||
) |
Definition at line 500 of file cub2h5m.cpp.
Referenced by generate_plots().
{ // find max and min double min = std::numeric_limits< double >::max(); double max = -std::numeric_limits< double >::max(); for( int i = 0; i < n_data; ++i ) { if( min > data[i] ) min = data[i]; if( max < data[i] ) max = data[i]; } // make bins for histogram double bin_width = ( max - min ) / n_bins; std::vector< int > bins( n_bins ); for( int i = 0; i < n_bins; ++i ) bins[i] = 0; // fill the bins for( int i = 0; i < n_data; ++i ) { double diff = data[i] - min; int bin = diff / bin_width; if( 9 < bin ) bin = 9; // cheap fix for numerical precision error if( 0 > bin ) bin = 0; // cheap fix for numerical precision error ++bins[bin]; } // create bars int max_bin = 0; for( int i = 0; i < n_bins; ++i ) if( max_bin < bins[i] ) max_bin = bins[i]; int bar_height; int max_bar_chars = 72; std::vector< std::string > bars( n_bins ); for( int i = 0; i < n_bins; ++i ) { bar_height = ( max_bar_chars * bins[i] ) / max_bin; for( int j = 0; j < bar_height; ++j ) bars[i] += "*"; } // print histogram header std::cout << std::endl; std::cout << " " << title << std::endl; // print results std::cout.width( 15 ); std::cout << min << std::endl; for( int i = 0; i < n_bins; ++i ) { std::cout.width( 15 ); std::cout << min + ( ( i + 1 ) * bin_width ); std::cout.width( max_bar_chars ); std::cout << bars[i] << bins[i] << std::endl; } // print histogram footer std::cout.width( 15 ); std::cout << y_axis_label; std::cout.width( max_bar_chars ); std::cout << " " << x_axis_label << std::endl; std::cout << std::endl; }
ErrorCode remove_empty_cgm_surfs_and_vols | ( | Interface * | MBI, |
const EntityHandle | cgm_file_set, | ||
const Tag | idTag, | ||
const Tag | dimTag, | ||
const bool | |||
) |
Definition at line 188 of file cub2h5m.cpp.
References moab::Range::begin(), moab::Interface::contains_entities(), moab::Interface::delete_entities(), moab::Range::end(), ErrorCode, moab::Interface::get_child_meshsets(), moab::Interface::get_entities_by_type(), moab::Interface::get_entities_by_type_and_tag(), moab::Interface::get_number_entities_by_type(), moab::Interface::get_parent_meshsets(), MB_SUCCESS, MBENTITYSET, MBTRI, moab::Interface::num_child_meshsets(), moab::Interface::remove_entities(), moab::Interface::remove_parent_child(), and moab::Interface::tag_get_data().
Referenced by main().
{ ErrorCode result; const int two = 2; const void* const two_val[] = { &two }; Range cgm_surfs; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, &dimTag, two_val, 1, cgm_surfs ); if( MB_SUCCESS != result ) return result; for( Range::iterator i = cgm_surfs.begin(); i != cgm_surfs.end(); ++i ) { int n_tris; result = MBI->get_number_entities_by_type( *i, MBTRI, n_tris ); if( MB_SUCCESS != result ) return result; if( 0 == n_tris ) { int surf_id; result = MBI->tag_get_data( idTag, &( *i ), 1, &surf_id ); assert( MB_SUCCESS == result ); Range parent_vols; result = MBI->get_parent_meshsets( *i, parent_vols ); assert( MB_SUCCESS == result ); for( Range::iterator j = parent_vols.begin(); j != parent_vols.end(); ++j ) { result = MBI->remove_parent_child( *j, *i ); assert( MB_SUCCESS == result ); } Range child_curves; result = MBI->get_child_meshsets( *i, child_curves ); assert( MB_SUCCESS == result ); for( Range::iterator j = child_curves.begin(); j != child_curves.end(); ++j ) { result = MBI->remove_parent_child( *i, *j ); assert( MB_SUCCESS == result ); } result = MBI->remove_entities( cgm_file_set, &( *i ), 1 ); assert( MB_SUCCESS == result ); // Is the set contained anywhere else? If the surface is in a CUBIT group, // such as "unmerged_surfs" it will cause write_mesh to fail. This should // be a MOAB bug. Range all_sets; result = MBI->get_entities_by_type( 0, MBENTITYSET, all_sets ); assert( MB_SUCCESS == result ); for( Range::iterator j = all_sets.begin(); j != all_sets.end(); ++j ) { if( MBI->contains_entities( *j, &( *i ), 1 ) ) { result = MBI->remove_entities( *j, &( *i ), 1 ); assert( MB_SUCCESS == result ); } } result = MBI->delete_entities( &( *i ), 1 ); assert( MB_SUCCESS == result ); std::cout << " Surface " << surf_id << ": removed because all of its mesh faces have been removed" << std::endl; } } // get all volumes const int three = 3; const void* const three_val[] = { &three }; Range cgm_vols; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, &dimTag, three_val, 1, cgm_vols ); if( MB_SUCCESS != result ) return result; for( Range::iterator i = cgm_vols.begin(); i != cgm_vols.end(); ++i ) { // get the volume's number of surfaces int n_surfs; result = MBI->num_child_meshsets( *i, &n_surfs ); assert( MB_SUCCESS == result ); // if no surfaces remain, remove the volume if( 0 == n_surfs ) { int vol_id; result = MBI->tag_get_data( idTag, &( *i ), 1, &vol_id ); assert( MB_SUCCESS == result ); // Is the set contained anywhere else? If the surface is in a CUBIT group, // such as "unmerged_surfs" it will cause write_mesh to fail. This should // be a MOAB bug. Range all_sets; result = MBI->get_entities_by_type( 0, MBENTITYSET, all_sets ); assert( MB_SUCCESS == result ); for( Range::iterator j = all_sets.begin(); j != all_sets.end(); ++j ) { if( MBI->contains_entities( *j, &( *i ), 1 ) ) { result = MBI->remove_entities( *j, &( *i ), 1 ); assert( MB_SUCCESS == result ); } } result = MBI->delete_entities( &( *i ), 1 ); assert( MB_SUCCESS == result ); std::cout << " Volume " << vol_id << ": removed because all of its surfaces have been removed" << std::endl; } } return MB_SUCCESS; }
ErrorCode replace_faceted_cgm_surfs | ( | Interface * | MBI, |
const EntityHandle | cgm_file_set, | ||
const EntityHandle | cub_file_set, | ||
const Tag | idTag, | ||
const Tag | dimTag, | ||
const bool | debug | ||
) |
Definition at line 791 of file cub2h5m.cpp.
References moab::Interface::add_entities(), moab::Range::begin(), moab::Interface::delete_entities(), moab::Range::end(), ErrorCode, moab::Range::front(), moab::Interface::get_entities_by_type(), moab::Interface::get_entities_by_type_and_tag(), idTag, make_tris_from_quads(), MB_SUCCESS, MBENTITYSET, MBQUAD, MBTRI, moab::Interface::remove_entities(), moab::Range::size(), and moab::Interface::tag_get_data().
Referenced by main().
{ ErrorCode result; const int two = 2; const void* const two_val[] = { &two }; Range cgm_surfs; result = MBI->get_entities_by_type_and_tag( cgm_file_set, MBENTITYSET, &dimTag, two_val, 1, cgm_surfs ); if( MB_SUCCESS != result ) return result; for( Range::iterator i = cgm_surfs.begin(); i != cgm_surfs.end(); ++i ) { int surf_id; result = MBI->tag_get_data( idTag, &( *i ), 1, &surf_id ); if( MB_SUCCESS != result ) return result; if( debug ) std::cout << "surf_id=" << surf_id << std::endl; // Find the meshed surface set with the same id Range cub_surf; const Tag tags[] = { idTag, dimTag }; const void* const tag_vals[] = { &surf_id, &two }; result = MBI->get_entities_by_type_and_tag( cub_file_set, MBENTITYSET, tags, tag_vals, 2, cub_surf ); if( MB_SUCCESS != result ) return result; if( 1 != cub_surf.size() ) { std::cout << " Surface " << surf_id << ": no meshed representation found, using CAD representation instead" << std::endl; continue; } // Get tris that represent the quads of the cub surf Range quads; result = MBI->get_entities_by_type( cub_surf.front(), MBQUAD, quads ); if( MB_SUCCESS != result ) return result; Range cub_tris; result = make_tris_from_quads( MBI, quads, cub_tris ); if( MB_SUCCESS != result ) return result; // Remove the tris from the cgm surf. Don't forget to remove them from the // cgm_file_set because it is not TRACKING. Range cgm_tris; result = MBI->get_entities_by_type( *i, MBTRI, cgm_tris ); if( MB_SUCCESS != result ) return result; result = MBI->remove_entities( *i, cgm_tris ); if( MB_SUCCESS != result ) return result; result = MBI->remove_entities( cgm_file_set, cgm_tris ); if( MB_SUCCESS != result ) return result; result = MBI->delete_entities( cgm_tris ); if( MB_SUCCESS != result ) return result; // Add the cub_tris to the cgm_surf result = MBI->add_entities( *i, cub_tris ); if( MB_SUCCESS != result ) return result; } return MB_SUCCESS; }
ErrorCode summarize_cell_volume_change | ( | Interface * | MBI, |
const EntityHandle | cgm_file_set, | ||
const Tag | categoryTag, | ||
const Tag | dimTag, | ||
const Tag | sizeTag, | ||
const Tag | nameTag, | ||
const Tag | idTag, | ||
const bool | conserve_mass, | ||
const bool | debug | ||
) |
Definition at line 59 of file cub2h5m.cpp.
References moab::Interface::add_entities(), moab::Range::begin(), CATEGORY_TAG_SIZE, moab::Interface::create_meshset(), moab::Range::end(), ErrorCode, moab::Interface::get_entities_by_type_and_tag(), get_group_names(), MB_SUCCESS, MBENTITYSET, moab::GeomQueryTool::measure_volume(), MESHSET_SET, moab::Interface::remove_entities(), moab::Interface::tag_get_data(), moab::Interface::tag_set_data(), and tokenize().
Referenced by main().
{ // get groups ErrorCode rval; const char group_category[CATEGORY_TAG_SIZE] = { "Group\0" }; const void* const group_val[] = { &group_category }; Range groups; rval = MBI->get_entities_by_type_and_tag( 0, MBENTITYSET, &categoryTag, group_val, 1, groups ); if( MB_SUCCESS != rval ) return rval; // Get the maximum group id. This is so that new groups do not have // duplicate ids. int max_grp_id = -1; for( Range::const_iterator i = groups.begin(); i != groups.end(); ++i ) { int grp_id; rval = MBI->tag_get_data( idTag, &( *i ), 1, &grp_id ); if( MB_SUCCESS != rval ) return rval; if( max_grp_id < grp_id ) max_grp_id = grp_id; } if( conserve_mass ) { std::cout << " Altering group densities to conserve mass for each volume." << std::endl; std::cout << " maximum group id=" << max_grp_id << std::endl; } for( Range::const_iterator i = groups.begin(); i != groups.end(); ++i ) { // get group names std::vector< std::string > grp_names; rval = get_group_names( MBI, *i, nameTag, grp_names ); if( MB_SUCCESS != rval ) return MB_FAILURE; // determine if it is a material group bool material_grp = false; int mat_id = -1; double rho = 0; for( std::vector< std::string >::const_iterator j = grp_names.begin(); j != grp_names.end(); ++j ) { if( std::string::npos != ( *j ).find( "mat" ) && std::string::npos != ( *j ).find( "rho" ) ) { material_grp = true; std::cout << " material group: " << *j << std::endl; // get the density and material id std::vector< std::string > tokens; tokenize( *j, tokens, "_" ); mat_id = atoi( tokens[1].c_str() ); rho = atof( tokens[3].c_str() ); } } if( !material_grp ) continue; // get the volume sets of the material group const int three = 3; const void* const three_val[] = { &three }; Range vols; rval = MBI->get_entities_by_type_and_tag( *i, MBENTITYSET, &dimTag, three_val, 1, vols ); if( MB_SUCCESS != rval ) return rval; // for each volume, sum predeformed and deformed volume double orig_grp_volume = 0, defo_grp_volume = 0; moab::GeomTopoTool gtt = moab::GeomTopoTool( MBI, false ); moab::GeomQueryTool gqt = moab::GeomQueryTool( >t ); for( Range::const_iterator j = vols.begin(); j != vols.end(); ++j ) { double defo_size = 0, orig_size = 0; rval = gqt.measure_volume( *j, defo_size ); if( MB_SUCCESS != rval ) return rval; defo_grp_volume += defo_size; rval = MBI->tag_get_data( sizeTag, &( *j ), 1, &orig_size ); if( MB_SUCCESS != rval ) return rval; orig_grp_volume += orig_size; // calculate a new density to conserve mass through the deformation if( !conserve_mass ) continue; double new_rho = rho * orig_size / defo_size; // create a group for the volume with modified density EntityHandle new_grp; rval = MBI->create_meshset( MESHSET_SET, new_grp ); if( MB_SUCCESS != rval ) return rval; std::stringstream new_name_ss; new_name_ss << "mat_" << mat_id << "_rho_" << new_rho << "\0"; std::string new_name; new_name_ss >> new_name; rval = MBI->tag_set_data( nameTag, &new_grp, 1, new_name.c_str() ); if( MB_SUCCESS != rval ) return rval; max_grp_id++; rval = MBI->tag_set_data( idTag, &new_grp, 1, &max_grp_id ); if( MB_SUCCESS != rval ) return rval; const char group_category2[CATEGORY_TAG_SIZE] = "Group\0"; rval = MBI->tag_set_data( categoryTag, &new_grp, 1, group_category2 ); if( MB_SUCCESS != rval ) return rval; // add the volume to the new group rval = MBI->add_entities( new_grp, &( *j ), 1 ); if( MB_SUCCESS != rval ) return rval; // add the new grp to the cgm_file_set rval = MBI->add_entities( cgm_file_set, &new_grp, 1 ); if( MB_SUCCESS != rval ) return rval; // remove the volume from the old group rval = MBI->remove_entities( *i, &( *j ), 1 ); if( MB_SUCCESS != rval ) return rval; if( debug ) std::cout << " new group: " << new_name << " id=" << max_grp_id << std::endl; } std::cout << " orig_volume=" << orig_grp_volume << " defo_volume=" << defo_grp_volume << " defo/orig=" << defo_grp_volume / orig_grp_volume << std::endl; } return MB_SUCCESS; }
static double tet_volume | ( | const CartVect & | v0, |
const CartVect & | v1, | ||
const CartVect & | v2, | ||
const CartVect & | v3 | ||
) | [inline, static] |
Definition at line 585 of file cub2h5m.cpp.
Referenced by measure().
{
return 1. / 6. * ( ( ( v1 - v0 ) * ( v2 - v0 ) ) % ( v3 - v0 ) );
}
void tokenize | ( | const std::string & | str, |
std::vector< std::string > & | tokens, | ||
const char * | delimiters | ||
) |
Definition at line 25 of file cub2h5m.cpp.
Referenced by summarize_cell_volume_change().
{ std::string::size_type last = str.find_first_not_of( delimiters, 0 ); std::string::size_type pos = str.find_first_of( delimiters, last ); if( std::string::npos == pos ) tokens.push_back( str ); else while( std::string::npos != pos && std::string::npos != last ) { tokens.push_back( str.substr( last, pos - last ) ); last = str.find_first_not_of( delimiters, pos ); pos = str.find_first_of( delimiters, last ); if( std::string::npos == pos ) pos = str.size(); } }
double DEFAULT_DISTANCE = 0.001 |
Definition at line 1103 of file cub2h5m.cpp.
double DEFAULT_LEN = 0.0 |
Definition at line 1104 of file cub2h5m.cpp.
int DEFAULT_NORM = 5 |
Definition at line 1105 of file cub2h5m.cpp.
Referenced by moab::ReadCGM::set_options().