Listing 5: The state machine
Listing 5: The state machine
/* Find the minimum of a function
* using a second-order interpolation
*/
double minimize(double (*f)(double), double X0,
double X1, long Steps, double eps){
// local variables
long iter = 60;
double xx, yy; // test point
double dx0, dx1; // space between points
int newmin = FALSE;
// qualify inputs
if(X1 <= X0){
cout << "minimize: Search range must be finite and positive\n";
return 0;
}
if(Steps < 2){
cout << "minimize: Must have at least two intervals\n";
return 0;
}
state = step(f, X0, X1, Steps);
// outer loop of state machine
while(--iter > 0){
switch(state){
case LEFT_EDGE:
return X0;
case RIGHT_EDGE:
return X1;
case HAVE_TRIPLET:
state = BISECT_LEFT;
break;
case BISECT_LEFT:
xx = bisect_left(index);
yy = f(xx);
if(yy < y[index]){
if(index == 3)
insert_left(xx, yy, 3);
else
insert_right(xx, yy, index);
}
else if(yy > y[index]){
if(index == 1){
insert_right(xx, yy, 1);
index = 2;
}
else
insert_left(xx, yy, index-1);
}
else{
x[index] = xx;
y[index] = yy;
}
state = BISECT_RIGHT;
break;
case BISECT_RIGHT:
xx = bisect_right(index);
yy = f(xx);
if(yy < y[index]){
if(index == 1)
insert_right(xx, yy, 1);
else
insert_left(xx, yy, index);
}
else if(yy > y[index]){
if(index == 3){
insert_left(xx, yy, 3);
index = 2;
}
else
insert_right(xx, yy, index+1);
}
else{
x[index] = xx;
y[index] = yy;
}
if((x[index+1]-x[index-1]) < eps)
state = DONE;
else
state = BISECT_LEFT;
break;
case PARA_FIT:
break;
case DONE:
return x[index];
}
}
cout << "minimize failed to converge" << endl;
return x[index];
}
|
