PowerPoint プレゼンテーション

微分方程式の数値解法
•
微分方程式
•
•
•
•
•
常微分方程式
偏微分方程式
初期値問題と境界値問題
微分方程式の解析解
微分方程式の数値解
•
•
•
•
差分法（常、偏）
オイラー法（常）
ルンゲ・クッタ法（常）
残差法（常、偏）
（選点法、最小２乗法、モーメント法、ガレルキン
法）
•
微分方程式
常微分方程式
独立変数＝１個
2
d y
dy
 p ( x, y )  g ( x, y )  0
2
dx
dx
偏微分方程式
独立変数＞１個
2 z 2 z
 2  f ( x, y, z )
2
x y
•
初期値問題と境界値問題
初期条件
y ( x0 )  y0
y( x0 )  y0
 y n1 ( x0 )  y0n1
n 階微分方程式において、独立変数 x の特定な一つの点で数値
で与えられているときには、その値を初期条件といい、初期条
件下で解を見出す問題を初期値問題と呼ぶ。
境界条件
y( x0 )  y0
y( x1 )  y1
 y( xn1 )  yn1
独立変数 x の二つ以上の複数の点において条件が与えられてい
るときには、その条件を境界条件といい、境界条件下で微分方
程式を解く問題を境界値問題と呼ぶ。
• 微分方程式の解析解
微分方程
式：
解析
解：
y
dy
x
dx
1 2
y  x c
2
初期値
y x 0  1
定数
x
Ｃ：任意定数
定
解
y x 0  0
c 1
c0
y  0.5 x 2  1
y x 0  1
c  1
y  0.5 x 2  1
y  0.5 x 2
例：落体問題
方程
式：
d2y
g
2
dt
積分
で：
解析
解：
dy
 gt  c1
dt
y t 0  y0
dy
 v0
dt t 0
c2  y0
c1  v0
t 0
yi
t  ti
y
y  0.5gt 2  c1t  c2
初期条件：
y0
v0
Ｃ1,Ｃ2 ：任意定数
定解：
1 2
y  gt  v0t  y0
2
例：落体問題
件）
（境界条
y0
t 0
方程
式：
d2y
g
2
dt
yi
t  ti
y1
t  t1
解析
解：
y  0.5gt  c1t  c2
2
y
境界条件：
y t 0  y0
c2  y0
y t t1  y1
v0  c1
c1  ( y1  y0  gt / 2) / t1
2
1
定解：
1 2
y  gt  v0t  y0
2
微分方程式の数値解
d2y
dy
 p ( x, y )  g ( x, y )  0
2
dx
dx
初期条件：
y( x0 )  y0
境界条件：
y( x0 )  y0
y (x)
y0
微分方程式：
y( x0 )  y0
y
y( x1 )  y1
x
x0 xi xi 1
x  xi  xi 1  h
微分方程式の数値解法
•
•
•
•
差分法（常、偏）
オイラー法（常）
ルンゲ・クッタ法（常）
残差法（常、偏）
（選点法、最小２乗法、
モーメント法、ガレルキン法）
微分と差分法
f (x)
前進差分：
１階差
分：
y
yi  yi 1  yi
yi
dy

 ( yi 1  yi ) / h
dx x  xi
x
２階差
分：
2 yi  yi 1  yi
 yi  2  yi 1  ( yi 1  yi )
 yi  2  2 yi 1  yi
2 yi
d2y
2


(
y

2
y

y
)
/
h
i2
i 1
i
2
2
dx x  xi
x 
yi 1
dy
dx
yi
xi
y
x
xi 1
yi  f ( xi )
x  xi 1  xi  h
x
差分法による微分方程式の数値解
後退差分：
１階差
分：
yi  yi  yi 1
yi
dy

 ( yi  yi 1 ) / h
dx x  xi
x
２階差
分：
f (x)
y
yi 1
yi 1 yi
y
x
dy
dx
y
x
xi 1 xi
 2 yi  yi  yi 1
 yi  yi 1  ( yi 1  yi 2 )
 yi  2 yi 1  yi 2
 2 yi
