Microsoft PowerPoint - 第5章_H27(MAC).ppt [互換モード]

講義予定

1. 第1回目 10/2 阿部 数値シミュレーションの手続き 2. 第2回目 10/9 阿部 差分方程式と差分スキーム 3. 第3回目 10/16 田中 方程式の代数化，連立1次方程式の解法 4. 第4回目 10/23 田中 並列計算法 5. 第5回目 10/30 阿部 MAC法など差分の計算方法 6. 第6回目 11/13 田中 有限要素法 7. 第7回目 11/20 田中 有限要素法 8. 第8回目 12/4 田中 有限要素法，N-Sプログラムによる実習 9. 第9回目 12/11 阿部 乱流・熱流体・多相流 10. 第10回目 12/18 田中 津波解析など事例紹介 11. 第11回目 12/25 (試験予定日)

流体運動の基礎方程式

質量保存式 ２次元非圧縮性流れのNavier-Stokes方程式                           2 2 2 2 1 y u x u x p y u v x u u t u                             2 2 2 2 1 y v x v y p y v v x v u t v _  密度 動粘性係数 圧力 : : : p   0       y v x u 方向流速 方向流速 _v:y x : u

渦度と流れ関数

渦度 y u x v        渦度輸送方程式                       2 2 2 2 y x y v x u t       流れ関数（コーシーリーマンの関係式） x v y u         ,     _      2 2 2 2 y x 流れ関数についてのPoisson方程式

渦度輸送方程式の導出

非保存系表示のNavier-Stokes方程式                         2 2 2 2 1 y u x u x p y u v x u u t u _                          2 2 2 2 1 y v x v y p y v v x u u t v _  ： 渦度輸送方程式 上式をyで微分した式から、下式をxで微分した式を差し引く。 ： 渦度                       2 2 2 2 y x y v x u t       y u x v        ここで

流れ関数についてのPoisson方程式の導出

を、以下の渦度の定義式に代入する y u x v        流れ関数の定義式（コーシーリーマンの関係式） x v y u        ,            2 2 2 2 y x 流れ関数についてのPoisson方程式 より 2 2 2 2 , x x v y y u              渦度移動方程式を１次元とし、 として書き換えると

渦度移動方程式

              2 2 x Re 1 x u t    Re / 1                         2 2 2 2 y x y v x u t       対流拡散方程式となる。 ２次元渦度移動方程式

圧力方程式

保存系表示のNavier-Stokes方程式                           2 2 2 2 2 ₁ y u x u x p y uv x u t u _                            2 2 2 2 2 1 y v x v y p y v x uv t v _  ： 圧力pについてのPoisson方程式 上式をxで、下式をyで微分して足し合わせる。 p S y p x p _      2 2 2 2                                                                2₂ 2₂ p y D x D t D x v y u y v x u 2 S   y v x u D       _{： 発散}

圧力についてのポアソン方程式

ここで  S y p x p       2 2 2 2 p 2 2 2 2 2 2 S y x y x 2 S                                          圧力pについてのPoisson方程式を流れ関数を用いて表すと、

計算フローチャート

速度場 コーシーリーマンの関係式 圧力 圧力に関するポアソン方程式 流れ関数表示のソース項 流れ関数 流れ関数に関するポアソン方程式 渦度 渦度移動方程式 流速分布                       2 2 2 2 y x y v x u t       y u x v           _      2 2 2 2 y x  p S y x y x S                                             2 2 2 2 2 2 2      S y p x p _      2 2 2 2 n n _v u , 1 1_,   n n _v u p

計算手順

i. t=0(n=0)での初期値から初めてt=ntまで計算が完了して いるとする。したがって、変数値un_{, v}n_{, }n_,n_{は既知である。} ii. 渦度移動方程式にしたがって、内点の新しい変数値n+1_を 計算する。 iii. 内点での新しいn+1_{を用い、流れ関数方程式にしたがって} 内点の新しいn+1_{を反復収束過程から計算する。} iv. 新しい値n+1_{を用いて新しい流速u, vを} から計算する。 i. 内点における新しいn+1_、n+1_{を用い、境界条件から新し} い境界値を計算する。 y / un1n1  x / vn1n1  下図に示したように境界B₁, B₂, …, B₆がある。 領域（ 、 ）は流れに置かれた 障害物である。障害物後方の流れを計算しなさい。

