対流項(移流項)を含む偏微分方程式の数値解法 [ネコ騙し数学]
対流項(移流項)を含む偏微分方程式の数値解法
問題 u(x,t)=f(x−ct)であるとき、
であることを示せ。ただし、cは定数とする。
【解】
ξ=x−ctとおき、合成関数の微分法を用いると、
したがって、
(解答終)
したがって、fを微分可能な任意の関数とすると、u(x,t)=f(x−ct)は
の(一般)解になる。
f(x)=sinxとおくと、
となり、これが(1)の(特殊)解になることから分かるように、(1)は波の方程式の一つである。
(1)式の左辺第1項を
と、時間微分に関しては前進差分を用いて近似し、(2)式の左辺の第2項を
と、空間微分に関しては中心差分を用いて近似すると、(1)式は
と近似することができ、これから
という漸化式を得ることができる。
とおくと、
となる。
(6)式の無次元数(物理の単位をもたない数のこと)をクーラン数という。
したがって、計算の起点となる、t=t₀におけるuの初期条件とx=x₀における初期条件が与えられていれば、(7)式を用い逐次的に(1)式の近似解を求めることができる。
ということで、初期条件を
さらに、c=0.2、Δx=1、Δt=1として、(1)を解いた結果は次のようになる。
赤い曲線はx=1、黄色の曲線はx=2、緑はx=3、茶色はx=4におけるuの数値解の時間tにおける変化を表している。
c=0.2だからx=1に高さ1の波が到達するのはt=5なのにも関わらず、t=Δt=1に既に波の一部(?)が到達しているだけでなく、t=5になっても波の高さは1にならない。しかも波の高さは、約t=20以降、正弦曲線のような曲線を描き、振動する。
x=2、x=3、x=4と下流にいけばいくほど、波の最大の高さは高くなり、より激しく振動するという、物理的にありえない、数値的に振動する解を(7)式は求めてしまう。
つまり、(1)式の時間微分に前進差分、空間部分に中心差分を用いるFTCS(Forward in Time Center in Space)は、数値的に不安定な解法で、この解法を(1)式に用いると物理的にありえない振動解を求めてしまうということ。
そこで、なぜ、FTCS法は不安定なのか、フォン・ノイマンの安定判定法を用いて調べることにする。
とし、(7)式に代入すると、
となり、
したがって、フォン・ノイマンの安定判定法によると、クーラン数に関わらず、FTCSは無条件不安定!!
ということで、FTCSは微分方程式(1)の数値解を求めることに使えないことが明らかになった。
そこで、(1)式の第2項の空間微分を
と風上化(上流化)すると、c>0のとき、(1)は
したがって、
となる。
(8)式を(13)式に代入すると、
よって、
したがって、
の時に安定になる。
風上(上流)差分を用いて、c=0.2、Δx=1、Δt=1、すなわち、クーラン数C=0.2として計算した結果は、次の通り。
高さ1の波の到達時刻は実際とは異なるけれど、今度は振動することなく計算できていることが分かる。
また、下図を見ると、この波の解は
でなければならないのに、風上差分を用いた(13)の数値解は、波の形が崩れていることが分かるだろう。
参考までに、同一の条件でFTCSを用いて計算した結果を下図に示す。