Salome HOME
Merge 'master' branch into 'V9_dev' branch.
[modules/smesh.git] / src / StdMeshers / StdMeshers_MEFISTO_2D.cxx
index 55091c7f9c598ba557cd8dcdf28890af4a7e16ba..bbc594b163c9e2a3de4670ff8c6e8aabd4783a9d 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2016  CEA/DEN, EDF R&D, OPEN CASCADE
 //
 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
@@ -6,7 +6,7 @@
 // This library is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
 // License as published by the Free Software Foundation; either
-// version 2.1 of the License.
+// version 2.1 of the License, or (at your option) any later version.
 //
 // This library is distributed in the hope that it will be useful,
 // but WITHOUT ANY WARRANTY; without even the implied warranty of
@@ -20,7 +20,7 @@
 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
 //
 
-//  SMESH SMESH : implementaion of SMESH idl descriptions
+//  SMESH SMESH : implementation of SMESH idl descriptions
 //  File   : StdMeshers_MEFISTO_2D.cxx
 //           Moved here from SMESH_MEFISTO_2D.cxx
 //  Author : Paul RASCLE, EDF
@@ -31,6 +31,7 @@
 #include "SMDS_EdgePosition.hxx"
 #include "SMDS_MeshElement.hxx"
 #include "SMDS_MeshNode.hxx"
+#include "SMESHDS_Mesh.hxx"
 #include "SMESH_Comment.hxx"
 #include "SMESH_Gen.hxx"
 #include "SMESH_Mesh.hxx"
@@ -78,10 +79,9 @@ using namespace std;
  */
 //=============================================================================
 
