#include "paraminhooke.h"
// ********************************************************
// functions for class ParaminHooke
// ********************************************************
ParaminHooke::ParaminHooke(NetInterface* netInt) : ParaminSearch(netInt) {
iters = 0;
returnID = -1;
lambda = 0;
rho = 0.5;
epsilon = 1e-4;
maxiterations = 1000;
type = OPTHOOKE;
// already set in optinfo
// converge = 0;
}
ParaminHooke::~ParaminHooke() {
}
void ParaminHooke::read(CommentStream& infile, char* text) {
int i = 0;
while (!infile.eof() && strcasecmp(text, "seed") && strcasecmp(text, "[hooke]") && strcasecmp(text, "[bfgs]")) {
infile >> ws;
if ((strcasecmp(text, "epsilon") == 0) || (strcasecmp(text, "hookeeps") == 0)) {
infile >> epsilon;
} else if ((strcasecmp(text, "maxiterations") == 0) || (strcasecmp(text, "hookeiter") == 0)) {
infile >> maxiterations;
} else if (strcasecmp(text, "rho") == 0) {
infile >> rho;
} else if (strcasecmp(text, "lambda") == 0) {
infile >> lambda;
} else if (strcasecmp(text, "bndcheck") == 0) {
//JMB - read and ignore bndcheck
infile >> text;
} else {
cerr << "Error while reading optinfo for Hooke - unknown option " << text << endl;
exit(EXIT_FAILURE);
}
infile >> text;
i++;
}
if (i == 0)
cerr << "Warning - no optinfo give for Hooke in file - using default parameters" << endl;
//check the values specified in the optinfo file ...
if ((rho < 0) || (rho > 1)) {
cerr << "\nError in value of rho - setting to default value of 0.5\n";
rho = 0.5;
}
if ((lambda < 0) || (lambda > 1)) {
cerr << "\nError in value of lambda - setting to default value of " << rho << endl;
lambda = rho;
}
}
void ParaminHooke::OptimiseLikelihood() {
double steplength;
int i, keep;
NumberOfHosts = net->getTotalNumProc();
int maxnumindatagroup = numvar*10;
// maybe not needed can do resize on the go????
// Vector tempV(maxnumindatagroup);
//previousf = tempV;
previousf.Reset();
previousf.resize(maxnumindatagroup, 0.0);
//par = new int[maxnumindatagroup];
par.resize(maxnumindatagroup, 0);
xbefore = net->getInitialX();
fbefore = net->getScore();
bestx = xbefore;
bestf = fbefore;
cout << "in opt hooke and Jeeves best x is " << endl;
for (i = 0; i < bestx.Size(); i++) {
cout << bestx[i] << " ";
};
cout << endl;
cout << "in opt hooke and jeeves best f is " << bestf << endl;
delta.Reset();
delta.resize(numvar, 0.0);
int numFromLineSeek;
// change = new int[numvar];
// param = new int[numvar];
change.Reset();
change.resize(numvar, 0);
param.Reset();
param.resize(numvar,0);
// The original definition of the delta array has not been changed, even
// though we're sometimes working with scaled x-values.
for (i = 0; i < numvar; i++) {
delta[i] = fabs(bestx[i] * rho);
if (delta[i] == 0.0)
delta[i] = rho;
param[i] = i;
}
steplength = ((lambda < verysmall) ? rho : lambda);
lineS = new LineSeeker();
while ((iters < maxiterations) && (steplength > epsilon)) {
/* randomize the order of the parameters once in a while, to avoid */
/* the order having an influence on which changes are accepted. */
randomOrder(param);
/* find best new point, one coord at a time */
bestNearby();
//JMB check for too many iterations here
/* if we made some improvements, pursue that direction */
keep = 1;
while ((bestf < fbefore) && (keep == 1) && (iters < maxiterations)) {
// cout << "some improvement from bestNearby" << endl;
for (i = 0; i < numvar; i++) {
/* firstly, arrange the sign of delta[] */
if (bestx[i] <= xbefore[i])
delta[i] = 0.0 - fabs(delta[i]);
else
delta[i] = fabs(delta[i]);
}
numFromLineSeek = lineS->doLineseek(xbefore, bestx, bestf, net);
iters += numFromLineSeek;
xbefore = bestx;
fbefore = bestf;
bestf = lineS->getBestF();
cout << "hooke and jeeves best f " << bestf << endl;
bestx = lineS->getBestX();
// if the further (optimistic) move was bad
if (bestf >= fbefore) {
bestf = fbefore;
bestx = xbefore;
} else {
// else, look around from that point
if (iters < maxiterations) {
xbefore = bestx;
fbefore = bestf;
bestNearby();
/* if the further (optimistic) move was bad */
if (bestf < fbefore) {
// bestf can only be equal or less than fbefore
/* make sure that the differences between the new and the old */
/* points are due to actual displacements - beware of roundoff */
/* errors that might cause newf < fbefore */
keep = 0;
for (i = 0; i < numvar; i++) {
keep = 1;
// AJ changing to be the same check as in gadget..
if (fabs(bestx[i] - xbefore[i]) > rathersmall)
break;
else
keep = 0;
}
}
}
}
}
if ((steplength >= epsilon) && (bestf >= fbefore)) {
steplength = steplength * rho;
for (i = 0; i < numvar; i++)
delta[i] *= rho;
}
}
// Must look into this cout.. things.. should be handle...
cout << "\nStopping Hooke and Jeeves\n\nThe optimisation stopped after " << iters
<< " iterations (max " << maxiterations << ")\nThe steplength was reduced to "
<< steplength << " (min " << epsilon << ")\n";
if (iters >= maxiterations)
cout << "The optimisation stopped because the maximum number of iterations" << "\nwas reached and NOT because an optimum was found for this run\n";
else {
cout << "The optimisation stopped because an optimum was found for this run\n";
converge = 1;
}
delete lineS;
// delete[] change;
// delete[] param;
// delete[] par;
net->setInitialScore(bestx, bestf);
score = bestf;
}
void ParaminHooke::bestNearby() {
double ftmp;
int i, j;
int withinbounds;
int newopt;
net->startNewDataGroup(10 * numvar);
for (i = 0; i < numvar; i++)
change[i] = 0;
numparamset = 0;
i = 0;
// sending as many points as there are processors in total.
// ************
// AJ if never within bound then will set numvar points!!!!
//
while ((i < NumberOfHosts) && (numparamset < numvar)) {
withinbounds = SetPoint(param[numparamset]);
if (withinbounds == 1) {
change[param[numparamset]]++;
i++;
}
numparamset++;
}
// send all available data to all free hosts and then receive them
// back one at a time, sending new points when possible. If point I was
// trying to set was not within bound then nothing has been set and
// nothing will be sent. This could lead to poor use of available hosts.
net->sendToAllIdleHosts();
while (net->getNumNotAns() > 0) {
ReceiveValue();
if (iters % 1000 == 0) {
cout << "\nAfter " << iters << " function evaluations, f(x) = " << bestf << " at\n";
for (j = 0; j < bestx.Size(); j++)
cout << bestx[j] << sep;
cout << endl;
}
if (iters < maxiterations) {
if (MyDataGroup()) {
// Update the optimum if received a better value
newopt = UpdateOpt();
if (!newopt) {
SetPointNearby();
};
}
if (net->allSent())
SetPointNextCoord();
if (!net->allSent())
net->sendToIdleHosts();
} else {
net->receiveAll();
// This takes time, ?? need to do this
// since I am receiving all, maybe should check if any
// of the return values are optimum.
}
}
net->stopUsingDataGroup();
}
int ParaminHooke::SetPoint(int n) {
// return 0 and do nothing if the changes goes out of bounds.
DoubleVector z(bestx);
double next = bestx[n] + delta[n];
int numset;
if (next < lowerbound[n]) {
return 0;
}
else if (next > upperbound[n]) {
return 0;
}
else {
numset = net->getNumDataItemsSet();
z[n] = next;
net->setX(z);
par[numset] = n;
previousf[numset] = bestf;
return 1;
}
}
int ParaminHooke::MyDataGroup() {
return (returnID >= 0);
}
void ParaminHooke::ReceiveValue() {
int receive = net->receiveOne();
if (receive == net->netError()) {
cerr << "Error in Hooke - failed to receive data in linesearch\n";
net->stopUsingDataGroup();
exit(EXIT_FAILURE);
}
returnID = net->getReceiveID();
if (MyDataGroup()) {
iters++;
freceive = net->getY(returnID);
}
}
int ParaminHooke::UpdateOpt() {
int newopt = 0;
if (freceive < bestf) {
bestf = freceive;
bestx = net->getX(returnID);
newopt = 1;
}
return newopt;
}
void ParaminHooke::SetPointNearby() {
int withinbounds = 0;
int returnparam;
returnparam = par[returnID];
if (change[returnparam] == 1) {
if (freceive < previousf[returnID])
withinbounds = SetPoint(returnparam);
else {
delta[returnparam] = 0.0 - delta[returnparam];
withinbounds = SetPoint(returnparam);
if (withinbounds)
change[returnparam]++;
}
}
}
void ParaminHooke::SetPointNextCoord() {
int withinbounds = 0;
while ((numparamset < numvar) && withinbounds == 0) {
withinbounds = SetPoint(param[numparamset]);
if (withinbounds)
change[param[numparamset]]++;
numparamset++;
}
}
void ParaminHooke::Print(ofstream& outfile, int prec) {
outfile << "; Hooke & Jeeves algorithm ran for " << iters
<< " function evaluations\n; and stopped when the likelihood value was "
<< setprecision(prec) << score;
if (converge == -1)
outfile << "\n; because an error occured during the optimisation\n";
else if (converge == 1)
outfile << "\n; because the convergence criteria were met\n";
else
outfile << "\n; because the maximum number of function evaluations was reached\n";
}