計算物理学:No.2

1 収束加速法： Aitken 法

1.1 Aitken 加速

与えられた数列に対し、 Aitken 加速された数列を返すサブルーチンを

subroutine Aitken(a,b,n)
c
implicit real*8(a-h,o-z)
dimension a(20),b(20)
.....
return
end

a: 配列。初めの数列を代入。 b: 配列。 Aitken 加速した結果を出力。 n: （最初の）数列の大きさ。

の形に作って見よ。このようにしておくと後で利用するのが（特に分割コンパイルを利用すると）便利になる。また、倍精度実数を使うこと。なお、関数副プログラム、サブルーチンの使い方、配列の定義と使い方については、 (Fortran,C,etc. の) 文法書を復習すること。

これは数値計算法概論で習った Richardson 加速のプログラムの応用であるといえる。若干注意しなければならないのは、最初の数列の大きさが n とすると、 Aitken 加速した後は n − 2 個の数列になる。

1.2 数値積分での台形公式

数値積分の台形公式で被積分関数が滑らかな関数の場合は、分割数 N に対し誤差が O(N^{−2}) となるので、分割数を N = 1, 2, 4, 8, 16, 32, 64, ⋯ と 2 の冪乗にとると 1 次収束となり、 Aitken 加速が使える。

区間の端で微係数が不連続な場合、上の誤差評価が成立しないが、端で被積分関数が constant+ x^{α}(0 < α < 1) と振舞う時、誤差は O(N^{−α−1}) で押えられる（解析的にこれを導いて見よ）。従ってこの場合も分割数を 2 の冪乗にとると Aitken 加速が使える。

1.3 交代級数

交代級数には Aitken 加速が有効である。交代級数として、次の例を上げる。

a_{n} = ∑_{k=1}^{n} \frac{(−1)^{k+1}}{2k − 1} ,  lim_{n→∞} a_{n} = arctan 1 = \frac{π}{4}  (1)

a_{n} = ∑_{k=1}^{n} \frac{(−1)^{k−1}}{k} ,  lim_{n→∞} a_{n} = log 2  (2)

はそれぞれ

a_{n} = a + (−1)^{n}(\frac{b_{1}}{n} + ⋯)  (3)

と振舞う。これらについて、 (多段の)Aitken 加速してみよ。

1.4 混合型

1 次収束と、 O(n^{−α}) の両者が混合したケース

a_{n} = a + λ^{n}(c_{1}n^{−α_{1}} + c_{2}n^{−α_{2}} + ⋯)  (4)

(ここで |λ| ≤ 1, λ ≠ 1, α_{1} < α_{2} < ⋯ である。 |λ| = 1 の場合はさらに 0 < α_{1} < α_{2} < ⋯ とする) の数列に対して No.1 の (17) 式のように Aitken 加速すると、

a_{n}^{(1)} ≈ a + \frac{c_{1}α_{1}}{(1 − λ)^{2}} λ^{n+2}n^{−α_{1}−2}  (5)

となる [1](ただし、 α_{2} > α_{1} + 2 を前提としている)。混合型を Aitken 加速した結果は再び混合型になっているので、さらに Aitken 加速することが可能である。交代数列は、 λ = −1 の場合に当たる。式 (5) の導出は面倒であるので付録で紹介する。

2 レポート

2.1 レポート課題

1. 以下の級数

a_{n} = ∑_{k=1}^{n}(−x)^{k−1},  lim_{n→∞} a_{n} = 1/(1 + x),  (|x| < 1)  (6)

を Aitken 加速してみよ（1 回で十分のはずである）。これは 1 次収束である。 No.1 の 2.2.2 の (16) 式を使い、収束性を調べてみよ。

2. 交代級数 (今回のプリントの (1),(2) その他岩波数学公式集 II p。 37 など参照) の内 1 つ以上の Aitken 加速をして見よ。 3. 数値積分の台形公式で区間の端で被積分関数が constant + x^{α}(0 < α < 1) と振舞う時、

誤差は O(N^{−α−1}) で押えられる（解析的にこれを導いて見よ）。従ってこの場合も分割数を N = 1, 2, 4, 8, 16, ⋯ と 2 の冪乗にとると Aitken 加速が使える。以前扱った関数の台形公式による積分で、

4 ∫_{0}^{1} \sqrt{1 − x^{2}}dx = π  (7)

