はじめに
目標
今回は以下のラプラス変換の数式を眺めながら実装してみた.
${\displaystyle F(s)=\int _{0}^{\infty }f(t)\mathrm {e} ^{-st}\mathrm {d} t}$
また, タイトルの通りpythonコードと数式を対応させながら実装していく.
※注意: 比較できる実装を見つけられず, 正しく実装出来たかの検証はしていない
準備
まず最初のおまじない.これで, πやeなどを使ったり, 求めた値をグラフ化できる.
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
関数fの定義
今回は以下の部分はsin波を用いる事にします.
${f(t)}$
ここは以下のように書くとfがsin関数になります.
f = np.sin
tを与えるとちゃんとsin(t)の値が出力される.
>>> t = np.pi / 2
>>> f(t)
1.0
上限とdtの設定
本来はそのまま以下の数式のように0から∞まで扱いたいのですが・・・
${\displaystyle \int _{0}^{\infty }f(t)\mathrm {e} ^{-st}\mathrm {d} t}$
さすがに無限は扱えないので, 上限とdtの値を設定します.（妥協・・・）
今回は上限を100, dt=0.001とします. そして以下で0から100まで0.001間隔で変化する数列Tが作れます.
dt = 0.001
max_t = 100
T = np.arange(int(max_t / dt)) * dt
この数列を関数fに与えてsin波を作ってみると以下のような感じ.
plt.plot(T, f(T))
plt.show()
結果 ちゃんとできてそうですね.
関数Fの定義
ここまで出来たらF(s)はlambda式を使ってさくっと書けます.
（数式まんまなのでわかりやすい・・・かな？）
F = lambda s: np.sum(f(T) * np.e**(-s * T) * dt)
そして以下の式が完成！ρ(・ω・、)
${\displaystyle F(s)=\int _{0}^{\infty }f(t)\mathrm {e} ^{-st}\mathrm {d} t}$
いぁ・・・∞じゃないからこうかな(´・ω・｀)
${\displaystyle F(s)=\int _{0}^{100}f(t)\mathrm {e} ^{-st}\mathrm {d} t}$
最後にsが1~100まで変化した時, F(s)がどのように変化するか見ておきます.
map使ってサクッと作った方が良いのですが・・・F(s)を順番に求めてるのが分かりやすいのであえてfor文で.
ついでにmap使った書き方はコメントで載せておく.
S = range(1, 100)
N = []
for s in S:
N.append(F(s))
plt.plot(S, N)
# plt.plot(S, list(map(F, S)))
plt.show()
まとめ
今回はラプラス変換の数式を眺めながら実装してみた.
未検証で正しく実装出来たかは分からないため, 今後検証は行いたい.
またwikiをみると
ラプラス変換は, 関数 f (t) にいったん e−σtθ(t) を乗じてからフーリエ変換する操作であると考えることができる
${\displaystyle F(s):=F(\sigma ,\omega )=\int _{-\infty }^{\infty }\theta (t)f(t)\mathrm {e} ^{-\sigma t}\mathrm {e} ^{-i\omega t}\,\mathrm {d} t{\overset {s=\sigma +i\omega }{=!=}}\int _{0}^{\infty }f(t)\mathrm {e} ^{-st}\,\mathrm {d} t}$
っとあるので, 今後は以前実装したフーリエ変換にかけてみたいρ(・ω・、)
（同時に検証もできるかな？）