複素周回積分による数値積分の理論誤差解析および応用

複素周回積分による数値積分の理論誤差解析および応用
緒方秀教 1 , 平山弘 2
1
電気通信大学大学院情報理工学研究科情報・通信工学専攻
神奈川工科大学創造工学部自動車システム開発工学科
e-mail : [email protected]
2
1
はじめに
台形公式は周期解析関数の 1 周期区間積分に
対して高精度な積分値を与える．[1] ではこれに
着目し，任意の有界閉区間上の実積分に対し，
これを複素周回積分に変換して台形公式を適用
した高精度な数値積分公式を提案している．本
研究ではこの数値積分公式に対し理論誤差評価
を与え，また，その応用，佐藤超関数論との関
係についても触れる．
[1] で示された数値積分公式の概要を，以下
に記す．区間 [a, b] ( −∞ < a < b < +∞ ) 上の
積分
∫
b
f (x)w(x)dx
I=
(1)
a
を考える．ここで，f (x) は与えられた関数 f (x)，
w(x) は重み関数である．C を複素平面内の閉曲
線で区間 [a, b] をその内部に含むものとし，f (x)
を複素関数に拡張した f (z) は C およびその内
部を含む領域で正則であるとする．このとき，
Cauchy の積分公式を用いると，積分 I は次の
ように書き換えられる．
{
}
∫ b
I
1
f (z)
I=
w(x)
dz dx
2πi C z − x
a
{∫ b
}
I
1
w(x)
=
dx dz.
f (z)
2π C
a z−x
ここで，
∫
Ψ(z) =
a
b
w(x)
dx
z−x
とおくと，積分 I は次のように表される．
I
1
f (z)Ψ(z)dz
I=
2πi C
∫ α
1
=
f (ϕ(u))Ψ(ϕ(u))ϕ0 (u)du, (2)
2πi 0
ここで，C のパラメータ表示を z = ϕ(u), 0 5
u 5 α （ϕ(u) は周期 α(> 0) の周期関数）と
した．関数 Ψ(z) の具体形は，例えば，(
(a, b) =)
z+1
(−1, 1)，w(x) = 1 の場合は Ψ(z) = log
z−1
となり，(a, b) = (0, 1)，w(x) = xα−1 (1−x)β−1
( α, β > 0 ) の場合は，
(
)
1
1
Ψ(z) = B(α, β) F α, 1; α + β,
z
z
となる．(2) 最右辺に台形公式を適用すること
により，次の数値積分公式を得る．
I ' IN
≡
N −1
h ∑
f (ϕ(kh))Ψ(ϕ(kh))ϕ0 (kh). (3)
2πi
k=0
後で述べるように，この数値積分法は佐藤超
関数の理論と密接な関係がある．そこで，(3)
で与えられる数値積分法を，ここでは「超関数
法」と呼ぶことにする．
2
理論誤差解析
[2]4.6.5 節の定理を用いて，数値積分公式 (3)
の誤差評価に対する下記の定理を得る．この定
理から，数値積分公式 (3) の誤差は標本点数 N
を大きくしていくと指数関数的に減衰すること
がわかる．
定理 1 積分路 C のパラメータ表示関数 ϕ(w)
は複素平面内の帯状領域 | Im w| < d0 ( d0 > 0 )
で正則，| Im w| 5 d0 で連続，そして，周期 α
の周期関数であるとする．−d0 5 δ 5 d0 なる
定数 δ に対し線分 0 5 Re w < α, Im w = δ の
関数 z = ϕ(w) による像である閉曲線を Cδ と
し，Cδ の内部である領域を Dδ とする．ただ
し，−d0 5 δ1 < δ2 5 d0 のとき Dδ2 $ Dδ1 で
あるとする．
関数 f (z) は D−d0 で正則，関数 Ψ(z) は D−d0 \
[a, b] で正則であるとする．このとき，0 < d <
d0 なる任意の d に対し，数値積分公式 (3) の誤
差に関して次の不等式が成立する．
|I − IN | 5 2α max |f (ϕ(w))Ψ(ϕ(w))ϕ0 (w)|
Im w=±d
×
exp(−(2πd/α)N )
.
1 − exp(−(2πd/α)N )
3
数値計算例
のことについて言及する．佐藤超関数とは，大
雑把に言えば，正則関数の境界値の差で表され
る関数である．すなわち，実軸上の開区間 (a, b)
を複素領域 D が含むとして，D \ [a, b] で正則
な関数 F (z) を用いて
積分
∫
1
I=
xα−1 (1 − x)β−1 ex dx
0
= B(α, β)F (α; α + β; 1)
= 3.71819 70362 84701 . . . × 104
( α = β = 10−4 )
に対し，超関数法と DE 公式で積分の近似値を
計算し，その誤差の標本点数 N に対する変化
の様子を図 1 に示した．ただし，本方法の方法
では，複素積分路 C として楕円
z = ϕ(u)
=
1 1
+
2 4
(
)
(
)
1
i
1
ρ+
cos u +
ρ−
sin u,
ρ
4
ρ
ρ = 10
をとっている．数値計算は C++でプログラミ
ングし，倍精度計算で行った．
図から，超関数法による数値積分の誤差は標
本点数 N に対し指数関数的に減衰しているが，
DE 公式ではまったく近似値が計算できていな
いことがわかる．この理由は次のように考えら
れる．DE 公式では，標本点を積分区間端点に
集積するようにとっているが，この数値例の場
合あまりにも特異性が強いので，十分な精度を
得ることができない．一方，超関数法では特異
点である積分区間端点を避け，複素平面上の被
積分関数が穏やかに変化する積分路上に標本点
を取っているので，被積分関数の特異性による
精度の低下を免れるのである．
2
log10(relative error)
0
a<x<b
と表される f (x) を（佐藤）超関数と呼ぶ．そ
して，正則関数 F (z) を超関数 f (x) の定義関数
と呼ぶ．超関数 f (x) の区間 (a, b) 上の積分を次
で定義する．
∫
a
0 5 u < 2π,
b
I
f (x)dx ≡ −
F (z)dz,
C
ここで，C は区間 (a, b) をその内部に含み，領
域 D に含まれる正の向きの閉積分路である．
本研究で扱っている数値積分公式に話を戻
すと，(1) の被積分関数 f (x)w(x) は F (z) =
−1
f (z)Ψ(z) を定義関数とする超関数と見な
2πi
すことができる．すると，元の積分 I の周回積
分表示 (2) は，超関数の意味での f (x)w(x) の
積分にほかならない．したがって，超関数法は
は考えている実積分 I を与える超関数積分を台
形公式で近似計算していることに相当する．
数値実験例で示したように，超関数法は端点
に極めて強い特異性のある積分に対して有効
であるので，境界要素法に現れるような特異
性をもつ積分に対して有用であると考えられ
る．さらに，超関数の理論では Cauchy の主値，
Hadamard の有限部分といった特異積分も通常
の積分と統一的に扱うことができるので，超関
数法はこうした特異積分の計算にも応用できる
と期待できる．
謝辞本研究は JSPS 科研費 25400196 の助成
を受けている．
-2
-4
-6
hyperfunction rule
DE rule
-8
参考文献
-10
-12
-14
-16
0
5
10
15
20
N
25
30
35
40
図 1. 超関数法と DE 公式の誤差
4
f (x) = F (x + i0) − F (x − i0),
佐藤超関数論との関係，応用
本研究で超関数法と呼んだ数値積分公式は，
その名の通り佐藤超関数論 [3] と関連がある．そ
[1] 平山弘：周回積分変換法による数値積分
法，第 44 回数値解析シンポジウム講演
予稿集，2015 年，21–24．
[2] P. J. Davis and P. Rabinowitz: Methods of Numerical Integration, 2nd. ed.,
Academic Press, San Diego, U.S.A.,
1984.
[3] 今井功：応用超関数論 I，II，サイエン
ス社，1981 年．

Download Report