演習問題（１）解答
演習問題（１）解答
棒の格子点2，3の温度を中心差分により求めなさい。
１）基本問題
u=100
d 2u
dx 2
d 2u
dx 2
u=0
1
3
2
4
 0 :
2
0
3
 0 :
棒の格子点2，3の温度を中心差分により求めなさい。
100
u3  2u2  u1
0
(x) 2
２）応用問題
u
0
x
u=100
u4  2u3  u2
0
(x) 2
1
 2 1  u2   100
 1  2 u    0 

 3  

2
3
4
注意：境界条件も離散化する
 200 
u
 2  1  2  1  100  3 


  

u3  3   1  2  0   100 
 3 
n 
j
u
c
x
n 
j
演習問題（２）解答
  1/ 2
0
Θ＝1/2、1としたときの格子点jに関する離散化式を示しなさい。
j-2
j-1
  1/ 2
j+1
j
u nj 1  u nj
u nj 1  u nj
j+2
t
n 1
 1   c
u nj1  u nj1
2x
 c
u nj11  u nj11
2x
n
n 1
1
1
1
1 c n
  u nj11  u nj 1  u nj11  u nj1  u nj  
u j 1
4
4
4
4 x
 1
u nj 1  u nj
t
c
u nj11  u nj11
2x
0
Crank-Nicolson法
（陰解法）
0
1 u j 1  u j 1 1 u j 1  u j 1
 c
 c
0
t
2
2x
2
2x
c n 1 n 1 1
c n 1 1
c n
c n
1
1
u j 1  u j  t
u j 1  t
u j 1  u nj  t
u j 1
 t
4 x
4 x
4 x
4 x
n
u3  2u2  u1
d u
0
 0 :
2 2
(x) 2
dx
d 2u
u4  2u3  u2
0
 0 :
2 3
(x) 2
dx
u
3
d 2u
u  2u4  u3

0
:

0
4
dx 2
(x) 2
中心差分
du
u3  u
 0 u  u
4  0 :

3
2x
dx
0  u2   100
 2 1
 1  2 1  u    0 


 3  



 0
2  2 u4   0 
演習問題（２）解答
u
t
100
2
1
1
c n 1
c n 1
 t
u j 1  u nj 1  t
u j 1  u nj
2 x
2 x
1
1
  u nj11  u nj 1  u nj11  u nj
2
2
 1
陰的Euler法
（陰解法）


 1
4
 


1
4
 4
・・





4





n 1
n


  u1 




 


 
・
 u j   RHS


   
1
 



4 u 


 n 

1

4



 1
2
 


1
2
 2
・・





2





n 1
n


  u1 




 


 
・
 u j   RHS


   
1
 


2  u 


 n 

1

2

連立一次方程式を解く必要がある→陰解法
Courant数が小さいと対角優位性が増す→収束性が増す（反復解法）
計算流体力学
（第3回）
物理法則に従った離散化
１）双曲型問題に対する風上化の必要性
u2
u
t
n
j
u
t

u
n 1
j
u
t
u
n 1
j
n
j
u
x
c

n
j
c
 u nj 1  u nj   (u nj1  u nj )
u1
u  u2 ( x)
ｘ
u
0
n
j
n 1
j
u
t
u u
n
j
 1 c 
 1
u ( x, t )  u ( x  ct )
差分法
誤った離散化
（例えば；時間に対して前進差分，空間に対して前進差分）
t
u1
t  0,  x  
u  u1 ( x  0)
初期条件
u
u2
u
u
c
0
x
t
物理法則に従った離散化
2014年10月10日
n
j 1
x
n
j
j-1
u
x
 O(t ) 時間：前進
n
j

x
j+1
j
u u
n
j
n
j 1
x
 O(x)
空間：後退
空間の1階微分：情報を時間と共にある方向に伝播する
対流（移流）を表す．
1階微分の係数が伝播する方向を表す．
＋：左から右へ伝播後退差分が適当
－：右から左へ伝播前進差分が適当
⇒風上差分法，風上（安定化）有限要素法
 u nj 1  u nj   (u nj  u nj1 )
x
 1 c 
 u nj 1  u nj1
t
0
 t 
 u  c (u nj  u nj1 )
 x 
n
j
厳密解と一致
Courant(クーラン)数
x
 u nj 1  2u nj  u nj1
t
厳密解と全く異なる
安定化手法
物理法則に従った離散化
２）放物型問題に対する中心差分，Galerkin法の有効性
u
 2u
 2  0
t
x
u
t
n
j
 2u
 2
