src/Extrema/Extrema_FuncExtCC.gxx

   1 // Copyright (c) 1995-1999 Matra Datavision
   2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   3 //
   4 // This file is part of Open CASCADE Technology software library.
   5 //
   6 // This library is free software; you can redistribute it and/or modify it under
   7 // the terms of the GNU Lesser General Public License version 2.1 as published
   8 // by the Free Software Foundation, with special exception defined in the file
   9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
  10 // distribution for complete text of the license and disclaimer of any warranty.
  11 //
  12 // Alternatively, this file may be used under the terms of Open CASCADE
  13 // commercial license or contractual agreement.
  14
  15 //Modified by MPS : 21-05-97   PRO 7598
  16 //                  Le nombre de solutions rendu est mySqDist.Length() et non
  17 //                  mySqDist.Length()/2
  18 //                  ajout des valeurs absolues dans le test d'orthogonalite de
  19 //                  GetStateNumber()
  20
  21 /*-----------------------------------------------------------------------------
  22 Fonctions permettant de rechercher une distance extremale entre 2 courbes
  23 C1 et C2 (en partant de points approches C1(u0) et C2(v0)).
  24 Cette classe herite de math_FunctionSetWithDerivatives et est utilisee par
  25 l'algorithme math_FunctionSetRoot.
  26 Si on note Du et Dv, les derivees en u et v, les 2 fonctions a annuler sont:
  27 { F1(u,v) = (C2(v)-C1(u)).Du(u)/||Du|| }
  28 { F2(u,v) = (C2(v)-C1(u)).Dv(v)||Dv|| }
  29 Si on note Duu et Dvv, les derivees secondes de C1 et C2, les derivees de F1
  30 et F2 sont egales a:
  31 { Duf1(u,v) = -||Du|| + C1C2.Duu/||Du||- F1(u,v)*Duu*Du/||Du||**2
  32 { Dvf1(u,v) = Dv.Du/||Du||
  33 { Duf2(u,v) = -Du.Dv/||Dv||
  34 { Dvf2(u,v) = ||Dv|| + C2C1.Dvv/||Dv||- F2(u,v)*Dv*Dvv/||Dv||**2
  35
  36 ----------------------------------------------------------------------------*/
  37
  38 #include <Precision.hxx>
  39
  40
  41 static const Standard_Real MinTol    = 1.e-20;
  42 static const Standard_Real TolFactor = 1.e-12;
  43 static const Standard_Real MinStep   = 1e-7;
  44
  45 static const Standard_Integer MaxOrder = 3;
  46
  47
  48
  49 //=============================================================================
  50 Standard_Real Extrema_FuncExtCC::SearchOfTolerance(const Standard_Address C)
  51 {
  52   const Standard_Integer NPoint = 10;
  53   Standard_Real aStartParam, anEndParam;
  54
  55   if(C==myC1)
  56   {
  57     aStartParam = myUinfium;
  58     anEndParam = myUsupremum;
  59   }
  60   else if(C==myC2)
  61   {
  62     aStartParam = myVinfium;
  63     anEndParam = myVsupremum;
  64   }
  65   else
  66   {
  67     //Warning: No curve for tolerance computing!
  68     return MinTol;
  69   }
  70
  71   const Standard_Real aStep = (anEndParam - aStartParam)/(Standard_Real)NPoint;
  72
  73   Standard_Integer aNum = 0;
  74   Standard_Real aMax = -Precision::Infinite();  //Maximum value of 1st derivative
  75   //(it is computed with using NPoint point)
  76
  77   do
  78   {
  79     Standard_Real u = aStartParam + aNum*aStep;    //parameter for every point
  80     if(u > anEndParam)
  81       u = anEndParam;
  82
  83     Pnt Ptemp;  //empty point (is not used below)
  84     Vec VDer;   // 1st derivative vector
  85     Tool1::D1(*((Curve1*)C), u, Ptemp, VDer);
  86     Standard_Real vm = VDer.Magnitude();
  87     if(vm > aMax)
  88       aMax = vm;
  89   }
  90   while(++aNum < NPoint+1);
  91
  92   return Max(aMax*TolFactor,MinTol);
  93 }
  94
  95 //=============================================================================
  96 Extrema_FuncExtCC::Extrema_FuncExtCC(const Standard_Real thetol) : myC1 (0), myC2 (0), myTol (thetol)
  97 {
  98   math_Vector V1(1,2), V2(1,2);
  99   V1(1) = 0.0;
 100   V2(1) = 0.0;
 101   V1(2) = 0.0;
 102   V2(2) = 0.0;
 103   SubIntervalInitialize(V1, V2);
 104   myMaxDerivOrderC1 = 0;
 105   myTolC1=MinTol;
 106   myMaxDerivOrderC2 = 0;
 107   myTolC2=MinTol;
 108 }
 109
 110 //=============================================================================
 111 Extrema_FuncExtCC::Extrema_FuncExtCC (const Curve1& C1,
 112                                       const Curve2& C2,
 113                                       const Standard_Real thetol) :
 114 myC1 ((Standard_Address)&C1), myC2 ((Standard_Address)&C2),
 115 myTol (thetol)
 116 {
 117   math_Vector V1(1,2), V2(1,2);
 118   V1(1) = Tool1::FirstParameter(*((Curve1*)myC1));
 119   V2(1) = Tool1::LastParameter(*((Curve1*)myC1));
 120   V1(2) = Tool2::FirstParameter(*((Curve2*)myC2));
 121   V2(2) = Tool2::LastParameter(*((Curve2*)myC2));
 122   SubIntervalInitialize(V1, V2);
 123
 124   switch(Tool1::GetType(*((Curve1*)myC1)))
 125   {
 126   case GeomAbs_BezierCurve:
 127   case GeomAbs_BSplineCurve:
 128   case GeomAbs_OffsetCurve:
 129   case GeomAbs_OtherCurve:
 130     myMaxDerivOrderC1 = MaxOrder;
 131     myTolC1 = SearchOfTolerance((Standard_Address)&C1);
 132     break;
 133   default:
 134     myMaxDerivOrderC1 = 0;
 135     myTolC1=MinTol;
 136     break;
 137   }
 138
 139   switch(Tool2::GetType(*((Curve2*)myC2)))
 140   {
 141   case GeomAbs_BezierCurve:
 142   case GeomAbs_BSplineCurve:
 143   case GeomAbs_OffsetCurve:
 144   case GeomAbs_OtherCurve:
 145     myMaxDerivOrderC2 = MaxOrder;
 146     myTolC2 = SearchOfTolerance((Standard_Address)&C2);
 147     break;
 148   default:
 149     myMaxDerivOrderC2 = 0;
 150     myTolC2=MinTol;
 151     break;
 152   }
 153 }
 154
 155 //=============================================================================
 156 void Extrema_FuncExtCC::SetCurve (const Standard_Integer theRank, const Curve1& C)
 157 {
 158   Standard_OutOfRange_Raise_if (theRank < 1 || theRank > 2, "Extrema_FuncExtCC::SetCurve()")
 159
 160     if (theRank == 1)
 161     {
 162       myC1 = (Standard_Address)&C;
 163       switch(/*Tool1::GetType(*((Curve1*)myC1))*/ C.GetType())
 164       {
 165       case GeomAbs_BezierCurve:
 166       case GeomAbs_BSplineCurve:
 167       case GeomAbs_OffsetCurve:
 168       case GeomAbs_OtherCurve:
 169         myMaxDerivOrderC1 = MaxOrder;
 170         myTolC1 = SearchOfTolerance((Standard_Address)&C);
 171         break;
 172       default:
 173         myMaxDerivOrderC1 = 0;
 174         myTolC1=MinTol;
 175         break;
 176       }
 177     }
 178     else if (theRank == 2)
 179     {
 180       myC2 = (Standard_Address)&C;
 181       switch(/*Tool2::GetType(*((Curve2*)myC2))*/C.GetType())
 182       {
 183       case GeomAbs_BezierCurve:
 184       case GeomAbs_BSplineCurve:
 185       case GeomAbs_OffsetCurve:
 186       case GeomAbs_OtherCurve:
 187         myMaxDerivOrderC2 = MaxOrder;
 188         myTolC2 = SearchOfTolerance((Standard_Address)&C);
 189         break;
 190       default:
 191         myMaxDerivOrderC2 = 0;
 192         myTolC2=MinTol;
 193         break;
 194       }
 195     }
 196 }
 197
 198 //=============================================================================
 199 Standard_Boolean Extrema_FuncExtCC::Value (const math_Vector& UV, math_Vector& F)
 200 {
 201   myU = UV(1);
 202   myV = UV(2);
 203   Tool1::D1(*((Curve1*)myC1), myU,myP1,myDu);
 204   Tool2::D1(*((Curve2*)myC2), myV,myP2,myDv);
 205
 206   Vec P1P2 (myP1,myP2);
 207
 208   Standard_Real Ndu = myDu.Magnitude();
 209
 210
 211   if(myMaxDerivOrderC1 != 0)
 212   {
 213     if (Ndu <= myTolC1)
 214     {
 215       //Derivative is approximated by Taylor-series
 216       const Standard_Real DivisionFactor = 1.e-3;
 217       Standard_Real du;
 218       if((myUsupremum >= RealLast()) || (myUinfium <= RealFirst()))
 219         du = 0.0;
 220       else
 221         du = myUsupremum-myUinfium;
 222
 223       const Standard_Real aDelta = Max(du*DivisionFactor,MinStep);
 224
 225       Standard_Integer n = 1; //Derivative order
 226       Vec V;
 227       Standard_Boolean IsDeriveFound;
 228
 229       do
 230       {
 231         V = Tool1::DN(*((Curve1*)myC1),myU,++n);
 232         Ndu = V.Magnitude();
 233         IsDeriveFound = (Ndu > myTolC1);
 234       }
 235       while(!IsDeriveFound && n < myMaxDerivOrderC1);
 236
 237       if(IsDeriveFound)
 238       {
 239         Standard_Real u;
 240
 241         if(myU-myUinfium < aDelta)
 242           u = myU+aDelta;
 243         else
 244           u = myU-aDelta;
 245
 246         Pnt P1, P2;
 247         Tool1::D0(*((Curve1*)myC1),Min(myU, u),P1);
 248         Tool1::D0(*((Curve1*)myC1),Max(myU, u),P2);
 249
 250         Vec V1(P1,P2);
 251         Standard_Real aDirFactor = V.Dot(V1);
 252
 253         if(aDirFactor < 0.0)
 254           myDu = -V;
 255         else
 256           myDu = V;
 257       }//if(IsDeriveFound)
 258       else
 259       {
 260         //Derivative is approximated by three points
 261
 262         Pnt Ptemp; //(0,0,0)-coordinate
 263         Pnt P1, P2, P3;
 264         Standard_Boolean IsParameterGrown;
 265
 266         if(myU-myUinfium < 2*aDelta)
 267         {
 268           Tool1::D0(*((Curve1*)myC1),myU,P1);
 269           Tool1::D0(*((Curve1*)myC1),myU+aDelta,P2);
 270           Tool1::D0(*((Curve1*)myC1),myU+2*aDelta,P3);
 271           IsParameterGrown = Standard_True;
 272         }
 273         else
 274         {
 275           Tool1::D0(*((Curve1*)myC1),myU-2*aDelta,P1);
 276           Tool1::D0(*((Curve1*)myC1),myU-aDelta,P2);
 277           Tool1::D0(*((Curve1*)myC1),myU,P3);
 278           IsParameterGrown = Standard_False;
 279         }
 280
 281         Vec V1(Ptemp,P1), V2(Ptemp,P2), V3(Ptemp,P3);
 282
 283         if(IsParameterGrown)
 284           myDu=-3*V1+4*V2-V3;
 285         else
 286           myDu=V1-4*V2+3*V3;
 287       }//else of if(IsDeriveFound)
 288       Ndu = myDu.Magnitude();
 289     }//if (Ndu <= myTolC1) condition
 290   }//if(myMaxDerivOrder != 0)
 291
 292   if (Ndu <= MinTol)
 293   {
 294     //Warning: 1st derivative of C1 is equal to zero!
 295     return Standard_False;
 296   }
 297
 298   Standard_Real Ndv = myDv.Magnitude();
 299
 300   if(myMaxDerivOrderC2 != 0)
 301   {
 302     if (Ndv <= myTolC2)
 303     {
 304       const Standard_Real DivisionFactor = 1.e-3;
 305       Standard_Real dv;
 306       if((myVsupremum >= RealLast()) || (myVinfium <= RealFirst()))
 307         dv = 0.0;
 308       else
 309         dv = myVsupremum-myVinfium;
 310
 311       const Standard_Real aDelta = Max(dv*DivisionFactor,MinStep);
 312
 313       //Derivative is approximated by Taylor-series
 314
 315       Standard_Integer n = 1; //Derivative order
 316       Vec V;
 317       Standard_Boolean IsDeriveFound;
 318
 319       do
 320       {
 321         V = Tool2::DN(*((Curve2*)myC2),myV,++n);
 322         Ndv = V.Magnitude();
 323         IsDeriveFound = (Ndv > myTolC2);
 324       }
 325       while(!IsDeriveFound && n < myMaxDerivOrderC2);
 326
 327       if(IsDeriveFound)
 328       {
 329         Standard_Real v;
 330
 331         if(myV-myVinfium < aDelta)
 332           v = myV+aDelta;
 333         else
 334           v = myV-aDelta;
 335
 336         Pnt P1, P2;
 337         Tool2::D0(*((Curve2*)myC2),Min(myV, v),P1);
 338         Tool2::D0(*((Curve2*)myC2),Max(myV, v),P2);
 339
 340         Vec V1(P1,P2);
 341         Standard_Real aDirFactor = V.Dot(V1);
 342
 343         if(aDirFactor < 0.0)
 344           myDv = -V;
 345         else
 346           myDv = V;
 347       }//if(IsDeriveFound)
 348       else
 349       {
 350         //Derivative is approximated by three points
 351
 352         Pnt Ptemp; //(0,0,0)-coordinate
 353         Pnt P1, P2, P3;
 354         Standard_Boolean IsParameterGrown;
 355
 356         if(myV-myVinfium < 2*aDelta)
 357         {
 358           Tool2::D0(*((Curve2*)myC2),myV,P1);
 359           Tool2::D0(*((Curve2*)myC2),myV+aDelta,P2);
 360           Tool2::D0(*((Curve2*)myC2),myV+2*aDelta,P3);
 361           IsParameterGrown = Standard_True;
 362         }
 363         else
 364         {
 365           Tool2::D0(*((Curve2*)myC2),myV-2*aDelta,P1);
 366           Tool2::D0(*((Curve2*)myC2),myV-aDelta,P2);
 367           Tool2::D0(*((Curve2*)myC2),myV,P3);
 368           IsParameterGrown = Standard_False;
 369         }
 370
 371         Vec V1(Ptemp,P1), V2(Ptemp,P2), V3(Ptemp,P3);
 372
 373         if(IsParameterGrown)
 374           myDv=-3*V1+4*V2-V3;
 375         else
 376           myDv=V1-4*V2+3*V3;
 377       }//else of if(IsDeriveFound)
 378
 379       Ndv = myDv.Magnitude();
 380     }//if (Ndv <= myTolC2)
 381   }//if(myMaxDerivOrder != 0)
 382
 383   if (Ndv <= MinTol)
 384   {
 385     //1st derivative of C2 is equal to zero!
 386     return Standard_False;
 387   }
 388
 389   F(1) = P1P2.Dot(myDu)/Ndu;
 390   F(2) = P1P2.Dot(myDv)/Ndv;
 391   return Standard_True;
 392 }
 393 //=============================================================================
 394
 395 Standard_Boolean Extrema_FuncExtCC::Derivatives (const math_Vector& UV,
 396                                                  math_Matrix& Df)
 397 {
 398   math_Vector F(1,2);
 399   return Values(UV,F,Df);
 400 }
 401 //=============================================================================
 402
 403 Standard_Boolean Extrema_FuncExtCC::Values (const math_Vector& UV,
 404                                             math_Vector& F,
 405                                             math_Matrix& Df)
 406 {
 407   myU = UV(1);
 408   myV = UV(2);
 409
 410   if(Value(UV, F) == Standard_False) //Computes F, myDu, myDv
 411   {
 412     //Warning: No function value found!
 413     return Standard_False;
 414   }//if(Value(UV, F) == Standard_False)
 415
 416   Vec Du, Dv, Duu, Dvv;
 417   Tool1::D2(*((Curve1*)myC1), myU,myP1,Du,Duu);
 418   Tool2::D2(*((Curve2*)myC2), myV,myP2,Dv,Dvv);
 419
 420   //Calling of "Value(...)" function change class member values.
 421   //After running it is necessary to return to previous values.
 422   const Standard_Real myU_old  = myU,    myV_old = myV;
 423   const Pnt           myP1_old = myP1,  myP2_old = myP2;
 424   const Vec           myDu_old = myDu,  myDv_old = myDv;
 425
 426
 427   //Attention:  aDelta value must be greater than same value for "Value(...)"
 428   //            function to avoid of points' collisions.
 429
 430   const Standard_Real DivisionFactor = 0.01;
 431
 432   Standard_Real du;
 433   if((myUsupremum >= RealLast()) || (myUinfium <= RealFirst()))
 434     du = 0.0;
 435   else
 436     du = myUsupremum-myUinfium;
 437
 438   const Standard_Real aDeltaU = Max(du*DivisionFactor,MinStep);
 439
 440   Standard_Real dv;
 441   if((myVsupremum >= RealLast()) || (myVinfium <= RealFirst()))
 442     dv = 0.0;
 443   else
 444     dv = myVsupremum-myVinfium;
 445
 446   const Standard_Real aDeltaV = Max(dv*DivisionFactor,MinStep);
 447
 448   Vec P1P2 (myP1,myP2);
 449
 450   if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1))
 451   {
 452     //Derivative is approximated by three points
 453
 454     math_Vector FF1(1,2), FF2(1,2), FF3(1,2);
 455     Standard_Real F1, F2, F3;
 456
 457     /////////////////////////// Search of DF1_u derivative (begin) ///////////////////
 458     if(myU-myUinfium < 2*aDeltaU)
 459     {
 460       F1=F(1);
 461       math_Vector UV2(1,2), UV3(1,2);
 462       UV2(1)=myU+aDeltaU;
 463       UV2(2)=myV;
 464       UV3(1)=myU+2*aDeltaU;
 465       UV3(2)=myV;
 466       if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
 467       {
 468         //There are many points close to singularity points and which have zero-derivative.
 469         //Try to decrease aDelta variable's value.
 470         return Standard_False;
 471       }
 472
 473       F2 = FF2(1);
 474       F3 = FF3(1);
 475
 476       Df(1,1) = (-3*F1+4*F2-F3)/(2.0*aDeltaU);
 477     }//if(myU-myUinfium < 2*aDeltaU)
 478     else
 479     {
 480       F3 = F(1);
 481       math_Vector UV2(1,2), UV1(1,2);
 482       UV2(1)=myU-aDeltaU;
 483       UV2(2)=myV;
 484       UV1(1)=myU-2*aDeltaU;
 485       UV1(2)=myV;
 486
 487       if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
 488       {
 489         //There are many points close to singularity points and which have zero-derivative.
 490         //Try to decrease aDelta variable's value.
 491         return Standard_False;
 492       }
 493
 494       F1 = FF1(1);
 495       F2 = FF2(1);
 496
 497       Df(1,1) = (F1-4*F2+3*F3)/(2.0*aDeltaU);
 498     }//else of if(myU-myUinfium < 2*aDeltaU) condition
 499     /////////////////////////// Search of DF1_u derivative (end) ///////////////////
 500
 501     //Return to previous values
 502     myU = myU_old;
 503     myV = myV_old;
 504
 505     /////////////////////////// Search of DF1_v derivative (begin) ///////////////////
 506     if(myV-myVinfium < 2*aDeltaV)
 507     {
 508       F1=F(1);
 509       math_Vector UV2(1,2), UV3(1,2);
 510       UV2(1)=myU;
 511       UV2(2)=myV+aDeltaV;
 512       UV3(1)=myU;
 513       UV3(2)=myV+2*aDeltaV;
 514
 515       if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
 516       {
 517         //There are many points close to singularity points and which have zero-derivative.
 518         //Try to decrease aDelta variable's value.
 519         return Standard_False;
 520       }
 521       F2 = FF2(1);
 522       F3 = FF3(1);
 523
 524       Df(1,2) = (-3*F1+4*F2-F3)/(2.0*aDeltaV);
 525     }//if(myV-myVinfium < 2*aDeltaV)
 526     else
 527     {
 528       F3 = F(1);
 529       math_Vector UV2(1,2), UV1(1,2);
 530       UV2(1)=myU;
 531       UV2(2)=myV-aDeltaV;
 532       UV1(1)=myU;
 533       UV1(2)=myV-2*aDeltaV;
 534       if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
 535       {
 536         //There are many points close to singularity points and which have zero-derivative.
 537         //Try to decrease aDelta variable's value.
 538         return Standard_False;
 539       }
 540
 541       F1 = FF1(1);
 542       F2 = FF2(1);
 543
 544       Df(1,2) = (F1-4*F2+3*F3)/(2.0*aDeltaV);
 545     }//else of if(myV-myVinfium < 2*aDeltaV)
 546     /////////////////////////// Search of DF1_v derivative (end) ///////////////////
 547
 548     //Return to previous values
 549     myU = myU_old;
 550     myV = myV_old;
 551     myP1 = myP1_old,  myP2 = myP2_old;
 552     myDu = myDu_old,  myDv = myDv_old;
 553   }//if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1))
 554   else
 555   {
 556     const Standard_Real Ndu = myDu.Magnitude();
 557     Df(1,1) = - Ndu + (P1P2.Dot(Duu)/Ndu) - F(1)*(myDu.Dot(Duu)/(Ndu*Ndu));
 558     Df(1,2) = myDv.Dot(myDu)/Ndu;
 559   }//else of if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1))
 560
 561   if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2))
 562   {
 563     //Derivative is approximated by three points
 564
 565     math_Vector FF1(1,2), FF2(1,2), FF3(1,2);
 566     Standard_Real F1, F2, F3;
 567
 568     /////////////////////////// Search of DF2_v derivative (begin) ///////////////////
 569     if(myV-myVinfium < 2*aDeltaV)
 570     {
 571       F1=F(2);
 572       math_Vector UV2(1,2), UV3(1,2);
 573       UV2(1)=myU;
 574       UV2(2)=myV+aDeltaV;
 575       UV3(1)=myU;
 576       UV3(2)=myV+2*aDeltaV;
 577
 578       if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
 579       {
 580         //There are many points close to singularity points and which have zero-derivative.
 581         //Try to decrease aDelta variable's value.
 582         return Standard_False;
 583       }
 584
 585       F2 = FF2(2);
 586       F3 = FF3(2);
 587
 588       Df(2,2) = (-3*F1+4*F2-F3)/(2.0*aDeltaV);
 589
 590     }//if(myV-myVinfium < 2*aDeltaV)
 591     else
 592     {
 593       F3 = F(2);
 594       math_Vector UV2(1,2), UV1(1,2);
 595       UV2(1)=myU;
 596       UV2(2)=myV-aDeltaV;
 597       UV1(1)=myU;
 598       UV1(2)=myV-2*aDeltaV;
 599
 600       if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
 601       {
 602         //There are many points close to singularity points and which have zero-derivative.
 603         //Try to decrease aDelta variable's value.
 604         return Standard_False;
 605       }
 606
 607       F1 = FF1(2);
 608       F2 = FF2(2);
 609
 610       Df(2,2) = (F1-4*F2+3*F3)/(2.0*aDeltaV);
 611     }//else of if(myV-myVinfium < 2*aDeltaV)
 612     /////////////////////////// Search of DF2_v derivative (end) ///////////////////
 613
 614     //Return to previous values
 615     myU = myU_old;
 616     myV = myV_old;
 617
 618     /////////////////////////// Search of DF2_u derivative (begin) ///////////////////
 619     if(myU-myUinfium < 2*aDeltaU)
 620     {
 621       F1=F(2);
 622       math_Vector UV2(1,2), UV3(1,2);
 623       UV2(1)=myU+aDeltaU;
 624       UV2(2)=myV;
 625       UV3(1)=myU+2*aDeltaU;
 626       UV3(2)=myV;
 627       if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
 628       {
 629         //There are many points close to singularity points and which have zero-derivative.
 630         //Try to decrease aDelta variable's value.
 631         return Standard_False;
 632       }
 633
 634       F2 = FF2(2);
 635       F3 = FF3(2);
 636
 637       Df(2,1) = (-3*F1+4*F2-F3)/(2.0*aDeltaU);
 638     }//if(myU-myUinfium < 2*aDelta)
 639     else
 640     {
 641       F3 = F(2);
 642       math_Vector UV2(1,2), UV1(1,2);
 643       UV2(1)=myU-aDeltaU;
 644       UV2(2)=myV;
 645       UV1(1)=myU-2*aDeltaU;
 646       UV1(2)=myV;
 647
 648       if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
 649       {
 650         //There are many points close to singularity points
 651         //and which have zero-derivative.
 652         //Try to decrease aDelta variable's value.
 653         return Standard_False;
 654       }
 655
 656       F1 = FF1(2);
 657       F2 = FF2(2);
 658
 659       Df(2,1) = (F1-4*F2+3*F3)/(2.0*aDeltaU);
 660     }//else of if(myU-myUinfium < 2*aDeltaU)
 661     /////////////////////////// Search of DF2_u derivative (end) ///////////////////
 662
 663     //Return to previous values
 664     myU = myU_old;
 665     myV = myV_old;
 666     myP1 = myP1_old;
 667     myP2 = myP2_old;
 668     myDu = myDu_old;
 669     myDv = myDv_old;
 670   }//if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2))
 671   else
 672   {
 673     Standard_Real Ndv = myDv.Magnitude();
 674     Df(2,2) = Ndv + (P1P2.Dot(Dvv)/Ndv) - F(2)*(myDv.Dot(Dvv)/(Ndv*Ndv));
 675     Df(2,1) = -myDu.Dot(myDv)/Ndv;
 676   }//else of if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2))
 677
 678   return Standard_True;
 679
 680 }//end of function
 681 //=============================================================================
 682
 683 Standard_Integer Extrema_FuncExtCC::GetStateNumber ()
 684 {
 685   Vec Du (myDu), Dv (myDv);
 686   Vec P1P2 (myP1, myP2);
 687
 688   Standard_Real mod = Du.Magnitude();
 689   if(mod > myTolC1) {
 690     Du /= mod;
 691   }
 692   mod = Dv.Magnitude();
 693   if(mod > myTolC2) {
 694     Dv /= mod;
 695   }
 696
 697   if (Abs(P1P2.Dot(Du)) <= myTol && Abs(P1P2.Dot(Dv)) <= myTol) {
 698     mySqDist.Append(myP1.SquareDistance(myP2));
 699     myPoints.Append(POnC(myU, myP1));
 700     myPoints.Append(POnC(myV, myP2));
 701   }
 702   return 0;
 703 }
 704 //=============================================================================
 705
 706 void Extrema_FuncExtCC::Points (const Standard_Integer N,
 707                                 POnC& P1,
 708                                 POnC& P2) const
 709 {
 710   P1 = myPoints.Value(2*N-1);
 711   P2 = myPoints.Value(2*N);
 712 }
 713 //=============================================================================
 714
 715 void Extrema_FuncExtCC::SubIntervalInitialize(const math_Vector& theInfBound,
 716                                               const math_Vector& theSupBound)
 717 {
 718   myUinfium = theInfBound(1);
 719   myUsupremum = theSupBound(1);
 720   myVinfium = theInfBound(2);
 721   myVsupremum = theSupBound(2);
 722 }
 723 //=============================================================================