数値解析で線形常微分方程式を解く
プログラミングの海 初の数学系記事です。
今回は「線形常微分方程式の数値解析」です。主に数値微分を使います。
数値微分とは
数値微分とは 数値微分では、具体的な計算を用いて関数の傾きを求めます。 解析的に解けない関数の傾きを求めるのに有効な手法で、勾配法のアルゴリズムなどにも用いられています。
数式で書くと$${y'(x)=\frac{y(x+d)-y(x)}{d}}$$を使います。（前進差分）
今回はオイラー法を使って解きます。
$${x≧0}$$について考えます。
また、微小値$${d=10^{-3}}$$とします。Scratchを使って解きます。
記法について
$${y=y(x)}$$
$${y}$$は$${x}$$の関数
$${d}$$は$${\frac{dy}{dx}}$$の$${d}$$ではなく、「ただ単に小さい値」
y(x)=y'(x)
$${y(x)=y'(x)}$$を解きます。$${y(x+d)}$$が出てくるまで移項します。
$${y=ce^x}$$になれば正解です。（$${c}$$は定数）
以上の条件を使って解きます。
微分の定義より
$${y(x)=\frac{y(x+d)-y(x)}{d}}$$
両辺に$${d}$$を掛けます
$${d*y(x)=y(x+d)-y(x)}$$
見にくなってきたので逆転させます
$${y(x+d)-y(x)=d*y(x)}$$
$${y(x)}$$を移項させます
$${y(x+d)=d*y(x)+y(x)}$$
$${y(x)}$$で括ってさらに簡単にします
$${y(x+d)=y(x)*(d+1)}$$
以上から、$${y(x+d)=y(x)*(d+1)}$$という関係が分かりました。
このことを使って作っていきます。TurboWarpの方が早いのでTurboWarpを使ってます。
y(0)を動かせます。
y=y'+xy+1
もうここまでくると、少なくとも私には解析的に求めるのが大変になってきます。見た目は簡単そうなのに不思議ですね。
上のものと同じように$${y(x+d)}$$の形にしていきます。
$${y(x)=y'(x)+x*y(x)+1}$$
$${y'(x)}$$を$${y(x)}$$の式で表します
$${y(x)=\frac{y(x+d)-y(x)}{d}+x*y(x)+1}$$
右辺左辺を交換します
$${\frac{y(x+d)-y(x)}{d}+x*y(x)+1=y(x)}$$
$${\frac{y(x+d)-y(x)}{d}}$$だけを残します
$${\frac{y(x+d)-y(x)}{d}=y(x)-x*y(x)-1}$$
両辺を$${d}$$倍します
$${y(x+d)-y(x)=d(y(x)-x*y(x)-1)}$$
右辺を展開して$${y(x)}$$を移項します
$${y(x+d)=d*y(x)-d*x*y(x)-d+y(x)}$$
以上から、$${y(x+d)=d*y(x)-d*x*y(x)-d+y(x)}$$という結果が出ました。これをscratchにプログラミングします。
大体Wolfram Alphaのものと一致してるのでは。。。？
まとめ
オイラー法を使うことで解くことの難しい線形微分方程式を数値的に解くことができた。
Scratchを使ってグラフを描画することができた。
Wolfram Alphaのすごさがよくわかった。