moab
moab::DagMC Class Reference

#include <DagMC.hpp>

List of all members.

Classes

class  RayHistory
 State object used in calls to ray_fire() More...

Public Member Functions

 ~DagMC ()
ErrorCode load_file (const char *cfile, const double facet_tolerance=0)
 Load a geometry description regardless of format.
ErrorCode load_existing_contents ()
ErrorCode init_OBBTree ()
 initialize the OBB tree structure for ray firing acceleration
ErrorCode ray_fire (const EntityHandle volume, const double ray_start[3], const double ray_dir[3], EntityHandle &next_surf, double &next_surf_dist, RayHistory *history=NULL, double dist_limit=0, int ray_orientation=1, OrientedBoxTreeTool::TrvStats *stats=NULL)
 find the next surface crossing from a given point in a given direction
ErrorCode point_in_volume (const EntityHandle volume, const double xyz[3], int &result, const double *uvw=NULL, const RayHistory *history=NULL)
 Test if a point is inside or outside a volume.
ErrorCode point_in_volume_slow (const EntityHandle volume, const double xyz[3], int &result)
 Robust test if a point is inside or outside a volume using unit sphere area method.
ErrorCode test_volume_boundary (const EntityHandle volume, const EntityHandle surface, const double xyz[3], const double uvw[3], int &result, const RayHistory *history=NULL)
 Given a ray starting at a surface of a volume, check whether the ray enters or exits the volume.
ErrorCode closest_to_location (EntityHandle volume, const double point[3], double &result)
 Find the distance to the point on the boundary of the volume closest to the test point.
ErrorCode measure_volume (EntityHandle volume, double &result)
ErrorCode measure_area (EntityHandle surface, double &result)
ErrorCode surface_sense (EntityHandle volume, int num_surfaces, const EntityHandle *surfaces, int *senses_out)
ErrorCode surface_sense (EntityHandle volume, EntityHandle surface, int &sense_out)
ErrorCode get_angle (EntityHandle surf, const double xyz[3], double angle[3], const RayHistory *history=NULL)
ErrorCode next_vol (EntityHandle surface, EntityHandle old_volume, EntityHandle &new_volume)
EntityHandle entity_by_index (int dimension, int index)
int id_by_index (int dimension, int index)
EntityHandle entity_by_id (int dimension, int id)
int index_by_handle (EntityHandle handle)
int get_entity_id (EntityHandle this_ent)
int num_entities (int dimension)
 get number of geometric sets corresponding to geometry of specified dimension
double overlap_thickness ()
double numerical_precision ()
double faceting_tolerance ()
bool use_CAD ()
void set_overlap_thickness (double new_overlap_thickness)
void set_numerical_precision (double new_precision)
void set_use_CAD (bool use_cad)
ErrorCode detect_available_props (std::vector< std::string > &keywords_out)
ErrorCode parse_properties (const std::vector< std::string > &keywords, const std::map< std::string, std::string > &synonyms=no_synonyms)
ErrorCode prop_value (EntityHandle eh, const std::string &prop, std::string &value)
ErrorCode prop_values (EntityHandle eh, const std::string &prop, std::vector< std::string > &value)
bool has_prop (EntityHandle eh, const std::string &prop)
ErrorCode get_all_prop_values (const std::string &prop, std::vector< std::string > &return_list)
ErrorCode entities_by_property (const std::string &prop, std::vector< EntityHandle > &return_list, int dimension=0, const std::string *value=NULL)
bool is_implicit_complement (EntityHandle volume)
Tag name_tag ()
Tag obb_tag ()
Tag geom_tag ()
Tag id_tag ()
Tag sense_tag ()
OrientedBoxTreeToolobb_tree ()
ErrorCode write_mesh (const char *ffile, const int flen)
ErrorCode getobb (EntityHandle volume, double minPt[3], double maxPt[3])
ErrorCode getobb (EntityHandle volume, double center[3], double axis1[3], double axis2[3], double axis3[3])
ErrorCode get_root (EntityHandle vol_or_surf, EntityHandle &root)
Interfacemoab_instance ()

Static Public Member Functions

static DagMCinstance (Interface *mb_impl=NULL)
static float version (std::string *version_string=NULL)
static unsigned int interface_revision ()

Private Types

typedef std::map< std::string,
std::string > 
prop_map

Private Member Functions

ErrorCode finish_loading ()
bool have_obb_tree ()
ErrorCode get_impl_compl ()
ErrorCode build_obbs (Range &surfs, Range &vols)
ErrorCode build_obb_impl_compl (Range &surfs)
ErrorCode CAD_ray_intersect (const double *point, const double *dir, const double huge_val, std::vector< double > &distances, std::vector< EntityHandle > &surfaces, double &len)
 pass the ray_intersection test to the solid modeling engine
ErrorCode boundary_case (EntityHandle volume, int &result, double u, double v, double w, EntityHandle facet, EntityHandle surface)
 determine the point membership when the point is effectively on the boundary
ErrorCode poly_solid_angle (EntityHandle face, const CartVect &point, double &area)
ErrorCode build_indices (Range &surfs, Range &vols)
void tokenize (const std::string &str, std::vector< std::string > &tokens, const char *delimiters) const
ErrorCode get_group_name (EntityHandle group_set, std::string &name)
ErrorCode parse_group_name (EntityHandle group_set, prop_map &result)
ErrorCode append_packed_string (Tag, EntityHandle, std::string &)
ErrorCode unpack_packed_string (Tag tag, EntityHandle eh, std::vector< std::string > &values)
std::vector< EntityHandle > & surf_handles ()
std::vector< EntityHandle > & vol_handles ()
std::vector< EntityHandle > & group_handles ()
Tag get_tag (const char *name, int size, TagType store, DataType type, const void *def_value=NULL, bool create_if_missing=true)
 DagMC (Interface *mb_impl)

Static Private Member Functions

static void create_instance (Interface *mb_impl=NULL)

Private Attributes

InterfacembImpl
OrientedBoxTreeTool obbTree
EntityHandle impl_compl_handle
Tag obbTag
Tag geomTag
Tag idTag
Tag nameTag
Tag senseTag
Tag facetingTolTag
std::vector< EntityHandleentHandles [5]
EntityHandle setOffset
std::vector< EntityHandlerootSets
std::vector< int > entIndices
std::vector< RefEntity * > geomEntities
std::map< std::string, Tagproperty_tagmap
char implComplName [NAME_TAG_SIZE]
double overlapThickness
double numericalPrecision
double facetingTolerance
double defaultFacetingTolerance
bool useCAD
bool have_cgm_geom
 true if user requested CAD-based ray firing
std::vector< double > distList
 true if CGM contains problem geometry; required for CAD-based ray firing.
std::vector< EntityHandleprevFacetList
std::vector< EntityHandlesurfList
std::vector< EntityHandlefacetList
std::vector< double > disList
std::vector< int > dirList
std::vector< EntityHandlesurList
std::vector< EntityHandlefacList
long long int n_pt_in_vol_calls
long long int n_ray_fire_calls

Static Private Attributes

static DagMCinstance_ = NULL
static const std::map
< std::string, std::string > 
no_synonyms

Detailed Description

Definition at line 32 of file DagMC.hpp.


Member Typedef Documentation

typedef std::map<std::string, std::string> moab::DagMC::prop_map [private]

Definition at line 469 of file DagMC.hpp.


Constructor & Destructor Documentation

moab::DagMC::DagMC ( Interface mb_impl) [private]

Definition at line 106 of file DagMC.cpp.

  : mbImpl(mb_impl), obbTree(mb_impl), have_cgm_geom(false),
    n_pt_in_vol_calls(0), n_ray_fire_calls(0)
{
    // This is the correct place to uniquely define default values for the dagmc settings
  overlapThickness = 0; // must be nonnegative 
  defaultFacetingTolerance = .001;
  numericalPrecision = .001; 
  useCAD = false;  
  
  memset( implComplName, 0, NAME_TAG_SIZE );
  strcpy( implComplName , "impl_complement" );

}

Member Function Documentation

ErrorCode moab::DagMC::append_packed_string ( Tag  tag,
EntityHandle  eh,
std::string &  new_string 
) [private]

Add a string value to a property tag for a given entity

Definition at line 1648 of file DagMC.cpp.

{
    // When properties have multiple values, the values are tagged in a single character array
    // with the different values separated by null characters
  ErrorCode rval;
  const void* p;
  const char* str;
  int len;
  rval = MBI->tag_get_by_ptr( tag, &eh, 1, &p, &len );
  if( rval == MB_TAG_NOT_FOUND ){
    // This is the first entry, and can be set directly
    p = new_string.c_str();
    return MBI->tag_clear_data( tag, &eh, 1, p, new_string.length()+1);
  }
  else if( rval != MB_SUCCESS ) return rval;
  else{ 
    str = static_cast<const char*>(p);
  }

  // append a new value for the property to the existing property string
  unsigned int tail_len = new_string.length() + 1;
  char* new_packed_string = new char[ len + tail_len ];
  memcpy( new_packed_string, str, len );
  memcpy( new_packed_string + len, new_string.c_str(), tail_len );

  int new_len = len + tail_len;
  p = new_packed_string;
  rval = MBI->tag_set_by_ptr( tag, &eh, 1, &p, &new_len );
  delete[] new_packed_string;
  return rval;
}
ErrorCode moab::DagMC::boundary_case ( EntityHandle  volume,
int &  result,
double  u,
double  v,
double  w,
EntityHandle  facet,
EntityHandle  surface 
) [private]

determine the point membership when the point is effectively on the boundary

Called by point_in_volume when the point is with tolerance of the boundary. Compares the ray direction with the surface normal to determine a volume membership.

Definition at line 1242 of file DagMC.cpp.

{
  ErrorCode rval;

  // test to see if uvx is provided
  if ( u <= 1.0 && v <= 1.0 && w <= 1.0 ) {

    const CartVect ray_vector(u, v, w);
    CartVect coords[3], normal(0.0);
    const EntityHandle *conn;
    int len, sense_out;
   
    rval = mbImpl->get_connectivity( facet, conn, len );
    assert( MB_SUCCESS == rval );
    if(MB_SUCCESS != rval) return rval;
    assert( 3 == len );
  
    rval = mbImpl->get_coords( conn, 3, coords[0].array() );
    assert(MB_SUCCESS == rval);
    if(MB_SUCCESS != rval) return rval;
   
    rval = surface_sense( volume, surface, sense_out );
    assert( MB_SUCCESS == rval);
    if(MB_SUCCESS != rval) return rval;

    coords[1] -= coords[0];
    coords[2] -= coords[0];
    normal = sense_out * (coords[1] * coords[2]);

    double sense = ray_vector % normal;

    if ( sense < 0.0 ) {
      result = 1;     // inside or entering
    } else  if ( sense > 0.0 ) {
      result = 0;     // outside or leaving
    } else  if ( sense == 0.0 ) {
      result = -1;    // tangent, therefore on boundary
    } else {
      result = -1;    // failure
      return MB_FAILURE;
    }
  
  // if uvw not provided, return on_boundary.
  } else {
    result = -1;      // on boundary
    return MB_SUCCESS;
  
  }

  return MB_SUCCESS;
}
ErrorCode moab::DagMC::build_indices ( Range surfs,
Range vols 
) [private]

