$30 off During Our Annual Pro Sale. View Details »

スレッド並列の基礎 for 数学と物理におけるJuliaの活用 2023-07-10

スレッド並列の基礎 for 数学と物理におけるJuliaの活用 2023-07-10

『数学と物理におけるJuliaの活用』研究会 チュートリアル講演(3) 講演資料スライド(2023/07/10)

GOTOH Shunsuke

July 11, 2023
Tweet

More Decks by GOTOH Shunsuke

Other Decks in Technology

Transcript

  1. 講義の目的と目標 • スレッド並列(の基本)の理解 ◦ 並列処理とは何か ◦ スレッドの基本的な概念 ◦ スレッドの利点と制限 ◦

    (スレッドのライフサイクル) • Juliaによるスレッド並列の理解 ◦ 最低限の事前準備 ◦ 各種API ◦ スレッドセーフ(とデッドロック)
  2. 実験コード • 本講演資料に記載のサンプルコードを 以下のURLで公開しています。 https://antimon2.github.io/julia_imi_workshop2023_tutorial/1st_examp le_nthreads.html • チュートリアル(1) で Pluto.jl

    を導入済の方は、上記URL→[Edit or Run] ボタン→「On your computer」の手順に従って手元の環境で実行できます。 • 実習課題も編集可能な状態で用意してありますので ご利用ください。
  3. 並列処理とは何か(3) アーキテクチャ的な話… • 一連の処理(タスク) を同時に走らせて処 理をすること • 目的は スループット 向上(パフォーマンス向上)

    • いくつか種類がある: ◦ スレッド並列 ◦ プロセス並列 julia> Threads.nthreads() 6 julia> fib(n) = n ≤ 1 ? n : fib(n - 2) + fib(n - 1) fib (generic function with 1 method) julia> @time fib(40) 0.510351 seconds 102334155 julia> @time for _=1:4 fib(40) end 2.073198 seconds (8.14 k allocations: 600.650 KiB, 1.68% compilation time) julia> @time Threads.@threads for _=1:4 fib(40) end 0.783171 seconds (34.33 k allocations: 2.571 MiB, 19.96% compilation time)
  4. 補足:『真のスレッド並列』 時間があれば(なければ飛ばすので読んでおいてください)… • 他言語にもスレッド(マルチスレッド)という機能はある • 一部言語のスレッドは、GIL(Global Interpreter Lock)という機構により 並列動作が制限されている •

    Julia にはその制限がない(『真のスレッド並列』) • 他言語で「並列処理するならマルチプロセス」として提案・提唱されている事案 でも、Julia ならスレッド並列で間に合う(むしろそちらの方が手っ取り早い)も のも多い • Julia では「オーバーヘッドによるデメリットよりも『スレッドアンセーフやデッド ロックの心配をしなくて済む』というメリットの方が勝る」ような場合には プロセ ス並列 も選択出来る(選択肢が多い)と建設的に考えることもできる
  5. 補足:スレッドのライフサイクル さらに概念的かつ内部仕様的な話(時間がなければ飛ばすので読んでおいてください)… • 一般的なスレッドの状態(ステージ):「1. 新規(New)」「2. 実行可能 (Runnable)」「3. 実行中(Running)」「4. ブロック(Blocked)」「5. 終了

    (Terminated)」 • 低レベルでスレッドを扱う必要のある環境では、これらの状態(ステージ)を適切 に把握・管理する必要がある • Julia ではあまり深く考える必要はない • 必要になったときに Julia が内部でスレッドを立ち上げ(「1. New」)タスクを 走らせて(「2. Runnable」「3. Running」)適宜協調動作(「4. Blocked」)し ながらタスクが終了したら終了(「5. Terminate」)してくれる • のでこの状態(ステージ)は、概念として参考までに
  6. 最低限の事前準備 Julia 起動時に最大スレッ ド数を指定する • julia -t 4 (絶対数指定) •

    julia -t auto (利用できるCPUコア最 大確保) • 環境変数による指定も可 能 • Pluto の場合デフォルト で「利用できる最大÷2」 が設定される > julia -t 4 # 最大スレッド数は4 > julia -t auto # 確保できる分だけスレッド数を確保 > env JULIA_NUM_THREADS=XXX julia # Bash/Zsh などで環境変数指定と同時にJulia起動 # Julia 内では ↓ で利用可能なスレッド数が取得できる julia> Threads.nthreads() 6
  7. 各種API(1) @threads • スレッド並列の基本 • for ~ end を簡単 にマルチスレッド化で

    きる • 本講ではモジュール プレフィックス付きで Threads.@threads という表記を利用(以 下同様) julia> Threads.nthreads() 6 julia> fib(n) = n ≤ 1 ? n : fib(n - 2) + fib(n - 1) fib (generic function with 1 method) julia> @time fib(40) 0.510351 seconds 102334155 julia> @time for _=1:4 fib(40) end 2.073198 seconds (8.14 k allocations: 600.650 KiB, 1.68% compilation time) julia> @time Threads.@threads for _=1:4 fib(40) end 0.783171 seconds (34.33 k allocations: 2.571 MiB, 19.96% compilation time)
  8. 各種API(2) @spawn • スレッドの起動 (その都度本当に新しいスレッド が起動するわけではない、詳細 割愛) • 続く式を1つのタスク としてまとめ、(空い

    ている)スレッドで処 理する • @sync マクロで同 期を取る必要あり julia> Threads.nthreads() 6 julia> @time Threads.@threads for _=1:4 fib(40) end 0.783171 seconds (34.33 k allocations: 2.571 MiB, 19.96% compilation time) julia> @time @sync for _=1:4 Threads.@spawn fib(40) end 0.715108 seconds (2.01 k allocations: 135.894 KiB, 8.50% compilation time)
  9. 各種API(3) @threads v.s. @spawn (1) • 大量の細かい(粒度が ほぼ揃っている)タス クを並列処理したい 場合は

    @threads の方が優秀 julia> let hist = zeros(Int, Threads.nthreads()) @time Threads.@threads for _=1:(10000*Threads.nthreads()) fib(15) # 毎回ほぼ同じ負荷の処理が実行される例 hist[Threads.threadid()] += 1 end hist end 0.134166 seconds (22.88 k allocations: 1.588 MiB, 169.92% compilation time) #> [10000, 10000, 10000, 10000, 10000, 10000] julia> let hist = zeros(Int, Threads.nthreads()) @time @sync for _=1:(10000*Threads.nthreads()) Threads.@spawn begin fib(15) # 毎回ほぼ同じ負荷の処理が実行される例 hist[Threads.threadid()] += 1 end end hist end 0.301360 seconds (569.85 k allocations: 35.267 MiB, 102.04% compilation time) #> [0, 12151, 12635, 11490, 11973, 11751]
  10. 各種API(4) @threads v.s. @spawn (2) • タスクの数がそこまで 大量でなく、粒度がバ ラバラの場合は、 @spawnの方が高パ

    フォーマンスとなるこ ともある (負荷分散) julia> let hist = zeros(Int, Threads.nthreads()) @time Threads.@threads for _=1:(250*Threads.nthreads()) fib(rand(25:35)) # 毎回わずかに異なる負荷の処理が実行 hist[Threads.threadid()] += 1 end hist end 4.370164 seconds (34.68 k allocations: 2.347 MiB, 6.02% compilation time) #> [250, 250, 250, 250, 250, 250] julia> let hist = zeros(Int, Threads.nthreads()) @time @sync for _=1:(250*Threads.nthreads()) Threads.@spawn begin fib(rand(25:35)) # 毎回わずかに異なる負荷の処理が実行 hist[Threads.threadid()] += 1 end end hist end 4.310160 seconds (19.18 k allocations: 1.329 MiB, 1.33% compilation time) #> [245, 244, 224, 283, 271, 233]
  11. 各種API(5) その他 • Channel(spawn=true) do ~ end ◦ 特定の処理だけをラップして別スレッドで処理できる ◦

    適宜値の受け渡しができる (元々Channelはタスク間で値を受け渡ししながら協調動作するためのもの) • Threads.foreach() ◦ Base.foreach() のマルチスレッド版 ◦ 処理(関数)とイテレータ(Channelでラップしたもの)を受け取りマルチスレッドで処理 ◦ 負荷分散に有用 • ※これらの詳細は割愛
  12. 実例(1) @threads の例(1) • 配列の異なるインデッ クスに同様の演算を 実施して値を格納す るパターン(このパ ターンが多い) •

    @threads は入れ 子にできる(Julia v1.8 以降) julia> function my_matmul(A::AbstractMatrix, B::AbstractMatrix) T = promote_type(eltype(A), eltype(B)) C = Matrix{T}(undef, (size(A, 1), size(B, 2))) Threads.@threads for x = axes(B, 2) Threads.@threads for y = axes(A, 1) C[y, x] = @view(A[y, :])' * @view(B[:, x]) end end C end my_matmul (generic function with 1 method) julia> let A=[1 2; 3 4; 5 6; 7 8], B=[1 2 3; 4 5 6] C = my_matmul(A, B) @assert C == A * B @show C; end C = [9 12 15; 19 26 33; 29 40 51; 39 54 69]
  13. 実例(2) @threads の例(2) • 前スライドで示したサ ンプルのベンチマーク • (※なお普通に A *

    B を計算した方がメモリ 使用量も少ないし断 然速い(あくまで参 考)) julia> function my_matmul_st(A::AbstractMatrix, B::AbstractMatrix) T = promote_type(eltype(A), eltype(B)) C = Matrix{T}(undef, (size(A, 1), size(B, 2))) for x = axes(B, 2) for y = axes(A, 1) C[y, x] = @view(A[y, :])' * @view(B[:, x]) end end C end # シングルスレッド版 my_matmul_st (generic function with 1 method) julia> using BenchmarkTools julia> @btime my_matmul_st(A, B) setup=(A=rand(100, 100); B=rand(100, 100)); 374.714 μs (2 allocations: 78.17 KiB) julia> @btime my_matmul(A, B) setup=(A=rand(100, 100); B=rand(100, 100)); 329.578 μs (3239 allocations: 419.11 KiB) julia> # 参考 @btime (A * B) setup=(A=rand(100, 100); B=rand(100, 100)); 48.012 μs (2 allocations: 78.17 KiB)
  14. 実例(3) @spawn の例 • 標準の map 関数を マルチスレッド化して みた •

    fetch() は、同期を 取って結果の値を取 得する関数 julia> function threaded_map(fn, array::AbstractArray) tasks = [Threads.@spawn(fn(v)) for v in array] [fetch(task) for task in tasks] end threaded_map (generic function with 1 method) julia> threaded_map(fib, 15:40) == map(fib, 15:40) true julia> using BenchmarkTools julia> @btime map(fib, 15:40); 1.419 s (1 allocation: 272 bytes) julia> @btime threaded_map(fib, 15:40); 652.908 ms (167 allocations: 13.77 KiB)
  15. 注意(1) スレッドセーフについて(1) • なんでもマルチスレッ ド化して期待通りに動 作するわけではない • 例:同じ変数(オブ ジェクト)の内容を更 新する場合など

    • ⇒スレッドアンセーフ julia> mutable struct UnsafeCounter count::Int UnsafeCounter() = new(0) end julia> begin counter1 = UnsafeCounter() for n=1:1000 counter1.count += 1 end counter1.count end 1000 julia> begin counter2 = UnsafeCounter() Threads.@threads for n=1:1000 counter2.count += 1 end counter2.count end 650 # 想定より少ない値になってしまう(ことがある)
  16. 注意(2) スレッドセーフについて(2) • 更新するフィールドに @atomic を付ける • そのフィールドを更新 するときも @atomic

    を付ける • ⇒スレッドセーフにな る julia> mutable struct AtomicCounter Threads.@atomic count::Int AtomicCounter() = new(0) end julia> begin counter3 = AtomicCounter() Threads.@threads for n=1:1000 Threads.@atomic counter3.count += 1 end counter3.count end 1000 # 期待通りの結果になる julia> begin counter4 = AtomicCounter() for n=1:1000 # ↓シングルスレッドでも `Threads.@atomic` 必須 Threads.@atomic counter4.count += 1 end counter4.count end 1000
  17. 注意(3) その他 • もっと広い範囲をスレッドセーフにしたい ⇒ロック機構(@lock 《lockオブジェクト》 begin ~ end)による 排他制御

    を導入する必要あり • デッドロック(竦みの状態になってどのスレッドもブロック状態(Blocked)に なってしまう)にも注意! • ※これらの詳細は割愛
  18. お題:N-Queen問題(2) N-Queen問題を解く(1) • Julia のコード例を示す • その1(前半): メイン部および再帰部 function nQueen(n::Int)

    counter = AtomicCounter() # 先ほど作ったカウンタを利用 for y=1:n Threads.@atomic counter.count += nQueen_sub(n, 1, [y]) end counter.count end function nQueen_sub(n::Int, k::Int, params::Vector{Int}) k == n && return 1 board = create_board(n, params) counter = AtomicCounter() for y=1:n if issafe(board, y, k + 1) Threads.@atomic counter.count += nQueen_sub(n, k+1, [params; y]) end end counter.count end
  19. お題:N-Queen問題(3) N-Queen問題を解く(2) • Julia のコード例を示す • その2(後半): 盤面を表す BitMatrix を生成する関数と、その盤

    面に指定位置にコマを置 いて安全かを判定する関 数 function create_board(n, params) board = falses(n, n) for (x, y) in pairs(params) board[y, x] = true end board end function issafe(board::BitMatrix, y, x) h, w = size(board) any(board[y, :]) && return false any(board[:, x]) && return false any(board[y-i, x-i] for i=1:min(y, x)-1) && return false any(board[y+i, x+i] for i=1:min(h-y, w-x)) && return false any(board[y-i, x+i] for i=1:min(y-1, w-x)) && return false any(board[y+i, x-i] for i=1:min(h-y, x-1)) && return false true end
  20. お題:N-Queen問題(4) N-Queen問題を解く(3) • Julia のコード例を示す • その3: 実行例(N=8~12まで) • ※算出した解は所謂

    バラエ ティ解(回転・反転で同一視 できる物も個々に列挙)で あり、基本解(回転・反転で 同一視できる物を1つと数 える方法)ではないことにも 注意 julia> using BenchmarkTools julia> @btime nQueen(8) 3.699 ms (75032 allocations: 2.95 MiB) 92 julia> @btime nQueen(9) 16.535 ms (326209 allocations: 13.11 MiB) 352 julia> @btime nQueen(10) 82.212 ms (1483298 allocations: 60.40 MiB) 724 julia> @btime nQueen(11) 436.700 ms (7371379 allocations: 303.26 MiB) 2680 julia> @btime nQueen(12) 2.611 s (39709620 allocations: 1.61 GiB) 14200
  21. お題:N-Queen問題(5) N-Queen解法の並列化(1) • 前半のコード (nQueen() 関数、 nQueen_sub() 関数) をマルチスレッド化(ス レッド並列化)してみよ

    う! • ヒント:今日覚えたAPI (@threads、@spawn など)を使えばOK # ↓マルチスレッド化してみよう! function nQueenMT(n::Int) counter = AtomicCounter() # 先ほど作ったカウンタを利用 for y=1:n Threads.@atomic counter.count += nQueen_subMT(n, 1, [y]) end counter.count end # ↓マルチスレッド化してみよう! # ↓適宜編集して最適な挙動にしてみよう! function nQueen_subMT(n::Int, k::Int, params::Vector{Int}) k == n && return 1 board = create_board(n, params) counter = AtomicCounter() for y=1:n if issafe(board, y, k + 1) Threads.@atomic counter.count += nQueen_sub(n, k+1, [params; y]) # Threads.@atomic counter.count += nQueen_subMT(n, k+1, [params; y]) end end counter.count end
  22. お題:N-Queen問題(6) N-Queen解法の並列化(2) • 実行(N=8~12まで)し てみて、パフォーマンス向 上しているか確認しよ う! (だいたい1/2~1/3くら いの処理時間になるはず) •

    ※うまく並列化できていな いと逆にパフォーマンス低 下(2~10倍の処理時間)に なることもあります、ご注 意! julia> using BenchmarkTools julia> @btime nQueenMT(8) ?.??? ms (????? allocations: ?.?? MiB) 92 julia> @btime nQueenMT(9) ?.??? ms (?????? allocations: ??.?? MiB) 352 julia> @btime nQueenMT(10) ??.??? ms (??????? allocations: ??.?? MiB) 724 julia> @btime nQueenMT(11) ???.??? ms (??????? allocations: ???.?? MiB) 2680 julia> @btime nQueenMT(12) ?.??? s (???????? allocations: ?.?? GiB) 14200
  23. 補足:アムダールの法則 • スレッド数 N で並列化しても、処理時間は 1/N とはならな い • 全体の

    90% が並列化できるとした場合、パフォーマンス は最大でも 1/((1-0.9)+0.9/N) 倍 (N→∞ で 10倍)に しかならない • 「いかに並列化不可能な部分を少なくするか」が課題 • 参考: Wikipedia[ja] アムダールの法則 等
  24. 講義の目的と目標 ✓ スレッド並列(の基本)の理解 ☑ 並列処理とは何か ☑ スレッドの基本的な概念 ☑ スレッドの利点と制限 ☑

    (スレッドのライフサイクル) ✓ Juliaによるスレッド並列の理解 ☑ 最低限の事前準備 ☑ 各種API ☑ スレッドセーフ(とデッドロック)