数値計算ツール

大石 進一 Shin’ichi Oishi

2000 年 8 月 30 日

1 数値計算ツール

精度保証付き数値計算はやはり数値計算である。したがって、まず、数値計算のツールとしてどのようなものがあるかを知ることが重要となる。まず、数値計算に関連するサイトを紹介しよう。よくある質問とその解答のサイトとして FAQ Numerical Analysis and Associated Fields Resource Guide がある。この中には数値計算パッケージや数値計算ソフトウエアへのリンクもあり、基本的な数値計算ツールを揃えるのに便利である。また、最適化に関連する数値計算のサイトとしては Bibliography for Optimization with Sensitivity Analysis がある。

1.1 数値線形代数の直接解法ライブラリ

まず、線形計算のためのツールとして、メモリ階層(キャッシュヒット率)を考えた、信頼性のあるパッケージとして LAPACK(レイパックと発音する) がある。これは、 Oak Ridge National Laboratory,Knoxville, Tennessee と AT&T Bell Laboratories, Murray Hill, NJ に置かれている Netlib サイトからダウンロードできる。日本のミラーサイトとしては phase がある。 LAPACK のユーザマニュアル第 3 版がある。第 2 版は日本語訳がある。 LAPACK はもともと FORTRAN で書かれているライブラリであるが、 CLAPACK はその C 版である。 LAPACK は CPU による違いを BLAS(Basic Linear ALgebraic Library) によって吸収してある。すなわち、線形代数の基本演算を BLAS の関数にし、これを標準化している。そして、 LAPACK 本体を BLAS を用いて書くことによってポータブル化を図っている。自分の使用する CPU に対して最適化されている BLAS(アセンブラ) を用いることによって、 C で書かれた BLAS を用いた場合の数十倍の高速化が達成されている。線形方程式の直接解法に関しては、最適化 BLAS と LAPACK の組み合わせは、最も優れた組み合わせの一つとなる。 Pentium II 上の Linux 用の最適化 BLAS がある。現在はまだ Pentium III 用の BLAS は開発されていないようである。 C のコンパイラの特性を考えてループ展開などの技法を駆使して、 C で書かれているにもかかわらず、コンパイル時にループ展開の数などのパラメータを最適化することによって最適化 BLAS に匹敵するオブジェクトコードを生成する試みもある。 ATLAS や PHIPACK などである。尚、 PHIPACK は dgemm という行列の乗算プログラムだけが対象にされている。原理的には、どのような CPU と OS であっても、 C コンパイラがあれば、 C でかかれた ATLAS や PHIPACK をコンパイルするだけで最適化 BLAS に匹敵するコードが得られる。線形代数計算の直接解法に関しては最も高速で、信頼性を求めるなら LAPACK と最適化 BLAS を実装するのが有力な解となる。 LAPACK++は C++ version であるが、 LAPACK の一部の関数しか実装されていない。また、分散マシン上の LAPACK として ScaLAPACK が開発されている。数値計算パッケージの紹介としてはつぎのようなサイトもある GAMS List of Packages。

注意

1. 数値計算のコードを載せた本として、しばしば利用される向きのある Numerical Recipes は数値計算の専門家でない著者による書であるが、つぎのような批判 Why not use Numerical Recipes 他がある。これらの批判はこの本に載せてあるコードが最適なものに比べて実行速度が数十倍遅く、また、信頼性も低い場合があるという点に集中している。 Netlib などの他の優良なソフトを利用すべきであろう。

2. 本稿では反復解法は扱わないが、数値線形代数の反復解法のサーベイがある。

1.2 数値計算パッケージ

LAPACK などの数値計算ライブラリは高速、高信頼であるが、 C や FORTRAN から呼び出して使うの、プログラム開発に時間がかかる。そこで、 LAPACK などのライブラリを基礎にして、使い易い GUI をもたせた数値計算用インタプリタパッケージが開発されている。 The MathWorks, Inc. から発売されている Matlab はその代表である。欧米では、研究者から、学生までこのようなパッケージを利用して数値計算をするのが一般的である。

