Salome HOME
0022098: EDF 2036 SMESH: Create groups from none conected parts of a mesh
[modules/smesh.git] / src / StdMeshers / StdMeshers_Prism_3D.cxx
index 6584b16652576dabb3e1e993186a55e764e716dc..149e35a0ebbfbc8f1b9b5ea872f9857e6ff69afe 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2013  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
@@ -68,6 +68,11 @@ using namespace std;
 // gp_Pnt p (xyz); \
 // cout << msg << " ("<< p.X() << "; " <<p.Y() << "; " <<p.Z() << ") " <<endl;\
 // }
+#ifdef _DEBUG_
+#define DBGOUT(msg) //cout << msg << endl;
+#else
+#define DBGOUT(msg)
+#endif
 
 namespace TAssocTool = StdMeshers_ProjectionUtils;
 
@@ -464,7 +469,8 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
   // look for meshed FACEs ("source" FACEs) that must be prism bottoms
   list< TopoDS_Face > meshedFaces, notQuadMeshedFaces, notQuadFaces;
   const bool meshHasQuads = ( theMesh.NbQuadrangles() > 0 );
-  for ( int iF = 1; iF < faceToSolids.Extent(); ++iF )
+  //StdMeshers_Quadrangle_2D* quadAlgo = TQuadrangleAlgo::instance( this );
+  for ( int iF = 1; iF <= faceToSolids.Extent(); ++iF )
   {
     const TopoDS_Face& face = TopoDS::Face( faceToSolids.FindKey( iF ));
     SMESH_subMesh*   faceSM = theMesh.GetSubMesh( face );
@@ -480,10 +486,13 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
       else
         meshedFaces.push_back( face );
     }
-    else if ( myHelper->Count( face, TopAbs_EDGE, /*ignoreSame=*/false ) != 4 )
-    {
-      notQuadFaces.push_back( face );
-    }
+    // not add not quadrilateral FACE as we can't compute it
+    // else if ( !quadAlgo->CheckNbEdges( theMesh, face ))
+    // // not add not quadrilateral FACE as it can be a prism side
+    // // else if ( myHelper->Count( face, TopAbs_EDGE, /*ignoreSame=*/false ) != 4 )
+    // {
+    //   notQuadFaces.push_back( face );
+    // }
   }
   // notQuadFaces are of medium priority, put them before ordinary meshed faces
   meshedFaces.splice( meshedFaces.begin(), notQuadFaces );
@@ -1043,6 +1052,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
 {
   SMESH_Mesh*     mesh = myHelper->GetMesh();
   SMESHDS_Mesh* meshDS = myHelper->GetMeshDS();
+  DBGOUT( endl << "COMPUTE Prism " << meshDS->ShapeToIndex( thePrism.myShape3D ));
 
   TProjction1dAlgo* projector1D = TProjction1dAlgo::instance( this );
   StdMeshers_Quadrangle_2D* quadAlgo = TQuadrangleAlgo::instance( this, myHelper );
@@ -1108,6 +1118,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
         const TopoDS_Edge& srcE = lftSide->Edge(i);
         SMESH_subMesh*    srcSM = mesh->GetSubMesh( srcE );
         if ( !srcSM->IsMeshComputed() ) {
+          DBGOUT( "COMPUTE V edge " << srcSM->GetId() );
           srcSM->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE );
           srcSM->ComputeStateEngine       ( SMESH_subMesh::COMPUTE );
           if ( !srcSM->IsMeshComputed() )
@@ -1142,7 +1153,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
         }
         continue;
       }
-      // Compute
+      // Compute 'vertical projection'
       if ( nbTgtMeshed == 0 )
       {
         // compute nodes on target VERTEXes
@@ -1161,8 +1172,9 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
         }
 
         // compute nodes on target EDGEs
+        DBGOUT( "COMPUTE V edge (proj) " << shapeID( lftSide->Edge(0)));
         rgtSide->Reverse(); // direct it same as the lftSide
