43 #ifndef O2SCL_ANNEAL_GSL_H
44 #define O2SCL_ANNEAL_GSL_H
51 #include <boost/numeric/ublas/vector.hpp>
52 #include <boost/numeric/ublas/matrix.hpp>
54 #include <o2scl/anneal.h>
55 #include <o2scl/multi_funct.h>
57 #ifndef DOXYGEN_NO_O2NS
141 class rng_t=int,
class rng_dist_t=rng_gsl>
164 virtual int mmin(
size_t nvar, vec_t &x0,
double &fmin,
168 O2SCL_ERR2(
"Tried to minimize over zero variables ",
176 double E, new_E, best_E, T, old_E;
180 for(j=0;j<nvar;j++) {
196 for(j=0;j<nvar;j++) old_x[j]=x[j];
201 for (i=0;i<this->
ntrial;++i) {
202 for (j=0;j<nvar;j++) new_x[j]=x[j];
205 new_E=func(nvar,new_x);
209 for(j=0;j<nvar;j++) best_x[j]=new_x[j];
216 for(j=0;j<nvar;j++) x[j]=new_x[j];
221 if (r < exp(-(new_E-E)/(
boltz*T))) {
222 for(j=0;j<nvar;j++) x[j]=new_x[j];
231 this->
print_iter(nvar,best_x,best_E,iter,T,
"anneal_gsl");
236 next(nvar,old_x,old_E,x,E,T,nmoves,done);
240 for(j=0;j<nvar;j++) x0[j]=best_x[j];
247 virtual const char *
type() {
return "anneal_gsl"; }
254 template<
class vec2_t>
int set_step(
size_t nv, vec2_t &stepv) {
257 for(
size_t i=0;i<nv;i++)
step_vec[i]=stepv[i];
313 #ifndef DOXYGEN_INTERNAL
319 ubvector x, new_x, best_x, old_x;
326 virtual int next(
size_t nvar, vec_t &x_old,
double min_old,
327 vec_t &x_new,
double min_new,
double &T,
328 size_t n_moves,
bool &finished) {
335 for(
size_t i=0;i<nvar;i++) {
346 virtual int start(
size_t nvar,
double &T) {
354 virtual int allocate(
size_t n,
double boltz_factor=1.0) {
365 virtual int step(vec_t &sx,
int nvar) {
367 for(
int i=0;i<nvar;i++) {
369 sx[i]=(2.0*u-1.0)*
step_vec[i%nstep]+sx[i];
378 #ifndef DOXYGEN_NO_O2NS
virtual int step(vec_t &sx, int nvar)
Make a step to a new attempted minimum.
double step_dec
Factor to decrease step size by (default 1.5)
invalid argument supplied by user
virtual int allocate(size_t n, double boltz_factor=1.0)
Allocate memory for a minimizer over n dimensions with stepsize step and Boltzmann factor boltz_facto...
int verbose
Output control.
virtual int next(size_t nvar, vec_t &x_old, double min_old, vec_t &x_new, double min_new, double &T, size_t n_moves, bool &finished)
Determine how to change the minimization for the next iteration.
virtual const char * type()
Return string denoting type ("anneal_gsl")
ubvector step_vec
Vector of step sizes.
virtual int start(size_t nvar, double &T)
Setup initial temperature and stepsize.
rng_dist_t rng_dist
The default random number distribution.
Simulated annealing base.
double T_dec
Factor to decrease temperature by (default 1.5)
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
rng_t rng
The default random number generator.
int set_step(size_t nv, vec2_t &stepv)
Set the step sizes.
double min_step_ratio
Ratio between minimum step size and tol_abs (default 100.0)
double T_start
Initial temperature (default 1.0)
double tol_abs
The independent variable tolerance.
Multidimensional minimization by simulated annealing (GSL)
anneal_base< func_t, vec_t, rng_t, rng_dist_t > & operator=(const anneal_base< func_t, vec_t, rng_t, rng_dist_t > &ab)
Copy constructor from operator=.
double boltz
Boltzmann factor (default 1.0).
std::function< double(size_t, const boost::numeric::ublas::vector< double > &)> multi_funct11
Multi-dimensional function typedef.
virtual int mmin(size_t nvar, vec_t &x0, double &fmin, func_t &func)
Calculate the minimum fmin of func w.r.t the array x0 of size nvar.
virtual int print_iter(size_t nv, vec_t &x, double y, int iter, double tptr, std::string comment)
Print out iteration information.
int ntrial
Maximum number of iterations.