日本ソフトウェア科学会第 37 回大会 (2020 年度) 講演論文集

ブラックボックステストによる軽量なリダクション
ループ自動並列化
森畑 明昌

佐藤 重幸

リダクションループは多数のデータを少数の結果へと集約する処理を表し、計数・平均値・最大値などの計算に広く用いられている。リダクションループは多くの場合効率的に並列計算できる。しかし、演算子の代数的性質を注意深く確認し計算の依存関係を寸断することが必要となるため、自動並列化は容易ではない。本発表では、リダクション
ループが半環上の線形式をなすかを、ブラックボックステストによって判定し、自動並列化を行う手法を提案する。
半環上の線形式に基づく既存の自動並列化手法は記号実行などの解析に基づいており、プログラムの末節の変更の影
響を受けやすく、大きなプログラムの扱いにも難があった。これに対し、提案手法では、リダクション変数の推測・
ループ分解・線形式の係数の推測を全てブラックボックステストとして行う。これにより、静的解析はほぼ不要となり、また複雑なプログラム自動合成も避けることができる。提案手法のプロトタイプ実装を用いた実験により、既存のリダクションループ自動並列化の論文で扱われていた例のほとんどが提案手法で容易に並列化可能であることが確認できた。

Reduction loops, which summarize huge data into a few values, are widely used for counting, averagecalculation, maximum-value finding, etc. While their eﬃcient parallel evaluation is commonly possible,
their automatic parallelization tends to be diﬃcult because of the necessity of using algebraic properties
for breaking value dependencies. In this presentation, we propose an automatic reduction parallelization
method based on black-box testing for recognizing the loop as a system of linear polynomials over semirings. There exist proposals for using linear polynomials for reduction parallelization. However, they relyon static analyses such as symbolic evaluations, and hence a small syntactic modification tends to aﬀect
their parallelizability. The current proposal instead uses black-box testing for every component, including
reduction variable inference, loop decomposition, and coeﬃcient inference. It is little susceptible to program
detail; moreover, it avoids complex program synthesis. Our proof-of-concept implementation could immediately parallelize nearly all test programs, which was exhaustively collected from the literature on automatic
reduction parallelization.

関係があるにもかかわらず効率的な並列計算が可能

1 はじめにであることがよく知られており、「並列リダクション」

リダクションループとは多数のデータを少数の結と呼ばれ重要な構成要素となっている。しかし、リダ

果へと集約する処理であり、データの個数のカウント

クションループの自動並列化は一般には容易でない。

や総和・平均値・最大値などの計算に広く用いられて並列化のためには計算の依存関係を寸断する必要がいる。特に並列計算の文脈では、計算に逐次的な依存あるが、これには演算子の代数的性質を注意深く確認しつつ計算の依存関係を認識することが求められる。

Lightweight Reduction-loop Parallelization by Blackbox Testing.
Akimasa Morihata, 東 京 大 学 大 学 院 総 合 文 化 研 究 科,
Graduate School of Arts and Sciences, the University of Tokyo.
Shigeyuki Sato, 東京大学大学院情報理工学系研究科,
Graduate School of Information Science and Technology, the University of Tokyo.

そのため、既存のフレームワーク（例えば OpenMP、

MPI、Inted TBB など）では、自明に自動並列化が可能なリダクションループ（例えば総和や最大値などの計算）のみを対象とするか、またはプログラマが明
示的に並列評価に必要な情報を与えることを求めている。

複雑なリダクションループや関連する構造（リダク

しの有無）などの影響を受けやすい。これに対し本研

ションを行う再帰関数など）の自動並列化については

究では、プログラムの詳細を全く無視し、プログラム

多くの既存研究がある [1–5, 8, 9, 11–16]. これらはいの入出力関係だけから、適切な並列プログラムを構

ずれも以下のような方針を採用している。

成する。このアプローチでは、静的解析ほど自由に情• 効率的な並列リダクションが可能なプログラムのパターンを与える。

• 並列化対象のリダクションループがそのパター報を得ることはできないものの、入出力関係だけはどんなプログラムでも均一に得ることができる。よって、入出力関係が十分な情報を含んでいれば、これはンに変形できるかを、静的解析（特に記号実行）

頑健な並列化手法となる。

や自動プログラム合成を用いて確認する。
これら 2 点の要素はトレードオフのような関係にあ例として、配列 array の総和を求める簡単なルー
プを考えよう †1 。

る。並列リダクションパターンの一般性の高いほど多s = 0

くのリダクションループが並列化できる可能性がある

for x in array:

が、そのパターンへの具体的な変換手順の構成が難し

s += x

くなり、複雑な自動プログラム合成などが必要となりこのループの本体は s += x という文である。この文やすい。一方、並列リダクションパターンとして一般

を「s と x を含む環境を更新する」プログラムだと性の低いものを採用すると、表現力は劣るが、そのパみなし、ブラックボックステストを行う。テストする

ターンへの変換に際しては特長を活かし効率的に行

内容は「プログラム終了後の s と x は、プログラムえる可能性がある。

開始前の s と x の線形式として表現できるか？」だ。

本研究では「半環上の線形式」を並列リダクション

テストによって、例えば以下のような結果が得られるが可能なプログラムのパターンとして採用する（詳細だろう。左側が実行前の環境、右側が実行後の環境では 2 節で解説する）。この方針は既存研究 [9,10,13,15]ある。

でも採用されているが、並列リダクションパターンと{s = 0, x = 0}→{s = 0, x = 0}

しての一般性が低いと思われること、また実際にその{s = 3, x = 0}→{s = 3, x = 0}{s = 0, x = 4}→

{s = 4, x = 4}

パターンに変換可能なリダクションループであったとしても、それが明示的に半環の演算子を用いて記述されていないと並列化できないことが問題だとされてきた。本研究はこれらの問題が杞憂であることを示す。
つまり、「半環上の線形式」は（少なくとも自動並列化が可能なリダクションループに対する抽象化としては）一般性が高く、また注意深く自動並列化アルゴ
リズムを構成すれば、プログラマが半環の演算子を使わなくても、またプログラムが複雑になろうとも、頑健に並列化が可能である、ということを示す。
本研究の主要な着想は、「半環上の線形式」を並列
リダクションのパターンとして採用すると、ブラック
ボックステストによって簡単に並列化可能性を推測できる、という点である。標準的なアプローチでは、適切な解析によってプログラム並列化に有益な情報を集め、それを用いることで自動並列化を行う。しかし、一般に静的解析は構文上の差異（例えば関数呼び出{s = 2, x = 4} → {s = 6, x = 4}
以上たった 4 回のテストの結果だけからであっても「プログラム終了後の s は、プログラム開始前の s とx の和、終了後の x は開始前の x と等しい」と推測するのは自然だろう。さらに多く、例えば 100 回ほどのテストを行えば、この推測がほぼ間違いないと確信することができるだろう（3 節では具体的なアルゴリズ
ムを示す）。
当然のことながら、上記テストは完全ではなく、推測は間違っているかも知れない。そもそも、ブラック
ボックステストは、原理的に完全にはなり得ない。しかし、以下の 2 つの理由から、このことは自動並列化そのものの本質的な困難に比べればマイナーな問題であると言える。
†1 本研究を通して Python 風の擬似コードををプログ
ラムの記述に用いる

まず、ほとんどの場合、ループ自動並列化は、プロ化する指針は既存研究で既に与えられていたが、多重グラマの「このループを並列化せよ」という指示に基

ループの並列化についてはほとんど議論がなかった。

づいて行われるものだ。ブラックボックステストによ

これに対し、本研究では多重ループに対しても同様にるアプローチが失敗するのは、典型的には特定の値を容易な並列化が可能であることを論じ、並列化アルゴ

閾値としてプログラムの挙動が本質的に変化するよ

リズムを提案した（4 節）。

うなケースだが、並列化対象として指示されるループ

また、複雑なリダクションループでは配列に一時的は通常そのような物ではないし、またプログラマ自身な結果を格納することがあるが、この場合に配列の添がそのようなループへの適用を避けることも難しく

字アクセスによって計算の依存関係が生じうるため、ない。また、プログラマが閾値を自覚しているなら、並列化が難しくなる。この問題に対し、ブラックボッ

その前後の値を両方システムにテストさせることも

クステストに基づき、ループ内で読み出される・書難しくない。

き込まれる配列の添字を推測する手法を提案する（5加えて、ブラックボックステスト完了後に 1 度だ
け、特定の具体的な線形式とプログラムが等価かどう節）。この手法は線形式係数の推測に基づいており、提案手法の応用ともなっている。

かを既存のプログラム検証手法で確認するのは、コス

さらに、提案手法に対しプロトタイプ実装を与え、トとして十分許容できると思われる。少なくとも、既既存のリダクションループ自動並列化の研究で扱わ存の少なくない自動並列化手法がプログラム検器（例

れていた例に対し網羅的な実験評価を行った（6 節）。

えば SAT や SMT のソルバー）を何度も呼び出すの76 例中 74 例については本実装によってごく短時間でに比べれば、大きな改善であるとも言えよう。

並列化が可能であった。並列化に失敗した 2 例は適切以上の洞察に基づき、実際にプログラム自動並列化

な半環を未実装だったために並列化に失敗したもので手法を構成した。これに際して、本研究では多数の技あった。この結果は、「半環上の線形式」が非常に広術的な貢献を行っている。

い範囲の、特に既存手法で自動並列化が可能であると

