効率的な数値的線形代数計算
浮動小数点データ - ハードウエア浮動小数点データと任意の精度のソフトウエア浮動小数点データの両方 - を含む大規模な行列とベクトルの計算は Maple において数値線形代数ルーチンの組み込みライブラリを利用することによって非常に効率的に実行することができます。 これらのルーチンの多くは Maplesoft と the Numerical Algorithms Group (NAG) の間の提携を通して提供されています。 さらに加えて、CLAPACK と最適化された ATLAS ライブラリの一部が統合されました。最近の改良に対する解説は Maple 7 の効率向上を見て下さい。
LinearAlgebra パッケージにおいて要求された操作を実行するサブルーチンの選択は、次のいくつかの要因に基づきます:
•
|
入力オブジェクトの shape, storage, order オプション
|
•
|
UseHardwareFloats 環境変数の設定
|
•
|
要求された操作を完了するための NAG ルーチンの有効性
|
これらの選択基準は以下のセクションで詳述します。
|
データの型
|
|
入力されたオブジェクト (行列、ベクトル、またはスカラー) が数字データ (すなわち、スカラー (もしあれば)、任意の行列やベクトルの入力においてその成分の型は complex (複素数)(extended_numeric) でなければならない) のみから成るならば、Maple は特別な LinearAlgebra ルーチンを実行するために NAG ルーチンを探すことだけを考えます。また、数値成分またはスカラーの少なくとも 1 つは浮動小数点数でなければなりません。これらの性質が行列またはベクトルの入力に対して満たされているかどうか判定するために、 Maple はまず、行列とベクトルの入力に関して datatype オプションをチェックします。もし、そのデータが浮動小数点 (複素浮動小数点を含む) のある形式であることをこのオプションが示しているときは、対応する成分をそれ以上に調査しません。もし、その datatype オプションがこの情報を提供しないときは、成分が調べられます。
すべての入力されたオブジェクトが数値データだけで構成され、少なくとも 1 つの値が浮動小数点であると判定されると、浮動小数点成分を指定するデータタイプを持たない全ての入力のコピーが作られます。コピーを作る際、最も近い (ソフトウェア浮動小数点の場合は Digits で決定される) 浮動小数点数に変換されます。この操作は高速ですが、明らかに入力されたオブジェクトが適切な (浮動小数点) datatype を持つことを保証することによって避けることができるオーバーヘッドです。これはオブジェクトが作られるとき、Matrix または Vector の datatype オプションにより指定されます。
コピー操作が必要とされるかどうかの判定の方法についての情報はセクション実効時間のフィードバックを参照して下さい。
|
|
Shape, Storage, Order オプションの指定
|
|
NAG ルーチンは入力された行列が Fortran (コラム メジャー) 命令で保存されていることを求めます。Maple においては行列に対してこのことをデフォルトとしていますが、行列を作成する際、それは order オプションを使用することによって無視することが可能です。必要ならば、入力された行列は NAG ルーチンが呼び出される前に、コピーされて Fortran 命令に送られます。
ある場合には、特定の LinearAlgebra 演算を実行するために、異なるいくつかの NAG ルーチンを利用することができます。ここでの違いは入力された行列またはベクトルオブジェクトを表すために使われるデータ構造に基づいています。行列またはベクトルの成分の大きな部分集合が常に 0 であると分かっているならば、これらの成分が保存される必要はなく、よりコンパクトなデータ構造が標準の長方構造の代わりに使用できるというのが、ここで意図するところです。
Matrix と Vector はオプションとして shape と storage を取ります。そしてそれは行列とベクトルを手元にある作業に対して最適の NAG ルーチンに会わせることができます。一般に、Maple は shape にどのような指定をしてもそれと互換性を持つ storage を選択するので、shape を指定すれば通常十分です。詳細は Matrix 、 Vector を参照して下さい。(MatrixOptions と VectorOptions もまた参照のこと。)
shape と storage が指定され、Matrix と Vector オブジェクトに対して要求された操作を実行する NAG ルーチンが利用できないが、より 一般の shape や storage を持つ Matrix または Vector オブジェクトに対する操作に適用できる NAG ルーチンが利用できれば、必要なだけコピーが作られ、より一般的な NAG ルーチンが演算を実行するために呼び出されます。 ( NAG ルーチンの Maple への組み込みは進行中のプロジェクトであるため、将来のリリースには入力の shapes と storage モードのより広い選択肢を扱う能力が含まれます。)
コピー操作が必要とされるかどうかの判定方法についての詳細はセクション実効時間のフィードバックを参照して下さい。
|
|
浮動小数点計算の環境
|
|
UseHardwareFloats 環境変数は、浮動小数点計算が実行される環境を制御する Boolean フラグです。最終的にそれがすべての浮動小数点計算に対してこの環境選択をコントロールする一方、ハードウエア精度の NAG ライブラリ (UseHardwareFloats = true) と任意精度の NAG ライブラリ (UseHardwareFloats = false) のどちらかを選んだり、datatype = float をつけて Matrix と Vector を作る際に datatype オプションを特別な値に設定するためだけに使われます。任意精度の NAG ライブラリは Maple のソフトウエア浮動小数点で直接働くことが出来る特別に作られたバージョンです。
(適切なコンストラクタへの呼び出しによって) Matrix または Vecotr がオプション datatype = float で作成されるならば、 UseHardwareFloats フラグは行列またはベクトルがハードウエア浮動小数点成分 (UseHardwareFloats = true) かソフトウエア浮動小数点成分 (UseHardwareFloats = false) かどうかを見るためにチェックされます。前者の場合、Matrix または Vector は datatype = float[8] で作成されます。これは、成分が 8 バイトのハードウェア浮動小数点数であることを意味します。一方、後者の場合、 datatype = sfloat ("ソフトウェア浮動小数点") で作成されることを意味します。
LinearAlgebra 演算への入力がすべて数値で、(データタイプのチェックを通るように) 少なくとも 1 つの値は浮動小数点数であり、NAG ルーチンがその作業を実行するために利用できると、UseHardwareFloats フラグは NAG ルーチンのハードウエア精度か任意精度のバージョンのどちらかが使われているか判断します。注意: 入力オブジェクトのデータタイプと結果オブジェクトのデータタイプはここでは考慮されません; 計算環境は UseHardwareFloats フラグだけによって制御されます。これは、たとえば、A と B が datatype = sfloat, Digits = 100, UseHardwareFloats = true の行列ならば、A+B の計算はハードウエア精度の NAG ルーチンで実行され、最高でも、ハードウエア精度 (およそ 15 進法の 10 桁) まで正確な結果しか返さないことを意味します。
上で述べたように、呼び出されているルーチンの outputoptions においてデータタイプの値を設定すると、UseHardwareFloats フラグを無効にします。しかしながら、この値が UseHardwareFloats フラグと一致しないならば、返されるオブジェクトを作るために最後にコピー操作が必要となります。
そのようなコピー操作が必要であるかどうかの決定方法に関する情報はセクション 実行時間のフィードバック を参照して下さい。
|
|
利用できる NAG ルーチン
|
|
NAG ライブラリのMapleへの編入は進行中の計画です。どの LinearAlgebra ルーチンに対してもそれを請け合う少なくとも1つの汎用 NAG ルーチンが利用できます。そして多くの LinearAlgebra ルーチンは、その操作の入力オブジェクトの特色によって、いくつかの NAG ルーチンのうちの 1 つに送ることができます。
NAG ルーチンが LinearAlgebra 計算を完了するために使用されるかどうかの決定方法に関する情報はセクション 実行時間のフィードバック を参照して下さい。
|
|
実行時間のフィードバック
|
|
LinearAlgebra ルーチンからユーザーフィードバックメッセージを求めるために、infolevel テーブルを用います。
infolevel[LinearAlgebra] が 1 以上ならば、LinearAlgebra 演算を実行すために NAG ルーチンを使用する度に userinfo メッセージが作成されます。呼び出される NAG ルーチンの名前の前には "hw_" または "sw_" がつきます。"hw_" はハードウエア精度の NAG ルーチンが呼び出されることを意味しています。"sw_" は、任意精度("ソフトウエア") NAG ルーチンが呼び出されることを意味しています。何らかの理由で、NAG ルーチンがこの計算を完了することができないときは、この結果に対する userinfo メッセージが表示されます。
infolevel[LinearAlgebra] が 2 以上ならば、入力または出力の Matrix または Vector の shape, storage, order, datatype のいずれかを変更するためにコピー操作が要求され、その度に userinfo メッセージが作成されます。
|
|
さらにShape と Storage について
|
|
一般に、LinearAlgebra の関数は、datatype=sfloat または datatype=float[8] を持ち、デフォルトの shape (none) と storage (rectangular) のみを持つ入力オブジェクトに対して計算を行うために NAG ルーチンを使います。いくつかの LinearAlgebra ルーチンは別の shape (例えば、triangular[upper]) および互換性のある storage (たとえば、shape=triangular[upper] で storage=triangular[upper] または storage=rectangular ならば、コピー操作の必要なしで NAG ルーチンを使うことができます)または疎な storage を持つ入力オブジェクトに対しても NAG ルーチンを使うことができます。
ある入力オブジェクトが認識されない shape または複数の shape を持つならば、LinearAlgebra のルーチンは使用できる NAG ルーチン (もしあれば) に対するフォーマットでオブジェクトをコピーします。これは、NAG ルーチンが根底にある物理的に格納されたデータだけを受け取り、オブジェクトの shape 情報を受け取らないからです。(shape 情報は indexing function と呼ばれる Maple ルーチンに格納され、この indexing function は NAG ルーチンに渡されません。)
複数の shape、または認識されない shape (例えば、ユーザー定義の shape) が入力されたオブジェクトに付加されるとき、Maple は根底にある物理的に格納されたデータがその shape 情報の補助無しに正しく解釈されることができないと判断します。
より広い範囲の入力オブジェクトに対して NAG サポートを加えて、LinearAlgebra の能力を拡張することは、進行中の計画です。
|
|
パックされたフォーマット
|
|
LinearAlgebra パッケージのルーチンの一部は、特別な "packed" フォーマットに格納されるオブジェクトを返すことができます。たとえば、 LinearAlgebra:-LUDecomposition(..) はベクトル (ピボット) と ("L" と "U" の要素を両方含む) 1 つの行列の形式で入力を分解された形式で返すことができます。これは基礎となる NAG ルーチンがその結果を返す形式となっています。したがって、浮動小数点での分解 (ただし、NAG ルーチンはサブルーチンとして呼び出される)、ここで分解が NAG ルーチンによって次の演算で使われると、この packed 形式を要求するのがもっとも効率的です。
たとえば、A は正方行列、 b は次数が合うベクトルと仮定します。連立方程式 A . x = b の解法は LinearSolve を 1 つ呼び出すことによって完了します。
>
|
with( LinearAlgebra ):
A := <<2.2,4.0,-6.1>|<0,1.8,3>|<-3.2,-5,7.4>>;
|
| (7.1) |
| (7.2) |
| (7.3) |
しかしながら、後の使用のために分解された形式を保存することが必要ならば、(例えば、同じ係数行列で右辺が違うような方程式を解くような場合)、方程式の解法は以下のステップにしたがって完了することができます。
>
|
p, lu := LUDecomposition( A, output='NAG' );
|
| (7.4) |
>
|
x := LinearSolve( [p, lu], b );
|
| (7.5) |
>
|
x2 := LinearSolve( [p, lu], <3.1,1.3,-.27> );
|
| (7.6) |
分解のパックされた形式を保存することによって、不必要なコピー操作と分解の繰り返しはこの計算の手順においては実行されません。
|
|
より高速なアクセス
|
|
数値線形代数計算のための LinearAlgebra パッケージを用いてユーザーが開発したプログラムがプログラミング・レイヤ・インタフェースを通して関連する LinearAlgebra ルーチンにアクセスすることが可能です (そしてそうすることが意図されています)。正しく検出されていない無効なデータにより起こりうる問題を犠牲にして、これはパラメータ確認とデフォルトの設定のオーバーヘッドを避けます。
そのようなプログラムを開発することへの手助けとして、プログラミング・レイヤ・ルーチンが判定を ("kernelopts( ASSERT=true )" または "kernelopts('assertlevel' = 1)" により) assertion を有効にすることでプログラム・レイヤ・ルーチンはパラメータに対し、いくつかの基本的な型チェックを行います。
プログラミング・レイヤ・インタフェースの詳細は LinearAlgebra パッケージでのプログラミング を参照して下さい。
|
|
Download Help Document
Was this information helpful?