Minimization

One-dimensional minimizers

One-dimensional minimization is performed by descendants of o2scl::min_base . There are two one-dimensional minimization algorithms, o2scl::min_cern and o2scl::min_brent_gsl, and they are both bracketing algorithms type where an interval and an initial guess must be provided. If only an initial guess and no bracket is given, these two classes will attempt to find a suitable bracket from the initial guess. While the o2scl::min_base base class is designed to allow future descendants to optionally use derivative information, this is not yet supported for any one-dimensional minimizers.

Multi-dimensional minimizers

Multi-dimensional minimization is performed by descendants of o2scl::mmin_base. O2scl includes a simplex minimizer (o2scl::mmin_simp2), traditional minimizers which use gradient information (o2scl::mmin_conp, o2scl::mmin_conf, and o2scl::mmin_bfgs2), and differential evolution minimizers (o2scl::diff_evo and o2scl::diff_evo_adapt). Minimization by simulated annealing is included and described in the Simulated Annealing section. Constrained minimization is also included and described in separately in Constrained Minimization.

See an example for the usage of the multi-dimensional minimizers in the Multidimensional minimizer example below.

Simplex minimizer

The class o2scl::mmin_simp2 minimizes a function using the Nelder-Mead simplex algorithm and does not require any information about the gradient. It is based on GSL and has been updated with the new "simplex2" method from GSL-1.12.

Restarts of the simplex minimizer are sometimes required to find the correct minimum, and o2scl::mmin_simp2 can get caught in infinite loops, especially for functions which have symmetries directly related to one or more of the parameters.

Traditional minimizers with gradient information

Classes o2scl::mmin_conp, o2scl::mmin_conf, and o2scl::mmin_bfgs2 are intended for use when the gradient of the function is available, but they can also automaticallly compute the gradient numerically. The standard way to provide the gradient is to use an object of type o2scl::grad_funct11 . The user may specify the automatic gradient object of type o2scl::gradient which is used by the minimizer to compute the gradient numerically when a function is not specified.

Generally, when a closed expression is available for the gradient, o2scl::mmin_bfgs2 is likely the best choice. However, the simplex methods can be more robust for sufficiently difficult functions.

It is important to note that not all of the minimization routines test the second derivative to ensure that it doesn't vanish to ensure that we have indeed found a true minimum.

Fixing Parameters

The class o2scl::mmin_fix_params provides a convenient way of fixing some of the parameters and minimizing over others, without requiring a the function interface to be rewritten. An example is given in the Minimizer fixing variables example below.

Multidimensional minimizer example

This example uses the O2scl minimizers based on GSL to minimize a rather complicated three-dimensional function which has constant level surfaces which look like springs oriented along the z-axis. This example function, originally created here for O2scl , was added later to the GSL library minimization test functions.

