NeuralNet.java のバグを修正しました

NeuralNet.java でバグが見つかり、本来

for (int y = 0; y < ysize; ++y) {
  d += w[x][y] * bderivative[y]
}

と書くべきところを x に関して和を取っていました。修正し、GitHub を更新しました。

修正の結果、収束がさらに速くなった気がします。また、ループを繰り返すと過学習してしまう傾向がみられます。

JCuda 記事を追加しました

JCuda の記事として、デバイスメモリの使い方〜カーネルの記述までを追加しました。

preparePtxFile() は JCuda というよりは JCudaDriverSample の機能なので、敢えて解説していません。

次はいよいよ JCuda を使ってニューラルネットを構築します。

JCudaSimple

JCuda ライブラリの利用方法を説明するためのプロジェクト JCudaSimple をGitHubにアップロードしました。

中身はほとんど JCudaDriverSample と同じですが、Java8 風に書き換えてあります。外部 JAR として jcuda, jcuda-natives を参照しているので、ここだけ各自の環境にあったファイル名に書き換えてください。

解説記事は一両日中に書きたいと思います。

scrayping.py にバグ

このホームページを拡張してサブフォルダを作ってみたり、他のサイトを参照したりした結果、scrayping.py がエラーを出力するようになりました。

このバグを修正した版を GitHub に上げておきました。

意外に修正が難しかった。

隠れ層のあるニューラルネットワーク

隠れ層のあるニューラルネットワークに対する逆誤差伝搬法について記述した記事をアプロードしました。

一応確認したつもりですが、式が間違っていたらごめんなさい…

最初  \bar{L} = \frac{1}{S} \sum^S_s L なので、Lの個別値を求めて満足していました。でもどうも収束が遅い。
サンプル全体について添字縮約しなければいけないことに気づいたのが午後10時30分。

なもので、まだ GitHub に上げたサンプルプログラムにはバグフィックスがかかっていません。そのうち直します。

自作ニューラルネットが動かない理由が判明

以前、Java で作成する AI みたいな本を買ってきてプログラムを打ち込んでいたが、MNIST の手書き認識精度が 10% 台とほとんど認識できない状態になっていた。

今回、改めて Java で組んでみてその理由が判明。

σ'(z) とのアダマール積を取るべきところを加算しておりました。
そりゃ、収束せんわなあ。

プログラムを組み直したら劇的に収束が改善され、10回ループくらいで90%の認識率を持つようになりました。(さすが最急降下法!!)
これ以上の認識精度は難しいようです。

まだコメントが一切ないコードなので、コメントを付けて GitHub に上げようと思います。

また、コードを JCuda 版に改造して劇的な性能改善を図る予定です。
(CPU による Java8 の parallelStream は Core i7 6900K では期待した性能が出ないようです。並列化のオーバーヘッドが大きい?)

今は AutoEncoder による Deep Learning を試みていますが、これがまたうまく行かない。

JCuda の記事を一部公開しました

開発が止まっている JCuda ですが、一部説明記事を公開しました。

デバイスメモリポインタの使い方とか、cu ファイルとの協調とかをそのうち追加していきたいと思います。まずはダウンロードからコンパイルまでを記述しました。

 

マルチコア vs SIMD vs CUDA

前回あまりにも Numba が遅いのでおかしいと思い、マルチコア、SIMD、CUDA を比較するためのテストプログラムを作ってみた。

結果は、以下の通り。シングル CPU の非SIMD演算が最も遅く、CUDA が最も速いという予想通りの結果になりました。

version: 201307
#20 processors, 20 threads
start parallel calculation.
<<single-nosimd>>
elapsed 1.115483 seconds

<<single-simd>>
elapsed 0.535289 seconds

<<multicore-nosimd>>
elapsed 0.189586 seconds

<<multicore-simd>>
elapsed 0.128129 seconds

<<cuda>>
elapsed 0.021330 seconds

1.234571 1.581539 1.109199 1.103452 0.831745 1.106268 0.878185 1.868425 1.353009 0.748571 
1.234571 1.581539 1.109199 1.103452 0.831745 1.106268 0.878185 1.868425 1.353009 0.748571 

ただし、cudaMallocManaged() による Unified Memory を使うと速度が極端に低下する現象が見られました。(20倍くらい遅くなる)

やはり PCIe バスの遅さが実行時間に影響していると考えてよいようです。Numba のように自動的にメモリ配置してしまうライブラリは却って速度低下の原因となることがわかりました。

というわけで、CUDA を使うときはなるべく Device Memory に閉じるように設計しましょう、という当たり前の結論となりました。

この下は今回使ったソースコードです。マルチコア化には OpenMP を、SIMDには AVX2 を使用しました。

