Salome HOME
Typo-fix by Kunda
[modules/smesh.git] / src / StdMeshers / StdMeshers_QuadToTriaAdaptor.cxx
index f133994fa1d1041d097d51ad43a7dfe833175ff5..95617d10a6979ea00e6ffe0f33946fd4e7fffb80 100644 (file)
@@ -174,7 +174,7 @@ namespace
                PrmJ->GetNodeIndex( otherFaceNode ) >= 0 ))
           continue; // f is a base quadrangle
 
-        // check projections of face direction (baOFN) to triange normals (nI and nJ)
+        // check projections of face direction (baOFN) to triangle normals (nI and nJ)
         gp_Vec baOFN( base2, SMESH_TNodeXYZ( otherFaceNode ));
         if ( nI * baOFN > 0 && nJ * baOFN > 0 &&
              baI* baOFN > 0 && baJ* baOFN > 0 ) // issue 0023212
@@ -296,12 +296,13 @@ namespace
       SMESH_ComputeErrorPtr& err = sm->GetComputeError();
       if ( !err || err->IsOK() )
       {
-        err = SMESH_ComputeError::New( COMPERR_BAD_INPUT_MESH, msg, sm->GetAlgo() );
-        err->myBadElements.push_back( face1 );
-        err->myBadElements.push_back( face2 );
+        SMESH_BadInputElements* badElems =
+          new SMESH_BadInputElements( mesh.GetMeshDS(),COMPERR_BAD_INPUT_MESH, msg, sm->GetAlgo() );
+        badElems->add( face1 );
+        badElems->add( face2 );
+        err.reset( badElems );
       }
     }
-    //throw SALOME_Exception( msg.c_str() );
 
     return false; // == "algo fails"
   }
@@ -485,9 +486,6 @@ static bool HasIntersection3(const gp_Pnt& P, const gp_Pnt& PC, gp_Pnt& Pint,
   gp_XYZ vert1 = P2.XYZ();
   gp_XYZ vert2 = P3.XYZ();
 
-  /* calculate distance from vert0 to ray origin */
-  gp_XYZ  tvec = orig - vert0;
-
   gp_XYZ edge1 = vert1 - vert0;
   gp_XYZ edge2 = vert2 - vert0;
 
@@ -497,9 +495,13 @@ static bool HasIntersection3(const gp_Pnt& P, const gp_Pnt& PC, gp_Pnt& Pint,
   /* if determinant is near zero, ray lies in plane of triangle */
   double det = edge1 * pvec;
 
-  if (det > -EPSILON && det < EPSILON)
+  const double ANGL_EPSILON = 1e-12;
+  if ( det > -ANGL_EPSILON && det < ANGL_EPSILON )
     return false;
 
+  /* calculate distance from vert0 to ray origin */
+  gp_XYZ  tvec = orig - vert0;
+
   /* calculate U parameter and test bounds */
   double u = ( tvec * pvec ) / det;
   //if (u < 0.0 || u > 1.0)
@@ -615,7 +617,7 @@ bool StdMeshers_QuadToTriaAdaptor::LimitHeight (gp_Pnt&
         }
       }
     }
-    if ( height < 1e-5 * idealHeight && intFace )
+    if ( height < 1e-2 * idealHeight && intFace )
       return overlapError( aMesh, NotCheckedFace, intFace, Shape );
   }
 
