This is George Bilchev's <G.Bilchev@plymouth.ac.uk>
inductive search algorithm (in C++) that won the first place at the 
<A HREF="http://iridia.ulb.ac.be/langerman/ICEO.html">
1st International Contest on Evolutionary Computation</A>
on a real-valued function testsuite. It builds a population of
points in dimensions i=1,2,...,n in turn by starting from a single point
population created from the best point of the previous dimension
a global one-dimensional search for the best new coordinate,
and then growing the population through a sequence of local searches.

************************************************************************
*** The location of this file may change soon.
*** Please refer at the moment to
*** http://solon.cma.univie.ac.at/~neum/glopt/test_results.html#bilchev
*** 
*** I'd appreciate if someone would remove the dependence of these
*** programs on routines not in the public domain, and adds a driver
*** routine that allows to optimize a function defined in a 
*** user-provided subroutine
*** 
************************************************************************


Dear Arnold,

Please find below exactly the version of my algorithm for the EC 
contest, i.e. an algorithm that is specialized to the test bed of the 
contest. However, this is the same algorithm which I used for your 
problems (excluding the toy protein-folding, which required the multi-
population generalization.) In order to be compiled the code requires 
the LEDA C++ library (public domain) and some routines from numerical
recipes in C (see Makefile).

Best wishes,
George

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
*** brent's univariate minimizer is also available publicly in 
http://www.netlib.no/netlib/c/brent.shar
http://www.netlib.no/netlib/c/serv.shar (needed for brent.shar)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

George Bilchev remarks elsewhere:
>>you can use your favourite local search and 1-D
global search (which I encourage you do since my choice
could have been better). For more complex cost functions
it is also better to maintain a population of "promising"
points at each dimension (the algorithm only explains for
one point).<<


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

--------------------
Makefile:
--------------------
all:
    gcc -c brent.c -ansi
    gcc -c nrutil.c -ansi
    g++ -o main main.cc nrutil.o brent.o -L/home/george/C++/LEDA-3.1c 
-I/home/george/C++/LEDA-3.1c/incl -lG -lL -lP -lm

------------------
main.cc
------------------
#include <Normal.h>
#include <Uniform.h>
#include <DiscUnif.h>
#include <ACG.h>
#include <sys/time.h>
#include <iostream.h>
#include <fstream.h>
#include <stdio.h>
#include <stdlib.h>
#include <LEDA/sortseq.h>
#include <LEDA/stream.h>
#include <math.h>
#include <iomanip.h>
#include "main.h"

/*
**  Include a test case here:
**
*/

#include "t1.c"



float func(float y)

/*
**
**   Auxiliary 1-D function 
*/

{
    iter_count++;   

    x[DIM-1] = y;
    float res =  f(x,DIM);
    return res;
}


float func_learning(float *y)

/*
**
**   Auxiliary DIM-dimensional function 
*/

{
    iter_count++;   
    float res = f(y,DIM);
    return res;
}


void oracle(int dim) 

/*
**
**   One possible deterministic implementation
**       of the non-deterministic ORACLE function 
*/

{
    float r,bren,xmin;  
    
    float ax = LOW;
    float cx = HIGH;
    float bx = (ax+cx)/2;

    DIM = dim;

    sortseq<float,Interval> Population;
    seq_item S;
    Interval I;
    I.low = ax;
    I.high = cx;
    Population.insert(cx-ax,I);

    int STOP = 0;
    int count = 0;
    float FMIN = 1e30, XMIN;

    while (STOP == 0) {

        S = Population.max();
        Population.del_item(S);
        I = Population.inf(S);
        
        ax = I.low;
        cx = I.high;
        bx = (ax+cx)/2;

        bren = brent(ax,bx,cx,func,TOL,&xmin);


        if (bren < FMIN) { FMIN = bren; XMIN = xmin; }
    

        if (bx < xmin) {
            I.low = ax; I.high = bx; Population.insert(I.high-I.low, 
I);
            I.low = bx; I.high = xmin; Population.insert(I.high-
I.low, I);
            I.low = xmin; I.high = cx; Population.insert(I.high-
I.low, I);
        } else {
            I.low = ax; I.high = xmin; Population.insert(I.high-
I.low, I);
            I.low = xmin; I.high = bx; Population.insert(I.high-
I.low, I);
            I.low = bx; I.high = cx; Population.insert(I.high-I.low, 
I);
        }
 
    count++;
    if (count > DELTA_N) break;
    }



//  Local learning

    x[DIM-1] = XMIN;

    if ( (LEARNING) && (DIM > 1) ) { 
        local_learn(func_learning, x, DIM); 
    }   
}


