use libcmaes instead of our copy of the CMA-ES C code

src/examples/wscript

@@ -7,19 +7,19 @@ def build(bld):
                       source='mono_dim.cpp',
                       includes='. .. ../../',
                       target='mono_dim',
-                      uselib='BOOST EIGEN TBB',
+                      uselib='BOOST EIGEN TBB LIBCMAES',
                       use='limbo')
     obj = bld.program(features='cxx',
                       source='obs_multi.cpp',
                       includes='. .. ../../',
                       target='obs_multi',
-                      uselib='BOOST EIGEN TBB',
+                      uselib='BOOST EIGEN TBB LIBCMAES',
                       use='limbo')
     obj = bld.program(features='cxx',
                       source='obs_multi_auto_mean.cpp',
                       includes='. .. ../../',
                       target='obs_multi_auto_mean',
-                      uselib='BOOST EIGEN TBB',
+                      uselib='BOOST EIGEN TBB LIBCMAES',
                       use='limbo')
  #   obj = bld.program(features = 'cxx',
  #                          source = 'parego.cpp',

src/limbo/opt/cmaes.hpp

@@ -3,17 +3,11 @@
 
 #include <vector>
 #include <iostream>
-#include <limits>
-#include <stdio.h>
-#include <stdlib.h> /* free() */
-#include <stddef.h> /* NULL */
-
 #include <Eigen/Core>
 
