Annotation of /trunk/gadget/optinfo.h

Parent Directory | Revision Log

---

Revision 12 - (view) (download)

1 :	agomez	1	#ifndef optinfo_h
2 :			#define optinfo_h
3 :
4 :			#include "maininfo.h"
5 :			#include "doublematrix.h"
6 :			#include "doublevector.h"
7 :			#include "intvector.h"
8 :	ulcessvp	11	#include "seq_optimize_template.h"
9 :	agomez	1
10 :			enum OptType { OPTHOOKE = 1, OPTSIMANN, OPTBFGS };
11 :
12 :			/**
13 :			* \class OptInfo
14 :			* \brief This is the base class used to perform the optimisation calculation for the model
15 :			*
16 :			* \note This will always be overridden by the derived classes that actually perform the optimisation calculation
17 :			*/
18 :			class OptInfo {
19 :			public:
20 :			/**
21 :			* \brief This is the default OptInfo constructor
22 :			*/
23 :			OptInfo() { converge = 0; iters = 0; score = 0.0; };
24 :			/**
25 :			* \brief This is the default OptInfo destructor
26 :			*/
27 :			~OptInfo() {};
28 :			/**
29 :			* \brief This is the function used to read in the optimisation parameters
30 :			* \param infile is the CommentStream to read the optimisation parameters from
31 :			* \param text is a text string used to compare parameter names
32 :			*/
33 :			virtual void read(CommentStream& infile, char* text) {};
34 :			/**
35 :			* \brief This function will print information from the optimisation algorithm
36 :			* \param outfile is the ofstream that the optimisation information gets sent to
37 :			* \param prec is the precision to use in the output file
38 :			*/
39 :			virtual void Print(ofstream& outfile, int prec) {};
40 :			/**
41 :			* \brief This is the function used to call the optimisation algorithms
42 :			*/
43 :			virtual void OptimiseLikelihood() {};
44 :	ulcessvp	12	/**
45 :			* \brief This is the function used to call the optimisation algorithms parallelized with OpenMP of the reproducible version
46 :			*/
47 :	ulcessvp	11	virtual void OptimiseLikelihoodOMP() {};
48 :	ulcessvp	12	/**
49 :			* \brief This function set the seeds used in SA
50 :			* \param val array of unsigned int with the seeds
51 :			*/
52 :	ulcessvp	11	void setSeed(unsigned* val) {seed = val[0]; seedM = val[1]; seedP = val[2];};
53 :	agomez	1	/**
54 :			* \brief This will return the type of optimisation class
55 :			* \return type
56 :			*/
57 :			OptType getType() const { return type; };
58 :			protected:
59 :			/**
60 :			* \brief This is the flag used to denote whether the optimisation converged or not
61 :			*/
62 :			int converge;
63 :			/**
64 :			* \brief This is the number of iterations that took place during the optimisation
65 :			*/
66 :			int iters;
67 :			/**
68 :			* \brief This is the value of the best likelihood score from the optimisation
69 :			*/
70 :			double score;
71 :			/**
72 :			* \brief This denotes what type of optimisation class has been created
73 :			*/
74 :			OptType type;
75 :	ulcessvp	12	/**
76 :			* \brief This is the seed used for the calculation of the new value of the parameters
77 :			*/
78 :	ulcessvp	11	unsigned seed;
79 :	ulcessvp	12	/**
80 :			* \brief This is the seed used for the acceptance of the metropolis
81 :			*/
82 :	ulcessvp	11	unsigned seedM;
83 :	ulcessvp	12	/**
84 :			* \brief This is the seed used to change the order of the parameters
85 :			*/
86 :	ulcessvp	11	unsigned seedP;
87 :	agomez	1	};
88 :
89 :			/**
90 :			* \class OptInfoHooke
91 :			* \brief This is the class used for the Hooke & Jeeves optimisation
92 :			*
93 :			* The Hooke & Jeeves optimisation is the default optimisation, and is a simple and fast optimising method, but somewhat unreliable, which is often described as a "hill climbing" technique. From the initial starting point the algorithm takes a step in various directions, and conducts a new model run. If the new likelihood score is better than the old one then the algorithm uses the new point as it's best guess. If it is worse then the algorithm retains the old point. The search proceeds in series of these steps, each step slightly smaller than the previous one. When the algorithm finds a point which it cannot improve on with a small step in any direction then it accepts this point as being the "solution", and exits. It is recommended that you re-run the optimisation, using the final point of one run as the start of the next.
94 :			*
95 :			* The Hooke & Jeeves algorithm used in Gadget is derived from that originally presented by R. Hooke and T. A. Jeeves, ''Direct Search Solution of Numerical and Statistical Problems'' in the April 1961 (Vol. 8, pp. 212-229) issue of the Journal of the ACM, with improvements presented by Arthur F Kaupe Jr., ''Algorithm 178: Direct Search'' in the June 1963 (Vol 6, pp.313-314) issue of the Communications of the ACM.
96 :
97 :			*/
98 :			class OptInfoHooke : public OptInfo {
99 :			public:
100 :			/**
101 :			* \brief This is the default OptInfoHooke constructor
102 :			*/
103 :			OptInfoHooke();
104 :			/**
105 :			* \brief This is the default OptInfoHooke destructor
106 :			*/
107 :			virtual ~OptInfoHooke() {};
108 :			/**
109 :			* \brief This is the function used to read in the Hooke & Jeeves parameters
110 :			* \param infile is the CommentStream to read the optimisation parameters from
111 :			* \param text is a text string used to compare parameter names
112 :			*/
113 :			virtual void read(CommentStream& infile, char* text);
114 :			/**
115 :			* \brief This function will print information from the optimisation algorithm
116 :			* \param outfile is the ofstream that the optimisation information gets sent to
117 :			* \param prec is the precision to use in the output file
118 :			*/
119 :			virtual void Print(ofstream& outfile, int prec);
120 :			/**
121 :			* \brief This is the function that will calculate the likelihood score using the Hooke & Jeeves optimiser
122 :			*/
123 :			virtual void OptimiseLikelihood();
124 :	ulcessvp	11	#ifdef GADGET_OPENMP
125 :	ulcessvp	12	/**
126 :			* \brief This is the function that will calculate the likelihood score using the Hooke & Jeeves optimiser parallelized with the reproducible version implemented OpenMP
127 :			*/
128 :	ulcessvp	11	virtual void OptimiseLikelihoodOMP();
129 :			#endif
130 :	agomez	1	private:
131 :			/**
132 :			* \brief This function will calculate the best point that can be found close to the current point
133 :			* \param delta is the DoubleVector of the steps to take when looking for the best point
134 :			* \param point is the DoubleVector that will contain the parameters corresponding to the best function value found from the search
135 :			* \param prevbest is the current best point value
136 :			* \param param is the IntVector containing the order that the parameters should be searched in
137 :			* \return the best function value found from the search
138 :			*/
139 :			double bestNearby(DoubleVector& delta, DoubleVector& point, double prevbest, IntVector& param);
140 :	ulcessvp	12	/**
141 :			* \brief This function implemented the reproducible version with OpenMP will calculate the best point that can be found close to the current point
142 :			* \param delta is the DoubleVector of the steps to take when looking for the best point
143 :			* \param point is the DoubleVector that will contain the parameters corresponding to the best function value found from the search
144 :			* \param prevbest is the current best point value
145 :			* \param param is the IntVector containing the order that the parameters should be searched in
146 :			* \return the best function value found from the search
147 :			*/
148 :	ulcessvp	11	double bestNearbyOMP(DoubleVector& delta, DoubleVector& point, double prevbest, IntVector& param);
149 :	ulcessvp	12	/**
150 :			* \brief This function implemented the speculative version with OpenMP will calculate the best point that can be found close to the current point
151 :			* \param delta is the DoubleVector of the steps to take when looking for the best point
152 :			* \param point is the DoubleVector that will contain the parameters corresponding to the best function value found from the search
153 :			* \param prevbest is the current best point value
154 :			* \param param is the IntVector containing the order that the parameters should be searched in
155 :			* \return the best function value found from the search
156 :			*/
157 :	ulcessvp	11	double bestNearbyOMP2(DoubleVector& delta, DoubleVector& point, double prevbest, IntVector& param);
158 :	agomez	1	/**
159 :			* \brief This is the maximum number of iterations for the Hooke & Jeeves optimisation
160 :			*/
161 :			int hookeiter;
162 :			/**
163 :			* \brief This is the reduction factor for the step length
164 :			*/
165 :			double rho;
166 :			/**
167 :			* \brief This is the initial step length
168 :			*/
169 :			double lambda;
170 :			/**
171 :			* \brief This is the minimum step length, use as the halt criteria for the optimisation process
172 :			*/
173 :			double hookeeps;
174 :			/**
175 :			* \brief This is the limit when checking if a parameter is stuck on the bound
176 :			*/
177 :			double bndcheck;
178 :			};
179 :
180 :			/**
181 :			* \class OptInfoSimann
182 :			* \brief This is the class used for the Simualted Annealing optimisation
183 :			*
184 :			* Simulated Annealing is a global optimisation method that distinguishes different local optima. Starting from an initial point, the algorithm takes a step and the function is evaluated. When minimizing a function, any downhill step is accepted and the process repeats from this new point. An uphill step may be accepted (thus, it can escape from local optima). This uphill decision is made by the Metropolis criteria. It uses a parameter known as "temperature" and the size of the uphill step in a probabilistic manner, and varying the temperature will affect the number of the uphill moves that are accepted. As the optimisation process proceeds, the length of the steps decline and the algorithm closes in on the global optimum.
185 :			*
186 :			* The Simulated Annealing algorithm used in Gadget is derived from that presented by Corana et al, ''Minimising Multimodal Functions of Continuous Variables with the 'Simulated Annealing' Algorithm'' in the September 1987 (Vol. 13, pp. 262-280) issue of the ACM Transactions on Mathematical Software and Goffe et al, ''Global Optimisation of Statistical Functions with Simulated Annealing'' in the January/February 1994 (Vol. 60, pp. 65-100) issue of the Journal of Econometrics.
187 :			*/
188 :			class OptInfoSimann : public OptInfo {
189 :			public:
190 :			/**
191 :			* \brief This is the default OptInfoSimann constructor
192 :			*/
193 :			OptInfoSimann();
194 :			/**
195 :			* \brief This is the default OptInfoSimann destructor
196 :			*/
197 :			virtual ~OptInfoSimann() {};
198 :			/**
199 :			* \brief This is the function used to read in the Simulated Annealing parameters
200 :			* \param infile is the CommentStream to read the optimisation parameters from
201 :			* \param text is a text string used to compare parameter names
202 :			*/
203 :			virtual void read(CommentStream& infile, char* text);
204 :			/**
205 :			* \brief This function will print information from the optimisation algorithm
206 :			* \param outfile is the ofstream that the optimisation information gets sent to
207 :			* \param prec is the precision to use in the output file
208 :			*/
209 :			virtual void Print(ofstream& outfile, int prec);
210 :			/**
211 :			* \brief This is the function that will calculate the likelihood score using the Simulated Annealing optimiser
212 :			*/
213 :			virtual void OptimiseLikelihood();
214 :	ulcessvp	11	#ifdef GADGET_OPENMP
215 :	ulcessvp	12	/**
216 :			* \brief This is the function that will calculate the likelihood score using the Simulated Annealing optimiser parallelized with the reproducible version implemented OpenMP
217 :			*/
218 :	ulcessvp	11	virtual void OptimiseLikelihoodOMP();
219 :			#endif
220 :	ulcessvp	12	/**
221 :			* \brief This function calculate a new valor for the parameter l
222 :			* \param nvars the number of variables to be optimised
223 :			* \param l the parameter to change
224 :			* \param param IntVector with the order of the parameters
225 :			* \param trialx DoubleVector that storage the values of the parameters to evaluate during this iteration
226 :			* \param x DoubleVector that storage the best values of the parameters
227 :			* \param lowerb DoubleVector with the lower bounds of the variables to be optimised
228 :			* \param upperb DoubleVector with the upper bounds of the variables to be optimised
229 :			* \param vm DoubleVector with the value for the maximum step length for each parameter
230 :			*/
231 :	ulcessvp	11	virtual void newValue(int nvars, int l, IntVector& param, DoubleVector& trialx,
232 :			DoubleVector& x, DoubleVector& lowerb, DoubleVector& upperb, DoubleVector& vm);
233 :	agomez	1	private:
234 :			/**
235 :			* \brief This is the temperature reduction factor
236 :			*/
237 :			double rt;
238 :			/**
239 :			* \brief This is the halt criteria for the Simulated Annealing algorithm
240 :			*/
241 :			double simanneps;
242 :			/**
243 :			* \brief This is the number of loops before the step length is adjusted
244 :			*/
245 :			int ns;
246 :			/**
247 :			* \brief This is the number of loops before the temperature is adjusted
248 :			*/
249 :			int nt;
250 :			/**
251 :			* \brief This is the "temperature" used for the Simulated Annealing algorithm
252 :			*/
253 :			double t;
254 :			/**
255 :			* \brief This is the factor used to adjust the step length
256 :			*/
257 :			double cs;
258 :			/**
259 :			* \brief This is the initial value for the maximum step length
260 :			*/
261 :			double vminit;
262 :			/**
263 :			* \brief This is the maximum number of function evaluations for the Simulated Annealing optimiation
264 :			*/
265 :			int simanniter;
266 :			/**
267 :			* \brief This is the upper bound when adjusting the step length
268 :			*/
269 :			double uratio;
270 :			/**
271 :			* \brief This is the lower bound when adjusting the step length
272 :			*/
273 :			double lratio;
274 :			/**
275 :			* \brief This is the number of temperature loops to check when testing for convergence
276 :			*/
277 :			int tempcheck;
278 :			/**
279 :			* \brief This is the flag to denote whether the parameters should be scaled or not (default 0, not scale)
280 :			*/
281 :			int scale;
282 :			};
283 :
284 :			/**
285 :			* \class OptInfoBFGS
286 :			* \brief This is the class used for the BFGS optimisation
287 :			*
288 :			* BFGS is a quasi-Newton global optimisation method that uses information about the gradient of the function at the current point to calculate the best direction to look in to find a better point. Using this information, the BFGS algorithm can iteratively calculate a better approximation to the inverse Hessian matrix, which will lead to a better approximation of the minimum value. From an initial starting point, the gradient of the function is calculated and then the algorithm uses this information to calculate the best direction to perform a linesearch for a point that is ''sufficiently better''. The linesearch that is used in Gadget to look for a better point in this direction is the ''Armijo'' linesearch. The algorithm will then adjust the current estimate of the inverse Hessian matrix, and restart from this new point. If a better point cannot be found, then the inverse Hessian matrix is reset and the algorithm restarts from the last accepted point.
289 :			*
290 :			* The BFGS algorithm used in Gadget is derived from that presented by Dimitri P Bertsekas, ''Nonlinear Programming'' (2nd edition, pp22-61) published by Athena Scientific.
291 :			*/
292 :			class OptInfoBFGS : public OptInfo {
293 :			public:
294 :			/**
295 :			* \brief This is the default OptInfoBFGS constructor
296 :			*/
297 :			OptInfoBFGS();
298 :			/**
299 :			* \brief This is the default OptInfoBFGS destructor
300 :			*/
301 :			~OptInfoBFGS() {};
302 :			/**
303 :			* \brief This is the function used to read in the BFGS parameters
304 :			* \param infile is the CommentStream to read the optimisation parameters from
305 :			* \param text is a text string used to compare parameter names
306 :			*/
307 :			virtual void read(CommentStream& infile, char* text);
308 :			/**
309 :			* \brief This function will print information from the optimisation algorithm
310 :			* \param outfile is the ofstream that the optimisation information gets sent to
311 :			* \param prec is the precision to use in the output file
312 :			*/
313 :			virtual void Print(ofstream& outfile, int prec);
314 :			/**
315 :			* \brief This is the function that will calculate the likelihood score using the BFGS optimiser
316 :			*/
317 :			virtual void OptimiseLikelihood();
318 :	ulcessvp	11	#ifdef GADGET_OPENMP
319 :	ulcessvp	12	/**
320 :			* \brief This function call the sequential function. BFGS isn't implemented with OpenMP
321 :			*/
322 :	ulcessvp	11	virtual void OptimiseLikelihoodOMP();
323 :			#endif
324 :	agomez	1	private:
325 :			/**
326 :			* \brief This function will numerically calculate the gradient of the function at the current point
327 :			* \param point is the DoubleVector that contains the parameters corresponding to the current function value
328 :			* \param pointvalue is the current function value
329 :			* \param newgrad is the DoubleVector that will contain the gradient vector for the current point
330 :			*/
331 :			void gradient(DoubleVector& point, double pointvalue, DoubleVector& newgrad);
332 :			/**
333 :			* \brief This function will calculate the smallest eigenvalue of the inverse Hessian matrix
334 :			* \param M is the DoubleMatrix containing the inverse Hessian matrix
335 :			* \return the smallest eigen value of the matrix
336 :			*/
337 :			double getSmallestEigenValue(DoubleMatrix M);
338 :			/**
339 :			* \brief This is the maximum number of function evaluations for the BFGS optimiation
340 :			*/
341 :			int bfgsiter;
342 :			/**
343 :			* \brief This is the halt criteria for the BFGS algorithm
344 :			*/
345 :			double bfgseps;
346 :			/**
347 :			* \brief This is the adjustment factor in the Armijo linesearch
348 :			*/
349 :			double beta;
350 :			/**
351 :			* \brief This is the halt criteria for the Armijo linesearch
352 :			*/
353 :			double sigma;
354 :			/**
355 :			* \brief This is the initial step size for the Armijo linesearch
356 :			*/
357 :			double step;
358 :			/**
359 :			* \brief This is the accuracy term used when calculating the gradient
360 :			*/
361 :			double gradacc;
362 :			/**
363 :			* \brief This is the factor used to adjust the gradient accuracy term
364 :			*/
365 :			double gradstep;
366 :			/**
367 :			* \brief This is the halt criteria for the gradient accuracy term
368 :			*/
369 :			double gradeps;
370 :			};
371 :
372 :			#endif

---

root@forge.cesga.es	ViewVC Help
Powered by ViewVC 1.0.0