cuda - CUDA 中的原子 Saxpy

Question

我在 CUDA 中遇到以下问题。

假设我们有一个索引列表，其中一些或所有索引可以出现多次:

inds = [1, 1, 1, 2, 2, 3, 4]

使用这些索引，我想对 float 组 x 执行原子 saxpy 操作(并行)。我不担心应用操作的顺序。也就是说，我想为 float a 和 k 这样做:

x[i] = x[i]*a + k;

如果 inds 中没有重复索引，这将是微不足道的。

我目前的解决方案(不起作用)是这样的:

// assume all values in adr are greater than or equal to 0.
// also assume a and k are strictly positive.

__device__ inline void atomicSaxpy(float *adr,
                                   const float a,
                                   const float k){

  float old = atomicExch(adr, -1.0f); // first exchange
  float new_;
  if (old <= -1.0f){
    new_ = -1.0f;
  } else {
    new_ = old*a + k;
  }

  while (true) {
    old = atomicExch(adr, new_); // second exchange
    if (old <= -1.0f){
      break;
    }
    new_ = old*a + k;
  }
}

在许多情况下，这似乎会返回正确答案。

以下是我认为当您没有得到正确答案时会发生的情况:

1. old 在第一次交换中得到 -1.0f 的值。 => new_ = -1.0f
2. old 在第二次交换中也得到 -1.0f 的值。
3. 函数退出时完全没有任何外部影响。

有点不同的方法是:

__device__ inline void atomicSaxpy(float *adr,
                                   const float ia,
                                   const float k){

  float val;

  while (true) {
    val = atomicExch(adr, -1.0f);
    if (val > 1.0f){
      break;
    }
    atomicExch(adr, val*ia + k);
  }
}

在我的机器上一直死锁。即使是非常简单的输入，例如上面的示例数据。

是否可以重写此函数以使其正常运行？

示例答案

假设k=0.1，a=0.95，args初值为0.5对于所有索引，结果应该是:

[0.5, 0.7139374999999998, 
 0.6462499999999999, 0.575, 0.575, ...]

我使用 Python 计算了这些值，它们在 CUDA 中看起来可能会有所不同。这是算法应该如何运行的示例，而不是遇到竞争条件问题的良好样本集。

引用

这是一个线程，他们使用 atomicExch 实现了 atomicAdd(此时已经存在 float ):

https://devtalk.nvidia.com/default/topic/458062/atomicadd-float-float-atomicmul-float-float-/

一个例子是这样的:

__device__ inline void atomicAdd(float* address, float value) {
  float old = value;  
  float new_old;

  do {
    new_old = atomicExch(address, 0.0f);
    new_old += old;
  }
  while ((old = atomicExch(address, new_old)) != 0.0f);
};

这似乎更容易一些，我不太明白如何适应它。

其他解决方案

能够以这种方式解决这个问题对我以后与内存 IO 相关的问题有几个好处。出于这个原因，我想知道这是否可能。

一种可能的不同方法是计算每个索引在 CPU 上出现的次数，然后在 GPU 上执行“常规”saxpy。我假设还有其他可能性，但我仍然对这个特定问题的答案感兴趣。

Answer 1

如果这是一个非并行问题，您只需这样做:

*adr = *adr * a + k;

由于adr上有多个线程操作，所以我们应该使用原子操作进行读写。

float adrValue = atomicExch(adr, -1.0f)
float newValue = adrValue * a + k
atomicExch(adr, newValue)

但是，我们必须意识到另一个线程在我们的读取步骤 (ln1) 和写入步骤 (ln3) 之间更新了 adr 的可能性。

所以我们这里的三步操作是非原子的。

为了使其成为原子的，我们应该使用比较和交换 (atomicCAS) 来确保我们只更新内存，如果它的值自从我们读取它以来没有改变。我们可以简单地重复我们的步骤，在每次迭代中使用 adr 中当时的当前值作为计算输入，直到 step3 返回预期的锁定值 -1.0f。

do {
    float adrValue = atomicExch(adr, -1.0f)
    float newValue = adrValue * a + k
    adrValue = __int_to_float(atomicCAS(adr, 
                                        __float_as_int(-1.0f),
                                        __float_as_int(newValue)))
} while (adrValue != -1.0f)

ps: 考虑上面的伪代码