-StdMeshers_MEFISTO_2D::StdMeshers_MEFISTO_2D(int hypId, int studyId, SMESH_Gen * gen):
-  SMESH_2D_Algo(hypId, studyId, gen)
+StdMeshers_MEFISTO_2D::StdMeshers_MEFISTO_2D(int hypId, SMESH_Gen * gen):
+  SMESH_2D_Algo(hypId, gen)
 {
-  MESSAGE("StdMeshers_MEFISTO_2D::StdMeshers_MEFISTO_2D");
   _name = "MEFISTO_2D";
   _shapeType = (1 << TopAbs_FACE);
   _compatibleHypothesis.push_back("MaxElementArea");
@@ -103,7 +103,6 @@ StdMeshers_MEFISTO_2D::StdMeshers_MEFISTO_2D(int hypId, int studyId, SMESH_Gen *
 
 StdMeshers_MEFISTO_2D::~StdMeshers_MEFISTO_2D()
 {
-  MESSAGE("StdMeshers_MEFISTO_2D::~StdMeshers_MEFISTO_2D");
 }
 
 //=============================================================================
@@ -122,6 +121,9 @@ bool StdMeshers_MEFISTO_2D::CheckHypothesis
   _edgeLength = 0;
   _maxElementArea = 0;
 
+  if ( !error( StdMeshers_ViscousLayers2D::CheckHypothesis( aMesh, aShape, aStatus )))
+    return false;
+
   list <const SMESHDS_Hypothesis * >::const_iterator itl;
   const SMESHDS_Hypothesis *theHyp;
 
@@ -185,8 +187,6 @@ bool StdMeshers_MEFISTO_2D::CheckHypothesis
 
 bool StdMeshers_MEFISTO_2D::Compute(SMESH_Mesh & aMesh, const TopoDS_Shape & aShape)
 {
-  MESSAGE("StdMeshers_MEFISTO_2D::Compute");
-
   TopoDS_Face F = TopoDS::Face(aShape.Oriented(TopAbs_FORWARD));
 
   // helper builds quadratic mesh if necessary
@@ -203,7 +203,7 @@ bool StdMeshers_MEFISTO_2D::Compute(SMESH_Mesh & aMesh, const TopoDS_Shape & aSh
   // get all edges of a face
   TError problem;
   TWireVector wires =
-    StdMeshers_FaceSide::GetFaceWires( F, aMesh, skipMediumNodes, problem, proxyMesh );
+    StdMeshers_FaceSide::GetFaceWires( F, aMesh, skipMediumNodes, problem, _helper, proxyMesh );
   int nbWires = wires.size();
   if ( problem && !problem->IsOK() ) return error( problem );
   if ( nbWires == 0 ) return error( "Problem in StdMeshers_FaceSide::GetFaceWires()");
@@ -280,8 +280,6 @@ bool StdMeshers_MEFISTO_2D::Compute(SMESH_Mesh & aMesh, const TopoDS_Shape & aSh
 
     if (ierr == 0)
     {
-      MESSAGE("... End Triangulation Generated Triangle Number " << nbt);
-      MESSAGE("                                    Node Number " << nbst);
       StoreResult(nbst, uvst, nbt, nust, mefistoToDS, scalex, scaley);
       isOk = true;
     }
@@ -309,8 +307,6 @@ bool StdMeshers_MEFISTO_2D::Evaluate(SMESH_Mesh & aMesh,
                                      const TopoDS_Shape & aShape,
                                      MapShapeNbElems& aResMap)
 {
-  MESSAGE("StdMeshers_MEFISTO_2D::Evaluate");
-
   TopoDS_Face F = TopoDS::Face(aShape.Oriented(TopAbs_FORWARD));
 
   double aLen = 0.0;
@@ -584,7 +580,7 @@ bool StdMeshers_MEFISTO_2D::LoadPoints(TWireVector &                 wires,
     F = TopoDS::Face( _helper->GetSubShape() );
     TopExp::MapShapesAndAncestors( F, TopAbs_VERTEX, TopAbs_WIRE, VWMap );
     int nbVertices = 0;
-    for ( int iW = 0; iW < wires.size(); ++iW )
+    for ( size_t iW = 0; iW < wires.size(); ++iW )
       nbVertices += wires[ iW ]->NbEdges();
     if ( nbVertices == VWMap.Extent() )
       VWMap.Clear(); // wires have no common vertices
@@ -592,10 +588,10 @@ bool StdMeshers_MEFISTO_2D::LoadPoints(TWireVector &                 wires,
 
   int m = 0;
 
-  for ( int iW = 0; iW < wires.size(); ++iW )
+  for ( size_t iW = 0; iW < wires.size(); ++iW )
   {
     const vector<UVPtStruct>& uvPtVec = wires[ iW ]->GetUVPtStruct();
-    if ( uvPtVec.size() != wires[ iW ]->NbPoints() ) {
+    if ((int) uvPtVec.size() != wires[ iW ]->NbPoints() ) {
       return error(COMPERR_BAD_INPUT_MESH,SMESH_Comment("Unexpected nb of points on wire ")
                    << iW << ": " << uvPtVec.size()<<" != "<<wires[ iW ]->NbPoints()
                    << ", probably because of invalid node parameters on geom edges");
@@ -693,7 +689,7 @@ void StdMeshers_MEFISTO_2D::ComputeScaleOnFace(SMESH_Mesh &        aMesh,
   double xmax = -1.e300;
   double ymin = 1.e300;
   double ymax = -1.e300;
-  int nbp = 23;
+  const int nbp = 23;
   scalex = 1;
   scaley = 1;
 
@@ -717,13 +713,8 @@ void StdMeshers_MEFISTO_2D::ComputeScaleOnFace(SMESH_Mesh &        aMesh,
         ymin = p.Y();
       if (p.Y() > ymax)
         ymax = p.Y();
-      //    MESSAGE(" "<< f<<" "<<l<<" "<<param<<" "<<xmin<<" "<<xmax<<" "<<ymin<<" "<<ymax);
     }
   }
-  //   SCRUTE(xmin);
-  //   SCRUTE(xmax);
-  //   SCRUTE(ymin);
-  //   SCRUTE(ymax);
   double xmoy = (xmax + xmin) / 2.;
   double ymoy = (ymax + ymin) / 2.;
   double xsize = xmax - xmin;
@@ -751,25 +742,33 @@ void StdMeshers_MEFISTO_2D::ComputeScaleOnFace(SMESH_Mesh &        aMesh,
   }
   scalex = length_x / xsize;
   scaley = length_y / ysize;
-//   SCRUTE(xsize);
-//   SCRUTE(ysize);
   double xyratio = xsize*scalex/(ysize*scaley);
   const double maxratio = 1.e2;
-  //SCRUTE(xyratio);
   if (xyratio > maxratio) {
-    SCRUTE( scaley );
     scaley *= xyratio / maxratio;
-    SCRUTE( scaley );
   }
   else if (xyratio < 1./maxratio) {
-    SCRUTE( scalex );
     scalex *= 1 / xyratio / maxratio;
-    SCRUTE( scalex );
   }
-  ASSERT(scalex);
-  ASSERT(scaley);
 }
 
+// namespace
+// {
+//   bool isDegenTria( const SMDS_MeshNode * nn[3] )
+//   {
+//     SMESH_TNodeXYZ p1( nn[0] );
+//     SMESH_TNodeXYZ p2( nn[1] );
+//     SMESH_TNodeXYZ p3( nn[2] );
+//     gp_XYZ vec1 = p2 - p1;
+//     gp_XYZ vec2 = p3 - p1;
+//     gp_XYZ cross = vec1 ^ vec2;
+//     const double eps = 1e-100;
+//     return ( fabs( cross.X() ) < eps &&
+//              fabs( cross.Y() ) < eps &&
+//              fabs( cross.Z() ) < eps );
+//   }
+// }
+
 //=============================================================================
 /*!
  *  
@@ -780,12 +779,13 @@ void StdMeshers_MEFISTO_2D::StoreResult(Z nbst, R2 * uvst, Z nbt, Z * nust,
                                         vector< const SMDS_MeshNode*>&mefistoToDS,
                                         double scalex, double scaley)
 {
-  SMESHDS_Mesh * meshDS = _helper->GetMeshDS();
-  int faceID = _helper->GetSubShapeID();
+  _helper->SetElementsOnShape( true );
 
   TopoDS_Face F = TopoDS::Face( _helper->GetSubShape() );
   Handle(Geom_Surface) S = BRep_Tool::Surface( F );
 
+  //const size_t nbInputNodes = mefistoToDS.size();
+
   Z n = mefistoToDS.size(); // nb input points
   mefistoToDS.resize( nbst );
   for ( ; n < nbst; n++)
@@ -796,12 +796,7 @@ void StdMeshers_MEFISTO_2D::StoreResult(Z nbst, R2 * uvst, Z nbt, Z * nust,
       double v = uvst[n][1] / scaley;
       gp_Pnt P = S->Value(u, v);
 
-      SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
-      meshDS->SetNodeOnFace(node, faceID, u, v);
-
-      //MESSAGE(P.X()<<" "<<P.Y()<<" "<<P.Z());
-      mefistoToDS[n] = node;
-      //MESSAGE("NEW: "<<n<<" "<<mefistoToDS[n+1]);
+      mefistoToDS[n] = _helper->AddNode( P.X(), P.Y(), P.Z(), 0, u, v );
     }
   }
 
@@ -810,26 +805,32 @@ void StdMeshers_MEFISTO_2D::StoreResult(Z nbst, R2 * uvst, Z nbt, Z * nust,
   // triangle points must be in trigonometric order if face is Forward
   // else they must be put clockwise
 
-  bool triangleIsWellOriented = ( F.Orientation() == TopAbs_FORWARD );
+  int i1 = 1, i2 = 2;
+  if ( F.Orientation() != TopAbs_FORWARD )
+    std::swap( i1, i2 );
 
+  const SMDS_MeshNode * nn[3];
   for (n = 1; n <= nbt; n++)
   {
-    const SMDS_MeshNode * n1 = mefistoToDS[ nust[m++] - 1 ];
-    const SMDS_MeshNode * n2 = mefistoToDS[ nust[m++] - 1 ];
-    const SMDS_MeshNode * n3 = mefistoToDS[ nust[m++] - 1 ];
+    // const bool allNodesAreOld = ( nust[m + 0] <= nbInputNodes &&
+    //                               nust[m + 1] <= nbInputNodes &&
+    //                               nust[m + 2] <= nbInputNodes );
+    nn[ 0 ] = mefistoToDS[ nust[m++] - 1 ];
+    nn[ 1 ] = mefistoToDS[ nust[m++] - 1 ];
+    nn[ 2 ] = mefistoToDS[ nust[m++] - 1 ];
+    m++;
 
     // avoid creating degenetrated faces
-    bool isDegen = ( _helper->HasDegeneratedEdges() && ( n1 == n2 || n1 == n3 || n2 == n3 ));
+    bool isDegen = ( _helper->HasDegeneratedEdges() &&
+                     ( nn[0] == nn[1] || nn[1] == nn[2] || nn[2] == nn[0] ));
+
+    // It was an attempt to fix a problem of a zero area face whose all nodes
+    // are on one staight EDGE. But omitting this face makes a hole in the mesh :(
+    // if ( !isDegen && allNodesAreOld )
+    //   isDegen = isDegenTria( nn );
+
     if ( !isDegen )
-    {
-      SMDS_MeshElement * elt;
-      if (triangleIsWellOriented)
-        elt = _helper->AddFace(n1, n2, n3);
-      else
-        elt = _helper->AddFace(n1, n3, n2);
-      meshDS->SetMeshElementOnShape(elt, faceID);
-    }
-    m++;
+      _helper->AddFace( nn[0], nn[i1], nn[i2] );
   }
 
   // remove bad elements built on vertices shared by wires
@@ -849,7 +850,7 @@ void StdMeshers_MEFISTO_2D::StoreResult(Z nbst, R2 * uvst, Z nbt, Z * nust,
           nbSame++;
       if (nbSame > 1) {
         MESSAGE( "RM bad element " << elem->GetID());
-        meshDS->RemoveElement( elem );
+        _helper->GetMeshDS()->RemoveElement( elem );
       }
     }
   }