0029858: Modeling Data - Regression in GeomAPI_ExtremaCurveCurve
[occt.git] / src / math / math_GlobOptMin.hxx
1 // Created on: 2014-01-20
2 // Created by: Alexaner Malyshev
3 // Copyright (c) 2014-2015 OPEN CASCADE SAS
4 //
5 // This file is part of Open CASCADE Technology software library.
6 //
7 // This library is free software; you can redistribute it and/or modify it under
8 // the terms of the GNU Lesser General Public License version 2.1 as published
9 // by the Free Software Foundation, with special exception defined in the file
10 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11 // distribution for complete text of the license and disclaimer of any warranty.
12 //
13 // Alternatively, this file may be used under the terms of Open CASCADE
14 // commercial license or contractual agreement.
15
16 #ifndef _math_GlobOptMin_HeaderFile
17 #define _math_GlobOptMin_HeaderFile
18
19 #include <gp_Pnt.hxx>
20 #include <gp_Pnt2d.hxx>
21 #include <NCollection_CellFilter.hxx>
22 #include <math_MultipleVarFunction.hxx>
23 #include <NCollection_Sequence.hxx>
24 #include <Standard_Type.hxx>
25
26 //! This class represents Evtushenko's algorithm of global optimization based on non-uniform mesh.
27 //! Article: Yu. Evtushenko. Numerical methods for finding global extreme (case of a non-uniform mesh).
28 //! U.S.S.R. Comput. Maths. Math. Phys., Vol. 11, N 6, pp. 38-54.
29 //!
30 //! This method performs search on non-uniform mesh. The search space is a box in R^n space.
31 //! The default behavior is to find all minimums in that box. Computation of maximums is not supported.
32 //!
33 //! The search box can be split into smaller boxes by discontinuity criteria.
34 //! This functionality is covered by SetGlobalParams and SetLocalParams API.
35 //!
36 //! It is possible to set continuity of the local boxes.
37 //! Such option can forcibly change local extrema search.
38 //! In other words if theFunc can be casted to the function with Hessian but, continuity is set to 1
39 //! Gradient based local optimization method will be used, not Hessian based method.
40 //! This functionality is covered by SetContinuity and GetContinuity API.
41 //!
42 //! It is possible to freeze Lipschitz const to avoid internal modifications on it.
43 //! This functionality is covered by SetLipConstState and GetLipConstState API.
44 //!
45 //! It is possible to perform single solution search.
46 //! This functionality is covered by first parameter in Perform method.
47 //!
48 //! It is possible to set / get minimal value of the functional.
49 //! It works well together with single solution search.
50 //! This functionality is covered by SetFunctionalMinimalValue and GetFunctionalMinimalValue API.
51 class math_GlobOptMin
52 {
53 public:
54
55   //! Constructor. Perform method is not called from it.
56   //! @param theFunc - objective functional.
57   //! @param theLowerBorder - lower corner of the search box.
58   //! @param theUpperBorder - upper corner of the search box.
59   //! @param theC - Lipschitz constant.
60   //! @param theDiscretizationTol - parameter space discretization tolerance.
61   //! @param theSameTol - functional value space indifference tolerance.
62   Standard_EXPORT math_GlobOptMin(math_MultipleVarFunction* theFunc,
63                                  const math_Vector& theLowerBorder,
64                                  const math_Vector& theUpperBorder,
65                                  const Standard_Real theC = 9,
66                                  const Standard_Real theDiscretizationTol = 1.0e-2,
67                                  const Standard_Real theSameTol = 1.0e-7);
68
69   //! @param theFunc - objective functional.
70   //! @param theLowerBorder - lower corner of the search box.
71   //! @param theUpperBorder - upper corner of the search box.
72   //! @param theC - Lipschitz constant.
73   //! @param theDiscretizationTol - parameter space discretization tolerance.
74   //! @param theSameTol - functional value space indifference tolerance.
75   Standard_EXPORT void SetGlobalParams(math_MultipleVarFunction* theFunc,
76                                        const math_Vector& theLowerBorder,
77                                        const math_Vector& theUpperBorder,
78                                        const Standard_Real theC = 9,
79                                        const Standard_Real theDiscretizationTol = 1.0e-2,
80                                        const Standard_Real theSameTol = 1.0e-7);
81
82   //! Method to reduce bounding box. Perform will use this box.
83   //! @param theLocalA - lower corner of the local box.
84   //! @param theLocalB - upper corner of the local box.
85   Standard_EXPORT void SetLocalParams(const math_Vector& theLocalA,
86                                       const math_Vector& theLocalB);
87
88   //! Method to set tolerances.
89   //! @param theDiscretizationTol - parameter space discretization tolerance.
90   //! @param theSameTol - functional value space indifference tolerance.
91   Standard_EXPORT void SetTol(const Standard_Real theDiscretizationTol,
92                               const Standard_Real theSameTol);
93
94   //! Method to get tolerances.
95   //! @param theDiscretizationTol - parameter space discretization tolerance.
96   //! @param theSameTol - functional value space indifference tolerance.
97   Standard_EXPORT void GetTol(Standard_Real& theDiscretizationTol,
98                               Standard_Real& theSameTol);
99
100   //! @param isFindSingleSolution - defines whether to find single solution or all solutions.
101   Standard_EXPORT void Perform(const Standard_Boolean isFindSingleSolution = Standard_False);
102
103   //! Return solution theIndex, 1 <= theIndex <= NbExtrema.
104   Standard_EXPORT void Points(const Standard_Integer theIndex, math_Vector& theSol);
105
106   //! Set / Get continuity of local borders splits (0 ~ C0, 1 ~ C1, 2 ~ C2).
107   inline void SetContinuity(const Standard_Integer theCont) { myCont = theCont; }
108   inline Standard_Integer GetContinuity() const {return myCont; }
109
110     //! Set / Get functional minimal value.
111   inline void SetFunctionalMinimalValue(const Standard_Real theMinimalValue)
112   { myFunctionalMinimalValue = theMinimalValue; }
113   inline Standard_Real GetFunctionalMinimalValue() const {return myFunctionalMinimalValue;}
114
115   //! Set / Get Lipchitz constant modification state. 
116   //! True means that the constant is locked and unlocked otherwise.
117   inline void SetLipConstState(const Standard_Boolean theFlag) {myIsConstLocked = theFlag; }
118   inline Standard_Boolean GetLipConstState() const { return myIsConstLocked; }
119
120   //! Return computation state of the algorithm.
121   inline Standard_Boolean isDone() const {return myDone; }
122
123   //! Get best functional value.
124   inline Standard_Real GetF() const {return myF;}
125
126   //! Return count of global extremas.
127   inline Standard_Integer NbExtrema() const {return mySolCount;}
128
129 private:
130
131   //! Class for duplicate fast search. For internal usage only.
132   class NCollection_CellFilter_Inspector
133   {
134   public:
135
136     //! Points and target type
137     typedef math_Vector Point;
138     typedef math_Vector Target;
139
140     NCollection_CellFilter_Inspector(const Standard_Integer theDim,
141                                      const Standard_Real theTol)
142       : myCurrent(1, theDim)
143     {
144       myTol = theTol * theTol;
145       myIsFind = Standard_False;
146       Dimension = theDim;
147     }
148
149     //! Access to coordinate.
150     static Standard_Real Coord (int i, const Point &thePnt)
151     {
152       return thePnt(i + 1);
153     }
154
155     //! Auxiliary method to shift point by each coordinate on given value;
156     //! useful for preparing a points range for Inspect with tolerance.
157     void Shift (const Point& thePnt,
158                 const NCollection_Array1<Standard_Real> &theTol,
159                 Point& theLowPnt,
160                 Point& theUppPnt) const
161     {
162       for(Standard_Integer anIdx = 1; anIdx <= Dimension; anIdx++)
163       {
164         theLowPnt(anIdx) = thePnt(anIdx) - theTol(anIdx - 1);
165         theUppPnt(anIdx) = thePnt(anIdx) + theTol(anIdx - 1);
166       }
167     }
168
169     void ClearFind()
170     {
171       myIsFind = Standard_False;
172     }
173
174     Standard_Boolean isFind()
175     {
176       return myIsFind;
177     }
178
179     //! Set current point to search for coincidence
180     void SetCurrent (const math_Vector& theCurPnt)
181     { 
182       myCurrent = theCurPnt;
183     }
184
185     //! Implementation of inspection method
186     NCollection_CellFilter_Action Inspect (const Target& theObject)
187     {
188       Standard_Real aSqDist = (myCurrent - theObject).Norm2();
189
190       if(aSqDist < myTol)
191       {
192         myIsFind = Standard_True;
193       }
194
195       return CellFilter_Keep;
196     }
197
198   private:
199     Standard_Real myTol;
200     math_Vector myCurrent;
201     Standard_Boolean myIsFind;
202     Standard_Integer Dimension;
203   };
204
205
206   // Compute cell size.
207   void initCellSize();
208
209   // Compute initial solution
210   void ComputeInitSol();
211
212   math_GlobOptMin & operator = (const math_GlobOptMin & theOther);
213
214   Standard_Boolean computeLocalExtremum(const math_Vector& thePnt, Standard_Real& theVal, math_Vector& theOutPnt);
215
216   void computeGlobalExtremum(Standard_Integer theIndex);
217
218   //! Check possibility to stop computations.
219   //! Find single solution + in neighborhood of best possible solution.
220   Standard_Boolean CheckFunctionalStopCriteria();
221
222   //! Computes starting value / approximation:
223   //! myF - initial best value.
224   //! myY - initial best point.
225   //! myC - approximation of Lipschitz constant.
226   //! to improve convergence speed.
227   void computeInitialValues();
228
229   //! Check that myA <= thePnt <= myB
230   Standard_Boolean isInside(const math_Vector& thePnt);
231
232   //! Check presence of thePnt in GlobOpt sequence.
233   Standard_Boolean isStored(const math_Vector &thePnt);
234
235   //! Check and add candidate point into solutions.
236   //! @param thePnt   Candidate point.
237   //! @param theValue Function value in the candidate point.
238   void checkAddCandidate(const math_Vector&  thePnt,
239                          const Standard_Real theValue);
240
241
242   // Input.
243   math_MultipleVarFunction* myFunc;
244   Standard_Integer myN;
245   math_Vector myA; // Left border on current C2 interval.
246   math_Vector myB; // Right border on current C2 interval.
247   math_Vector myGlobA; // Global left border.
248   math_Vector myGlobB; // Global right border.
249   Standard_Real myTol; // Discretization tolerance, default 1.0e-2.
250   Standard_Real mySameTol; // points with ||p1 - p2|| < mySameTol is equal,
251                            // function values |val1 - val2| * 0.01 < mySameTol is equal,
252                            // default value is 1.0e-7.
253   Standard_Real myC; //Lipchitz constant, default 9
254   Standard_Real myInitC; // Lipchitz constant initial value.
255   Standard_Boolean myIsFindSingleSolution; // Default value is false.
256   Standard_Real myFunctionalMinimalValue; // Default value is -Precision::Infinite
257   Standard_Boolean myIsConstLocked;  // Is constant locked for modifications.
258
259   // Output.
260   Standard_Boolean myDone;
261   NCollection_Sequence<Standard_Real> myY;// Current solutions.
262   Standard_Integer mySolCount; // Count of solutions.
263
264   // Algorithm data.
265   Standard_Real myZ;
266   Standard_Real myE1; // Border coefficient.
267   Standard_Real myE2; // Minimum step size.
268   Standard_Real myE3; // Local extrema starting parameter.
269
270   math_Vector myX; // Current modified solution.
271   math_Vector myTmp; // Current modified solution.
272   math_Vector myV; // Steps array.
273   math_Vector myMaxV; // Max Steps array.
274   Standard_Real myLastStep; // Last step.
275
276   NCollection_Array1<Standard_Real> myCellSize;
277   Standard_Integer myMinCellFilterSol;
278   Standard_Boolean isFirstCellFilterInvoke;
279   NCollection_CellFilter<NCollection_CellFilter_Inspector> myFilter;
280
281   // Continuity of local borders.
282   Standard_Integer myCont;
283
284   Standard_Real myF; // Current value of Global optimum.
285 };
286
287 #endif