CUDA 4.1の新機能(2): CURANDの新アルゴリズム

出典: トータル・ディスクロージャ・サイト(事実をありのままに)

モンテカルロ法などに効果を発揮する、CUDA Toolkit付属の乱数ライブラリ、CURANDCUDA 4.1では、このCURANDに、MRG32k3aMersenne Twisterの、2種類のアルゴリズムが新たに加わった。そこで、今回は、これらの新アルゴリズムを試してみた。

目次

テストの内容

モンテカルロ法を用いて、円周率の値を求めるプログラムを、CURANDの各種アルゴリズムで実装し、実行時間と精度の面から比較を行う。

CUDA Cプログラム

#include <cstdio>

#include <sys/time.h>

#include <curand.h>
#include <thrust/device_ptr.h>
#include <thrust/iterator/counting_iterator.h>
#include <thrust/iterator/transform_iterator.h>
#include <thrust/reduce.h>

struct func : thrust::unary_function<int, int> {
  double *d_a;
  func(double *d_a_) : d_a(d_a_)
  {}
  __device__ int
  operator()(int i) const
  {
    double x = d_a[i * 2 + 0];
    double y = d_a[i * 2 + 1];
    double r2 = x * x + y * y;
    return (r2 <= 1.0) ? 1 : 0;
  }
};

void
benchmark(const size_t &n, double *d_a, const curandRngType_t &rngtype,
          const char *name)
{
  struct timeval tv0, tv1;
  gettimeofday(&tv0, NULL);
  curandGenerator_t generator;
  curandCreateGenerator(&generator, rngtype);
  curandGenerateUniformDouble(generator, d_a, n * 2);
  thrust::counting_iterator<int> count(0);
  thrust::transform_iterator<func, thrust::counting_iterator<int> >
    trans(count, func(d_a));
  int sum = thrust::reduce(trans, trans + n);
  double pi = 4.0 * double(sum) / double(n);
  curandDestroyGenerator(generator);
  gettimeofday(&tv1, NULL);
  double error = std::abs(pi - M_PI) / M_PI;
  int sec = tv1.tv_sec - tv0.tv_sec;
  int usec = tv1.tv_usec - tv0.tv_usec;
  if (usec < 0) {
    --sec;
    usec += 1000000;
  }
  printf("%s (%d Points)\nPI:\t%.70f\nError:\t%.70f\nTime:\t%d.%06d\n",
         name, n, pi, error, sec, usec);
}

int
main()
{
  const size_t n = 300000000;

  double *d_a;
  cudaMalloc((void **)&d_a, n * 2 * sizeof(double));
  benchmark(n, d_a, CURAND_RNG_PSEUDO_XORWOW, "XORWOW");
  benchmark(n, d_a, CURAND_RNG_PSEUDO_MRG32K3A, "MRG32k3a");
  benchmark(n, d_a, CURAND_RNG_PSEUDO_MTGP32, "Mersenne Twister");
  benchmark(n, d_a, CURAND_RNG_QUASI_SOBOL64, "Sobol64");
  benchmark(n, d_a, CURAND_RNG_QUASI_SCRAMBLED_SOBOL64, "Scrambled Sobol64");
  cudaFree(d_a);
  return 0;
}

Thrustテンプレートを使っているため、見慣れないテンプレート関数が出てきているが、XとY座標に、それぞれ0 - 1の倍精度一様乱数を3億要素生成し、二乗和が1以下の場合に1、1を超えた場合に0を返すフィルタを通し、その結果の和を取って、4を掛けて要素数3億で割るプログラムである。後述するように、実行環境がメモリ6 GBのTesla M2090であるので、要素数を3億に設定しているが、メモリの少ないGPUによる実行の場合は、要素数nを少なくして計算することが望ましい。

実験環境

CPU: Six-Core Xeon X5690 3.46 GHz × 2
GPU: Tesla M2090 × 4
OS: CentOS 6.0
CUDA: CUDA 4.1 RC1

実行結果

XORWOW (300000000 Points)
PI:	3.1415300266666665862658192054368555545806884765625000000000000000000000
Error:	0.0000199347687724466886792767245095703287915966939181089401245117187500
Time:	0.110700
MRG32k3a (300000000 Points)
PI:	3.1416060933333334048711549257859587669372558593750000000000000000000000
Error:	0.0000042780032366486873097945371757777621724017080850899219512939453125
Time:	0.255586
Mersenne Twister (300000000 Points)
PI:	3.1414099733333333830387346097268164157867431640625000000000000000000000
Error:	0.0000581489316417232939178394435231211900827474892139434814453125000000
Time:	0.207856
Sobol64 (300000000 Points)
PI:	3.2915026266666664866988867288455367088317871093750000000000000000000000
Error:	0.0477178264679146910132168102336436277255415916442871093750000000000000
Time:	0.152308
Scrambled Sobol64 (300000000 Points)
PI:	3.9815186533333335461293245316483080387115478515625000000000000000000000
Error:	0.2673567493811729312014335846470203250646591186523437500000000000000000
Time:	0.152550

MRG32k3aは、従来手法に比べて、実行時間はかなり長いが、精度はかなり向上している。Mersenne Twisterに関しては、今回は実行時間でも精度でも、従来手法のXORWOWに負けているが、アプリケーションによっては有利な局面も出てくることも考えられる。結局速度面では従来手法のXORWOWが頭一つ抜けていて、誤差に関しても、今回の計算の標準偏差はsqrt((π/4×422)/(3×108))≒9×10-5であり、XORWOW、MRG32k3a、MTの3種類のアルゴリズムの誤差はこの範囲に収まっているため、残念ながら今回の試行回数では有意と考えられる差は見いだせなかった。

まとめ

CUDA 4.1に実装された、CURANDの新アルゴリズムを使うことによって、より精度の高いモンテカルロ法が、実現可能になると考えられる。今回の実験では結局XORWOWで十分であると考えられるが、より大規模な計算大規模なクラスタで実行する場合においては、2192-232周期のXORWOWの手に余ることは十分に考えられ、そういう用途のために選択肢を用意したのは興味深い点であろう。


この記事へのコメントをお寄せください

  • サイトへの書き込みに差し支えございましたら トータルディスクロージャーサイトサポート係へメールをお送りください
  • トータル・ディスクロージャ・サイトに投稿された文章と画像は、すべてその著作権がHPCシステムズ株式会社に帰属し、HPCシステムズ株式会社が著作権を所有することに同意してください。
  • あなたの文章が他人によって自由に編集、配布されることを望まない場合は、投稿を控えてください。
  • コメントを書き込む場合は名前にひらがなを織り交ぜてください。
  • あなたの投稿する文章と画像はあなた自身によって書かれたものであるか、パブリック・ドメインかそれに類する自由なリソースからの複製であることを約束してください。あなたが著作権を保持していない作品を許諾なしに投稿してはいけません!

<comments hideform="false" />


Comments

ノート:CUDA 4.1の新機能(2): CURANDの新アルゴリズム

個人用ツール