diff --git a/Geo/GModel.cpp b/Geo/GModel.cpp
index 4157e80cfda8c2bf2498598a9af7d9a295b1220c..95f0fa35ce593ee6b59cea4b8099acc16e6e946c 100644
--- a/Geo/GModel.cpp
+++ b/Geo/GModel.cpp
@@ -1268,6 +1268,14 @@ void GModel::_storeElementsInEntities(std::map< int, std::vector<MElement* > >&
}
}
+void GModel::_storeParentsInSubElements(std::map<int, std::vector<MElement*> > &map)
+{
+ std::map<int, std::vector<MElement*> >::const_iterator it;
+ for(it = map.begin(); it != map.end(); ++it)
+ for (int i=0; i<it->second.size(); ++i)
+ it->second[i]->updateParent(this);
+}
+
template<class T>
static void _associateEntityWithElementVertices(GEntity *ge, std::vector<T*> &elements,
bool force=false)
diff --git a/Geo/GModel.h b/Geo/GModel.h
index 9beb563afedc1cbc43fa1637d4408dddde11383b..96722dde8d249c3394d1bef3f0cbbf64b425bf92 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -119,6 +119,9 @@ class GModel
// the fly if needed
void _storeElementsInEntities(std::map<int, std::vector<MElement*> > &map);
+ // store the parent's pointer back into MSubElements (replacing numeric id)
+ void _storeParentsInSubElements(std::map< int, std::vector<MElement* > >& map);
+
// Discrete Entities have to have their mesh moved to a geometry container
void _createGeometryOfDiscreteEntities(bool force=false);
diff --git a/Geo/GModelIO_MSH.cpp b/Geo/GModelIO_MSH.cpp
index 4c986a60a48c28d16c47a79612f1e42f9e799b72..7cf835757f5ddbde1facf2e039cbec8e425bd308 100644
--- a/Geo/GModelIO_MSH.cpp
+++ b/Geo/GModelIO_MSH.cpp
@@ -354,6 +354,7 @@ int GModel::readMSH(const std::string &name)
_vertexMapCache.clear();
}
}
+
// $Elements section
else if(!strncmp(&str[1], "Elements", 8)) {
if(!fgets(str, sizeof(str), fp)){ fclose(fp); return 0; }
@@ -439,7 +440,7 @@ int GModel::readMSH(const std::string &name)
// entity does not exist, create a new (discrete) one.
for(int i = 0; i < (int)(sizeof(elements) / sizeof(elements[0])); i++)
_storeElementsInEntities(elements[i]);
-
+
// associate the correct geometrical entity with each mesh vertex
_associateEntityWithMeshVertices();
@@ -451,6 +452,9 @@ int GModel::readMSH(const std::string &name)
_createGeometryOfDiscreteEntities() ;
+ for(int i = 0; i < (int)(sizeof(elements) / sizeof(elements[0])); i++)
+ _storeParentsInSubElements(elements[i]);
+
fclose(fp);
return postpro ? 2 : 1;
diff --git a/Geo/GModelIO_MSH2.cpp b/Geo/GModelIO_MSH2.cpp
index 5304577b48bf91f2c9d193fcf06911b9ea8046e7..322c7520a380b118bde62c61e11cbe80cfc7ce65 100644
--- a/Geo/GModelIO_MSH2.cpp
+++ b/Geo/GModelIO_MSH2.cpp
@@ -73,7 +73,7 @@ static MElement *createElementMSH2(GModel *m, int num, int typeMSH, int physical
}
MElementFactory factory;
- MElement *e = factory.create(typeMSH, v, num, part, owner, parent, d1, d2);
+ MElement *e = factory.create(typeMSH, v, num, part, owner, 0, parent, d1, d2);
if(!e){
Msg::Error("Unknown type of element %d", typeMSH);
diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index 352467f3b776ba290211b805b7b37ab09f4ae7b4..e50141e44f52de1ca24223d9954e61eb15d40e16 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -1627,7 +1627,7 @@ MElement *MElement::copy(std::map<int, MVertex*> &vertexMap,
}
MElementFactory f;
- MElement *newEl = f.create(eType, vmv, getNum(), _partition, ownsParent(), parent);
+ MElement *newEl = f.create(eType, vmv, getNum(), _partition, ownsParent(), 0, parent);
for(int i = 0; i < 2; i++) {
MElement *dom = getDomain(i);
@@ -1646,7 +1646,8 @@ MElement *MElement::copy(std::map<int, MVertex*> &vertexMap,
}
MElement *MElementFactory::create(int type, std::vector<MVertex*> &v,
- int num, int part, bool owner, MElement *parent,
+ int num, int part, bool owner,
+ int parent, MElement* parent_ptr,
MElement *d1, MElement *d2)
{
switch (type) {
@@ -1662,7 +1663,7 @@ MElement *MElementFactory::create(int type, std::vector<MVertex*> &v,
case MSH_LIN_10: return new MLineN(v, num, part);
case MSH_LIN_11: return new MLineN(v, num, part);
case MSH_LIN_B: return new MLineBorder(v, num, part, d1, d2);
- case MSH_LIN_C: return new MLineChild(v, num, part, owner, parent);
+ case MSH_LIN_C: return new MLineChild(v, num, part, owner, parent_ptr);
case MSH_TRI_3: return new MTriangle(v, num, part);
case MSH_TRI_6: return new MTriangle6(v, num, part);
case MSH_TRI_10: return new MTriangleN(v, 3, num, part);
@@ -1701,7 +1702,7 @@ MElement *MElementFactory::create(int type, std::vector<MVertex*> &v,
case MSH_QUA_32: return new MQuadrangleN(v, 8, num, part);
case MSH_QUA_36I: return new MQuadrangleN(v, 9, num, part);
case MSH_QUA_40: return new MQuadrangleN(v,10, num, part);
- case MSH_POLYG_: return new MPolygon(v, num, part, owner, parent);
+ case MSH_POLYG_: return new MPolygon(v, num, part, owner, parent_ptr);
case MSH_POLYG_B: return new MPolygonBorder(v, num, part, d1, d2);
case MSH_TET_4: return new MTetrahedron(v, num, part);
case MSH_TET_10: return new MTetrahedron10(v, num, part);
@@ -1742,7 +1743,7 @@ MElement *MElementFactory::create(int type, std::vector<MVertex*> &v,
case MSH_TET_46: return new MTetrahedronN(v, 8, num, part);
case MSH_TET_52: return new MTetrahedronN(v, 9, num, part);
case MSH_TET_58: return new MTetrahedronN(v, 10, num, part);
- case MSH_POLYH_: return new MPolyhedron(v, num, part, owner, parent);
+ case MSH_POLYH_: return new MPolyhedron(v, num, part, owner, parent_ptr);
case MSH_HEX_32: return new MHexahedronN(v, 3, num, part);
case MSH_HEX_64: return new MHexahedronN(v, 3, num, part);
case MSH_HEX_125: return new MHexahedronN(v, 4, num, part);
@@ -1751,10 +1752,14 @@ MElement *MElementFactory::create(int type, std::vector<MVertex*> &v,
case MSH_HEX_512: return new MHexahedronN(v, 7, num, part);
case MSH_HEX_729: return new MHexahedronN(v, 8, num, part);
case MSH_HEX_1000:return new MHexahedronN(v, 9, num, part);
- case MSH_PNT_SUB: return new MSubPoint(v, num, part, owner, parent);
- case MSH_LIN_SUB: return new MSubLine(v, num, part, owner, parent);
- case MSH_TRI_SUB: return new MSubTriangle(v, num, part, owner, parent);
- case MSH_TET_SUB: return new MSubTetrahedron(v, num, part, owner, parent);
+ case MSH_PNT_SUB: return (parent_ptr) ? new MSubPoint(v, num, part, owner, parent_ptr)
+ : new MSubPoint(v, num, part, owner, parent);
+ case MSH_LIN_SUB: return (parent_ptr) ? new MSubLine(v, num, part, owner, parent_ptr)
+ : new MSubLine(v, num, part, owner, parent);
+ case MSH_TRI_SUB: return (parent_ptr) ? new MSubTriangle(v, num, part, owner, parent_ptr)
+ : new MSubTriangle(v, num, part, owner, parent);
+ case MSH_TET_SUB: return (parent_ptr) ? new MSubTetrahedron(v, num, part, owner, parent_ptr)
+ : new MSubTetrahedron(v, num, part, owner, parent);
case MSH_PYR_5: return new MPyramid(v, num, part);
case MSH_PYR_13: return new MPyramidN(v, 2, num, part);
case MSH_PYR_14: return new MPyramidN(v, 2, num, part);
@@ -1805,10 +1810,10 @@ MElement *MElementFactory::create(int num, int type, const std::vector<int> &dat
int part = 0;
int startPartitions = startVertices + numVertices;
- MElement *parent = 0;
+ int parent = 0;
if((type == MSH_PNT_SUB || type == MSH_LIN_SUB ||
type == MSH_TRI_SUB || type == MSH_TET_SUB)){
- parent = model->getMeshElementByTag(data[startPartitions]);
+ parent = data[startPartitions];
startPartitions += 1;
}
diff --git a/Geo/MElement.h b/Geo/MElement.h
index fda1ade945ebb1b1e218683571332657ebf71056..7aff9cc4dc4bff9cfa5d0cd8933194408deffde2 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -174,6 +174,7 @@ class MElement
// get and set parent and children for hierarchial grids
virtual MElement *getParent() const { return NULL; }
virtual void setParent(MElement *p, bool owner = false) {}
+ virtual void updateParent(GModel *gm) {} // used only for reading msh files
virtual int getNumChildren() const { return 0; }
virtual MElement *getChild(int i) const { return NULL; }
virtual bool ownsParent() const { return false; }
@@ -402,8 +403,8 @@ class MElement
class MElementFactory{
public:
- MElement *create(int type, std::vector<MVertex*> &v, int num=0, int part=0,
- bool owner=false, MElement *parent=0, MElement *d1=0, MElement *d2=0);
+ MElement *create(int type, std::vector< MVertex* > &v, int num=0, int part=0,
+ bool owner=false, int parent=0, MElement *parent_ptr=NULL, MElement *d1 = 0, MElement *d2 = 0);
MElement *create(int num, int type, const std::vector<int> &data, GModel *model);
};
diff --git a/Geo/MSubElement.cpp b/Geo/MSubElement.cpp
index 737251e3718488f4661eb1be841ecaa9250f318e..4526678a8e97e008d20bbc3b9ca73e5c33465166 100644
--- a/Geo/MSubElement.cpp
+++ b/Geo/MSubElement.cpp
@@ -8,6 +8,7 @@
#include "MSubElement.h"
#include "Numeric.h"
+#include "GModel.h"
// MSubTetrahedron
@@ -17,6 +18,11 @@ MSubTetrahedron::~MSubTetrahedron()
if(_base) delete _base;
}
+void MSubTetrahedron::updateParent(GModel *gm)
+{
+ _orig=gm->getMeshElementByTag(_orig_N);
+}
+
const nodalBasis* MSubTetrahedron::getFunctionSpace(int order, bool serendip) const
{
if(_orig) return _orig->getFunctionSpace(order, serendip);
@@ -209,6 +215,12 @@ MSubTriangle::~MSubTriangle()
if(_base) delete _base;
}
+void MSubTriangle::updateParent(GModel *gm)
+{
+ _orig=gm->getMeshElementByTag(_orig_N);
+}
+
+
const nodalBasis* MSubTriangle::getFunctionSpace(int order, bool serendip) const
{
if(_orig) return _orig->getFunctionSpace(order, serendip);
@@ -447,6 +459,11 @@ MSubLine::~MSubLine()
if(_base) delete _base;
}
+void MSubLine::updateParent(GModel *gm)
+{
+ _orig=gm->getMeshElementByTag(_orig_N);
+}
+
const nodalBasis* MSubLine::getFunctionSpace(int order, bool serendip) const
{
if(_orig) return _orig->getFunctionSpace(order, serendip);
@@ -676,6 +693,11 @@ MSubPoint::~MSubPoint()
if(_base) delete _base;
}
+void MSubPoint::updateParent(GModel *gm)
+{
+ _orig=gm->getMeshElementByTag(_orig_N);
+}
+
const nodalBasis* MSubPoint::getFunctionSpace(int order, bool serendip) const
{
if(_orig) return _orig->getFunctionSpace(order, serendip);
diff --git a/Geo/MSubElement.h b/Geo/MSubElement.h
index f54dae3f37aab01cf1bdc9de1c1f6e0f9b2398d6..82c17dc6e56df9850632948eda306aa26fa62b45 100644
--- a/Geo/MSubElement.h
+++ b/Geo/MSubElement.h
@@ -19,16 +19,26 @@
// A sub tetrahedron, contained in another element
class MSubTetrahedron : public MTetrahedron
{
+ friend class MElementFactory;
+ friend class GModel;
protected:
bool _owner;
- MElement* _orig;
- std::vector<MElement*> _parents;
+ union
+ {
+ MElement* _orig; // pointer to parent
+ int _orig_N ; // element number of parent (only used while reading file...)
+ };
+ std::vector<MElement*> _parents; // not used in msh 3 file format
mutable MElement* _base;
int _pOrder;
int _npts;
IntPt *_pts;
+ MSubTetrahedron(std::vector<MVertex*> v, int num, int part, bool owner, int orig)
+ : MTetrahedron(v, num, part), _owner(owner), _orig_N(orig), _base(0), _pOrder(-1), _npts(0), _pts(0) {}
+ virtual void updateParent(GModel *gm); // NEVER ever use this ! (except for reading msh files !)
+
public:
MSubTetrahedron(MVertex *v0, MVertex *v1, MVertex *v2, MVertex *v3, int num=0, int part=0,
bool owner=false, MElement* orig=NULL)
@@ -77,16 +87,25 @@ class MSubTetrahedron : public MTetrahedron
// A sub triangle, contained in another element
class MSubTriangle : public MTriangle
{
+ friend class MElementFactory;
+ friend class GModel;
protected:
bool _owner;
- MElement* _orig;
- std::vector<MElement*> _parents;
+ union
+ {
+ MElement* _orig; // pointer to parent
+ int _orig_N ; // element number of parent (only used while reading file...)
+ };
+ std::vector<MElement*> _parents; // not used in msh 3 file format
mutable MElement* _base;
int _pOrder;
int _npts;
IntPt *_pts;
-
+
+ MSubTriangle(std::vector<MVertex*> v, int num, int part, bool owner, int orig)
+ : MTriangle(v, num, part), _owner(owner), _orig_N(orig), _base(0), _pOrder(-1), _npts(0), _pts(0) {}
+ virtual void updateParent(GModel *gm); // NEVER ever use this ! (except for reading msh files !)
public:
MSubTriangle(MVertex *v0, MVertex *v1, MVertex *v2, int num=0, int part=0,
bool owner=false, MElement* orig=NULL)
@@ -135,16 +154,25 @@ class MSubTriangle : public MTriangle
// A sub line, contained in another element
class MSubLine : public MLine
{
+ friend class MElementFactory;
+ friend class GModel;
protected:
bool _owner;
- MElement* _orig;
- std::vector<MElement*> _parents;
+ union
+ {
+ MElement* _orig; // pointer to parent
+ int _orig_N ; // element number of parent (only used while reading file...)
+ };
+ std::vector<MElement*> _parents; // not used in msh 3 file format
mutable MElement* _base;
int _pOrder;
int _npts;
IntPt *_pts;
-
+
+ MSubLine(std::vector<MVertex*> v, int num, int part, bool owner, int orig)
+ : MLine(v, num, part), _owner(owner), _orig_N(orig), _base(0), _pOrder(-1), _npts(0), _pts(0) {}
+ virtual void updateParent(GModel *gm); // NEVER ever use this ! (except for reading msh files !)
public:
MSubLine(MVertex *v0, MVertex *v1, int num=0, int part=0,
bool owner=false, MElement* orig=NULL)
@@ -193,16 +221,25 @@ class MSubLine : public MLine
// A sub point, contained in another element
class MSubPoint : public MPoint
{
+ friend class MElementFactory;
+ friend class GModel;
protected:
bool _owner;
- MElement* _orig;
- std::vector<MElement*> _parents;
+ union
+ {
+ MElement* _orig; // pointer to parent
+ int _orig_N ; // element number of parent (only used while reading file...)
+ };
+ std::vector<MElement*> _parents; // not used in msh 3 file format
mutable MElement* _base;
int _pOrder;
int _npts;
IntPt *_pts;
-
+
+ MSubPoint(std::vector<MVertex*> v, int num, int part, bool owner, int orig)
+ : MPoint(v, num, part), _owner(owner), _orig_N(orig), _base(0), _pOrder(-1), _npts(0), _pts(0) {}
+ virtual void updateParent(GModel *gm); // NEVER ever use this ! (except for reading msh files !)
public:
MSubPoint(MVertex *v0, int num=0, int part=0,
bool owner=false, MElement* orig=NULL)