build internal index vectors that speed up handle-by-id, etc.

Definition at line 1409 of file DagMC.cpp.

{
  ErrorCode rval = MB_SUCCESS;

    // surf/vol offsets are just first handles
  setOffset = std::min( *surfs.begin(), *vols.begin() );

    // max
  EntityHandle tmp_offset = std::max( surfs.back(), vols.back() );

    // set size
  rootSets.resize(tmp_offset - setOffset + 1);
  entIndices.resize(rootSets.size());

    // store surf/vol handles lists (surf/vol by index) and
    // index by handle lists
  surf_handles().resize( surfs.size() + 1 );
  std::vector<EntityHandle>::iterator iter = surf_handles().begin();
  *(iter++) = 0;
  std::copy( surfs.begin(), surfs.end(), iter );
  int idx = 1;
  for (Range::iterator rit = surfs.begin(); rit != surfs.end(); rit++)
    entIndices[*rit-setOffset] = idx++;
  
  vol_handles().resize( vols.size() + 1 );
  iter = vol_handles().begin();
  *(iter++) = 0;
  std::copy( vols.begin(), vols.end(), iter );

  idx = 1;
  int max_id = -1;
  for (Range::iterator rit = vols.begin(); rit != vols.end(); rit++)    {
    entIndices[*rit-setOffset] = idx++;

    if( *rit != impl_compl_handle ){
      int result=0;
      MBI->tag_get_data( idTag, &*rit, 1, &result );
      max_id = std::max( max_id, result ); 
    }
  }
    // assign ID to implicit complement
    // for consistency with earlier versions of DagMC, make sure it always has the highest ID
  max_id++;
  MBI->tag_set_data(idTag, &impl_compl_handle, 1, &max_id);

#ifdef CGM
  if ( have_cgm_geom ) {
    // TODO: this block should only execute if the user has explicitly requested useCAD for ray firing.
    // however, this function curently executes before we know if useCAD will be specified, so do it every time.

    geomEntities.resize(rootSets.size());
      // get geometry entities by id and cache in this vector
    std::vector<int> ids;
    ids.resize(surfs.size());
    rval = MBI->tag_get_data(id_tag(), surfs, &ids[0]);
    if (MB_SUCCESS != rval) return MB_FAILURE;
    int i = 0;
    Range::iterator rit = surfs.begin();
    for (; rit != surfs.end(); rit++, i++) {
      RefEntity *this_surf = GeometryQueryTool::instance()->
        get_ref_face(ids[i]);
      assert(NULL != this_surf);
      geomEntities[*rit - setOffset] = this_surf;
    }
    ids.resize(vols.size());
    rval = MBI->tag_get_data(id_tag(), vols, &ids[0]);
    if (MB_SUCCESS != rval) return MB_FAILURE;
    i = 0;
    rit = vols.begin();
    for (; rit != vols.end(); rit++, i++) {
      if( is_implicit_complement( *rit ) ) continue;
      RefEntity *this_vol = GeometryQueryTool::instance()->
        get_ref_volume(ids[i]);
      assert(NULL != this_vol);
      geomEntities[*rit - setOffset] = this_vol;
    }
  }
#endif  

    // get group handles
  Tag category_tag = get_tag(CATEGORY_TAG_NAME, CATEGORY_TAG_LENGTH, 
                               MB_TAG_SPARSE, MB_TYPE_OPAQUE);
  char group_category[CATEGORY_TAG_SIZE];
  std::fill(group_category, group_category+CATEGORY_TAG_SIZE, '\0');
  sprintf(group_category, "%s", "Group");
  const void* const group_val[] = {&group_category};
  Range groups;
  rval = MBI->get_entities_by_type_and_tag(0, MBENTITYSET, &category_tag, 
                                           group_val, 1, groups);
  if (MB_SUCCESS != rval)
    return rval;
  group_handles().resize(groups.size()+1);
  group_handles()[0] = 0;
  std::copy(groups.begin(), groups.end(), &group_handles()[1]);

    // populate root sets vector
  std::vector<EntityHandle> rsets;
  rsets.resize(surfs.size());
  rval = MBI->tag_get_data(obb_tag(), surfs, &rsets[0]);
  if (MB_SUCCESS != rval) return MB_FAILURE;
  Range::iterator rit;
  int i;
  for (i = 0, rit = surfs.begin(); rit != surfs.end(); rit++, i++)
    rootSets[*rit-setOffset] = rsets[i];

  rsets.resize(vols.size());
  rval = MBI->tag_get_data(obb_tag(), vols, &rsets[0]);
  if (MB_SUCCESS != rval) return MB_FAILURE;
  for (i = 0, rit = vols.begin(); rit != vols.end(); rit++, i++)
    rootSets[*rit-setOffset] = rsets[i];

  return MB_SUCCESS;
}

build obb structure for the implicit complement

Definition at line 421 of file DagMC.cpp.

{
  EntityHandle comp_root, surf_obb_root;
  Range comp_tree;
  ErrorCode rval;
  std::vector<EntityHandle> parent_vols;
  
  int impl_compl_surf_count = 0;
  double impl_compl_surf_area = 0.0;

    // search through all surfaces  
  for (Range::iterator surf_i = surfs.begin(); surf_i != surfs.end(); ++surf_i) {
    
    parent_vols.clear();
      // get parents of each surface
    rval = MBI->get_parent_meshsets( *surf_i, parent_vols );
    if (MB_SUCCESS != rval)
      return rval;

      // if only one parent, get the OBB root for this surface
    if (parent_vols.size() == 1 ) {

      double a;
      measure_area( *surf_i, a );
      impl_compl_surf_count += 1;
      impl_compl_surf_area  += a;

      rval = MBI->tag_get_data( obbTag, &*surf_i, 1, &surf_obb_root );
      if (MB_SUCCESS != rval)
        return rval;
      if (!surf_obb_root)
        return MB_FAILURE;
      
        // add obb root to list of obb roots
      comp_tree.insert( surf_obb_root );

      // add this surf to the topology of the implicit complement volume
      rval = MBI->add_parent_child(impl_compl_handle,*surf_i);
      if (MB_SUCCESS != rval)
        return rval;

      // get the surface sense wrt original volume
      EntityHandle sense_data[2] = {0,0};
      rval = MBI->tag_get_data( sense_tag(), &(*surf_i), 1, sense_data );
      if (MB_SUCCESS != rval) return rval;
      
      // set the surface sense wrt implicit complement volume
      if(0==sense_data[0] && 0==sense_data[1]) return MB_FAILURE;
      if(0==sense_data[0])
        sense_data[0] = impl_compl_handle;
      else if(0==sense_data[1])
        sense_data[1] = impl_compl_handle;
      else
        return MB_FAILURE;
      rval = MBI->tag_set_data( sense_tag(), &(*surf_i), 1, sense_data );
      if (MB_SUCCESS != rval)  return rval;

    }  
  }
  // print info about the implicit complement if one was created
  if( impl_compl_surf_count ){
    bool one = (impl_compl_surf_count == 1);
    std::cout << "The implicit complement bounds " << impl_compl_surf_count
              << (one ? " surface" : " surfaces") << std::endl;
    std::cout << "The implicit complement's total surface area = " 
              << impl_compl_surf_area << std::endl;
  }
  
    // join surface trees to make OBB tree for implicit complement
  rval = obbTree.join_trees( comp_tree, comp_root );
  if (MB_SUCCESS != rval)
    return rval;
  
    // tag the implicit complement handle with the handle for its own OBB tree
  rval = MBI->tag_set_data( obbTag, &impl_compl_handle, 1, &comp_root );
  if (MB_SUCCESS != rval)
    return rval;
  
  // following ReadCGM, assign dimension and category tags
  int three = 3;
  rval = MBI->tag_set_data(geomTag, &impl_compl_handle, 1, &three );
  if (MB_SUCCESS != rval)
    return rval;

  Tag category_tag = get_tag(CATEGORY_TAG_NAME, CATEGORY_TAG_LENGTH, 
                             MB_TAG_SPARSE, MB_TYPE_OPAQUE);
  static const char volume_category[CATEGORY_TAG_SIZE] = "Volume\0";
  rval = MBI->tag_set_data(category_tag, &impl_compl_handle, 1, volume_category );
  if (MB_SUCCESS != rval)
    return rval;

  return MB_SUCCESS;

}
ErrorCode moab::DagMC::build_obbs ( Range surfs,
Range vols 
) [private]

build obb structure for each surface and volume

Definition at line 350 of file DagMC.cpp.

{
  ErrorCode rval = MB_SUCCESS;
  
  for (Range::iterator i = surfs.begin(); i != surfs.end(); ++i) {
    EntityHandle root;
    Range tris;
    rval = MBI->get_entities_by_dimension( *i, 2, tris );
    if (MB_SUCCESS != rval) 
      return rval;
    if (tris.empty()) 
      std::cerr << "WARNING: Surface " << get_entity_id(*i) << " has no facets." << std::endl;
    rval = obbTree.build( tris, root );
    if (MB_SUCCESS != rval) 
      return rval;
    rval = MBI->add_entities( root, &*i, 1 );
    if (MB_SUCCESS != rval)
      return rval;
    rval = MBI->tag_set_data( obbTag, &*i, 1, &root );
    if (MB_SUCCESS != rval)
      return rval;
  }
  
  for (Range::iterator i = vols.begin(); i != vols.end(); ++i) {
      // get all surfaces in volume
    Range tmp_surfs;
    rval = MBI->get_child_meshsets( *i, tmp_surfs );
    if (MB_SUCCESS != rval)
      return rval;
    
      // get OBB trees for each surface
    EntityHandle root;
    Range trees;
    for (Range::iterator j = tmp_surfs.begin();  j != tmp_surfs.end(); ++j) {
        // skip any surfaces that are non-manifold in the volume
        // because point containment code will get confused by them
      int sense;
      rval = surface_sense( *i, *j, sense );
      if (MB_SUCCESS != rval) {
        std::cerr << "Surface/Volume sense data missing." << std::endl;
        return rval;
      }
      if (!sense)
        continue;
      
      rval = MBI->tag_get_data( obbTag, &*j, 1, &root );
      if (MB_SUCCESS != rval || !root)  
        return MB_FAILURE;
      trees.insert( root );
    }
    
      // build OBB tree for volume
    rval = obbTree.join_trees( trees, root );
    if (MB_SUCCESS != rval)
      return rval;
    
    rval = MBI->tag_set_data( obbTag, &*i, 1, &root );
    if (MB_SUCCESS != rval)
      return rval;

  }

  rval = build_obb_impl_compl(surfs);
  if (MB_SUCCESS != rval) {
    std::cerr << "Unable to build OBB tree for implicit complement." << std::endl;
    return rval;
  }

  return MB_SUCCESS;
}
ErrorCode moab::DagMC::CAD_ray_intersect ( const double *  point,
const double *  dir,
const double  huge_val,
std::vector< double > &  distances,
std::vector< EntityHandle > &  surfaces,
double &  len 
) [private]

