src/StdMeshers/StdMeshers_Regular_1D.cxx

   1 //  SMESH SMESH : implementaion of SMESH idl descriptions
   2 //
   3 //  Copyright (C) 2003  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 //
  24 //  File   : StdMeshers_Regular_1D.cxx
  25 //           Moved here from SMESH_Regular_1D.cxx
  26 //  Author : Paul RASCLE, EDF
  27 //  Module : SMESH
  28 //  $Header$
  29
  30 #include "StdMeshers_Regular_1D.hxx"
  31 #include "StdMeshers_Distribution.hxx"
  32 #include "SMESH_Gen.hxx"
  33 #include "SMESH_Mesh.hxx"
  34 #include "SMESH_HypoFilter.hxx"
  35 #include "SMESH_subMesh.hxx"
  36
  37 #include "StdMeshers_LocalLength.hxx"
  38 #include "StdMeshers_NumberOfSegments.hxx"
  39 #include "StdMeshers_Arithmetic1D.hxx"
  40 #include "StdMeshers_StartEndLength.hxx"
  41 #include "StdMeshers_Deflection1D.hxx"
  42 #include "StdMeshers_AutomaticLength.hxx"
  43
  44 #include "SMDS_MeshElement.hxx"
  45 #include "SMDS_MeshNode.hxx"
  46 #include "SMDS_EdgePosition.hxx"
  47
  48 #include "Utils_SALOME_Exception.hxx"
  49 #include "utilities.h"
  50
  51 #include <BRep_Tool.hxx>
  52 #include <TopoDS_Edge.hxx>
  53 #include <TopoDS_Shape.hxx>
  54 #include <TopTools_ListIteratorOfListOfShape.hxx>
  55 #include <GeomAdaptor_Curve.hxx>
  56 #include <GCPnts_AbscissaPoint.hxx>
  57 #include <GCPnts_UniformAbscissa.hxx>
  58 #include <GCPnts_UniformDeflection.hxx>
  59 #include <Standard_ErrorHandler.hxx>
  60 #include <Precision.hxx>
  61 #include <Expr_GeneralExpression.hxx>
  62 #include <Expr_NamedUnknown.hxx>
  63 #include <Expr_Array1OfNamedUnknown.hxx>
  64 #include <TColStd_Array1OfReal.hxx>
  65 #include <ExprIntrp_GenExp.hxx>
  66 #include <OSD.hxx>
  67
  68 #include <string>
  69 #include <math.h>
  70
  71 using namespace std;
  72
  73 //=============================================================================
  74 /*!
  75  *
  76  */
  77 //=============================================================================
  78
  79 StdMeshers_Regular_1D::StdMeshers_Regular_1D(int hypId, int studyId,
  80         SMESH_Gen * gen):SMESH_1D_Algo(hypId, studyId, gen)
  81 {
  82         MESSAGE("StdMeshers_Regular_1D::StdMeshers_Regular_1D");
  83         _name = "Regular_1D";
  84         _shapeType = (1 << TopAbs_EDGE);
  85
  86         _compatibleHypothesis.push_back("LocalLength");
  87         _compatibleHypothesis.push_back("NumberOfSegments");
  88         _compatibleHypothesis.push_back("StartEndLength");
  89         _compatibleHypothesis.push_back("Deflection1D");
  90         _compatibleHypothesis.push_back("Arithmetic1D");
  91         _compatibleHypothesis.push_back("AutomaticLength");
  92
  93         _compatibleHypothesis.push_back("QuadraticMesh"); // auxiliary !!!
  94 }
  95
  96 //=============================================================================
  97 /*!
  98  *
  99  */
 100 //=============================================================================
 101
 102 StdMeshers_Regular_1D::~StdMeshers_Regular_1D()
 103 {
 104 }
 105
 106 //=============================================================================
 107 /*!
 108  *
 109  */
 110 //=============================================================================
 111
 112 bool StdMeshers_Regular_1D::CheckHypothesis
 113                          (SMESH_Mesh&                          aMesh,
 114                           const TopoDS_Shape&                  aShape,
 115                           SMESH_Hypothesis::Hypothesis_Status& aStatus)
 116 {
 117   _hypType = NONE;
 118   _quadraticMesh = false;
 119
 120   const bool ignoreAuxiliaryHyps = false;
 121   const list <const SMESHDS_Hypothesis * > & hyps =
 122     GetUsedHypothesis(aMesh, aShape, ignoreAuxiliaryHyps);
 123
 124   // find non-auxiliary hypothesis
 125   const SMESHDS_Hypothesis *theHyp = 0;
 126   list <const SMESHDS_Hypothesis * >::const_iterator h = hyps.begin();
 127   for ( ; h != hyps.end(); ++h ) {
 128     if ( static_cast<const SMESH_Hypothesis*>(*h)->IsAuxiliary() ) {
 129       if ( strcmp( "QuadraticMesh", (*h)->GetName() ) == 0 )
 130         _quadraticMesh = true;
 131     }
 132     else {
 133       if ( !theHyp )
 134         theHyp = *h; // use only the first non-auxiliary hypothesis
 135     }
 136   }
 137
 138   if ( !theHyp )
 139   {
 140     aStatus = SMESH_Hypothesis::HYP_MISSING;
 141     return false;  // can't work without a hypothesis
 142   }
 143
 144   string hypName = theHyp->GetName();
 145
 146   if (hypName == "LocalLength")
 147   {
 148     const StdMeshers_LocalLength * hyp =
 149       dynamic_cast <const StdMeshers_LocalLength * >(theHyp);
 150     ASSERT(hyp);
 151     _value[ BEG_LENGTH_IND ] = _value[ END_LENGTH_IND ] = hyp->GetLength();
 152     ASSERT( _value[ BEG_LENGTH_IND ] > 0 );
 153     _hypType = LOCAL_LENGTH;
 154     aStatus = SMESH_Hypothesis::HYP_OK;
 155   }
 156
 157   else if (hypName == "NumberOfSegments")
 158   {
 159     const StdMeshers_NumberOfSegments * hyp =
 160       dynamic_cast <const StdMeshers_NumberOfSegments * >(theHyp);
 161     ASSERT(hyp);
 162     _ivalue[ NB_SEGMENTS_IND  ] = hyp->GetNumberOfSegments();
 163     ASSERT( _ivalue[ NB_SEGMENTS_IND ] > 0 );
 164     _ivalue[ DISTR_TYPE_IND ] = (int) hyp->GetDistrType();
 165     switch (_ivalue[ DISTR_TYPE_IND ])
 166     {
 167     case StdMeshers_NumberOfSegments::DT_Scale:
 168       _value[ SCALE_FACTOR_IND ] = hyp->GetScaleFactor();
 169       break;
 170     case StdMeshers_NumberOfSegments::DT_TabFunc:
 171       _vvalue[ TAB_FUNC_IND ] = hyp->GetTableFunction();
 172       break;
 173     case StdMeshers_NumberOfSegments::DT_ExprFunc:
 174       _svalue[ EXPR_FUNC_IND ] = hyp->GetExpressionFunction();
 175       break;
 176     case StdMeshers_NumberOfSegments::DT_Regular:
 177       break;
 178     default:
 179       ASSERT(0);
 180       break;
 181     }
 182     if (_ivalue[ DISTR_TYPE_IND ] == StdMeshers_NumberOfSegments::DT_TabFunc ||
 183         _ivalue[ DISTR_TYPE_IND ] == StdMeshers_NumberOfSegments::DT_ExprFunc)
 184         _ivalue[ CONV_MODE_IND ] = hyp->ConversionMode();
 185     _hypType = NB_SEGMENTS;
 186     aStatus = SMESH_Hypothesis::HYP_OK;
 187   }
 188
 189   else if (hypName == "Arithmetic1D")
 190   {
 191     const StdMeshers_Arithmetic1D * hyp =
 192       dynamic_cast <const StdMeshers_Arithmetic1D * >(theHyp);
 193     ASSERT(hyp);
 194     _value[ BEG_LENGTH_IND ] = hyp->GetLength( true );
 195     _value[ END_LENGTH_IND ] = hyp->GetLength( false );
 196     ASSERT( _value[ BEG_LENGTH_IND ] > 0 && _value[ END_LENGTH_IND ] > 0 );
 197     _hypType = ARITHMETIC_1D;
 198     aStatus = SMESH_Hypothesis::HYP_OK;
 199   }
 200
 201   else if (hypName == "StartEndLength")
 202   {
 203     const StdMeshers_StartEndLength * hyp =
 204       dynamic_cast <const StdMeshers_StartEndLength * >(theHyp);
 205     ASSERT(hyp);
 206     _value[ BEG_LENGTH_IND ] = hyp->GetLength( true );
 207     _value[ END_LENGTH_IND ] = hyp->GetLength( false );
 208     ASSERT( _value[ BEG_LENGTH_IND ] > 0 && _value[ END_LENGTH_IND ] > 0 );
 209     _hypType = BEG_END_LENGTH;
 210     aStatus = SMESH_Hypothesis::HYP_OK;
 211   }
 212
 213   else if (hypName == "Deflection1D")
 214   {
 215     const StdMeshers_Deflection1D * hyp =
 216       dynamic_cast <const StdMeshers_Deflection1D * >(theHyp);
 217     ASSERT(hyp);
 218     _value[ DEFLECTION_IND ] = hyp->GetDeflection();
 219     ASSERT( _value[ DEFLECTION_IND ] > 0 );
 220     _hypType = DEFLECTION;
 221     aStatus = SMESH_Hypothesis::HYP_OK;
 222   }
 223
 224   else if (hypName == "AutomaticLength")
 225   {
 226     StdMeshers_AutomaticLength * hyp = const_cast<StdMeshers_AutomaticLength *>
 227       (dynamic_cast <const StdMeshers_AutomaticLength * >(theHyp));
 228     ASSERT(hyp);
 229     _value[ BEG_LENGTH_IND ] = _value[ END_LENGTH_IND ] = hyp->GetLength( &aMesh, aShape );
 230     ASSERT( _value[ BEG_LENGTH_IND ] > 0 );
 231     _hypType = LOCAL_LENGTH;
 232     aStatus = SMESH_Hypothesis::HYP_OK;
 233   }
 234   else
 235     aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
 236
 237   return ( _hypType != NONE );
 238 }
 239
 240 //=======================================================================
 241 //function : compensateError
 242 //purpose  : adjust theParams so that the last segment length == an
 243 //=======================================================================
 244
 245 static void compensateError(double a1, double an,
 246                             double U1, double Un,
 247                             double             length,
 248                             GeomAdaptor_Curve& C3d,
 249                             list<double> &     theParams)
 250 {
 251   int i, nPar = theParams.size();
 252   if ( a1 + an < length && nPar > 1 )
 253   {
 254     list<double>::reverse_iterator itU = theParams.rbegin();
 255     double Ul = *itU++;
 256     // dist from the last point to the edge end <Un>, it should be equal <an>
 257     double Ln = GCPnts_AbscissaPoint::Length( C3d, Ul, Un );
 258     double dLn = an - Ln; // error of <an>
 259     if ( Abs( dLn ) <= Precision::Confusion() )
 260       return;
 261     double dU = Abs( Ul - *itU ); // parametric length of the last but one segment
 262     double dUn = dLn * Abs( Un - U1 ) / length; // parametric error of <an>
 263     if ( dUn < 0.5 * dU ) { // last segment is a bit shorter than it should
 264       dUn = -dUn; // move the last parameter to the edge beginning
 265     }
 266     else {  // last segment is much shorter than it should -> remove the last param and
 267       theParams.pop_back(); nPar--; // move the rest points toward the edge end
 268       Ln = GCPnts_AbscissaPoint::Length( C3d, theParams.back(), Un );
 269       dUn = ( an - Ln ) * Abs( Un - U1 ) / length;
 270       if ( dUn < 0.5 * dU )
 271         dUn = -dUn;
 272     }
 273     if ( U1 > Un )
 274       dUn = -dUn;
 275     double q  = dUn / ( nPar - 1 );
 276     for ( itU = theParams.rbegin(), i = 1; i < nPar; itU++, i++ ) {
 277       (*itU) += dUn;
 278       dUn -= q;
 279     }
 280   }
 281 }
 282
 283 static bool computeParamByFunc(Adaptor3d_Curve& C3d, double first, double last,
 284                                double length, bool theReverse,
 285                                int nbSeg, Function& func,
 286                                list<double>& theParams)
 287 {
 288   // never do this way
 289   //OSD::SetSignal( true );
 290
 291   if( nbSeg<=0 )
 292     return false;
 293
 294   MESSAGE( "computeParamByFunc" );
 295
 296   int nbPnt = 1 + nbSeg;
 297   vector<double> x(nbPnt, 0.);
 298
 299   if( !buildDistribution( func, 0.0, 1.0, nbSeg, x, 1E-4 ) )
 300      return false;
 301
 302   MESSAGE( "Points:\n" );
 303   char buf[1024];
 304   for( int i=0; i<=nbSeg; i++ )
 305   {
 306     sprintf(  buf, "%f\n", float(x[i] ) );
 307     MESSAGE( buf );
 308   }
 309
 310
 311
 312   // apply parameters in range [0,1] to the space of the curve
 313   double prevU = first;
 314   double sign = 1.;
 315   if (theReverse)
 316   {
 317     prevU = last;
 318     sign = -1.;
 319   }
 320   for( int i = 1; i < nbSeg; i++ )
 321   {
 322     double curvLength = length * (x[i] - x[i-1]) * sign;
 323     GCPnts_AbscissaPoint Discret( C3d, curvLength, prevU );
 324     if ( !Discret.IsDone() )
 325       return false;
 326     double U = Discret.Parameter();
 327     if ( U > first && U < last )
 328       theParams.push_back( U );
 329     else
 330       return false;
 331     prevU = U;
 332   }
 333   return true;
 334 }
 335
 336 //=============================================================================
 337 /*!
 338  *
 339  */
 340 //=============================================================================
 341 bool StdMeshers_Regular_1D::computeInternalParameters(const TopoDS_Edge& theEdge,
 342                                                       list<double> &     theParams,
 343                                                       const bool         theReverse) const
 344 {
 345   theParams.clear();
 346
 347   double f, l;
 348   Handle(Geom_Curve) Curve = BRep_Tool::Curve(theEdge, f, l);
 349   GeomAdaptor_Curve C3d (Curve, f, l);
 350
 351   double length = EdgeLength(theEdge);
 352
 353   switch( _hypType )
 354   {
 355   case LOCAL_LENGTH:
 356   case NB_SEGMENTS: {
 357
 358     double eltSize = 1;
 359     if ( _hypType == LOCAL_LENGTH )
 360     {
 361       // Local Length hypothesis
 362       double nbseg = ceil(length / _value[ BEG_LENGTH_IND ]); // integer sup
 363       if (nbseg <= 0)
 364         nbseg = 1;                        // degenerated edge
 365       eltSize = length / nbseg;
 366     }
 367     else
 368     {
 369       // Number Of Segments hypothesis
 370       int NbSegm = _ivalue[ NB_SEGMENTS_IND ];
 371       if ( NbSegm < 1 )  return false;
 372       if ( NbSegm == 1 ) return true;
 373
 374       switch (_ivalue[ DISTR_TYPE_IND ])
 375       {
 376       case StdMeshers_NumberOfSegments::DT_Scale:
 377         {
 378           double scale = _value[ SCALE_FACTOR_IND ];
 379
 380           if (fabs(scale - 1.0) < Precision::Confusion()) {
 381             // special case to avoid division on zero
 382             for (int i = 1; i < NbSegm; i++) {
 383               double param = f + (l - f) * i / NbSegm;
 384               theParams.push_back( param );
 385             }
 386           } else {
 387             // general case of scale distribution
 388             if ( theReverse )
 389               scale = 1.0 / scale;
 390
 391             double alpha = pow(scale, 1.0 / (NbSegm - 1));
 392             double factor = (l - f) / (1.0 - pow(alpha, NbSegm));
 393
 394             for (int i = 1; i < NbSegm; i++) {
 395               double param = f + factor * (1.0 - pow(alpha, i));
 396               theParams.push_back( param );
 397             }
 398           }
 399           return true;
 400         }
 401         break;
 402       case StdMeshers_NumberOfSegments::DT_TabFunc:
 403         {
 404           FunctionTable func(_vvalue[ TAB_FUNC_IND ], _ivalue[ CONV_MODE_IND ]);
 405           return computeParamByFunc(C3d, f, l, length, theReverse,
 406                                     _ivalue[ NB_SEGMENTS_IND ], func,
 407                                     theParams);
 408         }
 409         break;
 410       case StdMeshers_NumberOfSegments::DT_ExprFunc:
 411         {
 412           FunctionExpr func(_svalue[ EXPR_FUNC_IND ].c_str(), _ivalue[ CONV_MODE_IND ]);
 413           return computeParamByFunc(C3d, f, l, length, theReverse,
 414                                     _ivalue[ NB_SEGMENTS_IND ], func,
 415                                     theParams);
 416         }
 417         break;
 418       case StdMeshers_NumberOfSegments::DT_Regular:
 419         eltSize = length / _ivalue[ NB_SEGMENTS_IND ];
 420         break;
 421       default:
 422         return false;
 423       }
 424     }
 425     GCPnts_UniformAbscissa Discret(C3d, eltSize, f, l);
 426     if ( !Discret.IsDone() )
 427       return false;
 428
 429     int NbPoints = Discret.NbPoints();
 430     for ( int i = 2; i < NbPoints; i++ )
 431     {
 432       double param = Discret.Parameter(i);
 433       theParams.push_back( param );
 434     }
 435     compensateError( eltSize, eltSize, f, l, length, C3d, theParams ); // for PAL9899
 436     return true;
 437   }
 438
 439   case BEG_END_LENGTH: {
 440
 441     // geometric progression: SUM(n) = ( a1 - an * q ) / ( 1 - q ) = length
 442
 443     double a1 = _value[ BEG_LENGTH_IND ];
 444     double an = _value[ END_LENGTH_IND ];
 445     double q  = ( length - a1 ) / ( length - an );
 446
 447     double U1 = theReverse ? l : f;
 448     double Un = theReverse ? f : l;
 449     double param = U1;
 450     double eltSize = theReverse ? -a1 : a1;
 451     while ( 1 ) {
 452       // computes a point on a curve <C3d> at the distance <eltSize>
 453       // from the point of parameter <param>.
 454       GCPnts_AbscissaPoint Discret( C3d, eltSize, param );
 455       if ( !Discret.IsDone() ) break;
 456       param = Discret.Parameter();
 457       if ( param > f && param < l )
 458         theParams.push_back( param );
 459       else
 460         break;
 461       eltSize *= q;
 462     }
 463     compensateError( a1, an, U1, Un, length, C3d, theParams );
 464     return true;
 465   }
 466
 467   case ARITHMETIC_1D: {
 468
 469     // arithmetic progression: SUM(n) = ( an - a1 + q ) * ( a1 + an ) / ( 2 * q ) = length
 470
 471     double a1 = _value[ BEG_LENGTH_IND ];
 472     double an = _value[ END_LENGTH_IND ];
 473
 474     double  q = ( an - a1 ) / ( 2 *length/( a1 + an ) - 1 );
 475     int     n = int( 1 + ( an - a1 ) / q );
 476
 477     double U1 = theReverse ? l : f;
 478     double Un = theReverse ? f : l;
 479     double param = U1;
 480     double eltSize = a1;
 481     if ( theReverse ) {
 482       eltSize = -eltSize;
 483       q = -q;
 484     }
 485     while ( n-- > 0 && eltSize * ( Un - U1 ) > 0 ) {
 486       // computes a point on a curve <C3d> at the distance <eltSize>
 487       // from the point of parameter <param>.
 488       GCPnts_AbscissaPoint Discret( C3d, eltSize, param );
 489       if ( !Discret.IsDone() ) break;
 490       param = Discret.Parameter();
 491       if ( param > f && param < l )
 492         theParams.push_back( param );
 493       else
 494         break;
 495       eltSize += q;
 496     }
 497     compensateError( a1, an, U1, Un, length, C3d, theParams );
 498
 499     return true;
 500   }
 501
 502   case DEFLECTION: {
 503
 504     GCPnts_UniformDeflection Discret(C3d, _value[ DEFLECTION_IND ], f, l, true);
 505     if ( !Discret.IsDone() )
 506       return false;
 507
 508     int NbPoints = Discret.NbPoints();
 509     for ( int i = 2; i < NbPoints; i++ )
 510     {
 511       double param = Discret.Parameter(i);
 512       theParams.push_back( param );
 513     }
 514     return true;
 515
 516   }
 517
 518   default:;
 519   }
 520
 521   return false;
 522 }
 523
 524 //=============================================================================
 525 /*!
 526  *
 527  */
 528 //=============================================================================
 529
 530 bool StdMeshers_Regular_1D::Compute(SMESH_Mesh & aMesh, const TopoDS_Shape & aShape)
 531 {
 532   MESSAGE("StdMeshers_Regular_1D::Compute");
 533
 534   if ( _hypType == NONE )
 535     return false;
 536
 537   SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
 538   aMesh.GetSubMesh(aShape);
 539
 540   const TopoDS_Edge & EE = TopoDS::Edge(aShape);
 541   TopoDS_Edge E = TopoDS::Edge(EE.Oriented(TopAbs_FORWARD));
 542   int shapeID = meshDS->ShapeToIndex( E );
 543
 544   double f, l;
 545   Handle(Geom_Curve) Curve = BRep_Tool::Curve(E, f, l);
 546
 547   TopoDS_Vertex VFirst, VLast;
 548   TopExp::Vertices(E, VFirst, VLast);   // Vfirst corresponds to f and Vlast to l
 549
 550   ASSERT(!VFirst.IsNull());
 551   SMDS_NodeIteratorPtr lid= aMesh.GetSubMesh(VFirst)->GetSubMeshDS()->GetNodes();
 552   if (!lid->more())
 553   {
 554     MESSAGE (" NO NODE BUILT ON VERTEX ");
 555     return false;
 556   }
 557   const SMDS_MeshNode * idFirst = lid->next();
 558
 559   ASSERT(!VLast.IsNull());
 560   lid=aMesh.GetSubMesh(VLast)->GetSubMeshDS()->GetNodes();
 561   if (!lid->more()) {
 562     MESSAGE (" NO NODE BUILT ON VERTEX ");
 563     return false;
 564   }
 565   const SMDS_MeshNode * idLast = lid->next();
 566
 567   if (!Curve.IsNull()) {
 568     list< double > params;
 569     bool reversed = false;
 570     if ( !_mainEdge.IsNull() )
 571       reversed = aMesh.IsReversedInChain( EE, _mainEdge );
 572     try {
 573       if ( ! computeInternalParameters( E, params, reversed )) {
 574         //cout << "computeInternalParameters() failed" <<endl;
 575         return false;
 576       }
 577     }
 578     catch ( Standard_Failure ) {
 579       //cout << "computeInternalParameters() failed, Standard_Failure" <<endl;
 580       return false;
 581     }
 582
 583     // edge extrema (indexes : 1 & NbPoints) already in SMDS (TopoDS_Vertex)
 584     // only internal nodes receive an edge position with param on curve
 585
 586     const SMDS_MeshNode * idPrev = idFirst;
 587     double parPrev = f;
 588     double parLast = l;
 589 //     if(reversed) {
 590 //       parPrev = l;
 591 //       parLast = f;
 592 //     }
 593
 594     for (list<double>::iterator itU = params.begin(); itU != params.end(); itU++) {
 595       double param = *itU;
 596       gp_Pnt P = Curve->Value(param);
 597
 598       //Add the Node in the DataStructure
 599       SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
 600       meshDS->SetNodeOnEdge(node, shapeID, param);
 601
 602       if(_quadraticMesh) {
 603         // create medium node
 604         double prm = ( parPrev + param )/2;
 605         gp_Pnt PM = Curve->Value(prm);
 606         SMDS_MeshNode * NM = meshDS->AddNode(PM.X(), PM.Y(), PM.Z());
 607         meshDS->SetNodeOnEdge(NM, shapeID, prm);
 608         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node, NM);
 609         meshDS->SetMeshElementOnShape(edge, shapeID);
 610       }
 611       else {
 612         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node);
 613         meshDS->SetMeshElementOnShape(edge, shapeID);
 614       }
 615
 616       idPrev = node;
 617       parPrev = param;
 618     }
 619     if(_quadraticMesh) {
 620       double prm = ( parPrev + parLast )/2;
 621       gp_Pnt PM = Curve->Value(prm);
 622       SMDS_MeshNode * NM = meshDS->AddNode(PM.X(), PM.Y(), PM.Z());
 623       meshDS->SetNodeOnEdge(NM, shapeID, prm);
 624       SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, idLast, NM);
 625       meshDS->SetMeshElementOnShape(edge, shapeID);
 626     }
 627     else {
 628       SMDS_MeshEdge* edge = meshDS->AddEdge(idPrev, idLast);
 629       meshDS->SetMeshElementOnShape(edge, shapeID);
 630     }
 631   }
 632   else {
 633     // Edge is a degenerated Edge : We put n = 5 points on the edge.
 634     const int NbPoints = 5;
 635     BRep_Tool::Range(E, f, l);
 636     double du = (l - f) / (NbPoints - 1);
 637     //MESSAGE("************* Degenerated edge! *****************");
 638
 639     TopoDS_Vertex V1, V2;
 640     TopExp::Vertices(E, V1, V2);
 641     gp_Pnt P = BRep_Tool::Pnt(V1);
 642
 643     const SMDS_MeshNode * idPrev = idFirst;
 644     for (int i = 2; i < NbPoints; i++) {
 645       double param = f + (i - 1) * du;
 646       SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
 647       if(_quadraticMesh) {
 648         // create medium node
 649         double prm = param - du/2.;
 650         SMDS_MeshNode * NM = meshDS->AddNode(P.X(), P.Y(), P.Z());
 651         meshDS->SetNodeOnEdge(NM, shapeID, prm);
 652         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node, NM);
 653         meshDS->SetMeshElementOnShape(edge, shapeID);
 654       }
 655       else {
 656         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node);
 657         meshDS->SetMeshElementOnShape(edge, shapeID);
 658       }
 659       meshDS->SetNodeOnEdge(node, shapeID, param);
 660       idPrev = node;
 661     }
 662     if(_quadraticMesh) {
 663       // create medium node
 664       double prm = l - du/2.;
 665       SMDS_MeshNode * NM = meshDS->AddNode(P.X(), P.Y(), P.Z());
 666       meshDS->SetNodeOnEdge(NM, shapeID, prm);
 667       SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, idLast, NM);
 668       meshDS->SetMeshElementOnShape(edge, shapeID);
 669     }
 670     else {
 671       SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, idLast);
 672       meshDS->SetMeshElementOnShape(edge, shapeID);
 673     }
 674   }
 675   return true;
 676 }
 677
 678 //=============================================================================
 679 /*!
 680  *  See comments in SMESH_Algo.cxx
 681  */
 682 //=============================================================================
 683
 684 const list <const SMESHDS_Hypothesis *> &
 685 StdMeshers_Regular_1D::GetUsedHypothesis(SMESH_Mesh &         aMesh,
 686                                          const TopoDS_Shape & aShape,
 687                                          const bool           ignoreAuxiliary)
 688 {
 689   _usedHypList.clear();
 690   _mainEdge.Nullify();
 691
 692   SMESH_HypoFilter auxiliaryFilter, compatibleFilter;
 693   auxiliaryFilter.Init( SMESH_HypoFilter::IsAuxiliary() );
 694   const bool ignoreAux = true;
 695   InitCompatibleHypoFilter( compatibleFilter, ignoreAux );
 696
 697   // get non-auxiliary assigned to aShape
 698   int nbHyp = aMesh.GetHypotheses( aShape, compatibleFilter, _usedHypList, false );
 699
 700   if (nbHyp == 0)
 701   {
 702     // Check, if propagated from some other edge
 703     if (aShape.ShapeType() == TopAbs_EDGE &&
 704         aMesh.IsPropagatedHypothesis(aShape, _mainEdge))
 705     {
 706       // Propagation of 1D hypothesis from <aMainEdge> on this edge;
 707       // get non-auxiliary assigned to _mainEdge
 708       nbHyp = aMesh.GetHypotheses( _mainEdge, compatibleFilter, _usedHypList, true );
 709     }
 710   }
 711
 712   if (nbHyp == 0) // nothing propagated nor assigned to aShape
 713   {
 714     SMESH_Algo::GetUsedHypothesis( aMesh, aShape, ignoreAuxiliary );
 715     nbHyp = _usedHypList.size();
 716   }
 717   else
 718   {
 719     // get auxiliary hyps from aShape
 720     aMesh.GetHypotheses( aShape, auxiliaryFilter, _usedHypList, true );
 721   }
 722   if ( nbHyp > 1 && ignoreAuxiliary )
 723     _usedHypList.clear(); //only one compatible non-auxiliary hypothesis allowed
 724
 725   return _usedHypList;
 726 }
 727
 728 //=============================================================================
 729 /*!
 730  *
 731  */
 732 //=============================================================================
 733
 734 ostream & StdMeshers_Regular_1D::SaveTo(ostream & save)
 735 {
 736   return save;
 737 }
 738
 739 //=============================================================================
 740 /*!
 741  *
 742  */
 743 //=============================================================================
 744
 745 istream & StdMeshers_Regular_1D::LoadFrom(istream & load)
 746 {
 747   return load;
 748 }
 749
 750 //=============================================================================
 751 /*!
 752  *
 753  */
 754 //=============================================================================
 755
 756 ostream & operator <<(ostream & save, StdMeshers_Regular_1D & hyp)
 757 {
 758   return hyp.SaveTo( save );
 759 }
 760
 761 //=============================================================================
 762 /*!
 763  *
 764  */
 765 //=============================================================================
 766
 767 istream & operator >>(istream & load, StdMeshers_Regular_1D & hyp)
 768 {
 769   return hyp.LoadFrom( load );
 770 }