狠狠撸

1
Halide for Memory
Halide勉強会＠フィックスターズ 2018/7/28 (7/29追記)

2
自己紹介
Twitter @koumeitomida
某社（製造業）勤務
HWエンジニア
計算機、ASIC、並列画像処理(multi-core, SIDM, GPU)、
最近はDeep Neural Network
管理職のお仕事も

3
自己紹介（続)
画像処理パイプライン
画像サイズ数千画素×数千画素×数色
有向グラフ DAG (Directed Acyclic Graph)
入力画像
画像
処理
中間画像
画像
処理
中間画像
画像
処理
出力画像
画像
処理
中間画像

5
Halideの概要
この勉強会できっと誰かが分かりやすく説明しているはず。
また、すでに優れた説明もあるので省略。
例えば
Fixstars 丸岡さん /fixstars/halide-82788728
名工大福嶋先生 https://qiita.com/fukushima1981/items/fa3537234e19baffc761
/FukushimaNorishige/ss-83087912

6
Halideの面白さ
JRKのDr.論文 http://people.csail.mit.edu/jrk/jrkthesis.pdf
Domain Specific Language (DSL)
Halide
特定領域にフォーカス
＋ High Performance
＋ Productivity
－ Turing Completeを放棄
完璧な計算機言語
... 残念ながら存在しない
「性能」と「生産性」が欲しい
場合、ドメインに応じて何かを
捨てる。
画像(テンソル)処理用DSL
Halide

7
Algorithm記述
画素値(テンソルの要素)を求める計算について記述する。
Algorithm記述 = 計算グラフ
f1 f2 f3
f1の出力
(中間画像)
f2の出力
(中間画像)
Func f3;
f3(x,y,z)
= f2(x,y,z) / 9;
Func f2;
f2(x,y,z)
= f1(x,y-1,z)
+ f1(x,y ,z)
+ f1(x,y+1,z);
Func f1;
f1(x,y,z)
= in(x-1,y,z)
+ in(x ,y,z)
+ in(x+1,y,z);
Buffer<T> in;

8
計算の依存性その1
要素に着目
f3(x,y,z)の計算のために、f2(x,y,z)が要る。
f2(x,y,z)の計算のために、f1(x,y-1,z), f1(x,y,z), f1(x,y+1,z)が要る。
f1(x,y,z)の計算のために、in(x-1,y,z), in(x,y,z), in(x+1,y,z)が要る。
f1 f2 f3
f1の出力
(中間画像)
f2の出力
(中間画像)

9
計算の依存性その2
テンソルに注目
f3(x,y,z)の計算のために、f2(x,y,z)が要る。
f2(x,y,z)の計算のために、f1(x,y-1,z), f1(x,y,z), f1(x,y+1,z)が要る。
f1(x,y,z)の計算のために、in(x-1,y,z), in(x,y,z), in(x+1,y,z)が要る。
f1 f2 f3
f1の出力
(中間画像)
f2の出力
(中間画像)
この赤い領域
を欲しい場合

10
Schedule
f3の出力画像(テンソル,直方体)を得るために、
f3(x,y,z)を、x,y,zでスキャン。スキャン順は?
f2の出力(直方体)を得るためのスキャン順は?
f1の出力(直方体)を得るためのスキャン順は?
更に、メモリにどう配置しようか?、etc.
f1 f2 f3
f1の出力
(中間画像)
f2の出力
(中間画像)
Var x,y,z;
Schedule
x
y
z

11
Schedule : reorder
HalideのFuncは副作用の無い関数なので
x,y,zのスキャン順は任意
ループの順番を変更 f3
x
y
z
//デフォルト
//疑似コード
for z:
for y:
for x:
f3[x,y,z] = ...
f3.reorder(z, x, y);
//疑似コード
for y:
for x:
for z:
f3[x,y,z] = ...
f3(x,y,z) = ...
zは外側
xは内側
C/C++とは逆順！

12
Schedule : parallel
Fork/Joinの並列処理
? Taskを生成しQueueに入れる
? Worker-threadが実行
? 全taskの実行完了を待つ
f3
x
y
z
//デフォルト (シリアル実行)
//疑似コード
for z:
for y:
for x:
f3[x,y,z] = ...
f3.parallel(z);
//疑似コード
parallel_for z:
for y:
for x:
f3[x,y,z] = ...
//ここでjoin
Task
Queue
Worker
Worker
Worker
Fork/Join

