輸送過程

図 1: 1D SHASTA の輸送過程(Boris and Book(1973),図1を基に作成)

SHASTA の輸送過程は以下のように与えられる. (Fig.1参照).

1. 隣接する格子点間で台形を作り(Fig.1(a)), その面積が 保存するように質量を輸送する.

   時刻 $t_{0} (\equiv n\Delta t)$ のときそれぞれ格子点上で密度は $\rho _ {i}^{n}, \rho _{i+1}^{n}$, 速度は $u_{i}^{n}, u_{i+1}^{n}$ とする. 格子間隔は $\Delta x$, 時間間隔は $\Delta t$ とする. 台形の面積を一定に保つので, その高さは底辺の長さに半 比例して与えられる. よって時刻 $t_{0}+\Delta t$ の台形の高さ (Fig.1(b))はそれぞれ,

   
   

   $$\rho _{1} = \frac{\rho _{i+1}^{n}\Delta x}{\Delta x + \Delta t (u_{i+1}^{n} - u_{i}^{n})},$$ (2')

   $$\rho _{2} = \frac{\rho _{i}^{n}\Delta x}{\Delta x + \Delta t (u_{i+1}^{n} - u_{i}^{n})},$$

   
   

   と与えられる. $\rho _{1}, \rho _{2}$ は $\rho _ {i}^{n}, \rho _{i+1}^{n}$ が正のとき及び

   
   

   $$\left\vert u \frac{\Delta t}{\Delta x}\right\vert < \frac{1}{2},$$ (3')

   
   

   の場合に常に正の値をとる. 以下ではこの条件が満たされるものとして 手続きを進める.

2. 輸送された台形を格子を含むセル上に分割し, それぞれのセルに含 まれる面積を格子点上に割り付ける. 具体的にはFig.1(c) の A は $i$に, B は $i+1$ に割り当てる. 割り当てた面積を格子間隔 で割った値が次ステップの値になる. セルの境界上の値は線形補間で求 める.

   Fig.1(c)A の面積は

   
   

   $$\frac{1}{2}\left[\rho _{1} + \left(\rho _{2}+\frac{\rho _{1... ...t)\right]\cdot (\frac{1}{2}\Delta x + u_{i+1}^{n}\Delta t),$$
   
   

   B の面積は

   
   

   $$\frac{1}{2}\left[\rho _{2} + \left(\rho _{2}+\frac{\rho _{1... ...ght)\right]\cdot (\frac{1}{2}\Delta x - u_{i}^{n}\Delta t),$$
   
   

   と与えられる. $i$ 格子に割り付けられる面積はFig.1(c)B と $i-1$ からの寄与

   
   

   $$\frac{1}{2}\left[\rho _{2} + \left(\rho _{3}+\frac{\rho _{2... ...ght)\right]\cdot (\frac{1}{2}\Delta x + u_{i}^{n}\Delta t),$$
   
   

   である. ただし,

   
   

   $$\rho _{3} = \frac{\rho _{i-1}^{n}\Delta x}{\Delta x + \Delta t (u_{i}^{n} - u_{i-1}^{n})},$$
   
   

   である. 割り付けられた面積を格子間隔で割って時刻 $t_{0}+\Delta t$ の$i$ 格子点での密度 $\rho _{i}^{n+1}$ を得る.

以上をまとめると,

$$\rho _{i}^{n+1} = \frac{1}{2}Q_{-}^{2}(\rho _{i-1}^{n}-\rho _{i}^{n}) + \frac{1}{2}Q_{+}^{2}(\rho _{i+1}^{n}-\rho _{i}^{n}) + (Q_{+}+Q_{-})\rho _{i}^{n},$$

$$Q_{\pm } = \left(\frac{1}{2} \mp u_{i}^{n}\frac{\Delta t}{\Delta x} \right)\... ...eft[1 \pm (u_{i\pm 1}^{n} - u_{i}^{n}) \frac{\Delta t}{\Delta x}\right.\right],$$ (4')

となる. 速度一定の場合は $\epsilon \equiv u_{i}^{n}\Delta t/\Delta x$ として,

$$\rho _{i}^{n+1} = \rho _{i}^{n} - \frac{\epsilon }{2} (\rh... ...{2}\right) (\rho _{i+1}^{n}-2\rho _{i}^{n}+\rho _{i-1}^{n}),$$ (5')

となる. これは中心差分で表現した移流に拡散が加わった形になっていて, 速度に依存しない拡散を除けば 2 段の Lax-Wendroff スキームと同じである ¹.