Commit | Line | Data |
---|---|---|
b311480e | 1 | // Copyright (c) 1995-1999 Matra Datavision |
2 | // Copyright (c) 1999-2012 OPEN CASCADE SAS | |
3 | // | |
4 | // The content of this file is subject to the Open CASCADE Technology Public | |
5 | // License Version 6.5 (the "License"). You may not use the content of this file | |
6 | // except in compliance with the License. Please obtain a copy of the License | |
7 | // at http://www.opencascade.org and read it completely before using this file. | |
8 | // | |
9 | // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its | |
10 | // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France. | |
11 | // | |
12 | // The Original Code and all software distributed under the License is | |
13 | // distributed on an "AS IS" basis, without warranty of any kind, and the | |
14 | // Initial Developer hereby disclaims all such warranties, including without | |
15 | // limitation, any warranties of merchantability, fitness for a particular | |
16 | // purpose or non-infringement. Please see the License for the specific terms | |
17 | // and conditions governing the rights and limitations under the License. | |
18 | ||
7fd59977 | 19 | //Modified by MPS : 21-05-97 PRO 7598 |
20 | // Le nombre de solutions rendu est mySqDist.Length() et non | |
21 | // mySqDist.Length()/2 | |
22 | // ajout des valeurs absolues dans le test d'orthogonalite de | |
23 | // GetStateNumber() | |
24 | ||
25 | /*----------------------------------------------------------------------------- | |
28cec2ba | 26 | Fonctions permettant de rechercher une distance extremale entre 2 courbes |
7fd59977 | 27 | C1 et C2 (en partant de points approches C1(u0) et C2(v0)). |
28cec2ba | 28 | Cette classe herite de math_FunctionSetWithDerivatives et est utilisee par |
7fd59977 | 29 | l'algorithme math_FunctionSetRoot. |
28cec2ba | 30 | Si on note Du et Dv, les derivees en u et v, les 2 fonctions a annuler sont: |
31 | { F1(u,v) = (C2(v)-C1(u)).Du(u)/||Du|| } | |
32 | { F2(u,v) = (C2(v)-C1(u)).Dv(v)||Dv|| } | |
33 | Si on note Duu et Dvv, les derivees secondes de C1 et C2, les derivees de F1 | |
7fd59977 | 34 | et F2 sont egales a: |
28cec2ba | 35 | { Duf1(u,v) = -||Du|| + C1C2.Duu/||Du||- F1(u,v)*Duu*Du/||Du||**2 |
36 | { Dvf1(u,v) = Dv.Du/||Du|| | |
37 | { Duf2(u,v) = -Du.Dv/||Dv|| | |
38 | { Dvf2(u,v) = ||Dv|| + C2C1.Dvv/||Dv||- F2(u,v)*Dv*Dvv/||Dv||**2 | |
7fd59977 | 39 | |
40 | ----------------------------------------------------------------------------*/ | |
7fd59977 | 41 | |
32ca7a51 | 42 | #include <Precision.hxx> |
43 | ||
44 | ||
45 | static const Standard_Real MinTol = 1.e-20; | |
46 | static const Standard_Real TolFactor = 1.e-12; | |
47 | static const Standard_Real MinStep = 1e-7; | |
48 | ||
49 | static const Standard_Integer MaxOrder = 3; | |
50 | ||
51 | ||
52 | ||
53 | //============================================================================= | |
54 | Standard_Real Extrema_FuncExtCC::SearchOfTolerance(const Standard_Address C) | |
28cec2ba | 55 | { |
32ca7a51 | 56 | const Standard_Integer NPoint = 10; |
57 | Standard_Real aStartParam, anEndParam; | |
28cec2ba | 58 | |
32ca7a51 | 59 | if(C==myC1) |
28cec2ba | 60 | { |
32ca7a51 | 61 | aStartParam = myUinfium; |
62 | anEndParam = myUsupremum; | |
28cec2ba | 63 | } |
32ca7a51 | 64 | else if(C==myC2) |
28cec2ba | 65 | { |
32ca7a51 | 66 | aStartParam = myVinfium; |
67 | anEndParam = myVsupremum; | |
28cec2ba | 68 | } |
32ca7a51 | 69 | else |
28cec2ba | 70 | { |
71 | //Warning: No curve for tolerance computing! | |
32ca7a51 | 72 | return MinTol; |
28cec2ba | 73 | } |
74 | ||
32ca7a51 | 75 | const Standard_Real aStep = (anEndParam - aStartParam)/(Standard_Real)NPoint; |
28cec2ba | 76 | |
32ca7a51 | 77 | Standard_Integer aNum = 0; |
78 | Standard_Real aMax = -Precision::Infinite(); //Maximum value of 1st derivative | |
28cec2ba | 79 | //(it is computed with using NPoint point) |
80 | ||
32ca7a51 | 81 | do |
28cec2ba | 82 | { |
32ca7a51 | 83 | Standard_Real u = aStartParam + aNum*aStep; //parameter for every point |
84 | if(u > anEndParam) | |
85 | u = anEndParam; | |
28cec2ba | 86 | |
32ca7a51 | 87 | Pnt Ptemp; //empty point (is not used below) |
88 | Vec VDer; // 1st derivative vector | |
89 | Tool1::D1(*((Curve1*)C), u, Ptemp, VDer); | |
90 | Standard_Real vm = VDer.Magnitude(); | |
91 | if(vm > aMax) | |
92 | aMax = vm; | |
28cec2ba | 93 | } |
32ca7a51 | 94 | while(++aNum < NPoint+1); |
28cec2ba | 95 | |
32ca7a51 | 96 | return Max(aMax*TolFactor,MinTol); |
28cec2ba | 97 | } |
7fd59977 | 98 | |
32ca7a51 | 99 | //============================================================================= |
7fd59977 | 100 | Extrema_FuncExtCC::Extrema_FuncExtCC(const Standard_Real thetol) : myC1 (0), myC2 (0), myTol (thetol) |
28cec2ba | 101 | { |
32ca7a51 | 102 | math_Vector V1(1,2), V2(1,2); |
103 | V1(1) = 0.0; | |
104 | V2(1) = 0.0; | |
105 | V1(2) = 0.0; | |
106 | V2(2) = 0.0; | |
107 | SubIntervalInitialize(V1, V2); | |
108 | myMaxDerivOrderC1 = 0; | |
109 | myTolC1=MinTol; | |
110 | myMaxDerivOrderC2 = 0; | |
111 | myTolC2=MinTol; | |
28cec2ba | 112 | } |
7fd59977 | 113 | |
32ca7a51 | 114 | //============================================================================= |
7fd59977 | 115 | Extrema_FuncExtCC::Extrema_FuncExtCC (const Curve1& C1, |
28cec2ba | 116 | const Curve2& C2, |
117 | const Standard_Real thetol) : | |
118 | myC1 ((Standard_Address)&C1), myC2 ((Standard_Address)&C2), | |
119 | myTol (thetol) | |
120 | { | |
32ca7a51 | 121 | math_Vector V1(1,2), V2(1,2); |
122 | V1(1) = Tool1::FirstParameter(*((Curve1*)myC1)); | |
123 | V2(1) = Tool1::LastParameter(*((Curve1*)myC1)); | |
124 | V1(2) = Tool2::FirstParameter(*((Curve2*)myC2)); | |
125 | V2(2) = Tool2::LastParameter(*((Curve2*)myC2)); | |
126 | SubIntervalInitialize(V1, V2); | |
28cec2ba | 127 | |
32ca7a51 | 128 | switch(Tool1::GetType(*((Curve1*)myC1))) |
28cec2ba | 129 | { |
130 | case GeomAbs_BezierCurve: | |
131 | case GeomAbs_BSplineCurve: | |
132 | case GeomAbs_OtherCurve: | |
133 | myMaxDerivOrderC1 = MaxOrder; | |
134 | myTolC1 = SearchOfTolerance((Standard_Address)&C1); | |
135 | break; | |
136 | default: | |
137 | myMaxDerivOrderC1 = 0; | |
138 | myTolC1=MinTol; | |
139 | break; | |
140 | } | |
141 | ||
32ca7a51 | 142 | switch(Tool2::GetType(*((Curve2*)myC2))) |
28cec2ba | 143 | { |
144 | case GeomAbs_BezierCurve: | |
145 | case GeomAbs_BSplineCurve: | |
146 | case GeomAbs_OtherCurve: | |
147 | myMaxDerivOrderC2 = MaxOrder; | |
148 | myTolC2 = SearchOfTolerance((Standard_Address)&C2); | |
149 | break; | |
150 | default: | |
151 | myMaxDerivOrderC2 = 0; | |
152 | myTolC2=MinTol; | |
153 | break; | |
32ca7a51 | 154 | } |
28cec2ba | 155 | } |
7fd59977 | 156 | |
32ca7a51 | 157 | //============================================================================= |
158 | void Extrema_FuncExtCC::SetCurve (const Standard_Integer theRank, const Curve1& C) | |
28cec2ba | 159 | { |
32ca7a51 | 160 | Standard_OutOfRange_Raise_if (theRank < 1 || theRank > 2, "Extrema_FuncExtCC::SetCurve()") |
28cec2ba | 161 | |
162 | if (theRank == 1) | |
32ca7a51 | 163 | { |
28cec2ba | 164 | myC1 = (Standard_Address)&C; |
165 | switch(/*Tool1::GetType(*((Curve1*)myC1))*/ C.GetType()) | |
32ca7a51 | 166 | { |
167 | case GeomAbs_BezierCurve: | |
168 | case GeomAbs_BSplineCurve: | |
169 | case GeomAbs_OtherCurve: | |
170 | myMaxDerivOrderC1 = MaxOrder; | |
171 | myTolC1 = SearchOfTolerance((Standard_Address)&C); | |
172 | break; | |
173 | default: | |
174 | myMaxDerivOrderC1 = 0; | |
175 | myTolC1=MinTol; | |
176 | break; | |
177 | } | |
178 | } | |
28cec2ba | 179 | else if (theRank == 2) |
32ca7a51 | 180 | { |
28cec2ba | 181 | myC2 = (Standard_Address)&C; |
182 | switch(/*Tool2::GetType(*((Curve2*)myC2))*/C.GetType()) | |
32ca7a51 | 183 | { |
184 | case GeomAbs_BezierCurve: | |
185 | case GeomAbs_BSplineCurve: | |
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 | } | |
28cec2ba | 196 | } |
32ca7a51 | 197 | |
198 | //============================================================================= | |
199 | Standard_Boolean Extrema_FuncExtCC::Value (const math_Vector& UV, math_Vector& F) | |
28cec2ba | 200 | { |
7fd59977 | 201 | myU = UV(1); |
202 | myV = UV(2); | |
0d167958 RL |
203 | Tool1::D1(*((Curve1*)myC1), myU,myP1,myDu); |
204 | Tool2::D1(*((Curve2*)myC2), myV,myP2,myDv); | |
32ca7a51 | 205 | |
7fd59977 | 206 | Vec P1P2 (myP1,myP2); |
207 | ||
0d167958 | 208 | Standard_Real Ndu = myDu.Magnitude(); |
32ca7a51 | 209 | |
210 | ||
211 | if(myMaxDerivOrderC1 != 0) | |
28cec2ba | 212 | { |
32ca7a51 | 213 | if (Ndu <= myTolC1) |
28cec2ba | 214 | { |
215 | //Derivative is approximated by Taylor-series | |
32ca7a51 | 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; | |
28cec2ba | 222 | |
32ca7a51 | 223 | const Standard_Real aDelta = Max(du*DivisionFactor,MinStep); |
28cec2ba | 224 | |
32ca7a51 | 225 | Standard_Integer n = 1; //Derivative order |
226 | Vec V; | |
227 | Standard_Boolean IsDeriveFound; | |
28cec2ba | 228 | |
32ca7a51 | 229 | do |
28cec2ba | 230 | { |
32ca7a51 | 231 | V = Tool1::DN(*((Curve1*)myC1),myU,++n); |
232 | Ndu = V.Magnitude(); | |
233 | IsDeriveFound = (Ndu > myTolC1); | |
28cec2ba | 234 | } |
32ca7a51 | 235 | while(!IsDeriveFound && n < myMaxDerivOrderC1); |
28cec2ba | 236 | |
32ca7a51 | 237 | if(IsDeriveFound) |
28cec2ba | 238 | { |
32ca7a51 | 239 | Standard_Real u; |
28cec2ba | 240 | |
32ca7a51 | 241 | if(myU-myUinfium < aDelta) |
242 | u = myU+aDelta; | |
243 | else | |
244 | u = myU-aDelta; | |
28cec2ba | 245 | |
32ca7a51 | 246 | Pnt P1, P2; |
247 | Tool1::D0(*((Curve1*)myC1),Min(myU, u),P1); | |
248 | Tool1::D0(*((Curve1*)myC1),Max(myU, u),P2); | |
28cec2ba | 249 | |
32ca7a51 | 250 | Vec V1(P1,P2); |
251 | Standard_Real aDirFactor = V.Dot(V1); | |
28cec2ba | 252 | |
32ca7a51 | 253 | if(aDirFactor < 0.0) |
254 | myDu = -V; | |
255 | else | |
256 | myDu = V; | |
28cec2ba | 257 | }//if(IsDeriveFound) |
32ca7a51 | 258 | else |
28cec2ba | 259 | { |
260 | //Derivative is approximated by three points | |
32ca7a51 | 261 | |
262 | Pnt Ptemp; //(0,0,0)-coordinate | |
263 | Pnt P1, P2, P3; | |
264 | Standard_Boolean IsParameterGrown; | |
28cec2ba | 265 | |
32ca7a51 | 266 | if(myU-myUinfium < 2*aDelta) |
28cec2ba | 267 | { |
32ca7a51 | 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; | |
28cec2ba | 272 | } |
32ca7a51 | 273 | else |
28cec2ba | 274 | { |
32ca7a51 | 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; | |
28cec2ba | 279 | } |
280 | ||
32ca7a51 | 281 | Vec V1(Ptemp,P1), V2(Ptemp,P2), V3(Ptemp,P3); |
28cec2ba | 282 | |
32ca7a51 | 283 | if(IsParameterGrown) |
284 | myDu=-3*V1+4*V2-V3; | |
285 | else | |
286 | myDu=V1-4*V2+3*V3; | |
28cec2ba | 287 | }//else of if(IsDeriveFound) |
32ca7a51 | 288 | Ndu = myDu.Magnitude(); |
28cec2ba | 289 | }//if (Ndu <= myTolC1) condition |
290 | }//if(myMaxDerivOrder != 0) | |
32ca7a51 | 291 | |
292 | if (Ndu <= MinTol) | |
28cec2ba | 293 | { |
294 | //Warning: 1st derivative of C1 is equal to zero! | |
32ca7a51 | 295 | return Standard_False; |
28cec2ba | 296 | } |
7fd59977 | 297 | |
0d167958 | 298 | Standard_Real Ndv = myDv.Magnitude(); |
32ca7a51 | 299 | |
300 | if(myMaxDerivOrderC2 != 0) | |
28cec2ba | 301 | { |
32ca7a51 | 302 | if (Ndv <= myTolC2) |
28cec2ba | 303 | { |
32ca7a51 | 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; | |
28cec2ba | 310 | |
32ca7a51 | 311 | const Standard_Real aDelta = Max(dv*DivisionFactor,MinStep); |
312 | ||
28cec2ba | 313 | //Derivative is approximated by Taylor-series |
314 | ||
32ca7a51 | 315 | Standard_Integer n = 1; //Derivative order |
316 | Vec V; | |
317 | Standard_Boolean IsDeriveFound; | |
28cec2ba | 318 | |
32ca7a51 | 319 | do |
28cec2ba | 320 | { |
32ca7a51 | 321 | V = Tool2::DN(*((Curve2*)myC2),myV,++n); |
322 | Ndv = V.Magnitude(); | |
323 | IsDeriveFound = (Ndv > myTolC2); | |
28cec2ba | 324 | } |
32ca7a51 | 325 | while(!IsDeriveFound && n < myMaxDerivOrderC2); |
28cec2ba | 326 | |
32ca7a51 | 327 | if(IsDeriveFound) |
28cec2ba | 328 | { |
32ca7a51 | 329 | Standard_Real v; |
28cec2ba | 330 | |
32ca7a51 | 331 | if(myV-myVinfium < aDelta) |
332 | v = myV+aDelta; | |
333 | else | |
334 | v = myV-aDelta; | |
28cec2ba | 335 | |
32ca7a51 | 336 | Pnt P1, P2; |
337 | Tool2::D0(*((Curve2*)myC2),Min(myV, v),P1); | |
338 | Tool2::D0(*((Curve2*)myC2),Max(myV, v),P2); | |
28cec2ba | 339 | |
32ca7a51 | 340 | Vec V1(P1,P2); |
341 | Standard_Real aDirFactor = V.Dot(V1); | |
28cec2ba | 342 | |
32ca7a51 | 343 | if(aDirFactor < 0.0) |
344 | myDv = -V; | |
345 | else | |
346 | myDv = V; | |
28cec2ba | 347 | }//if(IsDeriveFound) |
32ca7a51 | 348 | else |
28cec2ba | 349 | { |
350 | //Derivative is approximated by three points | |
32ca7a51 | 351 | |
352 | Pnt Ptemp; //(0,0,0)-coordinate | |
353 | Pnt P1, P2, P3; | |
354 | Standard_Boolean IsParameterGrown; | |
28cec2ba | 355 | |
32ca7a51 | 356 | if(myV-myVinfium < 2*aDelta) |
28cec2ba | 357 | { |
32ca7a51 | 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; | |
28cec2ba | 362 | } |
32ca7a51 | 363 | else |
28cec2ba | 364 | { |
32ca7a51 | 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; | |
28cec2ba | 369 | } |
370 | ||
32ca7a51 | 371 | Vec V1(Ptemp,P1), V2(Ptemp,P2), V3(Ptemp,P3); |
28cec2ba | 372 | |
32ca7a51 | 373 | if(IsParameterGrown) |
374 | myDv=-3*V1+4*V2-V3; | |
375 | else | |
376 | myDv=V1-4*V2+3*V3; | |
28cec2ba | 377 | }//else of if(IsDeriveFound) |
378 | ||
379 | Ndv = myDv.Magnitude(); | |
380 | }//if (Ndv <= myTolC2) | |
381 | }//if(myMaxDerivOrder != 0) | |
32ca7a51 | 382 | |
383 | if (Ndv <= MinTol) | |
28cec2ba | 384 | { |
385 | //1st derivative of C2 is equal to zero! | |
32ca7a51 | 386 | return Standard_False; |
28cec2ba | 387 | } |
7fd59977 | 388 | |
0d167958 RL |
389 | F(1) = P1P2.Dot(myDu)/Ndu; |
390 | F(2) = P1P2.Dot(myDv)/Ndv; | |
7fd59977 | 391 | return Standard_True; |
28cec2ba | 392 | } |
7fd59977 | 393 | //============================================================================= |
394 | ||
395 | Standard_Boolean Extrema_FuncExtCC::Derivatives (const math_Vector& UV, | |
28cec2ba | 396 | math_Matrix& Df) |
397 | { | |
7fd59977 | 398 | math_Vector F(1,2); |
399 | return Values(UV,F,Df); | |
28cec2ba | 400 | } |
7fd59977 | 401 | //============================================================================= |
32ca7a51 | 402 | |
7fd59977 | 403 | Standard_Boolean Extrema_FuncExtCC::Values (const math_Vector& UV, |
32ca7a51 | 404 | math_Vector& F, |
405 | math_Matrix& Df) | |
28cec2ba | 406 | { |
7fd59977 | 407 | myU = UV(1); |
408 | myV = UV(2); | |
32ca7a51 | 409 | |
410 | if(Value(UV, F) == Standard_False) //Computes F, myDu, myDv | |
28cec2ba | 411 | { |
412 | //Warning: No function value found! | |
32ca7a51 | 413 | return Standard_False; |
28cec2ba | 414 | }//if(Value(UV, F) == Standard_False) |
32ca7a51 | 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); | |
28cec2ba | 419 | |
420 | //Calling of "Value(...)" function change class member values. | |
421 | //After running it is necessary to return to previous values. | |
32ca7a51 | 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; | |
32ca7a51 | 425 | |
28cec2ba | 426 | |
427 | //Attention: aDelta value must be greater than same value for "Value(...)" | |
428 | // function to avoid of points' collisions. | |
32ca7a51 | 429 | |
430 | const Standard_Real DivisionFactor = 0.01; | |
28cec2ba | 431 | |
32ca7a51 | 432 | Standard_Real du; |
433 | if((myUsupremum >= RealLast()) || (myUinfium <= RealFirst())) | |
434 | du = 0.0; | |
435 | else | |
436 | du = myUsupremum-myUinfium; | |
28cec2ba | 437 | |
32ca7a51 | 438 | const Standard_Real aDeltaU = Max(du*DivisionFactor,MinStep); |
28cec2ba | 439 | |
32ca7a51 | 440 | Standard_Real dv; |
441 | if((myVsupremum >= RealLast()) || (myVinfium <= RealFirst())) | |
442 | dv = 0.0; | |
443 | else | |
444 | dv = myVsupremum-myVinfium; | |
28cec2ba | 445 | |
32ca7a51 | 446 | const Standard_Real aDeltaV = Max(dv*DivisionFactor,MinStep); |
7fd59977 | 447 | |
448 | Vec P1P2 (myP1,myP2); | |
28cec2ba | 449 | |
32ca7a51 | 450 | if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1)) |
28cec2ba | 451 | { |
452 | //Derivative is approximated by three points | |
7fd59977 | 453 | |
32ca7a51 | 454 | math_Vector FF1(1,2), FF2(1,2), FF3(1,2); |
455 | Standard_Real F1, F2, F3; | |
7fd59977 | 456 | |
28cec2ba | 457 | /////////////////////////// Search of DF1_u derivative (begin) /////////////////// |
32ca7a51 | 458 | if(myU-myUinfium < 2*aDeltaU) |
28cec2ba | 459 | { |
32ca7a51 | 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)))) | |
28cec2ba | 467 | { |
468 | //There are many points close to singularity points and which have zero-derivative. | |
469 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 470 | return Standard_False; |
28cec2ba | 471 | } |
472 | ||
32ca7a51 | 473 | F2 = FF2(1); |
474 | F3 = FF3(1); | |
28cec2ba | 475 | |
32ca7a51 | 476 | Df(1,1) = (-3*F1+4*F2-F3)/(2.0*aDeltaU); |
28cec2ba | 477 | }//if(myU-myUinfium < 2*aDeltaU) |
32ca7a51 | 478 | else |
28cec2ba | 479 | { |
32ca7a51 | 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)))) | |
28cec2ba | 488 | { |
489 | //There are many points close to singularity points and which have zero-derivative. | |
490 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 491 | return Standard_False; |
28cec2ba | 492 | } |
493 | ||
32ca7a51 | 494 | F1 = FF1(1); |
495 | F2 = FF2(1); | |
496 | ||
497 | Df(1,1) = (F1-4*F2+3*F3)/(2.0*aDeltaU); | |
28cec2ba | 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) /////////////////// | |
32ca7a51 | 506 | if(myV-myVinfium < 2*aDeltaV) |
28cec2ba | 507 | { |
32ca7a51 | 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)))) | |
28cec2ba | 516 | { |
517 | //There are many points close to singularity points and which have zero-derivative. | |
518 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 519 | return Standard_False; |
28cec2ba | 520 | } |
32ca7a51 | 521 | F2 = FF2(1); |
522 | F3 = FF3(1); | |
28cec2ba | 523 | |
32ca7a51 | 524 | Df(1,2) = (-3*F1+4*F2-F3)/(2.0*aDeltaV); |
28cec2ba | 525 | }//if(myV-myVinfium < 2*aDeltaV) |
32ca7a51 | 526 | else |
28cec2ba | 527 | { |
32ca7a51 | 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)))) | |
28cec2ba | 535 | { |
536 | //There are many points close to singularity points and which have zero-derivative. | |
537 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 538 | return Standard_False; |
28cec2ba | 539 | } |
540 | ||
32ca7a51 | 541 | F1 = FF1(1); |
542 | F2 = FF2(1); | |
543 | ||
544 | Df(1,2) = (F1-4*F2+3*F3)/(2.0*aDeltaV); | |
28cec2ba | 545 | }//else of if(myV-myVinfium < 2*aDeltaV) |
546 | /////////////////////////// Search of DF1_v derivative (end) /////////////////// | |
547 | ||
32ca7a51 | 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; | |
28cec2ba | 553 | }//if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1)) |
32ca7a51 | 554 | else |
28cec2ba | 555 | { |
32ca7a51 | 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; | |
28cec2ba | 559 | }//else of if((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1)) |
560 | ||
32ca7a51 | 561 | if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2)) |
28cec2ba | 562 | { |
563 | //Derivative is approximated by three points | |
32ca7a51 | 564 | |
565 | math_Vector FF1(1,2), FF2(1,2), FF3(1,2); | |
566 | Standard_Real F1, F2, F3; | |
567 | ||
28cec2ba | 568 | /////////////////////////// Search of DF2_v derivative (begin) /////////////////// |
32ca7a51 | 569 | if(myV-myVinfium < 2*aDeltaV) |
28cec2ba | 570 | { |
32ca7a51 | 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)))) | |
28cec2ba | 579 | { |
580 | //There are many points close to singularity points and which have zero-derivative. | |
581 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 582 | return Standard_False; |
28cec2ba | 583 | } |
584 | ||
32ca7a51 | 585 | F2 = FF2(2); |
586 | F3 = FF3(2); | |
28cec2ba | 587 | |
32ca7a51 | 588 | Df(2,2) = (-3*F1+4*F2-F3)/(2.0*aDeltaV); |
28cec2ba | 589 | |
590 | }//if(myV-myVinfium < 2*aDeltaV) | |
32ca7a51 | 591 | else |
28cec2ba | 592 | { |
32ca7a51 | 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)))) | |
28cec2ba | 601 | { |
602 | //There are many points close to singularity points and which have zero-derivative. | |
603 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 604 | return Standard_False; |
28cec2ba | 605 | } |
606 | ||
32ca7a51 | 607 | F1 = FF1(2); |
608 | F2 = FF2(2); | |
609 | ||
610 | Df(2,2) = (F1-4*F2+3*F3)/(2.0*aDeltaV); | |
28cec2ba | 611 | }//else of if(myV-myVinfium < 2*aDeltaV) |
612 | /////////////////////////// Search of DF2_v derivative (end) /////////////////// | |
32ca7a51 | 613 | |
28cec2ba | 614 | //Return to previous values |
615 | myU = myU_old; | |
616 | myV = myV_old; | |
32ca7a51 | 617 | |
28cec2ba | 618 | /////////////////////////// Search of DF2_u derivative (begin) /////////////////// |
32ca7a51 | 619 | if(myU-myUinfium < 2*aDeltaU) |
28cec2ba | 620 | { |
32ca7a51 | 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)))) | |
28cec2ba | 628 | { |
629 | //There are many points close to singularity points and which have zero-derivative. | |
630 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 631 | return Standard_False; |
28cec2ba | 632 | } |
633 | ||
32ca7a51 | 634 | F2 = FF2(2); |
635 | F3 = FF3(2); | |
28cec2ba | 636 | |
32ca7a51 | 637 | Df(2,1) = (-3*F1+4*F2-F3)/(2.0*aDeltaU); |
28cec2ba | 638 | }//if(myU-myUinfium < 2*aDelta) |
32ca7a51 | 639 | else |
28cec2ba | 640 | { |
32ca7a51 | 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; | |
28cec2ba | 647 | |
32ca7a51 | 648 | if(!((Value(UV2,FF2)) && (Value(UV1,FF1)))) |
28cec2ba | 649 | { |
650 | //There are many points close to singularity points | |
651 | //and which have zero-derivative. | |
652 | //Try to decrease aDelta variable's value. | |
32ca7a51 | 653 | return Standard_False; |
28cec2ba | 654 | } |
655 | ||
32ca7a51 | 656 | F1 = FF1(2); |
657 | F2 = FF2(2); | |
658 | ||
659 | Df(2,1) = (F1-4*F2+3*F3)/(2.0*aDeltaU); | |
28cec2ba | 660 | }//else of if(myU-myUinfium < 2*aDeltaU) |
661 | /////////////////////////// Search of DF2_u derivative (end) /////////////////// | |
662 | ||
32ca7a51 | 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; | |
28cec2ba | 670 | }//if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2)) |
32ca7a51 | 671 | else |
28cec2ba | 672 | { |
32ca7a51 | 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; | |
28cec2ba | 676 | }//else of if((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2)) |
677 | ||
7fd59977 | 678 | return Standard_True; |
679 | ||
28cec2ba | 680 | }//end of function |
7fd59977 | 681 | //============================================================================= |
682 | ||
683 | Standard_Integer Extrema_FuncExtCC::GetStateNumber () | |
684 | { | |
0d167958 RL |
685 | Vec Du (myDu), Dv (myDv); |
686 | Vec P1P2 (myP1, myP2); | |
687 | ||
7fd59977 | 688 | Standard_Real mod = Du.Magnitude(); |
32ca7a51 | 689 | if(mod > myTolC1) { |
7fd59977 | 690 | Du /= mod; |
691 | } | |
692 | mod = Dv.Magnitude(); | |
32ca7a51 | 693 | if(mod > myTolC2) { |
7fd59977 | 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, | |
28cec2ba | 707 | POnC& P1, |
708 | POnC& P2) const | |
7fd59977 | 709 | { |
710 | P1 = myPoints.Value(2*N-1); | |
711 | P2 = myPoints.Value(2*N); | |
712 | } | |
713 | //============================================================================= | |
714 | ||
32ca7a51 | 715 | void Extrema_FuncExtCC::SubIntervalInitialize(const math_Vector& theInfBound, |
716 | const math_Vector& theSupBound) | |
28cec2ba | 717 | { |
32ca7a51 | 718 | myUinfium = theInfBound(1); |
719 | myUsupremum = theSupBound(1); | |
720 | myVinfium = theInfBound(2); | |
721 | myVsupremum = theSupBound(2); | |
28cec2ba | 722 | } |
723 | //============================================================================= |