d2y
2


(
y

2
y

y
)
/
h
i
i 1
i 2
2
2
dx x  xi
x 
xi 1
yi  f ( xi )
x  xi  xi 1  h
x
差分法による微分方程式の数値解
中心差分：
１階差
分：
yi  yi 1  yi 1
yi yi 1  yi 1
dy


dx x  xi
x
2h
２階差
分：
f (x)
y
2 yi  yi 1  2 yi  yi 1
2 yi
d2y
2


(
y

2
y

y
)
/
h
i 1
i
i 1
2
2
dx x  xi
x 
yi 1
yi 1 yi
y
x
dy
dx
y
x
xi 1 xi
xi 1
yi  f ( xi )
x  xi  xi 1  h
x
差分法による微分方程式の数値解
中心差分：
１階差
分：
f (x)
y
yi  yi h / 2  yi h / 2
yi 1
yi 1 yi
y
x
dy
dx
y
x
２階差
分：
 yi  yi h / 2  yi h / 2
xi 1 xi
2 yi  yi 1  2 yi  yi 1
yi  f ( xi )
x  xi  xi 1  h
2
2 yi
d2y
2


(
y

2
y

y
)
/
h
i 1
i
i 1
2
2
dx x  xi
x 
xi 1
x
常微分方程式の差分化
y
中心差分を用いる
例：
y  2 x  3sin x
0.5
常微分方程式：
y  y  2 x  0
境界条件：
1
x0  0
y0  0
xn  1
yn  0.52
y
0.5
差分化：
yi 1  2 yi  yi 1
 yi  2 xi  0
2
h
h  1/ n
x
xi  ih
i  1,2,, n  1
xi
x
y
yi 1  (h  2) yi  yi 1  2ih  0
2
3
i  1 y0  (h  2) y1  y2  2h  0
2
3
0.5
i  2  y1  (h  2) y2  y3  4h  0
2
3
i  3  y2  (h  2) y3  y4  6h  0
2
3
     
i  n  1 yn2  (h 2  2) yn1  yn  2(n  1)h3  0
3
 2  h2 1
 y

2
h
 y0 

 1  

3
 1 2  h2 1
 y2  
4h










    

3
 yn  2 
2
2
(
n

2
)
h


1 2  h 1 





 2(n  1)h 3  y 
2  y n 1 

n
1 2  h 

x
3
 2  h2 1
 y

2
h
 y0 



 1


2


3
 1 2  h 1
 y2
4h
















 2( n  2 ) h 3 


2
y

 1 2  h  1 n 2  



 2(n  1)h 3  y 
2  y n 1 

n
1 2  h 

n5
h  1/ n  0.2
y0  0
y5  0.52
1.96  1 0 0  y1   0.016 

  

  1 1.96  1 0  y2   0.032 
 0  1 1.96  1 y    0.048

 3  

 0 0  1 1.96  y4   0.584 
y
0.5
0.2 0.4 0.6 0.8 1
x
偏微分方程式の差分化
ラプラス微分方程
式：
 2  2
 2 0
2
x
y
差分化：
y
j 1
j
j 1
i 1 i i  1
 2
1
 2 ( i 1, j  2 i , j   i 1, j )
2
x
x
 2
1
 2 ( i , j 1  2 i , j   i , j 1 )
2
y
y
2

x
2  2
y
差分方程式：
 i 1, j  2(1   2 ) i , j   i 1, j   2 i , j 1   2 i , j 1  0
x
偏微分方程式の差分化
ラプラス微分方程
式：
 2  2
 2 0
2
x
y
差分化：
2

x
2  2
y
 i 1, j  2(1   2 ) i , j   i 1, j   2 i , j 1   2 i , j 1  0
上図の例では、境界条件によって、境界上の点の
値は既知であり、未知数は16個で、線形方程式も
16個があるので、唯一の解が得られる。
差分方程式の解
j 1
j
j 1
ラプラス微分方程
式：
2
   2
 2 0
