使用scipy.optimize的fsolve,root函数求解非线性方程问题

from scipy.optimize import fsolve, root
import numpy as np

# 定义方程内容
def f(x, *arg):
    return arg[0] * 2 ** (1 - x) / (1 - x) + (1 - arg[0]) * 1.6 ** (1 - x) / (1 - x) - arg[0] * 3.85 ** (1 - x) / (
            1 - x) - (1 - arg[0]) * 0.1 ** (1 - x) / (1 - x)

# 参数p为一个超参数

results = [[p * 0.1, fsolve(f, x0=0, args=(p * 0.1))[0]] for p in range(2, 10)]

print(np.array(results), '\n')

results = [[p * 0.1, root(f, x0=0, args=(p * 0.1))['x'][0]] for p in range(2, 10)]

print(np.array(results))


# output results
[[ 0.2        -0.94683705]
 [ 0.3        -0.48657472]
 [ 0.4        -0.14263228]
 [ 0.5         0.14636333]
 [ 0.6         0.4114561 ]
 [ 0.7         0.67618002]
 [ 0.8         0.9705809 ]
 [ 0.9         1.3683912 ]] 

[[ 0.2        -0.94683705]
 [ 0.3        -0.48657472]
 [ 0.4        -0.14263228]
 [ 0.5         0.14636333]
 [ 0.6         0.4114561 ]
 [ 0.7         0.67618002]
 [ 0.8         0.9705809 ]
 [ 0.9         1.3683912 ]]

"""利用gekko求解非线性方程组"""
from gekko import GEKKO

m = GEKKO()
x = m.Var(value=0)  # 给定初值为0
y = m.Var(value=0)  # 给定初值为0

m.Equations([x ** 2 / 4 + y ** 2 == 1,
             (x - 0.2) ** 2 - y == 3])
m.solve(disp=False)
x, y = x.value, y.value
print(x, y)

from scipy.optimize import fsolve

def f(X):
    x = X[0]
    y = X[1]
    return [x ** 2 / 4 + y ** 2 - 1,
            (x - 0.2) ** 2 - y - 3]

X0 = [0, 0]
result = fsolve(f, X0)
print(result)

from scipy.optimize import fsolve, root

def f(X):
    x = X[0]
    y = X[1]
    return [x ** 2 / 4 + y ** 2 - 1,
            (x - 0.2) ** 2 - y - 3]

X0 = [10, 10]
result1 = fsolve(f, X0)
result2 = root(f, X0)

print(result2)

from scipy.optimize import leastsq

def f(X):
    x = X[0]
    y = X[1]
    return [x ** 2 / 4 + y ** 2 - 1,
            (x - 0.2) ** 2 - y - 3]

X0 = [10, 10]
h = leastsq(f, X0)
print(h)

from sympy import symbols, Eq, solve, nsolve

x, y = symbols('x y')
eqs = [Eq(x ** 2 / 4 + y ** 2, 1),
       Eq((x - 0.2) ** 2 - y, 3)]

print(solve(eqs, [x, y]))

from sympy import symbols, Eq, nsolve
x, y = symbols('x y')
eqs = [Eq(x ** 2 / 4 + y ** 2, 1),
       Eq((x - 0.2) ** 2 - y, 3)]

X0 = [3, 4]
print(nsolve(eqs, [x, y], X0))