/* Example: ex_mmin.cpp
-------------------------------------------------------------------
Example usage of the multidimensional minimizers with and without
gradients.
*/
#include <fstream>
#include <string>
#include <cmath>
#include <boost/numeric/ublas/vector.hpp>
#include <o2scl/test_mgr.h>
#include <o2scl/multi_funct.h>
#include <o2scl/constants.h>
#include <o2scl/mmin_simp2.h>
#include <o2scl/mmin_conf.h>
#include <o2scl/mmin_conp.h>
#include <o2scl/mmin_bfgs2.h>
using namespace std;
using namespace o2scl;
class cl {
public:
cl() {
param=30.0;
}
// To output function evaluations to a file
ofstream fout;
// Parameter of the quadratic
double param;
// Updated spring function
double spring_two(size_t nv, const ubvector &x) {
double theta=atan2(x[1],x[0]);
double r=sqrt(x[0]*x[0]+x[1]*x[1]);
double z=x[2];
double tmz=theta-z;
double fact=8.0-pow(sin(tmz+o2scl_const::pi/2.0)+1.0,3.0);
double rm1=r-1.0;
double ret=fact+exp(rm1*rm1)+z*z/param;
fout << x[0] << " " << x[1] << " " << x[2] << " " << ret << endl;
return ret;
}
// Gradient of the spring function
int sgrad(size_t nv, ubvector &x, ubvector &g) {
double theta=atan2(x[1],x[0]);
double r=sqrt(x[0]*x[0]+x[1]*x[1]);
double z=x[2];
double tmz=theta-z;
double rm1=r-1.0;
double fact=8.0-pow(sin(tmz+o2scl_const::pi/2.0)+1.0,3.0);
double dtdx=-x[1]/r/r;
double dtdy=x[0]/r/r;
double drdx=x[0]/r;
double drdy=x[1]/r;
double dfdt=-3.0*pow(sin(tmz+o2scl_const::pi/2.0)+1.0,2.0)*
cos(tmz+o2scl_const::pi/2.0);
double dfdz=2.0*z/param+3.0*pow(sin(tmz+o2scl_const::pi/2.0)+1.0,2.0)*
cos(tmz+o2scl_const::pi/2.0);
double dfdr=2.0*rm1*exp(rm1*rm1);
g[0]=dfdr*drdx+dfdt*dtdx;
g[1]=dfdr*drdy+dfdt*dtdy;
g[2]=dfdz;
return 0;
}
};
int main(void) {
cl acl;
ubvector x(3);
double fmin;
cout.setf(ios::scientific);
// Using a member function
#ifdef O2SCL_NO_CPP11
multi_funct_mfptr<cl> f1(&acl,&cl::spring_two);
grad_funct_mfptr<cl> f1g(&acl,&cl::sgrad);
#else
multi_funct11 f1=std::bind
(std::mem_fn<double(size_t,const ubvector &)>(&cl::spring_two),
&acl,std::placeholders::_1,std::placeholders::_2);
grad_funct11 f1g=std::bind
(std::mem_fn<int(size_t,ubvector &,ubvector &)>(&cl::sgrad),
&acl,std::placeholders::_1,std::placeholders::_2,
std::placeholders::_3);
#endif
// This function is difficult to minimize, so more trials
// are required.
gm1.ntrial*=10;
gm2.ntrial*=10;
gm3.ntrial*=10;
gm4.ntrial*=10;
// Simplex minimization
acl.fout.open("ex_mmin1.dat");
x[0]=1.0;
x[1]=1.0;
x[2]=7.0*o2scl_const::pi;
gm1.mmin(3,x,fmin,f1);
acl.fout.close();
cout << gm1.last_ntrial << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],1.0,1.0e-4,"1a");
t.test_rel(x[1],0.0,1.0e-4,"1b");
t.test_rel(x[2],0.0,1.0e-4,"1c");
// Fletcher-Reeves conjugate
acl.fout.open("ex_mmin2.dat");
x[0]=1.0;
x[1]=0.0;
x[2]=7.0*o2scl_const::pi;
gm2.mmin(3,x,fmin,f1);
acl.fout.close();
cout << gm2.last_ntrial << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],1.0,4.0e-3,"2a");
t.test_rel(x[1],0.0,4.0e-3,"2b");
t.test_rel(x[2],0.0,4.0e-3,"2c");
// Fletcher-Reeves conjugate with gradients
acl.fout.open("ex_mmin2g.dat");
x[0]=1.0;
x[1]=0.0;
x[2]=7.0*o2scl_const::pi;
gm2.mmin_de(3,x,fmin,f1,f1g);
acl.fout.close();
cout << gm2.last_ntrial << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],1.0,4.0e-3,"2a");
t.test_rel(x[1],0.0,4.0e-3,"2b");
t.test_rel(x[2],0.0,4.0e-3,"2c");
// Polak-Ribere conjugate
acl.fout.open("ex_mmin3.dat");
x[0]=1.0;
x[1]=0.0;
x[2]=7.0*o2scl_const::pi;
gm3.mmin(3,x,fmin,f1);
acl.fout.close();
cout << gm3.last_ntrial << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],1.0,4.0e-3,"3a");
t.test_rel(x[1],0.0,4.0e-3,"3b");
t.test_rel(x[2],0.0,4.0e-3,"3c");
// Polak-Ribere conjugate with gradients
acl.fout.open("ex_mmin3g.dat");
x[0]=1.0;
x[1]=0.0;
x[2]=7.0*o2scl_const::pi;
gm3.mmin_de(3,x,fmin,f1,f1g);
acl.fout.close();
cout << gm3.last_ntrial << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],1.0,4.0e-3,"3a");
t.test_rel(x[1],0.0,4.0e-3,"3b");
t.test_rel(x[2],0.0,4.0e-3,"3c");
// BFGS method
// BFGS has trouble converging (especially to zero, since the
// minimimum of x[0] is exactly at zero) if the derivative is not
// very accurate.
gm4.def_grad.epsrel=1.0e-8;
gm4.err_nonconv=false;
acl.fout.open("ex_mmin4.dat");
x[0]=1.0;
x[1]=0.0;
x[2]=7.0*o2scl_const::pi;
gm4.mmin(3,x,fmin,f1);
acl.fout.close();
cout << gm4.last_ntrial << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],1.0,1.0e-4,"4a");
t.test_rel(x[1],0.0,1.0e-4,"4b");
t.test_rel(x[2],0.0,1.0e-4,"4c");
t.report();
return 0;
}
// End of example