2
x
y
i 1 i i  1
差分方程式：
 i 1, j  2(1   2 ) i , j   i 1, j   2 i , j 1   2 i , j 1  0
i, j
1
2
2










 i , j 1 
i 1, j
i 1, j
i , j 1
2
2(1   )
SOR法による解

1
(k )
(k )
2 (k )
2 (k )
 i, j 








 i , j1
2
i 1, j
i 1, j
i , j 1
2(1   )
 i(,kj 1)   i(,kj )  SOR  *i , j   i(,kj )
*



偏微分方程式の差分化
j 1
j
j 1
ラプラス微分方程
式：
 2  2
 2 0
2
x
y
i 1 i i  1
差分法によるラプラス微分方程式の実用計算式：
 i, j
*

1
(k )
( k 1)
2 (k )
2 ( k 1)









 i , j1
2
i 1, j
i 1, j
i , j 1
2(1   )

 i(,kj 1)   i(,kj )  SOR  *i , j   i(,kj )


微分方程式の数値解：オイラー法
１階常微分方程
式：
境界条件：
y
dy
 f ( x, y )
dx
y( x0 )  y0
テイラー展開：
h2
y ( x0  h)  y ( x0 )  hy( x0 )  y( x0 )  
2
y (x)
y0
x0
xi
xi 1
y0  y( x0 )
x  xi  xi 1  h
h2 以降を無視して近似解を求
める方法がオイラー法である
y( x0  h)  y0  hy( x0 )
yi 1  yi  f ( xi , yi )h
y( x0 ) 
dy
 f ( x0 , y0 )
dx x  x0
x
微分方程式の数値解：オイラー法
１階常微分方程
式：
境界条件：
dy
 f ( x, y )
dx
y( x0 )  y0
一方、テイラーを利用して
y ( xi )  y ( xi 1  h)
 y ( xi 1 )  hy( xi 1 )
yi 1  yi  f ( xi 1 , yi 1 )h
y
y (x)
y0
x0
xi
xi 1
y0  y( x0 )
x  xi  xi 1  h
x
１階常微分方程
式：
境界条件：
dy
 f ( x, y )
dx
y( x0 )  y0
オイラー法による数値解
yi 1  yi  f ( xi , yi )h
yi1  yi  f ( xi1 , yi1 )h
y
y (x)
y0
x0
xi
xi 1
y0  y( x0 )
x  xi  xi 1  h
k1  f ( xi , yi )h
k2  f ( xi  h, yi  k1 )h
yi 1  yi  w1k1  w2 k2
ルンゲ―・クッタ法
x
k1  f ( xi , yi )h
k2  f ( xi  h, yi  k1 )h
yi 1  yi  w1k1  w2 k2
k2をテイラー展開
f ( xi  h, yi  k1 )  f i  f xh  f y k1  0.5 f xx (h) 2  f xyhk1  0.5 f yy ( k1 ) 2
2
f
f

f
f i  f ( xi , yi ) f x  x x f y 
f xx  2
x  xi
x y  yi

y
x
i
y  yi
x  xi
y  yi
2 f
2 f
f yy  2
f xy 
y
x

x
xy y  yi
i
x  xi
y  yi
yi 1  yi  w1 f i  w2 h( f i  f xh  f y k1  0.5 f xx (h) 2  f xyhk1  0.5 f yy ( k1 ) 2 )
yi 1  yi  h( w1  w2 ) f i  h 2 ( w2f x  w2 f i f y )
y(x)をテイラー展開
yi 1  yi  hyi  0.5h2 y
df
dy
yi  f i
yi 
 fx  fy
 f x  f y fi
dx
dx
yi 1  yi  h( w1  w2 ) f i  h 2 ( w2f x  w2 f i f y )
yi 1  yi  hyi  0.5h2 y
df
dy
yi 
 fx  fy
 f x  f y fi
