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

以前、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 上のメモリ転送がボトルネックを引き起こしているためと思われます。

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

unordered_map から unordered_set は作れないのか?

連想配列を使っていると、キーの集合だけ、値の集合だけを作成したいことがある。そんな時はどうするのがスマートなのだろうか。

1.for 文で地道に拾う

C++ っぽくないので却下。

2.for_each 式+ラムダ式で拾う

using namespace std;

void separate(unordered_set<MapKey> key_set, unordered_key<MapValue> value_set, const unordered_map<MapKey, MapValue> map) {
  for_each(begin(map), end(map), [&key_set, &key_value] (const pair<MapKey, MapValue> &i) { key_set.insert(i->first), value_set.insert(i->second) });
}

Java と違って標準の抜き出しメソッドを用意していないのは set とか他のクラスに実装する可能性を考えているからだろうか。

どっちが正解??

スクリプト処理系で配列変数を定義するにあたり、どちらの実装にするべきか悩んでいる。

  1. 演算子 “[” を定義し、第1引数にベース、第2引数にキーを設定する。
  2. 変数に添え字配列を用意する。

1.の利点は固定数の引数に対して演算子を定義すればよいことで、実装が簡単になる。その分、実行エンジンは遅くなりますが。

2.はイメージとしては

class Content {
  Content *c;
  vector <Content *> pc;
};

として、pc が empty でない場合は連想配列の引数とみなすようにする。

2.の利点は実装が多少煩雑になるが、多重配列を定義する場合の実行速度を上げることができる。(ツリーをたどらず、O(1) の配列要素が利用できるので。)

最初1、の実装で考えていたが、途中で考え直して2.の実装に切り替えた。(コード上はその残骸が残っていますが)ただし、まだ変数を左辺値に指定すると NullException が発生するというバグがある。

久々にスクリプトエンジンのコードを書く

半年ぶりにスクリプトの実行エンジンのコードを書いた。

lex, yacc に相当するところは既に書き終えており、式の評価以外の実行エンジンも書いてあるので、残るは数式の評価である。

しかし半年もたってしまうとどういう構造だったかとか、このメンバーが何を意図していたかとかが分からなくなってしまい、困った。

何とか +, 変数代入、ハッシュが書けたので、後は淡々と定義した演算子の実体を定義することになる。この作業、本当に地道で根気のいる作業なのです。

PC 見積もり代金 530万円也

hp のサイトで、以下のスペックで PC の見積もりを見てみた。

  • CPU: Xeon 2620v4 2.1G 8コア×2
  • RAM: 512GB
  • HDD: 2TB ×2
  • GPU: P6000 ×2
  • その他もろもろ(Blue-ray, GbE追加など)

結果は標題にある通り、530万円余りということになった。
(CPU には正直あまり期待していないので、必要最小限のオプションになっています。)

ちなみにさくらインターネットの高火力シリーズでP100×1を1年間の見積もりを計算すると363万円(ただしメモリ1TB)、マウスコンピュータのX299ベース、128GBのPCを4台見積もると440万円くらいになるので、電源事情を考えれば、2年以上使うのだったら hp から購入するのは悪くない選択肢に思える。
(Ubuntu のホームページを見ると Z840 筐体はサポートされているようだ)

納期が8週間かかりそうなのと、11月下旬に新しい筐体が出るようなので、今回は購入を見送ります。というか今使っている PC を取りあえずフルに使えと言うことなんですが。

N-gram 検索エンジンに and 機能を追加

N-gram 検索エンジンに and 検索機能を追加しました。

半角の空白文字を入力すると、空白で区切られた文字列を全て含むページを検索します。

実装がちょっとうまくないですが、おいおい直します。

表示についても統一化していないため、検索にヒットした文の数だけ結果が表示されます。例えば、「可逆 圧縮」と入力した場合、「可逆」で検索された後、「圧縮」で検索されるため、「可逆」で1行、「圧縮」で1行表示されます。

これは検索の問題というよりは表示の問題です。