diff --git a/Geo/MElement.h b/Geo/MElement.h index 0dd0b769be2887c184f1979aff11c2a58fdd99b0..38ba9eb16c1373ff3d4dbcabf350f2d4342d9894 100644 --- a/Geo/MElement.h +++ b/Geo/MElement.h @@ -192,6 +192,11 @@ class MElement // otherwise get the minimum radius of all the circles/spheres // tangent to the most boundaries of the element. virtual double getInnerRadius(){ return 0.; } + + // get the radius of the circumscribed circle/sphere if it exists, + // otherwise get the maximum radius of all the circles/spheres + // tangent to the most boundaries of the element. + virtual double getOuterRadius(){ return 0.; } // compute the barycenter virtual SPoint3 barycenter(); diff --git a/Geo/MQuadrangle.cpp b/Geo/MQuadrangle.cpp index bc9aa6f41b616dff4e4d3daed8f5d5d6c0f79c19..ff011fea99ac9e23b123df005424a3289d1c2b14 100644 --- a/Geo/MQuadrangle.cpp +++ b/Geo/MQuadrangle.cpp @@ -267,7 +267,11 @@ double MQuadrangle::angleShapeMeasure() return 1.; #endif } - +double MQuadrangle::getOuterRadius() +{ + // TO DO!!!!!!!!!!!!! (BRUNO SENY) + return 1.0; +} double MQuadrangle::getInnerRadius() { // get the coordinates (x, y, z) of the 4 points defining the Quad diff --git a/Geo/MQuadrangle.h b/Geo/MQuadrangle.h index 412224a8ff2043f9c0e3e21d7e7b26db358c9a41..8de5e45fcb021fc86e391da9581124751dafa48a 100644 --- a/Geo/MQuadrangle.h +++ b/Geo/MQuadrangle.h @@ -39,6 +39,7 @@ class MQuadrangle : public MElement { v[2] = _v[2]; v[3] = _v[3]; } +void projectInMeanPlane(double *xn, double *yn); public : MQuadrangle(MVertex *v0, MVertex *v1, MVertex *v2, MVertex *v3, int num=0, int part=0) : MElement(num, part) @@ -141,6 +142,7 @@ class MQuadrangle : public MElement { // planar, we compute the mean plane due to the least-square // criterion. virtual double getInnerRadius(); + virtual double getOuterRadius(); private: int edges_quad(const int edge, const int vert) const { diff --git a/Geo/MTriangle.cpp b/Geo/MTriangle.cpp index 92d0aab1fa3fd62b27f9f6b0f9c2880528226e47..21dc2f496be89131af49ad9b584716e7f8433ed3 100644 --- a/Geo/MTriangle.cpp +++ b/Geo/MTriangle.cpp @@ -44,7 +44,21 @@ double MTriangle::getInnerRadius() dist[i] = e.getVertex(0)->distance(e.getVertex(1)); k += 0.5 * dist[i]; } - return sqrt(k * (k - dist[0]) * (k - dist[1]) * (k - dist[2])) / k; + double area = sqrt(k*(k-dist[0])*(k-dist[1])*(k-dist[2])); + return area/k; +} + +double MTriangle::getOuterRadius() +{ + // radius of circle circumscribing a triangle + double dist[3], k = 0.; + for (int i = 0; i < 3; i++){ + MEdge e = getEdge(i); + dist[i] = e.getVertex(0)->distance(e.getVertex(1)); + k += 0.5 * dist[i]; + } + double area = sqrt(k*(k-dist[0])*(k-dist[1])*(k-dist[2])); + return dist[0]*dist[1]*dist[2]/(4*area); } double MTriangle::angleShapeMeasure() diff --git a/Geo/MTriangle.h b/Geo/MTriangle.h index 87150f25bc09891f0c6e82bd85f38c37cc9685d5..fdd769b25e054e316c6359270125a4b272acb4f4 100644 --- a/Geo/MTriangle.h +++ b/Geo/MTriangle.h @@ -54,6 +54,7 @@ class MTriangle : public MElement { virtual double gammaShapeMeasure(); virtual double distoShapeMeasure(); virtual double getInnerRadius(); + virtual double getOuterRadius(); virtual double angleShapeMeasure(); virtual int getNumVertices() const { return 3; } virtual MVertex *getVertex(int num){ return _v[num]; } diff --git a/Mesh/meshPartition.cpp b/Mesh/meshPartition.cpp index 7b1283f9132c638e930b751563c861c4ddc003ee..95ee43016aa382b807ed8193a0de35a2be81e58b 100644 --- a/Mesh/meshPartition.cpp +++ b/Mesh/meshPartition.cpp @@ -59,7 +59,7 @@ extern "C" void METIS_mCPartGraphRecursive (int *, int *, idxtype *, idxtype *, idxtype *, idxtype *, int *, int *, int *, int *, int *, idxtype *); extern "C" void METIS_mCPartGraphKway -(int *, int *, idxtype *, idxtype *, idxtype *, idxtype *, int *, int *, int *, int *, +(int *, int *, idxtype *, idxtype *, idxtype *, idxtype *, int *, int *, int *, float *, int *, int *, idxtype *); @@ -437,6 +437,7 @@ int PartitionGraph(Graph &graph, meshPartitionOptions &options) int numflag = 0; // if metisOptions[0]=0 then default options int metisOptions[5]; + float ubvec[options.ncon]; int edgeCut; const int iSec = 0; switch(options.algorithm) { @@ -507,11 +508,21 @@ int PartitionGraph(Graph &graph, meshPartitionOptions &options) metisOptions[2] = 1; metisOptions[3] = options.refine_algorithm; metisOptions[4] = 0; + printf("Tolerance for Constraints:["); + for(int u=0;u<options.ncon;u++){ + ubvec[u]=1.0; + if(options.tolerance[u]%options.num_partitions>0){ + ubvec[u] = (float) ceil((float)options.tolerance[u]/options.num_partitions)/((float)options.tolerance[u]/options.num_partitions); + } + printf(" %f", ubvec[u]); + } + printf("] \n"); + graph.fillWithMultipleWeights(options.ncon,options.getWeightMap()); if (options.num_partitions > 1) { METIS_mCPartGraphKway (&n,&options.ncon,&graph.xadj[graph.section[iSec]], &graph.adjncy[graph.section[iSec]], &graph.vwgts[graph.section[iSec]], NULL, &wgtflag, &numflag, - &options.num_partitions, metisOptions, &edgeCut, + &options.num_partitions,ubvec, metisOptions, &edgeCut, &graph.partition[graph.section[iSec]]); } break; diff --git a/Mesh/meshPartitionOptions.h b/Mesh/meshPartitionOptions.h index b90bdd9ecc2221b8093b6c06ac9ab8aba56ef0ae..549cb0af45bc26efe71e2250e35e8e41c913310c 100644 --- a/Mesh/meshPartitionOptions.h +++ b/Mesh/meshPartitionOptions.h @@ -96,6 +96,8 @@ class meshPartitionOptions std::vector<int> nodalWeights; std::map<int, std::vector<int> > weightMap; + + std::vector<int> tolerance; public: