[MOAB-dev] r1351 - MOAB/trunk
tautges at mcs.anl.gov
tautges at mcs.anl.gov
Fri Nov 2 09:51:59 CDT 2007
Author: tautges
Date: 2007-11-02 09:51:58 -0500 (Fri, 02 Nov 2007)
New Revision: 1351
Modified:
MOAB/trunk/GeomTopoTool.cpp
MOAB/trunk/MBCore.cpp
MOAB/trunk/MBCore.hpp
MOAB/trunk/MBInterface.hpp
MOAB/trunk/MBMeshSet.hpp
MOAB/trunk/Tqdcfr.cpp
MOAB/trunk/Tqdcfr.hpp
Log:
Optimizing reading of .cub files and restoring geometric topology from
cub files. Can now read a .cub file with ~5k geometric volumes in
under 10 mins on my laptop.
MBCore, MBInterface, MBMeshSet: implemented contains_entities, which
returns whether a vector of entities is contained in a set.
Tqdcfr:
- use maps of uid to set and gid to set, eliminates several expensive
calls to get_entities_by_type_and_tag
- fixed bug restoring groups which was breaking out after first group
- don't create geom topo sets for bodies or other entities with no
mesh
GeomTopoTool: revised algorithm for restoring geometric topology based
on owned entities and their adjacencies to other entities; new
algorithm is based on individual entities and their ownership in other
geometric topology sets; this algorithm won't scale as well for many
geom entities a few mesh entities as it does for few geom entities and
many mesh entities, but that's how it should be
Passes make check, except for dagmc test which was failing before.
Modified: MOAB/trunk/GeomTopoTool.cpp
===================================================================
--- MOAB/trunk/GeomTopoTool.cpp 2007-10-31 18:07:43 UTC (rev 1350)
+++ MOAB/trunk/GeomTopoTool.cpp 2007-11-02 14:51:58 UTC (rev 1351)
@@ -25,16 +25,15 @@
// look for geometric topology sets and restore parent/child links between them
// algorithm:
- // - for each entity
- // . get nodes inclusive in range and store as tag on entity
- // - for each entity of dimension, starting low & working upward:
- // . for each vertex: if part of boundary:
- // - get all connected entities of dim d-1
- // - look for intersections in ranges of (d-1)-entities with d-entity (whole intersection)
- // - if whole intersection, make parent/child link between (d-1)- and d-entities
+ // - for each entity of dimension d=D-1..0:
+ // . get d-dimensional entity in entity
+ // . get all (d+1)-dim adjs to that entity
+ // . for each geom entity if dim d+1, if it contains any of the ents,
+ // add it to list of parents
+ // . make parent/child links with parents
// get the geom topology tag
- MBTag geom_tag, verts_tag;
+ MBTag geom_tag;
MBErrorCode result = mdbImpl->tag_create(GEOM_DIMENSION_TAG_NAME, 4,
MB_TAG_SPARSE, geom_tag, NULL);
if (MB_SUCCESS != result && MB_ALREADY_ALLOCATED != result)
@@ -47,98 +46,45 @@
if (MB_SUCCESS != result || geom_sets.empty())
return result;
- // create the tag holding the vertex ranges
- MBRange *dum = NULL;
- result = mdbImpl->tag_create("__vertex_range", sizeof(MBRange*), MB_TAG_SPARSE, verts_tag, &dum);
- if (MB_SUCCESS != result)
- return result;
-
- result = construct_vertex_ranges(geom_sets, verts_tag);
- if (MB_SUCCESS != result)
- return result;
-
MBRange entities[4];
result = separate_by_dimension(geom_sets, entities, geom_tag);
if (MB_SUCCESS != result)
return result;
+ std::vector<MBEntityHandle> dp1ents;
+
// loop over dimensions
- std::vector<MBRange*> dm1_ranges, d_ranges, d0_ranges;
- MBRange::iterator d_it, v_it;
- std::vector<MBRange*>::iterator d_rit, v_rit;
- MBRange inters;
- std::vector<MBEntityHandle> dm1_ents;
- MBRange real_dm1_ents;
+ for (int dim = 2; dim >= 0; dim--) {
+ for (MBRange::iterator d_it = entities[dim].begin();
+ d_it != entities[dim].end(); d_it++) {
+ MBRange dents;
+ result = mdbImpl->get_entities_by_dimension(*d_it, dim, dents);
+ if (MB_SUCCESS != result) continue;
+ if (dents.empty()) continue;
+
+ // get (d+1)-dimensional adjs
+ dp1ents.clear();
+ result = mdbImpl->get_adjacencies(&(*dents.begin()), 1, dim+1,
+ false, dp1ents);
+ if (MB_SUCCESS != result || dp1ents.empty()) continue;
- // pre-load d_ranges with vertex ranges
- d_ranges.resize(entities[0].size());
- result = mdbImpl->tag_get_data(verts_tag, entities[0], &d_ranges[0]);
- if (MB_SUCCESS != result)
- return result;
- d0_ranges = d_ranges;
-
- for (int dim = 1; dim <= 3; dim++) {
+ // . for each geom entity if dim d+1, if it contains any of the ents,
+ // add it to list of parents
- // get the range * for these entities
- dm1_ranges.swap(d_ranges);
- d_ranges.resize(entities[dim].size());
- result = mdbImpl->tag_get_data(verts_tag, entities[dim], &d_ranges[0]);
- if (MB_SUCCESS != result)
- return result;
-
- for (d_it = entities[dim].begin(), d_rit = d_ranges.begin();
- d_it != entities[dim].end(); d_it++, d_rit++) {
-
- // iterate over vertices, finding any with intersected ranges
- //dm1_ents.clear();
- real_dm1_ents.clear();
- for (v_it = entities[dim-1].begin(), v_rit = dm1_ranges.begin();
- v_it != entities[dim-1].end(); v_it++, v_rit++) {
- inters = (*v_rit)->intersect(*(*d_rit));
- if (!inters.empty() && inters.size() == (*v_rit)->size()) {
- // non-zero intersection; get possible parent sets
-/*
- if (dim == 1) dm1_ents.push_back(*v_it);
- else {
- result = mdbImpl->get_parent_meshsets(*v_it, dm1_ents);
- if (MB_SUCCESS != result)
- return result;
- }
-*/
- real_dm1_ents.insert(*v_it);
- }
+ MBRange parents;
+ for (MBRange::iterator git = entities[dim+1].begin();
+ git != entities[dim+1].end(); git++) {
+ if (mdbImpl->contains_entities(*git, &dp1ents[0],
+ dp1ents.size()))
+ parents.insert(*git);
}
-
-/*
- // ok, we have possible children; check for real overlap, but only
- // if we're not doing edges
- if (dim == 1) std::copy(dm1_ents.begin(), dm1_ents.end(),
- mb_range_inserter(real_dm1_ents));
- else {
- // reuse v_it, v_rit here
- for (v_it = entities[dim-1].begin(), v_rit = dm1_ranges.begin();
- v_it != entities[dim-1].end(); v_it++, v_rit++) {
-
- // if this isn't one of the possible dm1 entities, go on
- if (std::find(dm1_ents.begin(), dm1_ents.end(), *v_it) == dm1_ents.end())
- continue;
-
- // check for real overlap
- inters = (*v_rit)->intersect(*(*d_rit));
- if (!inters.empty()) real_dm1_ents.insert(*v_it);
- }
- }
-
-*/
- // ok, we have the real children; add parent/child links; reuse v_it again
- for (v_it = real_dm1_ents.begin(); v_it != real_dm1_ents.end(); v_it++) {
- result = mdbImpl->add_parent_child(*d_it, *v_it);
- if (MB_SUCCESS != result)
- return result;
- }
-
- } // d_it
+ // . make parent/child links with parents
+ for (MBRange::iterator pit = parents.begin(); pit != parents.end(); pit++) {
+ result = mdbImpl->add_parent_child(*pit, *d_it);
+ }
+ }
+
} // dim
return result;
Modified: MOAB/trunk/MBCore.cpp
===================================================================
--- MOAB/trunk/MBCore.cpp 2007-10-31 18:07:43 UTC (rev 1350)
+++ MOAB/trunk/MBCore.cpp 2007-11-02 14:51:58 UTC (rev 1351)
@@ -2318,6 +2318,19 @@
return MB_ENTITY_NOT_FOUND;
}
+ //! return true if all entities are contained in set
+bool MBCore::contains_entities(MBEntityHandle meshset,
+ MBEntityHandle *entities,
+ int num_entities,
+ const int operation_type)
+{
+ MBMeshSet* set = get_mesh_set( sequence_manager(), meshset );
+ if (set)
+ return set->contains_entities(entities, num_entities, operation_type);
+ else
+ return false;
+}
+
// replace entities in a meshset; entities appear in pairs,
// old then new entity in each pair
bool MBCore::replace_entities(MBEntityHandle meshset,
Modified: MOAB/trunk/MBCore.hpp
===================================================================
--- MOAB/trunk/MBCore.hpp 2007-10-31 18:07:43 UTC (rev 1350)
+++ MOAB/trunk/MBCore.hpp 2007-11-02 14:51:58 UTC (rev 1351)
@@ -767,6 +767,13 @@
const MBEntityHandle *entities,
const int num_entities);
+ //! return true if all entities are contained in set
+ virtual bool contains_entities(MBEntityHandle meshset,
+ MBEntityHandle *entities,
+ int num_entities,
+ const int operation_type = MBInterface::INTERSECT);
+
+ //! replace entities
virtual bool replace_entities(MBEntityHandle meshset,
MBEntityHandle *entities,
int num_entities);
Modified: MOAB/trunk/MBInterface.hpp
===================================================================
--- MOAB/trunk/MBInterface.hpp 2007-10-31 18:07:43 UTC (rev 1350)
+++ MOAB/trunk/MBInterface.hpp 2007-11-02 14:51:58 UTC (rev 1351)
@@ -1272,6 +1272,19 @@
const MBEntityHandle *entities,
const int num_entities) = 0;
+ //! Return whether a set contains entities
+ /** Return whether a set contains entities. Returns true only
+ * if ALL entities are contained
+ * \param meshset Mesh set being queried
+ * \param entities Entities being queried
+ * \param num_entities Number of entities
+ * \return bool If true, all entities are contained in set
+ */
+ virtual bool contains_entities(MBEntityHandle meshset,
+ MBEntityHandle *entities,
+ int num_entities,
+ const int operation_type = MBInterface::INTERSECT) = 0;
+
//! Replace entities in a set with other entities
/** Replace entities in a set with other entities; entity list
* specified in pairs of (old, new)
Modified: MOAB/trunk/MBMeshSet.hpp
===================================================================
--- MOAB/trunk/MBMeshSet.hpp 2007-10-31 18:07:43 UTC (rev 1350)
+++ MOAB/trunk/MBMeshSet.hpp 2007-11-02 14:51:58 UTC (rev 1351)
@@ -34,6 +34,7 @@
#include "MBTypes.h"
#include "MBRange.hpp"
+#include "MBInterface.hpp"
#include "MBCN.hpp"
#include <vector>
@@ -165,7 +166,11 @@
MBRange& entity_list) const;
inline MBErrorCode get_non_set_entities( MBRange& range ) const;
-
+
+ inline bool contains_entities(MBEntityHandle *entities,
+ int num_entities,
+ const int operation_type) const;
+
//! subtract/intersect/unite meshset_2 from/with/into meshset_1; modifies meshset_1
inline MBErrorCode subtract(const MBMeshSet *meshset_2,
MBEntityHandle my_handle,
@@ -279,6 +284,10 @@
inline MBErrorCode get_entities_by_dimension(int dimension, \
MBRange& entity_list) const; \
\
+ inline bool contains_entities(MBEntityHandle *entities, \
+ int num_entities, \
+ const int operation_type) const; \
+ \
inline MBErrorCode get_non_set_entities( MBRange& range ) const; \
\
MBErrorCode subtract(const MBMeshSet*, \
@@ -460,6 +469,15 @@
reinterpret_cast<const MBMeshSet_MBRange*>(this)->get_non_set_entities( entities ) ;
}
+inline bool
+MBMeshSet::contains_entities(MBEntityHandle *entities, int num_entities,
+ const int operation_type) const
+{
+ return vector_based() ?
+ reinterpret_cast<const MBMeshSet_Vector* >(this)->contains_entities(entities, num_entities, operation_type) :
+ reinterpret_cast<const MBMeshSet_MBRange*>(this)->contains_entities(entities, num_entities, operation_type) ;
+}
+
//! subtract/intersect/unite meshset_2 from/with/into meshset_1; modifies meshset_1
inline MBErrorCode
MBMeshSet::subtract( const MBMeshSet *meshset_2,
@@ -626,6 +644,20 @@
return MB_SUCCESS;
}
+inline bool MBMeshSet_MBRange::contains_entities(MBEntityHandle *entities,
+ int num_entities,
+ const int op_type) const
+{
+ unsigned int num_in = 0;
+ for (int i = 0; i < num_entities; i++) {
+ if (mRange.find(entities[i]) != mRange.end()) num_in++;
+ else if (MBInterface::INTERSECT == op_type) break;
+ }
+
+ return ((MBInterface::INTERSECT == op_type && num_in == (unsigned int) num_entities) ||
+ MBInterface::UNION == op_type && 0 < num_in);
+}
+
inline unsigned int MBMeshSet_MBRange::num_entities() const
{
return mRange.size();
@@ -702,7 +734,21 @@
return MB_SUCCESS;
}
+inline bool MBMeshSet_Vector::contains_entities(MBEntityHandle *entities,
+ int num_entities,
+ const int op_type) const
+{
+ unsigned int num_in = 0;
+ for (int i = 0; i < num_entities; i++) {
+ if (std::find(mVector.begin(), mVector.end(), entities[i]) !=
+ mVector.end()) num_in++;
+ else if (MBInterface::INTERSECT == op_type) break;
+ }
+ return ((MBInterface::INTERSECT == op_type &&
+ num_in == (unsigned int) num_entities) ||
+ MBInterface::UNION == op_type && 0 < num_in);
+}
inline unsigned int MBMeshSet_Vector::num_entities() const
{
Modified: MOAB/trunk/Tqdcfr.cpp
===================================================================
--- MOAB/trunk/Tqdcfr.cpp 2007-10-31 18:07:43 UTC (rev 1350)
+++ MOAB/trunk/Tqdcfr.cpp 2007-11-02 14:51:58 UTC (rev 1351)
@@ -831,29 +831,10 @@
int *id_buf, const int id_buf_size,
std::vector<MBEntityHandle> &entities)
{
- MBErrorCode result = MB_SUCCESS;
- MBTag tags[2] = {geomTag, globalIdTag};
- int tag_vals[2];
- const void *tag_vals_ptr[2] = {tag_vals, tag_vals+1};
- if (this_type == GROUP || this_type == BODY) tag_vals[0] = 4;
- else tag_vals[0] = 5 - this_type;
+ for (int i = 0; i < id_buf_size; i++)
+ entities.push_back((gidSetMap[5-this_type])[id_buf[i]]);
- MBRange tmp_ents;
-
- for (int i = 0; i < id_buf_size; i++) {
- tmp_ents.clear();
- // set id tag value
- tag_vals[1] = id_buf[i];
-
- MBErrorCode tmp_result = mdbImpl->get_entities_by_type_and_tag(0, MBENTITYSET, tags,
- tag_vals_ptr, 2, tmp_ents);
- if (MB_SUCCESS != tmp_result && MB_TAG_NOT_FOUND != tmp_result) result = tmp_result;
-
- if (MB_SUCCESS == tmp_result && tmp_ents.size() == 1)
- entities.push_back(*tmp_ents.begin());
- }
-
- return result;
+ return MB_SUCCESS;
}
MBErrorCode Tqdcfr::get_mesh_entities(const int this_type,
@@ -1373,10 +1354,6 @@
for (int i = 0; i < info.numEntities; i++) {
- // create an entity set for this entity
- result = instance->create_set(geom_headers[i].setHandle);
- if (MB_SUCCESS != result) return result;
-
instance->FREADI(8);
geom_headers[i].nodeCt = instance->int_buf[0];
geom_headers[i].nodeOffset = instance->int_buf[1];
@@ -1386,6 +1363,14 @@
geom_headers[i].elemLength = instance->int_buf[5];
geom_headers[i].geomID = instance->int_buf[6];
+ // don't represent in MOAB if no mesh
+ if (geom_headers[i].nodeCt == 0 && geom_headers[i].elemCt == 0)
+ continue;
+
+ // create an entity set for this entity
+ result = instance->create_set(geom_headers[i].setHandle);
+ if (MB_SUCCESS != result) return result;
+
// set the dimension to -1; will have to reset later, after elements are read
dum_int = -1;
result = instance->mdbImpl->tag_set_data(instance->geomTag,
@@ -1397,6 +1382,9 @@
&(geom_headers[i].setHandle), 1,
&(geom_headers[i].geomID));
if (MB_SUCCESS != result) return result;
+
+ // put the set and uid into a map for later
+ instance->uidSetMap[geom_headers[i].geomID] = geom_headers[i].setHandle;
}
return MB_SUCCESS;
@@ -1444,8 +1432,8 @@
&(group_headers[i].setHandle), 1,
&(group_headers[i].grpID));
if (MB_SUCCESS != result) return result;
-
- break;
+
+ instance->gidSetMap[5][group_headers[i].grpID] = group_headers[i].setHandle;
}
return MB_SUCCESS;
@@ -2016,21 +2004,21 @@
// parsed the data; now put on mdb entities; first we need to find the entity
if (records[entity_rec_num].entity == 0) {
- // get the choices of entity
- MBRange entities;
- if (uid == -1) return MB_FAILURE;
- const void *dum_ptr = &uid;
- result = mdbImpl->get_entities_by_type_and_tag(0, MBENTITYSET, &uniqueIdTag,
- &dum_ptr, 1, entities);
- if (MB_SUCCESS != result) return result;
- if (entities.size() != 1) return MB_FAILURE;
- else records[entity_rec_num].entity = *entities.begin();
+ records[entity_rec_num].entity = uidSetMap[uid];
}
// set the id
if (id != -1) {
result = mdbImpl->tag_set_data(globalIdTag, &(records[entity_rec_num].entity), 1, &id);
if (MB_SUCCESS != result) return result;
+
+ int ent_dim = -1;
+ if (records[entity_rec_num].rec_type == aBODY) ent_dim = 4;
+ else if (records[entity_rec_num].rec_type == LUMP) ent_dim = 3;
+ else if (records[entity_rec_num].rec_type == FACE) ent_dim = 2;
+ else if (records[entity_rec_num].rec_type == aEDGE) ent_dim = 1;
+ else if (records[entity_rec_num].rec_type == aVERTEX) ent_dim = 0;
+ if (-1 != ent_dim) gidSetMap[ent_dim][id] = records[entity_rec_num].entity;
}
// set the name
Modified: MOAB/trunk/Tqdcfr.hpp
===================================================================
--- MOAB/trunk/Tqdcfr.hpp 2007-10-31 18:07:43 UTC (rev 1350)
+++ MOAB/trunk/Tqdcfr.hpp 2007-11-02 14:51:58 UTC (rev 1351)
@@ -28,6 +28,7 @@
#include <stdio.h>
#include <string>
#include <vector>
+#include <map>
class MBReadUtilIface;
class FEModelHeader;
@@ -276,6 +277,8 @@
int currElementIdOffset[MBMAXTYPE];
MBTag globalIdTag, cubIdTag, geomTag, uniqueIdTag, blockTag, nsTag, ssTag,
attribVectorTag, entityNameTag, categoryTag;
+ std::map<int, MBEntityHandle> uidSetMap;
+ std::map<int, MBEntityHandle> gidSetMap[6];
std::vector<int> int_buf;
std::vector<double> dbl_buf;
@@ -345,7 +348,7 @@
// map between cub ids and MOAB handles
std::vector<MBEntityHandle> *cubMOABVertexMap;
-
+
// enum used to identify element/entity type in groups
enum {GROUP = 0, BODY, VOLUME, SURFACE, CURVE, VERTEX, HEX, TET, PYRAMID, QUAD, TRI, EDGE, NODE};
static const MBEntityType group_type_to_mb_type[];
More information about the moab-dev
mailing list