-        myHelper->SetElementsOnShape( false );
+        myHelper->SetElementsOnShape( false ); // myHelper holds the prism shape
         TopoDS_Edge tgtEdge;
         for ( size_t iN = 1; iN < srcNodeStr.size()-1; ++iN ) // add nodes
         {
@@ -1173,25 +1185,24 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
         }
         for ( size_t iN = 1; iN < srcNodeStr.size(); ++iN ) // add segments
         {
-          SMDS_MeshElement* newEdge = myHelper->AddEdge( newNodes[ iN-1 ], newNodes[ iN ] );
+          // find an EDGE to set a new segment
           std::pair<int, TopAbs_ShapeEnum> id2type = 
             myHelper->GetMediumPos( newNodes[ iN-1 ], newNodes[ iN ] );
-          if ( id2type.second == TopAbs_EDGE )
-          {
-            meshDS->SetMeshElementOnShape( newEdge, id2type.first );
-          }
-          else // new nodes are on different EDGEs; put one of them on VERTEX
+          if ( id2type.second != TopAbs_EDGE )
           {
+            // new nodes are on different EDGEs; put one of them on VERTEX
             const int      edgeIndex = rgtSide->EdgeIndex( srcNodeStr[ iN-1 ].normParam );
             const double vertexParam = rgtSide->LastParameter( edgeIndex );
             const gp_Pnt           p = BRep_Tool::Pnt( rgtSide->LastVertex( edgeIndex ));
             const int         isPrev = ( Abs( srcNodeStr[ iN-1 ].normParam - vertexParam ) <
                                          Abs( srcNodeStr[ iN   ].normParam - vertexParam ));
-            meshDS->SetMeshElementOnShape( newEdge, newNodes[ iN-(1-isPrev) ]->getshapeId() );
             meshDS->UnSetNodeOnShape( newNodes[ iN-isPrev ] );
             meshDS->SetNodeOnVertex ( newNodes[ iN-isPrev ], rgtSide->LastVertex( edgeIndex ));
-            meshDS->MoveNode( newNodes[ iN-isPrev ], p.X(), p.Y(), p.Z() );
+            meshDS->MoveNode        ( newNodes[ iN-isPrev ], p.X(), p.Y(), p.Z() );
+            id2type.first = newNodes[ iN-(1-isPrev) ]->getshapeId();
           }
+          SMDS_MeshElement* newEdge = myHelper->AddEdge( newNodes[ iN-1 ], newNodes[ iN ] );
+          meshDS->SetMeshElementOnShape( newEdge, id2type.first );
         }
         myHelper->SetElementsOnShape( true );
         for ( int i = 0; i < rgtSide->NbEdges(); ++i ) // update state of sub-meshes
@@ -1228,28 +1239,43 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
       // Top EDGEs must be projections from the bottom ones
       // to compute stuctured quad mesh on wall FACEs
       // ---------------------------------------------------
