diff --git a/Mesh/meshGEdge.cpp b/Mesh/meshGEdge.cpp index 93e07a95fffd930f306dbe026aac3e04ef112e7c..297ea37e93d76ce8d1430c2170148397ac25bc82 100644 --- a/Mesh/meshGEdge.cpp +++ b/Mesh/meshGEdge.cpp @@ -1,4 +1,4 @@ -// $Id: meshGEdge.cpp,v 1.40 2007-09-04 13:47:02 remacle Exp $ +// $Id: meshGEdge.cpp,v 1.41 2007-09-12 14:28:29 geuzaine Exp $ // // Copyright (C) 1997-2007 C. Geuzaine, J.-F. Remacle // @@ -32,91 +32,82 @@ typedef struct{ double t, lc, p; }IntPoint; - struct xi2lc { - double xi,lc; - xi2lc (const double &_xi,const double _lc) - :xi(_xi),lc(_lc) - { } - bool operator < ( const xi2lc &other) + double xi, lc; + xi2lc(const double &_xi,const double _lc) + : xi(_xi), lc(_lc) + { + } + bool operator < (const xi2lc &other) { return xi < other.xi; } }; -std::vector<xi2lc> interpLc; +static std::vector<xi2lc> interpLc; -void smoothInterpLc (bool periodic, int nbSmooth) +void smoothInterpLc(bool periodic, int nbSmooth) { - if (periodic) - { - for (int i=0 ; i<interpLc.size()*nbSmooth; i++) - { - xi2lc &left = interpLc[(i-1)%interpLc.size()]; - xi2lc &mid = interpLc[i%interpLc.size()]; - xi2lc &right = interpLc[(i+1)%interpLc.size()]; - - if (1./mid.lc > 1.1 * 1./left.lc)mid.lc = left.lc/1.1; - if (1./mid.lc > 1.1 * 1./right.lc)mid.lc = right.lc/1.1; - } - } - else - { - for (int j=0 ; j<nbSmooth; j++) - { - for (int i=0 ; i<interpLc.size(); i++) - { - xi2lc &left = (i==0)?interpLc[0]:interpLc[i-1]; - xi2lc &mid = interpLc[i]; - xi2lc &right = (i==interpLc.size()-1)?interpLc[interpLc.size()-1]:interpLc[i+1]; - - if (1./mid.lc > 1.1 * 1./left.lc)mid.lc = left.lc/1.1; - if (1./mid.lc > 1.1 * 1./right.lc)mid.lc = right.lc/1.1; - } - } + if(periodic){ + for(int i = 0; i < interpLc.size() * nbSmooth; i++){ + xi2lc &left = interpLc[(i - 1) % interpLc.size()]; + xi2lc &mid = interpLc[i % interpLc.size()]; + xi2lc &right = interpLc[(i + 1) % interpLc.size()]; + if(1. / mid.lc > 1.1 * 1. / left.lc) mid.lc = left.lc/ 1.1; + if(1. / mid.lc > 1.1 * 1. / right.lc) mid.lc = right.lc/ 1.1; } + } + else{ + for(int j = 0; j < nbSmooth; j++){ + for(int i = 0 ; i < interpLc.size(); i++){ + xi2lc &left = (i == 0) ? interpLc[0] : interpLc[i - 1]; + xi2lc &mid = interpLc[i]; + xi2lc &right = (i == interpLc.size() - 1) ? + interpLc[interpLc.size() - 1] : interpLc[i+1]; + if(1. / mid.lc > 1.1 * 1. / left.lc) mid.lc = left.lc / 1.1; + if(1. / mid.lc > 1.1 * 1. / right.lc) mid.lc = right.lc / 1.1; + } + } + } } -void printInterpLc (const char *name) +void printInterpLc(const char *name) { - FILE *f = fopen (name,"w"); - for (int i=0 ; i<interpLc.size(); i++) - { - xi2lc &interp = interpLc[i]; - fprintf(f,"%12.5E %12.5E\n",interp.xi,1/interp.lc); - } + FILE *f = fopen(name,"w"); + for(int i = 0; i < interpLc.size(); i++){ + xi2lc &interp = interpLc[i]; + fprintf(f,"%12.5E %12.5E\n", interp.xi, 1 / interp.lc); + } fclose(f); } -void buildInterpLc (List_T *lcPoints) +void buildInterpLc(List_T *lcPoints) { IntPoint p; interpLc.clear(); - for (int i=0;i<List_Nbr(lcPoints);i++) - { - List_Read(lcPoints, i, &p); - interpLc.push_back(xi2lc ( p.t,p.lc)); - } - // printf("interpLc with %d points\n",interpLc.size()); + for(int i = 0; i < List_Nbr(lcPoints); i++){ + List_Read(lcPoints, i, &p); + interpLc.push_back(xi2lc( p.t, p.lc)); + } } double F_Lc_usingInterpLc(GEdge *ge, double t) { - std::vector<xi2lc>::iterator it = std::lower_bound (interpLc.begin(),interpLc.end(),xi2lc(t,0)); + std::vector<xi2lc>::iterator it = std::lower_bound(interpLc.begin(), + interpLc.end(), xi2lc(t, 0)); double t1 = it->xi; double l1 = it->lc; it++; SVector3 der = ge->firstDer(t); const double d = norm(der); - if (it == interpLc.end())return d*l1; + if(it == interpLc.end()) return d * l1; double t2 = it->xi; double l2 = it->lc; - double l = l1 + ((t-t1)/(t2-t1)) * (l2-l1); - return d*l; + double l = l1 + ((t - t1) / (t2 - t1)) * (l2 - l1); + return d * l; } - double F_Lc_usingInterpLcBis(GEdge *ge, double t) { GPoint p = ge->point(t); @@ -136,7 +127,6 @@ double F_Lc_usingInterpLcBis(GEdge *ge, double t) return 1 / lc_here; } - double F_Lc(GEdge *ge, double t) { GPoint p = ge->point(t); @@ -319,9 +309,6 @@ void meshGEdge::operator() (GEdge *ge) // first compute the length of the curve by integrating one double length = Integration(ge, t_begin, t_end, F_One, Points, 1.e-8); - - // printf("%d points for the length\n",List_Nbr(Points)); - ge->setLength(length); List_Reset(Points); @@ -334,29 +321,26 @@ void meshGEdge::operator() (GEdge *ge) N = ge->meshAttributes.nbPointsTransfinite; } else{ - if (CTX.mesh.lc_integration_precision > 1.e-8) - { - Integration(ge, t_begin, t_end, F_Lc_usingInterpLcBis, lcPoints, CTX.mesh.lc_integration_precision); - buildInterpLc (lcPoints); - printInterpLc ("toto1.dat"); - smoothInterpLc (ge->periodic(),20); - printInterpLc ("toto2.dat"); - a = Integration(ge, t_begin, t_end, F_Lc_usingInterpLc, Points, 1.e-8); - // printf("%d points for LC , %d points for the distribution interpLc %d\n",List_Nbr(lcPoints),List_Nbr(Points),interpLc.size()); - } - else - { - a = Integration(ge, t_begin, t_end, F_Lc, Points, 1.e-8); - } + if(CTX.mesh.lc_integration_precision > 1.e-8){ + Integration(ge, t_begin, t_end, F_Lc_usingInterpLcBis, lcPoints, + CTX.mesh.lc_integration_precision); + buildInterpLc(lcPoints); + printInterpLc("toto1.dat"); + smoothInterpLc(ge->periodic(), 20); + printInterpLc("toto2.dat"); + a = Integration(ge, t_begin, t_end, F_Lc_usingInterpLc, Points, 1.e-8); + } + else{ + a = Integration(ge, t_begin, t_end, F_Lc, Points, 1.e-8); + } N = std::max(ge->minimumMeshSegments() + 1, (int)(a + 1.)); } - const double b = a / (double)(N - 1); - - // if the curve is periodic and if the begin vertex is identical to the end vertex - // and if this vertex has only one model curve adjacent to it, then the vertex is - // not connecting any other curve. So, the mesh vertex and its associated geom vertex - // are not necessary at the same location + // if the curve is periodic and if the begin vertex is identical to + // the end vertex and if this vertex has only one model curve + // adjacent to it, then the vertex is not connecting any other + // curve. So, the mesh vertex and its associated geom vertex are not + // necessary at the same location GPoint beg_p, end_p; if(ge->getBeginVertex() == ge->getEndVertex() && ge->getBeginVertex()->edges().size() == 1){ @@ -369,44 +353,30 @@ void meshGEdge::operator() (GEdge *ge) end_p = GPoint(v1->x(), v1->y(), v1->z()); } - int count = 1, NUMP = 1, NUMP2 = 1; - IntPoint P1, P2; - // do not consider the first and the last vertex (those are not // classified on this mesh edge) - if(N > 1){ + const double b = a / (double)(N - 1); + int count = 1, NUMP = 1; + IntPoint P1, P2; ge->mesh_vertices.resize(N - 2); - GPoint last_p = beg_p; while(NUMP < N - 1) { List_Read(Points, count - 1, &P1); List_Read(Points, count, &P2); - const double d = (double)NUMP *b; + const double d = (double)NUMP * b; if((fabs(P2.p) >= fabs(d)) && (fabs(P1.p) < fabs(d))) { double dt = P2.t - P1.t; double dp = P2.p - P1.p; double t = P1.t + dt / dp * (d - P1.p); GPoint V = ge->point(t); - if(1 || ge->meshAttributes.Method == TRANSFINI){ - ge->mesh_vertices[NUMP2 - 1] = new MEdgeVertex(V.x(), V.y(), V.z(), ge, t); - NUMP2++; - } - else{ - double lc = BGM_MeshSize(ge, t, 0, V.x(), V.y(), V.z()); - if(V.distance(last_p) > 0.7 * lc && - !(NUMP == N - 2 && V.distance(end_p) < 0.7 * lc)){ - last_p = V; - ge->mesh_vertices[NUMP2 - 1] = new MEdgeVertex(V.x(), V.y(), V.z(), ge, t); - NUMP2++; - } - } + ge->mesh_vertices[NUMP - 1] = new MEdgeVertex(V.x(), V.y(), V.z(), ge, t); NUMP++; } else { count++; } } - ge->mesh_vertices.resize(NUMP2 - 1); + ge->mesh_vertices.resize(NUMP - 1); } List_Delete(Points); List_Delete(lcPoints); @@ -426,6 +396,4 @@ void meshGEdge::operator() (GEdge *ge) v0->y() = beg_p.y(); v0->z() = beg_p.z(); } - - } diff --git a/benchmarks/3d/core_coil.geo b/benchmarks/3d/core_coil.geo index c45163022fb1a555737cdb44141ffa568aa7b655..32bc06332eeb506883e9e92ae248d67bf76c49bd 100644 --- a/benchmarks/3d/core_coil.geo +++ b/benchmarks/3d/core_coil.geo @@ -1,5 +1,5 @@ -Mesh.CharacteristicLengthFactor = 0.7; // 32,000 tets +//Mesh.CharacteristicLengthFactor = 0.7; // 32,000 tets //Mesh.CharacteristicLengthFactor = 0.5; // 91,000 tets //Mesh.CharacteristicLengthFactor = 0.4; // 163,000 tets //Mesh.CharacteristicLengthFactor = 0.3; // 340,000 tets diff --git a/benchmarks/misc/faraway/coque.geo b/benchmarks/misc/faraway/coque.geo index fc0d07c5a191cd378e32cde1b858dd1381f7942c..3e7bb9f3e03c5550809baac167d3d9faa9410269 100644 --- a/benchmarks/misc/faraway/coque.geo +++ b/benchmarks/misc/faraway/coque.geo @@ -33,32 +33,32 @@ a=8.976-8.916436 ; b=(9.006-8.976) ; c = 0.139564 -a-b ; -Extrude Line {1308,{0,0,a}}{Recombine;Layers{{3},{7000},{1}};}; +Extrude Line {1308,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1309,{0,0,h}}{Recombine;Layers{{3,4,4},{7003,7004,7005},{a/h,(a+b)/h,1}};}; +Extrude Line {1309,{0,0,h}}{Recombine;Layers{{3,4,4},{a/h,(a+b)/h,1}};}; -Extrude Line {1310,{0,0,a}}{Recombine;Layers{{3},{7006},{1}};}; +Extrude Line {1310,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1311,{0,0,h}}{Recombine;Layers{{3,4,4},{7009,7010,7011},{a/h,(a+b)/h,1}};}; +Extrude Line {1311,{0,0,h}}{Recombine;Layers{{3,4,4},{a/h,(a+b)/h,1}};}; -Extrude Line {1312,{0,0,a}}{Recombine;Layers{{3},{7012},{1}};}; +Extrude Line {1312,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1313,{0,0,h}}{Recombine;Layers{{3,4,4},{7015,7016,7017},{a/h,(a+b)/h,1}};}; +Extrude Line {1313,{0,0,h}}{Recombine;Layers{{3,4,4},{a/h,(a+b)/h,1}};}; -Extrude Line {1314,{0,0,a}}{Recombine;Layers{{3},{7018},{1}};}; +Extrude Line {1314,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1315,{0,0,h}}{Recombine;Layers{{3,4,4},{7021,7022,7023},{a/h,(a+b)/h,1}};}; +Extrude Line {1315,{0,0,h}}{Recombine;Layers{{3,4,4},{a/h,(a+b)/h,1}};}; -Extrude Line {1316,{0,0,b}}{Recombine;Layers{{4},{7001},{1}};}; -Extrude Line {1324,{0,0,b}}{Recombine;Layers{{4},{7007},{1}};}; -Extrude Line {1332,{0,0,b}}{Recombine;Layers{{4},{7013},{1}};}; -Extrude Line {1340,{0,0,b}}{Recombine;Layers{{4},{7019},{1}};}; +Extrude Line {1316,{0,0,b}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1324,{0,0,b}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1332,{0,0,b}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1340,{0,0,b}}{Recombine;Layers{{4},{1}};}; -Extrude Line {1348,{0,0,c}}{Recombine;Layers{{4},{7002},{1}};}; -Extrude Line {1352,{0,0,c}}{Recombine;Layers{{4},{7008},{1}};}; -Extrude Line {1356,{0,0,c}}{Recombine;Layers{{4},{7014},{1}};}; -Extrude Line {1360,{0,0,c}}{Recombine;Layers{{4},{7020},{1}};}; +Extrude Line {1348,{0,0,c}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1352,{0,0,c}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1356,{0,0,c}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1360,{0,0,c}}{Recombine;Layers{{4},{1}};}; @@ -91,32 +91,32 @@ a=9.006-9.056 ; b=(8.976-9.006) ; c = -0.139564 -a-b ; -Extrude Line {1408,{0,0,a}}{Recombine;Layers{{3},{8000},{1}};}; +Extrude Line {1408,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1409,{0,0,-h}}{Recombine;Layers{{3,4,4},{8003,8004,8005},{-a/h,-(a+b)/h,1}};}; +Extrude Line {1409,{0,0,-h}}{Recombine;Layers{{3,4,4},{-a/h,-(a+b)/h,1}};}; -Extrude Line {1410,{0,0,a}}{Recombine;Layers{{3},{8006},{1}};}; +Extrude Line {1410,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1411,{0,0,-h}}{Recombine;Layers{{3,4,4},{8009,8010,8011},{-a/h,-(a+b)/h,1}};}; +Extrude Line {1411,{0,0,-h}}{Recombine;Layers{{3,4,4},{-a/h,-(a+b)/h,1}};}; -Extrude Line {1412,{0,0,a}}{Recombine;Layers{{3},{8012},{1}};}; +Extrude Line {1412,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1413,{0,0,-h}}{Recombine;Layers{{3,4,4},{8015,8016,8017},{-a/h,-(a+b)/h,1}};}; +Extrude Line {1413,{0,0,-h}}{Recombine;Layers{{3,4,4},{-a/h,-(a+b)/h,1}};}; -Extrude Line {1414,{0,0,a}}{Recombine;Layers{{3},{8018},{1}};}; +Extrude Line {1414,{0,0,a}}{Recombine;Layers{{3},{1}};}; -Extrude Line {1415,{0,0,-h}}{Recombine;Layers{{3,4,4},{8021,8022,8023},{-a/h,-(a+b)/h,1}};}; +Extrude Line {1415,{0,0,-h}}{Recombine;Layers{{3,4,4},{-a/h,-(a+b)/h,1}};}; -Extrude Line {1416,{0,0,b}}{Recombine;Layers{{4},{8001},{1}};}; -Extrude Line {1424,{0,0,b}}{Recombine;Layers{{4},{8007},{1}};}; -Extrude Line {1432,{0,0,b}}{Recombine;Layers{{4},{8013},{1}};}; -Extrude Line {1440,{0,0,b}}{Recombine;Layers{{4},{8019},{1}};}; +Extrude Line {1416,{0,0,b}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1424,{0,0,b}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1432,{0,0,b}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1440,{0,0,b}}{Recombine;Layers{{4},{1}};}; -Extrude Line {1448,{0,0,c}}{Recombine;Layers{{4},{8002},{1}};}; -Extrude Line {1452,{0,0,c}}{Recombine;Layers{{4},{8008},{1}};}; -Extrude Line {1456,{0,0,c}}{Recombine;Layers{{4},{8014},{1}};}; -Extrude Line {1460,{0,0,c}}{Recombine;Layers{{4},{8020},{1}};}; +Extrude Line {1448,{0,0,c}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1452,{0,0,c}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1456,{0,0,c}}{Recombine;Layers{{4},{1}};}; +Extrude Line {1460,{0,0,c}}{Recombine;Layers{{4},{1}};}; diff --git a/benchmarks/misc/faraway/faraway.geo b/benchmarks/misc/faraway/faraway.geo index 9601396c94d857103b7c95ed85decfcf77a8d2d4..0cfe9470792a8aa277f8e176596b23465de88c4f 100644 --- a/benchmarks/misc/faraway/faraway.geo +++ b/benchmarks/misc/faraway/faraway.geo @@ -110,11 +110,11 @@ Recombine Surface (2048) ; //5ème étape : extrusion -Extrude Surface {2038,{0,0,0.03}}{Layers{{4},{9000},{1}};Recombine;}; -Extrude Surface {2040,{0,0,0.03}}{Layers{{4},{9001},{1}};Recombine;}; -Extrude Surface {2042,{0,0,0.03}}{Layers{{4},{9002},{1}};Recombine;}; -Extrude Surface {2044,{0,0,0.03}}{Layers{{4},{9003},{1}};Recombine;}; -Extrude Surface {2048,{0,0,0.03}}{Layers{{4},{9004},{1}};Recombine;}; +Extrude Surface {2038,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {2040,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {2042,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {2044,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {2048,{0,0,0.03}}{Layers{4};Recombine;}; @@ -241,11 +241,11 @@ Ruled Surface(5048) = {5045,5047}; //5ème étape : extrusion -Extrude Surface {5038,{0,0,0.03}}{Layers{{4},{10000},{1}};Recombine;}; -Extrude Surface {5040,{0,0,0.03}}{Layers{{4},{10001},{1}};Recombine;}; -Extrude Surface {5042,{0,0,0.03}}{Layers{{4},{10002},{1}};Recombine;}; -Extrude Surface {5044,{0,0,0.03}}{Layers{{4},{10003},{1}};Recombine;}; -Extrude Surface {5048,{0,0,0.03}}{Layers{{4},{10004},{1}};Recombine;}; +Extrude Surface {5038,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {5040,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {5042,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {5044,{0,0,0.03}}{Layers{4};Recombine;}; +Extrude Surface {5048,{0,0,0.03}}{Layers{4};Recombine;};