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 170 virtual int mmin(
size_t nvar, vec_t &x0,
double &fmin,
174 O2SCL_ERR2(
"Tried to minimize over zero variables ",
184 double E, new_E, best_E, T, old_E;
188 for(j=0;j<nvar;j++) {
204 for(j=0;j<nvar;j++)
old_x[j]=
x[j];
209 for (i=0;i<this->
ntrial;++i) {
210 for (j=0;j<nvar;j++)
new_x[j]=
x[j];
213 new_E=func(nvar,
new_x);
224 for(j=0;j<nvar;j++)
x[j]=
new_x[j];
228 double r=this->
dist();
229 if (r < exp(-(new_E-E)/(
boltz*T))) {
230 for(j=0;j<nvar;j++)
x[j]=
new_x[j];
244 next(nvar,
old_x,old_E,
x,E,T,nmoves,best_E,done);
248 for(j=0;j<nvar;j++) x0[j]=
best_x[j];
255 virtual const char *
type() {
return "anneal_gsl"; }
262 template<
class vec2_t>
int set_step(
size_t nv, vec2_t &stepv) {
265 for(
size_t i=0;i<nv;i++)
step_vec[i]=stepv[i];
291 step_dec=ag.step_dec;
292 min_step_ratio=ag.min_step_ratio;
310 step_dec=ag.step_dec;
311 min_step_ratio=ag.min_step_ratio;
321 #ifndef DOXYGEN_INTERNAL 344 virtual int next(
size_t nvar, vec_t &x_old,
double min_old,
345 vec_t &x_new,
double min_new,
double &T,
346 size_t n_moves,
double best_E,
bool &finished) {
357 for(
size_t i=0;i<nvar;i++) {
367 virtual int start(
size_t nvar,
double &T) {
375 virtual int allocate(
size_t n,
double boltz_factor=1.0) {
386 virtual int step(vec_t &sx,
int nvar) {
387 size_t nstep=step_vec.size();
388 for(
int i=0;i<nvar;i++) {
389 double u=this->
dist();
392 double step_i=step_norm*step_vec[i%nstep];
394 if (step_i<this->
tol_abs*min_step_ratio) {
398 sx[i]=(2.0*u-1.0)*step_i+sx[i];
407 #ifndef DOXYGEN_NO_O2NS int set_step(size_t nv, vec2_t &stepv)
Set the step sizes.
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
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, double best_E, bool &finished)
Determine how to change the minimization for the next iteration.
ubvector best_x
Optimum point over all iterations.
ubvector step_vec
Vector of step sizes.
invalid argument supplied by user
o2scl::prob_dens_uniform dist
The random distribution object.
ubvector new_x
Proposed next point.
int verbose
Output control.
double T_dec
Factor to decrease temperature by (default 1.5)
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...
double boltz
Boltzmann factor (default 1.0).
Simulated annealing base.
double T_start
Initial temperature (default 1.0)
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
virtual int start(size_t nvar, double &T)
Setup initial temperature and stepsize.
double min_step_ratio
Ratio between minimum step size and tol_abs (default 100.0)
ubvector old_x
Last point from previous iteration.
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.
double tol_abs
The independent variable tolerance.
virtual int step(vec_t &sx, int nvar)
Make a step to a new attempted minimum.
Multidimensional minimization by simulated annealing (GSL)
Random number generator (GSL)
virtual int print_iter(size_t nv, vec_t &x, double y, int iter, double tptr, std::string comment)
Print out iteration information.
std::function< double(size_t, const boost::numeric::ublas::vector< double > &)> multi_funct11
Multi-dimensional function typedef.
virtual const char * type()
Return string denoting type ("anneal_gsl")
double step_dec
Factor to decrease step size by (default 1.5)
double step_norm
Normalization for step.
anneal_base< func_t, vec_t, rng_t > & operator=(const anneal_base< func_t, vec_t, rng_t > &ab)
Copy constructor from operator=.
int ntrial
Maximum number of iterations.