src/SMESH/SMESH_Algo.cxx

   1 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   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.
  10 //
  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.
  15 //
  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
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 //  SMESH SMESH : implementaion of SMESH idl descriptions
  24 //  File   : SMESH_Algo.cxx
  25 //  Author : Paul RASCLE, EDF
  26 //  Module : SMESH
  27
  28 #include "SMESH_Algo.hxx"
  29
  30 #include "SMDS_EdgePosition.hxx"
  31 #include "SMDS_FacePosition.hxx"
  32 #include "SMDS_MeshElement.hxx"
  33 #include "SMDS_MeshNode.hxx"
  34 #include "SMDS_VolumeTool.hxx"
  35 #include "SMESHDS_Mesh.hxx"
  36 #include "SMESHDS_SubMesh.hxx"
  37 #include "SMESH_Comment.hxx"
  38 #include "SMESH_Gen.hxx"
  39 #include "SMESH_HypoFilter.hxx"
  40 #include "SMESH_Mesh.hxx"
  41 #include "SMESH_TypeDefs.hxx"
  42 #include "SMESH_subMesh.hxx"
  43
  44 #include <Basics_OCCTVersion.hxx>
  45
  46 #include <BRepAdaptor_Curve.hxx>
  47 #include <BRepLProp.hxx>
  48 #include <BRep_Tool.hxx>
  49 #include <GCPnts_AbscissaPoint.hxx>
  50 #include <GeomAdaptor_Curve.hxx>
  51 #include <Geom_Surface.hxx>
  52 #include <TopExp.hxx>
  53 #include <TopExp_Explorer.hxx>
  54 #include <TopLoc_Location.hxx>
  55 #include <TopTools_ListIteratorOfListOfShape.hxx>
  56 #include <TopTools_ListOfShape.hxx>
  57 #include <TopoDS.hxx>
  58 #include <TopoDS_Edge.hxx>
  59 #include <TopoDS_Face.hxx>
  60 #include <TopoDS_Vertex.hxx>
  61 #include <TopoDS_Wire.hxx>
  62 #include <gp_Pnt.hxx>
  63 #include <gp_Pnt2d.hxx>
  64 #include <gp_Vec.hxx>
  65
  66 #include <Standard_ErrorHandler.hxx>
  67 #include <Standard_Failure.hxx>
  68
  69 #include "utilities.h"
  70
  71 #include <algorithm>
  72 #include <limits>
  73
  74 using namespace std;
  75
  76 //=============================================================================
  77 /*!
  78  *
  79  */
  80 //=============================================================================
  81
  82 SMESH_Algo::SMESH_Algo (int hypId, int studyId, SMESH_Gen * gen)
  83   : SMESH_Hypothesis(hypId, studyId, gen)
  84 {
  85   gen->_mapAlgo[hypId] = this;
  86
  87   _onlyUnaryInput = _requireDiscreteBoundary = _requireShape = true;
  88   _quadraticMesh = _supportSubmeshes = false;
  89   _error = COMPERR_OK;
  90   for ( int i = 0; i < 4; ++i )
  91     _neededLowerHyps[ i ] = false;
  92 }
  93
  94 //=============================================================================
  95 /*!
  96  *
  97  */
  98 //=============================================================================
  99
 100 SMESH_Algo::~SMESH_Algo()
 101 {
 102 }
 103
 104 //=============================================================================
 105 /*!
 106  *
 107  */
 108 //=============================================================================
 109
 110 SMESH_0D_Algo::SMESH_0D_Algo(int hypId, int studyId, SMESH_Gen* gen)
 111   : SMESH_Algo(hypId, studyId, gen)
 112 {
 113   _shapeType = (1 << TopAbs_VERTEX);
 114   _type = ALGO_0D;
 115   gen->_map0D_Algo[hypId] = this;
 116 }
 117 SMESH_1D_Algo::SMESH_1D_Algo(int hypId, int studyId, SMESH_Gen* gen)
 118   : SMESH_Algo(hypId, studyId, gen)
 119 {
 120   _shapeType = (1 << TopAbs_EDGE);
 121   _type = ALGO_1D;
 122   gen->_map1D_Algo[hypId] = this;
 123 }
 124 SMESH_2D_Algo::SMESH_2D_Algo(int hypId, int studyId, SMESH_Gen* gen)
 125   : SMESH_Algo(hypId, studyId, gen)
 126 {
 127   _shapeType = (1 << TopAbs_FACE);
 128   _type = ALGO_2D;
 129   gen->_map2D_Algo[hypId] = this;
 130 }
 131 SMESH_3D_Algo::SMESH_3D_Algo(int hypId, int studyId, SMESH_Gen* gen)
 132   : SMESH_Algo(hypId, studyId, gen)
 133 {
 134   _shapeType = (1 << TopAbs_SOLID);
 135   _type = ALGO_3D;
 136   gen->_map3D_Algo[hypId] = this;
 137 }
 138
 139 //=============================================================================
 140 /*!
 141  * Usually an algoritm has nothing to save
 142  */
 143 //=============================================================================
 144
 145 ostream & SMESH_Algo::SaveTo(ostream & save) { return save; }
 146 istream & SMESH_Algo::LoadFrom(istream & load) { return load; }
 147
 148 //=============================================================================
 149 /*!
 150  *
 151  */
 152 //=============================================================================
 153
 154 const vector < string > &SMESH_Algo::GetCompatibleHypothesis()
 155 {
 156   return _compatibleHypothesis;
 157 }
 158
 159 //=============================================================================
 160 /*!
 161  *  List the hypothesis used by the algorithm associated to the shape.
 162  *  Hypothesis associated to father shape -are- taken into account (see
 163  *  GetAppliedHypothesis). Relevant hypothesis have a name (type) listed in
 164  *  the algorithm. This method could be surcharged by specific algorithms, in
 165  *  case of several hypothesis simultaneously applicable.
 166  */
 167 //=============================================================================
 168
 169 const list <const SMESHDS_Hypothesis *> &
 170 SMESH_Algo::GetUsedHypothesis(SMESH_Mesh &         aMesh,
 171                               const TopoDS_Shape & aShape,
 172                               const bool           ignoreAuxiliary)
 173 {
 174   _usedHypList.clear();
 175   SMESH_HypoFilter filter;
 176   if ( InitCompatibleHypoFilter( filter, ignoreAuxiliary ))
 177   {
 178     aMesh.GetHypotheses( aShape, filter, _usedHypList, true );
 179     if ( ignoreAuxiliary && _usedHypList.size() > 1 )
 180       _usedHypList.clear(); //only one compatible hypothesis allowed
 181   }
 182   return _usedHypList;
 183 }
 184
 185 //=============================================================================
 186 /*!
 187  *  List the relevant hypothesis associated to the shape. Relevant hypothesis
 188  *  have a name (type) listed in the algorithm. Hypothesis associated to
 189  *  father shape -are not- taken into account (see GetUsedHypothesis)
 190  */
 191 //=============================================================================
 192
 193 const list<const SMESHDS_Hypothesis *> &
 194 SMESH_Algo::GetAppliedHypothesis(SMESH_Mesh &         aMesh,
 195                                  const TopoDS_Shape & aShape,
 196                                  const bool           ignoreAuxiliary)
 197 {
 198   _appliedHypList.clear();
 199   SMESH_HypoFilter filter;
 200   if ( InitCompatibleHypoFilter( filter, ignoreAuxiliary ))
 201     aMesh.GetHypotheses( aShape, filter, _appliedHypList, false );
 202
 203   return _appliedHypList;
 204 }
 205
 206 //=============================================================================
 207 /*!
 208  *  Compute length of an edge
 209  */
 210 //=============================================================================
 211
 212 double SMESH_Algo::EdgeLength(const TopoDS_Edge & E)
 213 {
 214   double UMin = 0, UMax = 0;
 215   if (BRep_Tool::Degenerated(E))
 216     return 0;
 217   TopLoc_Location L;
 218   Handle(Geom_Curve) C = BRep_Tool::Curve(E, L, UMin, UMax);
 219   GeomAdaptor_Curve AdaptCurve(C, UMin, UMax); //range is important for periodic curves
 220   double length = GCPnts_AbscissaPoint::Length(AdaptCurve, UMin, UMax);
 221   return length;
 222 }
 223
 224 //================================================================================
 225 /*!
 226  * \brief Calculate normal of a mesh face
 227  */
 228 //================================================================================
 229
 230 bool SMESH_Algo::FaceNormal(const SMDS_MeshElement* F, gp_XYZ& normal, bool normalized)
 231 {
 232   if ( !F || F->GetType() != SMDSAbs_Face )
 233     return false;
 234
 235   normal.SetCoord(0,0,0);
 236   int nbNodes = F->IsQuadratic() ? F->NbNodes()/2 : F->NbNodes();
 237   for ( int i = 0; i < nbNodes-2; ++i )
 238   {
 239     gp_XYZ p[3];
 240     for ( int n = 0; n < 3; ++n )
 241     {
 242       const SMDS_MeshNode* node = F->GetNode( i + n );
 243       p[n].SetCoord( node->X(), node->Y(), node->Z() );
 244     }
 245     normal += ( p[2] - p[1] ) ^ ( p[0] - p[1] );
 246   }
 247   double size2 = normal.SquareModulus();
 248   bool ok = ( size2 > numeric_limits<double>::min() * numeric_limits<double>::min());
 249   if ( normalized && ok )
 250     normal /= sqrt( size2 );
 251
 252   return ok;
 253 }
 254
 255 //================================================================================
 256 /*!
 257  * \brief Find out elements orientation on a geometrical face
 258  * \param theFace - The face correctly oriented in the shape being meshed
 259  * \param theMeshDS - The mesh data structure
 260  * \retval bool - true if the face normal and the normal of first element
 261  *                in the correspoding submesh point in different directions
 262  */
 263 //================================================================================
 264
 265 bool SMESH_Algo::IsReversedSubMesh (const TopoDS_Face&  theFace,
 266                                     SMESHDS_Mesh*       theMeshDS)
 267 {
 268   if ( theFace.IsNull() || !theMeshDS )
 269     return false;
 270
 271   // find out orientation of a meshed face
 272   int faceID = theMeshDS->ShapeToIndex( theFace );
 273   TopoDS_Shape aMeshedFace = theMeshDS->IndexToShape( faceID );
 274   bool isReversed = ( theFace.Orientation() != aMeshedFace.Orientation() );
 275
 276   const SMESHDS_SubMesh * aSubMeshDSFace = theMeshDS->MeshElements( faceID );
 277   if ( !aSubMeshDSFace )
 278     return isReversed;
 279
 280   // find element with node located on face and get its normal
 281   const SMDS_FacePosition* facePos = 0;
 282   int vertexID = 0;
 283   gp_Pnt nPnt[3];
 284   gp_Vec Ne;
 285   bool normalOK = false;
 286   SMDS_ElemIteratorPtr iteratorElem = aSubMeshDSFace->GetElements();
 287   while ( iteratorElem->more() ) // loop on elements on theFace
 288   {
 289     const SMDS_MeshElement* elem = iteratorElem->next();
 290     if ( elem && elem->NbNodes() > 2 ) {
 291       SMDS_ElemIteratorPtr nodesIt = elem->nodesIterator();
 292       const SMDS_FacePosition* fPos = 0;
 293       int i = 0, vID = 0;
 294       while ( nodesIt->more() ) { // loop on nodes
 295         const SMDS_MeshNode* node
 296           = static_cast<const SMDS_MeshNode *>(nodesIt->next());
 297         if ( i == 3 ) i = 2;
 298         nPnt[ i++ ].SetCoord( node->X(), node->Y(), node->Z() );
 299         // check position
 300         const SMDS_PositionPtr& pos = node->GetPosition();
 301         if ( !pos ) continue;
 302         if ( pos->GetTypeOfPosition() == SMDS_TOP_FACE ) {
 303           fPos = dynamic_cast< const SMDS_FacePosition* >( pos );
 304         }
 305         else if ( pos->GetTypeOfPosition() == SMDS_TOP_VERTEX ) {
 306           vID = node->getshapeId();
 307         }
 308       }
 309       if ( fPos || ( !normalOK && vID )) {
 310         // compute normal
 311         gp_Vec v01( nPnt[0], nPnt[1] ), v02( nPnt[0], nPnt[2] );
 312         if ( v01.SquareMagnitude() > RealSmall() &&
 313              v02.SquareMagnitude() > RealSmall() )
 314         {
 315           Ne = v01 ^ v02;
 316           normalOK = ( Ne.SquareMagnitude() > RealSmall() );
 317         }
 318         // we need position on theFace or at least on vertex
 319         if ( normalOK ) {
 320           vertexID = vID;
 321           if ((facePos = fPos))
 322             break;
 323         }
 324       }
 325     }
 326   }
 327   if ( !normalOK )
 328     return isReversed;
 329
 330   // node position on face
 331   double u,v;
 332   if ( facePos ) {
 333     u = facePos->GetUParameter();
 334     v = facePos->GetVParameter();
 335   }
 336   else if ( vertexID ) {
 337     TopoDS_Shape V = theMeshDS->IndexToShape( vertexID );
 338     if ( V.IsNull() || V.ShapeType() != TopAbs_VERTEX )
 339       return isReversed;
 340     gp_Pnt2d uv = BRep_Tool::Parameters( TopoDS::Vertex( V ), theFace );
 341     u = uv.X();
 342     v = uv.Y();
 343   }
 344   else
 345   {
 346     return isReversed;
 347   }
 348
 349   // face normal at node position
 350   TopLoc_Location loc;
 351   Handle(Geom_Surface) surf = BRep_Tool::Surface( theFace, loc );
 352   // if ( surf.IsNull() || surf->Continuity() < GeomAbs_C1 )
 353   // some surfaces not detected as GeomAbs_C1 are nevertheless correct for meshing
 354   if ( surf.IsNull() || surf->Continuity() < GeomAbs_C0 )
 355     {
 356       if (!surf.IsNull())
 357         MESSAGE("surf->Continuity() < GeomAbs_C1 " << (surf->Continuity() < GeomAbs_C1));
 358       return isReversed;
 359     }
 360   gp_Vec d1u, d1v;
 361   surf->D1( u, v, nPnt[0], d1u, d1v );
 362   gp_Vec Nf = (d1u ^ d1v).Transformed( loc );
 363
 364   if ( theFace.Orientation() == TopAbs_REVERSED )
 365     Nf.Reverse();
 366
 367   return Ne * Nf < 0.;
 368 }
 369
 370 //================================================================================
 371 /*!
 372  * \brief Just return false as the algorithm does not hold parameters values
 373  */
 374 //================================================================================
 375
 376 bool SMESH_Algo::SetParametersByMesh(const SMESH_Mesh* /*theMesh*/,
 377                                      const TopoDS_Shape& /*theShape*/)
 378 {
 379   return false;
 380 }
 381 bool SMESH_Algo::SetParametersByDefaults(const TDefaults& , const SMESH_Mesh*)
 382 {
 383   return false;
 384 }
 385 //================================================================================
 386 /*!
 387  * \brief Fill vector of node parameters on geometrical edge, including vertex nodes
 388  * \param theMesh - The mesh containing nodes
 389  * \param theEdge - The geometrical edge of interest
 390  * \param theParams - The resulting vector of sorted node parameters
 391  * \retval bool - false if not all parameters are OK
 392  */
 393 //================================================================================
 394
 395 bool SMESH_Algo::GetNodeParamOnEdge(const SMESHDS_Mesh* theMesh,
 396                                     const TopoDS_Edge&  theEdge,
 397                                     vector< double > &  theParams)
 398 {
 399   theParams.clear();
 400
 401   if ( !theMesh || theEdge.IsNull() )
 402     return false;
 403
 404   SMESHDS_SubMesh * eSubMesh = theMesh->MeshElements( theEdge );
 405   if ( !eSubMesh || !eSubMesh->GetElements()->more() )
 406     return false; // edge is not meshed
 407
 408   //int nbEdgeNodes = 0;
 409   set < double > paramSet;
 410   if ( eSubMesh )
 411   {
 412     // loop on nodes of an edge: sort them by param on edge
 413     SMDS_NodeIteratorPtr nIt = eSubMesh->GetNodes();
 414     while ( nIt->more() )
 415     {
 416       const SMDS_MeshNode* node = nIt->next();
 417       const SMDS_PositionPtr& pos = node->GetPosition();
 418       if ( pos->GetTypeOfPosition() != SMDS_TOP_EDGE )
 419         return false;
 420       const SMDS_EdgePosition* epos =
 421         static_cast<const SMDS_EdgePosition*>(node->GetPosition());
 422       if ( !paramSet.insert( epos->GetUParameter() ).second )
 423         return false; // equal parameters
 424     }
 425   }
 426   // add vertex nodes params
 427   TopoDS_Vertex V1,V2;
 428   TopExp::Vertices( theEdge, V1, V2);
 429   if ( VertexNode( V1, theMesh ) &&
 430        !paramSet.insert( BRep_Tool::Parameter(V1,theEdge) ).second )
 431     return false; // there are equal parameters
 432   if ( VertexNode( V2, theMesh ) &&
 433        !paramSet.insert( BRep_Tool::Parameter(V2,theEdge) ).second )
 434     return false; // there are equal parameters
 435
 436   // fill the vector
 437   theParams.resize( paramSet.size() );
 438   set < double >::iterator   par    = paramSet.begin();
 439   vector< double >::iterator vecPar = theParams.begin();
 440   for ( ; par != paramSet.end(); ++par, ++vecPar )
 441     *vecPar = *par;
 442
 443   return theParams.size() > 1;
 444 }
 445
 446 //================================================================================
 447 /*!
 448  * \brief Fill vector of node parameters on geometrical edge, including vertex nodes
 449  * \param theMesh - The mesh containing nodes
 450  * \param theEdge - The geometrical edge of interest
 451  * \param theParams - The resulting vector of sorted node parameters
 452  * \retval bool - false if not all parameters are OK
 453  */
 454 //================================================================================
 455
 456 bool SMESH_Algo::GetSortedNodesOnEdge(const SMESHDS_Mesh*                   theMesh,
 457                                       const TopoDS_Edge&                    theEdge,
 458                                       const bool                            ignoreMediumNodes,
 459                                       map< double, const SMDS_MeshNode* > & theNodes)
 460 {
 461   theNodes.clear();
 462
 463   if ( !theMesh || theEdge.IsNull() )
 464     return false;
 465
 466   SMESHDS_SubMesh * eSubMesh = theMesh->MeshElements( theEdge );
 467   if ( !eSubMesh || !eSubMesh->GetElements()->more() )
 468     return false; // edge is not meshed
 469
 470   int nbNodes = 0;
 471   set < double > paramSet;
 472   if ( eSubMesh )
 473   {
 474     // loop on nodes of an edge: sort them by param on edge
 475     SMDS_NodeIteratorPtr nIt = eSubMesh->GetNodes();
 476     while ( nIt->more() )
 477     {
 478       const SMDS_MeshNode* node = nIt->next();
 479       if ( ignoreMediumNodes ) {
 480         SMDS_ElemIteratorPtr elemIt = node->GetInverseElementIterator();
 481         if ( elemIt->more() && elemIt->next()->IsMediumNode( node ))
 482           continue;
 483       }
 484       const SMDS_PositionPtr& pos = node->GetPosition();
 485       if ( pos->GetTypeOfPosition() != SMDS_TOP_EDGE )
 486         return false;
 487       const SMDS_EdgePosition* epos =
 488         static_cast<const SMDS_EdgePosition*>(node->GetPosition());
 489       theNodes.insert( theNodes.end(), make_pair( epos->GetUParameter(), node ));
 490       //MESSAGE("U " << epos->GetUParameter() << " ID " << node->GetID());
 491       ++nbNodes;
 492     }
 493   }
 494   // add vertex nodes
 495   TopoDS_Vertex v1, v2;
 496   TopExp::Vertices(theEdge, v1, v2);
 497   const SMDS_MeshNode* n1 = VertexNode( v1, (SMESHDS_Mesh*) theMesh );
 498   const SMDS_MeshNode* n2 = VertexNode( v2, (SMESHDS_Mesh*) theMesh );
 499   //MESSAGE("Vertices ID " << n1->GetID() << " " << n2->GetID());
 500   Standard_Real f, l;
 501   BRep_Tool::Range(theEdge, f, l);
 502   if ( v1.Orientation() != TopAbs_FORWARD )
 503     std::swap( f, l );
 504   if ( n1 && ++nbNodes )
 505     theNodes.insert( make_pair( f, n1 ));
 506   if ( n2 && ++nbNodes )
 507     theNodes.insert( make_pair( l, n2 ));
 508
 509   return theNodes.size() == nbNodes;
 510 }
 511
 512 //================================================================================
 513 /*!
 514  * \brief Make filter recognize only compatible hypotheses
 515  * \param theFilter - the filter to initialize
 516  * \param ignoreAuxiliary - make filter ignore compatible auxiliary hypotheses
 517  */
 518 //================================================================================
 519
 520 bool SMESH_Algo::InitCompatibleHypoFilter( SMESH_HypoFilter & theFilter,
 521                                            const bool         ignoreAuxiliary) const
 522 {
 523   if ( !_compatibleHypothesis.empty() )
 524   {
 525     theFilter.Init( theFilter.HasName( _compatibleHypothesis[0] ));
 526     for ( int i = 1; i < _compatibleHypothesis.size(); ++i )
 527       theFilter.Or( theFilter.HasName( _compatibleHypothesis[ i ] ));
 528
 529     if ( ignoreAuxiliary )
 530       theFilter.AndNot( theFilter.IsAuxiliary() );
 531
 532     return true;
 533   }
 534   return false;
 535 }
 536
 537 //================================================================================
 538 /*!
 539  * \brief Return continuity of two edges
 540  * \param E1 - the 1st edge
 541  * \param E2 - the 2nd edge
 542  * \retval GeomAbs_Shape - regularity at the junction between E1 and E2
 543  */
 544 //================================================================================
 545
 546 GeomAbs_Shape SMESH_Algo::Continuity(TopoDS_Edge E1,
 547                                      TopoDS_Edge E2)
 548 {
 549   //E1.Orientation(TopAbs_FORWARD), E2.Orientation(TopAbs_FORWARD); // avoid pb with internal edges
 550   if (E1.Orientation() > TopAbs_REVERSED) // INTERNAL
 551     E1.Orientation( TopAbs_FORWARD );
 552   if (E2.Orientation() > TopAbs_REVERSED) // INTERNAL
 553     E2.Orientation( TopAbs_FORWARD );
 554
 555   TopoDS_Vertex V, VV1[2], VV2[2];
 556   TopExp::Vertices( E1, VV1[0], VV1[1], true );
 557   TopExp::Vertices( E2, VV2[0], VV2[1], true );
 558   if      ( VV1[1].IsSame( VV2[0] ))  { V = VV1[1]; }
 559   else if ( VV1[0].IsSame( VV2[1] ))  { V = VV1[0]; }
 560   else if ( VV1[1].IsSame( VV2[1] ))  { V = VV1[1]; E1.Reverse(); }
 561   else if ( VV1[0].IsSame( VV2[0] ))  { V = VV1[0]; E1.Reverse(); }
 562   else { return GeomAbs_C0; }
 563
 564   Standard_Real u1 = BRep_Tool::Parameter( V, E1 );
 565   Standard_Real u2 = BRep_Tool::Parameter( V, E2 );
 566   BRepAdaptor_Curve C1( E1 ), C2( E2 );
 567   Standard_Real tol = BRep_Tool::Tolerance( V );
 568   Standard_Real angTol = 2e-3;
 569   try {
 570 #if OCC_VERSION_LARGE > 0x06010000
 571     OCC_CATCH_SIGNALS;
 572 #endif
 573     return BRepLProp::Continuity(C1, C2, u1, u2, tol, angTol);
 574   }
 575   catch (Standard_Failure) {
 576   }
 577   return GeomAbs_C0;
 578 }
 579
 580 //================================================================================
 581 /*!
 582  * \brief Return the node built on a vertex
 583  * \param V - the vertex
 584  * \param meshDS - mesh
 585  * \retval const SMDS_MeshNode* - found node or NULL
 586  */
 587 //================================================================================
 588
 589 const SMDS_MeshNode* SMESH_Algo::VertexNode(const TopoDS_Vertex& V,
 590                                             const SMESHDS_Mesh*  meshDS)
 591 {
 592   if ( SMESHDS_SubMesh* sm = meshDS->MeshElements(V) ) {
 593     SMDS_NodeIteratorPtr nIt= sm->GetNodes();
 594     if (nIt->more())
 595       return nIt->next();
 596   }
 597   return 0;
 598 }
 599
 600 //=======================================================================
 601 //function : GetCommonNodes
 602 //purpose  : Return nodes common to two elements
 603 //=======================================================================
 604
 605 vector< const SMDS_MeshNode*> SMESH_Algo::GetCommonNodes(const SMDS_MeshElement* e1,
 606                                                          const SMDS_MeshElement* e2)
 607 {
 608   vector< const SMDS_MeshNode*> common;
 609   for ( int i = 0 ; i < e1->NbNodes(); ++i )
 610     if ( e2->GetNodeIndex( e1->GetNode( i )) >= 0 )
 611       common.push_back( e1->GetNode( i ));
 612   return common;
 613 }
 614
 615 //=======================================================================
 616 //function : GetMeshError
 617 //purpose  : Finds topological errors of a sub-mesh
 618 //WARNING  : 1D check is NOT implemented so far
 619 //=======================================================================
 620
 621 SMESH_Algo::EMeshError SMESH_Algo::GetMeshError(SMESH_subMesh* subMesh)
 622 {
 623   EMeshError err = MEr_OK;
 624
 625   SMESHDS_SubMesh* smDS = subMesh->GetSubMeshDS();
 626   if ( !smDS )
 627     return MEr_EMPTY;
 628
 629   switch ( subMesh->GetSubShape().ShapeType() )
 630   {
 631   case TopAbs_FACE: { // ====================== 2D =====================
 632
 633     SMDS_ElemIteratorPtr fIt = smDS->GetElements();
 634     if ( !fIt->more() )
 635       return MEr_EMPTY;
 636
 637     // We check that olny links on EDGEs encouter once, the rest links, twice
 638     set< SMESH_TLink > links;
 639     while ( fIt->more() )
 640     {
 641       const SMDS_MeshElement* f = fIt->next();
 642       int nbNodes = f->NbCornerNodes(); // ignore medium nodes
 643       for ( int i = 0; i < nbNodes; ++i )
 644       {
 645         const SMDS_MeshNode* n1 = f->GetNode( i );
 646         const SMDS_MeshNode* n2 = f->GetNode(( i+1 ) % nbNodes);
 647         std::pair< set< SMESH_TLink >::iterator, bool > it_added =
 648           links.insert( SMESH_TLink( n1, n2 ));
 649         if ( !it_added.second )
 650           // As we do NOT(!) check if mesh is manifold, we believe that a link can
 651           // encounter once or twice only (not three times), we erase a link as soon
 652           // as it encounters twice to speed up search in the <links> map.
 653           links.erase( it_added.first );
 654       }
 655     }
 656     // the links remaining in the <links> should all be on EDGE
 657     set< SMESH_TLink >::iterator linkIt = links.begin();
 658     for ( ; linkIt != links.end(); ++linkIt )
 659     {
 660       const SMESH_TLink& link = *linkIt;
 661       if ( link.node1()->GetPosition()->GetTypeOfPosition() > SMDS_TOP_EDGE ||
 662            link.node2()->GetPosition()->GetTypeOfPosition() > SMDS_TOP_EDGE )
 663         return MEr_HOLES;
 664     }
 665     // TODO: to check orientation
 666     break;
 667   }
 668   case TopAbs_SOLID: { // ====================== 3D =====================
 669
 670     SMDS_ElemIteratorPtr vIt = smDS->GetElements();
 671     if ( !vIt->more() )
 672       return MEr_EMPTY;
 673
 674     SMDS_VolumeTool vTool;
 675     while ( !vIt->more() )
 676     {
 677       if (!vTool.Set( vIt->next() ))
 678         continue; // strange
 679
 680       for ( int iF = 0; iF < vTool.NbFaces(); ++iF )
 681         if ( vTool.IsFreeFace( iF ))
 682         {
 683           int nbN = vTool.NbFaceNodes( iF );
 684           const SMDS_MeshNode** nodes =  vTool.GetFaceNodes( iF );
 685           for ( int i = 0; i < nbN; ++i )
 686             if ( nodes[i]->GetPosition()->GetTypeOfPosition() > SMDS_TOP_FACE )
 687               return MEr_HOLES;
 688         }
 689     }
 690     break;
 691   }
 692   default:;
 693   }
 694   return err;
 695 }
 696
 697 //================================================================================
 698 /*!
 699  * \brief Sets event listener to submeshes if necessary
 700  * \param subMesh - submesh where algo is set
 701  *
 702  * After being set, event listener is notified on each event of a submesh.
 703  * By default non listener is set
 704  */
 705 //================================================================================
 706
 707 void SMESH_Algo::SetEventListener(SMESH_subMesh* /*subMesh*/)
 708 {
 709 }
 710
 711 //================================================================================
 712 /*!
 713  * \brief Allow algo to do something after persistent restoration
 714  * \param subMesh - restored submesh
 715  *
 716  * This method is called only if a submesh has HYP_OK algo_state.
 717  */
 718 //================================================================================
 719
 720 void SMESH_Algo::SubmeshRestored(SMESH_subMesh* /*subMesh*/)
 721 {
 722 }
 723
 724 //================================================================================
 725 /*!
 726  * \brief Computes mesh without geometry
 727  * \param aMesh - the mesh
 728  * \param aHelper - helper that must be used for adding elements to \aaMesh
 729  * \retval bool - is a success
 730  */
 731 //================================================================================
 732
 733 bool SMESH_Algo::Compute(SMESH_Mesh & /*aMesh*/, SMESH_MesherHelper* /*aHelper*/)
 734 {
 735   return error( COMPERR_BAD_INPUT_MESH, "Mesh built on shape expected");
 736 }
 737
 738 //=======================================================================
 739 //function : CancelCompute
 740 //purpose  : Sets _computeCanceled to true. It's usage depends on
 741 //  *        implementation of a particular mesher.
 742 //=======================================================================
 743
 744 void SMESH_Algo::CancelCompute()
 745 {
 746   _computeCanceled = true;
 747   _error = COMPERR_CANCELED;
 748 }
 749
 750 //================================================================================
 751 /*!
 752  * \brief store error and comment and then return ( error == COMPERR_OK )
 753  */
 754 //================================================================================
 755
 756 bool SMESH_Algo::error(int error, const SMESH_Comment& comment)
 757 {
 758   _error   = error;
 759   _comment = comment;
 760   return ( error == COMPERR_OK );
 761 }
 762
 763 //================================================================================
 764 /*!
 765  * \brief store error and return ( error == COMPERR_OK )
 766  */
 767 //================================================================================
 768
 769 bool SMESH_Algo::error(SMESH_ComputeErrorPtr error)
 770 {
 771   if ( error ) {
 772     _error   = error->myName;
 773     _comment = error->myComment;
 774     _badInputElements = error->myBadElements;
 775     return error->IsOK();
 776   }
 777   return true;
 778 }
 779
 780 //================================================================================
 781 /*!
 782  * \brief return compute error
 783  */
 784 //================================================================================
 785
 786 SMESH_ComputeErrorPtr SMESH_Algo::GetComputeError() const
 787 {
 788   SMESH_ComputeErrorPtr err = SMESH_ComputeError::New( _error, _comment, this );
 789   // hope this method is called by only SMESH_subMesh after this->Compute()
 790   err->myBadElements.splice( err->myBadElements.end(),
 791                              (list<const SMDS_MeshElement*>&) _badInputElements );
 792   return err;
 793 }
 794
 795 //================================================================================
 796 /*!
 797  * \brief initialize compute error
 798  */
 799 //================================================================================
 800
 801 void SMESH_Algo::InitComputeError()
 802 {
 803   _error = COMPERR_OK;
 804   _comment.clear();
 805   list<const SMDS_MeshElement*>::iterator elem = _badInputElements.begin();
 806   for ( ; elem != _badInputElements.end(); ++elem )
 807     if ( (*elem)->GetID() < 1 )
 808       delete *elem;
 809   _badInputElements.clear();
 810
 811   _computeCanceled = false;
 812 }
 813
 814 //================================================================================
 815 /*!
 816  * \brief store a bad input element preventing computation,
 817  *        which may be a temporary one i.e. not residing the mesh,
 818  *        then it will be deleted by InitComputeError()
 819  */
 820 //================================================================================
 821
 822 void SMESH_Algo::addBadInputElement(const SMDS_MeshElement* elem)
 823 {
 824   if ( elem )
 825     _badInputElements.push_back( elem );
 826 }
 827
 828 //=======================================================================
 829 //function : addBadInputElements
 830 //purpose  : store a bad input elements or nodes preventing computation
 831 //=======================================================================
 832
 833 void SMESH_Algo::addBadInputElements(const SMESHDS_SubMesh* sm,
 834                                      const bool             addNodes)
 835 {
 836   if ( sm )
 837   {
 838     if ( addNodes )
 839     {
 840       SMDS_NodeIteratorPtr nIt = sm->GetNodes();
 841       while ( nIt->more() ) addBadInputElement( nIt->next() );
 842     }
 843     else
 844     {
 845       SMDS_ElemIteratorPtr eIt = sm->GetElements();
 846       while ( eIt->more() ) addBadInputElement( eIt->next() );
 847     }
 848   }
 849 }
 850
 851 //=============================================================================
 852 /*!
 853  *
 854  */
 855 //=============================================================================
 856
 857 // int SMESH_Algo::NumberOfWires(const TopoDS_Shape& S)
 858 // {
 859 //   int i = 0;
 860 //   for (TopExp_Explorer exp(S,TopAbs_WIRE); exp.More(); exp.Next())
 861 //     i++;
 862 //   return i;
 863 // }
 864
 865 //=============================================================================
 866 /*!
 867  *
 868  */
 869 //=============================================================================
 870
 871 int SMESH_Algo::NumberOfPoints(SMESH_Mesh& aMesh, const TopoDS_Wire& W)
 872 {
 873   int nbPoints = 0;
 874   for (TopExp_Explorer exp(W,TopAbs_EDGE); exp.More(); exp.Next()) {
 875     const TopoDS_Edge& E = TopoDS::Edge(exp.Current());
 876     int nb = aMesh.GetSubMesh(E)->GetSubMeshDS()->NbNodes();
 877     if(_quadraticMesh)
 878       nb = nb/2;
 879     nbPoints += nb + 1; // internal points plus 1 vertex of 2 (last point ?)
 880   }
 881   return nbPoints;
 882 }
 883
 884