まず、半環上の線形式によって並列化を行うには、報告されているもののほとんど全てに近いリダクショ

適切な半環とその線形式の係数を推測することが必

ンループを捉えられることを示している。もちろん、要になる。本研究では、これをブラックボックステス「半環上の線形式」が全ての自動並列化可能なリダク

トにより実現する方法を提案する（3. 3 節）。提案手

ションループを捉えられると主張するつもりはない。

法は数値の加算・乗算のなす半環だけでなく、真偽値

しかし、提案手法が軽量であることもあり、これはリに対する論理和・論理積のなす半環など、多くの半環

ダクション自動並列化の 1 つのベースラインとなりに対して係数推測方法を与えた。

うるだろう。

次に、リダクション並列化に際してはループ分解が効果的である一方、ループ分解は並列プログラム

2 半環とリダクションループ並列化の効率を大きく低下させる要因ともなるため闇雲に

2. 1 半環

行うのは現実的ではなく、適切な方針がよく分かって本研究では半環に基づいてリダクションループの並いなかった。これに対し、本研究ではブラックボック

列化を行う。半環は整数上の加算・乗算の関係を抽象

ステストよってループ分解を行う手法（3. 2 節）と推化したものである。形式的には半環は (S, ⊕, ⊗, 0, 1)

測された半環からループ再融合を行う手法を提案すからなる 5 つ組である。S は計算対象となる値の集る（3. 4 節）。提案手法は、並列化可能性の判定時に合である。⊗ と ⊕ は S に対する 2 項演算であり、そは on-the-fly で可能な限りループ分解を行い、そのれぞれ「加算」「乗算」と呼ばれる。0 と 1 は S の要

後に並列可能性を失わない範囲で可能な限りループ素である。これらは以下の性質を満たすことが要求さ

再融合を行う。

れる。

半環上の線形式によってリダクションループを並列

lm, gm = max(0, x + lm), max(gm, x + lm)

a⊕0=0⊕a=a

0 は ⊕ の単位元

a ⊕ (b ⊕ c) = (a ⊕ b) ⊕ c⊕ の結合性

この並列化を論じるためにいくつかの用語を準備

a⊕b=b⊕a⊕ の可換性する。リダクション変数とは、ループを繰り越す依存a⊗1=1⊗a=a

1 は ⊗ の単位元

関係のある値を保持している変数である。上記例でa ⊗ (b ⊗ c) = (a ⊗ b) ⊗ c⊗ の結合性は lm と gm がリダクション変数にあたる。一方、非a ⊗ (b ⊕ c) = (a ⊗ b) ⊕ (a ⊗ c)

⊗ の左分配性

リダクション変数とは、リダクション変数ではない変(b ⊕ c) ⊗ a = (b ⊗ a) ⊕ (c ⊗ a)

⊗ の右分配性数である。リダクション変数はループの実行に従って

a⊗0=0⊗a=0
0 は ⊗ の吸収元
実用上重要な半環の大部分は ⊗ が可換性をもつ。こ

計算を進めないと値が判明しないため、並列計算の妨れを可換半環と呼ぶ。⊗ が可換性をもつと議論が簡り返し毎にどの値を取るかがそれ以前の繰り返しを

単になることもあり、本稿では「半環」で単に「可換計算せずとも確定するため、並列計算の妨げにはなら半環」のことを指す。ただし、提案手法が本質的に ⊗ない。なお、ループカウンタの類（例えば上記ループの可換性を要求するわけではない。

の x）は、形式上はループを繰り越す依存関係を持つ半環の例としては以下のものが挙げられる。整数の加算と乗算のなす半環 (Z, +, ×, 0, 1)、真偽値のげになる。一方、非リダクション変数は、ループの繰

が、値がループを計算せずとも確定するため、非リダ
クション変数として扱う。

論理積と論理和のなす半環 ({0, 1}, ∨, ∧, 0, 1)、真偽

さて、上記ループの本体の計算は、リダクション変数

値の論理和と論理積のなす半環 ({0, 1}, ∧, ∨, 1, 0)、に対応する変数集合 {lm, gm} をとり、半環 (max, +)整数の二項最大値演算と加算のなす半環を考えることで以下の線形式系とみなすことができ(Z ∪ {−∞}, max, +, −∞, 0)、整 数 の 二 項 最 大る。なお、ここでは ↑ を二項最大値演算子として使う。
lm := 0 ↑ (x + lm) ↑ (−∞ + gm)

値 演 算 と 二 項 最 小 値 演 算 の な す 半 環 (Z ∪{∞, −∞}, min, max, ∞, −∞)、非負整数の二項最大値演算と乗算のなす半環 (Z+ , max, ×, 0, 1)、集合の和集合と積集合のなす半環 (2S , ∪, ∩, ∅, S)。これらいずれの例についても、上記性質は簡単に確認できる。
変数集合 {x1 , x2 , . . . , xn } が与えられたとき、半環(S, ⊕, ⊗, 0, 1) 上の線形式とは、係数 a0 , a1 , . . . , an ∈S で特定される以下の形式の式を指す。
a0 ⊕ (a1 ⊗ x1 ) ⊕ · · · ⊕ (an ⊗ xn )以 降 、誤 解 の な い 範 囲 で 半 環 (S, ⊕, ⊗, 0, 1) を(⊕, ⊗) と略記する。
2. 2 代数的単純化による分割統治計算
リダクションループが半環上の線形式で特定できる場合には並列化が可能であることが知られている [9, 10, 13, 15]。例として以下のリダクションループ
を考える。なお、これは最大部分列和問題 [4, 12] の計算に対応する。

gm = 0lm = 0
for x in array:gm := −∞ ↑ (x + lm) ↑ (0 + gm)
なお、左辺の lm と gm はプログラム中では次のルー
プでの lm と gm に対応する。
これをふまえこのプログラムを並列化する。並列プログラムでは、各プロセスはそれぞれループの連
続した繰り返しの処理を担当する。例として、あるプロセスがループの s 回目から t − 1 回目、つまりarray[s] から array[t − 1] までの処理を担当するとする。このプロセスは、本来ループの s − 1 回目で計算されたリダクション変数の値なしに計算を進めることができない。そこで、並列プログラムは代わりに線形式（正確にはその係数）を計算する。以降の繰
り返しについても、具体的な値ではなく線形式の形
で計算を続ける。このとき、半環の代数的性質によって線形式を単純化できる点が重要である。そのため、繰り返しを経ても式が肥大化することなく、係数を求め続ける形で計算を続けることができる。
例として再び最大部分列和問題を考える。いま、

a[s] の値が 3、a[s + 1] の値が −4 だとすると、s + 1
番目のループまで処理した時の線形式は以下のよう

に求められる。
lm := 0 ↑ (−4 + (0 ↑ (3 + lm)))

算するループと、count の値をふまえ valid を求めるループの 2 つに分解する。

=

0 ↑ (−1 + lm)

count = 0:=(−4 + (0 ↑ (3 + lm))) ↑ ((3 + lm) ↑ gm)

counts = [None] * len(array)= −4 ↑ (3 + lm) ↑ gm
上記の処理を繰り返すと、各プロセスが担当部分のfor i in len(1, len(array)):gm

計算を線形式の形まで集約できる。最後に、ループの計算の進行に従い初期値から順に、各プロセスが求めた線形式のリダクション変数へ値を代入し左辺値を求めることを繰り返すことで、最終的なリダクションの結果を求めることができる。
この並列計算の計算量は、プロセッサ数 p がループの回数 N に比べて非常に小さいとして、O(N/p) であり、漸近的に理想的な並列速度向上を示す。

2. 3 リダクション変数の依存関係によるループ
分解
前節で述べたアプローチだけでは様々な演算子が関係したループを並列化できない。例として括弧対
応 [1, 8] を考える。このプログラムは、正しく括弧が対応しているかどうかを判定する。

count = 0
valid = true
for i in len(1, len(array)):
if array[i] == ’(’:
count += 1
elif array[i] == ’)’:
valid = valid and count > 0
count -= 1
変数 count は現状の括弧のネストの深さを、変数

valid はそこまでの括弧対応が正しいかどうかを保持する。閉じ括弧があるときに count が 0 以下の場
合、つまり対応する開き括弧がない場合に valid は偽となる。このプログラムでは、整数に対する加算と真偽値に対する論理積演算が使われており、1 つの半環での計算には明らかに対応しない。
実は、変数 valid の値は変数 count の値次第で変化する（この状況を「valid が count に依存する」と呼ぶことにする）ものの、変数 count は変数 validには依存しない。そこで、count の値を事前に全て計

counts[i] = countif array[i] == ’(’:
count = count + 1
elif array[i] ==’)’:
count = count - 1
valid = true
for i in len(1, len(array)):
if array[i] == ’)’:
valid = valid and counts[i] > 0
前半のループでは、ループの各繰り返しごとの countの値を計算し、配列に記憶する。後半のループはその結果を用いて valid を求める。前半のループは整数の加算のみを、後半のループは真偽値の論理積のみを用いており、それぞれ半環上の線形式で捉えることができる。また、前半のループは単一の値ではな
く配列を計算するが、この計算パターンは並列 scan
（または並列 prefix-sum）と呼ばれ、並列リダクショ
ンと同様に並列計算できることが知られている。よって、このプログラム全体は効率的に並列計算できる。
以上の過程は以下のように一般化することができる。複数のリダクション変数がある場合には、変数の依存関係に基づきループ分解を行うことができる。その結果、相互に値が依存する変数に関する計算のみを（つまりそれ以外の変数の値は全て確定していると仮定した場合に）半環上の線形式によって捉えることができれば、プログラム全体が並列化できる。このこと
は先行研究、例えば [4, 11–13] などで既に広く利用されている。
ループの分解は並列化を簡単にするが、得られる並列プログラムには無視できないオーバーヘッドを課す。並列 scan は、多数の値を計算・記憶する必要がありアルゴリズムも複雑になるため、並列リダクショ
ンに比べ数倍程度遅い。また、メモリの読み書きの回数が増えるため、メモリ帯域により並列速度向上が律速されやすくなる。そのため、並列化に必須でないな