main() 
{
    for(DELTA_N=STOP_CRIT; DELTA_N<=STOP_CRIT; DELTA_N++) {
        iter_count = 0;
        for (int i=0; i<NDIM; i++) {
            oracle(i+1);
            cout << "Count: " << setw(5) << --iter_count;
            /* only a monitoring call - reduce iteration count */
            cout << "      dim: " << setw(3) << i+1 << "     f: " << 
setw(7) << f(x,i+1);
            cout << endl;
        }
    }
}

-------------------
main.h
-------------------
extern "C" float brent(float ax, float bx, float cx, float 
(*f)(float), float tol, float *xmin);
#define TOL 1.0e-4

unsigned long seed = time(NULL);
ACG generator(seed, 1000);

#include "rand.c"

#define NDIM 10

float   x[NDIM];
int DIM=0;
int DELTA_N;

static  int iter_count = 0;

#include "dhc.c"


struct Interval {
    float low,high;
};

void Print(Interval &I, ostream& out) {
    out << I.low << "   " << I.high;
}

void Read(Interval &I, istream& in) {
    in >> I.low >>  I.high;
}

void print(sortseq<float,Interval>& S)
{ seq_item it;
  forall_items(it,S) cout << S.key(it) << " ==> " << S.inf(it).low 
              << "   " << S.inf(it).high << endl;
  newline;
}

---------------
rand.c
----------------

float frandom(float l, float h) {
    float eps = 0.0000001;
    Uniform rnd(l-eps, h+eps, &generator);
    return rnd();
}

int irandom(int l, int h) {
    DiscreteUniform rnd(l, h, &generator);
    return (int)rnd();
}

int irandom(int l, int h, int prev) {
    DiscreteUniform rnd(l, h, &generator);
    int rand = (int)rnd();
    while ( rand == prev ) { rand = (int)rnd(); }
    return rand;
}

--------------------
t1.c
--------------------
#define STOP_CRIT 0
#define LEARNING 0

#define LOW -5.0
#define HIGH 5.0


float f(float *x,int n)
{
 
        register int    i;
        float       Sum;

        for (Sum = 0.0, i = 0; i < n; i++) {
            Sum += x[i]*x[i];   
        }

        return (Sum);
}


---------
t2.c
---------
#define STOP_CRIT 0
#define LEARNING 0

#define LOW -600.0
#define HIGH 600.0
#define D 4000.0

float f(float *x,int  n)

{
 
        register int    i;

        float       Val1,
            Val2,
            Sum;

    for (Val1 = 0.0, Val2 = 1.0, i = 0; i < n; i++) {
        Val1 += x[i] * x[i];
        Val2 *= cos(x[i] / sqrt((double) (i + 1)) );
    }

    Sum = Val1 / D - Val2 + 1.0;
        return (Sum);
}

---------
t3.c
---------
#include "shekel.h"

#define STOP_CRIT 1
#define LEARNING 1

#define LOW 0.0
#define HIGH 10.0

 

float   f(float *x,int n)
{
    register    i, j;
    float   sp, h, result = 0.0;

    for (i = 0; i < 30; i++) {
        sp = 0.0;
        for (j = 0; j < n; j++) {
            h   = x[j] - a[i][j];
            sp += h * h;
        }
        result += 1.0 / (sp + c[i]);
    }
    return( - result);
}

---------
t4.c
---------
#define STOP_CRIT  10  /* 3   if NDIM = 5 */
#define LEARNING 0


/* #define STOP_CRIT 10   if  NDIM = 10 */

#define LOW 0.0
#define HIGH 3.14

#define m  10.0

