EISPACK についての原作者講義

J. H. Wilkinson 理論、C. Moler FORTRAN 実装

Nick Higham のブログで、EISPACK についての原作者等による講義のビデオが紹介されていました。発掘の経緯についても書かれています。 The Argonne Tapesnickhigham.wordpress.com

３月に公開されて、現在まで 100 view 程度で草。

しかし、こういうリソースがあるのは白人様の文明の底力で恐れ入ります。ただで視聴乞食させていただきお靴を舐めたい所存であります。

あいまいみー４巻

短いし、FORTRAN の話なので、とりあえず Moler の分だけ見てみました。

Moler の講義では、実対称行列を Hausholder 法で三重対角化するサブルーチン tred2 と、それをさらに QL 法で完全に対角化する tql2 を基本例としてソースコードを OHP で写しながら解説しています。 tql の途中で終わっているのが残念です。tred + tql は、かつて一世を風靡した numerical recipes にも、後継改良 77 版がそのまま収録されています。

tred で扱うのは実対称行列なので下三角に変換行列を貯めてゆくことでメモリーを節約していますが、Moler も巧妙だと称賛しています。元々の ALGOL 版からそうなっていたようです。

www.youtube.com Youtube 自動字幕では ALGOL が alcohol になっていて草。

今も開発の続く LAPACK は、線形方程式ライブラリの LINPACK と固有値問題ライブラリの EISPACK を合体させ FORTRAN77 に書き換えたものです。EISPACK は元々はドイツ人（C. Reinsch ?）が ALGOL でコツコツ集めていたのを、イギリス人の Wilkinson が加わって集大成し、それをさらにアメリカ人が FORTRAN66 に移植して成ったものです。

ALGOL はコンピュータ言語としてのほか、アルゴリズムを共有するための表記法としての目的を持っていました。

Wilkinson は誤差解析で有名ですが、固有値問題の本も書いています。その Wilkinson による EISPACK ライブラリの理論的背景講義の後、FORTRAN 移植を行った Moler の FORTRAN 実装に関する講義がなされています。Moler は MATLAB の開発者として有名ですが、その他にも G. E. Forsythe, M. A. Malcolm と数値計算レシピ本の嚆矢のような本も出しています。

Rounding Errors in Algebraic Processes

基本的演算における丸め誤差解析 (1974年)

The Algebraic Eigenvalue Problem (MONOGRAPHS ON NUMERICAL ANALYSIS)

Computer Solution of Linear Algebraic Systems (Automatic Computation)

線形計算の基礎―計算機のための 連立1次方程式のプログラミング (1969年)

Computer Methods for Mathematical Computations

Numerical Recipes in Fortran 77: The Art of Scientific Computing

講義では Machine Epsilon や二次方程式の桁落ち予防変形についても触れています。 また Moler は、元の ALGOL プログラムで、同一名の変数をプログラム中で全く異なる目的で何度も使っていることにちょっとだけ苦言を呈しています。

講義冒頭で触れられている J. H. Wilkinson and C. Reinsch, Handbook for Automatic Computation の Moler 本人による紹介記事が MATLAB のサイトで読めます。 blogs.mathworks.com

Handbook for Automatic Computation: Volume II: Linear Algebra (Grundlehren der mathematischen Wissenschaften)

tred & tql

昔、これらのルーチンにはお世話になりました。あの頃こういうリソースがあれば良かったろうにと思います。

Hausholder 法での三重対角化 tred のアルゴリズムは有限回のステップで終わりますが、tql の三重対角行列の対角化は反復法なので適当な閾値になるま繰り返します。５次方程式以上の解の公式がない以上、５次以上の行列の固有値問題にも有限ステップの解法が無いであろうからもっともな気がします。

問題によっては、最初の三重対角化に Lanczos 法を使って、非対角項が十分小さくなった適当なところで打ち切ると、小さい行列の対角化になって（tql は O(N³) のアルゴリズムなので）計算量を大幅に減らせます。　

あいまいみー コミック 1-8巻 セット

あいまいみー　9 (バンブーコミックス WINセレクション)

あいまいみー [Blu-ray]