第6-4回 MUSCL法と流束制限関数[python]
今回は,第6-2回と第6-3回でそれぞれ説明したMUSCL法と流速制限関数を組み合わせます。
まず,なぜ流速制限関数を用いなければならないのか、考えます。
第6-2回で算出したMUSCL法は次の形でした。
$$u^L_{j+1/2}=u_j+\frac{1}{4}(1-k)(u_j-u_{j-1})+\frac{1}{4}(1+k)(u_{j+1}-u_j)$$
上記の式に対して、次のように変形を行います。
$$u^L_{j+1/2}-u_{j+1}=u_j-u_{j+1}+\frac{1}{4}(1-k)(u_j-u_{j-1})+\frac{1}{4}(1+k)(u_{j+1}-u_j)$$
$$⇔u^L_{j+1/2}-u_{j+1}=\frac{1}{4}(1-k)(u_j-u_{j-1})+\frac{1}{4}(-3+k)(u_{j+1}-u_j)$$
$$⇔u_{j+1}-u^L_{j+1/2} =\frac{1}{4}(3-k)(u_{j+1}-u_j) -\frac{1}{4}(1-k)(u_j-u_{j-1})$$
ここで,次のような条件を考えます。
$$\Delta ^+_{j}=u_{j+1}-u_j > 0 \ and \ k < 1 \ \ \ \ \dots ①$$
さらにこの条件下で,\(\Delta ^-_{j} \)が大きい,すなわち\(0 < \Delta ^+_{j} \ll \Delta ^-_{j}\)のとき
$$u_{j+1}-u^L_{j+1/2} = -\frac{1}{4}(1-k)(u_j-u_{j-1}) <0$$
$$⇔u_{j+1}<u^L_{j+1/2} \ \ \ \ \dots ②$$
①,②より
$$ u_{j}<u_{j+1}<u^L_{j+1/2}$$
このようにセル\(j\)と\(j+1\)の間の値を求めているにも関わらず,その値より大きい境界値が得られてしまうことがあります。
そこで,流速制限関数を用いてこのような余計な極大が生じないようにしなければなりません。
MUSCL法に対しては流束制限関数を用いて,次のように修正を行います.
$$u^L_{j+1/2}=u_j+\frac{1}{4}(1-k)\phi(1/r_j)(u_j-u_{j-1})+\frac{1}{4}(1+k)\phi(r_j)(u_{j+1}-u_j)$$
ここで,\(\phi\)は流束制限関数であり,\(r_j\)は次式で定義されます.
$$r_j=\frac{u_j-u_{j-1}}{u_{j+1}-u_j}$$
さらに,minmod関数は次のように一本の式で書きます.
$$minmod(r_j,1)=sgn(r_j) max[0,(min[| r_j |,sgn(r_j)×1])]$$
ここでsgn(x)はxの符号を表します。
minmod関数は次の特徴があり,これを用いて式変形を行います。
$$minmod(\frac{p}{q},1) q= minmod(p,q) $$
$$u^L_{j+1/2}=u_j+\frac{1}{4}(1-k)\phi(1/r_j)(u_j-u_{j-1})+\frac{1}{4}(1+k)\phi(r_j)(u_{j+1}-u_j)$$
$$⇔u^L_{j+1/2}=u_j+\frac{1}{4}(1-k)\phi(1/r_j)(u_j-u_{j-1})+\frac{1}{4}(1+k)\phi(r_j)(u_{j+1}-u_j)$$
$$⇔u^L_{j+1/2}=u_j+\frac{1}{4}(1-k)minmod(1/r_j)\dot (u_j-u_{j-1}) \\ + \frac{1}{4}(1+k)minmod(r_j)\dot (u_{j+1}-u_j)$$
$$⇔u^L_{j+1/2}=u_j+\frac{1}{4}(1-k)minmod(u_{j+1}-u_j,u_j-u_{j-1})$$ $$+\frac{1}{4}(1+k)minmod(u_j-u_{j-1},u_{j+1}-u_j)$$
さらに衝撃波や接触不連続面を鋭く捕らえられる圧縮パラメータb(>1)を導入します。 流束制限関数にminmodを用いる場合は(1\le b \le (3-k)/(1-k))として使用します。
まとめると次のような式になります。
$$u^L_{j+1/2}=u_j+\frac{1}{4}(1-k)minmod( \Delta^{-}_j ,b \Delta^{+}_j )+\frac{1}{4}(1+k)minmod( \Delta^{+}_j , b\Delta^{-}_j )$$
$$u^R_{j+1/2}=u_{j+1}-\frac{1}{4}(1-k)minmod( \Delta^{+}_{j+1} ,b\Delta^{-}_ {j+1} )-\frac{1}{4}(1+k)minmod( \Delta^{-}_ {j+1} ,b\Delta^{+}_ {j+1} )$$
$$\Delta^{+}_j= u_{j+1}-u_j , \Delta^{-}_j= u_j-u_{j-1} $$
次回は第6回のまとめとして,コードに落とし込みたいと思います。
*第6-5回にコードは書いてあります。
by hide
Tweet
コメントを残す