yi  f i
dx
dx
yi 1  yi  hf i  0.5h 2 ( f x  f y f i )
yi 1  yi  hf i  h 2 (0.5 f x  0.5 f y f i )
w1  w2  1
w2  0.5
w2   0.5
w1  w2  0.5
 1
 1
２次ルンゲ―・クッタ法
yi 1  yi  0.5k1  0.5k2
k1  f ( xi , yi )h
k2  f ( xi  h, yi  k1 )h
xi 1  xi  h
k1  f ( xi , yi )h
k2  f ( xi  h, yi  k1 )h
yi 1  yi  w1k1  w2 k2
w1  w2  0.5
 1
 1
４次ルンゲ―・クッタ法
yi 1  yi  (k1  2k2  2k3  k4 ) / 6
k1  f ( xi , yi )h
k2  f ( xi  0.5h, yi  0.5k1 )h
k3  f ( xi  0.5h, yi  0.5k2 )h
k4  f ( xi  h, yi  k3 )h
xi 1  xi  h
２階微分方程式
d2y
dy
 p ( x, y )  q ( x, y )  0
2
dx
dx
y( x0 )  y0
いま
dy
z
dx
y( x0 )  y0
（初期条件）
z0  y0
と置くと、２階微分方程式は１階微分方程式になる
dz
  p ( x, y ) z  q ( x, y )  g ( x, y , z )
dx
２階微分方程式
d2y
dy
 p ( x, y )  q ( x, y )  0
2
dx
dx
y( x0 )  y0
y( x0 )  y0
（初期条件）
連立１階微分方程式
dy
 f ( x, y , z )
dx
f ( x, y , z )  z
dz
 g ( x, y , z )
dx
x0 , y0 , z0
g ( x, y , z )   p ( x, y ) z  q ( x, y )
（初期条件）
連立１階微分方程式
dy
 f ( x, y , z )
dx
dz
 g ( x, y , z )
dx
x0 , y0 , z0
（初期条件）
オイラー法に数値解
yi 1  yi  hf ( xi , yi , zi )
zi 1  zi  hg ( xi , yi , zi )
xi 1  xi  h
x0 , y0 , z0
（初期条件）
２階微分方程式
d2y
dy

p
(
x
,
y
)
 q ( x, y )  0
2
dx
dx
y( x0 )  y0
y( x0 )  y0
（初期条
件）
４次ルンゲ･クッタ法の漸化式
yi 1  yi  (k1  2k2  2k3  k4 ) / 6
zi 1  zi  (l1  2l2  2l3  l4 ) / 6
k1  hzi
k  h( z  0.5l )
 2
i
1

k3  h( zi  0.5l2 )
k 4  h( zi  l3 )
x0 , y0 , z0  y0
l1  hg ( xi , yi , zi )
l  hg ( x  0.5h, y  0.5k , z  0.5l )
2
i
i
1
i
1

l3  hg ( xi  0.5h, yi  0.5k 2 , zi  0.5l2 )
l4  hg ( xi  h, yi  k3 , zi  l3 )
g ( xi , yi , zi )   p( xi , yi ) zi  q( xi , yi )
例：ベッセル微分方程式
d 2 y 1 dy

 y0
2
dx
x dx
初期条件：
y(0)  1
y(0)  0
h  0.1
y1  y0  (k1  2k2  2k3  k4 ) / 6
z1  z0  (l1  2l2  2l3  l4 ) / 6
k1  hzi
k  h( z  0.5l )
 2
i
1

k3  h( zi  0.5l2 )
k 4  h( zi  l3 )
x0 , y0 , z0  y0
y(0) / x  0.5
y (0.2) ?
y0  1
z0  y(0)  0
l1  hg ( xi , yi , zi )
l  hg ( x  0.5h, y  0.5k , z  0.5l )
2
i
i
1
i
1

l3  hg ( xi  0.5h, yi  0.5k 2 , zi  0.5l2 )
l4  hg ( xi  h, yi  k3 , zi  l3 )
g ( xi , yi , zi )   zi / xi  yi

Download Report