13
Schedule : split
例えばx軸を２つに分ける。
xo = x / 32;
xi = x % 32;
f3
x
y
z
//デフォルト
//疑似コード
for z:
for y:
for x:
f3[x,y,z] = ...
f3.split(z, zo, zi, 32);
//疑似コード
for zo:
for zi:
for y:
for x:
f3[x,y,zo*32+zi] = ...
z軸をzoとziにsplit

14
Schedule : fuse
Splitの逆、2つの軸を一つにする。
t = y * X + x
for loop (条件分岐)は遅いのでloopを減らしてみる価値あり。
並列処理(後述)のworkを増やすのに便利。
//デフォルト
//疑似コード
for z:
for y:
for x:
f3[x,y,z] = ...
f3.fuse(x, y, t);
//疑似コード
for z:
for t:
f3[t%X,t/X,z] = ...
この計算 vs. 条件分岐
どちらが低コスト？

15
Schedule : tile
split(y,yo,yi,yfactor)
split(x,xo,xi,xfactor)
reorder(xi,yi,xo,yo)
と等価 f.tile(x,y,xo,yo,xi,yi,32,32);
//疑似コード
for z:
for yo:
for xo:
for yi:
for xi:
f[xo*32+xi, yo*32+yi, z] = ...
f
x
y
z
Note:
着目点の周辺画素参照する場合
(例：フィルタ演算)、タイル化する
とキャッシュがhitしやすくなる。
Note:
xfactor, yfactorの値の
選択 ... 性能に影響す
るので試してみる。
アドレスインクリメント方
向にHW-prefetcherが
効く。

16
Schedule : Producer-Consumer
f3 の計算の前に、f3が必要とするデータをf2が用意する必要がある。
データの受け渡しが発生。
前段の計算タイミングのschedule記述
compute_inline
compute_root
compute_at
f2 f3
x
y
z
Producer Consumer

17
階層メモリシステム
ALU
Register
L1 Cache
L2 Cache
L3 Cache
DRAM
0.3 ns/cycle
0.3 ns
約1 ns
約百 ns ~
約十 ns
約数 ns
32 W
32KB
数GB ~
数MB
256KB
512×512×3 768KB
1K×1K×3 3MB
数K×数K×数ch 数十MB～
画像サイズとデータ量
（8bit/色の場合）
Coreごと
共有リソース
速い?小さい
遅い?大きい

18
Schedule : compute_inline
インライン展開 (デフォルト) ... レジスタ経由の受け渡し
f2 f3
x
y
z
f2.compute_inline();
//疑似コード
for z:
for y:
for x:
f3[x,y,z] = f2の展開式
Note :
? 速い
? メモリを消費しない
ただし、式が複雑すぎるとレジスタが
不足し、中間結果をメモリ（stack）に
Load/Storeする。? 遅くなる。

19
Schedule : compute_root
Func毎に画像?テンソル全体を計算
多くの画像処理ライブラリと同様
メモリ(cache, DRAM)経由の受け渡し
中間データの大きさに依存
(512×512×3ならL3に収まるかな)
f2 f3
x
y
z
f2.compute_root();
//疑似コード
バッファ f2;
for z:
for y:
for x:
f2[x,y,z] = ...
for z:
for y:
for x:
f3[x,y,z] = ...

20
Schedule : compute_at
下流側のFuncの軸のインクリメント時点で必要となるデータを計算
いくつかの例
f2.compute_at(f3,x)
f2.compute_at(f3,y)
f2.compute_at(f3,z)
f3.tile(x,y,xo,yo,xi,yi,32,32); f2.compute_at(f3,xo)
f3.tile(x,y,xo,yo,xi,yi,32,32); f2.compute_at(f3,yo)
f3.split(y,yo,yi,32); f2.compute_at(f3,yo)
f2 f3
x
y
z

