どうも、木村（@kimu3_slime）です。
今回は、Julia（SymPy）でラプラス変換、逆変換を求める方法を紹介します。
準備
SymPyを使うので、持っていなければインストールしておきましょう。
1 2 | using Pkg Pkg.add("SymPy") |
準備として、以下のコードを実行しておきます。
1 | using SymPy |
ラプラス変換、逆変換の求め方
まず、今回記号として使う変数を用意しましょう。
1 2 3 | @vars t positive=true @vars s @vars ω real=true |
\[ \begin{aligned}L(f)(s) := \int_0 ^\infty e^{-st}f(t)dt\end{aligned} \]
です。変数\(t\)から\(s\)へのラプラス変換は、「sympy.laplace_transform(関数 ,t,s)」によって求められます。試しに定数関数を変換してみましょう。
1 | sympy.laplace_transform(1 ,t,s) |
1 | (1/s, -oo, re(s) > 0) |
第1成分が変換後の関数\(1/s\)です。第2成分を\(a\)とすると、それは\(\mathrm{Re }(s)>a\)において収束することを意味しています。第3成分は追加で補助的な収束条件が表示されるようです。
変換結果だけ見るならば、
1 2 | sympy.laplace_transform(t ,t,s)[1] sympy.laplace_transform(t^2 ,t,s)[1] |
\[ \begin{aligned}\frac{1}{s^{2}}\end{aligned} \]
\[ \begin{aligned}\frac{2}{s^{3}}\end{aligned} \]
といったようにすれば良いです。
三角関数、指数関数のラプラス変換を求めてみましょう。
1 2 3 | sympy.laplace_transform(cos(ω *t) ,t,s)[1] sympy.laplace_transform(sin(ω *t) /ω ,t,s)[1] sympy.laplace_transform(exp(ω *t) ,t,s)[1] |
\[ \begin{aligned}\frac{s}{s^{2} + ω^{2}}\end{aligned} \]
\[ \begin{aligned}\frac{1}{s^{2} + ω^{2}}\end{aligned} \]
\[ \begin{aligned}\frac{1}{s – ω}\end{aligned} \]
変数\(s\)から\(t\)へのラプラス逆変換は、「sympy.inverse_laplace_transform(関数 ,s,t)」で求められます。
1 2 3 4 | sympy.inverse_laplace_transform(1/s ,s,t) sympy.inverse_laplace_transform(s/((s^2+1)*(s^2+4)) ,s,t) sympy.inverse_laplace_transform(s/(s(s+2)) ,s,t) sympy.inverse_laplace_transform(1/(s+1)^2 ,s,t) |
\[ \begin{aligned}1\end{aligned} \]
\[ \begin{aligned}\frac{\cos{\left(t \right)}}{3} – \frac{\cos{\left(2 t \right)}}{3}\end{aligned} \]
\[ \begin{aligned}- 2 e^{- 2 t}\end{aligned} \]
\[ \begin{aligned}t e^{- t}\end{aligned} \]
部分分数分解をしなくても、簡単に有理関数のラプラス逆変換が求められます。
\(\frac{e^{-s}}{s}\)を逆ラプラス変換してみましょう。
1 | sympy.inverse_laplace_transform(exp(-s)/s ,s,t) |
\[ \begin{aligned}\theta\left(t – 1\right)\end{aligned} \]
この\(\theta\)は、SymPyにおいてはヘビサイド関数（単位ステップ関数）を意味しています。
\[ \begin{aligned}\theta(x) = \begin{cases} 0 & \text{for}\; x < 0 \\ \frac{1}{2} &\text{for}\; x = 0 \\1 & \text{for}\; x > 0 \end{cases}\end{aligned} \]
ヘビサイド関数は「Heaviside(t)」として記述できます。これをラプラス変換すれば、さきほどの関数が得られます。
1 | sympy.laplace_transform(Heaviside(t-1) ,t,s)[1] |
\[ \begin{aligned}\frac{e^{- s}}{s}\end{aligned} \]
ディラックのデルタ関数も、「DiracDelta(t)」で扱えます。
1 2 | @vars t nonnegative=true sympy.laplace_transform(DiracDelta(t) ,t,s)[1] |
\[ \begin{aligned}1\end{aligned} \]
変数\(t\)の仮定を非負と指定しました。\(t\)が正であると仮定すると、デルタ関数の値が常に0になり、ラプラス変換が\(0\)になって意味がなくなってしまいます。
変数を少しずらして変換したり、逆変換してもきちんと結果が得られます。
1 2 | sympy.laplace_transform(DiracDelta(t-1) ,t,s)[1] sympy.inverse_laplace_transform(exp(-s) ,s,t) |
\[ \begin{aligned}e^{- s}\end{aligned} \]
\[ \begin{aligned}\delta\left(t – 1\right)\end{aligned} \]
以上、Julia（SymPy）でラプラス変換、逆変換を求める方法を紹介してきました。
ラプラス変換、逆変換の表や結果を見なくても、数式を打ち込めば的確な結果が表示されるので便利ですね。
木村すらいむ（@kimu3_slime）でした。ではでは。
コロナ社 (2020-03-26T00:00:01Z)
￥7,353 (コレクター商品)