-#include <cmaes/cmaes_interface.h>
-#include <cmaes/boundary_transformation.h>
-
-#include <limbo/tools/parallel.hpp>
+#ifdef USE_LIBCMAES
+# include <libcmaes/cmaes.h>
+#endif
 
 namespace limbo {
     namespace defaults {
@@ -30,99 +24,65 @@ namespace limbo {
             template <typename F>
             Eigen::VectorXd operator()(const F& f, double bounded) const
             {
-                // Currrently cmaes does not support unbounded search
-                assert(bounded);
-                int nrestarts = Params::cmaes::nrestarts();
                 size_t dim = f.param_size();
-                double incpopsize = 2;
-                cmaes_t evo;
-                double* const* pop;
-                double* fitvals;
-                double fbestever = 0, * xbestever = NULL;
-                double fmean;
-                int irun, lambda = 0, countevals = 0;
-                char const* stop;
-                boundary_transformation_t boundaries;
-                double lowerBounds[] = {0.0}; 
-                double upperBounds[] = {1.006309}; // Allows solution to be pretty close to 1
-                int nb_bounds = 1; /* numbers used from lower and upperBounds */
-
-                boundary_transformation_init(&boundaries, lowerBounds, upperBounds, nb_bounds);
-
-                double* x_in_bounds = cmaes_NewDouble(dim);
-                double init_point[dim];
-                for (int i = 0; i < dim; ++i)
-                    init_point[i] = f.init()(i);
-
-                for (irun = 0; irun < nrestarts + 1; ++irun) {
-
-                    fitvals = cmaes_init(&evo, dim, init_point, NULL, 0, lambda, NULL);
-
-                    evo.countevals = countevals;
-                    evo.sp.stopMaxFunEvals = Params::cmaes::max_fun_evals() < 0
-                        ? (900.0 * (dim + 3.0) * (dim + 3.0))
-                        : Params::cmaes::max_fun_evals();
-
-                    int pop_size = cmaes_Get(&evo, "popsize");
-                    double** all_x_in_bounds = new double* [pop_size];
-                    for (int i = 0; i < pop_size; ++i)
-                        all_x_in_bounds[i] = cmaes_NewDouble(dim);
-                    std::vector<Eigen::VectorXd> pop_eigen(pop_size, Eigen::VectorXd(dim));
-
-                    while (!(stop = cmaes_TestForTermination(&evo))) {
-                        pop = cmaes_SamplePopulation(&evo);
-                        tools::par::loop(0, pop_size, [&](int i) {
-                            // clang-format off
-                            boundary_transformation(&boundaries, pop[i], all_x_in_bounds[i], dim);
-                            for (int j = 0; j < dim; ++j)
-                              pop_eigen[i](j) = all_x_in_bounds[i][j];
-                            fitvals[i] = -f.utility(pop_eigen[i]);
-                            // clang-format on
-                        });
-                        cmaes_UpdateDistribution(&evo, fitvals);
-                    }
-
-                    for (int i = 0; i < pop_size; ++i)
-                        free(all_x_in_bounds[i]);
-
-                    lambda = incpopsize * cmaes_Get(&evo, "lambda");
-                    countevals = cmaes_Get(&evo, "eval");
-
-                    if (irun == 0 || cmaes_Get(&evo, "fbestever") < fbestever) {
-                        fbestever = cmaes_Get(&evo, "fbestever");
-                        xbestever = cmaes_GetInto(&evo, "xbestever",
-                            xbestever); /* alloc mem if needed */
-                    }
-                    const double* xmean = cmaes_GetPtr(&evo, "xmean");
-                    Eigen::VectorXd v(dim);
-                    for (int j = 0; j < v.size(); ++j)
-                        v(j) = xmean[j];
-
-                    if ((fmean = -f.utility(v)) < fbestever) {
-                        fbestever = fmean;
-                        xbestever = cmaes_GetInto(&evo, "xmean", xbestever);
-                    }
-
-                    cmaes_exit(&evo);
-
-                    if (stop) {
-                        if (strncmp(stop, "Fitness", 7) == 0 || strncmp(stop, "MaxFunEvals", 11) == 0) {
-                            //    printf("stop: %s", stop);
-                            break;
-                        }
-                    }
+#ifdef USE_LIBCMAES
+                using namespace libcmaes;
+
+                // create the parameter object
+                // boundary_transformation
+                double lbounds[dim],ubounds[dim]; // arrays for lower and upper parameter bounds, respectively
+                for (int i = 0; i < dim; i++) {
+                  lbounds[i] = 0.0;
+                  ubounds[i] = 1.005;
                 }
-                boundary_transformation(&boundaries, xbestever, x_in_bounds, dim);
-
-                Eigen::VectorXd result = Eigen::VectorXd::Zero(dim);
-                for (size_t i = 0; i < dim; ++i)
-                    result(i) = x_in_bounds[i];
-                free(xbestever);
-                boundary_transformation_exit(&boundaries);
-
-                free(x_in_bounds);
+                GenoPheno<pwqBoundStrategy> gp(lbounds, ubounds, dim);
+                // initial step-size, i.e. estimated initial parameter error.
+                // we suppose we are optimizing on [0, 1], but we have no idea where to start
+                double sigma = 0.5;
+                // initialize x0 as 0.50 in all 10 dimensions
+                // WARNING: we ignore the init() of the function to optimize!
+                // (because CMA-ES will perform its own random initialization)
+                std::vector<double> x0(dim, 0.5);
+                // -1 for automatically decided lambda, 0 is for random seeding of the internal generator.
+                CMAParameters<GenoPheno<pwqBoundStrategy>> cmaparams(dim, &x0.front(), sigma, -1 , 0, gp);
+                cmaparams.set_x0(0, 1.0);
+
+                // set multi-threading to true
+                cmaparams.set_mt_feval(true);
+                // aCMAES should be the best choice
+                // [see: https://github.com/beniz/libcmaes/wiki/Practical-hints ]
+                // but we want the restart -> aIPOP_CMAES
+                cmaparams.set_algo(aIPOP_CMAES);
+                cmaparams.set_restarts(Params::cmaes::nrestarts());
+                // if no max fun evals provided, we compute a recommended value
+                size_t max_evals = Params::cmaes::max_fun_evals() < 0
+                    ? (900.0 * (dim + 3.0) * (dim + 3.0))
+                    : Params::cmaes::max_fun_evals();
+                cmaparams.set_max_fevals(max_evals);
+                // max iteration is here only for security
+                cmaparams.set_max_iter(100000);
+                // we do not know if what is the actual maximum / minimum of the function
+                // therefore we deactivate this stopping criterion
+                cmaparams.set_stopping_criteria(FTARGET, false);
+
+                // wrap the function
+                FitFunc f_cmaes = [&](const double *x, const int n) {
+                  Eigen::Map<const Eigen::VectorXd> m(x, n);
+                  // remember that our optimizers maximize
+                  return -f.utility(m);
+                };
+
+                // the optimization itself
+                CMASolutions cmasols = cmaes<GenoPheno<pwqBoundStrategy>>(f_cmaes, cmaparams);
+                //cmasols.print(std::cout, 1, gp);
+                //to_f_representation
+                return gp.pheno(cmasols.get_best_seen_candidate().get_x_dvec());
+#else
+#warning NO libcmaes
+                assert(0);
+                return Eigen::VectorXd::Zero(dim);
+#endif
 
-                return result;
             }
         };
     }

src/tests/test_cmaes.cpp

+#define BOOST_TEST_DYN_LINK
+#define BOOST_TEST_MODULE nlopt_test
+
+#include <boost/test/unit_test.hpp>
+
+#include <limbo/limbo.hpp>
+#include <limbo/opt/cmaes.hpp>
+
+struct Params {
+    struct cmaes : public limbo::defaults::cmaes {
+    };
+};
+
+struct TestOpt {
+public:
+    double utility(const Eigen::VectorXd& params) const {
+        return -params(0) * params(0) - params(1) * params(1);
+    }
+
+    size_t param_size() const {
+          return 2;
+    }
+};
+
+BOOST_AUTO_TEST_CASE(test_cmaes)
+{
+    TestOpt util;
+    Eigen::VectorXd g = limbo::opt::Cmaes<Params>()(util, true);
+
+    BOOST_CHECK_SMALL(g(0), 0.00000001);
+    BOOST_CHECK_SMALL(g(1), 0.00000001);
+}

src/tests/wscript

@@ -9,7 +9,7 @@ def build(bld):
                 source='bo_functions.cpp',
                 includes='. .. ../../',
                 target='bo_functions',
-                uselib='BOOST EIGEN TBB',
+                uselib='BOOST EIGEN TBB LIBCMAES',
                 use='limbo')
     bld.program(features='cxx test',
                 source='test_gp.cpp',
@@ -27,7 +27,7 @@ def build(bld):
                 source='test_optimizers.cpp',
                 includes='. .. ../../',
                 target='test_optimizers',
-                uselib='BOOST EIGEN TBB',
+                uselib='BOOST EIGEN TBB LIBCMAES',
                 use='limbo')
     bld.program(features='cxx test',
                 source='test_macros.cpp',
@@ -48,6 +48,13 @@ def build(bld):
                     target='test_nlopt',
                     uselib='BOOST EIGEN TBB NLOPT',
                     use='limbo')
+    if bld.env.DEFINES_NLOPT:
+        bld.program(features='cxx test',
+                    source='test_cmaes.cpp',
+                    includes='. .. ../../',
+                    target='test_cmaes',
+                    uselib='BOOST EIGEN TBB LIBCMAES',
+                    use='limbo')
 
 
 def build_extensive_tests(ctx):

waf_tools/libcmaes.py

+#! /usr/bin/env python
+# encoding: utf-8
+# JB Mouret / Inria - 2015
+
+"""
+Quick n dirty libcmaes detection
+"""
+
+import os, glob, types
+from waflib.Configure import conf
+
+
+def options(opt):
+	opt.add_option('--libcmaes', type='string', help='path to libcmaes', dest='libcmaes')
+
+
+@conf
+def check_libcmaes(conf):
+	if conf.options.nlopt:
+		includes_check = [conf.options.nlopt + '/include']
+		libs_check = [conf.options.nlopt + '/lib']
+	else:
+		includes_check = ['/usr/local/include', '/usr/include']
+		libs_check = ['/usr/local/lib', '/usr/lib']
+
+	try:
+		conf.start_msg('Checking for libcmaes includes')
+		res = conf.find_file('libcmaes/cmaes.h', includes_check)
+		conf.end_msg('ok')
+	except:
+		conf.end_msg('Not found', 'RED')
+		return 1
+	conf.start_msg('Checking for libcmaes libs')
+	found = False
+	for lib in ['libcmaes.so', 'libcmaes.a', 'libcmaes.dylib']:
+		try:
+			found = found or conf.find_file(lib, libs_check)
+		except:
+			continue
+	if not found:
+		conf.end_msg('Not found', 'RED')
+		return 1
+	else:
+		conf.end_msg('ok')
+		conf.env.INCLUDES_LIBCMAES = includes_check
+		conf.env.LIBPATH_LIBCMAES = libs_check
+		conf.env.DEFINES_LIBCMAES = ['USE_LIBCMAES']
+		conf.env.LIB_LIBCMAES= ['cmaes']
+	return 1

wscript

@@ -26,7 +26,8 @@ def options(opt):
         opt.load('limbo')
         opt.load('openmp')
         opt.load('nlopt')
-        #opt.load('ode')
+        opt.load('libcmaes')
+
         opt.add_option('--exp', type='string', help='exp(s) to build, separate by comma', dest='exp')
         opt.add_option('--qsub', type='string', help='config file (json) to submit to torque', dest='qsub')
         opt.add_option('--oar', type='string', help='config file (json) to submit to oar', dest='oar')
@@ -45,6 +46,7 @@ def configure(conf):
         conf.load('mkl')
         conf.load('xcode')
         conf.load('nlopt')
+        conf.load('libcmaes')
 
         if conf.env.CXX_NAME in ["icc", "icpc"]:
             common_flags = "-Wall -std=c++11"
@@ -65,7 +67,7 @@ def configure(conf):
         conf.check_openmp()
         conf.check_mkl()
         conf.check_nlopt()
-        #conf.check_ode()
+        conf.check_libcmaes()
 
         if conf.env['CXXFLAGS_ODE']:
                 common_flags += ' ' + conf.env['CXXFLAGS_ODE']