15

私はスパース行列のために私自身の線形方程式ソルバを書く必要があります。 私はスパース行列に任意の種類のデータ構造を自由に使用できます。また、共役勾配を含むいくつかの解決策を実装する必要があります。クラスの高速スパース行列乗算

ベクトルとの乗算が比較的速いような疎行列を保存する有名な方法があるのだろうかと思っていました。

現在のところ、私のスパース行列は基本的にラップされたstd::map< std::pair<int, int>, double>で、データがあればそれを格納します。これは、各行列要素のルックアップを実行しなければならないので、ベクトルからO(n²)複雑度への行列の乗算をO(n²log(n))に変換します。 私はYale Sparse行列形式を調べました。要素の取得もO(log(n))にあるようですので、はるかに速いかどうかはわかりません。

参考のために、私は5000エントリで入力された800x800のマトリックスを持っています。共役勾配法を用いてこのようなシステムを解決するには、およそ450秒かかります。

別のデータ構造を使用するともっと早くできると思いますか?

ありがとうございました!

+0

まずはウィキペディアを読んでください。 http://en.wikipedia.org/wiki/Sparse_matrixには、効率的な操作を提供する一般的なストレージメソッドの良いリストがあります。 –

+0

@ソンワン:クラスの目的は、基本的には独自の有限要素法ソルバーをロールバックすることです – lezebulon

答えて

22

最も一般的な選択肢はCSC or CSR storageです。これらは両方とも行列 - ベクトル乗算に効率的である。あなた自身で行う必要がある場合は、それらの乗算ルーチンをコーディングすることも非常に簡単です。 Yaleストレージはまた、非常に効率的な行列 - ベクトル乗算をもたらすとも述べている。行列要素のルックアップを実行している場合は、その形式の使い方を誤解しています。行列ベクトルの乗算がどのように実装されているかを知るために、いくつかの標準的な疎なライブラリを調べることをお勧めします。

現在のストレージでも、O(n)の複雑さで行列乗算を実行できます。私が今まで見たスパース行列 - ベクトル乗算アルゴリズムはすべて同じステップにまで沸き立っています。例えばy = Axとする。

  1. 結果ベクトルをゼロにします。
  2. 行列の非零要素のイテレータを初期化する。
  3. 行列の次の非零要素を得る。A [i、j]は言う。 i、jのパターンは規則的なパターンに従わないことに注意してください。単にAの非ゼロ要素が格納される順序を反映しています。
  4. Y [i]が+ = A [I、J] * X [j]が
  5. Aの多くの要素がある場合は、後藤3.

私はあなたがループのための古典的なダブルを書いている疑いがあります高密度乗算コード:

for (i=0; i<N; i++) 
    for (j=0; j<N; j++) 
     y[i] += A[i,j]*x[j] 

これは、ルックアップの実行につながっているものです。

しかし、私はあなたがあなたのstd::mapのストレージをつけることを示唆していません。それは超効率的になることはありません。主にCSCが最も広く使われているので、私はCSCをお勧めしたいと思います。

+0

"行列要素のルックアップを実行する場合、形式を使用する方法を誤解しています" 私の元の問題。 CSRの場合は、同じルックアップの複雑さがありますが、行列のベクトル乗算を行うためにルックアップを行う必要はありません。 – lezebulon

+0

を編集して初めて配列をトラバースするだけです。これも同様に機能しますが、私のポイントはすでに行の最初の右にソートされているためです。 – lezebulon

+0

"アレイを一度トラバースする"。それはまさにそれです。あなたは今それを持っている!ちょっとした変化がありますが、そのように考えると、あなたはまっすぐに家にいます! –