x
n
j
双曲型問題に対する安定化
u
u
c
0
t
x
0
正しい離散化（時間に対して前進差分，空間に対して中心差分）
u
t

n
j

u nj 1  u nj
u nj 1  u nj
t
t

u nj 1  u nj  
u
x 2
2
 O(t )
u nj1  2u nj  u nj1
(x) 2
0
t
(u nj1  2u nj  u nj1 )
( x ) 2
n
j

u
n
j 1
 2u  u
n
j
n
j 1
(x) 2
1
t
k 
 
( x ) 2 2
u nj 1
 O((x) )
2
1
 (u nj1  u nj1 )
2
拡散現象⇒中心差分が適当
u
差分法
（陽解法） t
FTCS法
u
t
u
x
u
n 1
j
n
j
n
j


t
u nj1  u nj1
2x
t  0,  x  
c
u nj 1  u nj
u
x
n
j
 O(t )
0
時間：前進
j-1

 O(x)
空間：中心
u
j+1
j
n 1
j
u
t
n
j
c
u
n
j 1
u
2x
n
j 1
0
1  t 
 u  c (u nj1  u nj1 )
2  x 
Courant(クーラン)数
u
n 1
j
1
n 1
n
 u   (u nj1  u nj1 ) ⇒ u j  u j 発散の方向へ（不安定）
2
n
j
u  u2 ( x)
u ( x, t )  u ( x  ct )
  0.5
n
j
空間の２階微分：時間の経過と共に空間での値を緩やかにする拡散に対応する．
拡散現象：中心差分法，Galerkin有限要素法が適当
n
j
u  u1 ( x  0)
初期条件
x
安定化手法（Laxの方法）
安定化手法（1次精度風上差分）
u
t
n
j
u
t
u
x

c
n
j
n
j
u
x
n
j
u nj 1  u nj

 O(t )
t
n
u  u j 1

 O(x)
x
n
j
u nj 1  u nj
t
c
1
u nj 1  u nj   (u nj1  u nj1 )
2
u nj1  u nj1
n
(u j 
)
2
1
1
u nj 1  u nj   (u nj1  u nj1 )  (u nj1  2u nj  u nj1 )
2
2
0
u nj  u nj1
x
時間：前進
空間：後退
0
両辺を t で割り整理
u nj 1  u nj
 t 
u nj 1  u nj  c (u nj  u nj1 )
 x 
 u nj  (u nj  u nj1  u nj1  u nj1 )
t
u nj 1  u nj
1
1
 u nj   (u nj1  u nj1 )   (u nj1  2u nj  u nj1 )
2
2
u
n
j c
x
t
c
c
(u nj1  u nj1 )
2x
(u nj1  u nj1 )
2x
cx  2u
n
j 
2 x 2
n
j
u
t
0
n
j
c
u
x
n
j

n 1
 u j  t
n
u
t
n
j
1
 2u
 (t ) 2 2
2
t
n
j
u
x
n
j
 u
中心差分 
 t

n
j
n 1
 u j  ct
n
1
 2u
 c 2 (t ) 2 2
2
x

u nj1  u nj1
2x
n
j
 O ((t ) 3 )
 2u
, 2
x
1
1
 u nj 1  u nj   (u nj1  u nj1 )   2 (u nj1  2u nj  u nj1 )
2
2
u
t
u
n
j c
x

n
n
n
x 2 (u j 1  2u j  u j 1 )
2t
x 2
FTCS法
安定化項
c 2Δt  2u n
n
j 
j 0
2 x 2
人工拡散項
n
j
0
人工拡散項＝修正項の正体
安定化手法（まとめ）
いずれの方法もFTSC法による解を安定化させる修正項をもつ
 O((t ) 3 )
 u
 2u 
 
 
u
  u 

,    c   c   c 2 2 
 c x 
x  x 
x
t  t 
 t
uj
1
(u nj1  2u nj  u nj1 )
2t

x 2  2 u
2t x 2
安定化手法（Lax-Wendroffの方法）
uj
FTSC法による解を安定化させる修正項をもつ
安定化項
FTCS法
u
t
安定化項
FTCS法
n
j
u nj1  2u nj  u nj1 



x 2

