function [S,E,G]=golden(f,a,b,delta,epsilon)%Input  - f is the object function input as a string 'f'%       - a and b are the endpoints of the interval%       - delta is the tolerance for the abscissas%       - epsilon is the tolerance for the ordinates%Output - S=(p,yp) contains the abscissa p and%         the ordinate yp of the minimum%       - E=(dp,dy) contains the error bounds for p and yp%       - G is an n x 4 matrix: the kth row contains [ak ck dk bk];%         the values of a, c, d, and b at the kth iteration%  NUMERICAL METHODS: Matlab Programs% (c) 2004 by John H. Mathews and Kurtis D. Fink%  Complementary Software to accompany the textbook:%  NUMERICAL METHODS: Using Matlab, Fourth Edition%  ISBN: 0-13-065248-2%  Prentice-Hall Pub. Inc.%  One Lake Street%  Upper Saddle River, NJ 07458r1=(sqrt(5)-1)/2;r2=r1^2;h=b-a;ya=feval(f,a);yb=feval(f,b);c=a+r2*h;d=a+r1*h;yc=feval(f,c);yd=feval(f,d);k=1;A(k)=a; B(k)=b; C(k)=c; D(k)=d;while (abs(yb-ya)>epsilon)|(h>delta)   k=k+1;   if (yc<yd)      b=d;      yb=yd;      d=c;      yd=yc;      h=b-a;      c=a+r2*h;      yc=feval(f,c);   else      a=c;      ya=yc;      c=d;      yc=yd;      h=b-a;      d=a+r1*h;      yd=feval(f,d);   end   A(k)=a; B(k)=b; C(k)=c; D(k)=d;enddp=abs(b-a);dy=abs(yb-ya);p=a;yp=ya;if (yb<ya)   p=b;   yp=yb;endG=[A' C' D' B'];S=[p yp];E=[dp dy];