Minimization by simulated annealing is performed by descendants of sim_anneal (see o2scl::anneal_gsl).
Simulated annealing example
This example minimizes the function
over
where
is the Bessel function given in gsl_sf_bessel_J0
. The initial guess at
is far away from the global minimum.
The plot below plots the function above, the initial guess, and the minimum obtained by the example program.
Results of simulated annealing minimization
#include <iostream>
#include <cmath>
#include <boost/numeric/ublas/vector.hpp>
#include <gsl/gsl_sf_bessel.h>
#include <o2scl/multi_funct.h>
#include <o2scl/funct.h>
#include <o2scl/anneal_gsl.h>
#include <o2scl/test_mgr.h>
#include <o2scl/hdf_file.h>
#include <o2scl/hdf_io.h>
using namespace std;
void make_plot_data();
double bessel_fun(
size_t nvar,
const ubvector &x) {
double a, b;
a=(x[0]-2.0);
b=(x[1]+3.0);
if (fabs(x[0])>10.0 || fabs(x[1])>10.0) return 10.0;
return -gsl_sf_bessel_J0(a)*gsl_sf_bessel_J0(b);
}
int main(int argc, char *argv[]) {
cout.setf(ios::scientific);
double result;
double step[1]={10.0};
init[0]=6.0;
init[1]=7.0;
ga.
mmin(2,init,result,fx);
cout << "x: " << init[0] << " " << init[1]
<< ", minimum function value: " << result << endl;
cout << endl;
t.
test_rel(init[0],2.0,1.0e-3,
"another test - value");
t.
test_rel(init[1],-3.0,1.0e-3,
"another test - value 2");
t.
test_rel(result,-1.0,1.0e-3,
"another test - min");
make_plot_data();
return 0;
}