Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
GeomTopoTool.hpp
Go to the documentation of this file.
00001 /**
00002  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
00003  * storing and accessing finite element mesh data.
00004  *
00005  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
00006  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
00007  * retains certain rights in this software.
00008  *
00009  * This library is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU Lesser General Public
00011  * License as published by the Free Software Foundation; either
00012  * version 2.1 of the License, or (at your option) any later version.
00013  *
00014  */
00015 
00016 #ifndef MOAB_GEOM_TOPO_TOOL_HPP
00017 #define MOAB_GEOM_TOPO_TOOL_HPP
00018 
00019 #include "moab/Forward.hpp"
00020 #include "moab/Range.hpp"
00021 
00022 #include <map>
00023 #include <cassert>
00024 
00025 namespace moab
00026 {
00027 
00028 // forward declare this class to avoid the header leaking in here
00029 class OrientedBoxTreeTool;
00030 class GeomQueryTool;
00031 
00032 /** \class GeomTopoTool
00033  * \brief Tool for interpreting geometric topology sets in MOAB database
00034  * Tool for interpreting geometric topology sets in MOAB database; see MOAB metadata_info
00035  * document for information on how geometric topology sets are read and represented.
00036  */
00037 class GeomTopoTool
00038 {
00039   public:
00040     /** \brief Constructor (creates a GTT object)                       \
00041      *  Construct a GeomTopoTool object and search for geometric EntitySets if they
00042      *  exist in the provided moab instance.
00043      *  \param impl MOAB instance the GeomTopoTool will operate on.
00044      *  \param find_geoments if specified as True, geometric objects in the provided MOAB instance
00045                              will be searched for and added to the GTT.
00046         \param modelRootSet the GTT will operate only on geometric EntitySets contained by this
00047      EntitySet. If unprovided, the default value for the modelRootSet is the MOAB instance's root
00048      set, which contains everything in the instance. \param p_rootSets_vector determines the storage
00049      datastructure used to relate geometric EntitySets to their OrientedBoundingBox (OBB) Tree
00050      roots. If set to true (default) a vector will be used to store the root sets along with an
00051      EntityHandle offset for fast lookup of the root sets. If set to false, then a map will be used
00052      to link geometric EntitySets (keys) to the OBB Tree root sets (values). \param restore_rootSets
00053      determines whether or not to restore the internal index that links geomSets to their
00054      corresponding OBB Root.  Only relevant if find_geoments is true. (default = true)
00055      */
00056     GeomTopoTool( Interface* impl,
00057                   bool find_geoments        = false,
00058                   EntityHandle modelRootSet = 0,
00059                   bool p_rootSets_vector    = true,
00060                   bool restore_rootSets     = true );
00061 
00062     ~GeomTopoTool();
00063 
00064     //! Restore parent/child links between GEOM_TOPO mesh sets
00065     ErrorCode restore_topology_from_adjacency();
00066     //! Store sense of entity relative to wrt_entity.
00067     //!\return MB_MULTIPLE_ENTITIES_FOUND if surface already has a forward volume.
00068     //!        MB_SUCCESS if successful
00069     //!        otherwise whatever internal error code occured.
00070     ErrorCode set_sense( EntityHandle entity, EntityHandle wrt_entity, int sense );
00071     //! Get the sense of entity with respect to wrt_entity
00072     //! Returns MB_ENTITY_NOT_FOUND if no relationship found
00073     ErrorCode get_sense( EntityHandle entity, EntityHandle wrt_entity, int& sense );
00074     //! Get the sense of the surface(s) with respect to the volume
00075     ErrorCode get_surface_senses( EntityHandle volume, int num_surfs, const EntityHandle* surfs, int* senses_out );
00076     //! Get the senses of a surface with respect to its volumes
00077     ErrorCode get_surface_senses( EntityHandle surface_ent, EntityHandle& forward_vol, EntityHandle& reverse_vol );
00078 
00079     //! Set the senses of a surface with respect to its volumes
00080     ErrorCode set_surface_senses( EntityHandle surface_ent, EntityHandle forward_vol, EntityHandle reverse_vol );
00081     //! Get the senses of the lower dimension entity handle wrt the higher dimension entities
00082     ErrorCode get_senses( EntityHandle entity, std::vector< EntityHandle >& wrt_entities, std::vector< int >& senses );
00083     //! Set the senses of the entity wrt multiple higher dimension entities
00084     ErrorCode set_senses( EntityHandle entity, std::vector< EntityHandle >& wrt_entities, std::vector< int >& senses );
00085 
00086     /** \brief Get the volume on the other side of a surface
00087      *
00088      * @param A surface to query
00089      * @param old_volume A volume on one side of surface
00090      * @param new_volume Output parameter for volume on the other side of surface
00091      * @return MB_SUCCESS if new_volume was set successfully, error if not.
00092      */
00093     ErrorCode next_vol( EntityHandle surface, EntityHandle old_volume, EntityHandle& new_volume );
00094 
00095     //! Retrieve geometry sets of desired dimension from model set
00096     //  0 = verts, 1 = curves, 2 = surfs, 3 = vols
00097     ErrorCode get_gsets_by_dimension( int dim, Range& gset );
00098 
00099     /** \brief Build obb tree for the entity set given; entity can be surface or volume
00100      *
00101      * @param eh EntityHandle of the volume or surface to construct the OBB tree around
00102      */
00103     ErrorCode construct_obb_tree( EntityHandle eh );
00104 
00105     /** \brief Get the bouding points from a bounding box
00106      *
00107      * @param volume The volume for which the bounding coordinates are requested
00108      * @param minPt Location of the min xyz corner of the volume's axis-aligned bounding box
00109      * @param maxPt Location of the max xyz corner of the volume's axis-aligned bounding box
00110      */
00111     ErrorCode get_bounding_coords( EntityHandle volume, double minPt[3], double maxPt[3] );
00112 
00113     /** \brief Get the center point and three vectors for the OBB of a given volume
00114      *
00115      * @param volume The volume for which the OBB axes will be returned
00116      * @param center coordinates of the oriented bounding box's center point
00117      * @param axis1 scaled axis one of the oriented bounding box
00118      * @param axis2 scaled axis two of the oriented bounding box
00119      * @param axis3 scaled axis three of the oriented bounding box
00120      */
00121     ErrorCode get_obb( EntityHandle volume, double center[3], double axis1[3], double axis2[3], double axis3[3] );
00122 
00123     /** \brief Get the other (d-1)-dimensional entity bounding a set across a (d-2)-dimensional
00124      * entity
00125      *
00126      * Given a d-dimensional entity and one (d-1)-dimensional entity, return the (d-1) dimensional
00127      * entity across a specified (d-2)-dimensional entity.  For example, given a surface, edge, and
00128      * vertex, returns the other edge bounding the surface sharing the vertex.  In the case of
00129      * degenerate results, e.g. two loops bounding a surface and sharing a vertex, tries to step in
00130      * positively-oriented direction.  This won't always work; in those cases, will return
00131      * MB_MULTIPLE_ENTITIES_FOUND.
00132      *
00133      * In the special case where bounded is a curve, then not_this can be a vertex and across zero.
00134      * This function returns the other vertex on the curve.
00135      */
00136     ErrorCode other_entity( EntityHandle bounded, EntityHandle not_this, EntityHandle across, EntityHandle& other );
00137 
00138     /** \brief Return the dimension of the set, or -1 if it's not a geom_dimension set
00139      */
00140     int dimension( EntityHandle this_set );
00141 
00142     /** \brief Return the global ID of a given entity set
00143      *
00144      * @param this_set EntitySet for which the global ID will be returned
00145      */
00146     int global_id( EntityHandle this_set );
00147 
00148     //! Map from dimension & global ID to EntityHandle
00149     EntityHandle entity_by_id( int dimension, int id );
00150 
00151     ErrorCode find_geomsets( Range* ranges = NULL );
00152 
00153     //! Restore the internal cross-referencing of geometry sets and OBB roots
00154     //  The EntityHandle of an OBB Root can be tagged onto the geoemtry EntitySet
00155     //  that it represents so that this relationship can be recovered across
00156     //  write to/read from file.  Since finding the OBB Root for a given geomset
00157     //  is frequent, a faster lookup capability is enabled through data structures
00158     //  in GeomTopoTool (i.e. rootSets or mapRootSets).  This data structure
00159     //  needs to be populated upon file read.
00160     ErrorCode restore_obb_index();
00161 
00162     //! Build obb trees for all surfaces and volumes in model set.
00163     //  If make_one_vol true, joins trees from all surfaces in model into single
00164     //  volume obb tree.
00165     ErrorCode construct_obb_trees( bool make_one_vol = false );
00166 
00167     //! Delete the OBB tree of a volume or surface.
00168     //  If the passed entity is a volume, and the bool 'vol_only'
00169     //  is True, function will delete the volume OBB tree, but
00170     //  OBB trees of the surfaces that compose (are children of)
00171     //  the volume will remain in tact.  If the entity is a volume and
00172     //  'vol_only' is False, function will delete the volume OBB tree
00173     //  along with all child surface OBB trees.
00174     ErrorCode delete_obb_tree( EntityHandle gset, bool vol_only = false );
00175 
00176     ErrorCode delete_all_obb_trees();
00177 
00178     //! Delete the root of the obb tree from the set of all roots
00179     ErrorCode remove_root( EntityHandle vol_or_surf );
00180 
00181     //! Get the root of the obbtree for a given entity
00182     ErrorCode get_root( EntityHandle vol_or_surf, EntityHandle& root );
00183 
00184     //! If constructing one volume obb tree by joining all surface trees,
00185     //  get the root of that tree
00186     EntityHandle get_one_vol_root();
00187 
00188     //! Pointer to Oriented Box Tree Tool class
00189     OrientedBoxTreeTool* obb_tree()
00190     {
00191         return obbTree;
00192     }
00193 
00194     //! Adds a geometry set to the range of all geometry sets, the model set, and root set
00195     //  Make sure the set has the proper geometry dimension tag
00196     //  This could make the obb tree out of date
00197     ErrorCode add_geo_set( EntityHandle set, int dimension, int global_id = 0 );
00198 
00199     //! Will assume no geo sets are defined for this surface
00200     //  Will output a mesh_set that contains everything (all sets of interest), for proper output
00201     ErrorCode geometrize_surface_set( EntityHandle surface, EntityHandle& output );
00202 
00203     //! Checks to see if the entity is part of the model set
00204     ErrorCode is_owned_set( EntityHandle eh );
00205 
00206     //! This would be a deep copy, into a new geom topo tool
00207     //  sets will be duplicated, but entities not
00208     //  modelSet will be a new one;
00209     //  will take as input a pointer to a std::vector of gents (surfaces and volumes, usually),
00210     //  which will serve to filter the gents from modelSet (only dependents will be part of the new
00211     //  gtt) if the pointer is null, all gsets in the original modelSet are duplicated
00212     ErrorCode duplicate_model( GeomTopoTool*& duplicate, std::vector< EntityHandle >* pvGEnts = NULL );
00213 
00214     //! Return the model set handle (this is the full geometry)
00215     EntityHandle get_root_model_set()
00216     {
00217         return modelSet;
00218     }
00219 
00220     //! Checks that all geometric entities were created properly
00221     bool check_model();
00222 
00223     //! Should be used instead of keeping multiple ranges, for example in FBEngine
00224     const Range* geoRanges()
00225     {
00226         return geomRanges;
00227     }
00228 
00229     //! Return pointer to moab instance
00230     Interface* get_moab_instance()
00231     {
00232         return mdbImpl;
00233     }
00234 
00235     //! Returns the sense tag (sense2Tag) from check_face_sense_tag
00236     Tag get_sense_tag();
00237 
00238     //! Returns the global ID tag (gidTag) from check_gid_tag
00239     Tag get_gid_tag();
00240 
00241     //! Returns the geometry dimension tag (geomTag) from check_geom_tag
00242     Tag get_geom_tag();
00243 
00244     //! Returns true if obb trees have been added to the rootset
00245     bool have_obb_tree();
00246 
00247     //! returns the number of entities in the modelSet with specified geometric dimension
00248     int num_ents_of_dim( int dim );
00249 
00250     //! sets the implicit complement handle for this tool
00251     ErrorCode setup_implicit_complement();
00252 
00253     //! Get (or optionally, create) the implicit complement handle
00254     ErrorCode get_implicit_complement( EntityHandle& implicit_complement );
00255 
00256     //! detection method for the implicit complement
00257     bool is_implicit_complement( EntityHandle volume );
00258 
00259     /** \brief Discover and store the topological relationships among a set of volumes
00260      * This method may be used to discover the hierarchy that exists in a range of
00261      * volumes, that have no previous sense of hierarchy, and store it according
00262      * to the conventions of GeomTopoTool.
00263      * The following requirements about the range of flat_volumes must be met:
00264      * 1. Each volume must be represented by a single, closed surface
00265      *    a. The surface meshsets have triangles and vertices as members.
00266      *    b. For each "flat volume", there must be two meshsets: one for the
00267      *       volume and another for the surface that encloses it. These must be
00268      *     linked by a parent-child relationship.
00269      *  c. The SENSE_FORWARD tag on the surface meshset must be set to be
00270      *     the volume meshset it encloses.
00271      * 2. The surfaces must not touch or overlap
00272      *
00273      * After the hierarchy is established, the topological relationships between
00274      * surfaces and the volumes that enclose them are set.  This involves:
00275      * 1. Setting parent-child relationship between surfaces and the volumes that
00276      *    enclose them.
00277      * 2. Setting the SENSE_REVERSE tag on the surfaces to be the volume that
00278      *    encloses them.
00279      *
00280      */
00281     ErrorCode restore_topology_from_geometric_inclusion( const Range& flat_volumes );
00282 
00283   private:
00284     Interface* mdbImpl;
00285     Tag sense2Tag;
00286     Tag senseNEntsTag, senseNSensesTag;
00287     Tag geomTag;
00288     Tag gidTag;
00289     Tag nameTag;
00290     Tag obbRootTag;
00291     Tag obbGsetTag;
00292     // the model set encompasses a full topological model
00293     EntityHandle modelSet;
00294     // implicit complement handle cache
00295     EntityHandle impl_compl_handle;
00296 
00297     Range geomRanges[5];  // add one more dimension, for set of gentities; by default, they will
00298                           // have geom_dimension 4
00299     int maxGlobalId[5];   // one max global id for each dimension
00300     bool updated;
00301 
00302     OrientedBoxTreeTool* obbTree;
00303     EntityHandle setOffset;
00304     std::vector< EntityHandle > rootSets;
00305 
00306     bool m_rootSets_vector;
00307     std::map< EntityHandle, EntityHandle > mapRootSets;
00308     EntityHandle oneVolRootSet;
00309 
00310     //! Creates a volume for undefined space in the model
00311     // The implicit complement is composed of all surfaces that only
00312     // have one parent volume, i.e. surfaces that are in contact with the outside
00313     // world
00314     ErrorCode generate_implicit_complement( EntityHandle& implicit_complement_set );
00315 
00316     //! Compute vertices inclusive and put on tag on sets in geom_sets
00317     ErrorCode construct_vertex_ranges( const Range& geom_sets, const Tag verts_tag );
00318 
00319     //! Given a range of geom topology sets, separate by dimension
00320     ErrorCode separate_by_dimension( const Range& geom_sets );
00321 
00322     //! Verify global id tag
00323     ErrorCode check_gid_tag( bool create = false );
00324 
00325     //! Verify geometry tag
00326     ErrorCode check_geom_tag( bool create = false );
00327 
00328     //! Verify sense face tag
00329     ErrorCode check_face_sense_tag( bool create = false );
00330 
00331     //! Verify sense edge tags
00332     ErrorCode check_edge_sense_tags( bool create = false );
00333 
00334     ErrorCode resize_rootSets();
00335 
00336     ErrorCode set_root_set( EntityHandle vol_or_surf, EntityHandle root );
00337 
00338     //! Return a range of children of a desired geometric dimension
00339     Range get_ct_children_by_dimension( const EntityHandle parent, const int desired_dimension );
00340 
00341     //! Test if volume A is enclosed by volume B
00342     //  This will only produce the correct result if the conventions about
00343     //  volumes listed in the restore_topology_from_geometric_inclusion are
00344     //  upheld
00345     bool A_is_in_B( const EntityHandle volume_A, const EntityHandle volume_B, GeomQueryTool* GQT );
00346 
00347     //! Used by restore_topology_from_geometric_inclusion to generate the
00348     //  hierarchical tree of volumes
00349     ErrorCode insert_in_tree( const EntityHandle ct_root, const EntityHandle volume, GeomQueryTool* GQT );
00350 };
00351 
00352 inline int GeomTopoTool::num_ents_of_dim( int dim )
00353 {
00354     assert( 0 <= dim && 3 >= dim );
00355     return geomRanges[dim].size();
00356 }
00357 
00358 // get the root of the obbtree for a given entity
00359 inline ErrorCode GeomTopoTool::get_root( EntityHandle vol_or_surf, EntityHandle& root )
00360 {
00361     if( m_rootSets_vector )
00362     {
00363         unsigned int index = vol_or_surf - setOffset;
00364         root               = ( index < rootSets.size() ? rootSets[index] : 0 );
00365     }
00366     else
00367         root = mapRootSets[vol_or_surf];
00368     return ( root ? MB_SUCCESS : MB_INDEX_OUT_OF_RANGE );
00369 }
00370 
00371 inline EntityHandle GeomTopoTool::get_one_vol_root()
00372 {
00373     return oneVolRootSet;
00374 }
00375 
00376 inline Tag GeomTopoTool::get_sense_tag()
00377 {
00378     check_face_sense_tag( true );
00379     return sense2Tag;
00380 }
00381 
00382 inline Tag GeomTopoTool::get_gid_tag()
00383 {
00384     check_gid_tag( true );
00385     return gidTag;
00386 }
00387 
00388 inline Tag GeomTopoTool::get_geom_tag()
00389 {
00390     check_geom_tag( true );
00391     return geomTag;
00392 }
00393 
00394 inline bool GeomTopoTool::is_implicit_complement( EntityHandle volume )
00395 {
00396     return volume == impl_compl_handle;
00397 }
00398 
00399 }  // namespace moab
00400 
00401 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines