![]() |
Object-oriented Scientific Computing Library: Version 0.910
|
One-dimensional minimization is performed by descendants of minimize . There are two one-dimensional minimization algorithms, cern_minimize and gsl_min_brent, 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 minimize 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 minimization is performed by descendants of multi_min: gsl_mmin_simp2, gsl_mmin_conp, gsl_mmin_conf, and gsl_mmin_bfgs2. (The class gsl_mmin_simp2 has been updated with the new "simplex2" method from GSL-1.12. The older "simplex" method is also available in gsl_mmin_simp .) The classes gsl_mmin_simp and gsl_mmin_simp2 do not require any derivative information. The remaining minimization classes 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 a child of mm_funct Finally, the user may specify the automatic gradient object of type 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, the gsl_mmin_conp, gsl_mmin_conf, and gsl_mmin_bfgs2 are probably better. The simplex methods are somewhat robust for sufficiently difficult functions, though restarts are sometimes required to find the correct minimum.
See an example for the usage of the multi-dimensional minimizers in the Multidimensional minimizer example.
Simulated annealing methods are also provided for multi-dimensional minimization (see Simulated Annealing).
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.
The class multi_min_fix 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.
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 <o2scl/test_mgr.h> #include <o2scl/multi_funct.h> #include <o2scl/constants.h> #include <o2scl/gsl_mmin_simp2.h> #include <o2scl/gsl_mmin_conf.h> #include <o2scl/gsl_mmin_conp.h> #include <o2scl/gsl_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 ovector_base &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 << " " << ret << endl; return ret; } // Gradient of the spring function int sgrad(size_t nv, ovector_base &x, ovector_base &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; ovector x(3); double fmin; test_mgr t; t.set_output_level(1); cout.setf(ios::scientific); // Using a member function with \ref ovector objects multi_funct_mfptr<cl> f1(&acl,&cl::spring_two); grad_funct_mfptr<cl> f1g(&acl,&cl::sgrad); gsl_mmin_simp2<> gm1; gsl_mmin_conf<> gm2; gsl_mmin_conp<> gm3; gsl_mmin_bfgs2<> gm4; // 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 << 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 << 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 << 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 << 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 << 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 << 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
This example uses the multi_min_fix class to minimize the function
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 <o2scl/test_mgr.h> #include <o2scl/multi_funct.h> #include <o2scl/gsl_mmin_simp2.h> #include <o2scl/multi_min_fix.h> using namespace std; using namespace o2scl; class cl { public: double mfn(size_t nv, const ovector_base &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; ovector x(3); double fmin; test_mgr t; t.set_output_level(1); cout.setf(ios::scientific); /* Perform the minimization the standard way, with the simplex2 minimizer */ multi_funct_mfptr<cl> f1(&acl,&cl::mfn); gsl_mmin_simp2<> gm1; x[0]=0.5; x[1]=0.5; x[2]=0.5; gm1.mmin(3,x,fmin,f1); cout << gm1.last_ntrial << " iterations." << endl; cout << "Found minimum at: " << x << 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 multi_min_fix object multi_min_fix<bool[3]> gmf; // Create a base minimizer which can be used by the multi_min_fix // 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. gsl_mmin_simp2<multi_funct_mfptr<multi_min_fix<bool[3]> > > gm2; // 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,f1); cout << gmf.last_ntrial << " iterations." << endl; cout << "Found minimum at: " << x << 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,f1); cout << gmf.last_ntrial << " iterations." << endl; cout << "Found minimum at: " << x << 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).