/******************************************************* * Find minimum of a real function Y=F(X) * * using routine for Golden Section Search * * ---------------------------------------------------- * * REFERENCE: "Numerical Recipes, The Art of Scientific * * Computing By W.H. Press, B.P. Flannery, * * S.A. Teukolsky and W.T. Vetterling, * * Cambridge University Press, 1986" * * [BIBLI 08]. * * ---------------------------------------------------- * * Sample run: * * * * Find a minimum of function X*SIN(X)-2*COS(X) * * between X=4 and X=6. * * * * x1=4 x2=5 x3=6 tol=1e-6 * * * * Function minimum is -5.534529 * * for X= 5.232937 * * * * ---------------------------------------------------- * * C++ version by J-P Moreau. * * (www.jpmoreau.fr) * *******************************************************/ #include #include double x1,x2,x3,xmini,ymini,tol; //Function to be analyzed double F(double x) { return x*sin(x)-2*cos(x); } double GOLDEN(double ax,double bx,double cx,double tol,double *xmin) { /*------------------------------------------------------------------ Given a function F, and given a bracketing triplet of fabscissas ax, bx, cx (such that bx is between ax and cx, and F(bx) is less than both F(ax) and F(cx)), this routine performs a golden section search for the minimum, isolating it to a fractional precision of about tol. The fabscissa of the minimum is returned as XMIN, and the minimum function value is returned as GOLDEN, the returned function value. --------------------------------------------------------------------*/ double f0,f1,f2,f3,x0,x1,x2,x3; double R,C; R=0.61803399; C=1.0-R; // golden ratios x0=ax; //At any given time we will keep trace of 4 points: x0,x1,x2,x3 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=F(x1); f2=F(x2); //Initial function evaluations //main loop while (fabs(x3-x0) > tol*(fabs(x1)+fabs(x2))) { if (f2 < f1) { x0=x1; x1=x2; x2=R*x1+C*x3; f0=f1; f1=f2; f2=F(x2); } else { x3=x2; x2=x1; x1=R*x2+C*x0; f3=f2; f2=f1; f1=F(x1); } } if (f1 < f2) { *xmin=x1; return f1; } else { *xmin=x2; return f2; } } // Function Golden() void main() { x1=4.0; x2=5.0; x3=6.0; tol=1e-6; ymini=GOLDEN(x1,x2,x3,tol,&xmini); printf("\n Function minimum is %f\n\n",ymini); printf(" for X= %f\n\n",xmini); } //end of file Golden.cpp