定積分の誤差特性関数の評価
2009SE213 野間ちほ
指導教員：杉浦洋

はじめに

1

積分は解析学における最も基本的な演算である．しかし，積分値が解析的に与えられていることは稀で，その値を与えるために数値計算に頼らざるを得ないことが多い．したがって，数値積分による値の精度保証は極めて重要な課題である．ここで，精度保証とは近似積分値の絶対誤差の上界を計算することである．
今回は正則関数の定積分の精度保証付積分について研究する．数値の絶対誤差の上界を与える方法としては，被積分関数の高階導関数を用いるもの [2] と被積分関数の複素周回積分を用いるもの [1] がある．この論文では，後者の方法について研究を行った．

2複素区間計算2.1

実区間関数[定理 3.1] f (x) が [−1, 1] で正則のとき
I1
En f = Qf − Qn f =Φn (z)f (z)dz
2πi C∫ 1(4)

Φn (z) を数値積分公式の誤差の特性関数という．これは積分公式 Qn だけから定まり，被積分関数 f (x) にまったく依存しない．
補間型積分則については次の定理が成り立つ．
[定理 3.2] すべての標本点を零点とする n 次多項式をpn (x) = c(x − x1 )(x − x2 ) . . . (x − xn )，c ̸= 0(5)とすると，

f ([x]) = {f (x)|x ∈ [x]}が成立する．
複素区間解析

複素平面上の閉長方形領域 [z] = [x] + i[y] を複素区間という．実区間演算は，複素関数を実部と虚部に分けて計算すれば自然に複素区間演算に拡張される．

数値積分公式と誤差の特性関数
実軸上の有限区間 (-1, 1) 上の積分

f (x)dx1
pn (z)(1)を考える．この積分を数値的に近似計算するための公式をn∑

wi f (xi )pn (x)
dx−1 z − xn−1
∑Qn f =, z ̸= [−1, 1].//(6)

ωl f (ξl )(7)l=0を用いる．標本点は n 次 Chebyshev 多項式 Tn (x) の零点()π1ξl = cosl+(0 ≤ l < n)(8)n2である．Polya 則は補間型積分則で，その重み ωl は()⌊(n−1)/2⌋∑ cos 2πk l + 12
n2 
ωl =1−2(0 ≤ l < n)π4k 2 − 1
k=1(9)である．ここで，⌊x⌋ は実数 x を超えない最大正数である．Polya 則の誤差特性関数の絶対値の上界は以下の定理で与えられる [2]．

|Φn (z)| ≤ Fn (z),−1Qn f =∫ 1[定理 3.3] z ∈ C − [−1, 1] について，

∫ 1
Qf =

Φn (z) =具体的な積分則として n 点 Polya 則:

を f の区間関数という．Mathematica には四則演算と標準的な実数関数の区間関数が実装されている．
標準関数と四則演算による合成関数を f とする．f を構成する標準関数と四則演算をすべてその区間関数で置き換えたものを f の区間拡張といい [f ]([x]) で表す．このとき
f ([x]) ⊆ [f ]([x])

3n
∑
1wj
dx −．
//
z − xj−1 z − xj=1

Φn (z) =閉区間 [x, x] (x < x) を [x] = [x, x] と書く．実連続関数 f (x) に対し，2.2(3)

Fn (z) =(2)2n−n|w| − |w|√w = z + (z + 1)

j=1とする．xj ∈ [−1, 1] を積分公式 (2) の標本点，wj を重みという．この公式を n 点公式と呼ぶ．積分誤差は次のように複素周回積分で表せる [1]．

4log(|a| + 1) +(|a| − 1) +√√(|a| + 1)2 + b2(|a| − 1)2 + b2

z−1
, z = a + ib.//z+1

積分誤差の包囲定理 3.1，3.3 より，|En f | ≤12πI
Fn (z)|f (z)||dz|C(10),この式の最右辺がなるべく小さくなるような，積分路 C 5 数値実験を考える．G(z) = Fn (z)|f (x)| の鞍点にグラディエント
関数 f (x) = cos x について定積分方向をたどって登り，通過後は逆グラディエント方向を∫ 1
たどって降りれば式 (10) の最右辺の積分が小さくできる．
If =f (x)dx = 1.682941969615793 · · · ,これを鞍点法と呼ぶ．
−1
鞍点法を実現する積分路を J 角形 C で近似する．通過
f (x) = cos x
点を α0 , α1 , · · · , αJ = α0 とし，の精度保証付計算を行う．
C = C1 + C2 + · · · + CJ ,標本点数を n = 10 として Cn f を計算し，Cj : αj−1 →αj (1 ≤ j ≤ J)
(11)
C10 f = 1.682941969605210 · · ·とすると，

