1 // Created on: 1996-05-03
2 // Created by: Philippe MANGIN
3 // Copyright (c) 1996-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
6 // This file is part of Open CASCADE Technology software library.
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
18 #define No_Standard_RangeError
19 #define No_Standard_OutOfRange
20 #define No_Standard_DimensionError
23 #include <math_NewtonMinimum.ixx>
25 #include <math_Gauss.hxx>
26 #include <math_Jacobi.hxx>
28 //============================================================================
29 math_NewtonMinimum::math_NewtonMinimum(math_MultipleVarFunctionWithHessian& F,
30 const math_Vector& StartingPoint,
31 const Standard_Real Tolerance,
32 const Standard_Integer NbIterations,
33 const Standard_Real Convexity,
34 const Standard_Boolean WithSingularity)
35 //============================================================================
36 : TheLocation(1, F.NbVariables()),
37 TheGradient(1, F.NbVariables()),
38 TheStep(1, F.NbVariables(), 10*Tolerance),
39 TheHessian(1, F.NbVariables(), 1, F.NbVariables() )
43 Itermax = NbIterations;
44 NoConvexTreatement = WithSingularity;
45 Convex = Standard_True;
46 // Done = Standard_True;
47 // TheStatus = math_OK;
48 Perform ( F, StartingPoint);
51 //============================================================================
52 math_NewtonMinimum::math_NewtonMinimum(math_MultipleVarFunctionWithHessian& F,
53 const Standard_Real Tolerance,
54 const Standard_Integer NbIterations,
55 const Standard_Real Convexity,
56 const Standard_Boolean WithSingularity)
57 //============================================================================
58 : TheLocation(1, F.NbVariables()),
59 TheGradient(1, F.NbVariables()),
60 TheStep(1, F.NbVariables(), 10*Tolerance),
61 TheHessian(1, F.NbVariables(), 1, F.NbVariables() )
65 Itermax = NbIterations;
66 NoConvexTreatement = WithSingularity;
67 Convex = Standard_True;
68 Done = Standard_False;
69 TheStatus = math_NotBracketed;
71 //============================================================================
72 void math_NewtonMinimum::Delete()
74 //============================================================================
75 void math_NewtonMinimum::Perform(math_MultipleVarFunctionWithHessian& F,
76 const math_Vector& StartingPoint)
77 //============================================================================
79 math_Vector Point1 (1, F.NbVariables());
80 Point1 = StartingPoint;
81 math_Vector Point2(1, F.NbVariables());
82 math_Vector* precedent = &Point1;
83 math_Vector* suivant = &Point2;
84 math_Vector* auxiliaire = precedent;
86 Standard_Boolean Ok = Standard_True;
87 Standard_Integer NbConv = 0, ii, Nreduction;
88 Standard_Real VPrecedent, VItere;
94 while ( Ok && (NbConv < 2) ) {
99 Ok = F.Values(*precedent, VPrecedent, TheGradient, TheHessian);
101 Done = Standard_False;
102 TheStatus = math_FunctionError;
106 PreviousMinimum = VPrecedent;
107 TheMinimum = VPrecedent;
110 // Traitement de la non convexite
112 math_Jacobi CalculVP(TheHessian);
113 if ( !CalculVP.IsDone() ) {
114 Done = Standard_False;
115 TheStatus = math_FunctionError;
121 MinEigenValue = CalculVP.Values() ( CalculVP.Values().Min());
122 if ( MinEigenValue < CTol) {
123 Convex = Standard_False;
124 if (NoConvexTreatement) {
125 Standard_Real Delta = CTol+0.1*Abs(MinEigenValue) -MinEigenValue ;
126 for (ii=1; ii<=TheGradient.Length(); ii++) {
127 TheHessian(ii, ii) += Delta;
131 TheStatus = math_FunctionError;
138 math_Gauss LU(TheHessian, CTol/100);
140 Done = Standard_False;
141 TheStatus = math_DirectionSearchError;
145 LU.Solve(TheGradient, TheStep);
146 Standard_Boolean hasProblem = Standard_False;
149 *suivant = *precedent - TheStep;
151 // Gestion de la convergence
152 hasProblem = !(F.Value(*suivant, TheMinimum));
158 } while (hasProblem);
160 if (IsConverged()) { NbConv++; }
163 // Controle et corrections.
166 TheMinimum = PreviousMinimum;
168 while (VItere > VPrecedent && Nreduction < 10) {
170 *suivant = *precedent - TheStep;
171 F.Value(*suivant, VItere);
175 if (VItere <= VPrecedent) {
176 auxiliaire = precedent;
178 suivant = auxiliaire;
179 PreviousMinimum = VPrecedent;
181 Ok = (nbiter < Itermax);
182 if (!Ok && NbConv < 2) TheStatus = math_TooManyIterations;
186 TheStatus = math_DirectionSearchError;
189 TheLocation = *precedent;
192 //============================================================================
193 Standard_Boolean math_NewtonMinimum::IsConverged() const
194 //============================================================================
196 return ( (TheStep.Norm() <= XTol ) ||
197 ( Abs(TheMinimum-PreviousMinimum) <= XTol*Abs(PreviousMinimum) ));
200 //============================================================================
201 void math_NewtonMinimum::Dump(Standard_OStream& o) const
202 //============================================================================
204 o<< "math_Newton Optimisation: ";
205 o << " Done =" << Done << endl;
206 o << " Status = " << (Standard_Integer)TheStatus << endl;
207 o <<" Location Vector = " << Location() << endl;
208 o <<" Minimum value = "<< Minimum()<< endl;
209 o <<" Previous value = "<< PreviousMinimum << endl;
210 o <<" Number of iterations = " <<NbIterations() << endl;
211 o <<" Convexity = " << Convex << endl;
212 o <<" Eigen Value = " << MinEigenValue << endl;