3476

チュートリアル

科学技術計算の分野で始まった GPU（Graphic Processing Unit）による汎用計算は、製造業、医療、金融など様々な分野へ広まり、モバイル、組み込みの分野でも使われ始めています。
当初は専用言語を使う必要があり GPU による汎用計算は気軽に始められるものではありませんでしたが、用途が広がると同時に、従来の言語で GPU 向けにプログラミングする環境が整ってきました。OpenACC はその中の一つで、既存のプログラムにディレクティブと呼ばれるヒント情報を追加するだけで、コンパイラが自動で GPU コードを生成する開発手法です。今回のチュートリアルでは、全４回のシリーズとして、この OpenACC を使ってアプ
リケーションを効率良く GPU で加速する方法を解説します。

OpenACC で始める GPU コンピューティング :ループの並列化

っ た G P U ( G r a p h i c P r o c e s s i n g U n i t ) に よ る 汎 用 計 算 は 、製 造 業 、医 療 、り 、モ バ イ ル 、組 み 込 み の 分 野 で も 使 わ れ 始 め て い ま す 。当 初 は 専

U による汎用計算は気軽に始められるものではありませんでした、 従 来 の 言 語 で GPU 向 け に プ ロ グ ラ ミ ン グ す る 環 境 が 整 っ て き ま成瀬 彰

一つで、既存のプログラムにディレクティブと呼ばれるヒント情

パ イ ラ が 自 動 で GPU コ ー ド を 生 成 す る 開 発 手 法 で す 。 今 回 の チ ュ

シ リ ー ズ と し て 、 こ の OpenACC を 使 っ て ア プ リ ケ ー シ ョ ン を 効 率

解説します。

1

はじめに

前回は、OpenACC を使ってどのようにアプリケー
ションを並列化するのがよいのか、我々が推奨する方
U コンピューティング:
OpenACC 概要
法を紹介した。我々の推奨方法は、アプリケーション
computing in OpenACC: Summary of OpenACC
を一挙に並列化するのではなく、処理時間の長い部分
成瀬 （図1）
彰から一つずつ順番に並列化する方法である
。

は、多様化
パコンラン
プロセッサ
。その後、
、Intel はメ
米 CORAL
ロセッサがセッサも使して FPGA
ノードに、珍しくない。
様なプロセ
は朗報であ

の確保が問
様化してお
セッサ以外
もしくは、
ねない。後
だが、多様
性能ポータ
となってい
ロセッサ間
。
「並列性」
とである。
MD 幅を広ッサのコアになり[3]、て AVX-512
、粒度の異
方針を採っ
度のブロッ

ク並列性、細粒度な SIMT 並列性[5]、やはり粒度の異なる 2 つの並列性を高めることで高性能を実現するアーキテクチャである。近年の、そして、将来のプロセッサを有効に使うには、単にアプリケーションから大量の並列性を抽出するだけでなく、異なる粒度に合った並列性を抽出する必要があると言える。
もう一つの共通点は、メモリシステム階層の深化である。CPU のメモリ階層は、主記憶（典型的には DDRメモリ)、同一チップ内の CPU コア間でシェアされる
共有キャッシュ、CPU コア内のプライベートキャッシ
ュで構成される。一方、GPU のメモリ階層は、CPU の
主記憶、GPU の主記憶（DDR メモリよりスループッ
トの高い GDDR メモリが典型的）、複数階層のキャッシュ、スクラッチパッドメモリで構成される。
図1 OpenACC によるアプリケーション並列化のバグの修正が、いかに困難なものかを考慮すると、並列化に要する時間をトータルでは最も短くできる方法
であると我々は確信している。実際、コード行数が十万行を越える大規模アプリであっても、この方法に基づいて、アプリケーション全体の GPU 向け並列化を達成した事例も存在する。
今回は、このアプリケーション並列化サイクルの、
1）ホットスポットの特定と、2）ループの並列化を説明
してゆく。

2

ホットスポットの特定

コードの並列化を始める前に、プログラムのどの部分の処理に時間を要しているかを把握することは重要である。アプリケーション全体の実行時間に対して、相当の割合を占めるサブルーチンやループネストは
ホットスポットと呼ばれるが、並列化はこのホットス
ポットから始めるのが良い。

