Upgrade to Pro — share decks privately, control downloads, hide ads and more …

プログラムを高速化する話

 プログラムを高速化する話

プログラムを高速化するためのテクニックをまとめました。GPU編はこちら https://speakerdeck.com/primenumber/puroguramuwogao-su-hua-suruhua-ii-gpgpubian

prime number

July 10, 2023
Tweet

More Decks by prime number

Other Decks in Programming

Transcript

  1. プログラムを高速化する話
    @KMC 春合宿 2015
    KMC2 回生 prime

    View full-size slide

  2. 2
    自己紹介
    第 37 代 KMC 会計( 2014/02 〜 2015/01 )
    root ( 2013/11 〜)
    理学部数学系志望 ( もうすぐ系登録試験 )
    競技プログラミング練習会
    GR-SAKURA で遊ぼう
    コーディング大海
    KMC 学習発表会

    View full-size slide

  3. 3
    自己紹介
    過去の KMC での活動
    難解プログラミング言語勉強会
    オセロの AI を作ろう
    2048AI コンテスト
    RICOH and Java Developer Challange Plus 2013/2014
    KMC 流しそうめん
    KMC 静止そうめん
    KMC 鍋
    prime/ 免許合宿

    View full-size slide

  4. 4
    目次
    1.はじめに / 最適化について
    2.キャッシュを考慮した最適化
    3.ビット演算による高速化
    4.SIMD 命令による高速化

    View full-size slide

  5. 5
    はじめに
    現代ではそこまで頑張ってプログラムを高速化する必要
    は無くなってきている
    コンピュータの性能向上・省電力化
    既存の高速なライブラリの充実
    クラウドコンピューティングの普及

    View full-size slide

  6. 6
    はじめに
    ほとんどの場合、既存のライブラリを用いたり、アルゴ
    リズムを改良したりすることによって、必要なパフォー
    マンスを得ることが出来る
    それで不十分でもクラウドコンピューティングによって金
    の力で殴れば何とかなることが多い

    View full-size slide

  7. 7
    はじめに
    しかし、どうしてもプログラムを限界まで高速化したい
    場面も存在する
    利用できる計算資源が限られているとき
    消費電力を減らしたいとき
    時間制限のあるゲームの AI など、高速化するほど有利なと

    莫大な計算量がかかるとき

    View full-size slide

  8. 8
    はじめに
    今後は「最適化」と「高速化」を同じ意味で使います
    この講座では高速であることが正義です、最適です
    今回は特に高い効果の得やすい 3 つの高速化技法につい
    て扱います
    キャッシュを意識したプログラミング
    ビット演算
    SIMD 命令

    View full-size slide

  9. 9
    今回の目標
    コンパイラ・ライブラリ・ CPU 等の高速化技術を上手く
    活用して楽に高速化する
    アルゴリズムの改良をしたり既存ライブラリ等を使った
    りしても、必要なパフォーマンスが得られない時に、最
    後の手段として、手作業でプログラムを高速化するため
    に必要なテクニックを学ぶ

    View full-size slide

  10. 10
    最適化について
    手作業の泥臭い最適化をする前に、もっと高速なアルゴ
    リズムがないかを考えるべき
    オーダーレベルで計算量が減れば、圧倒的に速くなる事
    が多い
    競技プログラミング練習会に参加しましょう!!!
    以後アルゴリズム的には最適(もしくは最適に近い)も
    のを利用していると仮定することにします

    View full-size slide

  11. 11
    最適化について
    「細かい効率のことは忘れて、時間の 97% について考え
    よう。時期尚早な最適化は諸悪の根源だ。それでも残り
    3% についても機会を逃すべきではない」
    - Donald E. Knuth
    「プログラム最適化の第一法則 : 最適化するな。
    プログラム最適化の第二法則 ( 上級者限定 ): まだするな。

    - Michael A. Jackson

    View full-size slide

  12. 12
    最適化について
    最適化は、コードを複雑にすることが多いので、コード
    の変更やデバッグを困難にする
    そのうえ、パフォーマンスに重大な影響を与えるコード
    は全体のうちのほんの僅か
    パフォーマンスに大きな影響を与えないコードを最適化
    してもほとんど意味がない

    View full-size slide

  13. 13
    最適化の対象
    主に Intel の Haswell マイクロアーキテクチャ以降を対象
    多くのテクニックは他のプロセッサにも応用できます
    ベース マイクロアーキテクチャ プロセスルール 登場年
    Nehalem Nehalem 45nm 2008
    〃 Westmere 32nm 2010
    Sandy Bridge Sandy Bridge 32nm 2011
    〃 Ivy Bridge 22nm 2012
    Haswell Haswell 22nm 2013
    〃 Broadwell 14nm 2014

    View full-size slide

  14. 14
    開発環境
    C++ を使用することを仮定します
    低レベルな処理もサポート
    高性能なコンパイラが存在
    メモリ管理をプログラム側から制御可能
    C++ のコンパイラは GCC, Clang 等が有名

    View full-size slide

  15. 15
    Intel C/C++ Compiler
    Intel 謹製の C/C++ コンパイラ
    最適化能力が強いと言われている
    学生なら非営利目的に限り無料で入手可能
    Intel Parallel Studio XE についてくる
    最新版は gcc 4.9 互換なので C++ の最新規格のサポートも
    ぼちぼち

    View full-size slide

  16. 16
    最適化手法を学ぶには
    最適化手法を学ぶにはどうすればよいか?

    View full-size slide

  17. 17
    最適化手法を学ぶには
    「インテル ® 64 アーキテクチャー および IA-32 アーキ
    テクチャー 最適化リフ レンス マニュアル」 [0][1] を読め
    (完)

    View full-size slide

  18. 18
    最適化手法を学ぶには
    「インテル ® 64 アーキテクチャー および IA-32 アーキ
    テクチャー 最適化リフ レンス マニュアル」 [0][1] を読め
    (完)
    PDF ファイルをインテルの公式サイトからダウンロード
    できます
    日本語訳もあるけどちょっと古い
    700 ページ以上の超大作

    View full-size slide

  19. 19
    用語の説明
    CPU
    Central Processing Unit の略
    コンピュータなどにおいて中心的な処理装置として働く電
    子回路 (Wikipedia より )
    レジスタ
    CPU の中にある計算等に用いる容量の小さな記憶装置

    View full-size slide

  20. 20
    手作業による最適化
    キャッシュを意識したプログラム
    ビット演算の活用
    SIMD 命令の活用
    その他の最適化

    View full-size slide

  21. 21
    キャッシュを意識したプログラミング
    メモリアクセスの遅延とキャッシュ
    局所的でないメモリアクセスを避ける
    データ構造を SoA に
    ストリップマイニング
    ブロック化

    View full-size slide

  22. 22
    メモリアクセスの遅延
    メモリへのランダムアクセスは数十〜数百 cycle かかる
    レイテンシが大きい
    レジスタ上だけで計算が終わる場合に比べて単純計算で数
    十倍のレベルで遅くなる
    メモリアクセスが処理の律速になりやすい

    View full-size slide

  23. 23
    メモリアクセスの遅延
    転送速度も足りない
    メモリ側は DDR4-4266 でも 34.1GB/s
    CPU 側は最大で読み込み 64Bytes/cycle 、書き込み
    32Bytes/cycle なので、 CPU の動作周波数が 3GHz なら、
    読み込み 192GB/s 、書き込み 32Bytes/cycle
    線形にアクセスする場合でも速度が足りない!

    View full-size slide

  24. 24
    キャッシュメモリの導入
    メモリのうち、頻繁にアクセスする場所を高速でアクセ
    スできる場所に保持しておく
    頻度に応じて複数段階のキャッシュが設けられる
    最近では 1 次〜 3 次までの 3 段階
    4 次キャッシュを搭載したものもある

    View full-size slide

  25. 25
    キャッシュ階層
    1次キャッシュ
    2次キャッシュ
    3次キャッシュ
    メインメモリ
    4〜
    遅延(サイクル数)
    11〜
    20〜
    40
    高速
    大容量
    レジスタ
    数十

    View full-size slide

  26. 26
    局所的でないメモリアクセスを避ける
    局所的でないメモリアクセスをすると、アクセスする
    データがキャッシュに乗っている確率が低くなる
    例 ) 行列積
    for (int i = 0; i < ROWS; ++i)
    for (int j = 0; j < COLS; ++j)
    for (int k = 0; k < LEN; ++k)
    C[i][j] += A[i][k] * B[k][j];
    二次元配列 B に ROWS 要素飛びでアクセスしている

    View full-size slide

  27. 27
    局所的でないメモリアクセスを避ける
    局所的でないメモリアクセスをすると、アクセスする
    データがキャッシュに乗っている確率が低くなる
    例 ) 行列積
    for (int i = 0; i < ROWS; ++i)
    for (int k = 0; k < LEN; ++k)
    for (int j = 0; j < COLS; ++j)
    C[i][j] += A[i][k] * B[k][j];
    全ての配列に順番にアクセスするようになった
    入れ替えた

    View full-size slide

  28. 28
    データ構造を SoA に
    大量のデータを順番に処理するとき、
    AoS(Array of Structs; 構造体の配列 ) よりも、
    SoA(Struct of Arrays; 配列の構造体 ) の方が高速に動作
    する可能性がある

    View full-size slide

  29. 29
    データ構造を SoA に

    struct data {
    int a, b, c;
    double x, y, z;
    } d_ary[SIZE]; // AoS
    int a[SIZE], b[SIZE], c[SIZE];
    double x[SIZE], y[SIZE], z[SIZE]; //SoA

    View full-size slide

  30. 30
    データ構造を SoA に
    a0 b0 c0 x0 y0 z0 a1 b1 c1 x1 y1 z1 a2 b2 c2 …
    SoAだと順番にaにアクセスすると6要素ごとになる
    AoSだと順番にaにアクセスすると
    連続した領域にアクセスできる
    SIMD命令を使いやすくなる
    a0 a1 a2 …
    b0 b1 b2 …
    c0 c1 c2 …
    x0 x1 x2 …
    y0 y1 y2 …
    z0 z1 z2 …

    View full-size slide

  31. 31
    SoA のデメリットと対策
    多数の要素を読み取って計算しなければならない場合、
    キャッシュラインを使い尽くしてしまい、逆に遅くなる
    場合もある
    SoA と AoS の適切なハイブリッド構造にすることが必要
    頻繁に同時にアクセスする要素を一つの構造体に、など

    View full-size slide

  32. 32
    ストリップマイニング
    以下のコードを考える
    for (int i = 0; i < SIZE; ++i) {
    hoge(A[i]);
    }
    for (int i = 0; i < SIZE; ++i) {
    fuga(A[i]);
    }
    配列 A が十分長いとき、最初のループが終わった時点で、
    A の先頭はキャッシュから排出されている

    View full-size slide

  33. 33
    ストリップマイニング
    したがって、 A は 2 回メインメモリから読み込まれるこ
    とになり、効率が悪い
    一方で、 hoge や fuga は SIMD 命令で並列に処理できる
    ものとする
    ここで、ループをキャッシュに乗るサイズに分割すると、
    SIMD 命令を使いつつ、メインメモリへのアクセスを 1 回
    に減らせる

    View full-size slide

  34. 34
    ストリップマイニング
    for (int i = 0; i < SIZE; i += strip_size) {
    for (int j = i; j < min(SIZE, i+strip_size); ++j) {
    hoge(A[j]);
    }
    for (int j = i; j < min(SIZE, i+strip_size); ++j) {
    fuga(A[j]);
    }
    }

    View full-size slide

  35. 35
    ブロック化
    次のコードを考える
    for (int i = 0; i < SIZE; ++i) {
    for (int j = 0; j < SIZE; ++j) {
    A[i][j] += B[j][i];
    }
    }
    B には飛び飛びのアクセスをしているので SIZE が大きい
    と毎回キャッシュミスが発生して効率が悪い

    View full-size slide

  36. 36
    ブロック化
    ここで、ストリップマイニングの要領でループをキャッ
    シュに乗るサイズに分ける
    Aへのアクセス Bへのアクセス

    View full-size slide

  37. 37
    ブロック化
    ここで、ストリップマイニングの要領でループをキャッ
    シュに乗るサイズに分けると、キャッシュミスが減る
    Aへのアクセス Bへのアクセス

    View full-size slide

  38. 38
    ブロック化
    for (int i = 0; i < SIZE; i+=block_size) {
    for (int j = 0; j < SIZE; j+=block_size) {
    for (int ii = i; ii < i+block_size; ++ii) {
    for(int jj = j; jj < j+block_size; ++jj) {
    A[ii][jj] += B[jj][ii];
    }
    }
    }
    }

    View full-size slide

  39. 39
    ビット演算による最適化
    基礎知識
    ビット演算とは
    基本なビット演算
    ビット列の基本操作
    ビット演算テクニック集
    ビット演算関連の CPU 命令

    View full-size slide

  40. 40
    基礎知識
    コンピュータ内部ではデータは 2 進数で管理されている
    2 進数リテラルは C++ では 0b のあとに続けて 0/1 を書く
    ことで記述できる
    1 0 1 1 1 0 0 1
    1バイト = 8ビット
    0b10111001

    View full-size slide

  41. 41
    基礎知識
    符号なし整数は単純に 2 進数で表される
    1 バイトなら 0 〜 255=2 -1

    4 バイトなら 0 〜 4294967295=2³²-1
    n バイトなら 0 〜 2^8n-1
    0b10111001=185

    View full-size slide

  42. 42
    基礎知識
    符号付き整数は一番上の桁が符号を表す ( 負の数なら 1)
    n ビット符号付き整数は、正の数なら符号なしと同じ
    負の数なら、 mod 2^n で同じになる正の数の 2 進数表記
    0b10111001 を 1 バイト符号付き整数だとみなすと、
    0b10111001+0b01000111=2⁸ なので、
    0b10111001=-0b01000111=-71

    View full-size slide

  43. 43
    基礎知識
    符号付き整数は一番上の桁が符号を表す ( 負の数なら 1)
    n ビット符号付き整数は、正の数なら符号なしと同じ
    負の数なら、 mod 2^n で同じになる正の数の 2 進数表記
    1 バイトなら -128 〜 127
    4 バイトなら -2147483648 〜 2147483647
    n バイトなら -2^(n-1) 〜 2^(n-1)-1

    View full-size slide

  44. 44
    基礎知識
    実数は 1.x×2^i の形の有理数に丸めて表現する
    丸める精度に応じてバリエーションがある
    単精度浮動小数型
    4 バイト、 2 進数で 23 桁の精度、 ±10^±38 程度まで表現
    できる
    倍精度浮動小数型
    8 バイト、 2 進数で 53 桁の精度、 ±10^±308 程度まで表
    現できる

    View full-size slide

  45. 45
    ビット演算
    2 進数の 0/1 の列を操作するような演算の総称
    ビット論理和
    ビット論理積
    ビット排他的論理和
    ビット否定
    ビットシフト
    加減乗算などがビット列を操作するのに使われることも
    最近はビット列操作用の命令も追加されている

    View full-size slide

  46. 46
    ビット演算のメリット
    ビット演算自体が高速
    回路が単純なので、多くのプロセッサで高速に動作する
    一回の演算でビット幅分を一度に処理できるので高速
    64bit なら 64 個の 0/1 を一括で処理できる
    ビット単位にデータを詰めることで、メモリ使用量が減
    り、キャッシュヒット率が向上する

    View full-size slide

  47. 47
    ビット配列
    整数等の配列を 0/1 の配列として利用するテクニック
    符号なし 64bit 整数 (uint64_t) の配列が一番扱いやすい
    SIMD 命令とも相性がよく、高速化が期待できる

    View full-size slide

  48. 48
    基本的なビット演算
    ビット論理和 OR (C 言語 : |)
    片方でも 1 ならば 1 、そうでないなら 0
    0 0 1 0 1 1 0 0
    1 0 1 1 1 0 0 1
    1 0 1 1 1 1 0 1
    A
    B
    A OR B

    View full-size slide

  49. 49
    基本的なビット演算
    ビット論理積 AND (C 言語 : &)
    両方とも 1 ならば 1 、そうでないなら 0
    0 0 1 0 1 1 0 0
    1 0 1 1 1 0 0 1
    0 0 1 0 1 0 0 0
    A
    B
    A AND B

    View full-size slide

  50. 50
    基本的なビット演算
    ビット排他的論理和 XOR (C 言語 : ^)
    片方だけ 1 ならば 1 、そうでないなら 0
    0 0 1 0 1 1 0 0
    1 0 1 1 1 0 0 1
    1 0 0 1 0 1 0 1
    A
    B
    A XOR B

    View full-size slide

  51. 51
    基本的なビット演算
    ビット否定 NOT (C 言語 : ~)
    0/1 を反転する
    0 1 0 0 0 1 1 0
    1 0 1 1 1 0 0 1
    A
    NOT A

    View full-size slide

  52. 52
    基本的なビット演算
    ビットシフト (C 言語 : 左シフト <<, 右シフト >>)
    0/1 列を右や左にシフトする
    左シフトが上の桁の方にシフトする
    1 1 1 0 0 1 0 0
    1 0 1 1 1 0 0 1
    A
    A << 2
    シフトしたときに詰める数字によっていくつかバリ
    エーションが存在する

    View full-size slide

  53. 53
    C 言語でのビット演算
    // ビット演算は整数型に対してのみ使える
    A = B & C;
    A &= B;
    A = B << 4; // 下位ビットには 0 が詰められる
    A <<= 4;
    A = ~B;
    // 符号なし 64 ビットリテラルを扱うとき
    A = UINT64_C(0xCCCCCCCCCCCCCCCC);

    View full-size slide

  54. 54
    ビット列の基本的操作
    特定のビットを操作する
    ビットを立てる
    ビットを下ろす
    ビットを反転する
    ビットの値を取得する
    マスク

    View full-size slide

  55. 55
    特定のビットを操作する
    UINT64_C(1) << index で index ビット目のみ立った数を
    表せる
    //A の index ビット目を操作する
    A |= UINT64_C(1) << index; // ビットを立てる
    A &= ~(UINT64_C(1) << index); // ビットを下ろす
    A ^= UINT64_C(1) << index; // ビットを反転する
    result = (A >> index) & 1; // ビットの値を取得する

    View full-size slide

  56. 56
    マスク
    // 奇数ビット目をクリア
    // 0x5 = 0b0101
    A &= UINT64_C(0x5555555555555555);
    // 偶数ビットなら 0xAAA... = 0b1010...
    // 2 ビットごとに交互にクリアするなら
    0xCCC...=0b11001100...

    View full-size slide

  57. 57
    ビット演算テクニック集
    立っているビットの数を数える (popcount)
    ビット列のハミング距離
    立っている一番下のビットを求める
    立っているビット列を走査する
    立っている一番上のビットを求める
    ビット列の並びを反転する

    View full-size slide

  58. 58
    ビット演算テクニック集
    部分集合の列挙
    ビット列を一部だけスワップする
    ビット列の指定した場所を詰めて並べる

    View full-size slide

  59. 59
    立っているビットの数を数える
    (popcount)
    1 0 1 1 1 0 0 1
    1 1 1 0
    1 1 1 0
    A
    A&0xAA
    (A&0xAA) >> 1
    0 1 0 1
    A&0x55
    0 1 1 0 0 1 0 1
    ((A&0xAA) >> 1)
    + (A&0x55)
    2ビットごとの立っているビットの数の和

    View full-size slide

  60. 60
    立っているビットの数を数える
    (popcount)
    0 1 1 0 0 1 0 1
    0 1 0 1
    0 1 0 1
    A'
    A'&0xCC
    (A'&0xCC) >> 2
    1 0 0 1
    A'&0x33
    0 0 1 1 0 0 1 0
    ((A'&0xCC) >> 2)
    + (A'&0x33)
    4ビットごとの立っているビットの数の和

    View full-size slide

  61. 61
    立っているビットの数を数える
    (popcount)
    0 0 1 1 0 0 1 0
    0 0 1 1
    0 0 1 1
    A''
    A''&0xF0
    (A''&0xF0) >> 4
    0 0 1 0
    A''&0x0F
    0 0 0 0 0 1 0 1
    ((A''&0xF0) >> 4)
    + (A''&0x0F)
    8ビット全体の立っているビットの数の和

    View full-size slide

  62. 62
    ビット列のハミング距離を求める
    ハミング距離:同じ長さの配列の対応する位置にある、
    異なった値を持つ要素の数
    0 0 1 0 1 1 0 0
    1 0 1 1 1 0 0 1
    A
    B
    Hamming(A, B) = 4




    =
    =
    =
    =

    View full-size slide

  63. 63
    ビット列のハミング距離を求める
    // ビット列のハミング距離は XOR の popcount で求まる
    Hamming(A, B) = popcount(A ^ B);

    View full-size slide

  64. 64
    複数のビットから成るデータの配列の 
    ハミング距離を求める
    例えば、 2 ビットから成るデータの配列のハミング距離
    を求めるとき、
    for (int i = 0; i < ARY_SIZE; ++i) {
    uint8_t C = A[i] ^ B[i]; //8 ビットの場合
    C = ((C & 0xAA) >> 1) | (C & 0x55);
    result += popcount(C);
    }

    View full-size slide

  65. 65
    立っている一番下のビットを求める
    B = A & -A;
    1 0 1 1 1 0 0 0
    A
    0 1 0 0 1 0 0 0
    -A
    0 0 0 0 1 0 0 0
    A & -A
    -A = ~A + 1であることを利用

    View full-size slide

  66. 66
    立っている一番下のビットを中心に操作
    A & (A – 1); // 立っている一番下のビットをクリア
    A ^ -A; // 立っている一番下のビットより上の桁を 1 に
    A | -A; // さらに立っている一番下のビットも 1 に
    // 立っている一番下のビットより下の桁を 1 に
    A ^ (A – 1)

    View full-size slide

  67. 67
    立っているビット列を走査する
    // i &= i-1 で i の立っている一番下のビットをクリア
    for (uint64_t i = bits; i != 0; i &= i-1) {
    uint64_t rmb = i & -i;
    // 何らかの処理
    }
    立っているビットの数が少ない場合には、この方法でも
    立っているビットの数を高速に数えられる

    View full-size slide

  68. 68
    立っている一番上のビットを求める
    二分探索で求める
    0 0 1 1 1 0 0 1
    A
    0 0 1 1
    A&0xF0 != 0
    立っている一番上のビットは上4桁にある

    View full-size slide

  69. 69
    立っている一番上のビットを求める
    二分探索で求める
    0 0 1 1
    B=A&0xF0
    立っている一番上のビットは上から0,1,4,5桁目にはない
    B&0xCC 0 0 0 0 = 0

    View full-size slide

  70. 70
    立っている一番上のビットを求める
    二分探索で求める
    0 0 1 1
    C=B
    立っている一番上のビットは上から0,2,4,6桁目にある
    C&0xAA 0 1 0 0 != 0

    View full-size slide

  71. 71
    立っている一番上のビットを求める
    // 8 ビットの場合
    A = (A & 0xF0) ? (A & 0xF0) : A;
    A = (A & 0xCC) ? (A & 0xCC) : A;
    A = (A & 0xAA) ? (A & 0xAA) : A;

    View full-size slide

  72. 72
    ビット列の並びを反転する
    上位ビットを下に、下位ビットを上に持ってくる
    1 0 1 1 1 0 0 1
    A
    1 0 0 1 1 1 0 1
    Aの反転

    View full-size slide

  73. 73
    ビット列の並びを反転する
    分割統治法を用いる
    1 0 1 1 1 0 0 1
    A
    1 0 0 1 1 0 1 1
    (A>>4) | (A<<4)

    View full-size slide

  74. 74
    ビット列の並びを反転する
    分割統治法を用いる
    1 0 1 1 1 0 0 1
    A
    1 0 0 1 1 0 1 1
    B=(A>>4) | (A<<4)
    0 1 1 0 1 1 1 0
    ((B&0xCC)>>2)
    | ((B&0x33)<<2)

    View full-size slide

  75. 75
    ビット列の並びを反転する
    分割統治法を用いる
    1 0 1 1 1 0 0 1
    A
    1 0 0 1 1 0 1 1
    B=(A>>4) | (A<<4)
    0 1 1 0 1 1 1 0
    C=((B&0xCC)>>2)
    | ((B&0x33)<<2)
    ((C&0xAA)>>1)
    | ((C&0x55)<<1)
    1 0 0 1 1 1 0 1

    View full-size slide

  76. 76
    ビット列の並びを反転する
    分割統治法を用いる
    1 0 1 1 1 0 0 1
    A
    1 0 0 1 1 1 0 1
    Aの反転

    View full-size slide

  77. 77
    ビット列を一部だけスワップする
    Delta Swap という手法
    長さが等しく、重複のないビット列をスワップする
    * * * A B C * * * a b c * * * *
    この幅をdeltaとする
    0 0 0 0 0 0 0 0 0 1 1 1 0 0 0 0
    x
    mask
    ABCとabcをスワップする
    deltaとmaskを予め求めておく

    View full-size slide

  78. 78
    ビット列を一部だけスワップする
    Delta Swap という手法
    長さが等しく、重複のないビット列をスワップする
    * * * A B C * * * a b c * * * *
    0 0 0 0 0 0 0 0 0 P Q R 0 0 0 0
    x
    b := (x ^ (x >> delta)) & mask
    P = A^a, Q = B^b, R = C^cとなる
    b

    View full-size slide

  79. 79
    ビット列を一部だけスワップする
    Delta Swap という手法
    長さが等しく、重複のないビット列をスワップする
    * * * A B C * * * a b c * * * *
    0 0 0 P Q R 0 0 0 P Q R 0 0 0 0
    x
    ABCとabcが入れ替わった(他のビットは変化なし)
    c := b ^ (b << delta)
    c
    * * * a b c * * * A B C * * * *
    c ^ x

    View full-size slide

  80. 80
    ビット列の指定した場所を詰めて並べる
    64bit 整数で 8x8 の 2 次元データを表現するとする
    オセロやチェスなどの盤面の状態、など
    0 1 2 3 4 5 6 7
    8 9 10 11 12 13 14 15
    16 17 18 19 20 21 22 23
    24 25 26 27 28 29 30 31
    32 33 34 35 36 37 38 39
    40 41 42 43 44 45 46 47
    48 49 50 51 52 53 54 55
    56 57 58 59 60 61 62 63

    View full-size slide

  81. 81
    ビット列の指定した場所を詰めて並べる
    斜めの効きを調べたいときに下のようなビットを整列し
    て並べたい
    14
    17 21
    26 28
    42 44
    49 53
    17 26 44 53 49 42 28 21
    14

    View full-size slide

  82. 82
    ビット列の指定した場所を詰めて並べる
    delta swap などを使えば、一応実現可能だが、もっと高
    速な方法が存在する
    14
    17 21
    26 28
    42 44
    49 53
    17 26 44 53 49 42 28 21
    14

    View full-size slide

  83. 83
    magic bitboard
    実はビットの配置ごとに適切に選んだ整数 magic を掛け
    てやると 1 列に並ぶ
    magic bitboard として知られている
    14
    17 21
    26 28
    42 44
    49 53
    17 26 44 53 49 42 28 21
    14
    * magic ) >> some =
    (

    View full-size slide

  84. 84
    magic bitboard
    チェスのそれぞれの駒の効きの範囲に対する magic
    number は既に知られている
    そうでないようなビット配置に対しては、自分で一々
    magic number を求めなければならない
    あとで述べる pext 命令の追加により実質役割を終えた

    View full-size slide

  85. 85
    ビット演算関連の CPU 命令
    bts, btr, btc, bt
    blsi, blsmsk, blsr, tzcnt
    lzcnt
    bzhi
    bextr
    pext
    pdep

    View full-size slide

  86. 86
    bts, btr, btc, bt
    特定のビットを操作する命令
    bts
    unsigned char _bittestandset( __int32* a, __int32 b);
    unsigned char _bittestandset64( __int64* a, __int64 b);
    特定のビットを立て、立てる前のそのビットの状態を返す
    btr: ビットを下ろす
    btc: ビットを反転する
    bt: ビットを取得する

    View full-size slide

  87. 87
    blsi, blsmsk, blsr, tzcnt
    立っている一番下のビット関連の操作をする命令
    blsi: A & -A と同じ。立っている一番下のビットを求める
    blsmsk: A ^ -A と同じ。立っている一番下のビットより
    上のビットを 1 にした数を求める
    blsr: A & (A-1) と同じ。立っている一番下のビットを 0
    にした数を求める。
    tzcnt: 立っている一番下のビットの桁数を求める

    View full-size slide

  88. 88
    lzcnt, bzhi
    lzcnt
    立っている一番上のビットの上にある 0 の数を数える
    bzhi
    src の下 n 桁を dest に代入する

    View full-size slide

  89. 89
    bextr
    bextr
    start, len を指定して図のような操作をした結果を返す
    0 1 0 0 1 0 1 1 0 1 0 1 1 1 1 0
    start
    len
    0 0 0 0 0 0 0 0 0 0 0 0 1 1 0 1
    len

    View full-size slide

  90. 90
    pext
    mask で指定した位置にあるビットを下から詰めて並べる
    1 0 1 1 1 0 1 1
    mask
    0 1 0 1 0 0 1 1
    src
    0 0 0 0 1 0 1 1
    結果

    View full-size slide

  91. 91
    pdep
    mask で指定した位置にビットを入れていく
    pext のだいたい逆の操作ができる
    1 0 1 1 1 0 1 1
    mask
    0 1 0 1 0 0 1 1
    src
    0 0 1 0 0 0 1 1
    結果

    View full-size slide

  92. 92
    SIMD 命令による最適化
    SIMD 命令とは
    SIMD 命令が使える場合
    SIMD 命令の使い方
    SIMD プログラミングにおけるポイント
    メモリアラインメント
    マスクによる条件分岐の除去
    その他 SIMD 命令で可能なこと
    AVX-512 について

    View full-size slide

  93. 93
    SIMD 命令とは
    SIMD(Single Instruction Multiple Data の略 )
    一つの命令で複数のデータを処理する
    AVX 命令セットでは一度に 256bit(32Bytes) のデータを
    一括で処理できる
    単精度浮動小数 x8
    倍精度浮動小数 x4
    1/2/4/8Bytes 整数 x32/16/8/4

    View full-size slide

  94. 94
    SIMD 命令とは
    A[0] A[1] A[2] A[3]
    B[0] B[1] B[2] B[3]
    A[0]+B[0] A[1]+B[1] A[2]+B[2] A[3]+B[3]

    256ビット

    View full-size slide

  95. 95
    SIMD 命令が使える場合
    データがメモリ上で連続している
    そのデータに対しそれぞれ同じ、または似た操作をする
    という条件を満たさなければならない
    実際には、データが連続していなくても、近くにまとまって
    いれば、ギャザー命令を利用して一気にに計算可能だが、
    ギャザー命令は遅延が大きい

    View full-size slide

  96. 96
    SIMD 命令の使い方
    ●既に最適化されたライブラリを使う ( 割愛 )
    ●自動ベクトル化の機能のあるコンパイラでコンパイル
    ●AVX の C++ ベクトルクラス
    ●SIMD 組み込み関数の利用
    ●アセンブリ言語

    View full-size slide

  97. 97
    自動ベクトル化
    コンパイラによっては適切なコンパイルオプションを指
    定することによって、演算を SIMD 命令に変換してくれる
    icc なら -xAVX を指定すると AVX 命令を使ってくれる
    しかし、どんな場合も SIMD 命令に変換してくれるわけで
    はない
    そのときは、残り 3 通りのいずれかの方法を用いる

    View full-size slide

  98. 98
    AVX の C++ ベクトルクラス
    F64vec4/F32vec8 といったクラスが用意されている
    いい感じに算術演算子と添字演算子がオーバーロードさ
    れている

    View full-size slide

  99. 99
    SIMD 組み込み関数の利用
    #include “immintrin.h” で使えるようになる
    単精度浮動小数 x8 型 __m256
    倍精度浮動小数 x4 型 __m256d
    整数 x4/8/16/32 型 __m256i
    __m256d resv = _mm256_add_pd(a, b); のように使う
    今回は主に SIMD 組み込み関数を利用していきます

    View full-size slide

  100. 100
    アセンブリ言語
    コンパイラはなるべく一般的な最適化をしようとするの
    で、特殊な状況では最適でないプログラムを出力するこ
    とがある
    そのような場合の最終手段として直接アセンブリ言語を
    書いてプログラミングすることが考えられる
    ただし、多くの場合人間よりコンパイラのほうが賢い最
    適化をするので、プロファイラなどによりパフォーマン
    スが向上するとわかった時のみ行うべき

    View full-size slide

  101. 101
    SIMD 組み込み関数の利用例
    配列 A の各要素に配列 B の対応する要素を加算した結果
    を配列 C に格納する。
    配列の長さは 4 の倍数、 A,B,C はすべて倍精度浮動小数型
    for (int i = 0; i < SIZE/4; ++i) {
    __m256d va = _mm256_load_pd(A+4*i);
    __m256d vb = _mm256_load_pd(B+4*i);
    __m256d res = _mm256_add_pd(va, vb);
    _mm256_store_pd(C+4*i, res);
    }

    View full-size slide

  102. 102
    SIMD 組み込み関数の使い方の基本
    メモリ上からデータを取ってくるには _mm256_load*
    メモリ上にデータを書き込むには _mm256_store*
    _mm256_( 動作 )_( 型 ) となっていることが多い
    型には
    ps: 単精度浮動小数 x8
    pd: 倍精度浮動小数 x4
    epi8/16/13/64: 8/16/32/64 バイト整数 x32/16/8/4

    View full-size slide

  103. 103
    Intel Intrinsics Guide
    https://software.intel.com/sites/landingpage/Intrinsic
    sGuide/
    SIMD 組み込み関数を中心とした各種 CPU 命令の組み込
    み関数のリファレンス
    拡張の種類や命令のタイプで絞ったり検索したりできる
    便利!!!

    View full-size slide

  104. 104
    SIMD プログラミングにおけるポイント
    データを 32 バイトアラインメントに合わせる
    合っていないと実行速度が低下する
    同時にたくさんの処理が行える分、メモリアクセスが律
    速になりやすい
    キャッシュを考慮したプログラミング
    ループ中に条件分岐があるとベクトル化できない
    マスクを有効活用するなどして条件分岐を除去する

    View full-size slide

  105. 105
    メモリアラインメント 静的確保
    alignas(32) uint64_t ary[SIZE];
    ary の先頭は 32 バイト境界に合わせられる

    View full-size slide

  106. 106
    メモリアラインメント 動的確保
    メモリを動的確保する場合は __alignas は使えない
    void *_mm_malloc(size_t size, size_t align);
    size バイトの領域を確保。返されるポインタのアドレス
    は align の倍数になる
    _mm_malloc で確保したメモリは _mm_free で開放する
    他にも posix_memalign/aligned_alloc を使う方法や
    std::align を使う方法もある

    View full-size slide

  107. 107
    マスクによる条件分岐の除去
    __m256i _mm256_cmpeq_epi8(__m256i, __m256i);
    各バイトを比較して等しいなら FF 、異なるなら 00 を各
    バイトに入れて返す
    eq (等しい) /gt (大きい)、
    epi8/16/32/64 ( 1/2/4/8 バイト)のバリエーションが存
    在する

    View full-size slide

  108. 108
    マスクによる条件分岐の除去
    等しくないなら〜〜の場合は eq のビットを反転させれば
    良い
    <>≦≧の 4 種類の不等号もビットを反転させたり左右
    を入れ替えたりすれば実現可能
    浮動小数点型の場合は組み込み関数にどのような比較を
    するかを指定して渡せる

    View full-size slide

  109. 109
    マスクによる条件分岐の除去
    SIMD 比較命令によりマスクを作る
    演算対象と AND をとったりして正しい結果が得られるよ
    うにする

    View full-size slide

  110. 110
    マスクによる条件分岐の除去
    例 )
    for (int i = 0; i < SIZE; ++i) {
    if (a[i] > b[i]) {
    c[i] += a[i];
    }
    }
    条件分岐が入っているのでこのままでは並列に足し算で
    きない

    View full-size slide

  111. 111
    マスクによる条件分岐の除去
    例 )
    for (int i = 0; i < SIZE/32; ++i) {
    __m256i va = _mm256_load_si256((__m256i*)(A+32*i));
    __m256i vb = _mm256_load_si256((__m256i*)(B+32*i));
    __m256i mask = _mm256_cmpgt_epi8(va, vb);
    va = _mm256_and_si256(va, mask);
    __m256i vc = _mm256_load_si256((__m256i*)(C+32*i));
    vc = _mm256_add_epi8(vc, va);
    _mm256_store_si256((__m256i*)(C+32*i), vc);
    }

    View full-size slide

  112. 112
    マスクによる条件分岐の除去
    例 )
    for (int i = 0; i < SIZE/32; ++i) {
    __m256i va = _mm256_load_si256((__m256i*)(A+32*i));
    __m256i vb = _mm256_load_si256((__m256i*)(B+32*i));
    __m256i mask = _mm256_cmpgt_epi8(va, vb);
    va = _mm256_and_si256(va, mask);
    __m256i vc = _mm256_load_si256((__m256i*)(C+32*i));
    vc = _mm256_add_epi8(vc, va);
    _mm256_store_si256((__m256i*)(C+32*i), vc);
    }

    View full-size slide

  113. 113
    マスクによる条件分岐の除去
    例 )
    for (int i = 0; i < SIZE/32; ++i) {
    __m256i va = _mm256_load_si256((__m256i*)(A+32*i));
    __m256i vb = _mm256_load_si256((__m256i*)(B+32*i));
    __m256i mask = _mm256_cmpgt_epi8(va, vb);
    va = _mm256_and_si256(va, mask);
    __m256i vc = _mm256_load_si256((__m256i*)(C+32*i));
    vc = _mm256_add_epi8(vc, va);
    _mm256_store_si256((__m256i*)(C+32*i), vc);
    }

    View full-size slide

  114. 114
    マスクによる条件分岐の除去
    例 )
    for (int i = 0; i < SIZE/32; ++i) {
    __m256i va = _mm256_load_si256((__m256i*)(A+32*i));
    __m256i vb = _mm256_load_si256((__m256i*)(B+32*i));
    __m256i mask = _mm256_cmpgt_epi8(va, vb);
    va = _mm256_and_si256(va, mask);
    __m256i vc = _mm256_load_si256((__m256i*)(C+32*i));
    vc = _mm256_add_epi8(vc, va);
    _mm256_store_si256((__m256i*)(C+32*i), vc);
    }

    View full-size slide

  115. 115
    マスクによる条件分岐の除去
    例 )
    for (int i = 0; i < SIZE/32; ++i) {
    __m256i va = _mm256_load_si256((__m256i*)(A+32*i));
    __m256i vb = _mm256_load_si256((__m256i*)(B+32*i));
    __m256i mask = _mm256_cmpgt_epi8(va, vb);
    va = _mm256_and_si256(va, mask);
    __m256i vc = _mm256_load_si256((__m256i*)(C+32*i));
    vc = _mm256_add_epi8(vc, va);
    _mm256_store_si256((__m256i*)(C+32*i), vc);
    }

    View full-size slide

  116. 116
    その他 SIMD 命令で可能なこと
    複素数の計算
    水平加算(隣合う要素と足し算する)などを利用すると、
    愚直にやるより高速に計算できる
    バイト列の並び替え
    シャッフル命令によりバイト列を逆順に並び替えたり、よ
    り複雑な並び替えが高速に行える
    数学関連の関数の値の計算
    三角関数や指数関数などを並列に計算できる命令がある

    View full-size slide

  117. 117
    AVX-512 について
    Broadwell の次の Skylake で Xeon( サーバー向けプロ
    セッサ ) に AVX-512 拡張命令が追加される予定
    一度に 512 ビット扱えるだけでなく、ほぼすべての命令
    にマスクを掛けることが出来るようになる
    各種の便利命令がてんこ盛り
    Core i3/5/7 等には Skylake の次の Cannonlake で入る?

    View full-size slide

  118. 118
    今日触れられなかった内容
    並列実行時における最適化
    アウト・オブ・オーダーやスーパースカラーを意識した
    命令の選択
    プリフェッチ命令

    View full-size slide

  119. 119
    まとめ
    今回は 3 つの最適化手法について解説した
    キャッシュを意識したプログラミング
    ビット演算
    SIMD 命令
    いずれの手法も上手く使えれば数倍〜数十倍の高速化が
    見込める

    View full-size slide

  120. 120
    参考文献等
    [0] インテル ® 64 アーキテクチャーおよび IA-32 アーキテク
    チャー最適化リファレンス・マニュアル
    http://www.intel.co.jp/content/dam/www/public/ijkk/jp/ja/do
    cuments/developer/248966-024JA.pdf
    [1] 英語の最新版
    http://www.intel.co.jp/content/dam/www/public/us/en/docum
    ents/manuals/64-ia-32-architectures-optimization-manual.pdf
    [2] Intel Intrinsics Guide
    https://software.intel.com/sites/landingpage/Intrin
    sicsGuide/

    View full-size slide

  121. 121
    参考文献等
    [3] Intel® Architecture Instruction Set Extensions
    Programming Reference
    https://software.intel.com/sites/default/files/managed/0d/53/
    319433-022.pdf
    [4] Intel® 64 and IA-32 Architectures Software Developer’s
    Manual
    http://www.intel.co.jp/content/dam/www/public/us/en/docu
    ments/manuals/64-ia-32-architectures-software-developer-
    manual-325462.pdf
    [5] CPU – Wikipedia
    http://ja.wikipedia.org/wiki/CPU

    View full-size slide

  122. 122
    参考文献等
    [6] Chess Programming Wiki
    https://chessprogramming.wikispaces.com/

    View full-size slide