使用 Random123 的 OpenCL 随机数

Question

我一直在查看这个库 Random123 和相关引用：

一位神秘人来到我的展位，询问我对使用 OpenCL 生成随机数的了解。我告诉他有关 Mersenne Twister 的实现，但他没有留下深刻的印象。他告诉我一篇新的技术论文，该论文解释了如何通过结合整数计数器和分组密码在 GPU 上生成随机数。他恭敬地说，基于计数器的随机数生成器 (CBRNG) 生成的数字比 MT 具有更大的统计随机性和更快的速度。

我能够使用这个内核运行一个演示：

__kernel void counthits(unsigned n, __global uint2 *hitsp) {
    unsigned tid = get_global_id(0);
    unsigned hits = 0, tries = 0;
    threefry4x32_key_t k = {{tid, 0xdecafbad, 0xfacebead, 0x12345678}};
    threefry4x32_ctr_t c = {{0, 0xf00dcafe, 0xdeadbeef, 0xbeeff00d}};
    while (tries < n) {
        union {
            threefry4x32_ctr_t c;
            int4 i;
        } u;
        c.v[0]++;
        u.c = threefry4x32(c, k);
        long x1 = u.i.x, y1 = u.i.y;
        long x2 = u.i.z, y2 = u.i.w;
        if ((x1*x1 + y1*y1) < (1L<<62)) {
            hits++;
        }
        tries++;
        if ((x2*x2 + y2*y2) < (1L<<62)) {
            hits++;
        }
        tries++;
    }
    hitsp[tid].x = hits;
    hitsp[tid].y = tries;
}

我现在的问题是，这不会在每次运行时生成相同的随机数，随机数是基于全局 id 的吗？如何每次生成新的随机数。可以提供种子作为内核的参数，然后以某种方式使用它吗？

有谁一直在使用这个库，可以给我一些关于它的使用的更多见解吗？

Answer 1

是的。示例代码每次调用时都会生成相同的随机数序列。

要获得不同的随机数流，只需以不同方式初始化任何值 k[1..3] 和/或 c[1..3]。您可以从命令行参数、环境变量、时间、保存状态、/dev/urandom 或任何其他来源初始化它们。请注意：

a) 如果您在两次不同的运行中以完全相同的方式初始化它们，那么这两次运行将获得相同的随机数流

b) 如果您在两次不同的运行中以不同方式初始化它们，那么这两次运行将获得不同的随机数流。

有时您需要属性 a)。有时你想要属性 b)。花点时间想想你想要什么，并确保你正在做你想做的事。

更一般地说，库中的函数，例如threefry4x32，没有无状态。如果您更改输入中的任何位（即 c 或 k 的任何元素中的任何位），您将获得完全不同的随机、统计独立、均匀分布的输出。

附：我是该库和论文“Parallel Numbers: As Easy as 1, 2, 3”的作者之一： http://dl.acm.org/citation.cfm?id=2063405

如果您不是 ACM 数字图书馆的订阅者，上面的链接可能会遇到付费问题。或者，您可以通过点击此页面上的链接免费获取论文：

http://www.thesalmons.org/john/random123/index.html

Answer 2

库本身我无法帮助您，但我可以告诉您，在 OpenCL 中生成随机数的最常见方法是在调用内核之间保存一些状态。

随机数生成器通常使用一个状态，从中生成一个新状态和一个随机数。实际上，这并不复杂：您只需传递一个保存状态的额外数组。在我的代码中，我实现了如下随机数：

uint rand_uint(uint2* rvec) {  //Adapted from http://cas.ee.ic.ac.uk/people/dt10/research/rngs-gpu-mwc64x.html
    #define A 4294883355U
    uint x=rvec->x, c=rvec->y; //Unpack the state
    uint res = x ^ c;          //Calculate the result
    uint hi = mul_hi(x,A);     //Step the RNG
    x = x*A + c;
    c = hi + (x<c);
    *rvec = (uint2)(x,c);      //Pack the state back up
    return res;                //Return the next result
    #undef A
}
inline float rand_float(uint2* rvec) {
    return (float)(rand_uint(rvec)) / (float)(0xFFFFFFFF);
}
__kernel void my_kernel(/*more arguments*/ __global uint2* randoms) {
    int index = get_global_id(0);
    uint2 rvec = randoms[index];

    //Call rand_uint or rand_float a number of times with "rvec" as argument.
    //These calls update "rvec" with new state, and return a random number

    randoms[index] = rvec;
}

。 . .然后，您所做的就是传递一个额外的数组，将 RNG 的状态保存为随机数。在实践中，您需要为每个工作项以不同的方式播种此数组。

Answer 3

0xdecafbad、0xfacebead、0x12345678 和 0xf00dcafe、0xdeadbeef、0xbeeff00d 只是任意选择的数字，它们并不特殊。可以使用任何其他数字（甚至 0）代替它们——我将在示例代码中添加注释。

你可以用你传入的变量替换它们中的任何一个；避免在输出随机“流”中出现不必要的重复的唯一要求是避免重复 (c, k) 输入元组。示例代码使用线程 ID 和循环索引来确保唯一性，但您可以轻松添加更多变量以确保唯一性——例如计算主机代码中的内核调用并将该计数器传入，用它代替 k 或 c 的元素之一。

顺便说一句，尽管名称为“基于计数器的随机数生成器”，但并不要求输入 (c, k) 是“计数器”，只是计数器恰好是确保输入的最方便的习惯用法不要重复。