らループ分解は可能な限り避けたい。例えば、前小節する。

の最大部分列和問題でも、lm は gm に依存しないため

ii. E stmt {xj = vj′ } となる vj′ を求める。

ループ分解は可能であるが、ループ分解をせずとも並iii. vi ̸= wi となる値 wi をランダムに生列化できるため避けている。

成する。

iv. E[xi = wi ] stmt {xj = vj′′ } となる

3 ブラックボックステストによる並列化可能性判定

vj′′ を求める。
v. vj′ ̸= vj′′ ならば依存関係に xi 7→ xj を

3. 1 入力プログラム
提案手法の入力としては以下のようなループを想追加し、最内の繰り返しを終了する。

3. 依存関係の推移閉包をもとめ、最終的な依存関
係とする。

定する。

for a in array:

ここでは、ランダムな値を何度も生成しながら、xiの値のみを変更した場合に xj の値が変化しうるかど

stmt
つまり、ループの本体が stmt であるループである。

うかを確認している。xj の値が変化した場合、xj はstmt の実行は array の値を変更しないと仮定する。

xi の値に依存している。最終的な依存関係は、テス

この条件は、3. 2 節の解析時に同時に確認できる。こトによって求まった推移閉包に対応する。これは以下

れに加え、stmt が使用する変数の集合 X が既知であのようなケースを適切に扱うためである。

るものとする。変数集合を得るには簡単な静的解析

for a in array:を行っても良いが、プログラム全体をテスト実行し、

z = y + 1

このループ開始時・終了時の変数環境を確認しても

y = x + 1
このプログラムでは、上述のテストでは x を変更し良い。
各変数 xi ∈ X の値が vi の状況で stmt を実行するても z の値は変化しない。これは、上述のテストではと xj ∈ X の値が vj′ となる、という状況を Hoare 論

ループの 1 回分の計算しか行わないためである。し理にならい

かし、ループの複数回の計算を行えば、z の値は x の{x0 = v0 , x1 = v1 , . . .} stmt {xj = vj′ }と記述することにする。また、E = {x0 = v0 , x1 =値の変化の影響を受ける。そのため、これは「依存関係がある」と判定されるのが望ましい。

v1 , . . .} に 対 し 、E[xi1 = w1 , . . . , xim = wm ] は依存関係の解析はリダクション変数・非リダクショ

E 中 の xi1 , . . . , xim に つ い て の 等 式 を そ れ ぞ れ

ン変数の区別をせずに行っている。解析の結果とし

xi1 = w1 , . . . , xim = wm に 置 換 し た も の を 表 すて、自分自身に依存している変数がリダクション変とする。

数だと判明する。これ以外の変数は、本質的には非リダクション変数である。特に、どの変数にも依存し

3. 2 依存性解析とループ分解ていない変数は定数値である。なお、前述のとおり、提案手法では、2. 3 節に従い、まずリダクション変

変数 array は定数値でなければならない。他の変数数間の依存性解析を行い、ループの分解を行う。このに依存していても、自身に依存してない場合、ループために以下のブラックボックステストを行う。

をまたぐ依存関係はないため、自明に並列化できる。

手続き 3.1 (変数依存性の確認).

リダクション変数に依存していたとしても、ループ分

1. 依存関係を空集合で初期化する。
2. 各変数 xi , xj ∈ X に対し以下を行う。
(a) 十分な回数以下を確認する。
i. 各 xk ∈ X に対し値 vk をランダムに生成し、E = {x0 = v0 , x1 = v1 , . . .} と解を行えばループをまたぐ依存関係は消失する。

3. 3 半環・係数の推測と検証

m + 1 個推測することを考える。

3. 3. 1 半環・係数推測の概要

3. 3. 2 定数項の推測依存関係が解析できると、次は各リダクション変数

定数項 a0 の推測方法は半環によらず同じである。

xi について、その変数の計算がどの半環に対応する

E[xi1 = 0, . . . , xim = 0] stmt {xi = a0 }かを推測する。このためには以下のブラックボックスとなる値 a0 を求めればよい。これは、stmt の計算がテストを行う。以下は、変数 xi ∈ X の計算に対するもし線形式で表されるなら、各変数に 0 を代入すれ

半環 R = (S, ⊕, ⊗, 0, 1) のテストとなっており、必要ば他の係数の影響が消え定数項のみが残るという洞に応じて様々な変数・半環に対してテストを行う。ま

察に基づいている。

た、係数の推測は半環毎に手法が異なるため、後で詳しく論じる。

定数項以外の係数の推測方法は半環の性質毎に異なる。

手続き 3.2 (半環と係数の推測).

3. 3. 3 加算の逆元が存在する場合′

1. {xj | xj ∈ X, xi と xj は互いに依存 } = X =値 a−1 が値 a の加算 ⊕ の逆元であるとは、a−1 ⊕a ={xi1 , . . . , xim } とする。

0 であることを意味する。(+, ×) や行列の加算・乗算

2. 十分な回数以下を確認する:からなる半環では加算の逆元が存在する。この場合、(a) 各 xk ∈ X に対し値 vk をランダムに生成する。

(b) {x0 = v0 , x1 = v1 , . . .} stmt {xi = vi′ }となる vi′ を求める。

係数を簡単に推測できる。係数 aj を推測したい場合、

E[xi1 = 0, . . . , xij = 1, . . . , xim = 0] stmt {xi = vi′ }′となる vi′ を求め、aj = a−1
0 ⊕ vi とすればよい。

これは以下の洞察に基づいている。線形式 a0 ⊕(c) 半環 R について、v1 , v2 , . . . と vi′ をもとに

(a1 ⊗ xi1 ) ⊕ · · · ⊕ (am ⊗ xim ) に対し、xij のみに 1、X ′ を変数とする線形式の係数 a0 , a1 , . . . , am

他の変数に 0 を代入すると、a0 ⊕ aj となる。これがを推測する。

stmt の計算と等しいならば vi′ = a0 ⊕ aj であるから、(d) a0 ⊕ (a1 ⊗ vi1 ) ⊕ · · · ⊕ (am ⊗ vim ) ̸= vi′

aj を求めるには vi′ から a0 を引けばよい。

ならば半環 R は不適切である。手続きを終3. 3. 4 加算が乗算に分配する場合了する。

しばしば（可換）半環の加算と乗算を交換しても半3. xi の半環の候補に R を加える。

環をなす。具体的には (∨, ∧) と (∧, ∨)、(min, max)

ここでは、単純に stmt を実行した場合と、半環 Rと (max, min)、(∪, ∩) と (∩, ∪) などである。この場の線形式の係数を推測し、その線形式を計算した場合

合、加算の単位元が乗算の吸収元であることより、とで、結果が一致するかどうかを確認する。結果が一乗算の単位元 1 が加算の吸収元となる。すなわち、致しない場合、半環 R の線形式で計算を表すことは

1 ⊕ a = 1 を満たす。この場合も簡単に係数を求めできない。なお、相互に依存する変数は同時に計算することができる。具体的には、係数 aj を推測したい

る必要があるため、これらをリダクション変数（つま場合、り線形式の変数）として扱う。様々な変数・半環につE[xi1 = 0, . . . , xij = 1, . . . , xim = 0] stmt {xi = vi′ }いてこのテストをを行い、相互に依存する変数間で共となる vi′ を求め、aj = vi′ とすればよい。

通に利用できる半環が見つかった場合、それが実際に並列化に使われることになる。
こ の テ ス ト で は 、各 半 環 で ど の よ う に 線 形 式

a0 ⊕ (a1 ⊗ xi1 ) ⊕ · · · ⊕ (am ⊗ xim ) の係数を推測するかが重要となる。以下の説明では、手続き 3.2の通り、E = {x0 = v0 , x1 = v1 , . . .} を前提として、xi1 , . . . , xim を変数とする線形式のために係数を

これは以下の洞察に基づく。前述の通り vi′ = a0 ⊕ajである。しかし、以下に示すとおり、任意の b に対し、a0 ⊕ (aj ⊗ b) = a0 ⊕ ((a0 ⊕ aj ) ⊗ b) である。

こ の 手 法 は 、(max, +)、(min, +)、(max, ×)、a0 ⊕ ((a0 ⊕ aj ) ⊗ b)=(a0 ⊗ 1) ⊕ (a0 ⊗ b) ⊕ (aj ⊗ b)=(a0 ⊗ (1 ⊕ b)) ⊕ (aj ⊗ b)=(a0 ⊗ 1) ⊕ (aj ⊗ b)(min, ×) などの半環に有用である。例えば (max, +)であれば、ゼロ元 −∞ に十分近い値として、絶対値の非常に小さな数、例えば −230 などを用いること
で、高い精度で係数を推測することができる。

= a0 ⊕ (aj ⊗ b)
よって aj = vi′ として問題がない。

