I need to change this code into one that can solve complex roots

I have this Newton-Raphson method code, but I need to change it into a code that can actually solve polynomials and give as an answer complex roots.

Code:

function [tabla, raiz]=newtonraphsonMN(f,xa,errorD,imax)
tabla=[];
error=Inf;
i=0;
xr=NaN;
df=diff(f);

while error>errorD && i<=imax
    i=i+1;
    fxa=double(subs(f,xa));
    fpxa=double(subs(df,xa));
    xr=xa-fxa/fpxa;
    error=100*abs((xr-xa)/xr);
    tabla=[tabla; [xa fxa fpxa xr error]];
    xa=xr;
end
%raiz
raiz=xr;

Answer 1

[John D'Errico]

Interesting. It does not work? :) Gosh, You could have fooled me. I'll try an example. (I've attached the code you gave, so it will be used.)

syms X
f = X^2 + X + 1;

Does f has complex roots?

solve(f)

ans = 

vpa(ans)

ans = 

Of course. I'd not have used an example that lacks complex roots, since that is your question.

[tabla, raiz]=newtonraphsonMN(f,1 + i,1e-12,100);

Strange.

tabla

tabla = 
   1.0e+02 *

   0.0100 + 0.0100i   0.0200 + 0.0300i   0.0300 + 0.0200i   0.0008 + 0.0062i   1.6125 + 0.0000i
   0.0008 + 0.0062i   0.0070 + 0.0071i   0.0115 + 0.0123i  -0.0052 + 0.0063i   0.7267 + 0.0000i
  -0.0052 + 0.0063i   0.0035 - 0.0002i  -0.0003 + 0.0126i  -0.0049 + 0.0091i   0.2687 + 0.0000i
  -0.0049 + 0.0091i  -0.0008 + 0.0001i   0.0001 + 0.0182i  -0.0050 + 0.0087i   0.0424 + 0.0000i
  -0.0050 + 0.0087i  -0.0000 + 0.0000i   0.0000 + 0.0173i  -0.0050 + 0.0087i   0.0010 + 0.0000i
  -0.0050 + 0.0087i  -0.0000 + 0.0000i   0.0000 + 0.0173i  -0.0050 + 0.0087i   0.0000 + 0.0000i
  -0.0050 + 0.0087i  -0.0000 + 0.0000i   0.0000 + 0.0173i  -0.0050 + 0.0087i   0.0000 + 0.0000i
  -0.0050 + 0.0087i   0.0000 + 0.0000i   0.0000 + 0.0173i  -0.0050 + 0.0087i   0.0000 + 0.0000i

format long
raiz

raiz = 
 -0.500000000000000 + 0.866025403784439i

To me, it seems to have worked. But then, what do I know? :)

(Hint: Do you see what I did different?)

[Elena Cecilia Valero Garza]

Hello!

Your anwser works perfectly, thank you very much.

I'm really sorry that I wasn't clear enough :(, but I am supposed to use the polyval command in the code I appended (changing the subs command for the polyval). If it's not too much of a problem, I don't know if you can help me add it to our code. Thank you very much in advance.

[Torsten]

To solve for polynomial roots, use MATLAB's "roots" command.

Newton's method will only give you one solution per call. Which one you get depends on the initial guess.

Example:

To determine all roots of p(x) = x^2+x+1:

p = [1 1 1];
sol = roots(p)

sol = 
  -0.5000 + 0.8660i
  -0.5000 - 0.8660i

polyval(p,sol)

ans = 2×1

1.0e-15 * 0.3331 0.3331