2.1

サイクル

続きは Web で
日本計算工学会誌「計算工学(Vol.21, No.2)」HP:
この方法は、一見すると、回りくどい方法に見える
http://www.jsces.org/Issue/Journal/

プロファイリング・ツール

GNU gprof[1]、Intel VTune[2]、Vampir[3]、TAU[4] 等 々、アプリケーション情報を採取・収集する、いわゆるプロファイリングツールは世の中に多数存在する。各
ツールにはそれぞれ長所・短所があり、どのツールが良いとは一概には言えないが、より詳細なプロファイ
リングデータがあるほど、並列化の作業を効率良く進められる。

かもしれない。しかし、並列化作業の過程で混入した

筆者紹介

筆者紹介

なるせ
なるせあきらあきらAkira NARUSE
1996
年、名古屋大学大学院工学系研究科修士課
1996年、名古屋大学大学院工学系研究科修士課程修了（情報工学専攻）。同年、富士通研究所に入社、程修了(情報工学専攻)。同年、富士通研究所に入社、大規模サーバ開発、HPC システム開発な

大規模サーバ開発、HPC システム開発など、様々

ど、様々なプロジェクトに参加、関連するハーなプロジェクトに参加、関連するハードウェアやドウェアやソフトウェアの研究開発に従事。

ソフトウェアの研究開発に従事。2009年、山下記

続きは Web で
日本計算工学会誌「計算工学（Vol.21, No.3）」HP:
http://www.jsces.org/Issue/Journal/

2009 年、山下記念研究賞受賞。2010 年、念研究賞受賞。2010年、SACSIS2010最優秀論文賞

SACSIS2010 最優秀論文賞受賞。2013 年、NVIDIA に参加、シニアデベロッ受 賞。2013年、NVIDIA に 参 加、 シ ニ ア デ ベ ロ ッ

パーテクノロジーエンジニアとして様々なアルゴリズムとアプリケーショ

パーテクノロジーエンジニアとして様々なアルゴ

ンの GPU 向け並列化に従事、GPU コンピューティングの普及に努めている。

リズムとアプリケーションの GPU 向け並列化に従事、GPU コンピューティングの普及に努めている。

（26）

計算工学

チュートリアル

OpenACC で始める GPU コンピューティング : ループの並列化

例えば、コード行数が数百〜数千行のサブルーチンの場合、そのサブルーチンがホットスポットと分かるだ
けではあまり意味がない。具体的に、サブルーチン内のどの部分で時間がかかっているかを特定できないと、効率良く並列化を進められない。また、ホットス
ポットの性能律速要因が何なのかも重要である。もし
ボトルネックがディスク IO やプロセス間通信であれ
ば、GPU のようなアクセラレータに出番はない。演算性能やメモリアクセス性能が性能リミッターであれ
ば、GPU で加速可能だ。しかし、そこに十分なデータ
並列性がなければ、大きな性能向上は望めない。そのホットスポットがループネストだとすると、ループ反復回数は何回なのか、多いのか、少ないのか、このような情報も、並列化の効果が見込めるかどうかの観点
で、ホットスポットの重要度を決める上で有用である。
GPU 向 け 並 列 化 の 過 程 で は、CPU の プ ロ フ ァ イ ル
データだけでなく、GPU のプロファイルデータも採取できるのが望ましい。PGI のプロファイリングツール
pgprof[5] は、CPU 向けプロファイリング機能と GPU 向
けプロファイリング機能の両方をサポートしている。
現時点で OpenACC コンパイラとして普及しているのはPGI コンパイラであり、本チュートリアルでも PGI コン
パイラを使用することもあり、本稿ではプロファイリ
ング事例に pgprof を使用する。
2.2 正しいプロファイルデータ：
プロファイルデータはとても有用だが、正しいプロ
ファイルデータの採取に注意する必要がある。アプリ
ケーションには様々なパラメータやオプションが存在する。これらのちょっとした違いにより、アプリケー
ションは特性が変わることがある。例えば、実行の簡単な設定、例えば、小さな問題サイズの設定でアプリ
ケーションを実行した場合と、本来の用途の問題サイズで実行した場合とで、ホットスポットが全く変わることもある。適切ではないプロファイルに基づいて並列化を進めると、本当のホットスポットを見逃す危険性がある。本来の用途に近いパラメータ設定でアプリ
ケーションを実行したときのプロファイルデータを使うことに注意して欲しい。

