マルチコア・メニーコアプロセッサ向けの
行列計算アルゴリズム
－ 固有値・特異値計算を中心に －
2009年 4月 22日
研究会「アルゴリズムによる計算科学の融合と発
展」
筑波大学 計算科学研究センター
山本 有作 （名古屋大学）
はじめに： 本研究で扱う問題
• 固有値分解 A = XLX–1
– A： n×n 非対称密行列，L： 対角行列，X： 正則行列
– 応用分野
• 構造計算
• 量子化学
• 特異値分解 A = US VT
– A： n×n 密行列，S： 対角行列，U, V： 直交行列
– 応用分野
• 統計計算，信号処理，画像処理
• 電子状態計算 （filter diagonalization 法）
※ 対称密行列の固有値分解は特異値分解と共通点が多い
（1/48）
本研究の目的
• マルチコア／メニーコアプロセッサ上で高性能を達
成できる固有値・特異値分解ソルバの開発
• 背景
– 問題の大規模化
– マルチコア／メニーコア／GPGPUの普及
Quad Core Opteron
ClearSpeed CSX600
NVIDIA GeForce8800 GTX
（2/48）
目次
1. はじめに
2. マルチコアプロセッサでの特異値分解
3. GPUでの特異値分解
4. メニーコアプロセッサでの固有値分解
5. 終わりに
（3/48）
マルチコアプロセッサでの
特異値分解
（4/48）
標準的な特異値分解アルゴリズム
密行列 A
計算内容
計算手法
二重対角化
U0TAV0 = B
ハウスホルダー法
(U0, V0: 直交行列)
二重対角行列 B
二重対角行列の
特異値・特異ﾍﾞｸﾄﾙ計算
Bvi =σi xi
BTxi =σi yi
Bの特異値 {σi }，
特異ベクトル {xi }{yi }
vi = V0 yi
逆変換
ui = U 0 x i
QR法
分割統治法
MR3アルゴリズム
I-SVDアルゴリズム
逆変換
Aの特異ベクトル {ui }, {vi }
（5/48）
各部分の演算量と実行時間
演算量
密行列 A
実行時間(全特異ﾍﾞｸﾄﾙ)
n = 2500，Quad Core Xeon ×2
LAPACK での実行時間（秒）
(8/3) n3
二重対角化
二重対角行列の
特異値・特異ﾍﾞｸﾄﾙ計算
逆変換
O(n2)
～ O(n3)
実
行
時
間
（
秒
）
4mn2
（左右 m 本ずつの
特異ベクトル）
Aの特異ベクトル {ui }, {vi }
コア数
・二重対角化が実行時間の
大部分を占める
・速度向上率が低い
（6/48）
二重対角化の性能が出ない原因
• 二重対角化の演算パターン
– 左右からのハウスホルダー変換による行・列の消去
• A(k) := (I – aw wT ) A(k) = A(k) – aw (wTA(k))
行列ベクトル積
Rank-1更新
• 演算は level-2 BLAS（行列ベクトル積と rank-1 更新）
0
0
A(k)
0
右から
の変換
0
0
左から
の変換
0
非ゼロ要素
ゼロにしたい部分
影響を受ける部分
k
• 演算パターンに関する問題点
– Level-2 BLAS はデータ再利用性が低い。
キャッシュの有効利用が困難であり，単体性能が出にくい。
プロセッサ間のアクセス競合により，並列性能向上も困難
（7/48）
マルチコアプロセッサのメモリ階層
データ転送速度
非常に大
レジスタ
8～128 ワード
コア1
コア2
データ転送速度大
本当はキャッシュも
複数階層からなるが，
ここでは単純化
キャッシュ
数100Kバイト ～ 数Mバイト
帯域の
奪い合い
主メモリ
データ転送速度小
数100Mバイト ～ 数Gバイト
• 主メモリアクセスの場合，実効性能は実質的に転送速度で決まる
• Byte/Flop値: 1演算を行う時間で転送できるバイト数 ～
<1
• マルチコアでは複数のコアが帯域を奪い合うため，問題が更に深刻化
• 主メモリのデータ転送速度の遅さをカバーするには，キャッシュにデータ
がある間に演算をまとめて行う（データの再利用）ことが必要
（8/48）
BLASの利用による高性能化
• BLASとは
– Basic Linear Algebra Subprograms の略
– 個々のマシン向けにチューニングした基本行列演算のライブラリ
• BLASの種類
– Level-1 BLAS: ベクトルとベクトルの演算
• 内積 c := xT y
• AXPY演算 y: = ax + y など
– Level-2 BLAS: 行列とベクトルの演算
• 行列ベクトル積 y := Ax
• 行列のrank-1更新 A := A + xyT
– Level-3 BLAS: 行列と行列の演算
• 行列乗算 C := AB など
= A ×
A = A
＋
×
C = A × B
（9/48）
BLASにおけるデータ再利用性と並列粒度
• Level-1 BLAS
– 演算量： O(N)， データ量： O(N)
– データ再利用性： なし
– 並列粒度： O(N/p) （N： ベクトルの次元，p： プロセッサ台数）
• Level-2 BLAS
– 演算量： O(N2)， データ量： O(N2)
– ベクトルデータのみ再利用性あり
– 並列粒度： O(N2/p)
=
A ×
A = A
＋
×
• Level-3 BLAS
– 演算量： O(N3)， データ量： O(N2)
– O(N)回のデータ再利用が可能
C = A × B
• 行列乗算のサイズを大きくするほど，再利用性が向上
• Byte/Flop 値の小さいマシンでも性能を出せる。
– 並列粒度： O(N3/p)
演算をできる限り level-3 BLAS で行うことが高性能化のポイント
（10/48）
Level-3 BLAS に基づく特異値分解アルゴリズム
• ２段階の二重対角化アルゴリズム（Bischof et al., ’93）
– 密行列 A をまず帯幅 L の下三角帯行列 C に変換
– 次にこの帯行列を下二重対角行列 B に変換
次数 n
下三角
帯行列化
(8/3)n3
0
8n2L
0
A
村田法
C 帯幅 L
0
0
B
• 二重対角化を２段階で行うことの利点
– 前半の変換は，level-3 BLAS （行列乗算）のみを使って実行可能
キャッシュの有効利用が可能
– 後半の変換は level-2 BLAS が中心だが，演算量は O(n2L)
• 前半部に比べてずっと小さい。
（11/48）
下三角帯行列化のアルゴリズム
ブロック鏡像変換によるブロック列の消去
– ブロック鏡像変換 H = I – WaWT
ブロックベクトル
• H は直交行列
行列
• 与えられたブロックベクトルを上三角
行列（正確には右上三角部分のみ
非零でそれ以外が零の行列）に変形
第 K ステップでの処理
0
0
0
0
右からHKR
を乗算
非ゼロ要素
左からH
を乗算
0
左からHKL
を乗算
ゼロにしたい部分
・ 消去演算は行列乗算のみで行える
・ 行列乗算のサイズ ～ L
0
影響を受ける部分
（12/48）
Level-3 BLAS に基づく特異値分解アルゴリズム
• 特異ベクトル計算を含めたアルゴリズムの全体像
n
0
L
(8/3)n3
A
A の特異ﾍﾞｸﾄﾙ
{ui }{vi }
0
8n2L
0
C
4mn2
C の特異ﾍﾞｸﾄﾙ
{zi }{wi }
0
B
4mn2
B の特異ﾍﾞｸﾄﾙ
{xi }{yi }
特異値
{σi }
QR法
DC法
MR3
I-SVD
• 長所
– O(n3) の演算量を持つ部分はすべて level-3 BLAS で実行可能
• 短所
– 逆変換の演算量が 8mn2 （従来法の2倍）。ただし level-3 化可能
（13/48）
性能評価
• 評価条件
–
–
–
–
–
問題サイズ n : 250, 5000, 7500
帯幅 L
: 25, 50, 100
コア数
: 1, 2, 4, 8
Level-3 BLAS に基づく手法は Fortran で作成
従来法はLAPACKのルーチンを使用
• 計算機環境
–
–
–
–
Xeon X5355 (2.66GHz) Quad-core ×2ソケット
Red Hat Enterprise WS release 4
Intel Math Kernel Library
Intel Fortran Compiler 9
（14/48）
実験結果
• Xeon での実行時間
n=2500
n=5000
n=7500
実
行
時
間
（
秒
）
コア数
◆ LAPACK ■ Level-3 帯幅25 ▲ Level-3 帯幅50 × Level-3 帯幅100
• 結果（LAPACK(1CPU)に対する高速化率）
– LAPACK
– Level-3 BLAS
：8コアで最大 1.38倍
：8コアで最大 5.06倍
• 帯幅 L による性能差
– L = 100 が最も性能が良い
（15/48）
実験結果
• LAPACKとLevel-3 BLASの比較1
n = 2500，帯幅 L = 100
実
行
時
間
（
秒
）
1.38倍
4.01倍
コア数
LAPACK
Level-3 BLAS
• 結果(コア数による高速化率)
– LAPACK
– Level-3 BLAS
： 8コアで1.38倍
： 8コアで4.01倍
（16/48）
実験結果
• LAPACKとLevel-3 BLASの比較2
n = 7500，帯幅 L = 100
実
行
時
間
（
秒
）
1.31倍
4.36倍
コア数
LAPACK
Level-3 BLAS
• 実効時間の内訳
– 村田法（帯行列の二重対角化）の時間は L = 100 でも十分小さい
– 村田法の逆変換の実行時間が大きい
今後改良要
（17/48）
実験結果
• 対称行列の固有値分解の結果
n = 9000，帯幅 L = 100
1200
1200
帯化の逆変換
村田法の逆変換
分割統治法
村田法
帯行列化
三重対角化の逆変換
実
行
時
間
（
秒
）
1000
1000
分割統治法
三重対角化
800
800
600
600
400
400
200
200
0
0
コア数
1
2
4
8
LAPACK
1
2
4
8
Level-3 BLAS
• 結果
– 特異値分解とほぼ同じ傾向が見られる。
（18/48）
まとめ
• 特異値分解の標準的なアルゴリズムでは，二重対角化の部
分で level-2 BLAS を多用するため，マルチコアプロセッサ
で十分な性能を得られない。
• Level-3 BLAS を用いたアルゴリズムでは，演算量は増加す
るが，キャッシュの利用効率向上，アクセス競合の回避によ
り，従来法より高い性能が得られる。
• 対称行列の固有値分解についても，level-3 BLAS を用いた
アルゴリズムは効果的。
（19/48）
GPUでの特異値分解
（20/48）
研究背景
◆ GPGPU （General Purpose GPU）
一般の科学技術計算へのGPUの利用
• 単精度計算ではCPUの10倍以上の性能
• CPUを上回る性能向上ペース
• 開発環境の整備
NVIDIA GeForce8800 GTX
（単精度演算：345.6GFLOPS）
◆ CUDA （Compute Unified Device Architecture）
GPGPUのための統合開発環境
• 2006年にNVIDIA社が発表
• C/C++の拡張によりGPGPUのプログラミングが可能
• GPU向けにチューニング済みのライブラリ（BLAS，FFT）
CUDAの利用例：行列計算，重力多体計算，分子軌道計算，流体計算など
（Cf. http://www.nvidia.co.jp/object/cuda_home_jp.html）
（21/48）
CUDAによる行列計算の高速化
◆ GPU向けにプログラムの移植
• 拡張されたC言語でコーディングして，nvccでコンパイル
• 自由度の高いプログラミングが可能
• GPUの性能を十分に引き出すには様々なチューニングが必要
（スレッド数，条件分岐，メモリアクセスの連続性など）
◆ CUBLASの利用
• CUDAで提供されているBLAS（Basic Linear Algebra Subprograms）
• 標準のC言語のプログラム中で利用可能 （gccでコンパイル可能）
• 限られた基本演算（行列ベクトル積，行列乗算など）のみ
• GPU向けにチューニング済み
今回はできるだけCUBLASを使って高速化を行う
（22/48）
CUBLASの特徴
◆ 仕様の特徴
メインメモリ
• ユーザー自身がデータ転送を制御
PCI-Express
• ボードのデータはプログラム終了まで保持
(1)Set Data
(3)Get Data
GPU用メモリ
(2)SGEMM etc.
メモリ配置の工夫によるデータ転送コストの削減
GPU
グラフィックボード
◆ 性能
（GFLOPS）
※転送時間は含んでいない
140
120
100
80
60
40
20
0
SGEMV（行列ベクトル積，level-2 BLAS）
× ＝
SGEMM（行列乗算，level-3 BLAS）
×
0
1000
2000
3000
4000
GeForce8800 GTX & CUBLAS ver. 1.0
5000 （Size）
＝
Level-2 と 3 の性能差が
CPU以上に大きい
Level-3 BLAS活用の必要性
（23/48）
特異値分解の計算手順と特徴
(a) 二重対角化：
• U1, V1：直交行列， B：下二重対角行列
A
• 大半がBLASルーチンにより計算可能
• 演算量：O(n3)
B
(b) 二重対角行列の特異値分解：
• 様々なアルゴリズム（QR法，分割統治法，MR3，I-SVDなど）
• 丸め誤差の影響を受けやすい
• 演算量： O(n2) ～ O(n3)
(c) 逆変換：
• 大半がBLASルーチンにより計算可能
• 演算量：O(n3)
◆ 計算時間の内訳
100%
逆変換
Bの特異値分解（I-SVD）
二重対角化
80%
60%
40%
Core2 Duo (1.86GHz)
& Intel MKL ver. 8.1
20%
0%
1280
2560
3840
5120
(n)
（24/48）
高速化の方針
◆ GPUの使い方
• できるだけCUBLASを使う
行列乗算をできるだけ利用する
Bischof のアルゴリズムの利用
メインメモリとGPU用メモリ間のデータ通信コストを抑える
行列データは基本的にGPUに置き，必要な部分のみCPUに転送
• GPU向けのプログラムの移植は最小限にする
◆ 特異値分解の計算
• 二重対角化と逆変換の計算にはGPUを利用する
直交変換が中心のため，丸め誤差に強い
単精度計算でも安定性の問題がない
• Bの特異値分解の計算はCPUのみで行う
非線形方程式の求解あるいは漸化式の計算
安定性の面から倍精度計算が望ましい
（25/48）
Bischofのアルゴリズム（再掲）
◆ 二段階の二重対角化
(a-1) 下三角帯行列化：
• U11, V11：直交行列， C：半帯幅がLの下三角帯行列
• 大部分を行列乗算により計算可能
• 演算量： 8n3/ 3 （ただし n ≫ L ）
A
C
C
B
(a-2) 村田法：
• U12, V12：直交行列， B：下二重対角行列
• 演算量： 8n2L
• サイズの小さい level-2 BLAS のコールが多数発生
(b) 二重対角行列の特異値分解：
(c-1) 村田法の逆変換：
• 行列乗算により実行可能
• 演算量： 4n3
(c-2) 帯行列化の逆変換：
• 行列乗算により実行可能
• 演算量： 4n3
（26/48）
下三角帯行列化の実装法
CPU
GPU
＝
＝
行列乗算
ブロック鏡像変換の作成
（BLASのみでは実行不可）
（27/48）
GPUでの村田法の実装法
◆ 実装方針
• サイズの小さい BLAS2 のコールが多数発生
CUBLASでなくCUDAで実装
◆ 消去演算のパイプライン式並列化
• 第1列の二重対角化
・・・
• 第2列の二重対角化
更新範囲が重ならない
・・・
◆ GPU上での並列計算方法 （GeForce8800 GTX）
• 16個のMPで消去演算をパイプライン式に並列実行
• MP内の8個のSPで鏡像変換の作用を並列計算（MP内の共有メモリを利用）
（28/48）
特異値分解計算の全体の様子
CPU
GPU
A
A
ブロック鏡像変換の作成
ブロック鏡像変換の作用
C
C
CUBLAS
鏡像変換の作成・作用
B
H
B
H
特異値分解
U2
V2
U2
鏡像変換の合成
Q
V2
合成結果の作用
Q
U’
V’
ブロック鏡像変換の作用
U
V
U
CUBLAS
CUBLAS
V
（29/48）
性能評価
◆ 評価問題
• [-0.5 , 0.5]の乱数を要素とするn×nの正方行列の特異値分解
• n = 1280, 2560, 3840, 5120
◆ 手法
• 二重対角行列の特異値分解は全てI-SVDの倍精度計算で行う
• 二重対角化と逆変換の計算法は以下の通り（全て単精度で行う）
従来法（LAPACKのルーチン）をCPUのみで実行
Bischofの手法（自作プログラム）をCPUのみで実行
Bischofの手法（自作プログラム）をCPUとGPUで実行
※ Bischofの手法における半帯幅Lは32, 64, 128
◆ 実行環境
• CPU ： Intel Core2 Duo (1.86GHz), Intel MKL ver. 8.1, gcc オプション -O3
• GPU ： NVIDIA GeForce8800 GTX, CUBLAS ver. 1.0, nvcc ver. 1.0
（30/48）
二重対角化と逆変換の評価
◆ n=1280
実行時間 (sec)
10.00
逆変換（帯行列化）
逆変換（村田法）
逆変換
村田法
9.00
8.00
7.00
6.00
帯行列化
二重対角化
5.00
4.00
3.00
2.9倍の高速化
2.00
1.00
0.00
従来法
Bischof(32) Bischof(64) Bischof(128) Bischof(32) Bischof(64) Bischof(128)
CPU（2コア）
CPU（1コア） & GPU
（31/48）
特異分解計算全体の評価
◆ n=1280
実行時間 (sec)
10.00
逆変換（帯行列化）
逆変換（村田法）
逆変換
I-SVD
村田法
帯行列化
二重対角化
9.00
8.00
7.00
6.00
5.00
1.9倍の高速化
4.00
3.00
2.00
1.00
0.00
従来法
Bischof(32) Bischof(64) Bischof(128) Bischof(32) Bischof(64) Bischof(128)
CPU（2コア）
CPU（1コア） & GPU
（32/48）
二重対角化と逆変換の評価
◆ n=5120
実行時間 (sec)
400
逆変換（帯行列化）
逆変換（村田法）
逆変換
村田法
350
300
250
帯行列化
二重対角化
200
150
100
7.0倍の高速化
50
0
従来法
Bischof(32)
Bischof(64) Bischof(128) Bischof(32)
CPU（2コア）
Bischof(64) Bischof(128)
CPU（1コア） & GPU
（33/48）
特異分解計算全体の評価
◆ n=5120
実行時間 (sec)
400
逆変換（帯行列化）
逆変換（村田法）
逆変換
I-SVD
村田法
帯行列化
二重対角化
350
300
250
200
150
4.2倍の高速化
100
50
0
従来法
Bischof(32)
Bischof(64) Bischof(128) Bischof(32)
CPU（2コア）
Bischof(64) Bischof(128)
CPU（1コア） & GPU
（34/48）
Bischofの手法のステップごとの評価
◆ L=64の場合
Speedup
25.0
帯行列化
村田法
逆変換（村田法）
逆変換（帯行列化）
20.0
15.0
10.0
5.0
0.0
1280
2560
3840
5120
(n)
Speedup = CPU (2コア)での実行時間 / CPU (1コア) & GPUでの実行時間
（35/48）
精度評価
1.0E+00
CPU_LAPACK (double)
CPU_LAPACK (single)
CPU_Bischof (single)
CPU & GPU_Bischof (single)
1.0E-01
1.0E-02
1.0E-03
1.0E+00
1.0E-01
1.0E-02
1.0E-03
1.0E-04
1.0E-04
1.0E-05
1.0E-05
1.0E-06
1.0E-06
1.0E-07
1.0E-07
1.0E-08
1.0E-08
1.0E-09
1.0E-09
1.0E-10
1.0E-10
1.0E-11
1.0E-11
1.0E-12
1280
3840
2560
n
5120
1280
3840
2560
5120
n
（36/48）
まとめと今後の課題
• 本研究のまとめ
– CUBLAS/CUDAを用いて正方行列の特異値分解をGPU上で実装
– Bischof のアルゴリズムに基づき，演算量の大部分を行列乗算の形
で実行することで，高速な CUBLAS の SGEMM ルーチンを活用
– 行列を GPU メモリに置くことで，データ転送量を削減
– 特異値分解全体でデュアルコアCPUに比べ4倍程度の高速化を実現
• 今後の課題
– CPU と GPU の仕事の分担のさらなる効率化
– 新たな性能ネック部分の高速化
• 村田法（帯行列から二重対角行列への変換）
• 二重対角行列の特異値分解
– 最新の GPU や他のアクセラレータを用いた性能評価
（37/48）
メニーコアプロセッサでの
固有値分解
（38/48）
非対称行列の固有値計算
• ヘッセンベルグQR法
– 非対称行列を直交変換によりヘッセンベルグ行列に変換
– ヘッセンベルグ行列の固有値を QR 法により求める
高速化対象
対角成分 = 固有値
有限回演算
密行列
0
反復計算
ヘッセンベルグ行列
ハウスホルダー法
演算量 : （10/3）n3
0
右上三角行列
QR 法
演算量 : 10 n3 （経験値）
– 演算量の大部分を占める QR 法の高速化が必要
（39/48）
QR 法の高速化
• 在来研究
– マルチシフト QR 法 （K. Braman et al., 2002）
• 複数のシフトを導入
• 演算の大部分が行列乗算（Level-3 BLAS が利用可能）
– キャッシュの有効利用
– 並列性の向上
• 本研究
– 演算の大部分を占める行列乗算を高速化
• 浮動小数点アクセラレータの利用
– ClearSpeed社 CSX600
– 1 チップに96個の演算コアを集積
– 最適化された行列乗算ライブラリ
（40/48）
マルチシフトQR法をそのまま使った場合
12000
その他
行列乗算
実行時間（秒）
10000
8000
6000
• 乱数行列 （n = 6000）
• マルチシフト QR 法により
すべての固有値を計算
• 計算機環境
4000
– Xeon 3.2 GHz, Memory 8 GB
– ClearSpeed advance board
• CSX600 × 2
2000
0
CSX600 （無） / （有）
シフト数（ブロックサイズ）増加
• 問題点
– 実行時間の大部分を占める行列乗算（O(n3)）は，CSX600 により高速化
– ただし，CSX600 の性能を引き出すには，シフト数を大きく取る必要あり
• Byte/Flop が 0.011（通常のCPUの1/10以下）と極端に小さい
• そのため，性能を出すには行列乗算のサイズが最低 1500 程度必要
– その場合，行列乗算以外の部分（O(n2)）の演算量が増加してしまう
（41/48）
マルチシフトQR法の演算パターン
• バルジ追跡型の演算
– 直交変換により，（m/2）個のバルジ（出っ張り）
を対角線に沿って k 行移動させる
• 最初に 対角ブロック内でバルジを移動
– それに伴い，対角ブロック内部を更新
– 更新に用いた直交変換を1 個の行列に累積
BLAS 3
• 最後に 非対角ブロックを まとめて更新
– 累積した行列 と 非対角ブロックの行列乗算
k
0
バルジ（3×3）
逐次的に更新
まとめて更新
（BLAS 3 利用）
• シフト数 m が増えることの効果
– 行列乗算のサイズが増加
– 対角ブロック内部の更新の演算量増大
行列乗算部の性能向上
他の部分の実行時間増加
（42/48）
アルゴリズムの改良（1）
再帰型アルゴリズムへの変換
k/q
k
0
0
逐次的に更新
逐次的に更新
まとめて更新（BLAS 3 利用）
（例：再帰段数 d = 1）
(m / 2) / q 個の bulge を
k / q 行 chasing
まとめて更新（BLAS 3 利用）
Level-3 BLASの利用可能部分の増大
（43/48）
アルゴリズムの改良（2）
• 収束に伴う行列の分離
– 右下
小行列が分離 （
）
• 分離した小行列の固有値計算
– ダブルシフト QR 法を適用
– シフト数 m 増加により，小行列の次元増加
• 更新処理の分割
0
ブロック化
– 対角ブロックの更新 （右上三角化するまで）
• 更新に用いたハウスホルダー変換を
1 個の行列に累積
– 非対角ブロックの更新 （最後に 1 度だけ）
0
逐次的に更新
まとめて更新
（BLAS 3 利用）
• 累積した行列 と 非対角ブロックの行列乗算
・ Level-3 BLASの利用可能部分の増大
・ 演算量削減
（44/48）
数値実験
• テスト問題
– [0, 1] の成分を持つ
乱数行列
• ハウスホルダー法によりヘッセンベルグ化
– すべての固有値と固有ベクトルを計算
• 解法
– 従来のマルチシフトQR法
– 改良型のマルチシフトQR法
• 計算機環境
– Xeon 3.2 GHz , Memory 8 GB
– Fortran 77 倍精度演算
– ClearSpeed advance board
• CSX600 × 2
– Level-3 BLAS（行列乗算）
• CSXL Library (CSX600)
• Intel Math Kernel Library (CPU)
（45/48）
性能評価：アルゴリズムの効果
5000
実行時間（秒）
4000
3000
従来法，提案法ともに CSX600 を利用
その他
小行列の固有値計算
対角ブロックの更新
非対角ブロックの更新
2000
1000
0
従来法
提案法
従来法
提案法
（n = 3000, m = 160, q = 4） （n = 6000, m = 200, q = 5）
• 従来法に比べ，提案法は 約 1.4 倍高速
– 対角ブロックの更新：約 1.5 倍高速化（level-3 BLAS化）
– 小行列の固有値計算：10 倍以上高速化（演算量削減）
（46/48）
性能評価：CSX600 の効果
100000
実行時間（秒）
80000
60000
その他
小行列の固有値計算
対角ブロックの更新
非対角ブロックの更新
40000
n = 12000
20000
n = 6000
0
（無）
（有）
（有）
（無）
（有）
（有）
m = 100
q=5
m = 200
q=5
m = 100
q=5
m = 240
q=6
• CSX600 （+ 提案法） により
– n = 6000 のとき，約 3.5 倍高速化
– n = 12000 のとき，約 3.8 倍高速化
（47/48）
まとめと今後の課題
• 本研究のまとめ
– メニーコア型の浮動小数点アクセラレータ CSX600 を用いて非対称固有値
計算のためのマルチシフトQR法を実装
– アクセラレータの性能を引き出すには，ブロックサイズを大きく取る必要あり
– しかしその結果，level-3 BLAS 以外の(O(n2)の部分の)計算時間が増大
– 再帰型アルゴリズムへの変換などの手法を用いることで，O(n2)の部分につ
いても level-3 BLAS 化を行い，全体を高速化
• 今後の課題
– CPU とアクセラレータの仕事の分担のさらなる効率化
– 他のアクセラレータを用いた性能評価
（48/48）
© Copyright 2024 ExpyDoc
