/* Perform one step of a minimization by double
* bisection. Input required is three points, of
* which the middle one is strictly smaller.
* The points need not be evenly spaced.
*/
void bisect(double (*f)(double), double &x0, double &y0,
double &x2, double &y2,
double &x4, double &y4){
// Bisect the
first interval
double x1 = (x0+x2)/2;
double y1 = f(x1);
// if new point is lower,
// use as new middle
if(y1 
<
y2){
x4 = x2;
y4 = y2;
x2 = x1;
y2 = y1;
return;
}
// if it's higher, move in left limit
if(y1 > y2){
x0 = x1;
y0 = y1;
}
// Bisect the second interval
double x3 = (x2+x4)/2;
double y3 = f(x3);
// if new point is lower,
// use as new middle
if(y3 
<
y2){
x0 = x2;
y0 = y2;
x2 = x3;
y2 = y3;
return;
}
// if it's higher, move in right limit
if(y3 > y2){
x4 = x3;
y4 = y3;
}.
// if y1 was equal to y2, replace y2
if(y1 == y2){
x2 = x1;
y2 = y1;
}
}
/* Solve for minimum of function, using n successive bisections.
* Initial points _MUST_ represent legal configuration, i.e.,
* y0 > y1, y2 > y1.
* The points need not be evenly spaced.
*/
double minimize(double (*f)(double), double &x0, double &y0,
double &x1, double &y1,
double &x2, double &y2, long n){
for(long i=0; i<n; i++)
bisect(f, x0, y0, x1, y1, x2, y2);
return x1;
}
void main(void){
double x0 = 0;
double y0 = f(x0);
double x1 = 0.5;
double y1 = f(x1);
double x2 = 1;
double y2 = f(x2);
cout << minimize(f, x0, y0, x1, y1, x2, y2, 20);
}