1 /* Copyright (c) 2008-2013 Carlos Henrique da Silva Santos
3 * Permission is hereby granted, free of charge, to any person obtaining
4 * a copy of this software and associated documentation files (the
5 * "Software"), to deal in the Software without restriction, including
6 * without limitation the rights to use, copy, modify, merge, publish,
7 * distribute, sublicense, and/or sell copies of the Software, and to
8 * permit persons to whom the Software is furnished to do so, subject to
9 * the following conditions:
11 * The above copyright notice and this permission notice shall be
12 * included in all copies or substantial portions of the Software.
14 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
15 * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
16 * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
17 * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
18 * LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
19 * OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
20 * WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
27 /****************************************************************************/
28 /* Cauchy random number distribution */
29 static double randcauchy(const double params[7]) {
30 /* double min, double max, double mi, double t, double band */
31 double na_unif, cauchy_mit, limit_inf, limit_sup;
33 double min = params[1], max = params[2], mi = params[3],
34 t = params[4], band = params[5];
35 limit_inf = mi - (band*0.5);
36 limit_sup = mi + (band*0.5);
38 na_unif = nlopt_urand(0,1); /* ran2(0,1); */
39 cauchy_mit = t*tan((na_unif-0.5)*3.14159265358979323846) + mi;
40 } while ( (cauchy_mit<limit_inf) || (cauchy_mit>limit_sup) );
43 cauchy_mit = -cauchy_mit;
45 cauchy_mit = cauchy_mit + (band*0.5);
46 valor = cauchy_mit/band;
47 valor = min+(max-min)*valor;
51 /****************************************************************************/
53 /* main Individual representation type */
54 typedef struct IndividualStructure {
59 static int CompareIndividuals(void *unused, const void *a_, const void *b_) {
61 const Individual *a = (const Individual *) a_;
62 const Individual *b = (const Individual *) b_;
63 return a->fitness < b->fitness ? -1 : (a->fitness > b->fitness ? +1 : 0);
66 nlopt_result chevolutionarystrategy(
67 unsigned nparameters, /* Number of input parameters */
68 nlopt_func f, /* Recursive Objective Funtion Call */
69 void * data_f, /* Data to Objective Function */
70 const double* lb, /* Lower bound values */
71 const double* ub, /* Upper bound values */
72 double* x, /*in: initial guess, out: minimizer */
74 nlopt_stopping* stop, /* nlopt stop condition */
75 unsigned np, /* Number of Parents */
76 unsigned no) { /* Number of Offsprings */
78 /* variables from nlopt */
79 nlopt_result ret = NLOPT_SUCCESS;
83 unsigned crosspoint; /* crossover parameteres */
84 int contmutation, totalmutation; /* mutation parameters */
85 int idoffmutation, paramoffmutation; /* mutation parameters */
86 Individual * esparents; /* Parents population */
87 Individual * esoffsprings; /* Offsprings population */
88 Individual * estotal;/* copy containing Parents and Offsprings pops */
89 /* It is interesting to maintain the parents and offsprings
90 * populations stablished and sorted; when the final iterations
91 * is achieved, they are ranked and updated. */
93 /*********************************
94 * controling the population size
95 *********************************/
98 if ((np < 1)||(no<1)) {
99 nlopt_stop_msg(stop, "populations %d, %d are too small", np, no);
100 return NLOPT_INVALID_ARGS;
102 esparents = (Individual*) malloc(sizeof(Individual) * np);
103 esoffsprings = (Individual*) malloc(sizeof(Individual) * no);
104 estotal = (Individual*) malloc(sizeof(Individual) * (np+no));
105 if ((!esparents)||(!esoffsprings)||(!estotal)) {
106 free(esparents); free(esoffsprings); free(estotal);
107 return NLOPT_OUT_OF_MEMORY;
109 for (id=0; id < np; id++) esparents[id].parameters = NULL;
110 for (id=0; id < no; id++) esoffsprings[id].parameters = NULL;
111 /* From here the population is initialized */
112 /* we don't handle unbounded search regions;
113 this check is unnecessary since it is performed in nlopt_optimize.
114 for (j = 0; j < nparameters; ++j)
115 if (nlopt_isinf(low[j]) || nlopt_isinf(up[j]))
116 return NLOPT_INVALID_ARGS;
118 /* main vector of parameters to randcauchy */
119 vetor[0] = 4; /* ignored */
124 vetor[7] = 0; /* ignored */
125 /**************************************
126 * Initializing parents population
127 **************************************/
128 for (id=0; id < np; id++) {
129 esparents[id].parameters =
130 (double*) malloc(sizeof(double) * nparameters);
131 if (!esparents[id].parameters) {
132 ret = NLOPT_OUT_OF_MEMORY;
135 for (item=0; item<nparameters; item++) {
139 esparents[id].parameters[item] = randcauchy(vetor);
142 memcpy(esparents[0].parameters, x, nparameters * sizeof(double));
144 /**************************************
145 * Initializing offsprings population
146 **************************************/
147 for (id=0; id < no; id++) {
148 esoffsprings[id].parameters =
149 (double*) malloc(sizeof(double) * nparameters);
150 if (!esoffsprings[id].parameters) {
151 ret = NLOPT_OUT_OF_MEMORY;
154 for (item=0; item<nparameters; item++) {
158 esoffsprings[id].parameters[item] = randcauchy(vetor);
161 /**************************************
162 * Parents fitness evaluation
163 **************************************/
164 for (id=0; id < np; id++) {
165 esparents[id].fitness =
166 f(nparameters, esparents[id].parameters, NULL, data_f);
167 estotal[id].fitness = esparents[id].fitness;
169 if (*minf > esparents[id].fitness) {
170 *minf = esparents[id].fitness;
171 memcpy(x, esparents[id].parameters,
172 nparameters * sizeof(double));
174 if (nlopt_stop_forced(stop)) ret = NLOPT_FORCED_STOP;
175 else if (*minf < stop->minf_max) ret = NLOPT_MINF_MAX_REACHED;
176 else if (nlopt_stop_evals(stop)) ret = NLOPT_MAXEVAL_REACHED;
177 else if (nlopt_stop_time(stop)) ret = NLOPT_MAXTIME_REACHED;
178 if (ret != NLOPT_SUCCESS) goto done;
180 /**************************************
181 * Main Loop - Generations
182 **************************************/
184 /**************************************
186 **************************************/
187 for (id=0; id < no; id++)
189 parent1 = nlopt_iurand((int) np);
190 parent2 = nlopt_iurand((int) np);
191 crosspoint = (unsigned) nlopt_iurand((int) nparameters);
192 for (item=0; item < crosspoint; item++)
193 esoffsprings[id].parameters[item]
194 = esparents[parent1].parameters[item];
195 for (item=crosspoint; item < nparameters; item++)
196 esoffsprings[id].parameters[item]
197 = esparents[parent2].parameters[item];
199 /**************************************
201 **************************************/
202 totalmutation = (int) ((no * nparameters) / 10);
203 if (totalmutation < 1) totalmutation = 1;
204 for (contmutation=0; contmutation < totalmutation;
206 idoffmutation = nlopt_iurand((int) no);
207 paramoffmutation = nlopt_iurand((int) nparameters);
208 vetor[1] = lb[paramoffmutation];
209 vetor[2] = ub[paramoffmutation];
210 vetor[7] = vetor[7]+contmutation;
211 esoffsprings[idoffmutation].parameters[paramoffmutation]
214 /**************************************
215 * Offsprings fitness evaluation
216 **************************************/
217 for (id=0; id < no; id++){
218 /*esoffsprings[id].fitness = (double)fitness(esoffsprings[id].parameters, nparameters,fittype);*/
219 esoffsprings[id].fitness = f(nparameters, esoffsprings[id].parameters, NULL, data_f);
220 estotal[id+np].fitness = esoffsprings[id].fitness;
222 if (*minf > esoffsprings[id].fitness) {
223 *minf = esoffsprings[id].fitness;
224 memcpy(x, esoffsprings[id].parameters,
225 nparameters * sizeof(double));
227 if (nlopt_stop_forced(stop)) ret = NLOPT_FORCED_STOP;
228 else if (*minf < stop->minf_max)
229 ret = NLOPT_MINF_MAX_REACHED;
230 else if (nlopt_stop_evals(stop)) ret = NLOPT_MAXEVAL_REACHED;
231 else if (nlopt_stop_time(stop)) ret = NLOPT_MAXTIME_REACHED;
232 if (ret != NLOPT_SUCCESS) goto done;
234 /**************************************
235 * Individual selection
236 **************************************/
237 /* all the individuals are copied to one vector to easily identify best solutions */
238 for (i=0; i < np; i++)
239 estotal[i] = esparents[i];
240 for (i=0; i < no; i++)
241 estotal[np+i] = esoffsprings[i];
243 nlopt_qsort_r(estotal, no+np, sizeof(Individual), NULL,
245 /* copy after sorting: */
246 for (i=0; i < no+np; i++) {
248 esparents[i] = estotal[i];
250 esoffsprings[i-np] = estotal[i];
252 } /* generations loop */
255 for (id=0; id < np; id++) free(esparents[id].parameters);
256 for (id=0; id < no; id++) free(esoffsprings[id].parameters);
258 if (esparents) free(esparents);
259 if (esoffsprings) free(esoffsprings);
260 if (estotal) free(estotal);