【问题标题】:CUDA - Generating the Halton sequence in parallelCUDA - 并行生成 Halton 序列
【发布时间】:2015-04-10 03:59:23
【问题描述】:

我想在 CUDA 中编写一个内核,该内核将并行生成 the Halton sequence,每个线程生成并存储 1 个值。

查看序列,似乎生成序列中的每个后续值都涉及在生成前一个值时所做的工作。从头开始生成每个值将涉及冗余工作并导致线程执行时间之间存在很大差距。

有什么方法可以通过改进串行算法的并行内核来做到这一点?我对并行编程真的很陌生,所以如果答案是一些众所周知的模式,请原谅这个问题。

注意:我确实在教科书中找到了this link(它使用它但没有描述它是如何工作的)但是那里的文件链接已经失效。

【问题讨论】:

  • 您可能希望查看 CURAND 中的 Sobol 准随机生成器作为替代方案。如果 Halton 序列的递归可表示为变换矩阵和状态向量之间的矩阵向量乘法,则可以通过让线程 i 在 O(log(n)) 时间内计算 M**i 来并行化它,然后执行矩阵-向量乘法生成第 i 个状态向量。您还可以预先计算 M 的幂以减少冗余计算。

标签: c++ cuda parallel-processing montecarlo


【解决方案1】:

Halton 序列由以下方式生成:

  1. 获取 i 在 base-p 数字系统中的表示
  2. 反转位顺序

例如,base-2 Halton 序列:

index      binary     reversed     result
1             1           1           1 /   10 = 1 / 2
2            10          01          01 /  100 = 1 / 4
3            11          11          11 /  100 = 3 / 4
4           100         001         001 / 1000 = 1 / 8
5           101         101         101 / 1000 = 5 / 8
6           110         011         011 / 1000 = 3 / 8
7           111         111         111 / 1000 = 7 / 8

所以按位反转确实有很多重复的工作。我们可以做的第一件事就是重用以前的结果。

在计算base-p Halton序列中索引为i的元素时,我们首先确定i的前导位和base-p表示的剩余部分(这可以通过以base-p方式调度线程来完成) .然后我们有

out[i] = out[remaining_part] + leading_bit / p^(length_of_i_in_base_p_representation - 1)
//"^" is used for convenience

为避免不必要的全局内存读取,每个线程应处理具有相同“剩余部分”但不同“前导位”的所有元素。 如果我们在 p^n 和 p^(n+1) 之间生成 Halton 序列,则在概念上应该有 p^n 个并行任务。但是,如果我们为一个线程分配一组任务,这不会造成任何问题。

可以通过混合重新计算和从内存加载来进一步优化。

示例代码:

总线程数应为 p^m。

const int m = 3 //any value
__device__ void halton(float* out, int p, int N)
{
    const int tid = ... //globally unique and continuous thread id
    const int step = p^m; //you know what I mean
    int w = step; //w is the weight of the leading bit
    for(int n = m; n <= N; ++n) //n is the position of the leading bit
    {
        for(int r = tid; r < w; r += step) //r is the remaining part
        for(int l = 1; l < p; ++l) //l is the leading bit
            out[l*w + r] = out[r] + l/w;
        w *= p;
    }
}

注意:此示例不计算 Halton 序列中的前 p^m 个元素,但仍需要这些值。

【讨论】:

  • "这可以通过以 base-p 方式调度线程来完成" 你能解释一下吗?
  • p^(length_of_i_in_base_p_representation + 1) 不应该是 -1 吗?
  • 我已经更正了。 “以base-p方式调度线程”意味着线程数和循环周期是p的一些顺序。这个例子应该能说明这一点。
猜你喜欢
  • 2013-07-13
  • 2013-11-25
  • 2012-10-30
  • 2012-09-21
  • 2014-01-09
  • 2021-04-11
  • 1970-01-01
  • 2013-02-21
  • 1970-01-01
相关资源
最近更新 更多