問題設定

1 2 x₀ y₀ y x B₁ B₂ B₄ B₃ B₅ B₆ 0 0x x 0 0 yy

基礎方程式と基本関係式

バックステップ流れを渦度と流れ関数で記述すると基礎方 程式は                       2 2 2 2 y x y v x u t       x v y u         ,     _      2 2 2 2 y x である。 y u x v       

境界条件1

B₁は中心線とし、B₂とB₅はそれぞれ障害物の上面と底面とする。 B₃は上方境界、B₄は上流境界、B₆は流出境界と呼ばれる。 境界B₂-B₅-B₁は流線であるから=Constである。 記述の便利のために



B₂B₅B₁



0  と書く。上流境界B₄では一様な流速u₀で流れが体系に入る。 つまり次の条件である。















 





0, ,



0



1



1 , , 0 1 , , 0 0 0 0 0 0 0           y y t y y y y y u t y y y u t y u   u _y     y u x v        1 2 x₀ y₀ y x B₁ B₄ B₃ B₅ B₆ B₂ u₀

境界条件２

障害物の上面B₂ではすべりなしの壁とする。つまり、



x,y₀



0 















0





0



0 0 0 0 0 , 0 / , , y y y x y y x y x v y x          



x,y0

 

v x,y0



0



0 x x0



u     B₂のうえでは v/x0である。したがって



x,y₀



u



x,y₀



/y



0xx₀



 境界B₂の上では流線の条件式より 同じように障害物の底面B₅では 1 2 x₀ y₀ y x B₁ B₄ B₃ B₅ B₆ B₂ u₀

境界条件３

 

x,1 u

 

x,1 /y



0 x2



 境界B₁は中心線であるから、 u/y0, v0 したがって、 上方境界B₃もすべりなしの壁とした。この境界も流線であるか ら=Constである。B₄における流線の条件より

 

x,0 0



x0 x2



 B₄における速度の条件より

 

x,0 0



x₀  x2





 

x,1 u₀



1y₀





0x2



 また、B₂における条件より 1 2 x₀ y₀ y x B₁ B₄ B₃ B₅ B₆ B₂ u₀

計算フローチャート

速度場 コーシーリーマンの関係式 圧力 圧力に関するポアソン方程式 流れ関数表示のソース項 流れ関数 流れ関数に関するポアソン方程式 渦度 渦度移動方程式 流速分布                       2 2 2 2 y x y v x u t       y u x v           _      2 2 2 2 y x  p S y x y x S                                             2 2 2 2 2 2 2      S y p x p       2 2 2 2 n n _v u , 1 1 ,   n n v u p

渦度移動方程式の差分表現

渦度移動方程式のFTCS近似



fux fuy visx visy



t n j i n j i,1 ,      x 2 u fux n j , 1 i n j , 1 i j , i         y 2 v fuy n 1 j , i n 1 j , i j , i _   _  _   2 , 1 , , 1 2 x visx n j i n j i n j i         2 1 , , 1 , 2 y visy n j i n j i n j i        

SOR法によるPoisson方程式の数値解析

この式は連立一次方程式となる。これをSOR法で解くこととす る。つまり、第k回目の反復値を(k)_とすると    _      2 2 2 2 y x 流れ関数についてのPoisson方程式    







   



   







  n



j i k j i k j i k j i k j i k j i k j i k j i x2 _, , 2 1 1 , 1 1 , 2 1 , 1 1 , 1 2 , 1 , 1 2 1 2                       _ _ _ _  y x    / 

Poisson方程式の差分スキーム（1/6）

Poisson方程式 g y u x u 2 2 2 2        空間に中心差分近似を用いた差分方程式 j i j i j i j i j i j i j i g y u u u x u u u , 2 1 , , 1 , 2 , 1 , , 1 2 2 _            j i j i j i j i j i j i u u u u x g u , 2 1 , 1 , , , 1 , 1  4     _{   } x=yとする。

Poisson方程式の差分スキーム（2/6）

u1,1, u1,2,……, u3,3をu1, u2, ……, u9と番号を付け直す x y 1 2 3 4 5 6 7 8 9 11 12 13 14 15 16 17 18 19 21 22 23 24 25 20 10 x y (1,1) (1,2) (1,3) (2,1) (2,2) (2,3) (3,1) (3,2) (3,3) (0,0) (1,0) (2,0) (3,0) (4,0) (0,1) (4,1) (0,2) (4,2) (0,3) (4,3) (0,4) (1,4) (2,4) (3,4) (4,4) Natural ordering:

