11 数値積分

数値積分について練習してみましょう。 数値積分はOctaveのもっとも得意とするところです。

11.1 数値積分の例1

まずは、みなさんが良く分かっている例から計算してみましょう。

$$\int_0^2 (4 x^3 - 3 x^2 + 5) dx$$

手計算で算出してみてください。 では、Octaveで計算してみましょう。

############## Program ###############
function y = f(x)
  y = 4 * x^3 - 3 * x^2 + 5
endfunction

quad("f", 0, 2)

######################################

このプログラム中の ``quad'' と言うのが数値積分です。 内容は関数fに対して0から2の範囲で積分せよ、よいうことです。 続けて、実行例を見てみましょう。

############## 実行例 ################
[nisimiya]$ octave int_1.m
GNU Octave, version 2.0.16 (i386-vine-linux-gnu).
Copyright (C) 1996, 1997, 1998, 1999, 2000 John W. Eaton.
This is free software with ABSOLUTELY NO WARRANTY.
For details, type `warranty'.

y = 6
y = 4.9980
y = 24.075
y = 4.9552
y = 20.515
y = 4.8235
y = 15.485
y = 4.7645
y = 10.616
y = 5.2930
y = 7.1059
y = 4.9999
y = 24.844
y = 4.9867
y = 22.587
y = 4.8980
y = 18.076
y = 4.7608
y = 12.940
y = 4.9116
y = 8.6484
ans = 18

######################################

なにやら、数値がたくさん表示されました。一体これらは何を表しているので しょうか？ 最終行が答えになっているのは分かりますね。 他の行は、数値計算の反復値を示しております。

リストを以下のように書き換えてみましょう。

############## Program ###############
function y = f(x)
  y = 4 * x^3 - 3 * x^2 + 5;
endfunction

[a,b,c,d]=quad("f", 0, 2)

######################################

実行例を見てみましょう。

############## 実行例 ################
[nisimiya@localhost Renshu]$ octave int_1.m
GNU Octave, version 2.0.16 (i386-vine-linux-gnu).
Copyright (C) 1996, 1997, 1998, 1999, 2000 John W. Eaton.
This is free software with ABSOLUTELY NO WARRANTY.
For details, type `warranty'.

a = 18

b = 0

c = 21

d =  1.9984e-13
######################################

これらの数字は何を意味しているのでしょうか？ これは順番に、

[答え エラーの有無 反復回数 計算誤差]

を意味しています。 もちろん、記号a, b, c, dというのは、便宜上つけただけで 何でもかまいません。答えがベクトル形式で返ってくるという事です。 ここで重要なのは、１番目の答えと２番目のエラーの有無です。 ``0''が出力されれば、正常に終了したことになります。

11.2 数値積分の例2

もう少し具体例をやってみましょう。 こんな関数はどうでしょうか？

$$\int_0^{2 \pi} \frac{x}{1+\sin^2(x)} dx$$

############## Program ###############
function y = f(x)
  y = x / (1 + sin(x)^2);
endfunction

[a,b,c,d]=quad("f", 0, 2*pi)

######################################

以下は実行結果です。

[nisimiya@localhost Renshu]$ octave int_2.m
GNU Octave, version 2.0.16 (i386-vine-linux-gnu).
Copyright (C) 1996, 1997, 1998, 1999, 2000 John W. Eaton.
This is free software with ABSOLUTELY NO WARRANTY.
For details, type `warranty'.

a = 13.958

b = 0

c = 147

d =  1.7878e-10

11.3 数値積分の例3

以下の積分を求めてみよう。

$$\int_{-2}^2 e^{-x} \sin(x) dx$$

############## Program ###############
function y = f(x)
  y = exp(x).* sin(x);
endfunction

[a,b,c,d]=quad("f", -2, 2)

######################################

実行例を見てみましょう。乗算の記号に注意して下さい。 ``.*'' 各関数が行列である場合を仮定してこのような記号になってます。

############## 実行例 ################
[nisimiya@localhost Renshu]$ octave int_3.m
GNU Octave, version 2.0.16 (i386-vine-linux-gnu).
Copyright (C) 1996, 1997, 1998, 1999, 2000 John W. Eaton.
This is free software with ABSOLUTELY NO WARRANTY.
For details, type `warranty'.

a = 4.9303

b = 0

c = 21

d =  6.4933e-14

11.3.1 演習

1. $$\int_0^{2 \pi}\cos^2(x) dx$$を求めよ。
2. $\int_0^3 e^{-x} * \sqrt(1/x) dx$ を求めてみよう。