Writer: 近藤 鯛貴

岩手県立大学アドベントカレンダー18日目の記事です（大遅刻）

自己紹介

はじめまして近藤鯛貴です。

岩手県立大学ソフトウェア情報学研究科の博士課程に在学しつつ、Defios株式会社というベンチャー企業を経営しています。

大学ではラズパイとかJetsonなどの小型コンピュータを対象に高速化の研究を行っており、会社ではAI/IoT分野を中心に受託開発やコンサルティング業務を行ってます。

最近Defiosで技術ブログを立ち上げました。せっかくなのでこちらから投稿します。

1．はじめに

現代の大体のCPUではSIMD（Single Instruction Multiple Data)命令が実行可能です。SIMDとは一つの命令で複数のデータを処理する並列化形態です。上手く使うことでソフトウェアを高速化することが出来ます。
しかし、プロセッサのアーキテクチャや製品種類によって対応しているSIMDのbit幅や命令セットが異なり、様々なコンピュータ上でSIMDを活用するアプリケーションを実装するのは少し厄介です。

ふと、画像処理ライブラリであるOpenCVのSIMD対応がどうなっているのか気になって調べたところ以下の記事を見つけました。

どうやら、OpenCVではUniversal Intrinsicというライブラリを用いてSIMDアーキテクチャの差異を吸収しているようです。

本記事ではこのUniversal Intrinsic(以下UI)の使い方と簡単な実装例を紹介します。

2．使い方

普通にOpenCVをインストールして、<opencv2/core/simd_intrinsics.hpp>をincludeすればUIを使うことが出来ます。
SIMDが有効になっているか、どのbit幅が使えるかはOpenCVのサンプルのsimd_basic.cppを実行することで確認できます。
以下AVX2環境でのsimd_basic.cppの実行結果

==================  macro dump  ===================
CV_SIMD is defined: 1
CV_SIMD_WIDTH is defined: 32
CV_SIMD128 is defined: 1
CV_SIMD256 is defined: 1
CV_SIMD512 is defined: 0
CV_SIMD_64F is defined: 1
CV_SIMD_FP16 is defined: 0
=================  sizeof checks  =================
sizeof(v_uint8) = 32
sizeof(v_int32) = 32
sizeof(v_float32) = 32
==================  arithm check  =================
(vx_setall_u8(10) + vx_setall_u8(45)).get0() => 55
=====================  done  ======================

AVX2のSIMD幅は256bitです。CV_SIMD256が1となっているので有効になっていることが分かります。
ちなみに、私の環境ではコンパイルする際に-mavx2(gcc)オプションを入れないとAVX2が有効にならずにSSEで実行されました。

3．行列積を実装してみる

UIのドキュメントとしてチュートリアルや関数のリファレンスなどがあります。
これらを参考に行列積を簡単に書いてみました。
今回は使ってみることが主目的なので最適化はあまり考えてないです(許して..)。

int step=v_float32().nlanes;
#pragma omp parallel for
for(int i=0;i<SIZE;i++){
  for(int j=0;j<SIZE;j+=step){
    v_float32 c_vec = vx_setall_f32(0);
    for(int k=0;k<SIZE;k++){
      float *ptr = b + j + k * SIZE;
      v_float32 b_vec = vx_load(ptr);
      v_float32 a_vec = vx_setall_f32(a[(i*SIZE)+k]);
      c_vec = v_fma(a_vec,b_vec,c_vec);
    }
    v_store(c+i*SIZE+j,c_vec);
  }
}

行列のC=ABを計算します。最内ループにてAの１要素とBのベクトルを積和演算して計算してます。

v_float32().nlanesで32bit floatのレーン数を取得してループ回数を決定しています。
vx～系の命令は有効になっているSIMD幅で処理するものなので上記コードはどんなアーキテクチャでも対応していればSIMD並列化して計算してくれるコードになります（おそらく）。

ちなみに、今回スレッド並列化はOpenMPでやりましたが、OpenCVにはUIのように複数ある並列化APIを上手く使う_paralell_for_機能があります。

4．性能を見てみる

上記コードとSIMDなし実装の実行時間を比べてみました。 SIMDなし実装は以下になります。

#pragma omp parallel for
for(int i=0;i<SIZE;i++){
  for(int k=0;k<SIZE;k++){
    for(int j=0;j<SIZE;j++){
      c[((i*SIZE)+j)] += a[((i*SIZE)+k)] * b[((k*SIZE)+j)];
    }
  }
}

便宜上SIMDあり/なしと区別していますが、SIMD並列化以外にも様々な最適化でソフトウェアの実行速度は大きく変わります。今回の比較は純粋なSIMDの比較ではなく、あくまで参考程度に見ていただければと思います。

コンパイラによる最適化を防ぐために-O0でコンパイルしました。

30%ほど速くなりました。ちなみに、これくらいのコードだとSIMDなし実装にO3最適化を行う方が速くなりました。O3になると自動ベクトル化がつくので何もしなくても自動でAVX2が使われます。

5. 最後に

今回使ったコードはこちらに上げました。
最初vx~系の命令を知らずに256bitを指定して実装したコードも含まれてます(mm-256.cpp)

あまり自分のソフトウェアでUIを使う機会はないかと思いますが、簡単にSIMDプログラミングが出来て良いですね。

OpenCVのリポジトリを見た感じ、処理名.simd.hppにUIを使った実装が記述されているようです。
畳み込み計算のような基本的な処理は大体SIMD実装が用意されていますがORBやakazeのような特徴抽出はあまり対応されていないようです（単純にSIMD並列が難しいor他関数の組み合わせで事足りるだけかもしれないですが）。

対応していない関数をSIMD化してプルリク出してみるのも面白いかもですね。