float   f(float *x,int  n) /* Michalewitz */
{ 
float   PI=3.1415927;
float   u;
int     i;

    u=0; 

    for (i=0;i<n;i++) 
        u = u + sin(x[i])
          * pow(sin((i+1)*x[i]*x[i]/PI),2.0*m); 


    return(-u); 
}

-------------
t5.c
-------------

#include "langerman.h"

#define STOP_CRIT 2
#define LEARNING 1

#define LOW 0.0
#define HIGH 10.0

double min = 500;


float SqrDst(float *x1, float *x2, int n )
{
float dist = 0.0, d;
int i;

dist = 0;

for ( i = 0; i<n; i++ )
  {
    d =x1[i] - x2[i];
    dist += d*d;
  }

return (dist);
}



float   f(float *x, int n)  /* Langerman's function */

{
 
        register int    i;

        float       Sum,
                    PI,
                    dist;


    PI  = 3.141592653;
    
    Sum = 0;
    for ( i = 0; i < 5; i++ )
      {

        dist = SqrDst( x, a[i], n );
        Sum -= c[i] * (exp(-dist/PI) * cos( PI * dist ) );
      }

        return (Sum);
}

------
dhc.c  (dynamic hill climbing) (available from the authors)
------

#include <stdlib.h>
#include <stdio.h>
#include <math.h>
                                        

#define THRESHOLD 1E-4      /* Minimum step size */

#define INIT_SIZE 1.0       /* Initial step size */



/* You should really look at my thesis to understand the following 
code. */
/* This comes from D. Yuret, From GA's to efficient optimization, MSC
thesis MIT May 1994. */

void local_learn(float (*f)(float *x), float *x, int DIM)
{
  float u[DIM], v[DIM], xv[DIM];
  float fx,fxv;
  int vi,vvec;
  float vr;
  int i,iter,maxiter;

  for(i=0; i<DIM; i++)
    u[i] = v[i] = 0;
  vi = -1; vvec = 1;
  vr = -INIT_SIZE;
  fx = f(x);
  fxv = 1E10;

/*  printf("%d. %.7 <= ", iter_count, fx);
  for(i=0; i<DIM; i++) { printf("%.4f ", x[i]); }; printf("\n");
*/

  while(fabs(vr) >= THRESHOLD) {
    maxiter = ((fabs(vr) < 2*THRESHOLD) ? 2*DIM : 2);
    iter = 0;
    while((fxv >= fx) && (iter < maxiter)) {
      if(iter == 0) { for(i=0; i<DIM; i++) xv[i] = x[i]; }
      else xv[vi] -= vr;
      if(vvec) vvec = 0;
      vr = -vr;
      if(vr > 0) vi = ((vi+1) % DIM);
      xv[vi] += vr;
      fxv = f(xv);
      iter++;
    }
    if(fxv >= fx) {
      fxv = 1E10;
      vr /= 2;
    } else {
      fx = fxv;  /* printf("%d. %.7f <= ", iter_count, fx);*/
      for(i=0; i<DIM; i++) { x[i] = xv[i]; /* printf("%.4f ", x[i]); 
*/ }
      /* printf("\n");*/
      if(iter == 0) {
    if(vvec) {
      for(i=0; i<DIM; i++) {
        u[i] += v[i]; v[i] *= 2; xv[i] += v[i];
      }
      vr *= 2;
    } else {
      u[vi] += vr; vr *= 2; xv[vi] += vr;
    }
    fxv = f(xv);
      } else {
    for(i=0; i<DIM; i++) xv[i] += u[i];
    xv[vi] += vr;
    fxv = f(xv);
    if(fxv >= fx) {
      for(i=0; i<DIM; i++) { u[i] = 0; xv[i] = x[i]; }
      u[vi] = vr; vr *= 2;
      xv[vi] += vr; fxv = f(xv);
    } else {
      for(i=0; i<DIM; i++) x[i] = xv[i]; fx = fxv;
      u[vi] += vr;
      for(i=0; i<DIM; i++) v[i] = 2*u[i]; vvec = 1;
      for(i=0; i<DIM; i++) xv[i] += v[i]; fxv = f(xv);
      for(vr=0,i=0; i<DIM; i++) vr += v[i]*v[i];
      vr = sqrt(vr);
    }
      }
    }
  }
}