From c89adbc0173b6e3e13cd405e74cdb13958db5678 Mon Sep 17 00:00:00 2001 From: Thomas Toulorge <thomas.toulorge@mines-paristech.fr> Date: Tue, 8 Jul 2014 14:18:37 +0000 Subject: [PATCH] Changed tolerance w.r.t. level-set value in mesh splitting. Changed incompressible naca0012_ls test case. --- Geo/MElementCut.cpp | 6 +++--- Geo/gmshLevelset.cpp | 3 ++- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/Geo/MElementCut.cpp b/Geo/MElementCut.cpp index 7f087764ab..3d14181dce 100644 --- a/Geo/MElementCut.cpp +++ b/Geo/MElementCut.cpp @@ -606,7 +606,7 @@ static void elementSplitMesh(MElement *e, std::vector<gLevelset *> &RPN, int gLsTag = RPN[RPN.size() - 1]->getTag(); MElement *copy = e->copy(vertexMap, newParents, newDomains); - double eps = 1.e-6; + double eps = 1.e-10; bool splitElem = false; //split acording to center of gravity @@ -622,14 +622,14 @@ static void elementSplitMesh(MElement *e, std::vector<gLevelset *> &RPN, int ils = 0; for(int k = 1; k < e->getNumVertices(); k++){ int lsTag2 = (verticesLs(iLs, e->getVertex(k)->getIndex()) > eps) ? -1 : 1; - if (lsTag * lsTag2 < 0.0) { + if (lsTag * lsTag2 < 0) { lsTag = -1; splitElem = true; if(RPN.size() > 1){ for(; ils < verticesLs.size1() - 1; ils++){ int lsT1 = (verticesLs(ils, e->getVertex(0)->getIndex()) > eps) ? -1 : 1; int lsT2 = (verticesLs(ils, e->getVertex(k)->getIndex()) > eps) ? -1 : 1; - if(lsT1 * lsT2 < 0.0) break; + if(lsT1 * lsT2 < 0) break; } for(int i = 0; i <= ils; i++) if(!RPN[i]->isPrimitive()) ils++; diff --git a/Geo/gmshLevelset.cpp b/Geo/gmshLevelset.cpp index d695bb9cae..96599d5b7f 100644 --- a/Geo/gmshLevelset.cpp +++ b/Geo/gmshLevelset.cpp @@ -1250,7 +1250,7 @@ void gLevelsetNACA00::getClosestBndPoint(double x, double y, double z, { static const int maxIter = 100; - static const double tol = 1.e-8; + static const double tol = 1.e-10; const double tolr = tol/_c; // Tolerance (scaled bu chord) in = false; // Whether the point is inside the airfoil @@ -1283,6 +1283,7 @@ void gLevelsetNACA00::getClosestBndPoint(double x, double y, double z, if (fabs(mIncr) < tolr) break; else xtb -= mIncr; if (xtb < tolr) xtb = tolr; + if (xtb > _c-tolr) xtb = _c-tolr; } xb = _x0+xtb; yb = (y >= _y0) ? _y0+ytb : _y0-ytb; -- GitLab