MATCOM(Matlab to C++ Compiler は Matlab から C++へのコンパイラ、すなわち Matlab コードから MEX ファイルあるいは C++コードを生成するフリーなツールである。生成された C++コードを最適化された C++コンパイラでコンパイルすると Matlab より数段高速なコードとなる。 MAT<LIB> は Matlab コンパチな C++行列クラスライブラリである。同じサイトからダウンロードできる。 W. Kahan による Matlab’s Loss is Nobody’s Gain は Matlab に対する面白いコメントである。

また、つぎのようなフリーなパッケージソフトがある。 1. Octave は MATLAB に最も近い MATLAB クローンである。 Linux 用のソフ

トであるが、 Windows に実装することもできる。まだ、オブジェクトプログラミングが実装されていない点では MATLAB に譲るが、線型方程式の解や固有値を計算する組み込み関数の計算速度、信頼性では MATLAB に勝るとも劣らないと思われれる。グラフィックスは gnuplot を利用している。また、 ATLAS を使って Octave を高速化することもできる (Tuning up Octave by ATLAS)。 ATLAS による高速化の効果は大変大きく、 Alpha/EV56/533MHz で、 2 つの 1000 × 1000 行列の積の計算に reference BLAS(C で書かれた高速化していない BLAS) と Octave の組み合わせでは 54.7 秒かかっていたものが、 ATLAS を使うと 4.0 秒で計算できるようになったとの報告が書かれている。

2. Scilab は typed list によりオブジェクト指向プログラミングができる環境があり、 MATLAB クローンとしては一番完成している。グラフィックスも独自の優れたものが実装されている。 MAPLE とのリンクもできる。 Linux でも Windows でも実装できる。

3. Tela は大規模数値計算のプロトタイピング言語となることを目指して作られており、グラフィックス、線形代数、 FFT などを装備して偏微分方程式を解くことを念頭にしている。

4. RLaB は行列数学用パッケージ。

5. Euler は実数、複素数、区間などを扱える Matlab クローンな言語。 Linux 版はフリーで、 Windows 版はシェアウエアである。

6. Prophet は Unix 上のライフサイエンス用言語。

7. Yorick は型チェックをしない C ライクな Linux 上の数値計算インタプリタ言語。

8. SpYorick は Yorick のプラグインで疎行列処理をする。

9. PETSc は MPI を用いた偏微分方程式の並列数値計算用のツール。 unix 上の C や C++から呼び出して使う。

10. Ascend はフリーな強い型付けされたオブジェクト指向モデル記述言語。

11. Algae(Alki が前身) は疎系に対して Matlab より高速であることをねらった Matlab クローン的なインタプリタ数値計算言語。

著者は、 Matlab と Octave と Scilab を利用している。 Matlab 互換の Scilab による線形連立方程式の解法の例を次に示す。

-->A=[1 2;2 1]
A  =
!   1.    2. !
!   2.    1. !
-->b=[1;1]
 b  =
!   1. !
!   1. !
-->x=A\b
 x  =
!   0.3333333 !
!   0.3333333 !
-->A*x-b
 ans  =
!   0. !
!   0. !

このように非常に簡単に (かつ高速に) 線形方程式が解ける。詳細は著者による Matalab, Octave, Scilab の使い方の html ファイルなどを参照されたい (そのテキストも出版される [9])。

Demmel の著書 [10] の Matlab コードがあるのでよいプログラム例として参考にされたい。

参考文献

[1] 柏木雅英: “精度保証付き数値計算 [1]–区間解析”, シミュレーション 18 4, 261/267 (1999)

[2] 大石進一: “精度保証付き数値計算 [2]–線形方程式の高速精度保証とプログラミング技法”, シミュレーション 19 1, 39/45 (2000)

[3] 神澤雄智、相馬隆郎: “精度保証付き数値計算 [3]–常微分方程式の精度保証”, シミュレーション 19 2, 30/35 (2000)

[4] 中尾充宏: “精度保証付き数値計算 [4]–偏微分方程式の精度保証”, シミュレーション 19 3 (2000)

[5] 大石進一: 非線形解析入門、コロナ社 (2000)

[6] 中尾充宏、山本野人: 精度保証付き数値計算、日本評論社 (1998)

[7] 大石進一: 精度保証付き数値計算、コロナ社 (2000)

[8] 大石進一: 数値計算、裳華房 (1999)

[9] 大石進一: “Linux 数値計算ツール”, コロナ社 (2000.9)

[10] J. Demmel: “Applied Numerical Linear Algebra”, SIAM (1997)

[11] Shin’ichi Oishi and Siegfried M. Rump: “Fast verification of solutions of matrix equations”, submitted to Numer. Math..

[12] Nicholas J. Higham: “Accuracy and Stability of Numerical Algorithms”, SIAM (1996).

[13] Shin’ichi OISHI: “Fast Enclosure of Matrix Eigenvalues and Singular Values via Rounding Mode Controlled Computation”, to appear in Linear Algebra and its Applications