IG(z)|dz| =CJ I∑j=1Cj線分 Cj の m 等分点を αj−1 から近い順に

α0 = 10.61, α1 = 10.61 + 6.00i,
α2 = 8.77 + 8.22i, α3 = 5.63 + 10.20i,
α4 = 2.31 + 10.42i, α5 = 10.7978i(12)と定めた．式 (17) の Mi はとする．点 α, β を対角線上の頂点とする複素区間をR[α.β] = [min{ℜα, ℜβ}, max{ℜα, ℜβ}]+i [min{ℑα, ℑβ}, max{ℑα, ℑβ}](19)得る．被積分関数 cos x は実偶関数ゆえ，対称性を考慮し，第 1 象限のみを計算する．J 角形 C の第 1 象限部分の頂点を，G(z)|dz|．

αj−1 = β0 , β1 , · · · , βm = αj ,kβk = αj−1 + (αj − αj−1 ) (0 ≤ k ≤ m)
m(18)

(13)

M1 = 2π(1.51 · · · × 10−12 ), M2= 2π(1.18 · · · × 10−11 ),M3 = 2π(1.65 · · · × 10−10 ), M4= 2π(6.89 · · · × 10−10 ),M5 = 2π(8.97 · · · × 10−10 )と置くと，m∪Cj ⊂R[βk−1 , βk ](14)

k=1である．これにより，
Im I
11 ∑ βkG(z)|dz|
G(z)|dz| ≤2πi Cj2πβk−1≤12πk=1m∑|En f | ≤ M = 1.12 · · · × 10−9[G](R[βk−1 , βk ])|βk−1 − βk |k=1=: Mjとなった．2，3，4 象限は以上の繰り返しである．以上を
式 (18) に代入して(15)(20)が得られた．式 (18),(19),(20) より，|En f | = 1.05 · · · ×
10−11 ＜ M であり，精度保証は成功した．同じ問題に対する安田 [2] の結果 M = 6.05 · · · × 10−9 の 1/5 以下の上界が得られた．鞍点法の効果が実証された．

6

まとめ

今回は，複素周回積分による積分則の誤差評価式を使って，精度保証付きの数値積分を行った．安田 [2] が提案した Polya 則の誤差特性関数の絶対値の上界関数 Fn (z) を|En f | ≤ M = M1 + M2 + · · · + MJ(16) 改良した．F (z)|f (z)| の複素平面における等高線を用い，n鞍点法により積分路の決定を行った．結果，高精度の精を得る．
以上より，積分誤差の上界 M は次の手順で計算される． 度保証を得ることができた．Polya 則の誤差特性関数の絶対値の上界の改良と，周回積分の積分路を鞍点法によStep0: Polya 則による近似積分値 Cn f を計算．
Step1: G(z) = Fn (z)|f (z)| の等高線図を描き，鞍点をみ る改良が非常に有効であることがわかった．今後の課題
は，Fn (z) をさらに改良し，精度を向上させることと，人つける．
Step2: 鞍点を頂点に含む J 角形積分路 C を定め，その 間による等高線の観察によらず積分路を自動的に決定することである．
頂点を C にそって，α , α , · · · , α = α とする．
ここで，[G] は G の区間拡張である．
以上より，不等式

01

J0

Step3: 各辺 Cj : αj−1 →αj 上の積分の上界
I1G(z)|dz|
Mj ≥2π Cj7(17)

参考文献[1] 森正武:数値解析と複素関数論, 薩摩書房, 東京,1975.

[2] 安田理佳子：制度保証付き積分則, 南山大学数理情報学部情報システム数理学科卒業論文 2011 年度卒業論を式 (16) で計算する．辺の等分数 m は適切に定める (m文，2012.
が大きいほど Mj は小さくなるが，計算コストがかかる)．
Step4: M = M1 + M1 + · · · + MJ ．//