数学とか語学とか楽しいよね

数学とか語学とか楽しいよね。ドイツ語とかフランス語とか数値計算とか勉強したことをまとめます。右のカテゴリーから興味のある記事を探してください。最近はクラシックの名演も紹介しています。

【数値計算】二次精度風上差分の導出

二次精度風上差分の導出をやります。備忘録みたいなものです。点  i を中心に上流側二点の値を使うのが二次精度風上差分です。具体的には、流速が  V_i \geq 0 のときは  i i-1 i-2 の情報を使います。一方、流速が  V_i < 0 のときは  i i+1 i+2 の情報を使います。まず、 V_i \geq 0 のときから考えていきます。まず  u_{i-1} u_{i-2} の値を  u_i を中心としてTaylor展開します。


 \displaystyle u_{i-1} = u_{i} - \Delta x \frac{u'_{i}}{1!} + (\Delta x)^2 \frac{u''_{i}}{2!} + O[(\Delta x)^3]
 \displaystyle u_{i-2} = u_{i} - 2\Delta x \frac{u'_{i}}{1!} + (2\Delta x)^2 \frac{u''_{i}}{2!} + O[(\Delta x)^3]

二階微分の項を消すために  4u_{i-1} から  u_{i-2} を引きます。すると

 \displaystyle 4u_{i-1} - u_{i-2} = 3u_{i} - 2\Delta x u'_{i} + O[(\Delta x)^3]

となります。これを  u'_i について解くと

 \displaystyle u'_{i} = \frac{3u_{i} - 4u_{i-1} + u_{i-2}}{2\Delta x} + \frac{O[(\Delta x)^3]}{\Delta x} = \frac{3u_{i} - 4u_{i-1} + u_{i-2}}{2\Delta x} + O[(\Delta x)^2]

となります。つまり二次精度になっています。

同様にして次は、 V_i < 0 のときを考えていきます。まず  u_{i+1} u_{i+2} の値を  u_i を中心としてTaylor展開します。


 \displaystyle u_{i+1} = u_{i} + \Delta x \frac{u'_{i}}{1!} + (\Delta x)^2 \frac{u''_{i}}{2!} + O[(\Delta x)^3]
 \displaystyle u_{i+2} = u_{i} + 2\Delta x \frac{u'_{i}}{1!} + (2\Delta x)^2 \frac{u''_{i}}{2!} + O[(\Delta x)^3]

二階微分の項を消すために  4u_{i+1} から  u_{i+2} を引きます。すると

 \displaystyle 4u_{i+1} - u_{i+2} = 3u_{i} + 2\Delta x u'_{i} + O[(\Delta x)^3]

となります。これを  u'_i について解くと

 \displaystyle u'_{i} = \frac{-3u_{i} + 4u_{i+1} - u_{i+2}}{2\Delta x} + \frac{O[(\Delta x)^3]}{\Delta x} = \frac{-3u_{i} + 4u_{i+1} - u_{i+2}}{2\Delta x} + O[(\Delta x)^2]

となります。つまり二次精度になっています。

まとめると


 \displaystyle 
\left.V\frac{\partial u}{\partial x}\right|_i = 
\begin{cases} 
\displaystyle V_i\frac{3u_{i} - 4u_{i-1} + u_{i-2}}{2\Delta x} & (V_i \geq 0)\\ 
\displaystyle V_i\frac{-3u_{i} + 4u_{i+1} - u_{i+2}}{2\Delta x} & (V_i < 0)
\end{cases}

となります。これが二次精度の風上差分です。