1 // Copyright (C) 2007-2022 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 // SMESH SMESH : implementation of SMESH idl descriptions
24 // File : StdMeshers_Hexa_3D.cxx
25 // Moved here from SMESH_Hexa_3D.cxx
26 // Author : Paul RASCLE, EDF
29 #include "StdMeshers_Hexa_3D.hxx"
31 #include "SMDS_MeshNode.hxx"
32 #include "SMESH_Comment.hxx"
33 #include "SMESH_Gen.hxx"
34 #include "SMESH_Mesh.hxx"
35 #include "SMESH_MesherHelper.hxx"
36 #include "SMESH_subMesh.hxx"
37 #include "StdMeshers_BlockRenumber.hxx"
38 #include "StdMeshers_CompositeHexa_3D.hxx"
39 #include "StdMeshers_FaceSide.hxx"
40 #include "StdMeshers_HexaFromSkin_3D.hxx"
41 #include "StdMeshers_Penta_3D.hxx"
42 #include "StdMeshers_Prism_3D.hxx"
43 #include "StdMeshers_Quadrangle_2D.hxx"
44 #include "StdMeshers_ViscousLayers.hxx"
46 #include <BRep_Tool.hxx>
47 #include <Bnd_B3d.hxx>
49 #include <TopExp_Explorer.hxx>
50 #include <TopTools_MapOfShape.hxx>
51 #include <TopTools_SequenceOfShape.hxx>
54 #include "utilities.h"
55 #include "Utils_ExceptHandlers.hxx"
60 typedef SMESH_Comment TComm;
64 static SMESH_ComputeErrorPtr ComputePentahedralMesh(SMESH_Mesh &,
66 SMESH_ProxyMesh* proxyMesh=0);
68 static bool EvaluatePentahedralMesh(SMESH_Mesh &, const TopoDS_Shape &,
71 //=============================================================================
75 //=============================================================================
77 StdMeshers_Hexa_3D::StdMeshers_Hexa_3D(int hypId, SMESH_Gen * gen)
78 :SMESH_3D_Algo(hypId, gen)
81 _shapeType = (1 << TopAbs_SHELL) | (1 << TopAbs_SOLID); // 1 bit /shape type
82 _requireShape = false;
83 _compatibleHypothesis.push_back("ViscousLayers");
84 _compatibleHypothesis.push_back("BlockRenumber");
85 _quadAlgo = new StdMeshers_Quadrangle_2D( gen->GetANewId(), _gen );
88 //=============================================================================
92 //=============================================================================
94 StdMeshers_Hexa_3D::~StdMeshers_Hexa_3D()
100 //=============================================================================
102 * Retrieves defined hypotheses
104 //=============================================================================
106 bool StdMeshers_Hexa_3D::CheckHypothesis
108 const TopoDS_Shape& aShape,
109 SMESH_Hypothesis::Hypothesis_Status& aStatus)
111 // check nb of faces in the shape
113 aStatus = SMESH_Hypothesis::HYP_BAD_GEOMETRY;
115 for (TopExp_Explorer exp(aShape, TopAbs_FACE); exp.More(); exp.Next())
122 _viscousLayersHyp = nullptr;
123 _blockRenumberHyp = nullptr;
125 const list<const SMESHDS_Hypothesis*>& hyps =
126 GetUsedHypothesis(aMesh, aShape, /*ignoreAuxiliary=*/false);
127 list <const SMESHDS_Hypothesis* >::const_iterator h = hyps.begin();
128 if ( h == hyps.end())
130 aStatus = SMESH_Hypothesis::HYP_OK;
134 // only StdMeshers_ViscousLayers can be used
136 for ( ; h != hyps.end(); ++h )
138 if ( !_viscousLayersHyp &&
139 (_viscousLayersHyp = dynamic_cast< const StdMeshers_ViscousLayers*> ( *h )))
141 if ( !_blockRenumberHyp &&
142 (_blockRenumberHyp = dynamic_cast< const StdMeshers_BlockRenumber*> ( *h )))
146 if ((int) hyps.size() != (bool)_viscousLayersHyp + (bool)_blockRenumberHyp )
147 aStatus = HYP_INCOMPATIBLE;
150 if ( _viscousLayersHyp )
151 if ( !error( _viscousLayersHyp->CheckHypothesis( aMesh, aShape, aStatus )))
152 aStatus = HYP_BAD_PARAMETER;
154 if ( _blockRenumberHyp && aStatus == HYP_OK )
155 error( _blockRenumberHyp->CheckHypothesis( aMesh, aShape ));
158 return aStatus == HYP_OK;
163 //=============================================================================
165 typedef boost::shared_ptr< FaceQuadStruct > FaceQuadStructPtr;
166 typedef std::vector<gp_XYZ> TXYZColumn;
168 // symbolic names of box sides
169 enum EBoxSides{ B_BOTTOM=0, B_RIGHT, B_TOP, B_LEFT, B_FRONT, B_BACK, B_NB_SIDES };
171 // symbolic names of sides of quadrangle
172 enum EQuadSides{ Q_BOTTOM=0, Q_RIGHT, Q_TOP, Q_LEFT, Q_NB_SIDES };
174 enum EAxes{ COO_X=1, COO_Y, COO_Z };
176 //=============================================================================
178 * \brief Container of nodes of structured mesh on a qudrangular geom FACE
183 FaceQuadStructPtr _quad;
185 // map of (node parameter on EDGE) to (column (vector) of nodes)
186 TParam2ColumnMap _u2nodesMap;
188 // node column's taken from _u2nodesMap taking into account sub-shape orientation
189 vector<TNodeColumn> _columns;
191 // columns of normalized parameters of nodes within the unitary cube
192 vector<TXYZColumn> _ijkColumns;
194 // geometry of a cube side
197 const SMDS_MeshNode* GetNode(int iCol, int iRow) const
199 return _columns[iCol][iRow];
201 gp_XYZ GetXYZ(int iCol, int iRow) const
203 return SMESH_TNodeXYZ( GetNode( iCol, iRow ));
205 gp_XYZ& GetIJK(int iCol, int iRow)
207 return _ijkColumns[iCol][iRow];
211 //================================================================================
213 * \brief Converter of a pair of integers to a sole index
218 _Indexer( int xSize, int ySize ): _xSize(xSize), _ySize(ySize) {}
219 int size() const { return _xSize * _ySize; }
220 int operator()(const int x, const int y) const { return y * _xSize + x; }
223 //================================================================================
225 * \brief Finds FaceQuadStruct having a side equal to a given one and rearranges
226 * the found FaceQuadStruct::side to have the given side at a Q_BOTTOM place
228 FaceQuadStructPtr getQuadWithBottom( StdMeshers_FaceSidePtr side,
229 FaceQuadStructPtr quad[ 6 ])
231 FaceQuadStructPtr foundQuad;
232 for ( int i = 1; i < 6; ++i )
234 if ( !quad[i] ) continue;
235 for ( size_t iS = 0; iS < quad[i]->side.size(); ++iS )
237 const StdMeshers_FaceSidePtr side2 = quad[i]->side[iS];
238 if (( side->FirstVertex().IsSame( side2->FirstVertex() ) ||
239 side->FirstVertex().IsSame( side2->LastVertex() ))
241 ( side->LastVertex().IsSame( side2->FirstVertex() ) ||
242 side->LastVertex().IsSame( side2->LastVertex() ))
245 if ( iS != Q_BOTTOM )
247 vector< FaceQuadStruct::Side > newSides;
248 for ( size_t j = iS; j < quad[i]->side.size(); ++j )
249 newSides.push_back( quad[i]->side[j] );
250 for ( size_t j = 0; j < iS; ++j )
251 newSides.push_back( quad[i]->side[j] );
252 quad[i]->side.swap( newSides );
254 foundQuad.swap(quad[i]);
262 //================================================================================
264 * \brief Put quads to aCubeSide in the order of enum EBoxSides
266 //================================================================================
268 bool arrangeQuads( FaceQuadStructPtr quad[ 6 ], _FaceGrid aCubeSide[ 6 ], bool reverseBottom )
270 swap( aCubeSide[B_BOTTOM]._quad, quad[0] );
272 swap( aCubeSide[B_BOTTOM]._quad->side[ Q_RIGHT],// direct the bottom normal inside cube
273 aCubeSide[B_BOTTOM]._quad->side[ Q_LEFT ] );
275 aCubeSide[B_FRONT]._quad = getQuadWithBottom( aCubeSide[B_BOTTOM]._quad->side[Q_BOTTOM], quad );
276 aCubeSide[B_RIGHT]._quad = getQuadWithBottom( aCubeSide[B_BOTTOM]._quad->side[Q_RIGHT ], quad );
277 aCubeSide[B_BACK ]._quad = getQuadWithBottom( aCubeSide[B_BOTTOM]._quad->side[Q_TOP ], quad );
278 aCubeSide[B_LEFT ]._quad = getQuadWithBottom( aCubeSide[B_BOTTOM]._quad->side[Q_LEFT ], quad );
279 if ( aCubeSide[B_FRONT ]._quad )
280 aCubeSide[B_TOP]._quad = getQuadWithBottom( aCubeSide[B_FRONT ]._quad->side[Q_TOP ], quad );
282 for ( int i = 1; i < 6; ++i )
283 if ( !aCubeSide[i]._quad )
288 //================================================================================
290 * \brief Rearrange block sides according to StdMeshers_BlockRenumber hypothesis
292 //================================================================================
294 bool arrangeForRenumber( _FaceGrid blockSide[ 6 ],
296 TopoDS_Vertex& v001 )
298 std::swap( blockSide[B_BOTTOM]._quad->side[ Q_RIGHT],// restore after arrangeQuads()
299 blockSide[B_BOTTOM]._quad->side[ Q_LEFT ] );
302 TopTools_MapOfShape cornerVertices;
303 cornerVertices.Add( blockSide[B_BOTTOM]._quad->side[Q_BOTTOM].grid->LastVertex() );
304 cornerVertices.Add( blockSide[B_BOTTOM]._quad->side[Q_BOTTOM].grid->FirstVertex() );
305 cornerVertices.Add( blockSide[B_BOTTOM]._quad->side[Q_TOP ].grid->LastVertex() );
306 cornerVertices.Add( blockSide[B_BOTTOM]._quad->side[Q_TOP ].grid->FirstVertex() );
307 cornerVertices.Add( blockSide[B_TOP ]._quad->side[Q_BOTTOM].grid->FirstVertex() );
308 cornerVertices.Add( blockSide[B_TOP ]._quad->side[Q_BOTTOM].grid->LastVertex() );
309 cornerVertices.Add( blockSide[B_TOP ]._quad->side[Q_TOP ].grid->FirstVertex() );
310 cornerVertices.Add( blockSide[B_TOP ]._quad->side[Q_TOP ].grid->LastVertex() );
314 // block CS is not defined;
315 // renumber only if the block has an edge parallel to an axis of global CS
317 v000 = StdMeshers_RenumberHelper::GetVertex000( cornerVertices );
321 for ( auto it = cornerVertices.cbegin(); it != cornerVertices.cend(); ++it )
322 bbox.Add( BRep_Tool::Pnt( TopoDS::Vertex( *it )));
323 double tol = 1e-5 * Sqrt( bbox.SquareExtent() );
325 // get block edges starting at v000
327 std::vector< StdMeshers_FaceSidePtr > edgesAtV000;
328 std::vector< gp_Vec > edgeDir;
329 std::vector< int > iParallel; // 0 - none, 1 - X, 2 - Y, 3 - Z
330 TopTools_MapOfShape lastVertices;
331 for ( int iQ = 0; iQ < 6; ++iQ )
333 FaceQuadStructPtr quad = blockSide[iQ]._quad;
334 for ( size_t iS = 0; iS < quad->side.size() && edgesAtV000.size() < 3; ++iS )
336 StdMeshers_FaceSidePtr edge = quad->side[iS];
337 TopoDS_Vertex v1 = edge->FirstVertex(), v2 = edge->LastVertex();
338 if (( v1.IsSame( v000 ) && !lastVertices.Contains( v2 )) ||
339 ( v2.IsSame( v000 ) && !lastVertices.Contains( v1 )))
341 bool reverse = v2.IsSame( v000 );
344 lastVertices.Add( v2 );
346 edgesAtV000.push_back( edge );
348 gp_Pnt pf = BRep_Tool::Pnt( v1 );
349 gp_Pnt pl = BRep_Tool::Pnt( v2 );
350 gp_Vec vec( pf, pl );
351 edgeDir.push_back( vec );
353 iParallel.push_back( 0 );
354 if ( !v001.IsNull() )
356 if ( v001.IsSame( v2 ))
357 iParallel.back() = 3;
361 bool isStraight = true;
362 for ( int iE = 0; iE < edge->NbEdges() && isStraight; ++iE )
363 isStraight = SMESH_Algo::IsStraight( edge->Edge( iE ));
365 // is parallel to a GCS axis?
368 int nbDiff = (( Abs( vec.X() ) > tol ) +
369 ( Abs( vec.Y() ) > tol ) +
370 ( Abs( vec.Z() ) > tol ) );
372 iParallel.back() = ( Abs( vec.X() ) > tol ) ? 1 : ( Abs( vec.Y() ) > tol ) ? 2 : 3;
376 edgeDir.back() = gp_Vec( pf, edge->Value3d( reverse ? 0.99 : 0.01 ));
382 if ( std::accumulate( iParallel.begin(), iParallel.end(), 0 ) == 0 )
385 // find edge OZ and edge OX
386 StdMeshers_FaceSidePtr edgeOZ, edgeOX;
387 auto iZIt = std::find( iParallel.begin(), iParallel.end(), 3 );
388 if ( iZIt != iParallel.end() )
390 int i = std::distance( iParallel.begin(), iZIt );
391 edgeOZ = edgesAtV000[ i ];
392 int iE1 = SMESH_MesherHelper::WrapIndex( i + 1, edgesAtV000.size() );
393 int iE2 = SMESH_MesherHelper::WrapIndex( i + 2, edgesAtV000.size() );
394 if (( edgeDir[ iE1 ] ^ edgeDir[ iE2 ] ) * edgeDir[ i ] < 0 )
395 std::swap( iE1, iE2 );
396 edgeOX = edgesAtV000[ iE1 ];
400 for ( size_t i = 0; i < edgesAtV000.size(); ++i )
402 if ( !iParallel[ i ] )
404 int iE1 = SMESH_MesherHelper::WrapIndex( i + 1, edgesAtV000.size() );
405 int iE2 = SMESH_MesherHelper::WrapIndex( i + 2, edgesAtV000.size() );
406 if (( edgeDir[ iE1 ] ^ edgeDir[ iE2 ] ) * edgeDir[ i ] < 0 )
407 std::swap( iE1, iE2 );
408 edgeOZ = edgesAtV000[ iParallel[i] == 1 ? iE2 : iE1 ];
409 edgeOX = edgesAtV000[ iParallel[i] == 1 ? i : iE1 ];
414 if ( !edgeOZ || !edgeOX )
417 TopoDS_Vertex v100 = edgeOX->LastVertex();
418 if ( v100.IsSame( v000 ))
419 v100 = edgeOX->FirstVertex();
421 // Find the left quad, one including v000 but not v100
423 for ( int iQ = 0; iQ < 6; ++iQ )
425 FaceQuadStructPtr quad = blockSide[iQ]._quad;
426 bool hasV000 = false, hasV100 = false;
427 for ( size_t iS = 0; iS < quad->side.size(); ++iS )
429 StdMeshers_FaceSidePtr edge = quad->side[iS];
430 if ( edge->FirstVertex().IsSame( v000 ) || edge->LastVertex().IsSame( v000 ))
432 if ( edge->FirstVertex().IsSame( v100 ) || edge->LastVertex().IsSame( v100 ))
435 if ( hasV000 && !hasV100 )
437 // orient the left quad
438 for ( int i = 0; i < 4; ++i )
440 if ( quad->side[Q_BOTTOM].grid->Edge(0).IsSame( edgeOZ->Edge(0) ))
442 quad->shift( 1, true );
445 FaceQuadStructPtr quads[ 6 ];
446 quads[0].swap( blockSide[iQ]._quad );
447 for ( int i = 1, j = 0; i < 6; ++i, ++j )
448 if ( blockSide[ j ]._quad )
449 quads[ i ].swap( blockSide[ j ]._quad );
453 return arrangeQuads( quads, blockSide, false/* true*/ );
459 //================================================================================
461 * \brief Returns true if the 1st base node of sideGrid1 belongs to sideGrid2
463 //================================================================================
465 bool beginsAtSide( const _FaceGrid& sideGrid1,
466 const _FaceGrid& sideGrid2,
467 SMESH_ProxyMesh::Ptr proxymesh )
469 const TNodeColumn& col0 = sideGrid2._u2nodesMap.begin()->second;
470 const TNodeColumn& col1 = sideGrid2._u2nodesMap.rbegin()->second;
471 const SMDS_MeshNode* n00 = col0.front();
472 const SMDS_MeshNode* n01 = col0.back();
473 const SMDS_MeshNode* n10 = col1.front();
474 const SMDS_MeshNode* n11 = col1.back();
475 const SMDS_MeshNode* n = (sideGrid1._u2nodesMap.begin()->second)[0];
478 n00 = proxymesh->GetProxyNode( n00 );
479 n10 = proxymesh->GetProxyNode( n10 );
480 n01 = proxymesh->GetProxyNode( n01 );
481 n11 = proxymesh->GetProxyNode( n11 );
482 n = proxymesh->GetProxyNode( n );
484 return ( n == n00 || n == n01 || n == n10 || n == n11 );
487 //================================================================================
489 * \brief Fill in _FaceGrid::_ijkColumns
490 * \param [in,out] fg - a _FaceGrid
491 * \param [in] i1 - coordinate index along _columns
492 * \param [in] i2 - coordinate index along _columns[i]
493 * \param [in] v3 - value of the constant parameter
495 //================================================================================
497 void computeIJK( _FaceGrid& fg, int i1, int i2, double v3 )
499 gp_XYZ ijk( v3, v3, v3 );
500 const size_t nbCol = fg._columns.size();
501 const size_t nbRow = fg._columns[0].size();
503 fg._ijkColumns.resize( nbCol );
504 for ( size_t i = 0; i < nbCol; ++i )
505 fg._ijkColumns[ i ].resize( nbRow, ijk );
507 vector< double > len( nbRow );
509 for ( size_t i = 0; i < nbCol; ++i )
511 gp_Pnt pPrev = fg.GetXYZ( i, 0 );
512 for ( size_t j = 1; j < nbRow; ++j )
514 gp_Pnt p = fg.GetXYZ( i, j );
515 len[ j ] = len[ j-1 ] + p.Distance( pPrev );
518 for ( size_t j = 0; j < nbRow; ++j )
519 fg.GetIJK( i, j ).SetCoord( i2, len[ j ]/len.back() );
523 for ( size_t j = 0; j < nbRow; ++j )
525 gp_Pnt pPrev = fg.GetXYZ( 0, j );
526 for ( size_t i = 1; i < nbCol; ++i )
528 gp_Pnt p = fg.GetXYZ( i, j );
529 len[ i ] = len[ i-1 ] + p.Distance( pPrev );
532 for ( size_t i = 0; i < nbCol; ++i )
533 fg.GetIJK( i, j ).SetCoord( i1, len[ i ]/len.back() );
538 //=============================================================================
540 * Generates hexahedron mesh on hexaedron like form using algorithm from
541 * "Application de l'interpolation transfinie � la cr�ation de maillages
542 * C0 ou G1 continus sur des triangles, quadrangles, tetraedres, pentaedres
543 * et hexaedres d�form�s."
544 * Alain PERONNET - 8 janvier 1999
546 //=============================================================================
548 bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh & aMesh,
549 const TopoDS_Shape & aShape)
551 // PAL14921. Enable catching std::bad_alloc and Standard_OutOfMemory outside
552 //Unexpect aCatch(SalomeException);
553 SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
555 // Shape verification
556 // ----------------------
558 // shape must be a solid (or a shell) with 6 faces
559 TopExp_Explorer exp(aShape,TopAbs_SHELL);
561 return error(COMPERR_BAD_SHAPE, "No SHELL in the geometry");
562 if ( exp.Next(), exp.More() )
563 return error(COMPERR_BAD_SHAPE, "More than one SHELL in the geometry");
565 TopTools_IndexedMapOfShape FF, EE;
566 TopExp::MapShapes( aShape, TopAbs_FACE, FF);
567 if ( FF.Extent() != 6)
569 static StdMeshers_CompositeHexa_3D compositeHexa(_gen->GetANewId(), _gen);
570 compositeHexa.SetHypothesis( _blockRenumberHyp );
571 if ( !compositeHexa.Compute( aMesh, aShape ))
572 return error( compositeHexa.GetComputeError() );
574 return _blockRenumberHyp ? error( _blockRenumberHyp->CheckHypothesis( aMesh, aShape )) : true;
577 // Find sides of a cube
578 // ---------------------
580 // tool creating quadratic elements if needed
581 SMESH_MesherHelper helper (aMesh);
582 _quadraticMesh = helper.IsQuadraticSubMesh(aShape);
584 TopExp::MapShapes( aShape, TopAbs_EDGE, EE );
585 SMESH_MesherHelper* faceHelper = ( EE.Size() == 12 ) ? 0 : &helper;
587 FaceQuadStructPtr quad[ 6 ];
588 for ( int i = 0; i < 6; ++i )
591 faceHelper->SetSubShape( FF( i+1 ));
592 if ( !( quad[i] = FaceQuadStructPtr( _quadAlgo->CheckNbEdges( aMesh, FF( i+1 ),
593 /*considerMesh=*/true,
595 return error( _quadAlgo->GetComputeError() );
596 if ( quad[i]->side.size() != 4 )
597 return error( COMPERR_BAD_SHAPE, "Not a quadrangular box side" );
600 // put quads in a proper order
601 _FaceGrid aCubeSide[ 6 ];
602 if ( !arrangeQuads( quad, aCubeSide, true ))
603 return error( COMPERR_BAD_SHAPE );
606 // Make viscous layers
607 // --------------------
609 SMESH_ProxyMesh::Ptr proxymesh;
610 if ( _viscousLayersHyp )
612 proxymesh = _viscousLayersHyp->Compute( aMesh, aShape, /*makeN2NMap=*/ true );
617 // Check if there are triangles on cube sides
618 // -------------------------------------------
620 if ( aMesh.NbTriangles() > 0 )
622 for ( int i = 0; i < 6; ++i )
624 const TopoDS_Face& sideF = aCubeSide[i]._quad->face;
625 const SMESHDS_SubMesh* smDS =
626 proxymesh ? proxymesh->GetSubMesh( sideF ) : meshDS->MeshElements( sideF );
627 if ( !SMESH_MesherHelper::IsSameElemGeometry( smDS, SMDSGeom_QUADRANGLE,
628 /*nullSubMeshRes=*/false ))
630 SMESH_ComputeErrorPtr err = ComputePentahedralMesh(aMesh, aShape, proxymesh.get());
636 // Arrange sides according to _blockRenumberHyp
637 bool toRenumber = _blockRenumberHyp;
640 TopoDS_Vertex v000, v001;
641 _blockRenumberHyp->IsSolidIncluded( aMesh, aShape, v000, v001 );
643 toRenumber = arrangeForRenumber( aCubeSide, v000, v001 );
648 meshDS->CompactMesh(); // remove numbering holes
652 // Check presence of regular grid mesh on FACEs of the cube
653 // ------------------------------------------------------------
655 for ( int i = 0; i < 6; ++i )
657 const TopoDS_Face& F = aCubeSide[i]._quad->face;
658 StdMeshers_FaceSidePtr baseQuadSide = aCubeSide[i]._quad->side[ Q_BOTTOM ];
659 list<TopoDS_Edge> baseEdges( baseQuadSide->Edges().begin(), baseQuadSide->Edges().end() );
661 // assure correctness of node positions on baseE:
662 // helper.GetNodeU() will fix positions if they are wrong
663 helper.ToFixNodeParameters( true );
664 for ( int iE = 0; iE < baseQuadSide->NbEdges(); ++iE )
666 const TopoDS_Edge& baseE = baseQuadSide->Edge( iE );
667 if ( SMESHDS_SubMesh* smDS = meshDS->MeshElements( baseE ))
670 helper.SetSubShape( baseE );
671 SMDS_ElemIteratorPtr eIt = smDS->GetElements();
672 while ( eIt->more() )
674 const SMDS_MeshElement* e = eIt->next();
675 // expect problems on a composite side
676 try { helper.GetNodeU( baseE, e->GetNode(0), e->GetNode(1), &ok); }
678 try { helper.GetNodeU( baseE, e->GetNode(1), e->GetNode(0), &ok); }
686 helper.LoadNodeColumns( aCubeSide[i]._u2nodesMap, F, baseEdges, meshDS, proxymesh.get());
689 // check if the loaded grid corresponds to nb of quadrangles on the FACE
690 const SMESHDS_SubMesh* faceSubMesh =
691 proxymesh ? proxymesh->GetSubMesh( F ) : meshDS->MeshElements( F );
692 const smIdType nbQuads = faceSubMesh->NbElements();
693 const int nbHor = aCubeSide[i]._u2nodesMap.size() - 1;
694 const int nbVer = aCubeSide[i]._u2nodesMap.begin()->second.size() - 1;
695 ok = ( nbQuads == nbHor * nbVer );
699 SMESH_ComputeErrorPtr err = ComputePentahedralMesh(aMesh, aShape, proxymesh.get());
704 // Orient loaded grids of cube sides along axis of the unitary cube coord system
706 isReverse[B_BOTTOM] = beginsAtSide( aCubeSide[B_BOTTOM], aCubeSide[B_RIGHT ], proxymesh );
707 isReverse[B_TOP ] = beginsAtSide( aCubeSide[B_TOP ], aCubeSide[B_RIGHT ], proxymesh );
708 isReverse[B_FRONT ] = beginsAtSide( aCubeSide[B_FRONT ], aCubeSide[B_RIGHT ], proxymesh );
709 isReverse[B_BACK ] = beginsAtSide( aCubeSide[B_BACK ], aCubeSide[B_RIGHT ], proxymesh );
710 isReverse[B_LEFT ] = beginsAtSide( aCubeSide[B_LEFT ], aCubeSide[B_BACK ], proxymesh );
711 isReverse[B_RIGHT ] = beginsAtSide( aCubeSide[B_RIGHT ], aCubeSide[B_BACK ], proxymesh );
712 for ( int i = 0; i < 6; ++i )
714 aCubeSide[i]._columns.resize( aCubeSide[i]._u2nodesMap.size() );
716 size_t iFwd = 0, iRev = aCubeSide[i]._columns.size()-1;
717 size_t* pi = isReverse[i] ? &iRev : &iFwd;
718 TParam2ColumnMap::iterator u2nn = aCubeSide[i]._u2nodesMap.begin();
719 for ( ; iFwd < aCubeSide[i]._columns.size(); --iRev, ++iFwd, ++u2nn )
720 aCubeSide[i]._columns[ *pi ].swap( u2nn->second );
722 aCubeSide[i]._u2nodesMap.clear();
726 for ( int i = 0; i < 6; ++i )
727 for ( size_t j = 0; j < aCubeSide[i]._columns.size(); ++j)
728 for ( size_t k = 0; k < aCubeSide[i]._columns[j].size(); ++k)
730 const SMDS_MeshNode* & n = aCubeSide[i]._columns[j][k];
731 n = proxymesh->GetProxyNode( n );
734 // 4) Create internal nodes of the cube
735 // -------------------------------------
737 helper.SetSubShape( aShape );
738 helper.SetElementsOnShape(true);
740 // shortcuts to sides
741 _FaceGrid* fBottom = & aCubeSide[ B_BOTTOM ];
742 _FaceGrid* fRight = & aCubeSide[ B_RIGHT ];
743 _FaceGrid* fTop = & aCubeSide[ B_TOP ];
744 _FaceGrid* fLeft = & aCubeSide[ B_LEFT ];
745 _FaceGrid* fFront = & aCubeSide[ B_FRONT ];
746 _FaceGrid* fBack = & aCubeSide[ B_BACK ];
748 // cube size measured in nb of nodes
749 size_t x, xSize = fBottom->_columns.size() , X = xSize - 1;
750 size_t y, ySize = fLeft->_columns.size() , Y = ySize - 1;
751 size_t z, zSize = fLeft->_columns[0].size(), Z = zSize - 1;
753 // check sharing of FACEs (IPAL54417)
754 if ( fFront ->_columns.size() != xSize ||
755 fBack ->_columns.size() != xSize ||
756 fTop ->_columns.size() != xSize ||
758 fRight ->_columns.size() != ySize ||
759 fTop ->_columns[0].size() != ySize ||
760 fBottom->_columns[0].size() != ySize ||
762 fRight ->_columns[0].size() != zSize ||
763 fFront ->_columns[0].size() != zSize ||
764 fBack ->_columns[0].size() != zSize )
765 return error( COMPERR_BAD_SHAPE, "Not sewed faces" );
767 // columns of internal nodes "rising" from nodes of fBottom
768 _Indexer colIndex( xSize, ySize );
769 vector< vector< const SMDS_MeshNode* > > columns( colIndex.size() );
771 // fill node columns by front and back box sides
772 for ( x = 0; x < xSize; ++x ) {
773 vector< const SMDS_MeshNode* >& column0 = columns[ colIndex( x, 0 )];
774 vector< const SMDS_MeshNode* >& column1 = columns[ colIndex( x, Y )];
775 column0.resize( zSize );
776 column1.resize( zSize );
777 for ( z = 0; z < zSize; ++z ) {
778 column0[ z ] = fFront->GetNode( x, z );
779 column1[ z ] = fBack ->GetNode( x, z );
782 // fill node columns by left and right box sides
783 for ( y = 1; y < ySize-1; ++y ) {
784 vector< const SMDS_MeshNode* >& column0 = columns[ colIndex( 0, y )];
785 vector< const SMDS_MeshNode* >& column1 = columns[ colIndex( X, y )];
786 column0.resize( zSize );
787 column1.resize( zSize );
788 for ( z = 0; z < zSize; ++z ) {
789 column0[ z ] = fLeft ->GetNode( y, z );
790 column1[ z ] = fRight->GetNode( y, z );
793 // get nodes from top and bottom box sides
794 for ( x = 1; x < xSize-1; ++x ) {
795 for ( y = 1; y < ySize-1; ++y ) {
796 vector< const SMDS_MeshNode* >& column = columns[ colIndex( x, y )];
797 column.resize( zSize );
798 column.front() = fBottom->GetNode( x, y );
799 column.back() = fTop ->GetNode( x, y );
803 // compute normalized parameters of nodes on sides (PAL23189)
804 computeIJK( *fBottom, COO_X, COO_Y, /*z=*/0. );
805 computeIJK( *fRight, COO_Y, COO_Z, /*x=*/1. );
806 computeIJK( *fTop, COO_X, COO_Y, /*z=*/1. );
807 computeIJK( *fLeft, COO_Y, COO_Z, /*x=*/0. );
808 computeIJK( *fFront, COO_X, COO_Z, /*y=*/0. );
809 computeIJK( *fBack, COO_X, COO_Z, /*y=*/1. );
811 StdMeshers_RenumberHelper renumHelper( aMesh, _blockRenumberHyp );
813 // projection points of the internal node on cube sub-shapes by which
814 // coordinates of the internal node are computed
815 vector<gp_XYZ> pointsOnShapes( SMESH_Block::ID_Shell );
817 // projections on vertices are constant
818 pointsOnShapes[ SMESH_Block::ID_V000 ] = fBottom->GetXYZ( 0, 0 );
819 pointsOnShapes[ SMESH_Block::ID_V100 ] = fBottom->GetXYZ( X, 0 );
820 pointsOnShapes[ SMESH_Block::ID_V010 ] = fBottom->GetXYZ( 0, Y );
821 pointsOnShapes[ SMESH_Block::ID_V110 ] = fBottom->GetXYZ( X, Y );
822 pointsOnShapes[ SMESH_Block::ID_V001 ] = fTop->GetXYZ( 0, 0 );
823 pointsOnShapes[ SMESH_Block::ID_V101 ] = fTop->GetXYZ( X, 0 );
824 pointsOnShapes[ SMESH_Block::ID_V011 ] = fTop->GetXYZ( 0, Y );
825 pointsOnShapes[ SMESH_Block::ID_V111 ] = fTop->GetXYZ( X, Y );
827 gp_XYZ params; // normalized parameters of an internal node within the unit box
830 for ( y = 0; y < ySize; ++y )
832 vector< const SMDS_MeshNode* >& column0y = columns[ colIndex( 0, y )];
833 for ( z = 0; z < zSize; ++z )
834 renumHelper.AddReplacingNode( column0y[ z ] );
837 for ( x = 1; x < xSize-1; ++x )
841 vector< const SMDS_MeshNode* >& columnX0 = columns[ colIndex( x, 0 )];
842 for ( z = 0; z < zSize; ++z )
843 renumHelper.AddReplacingNode( columnX0[ z ] );
846 const double rX = x / double(X);
847 for ( y = 1; y < ySize-1; ++y )
849 const double rY = y / double(Y);
851 // a column to fill in during z loop
852 vector< const SMDS_MeshNode* >& column = columns[ colIndex( x, y )];
853 // projection points on horizontal edges
854 pointsOnShapes[ SMESH_Block::ID_Ex00 ] = fBottom->GetXYZ( x, 0 );
855 pointsOnShapes[ SMESH_Block::ID_Ex10 ] = fBottom->GetXYZ( x, Y );
856 pointsOnShapes[ SMESH_Block::ID_E0y0 ] = fBottom->GetXYZ( 0, y );
857 pointsOnShapes[ SMESH_Block::ID_E1y0 ] = fBottom->GetXYZ( X, y );
858 pointsOnShapes[ SMESH_Block::ID_Ex01 ] = fTop->GetXYZ( x, 0 );
859 pointsOnShapes[ SMESH_Block::ID_Ex11 ] = fTop->GetXYZ( x, Y );
860 pointsOnShapes[ SMESH_Block::ID_E0y1 ] = fTop->GetXYZ( 0, y );
861 pointsOnShapes[ SMESH_Block::ID_E1y1 ] = fTop->GetXYZ( X, y );
862 // projection points on horizontal faces
863 pointsOnShapes[ SMESH_Block::ID_Fxy0 ] = fBottom->GetXYZ( x, y );
864 pointsOnShapes[ SMESH_Block::ID_Fxy1 ] = fTop ->GetXYZ( x, y );
867 renumHelper.AddReplacingNode( column[ 0 ] );
869 for ( z = 1; z < zSize-1; ++z ) // z loop
871 const double rZ = z / double(Z);
873 const gp_XYZ& pBo = fBottom->GetIJK( x, y );
874 const gp_XYZ& pTo = fTop ->GetIJK( x, y );
875 const gp_XYZ& pFr = fFront ->GetIJK( x, z );
876 const gp_XYZ& pBa = fBack ->GetIJK( x, z );
877 const gp_XYZ& pLe = fLeft ->GetIJK( y, z );
878 const gp_XYZ& pRi = fRight ->GetIJK( y, z );
879 params.SetCoord( 1, 0.5 * ( pBo.X() * ( 1. - rZ ) + pTo.X() * rZ +
880 pFr.X() * ( 1. - rY ) + pBa.X() * rY ));
881 params.SetCoord( 2, 0.5 * ( pBo.Y() * ( 1. - rZ ) + pTo.Y() * rZ +
882 pLe.Y() * ( 1. - rX ) + pRi.Y() * rX ));
883 params.SetCoord( 3, 0.5 * ( pFr.Z() * ( 1. - rY ) + pBa.Z() * rY +
884 pLe.Z() * ( 1. - rX ) + pRi.Z() * rX ));
886 // projection points on vertical edges
887 pointsOnShapes[ SMESH_Block::ID_E00z ] = fFront->GetXYZ( 0, z );
888 pointsOnShapes[ SMESH_Block::ID_E10z ] = fFront->GetXYZ( X, z );
889 pointsOnShapes[ SMESH_Block::ID_E01z ] = fBack->GetXYZ( 0, z );
890 pointsOnShapes[ SMESH_Block::ID_E11z ] = fBack->GetXYZ( X, z );
891 // projection points on vertical faces
892 pointsOnShapes[ SMESH_Block::ID_Fx0z ] = fFront->GetXYZ( x, z );
893 pointsOnShapes[ SMESH_Block::ID_Fx1z ] = fBack ->GetXYZ( x, z );
894 pointsOnShapes[ SMESH_Block::ID_F0yz ] = fLeft ->GetXYZ( y, z );
895 pointsOnShapes[ SMESH_Block::ID_F1yz ] = fRight->GetXYZ( y, z );
897 // compute internal node coordinates
899 SMESH_Block::ShellPoint( params, pointsOnShapes, coords );
900 column[ z ] = helper.AddNode( coords.X(), coords.Y(), coords.Z() );
904 renumHelper.AddReplacingNode( column[ Z ] );
909 vector< const SMDS_MeshNode* >& columnX0 = columns[ colIndex( x, Y )];
910 for ( z = 0; z < zSize; ++z )
911 renumHelper.AddReplacingNode( columnX0[ z ] );
916 for ( y = 0; y < ySize; ++y )
918 vector< const SMDS_MeshNode* >& columnXy = columns[ colIndex( X, y )];
919 for ( z = 0; z < zSize; ++z )
920 renumHelper.AddReplacingNode( columnXy[ z ] );
923 // side data no more needed, free memory
924 for ( int i = 0; i < 6; ++i )
925 SMESHUtils::FreeVector( aCubeSide[i]._columns );
927 // 5) Create hexahedrons
928 // ---------------------
930 for ( x = 0; x < xSize-1; ++x ) {
931 for ( y = 0; y < ySize-1; ++y ) {
932 vector< const SMDS_MeshNode* >& col00 = columns[ colIndex( x, y )];
933 vector< const SMDS_MeshNode* >& col10 = columns[ colIndex( x+1, y )];
934 vector< const SMDS_MeshNode* >& col01 = columns[ colIndex( x, y+1 )];
935 vector< const SMDS_MeshNode* >& col11 = columns[ colIndex( x+1, y+1 )];
936 for ( z = 0; z < zSize-1; ++z )
938 // bottom face normal of a hexa mush point outside the volume
940 helper.AddVolume(col00[z], col01[z], col01[z+1], col00[z+1],
941 col10[z], col11[z], col11[z+1], col10[z+1]);
943 helper.AddVolume(col00[z], col01[z], col11[z], col10[z],
944 col00[z+1], col01[z+1], col11[z+1], col10[z+1]);
950 renumHelper.DoReplaceNodes();
953 if ( _blockRenumberHyp )
955 return error( _blockRenumberHyp->CheckHypothesis( aMesh, aShape ));
961 //=============================================================================
965 //=============================================================================
967 bool StdMeshers_Hexa_3D::Evaluate(SMESH_Mesh & aMesh,
968 const TopoDS_Shape & aShape,
969 MapShapeNbElems& aResMap)
971 vector < SMESH_subMesh * >meshFaces;
972 TopTools_SequenceOfShape aFaces;
973 for (TopExp_Explorer exp(aShape, TopAbs_FACE); exp.More(); exp.Next()) {
974 aFaces.Append(exp.Current());
975 SMESH_subMesh *aSubMesh = aMesh.GetSubMeshContaining(exp.Current());
977 meshFaces.push_back(aSubMesh);
979 if (meshFaces.size() != 6) {
980 //return error(COMPERR_BAD_SHAPE, TComm(meshFaces.size())<<" instead of 6 faces in a block");
981 static StdMeshers_CompositeHexa_3D compositeHexa(-10, aMesh.GetGen());
982 return compositeHexa.Evaluate(aMesh, aShape, aResMap);
987 //TopoDS_Shape aFace = meshFaces[i]->GetSubShape();
988 TopoDS_Shape aFace = aFaces.Value(i+1);
989 SMESH_Algo *algo = _gen->GetAlgo(aMesh, aFace);
991 std::vector<smIdType> aResVec(SMDSEntity_Last);
992 for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
993 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
994 aResMap.insert(std::make_pair(sm,aResVec));
995 SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
996 smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
999 string algoName = algo->GetName();
1000 bool isAllQuad = false;
1001 if (algoName == "Quadrangle_2D") {
1002 MapShapeNbElemsItr anIt = aResMap.find(meshFaces[i]);
1003 if( anIt == aResMap.end() ) continue;
1004 std::vector<smIdType> aVec = (*anIt).second;
1005 smIdType nbtri = std::max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
1009 if ( ! isAllQuad ) {
1010 return EvaluatePentahedralMesh(aMesh, aShape, aResMap);
1014 // find number of 1d elems for 1 face
1016 TopTools_MapOfShape Edges1;
1017 bool IsQuadratic = false;
1018 bool IsFirst = true;
1019 for (TopExp_Explorer exp(aFaces.Value(1), TopAbs_EDGE); exp.More(); exp.Next()) {
1020 Edges1.Add(exp.Current());
1021 SMESH_subMesh *sm = aMesh.GetSubMesh(exp.Current());
1023 MapShapeNbElemsItr anIt = aResMap.find(sm);
1024 if( anIt == aResMap.end() ) continue;
1025 std::vector<smIdType> aVec = (*anIt).second;
1026 nb1d += std::max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
1028 IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
1033 // find face opposite to 1 face
1035 for(i=2; i<=6; i++) {
1036 bool IsOpposite = true;
1037 for(TopExp_Explorer exp(aFaces.Value(i), TopAbs_EDGE); exp.More(); exp.Next()) {
1038 if( Edges1.Contains(exp.Current()) ) {
1048 // find number of 2d elems on side faces
1050 for(i=2; i<=6; i++) {
1051 if( i == OppNum ) continue;
1052 MapShapeNbElemsItr anIt = aResMap.find( meshFaces[i-1] );
1053 if( anIt == aResMap.end() ) continue;
1054 std::vector<smIdType> aVec = (*anIt).second;
1055 nb2d += std::max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
1058 MapShapeNbElemsItr anIt = aResMap.find( meshFaces[0] );
1059 std::vector<smIdType> aVec = (*anIt).second;
1060 smIdType nb2d_face0 = std::max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
1061 smIdType nb0d_face0 = aVec[SMDSEntity_Node];
1063 std::vector<smIdType> aResVec(SMDSEntity_Last);
1064 for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
1066 aResVec[SMDSEntity_Quad_Hexa] = nb2d_face0 * ( nb2d/nb1d );
1067 smIdType nb1d_face0_int = ( nb2d_face0*4 - nb1d ) / 2;
1068 aResVec[SMDSEntity_Node] = nb0d_face0 * ( 2*nb2d/nb1d - 1 ) - nb1d_face0_int * nb2d/nb1d;
1071 aResVec[SMDSEntity_Node] = nb0d_face0 * ( nb2d/nb1d - 1 );
1072 aResVec[SMDSEntity_Hexa] = nb2d_face0 * ( nb2d/nb1d );
1074 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
1075 aResMap.insert(std::make_pair(sm,aResVec));
1080 //================================================================================
1082 * \brief Computes hexahedral mesh from 2D mesh of block
1084 //================================================================================
1086 bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelper)
1088 static StdMeshers_HexaFromSkin_3D * algo = 0;
1090 SMESH_Gen* gen = aMesh.GetGen();
1091 algo = new StdMeshers_HexaFromSkin_3D( gen->GetANewId(), gen );
1093 algo->InitComputeError();
1094 algo->Compute( aMesh, aHelper );
1095 return error( algo->GetComputeError());
1098 //================================================================================
1100 * \brief Return true if the algorithm can mesh this shape
1101 * \param [in] aShape - shape to check
1102 * \param [in] toCheckAll - if true, this check returns OK if all shapes are OK,
1103 * else, returns OK if at least one shape is OK
1105 //================================================================================
1107 bool StdMeshers_Hexa_3D::IsApplicable( const TopoDS_Shape & aShape, bool toCheckAll )
1109 TopExp_Explorer exp0( aShape, TopAbs_SOLID );
1110 if ( !exp0.More() ) return false;
1112 for ( ; exp0.More(); exp0.Next() )
1114 int nbFoundShells = 0;
1115 TopExp_Explorer exp1( exp0.Current(), TopAbs_SHELL );
1116 for ( ; exp1.More(); exp1.Next(), ++nbFoundShells)
1117 if ( nbFoundShells == 2 ) break;
1118 if ( nbFoundShells != 1 ) {
1119 if ( toCheckAll ) return false;
1122 exp1.Init( exp0.Current(), TopAbs_FACE );
1123 int nbEdges = SMESH_MesherHelper::Count( exp1.Current(), TopAbs_EDGE, /*ignoreSame=*/true );
1124 bool ok = ( nbEdges > 3 );
1125 if ( toCheckAll && !ok ) return false;
1126 if ( !toCheckAll && ok ) return true;
1131 //=======================================================================
1132 //function : ComputePentahedralMesh
1134 //=======================================================================
1136 SMESH_ComputeErrorPtr ComputePentahedralMesh(SMESH_Mesh & aMesh,
1137 const TopoDS_Shape & aShape,
1138 SMESH_ProxyMesh* proxyMesh)
1140 SMESH_ComputeErrorPtr err = SMESH_ComputeError::New();
1143 StdMeshers_Penta_3D anAlgo;
1145 bOK=anAlgo.Compute(aMesh, aShape);
1147 err = anAlgo.GetComputeError();
1149 if ( !bOK && anAlgo.ErrorStatus() == 5 )
1151 static StdMeshers_Prism_3D * aPrism3D = 0;
1153 SMESH_Gen* gen = aMesh.GetGen();
1154 aPrism3D = new StdMeshers_Prism_3D( gen->GetANewId(), gen );
1156 SMESH_Hypothesis::Hypothesis_Status aStatus;
1157 if ( aPrism3D->CheckHypothesis( aMesh, aShape, aStatus ) ) {
1158 aPrism3D->InitComputeError();
1159 bOK = aPrism3D->Compute( aMesh, aShape );
1160 err = aPrism3D->GetComputeError();
1163 if ( !bOK && proxyMesh )
1165 // check if VL elements are present on block FACEs
1166 bool hasVLonFace = false;
1167 for ( TopExp_Explorer exp( aShape, TopAbs_FACE ); exp.More(); exp.Next() )
1169 const SMESHDS_SubMesh* sm1 = aMesh.GetSubMesh( exp.Current() )->GetSubMeshDS();
1170 const SMESHDS_SubMesh* sm2 = proxyMesh->GetSubMesh( exp.Current() );
1171 if (( hasVLonFace = ( sm2 && sm1->NbElements() != sm2->NbElements() )))
1176 err->myName = COMPERR_BAD_INPUT_MESH;
1177 err->myComment = "Can't build pentahedral mesh on viscous layers";
1185 //=======================================================================
1186 //function : EvaluatePentahedralMesh
1188 //=======================================================================
1190 bool EvaluatePentahedralMesh(SMESH_Mesh & aMesh,
1191 const TopoDS_Shape & aShape,
1192 MapShapeNbElems& aResMap)
1194 StdMeshers_Penta_3D anAlgo;
1195 bool bOK = anAlgo.Evaluate(aMesh, aShape, aResMap);
1197 //err = anAlgo.GetComputeError();
1198 //if ( !bOK && anAlgo.ErrorStatus() == 5 )
1200 static StdMeshers_Prism_3D * aPrism3D = 0;
1202 SMESH_Gen* gen = aMesh.GetGen();
1203 aPrism3D = new StdMeshers_Prism_3D( gen->GetANewId(), gen );
1205 SMESH_Hypothesis::Hypothesis_Status aStatus;
1206 if ( aPrism3D->CheckHypothesis( aMesh, aShape, aStatus ) ) {
1207 return aPrism3D->Evaluate(aMesh, aShape, aResMap);