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 2.1 (Fixed Point Iteration).
% Section 2.1, Iteration for Solving x = g(x), Page 51
echo on; clc; format long; hold off; clear
% This program implements fixed point iteration.

% Define and store g(x) in the M-file g.m
% function y = g(x)
% y = 1 + x - x.^2 ./4;

delete g.m
diary g.m; disp('function y = g(x)');...
disp('y = 1 + x - x.^2 ./4;');...
diary off;

% Remark. g.m and fixpt.m are used for Algorithm 2.1
g(0); % Test for file g.m
pause % Press any key to see the graph y = g(x).

clc; clg;
a = 0;
b = 5;
c = 0;
d = 5;
h = (b-a)/150;
X = a:h:b;
Y = g(X);
X1 = [-1000 1000];
Y1 = [-1000 1000];
axis([a b c d]);...
plot(X1,Y1,'-r',X,Y,'-g');...
hold on;...
plot([a b],[0 0],'b',[0 0],[c d],'b');...
xlabel('x');...
ylabel('y');...
title('The line y = x and the curve y = g(x).');...
grid;...
axis;...
hold off;...
shg; pause % Press any key to continue.

clc;
% Place the starting value in p0

% Place the number of iterations in max1

% Place the tolerance in delta

p0 = 4.0;
max1 = 100;
delta = 1e-9;

[pc,err,P] = fixpt('g',p0,delta,max1);

pause % Press any key for the list of iterations.

clc; clg;
max1 = length(P);
for j = 1:max1-1,
k1 = 2*j-1;
k2 = 2*j;
Vx(k1) = P(j);
Vy(k1) = P(j);
Vx(k2) = P(j);
Vy(k2) = P(j+1);
end
Vy(1) = 0;

c = 0;
d = 2;
Z0 = zeros(1,length(P));
axis([a b 0 2]);...
plot(X1,Y1,'-g',X,Y,'-g',Vx,Vy,'-r',P,Z0,'or');...
hold on;...
plot([a b],[0 0],'b',[0 0],[c d],'b');...
xlabel('x');...
ylabel('y');...
title('Graphical analysis for fixed point iteration.');...
grid;...
axis;...
hold off;...
shg; pause % Press any key to continue.

max1 = length(P);
J = 1:max1;
points = [J;P];
Mx1 = 'Computations for the fixed point iteration method.';
Mx2 = ' k p(k)';
Mx3 = 'The fixed point is g(p) = p = ';
Mx4 = 'The error estimate for p is ± ';
clc,echo off,diary output,...
disp(''),disp(Mx1),disp(''),disp(Mx2),disp(points'),...
disp(''),disp(Mx3),disp(pc),...
disp([Mx4,num2str(err)]),diary off,echo on