Poisson方程式の差分スキーム（3/6）

Natural ordering によって差分式は以下のように 表すことができる                                                                                                                                                              9 8 7 6 5 4 3 2 1 24 20 9 2 23 8 2 22 17 7 2 19 6 2 5 2 16 4 2 18 13 3 2 12 2 2 15 11 1 2 9 8 7 6 5 4 3 2 1 4 1 1 1 4 1 1 1 4 1 1 4 1 1 1 4 1 1 1 4 1 1 4 1 1 1 4 1 1 1 4 b b b b b b b b b u u g x u g x u u g x u g x g x u g x u u g x u g x u u g x u u u u u u u u u

Poisson方程式の差分スキーム（4/6）

以上の置き換えによって、行列式は、以下の形に書き換わる                                3 2 1 3 2 1 33 32 23 22 21 12 11 B B B U U U A A A A A A A        T  T  T T T T b b b B b b b B b b b B u u u U u u u U u u u U 9 8 7 3 6 5 4 2 3 2 1 1 9 8 7 3 6 5 4 2 3 2 1 1 , , , ,      







i j i j



i A A A ji ij ii                                  : 3 , 2 , 1 , 3 , 2 , 1 1 1 1 4 1 1 4 1 1 4 ここで

疎な行列(Sparse Matrix )

変数の数：格子分割の数：n×n=nn分割とすると2 係数行列の全要素数： n2_{× n}2 _{= n}4 非零の要素数の割合： (5n-4)/ n3 非零の要素数： (n+2(n-1))×n+n×(n-1)×2=(5n-4)n n=1 100 % n=2 75 % n=5 16.8 % n=10 4.6 % n=20 1.2 % n=50 0.1968 % n=100 0.0496 %                                3 2 1 3 2 1 33 32 23 22 21 12 11 B B B U U U A A A A A A A   i j i j A A i A ji ij ii                                  , 3 , 2 , 1 , 1 1 1 3 , 2 , 1 4 1 1 4 1 1 4

SOR法１

Gauss-Seidel法を書き、その左辺を仮の反復値 と見る。 すなわち、  1 ~k i u       i i i j k j j i i j k j j i k i i i u a u a u b a           1 1 , 1 1 1 , 1 ,~ と定義する。反復値 k1 は重みとして i u    



   k



i k i k i k i u u u u 1   ~ 1   1

_

_

  ~ 1 1   _{  } k i k i k i u u u   あるいは変形して は緩和係数と呼ばれる。

SOR法2

                        _{ }       k j i i i n i j k j j i i j k j j i k i i i k i i i u a u a u a u b a u a _, 1 , 1 1 1 , , 1 ,  ＞１を過緩和といい、 ＜１を不足緩和という。 上式から ~k1を消去すると、 i u この表現をSOR法という。行列Aを分離すると、この式のベク トル形式は



_D



_Euk1_

_

_

_

_

_DF

_

_u(k)_b 1 である。

SOR法3

と定義するとこの式は  



_L

 





_{I U}



_u 



_L



_{D b} uk1  1 1 1  k  1 1 1 となる。この式を整理して、 F D U E D L 1 _,  1    

_

_



L

 





I U



b D L u uk k                   1 1 1 1 1 1 1 L L が反復法の行列形式による表現である。 はSOR行列と呼 ばれ、 はGauss-Seidel行列を意味している。  L 1   L

SOR法4

よって、SOR法は残差修正法である。              _{   }           k j i i i n i j k j j i i j k j j i k i i i k i i i u a u a u a u b a u a , 1 , 1 1 1 , , 1 ,  この式中において



1   D R とし、右辺第2項の[ ]の中の1，2，3項は  k

u

A

b



~

とかける。こうすると上式は、 k  _u k _R



_{b Au} k



u 1    ~ Program SOR

境界条件の差分表現１

すべりなしの壁を考える。つまり、次の条件 , 0 , 0 _{, 0} 0 ,j  i j  i v u

が指定される壁である。(i,j0+1)の流れ関数_i,j0+1を点(i,j0)で Taylor展開すると             2 0 , 2 2 0 , 0 , 1 0 , 2 1 y y y y j i j i j i j i    

