日本免费高清视频-国产福利视频导航-黄色在线播放国产-天天操天天操天天操天天操|www.shdianci.com

學無先后，達者為師

網站首頁 編程語言 正文

使用scipy.optimize的fsolve,root函數求解非線性方程問題_python

作者：曲草 ? 更新時間： 2023-01-12 編程語言

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 ]]

python求解非線性方程組的幾種方式

問題

對此非線性方程組，根據方程圖像（如下圖）可知，應有4組不同的解。以下嘗試用不同的方法求出該方程組的解。


在這里插入圖片描述

1. 利用gekko的GEKKO求解

"""利用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)

輸出結果：

[-1.2961338938] [-0.7615833719]

換不同初值計算得到的結果如下：

[-1.6818042485] [0.54118722964] # ?給定初值為(-2,0)
[1.9760411678] [0.15432222765] # ?給定初值為(2,0)
[1.8018969861] [-0.43392604594] # ?給定初值為(2,-2)
[1.9760412095] [0.15432236862] # ?給定初值為(10,10)
[1.801896954] [-0.4339261545] # ?給定初值為(10,-10)

可知，用這種方法并不能得到方程組的全部解，并且最終得到的解為其解集中與給定的初值“距離”較近的一個。

2. 利用scipy.optimize的fsolve求解

optimize庫中的fsolve函數可以用來對非線性方程組進行求解。

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)

輸出結果：

[-1.29613389 -0.76158337]

換不同初值計算得到的結果如下：

[-1.68180425 ?0.54118723] # ?給定初值為(-2,0)
[1.97604116 0.15432219] # ?給定初值為(2,0)
[ 1.80189699 -0.43392605] # ?給定初值為(2,-2)
[1.97604116 0.15432219] # ?給定初值為(10,10)
[ 1.80189699 -0.43392605] # ?給定初值為(10,-10)

可知，用這種方法也不能得到方程組的全部解，并且最終得到的解與給定初值有關。

3. 利用scipy.optimize的root求解

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)

輸出結果：

? ?fjac: array([[-0.2547064 , -0.96701843],
? ? ? ?[ 0.96701843, -0.2547064 ]])
? ? ?fun: array([-3.34943184e-12, ?2.75734990e-12])
?message: 'The solution converged.'
? ? nfev: 22
? ? ?qtf: array([-1.65320424e-10, -2.73193431e-10])
? ? ? ?r: array([-3.70991104, ?0.8956477 , ?0.56891317])
? status: 1
?success: True
? ? ? ?x: array([1.97604116, 0.15432219])

結果與fsolve函數得到的結果相同。

4. 利用scipy.optimize的leastsq求解

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)

輸出結果：

(array([1.97604116, 0.15432219]), 2)

5. 利用sympy的solve和nsolve求解

5.1 利用solve求解所有精確解

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]))

輸出結果：

[[-1.68180424847377 + 1.56760579250585e-32*I
? 0.541187229573922 - 3.01196919624356e-31*I]
?[-1.29613389377477 + 1.95607066863502e-32*I
? -0.761583371898353 + 3.93313832308616e-31*I]
?[1.80189698634479 - 1.95607066863926e-32*I
? -0.433926045139482 - 8.10475677027422e-31*I]
?[1.97604115590375 - 1.56760579250161e-32*I
? 0.154322187463913 + 7.18358764343162e-31*I]]

可以看出，用這種方法能夠得到方程組的全部解，并且為精確解，缺點是求解時間較長。

5.2 利用nsolve求解數值解

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))

輸出結果：

Matrix([[1.97604115590375], [0.154322187463913]])

nsolve為數值求解，需要指定一個初始值，初始值會影響最終得到哪一個解（如果有多解的話），而且初始值設的不好，則可能找不到解。

scipy.optimize.root求解速度快，但只能得到靠近初始值的一個解。對形式簡單、有求根公式的方程，sympy.solve能夠得到所有嚴格解，但當方程組變量較多時，它求起來會很慢。而且對于不存在求根公式的復雜方程，sympy.solve無法求解。

總結

原文鏈接：https://blog.csdn.net/HelloWorld1108/article/details/104983983

欄目分類
最近更新