Scipy - All the Solutions of Non-linear Equations System

6.3k views Asked by At

I have a system of non-linear equations, where can be choosed any n, so length of vector x = (x1,...,xn) can be different. For example, system can be like that:

    f1(x1,...,xn) = sum( xi + xi^2 ) = 0, i={1,n}
    f2(x1,...,xn) = sum( e^xi + xi + sin(xi*pi) ) = 0, i={1,n}

According to this example, I use fsolve() of scipy library for solving such a NLE, but it returns only one solution for every single initial approximation of *x = x0. But as n can be large (for example, n = 100), and there can be a lot of solutions, so it's not very usefull to make initial conditions x = x0 for finding every solution.

So, please, could You give me an example, how to find All the solutions by fsolve() in such a situation? Or by any other simple method?

Additional For example, I have folloving simple system:

def equations(p):
    x, y = p
    return (x**2-1, x**3-1)

With different initial condiotions I have different solutions: x, y = fsolve(equations, (0, 0)) (0.0, 0.0)

x, y =  fsolve(equations, (1, 1))
(1.0, 1.0)

x, y =  fsolve(equations, (-1, 1))
(-0.47029706057873205, 0.41417128904566508)

Is it possible to use any Scipy-function, like fsolve(), to find All the solutions (roots), for example: x, y = some_scipy_solver(equations, (x0, y0)) 1. (1.0, 1.0) 2. (0.0, 0.0) 3. (-0.47029706057873205, 0.41417128904566508) ...

where (x0, y0) = any initial approximation:(0, 0),(1, 1),(-1, 1),(0.1, 10.0), etc, and where I determine just the constraints for x0, y0, like that: -1.0 <= x0 < 1.0, 0.0 <= x0 < 11.0.

2

There are 2 answers

3
rth On

It can be difficult (or impossible ) to find numerically all the solutions even for a single non-linear equation, let along a system. For instance, consider the equation,

sin(1/x) = 0

that has an infinity of solutions in the interval [0, 1]: you can't solve this with typical root-finding algorithms.

In particular, scipy.optimize.fsolve uses local optimization approaches to find one solution of the given equation. Conceptually, you can't use it (or anything else in scipy module, really) to find all the possible solutions. Of course, if you know that the system has a given number of solutions, you can just randomly set the initial conditions to fsolve (as you have done) until you find all of them. There is no general method that would work for this problem though.

Instead, you might have better luck in solving your system analytically with sympy, if it is simple enough.

0
nightcrawler On

with ipywidgets you can do. Root will change inetractively in range -100 .. 100

def root_finder(x):
    print('Root is %.3f' %fsolve(funct,x))
    
widgets.interact(root_finder,x=(-100,100,0.5));