MATLAB Code MAD 4401

Fixed Point Iteration

x = pi;
for k = 1:10
x = pi + 0.5*sin(x/2);
format long
disp(x);
end

Newton’s Method

x = 1.0;
nmax = 25;
eps = 1;
xvals = x;
n = 0;

while eps>= 1e-5 & n <= nmax
y=x – ((x-2)^2-log(x))/(2*(x-2)-1/x);
xvals = [xvals;y];
eps = abs(x-y);
x = y;
n = n+1;
end
format long

disp(xvals);

Secant Method

x1 = 1.0;
x2 = 1.5;
nmax = 25;
eps = 1;
xvals = x1;
xvals = [xvals;x2];
n = 0;

f = @(x) (x-2)^2 – log(x); % Function whose root is to be found

while eps>= 1e-5 & n <= nmax
y=x2 – (f(x2)*(x2 – x1))/(f(x2)-f(x1));
xvals = [xvals;y];
eps = abs(x2-y);
x1 = x2;
x2 = y;
n = n+1;
end
format long

disp(xvals);

Method of False Position I

% Uses the method of false position to find a root p of a function f
format long;

f = @(x) (x-2)^2 – log(x); % Function whose root is to be found

% Interval [a,b] on which the root is to be found f(a)f(b)<0;
a = 1; b = 2;

step_size = b – a; % Initializing step size
eps = 1e-5; % User specified error
n = 0; % Number of iterations counter

% Main loop
while (step_size >= eps)
c = b – (f(b)*(b-a)) / (f(b)-f(a));
if f(a)*f(c) < 0
% f(a) and f(c) have opposite signs, and so p is bracketed between a
% and c. We set b=c and shorten our interval to [a,c].
step_size = b – c;
b = c;
else
% f(a) and f(c) have the same sign, so r is bracketed between c and b.
% We set a=c and shorten our interval to [c,b].
step_size = c – a;
a = c;
end
n = n+1; % Add one to the number of iterations
disp(c);
end
fprintf(‘The approximate root after %d iterations is: %12.15f’,n,c);

Method of False Position II

x = 1.5;
fp = 1.0;
nmax = 25;
eps = 1;
xvals = x;
n = 0;

while eps>= 1e-5 & n <= nmax
y=x -(((x-2)^2-log(x))*(x -fp))/((x-2)^2-log(x) – (fp-2)^2 + log(fp));
xvals = [xvals;y];
eps = abs(x-y);
x = y;
n = n+1;
end
format long

disp(xvals);

Newton’s Divided Differences

%Newton’s Interpolating FORWARD Dividing Difference

%Formula
%Pn(x)= f[x0] + f[x0,x1](x-x0) f[x0,x1,x2](x-x0)(x-x1) +….+
%f[x0,…,xn](x-x0)…(x-xn-1)

%This example is from lecture on 9/19/12
degree = 4;
n = degree + 1;

matrix = zeros(n,n+1); %allocates matrix size
%Fill matrix values here
matrix(1,1) = 0;        %x0
matrix(1,2) = -6;       %f[x0]

matrix(2,1) = .1;       %x1
matrix(2,2) = -5.89483; %f[x1]

matrix(3,1) = .3;       %x2
matrix(3,2) = -5.65014; %f[x2]

matrix(4,1) = .6;       %x3
matrix(4,2) = -5.17788; %f[x3]

matrix(5,1) = 1;        %x4
matrix(5,2) = -4.28172; %f[x4]
%   :           :         :
%could continue this for as many elements as points given

for i = 2 : n;

for j = i : n;

matrix(j,i+1) =(matrix(j,i)-matrix(j-1,i))/((matrix(j,1)-matrix(j-i+1,1)));

end

end

disp(matrix);
%The interpolating polynomial can now be constructed using the farthest
%right non-zero term of each row