3. 3. 6 加算のみを用いる場合

3. 3. 5 乗算の逆元が存在する場合

依存関係によってループ分解をした後であれば、最後に加算ではなく乗算が逆元を持つ場合、すなわ

少なくないリダクションループでは 1 種類の可換

ち a ⊗ a−1 = a−1 ⊗ a = 1 の場合を考える。ただし、な演算子のみが計算に使われる。さらに、各変数が0 = 1 でない限り、0 には逆元が定義できない。なぜ

高々1 回しか使われないのであれば、線形式は必ず

なら、0

−1

⊗ 0 が 0 なのか 1 なのか定まらないからで

a0 ⊕ x1 ⊕ · · · ⊕ xm の形式となる。つまり、係数は 0か 1 のどちらかである。かつ、変数 xj の値が最終結ある。
乗算に逆元が存在するだけでは係数を求めること

果に直接依存する（つまり、手続き 3.1 で推移閉包をは難しいが、ここでゼロ元 0 に「十分近い」値 r ̸= 0

取る前の依存関係に含まれるなら）ならば、係数はの存在を仮定する。r に求めるのは、ほとんどの a

0 ではありえず、逆に直接依存しないなら 0 である。

および c について、(a ⊗ r) ⊕ c = c となることでよって、1 種類の演算子のみが使われ、変数が高々1ある。この時、以下の手順で係数 aj を推測できる。

回しか使われない場合には常に係数を適切に推測で′E = E[xi1 = 0, . . . , xij = r−1, . . . , xim = 0] とし、

E ′ stmt {xi = vi′ }きる。

3. 3. 7 その他の半環の場合となる vi′ を求め、aj = vi′ ⊗ r とする。

以上の手法で、実用上有用と思われるほとんどの半

これは以下の洞察に基づいている。線形式 a0 ⊕環を扱うことができる。しかし全てではない。例えば、(a1 ⊗ xi1 ) ⊕ · · · ⊕ (am ⊗ xim ) に対し、xij のみに

文字列集合に対する（可換でない）半環 (∪, ·, ∅, {ε})

−1

−1

)と

（X · Y = {vw | v ∈ X, w ∈ Y }、ε は空文字列）に対なる。つまり vi′ = a0 ⊕ (aj ⊗ r−1 ) であるはずであしては、上記いずれの条件も満たさないため、係数をる。これを以下のように式変形する。
a0 ⊕ (aj ⊗ r−1 )推測することはできない。これは今後の課題である。

