when optimizing 2D meshes: use the surface normal in computation of the jacobian

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;

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);

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);