b311480e |
1 | // Created on: 1996-05-03 |
2 | // Created by: Philippe MANGIN |
3 | // Copyright (c) 1996-1999 Matra Datavision |
4 | // Copyright (c) 1999-2012 OPEN CASCADE SAS |
5 | // |
6 | // The content of this file is subject to the Open CASCADE Technology Public |
7 | // License Version 6.5 (the "License"). You may not use the content of this file |
8 | // except in compliance with the License. Please obtain a copy of the License |
9 | // at http://www.opencascade.org and read it completely before using this file. |
10 | // |
11 | // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its |
12 | // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France. |
13 | // |
14 | // The Original Code and all software distributed under the License is |
15 | // distributed on an "AS IS" basis, without warranty of any kind, and the |
16 | // Initial Developer hereby disclaims all such warranties, including without |
17 | // limitation, any warranties of merchantability, fitness for a particular |
18 | // purpose or non-infringement. Please see the License for the specific terms |
19 | // and conditions governing the rights and limitations under the License. |
20 | |
7fd59977 |
21 | |
22 | //#ifndef DEB |
23 | #define No_Standard_RangeError |
24 | #define No_Standard_OutOfRange |
25 | #define No_Standard_DimensionError |
26 | //#endif |
27 | |
28 | #include <math_NewtonMinimum.ixx> |
29 | |
30 | #include <math_Gauss.hxx> |
31 | #include <math_Jacobi.hxx> |
32 | |
33 | //============================================================================ |
34 | math_NewtonMinimum::math_NewtonMinimum(math_MultipleVarFunctionWithHessian& F, |
35 | const math_Vector& StartingPoint, |
36 | const Standard_Real Tolerance, |
37 | const Standard_Integer NbIterations, |
38 | const Standard_Real Convexity, |
39 | const Standard_Boolean WithSingularity) |
40 | //============================================================================ |
41 | : TheLocation(1, F.NbVariables()), |
42 | TheGradient(1, F.NbVariables()), |
43 | TheStep(1, F.NbVariables(), 10*Tolerance), |
44 | TheHessian(1, F.NbVariables(), 1, F.NbVariables() ) |
45 | { |
46 | XTol = Tolerance; |
47 | CTol = Convexity; |
48 | Itermax = NbIterations; |
49 | NoConvexTreatement = WithSingularity; |
50 | Convex = Standard_True; |
51 | // Done = Standard_True; |
52 | // TheStatus = math_OK; |
53 | Perform ( F, StartingPoint); |
54 | } |
55 | |
56 | //============================================================================ |
57 | math_NewtonMinimum::math_NewtonMinimum(math_MultipleVarFunctionWithHessian& F, |
58 | const Standard_Real Tolerance, |
59 | const Standard_Integer NbIterations, |
60 | const Standard_Real Convexity, |
61 | const Standard_Boolean WithSingularity) |
62 | //============================================================================ |
63 | : TheLocation(1, F.NbVariables()), |
64 | TheGradient(1, F.NbVariables()), |
65 | TheStep(1, F.NbVariables(), 10*Tolerance), |
66 | TheHessian(1, F.NbVariables(), 1, F.NbVariables() ) |
67 | { |
68 | XTol = Tolerance; |
69 | CTol = Convexity; |
70 | Itermax = NbIterations; |
71 | NoConvexTreatement = WithSingularity; |
72 | Convex = Standard_True; |
73 | Done = Standard_False; |
74 | TheStatus = math_NotBracketed; |
75 | } |
76 | //============================================================================ |
77 | void math_NewtonMinimum::Delete() |
78 | {} |
79 | //============================================================================ |
80 | void math_NewtonMinimum::Perform(math_MultipleVarFunctionWithHessian& F, |
81 | const math_Vector& StartingPoint) |
82 | //============================================================================ |
83 | { |
84 | math_Vector Point1 (1, F.NbVariables()); |
85 | Point1 = StartingPoint; |
86 | math_Vector Point2(1, F.NbVariables()); |
87 | math_Vector* precedent = &Point1; |
88 | math_Vector* suivant = &Point2; |
89 | math_Vector* auxiliaire = precedent; |
90 | |
91 | Standard_Boolean Ok = Standard_True; |
92 | Standard_Integer NbConv = 0, ii, Nreduction; |
93 | Standard_Real VPrecedent, VItere; |
94 | |
95 | Done = Standard_True; |
96 | TheStatus = math_OK; |
97 | nbiter = 0; |
98 | |
99 | while ( Ok && (NbConv < 2) ) { |
100 | nbiter++; |
101 | |
102 | // Positionnement |
103 | |
104 | Ok = F.Values(*precedent, VPrecedent, TheGradient, TheHessian); |
105 | if (!Ok) { |
106 | Done = Standard_False; |
107 | TheStatus = math_FunctionError; |
108 | return; |
109 | } |
110 | if (nbiter==1) { |
111 | PreviousMinimum = VPrecedent; |
112 | TheMinimum = VPrecedent; |
113 | } |
114 | |
115 | // Traitement de la non convexite |
116 | |
117 | math_Jacobi CalculVP(TheHessian); |
118 | if ( !CalculVP.IsDone() ) { |
119 | Done = Standard_False; |
120 | TheStatus = math_FunctionError; |
121 | return; |
122 | } |
123 | |
124 | |
125 | |
126 | MinEigenValue = CalculVP.Values() ( CalculVP.Values().Min()); |
127 | if ( MinEigenValue < CTol) { |
128 | Convex = Standard_False; |
129 | if (NoConvexTreatement) { |
130 | Standard_Real Delta = CTol+0.1*Abs(MinEigenValue) -MinEigenValue ; |
131 | for (ii=1; ii<=TheGradient.Length(); ii++) { |
132 | TheHessian(ii, ii) += Delta; |
133 | } |
134 | } |
135 | else { |
136 | TheStatus = math_FunctionError; |
137 | return; |
138 | } |
139 | } |
140 | |
141 | // Schemas de Newton |
142 | |
143 | math_Gauss LU(TheHessian, CTol/100); |
144 | if ( !LU.IsDone()) { |
145 | Done = Standard_False; |
146 | TheStatus = math_DirectionSearchError; |
147 | return; |
148 | } |
149 | |
150 | LU.Solve(TheGradient, TheStep); |
151 | *suivant = *precedent - TheStep; |
152 | |
153 | |
154 | // Gestion de la convergence |
155 | |
156 | F.Value(*suivant, TheMinimum); |
157 | |
158 | if (IsConverged()) { NbConv++; } |
159 | else { NbConv=0; } |
160 | |
161 | // Controle et corrections. |
162 | |
163 | VItere = TheMinimum; |
164 | TheMinimum = PreviousMinimum; |
165 | Nreduction =0; |
166 | while (VItere > VPrecedent && Nreduction < 10) { |
167 | TheStep *= 0.4; |
168 | *suivant = *precedent - TheStep; |
169 | F.Value(*suivant, VItere); |
170 | Nreduction++; |
171 | } |
172 | |
173 | if (VItere <= VPrecedent) { |
174 | auxiliaire = precedent; |
175 | precedent = suivant; |
176 | suivant = auxiliaire; |
177 | PreviousMinimum = VPrecedent; |
178 | TheMinimum = VItere; |
179 | Ok = (nbiter < Itermax); |
180 | if (!Ok && NbConv < 2) TheStatus = math_TooManyIterations; |
181 | } |
182 | else { |
183 | Ok = Standard_False; |
184 | TheStatus = math_DirectionSearchError; |
185 | } |
186 | } |
187 | TheLocation = *precedent; |
188 | } |
189 | |
190 | //============================================================================ |
191 | Standard_Boolean math_NewtonMinimum::IsConverged() const |
192 | //============================================================================ |
193 | { |
194 | return ( (TheStep.Norm() <= XTol ) || |
195 | ( Abs(TheMinimum-PreviousMinimum) <= XTol*Abs(PreviousMinimum) )); |
196 | } |
197 | |
198 | //============================================================================ |
199 | void math_NewtonMinimum::Dump(Standard_OStream& o) const |
200 | //============================================================================ |
201 | { |
202 | o<< "math_Newton Optimisation: "; |
203 | o << " Done =" << Done << endl; |
204 | o << " Status = " << (Standard_Integer)TheStatus << endl; |
205 | o <<" Location Vector = " << Location() << endl; |
206 | o <<" Minimum value = "<< Minimum()<< endl; |
207 | o <<" Previous value = "<< PreviousMinimum << endl; |
208 | o <<" Number of iterations = " <<NbIterations() << endl; |
209 | o <<" Convexity = " << Convex << endl; |
210 | o <<" Eigen Value = " << MinEigenValue << endl; |
211 | } |
212 | |