[Nelem.h]

#define NELEM (400 * 1000 * 1000) # 400M elements

[sample_cuda.h]

#if !defined(SAMPLE_CUDA_H_)
#define SAMPLE_CUDA_H_

#if defined(__cplusplus)
extern "C"
#endif
void prepare(float *a, float *b, float *c);

#if defined(__cplusplus)
extern "C"
#endif
float *getResult(float *);

#if defined(__cplusplus)
extern "C"
#endif
void release();

#if defined(__cplusplus)
extern "C"
#endif
void wrapper_for_cuda(void);

#endif

[sample.c]

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <unistd.h>
#include <omp.h>
#include <time.h>
#include <x86intrin.h>
#include "Nelem.h"
#include "sample_cuda.h"

#define alignof(x) __alignof__(x)

float *a, *b, *c;

double elapsed(struct timespec start, struct timespec end) {
  time_t second = end.tv_sec - start.tv_sec;
  long nsecond = end.tv_nsec - start.tv_nsec;
  return (double)second + (double)nsecond / (1000.0 * 1000.0 * 1000.0);
}

typedef void (*CBFUNC)(void);

static inline void simd_add(int i) {
    __m256 ma = _mm256_load_ps(&a[i]);
    __m256 mb = _mm256_load_ps(&b[i]);
    __m256 mc = _mm256_add_ps(ma, mb);
    _mm256_store_ps(&c[i], mc);
}

static void single_nosimd_add(void) {
  int i;
  for (i = 0; i < NELEM; ++i) {
    c[i] = a[i] + b[i];
  }
}

static void single_simd_add(void) {
  int i;
  for (i = 0; i < NELEM; i += 8) {
    simd_add(i);
  }
}

static void multicore_nosimd_add(void) {
  int i;

  #pragma omp parallel for private(i)
  for (i = 0; i < NELEM; ++i) {
    c[i] = a[i] + b[i];
  }
}

static void multicore_simd_add(void) {
  int i;
  #pragma omp parallel for private(i)
  for (i = 0; i < NELEM; i += 8) {
    simd_add(i);
  }
}

char *cbtitle[] = { "single-nosimd", "single-simd",
                    "multicore-nosimd", "multicore-simd",
                    "cuda" };
CBFUNC cbfunc[] = {
  single_nosimd_add, single_simd_add,
  multicore_nosimd_add, multicore_simd_add,
  wrapper_for_cuda
};
#define NCBFUNC (sizeof(cbfunc)/sizeof(cbfunc[0]))

int main() {
  int i;
  float *f;
  struct timespec start, end;

#if defined(_OPENMP)
  {
    int np = omp_get_num_procs();
    int nth = omp_get_max_threads();
    printf("version: %d\n", _OPENMP);
    printf("#%d processors, %d threads\n", np, nth);
  }
#endif
 
  a = (float *)_mm_malloc(sizeof(float) * NELEM, 32);
  b = (float *)_mm_malloc(sizeof(float) * NELEM, 32);
  c = (float *)_mm_malloc(sizeof(float) * NELEM, 32);

  for (i = 0; i < NELEM; ++i) {
    a[i] = ((float)random() / (float)RAND_MAX);
    b[i] = ((float)random() / (float)RAND_MAX);
  }
  prepare(a, b, c);

  printf("start parallel calculation.\n");

  for (i = 0; i < NCBFUNC; ++i) {
    clock_gettime(CLOCK_REALTIME, &start);
    printf("<<%s>>\n", cbtitle[i]);
    (*cbfunc[i])();
    clock_gettime(CLOCK_REALTIME, &end);
    printf("elapsed %f seconds\n\n", elapsed(start, end));
  }
  

  for (i = 0; i < 10 && i < NELEM; ++i) {
    printf("%f ", c[i]);
  }
  fputs("\n", stdout);

  f = getResult(c);
  for (i = 0; i < 10 && i < NELEM; ++i) {
    printf("%f ", f[i]);
  }
  fputs("\n", stdout);
 
  _mm_free(a);
  _mm_free(b);
  _mm_free(c);
  release();
  return 0;
}

[sample_cuda.cu]

#include 
#include 
#include "Nelem.h"

#define NTHREAD (256)

float *pa, *pb, *pc;

__global__
void kernel(float *xc, float *a, float *b) {
  int i = blockIdx.x * blockDim.x + threadIdx.x;
  if (i < NELEM) {
    xc[i] = a[i] + b[i];
  }
  __syncthreads();
}

__host__
void invoke_kernel(float *c, float *a, float *b) {
  int block = NELEM / NTHREAD;
  kernel<<<block, NTHREAD>>>(c, a, b);
}

