ハイパフォーマンス コンピューティング: スパース行列等の取り扱い

ハイパフォーマンス
コンピューティング:
スパース行列等の取り扱い
担当: 荻田 武史
2
スパース行列
•  偏微分方程式を差分法や有限要素法で離散化するときなど、
応用問題では非ゼロ要素が極端に少ない行列(スパース行
列)が現れる。
•  たとえば、2次元Poisson方程式(長方形領域)を差分法で離
散化した場合は、以下のような行列が得られる。
青色が、非ゼロ要素
3
スパース行列のデータの持ち方
•  ほとんどの要素がゼロなので、データを圧縮して保持すること
ができる。
•  様々なデータ保持方法があるが、ここではCRSフォーマット
(Compressed Row Storage format)を用いる。
•  (参考)以下のような行列を帯行列と呼ぶ。これも、計算量やメ
モリ量を削減するための手法となる。
保持する要素数が、
n*nから、(帯幅)*n まで削減される。
4
例)
Row_Ptr(i): i行目が何番目の非ゼロ要素
から始まるかを格納。
最後の非ゼロ要素が何番目かを知るた
めに、行数をひとつ多く取っておく。
A(j): 非ゼロ要素のみを行方向
に連続に並べたもののうち、j番
目の要素。
Colmun_Index: A(j)の列番号。
5
スパース行列の行列・ベクトル積
•  行列・ベクトル積 y = Ax のプログラム(MATLAB)
for i=1:n
s = 0;
for j=Row_Ptr(i):Row_Ptr(i+1) - 1
jj = Column_Index(j);
s = s + A(j) * x(jj);
end
y(i) = s;
end
6
演習1
1.  C言語等でスパース行列の行列・ベクトル積を実装してみよ
う。
•  まずは4ページのような例題で試してみよう。
2.  密行列の場合と比べて、計算時間がどれくらい削減される
か、大きい行列に対して試してみよう。
•  たとえば、2ページのような行列で試してみよう。
3.  OpenMPによる並列化についても試してみよう。
7
帯行列のデータの持ち方
q
q
p
L
L
p
n
n
例)
n=5
p=2
q = 1
帯幅L = p + q + 1のとき、
保持する要素数が、n*nから、L*n まで削減。
a(1,2)
a(2,3)
a(3,4)
a(4,5)
a(1,1)
a(2,2)
a(3,3)
a(4,4)
a(5,5)
a(2,1)
a(3,2)
a(4,3)
a(5,4)
a(3,1)
a(4,2)
a(5,3)
8
帯行列のデータの変換
q
q
p
L
L
p
n
n
•  左側の行列Aを右のような配列Bに格納した場合、
A(i,j) = B(i-j+q+1,j)
となる。対角要素については A(i,i) = B(q+1,i) となる。
(実際に、上記が正しいか確かめてみよう)
9
演習2
•  帯行列のデータの持ち方の場合で、軸交換なしのガウスの消
去法を実行するプログラムを作成してみよう。
まずはMATLABで試そう。
2.  MATLABでうまくできたら、C言語等で実装してみよう。
3.  密行列として計算した場合と、計算時間を比較してみよう。
1.