% 30-May-2000
% Script File: ShowNmin
% Minimizes the objective function f using Newton.

close all
nSteps   = 15;      % Total number of steps.
con_size = 20;      % Size of the contour plot
xmin=-1;xmax=4;     % Korkeuskäyräpiirrosalueen
ymin=-5;ymax=5;     %               rajat
xc=[0;2]            % Alkupiste, jos et halua poimia hiirellä, editoi
                    % tähän ja kommentoi hiirihommat pois.

% Display contours of the f function.
% Määrittelemme f:n tässä vektoroidusti. Funktiota fungrad emme voi käyttää
% nätisti meshgrid:n kanssa, vaikka hajoittaisimme sen f- ja grad-osiin, koska
% haluamme yleisen, dimensiosta riippumattoman tavan. Mehgrid taas on raken-
% nettu 2 muuttfkt:lle f(x,y), jossa x ja y ovat (pitkiä) vektoreita.
% Siispä:

f=inline('(x-2).^4+(x-2*y).^2','x','y');

x=linspace(xmin,xmax,con_size);
y=linspace(ymin,ymax,con_size);
[X,Y]=meshgrid(x,y);
Z=f(X,Y);
figure
C = contour(x,y,Z,20);
%C=contour(x,y,Z,0:20);
%title('Enter Contour Labels (Optional). To quit, strike <return>.')
xlabel('x')
ylabel('y')
%clabel(C,'manual')               % Hiirikorkeusnumerot
clabel(C)
% Get ready for the iteration.            % Hiirialkup.
%title('Enter Starting Point for Newton') % Hiirialkup.
%xc = ginput(1);xc=xc'                    % Hiirialkup.
fc=objf(xc);
gc=objg(xc);

xvals = xc
fvals = fc;
gvals = norm(gc);

%hold on
%plot(xvals(1,end),xvals(2,end),'o')


title('')
% Newton iteration
clc
disp('Step       objf            x(1)         x(2)       norm(grad)')
disp('----------------------------------------------------------------')
disp(sprintf('%2.0f     %10.6f     %10.6f   %10.6f      %8.3e ',0,fc,xc,norm(gc)))

for step = 1:nSteps  
    xnew=Nminstep(xc);
    xvals = [xvals xnew];xc = xnew;fc=objf(xnew);gc=norm(objg(xnew));
    fvals = [fvals fc];       
    gvals = [gvals gc]; 
    disp(sprintf('%2.0f     %10.6f     %10.6f   %10.6f      %8.3e ',step,fc,xc,norm(gc)))
end

% Show the Newton path on the contour plot:

figure(1)
hold on
plot(xvals(1,:),xvals(2,:),'r',xvals(1,1),xvals(2,1),'bo')
plot(xvals(1,:),xvals(2,:),'kx')
hold off
title(sprintf('xmin = (%8.3f,%8.3f)   norm(gmin)= %8.4e',xvals(1,end),xvals(2,end),gvals(end)))

% Plot the descent of objf and its gradient:

figure(2)
subplot(2,1,1)
plot(fvals)
title('Value of objf')
xlabel('Iteration')
subplot(2,1,2)
semilogy(gvals)
title('Value of norm(objg).')
xlabel('Iteration')

figure(1);hold on
contour(x,y,Z,0:0.5:10);