プロファイルデータの採取
プロファイルデータの採取には、pgprof コマンドを使用する。オプション -o を付けて実行すると、このオプションで指定したファイルにプロファイルデータが記録される。以下の操作では、ファイル a.perf にデータが記録される。
$ pgprof –o a.prof ./a.out
プロファイルデータの表示
プロファイルデータの表示にも pgprof コマンドを使用する。表示の場合は、オプション -i でプロファイル
データが入ったファイルを指定する。
$ pgprof –i a.prof --cpu-proﬁling-mode ﬂat======== CPU proﬁling result (ﬂat):
Time(%) Time Name
100.00% 20.57s MAIN_
このプログラムには MAIN 以外にサブルーチンがないので、MIAN しか表示されないが、サブルーチンが多数ある場合は、各サブルーチンの実行時間と全実行時間に占める割合が、実行時間の長いものから表示される。なお、上記例のオプション --cpu-proﬁling-mode は
表 示 形 式 を 指 定 す る オ プ シ ョ ン で、bottom-up、topdown、ﬂat から選択できる。
プロファイルデータの表示（命令レベル）
このプログラムにはループネストが2つあるが、サブ
ルーチン単位のプロファイルからでは、どちらのルー
プネストがホットスポットか特定できない。pgprof は、
オプション --cpu-proﬁling-scope で表示する情報の粒度を、function、instruction から選択できる。デフォルトは function だが、このオプションに instruction を選択す
ると、命令レベルのプロファイルが表示される。

2.3 プリファイリング事例
CPU のプロファイルデータの採取例を説明する。サンプルコードは、本チュートリアルの第一回目の最後に示したヤコビ反復法の Fortran コードである。
プロファイリング実行ファイルの作成
PGI コンパイラでプロファイルデータを採取するには、コンパイル・オプションに -Mprof を追加して実行
ファイルを生成する。以下の操作で、プロファイリン
グ用の実行ファイル a.out が生成される。
$ pgf90 –Mprof jacobi.f90

計算工学$ pgprof --cpu-proﬁling-scope instruction --cpu-proﬁlingmode ﬂat –i a.prof======== CPU proﬁling result (ﬂat):
Time(%) Time Name
10.50% 2 . 16s MAIN_ (./jacobi.f90:42 0x2af)
5.49% 1 . 13s MAIN_ (./jacobi.f90:52 0x335)
4.67% 960ms MAIN_ (./jacobi.f90:42 0x24d)
4.38% 900ms MAIN_ (./jacobi.f90:43 0x2c0)
4.18% 860ms MAIN_ (./jacobi.f90:41 0x238)
4.18% 860ms MAIN_ (./jacobi.f90:52 0x329)
4.04% 830ms MAIN_ (./jacobi.f90:42 0x2a2)
3.99% 820ms MAIN_ (./jacobi.f90:42 0x269)
3.84% 790ms MAIN_ (./jacobi.f90:51 0x300)
3.84% 790ms MAIN_ (./jacobi.f90:42 0x275)...
各行の括弧内の「ファイル名 : 行番号」は、各命令に（26-2）

Vol.21, No.3 2016

チュートリアル

OpenACC で始める GPU コンピューティング : ループの並列化

対応するソースコード上の場所を示している。41-43行
は1番目のループネスト、51-52行は2番目のループネス
トに対応しており、上記の情報から1番目のループネス
トがよりホットスポットと分かる。
ここではテキストベースのプロファイリング事例を説明したが、pgprof は GUI 表示もサポートしている。
pgprof には他にも様々な機能があるが、詳細は [5] を参照して頂きたい。

-----------------------------[Fortran]!$acc kernels
do i=1,nx(i) = 0y(i) = ienddo

3

ループの並列化

アプリケーションのホットスポットが特定できた
ら、次はホットスポットの並列化である。アプリケー
ションによっては、ホットスポットが多数あるかもし
れないが、それらを一挙に並列化する必要はない。優先度の高いホットスポットから順番に OpenACC の並列化ディレクティブを挿入、一つずつループを並列化していけばよい。
GPU コンピューティングの識者は、CPU・GPU 間のデータ転送はどうなるのか、不安に思うかもしれな
い。ループを並列するディレクティブを入れただけでも、 コ ン パ イ ラ は ソ ー ス コ ー ド を 解 析 し て 必 要 な
CPU・GPU 間のデータ転送処理を自動で追加するので、この段階ではデータ移動を考慮する必要はない。
多くの場合、コンパイラが決定したデータ転送の方法は最適ではなく、CPU・GPU 間のデータ転送が必要以上に行われる。そのため、この段階ではアプリケー
ション性能が低下することも珍しくない。OpenACC にはデータ転送を制御するディレクティブがあるが、それを用いたデータ転送の最適化は、次のステップで取り組む。このステップではループの並列化に集中する。
OpenACC は並列化領域を指定する方法として2つのディレクティブを用意している。Kernels コンストラク
トと、Parallel コンストラクトである。以下、それぞれの特徴や使い方を説明する。
3.1 Kernels コンストラクト
Kernels コンストラクトは、ソースコード内の並列化
可能領域を指定するディレクティブの一つである。主
に、指定領域内の処理を具体的にどう並列化するかを、コンパイラに決定させるときに使用する。指定領
域内には、複数の処理ブロック、例えば、複数のルー
プネストが含まれていても良い。コンパイラはループ
ネスト毎に GPU カーネルを生成するかもしれないし、複数のループネストを一つにまとめて GPU カーネルを生成するかもしれない。また、各ループネストを並列化しても安全（＝正しい計算結果が得られる）か否かの判定はコンパイラに委ねられている。コンパイラが安
全でないと判断して、並列化が行われないこともある。Kernels コンストラクトは、「コンパイラおまかせ
方 式」と も 言 え る。 以 下、Fortran と C そ れ ぞ れ の、
Kernels コンストラクトの使用例を示す。

計算工学

do i=1,ny(i)=y(i)+2.0*x(i)enddo!$acc end kernels-----------------------------[C/C++]#pragma acc kernels{for (i=0;i<n;i++) {x[i] = 0.0;y[i] = i;}for (i=0;i<n;i++) {y[i] += 2.0*x[i];}}-----------------------------このコードサンプルは2つのループで構成されており、1番目のループで配列値を設定し、2番目のループ
ではその配列上で計算を行っている。どちらのループ
も、イタレーション間にデータ依存性が無く、並列化可能なので、コンパイラは2つのループのそれぞれにGPU カーネルを生成するだろう。Kernels コンストラク
トでは、内在する並列性を具体的にどう GPU リソースにマップするかはコンパイラに託されている。アクセ
ラレータの種類に依存しない、抽象度の高い記述が可能な方法であり、ソースコードのポータビリティを重視するときに有用である。
3.2 Parallel コンストラクト
Parallel コンストラクトも、並列化可能領域を指定するディレクティブである。Parallel コンストラクトは、後述する Loop コンストラクトと組み合わせた parallel
loop ディレクティブとして使われることが多い。以下、Fortran と C それぞれの、Parallel コンストラクトの使用例を示す。
-----------------------------[Fortran]!$acc parallel loop
do i=1,nx(i) = 0y(i) = i
enddo!$acc end parallel（26-3）

Vol.21, No.3 2016

チュートリアル

OpenACC で始める GPU コンピューティング : ループの並列化!$acc parallel loop
do i=1,ny(i)=y(i)+2.0*x(i)enddo!$acc end parallel-----------------------------[C/C++]#pragma acc parallel loopfor (i=0;i<n;i++) {x[i] = 0.0;y[i] = i;}#pragma acc parallel loopfor (i=0;i<n;i++) {y[i] += 2.0*x[i];}-----------------------------Parallel コンストラクトでは、Kernels コンストラクトとは違い、各処理ブロック（≒ループ or ループネスト）
に対してディレクティブを挿入する必要がある。上記の例では、2つのループの各々に Parallel loop ディレク
ティブが挿入されいている。1つの Parallel コンストラ
クトから、1つの GPU カーネルが生成されると捉えると分かりやすいだろう。これは、コード開発者が明示的に並列化方法を指示する場合に使うディレクティブとして、Parallel コンストラクトが設計されているからである。従って、GPU のタイプや種類に特化した並列化方法や最適化方法を指示する場合には、Parallel コン
ストラクトが適している。

