diff --git a/Geo/GModel.h b/Geo/GModel.h
index a71531f51e23b69f85380ec5afa7df476586476e..a6347c63442932e8373ae040da92d14e634b66fb 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -50,6 +50,10 @@ class GModel
std::vector<MElement*> _elementVectorCache;
std::map<int, MElement*> _elementMapCache;
+ // ghost cell information (stores partitions for each element acting
+ // as a ghost cell)
+ std::multimap<MElement*, short> _ghostCells;
+
// an octree for fast mesh element lookup
Octree *_octree;
@@ -314,6 +318,8 @@ class GModel
int getMinPartitionSize() const { return partitionSize[0]; }
int getMaxPartitionSize() const { return partitionSize[1]; }
+ std::multimap<MElement*, short> &getGhostCells(){ return _ghostCells; }
+
// perform various coherence tests on the mesh
void checkMeshCoherence(double tolerance);
diff --git a/Mesh/meshPartition.cpp b/Mesh/meshPartition.cpp
index ff9821f0a8436bf6bd9480cd82794b761bea9890..458b364a309e400ff2d0a61837131a426d7360ae 100644
--- a/Mesh/meshPartition.cpp
+++ b/Mesh/meshPartition.cpp
@@ -75,8 +75,16 @@ template <> struct LFaceTr<MFace>
template <unsigned DIM>
struct MakeGraphFromEntity;
+
template <unsigned DIM>
struct MatchBoElemToGrVertex;
+
+struct BoElemGr;
+typedef std::vector<BoElemGr> BoElemGrVec;
+
+int MakeGraph(GModel *const model, Graph &graph,
+ BoElemGrVec *const boElemGrVec = 0);
+
template <unsigned DIM, typename EntIter, typename EntIterBE>
void MakeGraphDIM(const EntIter begin, const EntIter end,
const EntIterBE beginBE, const EntIterBE endBE,
@@ -104,36 +112,35 @@ void MakeGraphDIM(const EntIter begin, const EntIter end,
-int PartitionMeshElements( std::vector<MElement*> &elements, meshPartitionOptions &options){
-
- GModel *tmp_model = new GModel();
- GFace *gf = new discreteFace(tmp_model, 1);
- std::set<MVertex *> setv;
- for (unsigned i=0;i<elements.size();++i)
- for (int j=0;j<elements[i]->getNumVertices();j++)
- setv.insert(elements[i]->getVertex(j));
-
- for (std::set<MVertex* >::iterator it = setv.begin(); it != setv.end(); it++)
- gf->mesh_vertices.push_back(*it);
-
- for (std::vector<MElement* >::iterator it = elements.begin(); it != elements.end(); it++){
- if ((*it)->getType() == TYPE_TRI)
- gf->triangles.push_back((MTriangle*)(*it));
- else if ((*it)->getType() == TYPE_QUA)
- gf->quadrangles.push_back((MQuadrangle*)(*it));
- }
- tmp_model->add(gf);
-
- PartitionMesh(tmp_model,options);
-
- tmp_model->remove(gf);
- delete tmp_model;
-
- return 1;
-
+int PartitionMeshElements(std::vector<MElement*> &elements, meshPartitionOptions &options)
+{
+ GModel *tmp_model = new GModel();
+ GFace *gf = new discreteFace(tmp_model, 1);
+ std::set<MVertex *> setv;
+ for (unsigned i=0;i<elements.size();++i)
+ for (int j=0;j<elements[i]->getNumVertices();j++)
+ setv.insert(elements[i]->getVertex(j));
+
+ for (std::set<MVertex* >::iterator it = setv.begin(); it != setv.end(); it++)
+ gf->mesh_vertices.push_back(*it);
+
+ for (std::vector<MElement* >::iterator it = elements.begin(); it != elements.end(); it++){
+ if ((*it)->getType() == TYPE_TRI)
+ gf->triangles.push_back((MTriangle*)(*it));
+ else if ((*it)->getType() == TYPE_QUA)
+ gf->quadrangles.push_back((MQuadrangle*)(*it));
+ }
+ tmp_model->add(gf);
+
+ PartitionMesh(tmp_model,options);
+
+ tmp_model->remove(gf);
+ delete tmp_model;
+
+ return 1;
}
-int PartitionMeshFace( std::list<GFace*> &cFaces, meshPartitionOptions &options)
+int PartitionMeshFace(std::list<GFace*> &cFaces, meshPartitionOptions &options)
{
GModel *tmp_model = new GModel();
for(std::list<GFace*>::iterator it = cFaces.begin(); it != cFaces.end(); it++)
@@ -187,7 +194,7 @@ int PartitionMesh(GModel *const model, meshPartitionOptions &options)
model->recomputeMeshPartitions();
if (options.createPartitionBoundaries)
- CreatePartitionBoundaries (model);
+ CreatePartitionBoundaries (model, true);
Msg::Info("Partitioning complete");
Msg::StatusBar(1, false, "Mesh");
@@ -725,7 +732,7 @@ struct MatchBoElemToGrVertex
template <class ITERATOR>
-void fillit_(std::multimap<MFace,MElement*,Less_Face> &faceToElement,
+void fillit_(std::multimap<MFace, MElement*, Less_Face> &faceToElement,
ITERATOR it_beg, ITERATOR it_end)
{
for (ITERATOR IT = it_beg; IT != it_end ; ++IT){
@@ -733,12 +740,12 @@ void fillit_(std::multimap<MFace,MElement*,Less_Face> &faceToElement,
for(int j=0;j < el->getNumFaces();j++) {
MFace e = el->getFace(j);
faceToElement.insert(std::make_pair(e,el));
- }// for every edge of the elem
- }// for every elem of the face
+ }
+ }
}
template <class ITERATOR>
-void fillit_(std::multimap<MEdge,MElement*,Less_Edge> &edgeToElement,
+void fillit_(std::multimap<MEdge, MElement*, Less_Edge> &edgeToElement,
ITERATOR it_beg, ITERATOR it_end)
{
for (ITERATOR IT = it_beg; IT != it_end ; ++IT){
@@ -746,12 +753,12 @@ void fillit_(std::multimap<MEdge,MElement*,Less_Edge> &edgeToElement,
for(int j=0;j < el->getNumEdges();j++) {
MEdge e = el->getEdge(j);
edgeToElement.insert(std::make_pair(e,el));
- }// for every edge of the elem
- }// for every elem of the face
+ }
+ }
}
template <class ITERATOR>
-void fillit_(std::multimap<MVertex*,MElement*> &vertexToElement,
+void fillit_(std::multimap<MVertex*, MElement*> &vertexToElement,
ITERATOR it_beg, ITERATOR it_end)
{
for (ITERATOR IT = it_beg; IT != it_end ; ++IT){
@@ -759,8 +766,8 @@ void fillit_(std::multimap<MVertex*,MElement*> &vertexToElement,
for(int j=0;j < el->getNumVertices();j++) {
MVertex* e = el->getVertex(j);
vertexToElement.insert(std::make_pair(e,el));
- }// for every edge of the elem
- }// for every elem of the face
+ }
+ }
}
void assignPartitionBoundary(GModel *model,
@@ -863,7 +870,7 @@ void splitBoundaryEdges(GModel *model, std::set<partitionEdge*, Less_partition
MVertex *v1 = (*it)->getVertex(0);
MVertex *v2 = (*it)->getVertex(1);
std::list<MLine*>::iterator itp;
- if ( v1 == vE ){
+ if (v1 == vE){
myLines.push_back(*it);
itp = it;
it++;
@@ -871,7 +878,7 @@ void splitBoundaryEdges(GModel *model, std::set<partitionEdge*, Less_partition
vE = v2;
i = -1;
}
- else if ( v2 == vE){
+ else if (v2 == vE){
myLines.push_back(*it);
itp = it;
it++;
@@ -948,110 +955,159 @@ void assignPartitionBoundary(GModel *model,
ppv->points.push_back(new MPoint (ve));
}
-int CreatePartitionBoundaries(GModel *model)
+static void addGhostCells(GEntity *ge,
+ std::multimap<MVertex*, MElement*> &vertexToElement,
+ std::multimap<MElement*, short> &ghosts)
+{
+ // get all the nodes on the partition boundary (there are not stored
+ // so we need to recompute them)
+ std::set<MVertex*> verts;
+ for(unsigned int i = 0; i < ge->getNumMeshElements(); i++){
+ MElement *e = ge->getMeshElement(i);
+ for(unsigned int j = 0; j < e->getNumVertices(); j++)
+ verts.insert(e->getVertex(j));
+ }
+
+ // get all the elements in touching the interface nodes
+ for(std::set<MVertex*>::iterator it = verts.begin(); it != verts.end(); it++){
+ MVertex *v = *it;
+ std::pair<std::multimap<MVertex*, MElement*>::iterator,
+ std::multimap<MVertex*, MElement*>::iterator> itp =
+ vertexToElement.equal_range(v);
+ // get all partitions sharing this node
+ std::set<short> parts;
+ for(std::multimap<MVertex*, MElement*>::iterator ite = itp.first;
+ ite != itp.second; ite++)
+ parts.insert(ite->second->getPartition());
+ // update partition info for each element touching the node
+ for(std::multimap<MVertex*, MElement*>::iterator ite = itp.first;
+ ite != itp.second; ite++){
+ MElement *e = ite->second;
+ for(std::set<short>::iterator its = parts.begin(); its != parts.end(); its++)
+ if(*its != e->getPartition())
+ ghosts.insert(std::pair<MElement*, short>(e, *its));
+ }
+ }
+#if 1
+ for(std::multimap<MElement*, short>::iterator it = ghosts.begin();
+ it != ghosts.end(); it++)
+ printf("ele %d (part %d) ghost %d\n", it->first->getNum(),
+ it->first->getPartition(), it->second);
+#endif
+}
+
+int CreatePartitionBoundaries(GModel *model, bool createGhostCells)
{
unsigned numElem[5];
const int meshDim = model->getNumMeshElements(numElem);
+ std::set<partitionFace*, Less_partitionFace> pfaces;
std::set<partitionEdge*, Less_partitionEdge> pedges;
std::set<partitionVertex*, Less_partitionVertex> pvertices;
- std::set<partitionFace*, Less_partitionFace> pfaces;
+ std::multimap<MFace, MElement*, Less_Face> faceToElement;
+ std::multimap<MEdge, MElement*, Less_Edge> edgeToElement;
+ std::multimap<MVertex*, MElement*> vertexToElement;
+
// assign partition faces
if (meshDim == 3){
- std::multimap<MFace,MElement*,Less_Face> faceToElement;
for(GModel::riter it = model->firstRegion(); it != model->lastRegion(); ++it){
- fillit_ ( faceToElement,(*it)->tetrahedra.begin(),(*it)->tetrahedra.end());
- fillit_ ( faceToElement,(*it)->hexahedra.begin(),(*it)->hexahedra.end());
- fillit_ ( faceToElement,(*it)->prisms.begin(),(*it)->prisms.end());
- fillit_ ( faceToElement,(*it)->pyramids.begin(),(*it)->pyramids.end());
- fillit_ ( faceToElement,(*it)->polyhedra.begin(),(*it)->polyhedra.end());
+ fillit_(faceToElement, (*it)->tetrahedra.begin(), (*it)->tetrahedra.end());
+ fillit_(faceToElement, (*it)->hexahedra.begin(), (*it)->hexahedra.end());
+ fillit_(faceToElement, (*it)->prisms.begin(), (*it)->prisms.end());
+ fillit_(faceToElement, (*it)->pyramids.begin(), (*it)->pyramids.end());
+ fillit_(faceToElement, (*it)->polyhedra.begin(), (*it)->polyhedra.end());
}
- {
- std::multimap<MFace,MElement*,Less_Face>::iterator it = faceToElement.begin();
- Equal_Face oper;
- while ( it != faceToElement.end() ){
- MFace e = it->first;
- std::vector<MElement*> voe;
- do {
- voe.push_back(it->second);
- ++it;
- if (it == faceToElement.end() ) break;
- }while (oper (e,it->first));
- assignPartitionBoundary (model,e,pfaces,voe);
- }
+ std::multimap<MFace, MElement*, Less_Face>::iterator it = faceToElement.begin();
+ Equal_Face oper;
+ while (it != faceToElement.end()){
+ MFace e = it->first;
+ std::vector<MElement*> voe;
+ do {
+ voe.push_back(it->second);
+ ++it;
+ if (it == faceToElement.end()) break;
+ } while (oper (e,it->first));
+ assignPartitionBoundary(model, e, pfaces, voe);
}
}
// assign partition edges
if (meshDim > 1){
- std::multimap<MEdge,MElement*,Less_Edge> edgeToElement;
if (meshDim == 2){
for(GModel::fiter it = model->firstFace(); it != model->lastFace(); ++it){
- fillit_ ( edgeToElement,(*it)->triangles.begin(),(*it)->triangles.end());
- fillit_ ( edgeToElement,(*it)->quadrangles.begin(),(*it)->quadrangles.end());
- fillit_ ( edgeToElement,(*it)->polygons.begin(),(*it)->polygons.end());
+ fillit_(edgeToElement, (*it)->triangles.begin(), (*it)->triangles.end());
+ fillit_(edgeToElement, (*it)->quadrangles.begin(), (*it)->quadrangles.end());
+ fillit_(edgeToElement, (*it)->polygons.begin(), (*it)->polygons.end());
}
}
else if (meshDim == 3){
for(GModel::riter it = model->firstRegion(); it != model->lastRegion(); ++it){
- fillit_ ( edgeToElement,(*it)->tetrahedra.begin(),(*it)->tetrahedra.end());
- fillit_ ( edgeToElement,(*it)->hexahedra.begin(),(*it)->hexahedra.end());
- fillit_ ( edgeToElement,(*it)->prisms.begin(),(*it)->prisms.end());
- fillit_ ( edgeToElement,(*it)->pyramids.begin(),(*it)->pyramids.end());
- fillit_ ( edgeToElement,(*it)->polyhedra.begin(),(*it)->polyhedra.end());
+ fillit_(edgeToElement, (*it)->tetrahedra.begin(), (*it)->tetrahedra.end());
+ fillit_(edgeToElement, (*it)->hexahedra.begin(), (*it)->hexahedra.end());
+ fillit_(edgeToElement, (*it)->prisms.begin(), (*it)->prisms.end());
+ fillit_(edgeToElement, (*it)->pyramids.begin(), (*it)->pyramids.end());
+ fillit_(edgeToElement, (*it)->polyhedra.begin(), (*it)->polyhedra.end());
}
}
- {
- std::multimap<MEdge,MElement*,Less_Edge>::iterator it = edgeToElement.begin();
- Equal_Edge oper;
- while ( it != edgeToElement.end() ){
- MEdge e = it->first;
- std::vector<MElement*> voe;
- do {
- voe.push_back(it->second);
- ++it;
- if (it == edgeToElement.end() ) break;
- }while (oper (e,it->first));
- assignPartitionBoundary (model,e,pedges,voe,pfaces);
- }
- //splitBoundaryEdges(model,pedges);
+ std::multimap<MEdge, MElement*, Less_Edge>::iterator it = edgeToElement.begin();
+ Equal_Edge oper;
+ while (it != edgeToElement.end()){
+ MEdge e = it->first;
+ std::vector<MElement*> voe;
+ do {
+ voe.push_back(it->second);
+ ++it;
+ if (it == edgeToElement.end()) break;
+ }while (oper (e,it->first));
+ assignPartitionBoundary(model, e, pedges, voe, pfaces);
}
+ //splitBoundaryEdges(model,pedges);
}
// make partition vertices
if (meshDim > 1){
- std::multimap<MVertex*,MElement*> vertexToElement;
if (meshDim == 2){
for(GModel::fiter it = model->firstFace(); it != model->lastFace(); ++it){
- fillit_ ( vertexToElement,(*it)->triangles.begin(),(*it)->triangles.end());
- fillit_ ( vertexToElement,(*it)->quadrangles.begin(),(*it)->quadrangles.end());
- fillit_ ( vertexToElement,(*it)->polygons.begin(),(*it)->polygons.end());
+ fillit_(vertexToElement, (*it)->triangles.begin(), (*it)->triangles.end());
+ fillit_(vertexToElement, (*it)->quadrangles.begin(), (*it)->quadrangles.end());
+ fillit_(vertexToElement, (*it)->polygons.begin(), (*it)->polygons.end());
}
}
else if (meshDim == 3){
for(GModel::riter it = model->firstRegion(); it != model->lastRegion(); ++it){
- fillit_ ( vertexToElement,(*it)->tetrahedra.begin(),(*it)->tetrahedra.end());
- fillit_ ( vertexToElement,(*it)->hexahedra.begin(),(*it)->hexahedra.end());
- fillit_ ( vertexToElement,(*it)->prisms.begin(),(*it)->prisms.end());
- fillit_ ( vertexToElement,(*it)->pyramids.begin(),(*it)->pyramids.end());
- fillit_ ( vertexToElement,(*it)->polyhedra.begin(),(*it)->polyhedra.end());
+ fillit_(vertexToElement, (*it)->tetrahedra.begin(), (*it)->tetrahedra.end());
+ fillit_(vertexToElement, (*it)->hexahedra.begin(), (*it)->hexahedra.end());
+ fillit_(vertexToElement, (*it)->prisms.begin(), (*it)->prisms.end());
+ fillit_(vertexToElement, (*it)->pyramids.begin(), (*it)->pyramids.end());
+ fillit_(vertexToElement, (*it)->polyhedra.begin(), (*it)->polyhedra.end());
}
}
- {
- std::multimap<MVertex*,MElement*>::iterator it = vertexToElement.begin();
- while ( it != vertexToElement.end() ){
- MVertex *v = it->first;
- std::vector<MElement*> voe;
- do {
- voe.push_back(it->second);
- ++it;
- if (it == vertexToElement.end() ) break;
- }while (v == it->first);
- assignPartitionBoundary (model,v,pvertices,voe,pedges,pfaces);
- }
+ std::multimap<MVertex*, MElement*>::iterator it = vertexToElement.begin();
+ while (it != vertexToElement.end()){
+ MVertex *v = it->first;
+ std::vector<MElement*> voe;
+ do {
+ voe.push_back(it->second);
+ ++it;
+ if (it == vertexToElement.end()) break;
+ }while (v == it->first);
+ assignPartitionBoundary(model, v, pvertices, voe, pedges, pfaces);
}
}
+ if(createGhostCells){
+ std::multimap<MElement*, short> &ghosts(model->getGhostCells());
+ ghosts.clear();
+ if(meshDim == 2)
+ for(std::set<partitionEdge*, Less_partitionEdge>::iterator it = pedges.begin();
+ it != pedges.end(); it++)
+ addGhostCells(*it, vertexToElement, ghosts);
+ else if(meshDim == 3)
+ for(std::set<partitionFace*, Less_partitionFace>::iterator it = pfaces.begin();
+ it != pfaces.end(); it++)
+ addGhostCells(*it, vertexToElement, ghosts);
+ }
+
return 1;
}
@@ -1063,7 +1119,7 @@ void createPartitionFaces(GModel *model, GFaceCompound *gf, int N,
//--------------------------------------------
std::vector<std::set<MVertex*> > allNodes;
int numMax = model->maxFaceNum() + 1;
- for( int i =0; i < N; i++){
+ for(int i = 0; i < N; i++){
//printf("*** Created discreteFace %d \n", numMax+i);
discreteFace *face = new discreteFace(model, numMax+i);
discreteFaces.push_back(face);
@@ -1075,7 +1131,7 @@ void createPartitionFaces(GModel *model, GFaceCompound *gf, int N,
std::list<GFace*> _compound = gf->getCompounds();
std::list<GFace*>::iterator it = _compound.begin();
- for( ; it != _compound.end() ; ++it){
+ for( ; it != _compound.end(); ++it){
for(unsigned int i = 0; i < (*it)->triangles.size(); ++i){
MTriangle *e = (*it)->triangles[i];
int part = e->getPartition();
@@ -1087,7 +1143,7 @@ void createPartitionFaces(GModel *model, GFaceCompound *gf, int N,
}
}
- for( int i =0; i < N; i++){
+ for(int i = 0; i < N; i++){
for (std::set<MVertex*>::iterator it = allNodes[i].begin(); it != allNodes[i].end(); it++){
discreteFaces[i]->mesh_vertices.push_back(*it);
}
diff --git a/Mesh/meshPartition.h b/Mesh/meshPartition.h
index 9e7780981a1dc0e13784c3ea064f27e2d0799c95..ac97f350f873233139c7702a1727b302a32ce2ad 100644
--- a/Mesh/meshPartition.h
+++ b/Mesh/meshPartition.h
@@ -3,22 +3,18 @@
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
-#ifndef _PARTITION_H_
-#define _PARTITION_H_
+#ifndef _MESH_PARTITION_H_
+#define _MESH_PARTITION_H_
#include <vector>
#include "partitionEdge.h"
#include "GFaceCompound.h"
#include "GFace.h"
-struct meshPartitionOptions;
-struct BoElemGr;
class GModel;
class GFace;
class Graph;
-
-typedef std::vector<BoElemGr> BoElemGrVec;
-typedef enum {LAPLACIAN= 0, MULTILEVEL=1} typeOfPartition;
+struct meshPartitionOptions;
/*******************************************************************************
*
@@ -26,14 +22,12 @@ typedef enum {LAPLACIAN= 0, MULTILEVEL=1} typeOfPartition;
*
******************************************************************************/
-int MakeGraph(GModel *const model, Graph &graph,
- BoElemGrVec *const boElemGrVec = 0);
int PartitionGraph(Graph &graph, meshPartitionOptions &options);
int PartitionMesh(GModel *const model, meshPartitionOptions &options);
int PartitionMeshFace(std::list<GFace*> &cFaces, meshPartitionOptions &options);
-int PartitionMeshElements( std::vector<MElement*> &elements, meshPartitionOptions &options);
+int PartitionMeshElements(std::vector<MElement*> &elements, meshPartitionOptions &options);
bool PartitionZeroGenus(std::list<GFace*> &cFaces, int &nbParts);
-int CreatePartitionBoundaries (GModel *model);
+int CreatePartitionBoundaries(GModel *model, bool createGhostCells);
void splitBoundaryEdges(GModel *model, std::set<partitionEdge*, Less_partitionEdge> &newEdges);
void createPartitionFaces(GModel *model, GFaceCompound * gf, int num_parts,
diff --git a/Mesh/multiscalePartition.cpp b/Mesh/multiscalePartition.cpp
index 7e7887400517339924449a395fb465dc9fb66e54..902c293bbc48a2d979a022e7bdcccc847f539a8d 100644
--- a/Mesh/multiscalePartition.cpp
+++ b/Mesh/multiscalePartition.cpp
@@ -7,11 +7,11 @@
#include "multiscaleLaplace.h"
#include "Context.h"
-static void recur_connect (MVertex *v,
- std::multimap<MVertex*,MEdge> &v2e,
- std::set<MEdge,Less_Edge> &group,
- std::set<MVertex*> &touched){
-
+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);
@@ -27,9 +27,8 @@ static void recur_connect (MVertex *v,
static int connected_bounds (std::vector<MEdge> &edges, std::vector<std::vector<MEdge> > &boundaries)
{
-
std::multimap<MVertex*,MEdge> v2e;
- for (unsigned i=0;i<edges.size();++i){
+ for (unsigned i = 0; i < edges.size(); ++i){
for (int j=0;j<edges[i].getNumVertices();j++){
v2e.insert(std::make_pair(edges[i].getVertex(j),edges[i]));
}
@@ -41,7 +40,7 @@ static int connected_bounds (std::vector<MEdge> &edges, std::vector<std::vector
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)
+ for (std::set<MVertex*>::iterator it = touched.begin() ; it != touched.end();++it)
v2e.erase(*it);
}
@@ -50,9 +49,10 @@ static int connected_bounds (std::vector<MEdge> &edges, std::vector<std::vector
static int getGenus (std::vector<MElement *> &elements,
- std::vector<std::vector<MEdge> > &boundaries){
+ std::vector<std::vector<MEdge> > &boundaries)
+{
- //We suppose MElements are simply connected
+ //We suppose MElements are simply connected
std::set<MEdge, Less_Edge> es;
std::set<MVertex*> vs;
@@ -91,10 +91,9 @@ static int getGenus (std::vector<MElement *> &elements,
}
-static int getAspectRatio (std::vector<MElement *> &elements,
- std::vector<std::vector<MEdge> > &boundaries)
+static int getAspectRatio(std::vector<MElement *> &elements,
+ std::vector<std::vector<MEdge> > &boundaries)
{
-
std::set<MVertex*> vs;
for(unsigned int i = 0; i < elements.size(); i++){
MElement *e = elements[i];
@@ -108,61 +107,56 @@ static int getAspectRatio (std::vector<MElement *> &elements,
for (std::set<MVertex* >::iterator it = vs.begin(); it != vs.end(); it++){
SPoint3 pt((*it)->x(),(*it)->y(), (*it)->z());
vertices.push_back(pt);
- bb +=pt;
+ bb += pt;
}
-
+
//obbox = SOrientedBoundingBox::buildOBB(vertices);
//double H = obbox.getMaxSize();
double H = norm(SVector3(bb.max(), bb.min()));
-
- double D = H;
- for (unsigned i=0; i< boundaries.size(); i++){
- std::set<MVertex*> vb;
- std::vector<MEdge> iBound = boundaries[i];
- for (unsigned j=0; j< iBound.size(); j++){
- MEdge e = iBound[j];
- vb.insert(e.getVertex(0));
- vb.insert(e.getVertex(1));
- }
- std::vector<SPoint3> vBounds;
- for (std::set<MVertex* >::iterator it = vb.begin(); it != vb.end(); it++){
- SPoint3 pt((*it)->x(),(*it)->y(), (*it)->z());
- vBounds.push_back(pt);
- }
- SOrientedBoundingBox obboxD = SOrientedBoundingBox::buildOBB(vBounds);
- D = std::min(D, obboxD.getMaxSize());
- }
- int AR = (int)ceil(H/D);
-
- return AR;
-
+
+ double D = H;
+ for (unsigned int i = 0; i < boundaries.size(); i++){
+ std::set<MVertex*> vb;
+ std::vector<MEdge> iBound = boundaries[i];
+ for (unsigned int j = 0; j < iBound.size(); j++){
+ MEdge e = iBound[j];
+ vb.insert(e.getVertex(0));
+ vb.insert(e.getVertex(1));
+ }
+ std::vector<SPoint3> vBounds;
+ for (std::set<MVertex* >::iterator it = vb.begin(); it != vb.end(); it++){
+ SPoint3 pt((*it)->x(),(*it)->y(), (*it)->z());
+ vBounds.push_back(pt);
+ }
+ SOrientedBoundingBox obboxD = SOrientedBoundingBox::buildOBB(vBounds);
+ D = std::min(D, obboxD.getMaxSize());
+ }
+ int AR = (int)ceil(H/D);
+
+ return AR;
}
-static void getGenusAndRatio(std::vector<MElement *> &elements, int & genus, int &AR){
+static void getGenusAndRatio(std::vector<MElement *> &elements, int & genus, int &AR)
+{
std::vector<std::vector<MEdge> > boundaries;
genus = getGenus(elements, boundaries);
AR = getAspectRatio(elements, boundaries);
-
- return;
-
}
-static void partitionRegions (std::vector<MElement*> &elements,
- std::vector<std::vector<MElement*> > ®ions){
-
- for (unsigned i=0;i<elements.size();++i){
+
+static void partitionRegions(std::vector<MElement*> &elements,
+ std::vector<std::vector<MElement*> > ®ions)
+{
+ for (unsigned int i = 0; i < elements.size(); ++i){
MElement *e = elements[i];
int part = e->getPartition();
regions[part-1].push_back(e);
}
-
- return;
-
}
-static void printLevel ( std::vector<MElement *> &elements, int recur, int region){
-
+static void printLevel(std::vector<MElement *> &elements, int recur, int region)
+{
char fn[256];
- sprintf(fn, "part_%d_%d.msh", recur, region );
+ sprintf(fn, "part_%d_%d.msh", recur, region);
double version = 2.0;
std::set<MVertex*> vs;
@@ -178,7 +172,7 @@ static void printLevel ( std::vector<MElement *> &elements, int recur, int regio
fprintf(fp, "%g %d %d\n", version, binary ? 1 : 0, (int)sizeof(double));
fprintf(fp, "$EndMeshFormat\n");
- fprintf(fp,"$Nodes\n%d\n",vs.size());
+ fprintf(fp,"$Nodes\n%d\n", (int)vs.size());
std::set<MVertex*> :: iterator it = vs.begin();
int index = 1;
for (; it != vs.end() ; ++it){
@@ -188,26 +182,25 @@ static void printLevel ( std::vector<MElement *> &elements, int recur, int regio
}
fprintf(fp,"$EndNodes\n",elements.size());
- fprintf(fp,"$Elements\n%d\n",elements.size());
+ fprintf(fp,"$Elements\n%d\n", (int)elements.size());
for (int i=0;i<elements.size();i++){
elements[i]->writeMSH(fp,version);
}
- fprintf(fp,"$EndElements\n%d\n",elements.size());
+ fprintf(fp,"$EndElements\n%d\n", (int)elements.size());
fclose(fp);
}
-multiscalePartition::multiscalePartition (std::vector<MElement *> &elements, int nbParts, typeOfPartition method)
+multiscalePartition::multiscalePartition(std::vector<MElement *> &elements, int nbParts, typeOfPartition method)
{
+ options = CTX::instance()->partitionOptions;
+ options.num_partitions = nbParts;
+ options.partitioner = 1; //1 CHACO, 2 METIS
+ if (options.partitioner == 1){
+ options.global_method = 2;// 1 Multilevel-KL, 2 Spectral
+ options.mesh_dims[0] = nbParts;
+ }
- options = CTX::instance()->partitionOptions;
- options.num_partitions = nbParts;
- options.partitioner = 1; //1 CHACO, 2 METIS
- if ( options.partitioner == 1){
- options.global_method = 2;// 1 Multilevel-KL, 2 Spectral
- options.mesh_dims[0] = nbParts;
- }
-
partitionLevel *level = new partitionLevel;
level->elements.insert(level->elements.begin(),elements.begin(),elements.end());
level->recur = 0;
@@ -218,18 +211,18 @@ multiscalePartition::multiscalePartition (std::vector<MElement *> &elements, int
partition(*level, nbParts, method);
totalParts = assembleAllPartitions();
-
}
-void multiscalePartition:: setNumberOfPartitions(int &nbParts)
+
+void multiscalePartition::setNumberOfPartitions(int &nbParts)
{
options.num_partitions = nbParts;
- if ( options.partitioner == 1){
+ if (options.partitioner == 1){
options.mesh_dims[0] = nbParts;
}
}
-void multiscalePartition::partition(partitionLevel & level, int nbParts, typeOfPartition method){
-
+void multiscalePartition::partition(partitionLevel & level, int nbParts, typeOfPartition method)
+{
#if defined(HAVE_METIS) || defined(HAVE_CHACO)
if (method == LAPLACIAN){
@@ -239,7 +232,7 @@ void multiscalePartition::partition(partitionLevel & level, int nbParts, typeOfP
setNumberOfPartitions(nbParts);
PartitionMeshElements(level.elements, options);
}
-
+
std::vector<std::vector<MElement*> > regions(nbParts);
partitionRegions(level.elements, regions);
level.elements.clear();
@@ -278,14 +271,14 @@ void multiscalePartition::partition(partitionLevel & level, int nbParts, typeOfP
Msg::Info("Multiscale part: level %d, region %d is ZERO-GENUS (AR=%d)",
nextLevel->recur,nextLevel->region, AR);
}
-
-}
-
+
+ }
+
#endif
}
-int multiscalePartition::assembleAllPartitions(){
-
+int multiscalePartition::assembleAllPartitions()
+{
int nbParts = 1;
for (unsigned i = 0; i< levels.size(); i++){
@@ -299,9 +292,6 @@ int multiscalePartition::assembleAllPartitions(){
nbParts++;
}
}
-
- return nbParts-1;
-
-
-
+
+ return nbParts - 1;
}
diff --git a/Mesh/multiscalePartition.h b/Mesh/multiscalePartition.h
index 3c071685426d444f2fd3eec523546b7e47a70b25..e21b8fb8d72d3678a5c9f5a4e8ef47e7ffcbf787 100644
--- a/Mesh/multiscalePartition.h
+++ b/Mesh/multiscalePartition.h
@@ -21,6 +21,8 @@ struct partitionLevel {
std::vector<MElement *> elements;
};
+typedef enum {LAPLACIAN= 0, MULTILEVEL=1} typeOfPartition;
+
class multiscalePartition{
private: