1 // Copyright (C) 2007-2019 CEA/DEN, EDF R&D, OPEN CASCADE
3 // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
23 // File : StdMeshers_PolyhedronPerSolid_3D.cxx
25 // Created : Fri Oct 20 11:37:07 2006
26 // Author : Edward AGAPOV (eap)
28 #include "StdMeshers_PolyhedronPerSolid_3D.hxx"
30 #include "SMESHDS_Mesh.hxx"
31 #include "SMESH_ControlsDef.hxx"
32 #include "SMESH_Gen.hxx"
33 #include "SMESH_Mesh.hxx"
34 #include "SMESH_MeshAlgos.hxx"
35 #include "SMESH_MeshEditor.hxx"
36 #include "SMESH_MesherHelper.hxx"
37 #include "SMESH_ProxyMesh.hxx"
38 #include "SMESH_subMesh.hxx"
39 #include "StdMeshers_PolygonPerFace_2D.hxx"
40 #include "StdMeshers_Regular_1D.hxx"
41 #include "StdMeshers_ViscousLayers.hxx"
43 #include <TopExp_Explorer.hxx>
50 struct _EdgeMesher : public StdMeshers_Regular_1D
52 _EdgeMesher( int hypId, SMESH_Gen* gen )
53 : StdMeshers_Regular_1D( hypId, gen )
55 _hypType = NB_SEGMENTS;
56 _ivalue[ NB_SEGMENTS_IND ] = 1;
60 //=======================================================================
63 //=======================================================================
65 const SMDS_MeshElement* addHexa( std::vector< const SMDS_MeshElement* >& faces,
66 const std::vector< int > & quantities,
67 SMESH_MesherHelper & helper )
69 const SMDS_MeshElement* newHexa = 0;
71 // check nb of nodes in faces
72 for ( size_t i = 0; i < quantities.size(); ++i )
73 if ( quantities[ i ] != 4 )
76 // look for a top face
77 const SMDS_MeshElement* topFace = 0;
78 const SMDS_MeshElement* botFace = faces[0];
79 std::vector< const SMDS_MeshNode* > nodes( 16 ); // last 8 is a working buffer
80 nodes.assign( botFace->begin_nodes(), botFace->end_nodes() );
81 for ( size_t iF = 1; iF < faces.size() && !topFace; ++iF )
83 bool hasCommonNode = false;
84 for ( int iN = 0; iN < quantities[ 0 ] && !hasCommonNode; ++iN )
85 hasCommonNode = ( faces[ iF ]->GetNodeIndex( nodes[ iN ]) >= 0 );
88 topFace = faces[ iF ];
91 nodes.resize( 8 ); // set top nodes after hexa nodes - [8-11]
92 nodes.insert( nodes.end(), topFace->begin_nodes(), topFace->end_nodes() );
94 nodes.insert( nodes.end(), nodes.begin() + 8, nodes.begin() + 12 );
96 // find corresponding nodes of top and bottom by finding a side face including 2 node of each
97 SMESHDS_Mesh* mesh = helper.GetMeshDS();
98 const SMDS_MeshElement* sideFace = 0;
100 for ( i = 8; i < nodes.size()-1 && !sideFace; ++i )
102 sideFace = mesh->FindFace( nodes[0], nodes[1], nodes[ i ], nodes[ i + 1 ]);
107 --i; // restore after ++i in the loop
108 bool botOriRight = SMESH_MeshAlgos::IsRightOrder( sideFace, nodes[ 0 ], nodes[ 1 ] );
109 bool topOriRight = SMESH_MeshAlgos::IsRightOrder( sideFace, nodes[ i ], nodes[ i + 1 ] );
110 if ( botOriRight == topOriRight )
112 nodes[ 4 ] = nodes[ i + 1 ];
113 nodes[ 5 ] = nodes[ i + 0 ];
114 nodes[ 6 ] = nodes[ i + 3 ];
115 nodes[ 7 ] = nodes[ i + 2 ];
119 nodes[ 4 ] = nodes[ i + 0 ];
120 nodes[ 5 ] = nodes[ i + 1 ];
121 nodes[ 6 ] = nodes[ i + 2 ];
122 nodes[ 7 ] = nodes[ i + 3 ];
125 newHexa = helper.AddVolume( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ], nodes[ 3 ],
126 nodes[ 4 ], nodes[ 5 ], nodes[ 6 ], nodes[ 7 ]);
131 //=======================================================================
132 //function : addTetra
134 //=======================================================================
136 const SMDS_MeshElement* addTetra( std::vector< const SMDS_MeshElement* >& faces,
137 const std::vector< int > & quantities,
138 SMESH_MesherHelper & helper )
140 const SMDS_MeshElement* newTetra = 0;
142 // check nb of nodes in faces
143 for ( size_t i = 0; i < quantities.size(); ++i )
144 if ( quantities[ i ] != 3 )
147 const SMDS_MeshElement* botFace = faces[0];
149 std::vector< const SMDS_MeshNode* > nodes( 6 );
150 nodes.assign( botFace->begin_nodes(), botFace->end_nodes() );
153 const SMDS_MeshNode* topNode = 0;
154 for ( size_t i = 0; i < 3 && !topNode; ++i )
156 topNode = faces[ 1 ]->GetNode( i );
157 if ( botFace->GetNodeIndex( topNode ) >= 0 )
163 nodes.push_back( topNode );
165 newTetra = helper.AddVolume( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ], nodes[ 3 ]);
170 //=======================================================================
171 //function : addPenta
173 //=======================================================================
175 const SMDS_MeshElement* addPenta( std::vector< const SMDS_MeshElement* >& faces,
176 const std::vector< int > & quantities,
177 SMESH_MesherHelper & helper )
179 const SMDS_MeshElement* newPenta = 0;
181 // check nb of nodes in faces and find triangle faces
182 int trias[2] = { -1, -1 };
183 for ( size_t i = 0; i < quantities.size(); ++i )
184 if ( quantities[ i ] != 4 )
186 if ( quantities[ i ] != 3 )
188 int iTria = ( trias[0] != -1 );
189 if ( trias[ iTria ] != -1 )
193 if ( trias[1] == -1 )
196 int iSide = trias[0] + 1;
197 if ( iSide == trias[1] )
200 const SMDS_MeshElement* botFace = faces[ trias[0]];
201 const SMDS_MeshElement* topFace = faces[ trias[1]];
202 const SMDS_MeshElement* sideFace = faces[ iSide ];
203 const SMDS_MeshNode* nodes[ 6 ] = { 0,0,0,0,0,0 };
204 for ( int i = 0 ; i < 3; ++i )
206 const SMDS_MeshNode* botNode = botFace->GetNode( i );
207 if ( sideFace->GetNodeIndex( botNode ) < 0 )
210 nodes[ bool( nodes[0] )] = botNode;
212 const SMDS_MeshNode* topNode = topFace->GetNode( i );
213 if ( sideFace->GetNodeIndex( topNode ) < 0 )
216 nodes[ 3 + bool( nodes[3]) ] = topNode;
218 bool botOriRight = SMESH_MeshAlgos::IsRightOrder( sideFace, nodes[ 0 ], nodes[ 1 ]);
219 bool topOriRight = SMESH_MeshAlgos::IsRightOrder( sideFace, nodes[ 3 ], nodes[ 4 ]);
220 if ( botOriRight == topOriRight )
221 std::swap( nodes[ 3 ], nodes[ 4 ]);
223 newPenta = helper.AddVolume( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ],
224 nodes[ 3 ], nodes[ 4 ], nodes[ 5 ]);
229 //=======================================================================
232 //=======================================================================
234 const SMDS_MeshElement* addPyra( std::vector< const SMDS_MeshElement* >& faces,
235 const std::vector< int > & quantities,
236 SMESH_MesherHelper & helper )
238 const SMDS_MeshElement* newPyra = 0;
240 // check nb of nodes in faces
242 for ( size_t i = 0; i < quantities.size(); ++i )
243 if ( quantities[ i ] != 3 )
245 if ( quantities[ i ] != 4 || iBot != -1 )
250 const SMDS_MeshElement* botFace = faces[ iBot ];
252 std::vector< const SMDS_MeshNode* > nodes( 8 );
253 nodes.assign( botFace->begin_nodes(), botFace->end_nodes() );
256 const SMDS_MeshNode* topNode = 0;
257 for ( size_t i = 0; i < 4 && !topNode; ++i )
259 topNode = faces[ 1 ]->GetNode( i );
260 if ( botFace->GetNodeIndex( topNode ) >= 0 )
266 nodes.push_back( topNode );
268 newPyra = helper.AddVolume( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ], nodes[ 3 ], nodes[4] );
273 //=======================================================================
274 //function : addHPrism
275 //purpose : add hexagonal prism
276 //=======================================================================
278 const SMDS_MeshElement* addHPrism( std::vector< const SMDS_MeshElement* >& faces,
279 const std::vector< int > & quantities,
280 SMESH_MesherHelper & helper )
282 const SMDS_MeshElement* newHexPrism = 0;
284 // check nb of nodes in faces and find hexagons
285 int hexa[2] = { -1, -1 };
286 for ( size_t i = 0; i < quantities.size(); ++i )
287 if ( quantities[ i ] != 4 )
289 if ( quantities[ i ] != 6 )
291 int iHex = ( hexa[0] != -1 );
292 if ( hexa[ iHex ] != -1 )
299 int iSide = hexa[0] + 1;
300 if ( iSide == hexa[1] )
303 const SMDS_MeshElement* botFace = faces[ hexa[ 0 ]];
304 const SMDS_MeshElement* topFace = faces[ hexa[ 1 ]];
305 std::vector< const SMDS_MeshNode* > nodes( 24 ); // last 12 is a working buffer
307 nodes.assign( botFace->begin_nodes(), botFace->end_nodes() );
308 nodes.resize( 12 ); // set top nodes after hexa nodes - [12-17]
309 nodes.insert( nodes.end(), topFace->begin_nodes(), topFace->end_nodes() );
311 nodes.insert( nodes.end(), nodes.begin() + 12, nodes.begin() + 18 );
313 // find corresponding nodes of top and bottom by finding a side face including 2 node of each
314 SMESHDS_Mesh* mesh = helper.GetMeshDS();
315 const SMDS_MeshElement* sideFace = 0;
317 for ( i = 12; i < nodes.size()-1 && !sideFace; ++i )
319 sideFace = mesh->FindFace( nodes[0], nodes[1], nodes[ i ], nodes[ i + 1 ]);
324 --i; // restore after ++i in the loop
325 bool botOriRight = SMESH_MeshAlgos::IsRightOrder( sideFace, nodes[ 0 ], nodes[ 1 ] );
326 bool topOriRight = SMESH_MeshAlgos::IsRightOrder( sideFace, nodes[ i ], nodes[ i + 1 ] );
327 if ( botOriRight == topOriRight )
329 nodes[ 6 ] = nodes[ i + 1 ];
330 nodes[ 7 ] = nodes[ i + 0 ];
331 nodes[ 8 ] = nodes[ i + 5 ];
332 nodes[ 9 ] = nodes[ i + 4 ];
333 nodes[ 10 ] = nodes[ i + 3 ];
334 nodes[ 11 ] = nodes[ i + 2 ];
338 nodes[ 6 ] = nodes[ i + 0 ];
339 nodes[ 7 ] = nodes[ i + 1 ];
340 nodes[ 8 ] = nodes[ i + 2 ];
341 nodes[ 9 ] = nodes[ i + 3 ];
342 nodes[ 10 ] = nodes[ i + 4 ];
343 nodes[ 11 ] = nodes[ i + 5 ];
346 newHexPrism = helper.AddVolume( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ],
347 nodes[ 3 ], nodes[ 4 ], nodes[ 5 ],
348 nodes[ 6 ], nodes[ 7 ], nodes[ 8 ],
349 nodes[ 9 ], nodes[10 ], nodes[11 ]);
354 //=======================================================================
357 //=======================================================================
359 const SMDS_MeshElement* addPoly( std::vector< const SMDS_MeshElement* >& faces,
360 const std::vector< int > & quantities,
361 SMESH_MesherHelper & helper )
363 const SMDS_MeshElement* newPoly = 0;
365 std::vector< const SMDS_MeshNode* > nodes;
366 for ( size_t iF = 0; iF < faces.size(); ++iF )
367 nodes.insert( nodes.end(), faces[iF]->begin_nodes(), faces[iF]->end_nodes() );
369 newPoly = helper.AddPolyhedralVolume( nodes, quantities );
376 //=======================================================================
377 //function : StdMeshers_PolyhedronPerSolid_3D
379 //=======================================================================
381 StdMeshers_PolyhedronPerSolid_3D::StdMeshers_PolyhedronPerSolid_3D(int hypId,
383 :SMESH_3D_Algo(hypId, gen),
384 myEdgeMesher( new _EdgeMesher( gen->GetANewId(), gen )),
385 myFaceMesher( new StdMeshers_PolygonPerFace_2D( gen->GetANewId(), gen ))
387 _name = "PolyhedronPerSolid_3D";
388 _requireDiscreteBoundary = false;
389 _supportSubmeshes = true;
390 _compatibleHypothesis.push_back("ViscousLayers");
391 _neededLowerHyps[0] = _neededLowerHyps[1] = _neededLowerHyps[2] = true;
394 //=======================================================================
395 //function : ~StdMeshers_PolyhedronPerSolid_3D
397 //=======================================================================
399 StdMeshers_PolyhedronPerSolid_3D::~StdMeshers_PolyhedronPerSolid_3D()
405 //=======================================================================
406 //function : CheckHypothesis
408 //=======================================================================
410 bool StdMeshers_PolyhedronPerSolid_3D::CheckHypothesis(SMESH_Mesh& theMesh,
411 const TopoDS_Shape& theShape,
412 Hypothesis_Status& theStatus)
414 myViscousLayersHyp = NULL;
416 const std::list<const SMESHDS_Hypothesis*>& hyps =
417 GetUsedHypothesis( theMesh, theShape, /*ignoreAuxiliary=*/false);
418 std::list <const SMESHDS_Hypothesis* >::const_iterator h = hyps.begin();
419 if ( h == hyps.end())
421 theStatus = SMESH_Hypothesis::HYP_OK;
425 // only StdMeshers_ViscousLayers can be used
427 for ( ; h != hyps.end(); ++h )
429 if ( !(myViscousLayersHyp = dynamic_cast< const StdMeshers_ViscousLayers*> ( *h )))
432 if ( !myViscousLayersHyp )
433 theStatus = HYP_INCOMPATIBLE;
435 error( myViscousLayersHyp->CheckHypothesis( theMesh, theShape, theStatus ));
437 return theStatus == HYP_OK;
440 //=======================================================================
443 //=======================================================================
445 bool StdMeshers_PolyhedronPerSolid_3D::Compute(SMESH_Mesh& theMesh,
446 const TopoDS_Shape& theShape)
448 const SMDS_MeshElement* newVolume = 0;
450 SMESH_subMesh* sm = theMesh.GetSubMesh( theShape );
451 SMESH_subMeshIteratorPtr smIt = sm->getDependsOnIterator( /*includeSelf=*/true,
452 /*complexFirst=*/false);
453 while ( smIt->more() )
456 if ( !sm->IsEmpty() )
459 const TopoDS_Shape & shape = sm->GetSubShape();
460 switch ( shape.ShapeType() )
463 sm->ComputeStateEngine( SMESH_subMesh::COMPUTE );
467 sm->ComputeStateEngine( SMESH_subMesh::COMPUTE );
469 myEdgeMesher->Compute( theMesh, shape );
473 sm->ComputeStateEngine( SMESH_subMesh::COMPUTE );
474 if ( sm->IsEmpty() && !myFaceMesher->Compute( theMesh, shape ))
476 sm->GetComputeError() = myFaceMesher->GetComputeError();
477 sm->GetComputeError()->myAlgo = myFaceMesher;
484 SMESH_MesherHelper helper( theMesh );
485 helper.SetElementsOnShape( true );
486 _quadraticMesh = helper.IsQuadraticSubMesh( shape );
488 SMESH_ProxyMesh::Ptr proxymesh( new SMESH_ProxyMesh( theMesh ));
489 if ( myViscousLayersHyp )
491 proxymesh = myViscousLayersHyp->Compute( theMesh, theShape );
496 std::vector< const SMDS_MeshElement* > faces;
499 for ( TopExp_Explorer faceEx( shape, TopAbs_FACE ); faceEx.More(); faceEx.Next() )
501 const SMESHDS_SubMesh* smDS = proxymesh->GetSubMesh( faceEx.Current() );
502 for ( SMDS_ElemIteratorPtr faceIt = smDS->GetElements(); faceIt->more(); )
503 faces.push_back( faceIt->next() );
506 bool useMediumNodes = false;
507 if ( !_quadraticMesh && theMesh.GetMeshDS()->GetMeshInfo().NbFaces( ORDER_QUADRATIC ))
508 for ( size_t i = 0; i < faces.size() && !useMediumNodes ; ++i )
509 useMediumNodes = faces[ i ]->IsQuadratic();
511 std::vector< int > quantities( faces.size() );
512 std::set< const SMDS_MeshNode* > nodes;
513 for ( size_t i = 0; i < faces.size(); ++i )
515 quantities[ i ] = useMediumNodes ? faces[ i ]->NbNodes() : faces[ i ]->NbCornerNodes();
516 for ( int iN = 0; iN < quantities[ i ]; ++iN )
517 nodes.insert( faces[ i ]->GetNode( iN ));
520 const size_t nbNodes = nodes.size(), nbFaces = faces.size();
521 if ( nbNodes == 8 && nbFaces == 6 ) newVolume = addHexa ( faces, quantities, helper );
522 else if ( nbNodes == 4 && nbFaces == 4 ) newVolume = addTetra ( faces, quantities, helper );
523 else if ( nbNodes == 6 && nbFaces == 5 ) newVolume = addPenta ( faces, quantities, helper );
524 else if ( nbNodes == 5 && nbFaces == 5 ) newVolume = addPyra ( faces, quantities, helper );
525 else if ( nbNodes == 12 && nbFaces == 8 ) newVolume = addHPrism( faces, quantities, helper );
527 newVolume = addPoly ( faces, quantities, helper );
531 SMESH::Controls::BadOrientedVolume checker;
532 checker.SetMesh( theMesh.GetMeshDS() );
533 if ( checker.IsSatisfy( newVolume->GetID() ))
535 SMESH_MeshEditor editor( &theMesh );
536 editor.Reorient( newVolume );
542 } // switch ( shape.ShapeType() )
543 } // loop on sub-meshes
548 //=======================================================================
549 //function : Evaluate
551 //=======================================================================
553 bool StdMeshers_PolyhedronPerSolid_3D::Evaluate(SMESH_Mesh& theMesh,
554 const TopoDS_Shape& theShape,
555 MapShapeNbElems& theResMap)
557 _quadraticMesh = false;
559 SMESH_subMesh* sm = theMesh.GetSubMesh( theShape );
560 SMESH_subMeshIteratorPtr smIt = sm->getDependsOnIterator( /*includeSelf=*/true,
561 /*complexFirst=*/false);
562 while ( smIt->more() )
566 MapShapeNbElems::iterator sm2vec = theResMap.find( sm );
567 if ( sm2vec != theResMap.end() && !sm2vec->second.empty() )
570 const TopoDS_Shape & shape = sm->GetSubShape();
571 switch ( shape.ShapeType() )
574 myEdgeMesher->Evaluate( theMesh, shape, theResMap );
579 myFaceMesher->Evaluate( theMesh, shape, theResMap );
580 std::vector<int> & quantities = theResMap[ sm ];
581 _quadraticMesh = ( !quantities.empty() &&
582 ( quantities[ SMDSEntity_Quad_Triangle ] +
583 quantities[ SMDSEntity_Quad_Quadrangle ] +
584 quantities[ SMDSEntity_Quad_Polygon ]));
590 std::vector<int> & quantities = theResMap[ sm ];
591 quantities.resize( SMDSEntity_Last, 0 );
593 SMESH_MesherHelper helper( theMesh );
594 const int nbNodes = helper.Count( shape, TopAbs_VERTEX, /*ignoreSame=*/true );
595 const int nbFaces = helper.Count( shape, TopAbs_FACE, /*ignoreSame=*/false );
597 if ( nbNodes == 8 && nbFaces == 6 )
598 quantities[ _quadraticMesh ? SMDSEntity_Quad_Hexa : SMDSEntity_Hexa ] = 1;
599 else if ( nbNodes == 4 && nbFaces == 4 )
600 quantities[ _quadraticMesh ? SMDSEntity_Quad_Tetra : SMDSEntity_Tetra ] = 1;
601 else if ( nbNodes == 6 && nbFaces == 5 )
602 quantities[ _quadraticMesh ? SMDSEntity_Quad_Penta : SMDSEntity_Penta ] = 1;
603 else if ( nbNodes == 5 && nbFaces == 5 )
604 quantities[ _quadraticMesh ? SMDSEntity_Quad_Pyramid : SMDSEntity_Pyramid ] = 1;
605 else if ( nbNodes == 12 && nbFaces == 8 )
606 quantities[ /*_quadraticMesh ? SMDSEntity_Quad_Pyramid :*/ SMDSEntity_Hexagonal_Prism ] = 1;
608 quantities[ /*_quadraticMesh ? SMDSEntity_Quad_Polyhedra : */SMDSEntity_Polyhedra ] = 1;
614 } // switch ( shape.ShapeType() )
615 } // loop on sub-meshes