CUDA rand (use MT19937)

For generator random number every time, we must put a seed from CPU. (eg : time(0))

cuda_kernel<<<grid,block>>>(time(0));

And in the parallel thread, we also use a offset value to generator difference random number in the same time.

c__global__ void cuda_kernel(int seed)
{
  int index= (blockIdx.y*blockDim.y+threadIdx.y)*(gridDim.x*blockDim.x)+(blockIdx.x*blockDim.x+threadIdx.x);
  MTRand mtrand(seed+index);  //seed + offset:index
  printf("%5.2f",mtrand.randf(0,10)");
}

MTRand Class code:

class MTRand
{
/*
  MT19937
  source: http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c
*/
  private:
#define __MTRand_N__ 624
#define __MTRand_M__ 397
#define __MTRand_MATRIX_A__ 0x9908b0dfUL
#define __MTRand_UPPER_MASK__ 0x80000000UL
#define __MTRand_LOWER_MASK__ 0x7fffffffUL
    unsigned long mt[__MTRand_N__];
    unsigned long mag01[2];
    int mti;
    __host__ __device__
    void init_genrand(unsigned long s)
    {
      mt[0]= s & 0xffffffffUL;
      for (mti=1; mti<__MTRand_N__; mti++) {
        mt[mti] = (1812433253UL * (mt[mti-1] ^ (mt[mti-1] >> 30)) + mti);
        mt[mti] &= 0xffffffffUL;
      }
    }
    __host__ __device__
    unsigned long genrand_int32(void)
    {
      unsigned long y;
      if (mti >= __MTRand_N__) {
        int kk;
        if (mti == __MTRand_N__+1)  init_genrand(5489UL);
        for (kk=0;kk<__MTRand_N__-__MTRand_M__;kk++) {
          y = (mt[kk]&__MTRand_UPPER_MASK__)|(mt[kk+1]&__MTRand_LOWER_MASK__);
          mt[kk] = mt[kk+__MTRand_M__] ^ (y >> 1) ^ mag01[y & 0x1UL];
        }
        for (;kk<__MTRand_N__-1;kk++) {
          y = (mt[kk]&__MTRand_UPPER_MASK__)|(mt[kk+1]&__MTRand_LOWER_MASK__);
          mt[kk] = mt[kk+(__MTRand_M__-__MTRand_N__)] ^ (y >> 1) ^ mag01[y & 0x1UL];
        }
        y = (mt[__MTRand_N__-1]&__MTRand_UPPER_MASK__)|(mt[0]&__MTRand_LOWER_MASK__);
        mt[__MTRand_N__-1] = mt[__MTRand_M__-1] ^ (y >> 1) ^ mag01[y & 0x1UL];
        mti = 0;
      }
      y = mt[mti++];
      y ^= (y >> 11);
      y ^= (y << 7) & 0x9d2c5680UL;
      y ^= (y << 15) & 0xefc60000UL;
      y ^= (y >> 18);
      return y;
    }
  public:
    __host__ __device__
    MTRand(int seed)
      :mti(__MTRand_N__+1)
    {
      mag01[0]=0x0UL;
      mag01[1]=__MTRand_MATRIX_A__;
      init_genrand(19650218UL+seed);
    }
    __host__ __device__
    unsigned long rand()
    {
      return genrand_int32();
    }
    __host__ __device__
    float randf()
    {
      return genrand_int32()*(1.0/4294967295.0);
    }
    __host__ __device__
    float randf(float min,float max)
    {
      if(max<min) { float t=min; min=max; max=t; }
      return randf()*(max-min)+min;
    }
};

發佈留言