ラマーの責任となる。従って、並列化できないループに対して誤って Parallel ディレクティブを適用すると、計算結果が不正となることもある。
Kernels ディレクティブと Parallel ディレクティブのどちらを使えば良いのか。それぞれ長所・短所があるので、その選択は開発者次第だが、プログラム並列化の経験の少ないプログラマーには Kernels ディレクティ
ブ、並列処理の経験豊富なプログラマーには Parallel
ディレクティブの使用が推奨される。これらディレク
ティブは併用可能なので、例えば、どう並列したらよ
いか分からない、そもそも並列化してよいかどうか分からない処理ブロックに対しては Kernels ディレクティ
ブを使用し、並列化の具体的なイメージのある部分には Parallel ディレクティブを使用するといった様に、同一プログラム内でも場所によって使い分けるても良い。

3.3 Kernels と Parallel の違い
並列領域を指定するディレクティブとして、何故、
Kernels と Parallel の2種類のディレクティブが OpenACCには存在するのか、戸惑う人は多いだろう。実際、これらはとても良く似ているが、微妙な違いが存在する。Kernels ディレクティブが並列化できそうなところ
をコンパイラに知らせる「ヒント」だとすると、Parallel
ディレクティブは並列化できるところをコンパイラに指示する一種の「命令」である。
Kernels ディレクティブの場合、コンパイラにはある
種の自由度が与えられる。あるループを並列化するか否か、ターゲットアクセラレータに対してどのように並列性をマッピングするのかは、コンパイラの責任で行われる。コンパイラが安全と判断すれば、ループは並列化されるが、コンパイラが安全と判断するに足る十分な情報がないと、ループは並列化されない。従って、Kernels ディレクティブでは、コンパイラを変えると、挙動や性能が大きく変わることがある。
Parallel ディレクティブはプログラマーからコンパイ
ラへの指示であり、コンパイラはその命令に従って、ループを並列化する。プログラマーからの指示が十分でない部分はコンパイラが補うが、原則としては、あるループを並列化しても安全か否かの判断は、プログ
計算工学

C/C++ での注意点：
C/C++ コ ー ド で は、 配 列 ポ イ ン タ ー が、 関 数 パ ラ
メータとして渡される場合、Parallel ディレクティブでは並列化できるのに、Kernels ディレクティブでは並列化できないといった現象がしばしば発生する。これは、ポインタ・エイリアス [6] が原因であることが多い。ポインタ・エイリアスとは、別々のポインタを経由する
メモリアクセスが、メモリ上の同じ場所をアクセスする問題であるが、ポインタ・エイリアスの可能性があるループの並列化は安全ではない。このようなループ
であっても Parallel ディレクティブを使うと、コンパイ
ラはプログラマーの指示に従ってそのループを並列化するが（Parallel ディレクティブは、コンパイラへの命令）、Kernels ディレクティブを使うと、コンパイラは安全ではないと判断してループを並列化しない。その配列ポインターでエイリアスが発生しないことが確かな場合は、C/C++ では restrict キーワードを使って、それをコンパイラに通知することができる。もし Kernels
ディレクティブでループが並列化されない現象が発生したら、配列ポインターに restrict キーワードを付けることを試すと良いだろう。
Fortran での注意点：
Fortran コードでは、例えば array（:,:,:) = 0 のような配列記述で、do ループを記述することなく、全部もしくは指定の配列要素に対して、同じ処理を適用できる。
この部分に Kernels ディレクティブを使うと、コンパイ
ラは自動で暗黙的に do ループを作成して並列化する
が、Parallel ディレクティブを使う場合には並列化が行われない。このような配列記述に Parallel ディレクティ
ブを適用するには、プログラマーが明示的にプログラ
ムを書き換えて do ループを入れる必要がある。

3.4

