1 /**********************************************************************
2   ga_gradient.c
3  **********************************************************************
4 
5   ga_gradient - Gradient methods for comparison and local search.
6   Copyright ©2002-2005, Stewart Adcock <stewart@linux-domain.com>
7   All rights reserved.
8 
9   The latest version of this program should be available at:
10   http://gaul.sourceforge.net/
11 
12   This program is free software; you can redistribute it and/or modify
13   it under the terms of the GNU General Public License as published by
14   the Free Software Foundation; either version 2 of the License, or
15   (at your option) any later version.  Alternatively, if your project
16   is incompatible with the GPL, I will probably agree to requests
17   for permission to use the terms of any other license.
18 
19   This program is distributed in the hope that it will be useful, but
20   WITHOUT ANY WARRANTY WHATSOEVER.
21 
22   A full copy of the GNU General Public License should be in the file
23   "COPYING" provided with this distribution; if not, see:
24   http://www.gnu.org/
25 
26  **********************************************************************
27 
28   Synopsis:     Gradient methods for comparison and local search.
29 
30 		Routines for local search and optimisation using
31 		non-evolutionary methods.  These methods are all
32 		first-order, that is, they require first derivatives.
33 
34 		Note that, these algorithms require that chromosomes
35 		may be reversibly mapped to arrays of double-precision
36 		floating-point array chromsomes.  If this is not
37 		possible then, hmmm, tough luck.
38 
39 		You might want to think carefully about your convergence
40 		criteria.
41 
42   References:
43 
44  **********************************************************************/
45 
46 #include "gaul/ga_gradient.h"
47 
48 /**********************************************************************
49   ga_population_set_gradient_parameters()
50   synopsis:     Sets the gradient-search parameters for a population.
51   parameters:	population *pop		Population to set parameters of.
52 		const GAto_double	Map chromosomal data to array of doubles.
53 		const GAfrom_double	Map array of doubles to chromosomal data.
54 		const int		Number of dimensions for double array (Needn't match dimensions of chromosome.)
55   return:	none
56   last updated: 19 Nov 2002
57  **********************************************************************/
58 
ga_population_set_gradient_parameters(population * pop,const GAto_double to_double,const GAfrom_double from_double,const GAgradient gradient,const int dimensions,const double step_size)59 void ga_population_set_gradient_parameters( population		*pop,
60                                         const GAto_double	to_double,
61                                         const GAfrom_double	from_double,
62                                         const GAgradient	gradient,
63 					const int		dimensions,
64 					const double		step_size)
65   {
66 
67   if ( !pop ) die("Null pointer to population structure passed.");
68 /*
69   if ( !to_double ) die("Null pointer to GAto_double callback passed.");
70   if ( !from_double ) die("Null pointer to GAfrom_double callback passed.");
71 */
72 
73   plog( LOG_VERBOSE, "Population's gradient methods parameters set" );
74 
75   if (pop->gradient_params == NULL)
76     pop->gradient_params = s_malloc(sizeof(ga_gradient_t));
77 
78   pop->gradient_params->to_double = to_double;
79   pop->gradient_params->from_double = from_double;
80   pop->gradient_params->gradient = gradient;
81   pop->gradient_params->step_size = step_size;
82   pop->gradient_params->dimensions = dimensions;
83   pop->gradient_params->alpha = 0.5;	/* Step-size scale-down factor. */
84   pop->gradient_params->beta = 1.2;	/* Step-size scale-up factor. */
85 
86   return;
87   }
88 
89 
90 /**********************************************************************
91   ga_steepestascent()
92   synopsis:	Performs optimisation on the passed entity by using a
93   		steepest ascents method (i.e. steepest descent, except
94 		maximising the fitness function).
95 		The passed entity will have its data overwritten.  The
96 		remainder of the population will be let untouched.
97 		Note that it is safe to pass a NULL initial structure,
98 		in which case a random starting structure wil be
99 		generated, however the final solution will not be
100 		available to the caller in any obvious way.
101   parameters:
102   return:
103   last updated:	18 Feb 2005
104  **********************************************************************/
105 
ga_steepestascent(population * pop,entity * current,const int max_iterations)106 int ga_steepestascent(	population	*pop,
107 			entity		*current,
108 			const int	max_iterations )
109   {
110   int		iteration=0;		/* Current iteration number. */
111   int		i;			/* Index into arrays. */
112   double	*current_d;		/* Current iteration solution array. */
113   double	*current_g;		/* Current iteration gradient array. */
114   entity	*putative;		/* New solution. */
115   double	*putative_d;		/* New solution array. */
116   entity	*tmpentity;		/* Used to swap working solutions. */
117   double	*tmpdoubleptr;		/* Used to swap working solutions. */
118   double	*buffer;		/* Storage for double arrays. */
119   double	step_size;		/* Current step size. */
120   double	grms;			/* Current RMS gradient. */
121   boolean	force_terminate=FALSE;	/* Force optimisation to terminate. */
122 
123 /*
124  * Checks.
125  */
126   if (!pop) die("NULL pointer to population structure passed.");
127   if (!pop->evaluate) die("Population's evaluation callback is undefined.");
128   if (!pop->gradient_params) die("ga_population_set_gradient_params(), or similar, must be used prior to ga_gradient().");
129   if (!pop->gradient_params->to_double) die("Population's genome to double callback is undefined.");
130   if (!pop->gradient_params->from_double) die("Population's genome from double callback is undefined.");
131   if (!pop->gradient_params->gradient) die("Population's first derivatives callback is undefined.");
132 
133 /*
134  * Prepare working entity and double arrays.
135  */
136   buffer = s_malloc(sizeof(double)*pop->gradient_params->dimensions*3);
137 
138   current_d = buffer;
139   putative_d = &(buffer[pop->gradient_params->dimensions]);
140   current_g = &(buffer[pop->gradient_params->dimensions*2]);
141 
142   putative = ga_get_free_entity(pop);
143 
144 /* Do we need to generate a random starting solution? */
145   if (current==NULL)
146     {
147     plog(LOG_VERBOSE, "Will perform gradient search with random starting solution.");
148 
149     current = ga_get_free_entity(pop);
150     ga_entity_seed(pop, current);
151     }
152   else
153     {
154     plog(LOG_VERBOSE, "Will perform gradient search with specified starting solution.");
155     }
156 
157 /*
158  * Get initial fitness and derivatives.
159  */
160   pop->evaluate(pop, current);
161   pop->gradient_params->to_double(pop, current, current_d);
162 
163   grms = pop->gradient_params->gradient(pop, current, current_d, current_g);
164 
165   plog( LOG_VERBOSE,
166         "Prior to the first iteration, the current solution has fitness score of %f and a RMS gradient of %f",
167          current->fitness, grms );
168 
169 /*
170  * Adjust step size based on gradient.
171  * This scales the step size according to the initial gradient so that the
172  * calculation doesn't blow-up completely.
173  */
174 /*  step_size=(pop->gradient_params->dimensions*pop->gradient_params->step_size)/grms;*/
175   step_size=pop->gradient_params->step_size;
176 
177 /*
178  * Do all the iterations:
179  *
180  * Stop when (a) max_iterations reached, or
181  *           (b) "pop->iteration_hook" returns FALSE.
182  * The iteration hook could evaluate the RMS gradient, or the maximum component
183  * of the gradient, or any other termination criteria that may be desirable.
184  */
185   while ( force_terminate==FALSE &&
186           (pop->iteration_hook?pop->iteration_hook(iteration, current):TRUE) &&
187           iteration<max_iterations )
188     {
189     iteration++;
190 
191     for( i=0; i<pop->gradient_params->dimensions; i++ )
192       putative_d[i]=current_d[i]+step_size*current_g[i];
193 
194     pop->gradient_params->from_double(pop, putative, putative_d);
195     pop->evaluate(pop, putative);
196 
197 #if GA_DEBUG>2
198     printf("DEBUG: current_d = %f %f %f %f\n", current_d[0], current_d[1], current_d[2], current_d[3]);
199     printf("DEBUG: current_g = %f %f %f %f grms = %f\n", current_g[0], current_g[1], current_g[2], current_g[3], grms);
200     printf("DEBUG: putative_d = %f %f %f %f fitness = %f\n", putative_d[0], putative_d[1], putative_d[2], putative_d[3], putative->fitness);
201 #endif
202 
203     if ( current->fitness > putative->fitness )
204       {	/* New solution is worse. */
205 
206       do
207         {
208         step_size *= pop->gradient_params->alpha;
209         /*printf("DEBUG: step_size = %e\n", step_size);*/
210 
211         for( i=0; i<pop->gradient_params->dimensions; i++ )
212           putative_d[i]=current_d[i]+step_size*current_g[i];
213 
214         pop->gradient_params->from_double(pop, putative, putative_d);
215         pop->evaluate(pop, putative);
216 
217 #if GA_DEBUG>2
218         printf("DEBUG: putative_d = %f %f %f %f fitness = %f\n", putative_d[0], putative_d[1], putative_d[2], putative_d[3], putative->fitness);
219 #endif
220         } while( current->fitness > putative->fitness && step_size > ApproxZero);
221 
222       if (step_size <= ApproxZero && grms <= ApproxZero) force_terminate=TRUE;
223       }
224     else
225       {	/* New solution is an improvement. */
226       step_size *= pop->gradient_params->beta;
227 #if GA_DEBUG>2
228       printf("DEBUG: step_size = %e\n", step_size);
229 #endif
230       }
231 
232 /* Store improved solution. */
233     tmpentity = current;
234     current = putative;
235     putative = tmpentity;
236 
237     tmpdoubleptr = current_d;
238     current_d = putative_d;
239     putative_d = tmpdoubleptr;
240 
241     grms = pop->gradient_params->gradient(pop, current, current_d, current_g);
242 
243 /*
244  * Use the iteration callback.
245  */
246     plog( LOG_VERBOSE,
247           "After iteration %d, the current solution has fitness score of %f and RMS gradient of %f (step_size = %f)",
248           iteration, current->fitness, grms, step_size );
249 
250     }	/* Iteration loop. */
251 
252 /*
253  * Cleanup.
254  */
255   ga_entity_dereference(pop, putative);
256 
257   s_free(buffer);
258 
259   return iteration;
260   }
261 
262 
263 /**********************************************************************
264   ga_steepestascent_double()
265   synopsis:	Performs optimisation on the passed entity by using a
266   		steepest ascents method (i.e. steepest descent, except
267 		maximising the fitness function).
268 		The passed entity will have its data overwritten.  The
269 		remainder of the population will be let untouched.
270 		Note that it is safe to pass a NULL initial structure,
271 		in which case a random starting structure wil be
272 		generated, however the final solution will not be
273 		available to the caller in any obvious way.
274 
275 		Only double chromosomes may be used in this optimised
276 		version of the algorithm.
277   parameters:
278   return:
279   last updated:	14 Apr 2005
280  **********************************************************************/
281 
ga_steepestascent_double(population * pop,entity * current,const int max_iterations)282 int ga_steepestascent_double(	population	*pop,
283 			entity		*current,
284 			const int	max_iterations )
285   {
286   int		iteration=0;		/* Current iteration number. */
287   int		i;			/* Index into arrays. */
288   double	*current_g;		/* Current iteration gradient array. */
289   entity	*putative;		/* New solution. */
290   entity	*tmpentity;		/* Used to swap working solutions. */
291   double	step_size;		/* Current step size. */
292   double	grms;			/* Current RMS gradient. */
293   boolean	force_terminate=FALSE;	/* Force optimisation to terminate. */
294 
295 /*
296  * Checks.
297  */
298   if (!pop) die("NULL pointer to population structure passed.");
299   if (!pop->evaluate) die("Population's evaluation callback is undefined.");
300   if (!pop->gradient_params) die("ga_population_set_gradient_params(), or similar, must be used prior to ga_gradient().");
301   if (!pop->gradient_params->gradient) die("Population's first derivatives callback is undefined.");
302 
303 /*
304  * Prepare working entity and gradient array.
305  */
306   current_g = s_malloc(sizeof(double)*pop->len_chromosomes);
307 
308   putative = ga_get_free_entity(pop);
309 
310 /* Do we need to generate a random starting solution? */
311   if (current==NULL)
312     {
313     plog(LOG_VERBOSE, "Will perform gradient search with random starting solution.");
314 
315     current = ga_get_free_entity(pop);
316     ga_entity_seed(pop, current);
317     }
318   else
319     {
320     plog(LOG_VERBOSE, "Will perform gradient search with specified starting solution.");
321     }
322 
323 /*
324  * Get initial fitness and derivatives.
325  */
326   pop->evaluate(pop, current);
327 
328   grms = pop->gradient_params->gradient(pop, current, (double *)current->chromosome[0], current_g);
329 
330   plog( LOG_VERBOSE,
331         "Prior to the first iteration, the current solution has fitness score of %f and a RMS gradient of %f",
332          current->fitness, grms );
333 
334 /*
335  * Adjust step size based on gradient.
336  * This scales the step size according to the initial gradient so that the
337  * calculation doesn't blow-up completely.
338  */
339 /*  step_size=(pop->len_chromosomes*pop->gradient_params->step_size)/grms;*/
340   step_size=pop->gradient_params->step_size;
341 
342 /*
343  * Do all the iterations:
344  *
345  * Stop when (a) max_iterations reached, or
346  *           (b) "pop->iteration_hook" returns FALSE.
347  * The iteration hook could evaluate the RMS gradient, or the maximum component
348  * of the gradient, or any other termination criteria that may be desirable.
349  */
350   while ( force_terminate==FALSE &&
351           (pop->iteration_hook?pop->iteration_hook(iteration, current):TRUE) &&
352           iteration<max_iterations )
353     {
354     iteration++;
355 
356     for( i=0; i<pop->len_chromosomes; i++ )
357       ((double *)putative->chromosome[0])[i]=((double *)current->chromosome[0])[i]+step_size*current_g[i];
358 
359     pop->evaluate(pop, putative);
360 
361 #if GA_DEBUG>2
362     printf("DEBUG: current_d = %f %f %f %f\n", current_d[0], current_d[1], current_d[2], current_d[3]);
363     printf("DEBUG: current_g = %f %f %f %f grms = %f\n", current_g[0], current_g[1], current_g[2], current_g[3], grms);
364     printf("DEBUG: putative_d = %f %f %f %f fitness = %f\n", putative_d[0], putative_d[1], putative_d[2], putative_d[3], putative->fitness);
365 #endif
366 
367     if ( current->fitness > putative->fitness )
368       {	/* New solution is worse. */
369 
370       do
371         {
372         step_size *= pop->gradient_params->alpha;
373         /*printf("DEBUG: step_size = %e\n", step_size);*/
374 
375         for( i=0; i<pop->len_chromosomes; i++ )
376           ((double *)putative->chromosome[0])[i]=((double *)current->chromosome[0])[i]+step_size*current_g[i];
377 
378         pop->evaluate(pop, putative);
379         } while( current->fitness > putative->fitness && step_size > ApproxZero);
380 
381       if (step_size <= ApproxZero && grms <= ApproxZero) force_terminate=TRUE;
382       }
383     else
384       {	/* New solution is an improvement. */
385       step_size *= pop->gradient_params->beta;
386 #if GA_DEBUG>2
387       printf("DEBUG: step_size = %e\n", step_size);
388 #endif
389       }
390 
391 /* Store improved solution. */
392     tmpentity = current;
393     current = putative;
394     putative = tmpentity;
395 
396     grms = pop->gradient_params->gradient(pop, current, (double *)current->chromosome[0], current_g);
397 
398 /*
399  * Use the iteration callback.
400  */
401     plog( LOG_VERBOSE,
402           "After iteration %d, the current solution has fitness score of %f and RMS gradient of %f (step_size = %f)",
403           iteration, current->fitness, grms, step_size );
404 
405     }	/* Iteration loop. */
406 
407 /*
408  * Cleanup.
409  */
410   ga_entity_dereference(pop, putative);
411 
412   return iteration;
413   }
414 
415