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