の収束性を調べ、 Aitken 加速してみよ。余裕があれば

∫_{0}^{1} （1 − x^{2}）^{α}dx  (8)

と一般化したものについても調べよ。

2.2 期限、提出先など

• 提出期限: 5 月 12 日

• 提出先:野村 (理学部 2 号館 6 階 2641 号室)

• e-mail:knomura@stat.phys.kyushu-u.ac.jp

e-mail でも受け付ける。提出者一覧をホームページ上に載せる予定であるので、各自提出後確認のこと。

何の数列に対して収束加速法を行なったのか明確に記し、プログラムと実行結果、収束性の考察、その他のコメントを付けて計算物理学レポート No.1 として提出せよ。講義名、名前、学籍番号、提出日時を表紙に明記すること。コピーレポートはオリジナルも含め認めない。

解析的な考察ができる場合は、その導出も書くこと。表を入れて分かり易すく工夫したものは評価する。参考文献などあれば明記すること。また、コピーレポートは認めないが、互いの議論によりレポートの質が向上した場合は評価するので、議論した部分・参考とした意見と相手の氏名を明記して出すこと (ホームページの\underline{レポート作法}参照)。 Fortran 以外の言語 (C,C++,Pascal) を使った場合は、記しておくこと。なお、日本語 L^{A}TEX によるレポートは歓迎する [2](コマンドは platex )。

これまで、数値計算法概論を取っていなくてこの講義を受講した人は、そのことをレポートに明記すること。また、その場合は数値計算法概論を取っている人から協力を得ることがあるので、協力者の名前も明記すること。

3 一般的注意

3.1 メールについて

電子メールによるレポート提出法は、例えば report2b と言うファイル (プログラム+解析結果+説明、日本語を含んでいても良い) を作成した後、

% nkf -j report2b > report2c
% mail knomura@stat.phys.kyushu-u.ac.jp < report2c

とすれば良い (nkf は日本語コード変換のため)。

さらに便利なメールを扱うツールとしては、 mule の中から立ち上げる mew を勧める ( M-x mew で起動。九大出身の人が作った世界的フリーソフト)。

A 付録

まず、

a ≈ a_{n}^{(1)} = a_{n} − \frac{(a_{n+1} − a_{n})^{2}}{a_{n+2} − 2a_{n+1} + a_{n}} = a_{n+2} − \frac{(a_{n+2} − a_{n+1})^{2}}{a_{n+2} − 2a_{n+1} + a_{n}}  (9)

は行列式を使った表記で

a_{n}^{(1)} = a + \frac{| \begin{matrix}a_{n+2} − a & a_{n+1} − a \\ a_{n+1} − a & a_{n} − a \end{matrix} |}{a_{n+2} − 2a_{n+1} + a_{n}}  (10)

と書き直せることに注意する。混合型の (4) を次のように書き直す。

a_{n} = a + λ^{n}f(n − 1),  f(n − 1) = c_{1}n^{−α_{1}}  (11)

これを上記の式に入れて計算すると、

a_{n}^{(1)} = a + λ^{n+2} \frac{| \begin{matrix}f(n + 1) & f(n) \\ f(n) & f(n − 1) \end{matrix} |}{λ^{2}f(n + 1) − 2λf(n) + f(n − 1)}  (12)

である (λ^{n+2} を括り出すには、行列式の多重線形性を用いた)。ここで、分母は、 c_{1}(1 − λ)^{2}n^{−α_{1}}, 分子の行列式は

| \begin{matrix}f(n) + f′(n) + \frac{1}{2}f′′(n) & f(n) \\ f(n) & f(n) − f′(n) + \frac{1}{2}f′′(n) \end{matrix} |
≈ f(n)f′′(n) − (f′(n))^{2} ≈ c_{1}^{2}α_{1}n^{−2α_{1}−2}  (13)

と近似されるので、

a_{n}^{(1)} ≈ a + \frac{c_{1}α_{1}}{(1 − λ)^{2}} λ^{n+2}n^{−α_{1}−2}  (14)

が導かれた。

参考文献

[1] 数値計算法の数理、杉山正顕、室田一雄共著、 (岩波書店) 第 12 章 4 節

[2] 文書処理システム “L^{A}TEX”, Leslie Lamport 著、倉沢良一監訳、アスキー出版社