Idea for Future:
Plot the spring function

Minimizer fixing variables

This example uses the o2scl::mmin_fix_params class to minimize the function

\[ y=\left(x_0-2\right)^2+\left(x_1-1\right)^2+x_2^2 \]

while fixing some of the parameters.

/* Example: ex_mmin_fix.cpp
-------------------------------------------------------------------
Example usage of the mmin_fix class, which fixes some of the
paramters for a multidimensional minimization.
*/
#include <cmath>
#include <boost/numeric/ublas/vector.hpp>
#include <o2scl/test_mgr.h>
#include <o2scl/multi_funct.h>
#include <o2scl/mmin_simp2.h>
#include <o2scl/mmin_fix.h>
using namespace std;
using namespace o2scl;
class cl {
public:
double mfn(size_t nv, const ubvector &x) {
return (x[0]-2.0)*(x[0]-2.0)+(x[1]-1.0)*(x[1]-1.0)+x[2]*x[2];
}
};
int main(void) {
cl acl;
ubvector x(3);
double fmin;
cout.setf(ios::scientific);
/*
Perform the minimization the standard way, with the
simplex2 minimizer
*/
std::bind(std::mem_fn<double(size_t,const ubvector &)>(&cl::mfn),
acl,std::placeholders::_1,std::placeholders::_2);
x[0]=0.5;
x[1]=0.5;
x[2]=0.5;
gm1.mmin(3,x,fmin,f1c11);
cout << gm1.last_ntrial << " iterations." << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],2.0,1.0e-4,"1a");
t.test_rel(x[1],1.0,1.0e-4,"1b");
t.test_rel(x[2],0.0,1.0e-4,"1c");
// Create a new mmin_fix_params object
// Create a base minimizer which can be used by the mmin_fix_params
// object. Note that we can't use 'gm1' here, because it has a
// different type than 'gm2', even though its functionality is
// effectively the same.
// Set the base minimizer
gmf.set_mmin(gm2);
/*
First perform the minimization as above.
*/
x[0]=0.5;
x[1]=0.5;
x[2]=0.5;
gmf.mmin(3,x,fmin,f1c11);
cout << gmf.last_ntrial << " iterations." << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],2.0,1.0e-4,"2a");
t.test_rel(x[1],1.0,1.0e-4,"2b");
t.test_rel(x[2],0.0,1.0e-4,"2c");
/*
Now fix the 2nd variable, and re-minimize.
*/
bool fix[3]={false,true,false};
x[0]=0.5;
x[1]=0.5;
x[2]=0.5;
gmf.mmin_fix(3,x,fmin,fix,f1c11);
cout << gmf.last_ntrial << " iterations." << endl;
cout << "Found minimum at: "
<< x[0] << " " << x[1] << " " << x[2] << endl;
t.test_rel(x[0],2.0,1.0e-4,"3a");
t.test_rel(x[1],0.5,1.0e-4,"3b");
t.test_rel(x[2],0.0,1.0e-4,"3c");
t.report();
return 0;
}
// End of example

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).