GPU 統合型 CPU を用いた科学技術計算の高速化

408
GPU 統合型 CPU を 用いた科学技術計算の高速化
Acceleration of Scientific Computing by using GPU-integrated CPU
○正
非
高橋 徹（名古屋大学） 非
延山 龍介（名古屋大学） 正
大澤 尚希（名古屋大学） 学
飯盛 浩司（名古屋大学） 正
浅野 朗（名古屋大学）
松本 敏郎（名古屋大学）
Toru TAKAHASHI, Naoki OSAWA, Hogara ASANO, Ryusuke NOBUYAMA, Hiroshi ISAKARI, Toshiro
MATSUMOTO, Nagoya University, Furo-cho, Chikusa-ku, Nagoya, Aichi
Key word : GPGPU, GPU-integrated CPU, Intel HD Graphics, FDTD, BEM, FEM
2)
研究背景と 目的
旧来はグラフィックス処理に特化した演算装置で あった
GPU（Graphical Processing Unit）は，その高い並列演算処
理能力が注目されて，汎用的な並列演算処理に応用すること
が 近 年 盛 ん に 行 わ れ て い る ．こ の 種 の 研 究 開 発 を GPGPU
（General Purpose computing on GPU）と 言 う ．GPGPU は
2000 年頃の試行錯誤的な段階 を経 て，Nvidia 社が販売した
ハードウェア（GeForce 8 シリーズ）と開発環境（CUDA）の
販売を契機として産学官に急速に普及した．今日では多種
多様なアプリケーション が GPU によ っ て加 速さ れて い る．
GPU の 主 流 な 形 態 は 外 部 装 置 と し て PC に 装 着 す る も の
で あ る（ 以 下 ，
「 外 付 け GPU」と 呼 ぶ ）．一 方 ，GPU の 別 の
形態として，CPU と同じダイ上に備え付けられた GPU（以
下，
「 内 蔵 型 GPU」と 呼 ぶ ）が 2010 年 頃 に 登 場 し た ．代 表
的な製品は，AMD 社の Accelerated Processing Unit（APU）
と Intel 社 の Intel HD Graphics で あ る ．こ れ ら 内 蔵 型 GPU
は ，CPU と メ モ リ を 物 理 的 に 共 有 し て い る た め ，外 付 け
GPU で は パ フ ォ ー マ ン ス 上 の ネ ッ ク と な り 得 る デ ー タ 転
送 が 不 要（ 論 理 的 に は 必 要 ）で あ る 点 が 特 徴 で あ る ．半
面 ，外 付 け GPU と 比 較 す る と ，内 蔵 GPU が 搭 載 し て い る
演算コア数は少ないために演算性能が小さい．また，GPU
が 単 独 で 利 用 可 能 な メ モ リ 容 量 が 小 さ い ．図 1 は ，Intel 社
の GPU 統 合 型 CPU で あ る Core i7 シ リ ー ズ に つ い て ，そ の
CPU 部 と GPU 部（Intel HD Graphics）の 演 算 性 能（ 単 精 度
浮 動 小 数 点 演 算 に 関 す る FLOPS 値 ）の 推 移 を 示 す 1) ．現
状 ，Intel HD Graphcs の 演 算 性 能 は ，外 付 け GPU で あ る
Nvidia 社 の GeForce あ る い は Tesla シ リ ー ズ の 10 分 の 1 程 度
で あ る も の の ，CPU 部 と ほ ぼ 同 等 の 性 能 を 有 す る GPU 部
をモニタ出力以外にも有効に活用できれば有意義である．
と こ ろ が ，内 蔵 型 GPU を 用 い て GPGPU を 実 践 す る に あ
たって，プログラミングのノウハウが十 分蓄 積し てい ると
は 言 い 難 い ．こ の よ う な 状 況 の 中 で ，著 者 ら は 試 行 錯 誤 の
結 果 ，OS を Windows7 と す る Core i5 を 搭 載 し た PC を 用 い
て ，OpenMP と OpenCL を 併 用 す る こ と に よ り ，Core i5 の
CPU 部 と GPU 部 を 併 用（ オ ー バ ラ ッ プ ）で き る こ と を 確
認 し た ．リ ス ト 1 は そ の 骨 組 み と な る 擬 似 コ ー ド で あ る ．
List 1 Parallelisation using OpenMP and OpenCL.
および有限要素法（FEM）3) の簡単な場合に対して適用
し ，そ の 有 用 性 を 議 論 し て い る 。本 研 究 の 目 的 は ，当 該 の
併 用 並 列 計 算 手 法 を 3 次 元 時 間 領 域 差 分 法（FDTD 法 ）に
適用することである．
104
K80
GTX TITAN Black
GTX TITAN
GTX680
M2090
103
GFLOPS
1
GTX280
C1060
GTX480
GTX285
8800GTX
8800GTX
Ultra
101
2006
L3406
W5590
Core i7
960
2007
Core2 Duo
E8600
2008
E7-4890 v2
E5-2678W
8800GTS
Core2 Duo
E6700K
HD Graphics
4000
9800GTX
102
HD Graphics
5200
GTX580
E5-2699 v3
Iris Pro Graphics
6200
2009
Core i7
4770K
E7-2870
Core i7 Core i7
2700K 3770K
HD Graphics
3000
HD Graphics
2010
2011
2012
2013
E5-2697 v2
Core i7
6700K
Intel CPU
Intel Xeon
Nvidia Geforce
Nvidia Tesla
Intel GPU
2014
2015
2016
Year
Fig. 1 Transition of FLOPS for single precision.
FDTD 法の概略と 並列化
2
2.1
概略
本 研 究 で 扱 う FDTD 法 は 標 準 的 な 定 式 化 に 従 う 4) ．
す な わ ち ，Maxwell 方 程 式 を Yee グ リ ッ ド（ 総 セ ル 数 は
Nx × Ny × Nz ）による空間分割と Leap flog 法による時間積
分を用いて離散化する．開領域における電磁場散乱問題を
解 析 対 象 と し ，Mur の 吸 収 境 界 条 件 を 用 い る ．
FDTD 法 に お け る 主 要 な 計 算 タ ス ク は ，初 期 化 に 続 き ，
(i) 内 部 セ ル の 電 磁 場 の 更 新（ 以 下 ，「 内 部 更 新 」），(ii) 境
界セルのそれ（「境界更新」），(iii) 電場・磁場ベクトルの出
力（「 ス ナ ッ プ シ ョ ッ ト 」）が 挙 げ ら れ る ．こ こ で ，(iii) は
所定の時間ステップ数毎に行う．
リ ス ト 2 は ，電 場 の 内 部 更 新 を す る た め の OpenMP 並 列
CPU コ ー ド で あ る ．こ こ で ，Ex，Ey，Ez は 電 場 を ，Hx，
Hy，Hz は 磁 場 を 格 納 す る float 型 配 列 で あ る ．ま た ，float
型 配 列 C，Cx，Cy，Cz は 諸 パ ラ メ ー タ（ 誘 電 率 ，透 磁 率 ，
電 気 伝 導 率 ，グ リ ッ ド サ イ ズ ∆x，∆y ，∆z ，お よ び ，時 間
ス テ ッ プ 長 さ ∆t）か ら 事 前 に 計 算 す る ．
List 2 Parallel CPU code to update the electric field.
#pragma omp parallel for
#pragma omp parallel {
for (int k = 1; k < Nz − 1; k ++) {
id = omp get thread num() // CPU スレッド番号の取得
for (int j = 1; j < Ny − 1; j ++) {
if (id == 0) {
for (int i = 0; i < Nx; i ++) {
CPU0 が OpenCL により任意の GPU コードを実行
int n = i + Nx ∗ (j + Ny ∗ k);
} else {
Ex[n] = C[n] ∗ Ex[n] + Cy[n] ∗ (Hz[n] − Hz[n − Nx]) −
CPU0 以外が任意の CPU コードを実行
Cz[n] ∗ (Hy[n] − Hy[n − Nx ∗ Ny]);
}}
}}}
日本機械学会東海支部第 65 期総会・講演会講演論文集(’16. 3. 17－18) No.163－1
// Similarly, Ey[n] and Ez[n] are updated here.
著者らは上記の併用並列計算手法を，境界要素法（BEM）
併用並列化
タ ス ク (iii) の ス ナ ッ プ シ ョ ッ ト は ，書 き 出 す べ き 時 間 ス
テ ッ プ の 電 磁 場 を 格 納 し た 配 列 Ex，Ey，Ez，Hx，Hy，Hz
が メ モ リ ま た は HDD 上 に 書 き 出 さ れ る 以 前 に 次 の 時 間 ス
テ ッ プ の 値 に よ っ て 更 新 さ れ て は な ら な い た め ，(iii) を (i)
あ る い は (ii) と オ ー バ ー ラ ッ プ さ せ る こ と は 難 し い ．
一方，タスク (i) と (ii) のそれぞれはセル単位で並列性があ
る．そこで，素朴に，境界更新は CPU 部，内部更新は GPU
部により行う．各 CPU コアの役割は，使用可能な CPU コア
数が３以上あることを前提（実際，4 コア搭載の PC を使用）
として，CPU0 が GPU コードの実行，CPU1 が電場の境界更
新 ，CPU2 が 磁 場 の 境 界 更 新 ，そ の 他 は 未 使 用 ，と 定 め る ．
以 上 の よ う な 役 割 分 担 は ，各 タ ス ク の 必 要 算 術 演 算 回
数 ，必 要 メ モ リ ア ク セ ス 回 数 ，各 プ ロ セ ッ サ の 演 算 性 能 等
を 勘 案 し て い な い の で ，必 ず し も 最 良 と は 言 え な い が ，初
期検討としては有意義であろう．
リ ス ト 3 は ，電 場 の 内 部 更 新 を 行 う た め の GPU カ ー ネ ル
（OpenCL プ ロ グ ラ ム ）で あ る ．二 つ の float 型 配 列 C と D は
事 前 に 計 算 す る ．な お ，CPU コ ー ド と は 異 な り ，四 つ の 配
列 C，Cx，Cy，Cz を カ ー ネ ル に 渡 す 場 合 よ り も 性 能 が 良
い ．ま た ，各 GPU コ ア が 一 度 に 64 byte の デ ー タ を 読 み 込
めるようにワークアイテムを構成している．
List 3 OpenCL kernel to update the electric field.
int i = get global id(0); // 3次元のワークアイテムを使用
int j = get global id(1);
int k = get global id(2);
int n = i + Nx ∗ (j + Ny ∗ k);
float D = D[n]; // D[n]:=1/(2∗Eps[n]−dt∗Sig[n]); Eps は誘電率
float C = C[n]; // C[n]:=2∗Eps[n]−dt∗Sig[n]; Sig は電気伝導率
float Cx = D ∗ ax; // ax:=2∗dt/dx; dt は時間ステップ長
float Cy = D ∗ ay; // ay:=2∗dt/dy; dx, dy, dz はグリッドサイズ
float Cz = D ∗ az; // az:=2∗dt/dz
Ex[n] = C ∗ Ex[n] + Cy ∗ Hz[n] − Cy ∗ Hz[n − Nx] − Cz ∗ Hy[
n] − Cz ∗ Hy[n − Nxy];
// Similarly, Ey[n] and Ez[n] are updated here.
1. セ ル 数 N が 大 き い 場 合 ，’HYBRID’ は’GPU’ と 同 等 で
ある．この理由は，境界更新の計算量が O(N 2 ) である
のに対して，内部更新は O(N 3 ) であるため，N が大き
く な る と 境 界 更 新 が 全 体 の 計 算 時 間 に寄 与し な くな る
た め と 考 え ら れ る ．こ れ よ り ，セ ル 数 が 大 き い 場 合 に
は ，内 部 更 新 に CPU と GPU を 併 用 す る の が 望 ま し い ．
少 な く と も ，使 用 し て い な い CPU3 を 内 部 計 算 の 一 部
に用いるとよいと考えられる。
2. 逆に，N が小さい場合には境界更新の寄 与は小 さく な
いので，本併用手法によって高い性能が得られている．
3. い ず れ の 手 法 の 性 能 も ，N ≈ 50 を ピ ー ク と し て 急
激 に 小 さ く な る ．こ の 理 由 と し て 考 え ら れ る の は ，
CPU 部 と GPU 部 が 共 に 使 用 し て い る LLC（Last Level
Cache；容 量 6MB；読 み 込 み 速 度 32B/cycle，書 き 込 み
速 度 64B/cycle）5) に ，CPU コ ー ド お よ び GPU カ ー ネ
ル で 用 い て い る 六 つ な い し は 八 つ の 配 列（Ex，Ey，
Ez，Hx，Hy，Hz，C，お よ び ，Cx，Cy，Cz，ま た は ，
D）を 格 納 で き る 限 界 が N ∼ 50 で あ る た め ，N > 50
で は 低 速 な メ モ リ（8B/cycle）へ の ア ク セ ス が 頻 発 す
るためと考えられる．
4. 上記 1 で言及した N が大きい場合について詳し く見 る
と，’HYBRID’ は’GPU’ よりもやや劣っている．この理
由 と し て は ，’HYBRID’ が 境 界 更 新 と 内 部 更 新 を 同 時
に 行 う た め に ，CPU 部 と GPU 部 が LLC を 食 い 合 う こ
とが考えられる．
400
HYBRID
GPU
CPU
350
mega cells per sec
2.2
300
250
200
150
100
50
3
3.1
数値実験
計算条件
グ リ ッ ド サ イ ズ は ∆(= ∆x = ∆y = ∆z ) = 6.11 × 10−3 m
と し て ，x，y ，z 方 向 に そ れ ぞ れ N (= Nx = Ny = Nz ) 個
のセルを設けた．誘電率および透磁率は一様に真空の値を
用 い た ．電 気 伝 導 率 は 零 と し た ．入 射 場 と し て ，計 算 領 域
の xy 断 面 の 中 央 を 通 り ， z 軸 に 平 行 な ダ イ ポ ー ル ア ン テ
ナ を 配 置 し た ．ア ン テ ナ の 駆 動 周 波 数 は 2.45 GHz と し た ．
ま た ，∆t = 4.08 × 10−12 s と し て ，時 間 ス テ ッ プ 数 は 1000
と し た ．な お ，ス ナ ッ プ シ ョ ッ ト は 行 わ な か っ た ．
3.3
0
計算環境
本研究では GPU 統合型 CPU である Intel Core i5-4440 を搭
載した PC を使用した．GPU 部は Intel HD Graphics 4600 で
ある．単精度浮動小数点数演算性能は，CPU 部が 422，GPU
部 が 352 GFLOPS で あ る ．な お ，Intel SDK for OpenCL
Application 2013 を OpenCL 開 発 キ ッ ト と し て 使 用 し た ．
3.2
0
結果
図 2 は，提案手法（’HYBRID’）の性能結果を，’GPU’（す
なわち，CPU 部による境界更新（2 スレッド使用）と，GPU
部 に よ る 内 部 更 新 を 逐 次 的 に 行 う 手 法 ），お よ び ，’CPU’
（ す な わ ち ，CPU 部 の み に よ る 境 界 更 新 と 内 部 更 新 を 逐 次
的 に 行 う 手 法；共 に 4 ス レ ッ ド 使 用 ）の 各 結 果 と 比 較 し て
い る ．次 の よ う に 考 察 し た：
50
100
150
N cells
200
250
Fig. 2 Comparison of performance.
4
結論
本 研 究 は ，GPU 統 合 型 CPU を 活 用 す る た め の OpenMP
と OpenCL に 基 づ く 併 用 並 列 計 算 手 法 を 3 次 元 FDTD 法 に
適 用 し ，数 値 実 験 に よ っ て そ の 有 効 性 を 確 認 し た ．
参考文献
1) Wikipedia（https://en.wikipedia.org/wiki）の 項 目
List_of_Nvidia_graphics_processing_units，Intel_
HD_Graphics，お よ び ，Xeon よ り デ ー タ を 採 取 し た ．
2) 大 澤 ほ か, GPU 統 合 型 CPU を 用 い た 二 次 元 ヘ ル ム ホ ル
ツ方程式に対する境界要 素法 の 高速 化, 日本機械学会計
算 力 学 講 演 会 論 文 集, Vol.27, 2014.
3) 延 山 龍 介, GPU 統 合 型 CPU を 用 い た ２ 次 元 有 限 要 素 法
の 並 列 計 算, 卒 業 論 文（ 名 古 屋 大 学 工 学 部 ）, 2015.
4) Allen Taflove. COMPUTATIONAL ELECTRODYNAMICS: The Finite-Difference Time-Domain Method. Artech
House: Boston 2005.
R Processor Graphics
5) The Compute Architecture of Intel
Gen7.5, Version 1.0, 2014.