Warning, /include/opencascade/Extrema_FuncExtCC.gxx is written in an unsupported language. File is not indexed.
0001 // Copyright (c) 1995-1999 Matra Datavision
0002 // Copyright (c) 1999-2014 OPEN CASCADE SAS
0003 //
0004 // This file is part of Open CASCADE Technology software library.
0005 //
0006 // This library is free software; you can redistribute it and/or modify it under
0007 // the terms of the GNU Lesser General Public License version 2.1 as published
0008 // by the Free Software Foundation, with special exception defined in the file
0009 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
0010 // distribution for complete text of the license and disclaimer of any warranty.
0011 //
0012 // Alternatively, this file may be used under the terms of Open CASCADE
0013 // commercial license or contractual agreement.
0014
0015 //Modified by MPS : 21-05-97 PRO 7598
0016 // Le nombre de solutions rendu est mySqDist.Length() et non
0017 // mySqDist.Length()/2
0018 // ajout des valeurs absolues dans le test d'orthogonalite de
0019 // GetStateNumber()
0020
0021 /*-----------------------------------------------------------------------------
0022 Fonctions permettant de rechercher une distance extremale entre 2 courbes
0023 C1 et C2 (en partant de points approches C1(u0) et C2(v0)).
0024 Cette classe herite de math_FunctionSetWithDerivatives et est utilisee par
0025 l'algorithme math_FunctionSetRoot.
0026 Si on note Du et Dv, les derivees en u et v, les 2 fonctions a annuler sont:
0027 { F1(u,v) = (C2(v)-C1(u)).Du(u)/||Du|| }
0028 { F2(u,v) = (C2(v)-C1(u)).Dv(v)||Dv|| }
0029 Si on note Duu et Dvv, les derivees secondes de C1 et C2, les derivees de F1
0030 et F2 sont egales a:
0031 { Duf1(u,v) = -||Du|| + C1C2.Duu/||Du||- F1(u,v)*Duu*Du/||Du||**2
0032 { Dvf1(u,v) = Dv.Du/||Du||
0033 { Duf2(u,v) = -Du.Dv/||Dv||
0034 { Dvf2(u,v) = ||Dv|| + C2C1.Dvv/||Dv||- F2(u,v)*Dv*Dvv/||Dv||**2
0035
0036 ----------------------------------------------------------------------------*/
0037
0038 #include <Precision.hxx>
0039
0040
0041 static const Standard_Real MinTol = 1.e-20;
0042 static const Standard_Real TolFactor = 1.e-12;
0043 static const Standard_Real MinStep = 1e-7;
0044
0045 static const Standard_Integer MaxOrder = 3;
0046
0047
0048
0049 //=============================================================================
0050 Standard_Real Extrema_FuncExtCC::SearchOfTolerance(const Standard_Address C)
0051 {
0052 const Standard_Integer NPoint = 10;
0053 Standard_Real aStartParam, anEndParam;
0054
0055 if(C==myC1)
0056 {
0057 aStartParam = myUinfium;
0058 anEndParam = myUsupremum;
0059 }
0060 else if(C==myC2)
0061 {
0062 aStartParam = myVinfium;
0063 anEndParam = myVsupremum;
0064 }
0065 else
0066 {
0067 //Warning: No curve for tolerance computing!
0068 return MinTol;
0069 }
0070
0071 const Standard_Real aStep = (anEndParam - aStartParam)/(Standard_Real)NPoint;
0072
0073 Standard_Integer aNum = 0;
0074 Standard_Real aMax = -Precision::Infinite(); //Maximum value of 1st derivative
0075 //(it is computed with using NPoint point)
0076
0077 do
0078 {
0079 Standard_Real u = aStartParam + aNum*aStep; //parameter for every point
0080 if(u > anEndParam)
0081 u = anEndParam;
0082
0083 Pnt Ptemp; //empty point (is not used below)
0084 Vec VDer; // 1st derivative vector
0085 Tool1::D1(*((Curve1*)C), u, Ptemp, VDer);
0086 Standard_Real vm = VDer.Magnitude();
0087 if(vm > aMax)
0088 aMax = vm;
0089 }
0090 while(++aNum < NPoint+1);
0091
0092 return Max(aMax*TolFactor,MinTol);
0093 }
0094
0095 //=============================================================================
0096 Extrema_FuncExtCC::Extrema_FuncExtCC(const Standard_Real thetol) : myC1 (0), myC2 (0), myTol (thetol)
0097 {
0098 math_Vector V1(1,2), V2(1,2);
0099 V1(1) = 0.0;
0100 V2(1) = 0.0;
0101 V1(2) = 0.0;
0102 V2(2) = 0.0;
0103 SubIntervalInitialize(V1, V2);
0104 myMaxDerivOrderC1 = 0;
0105 myTolC1=MinTol;
0106 myMaxDerivOrderC2 = 0;
0107 myTolC2=MinTol;
0108 }
0109
0110 //=============================================================================
0111 Extrema_FuncExtCC::Extrema_FuncExtCC (const Curve1& C1,
0112 const Curve2& C2,
0113 const Standard_Real thetol) :
0114 myC1 ((Standard_Address)&C1), myC2 ((Standard_Address)&C2),
0115 myTol (thetol)
0116 {
0117 math_Vector V1(1,2), V2(1,2);
0118 V1(1) = Tool1::FirstParameter(*((Curve1*)myC1));
0119 V2(1) = Tool1::LastParameter(*((Curve1*)myC1));
0120 V1(2) = Tool2::FirstParameter(*((Curve2*)myC2));
0121 V2(2) = Tool2::LastParameter(*((Curve2*)myC2));
0122 SubIntervalInitialize(V1, V2);
0123
0124 switch(Tool1::GetType(*((Curve1*)myC1)))
0125 {
0126 case GeomAbs_BezierCurve:
0127 case GeomAbs_BSplineCurve:
0128 case GeomAbs_OffsetCurve:
0129 case GeomAbs_OtherCurve:
0130 myMaxDerivOrderC1 = MaxOrder;
0131 myTolC1 = SearchOfTolerance((Standard_Address)&C1);
0132 break;
0133 default:
0134 myMaxDerivOrderC1 = 0;
0135 myTolC1=MinTol;
0136 break;
0137 }
0138
0139 switch(Tool2::GetType(*((Curve2*)myC2)))
0140 {
0141 case GeomAbs_BezierCurve:
0142 case GeomAbs_BSplineCurve:
0143 case GeomAbs_OffsetCurve:
0144 case GeomAbs_OtherCurve:
0145 myMaxDerivOrderC2 = MaxOrder;
0146 myTolC2 = SearchOfTolerance((Standard_Address)&C2);
0147 break;
0148 default:
0149 myMaxDerivOrderC2 = 0;
0150 myTolC2=MinTol;
0151 break;
0152 }
0153 }
0154
0155 //=============================================================================
0156 void Extrema_FuncExtCC::SetCurve (const Standard_Integer theRank, const Curve1& C)
0157 {
0158 Standard_OutOfRange_Raise_if (theRank < 1 || theRank > 2, "Extrema_FuncExtCC::SetCurve()")
0159
0160 if (theRank == 1)
0161 {
0162 myC1 = (Standard_Address)&C;
0163 switch(/*Tool1::GetType(*((Curve1*)myC1))*/ C.GetType())
0164 {
0165 case GeomAbs_BezierCurve:
0166 case GeomAbs_BSplineCurve:
0167 case GeomAbs_OffsetCurve:
0168 case GeomAbs_OtherCurve:
0169 myMaxDerivOrderC1 = MaxOrder;
0170 myTolC1 = SearchOfTolerance((Standard_Address)&C);
0171 break;
0172 default:
0173 myMaxDerivOrderC1 = 0;
0174 myTolC1=MinTol;
0175 break;
0176 }
0177 }
0178 else if (theRank == 2)
0179 {
0180 myC2 = (Standard_Address)&C;
0181 switch(/*Tool2::GetType(*((Curve2*)myC2))*/C.GetType())
0182 {
0183 case GeomAbs_BezierCurve:
0184 case GeomAbs_BSplineCurve:
0185 case GeomAbs_OffsetCurve:
0186 case GeomAbs_OtherCurve:
0187 myMaxDerivOrderC2 = MaxOrder;
0188 myTolC2 = SearchOfTolerance((Standard_Address)&C);
0189 break;
0190 default:
0191 myMaxDerivOrderC2 = 0;
0192 myTolC2=MinTol;
0193 break;
0194 }
0195 }
0196 }
0197
0198 //=============================================================================
0199 Standard_Boolean Extrema_FuncExtCC::Value (const math_Vector& UV, math_Vector& F)
0200 {
0201 myU = UV(1);
0202 myV = UV(2);
0203 Tool1::D1(*((Curve1*)myC1), myU,myP1,myDu);
0204 Tool2::D1(*((Curve2*)myC2), myV,myP2,myDv);
0205
0206 Vec P1P2 (myP1,myP2);
0207
0208 Standard_Real Ndu = myDu.Magnitude();
0209
0210
0211 if(myMaxDerivOrderC1 != 0)
0212 {
0213 if (Ndu <= myTolC1)
0214 {
0215 //Derivative is approximated by Taylor-series
0216 const Standard_Real DivisionFactor = 1.e-3;
0217 Standard_Real du;
0218 if((myUsupremum >= RealLast()) || (myUinfium <= RealFirst()))
0219 du = 0.0;
0220 else
0221 du = myUsupremum-myUinfium;
0222
0223 const Standard_Real aDelta = Max(du*DivisionFactor,MinStep);
0224
0225 Standard_Integer n = 1; //Derivative order
0226 Vec V;
0227 Standard_Boolean IsDeriveFound;
0228
0229 do
0230 {
0231 V = Tool1::DN(*((Curve1*)myC1),myU,++n);
0232 Ndu = V.Magnitude();
0233 IsDeriveFound = (Ndu > myTolC1);
0234 }
0235 while(!IsDeriveFound && n < myMaxDerivOrderC1);
0236
0237 if(IsDeriveFound)
0238 {
0239 Standard_Real u;
0240
0241 if(myU-myUinfium < aDelta)
0242 u = myU+aDelta;
0243 else
0244 u = myU-aDelta;
0245
0246 Pnt P1, P2;
0247 Tool1::D0(*((Curve1*)myC1),Min(myU, u),P1);
0248 Tool1::D0(*((Curve1*)myC1),Max(myU, u),P2);
0249
0250 Vec V1(P1,P2);
0251 Standard_Real aDirFactor = V.Dot(V1);
0252
0253 if(aDirFactor < 0.0)
0254 myDu = -V;
0255 else
0256 myDu = V;
0257 }//if(IsDeriveFound)
0258 else
0259 {
0260 //Derivative is approximated by three points
0261
0262 Pnt Ptemp; //(0,0,0)-coordinate
0263 Pnt P1, P2, P3;
0264 Standard_Boolean IsParameterGrown;
0265
0266 if(myU-myUinfium < 2*aDelta)
0267 {
0268 Tool1::D0(*((Curve1*)myC1),myU,P1);
0269 Tool1::D0(*((Curve1*)myC1),myU+aDelta,P2);
0270 Tool1::D0(*((Curve1*)myC1),myU+2*aDelta,P3);
0271 IsParameterGrown = Standard_True;
0272 }
0273 else
0274 {
0275 Tool1::D0(*((Curve1*)myC1),myU-2*aDelta,P1);
0276 Tool1::D0(*((Curve1*)myC1),myU-aDelta,P2);
0277 Tool1::D0(*((Curve1*)myC1),myU,P3);
0278 IsParameterGrown = Standard_False;
0279 }
0280
0281 Vec V1(Ptemp,P1), V2(Ptemp,P2), V3(Ptemp,P3);
0282
0283 if(IsParameterGrown)
0284 myDu=-3*V1+4*V2-V3;
0285 else
0286 myDu=V1-4*V2+3*V3;
0287 }//else of if(IsDeriveFound)
0288 Ndu = myDu.Magnitude();
0289 }//if (Ndu <= myTolC1) condition
0290 }//if(myMaxDerivOrder != 0)
0291
0292 if (Ndu <= MinTol)
0293 {
0294 //Warning: 1st derivative of C1 is equal to zero!
0295 return Standard_False;
0296 }
0297
0298 Standard_Real Ndv = myDv.Magnitude();
0299
0300 if(myMaxDerivOrderC2 != 0)
0301 {
0302 if (Ndv <= myTolC2)
0303 {
0304 const Standard_Real DivisionFactor = 1.e-3;
0305 Standard_Real dv;
0306 if((myVsupremum >= RealLast()) || (myVinfium <= RealFirst()))
0307 dv = 0.0;
0308 else
0309 dv = myVsupremum-myVinfium;
0310
0311 const Standard_Real aDelta = Max(dv*DivisionFactor,MinStep);
0312
0313 //Derivative is approximated by Taylor-series
0314
0315 Standard_Integer n = 1; //Derivative order
0316 Vec V;
0317 Standard_Boolean IsDeriveFound;
0318
0319 do
0320 {
0321 V = Tool2::DN(*((Curve2*)myC2),myV,++n);
0322 Ndv = V.Magnitude();
0323 IsDeriveFound = (Ndv > myTolC2);
0324 }
0325 while(!IsDeriveFound && n < myMaxDerivOrderC2);
0326
0327 if(IsDeriveFound)
0328 {
0329 Standard_Real v;
0330
0331 if(myV-myVinfium < aDelta)
0332 v = myV+aDelta;
0333 else
0334 v = myV-aDelta;
0335
0336 Pnt P1, P2;
0337 Tool2::D0(*((Curve2*)myC2),Min(myV, v),P1);
0338 Tool2::D0(*((Curve2*)myC2),Max(myV, v),P2);
0339
0340 Vec V1(P1,P2);
0341 Standard_Real aDirFactor = V.Dot(V1);
0342
0343 if(aDirFactor < 0.0)
0344 myDv = -V;
0345 else
0346 myDv = V;
0347 }//if(IsDeriveFound)
0348 else
0349 {
0350 //Derivative is approximated by three points
0351
0352 Pnt Ptemp; //(0,0,0)-coordinate
0353 Pnt P1, P2, P3;
0354 Standard_Boolean IsParameterGrown;
0355
0356 if(myV-myVinfium < 2*aDelta)
0357 {
0358 Tool2::D0(*((Curve2*)myC2),myV,P1);
0359 Tool2::D0(*((Curve2*)myC2),myV+aDelta,P2);
0360 Tool2::D0(*((Curve2*)myC2),myV+2*aDelta,P3);
0361 IsParameterGrown = Standard_True;
0362 }
0363 else
0364 {
0365 Tool2::D0(*((Curve2*)myC2),myV-2*aDelta,P1);
0366 Tool2::D0(*((Curve2*)myC2),myV-aDelta,P2);
0367 Tool2::D0(*((Curve2*)myC2),myV,P3);
0368 IsParameterGrown = Standard_False;
0369 }
0370
0371 Vec V1(Ptemp,P1), V2(Ptemp,P2), V3(Ptemp,P3);
0372
0373 if(IsParameterGrown)
0374 myDv=-3*V1+4*V2-V3;
0375 else
0376 myDv=V1-4*V2+3*V3;
0377 }//else of if(IsDeriveFound)
0378
0379 Ndv = myDv.Magnitude();
0380 }//if (Ndv <= myTolC2)
0381 }//if(myMaxDerivOrder != 0)
0382
0383 if (Ndv <= MinTol)
0384 {
0385 //1st derivative of C2 is equal to zero!
0386 return Standard_False;
0387 }
0388
0389 F(1) = P1P2.Dot(myDu)/Ndu;
0390 F(2) = P1P2.Dot(myDv)/Ndv;
0391 return Standard_True;
0392 }
0393 //=============================================================================
0394
0395 Standard_Boolean Extrema_FuncExtCC::Derivatives (const math_Vector& UV,
0396 math_Matrix& Df)
0397 {
0398 math_Vector F(1,2);
0399 return Values(UV,F,Df);
0400 }
0401 //=============================================================================
0402
0403 Standard_Boolean Extrema_FuncExtCC::Values (const math_Vector& UV,
0404 math_Vector& F,
0405 math_Matrix& Df)
0406 {
0407 myU = UV(1);
0408 myV = UV(2);
0409
0410 if(Value(UV, F) == Standard_False) //Computes F, myDu, myDv
0411 {
0412 //Warning: No function value found!
0413 return Standard_False;
0414 }//if(Value(UV, F) == Standard_False)
0415
0416 Vec Du, Dv, Duu, Dvv;
0417 Tool1::D2(*((Curve1*)myC1), myU,myP1,Du,Duu);
0418 Tool2::D2(*((Curve2*)myC2), myV,myP2,Dv,Dvv);
0419
0420 //Calling of "Value(...)" function change class member values.
0421 //After running it is necessary to return to previous values.
0422 const Standard_Real myU_old = myU, myV_old = myV;
0423 const Pnt myP1_old = myP1, myP2_old = myP2;
0424 const Vec myDu_old = myDu, myDv_old = myDv;
0425
0426
0427 //Attention: aDelta value must be greater than same value for "Value(...)"
0428 // function to avoid of points' collisions.
0429
0430 const Standard_Real DivisionFactor = 0.01;
0431
0432 Standard_Real du;
0433 if((myUsupremum >= RealLast()) || (myUinfium <= RealFirst()))
0434 du = 0.0;
0435 else
0436 du = myUsupremum-myUinfium;
0437
0438 const Standard_Real aDeltaU = Max(du*DivisionFactor,MinStep);
0439
0440 Standard_Real dv;
0441 if((myVsupremum >= RealLast()) || (myVinfium <= RealFirst()))
0442 dv = 0.0;
0443 else
0444 dv = myVsupremum-myVinfium;
0445
0446 const Standard_Real aDeltaV = Max(dv*DivisionFactor,MinStep);
0447
0448 Vec P1P2 (myP1,myP2);
0449
0450 if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1))
0451 {
0452 //Derivative is approximated by three points
0453
0454 math_Vector FF1(1,2), FF2(1,2), FF3(1,2);
0455 Standard_Real F1, F2, F3;
0456
0457 /////////////////////////// Search of DF1_u derivative (begin) ///////////////////
0458 if(myU-myUinfium < 2*aDeltaU)
0459 {
0460 F1=F(1);
0461 math_Vector UV2(1,2), UV3(1,2);
0462 UV2(1)=myU+aDeltaU;
0463 UV2(2)=myV;
0464 UV3(1)=myU+2*aDeltaU;
0465 UV3(2)=myV;
0466 if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
0467 {
0468 //There are many points close to singularity points and which have zero-derivative.
0469 //Try to decrease aDelta variable's value.
0470 return Standard_False;
0471 }
0472
0473 F2 = FF2(1);
0474 F3 = FF3(1);
0475
0476 Df(1,1) = (-3*F1+4*F2-F3)/(2.0*aDeltaU);
0477 }//if(myU-myUinfium < 2*aDeltaU)
0478 else
0479 {
0480 F3 = F(1);
0481 math_Vector UV2(1,2), UV1(1,2);
0482 UV2(1)=myU-aDeltaU;
0483 UV2(2)=myV;
0484 UV1(1)=myU-2*aDeltaU;
0485 UV1(2)=myV;
0486
0487 if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
0488 {
0489 //There are many points close to singularity points and which have zero-derivative.
0490 //Try to decrease aDelta variable's value.
0491 return Standard_False;
0492 }
0493
0494 F1 = FF1(1);
0495 F2 = FF2(1);
0496
0497 Df(1,1) = (F1-4*F2+3*F3)/(2.0*aDeltaU);
0498 }//else of if(myU-myUinfium < 2*aDeltaU) condition
0499 /////////////////////////// Search of DF1_u derivative (end) ///////////////////
0500
0501 //Return to previous values
0502 myU = myU_old;
0503 myV = myV_old;
0504
0505 /////////////////////////// Search of DF1_v derivative (begin) ///////////////////
0506 if(myV-myVinfium < 2*aDeltaV)
0507 {
0508 F1=F(1);
0509 math_Vector UV2(1,2), UV3(1,2);
0510 UV2(1)=myU;
0511 UV2(2)=myV+aDeltaV;
0512 UV3(1)=myU;
0513 UV3(2)=myV+2*aDeltaV;
0514
0515 if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
0516 {
0517 //There are many points close to singularity points and which have zero-derivative.
0518 //Try to decrease aDelta variable's value.
0519 return Standard_False;
0520 }
0521 F2 = FF2(1);
0522 F3 = FF3(1);
0523
0524 Df(1,2) = (-3*F1+4*F2-F3)/(2.0*aDeltaV);
0525 }//if(myV-myVinfium < 2*aDeltaV)
0526 else
0527 {
0528 F3 = F(1);
0529 math_Vector UV2(1,2), UV1(1,2);
0530 UV2(1)=myU;
0531 UV2(2)=myV-aDeltaV;
0532 UV1(1)=myU;
0533 UV1(2)=myV-2*aDeltaV;
0534 if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
0535 {
0536 //There are many points close to singularity points and which have zero-derivative.
0537 //Try to decrease aDelta variable's value.
0538 return Standard_False;
0539 }
0540
0541 F1 = FF1(1);
0542 F2 = FF2(1);
0543
0544 Df(1,2) = (F1-4*F2+3*F3)/(2.0*aDeltaV);
0545 }//else of if(myV-myVinfium < 2*aDeltaV)
0546 /////////////////////////// Search of DF1_v derivative (end) ///////////////////
0547
0548 //Return to previous values
0549 myU = myU_old;
0550 myV = myV_old;
0551 myP1 = myP1_old, myP2 = myP2_old;
0552 myDu = myDu_old, myDv = myDv_old;
0553 }//if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1))
0554 else
0555 {
0556 const Standard_Real Ndu = myDu.Magnitude();
0557 Df(1,1) = - Ndu + (P1P2.Dot(Duu)/Ndu) - F(1)*(myDu.Dot(Duu)/(Ndu*Ndu));
0558 Df(1,2) = myDv.Dot(myDu)/Ndu;
0559 }//else of if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1))
0560
0561 if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2))
0562 {
0563 //Derivative is approximated by three points
0564
0565 math_Vector FF1(1,2), FF2(1,2), FF3(1,2);
0566 Standard_Real F1, F2, F3;
0567
0568 /////////////////////////// Search of DF2_v derivative (begin) ///////////////////
0569 if(myV-myVinfium < 2*aDeltaV)
0570 {
0571 F1=F(2);
0572 math_Vector UV2(1,2), UV3(1,2);
0573 UV2(1)=myU;
0574 UV2(2)=myV+aDeltaV;
0575 UV3(1)=myU;
0576 UV3(2)=myV+2*aDeltaV;
0577
0578 if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
0579 {
0580 //There are many points close to singularity points and which have zero-derivative.
0581 //Try to decrease aDelta variable's value.
0582 return Standard_False;
0583 }
0584
0585 F2 = FF2(2);
0586 F3 = FF3(2);
0587
0588 Df(2,2) = (-3*F1+4*F2-F3)/(2.0*aDeltaV);
0589
0590 }//if(myV-myVinfium < 2*aDeltaV)
0591 else
0592 {
0593 F3 = F(2);
0594 math_Vector UV2(1,2), UV1(1,2);
0595 UV2(1)=myU;
0596 UV2(2)=myV-aDeltaV;
0597 UV1(1)=myU;
0598 UV1(2)=myV-2*aDeltaV;
0599
0600 if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
0601 {
0602 //There are many points close to singularity points and which have zero-derivative.
0603 //Try to decrease aDelta variable's value.
0604 return Standard_False;
0605 }
0606
0607 F1 = FF1(2);
0608 F2 = FF2(2);
0609
0610 Df(2,2) = (F1-4*F2+3*F3)/(2.0*aDeltaV);
0611 }//else of if(myV-myVinfium < 2*aDeltaV)
0612 /////////////////////////// Search of DF2_v derivative (end) ///////////////////
0613
0614 //Return to previous values
0615 myU = myU_old;
0616 myV = myV_old;
0617
0618 /////////////////////////// Search of DF2_u derivative (begin) ///////////////////
0619 if(myU-myUinfium < 2*aDeltaU)
0620 {
0621 F1=F(2);
0622 math_Vector UV2(1,2), UV3(1,2);
0623 UV2(1)=myU+aDeltaU;
0624 UV2(2)=myV;
0625 UV3(1)=myU+2*aDeltaU;
0626 UV3(2)=myV;
0627 if(!((Value(UV2,FF2)) && (Value(UV3,FF3))))
0628 {
0629 //There are many points close to singularity points and which have zero-derivative.
0630 //Try to decrease aDelta variable's value.
0631 return Standard_False;
0632 }
0633
0634 F2 = FF2(2);
0635 F3 = FF3(2);
0636
0637 Df(2,1) = (-3*F1+4*F2-F3)/(2.0*aDeltaU);
0638 }//if(myU-myUinfium < 2*aDelta)
0639 else
0640 {
0641 F3 = F(2);
0642 math_Vector UV2(1,2), UV1(1,2);
0643 UV2(1)=myU-aDeltaU;
0644 UV2(2)=myV;
0645 UV1(1)=myU-2*aDeltaU;
0646 UV1(2)=myV;
0647
0648 if(!((Value(UV2,FF2)) && (Value(UV1,FF1))))
0649 {
0650 //There are many points close to singularity points
0651 //and which have zero-derivative.
0652 //Try to decrease aDelta variable's value.
0653 return Standard_False;
0654 }
0655
0656 F1 = FF1(2);
0657 F2 = FF2(2);
0658
0659 Df(2,1) = (F1-4*F2+3*F3)/(2.0*aDeltaU);
0660 }//else of if(myU-myUinfium < 2*aDeltaU)
0661 /////////////////////////// Search of DF2_u derivative (end) ///////////////////
0662
0663 //Return to previous values
0664 myU = myU_old;
0665 myV = myV_old;
0666 myP1 = myP1_old;
0667 myP2 = myP2_old;
0668 myDu = myDu_old;
0669 myDv = myDv_old;
0670 }//if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2))
0671 else
0672 {
0673 Standard_Real Ndv = myDv.Magnitude();
0674 Df(2,2) = Ndv + (P1P2.Dot(Dvv)/Ndv) - F(2)*(myDv.Dot(Dvv)/(Ndv*Ndv));
0675 Df(2,1) = -myDu.Dot(myDv)/Ndv;
0676 }//else of if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2))
0677
0678 return Standard_True;
0679
0680 }//end of function
0681 //=============================================================================
0682
0683 Standard_Integer Extrema_FuncExtCC::GetStateNumber ()
0684 {
0685 Vec Du (myDu), Dv (myDv);
0686 Vec P1P2 (myP1, myP2);
0687
0688 Standard_Real mod = Du.Magnitude();
0689 if(mod > myTolC1) {
0690 Du /= mod;
0691 }
0692 mod = Dv.Magnitude();
0693 if(mod > myTolC2) {
0694 Dv /= mod;
0695 }
0696
0697 if (Abs(P1P2.Dot(Du)) <= myTol && Abs(P1P2.Dot(Dv)) <= myTol) {
0698 mySqDist.Append(myP1.SquareDistance(myP2));
0699 myPoints.Append(POnC(myU, myP1));
0700 myPoints.Append(POnC(myV, myP2));
0701 }
0702 return 0;
0703 }
0704 //=============================================================================
0705
0706 void Extrema_FuncExtCC::Points (const Standard_Integer N,
0707 POnC& P1,
0708 POnC& P2) const
0709 {
0710 P1 = myPoints.Value(2*N-1);
0711 P2 = myPoints.Value(2*N);
0712 }
0713 //=============================================================================
0714
0715 void Extrema_FuncExtCC::SubIntervalInitialize(const math_Vector& theInfBound,
0716 const math_Vector& theSupBound)
0717 {
0718 myUinfium = theInfBound(1);
0719 myUsupremum = theSupBound(1);
0720 myVinfium = theInfBound(2);
0721 myVsupremum = theSupBound(2);
0722 }
0723 //=============================================================================