で あ る 。 壁 面 (i,j0) は ひ と つ の 流 線 を 形 成 し 、_i,j0=const (i=1,2,･･･)である。B₂-B₅-B₁も流線であって、この場合

境界条件の差分表現２

すべりなしの壁では、 0 / 0 , 0 ,j   _{i j}  i y u  となる。しかし、y方向には流速uは分布をもつ。つまり 0 / / _,0 0 , 2 2 _{   }  _ij j i u y y  である。ところで、壁面での渦度 は、 であるから y u x v      / /  0 / _{, 0} v x_ij y u    / 

となる。すべりのない壁ではu_i,j0=v_i,j0=0であるが、「渦度の生 成」はある。

境界条件の差分表現３

これらの式を流れ関数をTaylor展開した式に代入すると、 2 1 0 , 0 ,j 2 i j / y i      となる。一般に壁面を点(i,w)とすると、この点の境界条件は、 である。ところで、壁面での渦度 は、 であるから y u x v      / /  0 / _,0 v x_ij y u    / 

となる。すべりのない壁ではu_i,j0=v_i,j0=0であるが、「渦度の生 成」はある。ここでn2_{は壁面に垂直にはなれた次の格子点ま} での距離である。





2 1 , 1 , ,w 2 iw iw / n i        2 0 , 1 0 , 2 1 y j i j i      あるいは

境界条件の差分表現４

x y A B i j j₀ i₀ j₀+1 i₀+1 鋭い角(i₀,j₀)について考える。鋭 い角の流れ関数 の評価 は心配ない。 鋭い角の渦度 の指定に は以下の式が用いられる。 0 , 0 j i  0 , 0 j i  • 不連続近似 2 0 , 1 0 2 1 0 , 0 / , 2 / 2 _{i j} y _{B i j} x A          • 対称近似 2 0 , 1 0 2 1 0 , 0 / 2 / 2 i j y i j x A         • 壁面平均値近似 2 0 , 1 0 2 1 0 , 0j / y i j / x i A      

境界条件の差分表現５

右図のように格子状に壁を設定 しないが点(i,j+1/2)は壁として表 現可能である。この境界をすべ りのない壁とする。つまり、 u_i,j+1/2=0である。図中のui,j+1に 対して境界点(i,ｊ)について、 1 , ,j  i j i u u j i j+1 i+1 1 j , i u  i-1 i+2 j-1 j+2 j , i u  とおけば点(i,j+1/2)ではu_i,j+1/2=0 となってすべりのない壁の近似 となる。

時間ステップ

tの制御法

Courant条件 1 /  t x c 拡散数の条件（：拡散係数） 2 / 1 / 2 t x  変係数の熱伝導方程式                x u u x t u _ FTCSスキームの安定条件

 

ut/x21/2  t₀の決定法







_{2 2}



1 0 0 / 2 /        t u x  x y

計算フローチャート

速度場 コーシーリーマンの関係式 圧力 圧力に関するポアソン方程式 流れ関数表示のソース項 流れ関数 流れ関数に関するポアソン方程式 渦度 渦度移動方程式 流速分布                       2 2 2 2 y x y v x u t       y u x v           _      2 2 2 2 y x  p S y x y x S                                             2 2 2 2 2 2 2      S y p x p _      2 2 2 2 n n _v u , 1 1_,   n n _v u p

圧力・流速場(u,v,p)を用いる解法

 Los Alamos国立研究所グループT-3(米国)

• MAC法(Marker and Cell method): (Hallow and Welch,1965; Welch et

al.,1966 )

• SMAC法(Simplified Marker and Cell method): (Amsden and

Hallow,1970)

• FS法(Frictional Step method または Two Step Projection method):

(Kothe et al., 1992)

• HSMAC法(Highly Simplified MAC method)またはSOLA法: (Hirt et al.,

1975)

• ICE法 (Hallow and Amsden, 1971)

• ALE法 (Hirt et al., 1974)

• VOF法 (Hirt and Nichols, 1981)

 Imperial Collage(英国) • SIMPLE法

基礎方程式：

連続の式とNavier-Stokes方程式