r、他の変数に 0 を代入すると、a0 ⊕ (aj ⊗ r=(a0 ⊗ r ⊗ r−1 ) ⊕ (aj ⊗ r−1 )= ((a0 ⊗ r) ⊕ aj ) ⊗ r−1
つまり、((a0 ⊗ r) ⊕ aj ) = aj であれば、vi′ ⊗ r =aj ⊗ r−1 ⊗ r = aj である。よって、ほとんどの a0 , ajに対して、vi′ ⊗ r = aj である。
この手法による推測は不完全である。実際の計算で用いられる a0 , aj では ((a0 ⊗ r) ⊕ aj ) = aj が満たされないかもしれない。つまり、テストをパスしたとしても、実際にはこのリダクションループはこの半環の線形式では表現できないかもしれない。しかし、この問題はこの係数推測手法にのみあるわけではない。
提案手法がブラックボックステストに基づいているた
め、全ての係数推測は完全ではない。よって、軽量に
高精度の推測を行えれば実用上十分である。

3. 4 ループ再融合
ここまでは可能な限りループを分解した上で各変数がどの半環に対応するかを推測してきた。しかし、

2. 3 節でも議論したとおり、ループ分解はオーバーヘッドが大きいため最小限にしたい。そこで、推測された半環から、並列化可能性を失わない範囲でループをできる限り再融合する。
ループを融合する際には、変数の依存関係に沿ったものでなければならない。加えて、複数の異なる半環
による線形式が依存関係のある変数の計算の表現に使われてはならない。異なる半環の線形式が混ざると係数を適切に求められなくなるためである。以上をまとめると、ループが融合できる・できないを判断する
条件は以下の通りである。

変数 xi と変数 xj を（xi = xj を含む）1 つのループで計算する場合、xi が xj に依存す

xi1 = 0.
.
.

るならば、xj に依存し xi に依存しない全て

xim = 0の変数は xi や xj と同じループで計算されstmtる必要があり、かつそれらは全て同じ半環で

a0 = xi

計算が表現できなければならない。
ここから導かれるいくつかの有用な結論を挙げる。

• 依存関係がある変数同士を 1 つのループで計算するときは同じ半環でなければならない。

• 依存関係が全くない変数同士は半環が異なっても 1 つのループで計算できる。

• 相互に依存する変数同士は 1 つのループで計算されなければならない。
以上の状況をもとに、融合できる範囲でループを融合してゆくことになる。
発見的にループ融合を行うのではなく、最適なループ融合を目指したいと考えるかもしれない。しかし、最適なループ融合の方針を与えるのは簡単ではない。
例えば、リダクション変数 x1 , x2 , x3 について、x1 と

x2 、x2 と x3 、の計算はそれぞれ 1 つの並列化可能なループにまとめることができるが、x1 と x3 （ひいては x1 , x2 , x3 の全て）を 1 つのループにまとめると半環で表現できなくなってしまうとする。この場合、x2を x1 と x3 のどちらと同じループで計算するのが、大域的にプログラムの効率を最善とするかの判断は自明ではない。詳しい議論は今後の課題である。

3. 5 コード生成
各リダクション変数の計算について、対応する半環が発見できた場合、コード生成を行う。コード生成の方針は、基本的には既存の線形式に基づく並列リダク
ション・並列 scan の生成方法 [10, 13, 15] に従う。提案手法で特徴的となるのは線形式の係数を算出する
コードの生成方法である。
線形式の各係数は、3. 3 節の方針に従って生成する。3. 3 節では stmt 文実行時の各変数の値が具体的に分かっている状況で係数の値を求めていたが、コード生成時には同じ方法で係数を求めるプログラムを出力することになる。例えば、xi の計算の定数項の生成に対応するコードは以下となる。

また、加算の逆元が存在する場合に、xi の計算での変数 xij に対応する係数 aj の生成に対応するコードは以下となる。

xi1 = 0.
.
.
xij = 1.
.
.
xim = 0stmt
aj = a−1⊕ xi0
以上の方針によって得られるコードには多数の stmtの計算が含まれることになり、非効率であると感じられるかもしれない。しかし、その大部分は不要ない
し共通である。そのため、既存のコンパイラ最適化
によって非効率性の大部分は取り除かれると期待できる。

3. 6 計算量提案手法に要するコストは、変数の数、試す半環の数、およびテストを行う回数に依存する。いま、変数の数を V 、半環の数を R、テストの実行回数を n とする。また、テストを 1 回行うのに要するコストは、
仮に変数の数 V に比例するものとする。
依存関係の解析に際しては、各変数の組に対してテストを行う。このコストは O(nV 3 ) である。推移閉
包計算のコストはアルゴリズムによるが、単純なアルゴリズムを用いるとして O(V 3 ) である。
半環の線形式係数の推測については、各変数の係数推測のコストが O(V ) であり、係数を V (V + 1) 回
推測することでテストを行う。そのため、コストはO(nRV 3 ) である。以上をまとめると、並列化可能性の推測に要するコストは O(nRV 3 ) である。
なお、実際にテストを行う回数は n よりかなり少
なくなることに注意してほしい。依存性解析であれば

依存性が発見された時点で、半環の推測であればそのランダム生成を難しくする。さらに、多重ループの場

半環が不適切であると判明した時点で、テストは終了

合、最外ループを並列化するのが最善の方針であるする。そのため、変数同士が複雑に関連し、複雑な計とは限らない。詳細は 4. 4 節で論じるが、内部ルー算を行っているようなプログラムの場合、テストを行

プを並列化するほうが良い効率となる場合もあるし、

う回数はむしろ少なくなる。逆に、変数は沢山あるがそもそも内部ループは並列化できるが最外ループは値を受け渡しているだけで、変数同士にはあまり関係並列化できない場合もある。

がなく、実質的な計算も行っていない（そのためどの以上の状況をふまえ、本節では、プログラムが多重半環にも対応しうる）ようなプログラムでは、多数のループ構造であることが分かっていることを前提とテストを行う必要が生じ、並列化可能性の推測にむしして、以下の方針に従って特化した並列化手法を与ろ時間がかかる。この点は提案手法の特徴である。

える。

• ループ処理を伴うプログラムに対するブラック

4 半環に基づく多重ループの並列化

ボックステストは避ける。上記プログラムであれ前節まででは 1 重のリダクションループに対する半ば、stmt 1・stmt 2・stmt 3 に対してブラックボッ環と線形式に基づく並列化を考えた。本節では前節のクステストを行い、stmt 2 を含むループに対するアプローチが多重ループにも適用できることを示す。

テストは避ける。

本節では以下の形式の 2 重ループを考える。

for a1 in array1:
stmt 1• 各ループに対してモジュラーな解析を行う。これによって、どのループを並列化するかの指針を与える。

for a2 in array2:
stmt 2stmt 3

4. 2 内部ループの独立計算並列化の基本的な方針は前節までの手法と同じでこの形式の 2 重ループを考えるのは本質的な制限である。つまり、各 i について、外側のループの i 回目はない。本節の手法は、複数の内部ループをもつルーの繰り返し処理が線形式の形で表現できればよい。こプ構造や 3 重以上のループ構造に対して自明に拡張のとき、「外側のループの i 回目の繰り返し処理」にできる。

は内側のループが含まれることに注意が必要である。
一般には、内側のループの計算は、それ以前の外側の4. 1 デザイン

ループ（やそれに伴う内側のループ）の計算結果に依

多重ループに対しても、3 節の方法をそのまま使っ存しうる。しかし、一重のループの場合と同様の観察て並列化を行うことはできる。そもそも 3 節の方法

により、内側のループの計算、具体的には stmt 2 のはブラックボックステストに基づいているため、ルー計算が半環上の線形式によって記述できていれば、そプ内部のプログラムがどのような構造をしているか

れ以前の処理とは独立に内側のループの「計算結果」に依存しない。しかし、以下の理由から、このアプを求めることができる。

ローチで常に十分であるとは言い難い。

さらに、外側のループの i 回目の繰り返し処理の結まず、多くの場合、多重ループを伴うプログラムの果と i + 1 回目の結果が適切に併合できるためには、テストには一重ループの場合に比べかなりの実行時間

内側のループの計算結果にあたる線形式が、stmt 1 おがかかる。そのため、ブラックボックステストを繰りよび stmt 3 の計算と適切に併合できる必要がある。こ返す手法は現実的ではないかもしれない。また、ルー

れは、stmt 2 と stmt 1 ・stmt 3 が共通のリダクション

プ内部のプログラムが複雑になると、それに伴いそこ

変数と半環をもつ線形式としてそれぞれ記述できてで処理される配列や添字等の間に一定の制約（不変いれば十分である。

条件）が課される可能性が高い。これは、入力環境の以上の洞察に基づき、2 重ループの並列化の方針は

以下となる。

1. stmt 2 に対して 3 節の方法を適用し、半環とその係数を推測する。

2. stmt 1・stmt 3 に対して 3 節の方法を適用し、半環とその係数を推測する。

ことをブラックボックステストで検証することは不可能ではない。s を変化させながらプログラムを実行
し、結果が s に依存するかどうかを確認すれば良い。
しかし、前提としてこのような解析は避ける方針である。また、リダクション変数を多めに見積もることに

3. 推測された半環・リダクション変数が同一であはデメリットが少ない。後の係数の推測等で多少コスれば並列化は成功であり、各ループについて、3

トはかかるが、最終的に係数が 0 になるだけである。

節の方針に従いコード生成を行う。

よって、各ループについてリダクション変数を推測

これに際してはいくつかの注意点がある。以下の小節ではこれらについて述べる。

し、その総和をとるアプローチの方が現実的である。
変数についての依存関係も同様である。いずれかのループで依存関係がある変数同士は、全てのループで

4. 3 変数依存性とリダクション変数

依存関係があると判断すべきである。やはり、これもまず、リダクション変数は、当該ループのみから判

保守的な判断である。そして、依存関係を保守的に見断されるリダクション変数だけでなく、他のループの積もることによって、本来可能であったループ分解がリダクション変数も考慮する必要がある。例として以不可能になり、ひいては並列化に失敗する可能性もあ下のループを考える。

る。現時点でこの問題が深刻となるような現実的な例

s = 0は見つかっていないが、場合によってはより詳細な解for i in range(0, n):

析を行うことも有意義かも知れない。

for j in range(0, m):s += a[i][j]4. 4 コード生成

このループでは、外側のループだけを見るとリダク

コード生成の方針は基本的には 3. 5 節のものに従ション変数は存在しない。しかし、内側のループではう。しかし、多重ループの場合、どの深さのループでs がリダクション変数であり、そのため、内側のルー

分割統治計算を行うかという選択肢がある。以下、プ

プの「計算結果」が s に対するループを繰り越す依

ロセッサ数を p、外側のループの繰返し回数を n、内存関係を導入する。そのため、s は外側のループでも側のループの繰返し回数を m とする。

リダクション変数とみなされなければならない。
上記例からも分かるとおり、原則として、いずれかのループでリダクション変数である変数は、どのルー
プでもリダクション変数であると判断せざるをえない。もちろん、これはやや保守的な判断ではある。例
として以下のループを考える。

k = 0
for i in range(0, n):s = 0
for j in range(0, m):s += a[i][j]k = max(k, s)

まず、ループ全体について並列化が可能だと判断されたと仮定しよう。この場合、以下の 3 通りの並列化の方針がある。

1. 外側のループを（例えば n/p 回の繰り返し毎に分割し）並列に処理する。

2. 外側のループの各繰り返しに伴う n 個の内側のループを各プロセッサに割当てて並列に計算し、その結果を外側のループで逐次的に集約する。

3. 外側のループの各繰り返しごとに、内側のループをそれぞれ（例えば m/p 回の繰り返し毎に分割し）並列計算しつつ、全体としては逐次的に計算を行う。

このプログラムでは、内側のループでは s がリダク

アルゴリズムは異なるが、n や m が p に比べて十ション変数となっているが、外側のループから見た場分大きければ、どの方法でも理論上の漸近計算量は

合 s はループを繰り越す依存関係をもたない。このO(nm/p) となり理想的な並列速度向上を達成する。

一方で、実際上の性能は多少異なる。一般に、(1) のp = m

方が並列計算タスクの生成回数が少なく（粒度が荒このプログラムは整数の配列 r を用いている。このい）、そのため並列計算に伴う仕事量のオーバーヘッ

配列は実質的には多数の変数の集合として機能してドが小さくなる可能性がある。一方、(3) の方が並列いる。

計算タスクを多く生成し（粒度が細かい）、そのた提案手法は半環上の線形式に基づいている。この

め仕事量は多くなるものの負荷分散に成功しやすい。

時、線形式の各変数に対応するのは配列ではなく値(2) は両者の中間的な性質となる。どのコード生成方（この場合は数値）であるし、また線形式の計算結果

針が適切かは具体的な処理内容やデータの大きさ・並も配列ではなく値である。そのため、配列に対して

列計算機環境等に依存するため、一概にはいえない。

は、配列そのものではなく配列の各インデックスをリ

また、 stmt 2 に対しては適切な半環と線形式が推

ダクション変数とみなし並列化を行う必要がある。し測できたが、stmt 1・ stmt 3 については上手くいかな

かし一方で、本研究ではブラックボックステストに基かった、というような場合もありうる。この場合、(1)

づいているため、具体的に配列のどのインデックスにの方針による並列化は不可能だが、 (2) ・ (3) の方針

注目すれば良いかを標準的な静的解析によって求めるによる並列化は可能である。

ことはできない。

提案手法は、各ループについてリダクション変数や以下では、この問題を解決する。すなわち、配列を半環を推測してゆき、それを最終的に集約して並列化含むプログラムに対して、ブラックボックステストに方針を判断するという方針を採用している。これに基づき、配列を用いたプログラムの並列化を本研究でよって、どのループまでであれば並列化が可能なのか扱えるものに帰着する方法を示す。

を解析結果から判断できる。さらに、半環上の線形式を用いたリダクション並列化ではリダクション変数の5. 1 配列変数に対するアクセスの推測数（つまり係数の数）に比例したオーバーヘッドを要

基本的な着想は、具体的に配列のどの変数に対してするため、並列化に伴うオーバーヘッドがどの程度で読み書きが行われたかを記録することである。3 節とあるのかも見積もることができる。そのため、提案手

同様の形式のループの並列化を考え、このループが配法では、これら見積もりをふまえフレキシブルにコー

列 x を扱っていると仮定しよう。この時、以下のよう

ド生成を行ったり、複数の方針でコードを生成し実験な手法によって、配列に対する読み出し・書き込みを判定できる †2 。

的に性能を比較したりといったことが可能である。

• stmt 文実行前の x[i] のみを変化させると stmt

5 中間データ構造として配列を用いる場合

文実行後のいずれかの変数値が変化する場合、

特に多重ループからなる並列リダクションでは、中x[i] に対する読み出しがあったと判断できる。

間的な値を配列に格納する場合が少なくない。例とし• stmt 文実行前の x[i] と実行後の x[i] の値が異て以下のプログラムを考える。これは最長共通部分なっている場合、x[i] に対する書き込みがあっ

列 [2, 13] の計算に対応する。

たと判断できる。

for i in range(0, n):p = 0
pp = 0for j in range(0, m):

これらにより、特定の実行における読み出し元・書き込み先を推測することができる。
しかし、特定の実行におけるアクセスを確認する
だけでは不十分である。並列化を行うためには、ルー

tmp = pp + 1 if a[i] == b[j] else 0m = max(tmp, r[j], p)
pp = r[j]r[j] = m†2 並列化対象の言語によっては、多少のメタプログラミ
ングの機構（例えば配列の添字アクセスや代入に対す
るオーバーロード）を使い、配列に対するアクセス・
書き込みを直接監視する方法もありうる。

プの何回目の繰り返し時には配列のどのインデック果が得られる。しかしこれは正しくない。x[i] の値

スにアクセスするかを知る必要がある。一般に、ほとは次のループでの x[i+1] の値になる。このように、んどの配列の添字アクセスは、他の変数値の（(+, ×)配列の場合、添字アクセスの変化による依存関係が発半環上の）線形式をなす。本研究でもこれを仮定し、生しうる。これを無視することは、正しくない並列化

3. 3 節の手法によってこの線形式を推測する。具体的を行うことに繋がりかねないには以下の手順による。
手続き 5.1 (配列インデックスアクセスの推測).

1. 3. 2 節の方法によって、配列の添字アクセスに影響する変数集合を得る。

上記問題を解消するには、配列の各添字に対する依存関係は同一視すればよい。つまり、いずれかの添字に対して依存関係があれば、どの添字に対しても依存関係があると判断するのである。これは保守的で

2. 添字アクセス先が (+, ×) 上の線形式であるとあり、実際には依存関係がない場合もあり得る。配列仮定し、係数を推測する。具体的な推測方法はの添字アクセスの形式とループカウンタの変化を注3. 3 節の通りである。

意深く確認すればより精密な解析も可能であろうが、

3. 十分な回数テストを行い、推測が正しいことをブラックボックステストに基づく本手法では、精密な確認する。推測が正しくなかった場合、並列化は解析を行う方法は少なくとも非自明であり、またそれ失敗する。

によってどの程度の利得があるかも不明瞭である。

前述の最長共通部分列のプログラムを例に説明を行う。まず、各変数値を変えつつ配列 r の読み出し
元・書き込み先を確認することで、変数 j の値が読み出し元・書き込み先に影響していることが分かる。

6 評価実験
本節では提案手法の有効性を調べるために行った評
価実験について報告する。

次に、j = 0 および j = 1 の場合の読み出し元・書き込み先から、線形式はどちらも 0 + 1 × j と推測でき

6. 1 実装る。最後に、十分な回数を行いテストを行い、推測が本稿で提案した並列化可能性判定手法を Python 言

正しいことを確認する。

語で実装した。なお、実装は 350 行ほどである。

以上によって、読み出し元・書き込み先が推測でき半環候補としては (+, ×)、(max, +)、(max, min)、

たら、これらを配列の代わりに変数集合に加え、3 節(min, max)、(∧, ∨)、(∨, ∧)、(max, ×) (非負値のみの方法に基づき並列化を行う。例えば最長共通部分列を扱う) を用意した。

の例であれば、配列 r ではなく r[j] を整数値変数とみなし、並列化を行う。

入力プログラムはループ内の計算を表す Python の関数とし、これを exec 関数を用いて実行する。さらに、各引数・返値の型を追加の情報として与えてい5. 2 配列に対する依存性解析

る。型としては整数・真偽値・整数値配列を許してい配列を中間的に用いている場合、変数を用いる場合る。浮動小数点数を入力とするプログラムについてとは異なる形で「変数」間の依存関係が生じる場合がは、入力のランダム生成と結果の同値判定が複雑となある。例として以下のプログラムを考える。

るため、原則として整数を扱うとして処理している。

for i in range(0,n):x[i] = x[i+1]このループに対し前小節の手法に基づき配列のアク実装はおおむね本稿の内容に従っているが、コード生成（3. 5 節）
・ループ再融合（3. 4 節）
・配列の添字の推測（5. 1 節）の機能は実装していない。

セスを推測すると、x[i] および x[i+1] という「変提案手法では入力プログラムの引数としてランダ数」を考える必要があることが分かる。これに続いてムな値を与えるが、プログラムによっては引数値に一変数依存関係の解析を行うと、x[i] は x[i+1] に依定の制約がある場合がある。今回の実装では、制約存するが x[i+1] は x[i] には依存しない、という結は assert 文としてプログラム中に記述することにし表1

実験結果: 数値の計算に関するリダクション

プログラム

ループ分割

演算子

並列化可能

実行時間（秒）総和 [1, 3, 4, 12, 14]

+

○

0.05

偶数和 [3]

+

○

0.04

正要素和 [3]

+

○

0.05

+, +

○

0.12

正要素数 [3]

+

○

0.04

特定要素数 [3]

+

○

0.07

平均 [1, 3]

○

+

○

0.09

max

○

0.06

max, max○

0.17

最大絶対値 [3]

max

○

0.08

最小値 [1, 3]

min, min

○

0.04

○

0.12

不連続部分数 [8]○

最大値 [1, 3, 4, 14]

2 番目最大値 [3]○2 番目最小値 [1, 8]○

min

最大最小差 [3]○max, min○0.06

最大要素数 [3]○max, +○

0.16

最小要素数 [3]○min, +○

0.11

内積 [14]

+

○

0.05

ハミング距離 [1]

+

○

0.05

(+, ×), ×△

0.09

複素数積 [14](+, ×)○

0.08

二重指数平滑化 [8](+, ×)△

0.09(+, ×)△

0.08(+, ×)△

0.13

多項式 [1, 8, 13]

三重対角行列 LU 分解 [13]○○差分法 [13]た。システムは、AssertionError を受け取った場合には入力を生成し直す †3 。

は任意の半環で表現できる。提案手法では不適切な半並列化可能性のテストに要する時間は主にいくつ

規定回テストする必要があり、よりコストがかかる。

の入力を試すかに依存する。今回の実験では各ステッ

値を受け渡すだけの変数は、適当な半環（実装上はプで 1000 通りの入力を試すようにした。

(+, ×)）の係数を推測する際に、追加で「0 でない係

環はテストですぐに発見される一方、適切な半環は数がたかだか 1 つであり、それが 1 であること」を

6. 2 実装上の効率化テストすることで確認できる。確認に成功した場合、本稿で述べた内容に加え、2 点実装上の効率化を

以降はテストを行わず「全ての半環で表現できる」と行った。

して処理する。

1 つは、実質的な計算を行わずただ値を受け渡すだ

もう 1 つは、多重ループ時のテストの順番である。

けのリダクション変数の検出である。このような計算

4 節では、各ループに対して独立に依存性解析と適切†3 ただし、係数推測のために単位元等を入力した際に、それが制約に反してしまい AssertionError となってしまう場合もある。この場合にはシステムは入力を再生成せず、その半環での並列化が失敗したものとして扱う。

な半環の推測を行う手法を提案した。しかし、多重ループ全体を並列化する前提であれば、それぞれの変数依存性や半環の推測において、各ループを巡回的にテストする方が効率が良い。例えば、2 重ループにお表2

実験結果: 最適化問題・判定問題・クエリに関するリダクション

プログラム

ループ分割

演算子

並列化可能

実行時間（秒）

最大連続 1 [1]○+, max○0.17

最大 1 間距離 [3, 8]+, max○0.13

最大 0 間要素和 [3]+, max○0.15

最大接頭辞和 [1, 8, 12]○+, max○0.23

最大接尾辞和 [8, 13]○(max, +)○

0.07

最大部分列和 [1, 3, 4, 8, 12, 13]○(max, +), max

○

0.52

最大部分列積 [8]○(max, ×), max

○

0.26

全て同じ [3]

○

∧

○

0.06

+

○

0.05

01 数一致 [3]括弧対応 [1, 8]

○

+, ∧

○

0.09

可視地点判定 [12]

○

max

○

0.11

dropwhile [1]○∨, max○0.07

要素検索 [3]○∨, max○

0.06

整列済 [1, 3]

○

∧

○

0.06

0*1* [1]

○

∧

○

0.06

(01)* [3]

○

∧

○

0.06

∧

○

0.02

0 が先頭以外になし [3]
0 が 1 の後にある [1]

○

∧,∧

○

0.06

1*のマッチ数 [3]

○

+, +

○

0.12

1*2 のマッチ数 [3]

○

+, +

○

0.14

10*2 のマッチ数 [3]○+, +, +○

0.28

1*2*3 のマッチ数 [3]○+, +, +○

0.31

10*20*3 のマッチ数 [3]○+, +, +, +, +

○

0.66

いて、外側のループでは変数に依存性がないが内側のLinux 上の Ubuntu 18.04、処理系は CPython 3.6.5

ループでは依存性がある場合、外側のループを何度テ（最適化オプション-OO 使用）である。

ストしても依存性は見つからないが、外側のループと並列化対象プログラムとしては、リダクションルー

内側のループを交互にテストすればすぐに依存性をプの自動並列化を行った先行研究 [1–4, 8, 12–14] で発見できる。
以上 2 つの効率化は、いずれも並列化に比較的時評価実験に用いられていたプログラムを全て用意した †4 。ただし、明らかに重複のある例については省間のかかるプログラムの効率を改善する。ただし、改いている。例えば、「最大部分列和」については、最善はせいぜい数倍程度であり、効率化なしでは並列化

大和だけでなくその部分列の開始インデックス・終了

が事実上不可能なプログラムの並列化が可能になる

インデックスも求めるプログラムのみを対象とし、最わけではない。

大和のみを求めるプログラムは省いている。

6. 3 実験環境と並列化対象プログラム
実験環境の CPU は Intel Core i7-7600U CPU (2.80

GHz)、メモリ 16GB、OS は Windows Subsystem†4 「行最小値の最大値」については、[2] のプログラムが明らかに正しくなかったため、少なからずプログラ
ムを修正している。なおこれにより条件分岐が増え、並列化はより難しくなっている。

表3

プログラム

実験結果: 多重ループのリダクション

ループ分割

2D 総和 [2]

演算子

並列化可能

実行時間（秒）

+

○

0.05

2D 整列済 [2]

○

∧

○

0.05

4D 最大値座標 [14]

○

min

○

0.63

上下整列済 [2]

∧

○

0.02

対角整列済 [2]

∧

○

0.03

行値範囲単調増加 [2]○max, ∧○

0.19

行値範囲 ovelap [2]○max, min, ∧

○

0.55

行値範囲単調縮小 [2]○max, min, ∧

○

0.55

∧

○

0.02

全行値範囲 overlap [2]行最小値の最大値 [2]○min, max○

0.16

列最小値の最大値 [2]○min, max○

0.23

鞍点 [2]○

min, max, min, max

○

0.51

2D 最大接頭辞和 [2]○2D 最大接尾辞和 [2]+, max○0.14(max, +)○

0.03

2D 最大部分列和 [2]○(max, +), max

○

0.11

最大左上部分配列和 [2]○+, +, max○

0.27

最大左下部分配列和 [2]○+, (max, +), max

○

0.24

最大右上部分配列和 [2]○+, (max, +), max

○

0.27

3D 最大接頭辞和 [2]○+, max○0.14(max, +)○

0.03

3D 最大接尾辞和 [2]3D 最大部分列和 [2]○(max, +), max

○

0.11

3D 最大左接頭辞和 [2]○+, +, +, max

○

0.43

括弧対応列数 [2]○+, ∧, +○

0.21

最頻出値 [2]○+, max○0.41

2 配列最大差 [2]

○

max

○

0.10

最長括弧対応 [2]○+, ∧, max○

0.33

(max, +)○

0.19

最長共通部分列 [2, 13]
独立要素数 [3]−2D ヒストグラム [14]−プログラムはできる限り自然な記述となるようつ6. 4 実験結果

とめた。特に、最大値・最小値に関するプログラムで実験の結果を表 1・表 2・表 3 に示す。「ループ分は、必要に応じて条件分岐も用いるなど、必ずしも半割」の列はループ分割が可能であると判定されたかど環に対応させることを意図しないプログラムとしてうかを表す。「演算子」の列は、どの演算子が使われいる。

ていると判断されたかを表す。ループ分割が行われた
場合、分割後の各ループについて使われた演算子がそれぞれ示されている。例えば「平均」の例では、2 つのループのそれぞれで + が使われていることを表す。

なお、提案手法は使われている半環を推測するが、演により、このゼロ除算を排除した。

算子が 1 種類しか使われいない場合、その演算子を

その他の例についても、提案手法での並列化を成功含む半環を列挙することになる。そのような場合にはさせる上で配慮が必要であったものがいくつかあった。

半環ではなく演算子のみを挙げている。
「並列化可能」• 「多項式（係数等に実数が現れる場合）」
・「二重の列は、最終的に並列化が可能であると判定された指数平滑化」
・「差分法」では、本質的に実数を扱場合に「○」、並列化可能と判定されるためにプログう。しかし、現状の実装では数値誤差を考慮してラムの書き方等に注意が必要だった場合に「△」、並いないため、推測した線形式と実際の計算の値が列化可能性の判定に失敗した場合に「−」をそれぞれ微妙に異なるため並列化に失敗する可能性があ

記載している。「実行時間」の列には、並列化可能性る。今回は数値誤差が生じない値（例えば 1/4）の判定に要した時間を記載している。システムが行っを用いて実験した。実用上は、並列化が数値誤差

た判定の正しさについては著者が目視で確認したが、にどう影響するかも含め、より踏み込んだ議論が誤りはなかった。

必要かもしれない。

実行時間は全ての例について 1 秒を大きく下回っ• 「最大部分列積」の計算では、最大値となる部

た。Python インタプリタでの実行であること、効率分列だけではなく最小値となる部分列も覚える

化等を全く考慮していないことを考えると、十分な性必要がある。これは、負数が現れたときに大小が能であると判断できる。

反転するためである。これをプログラムとして記並列化に失敗したのは「独立要素数」「2D ヒスト

述する際、負値を用いると (max, ×) の半環で扱グラム」の 2 例であった。理由は、「独立要素数」はえなくなるため、代わりにその絶対値を記憶する集合、「2D ヒストグラム」はベクトルを扱っている必要がある。

ためであった。これは今回の実装上の問題であり、適

しかし、配慮が必要だった例は多くなかった。プログ

切な型と半環を加えれば並列化が可能であると推測

ラム解析・変換ではプログラムの書かれ方によって利している。

用可能性が変化するのは一般的である。少なくとも、並列化が最も困難だったのは「三重対角行列 LU 分

他のリダクション自動並列化手法で必要な配慮と比解」であった。この計算では数値の加算・乗算に加え

べ、より注意深い配慮が必要であったとは考えていて除算を用いているため、半環で捉えることができない。

ない。Sato ら [13] はこのような問題もうまく式変形
することで半環上の線形式として定式化できることを論じているが、これはリダクション変数集合を変化させるものであり、提案手法の範囲で行うことはで

7 関連研究
リダクションループおよび関連する構造の自動並列化については少なくない研究がある。

きない。今回は手作業でこの変形を行った。さらに、

Fisher と Ghuloum [4] は複雑なリダクションルー

これを本実装で扱うには手作業でループ分解を行う

プの自動並列化について先駆的な研究を行っている。

必要があった。並列化には (+, ×) 半環の係数を推測

彼らの手法は、ループ本体のプログラムに対し、リする必要があり、そのためにはリダクション変数に 0

ダクション変数を「穴」と見なしたテンプレートを考を代入することになるのだが、これが除算の分母に使え、それらの複数回の合成が同じテンプレートに単純われてエラーとなってしまう。この除算は相互依存関化できるかをプログラムの等価変形を繰り返して確係のない変数の値を求める計算であり、係数推測には認するものである。後に Morihata と Matsuzaki [11]

関係がないのだが、ブラックボックステストである以は同様の手続きを限定記号除去を用いて定式化し直上、プログラムのどの部分でゼロ除算がなされているしている。また、Farzan と Nicolet [1] は構文主導合かわからず、しかもその部分の処理を取り除くことも成（Syntax Guided Synthesis）の枠組みに基づき、できない。そのため、手作業でループ分解を行うこと同様の自動並列化手法を提案している。これら手法は

一般性は高いものの、適切なテンプレートをどう得る

を再構成する。このアプローチは、複数回のテスト

かや合成時の等価変形をどう行うかの指針がないた

ケース実行の結果から係数を推測する本研究のアプ

めに、並列化手続きで多くの可能性を発見的に探索す

ローチによく似ている。しかし、両者には 2 点の重要る必要が生じやすく、複雑なプログラムの並列化は現な違いがある。まず、彼らのアプローチでは、ループ実的ではなかった。そのため、多くの研究では、並列の処理結果から要約結果を得る方針を静的解析によっ化の際に考えるプログラムのパターンをある程度限て決定している。これに対し本研究では、どの半環が定することで、現実的な自動並列化を目指している。

適切かのブラックボックステストによって、要約結果

Xu ら [15] は計算が半環（およびそれの拡張にあを得る方針を決定している。また、彼らの手法は数値たる代数構造）上の線形式であることを確認するの加算・乗算・最大値計算に特化した手法と、真偽値静的解析を提案し、これによってリストを集約するなどの有限値域の値に特化した手法から構成されて再帰関数が並列化できることを示した。Mastuzaki

おり、例えば和集合・積集合からなる半環やベクトルら [9] は木構造を集約する再帰関数について、Sato と計算などには適用できない。これに対し、本研究は半Iwasaki [13] はリダクションループについて、よく似環に注目することで、より広い範囲のプログラムを系た手法を提案している。1 節でも述べたとおり、これ統的に並列化することができる。

らの研究は静的解析に依存しており、明示的に半環のFarzan と Nicolet [2] は、多重ループに対する自動演算子を用いてプログラムが記述されていなければ適並列化手法を提案している。彼らの手法は、外側の用が難しかった。この問題を軽減するため、Sato とループの各繰り返しに対する内部ループの計算の集

Iwasaki [13] は条件分岐で記述されたプログラムから約結果を独立に求めることで、それら集約結果を集約

SMT ソルバーを使って max 演算子を抽出する手法

するネストのないリダクションループの並列化に帰着も与えている。本研究は、これらの研究の欠点をブする、というアプローチである。特に、内部ループのラックボックステストを用いることで改善している。

集約結果に注目し、それと外側のループとの関係を

これにより、SMT ソルバー等に依存することなく、確認する点で、我々の並列化手法との類似性がある。

明示的に半環の演算子を使っていないプログラムに対しかし、彼らの手法が「外側のループを並列化するた

しても並列化を行うことができる。

めに内部ループの構造を調べる」という、いわばトッ

Suriana ら [1] は Halide という高性能配列計算用

プダウンのものであるのに対し、我々の手法は「全領域限定言語のためのリダクション並列化手法についてのループの構造を独立に解析し、それらをもとにて論じている。彼らの手法はリダクション並列化が可

各ループの並列化可能性を調べる」という、いわばボ能な演算子群を事前に探索的に生成しておき、プログ

トムアップの構造をもつ点に違いがある。ボトムアッラムがその演算子群によるリダクションに対応するか

プにしたことにより、各ループに対する解析はよりモ

どうかを記号的に確認するものである。これは、半環ジュラーになり、またフレキシブルな並列化・コードの候補を事前に列挙しておき、その演算子でプログラ生成が可能になっている。このボトムアップな構造

ムが記述されているかどうか確認する提案手法と類は、半環という代数構造に着目したことの大きな利点

似している。しかし、提案手法では探索的な生成を避である。

け、代わりに半環に注目することでより系統的に演算

Fedyukovich ら [3] は反例主導帰納的合成（Counter

子群の特定を行っている。また提案手法では、明示的Example Guided Inductive Synthesis）の枠組みににその演算子でプログラムが記述されていなくても、基づき、SMT ソルバーを用いた自動並列化手法を提意味上等価なプログラムは同様に並列化できる。

案している。より具体的には、簡単な並列計算パター

Jiang ら [8] のリダクション自動並列化手法では、

ンから始め、並列化対象プログラムがそのパターン

実行時に複数の初期値からループの処理を行い、そので記述できるかを SMT ソルバーに問合わせ、不可能

結果を用いてループ処理の要約結果（例えば線形式）

なら徐々に複雑で表現力の高い並列計算パターンを

試す、という方針で並列プログラムを合成する。我々

象のプログラムを実行することで得る。特に、半環ののアプローチで得られる並列プログラムは彼らの考係数推測時には、半環の性質に基づき注意深く入力をえた並列計算パターンには直接は合致しないため、両与えている。これにより、プログラム合成を軽量かつ者の能力の比較は難しい。ただし、SMT ソルバーを効果的に実現できている。

何度も呼び出すという戦略上、彼らの手法が複雑で大きなプログラムにもスケールするという見通しは

8 まとめと今後の課題

乏しい。なお、単純な並列計算パターンから始め徐々

本稿ではブラックボックステストによってリダクに複雑なパターンを複雑化するアプローチは、本研究

ションループの並列化を行う手法を提案した。この手

でもより複雑なプログラムを並列化する際には有用

法は、ループでの計算が半環上の線形式に対応するかかもしれない。

どうかを推測することを要旨としている。また、ルーリダクションループの自動並列化については、これプの分解や配列へのアクセスの推測など、関連する問らの他にも第三リスト準同形定理 [6] に基づく手法が題についても同様にブラックボックステストによってある。Geser と Gorlatch [5] は第三リスト準同形定理行う手法を提案している。提案手法は多重ループに対と反単一化を組み合わせることで、Morita ら [12] はしても適用可能である。

第三リスト準同形定理とプログラム逆化を組み合わせ

コンセプト確認のための実装を用いた実験では、既ることで、森畑ら [16] は第三リスト準同形定理を並存のリダクションループの並列化で用いられたほとん

列プログラム候補の生成検査と組み合わせることで、ど全ての例を並列化可能であると正しく判定できた。

それぞれ自動並列化を行う手法を提案している。第三

このことは、半環上の線形式という抽象化が、リダク

リスト準同形定理に基づく場合、プログラマは通常のションループの自動並列化に際しては極めて有効な抽ループだけではなく、通常とは逆順に処理するループ象化であることの一つの証左であると考えている。

のプログラムも同時に与える必要がある。そのため現状の実装はコンセプト確認のためのものであり、プログラマへの負担が重い。一方で、少なくとも今回さらなる議論にはより緻密な実装が必要であろう。特の実験で確認できた範囲では、リダクションループのに、3. 5 節で論じたとおり、提案手法で得られる並列

自動並列化に際して 2 種類のループを与えることのプログラムにはかなりの冗長性が含まれるため、それ

恩恵はそれほど大きくないように見える。なお、森畑が本当に効率に影響を与えないかは確認する必要がら [16] の手法で用いている候補生成手法は、並列化あろう。また、適切なループ融合の方針についてもよ

対象プログラムの構文的な書き換えに基づくものであり詳細な検討が求められる。

り、半環に基づく本研究の手法よりはむしろ Farzanと Nicolet [1] の手法に近いものとなっている。
リダクションループの自動並列化は自動プログラム

謝辞 東北大学電気通信研究所共同プロジェクト研究集会での議論に感謝する。本研究は JSPS 科研費

18K18032 の助成を部分的に受けたものである。

合成 [7] の一種とみなすことができる。自動プログラ
ム合成は、特定の仕様（例えば与えられたプログラム

参 考 文 献と等価で、並列化可能なパターンに合致する）を満た[ 1 ] Farzan, A. and Nicolet, V.: Synthesis of divideand conquer parallelism for loops, Proceedings of the
38th ACM SIGPLAN Conference on Programming
Language Design and Implementation, PLDI 2017,
Barcelona, Spain, June 18-23, 2017, Cohen, A. andVechev, M. T.(eds.), ACM, 2017, pp. 540–555.
[ 2 ] Farzan, A. and Nicolet, V.: Modular divide-andconquer parallelization of nested loops, Proceedingsof the 40th ACM SIGPLAN Conference on Programming Language Design and Implementation,
PLDI 2019, Phoenix, AZ, USA, June 22-26, 2019,

すプログラムを自動的に構成する技法である。実際、少なくない自動並列化手法 [1–3, 11, 16] は明示的に自動プログラム合成に由来する手法を用いている。
自動プログラム合成の観点からみた場合、我々の手法は入出力例によるプログラミング（Programming

by Example）の一種だとみなすことができる。ただし、入出力例はユーザが与えるのではなく、並列化対

McKinley, K. S. and Fisher, K.(eds.), ACM, 2019,
pp. 610–624.
[ 3 ] Fedyukovich, G., Ahmad, M. B. S., and Bodı́k,R.: Gradual synthesis for static parallelization of
single-pass array-processing programs, Proceedingsof the 38th ACM SIGPLAN Conference on Programming Language Design and Implementation,
PLDI 2017, Barcelona, Spain, June 18-23, 2017,
Cohen, A. and Vechev, M. T.(eds.), ACM, 2017,
pp. 572–585.
[ 4 ] Fisher, A. L. and Ghuloum, A. M.: Parallelizing
Complex Scans and Reductions, Proceedings of the
ACM SIGPLAN’94 Conference on Programming
Language Design and Implementation (PLDI), Orlando, Florida, June 20-24, 1994, 1994, pp. 135–
146.
[ 5 ] Geser, A. and Gorlatch, S.: Parallelizing functional programs by generalization, Journal of Functional Programming, Vol. 9, No. 6(1999), pp. 649–
673.
[ 6 ] Gibbons, J.: The Third Homomorphism Theorem, Journal of Functional Programming, Vol. 6,
No. 4(1996), pp. 657–665.
[ 7 ] Gulwani, S., Polozov, A., and Singh, R.: Program Synthesis, NOW, 2017.
[ 8 ] Jiang, P., Chen, L., and Agrawal, G.: Revealing
parallel scans and reductions in recurrences throughfunction reconstruction, Proceedings of the 27th International Conference on Parallel Architectures
and Compilation Techniques, PACT 2018, Limassol, Cyprus, November 01-04, 2018, Evripidou, S.,
Stenström, P., and O’Boyle, M. F. P.(eds.), ACM,
2018, pp. 10:1–10:13.
[ 9 ] Matsuzaki, K., Hu, Z., and Takeichi, M.: Towards automatic parallelization of tree reductions
in dynamic programming, SPAA 2006: Proceedings of the 18th Annual ACM Symposium on Parallel Algorithms and Architectures, Cambridge, Massachusetts, USA, July 30 - August 2, 2006, Gibbons, P. B. and Vishkin, U.(eds.), ACM, 2006,
pp. 39–48.
[10] Morihata, A.: Lambda calculus with alge-

braic simplification for reduction parallelization byequational reasoning, Proc. ACM Program. Lang.,
Vol. 3, No. ICFP(2019), pp. 80:1–80:25.
[11] Morihata, A. and Matsuzaki, K.: Automatic
Parallelization of Recursive Functions Using Quantifier Elimination, Functional and Logic Programming, 10th International Symposium, FLOPS
2010, Sendai, Japan, April 19-21, 2010. Proceedings, Blume, M., Kobayashi, N., and Vidal, G.(eds.), Lecture Notes in Computer Science,Vol. 6009, Springer, 2010, pp. 321–336.
[12] Morita, K., Morihata, A., Matsuzaki, K., Hu, Z.,
and Takeichi, M.: Automatic inversion generates
divide-and-conquer parallel programs, Proceedingsof the ACM SIGPLAN 2007 Conference on Programming Language Design and Implementation,
San Diego, California, USA, June 10-13, 2007, Ferrante, J. and McKinley, K. S.(eds.), ACM, 2007,
pp. 146–155.
[13] Sato, S. and Iwasaki, H.: Automatic parallelization via matrix multiplication, Proceedings of the
32nd ACM SIGPLAN Conference on Programming
Language Design and Implementation, PLDI 2011,
San Jose, CA, USA, June 4-8, 2011, Hall, M. W.
and Padua, D. A.(eds.), ACM, 2011, pp. 470–479.
[14] Suriana, P., Adams, A., and Kamil, S.: Parallel associative reductions in Halide, Proceedings ofthe 2017 International Symposium on Code Generation and Optimization, CGO 2017, Austin, TX,
USA, February 4-8, 2017, Reddi, V. J., Smith, A.,
and Tang, L.(eds.), ACM, 2017, pp. 281–291.
[15] Xu, D. N., Khoo, S.-C., and Hu, Z.: PType
System: A Featherweight Parallelizability Detector, Programming Languages and Systems: Second Asian Symposium, APLAS 2004, Taipei, Taiwan, November 4-6, 2004. Proceedings, Chin,
W.-N.(ed.), Lecture Notes in Computer Science,Vol. 3302, Springer, 2004, pp. 197–212.
[16] 森畑明昌, 松崎公紀, 胡振江, 武市正人: 並列プロ
グラムの候補生成と適合性検査による並列化, 情報処
理学会論文誌：プログラミング, Vol. 2, No. 2(2009),
pp. 132–143.