0022922: Clean up warnings on uninitialized / unused variables
[occt.git] / src / math / math_BrentMinimum.cxx
CommitLineData
7fd59977 1#include <math_BrentMinimum.ixx>
2#include <math_Function.hxx>
3
4#define CGOLD 0.3819660
5#ifdef MAX
6#undef MAX
7#endif
8#define MAX(a,b) ((a) > (b) ? (a) : (b))
9#define SIGN(a,b) ((b) > 0.0 ? fabs(a) : -fabs(a))
10#define SHFT(a,b,c,d) (a)=(b);(b)=(c);(c)=(d)
11
12void math_BrentMinimum::Perform(math_Function& F,
13 const Standard_Real ax,
14 const Standard_Real bx,
15 const Standard_Real cx) {
16
17 Standard_Boolean OK;
18 Standard_Real etemp, fu, p, q, r;
19 Standard_Real tol1, tol2, u, v, w, xm;
20 Standard_Real e = 0.0;
21 Standard_Real d = RealLast();
22
23 a = ((ax < cx) ? ax : cx);
24 b = ((ax > cx) ? ax : cx);
25 x = w = v = bx;
26 if (!myF) {
27 OK = F.Value(x, fx);
28 if(!OK) return;
29 }
30 fw = fv = fx;
31 for(iter = 1; iter <= Itermax; iter++) {
32 xm = 0.5 * (a + b);
33 tol1 = XTol * fabs(x) + EPSZ;
34 tol2 = 2.0 * tol1;
35 if(IsSolutionReached(F)) {
36 Done = Standard_True;
37 return;
38 }
39 if(fabs(e) > tol1) {
40 r = (x - w) * (fx - fv);
41 q = (x - v) * (fx - fw);
42 p = (x - v) * q - (x - w) * r;
43 q = 2.0 * (q - r);
44 if(q > 0.0) p = -p;
45 q = fabs(q);
46 etemp = e;
47 e = d;
48 if(fabs(p) >= fabs(0.5 * q * etemp)
49 || p <= q * ( a - x) || p >= q * (b - x)) {
50 e = (x >= xm ? a - x : b - x);
51 d = CGOLD * e;
52 }
53 else {
54 d = p / q;
55 u = x + d;
56 if(u - a < tol2 || b - u < tol2) d = SIGN(tol1, xm - x);
57 }
58 }
59 else {
60 e = (x >= xm ? a - x : b - x);
61 d = CGOLD * e;
62 }
63 u = (fabs(d) >= tol1 ? x + d : x + SIGN(tol1, d));
64 OK = F.Value(u, fu);
65 if(!OK) return;
66 if(fu <= fx) {
67 if(u >= x) a = x; else b = x;
68 SHFT(v, w, x, u);
69 SHFT(fv, fw, fx, fu);
70 }
71 else {
72 if(u < x) a = u; else b = u;
73 if(fu <= fw || w == x) {
74 v = w;
75 w = u;
76 fv = fw;
77 fw = fu;
78 }
79 else if(fu <= fv || v == x || v == w) {
80 v = u;
81 fv = fu;
82 }
83 }
84 }
85 Done = Standard_False;
86 return;
87}
88
89
90math_BrentMinimum::math_BrentMinimum(math_Function& F,
91 const Standard_Real Ax,
92 const Standard_Real Bx,
93 const Standard_Real Cx,
94 const Standard_Real TolX,
95 const Standard_Integer NbIterations,
96 const Standard_Real ZEPS) {
97
98 XTol = TolX;
99 EPSZ = ZEPS;
100 Itermax = NbIterations;
101 myF = Standard_False;
102 Perform(F, Ax, Bx, Cx);
103}
104
105
106// Constructeur d'initialisation des champs.
107
108math_BrentMinimum::math_BrentMinimum(const Standard_Real TolX,
109 const Standard_Integer NbIterations,
110 const Standard_Real ZEPS) {
111 myF = Standard_False;
112 XTol = TolX;
113 EPSZ = ZEPS;
114 Itermax = NbIterations;
115}
116
117math_BrentMinimum::math_BrentMinimum(const Standard_Real TolX,
118 const Standard_Real Fbx,
119 const Standard_Integer NbIterations,
120 const Standard_Real ZEPS) {
121
122 fx = Fbx;
123 myF = Standard_True;
124 XTol = TolX;
125 EPSZ = ZEPS;
126 Itermax = NbIterations;
127}
128
129
130// Standard_Boolean math_BrentMinimum::IsSolutionReached(math_Function& F) {
131 Standard_Boolean math_BrentMinimum::IsSolutionReached(math_Function& ) {
132
6e6cd5d9 133// Standard_Real xm = 0.5 * (a + b);
7fd59977 134 // modified by NIZHNY-MKK Mon Oct 3 17:45:57 2005.BEGIN
135// Standard_Real tol = XTol * fabs(x) + EPSZ;
136// return fabs(x - xm) <= 2.0 * tol - 0.5 * (b - a);
137 Standard_Real TwoTol = 2.0 *(XTol * fabs(x) + EPSZ);
138 return ((x <= (TwoTol + a)) && (x >= (b - TwoTol)));
139 // modified by NIZHNY-MKK Mon Oct 3 17:46:00 2005.END
140 }
141
142
143
144 void math_BrentMinimum::Dump(Standard_OStream& o) const {
145 o << "math_BrentMinimum ";
146 if(Done) {
147 o << " Status = Done \n";
148 o << " Location value = " << x <<"\n";
149 o << " Minimum value = " << fx << "\n";
150 o << " Number of iterations = " << iter <<"\n";
151 }
152 else {
153 o << " Status = not Done \n";
154 }
155 }