// Example - how to use CMA-ES 
// package available at https://www.lri.fr/~hansen/cmaes_inmatlab.html#scilab

// test functions
function y=SE_f1(x) // sphere
//    n=size(x,2);
//    y=0;
//    for j=1:n
//        y=y+x(j)^2;
//    end
    y=sum(x(1:$).^2)
endfunction

function y=SE_f2(x)   // Griewank
    n=size(x,2);
    s=0;p=1;
    for j=1:n
        s=s+x(j)^2;
        p=p*cos(x(j)/sqrt(j));   
    end   
    y=s/4000-p+1;
endfunction

function y=SE_f3(x) // Rosenbrock
      y = 1e2*sum((x(1:$-1).^2 - x(2:$)).^2) + sum((x(1:$-1)-1).^2);
endfunction

// usage example (source: https://www.lri.fr/~hansen/scilabdoc/cma_new.htm)
// object creation
   p = cma_new();  // return (default) parameters
   disp(p);        // investigate default parameters
   p.readme        // display doc for parameters
   p.stop          // investigate the possible and default termination criteria
   clear p         // we do not need the defaults for initialization

   // Initialization: 
   // 1. setting of the start configuration (column vector with as many elements as is the problem size)
   // 2. setting the initial value of sigma
   // 3. setting of other parameters 
   p.x0 = rand(10,1);  // mandatory: initial solution vector
   p.sigma0 = 0.3;     // mandatory: initial standard deviation on x0 
                       //    (somewhat a search range) 
   // instantiation of a cmaes object with the specified parameters 
   es = cma_new(p);  // get an cma "object", constructor call
   // specification of the objective function
   fitfun = SE_f3;  // objective function needs to be defined

   // iteration loop 
   while ~ cma_stop(es)
     X = cma_ask(es);  // list of solution points (column vectors)
     y = [];            // just in case
     for i = 1:length(X)
       y(i) = fitfun(X(i));    // evaluate each point
       disp(X(i),"X(i)=");
       disp(y(i),"y(i)=");
     end
     es = cma_tell(es, X, y');  // update for the next iteration
   end

   // output
   [ybest, xbest] = cma_best(es);
   printf(" \n x(1) = %f   y = %f\n", xbest(1), ybest);

   for str = es.out.stopflags
     disp(' ' + str);  // display termination reasons
   end

   disp(es.out); // further investigate some useful output