21
Schedule : compute_at (基礎)
f2 f3
x
y
z f2.compute_at(f3,x)
f3の1画素(要素)ごとにf2の計算
を行い、結果をメモリに格納。
Cf. cpmpute_inline
f2 f3
x
y
z f2.compute_at(f3,y)
f3の1行ごとにf2の計算を行い、
結果(行)をメモリに格納。
f2 f3
x
y
z f2.compute_at(f3,z)
f3の1面ごとにf2の計算を行い、
結果(面)をメモリに格納。

22
Schedule : compute_at
f2.compute_at(f3,x)
//疑似コード
for z:
for y:
バッファ f2;
for x:
//producer
for z:
for y:
for x:
f2[x,y,z] = ...
//consumer
f3[x,y,z] = ...
f2.compute_at(f3,y)
//疑似コード
for z:
バッファ f2;
for y:
//producer
for z:
for y:
for x:
f2[x,y,z] = ...
//consumer
for x:
f3[x,y,z] = ...
f2.compute_at(f3,z)
//疑似コード
バッファ f2;
for z:
//producer
for z:
for y:
for x:
f2[x,y,z] = ...
//consumer
for y:
for x:
f3[x,y,z] = ...
f3の点に必要な
f2[x,y,z]を計算
f3の線に必要な
f2[x,y,z]を計算
f3の面に必要な
f2[x,y,z]を計算

23
Schedule : compute_at (応用)
f2 f3
x
y
z f3.tile(x,y,xo,yo,xi,yi,32,32)
f2.compute_at(f3,xo)
xoごとにf3のタイルに必要となるf2
の計算を行い、結果をメモリに格納。
f2 f3
x
y
z f3.split(y,yo,yi,32)
f2.compute_at(f3,yo)
f3の32行ごとに必要となるf2の計
算を行い、結果をメモリに格納。
f2 f3
x
y
z f3.tile(x,y,xo,yo,xi,yi,32,32)
f2.compute_at(f3,yo)
yoごとにf3のタイルに必要となるf2
の計算を行い、結果をメモリに格納。
f2 f3
x
y
z f3.split(x,xo,xi,32)
f3.reorder(xi,y,xo)
f2.compute_at(f3,xo)

24
Schedule : reorder_storage
Funcの出力を格納するバッファでの多次元配列の格納順の変更
Cacheは(例えば)32byte境界単位でDRAMをアクセス
Func f;
f(x,y,z) = ...
// デフォルト
// xが内側、zが外側になるよう配置
// Cふうに書くと
// f[Z][Y][X] = ...
Func f;
f(x,y,z) = ...
f.reorder_storage(z,x,y) = ...
// 内側からz, x, yの順番で配置
// Cふうに書くと
// f[Y][X][Z] = ...
Func f;
f(x,y,z) = ...
// 計算順を変えたらメモリ(バッファ)格納順も変えてみる
f.reorder(z,x,y).reorder_storage(z,x,y) = ...

25
Schedule : store_at, store_root
Loop-nestでのバッファの位置の指定
f2.compute_at(f3,x);
f2.store_at(f3, x); // デフォルト
//疑似コード
for z:
for y:
バッファf2; // f3のx
for x:
//producer
for z:
for y:
for x:
f2[x,y,z] = ...
//consumer
f3[x,y,z] = ...
f2.store_root();
//疑似コード
バッファf2; // 一番外側 (root)
for z:
for y:
for x:
//producer
for z:
for y:
for x:
f2[x,y,z] = ...
//consumer
f3[x,y,z] = ...
f2.store_at(f3, y);
//疑似コード
for z:
バッファf2; // f3のy
for y:
for x:
//producer
for z:
for y:
for x:
f2[x,y,z] = ...
//consumer
f3[x,y,z] = ...
Producerのstore位置はProducerのcompute位置と同じか外側

26
Schedule : storage関係
align_storage(Var dim, Expr alignment)
? 軸dimのstorage上のextentをalignmentの倍数にする。
? 注： alignmentの単位はByteではなく要素の個数
? SIMDのload/store時に有利(なハズ)
fold_storage(Var dim, Expr extent, Bool fold_forward=true)
? extentの大きさの循環バッファを使う。
? extentを2のべき乗にするのがよい。
? 同じ場所を再利用するのでcacheに優しい。

