From: stevenj Date: Fri, 8 Feb 2013 17:12:54 +0000 (-0500) Subject: added ESCH algorithm X-Git-Url: http://www.chiark.greenend.org.uk/ucgi/~ianmdlvl/git?a=commitdiff_plain;h=f40d6be74cae6fddbb148febe80fe2c8e4139f19;p=nlopt.git added ESCH algorithm Ignore-this: cc8b951db64798b0e9a6496591ef0896 darcs-hash:20130208171254-c8de0-39328682d256f5b2406e43e28c4c928e31c72718.gz --- diff --git a/Makefile.am b/Makefile.am index 3b22293..dc344e9 100644 --- a/Makefile.am +++ b/Makefile.am @@ -8,7 +8,7 @@ CXX_DIRS = stogo CXX_LIBS = stogo/libstogo.la endif -SUBDIRS = util direct cdirect $(CXX_DIRS) praxis luksan crs mlsl mma cobyla newuoa lbfgs neldermead auglag bobyqa isres slsqp api . octave test swig +SUBDIRS = util direct cdirect $(CXX_DIRS) praxis luksan crs mlsl mma cobyla newuoa lbfgs neldermead auglag bobyqa isres slsqp esch api . octave test swig EXTRA_DIST = autogen.sh nlopt.pc.in m4 if WITH_NOCEDAL @@ -21,7 +21,7 @@ cdirect/libcdirect.la $(CXX_LIBS) praxis/libpraxis.la $(NOCEDAL_LBFGS) \ luksan/libluksan.la crs/libcrs.la mlsl/libmlsl.la mma/libmma.la \ cobyla/libcobyla.la newuoa/libnewuoa.la neldermead/libneldermead.la \ auglag/libauglag.la bobyqa/libbobyqa.la isres/libisres.la \ -slsqp/libslsqp.la api/libapi.la util/libutil.la +slsqp/libslsqp.la esch/libesch.la api/libapi.la util/libutil.la if WITH_CXX libnlopt@NLOPT_SUFFIX@_la_LDFLAGS = -version-info @SHARED_VERSION_INFO@ diff --git a/api/Makefile.am b/api/Makefile.am index 4bfe52f..7605f8f 100644 --- a/api/Makefile.am +++ b/api/Makefile.am @@ -1,4 +1,4 @@ -AM_CPPFLAGS = -I$(top_srcdir)/cdirect -I$(top_srcdir)/direct -I$(top_srcdir)/stogo -I$(top_srcdir)/praxis -I$(top_srcdir)/lbfgs -I$(top_srcdir)/luksan -I$(top_srcdir)/crs -I$(top_srcdir)/mlsl -I$(top_srcdir)/mma -I$(top_srcdir)/cobyla -I$(top_srcdir)/newuoa -I$(top_srcdir)/neldermead -I$(top_srcdir)/auglag -I$(top_srcdir)/bobyqa -I$(top_srcdir)/isres -I$(top_srcdir)/slsqp -I$(top_srcdir)/util +AM_CPPFLAGS = -I$(top_srcdir)/cdirect -I$(top_srcdir)/direct -I$(top_srcdir)/stogo -I$(top_srcdir)/praxis -I$(top_srcdir)/lbfgs -I$(top_srcdir)/luksan -I$(top_srcdir)/crs -I$(top_srcdir)/mlsl -I$(top_srcdir)/mma -I$(top_srcdir)/cobyla -I$(top_srcdir)/newuoa -I$(top_srcdir)/neldermead -I$(top_srcdir)/auglag -I$(top_srcdir)/bobyqa -I$(top_srcdir)/isres -I$(top_srcdir)/slsqp -I$(top_srcdir)/esch -I$(top_srcdir)/util include_HEADERS = nlopt.h nlopt.f nlopt.hpp noinst_LTLIBRARIES = libapi.la diff --git a/api/general.c b/api/general.c index 8ce76d6..39287f6 100644 --- a/api/general.c +++ b/api/general.c @@ -97,6 +97,7 @@ static const char nlopt_algorithm_names[NLOPT_NUM_ALGORITHMS][256] = { "Multi-level single-linkage (MLSL), quasi-random (global, needs sub-algorithm)", "Sequential Quadratic Programming (SQP) (local, derivative)", "CCSA (Conservative Convex Separable Approximations) with simple quadratic approximations (local, derivative)", + "ESCH evolutionary strategy" }; const char * NLOPT_STDCALL nlopt_algorithm_name(nlopt_algorithm a) diff --git a/api/nlopt.h b/api/nlopt.h index e09877d..eecf231 100644 --- a/api/nlopt.h +++ b/api/nlopt.h @@ -150,6 +150,8 @@ typedef enum { NLOPT_LD_CCSAQ, + NLOPT_GN_ESCH, + NLOPT_NUM_ALGORITHMS /* not an algorithm, just the number of them */ } nlopt_algorithm; diff --git a/api/optimize.c b/api/optimize.c index 6522c51..329c5f9 100644 --- a/api/optimize.c +++ b/api/optimize.c @@ -60,6 +60,7 @@ static int my_isnan(double x) { return x != x; } #include "auglag.h" #include "bobyqa.h" #include "isres.h" +#include "esch.h" #include "slsqp.h" /*********************************************************************/ @@ -350,6 +351,7 @@ static int elimdim_wrapcheck(nlopt_opt opt) case NLOPT_LN_NELDERMEAD: case NLOPT_LN_SBPLX: case NLOPT_GN_ISRES: + case NLOPT_GN_ESCH: case NLOPT_GD_STOGO: case NLOPT_GD_STOGO_RAND: return 1; @@ -780,6 +782,13 @@ static nlopt_result nlopt_optimize_(nlopt_opt opt, double *x, double *minf) lb, ub, x, minf, &stop, (int) POP(0)); + case NLOPT_GN_ESCH: + if (!finite_domain(n, lb, ub)) return NLOPT_INVALID_ARGS; + return chevolutionarystrategy(n, f, f_data, + lb, ub, x, minf, &stop, + (unsigned) POP(0), + (unsigned) (POP(0)*1.5)); + case NLOPT_LD_SLSQP: return nlopt_slsqp(n, f, f_data, opt->m, opt->fc, diff --git a/configure.ac b/configure.ac index 3bc1c2d..69247c2 100644 --- a/configure.ac +++ b/configure.ac @@ -411,6 +411,7 @@ AC_CONFIG_FILES([ bobyqa/Makefile isres/Makefile slsqp/Makefile + esch/Makefile test/Makefile swig/Makefile swig/nlopt.scm diff --git a/esch/COPYRIGHT b/esch/COPYRIGHT new file mode 100644 index 0000000..e00ee0d --- /dev/null +++ b/esch/COPYRIGHT @@ -0,0 +1,22 @@ +/* Copyright (c) 2008-2013 Carlos Henrique da Silva Santos + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE + * LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION + * OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION + * WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + */ + diff --git a/esch/Makefile.am b/esch/Makefile.am new file mode 100644 index 0000000..377c00f --- /dev/null +++ b/esch/Makefile.am @@ -0,0 +1,6 @@ +AM_CPPFLAGS = -I$(top_srcdir)/util -I$(top_srcdir)/api + +noinst_LTLIBRARIES = libesch.la +libesch_la_SOURCES = esch.c esch.h + +EXTRA_DIST = README COPYRIGHT diff --git a/esch/README b/esch/README new file mode 100644 index 0000000..e34c8ce --- /dev/null +++ b/esch/README @@ -0,0 +1,20 @@ +This is Carlos Henrique da Silva Santos's implementation of +a modified version of the Evolutionary Algorithm described in +the following paper and Ph.D. thesis: + + C. H. da Silva Santos, M. S. Gonçalves, and H. E. Hernandez-Figueroa, + "Designing Novel Photonic Devices by Bio-Inspired Computing," + IEEE Photonics Technology Letters 22(15), pp. 1177-1179 (2010). + + C. H. da Silva Santos, "Parallel and Bio-Inspired Computing + Applied to Analyze Microwave and Photonic Metamaterial Strucutures," + University of Campinas, (2010) + http://www.bibliotecadigital.unicamp.br/document/?code=000767537&opt=4&lg=en_US + +It is distributed under the "MIT license" given in the attached +COPYRIGHT file (similar to the rest of NLopt), and was financially by +the São Paulo Science Foundation (FAPESP - Fundação de Amparo à +Pesquisa do Estado de São Paulo) under the grant 2012/14553-9. + +Carlos Henrique da Silva Santos +January 2013 diff --git a/esch/esch.c b/esch/esch.c new file mode 100644 index 0000000..19e4de9 --- /dev/null +++ b/esch/esch.c @@ -0,0 +1,258 @@ +/* Copyright (c) 2008-2013 Carlos Henrique da Silva Santos + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE + * LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION + * OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION + * WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + */ + +#include +#include "esch.h" + +// --------------------------------------------------------------------------- +// Cauchy random number distribution +static double randcauchy(const double params[7]) { + /* double min, double max, double mi, double t, double band */ + double na_unif, cauchy_mit, limit_inf, limit_sup; + double valor; + double min = params[1], max = params[2], mi = params[3], + t = params[4], band = params[5]; + limit_inf = mi - (band/2); + limit_sup = mi + (band/2); + do { + na_unif = nlopt_urand(0,1); // ran2(0,1); + cauchy_mit = t*tan((na_unif-(1/2))*M_PI) + mi; + } while ( (cauchy_mitlimit_sup) ); + + if (cauchy_mit < 0) + cauchy_mit = -cauchy_mit; + else + cauchy_mit = cauchy_mit + (band/2); + valor = (cauchy_mit*100/band)/100; + valor = min+(max-min)*valor; + return valor; +} + +// --------------------------------------------------------------------------- + +// main Individual representation type +typedef struct IndividualStructure { + double * parameters; + double fitness; +} Individual; + +static int CompareIndividuals(void *unused, const void *a_, const void *b_) { + (void) unused; + const Individual *a = (const Individual *) a_; + const Individual *b = (const Individual *) b_; + return a->fitness < b->fitness ? -1 : (a->fitness > b->fitness ? +1 : 0); +} + +nlopt_result chevolutionarystrategy( + unsigned nparameters, /* Number of input parameters */ + nlopt_func f, /* Recursive Objective Funtion Call */ + void * data_f, /* Data to Objective Function */ + const double* lb, /* Lower bound values */ + const double* ub, /* Upper bound values */ + double* x, /*in: initial guess, out: minimizer */ + double* minf, + nlopt_stopping* stop, /* nlopt stop condition */ + unsigned np, /* Number of Parents */ + unsigned no) { /* Number of Offsprings */ + + // variables from nlopt + nlopt_result ret = NLOPT_SUCCESS; + double vetor[8]; + unsigned i, id, item; + int parent1, parent2; + unsigned crosspoint; // crossover parameteres + int contmutation, totalmutation; // mutation parameters + int idoffmutation, paramoffmutation; // mutation parameters + Individual * esparents; // Parents population + Individual * esoffsprings; // Offsprings population + Individual * estotal;// copy containing Parents and Offsprings pops + /* It is interesting to maintain the parents and offsprings + * populations stablished and sorted; when the final iterations + * is achieved, they are ranked and updated. */ + + // ------------------------------- + // controling the population size + // ------------------------------- + if (!np) np = 40; + if (!no) no = 60; + if ((np < 1)||(no<1)) return NLOPT_INVALID_ARGS; + esparents = (Individual*) malloc(sizeof(Individual) * np); + esoffsprings = (Individual*) malloc(sizeof(Individual) * no); + estotal = (Individual*) malloc(sizeof(Individual) * (np+no)); + if ((!esparents)||(!esoffsprings)||(!estotal)) { + free(esparents); free(esoffsprings); free(estotal); + return NLOPT_OUT_OF_MEMORY; + } + for (id=0; id < np; id++) esparents[id].parameters = NULL; + for (id=0; id < no; id++) esoffsprings[id].parameters = NULL; + // From here the population is initialized + /* we don't handle unbounded search regions; + this check is unnecessary since it is performed in nlopt_optimize. + for (j = 0; j < nparameters; ++j) + if (nlopt_isinf(low[j]) || nlopt_isinf(up[j])) + return NLOPT_INVALID_ARGS; + */ + // main vector of parameters to randcauchy + vetor[0] = 4; // ignored + vetor[3] = 0; + vetor[4] = 1; + vetor[5] = 10; + vetor[6] = 1; + vetor[7] = 0; // ignored + // ------------------------------------ + // Initializing parents population + // ------------------------------------ + for (id=0; id < np; id++) { + esparents[id].parameters = + (double*) malloc(sizeof(double) * nparameters); + if (!esparents[id].parameters) { + ret = NLOPT_OUT_OF_MEMORY; + goto done; + } + for (item=0; itemnevals++; + if (*minf > esparents[id].fitness) { + *minf = esparents[id].fitness; + memcpy(x, esparents[id].parameters, + nparameters * sizeof(double)); + } + if (nlopt_stop_forced(stop)) ret = NLOPT_FORCED_STOP; + else if (*minf < stop->minf_max) ret = NLOPT_MINF_MAX_REACHED; + else if (nlopt_stop_evals(stop)) ret = NLOPT_MAXEVAL_REACHED; + else if (nlopt_stop_time(stop)) ret = NLOPT_MAXTIME_REACHED; + if (ret != NLOPT_SUCCESS) goto done; + } + // ------------------------------------ + // Main Loop - Generations + // ------------------------------------ + while (1) { + // ------------------------------------ + // Crossover + // ------------------------------------ + for (id=0; id < no; id++) + { + parent1 = nlopt_iurand((int) np); + parent2 = nlopt_iurand((int) np); + crosspoint = (unsigned) nlopt_iurand((int) nparameters); + for (item=0; item < crosspoint; item++) + esoffsprings[id].parameters[item] + = esparents[parent1].parameters[item]; + for (item=crosspoint; item < nparameters; item++) + esoffsprings[id].parameters[item] + = esparents[parent2].parameters[item]; + } + // ------------------------------------ + // Gaussian Mutation + // ------------------------------------ + totalmutation = (int) ((no * nparameters) / 10); + if (totalmutation < 1) totalmutation = 1; + for (contmutation=0; contmutation < totalmutation; + contmutation++) { + idoffmutation = nlopt_iurand((int) no); + paramoffmutation = nlopt_iurand((int) nparameters); + vetor[1] = lb[paramoffmutation]; + vetor[2] = ub[paramoffmutation]; + vetor[7] = vetor[7]+contmutation; + esoffsprings[idoffmutation].parameters[paramoffmutation] + = randcauchy(vetor); + } + // ------------------------------------ + // Offsprings fitness evaluation + // ------------------------------------ + for (id=0; id < no; id++){ + //esoffsprings[id].fitness = (double)fitness(esoffsprings[id].parameters, nparameters,fittype); + esoffsprings[id].fitness = f(nparameters, esoffsprings[id].parameters, NULL, data_f); + estotal[id+np].fitness = esoffsprings[id].fitness; + stop->nevals++; + if (*minf > esoffsprings[id].fitness) { + *minf = esoffsprings[id].fitness; + memcpy(x, esoffsprings[id].parameters, + nparameters * sizeof(double)); + } + if (nlopt_stop_forced(stop)) ret = NLOPT_FORCED_STOP; + else if (*minf < stop->minf_max) + ret = NLOPT_MINF_MAX_REACHED; + else if (nlopt_stop_evals(stop)) ret = NLOPT_MAXEVAL_REACHED; + else if (nlopt_stop_time(stop)) ret = NLOPT_MAXTIME_REACHED; + if (ret != NLOPT_SUCCESS) goto done; + } + // ------------------------------------ + // Individual selection + // ------------------------------------ + // all the individuals are copied to one vector to easily identify best solutions + for (i=0; i < np; i++) + estotal[i] = esparents[i]; + for (i=0; i < no; i++) + estotal[np+i] = esoffsprings[i]; + // Sorting + nlopt_qsort_r(estotal, no+np, sizeof(Individual), NULL, + CompareIndividuals); + // copy after sorting: + for (i=0; i < no+np; i++) { + if (i