CODEWITHSUNDEEP
matlabdebugging
user2466382
user2466382

Reputation: 117

Matlab: Fsolve giving incorrect roots

I'm trying to solve this system:

x = a + e(c - e*x/((x^2+y^2)^(3/2)))
y = b + c(d - e*y/((x^2+y^2)^(3/2)))

I'm using fsolve, but not matter what I put in as the starting points for the iteration, I get the answer that the starting points are the roots of the equation. 

close all, clear all, clc

a = 1;
b = 2;
c = 3;
d = 4;
e = 5;

fsolve(@u1FsolveFUNC, [1,2])

Function:

function outvar = u1FsolveFUNC(invar)

    global a b c d e

    outvar = [...
        -invar(1) + a + e*(c - e*(invar(1) / ((invar(1)^2 + invar(2)^2)^(3/2)))) ; 
        -invar(2) + b + e*(d - e*(invar(2) / ((invar(1)^2 + invar(2)^2)^(3/2))))]

end

I could try with [1,2] as invariables, and it will say that that is a root to the equation, alltough the correct answer for [1,2] is [12.76,15.52]

Ideas?

Upvotes: 1

Views: 388

Answers (2)

Rody Oldenhuis
Rody Oldenhuis

Reputation: 38042

If you write your script like this

a = 1;
b = 2;
c = 3;
d = 4;
e = 5;

f = @(XY) [...
        -XY(1) + a + e*(c - e*(XY(1) ./ (XY(1).^2 + XY(2).^2).^(3/2))) 
        -XY(2) + b + e*(d - e*(XY(2) ./ (XY(1).^2 + XY(2).^2).^(3/2)))];

fsolve(f, [1,2])

it is a lot clearer and cleaner. Moreover, it works :)

It works because you haven't declared a,b,c,d and e to be global before you assigned values to them. You then try to import them in your function, but at that time, they are still not defined as being global, so MATLAB thinks you just initialized a bunch of globals, setting their initial values to empty ([]).

And the solution to an empty equation is the initial value (I immediately admit, this is a bit counter-intuitive).

So this equation involves some inverse-square law...Gravity? Electrodynamics?

Note that, depending on the values of a-e, there may be multiple solutions; see this figure:

enter image description here

Solutions are those points [X,Y] where Z is simultaneously zero for both equations. for the values you give, there is a point like that at 

[X,Y] = [15.958213798693690  13.978039302961506]

but also at

[X,Y] = [0.553696225634946   0.789264790080377]

(there's possibly even more...)

Upvotes: 2

p8me
p8me

Reputation: 1860

When you are using global command you have to use the command with all the variables in each function (and main workspace).

eg.

Main script

global a b c d e % Note
a = 1; b = 2; c = 3; d = 4; e = 5;   
fsolve(@u1FsolveFUNC,[1,2])

Function

function outvar = u1FsolveFUNC(invar)
global a b c d e % Note
outvar = [-invar(1) + a + e*(c - e*(invar(1) / ((invar(1)^2 + invar(2)^2)^(3/2)))) ; -invar(2) + b + e*(d - e*(invar(2) / ((invar(1)^2 + invar(2)^2)^(3/2))))]

Upvotes: 1

Related Questions