pass the ray_intersection test to the solid modeling engine

The user has the options to specify that ray tracing should ultimately occur on the true CAD model rather than just on the faceted representation. This is called from within ray_fire if the user has selected that option

Definition at line 1166 of file DagMC.cpp.

{
#ifdef CGM
#ifndef HAVE_CGM_FIRE_RAY
  return MB_NOT_IMPLEMENTED;
#else
  std::vector<double>::iterator dit = distances.begin();
  std::vector<EntityHandle>::iterator sit = surfaces.begin();
  static DLIList<double> ray_params;
  
  for (; dit != distances.end(); dit++, sit++) {
      // get the RefFace
    RefEntity *this_face = geomEntities[*sit - setOffset];
      // get the ray distance to this face
    ray_params.clean_out();
    int result = GeometryQueryTool::instance()->
      fire_ray(dynamic_cast<RefFace*>(this_face), CubitVector(point),
               CubitVector(dir), ray_params);
    assert(CUBIT_SUCCESS == result);
    if(CUBIT_SUCCESS != result) return MB_FAILURE;
    if (ray_params.size() != 0) {
      ray_params.reset();
      *dit = ray_params.get();
    }
    else *dit = huge_val;
  }
  
    // now bubble sort list
  bool done = false;
  while (!done) {
    dit = distances.begin();
    sit = surfaces.begin();
    done = true;
    for (; dit != distances.end(); dit++, sit++) {
      if (dit+1 != distances.end() && *dit > *(dit+1)) {
        double tmp_dist = *dit;
        *dit = *(dit+1);
        *(dit+1) = tmp_dist;
        EntityHandle tmp_hand = *sit;
        *sit = *(sit+1);
        *(sit+1) = tmp_hand;
        done = false;
      }
    }
  }

  if (!distances.empty()) len = distances[0];
  
  return MB_SUCCESS;
#endif
#else
  return MB_FAILURE;
#endif
}
ErrorCode moab::DagMC::closest_to_location ( EntityHandle  volume,
const double  point[3],
double &  result 
)

Find the distance to the point on the boundary of the volume closest to the test point.

Parameters:
volumeVolume to query
pointCoordinates of test point
resultSet to the minimum distance from point to a surface in volume

Definition at line 911 of file DagMC.cpp.

{
    // Get OBB Tree for volume
  assert(volume - setOffset < rootSets.size());
  EntityHandle root = rootSets[volume - setOffset];

    // Get closest triangles in volume
  const CartVect point(coords);
  CartVect nearest;
  EntityHandle facet_out;
  ErrorCode rval = obbTree.closest_to_location( point.array(), root, nearest.array(), facet_out );
  if (MB_SUCCESS != rval) return rval;

  // calculate distance between point and nearest facet
  result = (point-nearest).length();
  
  return MB_SUCCESS;

}
void moab::DagMC::create_instance ( Interface mb_impl = NULL) [static, private]

Definition at line 82 of file DagMC.cpp.

{
  if (NULL == mb_impl) mb_impl = new Core();
  instance_ = new DagMC(mb_impl);
}
ErrorCode moab::DagMC::detect_available_props ( std::vector< std::string > &  keywords_out)

Detect all the property keywords that appear in the loaded geometry

Parameters:
keywords_outThe result list of keywords. This list could be validly passed to parse_properties().

Definition at line 1626 of file DagMC.cpp.

