Salome HOME
SALOME_TESTS/Grids/smesh/mesh_Projection_2D_00/A0
[modules/smesh.git] / src / StdMeshers / StdMeshers_Prism_3D.cxx
index c43a126c3e2823188a2992bec17ac34dbfbf1e8a..3b6cf64816c46f54fbdc81083d160747f42466c8 100644 (file)
@@ -157,6 +157,46 @@ namespace {
       return algo;
     }
   };
+  //=======================================================================
+  /*!
+   * \brief Returns already computed EDGEs
+   */
+  void getPrecomputedEdges( SMESH_MesherHelper&    theHelper,
+                            const TopoDS_Shape&    theShape,
+                            vector< TopoDS_Edge >& theEdges)
+  {
+    theEdges.clear();
+
+    SMESHDS_Mesh* meshDS = theHelper.GetMeshDS();
+    SMESHDS_SubMesh* sm;
+
+    TopTools_IndexedMapOfShape edges;
+    TopExp::MapShapes( theShape, TopAbs_EDGE, edges );
+    for ( int iE = 1; iE <= edges.Extent(); ++iE )
+    {
+      const TopoDS_Shape edge = edges( iE );
+      if (( ! ( sm = meshDS->MeshElements( edge ))) ||
+          ( sm->NbElements() == 0 ))
+        continue;
+
+      // there must not be FACEs meshed with triangles and sharing a computed EDGE
+      // as the precomputed EDGEs are used for propagation other to 'vertical' EDGEs
+      bool faceFound = false;
+      PShapeIteratorPtr faceIt =
+        theHelper.GetAncestors( edge, *theHelper.GetMesh(), TopAbs_FACE );
+      while ( const TopoDS_Shape* face = faceIt->next() )
+
+        if (( sm = meshDS->MeshElements( *face )) &&
+            ( sm->NbElements() > 0 ) &&
+            ( !theHelper.IsSameElemGeometry( sm, SMDSGeom_QUADRANGLE ) ))
+        {
+          faceFound;
+          break;
+        }
+      if ( !faceFound )
+        theEdges.push_back( TopoDS::Edge( edge ));
+    }
+  }
 
   //================================================================================
   /*!
@@ -172,6 +212,7 @@ namespace {
     quad->side[ QUAD_TOP_SIDE  ].grid->Reverse();
     quad->side[ QUAD_LEFT_SIDE ].grid->Reverse();
     int edgeIndex = 0;
+    bool isComposite = false;
     for ( size_t i = 0; i < quad->side.size(); ++i )
     {
       StdMeshers_FaceSidePtr quadSide = quad->side[i];
@@ -179,7 +220,7 @@ namespace {
         if ( botE.IsSame( quadSide->Edge( iE )))
         {
           if ( quadSide->NbEdges() > 1 )
-            return false;
+            isComposite = true; //return false;
           edgeIndex = i;
           i = quad->side.size(); // to quit from the outer loop
           break;
@@ -190,7 +231,7 @@ namespace {
 
     quad->face = TopoDS::Face( face );
 
-    return true;
+    return !isComposite;
   }
 
   //================================================================================
@@ -614,6 +655,7 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
   meshedFaces.splice( meshedFaces.begin(), notQuadMeshedFaces );
 
   Prism_3D::TPrismTopo prism;
+  myPropagChains = 0;
 
   if ( nbSolids == 1 )
   {
@@ -623,6 +665,21 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
              compute( prism ));
   }
 
+  // find propagation chains from already computed EDGEs
+  vector< TopoDS_Edge > computedEdges;
+  getPrecomputedEdges( helper, theShape, computedEdges );
+  myPropagChains = new TopTools_IndexedMapOfShape[ computedEdges.size() + 1 ];
+  SMESHUtils::ArrayDeleter< TopTools_IndexedMapOfShape > pcDel( myPropagChains );
+  for ( size_t i = 0, nb = 0; i < computedEdges.size(); ++i )
+  {
+    StdMeshers_ProjectionUtils::GetPropagationEdge( &theMesh, TopoDS_Edge(),
+                                                    computedEdges[i], myPropagChains + nb );
+    if ( myPropagChains[ nb ].Extent() < 2 ) // an empty map is a termination sign
+      myPropagChains[ nb ].Clear();
+    else
+      nb++;
+  }
+
   TopTools_MapOfShape meshedSolids;
   list< Prism_3D::TPrismTopo > meshedPrism;
   TopTools_ListIteratorOfListOfShape solidIt;
@@ -650,7 +707,11 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
                !compute( prism ))
             return false;
 
-          meshedFaces.push_front( prism.myTop );
+          SMESHDS_SubMesh* smDS = theMesh.GetMeshDS()->MeshElements( prism.myTop );
+          if ( !myHelper->IsSameElemGeometry( smDS, SMDSGeom_QUADRANGLE ))
+          {
+            meshedFaces.push_front( prism.myTop );
+          }
           meshedPrism.push_back( prism );
         }
       }
@@ -691,6 +752,7 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
               prism.myBottom  = candidateF;
               mySetErrorToSM = false;
               if ( !myHelper->IsSubShape( candidateF, prismIt->myShape3D ) &&
+                   myHelper->IsSubShape( candidateF, solid ) &&
                    !myHelper->GetMesh()->GetSubMesh( candidateF )->IsMeshComputed() &&
                    initPrism( prism, solid ) &&
                    project2dMesh( prismIt->myBottom, candidateF))
@@ -698,8 +760,12 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
                 mySetErrorToSM = true;
                 if ( !compute( prism ))
                   return false;
-                meshedFaces.push_front( prism.myTop );
-                meshedFaces.push_front( prism.myBottom );
+                SMESHDS_SubMesh* smDS = theMesh.GetMeshDS()->MeshElements( prism.myTop );
+                if ( !myHelper->IsSameElemGeometry( smDS, SMDSGeom_QUADRANGLE ))
+                {
+                  meshedFaces.push_front( prism.myTop );
+                  meshedFaces.push_front( prism.myBottom );
+                }
                 meshedPrism.push_back( prism );
                 meshedSolids.Add( solid );
               }
@@ -841,10 +907,17 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
           if ( !quadList.back() )
             return toSM( error(TCom("Side face #") << shapeID( face )
                                << " not meshable with quadrangles"));
-          if ( ! setBottomEdge( *edge, quadList.back(), face ))
-            return toSM( error(TCom("Composite 'horizontal' edges are not supported")));
-          thePrism.myWallQuads.push_back( quadList );
-          faceMap.Add( face );
+          bool isCompositeBase = ! setBottomEdge( *edge, quadList.back(), face );
+          if ( isCompositeBase )
+          {
+            // it's OK if all EDGEs of the bottom side belongs to the bottom FACE
+            StdMeshers_FaceSidePtr botSide = quadList.back()->side[ QUAD_BOTTOM_SIDE ];
+            for ( int iE = 0; iE < botSide->NbEdges(); ++iE )
+              if ( !myHelper->IsSubShape( botSide->Edge(iE), thePrism.myBottom ))
+                return toSM( error(TCom("Composite 'horizontal' edges are not supported")));
+          }
+          if ( faceMap.Add( face ))
+            thePrism.myWallQuads.push_back( quadList );
           break;
         }
       }
@@ -1012,7 +1085,7 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
   // Projections on the top and bottom faces are taken from nodes existing
   // on these faces; find correspondence between bottom and top nodes
   myBotToColumnMap.clear();
-  if ( !assocOrProjBottom2Top( bottomToTopTrsf ) ) // it also fills myBotToColumnMap
+  if ( !assocOrProjBottom2Top( bottomToTopTrsf, thePrism ) ) // it also fills myBotToColumnMap
     return false;
 
 
@@ -1187,7 +1260,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
   SMESHDS_Mesh* meshDS = myHelper->GetMeshDS();
   DBGOUT( endl << "COMPUTE Prism " << meshDS->ShapeToIndex( thePrism.myShape3D ));
 
-  TProjction1dAlgo* projector1D = TProjction1dAlgo::instance( this );
+  TProjction1dAlgo*      projector1D = TProjction1dAlgo::instance( this );
   StdMeshers_Quadrangle_2D* quadAlgo = TQuadrangleAlgo::instance( this, myHelper );
 
   SMESH_HypoFilter hyp1dFilter( SMESH_HypoFilter::IsAlgo(),/*not=*/true);
