#include <math.h>
#include <stdio.h>
double golden(double (*fname)(double), double ax, double bx, double cx, double tol) {
//GOLDEN Minimize function of one variable using golden section search
//
// [xmin, fmin] = golden(f, ax, bx, cx, tol) computes a local minimum
// of f. xmin is the computed local minimizer of f and fmin is
// f(xmin). xmin is computed to an relative accuracy of TOL.
//
// The parameters ax, bx and cx must satisfy the following conditions:
// ax < bx < cx, f(bx) < f(ax) and f(bx) < f(cx).
//
// xmin satisfies ax < xmin < cx. golden is guaranteed to succeed if f
// is continuous between ax and cx
//
// Roman Geus, ETH Zuerich, 9.12.97
double C, R, x0, x1, x2, x3, f1, f2, xmin, fmin;
long k;
C = (3-sqrt(5))/2;
R = 1-C;
x0 = ax;
x3 = cx;
if (fabs(cx-bx) < fabs(bx-ax)) {
x1 = bx;
x2 = bx + C*(cx-bx);
} else {
x2 = bx;
x1 = bx - C*(bx-ax);
}
f1 = (*fname)(x1);
f2 = (*fname)(x2);
k = 1;
while (fabs(x3-x0) > tol*(fabs(x1)+fabs(x2))) {
//fprintf(1,'k=%4d, |a-b|=%e\n', k, abs(x3-x0));
if (f2 < f1) {
x0 = x1;
x1 = x2;
x2 = R*x1 + C*x3; // x2 = x1+c*(x3-x1)
f1 = f2;
f2 = (*fname)(x2);
} else {
x3 = x2;
x2 = x1;
x1 = R*x2 + C*x0; // x1 = x2+c*(x0-x2)
f2 = f1;
f1 = (*fname)(x1);
}
k = k+1;
printf("%f %f %f\n", x1, x2, x3);
printf("%f %f %f\n", f1, f2, x3);
}
if (f1 < f2) {
xmin = x1;
fmin = f1;
} else {
xmin = x2;
fmin = f2;
}
return xmin;
}