{
  ErrorCode rval;
  std::set< std::string > keywords;
  for( std::vector<EntityHandle>::const_iterator grp=group_handles().begin();
       grp != group_handles().end(); ++grp )
  {
    std::map< std::string, std::string > properties;
    rval = parse_group_name( *grp, properties );
    if( rval == MB_TAG_NOT_FOUND ) continue;
    else if( rval != MB_SUCCESS ) return rval;

    for( prop_map::iterator i = properties.begin();
         i != properties.end(); ++i )
    {
      keywords.insert( (*i).first );
    }
  }
  keywords_list.assign( keywords.begin(), keywords.end() );
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::entities_by_property ( const std::string &  prop,
std::vector< EntityHandle > &  return_list,
int  dimension = 0,
const std::string *  value = NULL 
)

Get a list of all entities which have a given property

Parameters:
propThe canonical property name
return_listOutput param, a list of entity handles that have this property
dimensionIf nonzero, entities returned will be restricted to the given dimension, i.e. 2 for surfaces and 3 for volumes value If non-NULL, only entities for which the property takes on this value will be returned.
Returns:
MB_TAG_NOT_FOUND if prop is invalid. Otherwise return any errors from MOAB, or MB_SUCCESS if succesful

Definition at line 1846 of file DagMC.cpp.

{
  ErrorCode rval;
  std::map<std::string, Tag>::iterator it = property_tagmap.find(prop);
  if( it == property_tagmap.end() ){
    return MB_TAG_NOT_FOUND;
  }
  Tag proptag = (*it).second;
  Range all_ents;

  // Note that we cannot specify values for proptag here-- the passed value,
  // if it exists, may be only a subset of the packed string representation
  // of this tag.
  Tag tags[2] = {proptag, geomTag };
  void* vals[2] = {NULL, (dimension!=0) ? &dimension : NULL };
  rval = MBI->get_entities_by_type_and_tag( 0, MBENTITYSET, tags, vals, 2, all_ents );
  if( MB_SUCCESS != rval ) return rval;

  std::set<EntityHandle> handles;
  for( Range::iterator i = all_ents.begin(); i!=all_ents.end(); ++i){
    std::vector<std::string> values;
    rval = prop_values( *i, prop, values );
    if( MB_SUCCESS != rval ) return rval;
    if( value ){
      if( std::find(values.begin(), values.end(), *value) != values.end() ){
        handles.insert(*i);
      }
    }
    else{
      handles.insert(*i);
    }
  }

  return_list.assign( handles.begin(), handles.end() );
  return MB_SUCCESS;
}
EntityHandle moab::DagMC::entity_by_id ( int  dimension,
int  id 
)

map from dimension & global ID to EntityHandle

Definition at line 1362 of file DagMC.cpp.

{
  assert(0 <= dimension && 3 >= dimension);
  const Tag tags[] = { idTag, geomTag };
  const void* const vals[] = { &id, &dimension };
  ErrorCode rval;
  
  Range results;
  rval = MBI->get_entities_by_type_and_tag( 0, MBENTITYSET, tags, vals, 2, results );

  if ( MB_SUCCESS != rval ) 
      return 0;

  if ( results.empty() ){
    // old versions of dagmc did not set tags correctly on the implicit complement 'volume',
    // causing it to not be found by the call above.  This check allows this function to work
    // correctly, even on reloaded files from older versions.
    if( dimension == 3 && get_entity_id(impl_compl_handle) == id )
      return impl_compl_handle;
    else 
      return 0; 
  }
  
  return results.front();
}
EntityHandle moab::DagMC::entity_by_index ( int  dimension,
int  index 
) [inline]

map from dimension & base-1 ordinal index to EntityHandle

Definition at line 573 of file DagMC.hpp.

{
  assert(2 <= dimension && 3 >= dimension && (unsigned) index < entHandles[dimension].size());
  return entHandles[dimension][index];
}
double moab::DagMC::faceting_tolerance ( ) [inline]

retrieve faceting tolerance

Definition at line 359 of file DagMC.hpp.

loading code shared by load_file and load_existing_contents

Definition at line 190 of file DagMC.cpp.

{

  ErrorCode rval;


  nameTag = get_tag(NAME_TAG_NAME, NAME_TAG_SIZE, MB_TAG_SPARSE, MB_TYPE_OPAQUE, NULL, false);
  
  idTag = get_tag( GLOBAL_ID_TAG_NAME, 1, MB_TAG_DENSE, MB_TYPE_INTEGER );
  
  geomTag = get_tag( GEOM_DIMENSION_TAG_NAME, 1, MB_TAG_DENSE, MB_TYPE_INTEGER );

  obbTag = get_tag( MB_OBB_TREE_TAG_NAME, 1, MB_TAG_DENSE, MB_TYPE_HANDLE );

  facetingTolTag = get_tag(FACETING_TOL_TAG_NAME, 1, MB_TAG_SPARSE, MB_TYPE_DOUBLE );
  
    // get sense of surfaces wrt volumes
  senseTag = get_tag( "GEOM_SENSE_2", 2, MB_TAG_SPARSE, MB_TYPE_HANDLE );

  // search for a tag that has the faceting tolerance
  Range tagged_sets;
  double facet_tol_tagvalue = 0;
  bool other_set_tagged = false, root_tagged = false;
  
  // get list of entity sets that are tagged with faceting tolerance 
  // (possibly empty set)
  rval = MBI->get_entities_by_type_and_tag( 0, MBENTITYSET, &facetingTolTag,
                                            NULL, 1, tagged_sets );
  // if NOT empty set
  if (MB_SUCCESS == rval && !tagged_sets.empty()) {
    rval = MBI->tag_get_data( facetingTolTag, &(*tagged_sets.begin()), 1, &facet_tol_tagvalue );
    if (MB_SUCCESS != rval) return rval;
    other_set_tagged = true;
  }
  else if (MB_SUCCESS == rval) {
    // check to see if interface is tagged
    EntityHandle root = 0;
    rval = MBI->tag_get_data( facetingTolTag, &root, 1, &facet_tol_tagvalue );
    if (MB_SUCCESS == rval) root_tagged = true;
    else rval = MB_SUCCESS;
  }
  
  if ( (root_tagged || other_set_tagged) && facet_tol_tagvalue > 0) {
    facetingTolerance = facet_tol_tagvalue;
  }

  std::cout << "Using faceting tolerance: " << facetingTolerance << std::endl;

  return MB_SUCCESS;

}
Tag moab::DagMC::geom_tag ( ) [inline]

Definition at line 458 of file DagMC.hpp.

{return geomTag;}
ErrorCode moab::DagMC::get_all_prop_values ( const std::string &  prop,
std::vector< std::string > &  return_list 
)

Get a list of all unique values assigned to a named property on any entity

Parameters:
propThe canonical property name
return_listOutput param, a list of unique strings that are set as values for this property
Returns:
MB_TAG_NOT_FOUND if prop is invalid. Otherwise return any errors from MOAB, or MB_SUCCESS if succesful

Definition at line 1821 of file DagMC.cpp.

{
  ErrorCode rval;
  std::map<std::string, Tag>::iterator it = property_tagmap.find(prop);
  if( it == property_tagmap.end() ){
    return MB_TAG_NOT_FOUND;
  }
  Tag proptag = (*it).second;
  Range all_ents;

  rval = MBI->get_entities_by_type_and_tag( 0, MBENTITYSET, &proptag, NULL, 1, all_ents );
  if( MB_SUCCESS != rval ) return rval;

  std::set<std::string> unique_values;
  for( Range::iterator i = all_ents.begin(); i!= all_ents.end(); ++i){
    std::vector<std::string> values;
    rval = prop_values( *i, prop, values );
    if( MB_SUCCESS != rval ) return rval;
    unique_values.insert( values.begin(), values.end() );
  }

  return_list.assign( unique_values.begin(), unique_values.end() );
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::get_angle ( EntityHandle  surf,
const double  xyz[3],
double  angle[3],
const RayHistory history = NULL 
)

Get the normal to a given surface at the point on the surface closest to a given point

Parameters:
surfSurface on which to get normal
xyzPoint on surf
angleSet to coordinates of surface normal nearest xyz
historyOptional ray history from a previous call to ray_fire(). If present and non-empty, return the normal of the most recently intersected facet, ignoring xyz.

Definition at line 1099 of file DagMC.cpp.

{
  EntityHandle root = rootSets[surf - setOffset];
  ErrorCode rval;
  
  std::vector<EntityHandle> facets;

  // if no history or history empty, use nearby facets
  if( !history || (history->prev_facets.size() == 0) ){
    rval = obbTree.closest_to_location( in_pt, root, numericalPrecision, facets );
    assert(MB_SUCCESS == rval);
    if (MB_SUCCESS != rval) return rval;
  }
  // otherwise use most recent facet in history
  else{
    facets.push_back( history->prev_facets.back() );
  }

  CartVect coords[3], normal(0.0);
  const EntityHandle *conn;
  int len;
  for (unsigned i = 0; i < facets.size(); ++i) {
    rval = mbImpl->get_connectivity( facets[i], conn, len );
    assert( MB_SUCCESS == rval );
    assert( 3 == len );
  
    rval = mbImpl->get_coords( conn, 3, coords[0].array() );
    assert(MB_SUCCESS == rval);
    
    coords[1] -= coords[0];
    coords[2] -= coords[0];
    normal += coords[1] * coords[2];
  }
  
  normal.normalize();
  normal.get( angle );

  return MB_SUCCESS;
}

map from EntityHandle to global ID

Definition at line 1399 of file DagMC.cpp.

{
  int id = 0;
  ErrorCode result = MBI->tag_get_data(idTag, &this_ent, 1, &id);
  if (MB_TAG_NOT_FOUND == result)
    id = MBI->id_from_handle(this_ent);
    
  return id;
}
ErrorCode moab::DagMC::get_group_name ( EntityHandle  group_set,
std::string &  name 
) [private]

Store the name of a group in a string

Definition at line 1593 of file DagMC.cpp.

{
  ErrorCode rval;
  const void* v = NULL;
  int ignored;
  rval = MBI->tag_get_by_ptr(name_tag(), &group_set, 1, &v, &ignored);
  if( MB_SUCCESS != rval ) return rval;
  name = static_cast<const char*>(v);
  return MB_SUCCESS;
}

test for pre-existing implicit complement definition, or return a new one

Definition at line 307 of file DagMC.cpp.

{
  Range entities;
  const void* const tagdata[] = {implComplName};
  ErrorCode rval = mbImpl->get_entities_by_type_and_tag( 0, MBENTITYSET,
                                                           &nameTag, tagdata, 1,
                                                           entities );
  // query error
  if (MB_SUCCESS != rval) {
    std::cerr << "Unable to query for implicit complement." << std::endl;
    return rval;
  }

  // found too many
  if (entities.size() > 1) {
    std::cerr << "Too many implicit complement sets." << std::endl;
    return MB_MULTIPLE_ENTITIES_FOUND;
  }

  // found none
  if (entities.empty()) {
    rval= MBI->create_meshset(MESHSET_SET,impl_compl_handle);
    if (MB_SUCCESS != rval) {
      std::cerr << "Failed to create mesh set for implicit complement." << std::endl;
      return rval;
    }
      // tag this entity with name for implicit complement
    rval = MBI->tag_set_data(nameTag,&impl_compl_handle,1,&implComplName);
    if (MB_SUCCESS != rval) {
      std::cerr << "Failed to tag new entity as implicit complement." << std::endl;
    }

    return rval;

  } else {
    // found a single implicit complement
    impl_compl_handle = entities.front();
    return MB_SUCCESS;
  }
  
  
}
ErrorCode moab::DagMC::get_root ( EntityHandle  vol_or_surf,
EntityHandle root 
) [inline]

Definition at line 593 of file DagMC.hpp.

{
  unsigned int index = vol_or_surf - setOffset;
  root = (index < rootSets.size() ? rootSets[index] : 0);
  return (root ? MB_SUCCESS : MB_INDEX_OUT_OF_RANGE);
}
Tag moab::DagMC::get_tag ( const char *  name,
int  size,
TagType  store,
DataType  type,
const void *  def_value = NULL,
bool  create_if_missing = true 
) [private]

Definition at line 1907 of file DagMC.cpp.

{
  Tag retval = 0;
  unsigned flags = store|MB_TAG_CREAT;
  // NOTE: this function seens to be broken in that create_if_missing has
  // the opposite meaning from what its name implies.  However, changing the
  // behavior causes tests to fail, so I'm leaving the existing behavior
  // in place.  -- j.kraftcheck.
  if (!create_if_missing)  
    flags |= MB_TAG_EXCL;
  ErrorCode result = MBI->tag_get_handle(name, size, type, retval, flags, def_value);
  if (create_if_missing && MB_SUCCESS != result) 
    std::cerr << "Couldn't find nor create tag named " << name << std::endl;
  
  return retval;
}
ErrorCode moab::DagMC::getobb ( EntityHandle  volume,
double  minPt[3],
double  maxPt[3] 
)

Definition at line 1927 of file DagMC.cpp.

{
  double center[3], axis1[3], axis2[3], axis3[3];
 
    // get center point and vectors to OBB faces
  ErrorCode rval = getobb(volume, center, axis1, axis2, axis3);
  if (MB_SUCCESS != rval)
    return rval;
     
    // compute min and max verticies
  for (int i=0; i<3; i++) 
  {
    double sum = fabs(axis1[i]) + fabs(axis2[i]) + fabs(axis3[i]);
    minPt[i] = center[i] - sum;
    maxPt[i] = center[i] + sum;
  }
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::getobb ( EntityHandle  volume,
double  center[3],
double  axis1[3],
double  axis2[3],
double  axis3[3] 
)

Definition at line 1947 of file DagMC.cpp.

{
    //find EntityHandle node_set for use in box
  EntityHandle root = rootSets[volume - setOffset];
  
    // call box to get center and vectors to faces
  return obbTree.box(root, center, axis1, axis2, axis3);
  
}
std::vector<EntityHandle>& moab::DagMC::group_handles ( ) [inline, private]

Definition at line 483 of file DagMC.hpp.

{return entHandles[4];}
bool moab::DagMC::has_prop ( EntityHandle  eh,
const std::string &  prop 
)

Return true if a volume or surface has the named property set upon it

Parameters:
ehThe entity handle to query
propThe canonical property name True if the handle has the property set, or false if not. False is also returned if a MOAB error occurs.

Definition at line 1802 of file DagMC.cpp.

{
  ErrorCode rval;

  std::map<std::string, Tag>::iterator it = property_tagmap.find(prop);
  if( it == property_tagmap.end() ){
      return false;
  }

  Tag proptag = (*it).second;
  const void* data; 
  int ignored;

  rval = MBI->tag_get_by_ptr( proptag, &eh, 1, &data, &ignored );
  return ( rval == MB_SUCCESS );

}
bool moab::DagMC::have_obb_tree ( ) [private]

test for existing OBB Tree

Definition at line 297 of file DagMC.cpp.

{
  Range entities;
  ErrorCode rval = mbImpl->get_entities_by_type_and_tag( 0, MBENTITYSET,
                                                           &obbTag, 0, 1,
                                                           entities );
  return MB_SUCCESS == rval && !entities.empty();
}                                                    
int moab::DagMC::id_by_index ( int  dimension,
int  index 
)

map from dimension & base-1 ordinal index to global ID

Definition at line 1388 of file DagMC.cpp.

{
  EntityHandle h = entity_by_index( dimension, index );
  if (!h)
    return 0;
  
  int result = 0;
  MBI->tag_get_data( idTag, &h, 1, &result );
  return result;
}
Tag moab::DagMC::id_tag ( ) [inline]

Definition at line 459 of file DagMC.hpp.

{return idTag;}
int moab::DagMC::index_by_handle ( EntityHandle  handle) [inline]

PPHW: Missing dim & global ID ==> base-1 ordinal index map from EntityHandle to base-1 ordinal index

Definition at line 579 of file DagMC.hpp.

{
  assert(handle-setOffset < entIndices.size());
  return entIndices[handle-setOffset];
}

initialize the OBB tree structure for ray firing acceleration

This method generates an OBB tree from the faceted representation of the geometry. It also calls internal methods to generate the implicit complement and to build the cross-referencing indices.

Definition at line 242 of file DagMC.cpp.

{

  ErrorCode rval;

  Range surfs, vols;
  const int three = 3;
  const void* const three_val[] = {&three};
  rval = MBI->get_entities_by_type_and_tag( 0, MBENTITYSET, &geomTag, 
                                            three_val, 1, vols );
  if (MB_SUCCESS != rval)
    return rval;

  const int two = 2;
  const void* const two_val[] = {&two};
  rval = MBI->get_entities_by_type_and_tag( 0, MBENTITYSET, &geomTag, 
                                            two_val, 1, surfs );
  if (MB_SUCCESS != rval)
    return rval;

    // If it doesn't already exist, create implicit complement
    // Create data structures for implicit complement
  rval = get_impl_compl();
  if (MB_SUCCESS != rval) {
    std::cerr << "Failed to find or create implicit complement handle." << std::endl;
    return rval;
  }

  // Build OBB trees for everything, but only if we only read geometry
  // Changed to build obb tree if tree does not already exist. -- JK
  if (!have_obb_tree()) {
    rval = build_obbs(surfs, vols);
    if (MB_SUCCESS != rval) {
      std::cerr << "Failed to build obb." << std::endl;
      return rval;
    }
  }

  // build_indices expects the implicit complement to be in vols.
  if( vols.find(impl_compl_handle) == vols.end() ){
    vols.insert( vols.end(), impl_compl_handle );
  }

    // build the various index vectors used for efficiency
  rval = build_indices(surfs, vols);
  if (MB_SUCCESS != rval) {
    std::cerr << "Failed to build surface/volume indices." << std::endl;
    return rval;
  }
  
  return MB_SUCCESS;
}
DagMC * moab::DagMC::instance ( Interface mb_impl = NULL) [inline, static]

Definition at line 566 of file DagMC.hpp.

{
  if (NULL == instance_) create_instance(mb_impl);
  
  return instance_;
}
unsigned int moab::DagMC::interface_revision ( ) [static]

Get subversion revision of this file (DagMC.hpp)

Definition at line 95 of file DagMC.cpp.

{
  unsigned int result = 0;
  const char* interface_string = DAGMC_INTERFACE_REVISION; 
  if( strlen(interface_string) >= 5 ){
    // start looking for the revision number after "$Rev: " 
    result = strtol( interface_string+5, NULL, 10 ); 
  }
  return result;
}

Definition at line 1884 of file DagMC.cpp.

{
  return volume == impl_compl_handle;
}

Definition at line 185 of file DagMC.cpp.

                                        {
  
  return finish_loading();
}
ErrorCode moab::DagMC::load_file ( const char *  cfile,
const double  facet_tolerance = 0 
)

Load a geometry description regardless of format.

This method will load the geometry file with name cfile. In case this is a solid model geometry file, it will pass the facet_tolerance option as guidance for the faceting engine.

Parameters:
cfilethe file name to be loaded
facet_tolerancethe faceting tolerance guidance for the faceting engine
Returns:
- MB_SUCCESS if file loads correctly
  • other MB ErrorCodes returned from MOAB

Definition at line 126 of file DagMC.cpp.

{
  ErrorCode rval;

  std::cout << "Requested faceting tolerance: " << facet_tolerance << std::endl;
  
#ifdef CGM
  // cgm must be initialized so we can check it for CAD data after the load
  InitCGMA::initialize_cgma();
#endif

  facetingTolerance = defaultFacetingTolerance;
    // override default value of facetingTolerance with passed value
  if (facet_tolerance > 0 )
    facetingTolerance = facet_tolerance;

  char facetTolStr[16];

  sprintf(facetTolStr,"%g",facetingTolerance);

  char options[120] = "CGM_ATTRIBS=yes;FACET_DISTANCE_TOLERANCE=";
  strcat(options,facetTolStr);
  
  EntityHandle file_set;
  rval = MBI->create_meshset( MESHSET_SET, file_set );
  if (MB_SUCCESS != rval)
    return rval;

  rval = MBI->load_file(cfile, &file_set, options, NULL, 0, 0);
  
  if( MB_UNHANDLED_OPTION == rval ){
    // Some options were unhandled; this is common for loading h5m files.
    // Print a warning if an option was unhandled for a file that does not end in '.h5m'
    std::string filename(cfile);
    if( filename.length() < 4 || filename.substr(filename.length()-4) != ".h5m"){
      std::cerr << "DagMC warning: unhandled file loading options." << std::endl;
    }
  }
  else if (MB_SUCCESS != rval) {
    std::cerr << "DagMC Couldn't read file " << cfile << std::endl;
    std::string message;
    if (MB_SUCCESS == MBI->get_last_error(message) && !message.empty())
        std::cerr << "Error message: " << message << std::endl;
    
    return rval;
  }

#ifdef CGM  
  // check to see if CGM has data; if so, assume it corresponds to the data we loaded in.
  if( GeometryQueryTool::instance()->num_ref_volumes() > 0 ){
    have_cgm_geom = true;
  }
#endif

  return finish_loading();

}
ErrorCode moab::DagMC::measure_area ( EntityHandle  surface,
double &  result 
)

Calculate sum of area of triangles

Definition at line 1000 of file DagMC.cpp.

{
    // get triangles in surface
  Range triangles;
  ErrorCode rval = MBI->get_entities_by_dimension( surface, 2, triangles );
  if (MB_SUCCESS != rval) 
    return rval;
  if (!triangles.all_of_type(MBTRI)) {
    std::cout << "WARNING: Surface " << get_entity_id(surface)
              << " contains non-triangle elements. Area calculation may be incorrect." 
              << std::endl;
    triangles.clear();
    rval = MBI->get_entities_by_type( surface, MBTRI, triangles );
    if (MB_SUCCESS != rval) return rval;
  }

    // calculate sum of area of triangles
  result = 0.0;
  const EntityHandle *conn;
  int len;
  CartVect coords[3];
  for (Range::iterator j = triangles.begin(); j != triangles.end(); ++j) {
    rval = MBI->get_connectivity( *j, conn, len, true );
    if (MB_SUCCESS != rval) return rval;
    assert(3 == len);
    rval = MBI->get_coords( conn, 3, coords[0].array() );
    if (MB_SUCCESS != rval) return rval;

    coords[1] -= coords[0];
    coords[2] -= coords[0];
    coords[0] = coords[1] * coords[2];
    result += coords[0].length();
  }
  result *= 0.5;
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::measure_volume ( EntityHandle  volume,
double &  result 
)

Calculate the volume contained in a 'volume'

Definition at line 932 of file DagMC.cpp.

{
  ErrorCode rval;
  std::vector<EntityHandle> surfaces, surf_volumes;
  result = 0.0;
  
   // don't try to calculate volume of implicit complement
  if (volume == impl_compl_handle) {
    result = 1.0;
    return MB_SUCCESS;
  }

    // get surfaces from volume
  rval = MBI->get_child_meshsets( volume, surfaces );
  if (MB_SUCCESS != rval) return rval;
  
    // get surface senses
  std::vector<int> senses( surfaces.size() );
  rval = surface_sense( volume, surfaces.size(), &surfaces[0], &senses[0] );
  if (MB_SUCCESS != rval) {
    std::cerr << "ERROR: Surface-Volume relative sense not available. "
              << "Cannot calculate volume." << std::endl;
    return rval;
  }
  
  for (unsigned i = 0; i < surfaces.size(); ++i) {
      // skip non-manifold surfaces
    if (!senses[i])
      continue;
    
      // get triangles in surface
    Range triangles;
    rval = MBI->get_entities_by_dimension( surfaces[i], 2, triangles );
    if (MB_SUCCESS != rval) 
      return rval;
    if (!triangles.all_of_type(MBTRI)) {
      std::cout << "WARNING: Surface " << get_entity_id(surfaces[i])
                << " contains non-triangle elements. Volume calculation may be incorrect." 
                << std::endl;
      triangles.clear();
      rval = MBI->get_entities_by_type( surfaces[i], MBTRI, triangles );
      if (MB_SUCCESS != rval) return rval;
    }
    
      // calculate signed volume beneath surface (x 6.0)
    double surf_sum = 0.0;
    const EntityHandle *conn;
    int len;
    CartVect coords[3];
    for (Range::iterator j = triangles.begin(); j != triangles.end(); ++j) {
      rval = MBI->get_connectivity( *j, conn, len, true );
      if (MB_SUCCESS != rval) return rval;
      assert(3 == len);
      rval = MBI->get_coords( conn, 3, coords[0].array() );
      if (MB_SUCCESS != rval) return rval;
    
      coords[1] -= coords[0];
      coords[2] -= coords[0];
      surf_sum += (coords[0] % (coords[1] * coords[2]));
    }
    result += senses[i] * surf_sum;
  }
  
  result /= 6.0;
  return MB_SUCCESS;
}

Definition at line 506 of file DagMC.hpp.

{return mbImpl;}
Tag moab::DagMC::name_tag ( ) [inline]

get the tag for the "name" of a surface == global ID

Definition at line 454 of file DagMC.hpp.

{return nameTag;}
ErrorCode moab::DagMC::next_vol ( EntityHandle  surface,
EntityHandle  old_volume,
EntityHandle new_volume 
)

Get the volume on the other side of a surface

Parameters:
Asurface to query
old_volumeA volume on one side of surface
new_volumeOutput parameter for volume on the other side of surface
Returns:
MB_SUCCESS if new_volume was set successfully, error if not.

Definition at line 1139 of file DagMC.cpp.

{
  std::vector<EntityHandle> parents;
  ErrorCode rval = MBI->get_parent_meshsets( surface, parents );

  if (MB_SUCCESS == rval) {
    if (parents.size() != 2)
      rval = MB_FAILURE;
    else if (parents.front() == old_volume) 
      new_volume = parents.back();
    else if( parents.back() == old_volume )
      new_volume = parents.front();
    else 
      rval = MB_FAILURE;
  }

  if( rval != MB_SUCCESS ){
    std::cerr << "DAGMC: mesh error in next_vol for surf " << get_entity_id(surface) << std::endl;
  }

  return rval;

}
int moab::DagMC::num_entities ( int  dimension) [inline]

get number of geometric sets corresponding to geometry of specified dimension

For a given dimension (e.g. dimension=3 for volumes, dimension=2 for surfaces) return the number of entities of that dimension

Parameters:
dimensionthe dimensionality of the entities in question
Returns:
integer number of entities of that dimension

Definition at line 585 of file DagMC.hpp.

{
  assert(0 <= dimension && 3 >= dimension);
  
  return entHandles[dimension].size() - 1;
}
double moab::DagMC::numerical_precision ( ) [inline]

retrieve numerical precision

Definition at line 357 of file DagMC.hpp.

Tag moab::DagMC::obb_tag ( ) [inline]

Definition at line 457 of file DagMC.hpp.

{return obbTag;}

Definition at line 490 of file DagMC.hpp.

{return &obbTree;}
double moab::DagMC::overlap_thickness ( ) [inline]

retrieve overlap thickness

Definition at line 355 of file DagMC.hpp.

{return overlapThickness;}
ErrorCode moab::DagMC::parse_group_name ( EntityHandle  group_set,
prop_map result 
) [private]

Parse a group name into a set of key:value pairs

Definition at line 1604 of file DagMC.cpp.

{
  ErrorCode rval;
  std::string group_name;
  rval = get_group_name( group_set, group_name );
  if( rval != MB_SUCCESS ) return rval;

  std::vector< std::string > group_tokens;
  tokenize( group_name, group_tokens, "_" );

  // iterate over all the keyword positions 
  // keywords are even indices, their values (optional) are odd indices
  for( unsigned int i = 0; i < group_tokens.size(); i += 2 ){
    std::string groupkey = group_tokens[i];
    std::string groupval;
    if( i < group_tokens.size() - 1 )
      groupval = group_tokens[i+1];
    result[groupkey] = groupval;
  }
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::parse_properties ( const std::vector< std::string > &  keywords,
const std::map< std::string, std::string > &  synonyms = no_synonyms 
)

Parse properties from group names per metadata syntax standard

Parameters:
keywordsA list of keywords to parse. These are considered the canonical names of the properties, and constitute the valid inputs to has_prop() and prop_value().
synonymsAn optional mapping of synonym keywords to canonical keywords. This allows more than one group name keyword to take on the same meaning e.g. if synonyms["rest.of.world"] = "graveyard", then volumes in the "rest.of.world" group will behave as if they were in a group named "graveyard".

Definition at line 1700 of file DagMC.cpp.

{
  ErrorCode rval;

  // master keyword map, mapping user-set words in cubit to canonical property names
  std::map< std::string, std::string > keyword_map( keyword_synonyms );

  for( std::vector<std::string>::const_iterator i = keywords.begin();
       i != keywords.end(); ++i )
  {
    keyword_map[*i] = *i;
  }

  // the set of all canonical property names
  std::set< std::string > prop_names;
  for( prop_map::iterator i = keyword_map.begin();
       i != keyword_map.end(); ++i )
  {
    prop_names.insert((*i).second);
  }

  // set up DagMC's property tags based on what's been requested
  for( std::set<std::string>::iterator i = prop_names.begin();
       i != prop_names.end(); ++i )
  {
    std::string tagname("DAGMCPROP_");
    tagname += (*i);

    Tag new_tag;
    rval = MBI->tag_get_handle( tagname.c_str(), 0, MB_TYPE_OPAQUE, new_tag,
                                MB_TAG_SPARSE|MB_TAG_VARLEN|MB_TAG_CREAT );
    if( MB_SUCCESS != rval ) return rval;
    property_tagmap[(*i)] = new_tag;
  }

  // now that the keywords and tags are ready, iterate over all the actual geometry groups
  for( std::vector<EntityHandle>::iterator grp=group_handles().begin();
       grp != group_handles().end(); ++grp )
  {

    prop_map properties;
    rval = parse_group_name( *grp, properties );
    if( rval == MB_TAG_NOT_FOUND ) continue;
    else if( rval != MB_SUCCESS ) return rval;

    Range grp_sets;
    rval = MBI->get_entities_by_type( *grp, MBENTITYSET, grp_sets);
    if( MB_SUCCESS != rval ) return rval;
    if( grp_sets.size() == 0 ) continue;

    for( prop_map::iterator i = properties.begin(); 
         i != properties.end(); ++i )
    {
      std::string groupkey = (*i).first;
      std::string groupval = (*i).second;

      if( property_tagmap.find( groupkey ) != property_tagmap.end() ){
        Tag proptag = property_tagmap[groupkey];
        const unsigned int groupsize = grp_sets.size();
        for( unsigned int j = 0; j < groupsize; ++j){
            rval = append_packed_string( proptag, grp_sets[j], groupval );
        }
      }
    }
  }
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::point_in_volume ( const EntityHandle  volume,
const double  xyz[3],
int &  result,
const double *  uvw = NULL,
const RayHistory history = NULL 
)

Test if a point is inside or outside a volume.

This method finds the point on the boundary of the volume that is nearest the test point (x,y,z). If that point is "close" to a surface, a boundary test is performed based on the normal of the surface at that point and the optional ray direction (u,v,w).

Parameters:
volumeThe volume to test
xyzThe location to test for volume containment
resultSet to 0 if xyz it outside volume, 1 if inside, and -1 if on boundary.
Optionaldirection to use for underlying ray fire query. Used to ensure consistent results when a ray direction is known. If NULL or {0,0,0} is given, a random direction will be used.
historyOptional RayHistory object to pass to underlying ray fire query. The history is not modified by this call.

Definition at line 700 of file DagMC.cpp.

                                                            {
  // take some stats that are independent of nps
  if(counting) ++n_pt_in_vol_calls;

  // get OBB Tree for volume
  assert(volume - setOffset < rootSets.size());
  EntityHandle root = rootSets[volume - setOffset];

  // Don't recreate these every call. These cannot be the same as the ray_fire
  // vectors because both are used simultaneously.
  std::vector<double>       &dists = disList;
  std::vector<EntityHandle> &surfs = surList;
  std::vector<EntityHandle> &facets= facList;
  std::vector<int>          &dirs  = dirList;
  dists.clear();
  surfs.clear();
  facets.clear();
  dirs.clear();

  // if uvw is not given or is full of zeros, use a random direction
  double u = 0, v = 0, w = 0;

  if( uvw ){
    u = uvw[0]; v=uvw[1], w=uvw[2];
  }

  if( u == 0 && v == 0 && w == 0 )
  { 
    u = rand();
    v = rand();
    w = rand();
    const double magnitude = sqrt( u*u + v*v + w*w );
    u /= magnitude;
    v /= magnitude;
    w /= magnitude;
  }

  const double ray_direction[] = { u, v, w };
  
  // if overlaps, ray must be cast to infinity and all RTIs must be returned
  const double   large       = 1e15;
  const double   ray_length  = large;

  // If overlaps occur, the pt is inside if traveling along the ray from the
  // origin, there are ever more exits than entrances. In lieu of implementing
  // that, all intersections to infinity are required if overlaps occur (expensive)
  int min_tolerance_intersections;
  if(0 != overlapThickness) {
    min_tolerance_intersections = -1;
  // only the first intersection is needed if overlaps do not occur (cheap)
  } else {
    min_tolerance_intersections = 1;
  }

  // Get intersection(s) of forward and reverse orientation. Do not return 
  // glancing intersections or previous facets.
  ErrorCode rval = obbTree.ray_intersect_sets( dists, surfs, facets, root,
                                               numericalPrecision, 
                                               min_tolerance_intersections,
                                               xyz, ray_direction,
                                               &ray_length, NULL, NULL, &volume,
                                               &senseTag, NULL, 
                                               history ? &(history->prev_facets) : NULL );
  if(MB_SUCCESS != rval) return rval;

  // determine orientation of all intersections
  // 1 for entering, 0 for leaving, -1 for tangent
  // Tangent intersections are not returned from ray_tri_intersect.
  dirs.resize(dists.size());
  for(unsigned i=0; i<dists.size(); ++i) {
    rval = boundary_case( volume, dirs[i], u, v, w, facets[i], surfs[i] );
    if(MB_SUCCESS != rval) return rval;
  }

  // count all crossings
  if(0 != overlapThickness) {
    int sum = 0;
    for(unsigned i=0; i<dirs.size(); ++i) {
      if     ( 1==dirs[i]) sum+=1; // +1 for entering
      else if( 0==dirs[i]) sum-=1; // -1 for leaving
      else if(-1==dirs[i]) {       //  0 for tangent
        std::cout << "direction==tangent" << std::endl;
        sum+=0;
      } else {
        std::cout << "error: unknown direction" << std::endl;
        return MB_FAILURE;
      }
    }

    // inside/outside depends on the sum
    if(0<sum)                          result = 0; // pt is outside (for all vols)
    else if(0>sum)                     result = 1; // pt is inside  (for all vols)
    else if(impl_compl_handle==volume) result = 1; // pt is inside  (for impl_compl_vol)
    else                               result = 0; // pt is outside (for all other vols)

  // Only use the first crossing
  } else {
      if( dirs.empty() ) {
      result = 0; // pt is outside
    } else {
      int smallest = std::min_element( dists.begin(), dists.end() ) - dists.begin();
      if     ( 1==dirs[smallest] ) result = 0; // pt is outside
      else if( 0==dirs[smallest] ) result = 1; // pt is inside
      else if(-1==dirs[smallest] ) {
        // Should not be here because Plucker ray-triangle test does not 
        // return coplanar rays as intersections. 
        std::cout << "direction==tangent" << std::endl;
        result = -1;
      } else {
        std::cout << "error: unknown direction" << std::endl;
        return MB_FAILURE;
      }
    }
  }

  if(debug)
    std::cout << "pt_in_vol: result=" << result
              << " xyz=" << xyz[0] << " " << xyz[1] << " " << xyz[2] << " uvw=" << u << " " << v << " " << w
              << " vol_id=" << id_by_index(3, index_by_handle(volume)) << std::endl;
  
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::point_in_volume_slow ( const EntityHandle  volume,
const double  xyz[3],
int &  result 
)

Robust test if a point is inside or outside a volume using unit sphere area method.

This test may be more robust that the standard point_in_volume, but is much slower. It does not detect 'on boundary' situations as point_in_volume does.

Parameters:
volumeThe volume to test
xyzThe location to test for volume containment
resultSet to 0 if xyz it outside volume, 1 if inside.

Definition at line 865 of file DagMC.cpp.

{
  ErrorCode rval;
  Range faces;
  std::vector<EntityHandle> surfs;
  std::vector<int> senses;
  double sum = 0.0;
  const CartVect point(xyz);
  
  rval = MBI->get_child_meshsets( volume, surfs );
  if (MB_SUCCESS != rval)
    return rval;
  
  senses.resize( surfs.size() );
  rval = surface_sense( volume, surfs.size(), &surfs[0], &senses[0] );
  if (MB_SUCCESS != rval)
    return rval;
  
  for (unsigned i = 0; i < surfs.size(); ++i) {
    if (!senses[i])  // skip non-manifold surfaces
      continue;
     
    double surf_area = 0.0, face_area;
    faces.clear();
    rval = MBI->get_entities_by_dimension( surfs[i], 2, faces );
    if (MB_SUCCESS != rval)
      return rval;
    
    for (Range::iterator j = faces.begin(); j != faces.end(); ++j) {
      rval = poly_solid_angle( *j, point, face_area );
      if (MB_SUCCESS != rval)
        return rval;
      
      surf_area += face_area;
    }
    
    sum += senses[i] * surf_area;
  }
  
  result = fabs(sum) > 2.0*M_PI;
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::poly_solid_angle ( EntityHandle  face,
const CartVect point,
double &  area 
) [private]

get the solid angle projected by a facet on a unit sphere around a point

  • used by point_in_volume_slow

Definition at line 1305 of file DagMC.cpp.

{
  ErrorCode rval;
  
    // Get connectivity
  const EntityHandle* conn;
  int len;
  rval = MBI->get_connectivity( face, conn, len, true );
  if (MB_SUCCESS != rval)
    return rval;
  
    // Allocate space to store vertices
  CartVect coords_static[4];
  std::vector<CartVect> coords_dynamic;
  CartVect* coords = coords_static;
  if ((unsigned)len > (sizeof(coords_static)/sizeof(coords_static[0]))) {
    coords_dynamic.resize(len);
    coords = &coords_dynamic[0];
  }
  
    // get coordinates
  rval = MBI->get_coords( conn, len, coords->array() );
  if (MB_SUCCESS != rval)
    return rval;
  
    // calculate normal
  CartVect norm(0.0), v1, v0 = coords[1] - coords[0];
  for (int i = 2; i < len; ++i) {
    v1 = coords[i] - coords[0];
    norm += v0 * v1;
    v0 = v1;
  }
  
    // calculate area
  double s, ang;
  area = 0.0;
  CartVect r, n1, n2, b, a = coords[len-1] - coords[0];
  for (int i = 0; i < len; ++i) {
    r = coords[i] - point;
    b = coords[(i+1)%len] - coords[i];
    n1 = a * r; // = norm1 (magnitude is important)
    n2 = r * b; // = norm2 (magnitude is important)
    s = (n1 % n2) / (n1.length() * n2.length()); // = cos(angle between norm1,norm2)
    ang = s <= -1.0 ? M_PI : s >= 1.0 ? 0.0 : acos(s); // = acos(s)
    s = (b * a) % norm; // =orientation of triangle wrt point
    area += s > 0.0 ? M_PI - ang : M_PI + ang;
    a = -b;
  }
  
  area -= M_PI * (len - 2);
  if ((norm % r) > 0)
    area = -area;
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::prop_value ( EntityHandle  eh,
const std::string &  prop,
std::string &  value 
)

Get the value of a property on a volume or surface

Parameters:
ehThe entity handle to get a property value on
propThe canonical property name
valueOutput parameter, the value of the property. If no value was set on the handle, this will be the empty string.
Returns:
MB_TAG_NOT_FOUND if prop is invalid. Otherwise return any errors from MOAB, or MB_SUCCESS if successful

Definition at line 1769 of file DagMC.cpp.

{
  ErrorCode rval;

  std::map<std::string, Tag>::iterator it = property_tagmap.find(prop);
  if( it == property_tagmap.end() ){
      return MB_TAG_NOT_FOUND;
  }

  Tag proptag = (*it).second;
  const void* data; 
  int ignored;

  rval = MBI->tag_get_by_ptr( proptag, &eh, 1, &data, &ignored );
  if( rval != MB_SUCCESS ) return rval;
  value = static_cast<const char*>(data);
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::prop_values ( EntityHandle  eh,
const std::string &  prop,
std::vector< std::string > &  value 
)

Get the value of a property on a volume or surface

Parameters:
ehThe entity handle to get a property value on
propThe canonical property name
valuesOutput parameter, the values of the property will be appended to this list. If no value was set on the handle, no entries will be added.
Returns:
MB_TAG_NOT_FOUND if prop is invalid. Otherwise return any errors from MOAB, or MB_SUCCESS if successful

Definition at line 1788 of file DagMC.cpp.

{

  std::map<std::string, Tag>::iterator it = property_tagmap.find(prop);
  if( it == property_tagmap.end() ){
      return MB_TAG_NOT_FOUND;
  }
  Tag proptag = (*it).second;

  return unpack_packed_string( proptag, eh, values );

}
ErrorCode moab::DagMC::ray_fire ( const EntityHandle  volume,
const double  ray_start[3],
const double  ray_dir[3],
EntityHandle next_surf,
double &  next_surf_dist,
RayHistory history = NULL,
double  dist_limit = 0,
int  ray_orientation = 1,
OrientedBoxTreeTool::TrvStats stats = NULL 
)

find the next surface crossing from a given point in a given direction

This is the primary method of DagMC, enabling ray tracing through a geometry. Given a volume and a ray, it computes the surface ID and distance to the nearest intersection on that volume. The caller can compute the location of the intersection by adding the distance to the ray.

When a series of calls to this function are made along the same ray (e.g. for the purpose of tracking a ray through several volumes), the optional history argument should be given. The history prevents previously intersected facets from being intersected again. A single history should be used as long as a ray is proceeding forward without changing direction. This situation is sometimes referred to as "streaming."

If a ray changes direction at an intersection site, the caller should call reset_to_last_intersection() on the history object before the next ray fire.

Parameters:
volumeThe volume to fire the ray at.
ray_startAn array of x,y,z coordinates from which to start the ray.
ray_dirAn array of x,y,z coordinates indicating the direction of the ray. Must be of unit length.
next_surfOutput parameter indicating the next surface intersected by the ray. If no intersection is found, will be set to 0.
next_surf_distOutput parameter indicating distance to next_surf. If next_surf is 0, this value is undefined and should not be used.
historyOptional RayHistory object. If provided, the facets in the history are assumed to not intersect with the given ray. The facet intersected by this query will also be added to the history.
dist_limitOptional distance limit. If provided and > 0, no intersections at a distance further than this value will be returned.
ray_orientationOptional ray orientation. If provided determines intersections along the normal provided, e.g. if -1 allows intersections back along the the ray direction, Default is 1, i.e. exit intersections
statsOptional TrvStats object used to measure performance of underlying OBB ray-firing query. See OrientedBoxTreeTool.hpp for details.

Definition at line 535 of file DagMC.cpp.

                                                               { 

  // take some stats that are independent of nps
  if(counting) {
    ++n_ray_fire_calls;
    if(0==n_ray_fire_calls%10000000) {
      std::cout << "n_ray_fires="   << n_ray_fire_calls 
                << " n_pt_in_vols=" << n_pt_in_vol_calls << std::endl;
    }
  }

  if (debug) {
    std::cout << "ray_fire:" 
              << " xyz=" << point[0] << " " << point[1] << " " << point[2]
              << " uvw=" << dir[0] << " " << dir[1] << " " << dir[2]
              << " vol_id=" << id_by_index(3, index_by_handle(vol)) << std::endl;
    }

  const double huge_val = std::numeric_limits<double>::max();
  double dist_limit = huge_val;
  if( user_dist_limit > 0 )
    dist_limit = user_dist_limit;

  // don't recreate these every call
  std::vector<double>       &dists       = distList;
  std::vector<EntityHandle> &surfs       = surfList;
  std::vector<EntityHandle> &facets      = facetList;  
  dists.clear();
  surfs.clear();
  facets.clear();

  assert(vol - setOffset < rootSets.size());  
  const EntityHandle root = rootSets[vol - setOffset];
  ErrorCode rval;
 
  // check behind the ray origin for intersections
  double neg_ray_len;
  if(0 == overlapThickness) {
    neg_ray_len = -numericalPrecision;
  } else {
    neg_ray_len = -overlapThickness;
  }

  // optionally, limit the nonneg_ray_len with the distance to next collision.
  double nonneg_ray_len = dist_limit;

  // the nonneg_ray_len should not be less than -neg_ray_len, or an overlap 
  // may be missed due to optimization within ray_intersect_sets
  if(nonneg_ray_len < -neg_ray_len) nonneg_ray_len = -neg_ray_len;
  assert(0 <= nonneg_ray_len);
  assert(0 >     neg_ray_len);
  
  // min_tolerance_intersections is passed but not used in this call
  const int min_tolerance_intersections = 0;

  // numericalPrecision is used for box.intersect_ray and find triangles in the
  // neighborhood of edge/node intersections.
  rval = obbTree.ray_intersect_sets( dists, surfs, facets,
                                     root, numericalPrecision, 
                                     min_tolerance_intersections,
                                     point, dir, &nonneg_ray_len,
                                     stats, &neg_ray_len, &vol, &senseTag, 
                                     &ray_orientation, 
                                     history ? &(history->prev_facets) : NULL );
  assert( MB_SUCCESS == rval );
  if(MB_SUCCESS != rval) return rval;

  // if useCAD is true at this point, then we know we can call CGM's ray casting function.
  if (useCAD) {
    rval = CAD_ray_intersect( point, dir, huge_val, dists, surfs, nonneg_ray_len );
    if (MB_SUCCESS != rval) return rval;
  }
  
  // If no distances are returned, the particle is lost unless the physics limit
  // is being used. If the physics limit is being used, there is no way to tell
  // if the particle is lost. To avoid ambiguity, DO NOT use the distance limit 
  // unless you know lost particles do not occur.
  if( dists.empty() ) {
    next_surf = 0;
    if(debug) {
      std::cout << "          next_surf=0 dist=(undef)" << std::endl; 
    }
    return MB_SUCCESS;
  }

  // Assume that a (neg, nonneg) pair of RTIs could be returned,
  // however, only one or the other may exist. dists[] may be populated, but 
  // intersections are ONLY indicated by nonzero surfs[] and facets[].
  assert(2 == dists.size());
  assert(2 == facets.size());
  assert(0.0 >= dists[0]);
  assert(0.0 <= dists[1]);

  // If both negative and nonnegative RTIs are returned, the negative RTI must
  // closer to the origin.
  if(0!=facets[0] && 0!=facets[1]) {
    assert(-dists[0] <= dists[1]);
  }

  // If an RTI is found at negative distance, perform a PMT to see if the 
  // particle is inside an overlap.
  int exit_idx = -1;
  if(0!=facets[0]) {
    // get the next volume
    std::vector<EntityHandle> vols;
    EntityHandle nx_vol;
    rval = MBI->get_parent_meshsets( surfs[0], vols );
    if(MB_SUCCESS != rval) return rval;
    assert(2 == vols.size());
    if(vols.front() == vol) {
      nx_vol = vols.back();
    } else {
      nx_vol = vols.front();
    }
    // Check to see if the point is actually in the next volume.
    // The list of previous facets is used to topologically identify the 
    // "on_boundary" result of the PMT. This avoids a test that uses proximity 
    // (a tolerance).
    int result;
    rval = point_in_volume( nx_vol, point, result, dir, history );
    if(MB_SUCCESS != rval) return rval;
    if(1==result) exit_idx = 0;

  }

  // if the negative distance is not the exit, try the nonnegative distance
  if(-1==exit_idx && 0!=facets[1]) exit_idx = 1;
  
  // if the exit index is still unknown, the particle is lost
  if(-1 == exit_idx) {
    next_surf = 0;
    if (debug) {
      std::cout << "next surf hit = 0, dist = (undef)" << std::endl;
    }
    return MB_SUCCESS;
  }

  // return the intersection
  next_surf = surfs[exit_idx];
  next_surf_dist = ( 0>dists[exit_idx] ? 0 : dists[exit_idx]);

  if( history ){
    history->prev_facets.push_back( facets[exit_idx] );
  }

  if (debug) {
    if( 0 > dists[exit_idx] ){
      std::cout << "          OVERLAP track length=" << dists[exit_idx] << std::endl;
    }
    std::cout << "          next_surf = " <<  id_by_index(2, index_by_handle(next_surf)) 
              << ", dist = " << next_surf_dist << " new_pt=";
    for( int i = 0; i < 3; ++i ){
      std::cout << point[i]+dir[i]*next_surf_dist << " ";
    }
    std::cout << std::endl;
  }

  return MB_SUCCESS;
}

Definition at line 460 of file DagMC.hpp.

{ return senseTag; }
void moab::DagMC::set_numerical_precision ( double  new_precision)

Attempt to set a new numerical precision , first checking for sanity Use of this function is discouraged; see top of DagMC.cpp

Definition at line 1539 of file DagMC.cpp.

                                                         {

  if ( new_precision <= 0 || new_precision > 1) {
    std::cerr << "Invalid numerical_precision = " << numericalPrecision << std::endl;
  }
  else{
    numericalPrecision = new_precision;
  }

  std::cout << "Set numerical precision = " << numericalPrecision << std::endl;

}
void moab::DagMC::set_overlap_thickness ( double  new_overlap_thickness)

Attempt to set a new overlap thickness tolerance, first checking for sanity

Definition at line 1527 of file DagMC.cpp.

                                                       {

  if (new_thickness < 0 || new_thickness > 100) {
    std::cerr << "Invalid overlap_thickness = " << new_thickness << std::endl;
  }
  else{
    overlapThickness = new_thickness;
  }
  std::cout << "Set overlap thickness = " << overlapThickness << std::endl;

}
void moab::DagMC::set_use_CAD ( bool  use_cad)

attempt to set useCAD, first checking for availability

Definition at line 1552 of file DagMC.cpp.

                                     {
  useCAD = use_cad;
  if( useCAD ){
    if( !have_cgm_geom ){
      std::cerr << "Warning: CAD-based ray tracing not avaiable, because CGM has no data." << std::endl;
      std::cerr << "         your input file was probably not a CAD format." << std::endl;
      useCAD = false;
    }

#ifndef HAVE_CGM_FIRE_RAY
    {
      std::cerr << "Warning: use_cad = 1 not supported with this build of CGM/DagMC." << std:: endl;
      std::cerr << "         Required ray-fire query not available. (Cubit-based CGM?)" <<  std::endl;
      useCAD = false;
    }
#endif
  }

  std::cout << "Turned " << (useCAD?"ON":"OFF") << " ray firing on full CAD model." << std::endl;

}
std::vector<EntityHandle>& moab::DagMC::surf_handles ( ) [inline, private]

Definition at line 481 of file DagMC.hpp.

{return entHandles[2];}
ErrorCode moab::DagMC::surface_sense ( EntityHandle  volume,
int  num_surfaces,
const EntityHandle surfaces,
int *  senses_out 
)

Get the sense of surfaces wrt a volume. Sense values are: {-1 -> reversed, 0 -> both, 1 -> forward}

Definition at line 1038 of file DagMC.cpp.

{

  /* The sense tags do not reference the implicit complement handle.
     All surfaces that interact with the implicit complement should have
     a null handle in the direction of the implicit complement. */
  //if (volume == impl_compl_handle)
  //  volume = (EntityHandle) 0;

  std::vector<EntityHandle> surf_volumes( 2*num_surfaces );
  ErrorCode rval = MBI->tag_get_data( sense_tag(), surfaces, num_surfaces, &surf_volumes[0] );
  if (MB_SUCCESS != rval)  return rval;
  
  const EntityHandle* end = surfaces + num_surfaces;
  std::vector<EntityHandle>::const_iterator surf_vols = surf_volumes.begin();
  while (surfaces != end) {
    EntityHandle forward = *surf_vols; ++surf_vols;
    EntityHandle reverse = *surf_vols; ++surf_vols;
    if (volume == forward) 
      *senses_out = (volume != reverse); // zero if both, otherwise 1
    else if (volume == reverse)
      *senses_out = -1;
    else 
      return MB_ENTITY_NOT_FOUND;
    
    ++surfaces;
    ++senses_out;
  }
  
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::surface_sense ( EntityHandle  volume,
EntityHandle  surface,
int &  sense_out 
)

Get the sense of a single surface wrt a volume. Sense values are: {-1 -> reversed, 0 -> both, 1 -> forward}

Definition at line 1074 of file DagMC.cpp.

{
  /* The sense tags do not reference the implicit complement handle.
     All surfaces that interact with the implicit complement should have
     a null handle in the direction of the implicit complement. */
  //if (volume == impl_compl_handle)
  //  volume = (EntityHandle) 0;

    // get sense of surfaces wrt volumes
  EntityHandle surf_volumes[2];
  ErrorCode rval = MBI->tag_get_data( sense_tag(), &surface, 1, surf_volumes );
  if (MB_SUCCESS != rval)  return rval;
  
  if (surf_volumes[0] == volume)
    sense_out = (surf_volumes[1] != volume); // zero if both, otherwise 1
  else if (surf_volumes[1] == volume) 
    sense_out = -1;
  else
    return MB_ENTITY_NOT_FOUND;
  
  return MB_SUCCESS;
}
ErrorCode moab::DagMC::test_volume_boundary ( const EntityHandle  volume,
const EntityHandle  surface,
const double  xyz[3],
const double  uvw[3],
int &  result,
const RayHistory history = NULL 
)

Given a ray starting at a surface of a volume, check whether the ray enters or exits the volume.

This function is most useful for rays that change directions at a surface crossing. It can be used to check whether a direction change redirects the ray back into the originating volume.

Parameters:
volumeThe volume to test
surfaceA surface on volume
xyzA point location on surface
uvwA (unit) direction vector
resultSet to 1 if ray is entering volume, or 0 if it is leaving
historyOptional ray history object from a previous call to ray_fire. If present and non-empty, the history is used to look up the surface facet at which the ray begins. Absent a history, the facet nearest to xyz will be looked up. The history should always be provided if available, as it avoids the computational expense of a nearest-facet query.

Definition at line 827 of file DagMC.cpp.

{
  ErrorCode rval;
  int dir;

  if( history && history->prev_facets.size() ){
    // the current facet is already available
    rval = boundary_case( volume, dir, uvw[0], uvw[1], uvw[2], history->prev_facets.back(), surface );
    if (MB_SUCCESS != rval) return rval;
  }
  else{
    // look up nearest facet

    // Get OBB Tree for surface
    assert(volume - setOffset < rootSets.size());
    EntityHandle root = rootSets[volume - setOffset];
    
    // Get closest triangle on surface
    const CartVect point(xyz);
    CartVect nearest;
    EntityHandle facet_out;
    rval = obbTree.closest_to_location( point.array(), root, nearest.array(), facet_out );
    if (MB_SUCCESS != rval) return rval;

    rval = boundary_case( volume, dir, uvw[0], uvw[1], uvw[2], facet_out, surface );
    if (MB_SUCCESS != rval) return rval;

  }

  result = dir;

  return MB_SUCCESS;

}
void moab::DagMC::tokenize ( const std::string &  str,
std::vector< std::string > &  tokens,
const char *  delimiters 
) const [private]

tokenize the metadata stored in group names - basically borroed from ReadCGM.cpp

Definition at line 1889 of file DagMC.cpp.

{
  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();
    }
}
ErrorCode moab::DagMC::unpack_packed_string ( Tag  tag,
EntityHandle  eh,
std::vector< std::string > &  values 
) [private]

Convert a property tag's value on a handle to a list of strings

Definition at line 1681 of file DagMC.cpp.

{
  ErrorCode rval;
  const void* p;
  const char* str;
  int len;
  rval = MBI->tag_get_by_ptr( tag, &eh, 1, &p, &len );
  if( rval != MB_SUCCESS ) return rval;
  str = static_cast<const char*>(p);
  int idx = 0;
  while( idx < len ){
    std::string item(str + idx);
    values.push_back( item );
    idx += item.length() + 1;
  }
  return MB_SUCCESS;
}
bool moab::DagMC::use_CAD ( ) [inline]

retrieve use CAD toggle

Definition at line 361 of file DagMC.hpp.

{return useCAD;}
float moab::DagMC::version ( std::string *  version_string = NULL) [static]

Return the version of this library

Definition at line 88 of file DagMC.cpp.

{
  if (NULL != version_string)
    *version_string = std::string("DagMC version ") + std::string(DAGMC_VERSION_STRING);
  return DAGMC_VERSION;
}
std::vector<EntityHandle>& moab::DagMC::vol_handles ( ) [inline, private]

Definition at line 482 of file DagMC.hpp.

{return entHandles[3];}
ErrorCode moab::DagMC::write_mesh ( const char *  ffile,
const int  flen 
)

Definition at line 1574 of file DagMC.cpp.

{
  ErrorCode rval;
  
    // write out a mesh file if requested
  if (ffile && 0 < flen) {
    rval = MBI->write_mesh(ffile);
    if (MB_SUCCESS != rval) {
      std::cerr << "Failed to write mesh to " << ffile << "." << std::endl;
      return rval;
    }
  }
 
  return MB_SUCCESS;
}

Member Data Documentation

Definition at line 548 of file DagMC.hpp.

std::vector<int> moab::DagMC::dirList [private]

Definition at line 558 of file DagMC.hpp.

std::vector<double> moab::DagMC::disList [private]

Definition at line 557 of file DagMC.hpp.

std::vector<double> moab::DagMC::distList [private]

true if CGM contains problem geometry; required for CAD-based ray firing.

Definition at line 554 of file DagMC.hpp.

std::vector<EntityHandle> moab::DagMC::entHandles[5] [private]

Definition at line 524 of file DagMC.hpp.

std::vector<int> moab::DagMC::entIndices [private]

Definition at line 533 of file DagMC.hpp.

Definition at line 548 of file DagMC.hpp.

Definition at line 522 of file DagMC.hpp.

std::vector<EntityHandle> moab::DagMC::facetList [private]

Definition at line 555 of file DagMC.hpp.

std::vector<EntityHandle> moab::DagMC::facList [private]

Definition at line 559 of file DagMC.hpp.

std::vector<RefEntity *> moab::DagMC::geomEntities [private]

Definition at line 536 of file DagMC.hpp.

Definition at line 522 of file DagMC.hpp.

true if user requested CAD-based ray firing

Definition at line 550 of file DagMC.hpp.

Definition at line 522 of file DagMC.hpp.

Definition at line 521 of file DagMC.hpp.

Definition at line 544 of file DagMC.hpp.

DagMC * moab::DagMC::instance_ = NULL [static, private]

Definition at line 517 of file DagMC.hpp.

Definition at line 518 of file DagMC.hpp.

long long int moab::DagMC::n_pt_in_vol_calls [private]

Definition at line 562 of file DagMC.hpp.

long long int moab::DagMC::n_ray_fire_calls [private]

Definition at line 562 of file DagMC.hpp.

Definition at line 522 of file DagMC.hpp.

const std::map< std::string, std::string > moab::DagMC::no_synonyms [static, private]

Definition at line 540 of file DagMC.hpp.

Definition at line 547 of file DagMC.hpp.

Definition at line 522 of file DagMC.hpp.

Definition at line 520 of file DagMC.hpp.

Definition at line 546 of file DagMC.hpp.

std::vector<EntityHandle> moab::DagMC::prevFacetList [private]

Definition at line 555 of file DagMC.hpp.

std::map<std::string, Tag> moab::DagMC::property_tagmap [private]

Definition at line 542 of file DagMC.hpp.

std::vector<EntityHandle> moab::DagMC::rootSets [private]

Definition at line 531 of file DagMC.hpp.

Definition at line 522 of file DagMC.hpp.

Definition at line 528 of file DagMC.hpp.

std::vector<EntityHandle> moab::DagMC::surfList [private]

Definition at line 555 of file DagMC.hpp.

std::vector<EntityHandle> moab::DagMC::surList [private]

Definition at line 559 of file DagMC.hpp.

bool moab::DagMC::useCAD [private]

Definition at line 549 of file DagMC.hpp.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines