2018年3月14日水曜日

ClojureでOpenCLを利用して行列演算してみた


背景

ClojureとはJavaというプログラミング言語の上で動くLisp系のプログラミング言語です。
OpenCLとは、GPUを使って並列処理が行えるAPI(GPGPU)の一つです。

ClojureからOpenCLを呼び出すとどれくらいの早さで動くのか気になったので、確かめてみました。
実装する中で時間を取られた箇所があるので、それと合わせて内容を共有します。

全体像

  1. 使ったもの
  2. プログラムの説明
  3. 早さ比較
  4. まとめ

使ったもの


プログラムの説明

clojureclというライブラリを利用して、OpenCLを利用しました。

clojureclを利用するためのプロジェクトを作ります。
今回はclj-opencl-practiceという名前のプロジェクトを作りました。
lein new clj-opencl-practice

clojureclを呼び出すために、project.cljにclojureclを追加します。
project.clj
  :dependencies [..
                 [uncomplicate/clojurecl "0.8.0"]]

replを起動します。
replとは、プログラムを書いてエンターを押すとそれが実行される、対話型のプログラミング環境です。
emacsからciderを利用して呼び出すのが、補完などができるので楽ですが、プロジェクトのディレクトリで下記のコマンドを入力することでもreplを起動できます。
lein repl

requireでclojureclを呼んで、利用できる状態にします。
(require '[uncomplicate.clojurecl.core :refer :all]
         '[uncomplicate.clojurecl.info :refer :all]
         '[uncomplicate.commons.core :refer [release]])

OpenCLの実行に必要なGPUの情報を取得します。
(def first-platform (first (platforms)))
(def first-gpu (first (devices first-platform)))
(def ctx (context [first-gpu]))
(def queue (command-queue ctx first-gpu))

どのGPUを使おうとしているかは、 name-infoを利用して内容を確認できます。
(map name-info (platforms))
 => ("Intel(R) OpenCL")
(map name-info (devices first-platform))
 => ("Intel(R) Core(TM) i7-4600U CPU @ 2.10GHz")

実行したいOpenCLのカーネル(OpenCLの関数)を文字列として定義します。 行列Aと行列Bをかけ合わせて行列Resultを取得するカーネルを定義しました。
(def kernel-source "
  __kernel void matrix_dot_matrix(
    __global const float* A,
    __global const float* B,
    __global float* Result,
    const int wA,
    const int wB
  ) {
    const int x = get_global_id(0);
    const int y = get_global_id(1);
    float value = 0;
    for (int i = 0; i < wA; ++i) {
      int index_a = y * wA + i;
      int index_b = i * wB + x;
      float elementA = A[index_a];
      float elementB = B[index_b];
      if (elementA != 0.0 && elementB != 0.0) {
        value = value + elementA * elementB;
      }
    }
    Result[y * wB + x] = value;
  }")

文字列をOpenCLに登録して、定義したカーネルを呼び出します。
(def kernels (build-program! (program-with-source ctx [kernel-source])))
(def mul-matrixes-k (kernel kernels "matrix_dot_matrix"))

計算に使う高さと幅の高さと幅が10の行列を準備します。
(def matrix-len 10)
(def a-width matrix-len)
(def a-height matrix-len)
(def b-width matrix-len)
(def b-height a-width)
(def r-width b-width)
(def r-height a-height)
(def matrix-a (float-array (range (* a-width a-height))))
(def matrix-b (float-array (for [x (range b-width)
                                 y (range b-height)]
                             (if (not= x y)
                               0
                               (if (= x 1)
                                 2
                                 1)))))
(def matrix-result (float-array (* r-width r-height)))

行列AとBは、このような内容です。
行列Bは単位行列の座標が(左上を[0,0]とする場合の)[1, 1]の要素だけ2にしたものです。
AとBの掛け合わせにより、(左の列を0番目とする場合の)1番目の列だけ2倍になったAがResultとして出力されることを期待します。
(defn print-matrix [values w h]
  (if (> w 10)
    (prn (take 1 values))
    (doall
     (for [line (partition w values)]
       (prn line)))))

(print-matrix matrix-a a-width a-height)
(print-matrix matrix-b b-width b-height)

行列Aの内容
(0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0)
(10.0 11.0 12.0 13.0 14.0 15.0 16.0 17.0 18.0 19.0)
(20.0 21.0 22.0 23.0 24.0 25.0 26.0 27.0 28.0 29.0)
(30.0 31.0 32.0 33.0 34.0 35.0 36.0 37.0 38.0 39.0)
(40.0 41.0 42.0 43.0 44.0 45.0 46.0 47.0 48.0 49.0)
(50.0 51.0 52.0 53.0 54.0 55.0 56.0 57.0 58.0 59.0)
(60.0 61.0 62.0 63.0 64.0 65.0 66.0 67.0 68.0 69.0)
(70.0 71.0 72.0 73.0 74.0 75.0 76.0 77.0 78.0 79.0)
(80.0 81.0 82.0 83.0 84.0 85.0 86.0 87.0 88.0 89.0)
(90.0 91.0 92.0 93.0 94.0 95.0 96.0 97.0 98.0 99.0)

行列Bの内容
(1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0)
(0.0 2.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0)
(0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0)
(0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0)
(0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0)
(0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0)
(0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0)
(0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0)
(0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0)
(0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0)

行列A、B、ResultのOpenCLのバッファを確保し、行列A、Bの値を書き込みます。
(def sizeof-float 4)
(def matrix-a-buffer (cl-buffer ctx (* sizeof-float a-width a-height) :read-only))
(def matrix-b-buffer (cl-buffer ctx (* sizeof-float b-width b-height) :read-only))
(def matrix-result-buffer (cl-buffer ctx (* sizeof-float r-width r-height) :write-only))
(enq-write! queue matrix-a-buffer matrix-a)
(enq-write! queue matrix-b-buffer matrix-b)

上記のバッファ確保時は、バッファの長さに注意が必要です。
OpenCLのfloat型は4バイトなので、「行列の幅x行列の高さx4」のバッファを確保する必要があります。
そのため、上記の記述ではsizeof-floatを4として定義し、各行列のバッファ確保時に利用しています。
float以外の型を使う場合、OpenCLの仕様に合わせて、型のサイズを変える必要があると思います。

バッファのサイズが合っていないままカーネルの処理を行うと、下記のようなエラーが表示されます。
確保されてないメモリへのアクセスがこのエラーを引き起こしているのだと思いますが、エラーの内容からはそれを連想できなかったので原因が分かるまでに時間を取られました。

下記のエラーが表示された場合は、バッファのサイズを確認したり試しに十分大きく割り当ててみて、同じエラーが出るか、直るのか、確認してみても良いかもしれません。
#
# A fatal error has been detected by the Java Runtime Environment:
#
# SIGSEGV (0xb) at pc=0x00007f4aa246bad3, pid=12633, tid=0x00007f4a9b7fd700
#
# JRE version: OpenJDK Runtime Environment (8.0_151-b12) (build 1.8.0_151-8u151-b12-0ubuntu0.17.10.2-b12)
# Java VM: OpenJDK 64-Bit Server VM (25.151-b12 mixed mode linux-amd64 compressed oops)
# Problematic frame:
# C [libcpu_device.so+0x23ad3]
#
# Failed to write core dump. Core dumps have been disabled. To enable core dumping, try "ulimit -c unlimited" before starting Java again
#
# An error report file with more information is saved as:
# /home/asuki/gitprojects/clj/clj-opencl-practice/hs_err_pid12633.log
#
# If you would like to submit a bug report, please visit:
# http://bugreport.java.com/bugreport/crash.jsp
#

確保したバッファや変数をカーネルに割り当てます。
(set-arg! mul-matrixes-k 0 matrix-a-buffer)
(set-arg! mul-matrixes-k 1 matrix-b-buffer)
(set-arg! mul-matrixes-k 2 matrix-result-buffer)
(set-arg! mul-matrixes-k 3 (int-array [a-width]))
(set-arg! mul-matrixes-k 4 (int-array [b-width]))

上記のカーネルへの値割り当て時に注意すべきなのが、int型の変数を割り当てる際もint-array型に変換する必要があることです。
下記のように行列Aの幅をそのまま渡そうとすると、
(set-arg! mul-matrixes-k 3 a-width) ; Wrong

下記のようなエラーが表示されます。
#object[org.jocl.cl_kernel 0x73695b14 "cl_kernel[0x7fda7c001588]"]
clj-opencl-practice.core=> (set-arg! mul-matrixes-k 3 a-width)

ExceptionInfo OpenCL error: CL_INVALID_ARG_SIZE.  clojure.core/ex-info (core.clj:4739)

値の渡し方がどう悪いのか分からず時間を取られましたが、clojureclのテストプログラムに参考になる記述があり、それを真似たことで長さが1のintもint-arrayに変換する必要があると分かりました。
バッファを確保してない変数を渡す際は、ご注意ください。

変数の割り当てが終わったら、OpenCLのwork size(並列処理の大きさ)をResultベクトルの幅と高さにして、カーネルを実行します。
(enq-nd! queue mul-matrixes-k (work-size-2d r-width r-height))

Result行列を取り出して、表示します。
(enq-read! queue matrix-result-buffer matrix-result)
(print-matrix matrix-result r-width r-height)

(左が0番目としたときの)1番目の列が2倍された行列Aが取得できたので、期待通りに計算できたようです。
(0.0 2.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0)
(10.0 22.0 12.0 13.0 14.0 15.0 16.0 17.0 18.0 19.0)
(20.0 42.0 22.0 23.0 24.0 25.0 26.0 27.0 28.0 29.0)
(30.0 62.0 32.0 33.0 34.0 35.0 36.0 37.0 38.0 39.0)
(40.0 82.0 42.0 43.0 44.0 45.0 46.0 47.0 48.0 49.0)
(50.0 102.0 52.0 53.0 54.0 55.0 56.0 57.0 58.0 59.0)
(60.0 122.0 62.0 63.0 64.0 65.0 66.0 67.0 68.0 69.0)
(70.0 142.0 72.0 73.0 74.0 75.0 76.0 77.0 78.0 79.0)
(80.0 162.0 82.0 83.0 84.0 85.0 86.0 87.0 88.0 89.0)
(90.0 182.0 92.0 93.0 94.0 95.0 96.0 97.0 98.0 99.0)

計算が終わったら、リソースを開放します。
(release matrix-a)
(release matrix-b)
(release matrix-result)
(release kernels)
(release queue)
(release ctx)

今回作成したプログラムは、githubのリポジトリで公開しています。
asukiaaa/clj-opencl-practice

早さ比較

以前作ったcのプログラムclojureで並列無しで演算するプログラムとの実行速度を比較してみました。


Program10x10 matrix128x128 matrix1028x1028 matrix
Calc (sec)Total (sec) Calc (sec)Total (sec) Calc (sec)Total (sec)
Clojure naive 0.018 0.041 1.209 1.222 452.516 452.528
Clojurecl 0.001 0.104 0.002 0.120 0.624 1.334
C naive 0.000025 0.000721 0.0109 0.0118 10.979 11.025
C OpenCL 0.00019 0.106 0.004 0.108 2.912 3.069

驚くことに、自分の環境では、Cで書いたOpenCLのプログラムよりClojureで書いたOpenCLの方が早く演算できました。
OpenCLと言語間のデータ転送が、CよりClojureの方が早かったりするのでしょうか。

まとめ

ClojureでOpenCLを利用して行列演算できました。
自分の予想に反して、今回の例だとCで書くより処理速度が早いプログラムができたので、ClojureでOpenCLを使うアプリケーションを作るのは、実用的なのかもしれません。

「こうすると、Cの方が早くなる」など、アドバイスをいただける方は、コメントを残していただけると嬉しいです。

共有する情報は以上です。
何かの参考になれば嬉しいです。

参考

Interactive GPU Programming - Part 2 - Hello OpenCL

0 件のコメント :