1 // Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 // SMESH SMESH : implementaion of SMESH idl descriptions
21 // File : StdMeshers_RadialQuadrangle_1D2D.cxx
23 // Created : Fri Oct 20 11:37:07 2006
24 // Author : Edward AGAPOV (eap)
26 #include "StdMeshers_RadialQuadrangle_1D2D.hxx"
28 #include "StdMeshers_NumberOfLayers.hxx"
29 #include "StdMeshers_LayerDistribution.hxx"
30 #include "StdMeshers_Regular_1D.hxx"
31 #include "StdMeshers_NumberOfSegments.hxx"
33 #include "SMDS_MeshNode.hxx"
34 #include "SMESHDS_SubMesh.hxx"
35 #include "SMESH_Gen.hxx"
36 #include "SMESH_HypoFilter.hxx"
37 #include "SMESH_Mesh.hxx"
38 #include "SMESH_MesherHelper.hxx"
39 #include "SMESH_subMesh.hxx"
40 #include "SMESH_subMeshEventListener.hxx"
42 #include "utilities.h"
44 #include <BRepAdaptor_Curve.hxx>
45 #include <BRepBuilderAPI_MakeEdge.hxx>
46 #include <BRep_Tool.hxx>
47 #include <GeomAPI_ProjectPointOnSurf.hxx>
48 #include <Geom_Circle.hxx>
49 #include <Geom_Line.hxx>
50 #include <Geom_TrimmedCurve.hxx>
51 #include <TColgp_SequenceOfPnt.hxx>
52 #include <TColgp_SequenceOfPnt2d.hxx>
54 #include <TopExp_Explorer.hxx>
55 #include <TopTools_ListIteratorOfListOfShape.hxx>
61 #define RETURN_BAD_RESULT(msg) { MESSAGE(")-: Error: " << msg); return false; }
62 #define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
65 //=======================================================================
66 //function : StdMeshers_RadialQuadrangle_1D2D
68 //=======================================================================
70 StdMeshers_RadialQuadrangle_1D2D::StdMeshers_RadialQuadrangle_1D2D(int hypId,
73 :SMESH_2D_Algo(hypId, studyId, gen)
75 _name = "RadialQuadrangle_1D2D";
76 _shapeType = (1 << TopAbs_FACE); // 1 bit per shape type
78 _compatibleHypothesis.push_back("LayerDistribution2D");
79 _compatibleHypothesis.push_back("NumberOfLayers2D");
81 myDistributionHypo = 0;
82 _requireDescretBoundary = false;
83 _supportSubmeshes = true;
87 //================================================================================
91 //================================================================================
93 StdMeshers_RadialQuadrangle_1D2D::~StdMeshers_RadialQuadrangle_1D2D()
97 //=======================================================================
98 //function : CheckHypothesis
100 //=======================================================================
102 bool StdMeshers_RadialQuadrangle_1D2D::CheckHypothesis
104 const TopoDS_Shape& aShape,
105 SMESH_Hypothesis::Hypothesis_Status& aStatus)
109 myDistributionHypo = 0;
111 list <const SMESHDS_Hypothesis * >::const_iterator itl;
113 const list <const SMESHDS_Hypothesis * >&hyps = GetUsedHypothesis(aMesh, aShape);
114 if ( hyps.size() == 0 ) {
115 aStatus = SMESH_Hypothesis::HYP_OK;
116 return true; // can work with no hypothesis
119 if ( hyps.size() > 1 ) {
120 aStatus = SMESH_Hypothesis::HYP_ALREADY_EXIST;
124 const SMESHDS_Hypothesis *theHyp = hyps.front();
126 string hypName = theHyp->GetName();
128 if (hypName == "NumberOfLayers2D") {
129 myNbLayerHypo = static_cast<const StdMeshers_NumberOfLayers *>(theHyp);
130 aStatus = SMESH_Hypothesis::HYP_OK;
133 if (hypName == "LayerDistribution2D") {
134 myDistributionHypo = static_cast<const StdMeshers_LayerDistribution *>(theHyp);
135 aStatus = SMESH_Hypothesis::HYP_OK;
138 aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
144 // ------------------------------------------------------------------------------
146 * \brief Listener used to mark edges meshed by StdMeshers_RadialQuadrangle_1D2D
148 class TEdgeMarker : public SMESH_subMeshEventListener
150 TEdgeMarker(): SMESH_subMeshEventListener(/*isDeletable=*/false) {}
152 //!< Return static listener
153 static SMESH_subMeshEventListener* getListener()
155 static TEdgeMarker theEdgeMarker;
156 return &theEdgeMarker;
158 //! Clear face sumbesh if something happens on edges
159 void ProcessEvent(const int event,
161 SMESH_subMesh* edgeSubMesh,
162 EventListenerData* data,
163 const SMESH_Hypothesis* /*hyp*/)
165 if ( data && !data->mySubMeshes.empty() && eventType == SMESH_subMesh::ALGO_EVENT)
167 ASSERT( data->mySubMeshes.front() != edgeSubMesh );
168 SMESH_subMesh* faceSubMesh = data->mySubMeshes.front();
169 faceSubMesh->ComputeStateEngine( SMESH_subMesh::CLEAN );
174 // ------------------------------------------------------------------------------
176 * \brief Mark an edge as computed by StdMeshers_RadialQuadrangle_1D2D
178 void markEdgeAsComputedByMe(const TopoDS_Edge& edge, SMESH_subMesh* faceSubMesh)
180 if ( SMESH_subMesh* edgeSM = faceSubMesh->GetFather()->GetSubMeshContaining( edge ))
182 if ( !edgeSM->GetEventListenerData( TEdgeMarker::getListener() ))
183 faceSubMesh->SetEventListener( TEdgeMarker::getListener(),
184 SMESH_subMeshEventListenerData::MakeData(faceSubMesh),
188 // ------------------------------------------------------------------------------
190 * \brief Return true if a radial edge was meshed with StdMeshers_RadialQuadrangle_1D2D with
191 * the same radial distribution
193 // bool isEdgeCompatiballyMeshed(const TopoDS_Edge& edge, SMESH_subMesh* faceSubMesh)
195 // if ( SMESH_subMesh* edgeSM = faceSubMesh->GetFather()->GetSubMeshContaining( edge ))
197 // if ( SMESH_subMeshEventListenerData* otherFaceData =
198 // edgeSM->GetEventListenerData( TEdgeMarker::getListener() ))
200 // // compare hypothesis aplied to two disk faces sharing radial edges
201 // SMESH_Mesh& mesh = *faceSubMesh->GetFather();
202 // SMESH_Algo* radialQuadAlgo = mesh.GetGen()->GetAlgo(mesh, faceSubMesh->GetSubShape() );
203 // SMESH_subMesh* otherFaceSubMesh = otherFaceData->mySubMeshes.front();
204 // list <const SMESHDS_Hypothesis *> hyps1 =
205 // radialQuadAlgo->GetUsedHypothesis( mesh, faceSubMesh->GetSubShape());
206 // list <const SMESHDS_Hypothesis *> hyps2 =
207 // radialQuadAlgo->GetUsedHypothesis( mesh, otherFaceSubMesh->GetSubShape());
208 // if( hyps1.empty() && hyps2.empty() )
209 // return true; // defaul hyps
210 // if ( hyps1.size() != hyps2.size() )
212 // return *hyps1.front() == *hyps2.front();
218 //================================================================================
220 * \brief Return base curve of the edge and extremum parameters
222 //================================================================================
224 Handle(Geom_Curve) getCurve(const TopoDS_Edge& edge, double* f=0, double* l=0)
226 Handle(Geom_Curve) C;
227 if ( !edge.IsNull() )
229 double first = 0., last = 0.;
230 C = BRep_Tool::Curve(edge, first, last);
233 Handle(Geom_TrimmedCurve) tc = Handle(Geom_TrimmedCurve)::DownCast(C);
234 while( !tc.IsNull() ) {
235 C = tc->BasisCurve();
236 tc = Handle(Geom_TrimmedCurve)::DownCast(C);
245 //================================================================================
247 * \brief Return edges of the face
248 * \retval int - nb of edges
250 //================================================================================
252 int analyseFace(const TopoDS_Shape& face,
253 TopoDS_Edge& CircEdge,
254 TopoDS_Edge& LinEdge1,
255 TopoDS_Edge& LinEdge2)
257 CircEdge.Nullify(); LinEdge1.Nullify(); LinEdge2.Nullify();
260 for ( TopExp_Explorer exp( face, TopAbs_EDGE ); exp.More(); exp.Next(), ++nbe )
262 const TopoDS_Edge& E = TopoDS::Edge( exp.Current() );
264 Handle(Geom_Curve) C = getCurve(E,&f,&l);
267 if ( C->IsKind( STANDARD_TYPE(Geom_Circle)))
269 if ( CircEdge.IsNull() )
274 else if ( LinEdge1.IsNull() )
283 //================================================================================
284 //================================================================================
286 * \brief Class computing layers distribution using data of
287 * StdMeshers_LayerDistribution hypothesis
289 //================================================================================
290 //================================================================================
292 class TNodeDistributor: public StdMeshers_Regular_1D
294 list <const SMESHDS_Hypothesis *> myUsedHyps;
296 // -----------------------------------------------------------------------------
297 static TNodeDistributor* GetDistributor(SMESH_Mesh& aMesh)
299 const int myID = -1000;
300 map < int, SMESH_1D_Algo * > & algoMap = aMesh.GetGen()->_map1D_Algo;
301 map < int, SMESH_1D_Algo * >::iterator id_algo = algoMap.find( myID );
302 if ( id_algo == algoMap.end() )
303 return new TNodeDistributor( myID, 0, aMesh.GetGen() );
304 return static_cast< TNodeDistributor* >( id_algo->second );
306 // -----------------------------------------------------------------------------
307 //! Computes distribution of nodes on a straight line ending at pIn and pOut
308 bool Compute( vector< double > & positions,
312 const SMESH_Hypothesis* hyp1d)
314 if ( !hyp1d ) return error( "Invalid LayerDistribution hypothesis");
316 double len = pIn.Distance( pOut );
317 if ( len <= DBL_MIN ) return error("Too close points of inner and outer shells");
320 myUsedHyps.push_back( hyp1d );
322 TopoDS_Edge edge = BRepBuilderAPI_MakeEdge( pIn, pOut );
323 SMESH_Hypothesis::Hypothesis_Status aStatus;
324 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, edge, aStatus ))
325 return error( "StdMeshers_Regular_1D::CheckHypothesis() failed "
326 "with LayerDistribution hypothesis");
328 BRepAdaptor_Curve C3D(edge);
329 double f = C3D.FirstParameter(), l = C3D.LastParameter();
330 list< double > params;
331 if ( !StdMeshers_Regular_1D::computeInternalParameters( aMesh, C3D, len, f, l, params, false ))
332 return error("StdMeshers_Regular_1D failed to compute layers distribution");
335 positions.reserve( params.size() );
336 for (list<double>::iterator itU = params.begin(); itU != params.end(); itU++)
337 positions.push_back( *itU / len );
340 // -----------------------------------------------------------------------------
341 //! Make mesh on an adge using assigned 1d hyp or defaut nb of segments
342 bool ComputeCircularEdge(SMESH_Mesh& aMesh,
343 const TopoDS_Edge& anEdge)
345 _gen->Compute( aMesh, anEdge);
346 SMESH_subMesh *sm = aMesh.GetSubMesh(anEdge);
347 if ( sm->GetComputeState() != SMESH_subMesh::COMPUTE_OK)
349 // find any 1d hyp assigned (there can be a hyp w/o algo)
350 myUsedHyps = SMESH_Algo::GetUsedHypothesis(aMesh, anEdge, /*ignoreAux=*/true);
351 Hypothesis_Status aStatus;
352 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, anEdge, aStatus ))
354 // no valid 1d hyp assigned, use default nb of segments
355 _hypType = NB_SEGMENTS;
356 _ivalue[ DISTR_TYPE_IND ] = StdMeshers_NumberOfSegments::DT_Regular;
357 _ivalue[ NB_SEGMENTS_IND ] = _gen->GetDefaultNbSegments();
359 return StdMeshers_Regular_1D::Compute( aMesh, anEdge );
363 // -----------------------------------------------------------------------------
364 //! Make mesh on an adge using assigned 1d hyp or defaut nb of segments
365 bool EvaluateCircularEdge(SMESH_Mesh& aMesh,
366 const TopoDS_Edge& anEdge,
367 MapShapeNbElems& aResMap)
369 _gen->Evaluate( aMesh, anEdge, aResMap );
370 if ( aResMap.count( aMesh.GetSubMesh( anEdge )))
373 // find any 1d hyp assigned
374 myUsedHyps = SMESH_Algo::GetUsedHypothesis(aMesh, anEdge, /*ignoreAux=*/true);
375 Hypothesis_Status aStatus;
376 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, anEdge, aStatus ))
378 // no valid 1d hyp assigned, use default nb of segments
379 _hypType = NB_SEGMENTS;
380 _ivalue[ DISTR_TYPE_IND ] = StdMeshers_NumberOfSegments::DT_Regular;
381 _ivalue[ NB_SEGMENTS_IND ] = _gen->GetDefaultNbSegments();
383 return StdMeshers_Regular_1D::Evaluate( aMesh, anEdge, aResMap );
386 // -----------------------------------------------------------------------------
387 TNodeDistributor( int hypId, int studyId, SMESH_Gen* gen)
388 : StdMeshers_Regular_1D( hypId, studyId, gen)
391 // -----------------------------------------------------------------------------
392 virtual const list <const SMESHDS_Hypothesis *> &
393 GetUsedHypothesis(SMESH_Mesh &, const TopoDS_Shape &, const bool)
397 // -----------------------------------------------------------------------------
401 //=======================================================================
403 * \brief Allow algo to do something after persistent restoration
404 * \param subMesh - restored submesh
406 * call markEdgeAsComputedByMe()
408 //=======================================================================
410 void StdMeshers_RadialQuadrangle_1D2D::SubmeshRestored(SMESH_subMesh* faceSubMesh)
412 if ( !faceSubMesh->IsEmpty() )
414 TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
415 analyseFace( faceSubMesh->GetSubShape(), CircEdge, LinEdge1, LinEdge2 );
416 if ( !CircEdge.IsNull() ) markEdgeAsComputedByMe( CircEdge, faceSubMesh );
417 if ( !LinEdge1.IsNull() ) markEdgeAsComputedByMe( LinEdge1, faceSubMesh );
418 if ( !LinEdge2.IsNull() ) markEdgeAsComputedByMe( LinEdge2, faceSubMesh );
422 //=======================================================================
425 //=======================================================================
427 bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh& aMesh,
428 const TopoDS_Shape& aShape)
431 SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
433 myHelper = new SMESH_MesherHelper( aMesh );
434 myHelper->IsQuadraticSubMesh( aShape );
435 // to delete helper at exit from Compute()
436 auto_ptr<SMESH_MesherHelper> helperDeleter( myHelper );
438 TNodeDistributor* algo1d = TNodeDistributor::GetDistributor(aMesh);
440 TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
441 int nbe = analyseFace( aShape, CircEdge, LinEdge1, LinEdge2 );
442 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
443 if( nbe>3 || nbe < 1 || aCirc.IsNull() )
444 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)");
447 // points for rotation
448 TColgp_SequenceOfPnt Points;
449 // angles for rotation
450 TColStd_SequenceOfReal Angles;
451 // Nodes1 and Nodes2 - nodes along radiuses
452 // CNodes - nodes on circle edge
453 vector< const SMDS_MeshNode* > Nodes1, Nodes2, CNodes;
455 // parameters edge nodes on face
456 TColgp_SequenceOfPnt2d Pnts2d1;
459 int faceID = meshDS->ShapeToIndex(aShape);
460 TopoDS_Face F = TopoDS::Face(aShape);
461 Handle(Geom_Surface) S = BRep_Tool::Surface(F);
466 if (!algo1d->ComputeCircularEdge( aMesh, CircEdge ))
467 return error( algo1d->GetComputeError() );
468 map< double, const SMDS_MeshNode* > theNodes;
469 if ( !GetSortedNodesOnEdge(aMesh.GetMeshDS(),CircEdge,true,theNodes))
470 return error("Circular edge is incorrectly meshed");
473 map< double, const SMDS_MeshNode* >::iterator itn = theNodes.begin();
474 const SMDS_MeshNode* NF = (*itn).second;
475 CNodes.push_back( (*itn).second );
476 double fang = (*itn).first;
477 if ( itn != theNodes.end() ) {
479 for(; itn != theNodes.end(); itn++ ) {
480 CNodes.push_back( (*itn).second );
481 double ang = (*itn).first - fang;
482 if( ang>PI ) ang = ang - 2*PI;
483 if( ang<-PI ) ang = ang + 2*PI;
484 Angles.Append( ang );
487 P1 = gp_Pnt( NF->X(), NF->Y(), NF->Z() );
488 P0 = aCirc->Location();
490 if ( !computeLayerPositions(P0,P1))
493 exp.Init( CircEdge, TopAbs_VERTEX );
494 TopoDS_Vertex V1 = TopoDS::Vertex( exp.Current() );
495 gp_Pnt2d p2dV = BRep_Tool::Parameters( V1, TopoDS::Face(aShape) );
497 NC = meshDS->AddNode(P0.X(), P0.Y(), P0.Z());
498 GeomAPI_ProjectPointOnSurf PPS(P0,S);
500 PPS.Parameters(1,U0,V0);
501 meshDS->SetNodeOnFace(NC, faceID, U0, V0);
502 PC = gp_Pnt2d(U0,V0);
505 gp_Vec2d aVec2d(PC,p2dV);
506 Nodes1.resize( myLayerPositions.size()+1 );
507 Nodes2.resize( myLayerPositions.size()+1 );
509 for(; i<myLayerPositions.size(); i++) {
510 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
511 P0.Y() + aVec.Y()*myLayerPositions[i],
512 P0.Z() + aVec.Z()*myLayerPositions[i] );
514 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
517 double U = PC.X() + aVec2d.X()*myLayerPositions[i];
518 double V = PC.Y() + aVec2d.Y()*myLayerPositions[i];
519 meshDS->SetNodeOnFace( node, faceID, U, V );
520 Pnts2d1.Append(gp_Pnt2d(U,V));
522 Nodes1[Nodes1.size()-1] = NF;
523 Nodes2[Nodes1.size()-1] = NF;
525 else if(nbe==2 && LinEdge1.Orientation() != TopAbs_INTERNAL )
527 // one curve must be a half of circle and other curve must be
530 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge, &fp, &lp ));
531 if( fabs(fabs(lp-fp)-PI) > Precision::Confusion() ) {
532 // not half of circle
533 return error(COMPERR_BAD_SHAPE);
535 Handle(Geom_Line) aLine = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
536 if( aLine.IsNull() ) {
537 // other curve not line
538 return error(COMPERR_BAD_SHAPE);
541 if ( !algo1d->ComputeCircularEdge( aMesh, CircEdge ))
542 return error( algo1d->GetComputeError() );
543 map< double, const SMDS_MeshNode* > theNodes;
544 if ( !GetSortedNodesOnEdge(aMesh.GetMeshDS(),CircEdge,true,theNodes) )
545 return error("Circular edge is incorrectly meshed");
547 map< double, const SMDS_MeshNode* >::iterator itn = theNodes.begin();
549 CNodes.push_back( itn->second );
550 double fang = (*itn).first;
552 for(; itn != theNodes.end(); itn++ ) {
553 CNodes.push_back( (*itn).second );
554 double ang = (*itn).first - fang;
555 if( ang>PI ) ang = ang - 2*PI;
556 if( ang<-PI ) ang = ang + 2*PI;
557 Angles.Append( ang );
559 const SMDS_MeshNode* NF = theNodes.begin()->second;
560 const SMDS_MeshNode* NL = theNodes.rbegin()->second;
561 P1 = gp_Pnt( NF->X(), NF->Y(), NF->Z() );
562 gp_Pnt P2( NL->X(), NL->Y(), NL->Z() );
563 P0 = aCirc->Location();
565 bool linEdgeComputed;
566 if ( !computeLayerPositions(P0,P1,LinEdge1,&linEdgeComputed))
569 if ( linEdgeComputed )
571 if (!GetSortedNodesOnEdge(aMesh.GetMeshDS(),LinEdge1,true,theNodes))
572 return error("Invalid mesh on a straight edge");
574 Nodes1.resize( myLayerPositions.size()+1 );
575 Nodes2.resize( myLayerPositions.size()+1 );
576 vector< const SMDS_MeshNode* > *pNodes1 = &Nodes1, *pNodes2 = &Nodes2;
577 bool nodesFromP0ToP1 = ( theNodes.rbegin()->second == NF );
578 if ( !nodesFromP0ToP1 ) std::swap( pNodes1, pNodes2 );
580 map< double, const SMDS_MeshNode* >::reverse_iterator ritn = theNodes.rbegin();
581 itn = theNodes.begin();
582 for ( int i = Nodes1.size()-1; i > -1; ++itn, ++ritn, --i )
584 (*pNodes1)[i] = ritn->second;
585 (*pNodes2)[i] = itn->second;
586 Points.Prepend( gpXYZ( Nodes1[i]));
587 Pnts2d1.Prepend( myHelper->GetNodeUV( F, Nodes1[i]));
589 NC = const_cast<SMDS_MeshNode*>( itn->second );
590 Points.Remove( Nodes1.size() );
595 int edgeID = meshDS->ShapeToIndex(LinEdge1);
597 Handle(Geom_Curve) Crv = BRep_Tool::Curve(LinEdge1,fp,lp);
601 if( P1.Distance(Ptmp) > Precision::Confusion() )
603 // get UV points for edge
605 BRep_Tool::UVPoints( LinEdge1, TopoDS::Face(aShape), PF, PL );
606 PC = gp_Pnt2d( (PF.X()+PL.X())/2, (PF.Y()+PL.Y())/2 );
608 if(ori) V2d = gp_Vec2d(PC,PF);
609 else V2d = gp_Vec2d(PC,PL);
611 double cp = (fp+lp)/2;
612 double dp2 = (lp-fp)/2;
613 NC = meshDS->AddNode(P0.X(), P0.Y(), P0.Z());
614 meshDS->SetNodeOnEdge(NC, edgeID, cp);
615 Nodes1.resize( myLayerPositions.size()+1 );
616 Nodes2.resize( myLayerPositions.size()+1 );
618 for(; i<myLayerPositions.size(); i++) {
619 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
620 P0.Y() + aVec.Y()*myLayerPositions[i],
621 P0.Z() + aVec.Z()*myLayerPositions[i] );
623 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
627 param = fp + dp2*(1-myLayerPositions[i]);
629 param = cp + dp2*myLayerPositions[i];
630 meshDS->SetNodeOnEdge(node, edgeID, param);
631 P = gp_Pnt( P0.X() - aVec.X()*myLayerPositions[i],
632 P0.Y() - aVec.Y()*myLayerPositions[i],
633 P0.Z() - aVec.Z()*myLayerPositions[i] );
634 node = meshDS->AddNode(P.X(), P.Y(), P.Z());
637 param = fp + dp2*(1-myLayerPositions[i]);
639 param = cp + dp2*myLayerPositions[i];
640 meshDS->SetNodeOnEdge(node, edgeID, param);
641 // parameters on face
642 gp_Pnt2d P2d( PC.X() + V2d.X()*myLayerPositions[i],
643 PC.Y() + V2d.Y()*myLayerPositions[i] );
646 Nodes1[ myLayerPositions.size() ] = NF;
647 Nodes2[ myLayerPositions.size() ] = NL;
648 // create 1D elements on edge
649 vector< const SMDS_MeshNode* > tmpNodes;
650 tmpNodes.resize(2*Nodes1.size()+1);
651 for(i=0; i<Nodes2.size(); i++)
652 tmpNodes[Nodes2.size()-i-1] = Nodes2[i];
653 tmpNodes[Nodes2.size()] = NC;
654 for(i=0; i<Nodes1.size(); i++)
655 tmpNodes[Nodes2.size()+1+i] = Nodes1[i];
656 for(i=1; i<tmpNodes.size(); i++) {
657 SMDS_MeshEdge* ME = myHelper->AddEdge( tmpNodes[i-1], tmpNodes[i] );
658 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
660 markEdgeAsComputedByMe( LinEdge1, aMesh.GetSubMesh( F ));
663 else // nbe==3 or ( nbe==2 && linEdge is INTERNAL )
665 if (nbe==2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
668 // one curve must be a part of circle and other curves must be
671 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
672 Handle(Geom_Line) aLine1 = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
673 Handle(Geom_Line) aLine2 = Handle(Geom_Line)::DownCast( getCurve( LinEdge2 ));
674 if( aCirc.IsNull() || aLine1.IsNull() || aLine2.IsNull() )
675 return error(COMPERR_BAD_SHAPE);
677 if ( !algo1d->ComputeCircularEdge( aMesh, CircEdge ))
678 return error( algo1d->GetComputeError() );
679 map< double, const SMDS_MeshNode* > theNodes;
680 if ( !GetSortedNodesOnEdge(aMesh.GetMeshDS(),CircEdge,true,theNodes))
681 return error("Circular edge is incorrectly meshed");
683 const SMDS_MeshNode* NF = theNodes.begin()->second;
684 const SMDS_MeshNode* NL = theNodes.rbegin()->second;
686 CNodes.push_back( NF );
687 map< double, const SMDS_MeshNode* >::iterator itn = theNodes.begin();
688 double fang = (*itn).first;
690 for(; itn != theNodes.end(); itn++ ) {
691 CNodes.push_back( (*itn).second );
692 double ang = (*itn).first - fang;
693 if( ang>PI ) ang = ang - 2*PI;
694 if( ang<-PI ) ang = ang + 2*PI;
695 Angles.Append( ang );
697 P1 = gp_Pnt( NF->X(), NF->Y(), NF->Z() );
698 gp_Pnt P2( NL->X(), NL->Y(), NL->Z() );
699 P0 = aCirc->Location();
701 bool linEdge1Computed, linEdge2Computed;
702 if ( !computeLayerPositions(P0,P1,LinEdge1,&linEdge1Computed))
705 Nodes1.resize( myLayerPositions.size()+1 );
706 Nodes2.resize( myLayerPositions.size()+1 );
708 // check that both linear edges have same hypotheses
709 if ( !computeLayerPositions(P0,P1,LinEdge2, &linEdge2Computed))
711 if ( Nodes1.size() != myLayerPositions.size()+1 )
712 return error("Different hypotheses apply to radial edges");
714 exp.Init( LinEdge1, TopAbs_VERTEX );
715 TopoDS_Vertex V1 = TopoDS::Vertex( exp.Current() );
717 TopoDS_Vertex V2 = TopoDS::Vertex( exp.Current() );
718 gp_Pnt PE1 = BRep_Tool::Pnt(V1);
719 gp_Pnt PE2 = BRep_Tool::Pnt(V2);
720 if( ( P1.Distance(PE1) > Precision::Confusion() ) &&
721 ( P1.Distance(PE2) > Precision::Confusion() ) )
723 std::swap( LinEdge1, LinEdge2 );
724 std::swap( linEdge1Computed, linEdge2Computed );
726 TopoDS_Vertex VC = V2;
727 if( ( P1.Distance(PE1) > Precision::Confusion() ) &&
728 ( P2.Distance(PE1) > Precision::Confusion() ) )
730 int vertID = meshDS->ShapeToIndex(VC);
733 if ( linEdge1Computed )
735 if (!GetSortedNodesOnEdge(aMesh.GetMeshDS(),LinEdge1,true,theNodes))
736 return error("Invalid mesh on a straight edge");
738 bool nodesFromP0ToP1 = ( theNodes.rbegin()->second == NF );
739 NC = const_cast<SMDS_MeshNode*>
740 ( nodesFromP0ToP1 ? theNodes.begin()->second : theNodes.rbegin()->second );
741 int i = 0, ir = Nodes1.size()-1;
742 int * pi = nodesFromP0ToP1 ? &i : &ir;
743 itn = theNodes.begin();
744 if ( nodesFromP0ToP1 ) ++itn;
745 for ( ; i < Nodes1.size(); ++i, --ir, ++itn )
747 Nodes1[*pi] = itn->second;
749 for ( i = 0; i < Nodes1.size()-1; ++i )
751 Points.Append( gpXYZ( Nodes1[i]));
752 Pnts2d1.Append( myHelper->GetNodeUV( F, Nodes1[i]));
757 int edgeID = meshDS->ShapeToIndex(LinEdge1);
760 Handle(Geom_Curve) Crv = BRep_Tool::Curve(LinEdge1,fp,lp);
761 gp_Pnt Ptmp = Crv->Value(fp);
763 if( P1.Distance(Ptmp) > Precision::Confusion() )
765 // get UV points for edge
767 BRep_Tool::UVPoints( LinEdge1, TopoDS::Face(aShape), PF, PL );
770 V2d = gp_Vec2d(PF,PL);
774 V2d = gp_Vec2d(PL,PF);
777 NC = const_cast<SMDS_MeshNode*>( VertexNode( VC, meshDS ));
780 NC = meshDS->AddNode(P0.X(), P0.Y(), P0.Z());
781 meshDS->SetNodeOnVertex(NC, vertID);
785 for(; i<myLayerPositions.size(); i++) {
786 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
787 P0.Y() + aVec.Y()*myLayerPositions[i],
788 P0.Z() + aVec.Z()*myLayerPositions[i] );
790 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
794 param = fp + dp*(1-myLayerPositions[i]);
796 param = fp + dp*myLayerPositions[i];
797 meshDS->SetNodeOnEdge(node, edgeID, param);
798 // parameters on face
799 gp_Pnt2d P2d( PC.X() + V2d.X()*myLayerPositions[i],
800 PC.Y() + V2d.Y()*myLayerPositions[i] );
803 Nodes1[ myLayerPositions.size() ] = NF;
804 // create 1D elements on edge
805 SMDS_MeshEdge* ME = myHelper->AddEdge( NC, Nodes1[0] );
806 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
807 for(i=1; i<Nodes1.size(); i++) {
808 ME = myHelper->AddEdge( Nodes1[i-1], Nodes1[i] );
809 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
811 if (nbe==2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
814 markEdgeAsComputedByMe( LinEdge1, aMesh.GetSubMesh( F ));
817 if ( linEdge2Computed )
819 if (!GetSortedNodesOnEdge(aMesh.GetMeshDS(),LinEdge2,true,theNodes))
820 return error("Invalid mesh on a straight edge");
822 bool nodesFromP0ToP2 = ( theNodes.rbegin()->second == NL );
823 int i = 0, ir = Nodes1.size()-1;
824 int * pi = nodesFromP0ToP2 ? &i : &ir;
825 itn = theNodes.begin();
826 if ( nodesFromP0ToP2 ) ++itn;
827 for ( ; i < Nodes2.size(); ++i, --ir, ++itn )
828 Nodes2[*pi] = itn->second;
832 int edgeID = meshDS->ShapeToIndex(LinEdge2);
833 gp_Vec aVec = gp_Vec(P0,P2);
835 Handle(Geom_Curve) Crv = BRep_Tool::Curve(LinEdge2,fp,lp);
836 gp_Pnt Ptmp = Crv->Value(fp);
838 if( P2.Distance(Ptmp) > Precision::Confusion() )
840 // get UV points for edge
842 BRep_Tool::UVPoints( LinEdge2, TopoDS::Face(aShape), PF, PL );
845 V2d = gp_Vec2d(PF,PL);
849 V2d = gp_Vec2d(PL,PF);
853 for(int i=0; i<myLayerPositions.size(); i++) {
854 gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
855 P0.Y() + aVec.Y()*myLayerPositions[i],
856 P0.Z() + aVec.Z()*myLayerPositions[i] );
857 SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
861 param = fp + dp*(1-myLayerPositions[i]);
863 param = fp + dp*myLayerPositions[i];
864 meshDS->SetNodeOnEdge(node, edgeID, param);
865 // parameters on face
866 gp_Pnt2d P2d( PC.X() + V2d.X()*myLayerPositions[i],
867 PC.Y() + V2d.Y()*myLayerPositions[i] );
869 Nodes2[ myLayerPositions.size() ] = NL;
870 // create 1D elements on edge
871 SMDS_MeshEdge* ME = myHelper->AddEdge( NC, Nodes2[0] );
872 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
873 for(int i=1; i<Nodes2.size(); i++) {
874 ME = myHelper->AddEdge( Nodes2[i-1], Nodes2[i] );
875 if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
878 markEdgeAsComputedByMe( LinEdge2, aMesh.GetSubMesh( F ));
880 markEdgeAsComputedByMe( CircEdge, aMesh.GetSubMesh( F ));
883 bool IsForward = ( CircEdge.Orientation()==TopAbs_FORWARD );
885 // create nodes and mesh elements on face
886 // find axis of rotation
887 gp_Pnt P2 = gp_Pnt( CNodes[1]->X(), CNodes[1]->Y(), CNodes[1]->Z() );
890 gp_Vec Axis = Vec1.Crossed(Vec2);
893 //cout<<"Angles.Length() = "<<Angles.Length()<<" Points.Length() = "<<Points.Length()<<endl;
894 //cout<<"Nodes1.size() = "<<Nodes1.size()<<" Pnts2d1.Length() = "<<Pnts2d1.Length()<<endl;
895 for(; i<Angles.Length(); i++) {
896 vector< const SMDS_MeshNode* > tmpNodes;
897 tmpNodes.reserve(Nodes1.size());
899 gp_Ax1 theAxis(P0,gp_Dir(Axis));
900 aTrsf.SetRotation( theAxis, Angles.Value(i) );
902 aTrsf2d.SetRotation( PC, Angles.Value(i) );
905 for(; j<=Points.Length(); j++) {
907 Points.Value(j).Coord( cx, cy, cz );
908 aTrsf.Transforms( cx, cy, cz );
909 SMDS_MeshNode* node = myHelper->AddNode( cx, cy, cz );
910 // find parameters on face
911 Pnts2d1.Value(j).Coord( cx, cy );
912 aTrsf2d.Transforms( cx, cy );
914 meshDS->SetNodeOnFace( node, faceID, cx, cy );
915 tmpNodes[j-1] = node;
918 tmpNodes[Points.Length()] = CNodes[i];
920 for(j=0; j<Nodes1.size()-1; j++) {
923 MF = myHelper->AddFace( tmpNodes[j], Nodes1[j],
924 Nodes1[j+1], tmpNodes[j+1] );
926 MF = myHelper->AddFace( tmpNodes[j], tmpNodes[j+1],
927 Nodes1[j+1], Nodes1[j] );
928 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
933 MF = myHelper->AddFace( NC, Nodes1[0], tmpNodes[0] );
935 MF = myHelper->AddFace( NC, tmpNodes[0], Nodes1[0] );
936 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
937 for(j=0; j<Nodes1.size(); j++) {
938 Nodes1[j] = tmpNodes[j];
943 for(i=0; i<Nodes1.size()-1; i++) {
946 MF = myHelper->AddFace( Nodes2[i], Nodes1[i],
947 Nodes1[i+1], Nodes2[i+1] );
949 MF = myHelper->AddFace( Nodes2[i], Nodes2[i+1],
950 Nodes1[i+1], Nodes1[i] );
951 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
956 MF = myHelper->AddFace( NC, Nodes1[0], Nodes2[0] );
958 MF = myHelper->AddFace( NC, Nodes2[0], Nodes1[0] );
959 if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
964 //================================================================================
966 * \brief Compute positions of nodes on the radial edge
967 * \retval bool - is a success
969 //================================================================================
971 bool StdMeshers_RadialQuadrangle_1D2D::computeLayerPositions(const gp_Pnt& p1,
973 const TopoDS_Edge& linEdge,
974 bool* linEdgeComputed)
976 // First, try to compute positions of layers
978 myLayerPositions.clear();
980 SMESH_Mesh * mesh = myHelper->GetMesh();
982 const SMESH_Hypothesis* hyp1D = myDistributionHypo ? myDistributionHypo->GetLayerDistribution() : 0;
983 int nbLayers = myNbLayerHypo ? myNbLayerHypo->GetNumberOfLayers() : 0;
985 if ( !hyp1D && !nbLayers )
987 // No own algo hypotheses assigned, so first try to find any 1D hypothesis.
989 TopoDS_Shape edge = linEdge;
990 if ( edge.IsNull() && !myHelper->GetSubShape().IsNull())
991 for ( TopExp_Explorer e(myHelper->GetSubShape(), TopAbs_EDGE); e.More(); e.Next())
993 if ( !edge.IsNull() )
995 // find a hyp usable by TNodeDistributor
996 SMESH_HypoFilter hypKind;
997 TNodeDistributor::GetDistributor(*mesh)->InitCompatibleHypoFilter(hypKind,/*ignoreAux=*/1);
998 hyp1D = mesh->GetHypothesis( edge, hypKind, /*fromAncestors=*/true);
1001 if ( hyp1D ) // try to compute with hyp1D
1003 if ( !TNodeDistributor::GetDistributor(*mesh)->Compute( myLayerPositions,p1,p2,*mesh,hyp1D )) {
1004 if ( myDistributionHypo ) { // bad hyp assigned
1005 return error( TNodeDistributor::GetDistributor(*mesh)->GetComputeError() );
1008 // bad hyp found, its Ok, lets try with default nb of segnents
1013 if ( myLayerPositions.empty() ) // try to use nb of layers
1016 nbLayers = _gen->GetDefaultNbSegments();
1020 myLayerPositions.resize( nbLayers - 1 );
1021 for ( int z = 1; z < nbLayers; ++z )
1022 myLayerPositions[ z - 1 ] = double( z )/ double( nbLayers );
1026 // Second, check presence of a mesh built by other algo on linEdge
1027 // and mesh conformity to my hypothesis
1029 bool meshComputed = (!linEdge.IsNull() && !mesh->GetSubMesh(linEdge)->IsEmpty() );
1030 if ( linEdgeComputed ) *linEdgeComputed = meshComputed;
1034 vector< double > nodeParams;
1035 GetNodeParamOnEdge( mesh->GetMeshDS(), linEdge, nodeParams );
1037 // nb of present nodes must be different in cases of 1 and 2 straight edges
1039 TopoDS_Vertex VV[2];
1040 TopExp::Vertices( linEdge, VV[0], VV[1]);
1041 const gp_Pnt* points[] = { &p1, &p2 };
1042 gp_Pnt vPoints[] = { BRep_Tool::Pnt(VV[0]), BRep_Tool::Pnt(VV[1]) };
1043 const double tol[] = { BRep_Tool::Tolerance(VV[0]), BRep_Tool::Tolerance(VV[1]) };
1044 bool pointsAreOnVertices = true;
1045 for ( int iP = 0; iP < 2 && pointsAreOnVertices; ++iP )
1046 pointsAreOnVertices = ( points[iP]->Distance( vPoints[0] ) < tol[0] ||
1047 points[iP]->Distance( vPoints[1] ) < tol[1] );
1049 int nbNodes = nodeParams.size() - 2; // 2 straight edges
1050 if ( !pointsAreOnVertices )
1051 nbNodes = ( nodeParams.size() - 3 ) / 2; // 1 straight edge
1053 if ( myLayerPositions.empty() )
1055 myLayerPositions.resize( nbNodes );
1057 else if ( myDistributionHypo || myNbLayerHypo )
1059 // linEdge is computed by other algo. Check if there is a meshed face
1060 // using nodes on linEdge
1061 bool nodesAreUsed = false;
1062 TopTools_ListIteratorOfListOfShape ancestIt = mesh->GetAncestors( linEdge );
1063 for ( ; ancestIt.More() && !nodesAreUsed; ancestIt.Next() )
1064 if ( ancestIt.Value().ShapeType() == TopAbs_FACE )
1065 nodesAreUsed = (!mesh->GetSubMesh( ancestIt.Value() )->IsEmpty());
1066 if ( !nodesAreUsed ) {
1068 mesh->GetSubMesh( linEdge )->ComputeStateEngine( SMESH_subMesh::CLEAN );
1069 if ( linEdgeComputed ) *linEdgeComputed = false;
1073 if ( myLayerPositions.size() != nbNodes )
1074 return error("Radial edge is meshed by other algorithm");
1079 return !myLayerPositions.empty();
1083 //=======================================================================
1084 //function : Evaluate
1086 //=======================================================================
1088 bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh,
1089 const TopoDS_Shape& aShape,
1090 MapShapeNbElems& aResMap)
1092 if( aShape.ShapeType() != TopAbs_FACE ) {
1095 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
1096 if( aResMap.count(sm) )
1099 vector<int>& aResVec =
1100 aResMap.insert( make_pair(sm, vector<int>(SMDSEntity_Last,0))).first->second;
1102 myHelper = new SMESH_MesherHelper( aMesh );
1103 myHelper->SetSubShape( aShape );
1104 auto_ptr<SMESH_MesherHelper> helperDeleter( myHelper );
1106 TNodeDistributor* algo1d = TNodeDistributor::GetDistributor(aMesh);
1108 TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
1109 int nbe = analyseFace( aShape, CircEdge, LinEdge1, LinEdge2 );
1110 if( nbe>3 || nbe < 1 || CircEdge.IsNull() )
1113 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
1114 if( aCirc.IsNull() )
1115 return error(COMPERR_BAD_SHAPE);
1117 gp_Pnt P0 = aCirc->Location();
1118 gp_Pnt P1 = aCirc->Value(0.);
1119 computeLayerPositions( P0, P1, LinEdge1 );
1121 int nb0d=0, nb2d_tria=0, nb2d_quad=0;
1122 bool isQuadratic = false, ok = true;
1125 // C1 must be a circle
1126 ok = algo1d->EvaluateCircularEdge( aMesh, CircEdge, aResMap );
1128 const vector<int>& aVec = aResMap[aMesh.GetSubMesh(CircEdge)];
1129 isQuadratic = aVec[SMDSEntity_Quad_Edge]>aVec[SMDSEntity_Edge];
1132 nb0d = (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1133 // radial medium nodes
1134 nb0d += (aVec[SMDSEntity_Node]+1) * (myLayerPositions.size()+1);
1135 // other medium nodes
1136 nb0d += (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1139 nb0d = (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1141 nb2d_tria = aVec[SMDSEntity_Node] + 1;
1145 else if(nbe==2 && LinEdge1.Orientation() != TopAbs_INTERNAL)
1147 // one curve must be a half of circle and other curve must be
1148 // a segment of line
1150 Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge, &fp, &lp ));
1151 if( fabs(fabs(lp-fp)-PI) > Precision::Confusion() ) {
1152 // not half of circle
1153 return error(COMPERR_BAD_SHAPE);
1155 Handle(Geom_Line) aLine = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
1156 if( aLine.IsNull() ) {
1157 // other curve not line
1158 return error(COMPERR_BAD_SHAPE);
1160 ok = !aResMap.count( aMesh.GetSubMesh(LinEdge1) );
1162 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge1) ];
1163 ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
1166 ok = algo1d->EvaluateCircularEdge( aMesh, CircEdge, aResMap );
1169 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(CircEdge) ];
1170 isQuadratic = aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge];
1173 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1174 // radial medium nodes
1175 nb0d += aVec[SMDSEntity_Node] * (myLayerPositions.size()+1);
1176 // other medium nodes
1177 nb0d += (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1180 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1182 nb2d_tria = aVec[SMDSEntity_Node] + 1;
1183 nb2d_quad = nb2d_tria * myLayerPositions.size();
1184 // add evaluation for edges
1185 vector<int> aResVec(SMDSEntity_Last,0);
1187 aResVec[SMDSEntity_Node] = 4*myLayerPositions.size() + 3;
1188 aResVec[SMDSEntity_Quad_Edge] = 2*myLayerPositions.size() + 2;
1191 aResVec[SMDSEntity_Node] = 2*myLayerPositions.size() + 1;
1192 aResVec[SMDSEntity_Edge] = 2*myLayerPositions.size() + 2;
1194 aResMap[ aMesh.GetSubMesh(LinEdge1) ] = aResVec;
1197 else // nbe==3 or ( nbe==2 && linEdge is INTERNAL )
1199 if (nbe==2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
1200 LinEdge2 = LinEdge1;
1202 // one curve must be a part of circle and other curves must be
1204 Handle(Geom_Line) aLine1 = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
1205 Handle(Geom_Line) aLine2 = Handle(Geom_Line)::DownCast( getCurve( LinEdge2 ));
1206 if( aLine1.IsNull() || aLine2.IsNull() ) {
1207 // other curve not line
1208 return error(COMPERR_BAD_SHAPE);
1210 int nbLayers = myLayerPositions.size();
1211 computeLayerPositions( P0, P1, LinEdge2 );
1212 if ( nbLayers != myLayerPositions.size() )
1213 return error("Different hypotheses apply to radial edges");
1215 bool ok = !aResMap.count( aMesh.GetSubMesh(LinEdge1));
1217 if ( myDistributionHypo || myNbLayerHypo )
1218 ok = true; // override other 1d hyps
1220 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge1) ];
1221 ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
1224 if( ok && aResMap.count( aMesh.GetSubMesh(LinEdge2) )) {
1225 if ( myDistributionHypo || myNbLayerHypo )
1226 ok = true; // override other 1d hyps
1228 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge2) ];
1229 ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
1233 ok = algo1d->EvaluateCircularEdge( aMesh, CircEdge, aResMap );
1236 const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(CircEdge) ];
1237 isQuadratic = aVec[SMDSEntity_Quad_Edge]>aVec[SMDSEntity_Edge];
1240 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1241 // radial medium nodes
1242 nb0d += aVec[SMDSEntity_Node] * (myLayerPositions.size()+1);
1243 // other medium nodes
1244 nb0d += (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
1247 nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
1249 nb2d_tria = aVec[SMDSEntity_Node] + 1;
1250 nb2d_quad = nb2d_tria * myLayerPositions.size();
1251 // add evaluation for edges
1252 vector<int> aResVec(SMDSEntity_Last, 0);
1254 aResVec[SMDSEntity_Node] = 2*myLayerPositions.size() + 1;
1255 aResVec[SMDSEntity_Quad_Edge] = myLayerPositions.size() + 1;
1258 aResVec[SMDSEntity_Node] = myLayerPositions.size();
1259 aResVec[SMDSEntity_Edge] = myLayerPositions.size() + 1;
1261 sm = aMesh.GetSubMesh(LinEdge1);
1262 aResMap[sm] = aResVec;
1263 sm = aMesh.GetSubMesh(LinEdge2);
1264 aResMap[sm] = aResVec;
1271 aResVec[SMDSEntity_Quad_Triangle] = nb2d_tria;
1272 aResVec[SMDSEntity_Quad_Quadrangle] = nb2d_quad;
1275 aResVec[SMDSEntity_Triangle] = nb2d_tria;
1276 aResVec[SMDSEntity_Quadrangle] = nb2d_quad;
1282 sm = aMesh.GetSubMesh(aShape);
1283 SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
1284 smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,
1285 "Submesh can not be evaluated",this));