![]() |
Object-oriented Scientific Computing Library: Version 0.910
|
Minimization by simulated annealing is performed by descendants of sim_anneal (see gsl_anneal). Because simulated annealing is particularly well-suited to parallelization, a multi-threaded minimizer analogous to gsl_anneal is given in anneal_mt. This header-only class uses the Boost libraries and requires it for use.
This example minimizes the function
over where
is the Bessel function given in
gsl_sf_bessel_J0
. The initial guess at is several local minima away from the global minimum.
The plot below plots the function above, the initial guess, and the minimum obtained by the example program.
/* Example: ex_anneal.cpp ------------------------------------------------------------------- An example to demonstrate minimization by simulated annealing */ #include <iostream> #include <cmath> #include <gsl/gsl_sf_bessel.h> #include <o2scl/ovector_tlate.h> #include <o2scl/multi_funct.h> #include <o2scl/funct.h> #include <o2scl/gsl_anneal.h> #include <o2scl/test_mgr.h> using namespace std; using namespace o2scl; // A simple function with many local minima. A "greedy" minimizer // would likely fail to find the correct minimum. double function(size_t nvar, const ovector_base &x) { double a, b; a=(x[0]-2.0); b=(x[1]+3.0); // This is important to prevent the annealing algorithm to // random walk to infinity since the product of Bessel // functions is flat far from the origin 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[]) { test_mgr t; t.set_output_level(1); cout.setf(ios::scientific); gsl_anneal<> ga; double result; ovector init(2); multi_funct_fptr<> fx(function); ga.ntrial=4000; ga.verbose=1; ga.tol_abs=1.0e-7; ga.T_dec=1.1; // Set a large initial step size double step[1]={10.0}; ga.set_step(1,step); // Choose an initial point at a local minimum away from // the global minimum init[0]=6.0; init[1]=7.0; // Perform the minimization ga.mmin(2,init,result,fx); cout << "x: " << init[0] << " " << init[1] << ", minimum function value: " << result << endl; cout << endl; // Test that it found the global minimum 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"); t.report(); return 0; } // End of example
Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).