@@ -1252,8 +1325,17 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
         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 );
+          TopoDS_Edge prpgSrcE = findPropagationSource( srcE );
+          if ( !prpgSrcE.IsNull() ) {
+            srcSM->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE );
+            projector1D->myHyp.SetSourceEdge( prpgSrcE );
+            projector1D->Compute( *mesh, srcE );
+            srcSM->ComputeStateEngine( SMESH_subMesh::CHECK_COMPUTE_STATE );
+          }
+          else {
+            srcSM->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE );
+            srcSM->ComputeStateEngine       ( SMESH_subMesh::COMPUTE );
+          }
           if ( !srcSM->IsMeshComputed() )
             return toSM( error( "Can't compute 1D mesh" ));
         }
@@ -1369,10 +1451,13 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
     Prism_3D::TQuadList::const_iterator quad = thePrism.myWallQuads[iW].begin();
     for ( ; quad != thePrism.myWallQuads[iW].end(); ++quad )
     {
-      // Top EDGEs must be projections from the bottom ones
-      // to compute stuctured quad mesh on wall FACEs
-      // ---------------------------------------------------
+      const TopoDS_Face& face = (*quad)->face;
+      SMESH_subMesh* fSM = mesh->GetSubMesh( face );
+      if ( ! fSM->IsMeshComputed() )
       {
+        // 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 ].grid->Edge(0);
         const TopoDS_Edge& topE = (*quad)->side[ QUAD_TOP_SIDE    ].grid->Edge(0);
         SMESH_subMesh*    botSM = mesh->GetSubMesh( botE );
@@ -1414,7 +1499,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
           // compute nodes on VERTEXes
           SMESH_subMeshIteratorPtr smIt = tgtSM->getDependsOnIterator(/*includeSelf=*/false);
           while ( smIt->more() )
-            smIt->next()->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE );
+            smIt->next()->ComputeStateEngine( SMESH_subMesh::COMPUTE );
           // project segments
           DBGOUT( "COMPUTE H edge (proj) " << tgtSM->GetId());
           projector1D->myHyp.SetSourceEdge( TopoDS::Edge( srcSM->GetSubShape() ));
@@ -1429,14 +1514,11 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
           }
         }
         tgtSM->ComputeStateEngine( SMESH_subMesh::CHECK_COMPUTE_STATE );
