user5513029
user5513029

Reputation: 11

Newton’s method pth root- MATLAB

Newton’s method can be used to calculate the p-th root of a positive real number. Write a Matlab function myrootp.m which calculates a^(1/p) for any a> 0 and a positive integer p. The function should stop when abs((Xk^(p)/a)-1)<10^(-15).

I used the following code but I cannot find the error:

   function root1 = root1( a,X0,p)

    Xk=X0;
    fd= p*Xk^(p-1);
    f=(Xk^p)-a;
    if a<0 
then disp('Error:a must be positive and real')
    else
for k=1:p
    if abs(Xk^p/a-1)<10^(-15)
        then return disp('Error');
    else            
    Xk1=Xk-f/fd;
    Xk=Xk1;
    fd= p*Xk^(p-1);
    f=(Xk^p)-a;
    end
end

  end

  root1 = Xk;
  return
  end

Upvotes: 0

Views: 355

Answers (1)

rayryeng
rayryeng

Reputation: 104503

I'd like to point out for those who are unfamiliar with this algorithm that this is finding the nth root of a number using Newton's method.


Back to your problem, there are two errors:

  1. There is no such keyword as then in MATLAB. Get rid of it.
  2. You'll want to iterate until the error is as you specified. You are looping for as many times as there is the power p, and it's highly possible that the algorithm can converge before then, or the worst case require more iterations than p to converge.

Making these two changes:

function root1 = root1(a,X0,p)

Xk=X0;
fd= p*Xk^(p-1);
f=(Xk^p)-a;
if a<0 
    disp('Error:a must be positive and real')
else
    while true
        if abs(Xk^p/a-1)<10^(-15)
            break;
        else            
            Xk1=Xk-f/fd;
            Xk=Xk1;
            fd= p*Xk^(p-1);
            f=(Xk^p)-a;
        end
    end
end

root1 = Xk;

To test this out, using an initial guess of 5 or X0=5, and we want to find the square root of 10 (a = 10, p = 2), we get:

>> format long g;
>> root1(10, 5, 2)

ans =

      3.16227766016838

Comparing this with MATLAB's sqrt function, we get:

>> format long g;
>> sqrt(10)

ans =

      3.16227766016838

Upvotes: 1

Related Questions