1
1
u nj 1  u nj   (u nj1  u nj1 )  (u nj1  2u nj  u nj1 )
2
2
Laxの方法：
Lax-Wendroffの方法：
1
1
u nj 1  u nj   (u nj1  u nj1 )   2 (u nj1  2u nj  u nj1 )
2
2
1次精度風上差分
1
1
u nj 1  u nj   (u nj1  u nj1 )   (u nj1  2u nj  u nj1 )
2
2
上記の手法の中では、Lax-Wendroffの方法が最も人工拡散が弱い
（陽解法の場合ν≦１なので）
安定条件（これから安定条件について考える）
安定性解析
安定性解析（陽解法:FTCS法）
von Neumannの解析法
適用条件：等間隔格子、線形連立方程式、初期値問題
g n : 複素振幅率
i : 虚数
j : 格子点位置
   0 : 空間に対して一定　
：
振動成分 
   : jの添え字１つおきに符号を変える解

Neumannの安定性基準
 

2

3
g : 増幅率＞１：不安定
0
e i  e  i
2
e i  e i
sin  
2i
cos  
 g  1  i sin 
2
⇒
g  1   2 sin 2   1
宿題
Lax法およびLax-Wendroff法の安定性を調べよ。
g n 1eij  g n eij  ( g n e ij  g n ei ( j 1) ) 両辺を g n eij で割る
e i  e  i
2
e i  e i
sin  
2i
cos  
2
g  (1    cos  ) 2  ( sin  ) 2
 1  2 (1  )(1  cos  )

g n 1
 1  ( e i  e  i ) gn
2
 1  i sin 
 1   2 sin 2 
u nj  g n expi j 　　(1) :
 1  (1  cos  )  i sin 
両辺を g n eij で割る
＜１：安定（減衰）
u nj 1  u nj  (u nj  u nj1 ) (1)
g
 1  (1  e i ) gn
 1  (1  cos   i sin  )

g n 1e ij  g n e ij  ( g n e i ( j 1)  g n e i ( j 1) ) 2
g  12  ( sin  ) 2
0
n 1
n
＝１：安定（中立）
安定性解析（陽解法:1次精度風上差分法）
n
j
n
j
u  g expi j 　　(1) :
解の重ね合わせ⇒実際の解

u
x
u
x
n
j

c
c
2
u nj  g n expi  j 
n
j
n
j

(1)
FTCS法 u nj 1  u nj  (u nj1  u nj1 ) 差分式の解
u
t
u
t
条件付安定
⇒ g  1  2 (1  )(1  cos  )  1 となるためには 1   0    1
無条件に不安定
安定性解析（陰解法：完全陰解法）
安定性解析（陰解法：Crank-Nicolson法）
u
t
n 
c
j
u
x
n 
0
u nj 1  u nj

t
j
Crank-Nicolson法：   1
2
 n 1 n 1  n 1  n
 1   c
u nj1  u nj1
2x
 c
u nj11  u nj11
2x
0

4

u nj  g n expi j 　　(1) :
g n 1e i ( j 1)  g n 1e ij 

4
j
u
x
n 
0

u nj 1  u nj
t
j

g n 1e i ( j 1) 

4
g n e i ( j 1)  g n e ij 

4
g n e i ( j 1) 両辺を g n eij で割る
 i   ij
 i
g   ij
  e 1 e   e 1 e gn  4
4  4
4
1
1
1  (  sin  ) 2
1  i sin 
g n 1
2
2
2
g 
 1 ⇒  g  1 安定（中立）

1
g n 1  1 i sin 
1  (  sin  ) 2
2
2
安定性解析（陰解法）
完全陰解法
c
 1   c
u nj1  u nj1
2x
 c
u nj11  u nj11
2x
0

 u nj11  u nj 1  u nj11  u nj (1)
2
2
n 1
Crank-Nicolson法
n 
完全陰解法：   1
 u j 1  u j  u j 1  u j 1  u nj  u nj1　(1)
4
4
4
4

u
t


 1
4
 


1
4
 4
・・





4





n 1
n


  u1 


  


 





u
・
RHS
 j


   
1
 



4 u 


 n 

1

4



 1
2
 


1
2
 2
・・





2




a ii 
a ij

n 1
n


  u1 
i j


  
 であれば収束する

 
・
 u j   RHS


   
1
 


2  u 



n




1

2

例：Gauss-Seidel法

Courant数が下がると対角優位性が増す⇒収束性が増す（反復解法）


u nj  g n expi j 　　(1) :
2
g n 1e i ( j 1)  g n 1e ij 

2
g n 1e i ( j 1)  g n e ij 両辺を g n eij で割る
 ij 
g n 1   ij
  e  1  e   1　gn  2
2

2
g 
g n 1
1  i sin    1　gn
1
1
1
⇒ g 
2
1   sin 
1   2 sin 2 
2
無条件安定（減衰）

Download Report