echo off;
% NUMERICAL METHODS: MATLAB Programs, (c) John H. Mathews 1994
% To accompany the text:
% NUMERICAL METHODS for Mathematics, Science and Engineering, 2nd Ed, 1992
% Prentice Hall, Englewood Cliffs, New Jersey, 07632, U.S.A.
% This free software is complements of the author.

% Algorithm 8.2 (Nelder-Mead's Minimization Method).
% Section 8.1, Minimization of a Function, Page 414
echo on; clc; format short; hold off; clear

% This program finds the minimum of a function Fn(X),
% where X is an n-dimensional vector.
% The technique is the Nelder-Mead method for n-dimensions.
% The function Fn(X) is to be placed in two M-files
% The function f.m is for graphing and Fn.m is for computing.

% function z = f(x,y)
% z = x.^2 - 4.*x + y.^2 - y - x.*y;

% function z = Fn(X)
% x = X(1); y = X(2);
% z = x^2 - 4*x + y^2 - y - x*y;

delete f.m
diary f.m; disp('function z = f(x,y)');...
disp('z = x.^2 - 4.*x + y.^2 - y - x.*y;');...
diary off;

delete Fn.m
diary Fn.m; disp('function z = Fn(X)');...
disp('x = X(1); y = X(2);');...
disp('z = x^2 - 4*x + y^2 - y - x*y;');...
diary off;

% Remark. f.m and Fn.m nelder.m are used for Algorithm 8.2
% Remark. f(x,y) is used only to graph 2-dimensional functions.
f(0,0); Fn([0 0]);
pause % Press any key to continue.

clc;
% The example is 2-dimensional, hence a surface z = f(x,y)

% can be graphed. (Omit this section for higher dimensions.)

% The endpoints of the rectangle [a,b] x [c,d] are a,b,c,d.

% Place the minimum and maximum number of iterations in

% min1 and max1, respectively.

% Place the tolerance in epsilon.

a = 0;
b = 4;
c = 0;
d = 3;
min1 = 8;
max1 = 80;
epsilon = 1e-6;

pause % Press any key to graph the 2-dim function.

clc; clg;
hs = (b-a)/16;
ks = (d-c)/16;
[Xs,Ys] = meshdom(a:hs:b, c:ks:d);
Zs = f(Xs,Ys);
mesh(Zs);...
title('To search for the minimum of z = f(x,y).');...
shg; pause % Press any key to enter the starting vertices.

clc;

% The Nelder-Mead simplex method or "polytope method is used to

% find the minimum of the n-dim function Fn(X), for convenience,

% Fn(X) can be expressed using X = (x1,x2,x3,x4,x5,x6) and the

% variables x = x1 , y = x2 , z = x3 , u = x4 , v = x5 , w = x6

% You must supply n+1 linearly independent starting points:

% For n=2, supply: V(1,1:n), V(2,1:n), V(3,1:n)

% For n=3, supply: V(1,1:n), V(2,1:n), V(3,1:n), V(4,1:n)

n = 2;
V(1,1:n) = [0.0 0.0];
V(2,1:n) = [1.2 0.0];
V(3,1:n) = [0.0 0.8];

pause % Press any key to continue.

clc; clg;
a0 = 0;
b0 = 4;
c0 = 0;
d0 = 3;
hold off;
XS = V(1:n+1,1)'; XSL = [XS,XS(1)];
YS = V(1:n+1,2)'; YSL = [YS,YS(1)];
axis([a0 b0 c0 d0]);...
plot(XS,YS,'or',XSL,YSL,'-g');...
hold on;...
plot([a0 b0],[0 0],'b',[0 0],[c0 d0],'b');...
xlabel('x');...
ylabel('y');...
title('The simplex vertices for Nelder-Mead`s method.');...
grid;...
shg; pause % Press any key to continue.

[V0,y0,dV,dy,P,Q] = nelder('Fn',V,min1,max1,epsilon,1);

pause % Press any key to find the minimum of Fn(X).

clc;
format long;
Mx1 = 'The coordinates of the point P are:';
Mx2 = 'Error bound for the coordinates of P:';
Mx3 = 'The minimum value F(P) is:';
Mx4 = 'Error bound for F(P) is:';
clc,echo off,diary output,...
disp(''),disp(Mx1),disp(V0),disp(Mx2),...
disp(['± ',num2str(dV),' ± ',num2str(dV)]),disp(''),...
disp(Mx3),disp(y0),disp(Mx4),disp(['± ',num2str(dy)]),...
diary off,echo on