1 // Copyright (C) 2007-2008 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.
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
22 // SMESH SMESH : implementaion of SMESH idl descriptions
23 // File : StdMeshers_RadialQuadrangle_1D2D.cxx
25 // Created : Fri Oct 20 11:37:07 2006
26 // Author : Edward AGAPOV (eap)
28 #include "StdMeshers_RadialQuadrangle_1D2D.hxx"
30 #include "StdMeshers_NumberOfLayers.hxx"
31 #include "StdMeshers_LayerDistribution.hxx"
32 #include "StdMeshers_Regular_1D.hxx"
33 #include "StdMeshers_NumberOfSegments.hxx"
35 #include "SMDS_MeshNode.hxx"
36 #include "SMESHDS_SubMesh.hxx"
37 #include "SMESH_Gen.hxx"
38 #include "SMESH_HypoFilter.hxx"
39 #include "SMESH_Mesh.hxx"
40 #include "SMESH_MesherHelper.hxx"
41 #include "SMESH_subMesh.hxx"
42 #include "SMESH_subMeshEventListener.hxx"
44 #include "utilities.h"
46 #include <BRepAdaptor_Curve.hxx>
47 #include <BRepBuilderAPI_MakeEdge.hxx>
48 #include <BRep_Tool.hxx>
49 #include <GeomAPI_ProjectPointOnSurf.hxx>
50 #include <Geom_Circle.hxx>
51 #include <Geom_Line.hxx>
52 #include <Geom_TrimmedCurve.hxx>
53 #include <TColgp_SequenceOfPnt.hxx>
54 #include <TColgp_SequenceOfPnt2d.hxx>
55 #include <TopExp_Explorer.hxx>
56 #include <TopTools_ListIteratorOfListOfShape.hxx>
62 #define RETURN_BAD_RESULT(msg) { MESSAGE(")-: Error: " << msg); return false; }
63 #define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
66 //=======================================================================
67 //function : StdMeshers_RadialQuadrangle_1D2D
69 //=======================================================================
71 StdMeshers_RadialQuadrangle_1D2D::StdMeshers_RadialQuadrangle_1D2D(int hypId,
74 :SMESH_2D_Algo(hypId, studyId, gen)
76 _name = "RadialQuadrangle_1D2D";
77 _shapeType = (1 << TopAbs_FACE); // 1 bit per shape type
79 _compatibleHypothesis.push_back("LayerDistribution2D");
80 _compatibleHypothesis.push_back("NumberOfLayers2D");
82 myDistributionHypo = 0;
83 _requireDescretBoundary = false;
84 _supportSubmeshes = true;
88 //================================================================================
92 //================================================================================
94 StdMeshers_RadialQuadrangle_1D2D::~StdMeshers_RadialQuadrangle_1D2D()
98 //=======================================================================
99 //function : CheckHypothesis
101 //=======================================================================
103 bool StdMeshers_RadialQuadrangle_1D2D::CheckHypothesis
105 const TopoDS_Shape& aShape,
106 SMESH_Hypothesis::Hypothesis_Status& aStatus)
110 myDistributionHypo = 0;
112 list <const SMESHDS_Hypothesis * >::const_iterator itl;
114 const list <const SMESHDS_Hypothesis * >&hyps = GetUsedHypothesis(aMesh, aShape);
115 if ( hyps.size() == 0 ) {
116 aStatus = SMESH_Hypothesis::HYP_OK;
117 return true; // can work with no hypothesis
120 if ( hyps.size() > 1 ) {
121 aStatus = SMESH_Hypothesis::HYP_ALREADY_EXIST;
125 const SMESHDS_Hypothesis *theHyp = hyps.front();
127 string hypName = theHyp->GetName();
129 if (hypName == "NumberOfLayers2D") {
130 myNbLayerHypo = static_cast<const StdMeshers_NumberOfLayers *>(theHyp);
131 aStatus = SMESH_Hypothesis::HYP_OK;
134 if (hypName == "LayerDistribution2D") {
135 myDistributionHypo = static_cast<const StdMeshers_LayerDistribution *>(theHyp);
136 aStatus = SMESH_Hypothesis::HYP_OK;
139 aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
145 // ------------------------------------------------------------------------------
147 * \brief Listener used to mark edges meshed by StdMeshers_RadialQuadrangle_1D2D
149 class TEdgeMarker : public SMESH_subMeshEventListener
151 TEdgeMarker(): SMESH_subMeshEventListener(/*isDeletable=*/false) {}
153 //!< Return static listener
154 static SMESH_subMeshEventListener* getListener()
156 static TEdgeMarker theEdgeMarker;
157 return &theEdgeMarker;
159 //! Clear face sumbesh if something happens on edges
160 void ProcessEvent(const int event,
162 SMESH_subMesh* edgeSubMesh,
163 EventListenerData* data,
164 const SMESH_Hypothesis* /*hyp*/)
166 if ( data && !data->mySubMeshes.empty() && eventType == SMESH_subMesh::ALGO_EVENT)
168 ASSERT( data->mySubMeshes.front() != edgeSubMesh );
169 SMESH_subMesh* faceSubMesh = data->mySubMeshes.front();
170 faceSubMesh->ComputeStateEngine( SMESH_subMesh::CLEAN );
175 // ------------------------------------------------------------------------------
177 * \brief Mark an edge as computed by StdMeshers_RadialQuadrangle_1D2D
179 void markEdgeAsComputedByMe(const TopoDS_Edge& edge, SMESH_subMesh* faceSubMesh)
181 if ( SMESH_subMesh* edgeSM = faceSubMesh->GetFather()->GetSubMeshContaining( edge ))
183 if ( !edgeSM->GetEventListenerData( TEdgeMarker::getListener() ))
184 faceSubMesh->SetEventListener( TEdgeMarker::getListener(),
185 SMESH_subMeshEventListenerData::MakeData(faceSubMesh),
189 // ------------------------------------------------------------------------------
191 * \brief Return true if a radial edge was meshed with StdMeshers_RadialQuadrangle_1D2D with
192 * the same radial distribution
194 // bool isEdgeCompatiballyMeshed(const TopoDS_Edge& edge, SMESH_subMesh* faceSubMesh)
196 // if ( SMESH_subMesh* edgeSM = faceSubMesh->GetFather()->GetSubMeshContaining( edge ))
198 // if ( SMESH_subMeshEventListenerData* otherFaceData =
199 // edgeSM->GetEventListenerData( TEdgeMarker::getListener() ))
201 // // compare hypothesis aplied to two disk faces sharing radial edges
202 // SMESH_Mesh& mesh = *faceSubMesh->GetFather();
203 // SMESH_Algo* radialQuadAlgo = mesh.GetGen()->GetAlgo(mesh, faceSubMesh->GetSubShape() );
204 // SMESH_subMesh* otherFaceSubMesh = otherFaceData->mySubMeshes.front();
205 // list <const SMESHDS_Hypothesis *> hyps1 =
206 // radialQuadAlgo->GetUsedHypothesis( mesh, faceSubMesh->GetSubShape());
207 // list <const SMESHDS_Hypothesis *> hyps2 =
208 // radialQuadAlgo->GetUsedHypothesis( mesh, otherFaceSubMesh->GetSubShape());
209 // if( hyps1.empty() && hyps2.empty() )
210 // return true; // defaul hyps
211 // if ( hyps1.size() != hyps2.size() )
213 // return *hyps1.front() == *hyps2.front();
219 //================================================================================
221 * \brief Return base curve of the edge and extremum parameters
223 //================================================================================
225 Handle(Geom_Curve) getCurve(const TopoDS_Edge& edge, double* f=0, double* l=0)
227 Handle(Geom_Curve) C;
228 if ( !edge.IsNull() )
230 double first = 0., last = 0.;
231 C = BRep_Tool::Curve(edge, first, last);
234 Handle(Geom_TrimmedCurve) tc = Handle(Geom_TrimmedCurve)::DownCast(C);
235 while( !tc.IsNull() ) {
236 C = tc->BasisCurve();
237 tc = Handle(Geom_TrimmedCurve)::DownCast(C);
246 //================================================================================
248 * \brief Return edges of the face
249 * \retval int - nb of edges
251 //================================================================================
253 int analyseFace(const TopoDS_Shape& face,
254 TopoDS_Edge& CircEdge,
255 TopoDS_Edge& LinEdge1,
256 TopoDS_Edge& LinEdge2)
258 CircEdge.Nullify(); LinEdge1.Nullify(); LinEdge2.Nullify();
261 for ( TopExp_Explorer exp( face, TopAbs_EDGE ); exp.More(); exp.Next(), ++nbe )
263 const TopoDS_Edge& E = TopoDS::Edge( exp.Current() );
265 Handle(Geom_Curve) C = getCurve(E,&f,&l);
268 if ( C->IsKind( STANDARD_TYPE(Geom_Circle)))
270 if ( CircEdge.IsNull() )
275 else if ( LinEdge1.IsNull() )
284 //================================================================================
285 //================================================================================
287 * \brief Class computing layers distribution using data of
288 * StdMeshers_LayerDistribution hypothesis
290 //================================================================================
291 //================================================================================
293 class TNodeDistributor: public StdMeshers_Regular_1D
295 list <const SMESHDS_Hypothesis *> myUsedHyps;
297 // -----------------------------------------------------------------------------
298 static TNodeDistributor* GetDistributor(SMESH_Mesh& aMesh)
300 const int myID = -1000;
301 map < int, SMESH_1D_Algo * > & algoMap = aMesh.GetGen()->_map1D_Algo;
302 map < int, SMESH_1D_Algo * >::iterator id_algo = algoMap.find( myID );
303 if ( id_algo == algoMap.end() )
304 return new TNodeDistributor( myID, 0, aMesh.GetGen() );
305 return static_cast< TNodeDistributor* >( id_algo->second );
307 // -----------------------------------------------------------------------------
308 //! Computes distribution of nodes on a straight line ending at pIn and pOut
309 bool Compute( vector< double > & positions,
313 const SMESH_Hypothesis* hyp1d)
315 if ( !hyp1d ) return error( "Invalid LayerDistribution hypothesis");
317 double len = pIn.Distance( pOut );
318 if ( len <= DBL_MIN ) return error("Too close points of inner and outer shells");
321 myUsedHyps.push_back( hyp1d );
323 TopoDS_Edge edge = BRepBuilderAPI_MakeEdge( pIn, pOut );
324 SMESH_Hypothesis::Hypothesis_Status aStatus;
325 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, edge, aStatus ))
326 return error( "StdMeshers_Regular_1D::CheckHypothesis() failed "
327 "with LayerDistribution hypothesis");
329 BRepAdaptor_Curve C3D(edge);
330 double f = C3D.FirstParameter(), l = C3D.LastParameter();
331 list< double > params;
332 if ( !StdMeshers_Regular_1D::computeInternalParameters( aMesh, C3D, len, f, l, params, false ))
333 return error("StdMeshers_Regular_1D failed to compute layers distribution");
336 positions.reserve( params.size() );
337 for (list<double>::iterator itU = params.begin(); itU != params.end(); itU++)
338 positions.push_back( *itU / len );
341 // -----------------------------------------------------------------------------
342 //! Make mesh on an adge using assigned 1d hyp or defaut nb of segments
343 bool ComputeCircularEdge(SMESH_Mesh& aMesh,
344 const TopoDS_Edge& anEdge)
346 _gen->Compute( aMesh, anEdge);
347 SMESH_subMesh *sm = aMesh.GetSubMesh(anEdge);
348 if ( sm->GetComputeState() != SMESH_subMesh::COMPUTE_OK)
350 // find any 1d hyp assigned (there can be a hyp w/o algo)
351 myUsedHyps = SMESH_Algo::GetUsedHypothesis(aMesh, anEdge, /*ignoreAux=*/true);
352 Hypothesis_Status aStatus;
353 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, anEdge, aStatus ))
355 // no valid 1d hyp assigned, use default nb of segments
356 _hypType = NB_SEGMENTS;
357 _ivalue[ DISTR_TYPE_IND ] = StdMeshers_NumberOfSegments::DT_Regular;
358 _ivalue[ NB_SEGMENTS_IND ] = _gen->GetDefaultNbSegments();
360 return StdMeshers_Regular_1D::Compute( aMesh, anEdge );
364 // -----------------------------------------------------------------------------
365 //! Make mesh on an adge using assigned 1d hyp or defaut nb of segments
366 bool EvaluateCircularEdge(SMESH_Mesh& aMesh,
367 const TopoDS_Edge& anEdge,
368 MapShapeNbElems& aResMap)
370 _gen->Evaluate( aMesh, anEdge, aResMap );
371 if ( aResMap.count( aMesh.GetSubMesh( anEdge )))
374 // find any 1d hyp assigned
375 myUsedHyps = SMESH_Algo::GetUsedHypothesis(aMesh, anEdge, /*ignoreAux=*/true);
376 Hypothesis_Status aStatus;
377 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, anEdge, aStatus ))
379 // no valid 1d hyp assigned, use default nb of segments
380 _hypType = NB_SEGMENTS;
381 _ivalue[ DISTR_TYPE_IND ] = StdMeshers_NumberOfSegments::DT_Regular;
382 _ivalue[ NB_SEGMENTS_IND ] = _gen->GetDefaultNbSegments();
384 return StdMeshers_Regular_1D::Evaluate( aMesh, anEdge, aResMap );
387 // -----------------------------------------------------------------------------
388 TNodeDistributor( int hypId, int studyId, SMESH_Gen* gen)
389 : StdMeshers_Regular_1D( hypId, studyId, gen)
392 // -----------------------------------------------------------------------------
393 virtual const list <const SMESHDS_Hypothesis *> &
394 GetUsedHypothesis(SMESH_Mesh &, const TopoDS_Shape &, const bool)
398 // -----------------------------------------------------------------------------
402 //=======================================================================
404 * \brief Allow algo to do something after persistent restoration
405 * \param subMesh - restored submesh
407 * call markEdgeAsComputedByMe()
409 //=======================================================================
411 void StdMeshers_RadialQuadrangle_1D2D::SubmeshRestored(SMESH_subMesh* faceSubMesh)
413 if ( !faceSubMesh->IsEmpty() )
415 TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
416 analyseFace( faceSubMesh->GetSubShape(), CircEdge, LinEdge1, LinEdge2 );
417 if ( !CircEdge.IsNull() ) markEdgeAsComputedByMe( CircEdge, faceSubMesh );
418 if ( !LinEdge1.IsNull() ) markEdgeAsComputedByMe( LinEdge1, faceSubMesh );
419 if ( !LinEdge2.IsNull() ) markEdgeAsComputedByMe( LinEdge2, faceSubMesh );
423 //=======================================================================
426 //=======================================================================
428 bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh& aMesh,
429 const TopoDS_Shape& aShape)
432 SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
434 myHelper = new SMESH_MesherHelper( aMesh );
435 myHelper->IsQuadraticSubMesh( aShape );
436 // to delete helper at exit from Compute()
437 auto_ptr<SMESH_MesherHelper> helperDeleter( myHelper );
439 TNodeDistributor* algo1d = TNodeDistributor::GetDistributor(aMesh);
441 TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
442 int nbe = analyseFace( aShape, CircEdge, LinEdge1, LinEdge2 );
443 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
444 if( nbe>3 || nbe < 1 || aCirc.IsNull() )
445 return error("The face must be a full circle or a part of circle (i.e. the number of edges is less or equal to 3 and one of them is a circle curve)");
448 // points for rotation
449 TColgp_SequenceOfPnt Points;
450 // angles for rotation
451 TColStd_SequenceOfReal Angles;
452 // Nodes1 and Nodes2 - nodes along radiuses
453 // CNodes - nodes on circle edge
454 vector< const SMDS_MeshNode* > Nodes1, Nodes2, CNodes;
456 // parameters edge nodes on face
457 TColgp_SequenceOfPnt2d Pnts2d1;
460 int faceID = meshDS->ShapeToIndex(aShape);
461 TopoDS_Face F = TopoDS::Face(aShape);
462 Handle(Geom_Surface) S = BRep_Tool::Surface(F);
467 if (!algo1d->ComputeCircularEdge( aMesh, CircEdge ))
468 return error( algo1d->GetComputeError() );
469 map< double, const SMDS_MeshNode* > theNodes;
470 if ( !GetSortedNodesOnEdge(aMesh.GetMeshDS(),CircEdge,true,theNodes))
471 return error("Circular edge is incorrectly meshed");
474 map< double, const SMDS_MeshNode* >::iterator itn = theNodes.begin();
475 const SMDS_MeshNode* NF = (*itn).second;
476 CNodes.push_back( (*itn).second );
477 double fang = (*itn).first;
478 if ( itn != theNodes.end() ) {
480 for(; itn != theNodes.end(); itn++ ) {
481 CNodes.push_back( (*itn).second );
482 double ang = (*itn).first - fang;
483 if( ang>PI ) ang = ang - 2*PI;
484 if( ang<-PI ) ang = ang + 2*PI;
485 Angles.Append( ang );
488 P1 = gp_Pnt( NF->X(), NF->Y(), NF->Z() );
489 P0 = aCirc->Location();
491 if ( !computeLayerPositions(P0,P1))
494 exp.Init( CircEdge, TopAbs_VERTEX );
495 TopoDS_Vertex V1 = TopoDS::Vertex( exp.Current() );
496 gp_Pnt2d p2dV = BRep_Tool::Parameters( V1, TopoDS::Face(aShape) );
498 NC = meshDS->AddNode(P0.X(), P0.Y(), P0.Z());
499 GeomAPI_ProjectPointOnSurf PPS(P0,S);
501 PPS.Parameters(1,U0,V0);
502 meshDS->SetNodeOnFace(NC, faceID, U0, V0);
503 PC = gp_Pnt2d(U0,V0);
506 gp_Vec2d aVec2d(PC,p2dV);
507 Nodes1.resize( myLayerPositions.size()+1 );
508 Nodes2.resize( myLayerPositions.size()+1 );
510 for(; i<myLayerPositions.size(); i++) {
511 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
512 P0.Y() + aVec.Y()*myLayerPositions[i],
513 P0.Z() + aVec.Z()*myLayerPositions[i] );
515 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
518 double U = PC.X() + aVec2d.X()*myLayerPositions[i];
519 double V = PC.Y() + aVec2d.Y()*myLayerPositions[i];
520 meshDS->SetNodeOnFace( node, faceID, U, V );
521 Pnts2d1.Append(gp_Pnt2d(U,V));
523 Nodes1[Nodes1.size()-1] = NF;
524 Nodes2[Nodes1.size()-1] = NF;
526 else if(nbe==2 && LinEdge1.Orientation() != TopAbs_INTERNAL )
528 // one curve must be a half of circle and other curve must be
531 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge, &fp, &lp ));
532 if( fabs(fabs(lp-fp)-PI) > Precision::Confusion() ) {
533 // not half of circle
534 return error(COMPERR_BAD_SHAPE);
536 Handle(Geom_Line) aLine = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
537 if( aLine.IsNull() ) {
538 // other curve not line
539 return error(COMPERR_BAD_SHAPE);
542 if ( !algo1d->ComputeCircularEdge( aMesh, CircEdge ))
543 return error( algo1d->GetComputeError() );
544 map< double, const SMDS_MeshNode* > theNodes;
545 if ( !GetSortedNodesOnEdge(aMesh.GetMeshDS(),CircEdge,true,theNodes) )
546 return error("Circular edge is incorrectly meshed");
548 map< double, const SMDS_MeshNode* >::iterator itn = theNodes.begin();
550 CNodes.push_back( itn->second );
551 double fang = (*itn).first;
553 for(; itn != theNodes.end(); itn++ ) {
554 CNodes.push_back( (*itn).second );
555 double ang = (*itn).first - fang;
556 if( ang>PI ) ang = ang - 2*PI;
557 if( ang<-PI ) ang = ang + 2*PI;
558 Angles.Append( ang );
560 const SMDS_MeshNode* NF = theNodes.begin()->second;
561 const SMDS_MeshNode* NL = theNodes.rbegin()->second;
562 P1 = gp_Pnt( NF->X(), NF->Y(), NF->Z() );
563 gp_Pnt P2( NL->X(), NL->Y(), NL->Z() );
564 P0 = aCirc->Location();
566 bool linEdgeComputed;
567 if ( !computeLayerPositions(P0,P1,LinEdge1,&linEdgeComputed))
570 if ( linEdgeComputed )
572 if (!GetSortedNodesOnEdge(aMesh.GetMeshDS(),LinEdge1,true,theNodes))
573 return error("Invalid mesh on a straight edge");
575 vector< const SMDS_MeshNode* > *pNodes1 = &Nodes1, *pNodes2 = &Nodes2;
576 bool nodesFromP0ToP1 = ( theNodes.rbegin()->second == NF );
577 if ( !nodesFromP0ToP1 ) std::swap( pNodes1, pNodes2 );
579 map< double, const SMDS_MeshNode* >::reverse_iterator ritn = theNodes.rbegin();
580 itn = theNodes.begin();
581 for ( int i = Nodes1.size()-1; i > -1; ++itn, ++ritn, --i )
583 (*pNodes1)[i] = ritn->second;
584 (*pNodes2)[i] = itn->second;
585 Points.Append( gpXYZ( Nodes1[i]));
586 Pnts2d1.Append( myHelper->GetNodeUV( F, Nodes1[i]));
588 NC = const_cast<SMDS_MeshNode*>( itn->second );
589 Points.Remove( Nodes1.size() );
594 int edgeID = meshDS->ShapeToIndex(LinEdge1);
596 Handle(Geom_Curve) Crv = BRep_Tool::Curve(LinEdge1,fp,lp);
600 if( P1.Distance(Ptmp) > Precision::Confusion() )
602 // get UV points for edge
604 BRep_Tool::UVPoints( LinEdge1, TopoDS::Face(aShape), PF, PL );
605 PC = gp_Pnt2d( (PF.X()+PL.X())/2, (PF.Y()+PL.Y())/2 );
607 if(ori) V2d = gp_Vec2d(PC,PF);
608 else V2d = gp_Vec2d(PC,PL);
610 double cp = (fp+lp)/2;
611 double dp2 = (lp-fp)/2;
612 NC = meshDS->AddNode(P0.X(), P0.Y(), P0.Z());
613 meshDS->SetNodeOnEdge(NC, edgeID, cp);
614 Nodes1.resize( myLayerPositions.size()+1 );
615 Nodes2.resize( myLayerPositions.size()+1 );
617 for(; i<myLayerPositions.size(); i++) {
618 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
619 P0.Y() + aVec.Y()*myLayerPositions[i],
620 P0.Z() + aVec.Z()*myLayerPositions[i] );
622 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
626 param = fp + dp2*(1-myLayerPositions[i]);
628 param = cp + dp2*myLayerPositions[i];
629 meshDS->SetNodeOnEdge(node, edgeID, param);
630 P = gp_Pnt( P0.X() - aVec.X()*myLayerPositions[i],
631 P0.Y() - aVec.Y()*myLayerPositions[i],
632 P0.Z() - aVec.Z()*myLayerPositions[i] );
633 node = meshDS->AddNode(P.X(), P.Y(), P.Z());
636 param = fp + dp2*(1-myLayerPositions[i]);
638 param = cp + dp2*myLayerPositions[i];
639 meshDS->SetNodeOnEdge(node, edgeID, param);
640 // parameters on face
641 gp_Pnt2d P2d( PC.X() + V2d.X()*myLayerPositions[i],
642 PC.Y() + V2d.Y()*myLayerPositions[i] );
645 Nodes1[ myLayerPositions.size() ] = NF;
646 Nodes2[ myLayerPositions.size() ] = NL;
647 // create 1D elements on edge
648 vector< const SMDS_MeshNode* > tmpNodes;
649 tmpNodes.resize(2*Nodes1.size()+1);
650 for(i=0; i<Nodes2.size(); i++)
651 tmpNodes[Nodes2.size()-i-1] = Nodes2[i];
652 tmpNodes[Nodes2.size()] = NC;
653 for(i=0; i<Nodes1.size(); i++)
654 tmpNodes[Nodes2.size()+1+i] = Nodes1[i];
655 for(i=1; i<tmpNodes.size(); i++) {
656 SMDS_MeshEdge* ME = myHelper->AddEdge( tmpNodes[i-1], tmpNodes[i] );
657 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
659 markEdgeAsComputedByMe( LinEdge1, aMesh.GetSubMesh( F ));
662 else // nbe==3 or ( nbe==2 && linEdge is INTERNAL )
664 if (nbe==2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
667 // one curve must be a part of circle and other curves must be
670 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
671 Handle(Geom_Line) aLine1 = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
672 Handle(Geom_Line) aLine2 = Handle(Geom_Line)::DownCast( getCurve( LinEdge2 ));
673 if( aCirc.IsNull() || aLine1.IsNull() || aLine2.IsNull() )
674 return error(COMPERR_BAD_SHAPE);
676 if ( !algo1d->ComputeCircularEdge( aMesh, CircEdge ))
677 return error( algo1d->GetComputeError() );
678 map< double, const SMDS_MeshNode* > theNodes;
679 if ( !GetSortedNodesOnEdge(aMesh.GetMeshDS(),CircEdge,true,theNodes))
680 return error("Circular edge is incorrectly meshed");
682 const SMDS_MeshNode* NF = theNodes.begin()->second;
683 const SMDS_MeshNode* NL = theNodes.rbegin()->second;
685 CNodes.push_back( NF );
686 map< double, const SMDS_MeshNode* >::iterator itn = theNodes.begin();
687 double fang = (*itn).first;
689 for(; itn != theNodes.end(); itn++ ) {
690 CNodes.push_back( (*itn).second );
691 double ang = (*itn).first - fang;
692 if( ang>PI ) ang = ang - 2*PI;
693 if( ang<-PI ) ang = ang + 2*PI;
694 Angles.Append( ang );
696 P1 = gp_Pnt( NF->X(), NF->Y(), NF->Z() );
697 gp_Pnt P2( NL->X(), NL->Y(), NL->Z() );
698 P0 = aCirc->Location();
700 bool linEdge1Computed, linEdge2Computed;
701 if ( !computeLayerPositions(P0,P1,LinEdge1,&linEdge1Computed))
704 Nodes1.resize( myLayerPositions.size()+1 );
705 Nodes2.resize( myLayerPositions.size()+1 );
707 // check that both linear edges have same hypotheses
708 if ( !computeLayerPositions(P0,P1,LinEdge2, &linEdge2Computed))
710 if ( Nodes1.size() != myLayerPositions.size()+1 )
711 return error("Different hypotheses apply to radial edges");
713 exp.Init( LinEdge1, TopAbs_VERTEX );
714 TopoDS_Vertex V1 = TopoDS::Vertex( exp.Current() );
716 TopoDS_Vertex V2 = TopoDS::Vertex( exp.Current() );
717 gp_Pnt PE1 = BRep_Tool::Pnt(V1);
718 gp_Pnt PE2 = BRep_Tool::Pnt(V2);
719 if( ( P1.Distance(PE1) > Precision::Confusion() ) &&
720 ( P1.Distance(PE2) > Precision::Confusion() ) )
722 std::swap( LinEdge1, LinEdge2 );
723 std::swap( linEdge1Computed, linEdge2Computed );
725 TopoDS_Vertex VC = V2;
726 if( ( P1.Distance(PE1) > Precision::Confusion() ) &&
727 ( P2.Distance(PE1) > Precision::Confusion() ) )
729 int vertID = meshDS->ShapeToIndex(VC);
732 if ( linEdge1Computed )
734 if (!GetSortedNodesOnEdge(aMesh.GetMeshDS(),LinEdge1,true,theNodes))
735 return error("Invalid mesh on a straight edge");
737 bool nodesFromP0ToP1 = ( theNodes.rbegin()->second == NF );
738 NC = const_cast<SMDS_MeshNode*>
739 ( nodesFromP0ToP1 ? theNodes.begin()->second : theNodes.rbegin()->second );
740 int i = 0, ir = Nodes1.size()-1;
741 int * pi = nodesFromP0ToP1 ? &i : &ir;
742 itn = theNodes.begin();
743 if ( nodesFromP0ToP1 ) ++itn;
744 for ( ; i < Nodes1.size(); ++i, --ir, ++itn )
746 Nodes1[*pi] = itn->second;
748 for ( i = 0; i < Nodes1.size()-1; ++i )
750 Points.Append( gpXYZ( Nodes1[i]));
751 Pnts2d1.Append( myHelper->GetNodeUV( F, Nodes1[i]));
756 int edgeID = meshDS->ShapeToIndex(LinEdge1);
759 Handle(Geom_Curve) Crv = BRep_Tool::Curve(LinEdge1,fp,lp);
760 gp_Pnt Ptmp = Crv->Value(fp);
762 if( P1.Distance(Ptmp) > Precision::Confusion() )
764 // get UV points for edge
766 BRep_Tool::UVPoints( LinEdge1, TopoDS::Face(aShape), PF, PL );
769 V2d = gp_Vec2d(PF,PL);
773 V2d = gp_Vec2d(PL,PF);
776 NC = const_cast<SMDS_MeshNode*>( VertexNode( VC, meshDS ));
779 NC = meshDS->AddNode(P0.X(), P0.Y(), P0.Z());
780 meshDS->SetNodeOnVertex(NC, vertID);
784 for(; i<myLayerPositions.size(); i++) {
785 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
786 P0.Y() + aVec.Y()*myLayerPositions[i],
787 P0.Z() + aVec.Z()*myLayerPositions[i] );
789 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
793 param = fp + dp*(1-myLayerPositions[i]);
795 param = fp + dp*myLayerPositions[i];
796 meshDS->SetNodeOnEdge(node, edgeID, param);
797 // parameters on face
798 gp_Pnt2d P2d( PC.X() + V2d.X()*myLayerPositions[i],
799 PC.Y() + V2d.Y()*myLayerPositions[i] );
802 Nodes1[ myLayerPositions.size() ] = NF;
803 // create 1D elements on edge
804 SMDS_MeshEdge* ME = myHelper->AddEdge( NC, Nodes1[0] );
805 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
806 for(i=1; i<Nodes1.size(); i++) {
807 ME = myHelper->AddEdge( Nodes1[i-1], Nodes1[i] );
808 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
810 if (nbe==2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
813 markEdgeAsComputedByMe( LinEdge1, aMesh.GetSubMesh( F ));
816 if ( linEdge2Computed )
818 if (!GetSortedNodesOnEdge(aMesh.GetMeshDS(),LinEdge2,true,theNodes))
819 return error("Invalid mesh on a straight edge");
821 bool nodesFromP0ToP2 = ( theNodes.rbegin()->second == NL );
822 int i = 0, ir = Nodes1.size()-1;
823 int * pi = nodesFromP0ToP2 ? &i : &ir;
824 itn = theNodes.begin();
825 if ( nodesFromP0ToP2 ) ++itn;
826 for ( ; i < Nodes2.size(); ++i, --ir, ++itn )
827 Nodes2[*pi] = itn->second;
831 int edgeID = meshDS->ShapeToIndex(LinEdge2);
832 gp_Vec aVec = gp_Vec(P0,P2);
834 Handle(Geom_Curve) Crv = BRep_Tool::Curve(LinEdge2,fp,lp);
835 gp_Pnt Ptmp = Crv->Value(fp);
837 if( P2.Distance(Ptmp) > Precision::Confusion() )
839 // get UV points for edge
841 BRep_Tool::UVPoints( LinEdge2, TopoDS::Face(aShape), PF, PL );
844 V2d = gp_Vec2d(PF,PL);
848 V2d = gp_Vec2d(PL,PF);
852 for(int i=0; i<myLayerPositions.size(); i++) {
853 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
854 P0.Y() + aVec.Y()*myLayerPositions[i],
855 P0.Z() + aVec.Z()*myLayerPositions[i] );
856 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
860 param = fp + dp*(1-myLayerPositions[i]);
862 param = fp + dp*myLayerPositions[i];
863 meshDS->SetNodeOnEdge(node, edgeID, param);
864 // parameters on face
865 gp_Pnt2d P2d( PC.X() + V2d.X()*myLayerPositions[i],
866 PC.Y() + V2d.Y()*myLayerPositions[i] );
868 Nodes2[ myLayerPositions.size() ] = NL;
869 // create 1D elements on edge
870 SMDS_MeshEdge* ME = myHelper->AddEdge( NC, Nodes2[0] );
871 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
872 for(int i=1; i<Nodes2.size(); i++) {
873 ME = myHelper->AddEdge( Nodes2[i-1], Nodes2[i] );
874 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
877 markEdgeAsComputedByMe( LinEdge2, aMesh.GetSubMesh( F ));
879 markEdgeAsComputedByMe( CircEdge, aMesh.GetSubMesh( F ));
882 bool IsForward = ( CircEdge.Orientation()==TopAbs_FORWARD );
884 // create nodes and mesh elements on face
885 // find axis of rotation
886 gp_Pnt P2 = gp_Pnt( CNodes[1]->X(), CNodes[1]->Y(), CNodes[1]->Z() );
889 gp_Vec Axis = Vec1.Crossed(Vec2);
892 //cout<<"Angles.Length() = "<<Angles.Length()<<" Points.Length() = "<<Points.Length()<<endl;
893 //cout<<"Nodes1.size() = "<<Nodes1.size()<<" Pnts2d1.Length() = "<<Pnts2d1.Length()<<endl;
894 for(; i<Angles.Length(); i++) {
895 vector< const SMDS_MeshNode* > tmpNodes;
896 tmpNodes.reserve(Nodes1.size());
898 gp_Ax1 theAxis(P0,gp_Dir(Axis));
899 aTrsf.SetRotation( theAxis, Angles.Value(i) );
901 aTrsf2d.SetRotation( PC, Angles.Value(i) );
904 for(; j<=Points.Length(); j++) {
906 Points.Value(j).Coord( cx, cy, cz );
907 aTrsf.Transforms( cx, cy, cz );
908 SMDS_MeshNode* node = myHelper->AddNode( cx, cy, cz );
909 // find parameters on face
910 Pnts2d1.Value(j).Coord( cx, cy );
911 aTrsf2d.Transforms( cx, cy );
913 meshDS->SetNodeOnFace( node, faceID, cx, cy );
914 tmpNodes[j-1] = node;
917 tmpNodes[Points.Length()] = CNodes[i];
919 for(j=0; j<Nodes1.size()-1; j++) {
922 MF = myHelper->AddFace( tmpNodes[j], Nodes1[j],
923 Nodes1[j+1], tmpNodes[j+1] );
925 MF = myHelper->AddFace( tmpNodes[j], tmpNodes[j+1],
926 Nodes1[j+1], Nodes1[j] );
927 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
932 MF = myHelper->AddFace( NC, Nodes1[0], tmpNodes[0] );
934 MF = myHelper->AddFace( NC, tmpNodes[0], Nodes1[0] );
935 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
936 for(j=0; j<Nodes1.size(); j++) {
937 Nodes1[j] = tmpNodes[j];
942 for(i=0; i<Nodes1.size()-1; i++) {
945 MF = myHelper->AddFace( Nodes2[i], Nodes1[i],
946 Nodes1[i+1], Nodes2[i+1] );
948 MF = myHelper->AddFace( Nodes2[i], Nodes2[i+1],
949 Nodes1[i+1], Nodes1[i] );
950 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
955 MF = myHelper->AddFace( NC, Nodes1[0], Nodes2[0] );
957 MF = myHelper->AddFace( NC, Nodes2[0], Nodes1[0] );
958 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
963 //================================================================================
965 * \brief Compute positions of nodes on the radial edge
966 * \retval bool - is a success
968 //================================================================================
970 bool StdMeshers_RadialQuadrangle_1D2D::computeLayerPositions(const gp_Pnt& p1,
972 const TopoDS_Edge& linEdge,
973 bool* linEdgeComputed)
975 // First, try to compute positions of layers
977 myLayerPositions.clear();
979 SMESH_Mesh * mesh = myHelper->GetMesh();
981 const SMESH_Hypothesis* hyp1D = myDistributionHypo ? myDistributionHypo->GetLayerDistribution() : 0;
982 int nbLayers = myNbLayerHypo ? myNbLayerHypo->GetNumberOfLayers() : 0;
984 if ( !hyp1D && !nbLayers )
986 // No own algo hypotheses assigned, so first try to find any 1D hypothesis.
988 TopoDS_Shape edge = linEdge;
989 if ( edge.IsNull() && !myHelper->GetSubShape().IsNull())
990 for ( TopExp_Explorer e(myHelper->GetSubShape(), TopAbs_EDGE); e.More(); e.Next())
992 if ( !edge.IsNull() )
994 // find a hyp usable by TNodeDistributor
995 SMESH_HypoFilter hypKind;
996 TNodeDistributor::GetDistributor(*mesh)->InitCompatibleHypoFilter(hypKind,/*ignoreAux=*/1);
997 hyp1D = mesh->GetHypothesis( edge, hypKind, /*fromAncestors=*/true);
1000 if ( hyp1D ) // try to compute with hyp1D
1002 if ( !TNodeDistributor::GetDistributor(*mesh)->Compute( myLayerPositions,p1,p2,*mesh,hyp1D ))
1003 if ( myDistributionHypo ) { // bad hyp assigned
1004 return error( TNodeDistributor::GetDistributor(*mesh)->GetComputeError() );
1007 // bad hyp found, its Ok, lets try with default nb of segnents
1011 if ( myLayerPositions.empty() ) // try to use nb of layers
1014 nbLayers = _gen->GetDefaultNbSegments();
1018 myLayerPositions.resize( nbLayers - 1 );
1019 for ( int z = 1; z < nbLayers; ++z )
1020 myLayerPositions[ z - 1 ] = double( z )/ double( nbLayers );
1024 // Second, check presence of a mesh built by other algo on linEdge
1025 // and mesh conformity to my hypothesis
1027 bool meshComputed = (!linEdge.IsNull() && !mesh->GetSubMesh(linEdge)->IsEmpty() );
1028 if ( linEdgeComputed ) *linEdgeComputed = meshComputed;
1032 vector< double > nodeParams;
1033 GetNodeParamOnEdge( mesh->GetMeshDS(), linEdge, nodeParams );
1035 if ( myLayerPositions.empty() )
1037 myLayerPositions.resize( nodeParams.size() - 2 );
1039 else if ( myDistributionHypo || myNbLayerHypo )
1041 // linEdge is computed by other algo. Check if there is a meshed face
1042 // using nodes on linEdge
1043 bool nodesAreUsed = false;
1044 TopTools_ListIteratorOfListOfShape ancestIt = mesh->GetAncestors( linEdge );
1045 for ( ; ancestIt.More() && !nodesAreUsed; ancestIt.Next() )
1046 if ( ancestIt.Value().ShapeType() == TopAbs_FACE )
1047 nodesAreUsed = (!mesh->GetSubMesh( ancestIt.Value() )->IsEmpty());
1048 if ( !nodesAreUsed ) {
1050 mesh->GetSubMesh( linEdge )->ComputeStateEngine( SMESH_subMesh::CLEAN );
1051 if ( linEdgeComputed ) *linEdgeComputed = false;
1053 else if ( myLayerPositions.size() != nodeParams.size()-2 ) {
1054 return error("Radial edge is meshed by other algorithm");
1059 return !myLayerPositions.empty();
1063 //=======================================================================
1064 //function : Evaluate
1066 //=======================================================================
1068 bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh,
1069 const TopoDS_Shape& aShape,
1070 MapShapeNbElems& aResMap)
1072 if( aShape.ShapeType() != TopAbs_FACE ) {
1075 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
1076 if( aResMap.count(sm) )
1079 vector<int>& aResVec =
1080 aResMap.insert( make_pair(sm, vector<int>(SMDSEntity_Last,0))).first->second;
1082 myHelper = new SMESH_MesherHelper( aMesh );
1083 myHelper->SetSubShape( aShape );
1084 auto_ptr<SMESH_MesherHelper> helperDeleter( myHelper );
1086 TNodeDistributor* algo1d = TNodeDistributor::GetDistributor(aMesh);
1088 TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
1089 int nbe = analyseFace( aShape, CircEdge, LinEdge1, LinEdge2 );
1090 if( nbe>3 || nbe < 1 || CircEdge.IsNull() )
1093 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
1094 if( aCirc.IsNull() )
1095 return error(COMPERR_BAD_SHAPE);
1097 gp_Pnt P0 = aCirc->Location();
1098 gp_Pnt P1 = aCirc->Value(0.);
1099 computeLayerPositions( P0, P1, LinEdge1 );
1101 int nb0d=0, nb2d_tria=0, nb2d_quad=0;
1102 bool isQuadratic = false, ok = true;
1105 // C1 must be a circle
1106 ok = algo1d->EvaluateCircularEdge( aMesh, CircEdge, aResMap );
1108 const vector<int>& aVec = aResMap[aMesh.GetSubMesh(CircEdge)];
1109 isQuadratic = aVec[SMDSEntity_Quad_Edge]>aVec[SMDSEntity_Edge];
1112 nb0d = (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1113 // radial medium nodes
1114 nb0d += (aVec[SMDSEntity_Node]+1) * (myLayerPositions.size()+1);
1115 // other medium nodes
1116 nb0d += (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1119 nb0d = (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1121 nb2d_tria = aVec[SMDSEntity_Node] + 1;
1125 else if(nbe==2 && LinEdge1.Orientation() != TopAbs_INTERNAL)
1127 // one curve must be a half of circle and other curve must be
1128 // a segment of line
1130 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge, &fp, &lp ));
1131 if( fabs(fabs(lp-fp)-PI) > Precision::Confusion() ) {
1132 // not half of circle
1133 return error(COMPERR_BAD_SHAPE);
1135 Handle(Geom_Line) aLine = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
1136 if( aLine.IsNull() ) {
1137 // other curve not line
1138 return error(COMPERR_BAD_SHAPE);
1140 ok = !aResMap.count( aMesh.GetSubMesh(LinEdge1) );
1142 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge1) ];
1143 ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
1146 ok = algo1d->EvaluateCircularEdge( aMesh, CircEdge, aResMap );
1149 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(CircEdge) ];
1150 isQuadratic = aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge];
1153 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1154 // radial medium nodes
1155 nb0d += aVec[SMDSEntity_Node] * (myLayerPositions.size()+1);
1156 // other medium nodes
1157 nb0d += (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1160 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1162 nb2d_tria = aVec[SMDSEntity_Node] + 1;
1163 nb2d_quad = nb2d_tria * myLayerPositions.size();
1164 // add evaluation for edges
1165 vector<int> aResVec(SMDSEntity_Last,0);
1167 aResVec[SMDSEntity_Node] = 4*myLayerPositions.size() + 3;
1168 aResVec[SMDSEntity_Quad_Edge] = 2*myLayerPositions.size() + 2;
1171 aResVec[SMDSEntity_Node] = 2*myLayerPositions.size() + 1;
1172 aResVec[SMDSEntity_Edge] = 2*myLayerPositions.size() + 2;
1174 aResMap[ aMesh.GetSubMesh(LinEdge1) ] = aResVec;
1177 else // nbe==3 or ( nbe==2 && linEdge is INTERNAL )
1179 if (nbe==2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
1180 LinEdge2 = LinEdge1;
1182 // one curve must be a part of circle and other curves must be
1184 Handle(Geom_Line) aLine1 = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
1185 Handle(Geom_Line) aLine2 = Handle(Geom_Line)::DownCast( getCurve( LinEdge2 ));
1186 if( aLine1.IsNull() || aLine2.IsNull() ) {
1187 // other curve not line
1188 return error(COMPERR_BAD_SHAPE);
1190 int nbLayers = myLayerPositions.size();
1191 computeLayerPositions( P0, P1, LinEdge2 );
1192 if ( nbLayers != myLayerPositions.size() )
1193 return error("Different hypotheses apply to radial edges");
1195 bool ok = !aResMap.count( aMesh.GetSubMesh(LinEdge1));
1197 if ( myDistributionHypo || myNbLayerHypo )
1198 ok = true; // override other 1d hyps
1200 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge1) ];
1201 ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
1204 if( ok && aResMap.count( aMesh.GetSubMesh(LinEdge2) )) {
1205 if ( myDistributionHypo || myNbLayerHypo )
1206 ok = true; // override other 1d hyps
1208 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge2) ];
1209 ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
1213 ok = algo1d->EvaluateCircularEdge( aMesh, CircEdge, aResMap );
1216 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(CircEdge) ];
1217 isQuadratic = aVec[SMDSEntity_Quad_Edge]>aVec[SMDSEntity_Edge];
1220 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1221 // radial medium nodes
1222 nb0d += aVec[SMDSEntity_Node] * (myLayerPositions.size()+1);
1223 // other medium nodes
1224 nb0d += (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1227 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1229 nb2d_tria = aVec[SMDSEntity_Node] + 1;
1230 nb2d_quad = nb2d_tria * myLayerPositions.size();
1231 // add evaluation for edges
1232 vector<int> aResVec(SMDSEntity_Last, 0);
1234 aResVec[SMDSEntity_Node] = 2*myLayerPositions.size() + 1;
1235 aResVec[SMDSEntity_Quad_Edge] = myLayerPositions.size() + 1;
1238 aResVec[SMDSEntity_Node] = myLayerPositions.size();
1239 aResVec[SMDSEntity_Edge] = myLayerPositions.size() + 1;
1241 sm = aMesh.GetSubMesh(LinEdge1);
1242 aResMap[sm] = aResVec;
1243 sm = aMesh.GetSubMesh(LinEdge2);
1244 aResMap[sm] = aResVec;
1251 aResVec[SMDSEntity_Quad_Triangle] = nb2d_tria;
1252 aResVec[SMDSEntity_Quad_Quadrangle] = nb2d_quad;
1255 aResVec[SMDSEntity_Triangle] = nb2d_tria;
1256 aResVec[SMDSEntity_Quadrangle] = nb2d_quad;
1262 sm = aMesh.GetSubMesh(aShape);
1263 SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
1264 smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,
1265 "Submesh can not be evaluated",this));