Loop コンストラクト

Loop コンストラクトは、それが挿入された次の行のループに関するヒント・指示を、コンパイラに渡すときに使用する。Loop コンストラクトは、例えば parallel（26-4）

Vol.21, No.3 2016

チュートリアル

OpenACC で始める GPU コンピューティング : ループの並列化

loop ディレクティブのように、parallel コンストラクトと組み合わせて使用できるが、kernels コンストラクトで指定した並列領域内でも使うことができる。
Loop ディレクティブを使うと、様々なヒントや指示をコンパイラに渡すことができる。その指示の種類は、性能を最適化するため指示と、正しい結果を得るための指示に、大別される。以下では、正しい結果を得る
ための Loop ディレクティブの使い方を一部説明する。

do i=1,nx(i) = ...
sum = sum+x(i)enddo!$acc end kernels-----------------------------[C/C++]#pragma acc parallel loop reduction(+:sum)for (i=0;i<n;i++) {x[i] = ...;
sum += x[i];}-----------------------------上の例は総和を求めているので、縮約操作として +を指定している。縮約計算の内容が、総乗、最大値、最小値の場合は、それぞれ、*、max、min を縮約操作として指定すればよい。なお、縮約対象として指定できるのは変数のみである。Reduction 節では、配列を縮約対象に指定できない。

3.4.1 Private 節
Private 節は、ループイタレーション毎に独立して使用すべき変数・配列（private 変数・配列）をコンパイラに伝えるのに使われる。例えば、ループ内で、配列tmp が計算過程の値を保持する一時配列として使われていたとする。もし、この配列 tmp をイタレーション間で共用する状態で、全イタレーションを並列実行す
ると、配列 tmp 上でデータ競合が発生し、正しい結果が得られない。多くの場合、private 化すべき変数・配列は、コンパイラが自動で識別するので、常に private
節を使用せねばならない、ということはない。しかし、コンパイラが正しく識別できないケースもある。
そ の 場 合 は、 以 下 の よ う に、Private 節 を 使 っ て、
private な変数・配列を指示する必要がある。
-----------------------------[Fortran]!$acc parallel loop private(tmp)
do i=1,ntmp(1) = ...
tmp(2) = ...
x(i) = x(i)+tmp(1)*tmp(2)enddo!$acc end parallel-----------------------------[C/C++]#pragma acc parallel loop private(tmp)for (i=0;i<n;i++) {tmp[0] = ...;
tmp[1] = ...;
x[i] += tmp[0]*tmp[1];}------------------------------

3.4.2

3.5 ヤコビ反復法の並列化
ここまで説明したディレクティブを使って、ヤコビ
反復法のループを並列化してみよう。ヤコビ反復法のサンプルコードのメインループには、2つの処理ブロック（ループネスト）が含まれている。1番目のループネ
ストは、配列 A に記録されている現時点の配列値から、各配列要素の次の値を計算し、それを配列 Anew に記録する処理である。2番目のループネストでは、配列
Anew の値が配列 A にコピーバックされる。どちらのループネストも、基本的には、全ての配列要素を独立して並列に処理可能だが、1番目のループネストでは、各配列要素の現在の値と次の値の差分の絶対値の最大値の計算（変数 error に関わる計算 )、いわゆる縮約計算が行われていることに注意が必要である。
3.5.1 Parallel コンストラクト使用
Parallel コンストラクトで並列化する場合、ディレク
ティブの入れ方は以下の通りとなる。

Reduction 節

Reduction 節は、各イタレーションの計算結果をまとめて一つの値を計算する処理、いわゆる縮約計算に使われる。縮約計算とは、例えば、総和、総乗、最大値、最小値を計算する処理である。例えば、イタレー
ション毎に計算した結果の合計（sum）を求める処理では、以下の様に Reduction 節を追加することで、ループを並列化できる。
-----------------------------[Fortran]!$acc parallel loop reduction(+:sum)
計算工学（26-5）

Vol.21, No.3 2016

チュートリアル

OpenACC で始める GPU コンピューティング : ループの並列化-----------------------------[Fortran]
do while ( error > tol .and. iter < iter_max )error = 0.0!$acc parallel loop reduction(max:error)
do j=2, N-1!$acc loop reduction(max:error)
do i=2, M-1
Anew(i,j) = 0.25 * &(A(i-1,j) + A(i+1,j) + A(i,j-1) + A(i,j+1))
error = max(error, abs(Anew(i,j) – A(i,j) ))end doend do!$acc end parallel!$acc parallel loop
do j=2, N-1!$acc loop
do i=2, M-1
A(i,j) = Anew(i,j)
end doend do!$acc end parallel
if (mod(iter,100)==0) write(*,*) iter, error
iter = iter + 1
end do-----------------------------[C/C++]
while ( error > tol && iter < iter_max ) {error = 0.0;
#pragma acc parallel loop reduction(max:error)for (int j = 1; j < N-1; j++ ) {#pragma acc loop reduction(max:error)for (int i = 1; i < M-1; i++ ) {Anew[j][i] = 0.25 *(A[j][i-1] + A[j][i+1] + A[j-1][i] + A[j+1][i]);error = fmax(error, fabs(Anew[j][i] – A[j][i] ));}}#pragma acc parallel loopfor (int j = 1; j < N-1; j++ ) {#pragma acc loopfor (int i = 1; i < M-1; i++ ) {A[j][i] = Anew[j][i];}}if (iter % 100 == 0) printf("%d, %f\n", iter, error)
iter++;}-----------------------------C/C++ コード、Fortran コード、どちらもループネス
トの外側のループに parallel loop ディレクティブを挿入
し、内側のループの前に loop ディレクティブを挿入する。多重ループの場合、内側のループに対しても loop

ディレクティブの挿入が必要なことに注意されたい。
Parallel コンストラクトでは、loop コンストラクトが付けられたループは、並列化可能であることを意味する。従って、これでループネスト内の全てのループが並列化可能であることをコンパイラに指示したことになる。
1番目のループネストでは、変数 error に対して縮約計算（最大値の算出）が行われているので、reduction 節を使って、これを指示する。なお、コンパイラによっては、reduction 節がなくても、自動で縮約計算処理を認識して、暗黙的に縮約計算を行うコードを生成することがある。コンパイラの自動認識機能に依存してもよいが、コードのポータビリティを重視するなら明示的に reduction 節を入れるのがよいだろう。
OpenACC コードから、PGI コンパイラを使って GPU
用 実 行 フ ァ イ ル を 生 成 す る に は、 オ プ シ ョ ン -acc –ta=tesla を付ける。追加でオプション -Minfo=accel を付
けると、コンパイラのコード解析結果に関する情報が表示される。特に開発過程では、このオプションの使用 を 推 奨 す る。 以 下 に、 こ の サ ン プ ル コ ード を PGI
Fortran コンパイラ（v16.4）でビルドしたときの出力を示す。コンパイラメッセージの詳細は、ここでは説明しないが、このメッセージを見ると、どのループが並列化されたのか、どのように並列化されたのか、縮約計算は
指示通り行われているか、計算に必要なデータは正し
く転送されているか、などを確認することができる。
-----------------------------$ pgf90 -acc -ta=tesla -Minfo=accel jacobi.f90
jacobi:37, Accelerator kernel generated
Generating Tesla code
37, Generating reduction(max:error)
38, !$acc loop gang ! blockidx%x
40, !$acc loop vector(128) ! threadidx%x
Generating reduction(max:error)
37, Generating copyout(anew(2:4095,2:4095))
Generating copyin(a(1:4096,1:4096))
40, Loop is parallelizable
47, Accelerator kernel generated
Generating Tesla code
48, !$acc loop gang ! blockidx%x
50, !$acc loop vector(128) ! threadidx%x
47, Generating copyin(anew(2:4095,2:4095))
Generating copyout(a(2:4095,2:4095))
50, Loop is parallelizable------------------------------

計算工学

3.5.2

Kernels コンストラクト使用

次に、Kernels コンストラクトを使用して OpenACC
化する。Kernels コンストラクトの基本は、コンパイラにおまかせ、である。その考えに従うと、kernels ディ
レクティブの入れ方は以下の通りとなる。
------------------------------（26-6）

Vol.21, No.3 2016

