moab
SimpleMoab.h
Go to the documentation of this file.
00001 
00002 #ifndef __smoab_SimpleMoab_h
00003 #define __smoab_SimpleMoab_h
00004 
00005 #include "moab/Core.hpp"
00006 #include "moab/Interface.hpp"
00007 #include "moab/Range.hpp"
00008 #include "MBTagConventions.hpp"
00009 
00010 #include <iostream>
00011 
00012 namespace smoab
00013 {
00014 //adjacency intersect / union named enum to match
00015 //the types from moab
00016 enum adjacency_type{INTERSECT=moab::Interface::INTERSECT,
00017                     UNION=moab::Interface::UNION};
00018 
00019 //make our range equal to moabs range
00020 typedef moab::Range Range;
00021 typedef moab::EntityHandle EntityHandle;
00022 typedef moab::EntityType EntityType;
00023 
00024 //bring in range functions
00025 using moab::intersect;
00026 using moab::subtract;
00027 using moab::unite;
00028 
00029 //forward declare this->Moab for Tag
00030 struct Interface;
00031 
00032 //forward declar the DataSetConverter so it can be a friend of Interface
00033 class DataSetConverter;
00034 
00035 class Tag
00036 {
00037   const std::string Name_;
00038 public:
00039   Tag(std::string const& n):Name_(n)
00040     {
00041     }
00042 
00043   virtual ~Tag()
00044     {
00045     }
00046 
00047   const char* name() const { return this->Name_.c_str(); }
00048   moab::DataType virtual dataType() const { return moab::MB_TYPE_INTEGER; }
00049   virtual bool isComparable() const { return false; }
00050   virtual int value() const { return int(); }
00051 };
00052 
00053 //lightweight structs to wrap set names, so we detected
00054 //incorrect names at compile time. In the future I expect material and
00055 //boundary conditions to be comparable
00056 class MaterialTag : public Tag{ public: MaterialTag():Tag("MATERIAL_SET"){}};
00057 class DirichletTag : public Tag{ public: DirichletTag():Tag("DIRICHLET_SET"){}};
00058 class NeumannTag: public Tag{public:  NeumannTag():Tag("NEUMANN_SET"){}};
00059 class GroupTag: public Tag{ public: GroupTag():Tag("GROUP"){}};
00060 
00061 //geom is the only comparable tag, since it can have a dimension.
00062 class GeomTag: public Tag
00063   {
00064   int dim;
00065 public:
00066   GeomTag(int d):Tag("GEOM_DIMENSION"),dim(d){}
00067   GeomTag():Tag("GEOM_DIMENSION"), dim(0){}
00068 
00069   bool isComparable() const { return dim > 0; }
00070   int value() const { return dim; }
00071   };
00072 
00073 //light weight wrapper on a moab this->Moab that exposes only the reduced class
00074 //that we need
00075 class Interface
00076 {
00077 public:
00078   Interface(const std::string &file)
00079     {
00080     this->Moab = new moab::Core();
00081     this->Moab->load_file(file.c_str());
00082     }
00083 
00084   ~Interface()
00085     {
00086     if(this->Moab)
00087       {
00088       delete this->Moab;
00089       this->Moab = NULL;
00090       }
00091     }
00092 
00093   //----------------------------------------------------------------------------
00094   moab::Tag getMoabTag(const smoab::Tag& simpleTag) const
00095     {
00096     moab::Tag tag;
00097     this->Moab->tag_get_handle(simpleTag.name(),
00098                                1,
00099                                simpleTag.dataType(),
00100                                tag);
00101     return tag;
00102     }
00103 
00104   //----------------------------------------------------------------------------
00105   //returns the moab name for the given entity handle if it has a sparse Name tag
00106   std::string name(const smoab::EntityHandle& entity) const
00107     {
00108     moab::Tag nameTag;
00109     moab::ErrorCode rval = this->Moab->tag_get_handle(NAME_TAG_NAME,
00110                                                       NAME_TAG_SIZE,
00111                                                       moab::MB_TYPE_OPAQUE,
00112                                                       nameTag);
00113     if(rval != moab::MB_SUCCESS) { return std::string(); }
00114 
00115     char name[NAME_TAG_SIZE];
00116     rval = this->Moab->tag_get_data(nameTag,&entity,1,&name);
00117     if(rval != moab::MB_SUCCESS) { return std::string(); }
00118 
00119     return std::string(name);
00120     }
00121 
00122 
00123   //----------------------------------------------------------------------------
00124   smoab::EntityHandle getRoot() const { return this->Moab->get_root_set(); }
00125 
00126   //----------------------------------------------------------------------------
00127   smoab::Range findEntities(const smoab::EntityHandle root, moab::EntityType type) const
00128     {
00129     smoab::Range result;
00130     // get all the sets of that type in the mesh
00131     this->Moab->get_entities_by_type(root, type, result);
00132     return result;
00133     }
00134 
00135   //----------------------------------------------------------------------------
00136   //Find all entities with a given tag. We don't use geom as a tag as that
00137   //isn't a fast operation. Yes finding the intersection of geom entities and
00138   //a material / boundary tag will be more work, but it is rarely done currently
00139   //Returns the found group of entities
00140   smoab::Range findEntitiesWithTag (const smoab::Tag& tag, smoab::EntityHandle root,
00141                                     moab::EntityType type = moab::MBENTITYSET) const
00142     {
00143     smoab::Range result;
00144 
00145     moab::Tag t = this->getMoabTag(tag);
00146 
00147     // get all the entities of that type in the mesh
00148     this->Moab->get_entities_by_type_and_tag(root, type, &t, NULL, 1,result);
00149 
00150 
00151     if(tag.isComparable())
00152       {
00153       int value=0;
00154       //now we have to remove any that doesn't match the tag value
00155       smoab::Range resultMatchingTag;
00156       typedef moab::Range::const_iterator iterator;
00157       for(iterator i=result.begin();
00158           i != result.end();
00159           ++i)
00160         {
00161         value = 0;
00162         moab::EntityHandle handle = *i;
00163         this->Moab->tag_get_data(t, &handle, 1, &value);
00164         if(value == tag.value())
00165           {
00166           resultMatchingTag.insert(*i);
00167           }
00168         }
00169 
00170       return resultMatchingTag;
00171       }
00172     else
00173       {
00174       //we return all the items we found
00175       return result;
00176       }
00177     }
00178 
00179   //----------------------------------------------------------------------------
00180   //Find all entities from a given root of a given dimensionality
00181   smoab::Range findEntitiesWithDimension(const smoab::EntityHandle root,
00182                                          int dimension) const
00183     {
00184     typedef smoab::Range::const_iterator iterator;
00185 
00186     smoab::Range result;
00187     this->Moab->get_entities_by_dimension(root,dimension,result);
00188 
00189 
00190     smoab::Range children;
00191     this->Moab->get_child_meshsets(root,children,0);
00192     for(iterator i=children.begin(); i !=children.end();++i)
00193       {
00194       this->Moab->get_entities_by_dimension(*i,dimension,result);
00195       }
00196     return result;
00197     }
00198 
00199   //----------------------------------------------------------------------------
00200   smoab::Range findAdjacentEntities(const smoab::EntityHandle& entity,
00201                                     int dimension) const
00202     {
00203     const int adjType = static_cast<int>(smoab::INTERSECT);
00204     smoab::Range result;
00205     const bool create_if_missing = false;
00206     this->Moab->get_adjacencies(&entity,
00207                                 1,
00208                                 dimension,
00209                                 create_if_missing,
00210                                 result,
00211                                 adjType);
00212     return result;
00213     }
00214 
00215     //----------------------------------------------------------------------------
00216   smoab::Range findAdjacentEntities(const smoab::Range& range,
00217                                     int dimension,
00218                                     const smoab::adjacency_type type = smoab::UNION) const
00219     {
00220     //the smoab and moab adjacent intersection enums are in the same order
00221     const int adjType = static_cast<int>(type);
00222     smoab::Range result;
00223     const bool create_if_missing = false;
00224     this->Moab->get_adjacencies(range,dimension,
00225                                 create_if_missing,
00226                                 result,
00227                                 adjType);
00228 
00229     return result;
00230     }
00231 
00232   //----------------------------------------------------------------------------
00233   //Find all elements in the database that have children and zero parents.
00234   //this doesn't find
00235   smoab::Range findEntityRootParents(smoab::EntityHandle const& root) const
00236     {
00237     smoab::Range parents;
00238 
00239     typedef moab::Range::const_iterator iterator;
00240     moab::Range sets;
00241 
00242     this->Moab->get_entities_by_type(root, moab::MBENTITYSET, sets);
00243     for(iterator i=sets.begin(); i!=sets.end();++i)
00244       {
00245       int numParents=0,numChildren=0;
00246       this->Moab->num_parent_meshsets(*i,&numParents);
00247       if(numParents==0)
00248         {
00249         this->Moab->num_child_meshsets(*i,&numChildren);
00250         if(numChildren>=0)
00251           {
00252           parents.insert(*i);
00253           }
00254         }
00255       }
00256     return parents;
00257     }
00258 
00259   //----------------------------------------------------------------------------
00260   //finds entities that have zero children and zero parents
00261   smoab::Range findDetachedEntities(moab::EntityHandle const& root) const
00262     {
00263     smoab::Range detached;
00264 
00265     typedef moab::Range::const_iterator iterator;
00266     moab::Range sets;
00267 
00268     this->Moab->get_entities_by_type(root, moab::MBENTITYSET, sets);
00269     for(iterator i=sets.begin(); i!=sets.end();++i)
00270       {
00271       int numParents=0,numChildren=0;
00272       this->Moab->num_parent_meshsets(*i,&numParents);
00273       if(numParents==0)
00274         {
00275         this->Moab->num_child_meshsets(*i,&numChildren);
00276         if(numChildren==0)
00277           {
00278           detached.insert(*i);
00279           }
00280         }
00281       }
00282     return detached;
00283     }
00284 
00285   //----------------------------------------------------------------------------
00286   //find all children of the entity passed in that has multiple parents
00287   smoab::Range findEntitiesWithMultipleParents(smoab::EntityHandle const& root)
00288     {
00289     smoab::Range multipleParents;
00290     typedef moab::Range::const_iterator iterator;
00291 
00292     //for all the elements in the range, find all items with multiple parents
00293     moab::Range children;
00294     this->Moab->get_child_meshsets(root,children,0);
00295     for(iterator i=children.begin(); i!=children.end();++i)
00296       {
00297       int numParents=0;
00298       this->Moab->num_parent_meshsets(*i,&numParents);
00299       if(numParents>1)
00300         {
00301         multipleParents.insert(*i);
00302         }
00303       }
00304     return multipleParents;
00305     }
00306 
00307   //----------------------------------------------------------------------------
00308   //prints all elements in a range objects
00309   void printRange(smoab::Range const& range)
00310     {
00311     typedef Range::const_iterator iterator;
00312     for(iterator i=range.begin(); i!=range.end(); ++i)
00313       {
00314       std::cout << "entity id: " << *i << std::endl;
00315       this->Moab->list_entity(*i);
00316       }
00317     }
00318 
00319   friend class smoab::DataSetConverter;
00320 private:
00321   moab::Interface* Moab;
00322 };
00323 
00324 }
00325 
00326 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines