Random number generator on GPU

-
Pseudo-random number generation on CPU is simple - we need to store previous state and use it to generate the next number.
On a GPU, it's need more attention, we need to generate random numbers in each instance independently.

That is, we need to keep for each instance a state, randomly initialized on CPU, and then update the state, using corresponding functions.

The good explanation and particular implementation of good random generator is discussed here:

It's about method explained in NVidia's book "GPU Gems", Chapter 37. Efficient Random Number Generation and Application Using CUDA
https://developer.nvidia.com/gpugems/GPUGems3/gpugems3_ch37.html

So, here is my code:

//--------------------------------------------------------------------------
//Random generator for GPU
//Based on
//https://math.stackexchange.com/questions/337782/pseudo-random-number-generation-on-the-gpu#answers-header
//Note that the initial seed values for state should be larger than 128!

//typedef unsigned int uint; //32-bit unsigned integer
//struct uvec4 {
// uint x = 0;
// uint y = 0;
// uint z = 0;
// uint w = 0;
// uvec4() {}
// uvec4(uint x0, uint y0, uint z0, uint w0) {
// x = x0; y = y0; z = z0; w = w0;
// }
//};

//Usage:
/*

//one-time generation, such as particles starting position
//(Note, it works good for N=1..100 particles, for N>100000 - use two-level generation,
//that is create N/1000 seeds and then use them for create per 1000 particles on each, see example below

RandomState rnd;
rnd.state = uvec4(1884,2991,222,334);

for (int i=0; i<N; i++) {
  random2(rnd);

  vec2 pos = rnd.value.xy;
  pos.x *= w;
  pos.y *= h;

  worm[i].pos = vec4(pos,0,1);
}


//Generate many seeds
int N = 100;

Buffer<uvec4> states(N);
...

//Initialization
RandomState rnd;
rnd.state = uvec4(1884,2991,222,334);

//states
for (int i=0; i<N; i++) {
  uvec4 s;
  random4(rnd);
  s.x = int(rnd.value.x*10000)+150;
  s.y = int(rnd.value.y*10000)+150;
  s.z = int(rnd.value.z*10000)+150;
  s.w = int(rnd.value.w*10000)+150;
  states[i] = s;
}

//Using seeds for generating new values
//1) get seed
RandomState rnd;
rnd.state = states[i];

//2) generate values
random2(rnd);
float x = rnd.value.x * w;
float y = rnd.value.y * h;

//3) store updated seed
states[i] = rnd.state;

*/


//--------------------------------------------------------------


struct RandomState
{
  uvec4 state; //16 bytes; Note that the initial seed values for state should be larger than 128!
  vec4 value; //16 bytes
};

//--------------------------------------------------------------------------
//Random generator for GPU

uint TausStep(uint z, int S1, int S2, int S3, uint M)
{
  uint b = (((z << S1) ^ z) >> S2);
  return (((z & M) << S3) ^ b);
}

uint LCGStep(uint z, uint A, uint C)
{
  return (A * z + C);
}

//------------------------------------------------
//Generate 1,2,3,4 random values, stored in value.x,y,z,w
void random1(inout RandomState rnd)
{
  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.x = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

}

void random2(inout RandomState rnd)
{
  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.x = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.y = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

}

void random3(inout RandomState rnd)
{
  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.x = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.y = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.z = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

}

void random4(inout RandomState rnd)
{
  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.x = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.y = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.z = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));

  rnd.state.x = TausStep(rnd.state.x, 13, 19, 12, 4294967294);
  rnd.state.y = TausStep(rnd.state.y, 2, 25, 4, 4294967288);
  rnd.state.z = TausStep(rnd.state.z, 3, 11, 17, 4294967280);
  rnd.state.w = LCGStep(rnd.state.w, 1664525, 1013904223);

  rnd.value.w = float(2.3283064365387e-10 * double(rnd.state.x ^ rnd.state.y ^ rnd.state.z ^ rnd.state.w));
}


Comments

  1. Denis Perevalov23 October 2019 at 23:56

    Actually, seeds initialization can be made in GPU too :)

    ReplyDelete

Post a Comment

Popular posts from this blog

Computing ray origin and direction from Model View Projection matrices for raymarching

-
When performing raymarching (as well a raycasting and raytracing) using fragment shaders in OpenGL, it's required to compute ray's origin and direction. Here we will give solution for finding it from the Model View Projection matrices in OpenGL, implemented in openFrameworks. Using it, you can combine raymarching/raycasting/raytracing with forward OpenGL rendering, so obtain hybrid rendering modes. Also, you can use it for creating VE applications based on raymarching or hybrid rendering. (Actually, we checked what this approach is corrent by creating VR rendering app in HTC Vive). The algorithm can be used on any OpenGL platform with any language supporting matrix inversion computation (in our case it's C++ with GLM library, included in openFrameworks 0.10.1). The approach is based on this hint by  GClements in this discussion on the topic: https://community.khronos.org/t/ray-origin-through-view-and-projection-matrices/7257...
Read more

Create Blueprint Library, print to log and using windows.h in a Unreal Engine C++ project

-
When you need to use some specific Windows API functions, such as LoadLibrary for loading DLL, or CreateFileMapping for using shared memory, you need to use #include <windows.h> in your project's C++ code. To be sure all building successfully, please insert this with "WIN32_LEAN_AND_MEAN" define, and use it as early as possible (I use it as the second #include in CPP code): #include "MyBlueprintFunctionLibrary.h" #define WIN32_LEAN_AND_MEAN #include <windows.h> For advanced capabilites, such as using GEngine class, it's required to include "Engine.h" . It's required to do before including <windows.h> : #include "MyBlueprintFunctionLibrary.h" #include "Engine.h" #define WIN32_LEAN_AND_MEAN #include <windows.h> Here is example how to load DLL using Windows API in Blueprint Function Library C++ class. Also, it demonstrates print to screen and to log in Unreal Engine C+...
Read more

Forward, Deferred and Raytracing rendering in openFrameworks and web

-
Image
Currently, there are three basic approaches for rendering computer graphics scenes: forward rendering, deferred rendering, raytracing. 1. Forward rendering It's that OpenGL does, you create vertices, then make faces, and rasterize them, all passing via shaders. The rendering process is straightworward and fast, but rendering shadows, making degree-of-field effect and drawing clouds is hard here. Most of openFrameworks examples use forward rendering. 2. Deferred rendering It's the same as forward rendering, but additionally some extra data is rendered for computing shadows. The lighting pass is performed after the scene geometry is drawn. See detailed explanation here: https://learnopengl.com/Advanced-Lighting/Deferred-Shading Shadows, degree-of-field effect and some kind of reflections are possible with deferred rendering. Unreal Engine 4 and Unity currently (2019) are based on deferred rendering. openFrameworks: https://github.com/perevalovds...
Read more