■ 問題提出

在博文 逆Laplace數(shù)值逆變換^[1] 給出了數(shù)值計(jì)算Laplace逆變換的簡易程序。其中存在以下幾個(gè)問題需要討論：

• 問題1： 程序?qū)崿F(xiàn)過程原理以及優(yōu)化；
• 問題2： 運(yùn)算參數(shù)： ， ， 對于積分?jǐn)?shù)值的影響。從上篇博文中明顯看到有些計(jì)算出的結(jié)果大大偏離的實(shí)際值。比如其中的函數(shù) 所計(jì)算出來的幅值超過了1.

01程序?qū)崿F(xiàn)原理以及優(yōu)化

1.Laplace逆運(yùn)算

所以，算法的核心是對 進(jìn)行傅里葉反變換，然后在乘以 。

由于確認(rèn)變換后的函數(shù) 是實(shí)函數(shù)，因此，為了節(jié)省計(jì)算時(shí)間，只對傅里葉反變換的積分，進(jìn)行正半軸的積分，同時(shí)積分的上限由參數(shù) 決定。對積分的數(shù)值取齊實(shí)部，再乘以2便可以得到 。

2.Python積分程序?qū)崿F(xiàn)優(yōu)化

使用梯形積分來實(shí)現(xiàn)函數(shù)的積分，可以獲得更精確的積分值。理論分析可知：

其中 以及 。那么積分誤差上限為：

其中， ， 。

計(jì)算 可以有兩種方式：

方式1： 對左、右黎曼積分加權(quán)平均：

方式2： 利用如下的公式計(jì)算：

最后實(shí)現(xiàn)的代碼為：

def trapz(f, a, b, N=50): x = linspace(a, b, N+1) y = f(x) y_right = y[1:] y_left = y[:-1] dx = (b-a) / N T = dx/2 * sum(y_right + y_left) return T

利用上述公式，對 進(jìn)行積分測試：

printf(trapz(sin, 0, pi/2, 1000))

所得到結(jié)果為：0.9999997943832332?？梢钥吹绞褂胣=1000對應(yīng)的結(jié)果后面的積分精度達(dá)到了小數(shù)點(diǎn)后面6位9的小數(shù)點(diǎn)的位數(shù)。

3.實(shí)現(xiàn)Laplace逆運(yùn)算

通過上面梯形積分方法，實(shí)現(xiàn)Laplace數(shù)值逆變換，具體的子程序如下面所示。

#------------------------------------------------------------def invlt(t, fs, sigma, omiga, nint): omigadim = linspace(0, omiga, nint+1, endpoint=True) y = [(exp(1j*o*t) * fs(sigma+1j*o)).real for o in omigadim] y_left = y[:-1] y_right = y[0:] T = sum(y_right + y_left) * omiga/nint return exp(sigma*t) * T/ pi / 2
#------------------------------------------------------------def fs(s): return 1/(s*s+1)
#------------------------------------------------------------sigma = 0.2omiga=200nint=omiga*50
tdim = linspace(0, 2*pi* 3, 200)ft = [invlt(t, fs, sigma, omiga, nint) for t in tdim]

02一些基本函數(shù)的實(shí)驗(yàn)

下面通過對一些基本常見函數(shù)的laplace變換，來測試一下上述程序的性能。

Ⅰ.sin(t)

sigma=0.2, omiga=200, nint=omiga*50

Ⅱ.exp(-t)

sigam=-1+0.1, omiga=200, nint=omiga*50

Ⅲ.u(t)

Ⅳ.u(t-1)

Ⅴ.周期化脈沖信號(hào)

※ 結(jié)論

通過原理分析，可以獲得建議的Laplace數(shù)值逆運(yùn)算的正確的PYTHON算法程序。

這個(gè)程序是直接對Laplace反變換公式利用梯形積分方法獲得計(jì)算結(jié)果。通過對幾種常見的信號(hào)Laplace的反變換，驗(yàn)證了這個(gè)算法的正確性。

通過在此過程，可以看到，對于參數(shù)sigma, omiga, nint對于計(jì)算結(jié)果還是有很大的影響。另外，對于時(shí)間t，只能在比較小的范圍內(nèi)有效，當(dāng)t超過一定長度，前面所計(jì)算的結(jié)果都會(huì)出現(xiàn)比較大的誤差。

參考資料