pythonscipyscipy-optimizefsolve

scipy fsolve() method throws different first value when the second value changes


I am trying to get the roots of f(x) = x**3-100*x**2-x+100

I know two of the roots are 1, 100. So set guesses close to this should work.

This is the code (with two different guesses):

from scipy.optimize import fsolve

def f(x): return x**3-100*x**2-x+100
s0 = fsolve(f, [2,80])
s1 = fsolve(f, [2,75])  # <-- changed from 80 to 75

print(s0)
print(s1)

result:

[  1. 100.]
[ -1. 100.]   # <-- first argument is -1

I was expecting to get the same result twice.

So my question is why when I changed the second input parameter in the list, did the first value in the result change from +1 to -1 ?

additional notes (based on the comments: I run this extra piece of code to see what is returned.

from scipy.optimize import fsolve

def f(x): return x**3-100*x**2-x+100

for i in range(65, 85):
    x = fsolve(f, [2,i])
    print(i, x)

and get this back:

65 [-1.  1.]
66 [ 1.94287787 65.87673804]
67 [ 2. 67.]
68 [  1. 100.]
69 [  1. 100.]
70 [  1. 100.]
71 [  1. 100.]
72 [ -1.00000028 100.        ]
73 [  1. 100.]
74 [  1.00000006 100.        ]
75 [ -1. 100.]
76 [  1. 100.]
77 [  1. 100.]
78 [  1. 100.]
79 [  1. 100.]
80 [  1. 100.]
81 [  1. 100.]
82 [  1. 100.]
83 [  1. 100.]
84 [  1. 100.]

I appreciate that -1 is the third root, but it was never the closest root given the first element of the guess is always 2.

And also, we see that fsolve(f, 2) in isolation returns array([1.])...


Solution

  • fsolve solves does not accept guesses. It takes only one guess for a multidimensional system

    f1(x1, x2, x3...) = 0
    f2(x1, x2, x3...) = 0
    f3(x1, x2, x3...) = 0
    …
    

    When you run fsolve(f, [2,80]), you attempt to solve a 2D system. When you run fsolve(f, [2]), you attempt to solve a 1D system. You may try to solve a 3D system with fsolve(f, [2,80, 12]) :)

    fsolve passes numpy.ndarray to your callback f and it’s because NumPy magic f(x) returns multidimensional output in your case. Most people have to write different codes for each function in the system to get multidimensional output from a callback they pass to fsolve. And in general, such functions are supposed to be different. fsolve does not assume that they can be the same, but it seems fsolve can handle the case with the same functions quite efficiently.

    The Powell method (or its MINPACK modifications used in scipy.fsolve) changes elements of multidimensional vector x = (x1, x2, x3, ...) simultaneously. I should admit that I cannot understand MINPACK’s Fotran code. But it should be some kind of linesearch. And when you change the 2D initial point you change the starting point for the first line in 2D space that is constructed inside the root-searching procedure. That's why changing "the second guess" influences on "the first result".

    If you do not want such dependence run fsolve independently for each guess.