-      const TopoDS_Edge& botE = (*quad)->side[ QUAD_BOTTOM_SIDE ]->Edge(0);
-      const TopoDS_Edge& topE = (*quad)->side[ QUAD_TOP_SIDE    ]->Edge(0);
-
-      projector1D->myHyp.SetSourceEdge( botE );
-
-      SMESH_subMesh* tgtEdgeSm = mesh->GetSubMesh( topE );
-      if ( !tgtEdgeSm->IsMeshComputed() )
       {
-        // compute nodes on VERTEXes
-        tgtEdgeSm->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE );
-        // project segments
-        projector1D->InitComputeError();
-        bool ok = projector1D->Compute( *mesh, topE );
-        if ( !ok )
+        const TopoDS_Edge& botE = (*quad)->side[ QUAD_BOTTOM_SIDE ]->Edge(0);
+        const TopoDS_Edge& topE = (*quad)->side[ QUAD_TOP_SIDE    ]->Edge(0);
+        SMESH_subMesh*    botSM = mesh->GetSubMesh( botE );
+        SMESH_subMesh*    topSM = mesh->GetSubMesh( topE );
+        SMESH_subMesh*    srcSM = botSM;
+        SMESH_subMesh*    tgtSM = topSM;
+        if ( !srcSM->IsMeshComputed() && topSM->IsMeshComputed() )
+          std::swap( srcSM, tgtSM );
+
+        if ( !srcSM->IsMeshComputed() )
         {
-          SMESH_ComputeErrorPtr err = projector1D->GetComputeError();
-          if ( err->IsOK() ) err->myName = COMPERR_ALGO_FAILED;
-          tgtEdgeSm->GetComputeError() = err;
-          return false;
+          DBGOUT( "COMPUTE H edge " << srcSM->GetId());
+          srcSM->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE ); // nodes on VERTEXes
+          srcSM->ComputeStateEngine( SMESH_subMesh::COMPUTE );        // segments on the EDGE
+        }
+        srcSM->ComputeStateEngine( SMESH_subMesh::CHECK_COMPUTE_STATE );
+
+        if ( !tgtSM->IsMeshComputed() )
+        {
+          // compute nodes on VERTEXes
+          tgtSM->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE );
+          // project segments
+          DBGOUT( "COMPUTE H edge (proj) " << tgtSM->GetId());
+          projector1D->myHyp.SetSourceEdge( TopoDS::Edge( srcSM->GetSubShape() ));
+          projector1D->InitComputeError();
+          bool ok = projector1D->Compute( *mesh, tgtSM->GetSubShape() );
+          if ( !ok )
+          {
+            SMESH_ComputeErrorPtr err = projector1D->GetComputeError();
+            if ( err->IsOK() ) err->myName = COMPERR_ALGO_FAILED;
+            tgtSM->GetComputeError() = err;
+            return false;
+          }
         }
+        tgtSM->ComputeStateEngine( SMESH_subMesh::CHECK_COMPUTE_STATE );
       }
-      tgtEdgeSm->ComputeStateEngine( SMESH_subMesh::CHECK_COMPUTE_STATE );
 
       // Compute quad mesh on wall FACEs
       // -------------------------------
@@ -1264,6 +1290,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
                               "Not all edges have valid algorithm and hypothesis"));
         // mesh the <face>
         quadAlgo->InitComputeError();
+        DBGOUT( "COMPUTE Quad face " << fSM->GetId());
         bool ok = quadAlgo->Compute( *mesh, face );
         fSM->GetComputeError() = quadAlgo->GetComputeError();
         if ( !ok )
@@ -1545,7 +1572,7 @@ bool StdMeshers_Prism_3D::assocOrProjBottom2Top()
     _gen->Compute( *myHelper->GetMesh(), botSM->GetSubShape() );
     botSMDS = botSM->GetSubMeshDS();
     if ( !botSMDS || botSMDS->NbElements() == 0 )
-      return toSM( error(TCom("No elememts on face #") << botSM->GetId() ));
+      return toSM( error(TCom("No elements on face #") << botSM->GetId() ));
   }
 
   bool needProject = !topSM->IsMeshComputed();
@@ -2195,8 +2222,10 @@ bool StdMeshers_PrismAsBlock::Init(SMESH_MesherHelper*         helper,
     {
       if ( len2edgeMap.size() != nbEdges )
         RETURN_BAD_RESULT("Uniqueness of edge lengths not assured");
-      map< double, int >::reverse_iterator maxLen_i = len2edgeMap.rbegin();
-      map< double, int >::reverse_iterator midLen_i = ++len2edgeMap.rbegin();
+
+      multimap< double, int >::reverse_iterator maxLen_i = len2edgeMap.rbegin();
+      multimap< double, int >::reverse_iterator midLen_i = ++len2edgeMap.rbegin();
+
       double maxLen = maxLen_i->first;
       double midLen = ( len2edgeMap.size() == 1 ) ? 0 : midLen_i->first;
       switch ( nbEdges ) {