27
Schedule 例
f3のタイルを計算する際に
f1でタイルを計算 ? f2でタイルを計算 ? f3でタイルを計算
タイルのサイズを選べば画像が大きくてもcache hitを期待できる。
f1 f2 f3
x
y
z
// zを内側にしてアンロール
f3.reorder(z,x,y);
f3.unroll(z);
// タイル化
f3.tile(x,y,xo,yo,xi,yi,64,64);
f2.compute_at(f3,xo);
f1.compute_at(f3,xo);
// 更にSIMD化したい場合
f1.vectorize(xi);
f2.vectorize(xi);
f3.vectorize(xi);
ループ(条件分
岐)は遅いので
繰り返しが少な
いなら展開する。

28
並列処理の時代
“The Free Lunch is Over.”
by Herb Sutter
http://www.gotw.ca/publications/concurrency-ddj.htm
並列化で頑張る時代
Halideの場合
並列化
f.parallel(軸,factor)
factor単位でtaskを作る。
SIMD化
f.vectorize(軸)
並列化

29
アムダールの法則
並列化は簡単だが、N coreでN倍の性能UPのためには努力が要る。
プログラムにはパラレルにできる部分とできない部分がある。
1 coreでの実行時間を 1
並列化不可能部分 α
並列化可能部分 1 – α
N core利用時
Speed Up = N / (1+α(N-1)) 1
2
4
8
16
32
1 2 4 8 16 32
Core数
SpeedUp
α= 0
α= 0.01
α= 0.1

30
共有リソース
ALU
Register
L1 Cache
L2 Cache
L3 Cache
DRAM
0.3 ns
0.3 ns
約1 ns
約百 ns ~
約十 ns
約数 ns
32 W
32KB
数GB ~
数MB
256KB
Coreごと
共有リソース
bus
L1, L2でmissが発生
missが多発すると...
bus, L3, DRAMでア
クセス競合 ? 待ち。
残念な結果

31
残念な結果になるケース
Halideのparallel(軸,factor)
f1.compute_root().parallel(y,32);
f1 f2 f3
? y軸方向に32行ごとの(Y+1)/32個のworkを作成
? workを環境変数HALIDE_NUM_THREADS個
のworker threadで実行
? 生成する中間データが大きい ? cacheからあふれる。
? cacheに収まっても、f1のある領域を担当したcoreが
f2でその領域を担当するとは限らない ? core間の
local cacheデータ移動(内部bus経由)の発生。
並列実行
Fork/Join
並列実行
Fork/Join
並列実行
Fork/Join

32
共有リソース
演算器
Register
L1 Cache
L2 Cache
L3 Cache
DRAM
演算器
Register
L1 Cache
L2 Cache
演算器
Register
L1 Cache
L2 Cache
演算器
Register
L1 Cache
L2 Cache
共有
リソース
共有
リソース
Core
パラレル実行
アクセス競合すると
シリアル実行
数年前...
- Core i7 約15GB/sec
- Xeron 約25GB/sec
1 coreあたりMAX約
10GB/sec DRAMバ
ンド幅を使う。

33
並列schedule for cache 例
並列処理でも
適切なタイルサイズで
Local cache経由データ受け渡し
f1 f2 f3
x
y
z
// タイル化
f3.tile(x,y,xo,yo,xi,yi,64,64);
f3.fuse(xo,yo,t);
f2.compute_at(f3,t);
f1.compute_at(f3,t);
f3.parallel(t);
//疑似コード
for z:
parallel_for t:
for yi:
for xi:
for z:
for y:
for x:
f1[x,y,z] = ...
for z:
for y:
for x:
f2[x,y,z] = ...
f3[x,y,z] = ...
Task (同一thread)
Task数は多い方
が自動ロードバラ
ンスが効きやすい。
fuseでworkを増や
す。
Thread(論理CPU)
≒ 物理CPU core

35
Auto Scheduler
CMUで開発。今年、Halide本流にマージ済み。
http://graphics.cs.cmu.edu/projects/halidesched/
https://github.com/halide/Halide/tree/master/test/auto_schedule
万能ではないが、面倒な時とか、とりあ
えず試してみる価値アリ。
Func f;
f(x) = ... ;
Target target = ... ;
f.estimate(x, 0, SIZE);
f.auto_schedule(target);
サンプル
コード
auto schedulerに
大きさを伝える
ターゲット
デバイス

