MeshKit  1.0
ModelEnt.cpp
Go to the documentation of this file.
00001 #include "meshkit/iGeom.hpp"
00002 #include "meshkit/ModelEnt.hpp"
00003 #include "meshkit/SizingFunction.hpp"
00004 #include "meshkit/Error.hpp"
00005 #include "meshkit/MKCore.hpp"
00006 #include "meshkit/Types.hpp"
00007 #include "meshkit/VecUtil.hpp"
00008 #include "moab/GeomTopoTool.hpp"
00009 #include "moab/CN.hpp"
00010 #include "RefEntity.hpp"
00011 
00012 namespace MeshKit
00013 {
00014 
00015 ModelEnt::ModelEnt(MKCore *mk,
00016                    iGeom::EntityHandle geom_ent,
00017                    int geom_index,
00018                    moab::EntityHandle mesh_ent,
00019                    int mesh_index,
00020                    int irel_index,
00021                    int sizing_index)
00022         : mkCore(mk), igeomIndex(geom_index), meshIndex(mesh_index), irelIndex(irel_index),
00023           iGeomEnt(geom_ent), iGeomSet(NULL),
00024           moabEntSet(mesh_ent), sizingFunctionIndex(sizing_index),
00025           meshIntervals(-1), intervalFirmness(DEFAULT), constrainEven(false),
00026           meshedState(NO_MESH), iaVariable(NULL)
00027 {
00028     // mark the geometry entity with this modelEnt
00029   if (-1 != igeomIndex && geom_ent) {
00030     ModelEnt *dum_this = this;
00031     iGeom::Error err = igeom_instance()->setData(geom_ent, mkCore->igeom_model_tag(igeomIndex), &dum_this);
00032     IBERRCHK(err, "Failed to set iGeom ModelEnt tag.");
00033   }
00034   
00035   if (!mesh_ent && geom_ent && -1 != meshIndex && -1 != irelIndex) {
00036       // get a corresponding mesh set handle, if any; create one if missing and have a meshIndex
00037     iBase_EntitySetHandle h;
00038     iRel::Error err = mkCore->irel_pair(irelIndex)->getEntSetRelation(geom_ent, false, h);
00039     if (iBase_SUCCESS == err)
00040       moabEntSet = reinterpret_cast<moab::EntityHandle>(h);
00041     else moabEntSet = 0;
00042   }
00043 
00044     // if still no mesh entity and non-default mesh index, create one
00045   if (-1 != meshIndex && !moabEntSet) {
00046       // should have a geometry entity here
00047     assert(iGeomEnt || iGeomSet);
00048     create_mesh_set();
00049   }
00050   
00051     // if there's a mesh entity, set it to point here too
00052   if (moabEntSet && -1 != meshIndex) {
00053     ModelEnt *dum_this = this;
00054     moab::ErrorCode err = mkCore->moab_instance(meshIndex)->tag_set_data(mkCore->moab_model_tag(meshIndex), 
00055                                                                           &moabEntSet, 1, &dum_this);
00056     IBERRCHK(err, "Failed to set moab ModelEnt tag.");
00057   }
00058 }
00059 
00060 ModelEnt::ModelEnt(MKCore *mk,
00061                    iGeom::EntitySetHandle geom_ent,
00062                    int geom_index,
00063                    moab::EntityHandle mesh_ent,
00064                    int mesh_index,
00065                    int irel_index,
00066                    int sizing_index,
00067                    IAVariable *ia_var)
00068         : mkCore(mk), igeomIndex(geom_index), meshIndex(mesh_index), irelIndex(irel_index),
00069           iGeomEnt(NULL), iGeomSet(geom_ent),
00070           moabEntSet(mesh_ent), sizingFunctionIndex(sizing_index),
00071           meshIntervals(-1), intervalFirmness(DEFAULT), constrainEven(false),
00072           meshedState(NO_MESH), iaVariable(ia_var)
00073 {
00074     // mark the geometry entity with this modelEnt
00075   if (-1 != igeomIndex && geom_ent) {
00076     ModelEnt *dum_this = this;
00077     iGeom::Error err = igeom_instance()->setEntSetData(geom_ent, 
00078     mkCore->igeom_model_tag(igeomIndex),
00079     &dum_this);
00080     IBERRCHK(err, "Failed to set iGeom ModelSet tag.");
00081   }
00082   
00083   if (!mesh_ent && geom_ent && -1 != meshIndex && -1 != irelIndex) {
00084       // get a corresponding mesh set handle, if any; create one if missing and have a meshIndex
00085     iBase_EntitySetHandle h;
00086     iRel::Error err = mkCore->group_set_pair(irelIndex)->getSetSetRelation(geom_ent, false, h);
00087     IBERRCHK(err, "Failed to get mesh set handle for geometry entity set.");
00088     moabEntSet = reinterpret_cast<moab::EntityHandle>(h);
00089   }
00090 
00091     // if still no mesh entity and non-default mesh index, create one
00092   if (-1 != meshIndex && !moabEntSet) {
00093       // should have a geometry entity here
00094     assert(iGeomEnt || iGeomSet);
00095     create_mesh_set();
00096   }
00097   
00098     // if there's a mesh entity, set it to point here too
00099   if (moabEntSet && -1 != meshIndex) {
00100     ModelEnt *dum_this = this;
00101     moab::ErrorCode err = mkCore->moab_instance(meshIndex)->tag_set_data(mkCore->moab_model_tag(meshIndex), 
00102                                                                           &moabEntSet, 1, &dum_this);
00103     IBERRCHK(err, "Failed to set moab ModelEnt tag.");
00104   }
00105 }
00106     
00107 ModelEnt::~ModelEnt() 
00108 {
00109   // if we delete the Model Entity, remove the tag on the gent that
00110   // was pointing to it (or set it to NULL)
00111   /*if (-1 != igeomIndex && iGeomEnt) {
00112     iGeom::Error err = igeom_instance()->rmvTag(iGeomEnt,
00113         mkCore->igeom_model_tag(igeomIndex));
00114     IBERRCHK(err, "Failed to set iGeom ModelEnt tag.");
00115   }
00116 
00117   if (-1 != igeomIndex && iGeomSet) {
00118     iGeom::Error err = igeom_instance()->rmvEntSetTag(iGeomSet,
00119         mkCore->igeom_model_tag(igeomIndex));
00120     IBERRCHK(err, "Failed to set iGeom ModelSet tag.");
00121   }*/
00122   // if we delete the Model Entity, remove the tag on the moabEntSet
00123   // that was pointing to it (set the value to NULL)
00124   if (moabEntSet && -1 != meshIndex) {
00125     moab::ErrorCode err = mkCore->moab_instance(meshIndex)->tag_delete_data(
00126           mkCore->moab_model_tag(meshIndex), &moabEntSet, 1);
00127     IBERRCHK(err, "Failed to set moab ModelEnt tag.");
00128   }
00129 }
00130 
00131 void ModelEnt::sizing_function_index(int ind, bool children_too)
00132 {
00133   sizingFunctionIndex = ind;
00134 
00135     // set on children too, if requested
00136   if (children_too && dimension() > 1) {
00137     MEntVector childrn;
00138     get_adjacencies(dimension()-1, childrn);
00139     for (MEntVector::iterator vit = childrn.begin(); vit != childrn.end(); vit++) 
00140       (*vit)->sizing_function_index(ind, children_too);
00141   }
00142 }
00143     
00144 void ModelEnt::create_mesh_set(int ordered_flag)
00145 {
00146   if (moabEntSet) 
00147     throw Error(MK_FAILURE, "Tried to create meshset for an entity that already had one.");
00148 
00149   else if (!iGeomEnt && !iGeomSet) 
00150     throw Error(MK_FAILURE, "Tried to create ModelEnt missing a geometry entity or set.");
00151   
00152   moab::EntitySetProperty ordered = moab::MESHSET_SET;
00153 
00154     // need type for later
00155   iBase_EntityType this_tp = iBase_ALL_TYPES;
00156   iGeom::Error err;
00157   if (-1 > ordered_flag || 1 < ordered_flag) throw Error(MK_FAILURE, "Invalid value for ordered flag.");
00158   else if (-1 == ordered && iGeomSet) ordered = moab::MESHSET_ORDERED;
00159   else if (-1 == ordered_flag && iGeomEnt) {
00160     err = igeom_instance()->getEntType(iGeomEnt, this_tp);
00161     IBERRCHK(err, "Trouble getting entity type.");
00162     if (iBase_EDGE == this_tp) 
00163       ordered = moab::MESHSET_ORDERED;
00164   }
00165   
00166     // create the set
00167   moab::ErrorCode rval = mkCore->moab_instance(meshIndex)->create_meshset(ordered, moabEntSet);
00168   MBERRCHK(rval, mkCore->moab_instance());
00169 
00170     // set dimension tag and global id tag
00171   if (iBase_ALL_TYPES != this_tp) {
00172     rval = mkCore->moab_instance(meshIndex)->tag_set_data(mkCore->moab_geom_dim_tag(), &moabEntSet, 1, &this_tp);
00173     MBERRCHK(rval, mkCore->moab_instance());
00174     // set the id tag; get it with id(); if it returns -1, set a value 0
00175     int id_from_geometry = id();
00176     if (id_from_geometry==-1)
00177       id_from_geometry = 0;
00178     rval = mkCore->moab_instance(meshIndex)->tag_set_data(mkCore->moab_global_id_tag(), &moabEntSet, 1, &id_from_geometry);
00179     MBERRCHK(rval, mkCore->moab_instance());
00180   }
00181   
00182     // tag it with this ModelEnt
00183   ModelEnt *this_ptr = this;
00184   rval = mkCore->moab_instance()->tag_set_data(mkCore->moab_model_tag(), &moabEntSet, 1, &this_ptr);
00185   MBERRCHK(rval, mkCore->moab_instance());
00186 
00187   // tag it with geometry dimension
00188    int geom_dim = this_tp;
00189    rval = mkCore->moab_instance()->tag_set_data(mkCore->moab_geom_dim_tag(), &moabEntSet, 1, &geom_dim);
00190    MBERRCHK(rval, mkCore->moab_instance());
00191 
00192     // relate the mesh to the geom, only if iRelFlag is true
00193   if (-1 != irelIndex)
00194   {
00195     iRel::Error merr;
00196     if (iGeomEnt) {
00197       merr = mkCore->irel_pair(irelIndex)->setEntSetRelation(iGeomEnt, IBSH(moabEntSet));
00198     }
00199     else {
00200       merr = mkCore->group_set_pair()->setSetSetRelation(iGeomSet, IBSH(moabEntSet));
00201     }
00202     IBERRCHK(merr, "Failed to set iRel relation for a mesh set.");
00203   }
00204 
00205   if (iGeomEnt) init_parents_children();
00206   else if (iGeomSet) init_group_contents();
00207   
00208     // now do senses wrt lower-dimensional entities, which should be initialized by now
00209   if (iGeomEnt)
00210     set_downward_senses();
00211 }
00212 
00213 void ModelEnt::init_parents_children() 
00214 {
00215     // get parents and children, and link to corresponding sets; don't do this for any missing mesh sets, 
00216     // will get done when those sets get created
00217   assert(iGeomEnt);
00218   std::vector<iGeom::EntityHandle> geom_adjs;
00219   std::vector<iGeom::EntityHandle>::iterator vit;
00220   moab::EntityHandle mseth;
00221   //moab::GeomTopoTool gt(mkCore->moab_instance(meshIndex), false);
00222   iBase_EntityType this_tp;
00223   moab::ErrorCode rval;
00224   iGeom::Error err = igeom_instance()->getEntType(iGeomEnt, this_tp);
00225   IBERRCHK(err, "Trouble getting entity type.");
00226   if (this_tp < iBase_REGION) {
00227     err = igeom_instance()->getEntAdj(iGeomEnt, (iBase_EntityType)(this_tp+1), geom_adjs);
00228     IBERRCHK(err, "Trouble finding parent entities.");
00229     for (vit = geom_adjs.begin(); vit != geom_adjs.end(); vit++) {
00230       try {mseth = mesh_handle(*vit);
00231       }
00232       catch (Error err) {
00233           // just continue, will get this when parent gets created later
00234         continue;
00235       }
00236       if (mseth) {
00237         rval = mkCore->moab_instance(meshIndex)->add_parent_child(mseth, moabEntSet);
00238         MBERRCHK(rval, mkCore->moab_instance(meshIndex));
00239       }
00240     }
00241   }
00242   
00243   if (this_tp > iBase_VERTEX) {
00244     geom_adjs.clear();
00245     err = igeom_instance()->getEntAdj(iGeomEnt, (iBase_EntityType)(this_tp-1), geom_adjs);
00246     IBERRCHK(err, "Trouble finding child entities.");
00247     for (vit = geom_adjs.begin(); vit != geom_adjs.end(); vit++) {
00248         // should definitely have a mesh handle here, error if not
00249       mseth = mesh_handle(*vit);
00250       if (mseth) {
00251         rval = mkCore->moab_instance(meshIndex)->add_parent_child(moabEntSet, mseth);
00252         MBERRCHK(rval, mkCore->moab_instance(meshIndex));
00253       }
00254     }
00255   }
00256 }
00257 
00258 void ModelEnt::init_group_contents() 
00259 {
00260     // should only be calling this function if we have a geometry set and a mesh index
00261   assert(-1 != meshIndex && iGeomSet);
00262   
00263     // get the geometry entities in this group and add the corresponding mesh entity sets
00264   std::vector<iGeom::EntityHandle> geom_adjs;
00265   std::vector<iGeom::EntityHandle>::iterator vit;
00266   moab::EntityHandle mseth;
00267   moab::ErrorCode rval;
00268   iGeom::Error err = igeom_instance()->getEntities(iGeomSet, iBase_ALL_TYPES, geom_adjs);
00269   IBERRCHK(err, "Trouble finding parent entities.");
00270   for (vit = geom_adjs.begin(); vit != geom_adjs.end(); vit++) {
00271     mseth = 0;
00272     try {mseth = mesh_handle(*vit);
00273     }
00274     catch (Error) {
00275         // just continue, will get this when parent gets created later
00276       continue;
00277     }
00278     if (mseth) {
00279       rval = mkCore->moab_instance(meshIndex)->add_entities(moabEntSet, &mseth, 1);
00280       MBERRCHK(rval, mkCore->moab_instance(meshIndex));
00281     }
00282   }
00283 
00284     // same for set members
00285   std::vector<iGeom::EntitySetHandle> geom_sadjs;
00286   std::vector<iGeom::EntitySetHandle>::iterator svit;
00287   err = igeom_instance()->getEntSets(iGeomSet, -1, geom_sadjs);
00288   IBERRCHK(err, "Trouble finding parent entities.");
00289   for (svit = geom_sadjs.begin(); svit != geom_sadjs.end(); svit++) {
00290     mseth = 0;
00291     try {mseth = mesh_handle(*svit);
00292     }
00293     catch (Error) {
00294         // just continue, will get this when parent gets created later
00295       continue;
00296     }
00297     if (mseth) {
00298       rval = mkCore->moab_instance(meshIndex)->add_entities(moabEntSet, &mseth, 1);
00299       MBERRCHK(rval, mkCore->moab_instance(meshIndex));
00300     }
00301   }
00302 }
00303   
00304 void ModelEnt::set_downward_senses() 
00305 {
00306     // should have both mesh and geometry
00307   assert(iGeomEnt && moabEntSet && igeomIndex != -1 && meshIndex != -1);
00308 
00309     // skip if vertex or edge
00310   if (dimension() <= 1) return;
00311   
00312   int dim = dimension();
00313   iGeom::Error err;
00314   MEntVector adjs;
00315   get_adjacencies(dim-1, adjs);    
00316 
00317   std::vector<int> senses;
00318   int dum_sense;
00319   std::vector<moab::EntityHandle> ments;
00320   moab::GeomTopoTool gt(mkCore->moab_instance(meshIndex));
00321   for (MEntVector::iterator vit = adjs.begin(); vit != adjs.end(); vit++) {
00322     err = igeom_instance()->getSense((*vit)->geom_handle(), geom_handle(), dum_sense);
00323     IBERRCHK(err, "Problem getting senses.");
00324     moab::ErrorCode rval = gt.set_sense((*vit)->mesh_handle(), mesh_handle(), dum_sense);
00325     MBERRCHK(rval, mkCore->moab_instance(meshIndex));
00326   }
00327 }
00328 
00337 void ModelEnt::commit_mesh(moab::Range &mesh_ents, MeshedState mstate) 
00338 {
00339   std::vector<moab::EntityHandle> ent_vec;
00340   std::copy(mesh_ents.begin(), mesh_ents.end(), std::back_inserter(ent_vec));
00341   commit_mesh(&ent_vec[0], ent_vec.size(), mstate);
00342 }
00343 
00350 void ModelEnt::commit_mesh(moab::EntityHandle *mesh_ents,
00351                            int num_ents,
00352                            MeshedState mstate) 
00353 {
00354   assert(-1 != meshIndex);
00355   
00356     // add them to the set
00357   moab::ErrorCode rval = mkCore->moab_instance(meshIndex)->add_entities(moabEntSet, mesh_ents, num_ents);
00358   MBERRCHK(rval, mkCore->moab_instance(meshIndex));
00359   
00360   meshedState = mstate;
00361 }
00362   
00363 void ModelEnt::get_adjacencies(int dim, MEntVector &adjs) const 
00364 {
00365   std::vector<iGeom::EntityHandle> gents;
00366   // this call assumes that we have geometry. Is that always true?
00367   // it would be possible to have ModelEnt without geometry  (geometry handle null),
00368   // so then we should get adjacencies from mesh based geometry model.
00369   //we should really check if (igeomIndex>=0);
00370   iGeom::Error err = igeom_instance()->getEntAdj(geom_handle(), (iBase_EntityType)dim, gents);
00371   IBERRCHK(err, "Trouble getting geom adjacencies.");
00372   adjs.resize(gents.size());
00373   iGeom::TagHandle mkmodeltag;
00374   err = igeom_instance()->getTagHandle("__MKModelEntityGeo", mkmodeltag);
00375   IBERRCHK(err, "Failed to get tag handle for model entity.");
00376   err = igeom_instance()->getArrData(&gents[0], adjs.size(), mkmodeltag,
00377                                              &adjs[0]);
00378   IBERRCHK(err, "Trouble getting ModelEnts for geom adjacencies.");
00379 }
00380       
00381 int ModelEnt::dimension() const 
00382 {
00383   if (iGeomEnt) {
00384     iBase_EntityType tp;
00385     iGeom::Error err = igeom_instance()->getEntType(iGeomEnt, tp);
00386     IBERRCHK(err, "Couldn't get geom entity type.");
00387     return (tp - iBase_VERTEX);
00388   }
00389   else if (moabEntSet) {
00390     int dim;
00391     moab::ErrorCode rval = moab_instance()->tag_get_data(mkCore->moab_geom_dim_tag(), &moabEntSet, 1, &dim);
00392     MBERRCHK(rval, moab_instance());
00393     return dim;
00394   }
00395   else throw Error(MK_BAD_INPUT, "Couldn't get dimension for this ModelEnt.");
00396 }
00397   
00398 double ModelEnt::measure() const
00399 {
00400   double meas;
00401   iGeom::Error err = igeom_instance()->measure(&iGeomEnt, 1, &meas);
00402   IBERRCHK(err, "Couldn't get geom entity measure.");
00403   return meas;
00404 }
00405     
00406 double ModelEnt::measure_discrete() const
00407 {
00408   return measure();
00409 }
00410 
00411 void ModelEnt::evaluate(double x, double y, double z, 
00412                         double *close,
00413                         double *direction,
00414                         double *curvature1,
00415                         double *curvature2) const
00416 {
00417   iGeom::Error err = iBase_SUCCESS;
00418   if (0 == dimension() || 3 == dimension()) {
00419     if (direction || curvature1 || curvature2) {
00420       MKERRCHK(Error(MK_BAD_INPUT), "Direction or curvature not available for entities of this type.");
00421     }
00422     else if (!close) {
00423       MKERRCHK(Error(MK_BAD_INPUT), "Must pass closest point pointer for output.");
00424     }
00425   }
00426 
00427   if (0 == dimension()) {
00428     err = igeom_instance()->getVtxCoord(geom_handle(), close[0], close[1], close[2]);
00429     IBERRCHK(err, "Problem getting vertex coordinates.");
00430     return;
00431   }
00432   else if (3 == dimension()) {
00433     err = igeom_instance()->getEntClosestPt(geom_handle(), x, y, z,
00434                                                        close[0], close[1], close[2]);
00435     IBERRCHK(err, "Problem getting closest point for this model entity.");
00436   }
00437   else {
00438     double cls[3], dir[3], curve1[3], curve2[3];
00439     if (!close) close = cls;
00440     if (!direction) direction = dir;
00441     if (!curvature1) curvature1 = curve1;
00442     if (1 == dimension()) 
00443       err = igeom_instance()->getEgEvalXYZ(geom_handle(), x, y, z,
00444           close[0], close[1], close[2],
00445           direction[0], direction[1], direction[2],
00446           curvature1[0], curvature1[1], curvature1[2]);
00447     else if (2 == dimension()) {
00448       if (!curvature2) curvature2 = curve2;
00449       err = igeom_instance()->getFcEvalXYZ(geom_handle(), x, y, z,
00450         close[0], close[1], close[2],
00451         direction[0], direction[1], direction[2],
00452         curvature1[0], curvature1[1], curvature1[2],
00453         curvature2[0], curvature2[1], curvature2[2]);
00454     }
00455 
00456     IBERRCHK(err, "Couldn't evaluate model entity.");
00457   }
00458 }
00459 
00460 int ModelEnt::id() const
00461 {
00462     // get a global id for this entity, if there is one
00463   if (geom_handle()) {
00464     iGeom::TagHandle gid_tag;
00465     iGeom::Error err = igeom_instance()->getTagHandle("GLOBAL_ID", gid_tag);
00466     if (iBase_SUCCESS == err) {
00467       int gid;
00468       err = igeom_instance()->getIntData(geom_handle(), gid_tag, gid);
00469       if (iBase_SUCCESS == err) return gid;
00470     }
00471   }
00472   
00473   if (mesh_handle()) {
00474     moab::Tag gid_tag;
00475 
00476     moab::ErrorCode err = mk_core()->moab_instance()->tag_get_handle("GLOBAL_ID", 
00477                                                                      1, moab::MB_TYPE_INTEGER, 
00478                                                                      gid_tag, moab::MB_TAG_SPARSE);
00479     if (moab::MB_SUCCESS == err) {
00480       int gid;
00481       moab::EntityHandle this_mesh = mesh_handle();
00482       err = mk_core()->moab_instance()->tag_get_data(gid_tag, &this_mesh, 1, &gid);
00483       if (moab::MB_SUCCESS == err) return gid;
00484     }
00485   }
00486   
00487   return -1;
00488 }
00489 
00490 void ModelEnt::get_mesh(int dim,
00491                         std::vector<moab::EntityHandle> &ments,
00492                         bool bdy_too)
00493 {
00494   if (dim > dimension()) throw Error(MK_BAD_INPUT, "Called get_mesh for dimension %d on a %d-dimensional entity.", 
00495                                      dim, dimension());
00496 
00497   if (!bdy_too || dim == dimension()) {
00498       // just owned entities, which will be in the set
00499     moab::ErrorCode rval = mk_core()->moab_instance()->get_entities_by_dimension(moabEntSet, dimension(), ments);
00500     MBERRCHK(rval, mkCore->moab_instance());
00501     return;
00502   }
00503 
00504     // filter out the cases where order matters and boundary is requested too
00505   else if (1 == dimension() && 0 == dim) {
00506     std::vector<moab::EntityHandle> tmp_edgevs, tmp_vvs;
00507     moab::ErrorCode rval = mk_core()->moab_instance()->get_entities_by_dimension(moabEntSet, 0, tmp_edgevs);
00508     MBERRCHK(rval, mkCore->moab_instance());
00509     boundary(0, tmp_vvs);
00510     if (tmp_vvs.empty()) throw Error(MK_NOT_FOUND, "No mesh on bounding entity.");
00511     ments.push_back(tmp_vvs[0]);
00512     std::copy(tmp_edgevs.begin(), tmp_edgevs.end(), std::back_inserter(ments));
00513     if (2 == tmp_vvs.size()) ments.push_back(tmp_vvs[1]);
00514     else if (1 == tmp_vvs.size()) ments.push_back(tmp_vvs[0]);
00515     return;
00516   }
00517   else {
00518       // remaining case is bdy_too=true and dim < dimension()
00519       // first, get dimension()-dimensional entities, then get dim-dimensional neighbors of those
00520     std::vector<moab::EntityHandle> tmp_ments;
00521     moab::ErrorCode rval = mk_core()->moab_instance()->get_entities_by_dimension(moabEntSet, dimension(), tmp_ments);
00522     MBERRCHK(rval, mkCore->moab_instance());
00523     rval = mk_core()->moab_instance()->get_adjacencies(&tmp_ments[0], tmp_ments.size(), dim, false, ments, 
00524                                                        moab::Interface::UNION);
00525     MBERRCHK(rval, mkCore->moab_instance());
00526   }
00527 }
00528     
00529 void ModelEnt::boundary(int dim,
00530                         std::vector<moab::EntityHandle> &bdy,
00531                         std::vector<int> *senses,
00532                         std::vector<int> *group_sizes) 
00533 {
00534   MEntVector me_loop;
00535   std::vector<int> me_senses, me_group_sizes;
00536   boundary(dim, me_loop, &me_senses, &me_group_sizes);
00537   
00538   MEntVector::iterator vit = me_loop.begin();
00539   std::vector<int>::iterator sense_it = me_senses.begin(), grpsize_it = me_group_sizes.begin();
00540   int gents_ctr = 0, first_ment = 0;
00541 
00542     // packing into a single list, but need to keep track of loop/shell sizes; don't increment
00543     // grpsize_it here, just when we've done all the mes in this group
00544   for (vit = me_loop.begin(); vit != me_loop.end(); vit++, sense_it++) {
00545     std::vector<moab::EntityHandle> tmp_ments;
00546     (*vit)->get_mesh(dim, tmp_ments, true);
00547     if (*sense_it == SENSE_REVERSE) 
00548       std::reverse(tmp_ments.begin(), tmp_ments.end());
00549     if (2 == dimension() && 0 == dim)
00550         // assembling vertices on loops; don't do last, since that'll be the first of the
00551         // next edge
00552       std::copy(tmp_ments.begin(), tmp_ments.end()-1, std::back_inserter(bdy));
00553     else
00554       std::copy(tmp_ments.begin(), tmp_ments.end(), std::back_inserter(bdy));
00555     if (senses) {
00556       for (unsigned int i = 0; i < tmp_ments.size(); i++) senses->push_back(*sense_it);
00557     }
00558     if (group_sizes) {
00559       gents_ctr++;
00560       if (gents_ctr == *grpsize_it) {
00561         group_sizes->push_back(bdy.size()-first_ment);
00562           // first_ment will be the first element in the next group
00563         first_ment = bdy.size();
00564         gents_ctr = 0;
00565         grpsize_it++;
00566       } // block at end of a given loop
00567     } // block updating group size
00568   } // over loop
00569   
00570 }
00571     
00572 void ModelEnt::boundary(int dim,
00573                         MEntVector &entities,
00574                         std::vector<int> *senses,
00575                         std::vector<int> *group_sizes)
00576 {
00577     // for a given entity, return the bounding edges in the form of edge groups,
00578     // oriented ccw around surface
00579   int this_dim = dimension();
00580 
00581     // shouldn't be calling this function if we're a vertex
00582   if (this_dim == iBase_VERTEX) throw Error(MK_BAD_INPUT, "Shouldn't call this for a vertex.");
00583   else if (dim >= this_dim) throw Error(MK_BAD_INPUT, "Calling for boundary entities of equal or greater dimension.");;
00584   
00585     // get (d-1)-adj ents & senses
00586   MEntVector tmp_ents;
00587   get_adjacencies(this_dim-1, tmp_ents);
00588 
00589     // get out here if we're a curve
00590   if (1 == this_dim) {
00591     // it is important to order the vertices, especially for boundary
00592     // cases get_adjacency does not say anything about order
00593     if (tmp_ents.size() <= 1) // no worry
00594     {
00595       for (unsigned int i = 0; i < tmp_ents.size(); i++) {
00596         entities.push_back(tmp_ents[i]);
00597         if (senses) senses->push_back(SENSE_FORWARD);
00598       }
00599     }
00600     else
00601     {
00602       // we have at least 2 vertices; if more than 2, this is an error
00603       if (tmp_ents.size() > 2)
00604         throw Error(MK_FAILURE, " edge with too many vertices ");
00605       // we have now exactly 2 vertices; order them according to the
00606       // edge they come from
00607       int sense = 0;
00608       iGeom::Error rg = igeom_instance()->getEgVtxSense( geom_handle(),tmp_ents[0]->geom_handle(),
00609           tmp_ents[1]->geom_handle(),  sense );
00610       IBERRCHK(rg, "Trouble getting edge sense");
00611       if (-1==sense)
00612       {
00613         // the vertices are reversed in adjacency list
00614         entities.push_back(tmp_ents[1]);
00615         entities.push_back(tmp_ents[0]);
00616       }
00617       else
00618       {
00619         // vertices are fine
00620         entities.push_back(tmp_ents[0]);
00621         entities.push_back(tmp_ents[1]);
00622       }
00623       if (senses) {
00624         senses->push_back(SENSE_FORWARD);
00625         senses->push_back(SENSE_FORWARD);
00626       }
00627     }
00628     return; // we need to get out, we treated a dim 1 entity
00629   }
00630   
00631     // get adjacent entities into a sorted, mutable list
00632   MEntVector b_ents = tmp_ents;
00633   MEntSet dbl_curves, sgl_curves;
00634   
00635   std::sort(b_ents.begin(), b_ents.end());
00636 
00637   while (!b_ents.empty()) {
00638       // get 1st in group
00639     
00640     MEntVector group_stack, this_group;
00641     group_stack.push_back(b_ents.front());
00642     
00643     int current_senses_size = 0;
00644     if (senses)
00645       current_senses_size= senses->size();
00646       // while there are still entities on the stack
00647     while (!group_stack.empty()) {
00648       
00649         // pop one off & get its sense
00650       ModelEnt *this_entity = group_stack.back(); group_stack.pop_back();
00651       
00652       int this_sense;
00653       iGeom::Error err = 
00654           igeom_instance()->getSense(this_entity->geom_handle(), geom_handle(), this_sense);
00655       IBERRCHK(err,"Error getting geometry sense");
00656 
00657         // if we already have this one, continue
00658       if ((0 != this_sense && 
00659            std::find(this_group.begin(), this_group.end(), this_entity) != this_group.end()) ||
00660           std::find(dbl_curves.begin(), dbl_curves.end(), this_entity) != dbl_curves.end())
00661         continue;
00662 
00663         // either way we need the d-2 entities
00664       MEntVector bridges;
00665       this_entity->get_adjacencies(dimension()-2, bridges);
00666       if (dimension()-2==0 && bridges.size()==2)// edge case; this_entity is an edge
00667       {
00668         // if the vertices are reversed, change them back
00669         int sense_vertices = 0;
00670         iGeom::Error rg = igeom_instance()->getEgVtxSense(this_entity->geom_handle(),
00671             bridges[0]->geom_handle(),
00672             bridges[1]->geom_handle(),  sense_vertices );
00673         IBERRCHK(rg, "Trouble getting edge sense");
00674         if (-1==sense_vertices)
00675         {
00676           // reverse the order of bridges;
00677           // a better fix would be to have the bridges from get_adjacencies in order
00678           // this would mean that children of sets should be returned in the order they
00679           // were inserted as children (maybe too hard to control that)
00680           ModelEnt * tmp = bridges[0];
00681           bridges[0] = bridges[1];
00682           bridges[1] = tmp;
00683         }
00684 
00685       }
00686       
00687         // only remove from the list of candidates if it's not dual-sensed
00688       if (0 != this_sense) 
00689         b_ents.erase(std::remove(b_ents.begin(), b_ents.end(), this_entity), b_ents.end());
00690 
00691         // if it's double-sensed and this is the first time we're seeing it, find the right sense
00692       else {
00693         assert(dimension() == 2);
00694         if (std::find(sgl_curves.begin(), sgl_curves.end(), this_entity) == sgl_curves.end()) {
00695           sgl_curves.insert(this_entity);
00696           if (!this_group.empty()) {
00697             ModelEnt *common_v = this_entity->shared_entity(this_group.back(), 0);
00698             if (common_v == bridges[0]) this_sense = 1;
00699             else if (common_v == bridges[1]) this_sense = -1;
00700             else assert(false);
00701           }
00702             // else, if this is the first one in the loop, just choose a sense
00703           else {
00704             this_sense = 1;
00705           }
00706         }
00707         else {
00708             // else this is the second time we're seeing it, move it to the dbl_curves list and remove
00709             // from the candidates list
00710           dbl_curves.insert(this_entity);
00711           sgl_curves.erase(this_entity);
00712           b_ents.erase(std::remove(b_ents.begin(), b_ents.end(), this_entity), b_ents.end());
00713           if (senses) {
00714               // need to get the sense of the first appearance in list
00715             MEntVector::iterator vit = std::find(this_group.begin(), this_group.end(), this_entity);
00716             assert(vit != this_group.end());
00717             int other_sense = (*senses)[current_senses_size + vit-this_group.begin()];
00718             this_sense = SENSE_REVERSE*other_sense;
00719           }
00720         }
00721       }
00722           
00723         // it's in the group; put on group & remove from untreated ones
00724       this_group.push_back(this_entity);
00725       if (senses) senses->push_back(this_sense);
00726       MEntVector tmp_from, tmp_adjs;
00727 
00728         // if we're on a face and we're the first in a group, check sense of this first
00729         // edge; make sure "next" in loop sense is last on list
00730       if (2 == dimension() && this_group.size() == 1) {
00731 
00732           // get vertex which we know is shared by the "right" next edge; first get the vertices
00733           // if sense of current edge is forward, it's the 2nd vertex we want,
00734           // otherwise the first
00735         if (1 == this_sense && bridges.size() > 1) tmp_from.push_back(bridges[1]);
00736         else tmp_from.push_back(bridges[0]);
00737         tmp_adjs = b_ents;
00738         get_adjs_bool(tmp_from, dimension()-1, tmp_adjs, INTERSECT);
00739       }
00740       else {
00741         std::copy(bridges.begin(), bridges.end(), std::back_inserter(tmp_from));
00742         MEntVector tmp_adjs2;
00743         get_adjs_bool(tmp_from, dimension()-1, tmp_adjs2, UNION);
00744         std::sort(tmp_adjs2.begin(), tmp_adjs2.end());
00745         tmp_adjs.resize(tmp_adjs2.size());
00746         tmp_adjs.erase(std::set_intersection(b_ents.begin(), b_ents.end(),
00747                                              tmp_adjs2.begin(), tmp_adjs2.end(),
00748                                              tmp_adjs.begin()), tmp_adjs.end());
00749       }
00750 
00751       if (2 == dimension() && tmp_adjs.size() > 1) {
00752           // more than one adjacent edge - need to evaluate winding angle to find right one
00753         ModelEnt *next_ent = next_winding(this_entity, this_sense, tmp_adjs);
00754         if (NULL == next_ent) return;
00755         group_stack.push_back(next_ent);
00756       }
00757       else if (!tmp_adjs.empty()) 
00758         std::copy(tmp_adjs.begin(), tmp_adjs.end(), std::back_inserter(group_stack));
00759     }
00760     
00761     // put group in group list
00762     std::copy(this_group.begin(), this_group.end(), std::back_inserter(entities));
00763     if (group_sizes) group_sizes->push_back(this_group.size());
00764   }
00765 }
00766 
00767 ModelEnt *ModelEnt::shared_entity(ModelEnt *ent2, int to_dim) 
00768 {
00769     // find the shared entity between the two entities of the prescribed dimension
00770   MEntVector from_ents(2), to_ents;
00771   from_ents[0] = this;
00772   from_ents[1] = ent2;
00773   try {
00774     get_adjs_bool(from_ents, to_dim, to_ents, INTERSECT, false);
00775   }
00776   catch (Error err) {
00777     if (err.error_code() != MK_SUCCESS || to_ents.empty()) return NULL;
00778   }
00779 
00780   if (to_ents.size() > 1) throw Error(MK_MULTIPLE_FOUND, "Multiple shared entities found.");
00781   
00782   return to_ents[0];
00783 }
00784 
00785 void ModelEnt::get_adjs_bool(MEntVector &from_ents,
00786                              int to_dim,
00787                              MEntVector &to_ents,
00788                              BooleanType op_type,
00789                              bool only_to_ents) 
00790 {
00791   if (from_ents.empty()) {
00792     to_ents.clear();
00793     return;
00794   }
00795 
00796   MEntVector bridges;
00797 
00798   MEntVector::iterator from_it = from_ents.begin();
00799   if (to_ents.empty() && !only_to_ents && op_type == INTERSECT) {
00800     from_ents.front()->get_adjacencies(to_dim, bridges);
00801     std::copy(bridges.begin(), bridges.end(), std::back_inserter(to_ents));
00802     from_it++;
00803   }
00804 
00805   std::sort(to_ents.begin(), to_ents.end());
00806   MEntVector result_ents(to_ents.size());
00807   for (; from_it != from_ents.end(); from_it++) {
00808     bridges.clear();
00809     (*from_it)->get_adjacencies(to_dim, bridges);
00810     if (op_type == INTERSECT) {
00811       std::sort(bridges.begin(), bridges.end());
00812       result_ents.erase(std::set_intersection(to_ents.begin(), to_ents.end(),
00813                                               bridges.begin(), bridges.end(),
00814                                               result_ents.begin()), result_ents.end());
00815     }
00816     else {
00817       std::copy(bridges.begin(), bridges.end(), std::back_inserter(result_ents));
00818     }
00819     
00820     to_ents = result_ents;
00821   }
00822   
00823   if (op_type == UNION)
00824     std::sort(to_ents.begin(), to_ents.end());
00825   
00826   to_ents.erase(std::unique(to_ents.begin(), to_ents.end()), to_ents.end());
00827 }
00828 
00829 ModelEnt *ModelEnt::next_winding(ModelEnt *this_edge, 
00830                                  int this_sense, 
00831                                  MEntVector &tmp_adjs) 
00832 {
00833     // given this_entity, a d-1 entity bounding gentity, and optional "next"
00834     // entities in tmp_adjs, find the next one based on windings around the shared
00835     // vertex
00836   
00837     // first, get the shared vertex 
00838   MEntVector verts;
00839   this_edge->get_adjacencies(0, verts);
00840 
00841   ModelEnt *shared_vert = verts[0];
00842   if (this_sense == 1 && verts.size() > 1) shared_vert = verts[1];
00843 
00844     // get locations just before the vertex, at the vertex, and just after the vertex
00845   double v1[3], v2[3], v3[3];
00846   double umin, umax;
00847   igeom_instance()->getEntURange(this_edge->geom_handle(), umin, umax);
00848   double utgt;
00849   if (1 == this_sense) utgt = umin + 0.9 * (umax - umin);
00850   else utgt = umin + 0.1 * (umax - umin);
00851   igeom_instance()->getEntUtoXYZ(this_edge->geom_handle(), utgt, v1[0], v1[1], v1[2]);
00852   igeom_instance()->getVtxCoord(shared_vert->geom_handle(), v2[0], v2[1], v2[2]);
00853   double v21[] = {v1[0]-v2[0], v1[1]-v2[1], v1[2]-v2[2]};
00854   VecUtil::normalize(v21);
00855 
00856     // get the normal vector at the vertex
00857   double normal[3];
00858   igeom_instance()->getEntNrmlXYZ(geom_handle(), v2[0], v2[1], v2[2],
00859                                           normal[0], normal[1], normal[2]);
00860   
00861     // now loop over candidates, finding magnitude of swept angle
00862   ModelEnt *other = NULL;
00863   double angle = VecUtil::TWO_PI;
00864   bool is_adj;
00865   for (MEntVector::iterator vit = tmp_adjs.begin(); vit != tmp_adjs.end(); vit++) {
00866       // if we're here, we have multiple candidates, therefore don't choose the same one
00867     if (*vit == this_edge)
00868       continue;
00869     igeom_instance()->isEntAdj((*vit)->geom_handle(), shared_vert->geom_handle(), is_adj);
00870     if (!is_adj)
00871       continue;
00872 
00873       // get param range
00874     igeom_instance()->getEntURange((*vit)->geom_handle(), umin, umax);
00875       // get sense
00876     int tmp_sense;
00877     igeom_instance()->getSense((*vit)->geom_handle(), geom_handle(), tmp_sense);
00878     if (1 == tmp_sense) utgt = umin + 0.1 * (umax - umin);
00879     else utgt = umin + 0.9 * (umax - umin);
00880     igeom_instance()->getEntUtoXYZ((*vit)->geom_handle(), utgt, v3[0], v3[1], v3[2]);
00881     double v23[] = {v3[0]-v2[0], v3[1]-v2[1], v3[2]-v2[2]};
00882     VecUtil::normalize(v23);
00883 
00884     VecUtil::cross(normal, v23, v1);
00885     
00886     VecUtil::cross(v1, normal, v3);
00887 
00888     double x = VecUtil::dot(v21, v3);
00889     double y = VecUtil::dot(v21, v1);
00890 
00891     assert(x != 0.0 || y != 0.0);
00892     double this_angle = atan2(y, x);
00893     if (this_angle < 0.0)
00894     {
00895       this_angle += VecUtil::TWO_PI;
00896     }
00897     if (this_angle < angle) {
00898       other = *vit;
00899       angle = this_angle;
00900     }
00901   }
00902   
00903   return other;
00904 }
00905 
00906 void ModelEnt::boundary(int dim, 
00907                         moab::Range &ents) const
00908 {
00909   throw Error(MK_NOT_IMPLEMENTED, "Not implemented yet.");
00910 }
00911 
00912 void ModelEnt::get_indexed_connect_coords(std::vector<moab::EntityHandle> &ents,
00913                                           std::vector<int> *senses,
00914                                           moab::Tag tagh,
00915                                           std::vector<int> &ents_as_ids,
00916                                           std::vector<double> &coords,
00917                                           moab::Range *verts_range,
00918                                           int start_index) 
00919 {
00920     // if we need to, make a tag
00921   moab::ErrorCode rval;
00922   bool i_created_tag = false;
00923   if (0 == tagh) {
00924     //rval = mk_core()->moab_instance()->tag_create("__ModelEntidtag", sizeof(int), moab::MB_TAG_DENSE, 
00925     //                                              moab::MB_TYPE_INTEGER, tagh, NULL);
00926     rval = mk_core()->moab_instance()->tag_get_handle("__ModelEntidtag", 1, moab::MB_TYPE_INTEGER, tagh, 
00927                                                       moab::MB_TAG_DENSE|moab::MB_TAG_CREAT);
00928     MBERRCHK(rval, mkCore->moab_instance());
00929     i_created_tag = true;
00930   }
00931   
00932     // put entities into range, after clearing it
00933   moab::Range tmp_range, ents_range;
00934   if (!verts_range) verts_range = &tmp_range;
00935   verts_range->clear();
00936 
00937   std::copy(ents.begin(), ents.end(), moab::range_inserter(ents_range));
00938   bool all_verts = (mkCore->moab_instance()->type_from_handle(*ents_range.begin()) ==
00939                     mkCore->moab_instance()->type_from_handle(*ents_range.begin()) &&
00940                     mkCore->moab_instance()->type_from_handle(*ents_range.begin()) == moab::MBVERTEX);
00941 
00942 
00943     // get connectivity of all ents and store in range
00944   if (!all_verts) {
00945     rval = mkCore->moab_instance()->get_connectivity(ents_range, *verts_range);
00946     MBERRCHK(rval, mkCore->moab_instance());
00947   }
00948   else *verts_range = ents_range;
00949 
00950     // resize index array, to max number of connectivity entries
00951   int max_numconnect = moab::CN::VerticesPerEntity(mkCore->moab_instance()->type_from_handle(*(ents_range.rbegin())));
00952   ents_as_ids.resize(ents.size()*max_numconnect);
00953   
00954     // temporarily store ids in ents_as_ids array, and set id tag 
00955   assert(ents_as_ids.size() >= verts_range->size());
00956   for (unsigned int i = 0; i < verts_range->size(); i++) ents_as_ids[i] = i+start_index;
00957   rval = mk_core()->moab_instance()->tag_set_data(tagh, *verts_range, &ents_as_ids[0]);
00958   MBERRCHK(rval, mkCore->moab_instance());
00959 
00960     // get ids into ids vector in connectivity or ents order
00961   if (all_verts) {
00962     rval = mk_core()->moab_instance()->tag_get_data(tagh, &ents[0], ents.size(), &ents_as_ids[0]);
00963     MBERRCHK(rval, mkCore->moab_instance());
00964     ents_as_ids.resize(ents.size());
00965   }
00966   else {
00967       // loop over entities
00968     std::vector<moab::EntityHandle> storage;
00969     const moab::EntityHandle *connect;
00970     int num_connect;
00971     unsigned int last = 0;
00972     for (unsigned int i = 0; i < ents.size(); i++) {
00973         // first, get connect vector
00974       mk_core()->moab_instance()->get_connectivity(ents[i], connect, num_connect, true, &storage);
00975       MBERRCHK(rval, mk_core()->moab_instance());
00976       assert(ents_as_ids.size() >= last + num_connect);
00977         // next, get ids and put into ids vector
00978       rval = mk_core()->moab_instance()->tag_get_data(tagh, connect, num_connect, &ents_as_ids[last]);
00979       MBERRCHK(rval, mk_core()->moab_instance());
00980         // reverse if necessary
00981       assert(!senses || (*senses)[i] != SENSE_BOTH);
00982       if (senses && (*senses)[i] == SENSE_REVERSE) std::reverse(&ents_as_ids[last], &ents_as_ids[last+num_connect]);
00983         // update last
00984       last += num_connect;
00985     }
00986   }
00987 
00988     // get coords for range-ordered vertices
00989   coords.resize(3*verts_range->size());
00990   rval = mk_core()->moab_instance()->get_coords(*verts_range, &coords[0]);
00991   MBERRCHK(rval, mk_core()->moab_instance());
00992 
00993     // delete the tag if I created it
00994   if (i_created_tag) {
00995     rval = mk_core()->moab_instance()->tag_delete(tagh);
00996     MBERRCHK(rval, mkCore->moab_instance());
00997   }
00998 }
00999 
01000 iGeom::EntityHandle ModelEnt::geom_handle(moab::EntityHandle ment) const 
01001 {
01002   // use iRel to get this information or use model entity tag
01003   iBase_EntityHandle gent = 0;
01004   if (-1 != irelIndex)
01005   {
01006     Error err = mkCore->irel_pair(irelIndex)->getSetEntRelation(IBSH(ment), true, gent);
01007     MKERRCHK(err, "Failed to get geometry handle for mesh set.");
01008     return gent;
01009   }
01010   else
01011   {
01012     // use the model entity tag from moab instance... assume only one here
01013     //moab::Tag moabModelTag = mkCore->moab_model_tag();
01014     ModelEnt *modelEnt = NULL;
01015     moab::ErrorCode rval = mkCore->moab_instance()->tag_get_data(mkCore->moab_model_tag(), &ment, 1, &modelEnt);
01016     MBERRCHK(rval, mkCore->moab_instance());
01017     if (NULL == modelEnt)
01018        return gent;// still null
01019     return modelEnt->geom_handle();
01020   }
01021 }
01022 
01023 iGeom *ModelEnt::igeom_instance() const 
01024 {
01025   if (-1 == igeomIndex) throw Error(MK_FAILURE, "No geometry instance associated with this ModelEnt.");
01026   return mkCore->igeom_instance(igeomIndex);
01027 }
01028 
01029 moab::Interface *ModelEnt::moab_instance() const 
01030 {
01031   if (-1 == meshIndex) throw Error(MK_FAILURE, "No moab instance associated with this ModelEnt.");
01032   return mkCore->moab_instance(meshIndex);
01033 }
01034 
01035 moab::EntityHandle ModelEnt::mesh_handle(iGeom::EntityHandle gent) const 
01036 {
01037     // use iRel to get this information or use model entity tag
01038   moab::EntityHandle ment = 0;
01039   if (-1 != irelIndex)
01040   {
01041     iBase_EntitySetHandle h;
01042     iRel::Error err = mkCore->irel_pair(irelIndex)->getEntSetRelation(gent, false, h);
01043     IBERRCHK(err, "Failed to get mesh set handle for geometry entity.");
01044     return reinterpret_cast<moab::EntityHandle>(h);
01045   }
01046   else
01047   {
01048     // if no irel, use the model entity tag from the igeom instance
01049     //
01050     if (!igeom_instance())
01051       return ment;// NULL so far
01052     iGeom::TagHandle mkmodeltag;
01053     iBase_ErrorType err = igeom_instance()->getTagHandle("__MKModelEntityGeo", mkmodeltag);
01054     IBERRCHK(err, "Failed to get tag handle for model entity.");
01055     // now get the model entity
01056     ModelEnt * modelEnt = NULL;
01057     err = igeom_instance()->getData(gent, mkmodeltag, &modelEnt);
01058     if (NULL == modelEnt || iBase_TAG_NOT_FOUND == err)
01059        return ment;// still null
01060     return modelEnt->mesh_handle();
01061   }
01062 
01063 }
01064 
01065 moab::EntityHandle ModelEnt::mesh_handle(iGeom::EntitySetHandle gent) const 
01066 {
01067     // use iRel to get this information or use model entity tag
01068   moab::EntityHandle ment = 0;
01069   if (-1 != irelIndex)
01070   {
01071     iBase_EntitySetHandle h;
01072     iRel::Error err = mkCore->group_set_pair(irelIndex)->getSetSetRelation(gent, false, h);
01073     IBERRCHK(err, "Failed to get mesh set handle for geometry entity.");
01074     return reinterpret_cast<moab::EntityHandle>(h);
01075   }
01076   else
01077   {
01078     // if no irel, use the model entity tag from the igeom instance
01079     //
01080     if (-1 == igeomIndex)
01081       return ment;// NULL so far
01082     iGeom::TagHandle mkmodeltag;
01083     iBase_ErrorType err = igeom_instance()->getTagHandle("__MKModelEntityGeo", mkmodeltag);
01084     IBERRCHK(err, "Failed to get tag handle for model entity.");
01085     // now get the model entity
01086     ModelEnt * modelEnt = NULL;
01087     err = igeom_instance()->getEntSetData(gent, mkmodeltag, &modelEnt);
01088     if (NULL == modelEnt || iBase_TAG_NOT_FOUND == err)
01089        return ment;// still null
01090     return modelEnt->mesh_handle();
01091   }
01092 
01093 }
01094 
01099 double ModelEnt::mesh_interval_size() const 
01100 {
01101   if (-1 != sizing_function_index() && mk_core()->sizing_function(sizing_function_index()) &&
01102       -1 != mk_core()->sizing_function(sizing_function_index())->size())
01103     return mk_core()->sizing_function(sizing_function_index())->size();
01104   
01105   else if (1 == dimension() && -1 != mesh_intervals() && 0 != mesh_intervals())
01106     return measure() / mesh_intervals();
01107     
01108   else return -1;
01109 }
01110 
01111 } // namespace MeshKit
01112 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines