5.1 疎行列の格納形式

5.1.1 疎行列用の線形計算プログラムの構成

疎行列用の格納形式には使用目的に適した多くの種類があり, 絶対的な標準形式といったものはない.

そこで, 疎行列用の線形計算プログラムでは疎行列の格納形式に依存しないような作りになっているものが多い.

具体的には, リバースコミュニケーションやサブルーチン形式により疎行列とベクトルの積の計算を外付けにすることにより, アルゴリズムの実装から格納形式に依存する部分を切り離すといったものである. これにより, 線形計算プログラムを変更することなく, 数多くある格納形式から好みのものを選択して使用することができるようになる.

なお, 最近のプログラムでは, 格納形式の詳細は隠してオブジェクトとして扱うものもある. この場合は格納データの入出力方法 (メソッド) だけが公開される.

5.3 で説明するクリロフ部分空間法では, 次の 2 つの操作を定義すればよいようにアルゴリズムが記述されている. なお, 係数行列を A, 前処理行列を M, ベクトルを x, y で表す.

  • Matvec演算: \(Ax\) を求める
  • Psolve演算: \(x = M^{-1}y\) (\(Mx = y\) の解 \(x\)) を求める (前処理なしの場合には不要)

一部の解法 (BICG 法と QMR 法) ではさらに次の 2 つが必要である. 説明の中で「\(A^T\) に関する計算」と表記される場合にはこれらを指す.

  • MatvecTrans演算: \(A^Tx\) を求める
  • PsolveTrans演算: \(x = (M^T)^{-1}y\) (\(M^Tx = y\) の解 \(x\)) を求める (前処理なしの場合には不要)

XLPack のクリロフ部分空間法では 1 つのアルゴリズムに対して 3 種類のプログラムを用意している. すなわち, ① リバースコミュニケーション版, ② サブルーチン版, ③ CSR / CSC 版である. ③は格納形式を CSR / CSC 形式 (下記参照) に限定する代わりにプログラムを容易に使用できるようにしたものである.

5.1.2 疎行列の格納形式

疎行列用の線形計算プログラムが任意の疎行列格納形式を使うことができるように作られているとはいえ, できるだけメジャーな形式を使用しておいた方がよい.

ここでは汎用的に使える次の代表的な 3 種類について説明する.

(1) 圧縮行格納 (Compressed Sparse Row (CSR)) 形式
(2) 圧縮列格納 (Compressed Sparse Column (CSC)) 形式
(3) 座標 (Coordinate (COO)) 形式

これらは XLPack で使用されている形式である. (1) と (2) は主に計算に使われ, (3) は主に入出力に使われる.

5.1.2.1 例題疎行列

以下, 格納形式の説明で使用する例題として次の 6 x 6 行列を考える.

10 0 0 0 -2 0
3 9 0 0 0 3
0 7 8 7 0 0
3 0 8 7 5 0
0 8 0 9 0 13
0 4 0 0 2 -1

5.1.2.2 圧縮行格納 (Compressed Sparse Row (CSR)) 形式

CRS (Compressed Row Storage) 形式とよばれることもある. 疎行列の構造について何も仮定しない一般的な形式であり広く用いられている. 特定の構造や計算に特化したような形式に比べると効率は落ちるものの汎用性を特長とする.

1 本の浮動小数配列 val と 2 本の整数配列 rowptr および colind に非ゼロ要素だけを格納する. val には非ゼロ要素の値を行方向順に格納する. rowptr は各行の先頭要素の val 配列中の位置を表す. colind は対応する val の要素の列番号である. rowptr の最後には val 配列の最後 + 1 を入れておく.

行番号と列番号の振り方には 2 種類あるので注意が必要である, すなわち, FORTRAN プログラムの場合には 1 から, C プログラムの場合には 0 から始まる.

上の例題行列を CSR 形式で表すと次のようになる.

(FORTRAN 形式)

i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
val(i) 10 -2 3 9 3 7 8 7 3 8 7 5 8 9 13 4 2 -1
colind(i) 1 5 1 2 6 2 3 4 1 3 4 5 2 4 6 2 5 6
rowptr(i) 1 3 6 9 13 16 19

(C 形式)

i 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
val[i] 10 -2 3 9 3 7 8 7 3 8 7 5 8 9 13 4 2 -1
colind[i] 0 4 0 1 5 1 2 3 0 2 3 4 1 3 5 1 4 5
rowptr[i] 0 2 5 8 12 15 18

対称行列の場合には, 対角要素とその左側のみ, あるいは, 対角要素とその右側のみを格納して使用メモリーを減らすことができる. XLPack ではこれを独自に SSR 形式と表示している.

上の例では各行の中で列番号の昇順に並んでいるが, プログラムによっては行内の並び順は任意でよいことがある. XLPack では基本的には昇順に並んでいることを仮定して計算効率を上げるようなプログラムになっている.

CSR 形式は実は脆弱なデータ構造で, rowptr および colind は整数値ではあるが val を指すポインターみたいなものであるから, 誤った値を設定すると簡単にメモリーアクセスエラーを引き起こしてしまうので注意が必要である.

5.1.2.3 圧縮列格納 (Compressed Sparse Column (CSC)) 形式

CCS (Compressed Column Storage) 形式とよばれることもある. CSR 形式とは転置の関係にある. すなわち, CSR 形式で格納した行列 A を CSC 形式とみなして読みだすと転置行列 AT になる.

1 本の浮動小数配列 val と 2 本の整数配列 colptr および rowind に非ゼロ要素だけを格納する. val には非ゼロ要素の値を列方向順に格納する. colptr は各列の先頭要素の val 配列中の位置を表す. rowind は対応する val の要素の行番号である. colptr の最後には val 配列の最後 + 1 を入れておく.

行番号と列番号の振り方が 2 種類あるのは CSR 形式と同様である, すなわち, FORTRAN プログラムの場合には 1 から, C プログラムの場合には 0 から始まる.

上の例題行列を CSC 形式で表すと次のようになる.

(FORTRAN 形式)

i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
val(i) 10 3 3 9 7 8 4 8 8 7 7 9 -2 5 2 3 13 -1
rowind(i) 1 2 3 2 3 5 6 3 4 3 4 5 1 4 6 2 5 6
colptr(i) 1 4 8 10 13 16 19

(C 形式)

i 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
val[i] 10 3 3 9 7 8 4 8 8 7 7 9 -2 5 2 3 13 -1
rowind[i] 0 1 2 1 2 4 5 2 3 2 3 4 0 3 5 1 4 5
colptr[i] 0 3 7 9 12 15 18

対称行列の場合には, 対角要素とその左側のみ, あるいは, 対角要素とその右側のみを格納して使用メモリーを減らすことができる. XLPack ではこれを独自に SSC 形式と表示している.

上の例では各列の中で行番号の昇順に並んでいるが, プログラムによっては列内の並び順は任意でよいことがある. XLPack では基本的には昇順に並んでいることを仮定して計算効率を上げるようなプログラムになっている.

CSC 形式は CSR 形式と同様に colptr および rowind に誤った値を設定すると簡単にメモリーアクセスエラーを引き起こしてしまうので注意が必要である.

5.1.2.4 座標 (Coordinate (COO)) 形式

1 本の浮動小数配列 val と 2 本の整数配列 rowind および colind に非ゼロ要素だけを格納する. val には非ゼロ要素の値を格納する. rowind, colind には対応する行番号, 列番号を入れる.

行番号と列番号の振り方には 2 種類あるので注意が必要である, すなわち, FORTRAN プログラムの場合には 1 から, C プログラムの場合には 0 から始まる.

上の例題行列を COO 形式で表すと次のようになる.

(FORTRAN 形式)

i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
val(i) 10 3 3 9 7 8 4 8 8 7 7 9 -2 5 2 3 13 -1
rowind(i) 1 2 4 2 3 5 6 3 4 3 4 5 1 4 6 2 5 6
colind(i) 1 1 1 2 2 2 2 3 3 4 4 4 5 5 5 6 6 6

(C 形式)

i 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
val[i] 10 3 3 9 7 8 4 8 8 7 7 9 -2 5 2 3 13 -1
rowind[i] 0 1 3 1 2 4 5 2 3 2 3 4 0 3 5 1 4 5
colind[i] 0 0 0 1 1 1 1 2 2 3 3 3 4 4 4 5 5 5

上の例では列番号, 行番号の昇順に並んでいるが, プログラムによって並び順は任意でよいことがある.

COO 形式は入出力に向いた形式であり, そのまま計算に使うと計算効率が悪い. そのため, COO 形式を入力用に使い, 計算時には CSR など他の形式に変換して使うことが多い.

5.1.2.5 基本演算プログラムと変換プログラム

それぞれの疎行列格納形式に対しては, 基本的な演算を行う基本演算プログラムが用意されていることが多い. 標準的な基本演算セットもいくつか提案されていが決定版はないようである.

XLPack では, CSR 形式の行列に対して次のような基本演算プログラムが用意されている.

 
  CsrDusmv  y <- αAx + βy または y <- αATx + βy
  CsrDussv  Ax = b または ATx = b の解 (三角行列)
  CsrDusmm  C <- αAB + βC または C <- αATB + βC
  CsrDussm  AX = B または ATX = B の解 (三角行列)
  SsrDusmv  y <- αAx + βy (対称行列)
  CsrTrans  行列の転置
  CsxDiag  行列の対角要素
  CsxDiagInd  行列の対角要素のインデックス
  CsxSort  行列の要素のソート

 
また, 次の変換プログラムが用意されている.

 
  CooCsr    COO -> CSR
  CsrCoo    CSR -> COO
  SsrCsr    SSR (CSR 対称行列形式) -> CSR (対称なフル行列)
  CsrSsr    CSR (対称なフル行列) -> SSR (CSR 対称行列形式)
  CsrDense  CSR -> 密行列
  DenseCsr  密行列 -> CSR

 
同様の基本演算プログラムと変換プログラムが CSC 形式にも用意されている. また, 複素数型の同様のプログラムも用意されている.