Skip to content
Snippets Groups Projects
Commit 660bc529 authored by Christophe Geuzaine's avatar Christophe Geuzaine
Browse files

work on createTopology

parent eeab1978
Branches
Tags
No related merge requests found
Show whitespace changes
Inline Side-by-side
Geo/GModel.cpp
+ 228
− 237
View file @ 660bc529
...@@ -46,88 +46,6 @@ ...@@ -46,88 +46,6 @@
std::vector<GModel*> GModel::list; std::vector<GModel*> GModel::list;
int GModel::_current = -1; int GModel::_current = -1;
static void recur_connect(MVertex *v,
std::multimap<MVertex*,MEdge> &v2e,
std::set<MEdge,Less_Edge> &group,
std::set<MVertex*> &touched)
{
if (touched.find(v) != touched.end()) return;
touched.insert(v);
for (std::multimap <MVertex*,MEdge>::iterator it = v2e.lower_bound(v);
it != v2e.upper_bound(v) ; ++it){
group.insert(it->second);
for (int i = 0; i < it->second.getNumVertices(); ++i){
recur_connect (it->second.getVertex(i), v2e, group, touched);
}
}
}
// starting form a list of elements, returns lists of lists that are
// all simply connected
static void recur_connect_e(const MEdge &e,
std::multimap<MEdge,MElement*, Less_Edge> &e2e,
std::set<MElement*> &group,
std::set<MEdge, Less_Edge> &touched)
{
if (touched.find(e) != touched.end())return;
touched.insert(e);
for (std::multimap <MEdge,MElement*,Less_Edge>::iterator it = e2e.lower_bound(e);
it != e2e.upper_bound(e) ; ++it){
group.insert(it->second);
for (int i = 0; i < it->second->getNumEdges(); ++i){
recur_connect_e(it->second->getEdge(i), e2e, group, touched);
}
}
}
static int connected_bounds(std::set<MEdge, Less_Edge> &edges,
std::vector<std::vector<MEdge> > &boundaries)
{
std::multimap<MVertex*,MEdge> v2e;
for(std::set<MEdge, Less_Edge>::iterator it = edges.begin(); it != edges.end(); it++){
for (int j=0;j<it->getNumVertices();j++){
v2e.insert(std::make_pair(it->getVertex(j),*it));
}
}
while (!v2e.empty()){
std::set<MEdge, Less_Edge> group;
std::set<MVertex*> touched;
recur_connect(v2e.begin()->first, v2e, group, touched);
std::vector<MEdge> temp;
temp.insert(temp.begin(), group.begin(), group.end());
boundaries.push_back(temp);
for (std::set<MVertex*>::iterator it = touched.begin() ; it != touched.end();++it)
v2e.erase(*it);
}
return boundaries.size();
}
static int connectedRegions(std::vector<MElement*> &elements,
std::vector<std::vector<MElement*> > &regions)
{
std::multimap<MEdge,MElement*,Less_Edge> e2e;
for (unsigned int i = 0; i < elements.size(); ++i){
for (int j = 0; j < elements[i]->getNumEdges(); j++){
e2e.insert(std::make_pair(elements[i]->getEdge(j),elements[i]));
}
}
while (!e2e.empty()){
std::set<MElement*> group;
std::set<MEdge,Less_Edge> touched;
recur_connect_e (e2e.begin()->first,e2e,group,touched);
std::vector<MElement*> temp;
temp.insert(temp.begin(), group.begin(), group.end());
regions.push_back(temp);
for ( std::set<MEdge,Less_Edge>::iterator it = touched.begin() ; it != touched.end();++it)
e2e.erase(*it);
}
return regions.size();
}
GModel::GModel(std::string name) GModel::GModel(std::string name)
: _name(name), _visible(1), _octree(0), : _name(name), _visible(1), _octree(0),
_geo_internals(0), _occ_internals(0), _acis_internals(0), _fm_internals(0), _geo_internals(0), _occ_internals(0), _acis_internals(0), _fm_internals(0),
...@@ -1133,37 +1051,92 @@ int GModel::removeDuplicateMeshVertices(double tolerance) ...@@ -1133,37 +1051,92 @@ int GModel::removeDuplicateMeshVertices(double tolerance)
return num; return num;
} }
void GModel::createTopologyFromMesh() static void recurConnectByEdge(const MEdge &e,
std::multimap<MEdge,MElement*, Less_Edge> &e2e,
std::set<MElement*> &group,
std::set<MEdge, Less_Edge> &touched)
{ {
if (touched.find(e) != touched.end())return;
touched.insert(e);
for (std::multimap <MEdge,MElement*,Less_Edge>::iterator it = e2e.lower_bound(e);
it != e2e.upper_bound(e) ; ++it){
group.insert(it->second);
for (int i = 0; i < it->second->getNumEdges(); ++i){
recurConnectByEdge(it->second->getEdge(i), e2e, group, touched);
}
}
}
Msg::StatusBar(2, true, "Creating topology from mesh..."); static int connectedSurfaces(std::vector<MElement*> &elements,
double t1 = Cpu(); std::vector<std::vector<MElement*> > &regions)
{
std::multimap<MEdge, MElement*, Less_Edge> e2e;
for(unsigned int i = 0; i < elements.size(); ++i){
for(int j = 0; j < elements[i]->getNumEdges(); j++){
e2e.insert(std::make_pair(elements[i]->getEdge(j), elements[i]));
}
}
while(!e2e.empty()){
std::set<MElement*> group;
std::set<MEdge, Less_Edge> touched;
recurConnectByEdge(e2e.begin()->first, e2e, group, touched);
std::vector<MElement*> temp;
temp.insert(temp.begin(), group.begin(), group.end());
regions.push_back(temp);
for(std::set<MEdge, Less_Edge>::iterator it = touched.begin() ; it != touched.end(); ++it)
e2e.erase(*it);
}
return regions.size();
}
removeDuplicateMeshVertices(CTX::instance()->geom.tolerance); static void recurConnectByVertex(MVertex *v,
std::multimap<MVertex*,MEdge> &v2e,
std::set<MEdge,Less_Edge> &group,
std::set<MVertex*> &touched)
{
if (touched.find(v) != touched.end()) return;
// for each discreteRegion, create topology touched.insert(v);
std::vector<discreteRegion*> discRegions; for (std::multimap <MVertex*,MEdge>::iterator it = v2e.lower_bound(v);
for(riter it = firstRegion(); it != lastRegion(); it++) it != v2e.upper_bound(v) ; ++it){
if((*it)->geomType() == GEntity::DiscreteVolume) group.insert(it->second);
discRegions.push_back((discreteRegion*) *it); for (int i = 0; i < it->second.getNumVertices(); ++i){
recurConnectByVertex(it->second.getVertex(i), v2e, group, touched);
}
}
}
//add mesh vertices static int connectedSurfaceBoundaries(std::set<MEdge, Less_Edge> &edges,
for (std::vector<discreteRegion*>::iterator it = discRegions.begin(); std::vector<std::vector<MEdge> > &boundaries)
it != discRegions.end(); it++){ {
for( std::vector<MVertex*>::const_iterator itv = (*it)->mesh_vertices.begin(); std::multimap<MVertex*,MEdge> v2e;
itv != (*it)->mesh_vertices.end(); itv++) for(std::set<MEdge, Less_Edge>::iterator it = edges.begin(); it != edges.end(); it++){
(*itv)->setEntity(*it); for (int j = 0; j < it->getNumVertices(); j++){
v2e.insert(std::make_pair(it->getVertex(j), *it));
}
}
while (!v2e.empty()){
std::set<MEdge, Less_Edge> group;
std::set<MVertex*> touched;
recurConnectByVertex(v2e.begin()->first, v2e, group, touched);
std::vector<MEdge> temp;
temp.insert(temp.begin(), group.begin(), group.end());
boundaries.push_back(temp);
for (std::set<MVertex*>::iterator it = touched.begin() ; it != touched.end();++it)
v2e.erase(*it);
} }
//for each discreteFace, createTopology return boundaries.size();
}
void GModel::makeDiscreteFacesSimplyConnected()
{
std::vector<discreteFace*> discFaces; std::vector<discreteFace*> discFaces;
for(fiter it = firstFace(); it != lastFace(); it++) for(fiter it = firstFace(); it != lastFace(); it++)
if((*it)->geomType() == GEntity::DiscreteSurface) if((*it)->geomType() == GEntity::DiscreteSurface)
discFaces.push_back((discreteFace*) *it); discFaces.push_back((discreteFace*) *it);
//--------
//for each discrete face, compute if it is made of connected faces
std::vector<discreteFace*> newDiscFaces;
for (std::vector<discreteFace*>::iterator itF = discFaces.begin(); for (std::vector<discreteFace*>::iterator itF = discFaces.begin();
itF != discFaces.end(); itF++){ itF != discFaces.end(); itF++){
...@@ -1171,59 +1144,87 @@ void GModel::createTopologyFromMesh() ...@@ -1171,59 +1144,87 @@ void GModel::createTopologyFromMesh()
for(unsigned int i = 0; i < (*itF)->getNumMeshElements(); i++) for(unsigned int i = 0; i < (*itF)->getNumMeshElements(); i++)
allElements.push_back((*itF)->getMeshElement(i)); allElements.push_back((*itF)->getMeshElement(i));
std::vector<std::vector<MElement*> > conRegions; std::vector<std::vector<MElement*> > conFaces;
int nbFaces = connectedRegions (allElements, conRegions); int nbFaces = connectedSurfaces(allElements, conFaces);
remove(*itF); remove(*itF);
for (int ifa = 0; ifa < nbFaces; ifa++){ for (int ifa = 0; ifa < nbFaces; ifa++){
std::vector<MElement*> myElements = conRegions[ifa]; std::vector<MElement*> myElements = conFaces[ifa];
int numF; int numF;
if (nbFaces == 1) numF = (*itF)->tag(); if (nbFaces == 1) numF = (*itF)->tag();
else numF = getMaxElementaryNumber(2) + 1; else numF = getMaxElementaryNumber(2) + 1;
discreteFace *f = new discreteFace(this, numF); discreteFace *f = new discreteFace(this, numF);
//printf("*** Created discreteFace %d \n", numF);
add(f); add(f);
newDiscFaces.push_back(f);
//fill new face with mesh vertices
std::set<MVertex*> allV; std::set<MVertex*> allV;
for(unsigned int i = 0; i < myElements.size(); i++) { for(unsigned int i = 0; i < myElements.size(); i++) {
MVertex *v0 = myElements[i]->getVertex(0); MElement *e = myElements[i];
MVertex *v1 = myElements[i]->getVertex(1); std::vector<MVertex*> verts;
MVertex *v2 = myElements[i]->getVertex(2); e->getVertices(verts);
f->triangles.push_back(new MTriangle( v0, v1, v2)); for(unsigned int k = 0; k < verts.size(); k++){
allV.insert(v0); allV.insert(v1); allV.insert(v2); allV.insert(verts[k]);
v0->setEntity(f); v1->setEntity(f); v2->setEntity(f); verts[k]->setEntity(f);
}
MElementFactory factory;
MElement *e2 = factory.create(e->getTypeForMSH(), verts, e->getNum(),
e->getPartition());
if(e2->getType() == TYPE_TRI)
f->triangles.push_back((MTriangle*)e2);
else
f->quadrangles.push_back((MQuadrangle*)e2);
} }
f->mesh_vertices.insert(f->mesh_vertices.begin(), allV.begin(), allV.end()); f->mesh_vertices.insert(f->mesh_vertices.begin(), allV.begin(), allV.end());
//delete new mesh vertices of face from adjacent regions
for (std::vector<discreteRegion*>::iterator itR = discRegions.begin();
itR != discRegions.end(); itR++){
for (std::set<MVertex*>::iterator itv = allV.begin();itv != allV.end(); itv++) {
std::vector<MVertex*>::iterator itve =
std::find((*itR)->mesh_vertices.begin(), (*itR)->mesh_vertices.end(), *itv);
if (itve != (*itR)->mesh_vertices.end()) (*itR)->mesh_vertices.erase(itve);
} }
} }
} }
} void GModel::createTopologyFromMesh()
discFaces = newDiscFaces; {
//------ Msg::StatusBar(2, true, "Creating topology from mesh...");
double t1 = Cpu();
//set boundary faces for each volume removeDuplicateMeshVertices(CTX::instance()->geom.tolerance);
makeDiscreteFacesSimplyConnected();
// create topology for all discrete regions
std::vector<discreteRegion*> discRegions;
for(riter it = firstRegion(); it != lastRegion(); it++)
if((*it)->geomType() == GEntity::DiscreteVolume)
discRegions.push_back((discreteRegion*) *it);
createTopologyFromRegions(discRegions);
// create topology for all discrete faces
std::vector<discreteFace*> discFaces;
for(fiter it = firstFace(); it != lastFace(); it++)
if((*it)->geomType() == GEntity::DiscreteSurface)
discFaces.push_back((discreteFace*) *it);
createTopologyFromFaces(discFaces);
// create old format (necessary for boundary layers)
exportDiscreteGEOInternals();
double t2 = Cpu();
Msg::StatusBar(2, true, "Done creating topology from mesh (%g s)", t2 - t1);
}
void GModel::createTopologyFromRegions(std::vector<discreteRegion*> &discRegions)
{
// find boundary faces of each region and put them in a map_faces that // find boundary faces of each region and put them in a map_faces that
// associates the MElements with the tags of the adjacent regions // associates the faces with the tags of the adjacent regions
std::map<MFace, std::vector<int>, Less_Face > map_faces; std::map<MFace, std::vector<int>, Less_Face > map_faces;
for (std::vector<discreteRegion*>::iterator it = discRegions.begin(); for (std::vector<discreteRegion*>::iterator it = discRegions.begin();
it != discRegions.end(); it++){ it != discRegions.end(); it++){
(*it)->findFaces(map_faces); (*it)->findFaces(map_faces);
} }
// get all discrete faces
std::vector<discreteFace*> discFaces;
for(fiter it = firstFace(); it != lastFace(); it++)
if((*it)->geomType() == GEntity::DiscreteSurface)
discFaces.push_back((discreteFace*) *it);
// create reverse map, for each region find set of MFaces that are // create reverse map, for each region find set of MFaces that are
// candidate for new discrete Face // candidate for new discrete Face
std::map<int, std::set<int> > region2Faces; std::map<int, std::set<int> > region2Faces;
...@@ -1235,10 +1236,12 @@ void GModel::createTopologyFromMesh() ...@@ -1235,10 +1236,12 @@ void GModel::createTopologyFromMesh()
std::map<MFace, std::vector<int>, Less_Face >::iterator itset = map_faces.find(mf); std::map<MFace, std::vector<int>, Less_Face >::iterator itset = map_faces.find(mf);
if (itset != map_faces.end()){ if (itset != map_faces.end()){
std::vector<int> tagRegions = itset->second; std::vector<int> tagRegions = itset->second;
for (std::vector<int>::iterator itR =tagRegions.begin(); itR != tagRegions.end(); itR++){ for (std::vector<int>::iterator itR = tagRegions.begin(); itR != tagRegions.end();
itR++){
std::map<int, std::set<int> >::iterator itmap = region2Faces.find(*itR); std::map<int, std::set<int> >::iterator itmap = region2Faces.find(*itR);
if (itmap == region2Faces.end()){ if (itmap == region2Faces.end()){
std::set<int> oneFace; oneFace.insert((*itF)->tag()); std::set<int> oneFace;
oneFace.insert((*itF)->tag());
region2Faces.insert(std::make_pair(*itR, oneFace)); region2Faces.insert(std::make_pair(*itR, oneFace));
} }
else{ else{
...@@ -1251,22 +1254,14 @@ void GModel::createTopologyFromMesh() ...@@ -1251,22 +1254,14 @@ void GModel::createTopologyFromMesh()
} }
} }
for (std::vector<discreteRegion*>::iterator it = discRegions.begin(); it != discRegions.end(); it++){ for (std::vector<discreteRegion*>::iterator it = discRegions.begin();
it != discRegions.end(); it++){
std::map<int, std::set<int> >::iterator itr = region2Faces.find((*it)->tag()); std::map<int, std::set<int> >::iterator itr = region2Faces.find((*it)->tag());
if (itr != region2Faces.end()){ if (itr != region2Faces.end()){
std::set<int> bcFaces = itr->second; std::set<int> bcFaces = itr->second;
(*it)->setBoundFaces(bcFaces); (*it)->setBoundFaces(bcFaces);
} }
} }
//create Topo for faces
createTopologyFromFaces(discFaces);
double t2 = Cpu();
Msg::StatusBar(2, true, "Done creating topology from mesh (%g s)", t2-t1);
//create old format (necessary for boundary layers)
exportDiscreteGEOInternals();
} }
void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
...@@ -1294,7 +1289,6 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1294,7 +1289,6 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
while (!map_edges.empty()){ while (!map_edges.empty()){
std::set<MEdge, Less_Edge> myEdges; std::set<MEdge, Less_Edge> myEdges;
myEdges.clear();
std::vector<int> tagFaces = map_edges.begin()->second; std::vector<int> tagFaces = map_edges.begin()->second;
myEdges.insert(map_edges.begin()->first); myEdges.insert(map_edges.begin()->first);
map_edges.erase(map_edges.begin()); map_edges.erase(map_edges.begin());
...@@ -1310,10 +1304,10 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1310,10 +1304,10 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
itmap++; itmap++;
} }
// if the loaded mesh already contains discrete Edges, check if // if the loaded mesh already contains discrete edges, check if
// the candidate discrete Edge does contain any of those; if not, // the candidate discrete edge does contain any of those; if not,
// create discreteEdges and create a map face2Edges that associate // create discreteEdges and create a map face2Edges that associate
// for each face the boundary discrete Edges // for each face the boundary discrete edges
for (std::vector<discreteEdge*>::iterator itE = discEdges.begin(); for (std::vector<discreteEdge*>::iterator itE = discEdges.begin();
itE != discEdges.end(); itE++){ itE != discEdges.end(); itE++){
...@@ -1336,45 +1330,42 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1336,45 +1330,42 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
std::set<MEdge, Less_Edge >::iterator itset = myEdges.find(me); std::set<MEdge, Less_Edge >::iterator itset = myEdges.find(me);
myEdges.erase(itset); myEdges.erase(itset);
} }
for (std::vector<int>::iterator itFace = tagFaces.begin(); for (std::vector<int>::iterator itFace = tagFaces.begin();
itFace != tagFaces.end(); itFace++) { itFace != tagFaces.end(); itFace++) {
std::map<int, std::vector<int> >::iterator it = face2Edges.find(*itFace); std::map<int, std::vector<int> >::iterator it = face2Edges.find(*itFace);
if (it == face2Edges.end()) face2Edges.insert(std::make_pair(*itFace, tagEdges)); if (it == face2Edges.end())
face2Edges.insert(std::make_pair(*itFace, tagEdges));
else{ else{
std::vector<int> allEdges = it->second; std::vector<int> allEdges = it->second;
allEdges.insert(allEdges.begin(), tagEdges.begin(), tagEdges.end()); allEdges.insert(allEdges.begin(), tagEdges.begin(), tagEdges.end());
it->second = allEdges; it->second = allEdges;
} }
// delete new mesh vertices of edge from adjacent faces // delete new mesh vertices of edge from adjacent faces
GFace *gf = getFaceByTag(*itFace); GFace *gf = getFaceByTag(*itFace);
for (std::vector<int>::iterator itEdge = tagEdges.begin(); for (std::vector<int>::iterator itEdge = tagEdges.begin();
itEdge != tagEdges.end(); itEdge++) { itEdge != tagEdges.end(); itEdge++) {
int count = 0;
GEdge *ge = getEdgeByTag(*itEdge); GEdge *ge = getEdgeByTag(*itEdge);
for(std::vector<MVertex*>::const_iterator itv = ge->mesh_vertices.begin(); for(std::vector<MVertex*>::const_iterator itv = ge->mesh_vertices.begin();
itv != ge->mesh_vertices.end(); itv++){ itv != ge->mesh_vertices.end(); itv++){
std::vector<MVertex*>::iterator itve = std::find std::vector<MVertex*>::iterator itve = std::find
(gf->mesh_vertices.begin(), gf->mesh_vertices.end(), *itv); (gf->mesh_vertices.begin(), gf->mesh_vertices.end(), *itv);
if (itve != gf->mesh_vertices.end()){ count++; gf->mesh_vertices.erase(itve);} if (itve != gf->mesh_vertices.end()){
gf->mesh_vertices.erase(itve);
}
} }
} }
} }
} }
} }
std::vector<std::vector<MEdge> > boundaries; std::vector<std::vector<MEdge> > boundaries;
int nbBounds = connected_bounds(myEdges, boundaries); int nbBounds = connectedSurfaceBoundaries(myEdges, boundaries);
for (int ib = 0; ib < nbBounds; ib++){ for (int ib = 0; ib < nbBounds; ib++){
std::vector<MEdge> myLines = boundaries[ib]; std::vector<MEdge> myLines = boundaries[ib];
int numE = getMaxElementaryNumber(1) + 1; int numE = getMaxElementaryNumber(1) + 1;
discreteEdge *e = new discreteEdge(this, numE, 0, 0); discreteEdge *e = new discreteEdge(this, numE, 0, 0);
//printf("*** Created discreteEdge %d \n", numE);
add(e); add(e);
discEdges.push_back(e); discEdges.push_back(e);
...@@ -1384,7 +1375,7 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1384,7 +1375,7 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
MVertex *v0 = myLines[i].getVertex(0); MVertex *v0 = myLines[i].getVertex(0);
MVertex *v1 = myLines[i].getVertex(1); MVertex *v1 = myLines[i].getVertex(1);
e->lines.push_back(new MLine(v0, v1)); e->lines.push_back(new MLine(v0, v1));
allV.insert(v0);//before it was std::vector with find ?? allV.insert(v0);
allV.insert(v1); allV.insert(v1);
v0->setEntity(e); v0->setEntity(e);
v1->setEntity(e); v1->setEntity(e);
...@@ -1394,7 +1385,7 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1394,7 +1385,7 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
for (std::vector<int>::iterator itFace = tagFaces.begin(); for (std::vector<int>::iterator itFace = tagFaces.begin();
itFace != tagFaces.end(); itFace++) { itFace != tagFaces.end(); itFace++) {
//delete new mesh vertices of edge from adjacent faces // delete mesh vertices of new edge from adjacent faces
GFace *dFace = getFaceByTag(*itFace); GFace *dFace = getFaceByTag(*itFace);
for (std::set<MVertex*>::iterator itv = allV.begin();itv != allV.end(); itv++) { for (std::set<MVertex*>::iterator itv = allV.begin();itv != allV.end(); itv++) {
std::vector<MVertex*>::iterator itve = std::vector<MVertex*>::iterator itve =
...@@ -1418,7 +1409,6 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1418,7 +1409,6 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
} }
} }
if (map_edges.empty()) break;
} }
// set boundary edges for each face // set boundary edges for each face
...@@ -1431,12 +1421,13 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1431,12 +1421,13 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
} }
} }
// for each discreteEdge, create Topology // for each discreteEdge, create topology
std::map<GFace*, std::map<MVertex*, MVertex*, std::less<MVertex*> > > face2Vert; std::map<GFace*, std::map<MVertex*, MVertex*, std::less<MVertex*> > > face2Vert;
std::map<GRegion*, std::map<MVertex*, MVertex*, std::less<MVertex*> > > region2Vert; std::map<GRegion*, std::map<MVertex*, MVertex*, std::less<MVertex*> > > region2Vert;
face2Vert.clear(); face2Vert.clear();
region2Vert.clear(); region2Vert.clear();
for (std::vector<discreteEdge*>::iterator it = discEdges.begin(); it != discEdges.end(); it++){ for (std::vector<discreteEdge*>::iterator it = discEdges.begin();
it != discEdges.end(); it++){
(*it)->createTopo(); (*it)->createTopo();
(*it)->parametrize(face2Vert, region2Vert); (*it)->parametrize(face2Vert, region2Vert);
} }
...@@ -1451,9 +1442,8 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1451,9 +1442,8 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
// edgeLoops.push_back(el); // edgeLoops.push_back(el);
// } // }
// we need to recreate lines, triangles and tets that contain those
//we need to recreate lines, triangles and tets // new MEdgeVertices
//that contain those new MEdgeVertices
for (std::map<GFace*, std::map<MVertex*, MVertex*, std::less<MVertex*> > >::iterator for (std::map<GFace*, std::map<MVertex*, MVertex*, std::less<MVertex*> > >::iterator
iFace = face2Vert.begin(); iFace != face2Vert.end(); iFace++){ iFace = face2Vert.begin(); iFace != face2Vert.end(); iFace++){
std::map<MVertex*, MVertex*, std::less<MVertex*> > old2new = iFace->second; std::map<MVertex*, MVertex*, std::less<MVertex*> > old2new = iFace->second;
...@@ -1462,25 +1452,26 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1462,25 +1452,26 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
std::vector<MQuadrangle*> newQuadrangles; std::vector<MQuadrangle*> newQuadrangles;
for (unsigned int i = 0; i < gf->getNumMeshElements(); ++i){ for (unsigned int i = 0; i < gf->getNumMeshElements(); ++i){
MElement *e = gf->getMeshElement(i); MElement *e = gf->getMeshElement(i);
int N = e->getNumVertices(); std::vector<MVertex *> v;
std::vector<MVertex *> v(N); e->getVertices(v);
for(int j = 0; j < N; j++) v[j] = e->getVertex(j); for (unsigned int j = 0; j < v.size(); j++){
for (int j = 0; j < N; j++){ std::map<MVertex*, MVertex*, std::less<MVertex*> >::iterator
std::map<MVertex*, MVertex*, std::less<MVertex*> >::iterator itmap = old2new.find(v[j]); itmap = old2new.find(v[j]);
MVertex *vNEW; if (itmap != old2new.end())
if (itmap != old2new.end()) { v[j] = itmap->second;
vNEW = itmap->second;
v[j]=vNEW;
} }
MElementFactory factory;
MElement *e2 = factory.create(e->getTypeForMSH(), v, e->getNum(),
e->getPartition());
switch(e2->getType()){
case TYPE_TRI: newTriangles.push_back((MTriangle*)e2); break;
case TYPE_QUA: newQuadrangles.push_back((MQuadrangle*)e2); break;
} }
if (N == 3) newTriangles.push_back(new MTriangle(v[0], v[1], v[2]));
else if ( N == 4) newQuadrangles.push_back(new MQuadrangle(v[0], v[1], v[2], v[3]));
} }
gf->deleteVertexArrays(); gf->deleteVertexArrays();
gf->triangles.clear(); for(unsigned int i = 0; i < gf->triangles.size(); i++) delete gf->triangles[i];
for(unsigned int i = 0; i < gf->quadrangles.size(); i++) delete gf->quadrangles[i];
gf->triangles = newTriangles; gf->triangles = newTriangles;
gf->quadrangles.clear();
gf->quadrangles = newQuadrangles; gf->quadrangles = newQuadrangles;
} }
...@@ -1494,34 +1485,34 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces) ...@@ -1494,34 +1485,34 @@ void GModel::createTopologyFromFaces(std::vector<discreteFace*> &discFaces)
std::vector<MPyramid*> newPyramids; std::vector<MPyramid*> newPyramids;
for (unsigned int i = 0; i < gr->getNumMeshElements(); ++i){ for (unsigned int i = 0; i < gr->getNumMeshElements(); ++i){
MElement *e = gr->getMeshElement(i); MElement *e = gr->getMeshElement(i);
int N = e->getNumVertices(); std::vector<MVertex *> v;
std::vector<MVertex *> v(N); e->getVertices(v);
for(int j = 0; j < N; j++) v[j] = e->getVertex(j); for (unsigned int j = 0; j < v.size(); j++){
for (int j = 0; j < N; j++){ std::map<MVertex*, MVertex*, std::less<MVertex*> >::iterator
std::map<MVertex*, MVertex*, std::less<MVertex*> >::iterator itmap = old2new.find(v[j]); itmap = old2new.find(v[j]);
MVertex *vNEW; if (itmap != old2new.end())
if (itmap != old2new.end()) { v[j] = itmap->second;
vNEW = itmap->second; }
v[j]=vNEW; MElementFactory factory;
} MElement *e2 = factory.create(e->getTypeForMSH(), v, e->getNum(),
} e->getPartition());
if (N == 4) newTetrahedra.push_back(new MTetrahedron(v[0], v[1], v[2], v[3])); switch(e2->getType()){
else if (N == 5) newPyramids.push_back(new MPyramid(v[0], v[1], v[2], v[3], v[4])); case TYPE_TET: newTetrahedra.push_back((MTetrahedron*)e2); break;
else if (N == 6) newPrisms.push_back(new MPrism(v[0], v[1], v[2], v[3], v[4], v[5])); case TYPE_HEX: newHexahedra.push_back((MHexahedron*)e2); break;
else if (N == 8) newHexahedra.push_back(new MHexahedron(v[0], v[1], v[2], v[3], case TYPE_PRI: newPrisms.push_back((MPrism*)e2); break;
v[4], v[5], v[6], v[7])); case TYPE_PYR: newPyramids.push_back((MPyramid*)e2); break;
}
} }
gr->deleteVertexArrays(); gr->deleteVertexArrays();
gr->tetrahedra.clear(); for(unsigned int i = 0; i < gr->tetrahedra.size(); i++) delete gr->tetrahedra[i];
for(unsigned int i = 0; i < gr->hexahedra.size(); i++) delete gr->hexahedra[i];
for(unsigned int i = 0; i < gr->prisms.size(); i++) delete gr->prisms[i];
for(unsigned int i = 0; i < gr->pyramids.size(); i++) delete gr->pyramids[i];
gr->tetrahedra = newTetrahedra; gr->tetrahedra = newTetrahedra;
gr->pyramids.clear();
gr->pyramids = newPyramids;
gr->prisms.clear();
gr->prisms = newPrisms;
gr->hexahedra.clear();
gr->hexahedra = newHexahedra; gr->hexahedra = newHexahedra;
gr->prisms = newPrisms;
gr->pyramids = newPyramids;
} }
} }
GModel *GModel::buildCutGModel(gLevelset *ls, bool cutElem) GModel *GModel::buildCutGModel(gLevelset *ls, bool cutElem)
......
Geo/GModel.h
+ 3
− 0
View file @ 660bc529
...@@ -27,6 +27,7 @@ class FieldManager; ...@@ -27,6 +27,7 @@ class FieldManager;
class CGNSOptions; class CGNSOptions;
class gLevelset; class gLevelset;
class discreteFace; class discreteFace;
class discreteRegion;
class binding; class binding;
class MElementOctree; class MElementOctree;
class GModelFactory; class GModelFactory;
...@@ -352,7 +353,9 @@ class GModel ...@@ -352,7 +353,9 @@ class GModel
// create topology from mesh // create topology from mesh
void createTopologyFromMesh(); void createTopologyFromMesh();
void createTopologyFromRegions(std::vector<discreteRegion*> &discRegions);
void createTopologyFromFaces(std::vector<discreteFace*> &pFaces); void createTopologyFromFaces(std::vector<discreteFace*> &pFaces);
void makeDiscreteFacesSimplyConnected();
// a container for smooth normals // a container for smooth normals
smooth_normals *normals; smooth_normals *normals;
......
Geo/MElement.h
+ 11
− 2
View file @ 660bc529
...@@ -75,7 +75,16 @@ class MElement ...@@ -75,7 +75,16 @@ class MElement
// get & set the vertices // get & set the vertices
virtual int getNumVertices() const = 0; virtual int getNumVertices() const = 0;
virtual MVertex *getVertex(int num) = 0; virtual MVertex *getVertex(int num) = 0;
virtual void setVertex(int num, MVertex *v) {throw;} void getVertices(std::vector<MVertex*> &verts)
{
int N = getNumVertices();
verts.resize(N);
for(int i = 0; i < N; i++) verts[i] = getVertex(i);
}
virtual void setVertex(int num, MVertex *v)
{
Msg::Error("Vertex set not supported for this element");
}
// give an MVertex as input and get its local number // give an MVertex as input and get its local number
virtual void getVertexInfo(const MVertex *vertex, int &ithVertex) const virtual void getVertexInfo(const MVertex *vertex, int &ithVertex) const
......
Geo/discreteFace.cpp
+ 17
− 16
View file @ 660bc529
...@@ -21,6 +21,16 @@ discreteFace::discreteFace(GModel *model, int num) : GFace(model, num) ...@@ -21,6 +21,16 @@ discreteFace::discreteFace(GModel *model, int num) : GFace(model, num)
meshStatistics.status = GFace::DONE; meshStatistics.status = GFace::DONE;
} }
void discreteFace::setBoundEdges(std::vector<int> tagEdges)
{
for (std::vector<int>::iterator it = tagEdges.begin(); it != tagEdges.end(); it++){
GEdge *ge = GModel::current()->getEdgeByTag(*it);
l_edges.push_back(ge);
l_dirs.push_back(1);
ge->addFace(this);
}
}
void discreteFace::findEdges(std::map<MEdge, std::vector<int>, Less_Edge> &map_edges) void discreteFace::findEdges(std::map<MEdge, std::vector<int>, Less_Edge> &map_edges)
{ {
std::set<MEdge, Less_Edge> bound_edges; std::set<MEdge, Less_Edge> bound_edges;
...@@ -29,12 +39,14 @@ void discreteFace::findEdges(std::map<MEdge, std::vector<int>, Less_Edge> &map_e ...@@ -29,12 +39,14 @@ void discreteFace::findEdges(std::map<MEdge, std::vector<int>, Less_Edge> &map_e
for (int iEdge = 0; iEdge < e->getNumEdges(); iEdge++) { for (int iEdge = 0; iEdge < e->getNumEdges(); iEdge++) {
MEdge tmp_edge = e->getEdge(iEdge); MEdge tmp_edge = e->getEdge(iEdge);
std::set<MEdge, Less_Edge >::iterator itset = bound_edges.find(tmp_edge); std::set<MEdge, Less_Edge >::iterator itset = bound_edges.find(tmp_edge);
if (itset == bound_edges.end()) bound_edges.insert(tmp_edge); if(itset == bound_edges.end())
else bound_edges.erase(itset); bound_edges.insert(tmp_edge);
else
bound_edges.erase(itset);
} }
} }
// for the boundary edges, associate the tag of the current discrete face // for the boundary edges, associate the tag of the discrete face
for (std::set<MEdge, Less_Edge>::iterator itv = bound_edges.begin(); for (std::set<MEdge, Less_Edge>::iterator itv = bound_edges.begin();
itv != bound_edges.end(); ++itv){ itv != bound_edges.end(); ++itv){
std::map<MEdge, std::vector<int>, Less_Edge >::iterator itmap = map_edges.find(*itv); std::map<MEdge, std::vector<int>, Less_Edge >::iterator itmap = map_edges.find(*itv);
...@@ -51,16 +63,6 @@ void discreteFace::findEdges(std::map<MEdge, std::vector<int>, Less_Edge> &map_e ...@@ -51,16 +63,6 @@ void discreteFace::findEdges(std::map<MEdge, std::vector<int>, Less_Edge> &map_e
} }
} }
void discreteFace::setBoundEdges(std::vector<int> tagEdges)
{
for (std::vector<int>::iterator it = tagEdges.begin(); it != tagEdges.end(); it++){
GEdge *ge = GModel::current()->getEdgeByTag(*it);
l_edges.push_back(ge);
l_dirs.push_back(1);
ge->addFace(this);
}
}
GPoint discreteFace::point(double par1, double par2) const GPoint discreteFace::point(double par1, double par2) const
{ {
Msg::Error("Cannot evaluate point on discrete face"); Msg::Error("Cannot evaluate point on discrete face");
...@@ -134,5 +136,4 @@ void discreteFace::writeGEO(FILE *fp) ...@@ -134,5 +136,4 @@ void discreteFace::writeGEO(FILE *fp)
count ++; count ++;
} }
fprintf(fp, "};\n"); fprintf(fp, "};\n");
} }
Geo/discreteRegion.cpp
+ 2
− 11
View file @ 660bc529
...@@ -16,24 +16,15 @@ discreteRegion::discreteRegion(GModel *model, int num) : GRegion(model, num) ...@@ -16,24 +16,15 @@ discreteRegion::discreteRegion(GModel *model, int num) : GRegion(model, num)
void discreteRegion::setBoundFaces(std::set<int> tagFaces) void discreteRegion::setBoundFaces(std::set<int> tagFaces)
{ {
for (std::set<int>::iterator it = tagFaces.begin() ; it != tagFaces.end();++it){ for (std::set<int>::iterator it = tagFaces.begin() ; it != tagFaces.end();++it){
GFace *face = model()->getFaceByTag(*it); GFace *face = model()->getFaceByTag(*it);
l_faces.push_back(face); l_faces.push_back(face);
face->addRegion(this); face->addRegion(this);
} }
// in case discrete region already exist
// to modify to take into account appropriate faces
// for(GModel::fiter face = model()->firstFace(); face != model()->lastFace(); face++){
// l_faces.push_back(*face);
// (*face)->addRegion(this);
// }
} }
void discreteRegion::findFaces(std::map<MFace, std::vector<int>, Less_Face> &map_faces) void discreteRegion::findFaces(std::map<MFace, std::vector<int>, Less_Face> &map_faces)
{ {
std::set<MFace, Less_Face> bound_faces; std::set<MFace, Less_Face> bound_faces;
for (unsigned int elem = 0; elem < getNumMeshElements() ; elem++) { for (unsigned int elem = 0; elem < getNumMeshElements() ; elem++) {
MElement *e = getMeshElement(elem); MElement *e = getMeshElement(elem);
...@@ -45,7 +36,7 @@ void discreteRegion::findFaces(std::map<MFace, std::vector<int>, Less_Face> &map ...@@ -45,7 +36,7 @@ void discreteRegion::findFaces(std::map<MFace, std::vector<int>, Less_Face> &map
} }
} }
// for the boundary faces, associate the tag of the current discrete face // for the boundary faces, associate the tag of the discrete face
for (std::set<MFace, Less_Face>::iterator itv = bound_faces.begin(); for (std::set<MFace, Less_Face>::iterator itv = bound_faces.begin();
itv != bound_faces.end(); ++itv){ itv != bound_faces.end(); ++itv){
std::map<MFace, std::vector<int>, Less_Face >::iterator itmap = map_faces.find(*itv); std::map<MFace, std::vector<int>, Less_Face >::iterator itmap = map_faces.find(*itv);
...@@ -60,5 +51,5 @@ void discreteRegion::findFaces(std::map<MFace, std::vector<int>, Less_Face> &map ...@@ -60,5 +51,5 @@ void discreteRegion::findFaces(std::map<MFace, std::vector<int>, Less_Face> &map
itmap->second = tagRegions; itmap->second = tagRegions;
} }
} }
} }
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment