Sympy是python中非常強大的符號運算庫,可以以書寫習慣表示數學表達式。下面介紹用Sympy求方程數值解的方法。
下面代碼全部在
from sympy import *
init_printing(use_unicode=True) # 按書寫習慣輸出
下運行。
數學表達式的輸入
首先聲明符號:
x = symbols('x')
即計算機中的變量x代表數學表達式中的x。在后文輸出中所有的x會顯示為x。如果x=symbols('x0'),則輸入的方程中所有x將在輸出中以x0表示。
如果需要希臘字母
l, r = symbol('lambda rho')
l, r將分別以\(\lambda, \rho\)表示。可以在一個表達式中同時聲明多個符號。
或者使用var()聲明:
var('x')
與上面等效。
聲明表達式:
f = (5/x)*(exp(x)-1)-exp(x)
此時若輸出f可以看到書寫習慣的表達式。由於表達式在markdown下顯示不正常,在此不放置示例。注意f的類型是class 'sympy.core.add.Add'
求\(f(x)=0\)數值解
因為有的函數零點不止一個,因此在Sympy中解的輸出為一個list。使用solve(表達式,自變量符號)可以解析地解方程:
s, = solve(f, x)
這里根據上面f的賦值,得到s為
LambertW(-5e**-5)+5
其中用了特殊函數表達。
我們需要求這個結果的數值近似,則輸出
s.evalf()
得到輸出
4.96511423174428
就是方程\(f(x)=0\)的數值解。
更多solve()的功能見Sympy documentation: http://docs.sympy.org/latest/tutorial/solvers.html
求給定自變量\(x\)值時函數\(f(x)\)的值 | 將表達式轉化為函數
f.evalf(subs = {x:4.96})
得到\(f(4.96)\)的數值
0.141885450782171
如果需要以計算機函數的形式定義函數f(x),則可以使用lambdify()進行轉化:
f_func = lambdify(x, f)
之后可以調用
f_func(4.96)
輸出
0.141885450782
利用這個方法可以測試方程的數值算法,如使用sympy接口寫牛頓法等。
