このピュータ基礎演習 TUAT_E2015

Excel による数値計算 (1) : 数値微分と数値積分
「コンピュータ基礎演習」資料
E1:桝田, E2:竹村
1
目的
数値計算法とは，解析的な解法を用いるのではなく，数値データを処理して近似的な計算を行
う手法である。今回は，数値微分と数値積分の方法について学び，Excel を使って実際にそれらの
計算を行う。
2
数値微分
(微分) f (x) = lim
Δx→0
f (x + Δx) − f (x)
Δx
(1)
数値微分は差分で代用する。
(前方差分)
f (x) ∼
=
f (x + Δx) − f (x)
Δx
(2)
(後方差分)
f (x) ∼
=
f (x) − f (x − Δx)
Δx
(3)
(中心差分)
f (x + Δx) − f (x − Δx)
(前方差分) + (後方差分)
=
f (x) ∼
=
2
2Δx
(4)
式 (2) の 1 階前方差分を，もう一度前方差分すると 2 階微分 f (x) を求めることができる。
f (x + 2Δx) − 2f (x + Δx) + f (x)
f (x + Δx) − f (x)
=
f (x) ∼
=
Δx
(Δx)2
(5)
一般には，次の式がよく用いられる。
f (x + Δx) − 2f (x) + f (x − Δx)
f (x) ∼
=
(Δx)2
3
数値積分
(定積分) I =
b
a
(6)
f (x)dx
(7)
f (x) が連続であるとき I は，積分の範囲 [a, b] で x 軸と f (x) で挟まれた部分の面積となる。
1
f (x)
f (x )
x
x i−1 x i x i +1
0
0
x
(a) 台形法
x
x i−1 x i x i +1
x
(b) シンプソン法
図 1: 数値積分の考え方
ここで，積分範囲 [a, b] を n 等分して
x1 = a,
xn = b,
Δx =
b−a
n−1
(8)
とおく。これにより，積分範囲 [a, b] における x を Δx の幅で刻んで，サンプル点 x0 , x1 , · · · , xn
を得ることができ，隣り合うサンプル点間には次式の関係がある。
xi+1 = xi + Δx,
(i = 1, 2, · · · , n − 1)
(9)
式 (7) の近似値を，サンプル区間ごとの長方形の面積の和 S長方形 とすると，次式が成り立つ。
I∼
= S長方形 =
n−1
f (xi )Δx
(10)
i=1
3.1
台形法
図 1 (a) のように，隣り合う 2 点間の f (x) を直線で近似すると，台形の面積の和を求めること
によって積分値 I の近似計算ができる。これを式で表すと，式 (11) のようになる。
S台形 =
n−1
i=1
=
Δx f (xi ) + f (xi+1 )
2
Δx f (x1 ) + 2f (x2 ) + 2f (x3 ) + · · · · · · + 2f (xn−1 ) + f (xn )
2
一般に，この方法 (台形法) を用いると，式 (10) よりも精度の良い積分値が得られる。
2
(11)
3.2
シンプソン (Simpson) 法
台形法は 隣接する 2 点間を直線で近似したが，次に説明するシンプソン法では，図 1 (b) のよ
うに隣接する 3 点間を 2 次関数 (放物線) で近似するため，台形法よりも精度を良くすることがで
きる。
図 1 (b) のように区間 [xi−1 , xi+1 ] の 3 つのサンプル点
( xi−1 , f (xi−1 ) ),
( xi , f (xi ) ),
( xi+1 , f (xi+1 ) )
を通る放物線と x 軸とで挟まれた部分の面積 Si は，
Si =
Δx f (xi−1 ) + 4f (xi ) + f (xi+1 )
3
(12)
となる。なお，この式の導出は付録を参照のこと。これより，求める積分値 I の近似値 SSimpson
は，次式により得ることができる。
SSimpson =
Δx f (x1 ) + 4f (x2 ) + 2f (x3 ) + 4f (x4 ) + 2f (x5 )
3
+ · · · · · · + 4f (xn−3 ) + 2f (xn−2 ) + 4f (xn−1 ) + f (xn )
(13)
(注) 式 (13) を見ればわかるとおり，シンプソン法では，サンプル点数 n は奇数でなければなら
ない。
4
課題
(4.1) 式 (14) の定積分を，台形法及びシンプソン法により数値積分する Excel テンプレートを
作成しなさい。ただし，x のサンプル間隔 (刻み幅) Δx の設定は，授業中の指示に従うこと。
1
0
dx
+4
(14)
x2
(4.2) 課題 (4.1) で求めた各数値積分値と真値（手計算によって解析的に求めた解）より，式 (15)
の相対誤差 [%] を，Excel ワークシート上で計算しなさい。
真値 − 数値積分値 × 100 %
相対誤差 = (15)
真値
(4.3) Excel ワークシート上にテキストボックスを配置し，その中に各数値計算法の精度につい
ての考察を述べよ。
(4.4) 式 (14) とは異なる定積分 (授業中で指示する) についても，(4.1) と同様に台形法及びシ
ンプソン法により数値計算しなさい。また，この場合においても (4.2) 及び (4.3) と同様に計算
精度の評価を行いなさい。
課題の提出
作成した Excel ファイルを提出。提出方法・期限は，授業中の指示に従うこと。
3
5
付録 : 式 (12) の導出
閉区間 [xi−1 , xi+1 ] にある３つの点を，(xi−1 , yi−1 ), (xi−1 , yi ), (xi−1 , yi+1 ) と表す。これら 3 点
を通る放物線は，係数を a 及び b とおくと
y = a (x − xi )2 + b (x − xi ) + yi
(16)
と表すことができる。また，x 軸上の刻み幅 Δx は，
Δx = xi+1 − xi = xi − xi−1 =
xi+1 − xi−1
2
となる。式 (16) が (xi−1 , yi−1 ) 及び (xi+1 , yi+1 ) を通ることから，
yi−1 = a Δx2 − b Δx + yi
yi+1 = a Δx2 + b Δx + yi
が成り立つ。従って，係数 a, b は
a = (yi−1 − 2yi + yi+1 ) / (2 Δx2 )
(17)
b = (yi+1 − yi−1 ) / (2 Δx)
(18)
となる。
放物線と x 軸上の区間 [xi−1 , xi+1 ] とで囲まれた部分の面積 Si は，式 (16) を積分することに
よって，次のようになる。
Si =
=
xi+1
xi−1
y dx
2
a Δx3 + 2 Δx yi
3
(19)
式 (19) の a に式 (17) を代入すると，
Si =
=
Δx
(yi−1 + 4yi + yi+1 )
3
Δx
(f (xi−1 ) + 4f (xi ) + f (xi+1 ))
3
(20)
となり，式 (12) が導出された。
4