36
rfactor
“Parallel Associative Reductions in Halide”
https://andrew.adams.pub/par.pdf
HalideでのReduction計算
RDom r(0,SIZE);
Func f;
f() = 0; // スカラー、0-D テンソル
f() += in(r); // Σ計算、Reduction
// fへの初期値代入部分のschedule
f.some_schedule();
// f() += 部分のschedule
f.update(0).some_schedule();
出典：上記
並列化ダメ
f.split(r, rx, ry, 4);
// 並列化する場合
Var v;
f.rfactor(ry, v).parallel(v);
// ベクトル化する場合
Var u;
f.rfactor(rx, u).parallel(u);
内側

37
自動微分
“Differentiable Image Processing in Halide”
https://people.csail.mit.edu/tzumao/gradient_halide/
https://people.csail.mit.edu/tzumao/gradient_halide/gradient_halide.pdf
自動微分の導入とbackwardプログラムの自動生成
出典：上記

39
Convolution計算1
出典 http://www.rle.mit.edu/eems/wp-content/uploads/2018/06/2018_vlsi_forum.pdf
画像処理のConvolutionと同じ

40
Convolution計算2
入力のchannel数がCの場合
Filterのchannel数もC 更に、出力のchannel数をMとする
FilterをMセット用意

41
Convolution計算3
N batch size 1~
C, M 3, 16, 32, 64, 128, 256, ...
H over 200, ..., 1
R 1, 3, 5, ...
s stride, 1, 2, ...
outputの計算
x,y,m,n,u,v,c の7重ループ
Batch処理
入力をNセット
出力もNセット

42
Convolution計算別のアルゴリズム
? im2col データを並び替えてGEMM(行列乗算)
? Winograd 2倍ぐらい早くなる
? etc.
Idein 中村さん
「Convolutionの数理とアルゴリズム」
https://speakerdeck.com/nineties/convolutionfalseshu-li-toarugorizumu

43
Deep Neural Network
ResNet
出典 https://medium.com/@siddharthdas_32104/cnns-architectures-lenet-alexnet-vgg-googlenet-resnet-and-more-666091488df5
- Convolution各種 (積和算)
- 要素ごとの加算
- ReLU ( max(x,0) )
- Full Combination (積和算)
- etc.
各中間data ～数百K語
各conv.係数～数M語
参考 https://dgschwend.github.io/netscope/quickstart.html
ScheduleでCacheをうまく使いたい
Convolution等の深いパイプライン

45
HalideとDeep Neural Network
Halide MIT
CMU
Auto Scheduling
Google
rfactor
TVM
Tensor Comprehensions
Tiramisu
Facebook
MIT
Univ. of Washington
Differentiable
Programming
Reduction計算部
(例:Σ)の並列化
自動微分の導入
Deep Neural
Network等の計算
Auto Tuning
進化的
アルゴリズム
Polyhedral
Optimization

46
Tensor Comprehensions
https://research.fb.com/announcing-tensor-comprehensions/
Fixstars 飯塚さん
TVM
https://tvm.ai/
Halide-IR(中間表現)の利用
Tiramisu
https://github.com/Tiramisu-Compiler/tiramisu
@Vengieerさん

47
Domain Specific XXXの時代
Mooreの法則はもうおしまい。(かな?)
Single-core CPU Multi-core CPU
GP-GPU
専用回路
Domain Specific Architecture
? Google TPU
? NVIDIA TensroRT, NVDLA
ISA互換、
HWで頑張る
でも、最後は
結局メモリ
小さいcoreと
threadを沢山。
CUDAでSW
をがんばろう。
HW的詳細を
知らないと最
適化困難かも
Domain
Specific
Language
コンパイラをいじれることが重要になりそう
でも、最後は
結局メモリ

49
まとめ
? HalideのAlgorithm記述 = 計算グラフ
? データ?計算の依存性
? Scheduleでテンソル計算順?メモリ格納順を操作
? Func間のデータの渡し方、Cache/Memoryをどう使うか。
? Multi-core時代、Cacheの利用
? 新機能 Auto Scheduler, Reduction内の並列化, 自動微分
? Deep Neural Network
? Halide-IRとpost Moore時代のDomain Specific XXX
? コンパイラ

狠狠撸

Halide for Memory

More Related Content

Halide for Memory