#ifndef SEQ_OPTIMIZE_TEMPLATE_H
#define SEQ_OPTIMIZE_TEMPLATE_H
#include <iostream>
#include <cstdlib>
#include <sys/time.h>
#ifdef _OPENMP
#include <omp.h>
#endif
class Siman {
int _nvars;
unsigned _seed;
unsigned _seedM;
unsigned _seedP;
int _l;
int _nt;
int _NT;
int _ns;
int _NS;
IntVector _param;
DoubleVector* _x;
DoubleVector* _lowerb;
DoubleVector* _upperb;
DoubleVector _vm;
double _t;
double _rt;
IntVector _nacp;
int _nacc; //The number of accepted function evaluations
int _nrej; //The number of rejected function evaluations
int _naccmet; //The number of metropolis accepted function evaluations
double _nsdiv;
int _tempcheck;
double _simanneps;
DoubleVector _fstar;
double _lratio;
double _uratio;
double _cs;
DoubleVector* _bestx;
int _scale;
int* _converge; //NN
double* _score; //NN
public:
//Required Default constructor : Params()
Siman(){
_nvars = 0;
_seed = 1234;
_seedM = 1234;
_seedP = 1234;
_l = 0;
_nt = 0;
_NT = 0;
_ns = 0;
_NS = 0;
_param = IntVector();
_x = NULL;
_lowerb = NULL;
_upperb = NULL;
_vm = DoubleVector();
_t = 0.;
_rt = 0.;
_nacp = IntVector();
_nacc = 0; //The number of accepted function evaluations
_nrej = 0; //The number of rejected function evaluations
_naccmet = 0; //The number of metropolis accepted function evaluations
_nsdiv = 0;
_tempcheck = 0;
_simanneps = 0;
_fstar = DoubleVector();
_lratio = 0.;
_uratio = 0.;
_cs = 0.;
_bestx = NULL;
_scale = 0;
_converge =NULL;
_score = 0;
}
Siman(unsigned &seed, unsigned &seedM, unsigned &seedP, int &nvars, int &nt, int &ns, IntVector param, DoubleVector* x,
DoubleVector* lowerb, DoubleVector* upperb, DoubleVector vm, double t, double &rt,
double nsdiv, int &tempcheck, double &simanneps, DoubleVector fstar,
double &lratio, double &uratio, double &cs, DoubleVector* bestx, int &scale, int* converge,
double* score):
_seed(seed), _seedM(seedM), _seedP(seedP), _l(0), _nvars(nvars), _nt(nt), _NT(0), _ns(ns), _NS(0), _t(t),
_rt(rt), _nsdiv(nsdiv), _nacc(1), _tempcheck(tempcheck), _simanneps(simanneps), _lratio(lratio),
_uratio(uratio), _cs(cs), _scale(scale){
_param = IntVector(param);
_x = x;
_lowerb = lowerb;
_upperb = upperb;
_vm = DoubleVector(vm);
_nacp = IntVector(_nvars, 0);
_nrej = 0;
_naccmet = 0;
randomizeParams();
_fstar = DoubleVector(fstar);
_bestx = bestx;
_converge = converge;
_score = score;
}
//Required Copy constructor : Params(const Params& other)
Siman(const Siman& other): _nvars(other._nvars), _seed(other._seed), _seedM(other._seedM), _seedP(other._seedP),
_l(other._l), _nt(other._nt), _NT(other._NT), _ns(other._ns), _NS(other._NS), _param(IntVector(other._param)),
_x(other._x), _lowerb(other._lowerb),_upperb(other._upperb), _vm(DoubleVector(other._vm)), _t(other._t),
_rt(other._rt), _nacp(IntVector(other._nacp)), _nacc(other._nacc), _nrej(other._nrej),
_naccmet(other._naccmet), _nsdiv(other._nsdiv), _tempcheck(other._tempcheck), _simanneps(other._simanneps),
_fstar(DoubleVector(other._fstar)), _lratio(other._lratio), _uratio(other._uratio), _cs(other._cs),
_bestx(other._bestx), _scale(other._scale), _converge(other._converge), _score(other._score){}
//Required operator= : Params& operator= (const Params& other)
Siman& operator= (const Siman& other){
_nvars = other._nvars;
_seed = other._seed;
_seedM = other._seedM;
_seedP = other._seedP;
_l = other._l;
_nt = other._nt;
_NT = other._NT;
_ns = other._ns;
_NS = other._NS;
_param = IntVector(other._param);
_x = other._x;
_lowerb = other._lowerb;
_upperb = other._upperb;
_vm = DoubleVector(other._vm);
_t = other._t;
_rt = other._rt;
_nacp = IntVector(other._nacp);
_nacc = other._nacc; //The number of accepted function evaluations
_nrej = other._nrej; //The number of rejected function evaluations
_naccmet = other._naccmet; //The number of metropolis accepted function evaluations
_nsdiv = other._nsdiv;
_tempcheck = other._tempcheck;
_simanneps = other._simanneps;
_fstar = DoubleVector(other._fstar);
_lratio = other._lratio;
_uratio = other._uratio;
_cs = other._cs;
_bestx = other._bestx;
_scale = other._scale;
_converge =other._converge;
_score = other._score;
return *this;
}
//Destructor : ~Params()
~Siman(){
}
unsigned* getSeed(){return &_seed;}
unsigned* getSeedM(){return &_seedM;}
unsigned getSeedM_(){return _seedM;}
void setSeedM(unsigned val){_seedM = val;}
int getNvars(){return _nvars;}
int getL(){return _l;}
void incrementL() { _l = (_l +1) % _nvars;}
IntVector& getParam(){return _param;}
DoubleVector& getX(){return *_x;}
void setX(DoubleVector* point) {_x=point;}
DoubleVector& getLowerb(){return *_lowerb;}
DoubleVector& getUpperb(){return *_upperb;}
DoubleVector& getVm(){return _vm;}
void setVm(DoubleVector vm){_vm = vm;}
DoubleVector& getFstar(){return _fstar;}
double getT(){return _t;}
void setT(double val){_t = val;}
void incrementNacc() {_nacc++;}
void incrementNrej(){_nrej++;}
void incrementNaccmet(){_naccmet++;}
void incrementNacp(int i){_nacp[i]++;}
void randomizeParams(){
int rnumber, rcheck, i, rchange = 0;
while (rchange < _nvars) {
rnumber = rand_r(&_seedP) % _nvars;
rcheck = 1;
for (i = 0; i < rchange; i++)
if ((_param)[i] == rnumber)
rcheck = 0;
if (rcheck) {
(_param)[rchange] = rnumber;
rchange++;
}
}
}
int getNs() const {return _ns;}
int getNS() const {return _NS;}
void setNS(int ns = 0) {_NS = ns;}
void incrementNS() {_NS++;}
int getNt() const {return _nt;}
int getNT() const {return _NT;}
void setNT(int nt) {_NT = nt;}
void incrementNT() {_NT++;};
double getSimanneps() const {return _simanneps;}
int getTempcheck() const {return _tempcheck;}
const int getNacp(int i) const {return _nacp[i];}
void setNacp(int i, int val) {_nacp[i] = val;}
double getCs() const {return _cs;}
double getLratio() const {return _lratio;}
double getUratio() const {return _uratio;}
void setConverge(int val) {(*_converge) = val;}
DoubleVector* getBestx() {return _bestx;}
int getScale() const {return _scale;}
double getRt() const {return _rt;}
double getNsdiv() const {return _nsdiv;}
int* getConverge() const {return _converge;}
int getNacc() const {return _nacc;}
int getNaccmet() const {return _naccmet;}
int getNrej() const {return _nrej;}
void setNrej(int nrej) {_nrej = nrej;};
double* getScore() const {return _score;}
};
/* IMPROVEMENTS:
- OpenMP parallel implementation added
- Can be compiled with OpenMP using macro _OPENMP
- Memory for items is dynamically reused in the parallel searches
- Items whose evaluation has not begun when a new optimum is found are not evaluated (only possible in TBB search)
*/
/// Performs a reproducible optimization process, sequentially in seq_opt or in parallel in paral_opt
template<class SEED_T, class PARAMS_T, class CONTROL_T, double Evaluator(
const PARAMS_T&), void BuildNewParams(SEED_T&, PARAMS_T&)>
class ReproducibleSearch {
const SEED_T initialSeed_;
const PARAMS_T initialParams_;
const int maxIters_;
SEED_T seed_;
PARAMS_T p_;
int iters_;
double optimumValue_;
double funcval_;
double searchTime_;
/// Stores an individual item to evaluate and the data to restore the search state
/// if it is chosen as new optimum item
struct ParalEvalItem_t {
int order_;
PARAMS_T p_;
SEED_T seed_; ///< seed AFTER generation of p_
double value_;
ParalEvalItem_t(int order, const PARAMS_T& p, const SEED_T& seed,
double value = 0.) :
order_(order), p_(p), seed_(seed), value_(value) {
}
ParalEvalItem_t() :
value_(0.) {
}
void evaluateItem() {
value_ = Evaluator(p_);
}
void print() const {
std::cout << order_ << ' ' << p_ << " v= " << value_ << std::endl;
}
};
/// Internal storage for search items with a padding field to avoid false sharing
struct ParalEvalItemStorage_t {
ParalEvalItem_t pei_;
char padding_[128];
};
public:
ReproducibleSearch(const SEED_T& seed, const PARAMS_T& params, int max_iters = 200) :
initialSeed_(seed), initialParams_(params), maxIters_(max_iters) {}
// Set the problem in its initial state
void reset() {
seed_ = initialSeed_;
p_ = initialParams_;
std::cout << ": " << p_ << " v= " << Evaluator(p_) << std::endl;
}
//########################################################
void seq_opt(double funcval) {
struct timeval t0, t1, t;
reset();
CONTROL_T control_object;
gettimeofday(&t0, NULL);
seed_.setX(&p_);
funcval_ = funcval;
optimumValue_ = funcval_;
bool quit = false;
for (iters_ = 1; iters_ <= maxIters_; ++iters_) {
PARAMS_T p = p_;
BuildNewParams(seed_, p);
const double prev_val = funcval_;
const double val = Evaluator(p);
if (control_object.mustAccept(prev_val, val, seed_, iters_)) {
int i;
for (i = 0; i < seed_.getNvars(); i++){
p_[i] = p[i];
}
funcval_ = val;
} else if ((*seed_.getConverge()) == -1) return;
control_object.optimum(val, optimumValue_, iters_, p_, initialParams_, seed_);
if ((quit=control_object.mustTerminate(optimumValue_, funcval_, seed_, iters_))) {
break;
}
}
control_object.printResult(quit, seed_, iters_);
gettimeofday(&t1, NULL);
timersub(&t1, &t0, &t);
searchTime_ = t.tv_sec + t.tv_usec / 1000000.0;
}
void paral_opt_omp(double funcval, int nthreads, int paral_tokens = -1) {
struct timeval t0, t1, t;
reset();
CONTROL_T control_object;
gettimeofday(&t0, NULL);
#if defined(_OPENMP)
omp_set_num_threads(nthreads);
#endif
if (paral_tokens < 0) {
paral_tokens = nthreads;
}
seed_.setX(&p_);
funcval_ = funcval;
optimumValue_ = funcval_;
bool quit = false;
ParalEvalItemStorage_t * const storage =
new ParalEvalItemStorage_t[paral_tokens]; //Guarda los datos de los diferentes threads
int old_paral_tokens = paral_tokens;
for (iters_ = 1; iters_ <= maxIters_;) {
if ((seed_.getL() + paral_tokens -1) >= seed_.getNvars())
paral_tokens = seed_.getNvars() - seed_.getL();
if ((iters_ + paral_tokens - 1) > maxIters_) {
paral_tokens = maxIters_ - (iters_ - 1);
}
for (int i = 0; i < paral_tokens; i++) {
storage[i].pei_.p_ = p_;
BuildNewParams(seed_, storage[i].pei_.p_);
storage[i].pei_.seed_ = seed_;
seed_.incrementL();
}
const double prev_val = funcval_;
#pragma omp parallel for
for (int i = 0; i < paral_tokens; i++) {
storage[i].pei_.evaluateItem();
}
int j, i=0; //Notice that we need it after the loop
if (control_object.mustAccept(prev_val, storage[i].pei_.value_, storage[i].pei_.seed_, iters_+i)) {
for (j = 0; j < seed_.getNvars(); j++)
p_[j] = storage[i].pei_.p_[j];
funcval_ = storage[i].pei_.value_;
} else{
if ((*storage[i].pei_.seed_.getConverge()) == -1) return;
for (i = 1; i < paral_tokens; i++) {
storage[i].pei_.seed_.setSeedM(storage[i-1].pei_.seed_.getSeedM_());
if (control_object.mustAccept(prev_val, storage[i].pei_.value_, storage[i].pei_.seed_, iters_+i)) {
for (j = 0; j < seed_.getNvars(); j++)
p_[j] = storage[i].pei_.p_[j];
funcval_ = storage[i].pei_.value_;
break;
} else if ((*storage[i].pei_.seed_.getConverge()) == -1) return;
}
}
i = (i == paral_tokens) ? (i - 1) : i;
seed_ = storage[i].pei_.seed_;
seed_.setNrej(seed_.getNrej()+i);
control_object.optimum(storage[i].pei_.value_, optimumValue_, iters_,
storage[i].pei_.p_, initialParams_, seed_);
iters_ += i+1;
if ((quit=control_object.mustTerminate(optimumValue_, funcval_, seed_, iters_))) {
break;
}
paral_tokens = old_paral_tokens;
}
control_object.printResult(quit, seed_ ,iters_);
gettimeofday(&t1, NULL);
timersub(&t1, &t0, &t);
searchTime_ = t.tv_sec + t.tv_usec / 1000000.0;
delete[] storage;
}
/// Only used by sequential search when verbose
void print(int iter, double val) const {
ParalEvalItem_t v(iter, p_, seed_, val);
v.print();
}
const PARAMS_T& getResult() const {
return p_;
}
double getOptimumValue() const {
return funcval_;
}
double getSearchTime() const {
return searchTime_;
}
int iterations() const {
return iters_;
}
};
#endif