0    V --- ① --- ② --- ③ --- ④     F g p V V t V   _               1 1 ) ( 差分式(流速と圧力のみを陰的に取り扱う) 基礎方程式 0 1_   n V n n n n n n n n n _F g p V V t V V                         _{1 1} ) ( 1 1

MAC法

0 1_   _Vn ) 1 ( 1 ) ( 2 1 2 1 n n n n n n n n n _F g p V V t V V                               差分式④の両辺の発散(divergence,  )をとる 連続の式③を満足する圧力を _pn1 とすると、このとき であるから、 ) 1 ( ) ( 2 1 2 n n n n n n n n n n F g V V V t p       _ _{ } _ _{   } _      時刻nにおいて、右辺は全て既知であるから、上式は についてのポアッソン方程式となる。一旦 が求まると、 1  n p 1  n p ] 1 1 ) ( [ 1 1 n n n n n n n n n F g p V V t V V                          より、_Vn1 が決定できる。

MAC法の解法手順

初期値：_Vn ポアソン方程式より _pn1を求める 速度_Vn1を求める ) 1 ( ) ( 2 1 2 n n n n n n n n n _{V VV g} F t p                            ] 1 1 ) ( [ 1 1 n n n n n n n n n _{V t VV p g} F V           _{    }  _{   }

SMAC法

n n n n n n n n n n F g p p V V t V V                            1 ) ( 1 1 ) ( 1 上式の両辺の発散(divergence,  )をとる 連続の式_ n1_0 が満足される圧力を とすると、 V pn1 P p pn1 n 求めるべき圧力場を と分解すると、 上式を、陽的部分と陰的部分に分解する。 n n n n n n n n F g p V V t V V                      1 1 ) ( ~ ) ( 1 ~ 1 p t V V n n                             ) ( 1 ~ 1 p t V V n n    V~ t ) p ( n 2 _{ }     (陰的部分) (陽的部分)

SMAC法の解法手順

を求める p p pn1 n 速度_Vn1を求める 初期値：_Vn pn を求める V~ ポアソン方程式より p を求める n n n n n n n n F g p V V t V V                      1 1 ) ( ~ 2 V~ t ) p ( n 2 _{ }     ) ( 1 ~ 1 p t V V n n         

FS法 (Two Step Projection法)

上式の両辺の発散(divergence,  )をとると、 連続の式が満足されているとするとVn10であるから このポアッソン方程式を解くことによって が求まり、 最終的に、速度場が求められる。 1  n p を、直接陽的部分と陰的部分に分解する。 n n n n n n _F g V V t V V                   1 ) ( ~ ) ( 1 ~ 1 1        n n n p t V V                  ) ( 1 ~ ₁ 1 n n n _{V t p} V   V~ t ) p ( n 1 n 2 _{ }     (陰的部分) (陽的部分) n n n n n n n n n _F g p V V t V V                         _{1 1} ) ( 1 1

FS法の解法手順

速度_Vn1を求める 初期値：_Vn pn を求める V~ ポアソン方程式より _pn1を求める n n n n n n _F g V V t V V                   1 ) ( ~ V~ t ) p ( n 1 n 2 _{ }                    ) ( 1 ~ 1 1 n n n p t V V   これは、変数 が質量保存式を満足していない程 度を表す指標となりうる。 D<0ならば、このセルに向かって正味で質量の流入が生じる ことに対応し、このセルにおける圧力を低下させる必要がある。 逆に、D>0ならば、このセルより正味で質量の流出が生じる ことに対応し、このセルにおける圧力を増加させる必要がある。 このように、圧力を調整してDを0に調節すれば、最終的な速 度場が得られることがわかる。

HSMAC法(SOLA法)

セル(i,j)における発散D_i,j を定義する



 

1



1 , 1 , 1 , 1 1 , 1 , 1 1        _      n j i n j i n j i n j i n j i v v y u u x D

   

1 , 1 , , inj n j i v u ) ( 1 1 p t V V n k k          p t V _{ }    1  差分式(流速と圧力のみを陰的に取り扱う)     F g p V V t V                  1 1 ) (

HSMAC法(SOLA法)アルゴリズム

圧力変化分 pが求められれば、速度場_Vk1 が求められる n n n n n n n n F g p V V t V V                      1 1 ) ( ~ ) ( 1 ~ 1 1        n n n p t V V  

流速の更新

) ( 1 1 p t V V n k k                              _v  _{t p}  _y v y p t v v x p t u u x p t u u k k j i k j i k k j i k j i k k j i k j i k k j i k j i                         / / / / 1 , 1 1 , , 1 , , 1 1 , 1 , 1 ,     p(k)_{が計算されると、速度場は、} より p_n p とすると、   _

流速の更新







1



1 , 1 , 1 , 1 1 , 1 , 1 1        _{   } k j i k j i k j i k j i k j i v v y u u x D                       _v  _{t p}  _y v y p t v v x p t u u x p t u u k k j i k j i k k j i k j i k k j i k j i k k j i k j i                         / / / / 1 , 1 1 , , 1 , , 1 1 , 1 , 1 ,     速度場の評価値、 を発散 の式に代入すると                                                  2 2 1 1 , 1 , 1 , 1 1 , 1 1 2 2 1 2 1 1 1 y x t y t y x t x p v p v y p u p u x p D ikj k j i k j i k j i           

圧力の補正

0 ) 1 (k _ D

いま、セル(i,j)における発散D_i,jをそのセルの圧力p_i,jの関数と みなすことにする。つまり、D_i,j=D(p_i,j)である。D=0の根を求める ためにNewton法を用いると、        k j i j i k j i k j i k j i j i p D D p p p _              , , , , 1 , , /  ここで、kは第k番目の反復を表す。まず質量保存式をDとおき、 を求めてから同式に代入する。求め られたD_i,jをp_i,jで偏微分すれば次式のようになる。 1 1 , 1 , 1 , 1 1 , , in j, inj , inj n j i u u u u                       2 1₂ 1₂ y x t D p  ) ( ( 1) () ) ( ) 1 (k k k k D D D p p p             

流速の更新

ここで加速係数を取り入れて、早く収束させる。                     _v  _{t p}  _y v y p t v v x p t u u x p t u u k k j i k j i k k j i k j i k k j i k j i k k j i k j i                         / / / / 1 , 1 1 , , 1 , , 1 1 , 1 , 1 ,     この式にD(k)_{を代入してp}(k)_{が計算されると、発散を0とすべく} 流速が更新される。     1 1 _{1 2} 2 _{2 2} 1 _{ }                        _{ }  y x t D pk k

対流項の差分

対流項 fuxn, fuyn, fvxn, fvyn



















n



j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n uy u u v v u u v v v v v v u u v v y f , 1 , 1 , 1 1 , , 1 , 1 , 1 1 , 1 , , , 1 , 1 , , , 1 , 4 1         _{    }                 











_



_{ }



_



_    _{   }          n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n ux u u u u u u u u u u u u x f , , 1 , , 1 2 , , 1 , 1 , , 1 , 2 , 1 , 4 1  

対流項の差分



















n



j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n vx v v u u v v u u v v u u v v u u x f , 1 , 1 , 1 1 , , 1 , 1 , 1 1 , 1 , , , 1 , , 1 , 1 , , 4 1         _{    }                 











_



_{ }



_



_    _{   }          n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n j i n vy u u u u u u u u v v v v y f , 1 , , 1 , 2 , 1 , 1 , , 1 , , 2 1 , , 4 1   ここで重みaは、=0のときは中心差分となり差分近似の精度 が高くなり、=1の時には風上差分となり安定性が確保されや すくなるという効果をもつ。

ドナー・セル差分(風上差分)

変化の激しい現象を扱う場合、 , 0 , 0 0 , 0 ₁1_{/2 1/2 1/2} 2 / 1         in in in n i u u u u または であった場合、数値的不安定性が発生しやすい。たとえば、















0



0 2 / 1 2 / 1 1 1 2 / 1 2 / 1 1 1                     u u u x t u u u u x t u n i n i n i n i n i n i n i n i n i n i n i n i         とする。つまり流速の向きを考えてドナーが風上にいるように する。上式では(i-1)がドナーで(i)がアクセプターとなる。この差 分をドナー・セル差分という。

粘性項の差分









                in j inj in j inj inj inj n visx u u u y u u u x f _{, 1 , , 1} 2 , 1 , , 1 2 2 1 2 1  粘性項 n visy n visx f f ,









                _{  } n _ j i n j i n j i n j i n j i n j i n visy v v v y v v v x f _{, 1 , , 1} 2 , 1 , , 1 2 2 1 2 1 