From 97ce1460f3ee38a86fc316fb18107d47a428618e Mon Sep 17 00:00:00 2001
From: Amaury Johnan <amjohnen@gmail.com>
Date: Fri, 19 Aug 2016 10:53:43 +0000
Subject: [PATCH] when optimizing 2D meshes: use the surface normal in
computation of the jacobian
---
Plugin/AnalyseCurvedMesh.cpp | 1 -
contrib/HighOrderMeshOptimizer/OptHomMesh.cpp | 13 ++++++++++++-
contrib/HighOrderMeshOptimizer/OptHomRun.cpp | 5 +++++
3 files changed, 17 insertions(+), 2 deletions(-)
diff --git a/Plugin/AnalyseCurvedMesh.cpp b/Plugin/AnalyseCurvedMesh.cpp
index b2c14347c8..918a4fbabd 100644
--- a/Plugin/AnalyseCurvedMesh.cpp
+++ b/Plugin/AnalyseCurvedMesh.cpp
@@ -275,7 +275,6 @@ void GMSH_AnalyseCurvedMeshPlugin::_computeMinMaxJandValidity(int dim)
normals->set(0, 0, 0);
normals->set(0, 1, 0);
normals->set(0, 2, 1);
- Msg::Error("Discrete face is 2D!");
}
}
break;
diff --git a/contrib/HighOrderMeshOptimizer/OptHomMesh.cpp b/contrib/HighOrderMeshOptimizer/OptHomMesh.cpp
index 1121b97601..d91c48c9c1 100644
--- a/contrib/HighOrderMeshOptimizer/OptHomMesh.cpp
+++ b/contrib/HighOrderMeshOptimizer/OptHomMesh.cpp
@@ -128,7 +128,7 @@ void Mesh::calcScaledNormalEl2D(const std::map<MElement*,GEntity*> &element2enti
SVector3 geoNorm(0.,0.,0.);
std::map<MElement*,GEntity*>::const_iterator itEl2ent = element2entity.find(_el[iEl]);
GEntity *ge = (itEl2ent == element2entity.end()) ? 0 : itEl2ent->second;
- const bool hasGeoNorm = ge && (ge->dim() == 2) && ge->haveParametrization();
+ bool hasGeoNorm = ge && (ge->dim() == 2) && ge->haveParametrization();
for (int i=0; i<jac->getNumPrimMapNodes(); i++) {
const int &iV = _el2V[iEl][i];
primNodesXYZ(i,0) = _xyz[iV].x();
@@ -145,6 +145,17 @@ void Mesh::calcScaledNormalEl2D(const std::map<MElement*,GEntity*> &element2enti
SPoint2 param = ((GFace*)ge)->parFromPoint(_el[iEl]->barycenter(true),false);
geoNorm = ((GFace*)ge)->normal(param);
}
+ if (!hasGeoNorm && ge && ge->geomType() == GEntity::DiscreteSurface) {
+ SBoundingBox3d bb = ge->bounds();
+ // If we don't have the CAD, check if the mesh is 2D:
+ if (!bb.empty() && bb.max().z() - bb.min().z() == .0) {
+ hasGeoNorm = true;
+ geoNorm = SVector3(0, 0, 1);
+ }
+ }
+ else if (!ge) {
+ Msg::Warning("don't have the entity");
+ }
fullMatrix<double> &elNorm = _scaledNormEl[iEl];
elNorm.resize(1,3);
diff --git a/contrib/HighOrderMeshOptimizer/OptHomRun.cpp b/contrib/HighOrderMeshOptimizer/OptHomRun.cpp
index 83ff11f5b0..3d51a4f3ae 100644
--- a/contrib/HighOrderMeshOptimizer/OptHomRun.cpp
+++ b/contrib/HighOrderMeshOptimizer/OptHomRun.cpp
@@ -644,7 +644,12 @@ void HighOrderMeshOptimizer(GModel *gm, OptHomParameters &p)
Msg::Info("Computing connectivity and bad elements for entity %d...",
entity->tag());
calcVertex2Elements(p.dim,entity,vertex2elements);
+
if (p.optPrimSurfMesh) calcElement2Entity(entity,element2entity);
+ // Otherwise, still compute element2entity: Hack for using the geometry
+ // normal in the computation of the Jacobian when optimizing surface meshes.
+ // Warning: Accurate for planar surface but not really for curved surface...
+ else calcElement2Entity(entity,element2entity);
}
OptHomPeriodicity periodicity = OptHomPeriodicity(entities);
--
GitLab