@@ -665,7 +667,7 @@ bool StdMeshers_QuadToTriaAdaptor::LimitHeight (gp_Pnt&
 
 //================================================================================
 /*!
- * \brief Prepare data for the given face
+ * \brief Retrieve data of the given face
  *  \param PN - coordinates of face nodes
  *  \param VN - cross products of vectors (PC-PN(i)) ^ (PC-PN(i+1))
  *  \param FNodes - face nodes
@@ -704,7 +706,8 @@ int StdMeshers_QuadToTriaAdaptor::Preparation(const SMDS_MeshElement*       face
   int nbp = 4;
 
   int j = 0;
-  for(i=1; i<4; i++) {
+  for ( i = 1; i < 4; i++ )
+  {
     j = i+1;
     for(; j<=4; j++) {
       if( PN(i).Distance(PN(j)) < 1.e-6 )
@@ -712,11 +715,10 @@ int StdMeshers_QuadToTriaAdaptor::Preparation(const SMDS_MeshElement*       face
     }
     if(j<=4) break;
   }
-  //int deg_num = IsDegenarate(PN);
-  //if(deg_num>0) {
+
   bool hasdeg = false;
-  if(i<4) {
-    //cout<<"find degeneration"<<endl;
+  if ( i < 4 )
+  {
     hasdeg = true;
     gp_Pnt Pdeg = PN(i);
 
@@ -727,7 +729,6 @@ int StdMeshers_QuadToTriaAdaptor::Preparation(const SMDS_MeshElement*       face
       gp_Pnt Ptmp(N->X(),N->Y(),N->Z());
       if(Pdeg.Distance(Ptmp)<1.e-6) {
         DegNode = N;
-        //DegNode = const_cast<SMDS_MeshNode*>(N);
         break;
       }
     }
@@ -747,6 +748,7 @@ int StdMeshers_QuadToTriaAdaptor::Preparation(const SMDS_MeshElement*       face
 
   PN.SetValue(nbp+1,PN(1));
   FNodes[nbp] = FNodes[0];
+
   // find normal direction
   gp_Vec V1(PC,PN(nbp));
   gp_Vec V2(PC,PN(1));
@@ -788,7 +790,6 @@ int StdMeshers_QuadToTriaAdaptor::Preparation(const SMDS_MeshElement*       face
     }
   }
 
-  //cout<<"  VNorm("<<VNorm.X()<<","<<VNorm.Y()<<","<<VNorm.Z()<<")"<<endl;
   return hasdeg ? DEGEN_QUAD : QUAD;
 }
 
@@ -1012,6 +1013,8 @@ bool StdMeshers_QuadToTriaAdaptor::Compute(SMESH_Mesh& aMesh)
       groupDS = 0;
   }
 
+  const bool toFindVolumes = aMesh.NbVolumes() > 0;
+
   vector<const SMDS_MeshElement*> myPyramids;
   SMESH_MesherHelper helper(aMesh);
   helper.IsQuadraticSubMesh(aMesh.GetShapeToMesh());
@@ -1023,13 +1026,16 @@ bool StdMeshers_QuadToTriaAdaptor::Compute(SMESH_Mesh& aMesh)
   if ( !myElemSearcher )
     myElemSearcher = SMESH_MeshAlgos::GetElementSearcher( *meshDS );
   SMESH_ElementSearcher* searcher = const_cast<SMESH_ElementSearcher*>(myElemSearcher);
+  SMESHUtils::Deleter<SMESH_ElementSearcher>
+    volSearcher( SMESH_MeshAlgos::GetElementSearcher( *meshDS ));
+  vector< const SMDS_MeshElement* > suspectFaces, foundVolumes;
 
   TColgp_Array1OfPnt PN(1,5);
   TColgp_Array1OfVec VN(1,4);
   vector<const SMDS_MeshNode*> FNodes(5);
   TColgp_SequenceOfPnt aContour;
 
-  SMDS_FaceIteratorPtr fIt = meshDS->facesIterator(/*idInceasingOrder=*/true);
+  SMDS_FaceIteratorPtr fIt = meshDS->facesIterator();
   while( fIt->more())
   {
     const SMDS_MeshElement* face = fIt->next();
@@ -1042,7 +1048,7 @@ bool StdMeshers_QuadToTriaAdaptor::Compute(SMESH_Mesh& aMesh)
     if ( what == NOT_QUAD )
       continue;
     if ( volumes[0] && volumes[1] )
-      continue; // face is shared by two volumes - no space for a pyramid
+      continue; // face is shared by two volumes - no room for a pyramid
 
     if ( what == DEGEN_QUAD )
     {
@@ -1143,6 +1149,7 @@ bool StdMeshers_QuadToTriaAdaptor::Compute(SMESH_Mesh& aMesh)
     }
 
     // Restrict pyramid height by intersection with other faces
+
     gp_Vec tmpDir(PC,PCbest); tmpDir.Normalize();
     double tmp = PN(1).Distance(PN(3)) + PN(2).Distance(PN(4));
     // far points: in (PC, PCbest) direction and vice-versa
@@ -1154,8 +1161,24 @@ bool StdMeshers_QuadToTriaAdaptor::Compute(SMESH_Mesh& aMesh)
     gp_Pnt intPnt     [2];
     int    intFaceInd [2] = { 0, 0 };
 
+    if ( toFindVolumes && 0 ) // non-conformal mesh is not suitable for any mesher so far
+    {
+      // there are volumes in the mesh, in a non-conformal mesh an neighbor
+      // volume can be not found yet
+      for ( int isRev = 0; isRev < 2; ++isRev )
+      {
+        if ( volumes[isRev] ) continue;
+        gp_Pnt testPnt = PC.XYZ() + tmpDir.XYZ() * height * ( isRev ? -0.1: 0.1 );
+        foundVolumes.clear();
+        if ( volSearcher->FindElementsByPoint( testPnt, SMDSAbs_Volume, foundVolumes ))
+          volumes[isRev] = foundVolumes[0];
+      }
+      if ( volumes[0] && volumes[1] )
+        continue; // no room for a pyramid
+    }
+
     gp_Ax1 line( PC, tmpDir );
-    vector< const SMDS_MeshElement* > suspectFaces;
+    suspectFaces.clear();
     searcher->GetElementsNearLine( line, SMDSAbs_Face, suspectFaces);
 
     for ( size_t iF = 0; iF < suspectFaces.size(); ++iF )
@@ -1208,7 +1231,7 @@ bool StdMeshers_QuadToTriaAdaptor::Compute(SMESH_Mesh& aMesh)
       if ( !intersected[isRev] ) continue;
       double pyramidH = Min( height, dist2int[isRev]/3. );
       gp_Pnt    Papex = PC.XYZ() + tmpDir.XYZ() * (isRev ? -pyramidH : pyramidH);
-      if ( pyramidH < 1e-5 * height )
+      if ( pyramidH < 1e-2 * height )
         return overlapError( aMesh, face, suspectFaces[ intFaceInd[isRev] ] );
 
       if ( !LimitHeight( Papex, PC, PN, FNodes, aMesh, face, /*UseApexRay=*/false ))