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
  33 #include "StdMeshers_LocalLength.hxx"
  34 #include "StdMeshers_NumberOfSegments.hxx"
  35 #include "StdMeshers_Arithmetic1D.hxx"
  36 #include "StdMeshers_StartEndLength.hxx"
  37 #include "StdMeshers_Deflection1D.hxx"
  38 #include "StdMeshers_AutomaticLength.hxx"
  39 #include "StdMeshers_SegmentLengthAroundVertex.hxx"
  40
  41 #include "SMESH_Gen.hxx"
  42 #include "SMESH_Mesh.hxx"
  43 #include "SMESH_HypoFilter.hxx"
  44 #include "SMESH_subMesh.hxx"
  45 #include "SMESH_subMeshEventListener.hxx"
  46
  47 #include "SMDS_MeshElement.hxx"
  48 #include "SMDS_MeshNode.hxx"
  49
  50 #include "Utils_SALOME_Exception.hxx"
  51 #include "utilities.h"
  52
  53 #include <BRepAdaptor_Curve.hxx>
  54 #include <BRep_Tool.hxx>
  55 #include <TopoDS_Edge.hxx>
  56 #include <TopExp_Explorer.hxx>
  57 #include <GCPnts_AbscissaPoint.hxx>
  58 #include <GCPnts_UniformAbscissa.hxx>
  59 #include <GCPnts_UniformDeflection.hxx>
  60 #include <Precision.hxx>
  61
  62 #include <Standard_ErrorHandler.hxx>
  63 #include <Standard_Failure.hxx>
  64
  65 #include <string>
  66 //#include <math.h>
  67
  68 using namespace std;
  69
  70 //=============================================================================
  71 /*!
  72  *
  73  */
  74 //=============================================================================
  75
  76 StdMeshers_Regular_1D::StdMeshers_Regular_1D(int hypId, int studyId,
  77         SMESH_Gen * gen):SMESH_1D_Algo(hypId, studyId, gen)
  78 {
  79         MESSAGE("StdMeshers_Regular_1D::StdMeshers_Regular_1D");
  80         _name = "Regular_1D";
  81         _shapeType = (1 << TopAbs_EDGE);
  82
  83         _compatibleHypothesis.push_back("LocalLength");
  84         _compatibleHypothesis.push_back("NumberOfSegments");
  85         _compatibleHypothesis.push_back("StartEndLength");
  86         _compatibleHypothesis.push_back("Deflection1D");
  87         _compatibleHypothesis.push_back("Arithmetic1D");
  88         _compatibleHypothesis.push_back("AutomaticLength");
  89
  90         _compatibleHypothesis.push_back("QuadraticMesh"); // auxiliary !!!
  91 }
  92
  93 //=============================================================================
  94 /*!
  95  *
  96  */
  97 //=============================================================================
  98
  99 StdMeshers_Regular_1D::~StdMeshers_Regular_1D()
 100 {
 101 }
 102
 103 //=============================================================================
 104 /*!
 105  *
 106  */
 107 //=============================================================================
 108
 109 bool StdMeshers_Regular_1D::CheckHypothesis
 110                          (SMESH_Mesh&                          aMesh,
 111                           const TopoDS_Shape&                  aShape,
 112                           SMESH_Hypothesis::Hypothesis_Status& aStatus)
 113 {
 114   _hypType = NONE;
 115   _quadraticMesh = false;
 116
 117   const bool ignoreAuxiliaryHyps = false;
 118   const list <const SMESHDS_Hypothesis * > & hyps =
 119     GetUsedHypothesis(aMesh, aShape, ignoreAuxiliaryHyps);
 120
 121   // find non-auxiliary hypothesis
 122   const SMESHDS_Hypothesis *theHyp = 0;
 123   list <const SMESHDS_Hypothesis * >::const_iterator h = hyps.begin();
 124   for ( ; h != hyps.end(); ++h ) {
 125     if ( static_cast<const SMESH_Hypothesis*>(*h)->IsAuxiliary() ) {
 126       if ( strcmp( "QuadraticMesh", (*h)->GetName() ) == 0 )
 127         _quadraticMesh = true;
 128     }
 129     else {
 130       if ( !theHyp )
 131         theHyp = *h; // use only the first non-auxiliary hypothesis
 132     }
 133   }
 134
 135   if ( !theHyp )
 136   {
 137     aStatus = SMESH_Hypothesis::HYP_MISSING;
 138     return false;  // can't work without a hypothesis
 139   }
 140
 141   string hypName = theHyp->GetName();
 142
 143   if (hypName == "LocalLength")
 144   {
 145     const StdMeshers_LocalLength * hyp =
 146       dynamic_cast <const StdMeshers_LocalLength * >(theHyp);
 147     ASSERT(hyp);
 148     _value[ BEG_LENGTH_IND ] = _value[ END_LENGTH_IND ] = hyp->GetLength();
 149     ASSERT( _value[ BEG_LENGTH_IND ] > 0 );
 150     _hypType = LOCAL_LENGTH;
 151     aStatus = SMESH_Hypothesis::HYP_OK;
 152   }
 153
 154   else if (hypName == "NumberOfSegments")
 155   {
 156     const StdMeshers_NumberOfSegments * hyp =
 157       dynamic_cast <const StdMeshers_NumberOfSegments * >(theHyp);
 158     ASSERT(hyp);
 159     _ivalue[ NB_SEGMENTS_IND  ] = hyp->GetNumberOfSegments();
 160     ASSERT( _ivalue[ NB_SEGMENTS_IND ] > 0 );
 161     _ivalue[ DISTR_TYPE_IND ] = (int) hyp->GetDistrType();
 162     switch (_ivalue[ DISTR_TYPE_IND ])
 163     {
 164     case StdMeshers_NumberOfSegments::DT_Scale:
 165       _value[ SCALE_FACTOR_IND ] = hyp->GetScaleFactor();
 166       break;
 167     case StdMeshers_NumberOfSegments::DT_TabFunc:
 168       _vvalue[ TAB_FUNC_IND ] = hyp->GetTableFunction();
 169       break;
 170     case StdMeshers_NumberOfSegments::DT_ExprFunc:
 171       _svalue[ EXPR_FUNC_IND ] = hyp->GetExpressionFunction();
 172       break;
 173     case StdMeshers_NumberOfSegments::DT_Regular:
 174       break;
 175     default:
 176       ASSERT(0);
 177       break;
 178     }
 179     if (_ivalue[ DISTR_TYPE_IND ] == StdMeshers_NumberOfSegments::DT_TabFunc ||
 180         _ivalue[ DISTR_TYPE_IND ] == StdMeshers_NumberOfSegments::DT_ExprFunc)
 181         _ivalue[ CONV_MODE_IND ] = hyp->ConversionMode();
 182     _hypType = NB_SEGMENTS;
 183     aStatus = SMESH_Hypothesis::HYP_OK;
 184   }
 185
 186   else if (hypName == "Arithmetic1D")
 187   {
 188     const StdMeshers_Arithmetic1D * hyp =
 189       dynamic_cast <const StdMeshers_Arithmetic1D * >(theHyp);
 190     ASSERT(hyp);
 191     _value[ BEG_LENGTH_IND ] = hyp->GetLength( true );
 192     _value[ END_LENGTH_IND ] = hyp->GetLength( false );
 193     ASSERT( _value[ BEG_LENGTH_IND ] > 0 && _value[ END_LENGTH_IND ] > 0 );
 194     _hypType = ARITHMETIC_1D;
 195     aStatus = SMESH_Hypothesis::HYP_OK;
 196   }
 197
 198   else if (hypName == "StartEndLength")
 199   {
 200     const StdMeshers_StartEndLength * hyp =
 201       dynamic_cast <const StdMeshers_StartEndLength * >(theHyp);
 202     ASSERT(hyp);
 203     _value[ BEG_LENGTH_IND ] = hyp->GetLength( true );
 204     _value[ END_LENGTH_IND ] = hyp->GetLength( false );
 205     ASSERT( _value[ BEG_LENGTH_IND ] > 0 && _value[ END_LENGTH_IND ] > 0 );
 206     _hypType = BEG_END_LENGTH;
 207     aStatus = SMESH_Hypothesis::HYP_OK;
 208   }
 209
 210   else if (hypName == "Deflection1D")
 211   {
 212     const StdMeshers_Deflection1D * hyp =
 213       dynamic_cast <const StdMeshers_Deflection1D * >(theHyp);
 214     ASSERT(hyp);
 215     _value[ DEFLECTION_IND ] = hyp->GetDeflection();
 216     ASSERT( _value[ DEFLECTION_IND ] > 0 );
 217     _hypType = DEFLECTION;
 218     aStatus = SMESH_Hypothesis::HYP_OK;
 219   }
 220
 221   else if (hypName == "AutomaticLength")
 222   {
 223     StdMeshers_AutomaticLength * hyp = const_cast<StdMeshers_AutomaticLength *>
 224       (dynamic_cast <const StdMeshers_AutomaticLength * >(theHyp));
 225     ASSERT(hyp);
 226     _value[ BEG_LENGTH_IND ] = _value[ END_LENGTH_IND ] = hyp->GetLength( &aMesh, aShape );
 227     ASSERT( _value[ BEG_LENGTH_IND ] > 0 );
 228     _hypType = LOCAL_LENGTH;
 229     aStatus = SMESH_Hypothesis::HYP_OK;
 230   }
 231   else
 232     aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
 233
 234   return ( _hypType != NONE );
 235 }
 236
 237 static bool computeParamByFunc(Adaptor3d_Curve& C3d, double first, double last,
 238                                double length, bool theReverse,
 239                                int nbSeg, Function& func,
 240                                list<double>& theParams)
 241 {
 242   // never do this way
 243   //OSD::SetSignal( true );
 244
 245   if( nbSeg<=0 )
 246     return false;
 247
 248   MESSAGE( "computeParamByFunc" );
 249
 250   int nbPnt = 1 + nbSeg;
 251   vector<double> x(nbPnt, 0.);
 252
 253   if( !buildDistribution( func, 0.0, 1.0, nbSeg, x, 1E-4 ) )
 254      return false;
 255
 256   MESSAGE( "Points:\n" );
 257   char buf[1024];
 258   for( int i=0; i<=nbSeg; i++ )
 259   {
 260     sprintf(  buf, "%f\n", float(x[i] ) );
 261     MESSAGE( buf );
 262   }
 263
 264
 265
 266   // apply parameters in range [0,1] to the space of the curve
 267   double prevU = first;
 268   double sign = 1.;
 269   if (theReverse)
 270   {
 271     prevU = last;
 272     sign = -1.;
 273   }
 274   for( int i = 1; i < nbSeg; i++ )
 275   {
 276     double curvLength = length * (x[i] - x[i-1]) * sign;
 277     GCPnts_AbscissaPoint Discret( C3d, curvLength, prevU );
 278     if ( !Discret.IsDone() )
 279       return false;
 280     double U = Discret.Parameter();
 281     if ( U > first && U < last )
 282       theParams.push_back( U );
 283     else
 284       return false;
 285     prevU = U;
 286   }
 287   return true;
 288 }
 289
 290
 291 //================================================================================
 292 /*!
 293  * \brief adjust internal node parameters so that the last segment length == an
 294   * \param a1 - the first segment length
 295   * \param an - the last segment length
 296   * \param U1 - the first edge parameter
 297   * \param Un - the last edge parameter
 298   * \param length - the edge length
 299   * \param C3d - the edge curve
 300   * \param theParams - internal node parameters to adjust
 301   * \param adjustNeighbors2an - to adjust length of segments next to the last one
 302   *  and not to remove parameters
 303  */
 304 //================================================================================
 305
 306 static void compensateError(double a1, double an,
 307                             double U1, double Un,
 308                             double            length,
 309                             Adaptor3d_Curve&  C3d,
 310                             list<double> &    theParams,
 311                             bool              adjustNeighbors2an = false)
 312 {
 313   int i, nPar = theParams.size();
 314   if ( a1 + an < length && nPar > 1 )
 315   {
 316     bool reverse = ( U1 > Un );
 317     GCPnts_AbscissaPoint Discret(C3d, reverse ? an : -an, Un);
 318     if ( !Discret.IsDone() )
 319       return;
 320     double Utgt = Discret.Parameter(); // target value of the last parameter
 321     list<double>::reverse_iterator itU = theParams.rbegin();
 322     double Ul = *itU++; // real value of the last parameter
 323     double dUn = Utgt - Ul; // parametric error of <an>
 324     if ( Abs(dUn) <= Precision::Confusion() )
 325       return;
 326     double dU = Abs( Ul - *itU ); // parametric length of the last but one segment
 327     if ( adjustNeighbors2an || Abs(dUn) < 0.5 * dU ) { // last segment is a bit shorter than it should
 328       // move the last parameter to the edge beginning
 329     }
 330     else {  // last segment is much shorter than it should -> remove the last param and
 331       theParams.pop_back(); nPar--; // move the rest points toward the edge end
 332       dUn = Utgt - theParams.back();
 333     }
 334
 335     double q  = dUn / ( nPar - 1 );
 336     if ( !adjustNeighbors2an ) {
 337       for ( itU = theParams.rbegin(), i = 1; i < nPar; itU++, i++ ) {
 338         (*itU) += dUn;
 339         dUn -= q;
 340       }
 341     }
 342     else {
 343       theParams.back() += dUn;
 344       double sign = reverse ? -1 : 1;
 345       double prevU = theParams.back();
 346       itU = theParams.rbegin();
 347       for ( ++itU, i = 2; i < nPar; ++itU, i++ ) {
 348         double newU = *itU + dUn;
 349         if ( newU*sign < prevU*sign ) {
 350           prevU = *itU = newU;
 351           dUn -= q;
 352         }
 353         else { // set U between prevU and next valid param
 354           list<double>::reverse_iterator itU2 = itU;
 355           ++itU2;
 356           int nb = 2;
 357           while ( (*itU2)*sign > prevU*sign ) {
 358             ++itU2; ++nb;
 359           }
 360           dU = ( *itU2 - prevU ) / nb;
 361           while ( itU != itU2 ) {
 362             *itU += dU; ++itU;
 363           }
 364           break;
 365         }
 366       }
 367     }
 368   }
 369 }
 370
 371 //================================================================================
 372 /*!
 373  * \brief Class used to clean mesh on edges when 0D hyp modified.
 374  * Common approach doesn't work when 0D algo is missing because the 0D hyp is
 375  * considered as not participating in computation whereas it is used by 1D algo.
 376  */
 377 //================================================================================
 378
 379 struct VertexEventListener : public SMESH_subMeshEventListener
 380 {
 381   VertexEventListener():SMESH_subMeshEventListener(0) // won't be deleted by submesh
 382   {}
 383   /*!
 384    * \brief Clean mesh on edges
 385    * \param event - algo_event or compute_event itself (of SMESH_subMesh)
 386    * \param eventType - ALGO_EVENT or COMPUTE_EVENT (of SMESH_subMesh)
 387    * \param subMesh - the submesh where the event occures
 388    */
 389   void ProcessEvent(const int event, const int eventType, SMESH_subMesh* subMesh,
 390                     EventListenerData*, SMESH_Hypothesis*)
 391   {
 392     if ( eventType == SMESH_subMesh::ALGO_EVENT) // all algo events
 393     {
 394       subMesh->ComputeStateEngine( SMESH_subMesh::MODIF_ALGO_STATE );
 395     }
 396   }
 397 }; // struct VertexEventListener
 398
 399 //=============================================================================
 400 /*!
 401  * \brief Sets event listener to vertex submeshes
 402  * \param subMesh - submesh where algo is set
 403  *
 404  * This method is called when a submesh gets HYP_OK algo_state.
 405  * After being set, event listener is notified on each event of a submesh.
 406  */
 407 //=============================================================================
 408
 409 void StdMeshers_Regular_1D::SetEventListener(SMESH_subMesh* subMesh)
 410 {
 411   static VertexEventListener listener;
 412   const map < int, SMESH_subMesh * >& vSMs = subMesh->DependsOn();
 413   map < int, SMESH_subMesh * >::const_iterator itsub;
 414   for (itsub = vSMs.begin(); itsub != vSMs.end(); itsub++)
 415   {
 416     subMesh->SetEventListener( &listener, 0, itsub->second );
 417   }
 418 }
 419
 420 //=============================================================================
 421 /*!
 422  * \brief Do nothing
 423  * \param subMesh - restored submesh
 424  *
 425  * This method is called only if a submesh has HYP_OK algo_state.
 426  */
 427 //=============================================================================
 428
 429 void StdMeshers_Regular_1D::SubmeshRestored(SMESH_subMesh* subMesh)
 430 {
 431 }
 432
 433 //=============================================================================
 434 /*!
 435  * \brief Return StdMeshers_SegmentLengthAroundVertex assigned to vertex
 436  */
 437 //=============================================================================
 438
 439 const StdMeshers_SegmentLengthAroundVertex*
 440 StdMeshers_Regular_1D::getVertexHyp(SMESH_Mesh &          theMesh,
 441                                     const TopoDS_Vertex & theV)
 442 {
 443   static SMESH_HypoFilter filter( SMESH_HypoFilter::HasName("SegmentLengthAroundVertex"));
 444   const SMESH_Hypothesis * hyp = theMesh.GetHypothesis( theV, filter, true );
 445   return static_cast<const StdMeshers_SegmentLengthAroundVertex*>( hyp );
 446 }
 447
 448 //================================================================================
 449 /*!
 450  * \brief Tune parameters to fit "SegmentLengthAroundVertex" hypothesis
 451   * \param theC3d - wire curve
 452   * \param theLength - curve length
 453   * \param theParameters - internal nodes parameters to modify
 454   * \param theVf - 1st vertex
 455   * \param theVl - 2nd vertex
 456  */
 457 //================================================================================
 458
 459 void StdMeshers_Regular_1D::redistributeNearVertices (SMESH_Mesh &          theMesh,
 460                                                       Adaptor3d_Curve &     theC3d,
 461                                                       double                theLength,
 462                                                       std::list< double > & theParameters,
 463                                                       const TopoDS_Vertex & theVf,
 464                                                       const TopoDS_Vertex & theVl) const
 465 {
 466   double f = theC3d.FirstParameter(), l = theC3d.LastParameter();
 467   int nPar = theParameters.size();
 468   for ( int isEnd1 = 0; isEnd1 < 2; ++isEnd1 )
 469   {
 470     const TopoDS_Vertex & V = isEnd1 ? theVf : theVl;
 471     const StdMeshers_SegmentLengthAroundVertex* hyp = getVertexHyp (theMesh, V );
 472     if ( hyp ) {
 473       double vertexLength = hyp->GetLength();
 474       if ( vertexLength > theLength / 2.0 )
 475         continue;
 476       if ( isEnd1 ) { // to have a segment of interest at end of theParameters
 477         theParameters.reverse();
 478         std::swap( f, l );
 479       }
 480       if ( _hypType == NB_SEGMENTS )
 481       {
 482         compensateError(0, vertexLength, f, l, theLength, theC3d, theParameters, true );
 483       }
 484       else if ( nPar <= 3 )
 485       {
 486         if ( !isEnd1 )
 487           vertexLength = -vertexLength;
 488         GCPnts_AbscissaPoint Discret(theC3d, vertexLength, l);
 489         if ( Discret.IsDone() ) {
 490           if ( nPar == 0 )
 491             theParameters.push_back( Discret.Parameter());
 492           else {
 493             double L = GCPnts_AbscissaPoint::Length( theC3d, theParameters.back(), l);
 494             if ( vertexLength < L / 2.0 )
 495               theParameters.push_back( Discret.Parameter());
 496             else
 497               compensateError(0, vertexLength, f, l, theLength, theC3d, theParameters, true );
 498           }
 499         }
 500       }
 501       else
 502       {
 503         // recompute params between the last segment and a middle one.
 504         // find size of a middle segment
 505         int nHalf = ( nPar-1 ) / 2;
 506         list< double >::reverse_iterator itU = theParameters.rbegin();
 507         std::advance( itU, nHalf );
 508         double Um = *itU++;
 509         double Lm = GCPnts_AbscissaPoint::Length( theC3d, Um, *itU);
 510         double L = GCPnts_AbscissaPoint::Length( theC3d, *itU, l);
 511         StdMeshers_Regular_1D algo( *this );
 512         algo._hypType = BEG_END_LENGTH;
 513         algo._value[ BEG_LENGTH_IND ] = Lm;
 514         algo._value[ END_LENGTH_IND ] = vertexLength;
 515         double from = *itU, to = l;
 516         if ( isEnd1 ) {
 517           std::swap( from, to );
 518           std::swap( algo._value[ BEG_LENGTH_IND ], algo._value[ END_LENGTH_IND ]);
 519         }
 520         list<double> params;
 521         if ( algo.computeInternalParameters( theC3d, L, from, to, params, false ))
 522         {
 523           if ( isEnd1 ) params.reverse();
 524           while ( 1 + nHalf-- )
 525             theParameters.pop_back();
 526           theParameters.splice( theParameters.end(), params );
 527         }
 528         else
 529         {
 530           compensateError(0, vertexLength, f, l, theLength, theC3d, theParameters, true );
 531         }
 532       }
 533       if ( isEnd1 )
 534         theParameters.reverse();
 535     }
 536   }
 537 }
 538
 539 //=============================================================================
 540 /*!
 541  *
 542  */
 543 //=============================================================================
 544 bool StdMeshers_Regular_1D::computeInternalParameters(Adaptor3d_Curve& theC3d,
 545                                                       double           theLength,
 546                                                       double           theFirstU,
 547                                                       double           theLastU,
 548                                                       list<double> &   theParams,
 549                                                       const bool       theReverse) const
 550 {
 551   theParams.clear();
 552
 553   double f = theFirstU, l = theLastU;
 554
 555   switch( _hypType )
 556   {
 557   case LOCAL_LENGTH:
 558   case NB_SEGMENTS: {
 559
 560     double eltSize = 1;
 561     if ( _hypType == LOCAL_LENGTH )
 562     {
 563       // Local Length hypothesis
 564       double nbseg = ceil(theLength / _value[ BEG_LENGTH_IND ]); // integer sup
 565       if (nbseg <= 0)
 566         nbseg = 1;                        // degenerated edge
 567       eltSize = theLength / nbseg;
 568     }
 569     else
 570     {
 571       // Number Of Segments hypothesis
 572       int NbSegm = _ivalue[ NB_SEGMENTS_IND ];
 573       if ( NbSegm < 1 )  return false;
 574       if ( NbSegm == 1 ) return true;
 575
 576       switch (_ivalue[ DISTR_TYPE_IND ])
 577       {
 578       case StdMeshers_NumberOfSegments::DT_Scale:
 579         {
 580           double scale = _value[ SCALE_FACTOR_IND ];
 581
 582           if (fabs(scale - 1.0) < Precision::Confusion()) {
 583             // special case to avoid division on zero
 584             for (int i = 1; i < NbSegm; i++) {
 585               double param = f + (l - f) * i / NbSegm;
 586               theParams.push_back( param );
 587             }
 588           } else {
 589             // general case of scale distribution
 590             if ( theReverse )
 591               scale = 1.0 / scale;
 592
 593             double alpha = pow(scale, 1.0 / (NbSegm - 1));
 594             double factor = (l - f) / (1.0 - pow(alpha, NbSegm));
 595
 596             for (int i = 1; i < NbSegm; i++) {
 597               double param = f + factor * (1.0 - pow(alpha, i));
 598               theParams.push_back( param );
 599             }
 600           }
 601           return true;
 602         }
 603         break;
 604       case StdMeshers_NumberOfSegments::DT_TabFunc:
 605         {
 606           FunctionTable func(_vvalue[ TAB_FUNC_IND ], _ivalue[ CONV_MODE_IND ]);
 607           return computeParamByFunc(theC3d, f, l, theLength, theReverse,
 608                                     _ivalue[ NB_SEGMENTS_IND ], func,
 609                                     theParams);
 610         }
 611         break;
 612       case StdMeshers_NumberOfSegments::DT_ExprFunc:
 613         {
 614           FunctionExpr func(_svalue[ EXPR_FUNC_IND ].c_str(), _ivalue[ CONV_MODE_IND ]);
 615           return computeParamByFunc(theC3d, f, l, theLength, theReverse,
 616                                     _ivalue[ NB_SEGMENTS_IND ], func,
 617                                     theParams);
 618         }
 619         break;
 620       case StdMeshers_NumberOfSegments::DT_Regular:
 621         eltSize = theLength / _ivalue[ NB_SEGMENTS_IND ];
 622         break;
 623       default:
 624         return false;
 625       }
 626     }
 627     GCPnts_UniformAbscissa Discret(theC3d, eltSize, f, l);
 628     if ( !Discret.IsDone() )
 629       return false;
 630
 631     int NbPoints = Discret.NbPoints();
 632     for ( int i = 2; i < NbPoints; i++ )
 633     {
 634       double param = Discret.Parameter(i);
 635       theParams.push_back( param );
 636     }
 637     compensateError( eltSize, eltSize, f, l, theLength, theC3d, theParams ); // for PAL9899
 638     return true;
 639   }
 640
 641   case BEG_END_LENGTH: {
 642
 643     // geometric progression: SUM(n) = ( a1 - an * q ) / ( 1 - q ) = theLength
 644
 645     double a1 = _value[ BEG_LENGTH_IND ];
 646     double an = _value[ END_LENGTH_IND ];
 647     double q  = ( theLength - a1 ) / ( theLength - an );
 648
 649     double U1 = theReverse ? l : f;
 650     double Un = theReverse ? f : l;
 651     double param = U1;
 652     double eltSize = theReverse ? -a1 : a1;
 653     while ( 1 ) {
 654       // computes a point on a curve <theC3d> at the distance <eltSize>
 655       // from the point of parameter <param>.
 656       GCPnts_AbscissaPoint Discret( theC3d, eltSize, param );
 657       if ( !Discret.IsDone() ) break;
 658       param = Discret.Parameter();
 659       if ( param > f && param < l )
 660         theParams.push_back( param );
 661       else
 662         break;
 663       eltSize *= q;
 664     }
 665     compensateError( a1, an, U1, Un, theLength, theC3d, theParams );
 666     return true;
 667   }
 668
 669   case ARITHMETIC_1D: {
 670
 671     // arithmetic progression: SUM(n) = ( an - a1 + q ) * ( a1 + an ) / ( 2 * q ) = theLength
 672
 673     double a1 = _value[ BEG_LENGTH_IND ];
 674     double an = _value[ END_LENGTH_IND ];
 675
 676     double  q = ( an - a1 ) / ( 2 *theLength/( a1 + an ) - 1 );
 677     int     n = int( 1 + ( an - a1 ) / q );
 678
 679     double U1 = theReverse ? l : f;
 680     double Un = theReverse ? f : l;
 681     double param = U1;
 682     double eltSize = a1;
 683     if ( theReverse ) {
 684       eltSize = -eltSize;
 685       q = -q;
 686     }
 687     while ( n-- > 0 && eltSize * ( Un - U1 ) > 0 ) {
 688       // computes a point on a curve <theC3d> at the distance <eltSize>
 689       // from the point of parameter <param>.
 690       GCPnts_AbscissaPoint Discret( theC3d, eltSize, param );
 691       if ( !Discret.IsDone() ) break;
 692       param = Discret.Parameter();
 693       if ( param > f && param < l )
 694         theParams.push_back( param );
 695       else
 696         break;
 697       eltSize += q;
 698     }
 699     compensateError( a1, an, U1, Un, theLength, theC3d, theParams );
 700
 701     return true;
 702   }
 703
 704   case DEFLECTION: {
 705
 706     GCPnts_UniformDeflection Discret(theC3d, _value[ DEFLECTION_IND ], f, l, true);
 707     if ( !Discret.IsDone() )
 708       return false;
 709
 710     int NbPoints = Discret.NbPoints();
 711     for ( int i = 2; i < NbPoints; i++ )
 712     {
 713       double param = Discret.Parameter(i);
 714       theParams.push_back( param );
 715     }
 716     return true;
 717
 718   }
 719
 720   default:;
 721   }
 722
 723   return false;
 724 }
 725
 726 //=============================================================================
 727 /*!
 728  *
 729  */
 730 //=============================================================================
 731
 732 bool StdMeshers_Regular_1D::Compute(SMESH_Mesh & aMesh, const TopoDS_Shape & aShape)
 733 {
 734   //MESSAGE("StdMeshers_Regular_1D::Compute");
 735
 736   if ( _hypType == NONE )
 737     return false;
 738
 739   SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
 740
 741   const TopoDS_Edge & EE = TopoDS::Edge(aShape);
 742   TopoDS_Edge E = TopoDS::Edge(EE.Oriented(TopAbs_FORWARD));
 743   int shapeID = meshDS->ShapeToIndex( E );
 744
 745   double f, l;
 746   Handle(Geom_Curve) Curve = BRep_Tool::Curve(E, f, l);
 747
 748   TopoDS_Vertex VFirst, VLast;
 749   TopExp::Vertices(E, VFirst, VLast);   // Vfirst corresponds to f and Vlast to l
 750
 751   ASSERT(!VFirst.IsNull());
 752   const SMDS_MeshNode * idFirst = SMESH_Algo::VertexNode( VFirst, meshDS );
 753   if (!idFirst) {
 754     MESSAGE (" NO NODE BUILT ON VERTEX ");
 755     return false;
 756   }
 757
 758   ASSERT(!VLast.IsNull());
 759   const SMDS_MeshNode * idLast = SMESH_Algo::VertexNode( VLast, meshDS );
 760   if (!idLast) {
 761     MESSAGE (" NO NODE BUILT ON VERTEX ");
 762     return false;
 763   }
 764
 765   if (!Curve.IsNull()) {
 766     list< double > params;
 767     bool reversed = false;
 768     if ( !_mainEdge.IsNull() )
 769       reversed = aMesh.IsReversedInChain( EE, _mainEdge );
 770     BRepAdaptor_Curve C3d( E );
 771     double length = EdgeLength( E );
 772     try {
 773 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 774       OCC_CATCH_SIGNALS;
 775 #endif
 776       if ( ! computeInternalParameters( C3d, length, f, l, params, reversed )) {
 777         return false;
 778       }
 779       redistributeNearVertices( aMesh, C3d, length, params, VFirst, VLast );
 780     }
 781     catch ( Standard_Failure ) {
 782       return false;
 783     }
 784
 785     // edge extrema (indexes : 1 & NbPoints) already in SMDS (TopoDS_Vertex)
 786     // only internal nodes receive an edge position with param on curve
 787
 788     const SMDS_MeshNode * idPrev = idFirst;
 789     double parPrev = f;
 790     double parLast = l;
 791
 792     for (list<double>::iterator itU = params.begin(); itU != params.end(); itU++) {
 793       double param = *itU;
 794       gp_Pnt P = Curve->Value(param);
 795
 796       //Add the Node in the DataStructure
 797       SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
 798       meshDS->SetNodeOnEdge(node, shapeID, param);
 799
 800       if(_quadraticMesh) {
 801         // create medium node
 802         double prm = ( parPrev + param )/2;
 803         gp_Pnt PM = Curve->Value(prm);
 804         SMDS_MeshNode * NM = meshDS->AddNode(PM.X(), PM.Y(), PM.Z());
 805         meshDS->SetNodeOnEdge(NM, shapeID, prm);
 806         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node, NM);
 807         meshDS->SetMeshElementOnShape(edge, shapeID);
 808       }
 809       else {
 810         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node);
 811         meshDS->SetMeshElementOnShape(edge, shapeID);
 812       }
 813
 814       idPrev = node;
 815       parPrev = param;
 816     }
 817     if(_quadraticMesh) {
 818       double prm = ( parPrev + parLast )/2;
 819       gp_Pnt PM = Curve->Value(prm);
 820       SMDS_MeshNode * NM = meshDS->AddNode(PM.X(), PM.Y(), PM.Z());
 821       meshDS->SetNodeOnEdge(NM, shapeID, prm);
 822       SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, idLast, NM);
 823       meshDS->SetMeshElementOnShape(edge, shapeID);
 824     }
 825     else {
 826       SMDS_MeshEdge* edge = meshDS->AddEdge(idPrev, idLast);
 827       meshDS->SetMeshElementOnShape(edge, shapeID);
 828     }
 829   }
 830   else {
 831     // Edge is a degenerated Edge : We put n = 5 points on the edge.
 832     const int NbPoints = 5;
 833     BRep_Tool::Range( E, f, l ); // PAL15185
 834     double du = (l - f) / (NbPoints - 1);
 835     //MESSAGE("************* Degenerated edge! *****************");
 836
 837     gp_Pnt P = BRep_Tool::Pnt(VFirst);
 838
 839     const SMDS_MeshNode * idPrev = idFirst;
 840     for (int i = 2; i < NbPoints; i++) {
 841       double param = f + (i - 1) * du;
 842       SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
 843       if(_quadraticMesh) {
 844         // create medium node
 845         double prm = param - du/2.;
 846         SMDS_MeshNode * NM = meshDS->AddNode(P.X(), P.Y(), P.Z());
 847         meshDS->SetNodeOnEdge(NM, shapeID, prm);
 848         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node, NM);
 849         meshDS->SetMeshElementOnShape(edge, shapeID);
 850       }
 851       else {
 852         SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, node);
 853         meshDS->SetMeshElementOnShape(edge, shapeID);
 854       }
 855       meshDS->SetNodeOnEdge(node, shapeID, param);
 856       idPrev = node;
 857     }
 858     if(_quadraticMesh) {
 859       // create medium node
 860       double prm = l - du/2.;
 861       SMDS_MeshNode * NM = meshDS->AddNode(P.X(), P.Y(), P.Z());
 862       meshDS->SetNodeOnEdge(NM, shapeID, prm);
 863       SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, idLast, NM);
 864       meshDS->SetMeshElementOnShape(edge, shapeID);
 865     }
 866     else {
 867       SMDS_MeshEdge * edge = meshDS->AddEdge(idPrev, idLast);
 868       meshDS->SetMeshElementOnShape(edge, shapeID);
 869     }
 870   }
 871   return true;
 872 }
 873
 874 //=============================================================================
 875 /*!
 876  *  See comments in SMESH_Algo.cxx
 877  */
 878 //=============================================================================
 879
 880 const list <const SMESHDS_Hypothesis *> &
 881 StdMeshers_Regular_1D::GetUsedHypothesis(SMESH_Mesh &         aMesh,
 882                                          const TopoDS_Shape & aShape,
 883                                          const bool           ignoreAuxiliary)
 884 {
 885   _usedHypList.clear();
 886   _mainEdge.Nullify();
 887
 888   SMESH_HypoFilter auxiliaryFilter, compatibleFilter;
 889   auxiliaryFilter.Init( SMESH_HypoFilter::IsAuxiliary() );
 890   const bool ignoreAux = true;
 891   InitCompatibleHypoFilter( compatibleFilter, ignoreAux );
 892
 893   // get non-auxiliary assigned to aShape
 894   int nbHyp = aMesh.GetHypotheses( aShape, compatibleFilter, _usedHypList, false );
 895
 896   if (nbHyp == 0)
 897   {
 898     // Check, if propagated from some other edge
 899     if (aShape.ShapeType() == TopAbs_EDGE &&
 900         aMesh.IsPropagatedHypothesis(aShape, _mainEdge))
 901     {
 902       // Propagation of 1D hypothesis from <aMainEdge> on this edge;
 903       // get non-auxiliary assigned to _mainEdge
 904       nbHyp = aMesh.GetHypotheses( _mainEdge, compatibleFilter, _usedHypList, true );
 905     }
 906   }
 907
 908   if (nbHyp == 0) // nothing propagated nor assigned to aShape
 909   {
 910     SMESH_Algo::GetUsedHypothesis( aMesh, aShape, ignoreAuxiliary );
 911     nbHyp = _usedHypList.size();
 912   }
 913   else
 914   {
 915     // get auxiliary hyps from aShape
 916     aMesh.GetHypotheses( aShape, auxiliaryFilter, _usedHypList, true );
 917   }
 918   if ( nbHyp > 1 && ignoreAuxiliary )
 919     _usedHypList.clear(); //only one compatible non-auxiliary hypothesis allowed
 920
 921   return _usedHypList;
 922 }