吾未知足 唯修身爾

一日一哩を目標としています

MENU
トップ > 物理数学 > 非粘性Burgers方程式のRによる計算

非粘性Burgers方程式のRによる計算

物理数学

先日、Burgers方程式をRで計算して解の挙動を確認したが、折角なので、Burgers方程式で動粘性係数\nu=0 とした場合の解を確認することにする。

非粘性Burgers方程式とは

非粘性Burgers方程式とは、次の形の非線形偏微分方程式である。

\dfrac{\partial u}{\partial t} + u\dfrac{\partial u}{\partial x} = 0

Burgers方程式の右辺は拡散項\nu \frac{\partial^{2} u}{\partial x^{2}} だったが、非粘性Burgers方程式は動粘性係数を\nu=0 としたものであるため、右辺は0 となる。粘性がゼロのため、波高は減衰しないことを意味するが、果たして...

非粘性Burgers方程式の差分方程式

それでは、数値計算用に非粘性Burgers方程式を差分方程式に書き換える。Burgers方程式の計算の際と同様に、時間微分に1次の前進差分、空間微分に2次の中心差分を用いることとすると、

\dfrac{u_{m}^{n+1}-u_{m}^{n}}{\Delta t}+u_{m}^{n} \dfrac{u_{m+1}^{n}-u_{m-1}^{n}}{2\Delta x} = 0

となるので、これを整理すると、

u_{m}^{n+1} = - B u_{m+1}^{n} + u_{m}^{n}+ B u_{m-1}^{n}

を得る。ここで、B=u_{m}^{n} \frac{ \Delta t}{2 \Delta x}である。（Burgers方程式の差分方程式でA=0としたものと同じ）

計算結果

これまで同様、初期条件 (1.0+\cos \pi x)および周期境界条件 u(0,t)=u(2,t) のもと上記の差分方程式に従って計算する。ここで、 dx=0.01、 dt=5.0\times10^{-4}とする。（これまでの計算と全く同じ初期条件・境界条件、メッシュサイズ）

解の挙動は以下の通り。今回も一時停止や巻戻しが出来るようyoutube経由。


www.youtube.com

計算結果の考察

コサインカーブから始まるが、波高が高いほど流速が速いため、時間の経過と共に波が突き立つ（衝撃波が生じる）。しかしながら、Burgers方程式とは異なり波高が減衰しないため、波面がほぼ垂直になった後は解が多値となって計算不能になる。

計算コード（Rコード）は省略

Burgers方程式の際のコードで\nu=0 とするだけなので省略。なお、面グラフの色は"purple"としている。

*

非粘性Burgers方程式では解がすぐに計算不能になるため、あまり面白くないかと思って当初記事にしていませんでしたが、後々振り返ってみたときに、なぜこのパターンは記事にしなかったのだろうと思いそうなため、あらためて記事にしておきました。間違い等ありましたらお手数ですがコメント等いただけますと幸甚です。