How to add the count of iterations and relative error to the code?

6 views (last 30 days)
I have this Newton-Rhapson and bisection method functions, but they cannot track how mant literations operated and what the relative error is. I am wondering where I can add lines of code to let them able to achieve those?
function root = newton(fun,x0,fprime,tol,maxiter)
% newton function computes the root of the function
% fun using Newton-Raphson method
% Check whether the initial guess is a root
% Error occurs if root is not determined within the
% desired tolerance in the max number iterations
root = x0;
iter = 0;
if ~isa(fprime, 'function_handle')
fprime = @(x) central_difference(fun,x,1e-4);
end
while iter < maxiter
if fprime(root) == x0
disp('The derivative is equal to 0, the current root is %.2f\n',root)
break
end
dx = -fun(root)/fprime(root);
root = root + dx;
if abs(dx) <= tol
break
end
iter = iter + 1;
end
if iter == maxiter && abs(dx) > tol
error('root not found in given iterations')
end
end
function [c] = bisection(fun,a,b,tol,maxlits,plot_output)
% bisection function computes the root, c, of a function, fun
% initial intercal
c_old = 0;
iters = 0;
e_a_vec = [];
% While loop
while iters < maxlits
iters = iters + 1;
c = (a + b)/2;
if plot_output == 1
x = linspace(a,b);
y = [];
for i = 1:length(x)
y(i) = fun(x(i));
end
nexttile
hold on
plot(a,fun(a),'-o')
plot(b,fun(b),'-o')
plot(c,fun(c),'-o')
plot(x,y,'k-')
yline(0,'--')
end
if iters >= 1
e_a = (c - c_old)/c;
e_a_vec(iters) = e_a;
if abs(b-a)/2 < tol | abs(fun(c)) <= tol | abs(e_a) <= tol
break
end
else abs(fun(c)) > tol
error('The root is not within the desired tolerance')
end
if fun(a) * fun(b) < 0
if fun(a) * fun(c) > 0
a = c;
else
b = c;
end
c_old = c_old + c;
else
error('The limites do not bracket a root')
end
end

Answers (1)

Sulaymon Eshkabilov
Sulaymon Eshkabilov on 21 Oct 2021
Here are how you can obtain the final iteration and rel./abs error values when the roots are found.
function root = newton(fun,x0,fprime,tol,maxiter)
% newton function computes the root of the function
% fun using Newton-Raphson method
% Check whether the initial guess is a root
% Error occurs if root is not determined within the
% desired tolerance in the max number iterations
root = x0;
iter = 0;
if ~isa(fprime, 'function_handle')
fprime = @(x) central_difference(fun,x,1e-4);
end
while iter < maxiter
if fprime(root) == x0
disp('The derivative is equal to 0, the current root is %.2f\n',root)
break
end
dx = -fun(root)/fprime(root);
root = root + dx;
if abs(dx) <= tol
fprintf('Simulation is halted after %d iterations \n', iter) % Number of iterations
fprintf('Simulation is halted after %f of Rel. Error (Derivative Value) \n', dx) % Rel. Error
break
end
iter = iter + 1;
end
if iter == maxiter && abs(dx) > tol
error('root not found in given iterations')
end
end
function [c] = bisection(fun,a,b,tol,maxlits,plot_output)
% bisection function computes the root, c, of a function, fun
% initial intercal
c_old = 0;
iters = 0;
e_a_vec = [];
% While loop
while iters < maxlits
iters = iters + 1;
c = (a + b)/2;
if plot_output == 1
x = linspace(a,b);
y = [];
for i = 1:length(x)
y(i) = fun(x(i));
end
nexttile
hold on
plot(a,fun(a),'-o')
plot(b,fun(b),'-o')
plot(c,fun(c),'-o')
plot(x,y,'k-')
yline(0,'--')
end
if iters >= 1
e_a = (c - c_old)/c;
e_a_vec(iters) = e_a;
if abs(b-a)/2 < tol | abs(fun(c)) <= tol | abs(e_a) <= tol
fprintf('Simulation is halted after %d iterations \n', iters) % Number of iterations
fprintf('Simulation is halted after %f of Rel Diff \n', abs(b-a)/2) % Rel. Difference
fprintf('Simulation is halted after %f of Abs. Fcn. Value \n', abs(fun(c))) % Abs. Value
fprintf('Simulation is halted after %f of Abs. Error \n', abs(e_a)) % Abs. Error
break
end
else abs(fun(c)) > tol
error('The root is not within the desired tolerance')
end
if fun(a) * fun(b) < 0
if fun(a) * fun(c) > 0
a = c;
else
b = c;
end
c_old = c_old + c;
else
error('The limites do not bracket a root')
end
end

Categories

Find more on Symbolic Math Toolbox in Help Center and File Exchange

Products


Release

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!