function [Xpt,h1,numFcnEval] = SteepestDDescentF(fcn,grad,Xpt,maxIt)
% SteepestDescentF.m
% Steepest descent minimization of a scalar function
% of n variables.
% fcn - string containing name of the scalar function
% grad - string containing the name of the gradient function
% Xpt - initial guess of the minimum (input)
% - approximation of the minimum (output)
% maxIt - maximum number of iterations
% h1 - function value at the minimum
% numFcnEval - total number of function evaluations
% Usage: [X,f,n]=steepestDescentF('CP6p7F','CP6p7DF',[3;1],19)
% Usage: [X,f,n]=steepestDescentF('RosenbrockF','RosenbrockDF',[3;1],1000)
% K. Ming Leung, 09/26/03

format short;
numFcnEval = 0; % number of function evaluations.
h1 = feval(fcn,Xpt);
numFcnEval = numFcnEval+1;

disp(' Iteration x1 x2 f(x1,x2)');
disp([0 Xpt.' h1]);
alpha3 = 1;

for k = 1:maxIt
Step = -feval(grad,Xpt); % Step along negative gradient direction

alpha3 = 0.5; % much better convergence & fewer function calls (37 rather than 98)!
h3 = feval(fcn,Xpt+alpha3*Step);
numFcnEval = numFcnEval+1;

while (h3 >= h1)
alpha3 = 0.5*alpha3;
h3 = feval(fcn,Xpt+alpha3*Step);
numFcnEval = numFcnEval+1;
end
alpha2 = 0.5*alpha3;
h2 = feval(fcn,Xpt+alpha2*Step);
numFcnEval = numFcnEval+1;

% quadratic interpolation
c2 = (h2-h1)/alpha2;
c22 = (h3-h1)/alpha3;
c3 = (c22-c2)/(alpha3-alpha2);
alphaTemp = 0.5*(alpha2-c2/c3);

% find min. of quadratic in [alpha1, alpha3]
if alphaTemp > 0
alphaMin = alphaTemp;
else
alphaMin = alpha3;
end

% use the min. to take a step
Xpt = Xpt+alphaMin*Step;
h1 = feval(fcn,Xpt);
numFcnEval = numFcnEval+1;

disp([k Xpt.' h1]);

end