数千の小さな行列 (前に書いたように 4x3 ではなく 8x9) の nullspace を並列に計算する必要があります (CUDA)。すべての参照はSVDを指していますが、数値レシピのアルゴリズムは非常に高価に見え、実際には必要のないnullスペース以外に多くのものを与えてくれます. ガウスの消去法は本当にオプションではないのですか? 他に一般的に使用される方法はありますか?
7 に答える
あなたの質問に直接答えるために...はい!QR分解!
Aをランクnのm行n列の行列とします。QR分解は、A=QRとなる正規直交m行m列行列Qと上三角m行n列行列Rを見つけます。Q = [Q1 Q2]を定義すると、Q1はm x n、Q2はm x(mn)であり、Q2の列はA^Tの零空間を形成します。
QR分解は、グラムシュミット、ギブンス回転、またはハウスホルダー反射のいずれかによって計算されます。それらは異なる安定性特性と操作カウントを持っています。
あなたは正しいです:SVDは高価です!最先端のものが何を使っているのかは言えませんが、「ヌル空間を計算する」(編集:わかりやすい方法で)と聞くとQRだと思います。
上記の提案された方法が常にヌルスペース全体を与えるとは思いません。要約すると、「A = QR、ここで Q = [Q1 Q2]、Q1 は m x n、Q2 は m x (mn) です。その後、Q2 の列は A^T のヌル空間を形成します。」
実際、これはヌル空間の部分空間しか与えないかもしれません。単純な反例は、A=0 の場合です。この場合、A^T のヌル スペースは R^m 全体です。
したがって、R もチェックする必要があります。Matlab での私の経験に基づいて、R の行がストレート 0 の場合、Q の対応する列も A^T のヌル空間の基礎になるはずです。明らかに、この観察結果はヒューリスティックであり、QR 分解に使用される特定のアルゴリズムに依存しています。
ガウス消去法は、4x3 行列に対して十分高速です。IIRC 並列処理なしで Java を使用して、1 秒あたり約 500 万回実行しました。このような小さな問題の場合、最善の策は、ルーチン (行削減など) を自分でコーディングすることです。そうしないと、データを外部ルーチンに適した形式に変換するのにほとんどの時間が無駄になります。
行列が単なるランダムではなく関連しているかどうか疑問に思ったので、探しているヌル空間は N 空間 (N = 9) の曲線に対する 1 次元の接線のように見なすことができます。その場合、ニュートン法を使用して、前のヌル空間ベクトルから開始して、二次方程式 Ax = 0、|x|^2 = 1 のシステムの連続するインスタンスを解くことにより、処理を高速化できる場合があります。ニュートン法は一次導関数を使用して解に収束するため、ガウス消去法を使用して 9x9 システムを解決します。この手法を使用すると、パラメーターを変更するなどして、マトリックスからマトリックスへの小さなステップを作成できる必要があります。
したがって、最初のマトリックスで SVD を使用して初期化しますが、その後、次の反復の開始点として 1 つのヌル空間ベクトルを使用して、マトリックスからマトリックスへとステップします。収束するには、1 回または 2 回の反復が必要です。収束しない場合は、SVD を使用して再起動します。このような状況である場合、各マトリックスで新たに開始するよりもはるかに高速です。
私はずっと前にこれを使用して、電力システムの動作に関連する 50 x 50 の二次方程式のセットの解の等高線をマッピングしました。
CUDA にとって最も重要なことは、条件分岐 (グラフィックス ハードウェアでは非常に遅い) に依存しないアルゴリズムを見つけることだと思います。条件付き代入に最適化できる単純な if ステートメントの方がはるかに優れています (または ?: 演算子を使用できます)。
必要に応じて、条件付き割り当てを使用して何らかの形式のピボットを実行できる必要があります。結果を格納する方法を決定するのは実際には難しいかもしれません: 行列が階数不足の場合、CUDA プログラムにそれに対して何をしてもらいたいでしょうか?
4x3 行列が実際にはランク不足ではないと仮定すると、(単一の) ヌル空間ベクトルを条件なしで見つけることができます。行列は十分に小さいため、Cramer の規則を効率的に使用できます。
実際には、null ベクトルのスケールは実際には気にしないので、行列式で割る必要はありません。マイナーの行列式を取ることができます。
x1 x2 x3
M = y1 y2 y3
z1 z2 z3
w1 w2 w3
|y1 y2 y3| |x1 x2 x3| |x1 x2 x3| |x1 x2 x3|
-> x0 = |z1 z2 z3| y0 = -|z1 z2 z3| z0 = |y1 y2 y3| w0 = -|y1 y2 y3|
|w1 w2 w3| |w1 w2 w3| |w1 w2 w3| |z1 z2 z3|
これらの 3x3 行列式は単なる三重積であることに注意してください。外積を再利用することで計算を節約できます。