How to use Newton-Raphson for numerical solution of two variables?

Hello community,

I have a non-linear function f(x,y), which I would like to find the roots of with the Newton-Raphson method. However, I haven't yet found a simple code on the internet in case of two variables, which just lets me enter my function and get the result. If I find one, the solutions on the internet always need two equations, but I only have one?

[Matt J]

Newton-Raphson only applies to N equations in N unknowns. If you have fewer equations than unknowns, there will normally be an infinite continuum of roots. For example,

 f(x,y)=x+y

has roots at all x=-y.

[Jonas]

What if the function is complex? Can i divide it into a real and an imaginary equation - in that case I would have two equations? If, how?

Answer 1

[Matt J]

Can i divide it into a real and an imaginary equation

You could just rewrite your function to return a 2x1 vector containing the real and imaginary part respectively.

If f() is your existing function, you could also do

 f_new=[real(f(x,y)); imag(f(x,y))]

[Jonas]

... unfortunately it seems like that the equation is too huge for matlab beeing able to divide the equation into a real and imaginary part :(

[Jonas]

fname=@(x,y)x*(- 4157552041235969998848 + 14994143880568163532800*i) + ...
    y*(- 14994143880568163532800 - 4157552041235969998848*i) + ...
    (x + y*i)^2*(29558214356309966848 - 131338673823403081728*i) + ...
    (x + y*i)^3*(100369407825487824 + 178402860122995232*i) + ...
    (x + y*i)^4*(- 3688268389229367/16 + (578208319202799*i)/16) + ...
    (x + y*i)^5*(3404177150896959/32768 - (8024424565434255*i)/65536) + ...
    (x + y*i)^6*(- 3003839359850141/268435456 + (2314057354380575*i)/33554432) + ...
    (x + y*i)^7*(- 4711673569865107/549755813888 - ...
    (7173966121489549*i)/549755813888) + (x + y*i)^8 - ...
    200770287297318937952256 + 213383934120794144112640*i;

[Matt J]

I'm having no difficulty with it.

    >> f=@(x,y) [real(fname(x,y)) ; imag(fname(x,y)) ]

    >> fname(.1,.2)

    ans =

      -2.0418e+23 + 2.1406e+23i

    >> f(.1,.2)

    ans =

       1.0e+23 *

       -2.0418
        2.1406

[Jonas]

Ah sorry, my mistake :). But I'm still stuck here - to get the final values. Do you have an idea?

[Matt J]

Now that you have a system of 2 equations in 2 unknowns, the Newton-Raphson code that you've been using should work fine. There is also FZERO, if you have it.

[Jonas]

I have FZERO and found this code <http://www.mathworks.com/matlabcentral/fileexchange/27120-newton-raphson-method-for-2-variables/content/newton2v2.m> but I'm stuck with both of it. Though I would prefer FZERO :)

[Matt J]

Sorry, I meant to say FSOLVE, not FZERO.

Answer 2

[Jonas]

I went one step back and did following: I got following equation with pt as a variable instead of pt = x + i*y.

eq = pt*(- 4157552041235969998848 + 14994143880568163532800*i) + pt^2*(29558214356309966848 - 131338673823403081728*i) + pt^3*(100369407825487824 + 178402860122995232*i) + pt^4*(- 3688268389229367/16 + (578208319202799*i)/16) + pt^5*(3404177150896959/32768 - (8024424565434255*i)/65536) + pt^6*(- 3003839359850141/268435456 + (2314057354380575*i)/33554432) + pt^7*(- 4711673569865107/549755813888 - (7173966121489549*i)/549755813888) + pt^8 - 200770287297318937952256 + 213383934120794144112640*i

Here it is possible to say e.g. solve(eq == 0) and I get some nice results. However, if I check one of the results and insert one of the results (e.g. pt = -15-7*i) the final result is not equal 0. How does that make sense?

[Matt J]

Floating point errors probably. You have very large coefficients and also a rather large polynomial order. ROOTS might give you a more exact result (and is more appropriate anyway seeing as your function is a polynomial), but it's hard to say.

[Jonas]

Alright - and thank you very much for all your help!