178 |
#include "ecosystem.h" |
#include "ecosystem.h" |
179 |
#include "global.h" |
#include "global.h" |
180 |
#include "seq_optimize_template.h" |
#include "seq_optimize_template.h" |
181 |
#ifdef GADGET_OPENMP |
#ifdef SPECULATIVE |
182 |
#include <omp.h> |
#include <omp.h> |
183 |
#endif |
#endif |
184 |
|
|
304 |
//// trialx[i] = x[i]; |
//// trialx[i] = x[i]; |
305 |
//// } |
//// } |
306 |
// |
// |
|
//// cout << "param:" << param[l] << "-" << trialx[param[l]] << endl; |
|
|
// |
|
307 |
// //Evaluate the function with the trial point trialx and return as -trialf |
// //Evaluate the function with the trial point trialx and return as -trialf |
308 |
// trialf = EcoSystem->SimulateAndUpdate(trialx); |
// trialf = EcoSystem->SimulateAndUpdate(trialx); |
309 |
// trialf = -trialf; |
// trialf = -trialf; |
324 |
// handle.logMessage(LOGINFO, "\nSimulated Annealing finished with a likelihood score of", score); |
// handle.logMessage(LOGINFO, "\nSimulated Annealing finished with a likelihood score of", score); |
325 |
// return; |
// return; |
326 |
// } |
// } |
|
//// cout << "f:" << trialf << endl; |
|
327 |
// //Accept the new point if the new function value better |
// //Accept the new point if the new function value better |
|
//// cout << "mustAccept:" << trialf << "|" << funcval<< endl; |
|
328 |
// if ((trialf - funcval) > verysmall) { |
// if ((trialf - funcval) > verysmall) { |
329 |
// for (i = 0; i < nvars; i++) |
// for (i = 0; i < nvars; i++) |
330 |
// x[i] = trialx[i]; |
// x[i] = trialx[i]; |
348 |
// nrej++; |
// nrej++; |
349 |
// } |
// } |
350 |
// } |
// } |
|
//// cout << "goog VALUE:" << funcval << endl; |
|
351 |
// // JMB added check for really silly values |
// // JMB added check for really silly values |
352 |
// if (isZero(trialf)) { |
// if (isZero(trialf)) { |
353 |
// handle.logMessage(LOGINFO, "Error in Simulated Annealing optimisation after", iters, "function evaluations, f(x) = 0"); |
// handle.logMessage(LOGINFO, "Error in Simulated Annealing optimisation after", iters, "function evaluations, f(x) = 0"); |
436 |
//} |
//} |
437 |
|
|
438 |
|
|
439 |
#ifdef GADGET_OPENMP |
#ifdef SPECULATIVE |
440 |
void OptInfoSimann::OptimiseLikelihoodOMP() { |
void OptInfoSimann::OptimiseLikelihoodOMP() { |
441 |
|
|
442 |
//set initial values |
//set initial values |
446 |
|
|
447 |
double tmp, p, pp, ratio, nsdiv; |
double tmp, p, pp, ratio, nsdiv; |
448 |
double fopt, funcval, trialf; |
double fopt, funcval, trialf; |
449 |
int a, i, j, k, l, offset, quit; |
int a, i, j, k, l, quit; |
450 |
int rchange, rcheck, rnumber; //Used to randomise the order of the parameters |
int rchange, rcheck, rnumber; //Used to randomise the order of the parameters |
451 |
|
|
452 |
|
// store the info of the different threads |
453 |
struct Storage { |
struct Storage { |
454 |
DoubleVector trialx; |
DoubleVector trialx; |
455 |
double newLikelihood; |
double newLikelihood; |
470 |
IntVector nacp(nvars, 0); |
IntVector nacp(nvars, 0); |
471 |
|
|
472 |
EcoSystem->resetVariables(); //JMB need to reset variables in case they have been scaled |
EcoSystem->resetVariables(); //JMB need to reset variables in case they have been scaled |
473 |
if (scale) |
if (scale) { |
474 |
EcoSystem->scaleVariables(); |
EcoSystem->scaleVariables(); |
475 |
|
int numThr = omp_get_max_threads ( ); |
476 |
|
for(i = 0; i < numThr; i++) // scale the variables for the ecosystem of every thread |
477 |
|
EcoSystems[i]->scaleVariables(); |
478 |
|
} |
479 |
EcoSystem->getOptScaledValues(x); |
EcoSystem->getOptScaledValues(x); |
480 |
EcoSystem->getOptLowerBounds(lowerb); |
EcoSystem->getOptLowerBounds(lowerb); |
481 |
EcoSystem->getOptUpperBounds(upperb); |
EcoSystem->getOptUpperBounds(upperb); |
511 |
|
|
512 |
//the function is to be minimised so switch the sign of funcval (and trialf) |
//the function is to be minimised so switch the sign of funcval (and trialf) |
513 |
funcval = -funcval; |
funcval = -funcval; |
|
offset = EcoSystem->getFuncEval(); //number of function evaluations done before loop |
|
514 |
nacc++; |
nacc++; |
515 |
cs /= lratio; //JMB save processing time |
cs /= lratio; //JMB save processing time |
516 |
nsdiv = 1.0 / ns; |
nsdiv = 1.0 / ns; |
529 |
for (i=0; i<numThr; ++i) |
for (i=0; i<numThr; ++i) |
530 |
storage[i].trialx = trialx; |
storage[i].trialx = trialx; |
531 |
|
|
|
//Start the main loop. Note that it terminates if |
|
|
//(i) the algorithm succesfully optimises the function or |
|
|
//(ii) there are too many function evaluations |
|
|
int aux; |
|
532 |
|
|
533 |
DoubleVector vns(nvars, 0); |
int aux; //store the number of evaluations that are not useful for the algorithm |
534 |
|
|
535 |
|
DoubleVector vns(nvars, 0); //vector of ns |
536 |
int ns_ = ceil(numThr/2.); |
int ns_ = ceil(numThr/2.); |
537 |
double res; |
double res; |
538 |
aux=0; |
aux=0; |
539 |
|
//Start the main loop. Note that it terminates if |
540 |
|
//(i) the algorithm succesfully optimises the function or |
541 |
|
//(ii) there are too many function evaluations |
542 |
while (1) { |
while (1) { |
543 |
for (a = 0; a < nt; ++a) { |
for (a = 0; a < nt; ++a) { |
544 |
//Randomize the order of the parameters once in a while, to avoid |
//Randomize the order of the parameters once in a while, to avoid |
740 |
} |
} |
741 |
#endif |
#endif |
742 |
|
|
743 |
|
// calcule a new point |
744 |
void OptInfoSimann::newValue(int nvars, int l, IntVector& param, DoubleVector& trialx, |
void OptInfoSimann::newValue(int nvars, int l, IntVector& param, DoubleVector& trialx, |
745 |
DoubleVector& x, DoubleVector& lowerb, DoubleVector& upperb, DoubleVector& vm) |
DoubleVector& x, DoubleVector& lowerb, DoubleVector& upperb, DoubleVector& vm) |
746 |
{ |
{ |
747 |
int i, k; |
int i, k; |
|
// double old = trialx[param[l]]; |
|
|
// cout << "[" << endl; |
|
748 |
for (i = 0; i < nvars; ++i) { |
for (i = 0; i < nvars; ++i) { |
749 |
if (i == param[l]) { |
if (i == param[l]) { |
750 |
trialx[i] = x[i] + ((randomNumber(&seed) * 2.0) - 1.0) * vm[i]; |
trialx[i] = x[i] + ((randomNumber(&seed) * 2.0) - 1.0) * vm[i]; |
762 |
} |
} |
763 |
} else |
} else |
764 |
trialx[i] = x[i]; |
trialx[i] = x[i]; |
|
// cout <<" " << trialx[i]; |
|
765 |
} |
} |
|
|
|
|
// cout << endl<< "old[" <<param[l] << "]:" << old <<"-" << trialx[param[l]]<< " - " << vm[param[l]]<<endl; |
|
766 |
} |
} |
767 |
|
|
768 |
|
/*################################################################################ |
769 |
|
* code use in the template (seq_optimize_template.h) |
770 |
|
################################################################################*/ |
771 |
|
|
772 |
|
|
773 |
//Generate trialx, the trial value of x |
//Generate trialx, the trial value of x |
774 |
void newValue(int nvars, int l, IntVector& param, DoubleVector& trialx, |
void newValue(int nvars, int l, IntVector& param, DoubleVector& trialx, |
1037 |
EcoSystem->scaleVariables(); |
EcoSystem->scaleVariables(); |
1038 |
#ifndef NO_OPENMP |
#ifndef NO_OPENMP |
1039 |
int numThr = omp_get_max_threads ( ); |
int numThr = omp_get_max_threads ( ); |
1040 |
for(i = 0; i < numThr; i++) |
for(i = 0; i < numThr; i++) // scale the variables for the ecosystem of every thread |
1041 |
EcoSystems[i]->scaleVariables(); |
EcoSystems[i]->scaleVariables(); |
1042 |
#endif |
#endif |
1043 |
} |
} |