extern "C" void prepare(float *a, float *b, float *c) {
  cudaMalloc((void**)&pa, sizeof(float) * NELEM);
  cudaMalloc((void**)&pb, sizeof(float) * NELEM);
  cudaMalloc((void**)&pc, sizeof(float) * NELEM);
  cudaMemcpy(pa, a, sizeof(float) * NELEM, cudaMemcpyHostToDevice);
  cudaMemcpy(pb, b, sizeof(float) * NELEM, cudaMemcpyHostToDevice);
}

extern "C" float *getResult(float *c) {
  cudaMemcpy(c, pc, sizeof(float) * NELEM, cudaMemcpyDeviceToHost);
  return c;
}

extern "C" void release(void) {
  cudaFree(pa);
  cudaFree(pb);
  cudaFree(pc);
}

extern "C" void wrapper_for_cuda(void) {
  invoke_kernel(pc, pa, pb);
  cudaDeviceSynchronize();
}

[Makefile]

SRC=sample.c
OBJ=$(SRC:.c=.o)
CUDASRC=sample_cuda.cu
CUDAOBJ=$(CUDASRC:.cu=.o)
TARGET_OPENMP=sample

all: $(TARGET_OPENMP)

$(TARGET_OPENMP): $(OBJ) $(CUDAOBJ)
        nvcc -g -O0 -o $(TARGET_OPENMP) --gpu-architecture=compute_61 --gpu-code=sm_61 --compiler-options="-fopenmp -DAVX2 -march=native" $(OBJ) $(CUDAOBJ)

$(CUDAOBJ): $(CUDASRC)
        nvcc --gpu-architecture=compute_61 --gpu-code=sm_61 -c $(CUDASRC)

.c.o:
        gcc -c -g -O0 -fopenmp -DAVX2 -march=native $<

clean:
        -rm $(TARGET_OPENMP) $(OBJ) $(CUDAOBJ)

 

Anaconda の Numba を CUDA 9.0 で使う

Volta 対応の CUDA 9.0 が出ましたが、Anaconda の cudatoolkit は 8.0 までしか対応していないため、そのままでは Python  の Numba パッケージが使えません。

というわけで、悪戦苦闘した結果を書いておきます。

環境変数の設定

Numba を利用するためには以下の環境変数を設定する必要があります。

  • NUMBAPRO_NVVM
  • NUMBAPRO_LIBDEVICE

CUDA 9.0 をインストールすると /usr/local/cuda/nvvm/lib64/libnvvm.so.* というライブラリがが作成されますので、このパスを NUMBAPRO_NVVM に設定します。(ファイル名まで含む。)

また、/usr/local/cuda/nvvm/libdevice の下に libdevice.10.bc というファイルがありますので、シンボリックリンクで
libdevice.compute_XX.10.bc というリンクを貼っておきます。(XX は使用している GPU の Compute Capability を指定します。例えば GTX 1080 なら compute_61 となります。)NUMBAPRO_LIBDEVICE には /usr/local/cuda/nvvm/libdevice を指定します。

compute_XX の追加

Anaconda 5.0.0 の Numba には Compute Capability が 5.0 までしか登録されていないため、GPU によっては実行エラーが出てしまいます。ちょっと汚いやり方ですが、Numba パッケージに新しい Compute Capability を追加する方法で回避しました。

GTX 1080 用なら、具体的には ~/anaconda3/lib/python3.6/site-packages/numba/cuda/cudadrv/nvvm.py に以下の2点を追加します。

変数やクラス 修正内容
SUPPORTED_CC リストの最後に (6, 1) を追加
LibDevice クラスの _known_arch 変数 “compute_61” をリストの最後に追加

nVidia のデモプログラムへの修正

NVidia のサイトで説明しているサンプルは numbapro を使ったものなので、numba 用に修正が必要となります。

import numpy as np
from numpy.random import rand
from numba import vectorize

@vectorize(['float32(float32, float32)'], target='cuda')
def VectorAdd(a, b):
    return a + b

def main():
    N = 100000000 # 100M
    print("wait for initialization...")
    A = rand(N).astype(np.float32)
    B = rand(N).astype(np.float32)
    C = np.zeros(N, dtype=np.float32)
    print("OK. Let's start!")

    C = VectorAdd(A, B)
    print("C[:5] = " + str(C[:5]))

if __name__ == "__main__":
    main()

ちなみに、Anaconda3 は Intel の MKL を抱えているっぽいので、この程度の演算だと GPU の方が CPU より遅いという結果となりました。PCIe 上のメモリ転送がボトルネックを引き起こしているためと思われます。

デバイスメモリに結果をキャッシュするような例なら速くなるはず。