チュートリアル

OpenACC で始める GPU コンピューティング : ループの並列化[Fortran]
do while ( error > tol .and. iter < iter_max )error = 0.0!$acc kernels
do j=2, N-1
do i=2, M-1
Anew(i,j) = 0.25 * &(A(i-1,j) + A(i+1,j) + A(i,j-1) + A(i,j+1))
error = max(error, abs(Anew(i,j) – A(i,j) ))end doend dodo j=2, N-1
do i=2, M-1
A(i,j) = Anew(i,j)
end doend do!$acc end kernels
if (mod(iter,100)==0) write(*,*) iter, error
iter = iter + 1
end do-----------------------------[C/C++]
while ( error > tol && iter < iter_max ) {error = 0.0;
#pragma acc kernels{for (int j = 1; j < N-1; j++ ) {for (int i = 1; i < M-1; i++ ) {Anew[j][i] = 0.25 *(A[j][i-1] + A[j][i+1] + A[j-1][i] + A[j+1][i]);error = fmax(error, fabs(Anew[j][i] – A[j][i] ));}}for (int j = 1; j < N-1; j++ ) {for (int i = 1; i < M-1; i++ ) {A[j][i] = Anew[j][i];}}}if (iter % 100 == 0) printf("%d, %f\n", iter, error)
iter++;}-----------------------------Kernels デ ィ レ ク テ ィ ブ を 使 う と、 上 記 の よ う に、
ディレクティブ挿入量を少なくすることができる。
ソースコードをポータビリティの高い状態にできると言ってよいだろう。一方で、コンパイラがループをどう認識したか、ループを正しく並列化できているか、並列化可能なところで逐次コードが生成されていないか、その確認は開発者の責任となる。以下に kernels
ディレクティブで OpenACC 化した Fortran コードを PGI
コンパイラ（v16.4）でビルドしたときのコンパイラ出力を示す。このメッセージを見ると、全てのループが並列化されており、集約計算も正しく対処されていることが分かる。
-----------------------------$ pgf90 -acc –ta=tesla -Minfo=accel jacobi-k.f90
jacobi:37, Generating copyout(anew(2:4095,2:4095))
Generating copyin(a(1:4096,1:4096))
Generating copyout(a(2:4095,2:4095))
38, Loop is parallelizable
39, Loop is parallelizable
Accelerator kernel generated
Generating Tesla code
38, !$acc loop gang, vector(4) ! blockidx%ythreadidx%y
39, !$acc loop gang, vector(32) ! blockidx%x
threadidx%x
41, Generating implicit reduction(max:error)
45, Loop is parallelizable
46, Loop is parallelizable
Accelerator kernel generated
Generating Tesla code
45, !$acc loop gang, vector(4) ! blockidx%ythreadidx%y
46, !$acc loop gang, vector(32) ! blockidx%x
threadidx%x------------------------------

計算工学

4

おわりに

今回は、我々が推奨する OpenACC によるアプリケー
ション並列化サイクルの前半、1）ホットスポットの特定、2）ループの並列化、この2ステップを具体的な事例を用いて説明した。サンプルコードに用いたヤコビ反復法で時間のかかる処理は、この2ステップを適用する
ことで、既に GPU 上で実行できる状態になっている。
しかし、実はこの状態では、CPU 実行よりも GPU 実行の方が遅い。その理由は、CPU と GPU 間のデータ転送が最適化されておらず、データ転送が必要以上に行われているためである。
次回は、OpenACC によるアプリケーション並列化サ
イクルの後半、データ転送を最適化する方法、そして、
ループをより良く並列化する方法について説明する。
■参考文献
[1]
[2]
[3]
[4]
[5]
[6]

（26-7）

GNU gprof: https://sourceware.org/binutils/docs/gprof/
Intel VTune: https://software.intel.com/en-us/intel-vtuneampliﬁer-xe
Vampir: https://www.vampir.eu/
TAU: https://www.cs.uoregon.edu/research/tau/home.php
PGI Profiler User's Guide: http://www.pgroup.com/doc/
pgprofug.pdf
Pointer aliasing: https://en.wikipedia.org/wiki/Pointer_aliasing

Vol.21, No.3 2016