-      }
 
-      // Compute quad mesh on wall FACEs
-      // -------------------------------
-      const TopoDS_Face& face = (*quad)->face;
-      SMESH_subMesh* fSM = mesh->GetSubMesh( face );
-      if ( ! fSM->IsMeshComputed() )
-      {
+
+        // Compute quad mesh on wall FACEs
+        // -------------------------------
+
         // make all EDGES meshed
         fSM->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE );
         if ( !fSM->SubMeshesComputed() )
@@ -1464,6 +1546,22 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
   return true;
 }
 
+//=======================================================================
+/*!
+ * \brief Returns a source EDGE of propagation to a given EDGE
+ */
+//=======================================================================
+
+TopoDS_Edge StdMeshers_Prism_3D::findPropagationSource( const TopoDS_Edge& E )
+{
+  if ( myPropagChains )
+    for ( size_t i = 0; !myPropagChains[i].IsEmpty(); ++i )
+      if ( myPropagChains[i].Contains( E ))
+        return TopoDS::Edge( myPropagChains[i].FindKey( 1 ));
+
+  return TopoDS_Edge();
+}
+
 //=======================================================================
 //function : Evaluate
 //purpose  : 
@@ -1713,7 +1811,8 @@ void StdMeshers_Prism_3D::AddPrisms( vector<const TNodeColumn*> & columns,
  */
 //================================================================================
 
-bool StdMeshers_Prism_3D::assocOrProjBottom2Top( const gp_Trsf & bottomToTopTrsf )
+bool StdMeshers_Prism_3D::assocOrProjBottom2Top( const gp_Trsf & bottomToTopTrsf,
+                                                 const Prism_3D::TPrismTopo& thePrism)
 {
   SMESH_subMesh * botSM = myBlock.SubMesh( ID_BOT_FACE );
   SMESH_subMesh * topSM = myBlock.SubMesh( ID_TOP_FACE );
@@ -1730,7 +1829,7 @@ bool StdMeshers_Prism_3D::assocOrProjBottom2Top( const gp_Trsf & bottomToTopTrsf
   }
 
   bool needProject = !topSM->IsMeshComputed();
-  if ( !needProject && 
+  if ( !needProject &&
        (botSMDS->NbElements() != topSMDS->NbElements() ||
         botSMDS->NbNodes()    != topSMDS->NbNodes()))
   {
@@ -1756,19 +1855,70 @@ bool StdMeshers_Prism_3D::assocOrProjBottom2Top( const gp_Trsf & bottomToTopTrsf
   TopoDS_Face topFace = TopoDS::Face( myBlock.Shape( ID_TOP_FACE ));
   // associate top and bottom faces
   TAssocTool::TShapeShapeMap shape2ShapeMap;
-  if ( !TAssocTool::FindSubShapeAssociation( botFace, myBlock.Mesh(),
-                                             topFace, myBlock.Mesh(),
-                                             shape2ShapeMap) )
-    return toSM( error(TCom("Topology of faces #") << botSM->GetId()
-                       <<" and #"<< topSM->GetId() << " seems different" ));
+  const bool sameTopo = 
+    TAssocTool::FindSubShapeAssociation( botFace, myBlock.Mesh(),
+                                         topFace, myBlock.Mesh(),
+                                         shape2ShapeMap);
+  if ( !sameTopo )
+    for ( size_t iQ = 0; iQ < thePrism.myWallQuads.size(); ++iQ )
+    {
+      const Prism_3D::TQuadList& quadList = thePrism.myWallQuads[iQ];
+      StdMeshers_FaceSidePtr      botSide = quadList.front()->side[ QUAD_BOTTOM_SIDE ];
+      StdMeshers_FaceSidePtr      topSide = quadList.back ()->side[ QUAD_TOP_SIDE ];
+      if ( botSide->NbEdges() == topSide->NbEdges() )
+      {
+        for ( int iE = 0; iE < botSide->NbEdges(); ++iE )
+        {
+          TAssocTool::InsertAssociation( botSide->Edge( iE ),
+                                         topSide->Edge( iE ), shape2ShapeMap );
+          TAssocTool::InsertAssociation( myHelper->IthVertex( 0, botSide->Edge( iE )),
+                                         myHelper->IthVertex( 0, topSide->Edge( iE )),
+                                         shape2ShapeMap );
+        }
+      }
+      else
+      {
+        TopoDS_Vertex vb, vt;
+        StdMeshers_FaceSidePtr sideB, sideT;
+        vb = myHelper->IthVertex( 0, botSide->Edge( 0 ));
+        vt = myHelper->IthVertex( 0, topSide->Edge( 0 ));
+        sideB = quadList.front()->side[ QUAD_LEFT_SIDE ];
+        sideT = quadList.back ()->side[ QUAD_LEFT_SIDE ];
+        if ( vb.IsSame( sideB->FirstVertex() ) &&
+             vt.IsSame( sideT->LastVertex() ))
+        {
+          TAssocTool::InsertAssociation( botSide->Edge( 0 ),
+                                         topSide->Edge( 0 ), shape2ShapeMap );
+          TAssocTool::InsertAssociation( vb, vt, shape2ShapeMap );
+        }
+        vb = myHelper->IthVertex( 1, botSide->Edge( botSide->NbEdges()-1 ));
+        vt = myHelper->IthVertex( 1, topSide->Edge( topSide->NbEdges()-1 ));
+        sideB = quadList.front()->side[ QUAD_RIGHT_SIDE ];
+        sideT = quadList.back ()->side[ QUAD_RIGHT_SIDE ];
+        if ( vb.IsSame( sideB->FirstVertex() ) &&
+             vt.IsSame( sideT->LastVertex() ))
+        {
+          TAssocTool::InsertAssociation( botSide->Edge( botSide->NbEdges()-1 ),
+                                         topSide->Edge( topSide->NbEdges()-1 ),
+                                         shape2ShapeMap );
+          TAssocTool::InsertAssociation( vb, vt, shape2ShapeMap );
+        }
+      }
+    }
 
   // Find matching nodes of top and bottom faces
   TNodeNodeMap n2nMap;
   if ( ! TAssocTool::FindMatchingNodesOnFaces( botFace, myBlock.Mesh(),
                                                topFace, myBlock.Mesh(),
                                                shape2ShapeMap, n2nMap ))
-    return toSM( error(TCom("Mesh on faces #") << botSM->GetId()
-                       <<" and #"<< topSM->GetId() << " seems different" ));
+  {
+    if ( sameTopo )
+      return toSM( error(TCom("Mesh on faces #") << botSM->GetId()
+                         <<" and #"<< topSM->GetId() << " seems different" ));
+    else
+      return toSM( error(TCom("Topology of faces #") << botSM->GetId()
+                         <<" and #"<< topSM->GetId() << " seems different" ));
+  }
 
   // Fill myBotToColumnMap
 
@@ -2408,6 +2558,29 @@ namespace Prism_3D
     myWallQuads.clear();
   }
 
+  //================================================================================
+  /*!
+   * \brief Set upside-down
+   */
+  //================================================================================
+
+  void TPrismTopo::SetUpsideDown()
+  {
+    std::swap( myBottom, myTop );
+    myBottomEdges.clear();
+    std::reverse( myBottomEdges.begin(), myBottomEdges.end() );
+    for ( size_t i = 0; i < myWallQuads.size(); ++i )
+    {
+      myWallQuads[i].reverse();
+      TQuadList::iterator q = myWallQuads[i].begin();
+      for ( ; q != myWallQuads[i].end(); ++q )
+      {
+        (*q)->shift( 2, /*keepUnitOri=*/true );
+      }
+      myBottomEdges.push_back( myWallQuads[i].front()->side[ QUAD_BOTTOM_SIDE ].grid->Edge(0) );
+    }
+  }
+
 } // namespace Prism_3D
 
 //================================================================================
@@ -2514,7 +2687,7 @@ bool StdMeshers_Prism_3D::initPrism(Prism_3D::TPrismTopo& thePrism,
   SMESH_subMesh * botSM = 0;
   SMESH_subMesh * topSM = 0;
 
-  if ( hasNotQuad ) // can chose a bottom FACE
+  if ( hasNotQuad ) // can choose a bottom FACE
   {
     if ( nbNotQuadMeshed > 0 ) botSM = notQuadElemSubMesh.front();
     else                       botSM = notQuadGeomSubMesh.front();
@@ -2622,6 +2795,12 @@ bool StdMeshers_Prism_3D::initPrism(Prism_3D::TPrismTopo& thePrism,
                       "Non-quadrilateral faces are not opposite"));
   }
 
+  if ( thePrism.myBottomEdges.size() > thePrism.myWallQuads.size() )
+  {
+    // composite bottom sides => set thePrism upside-down
+    thePrism.SetUpsideDown();
+  }
+
   return true;
 }