如何在 CUDA 中为多个 Thrust 对象成员函数调用内核函数？

Question

请注意，我是 CUDA 的绝对初学者，下面的所有内容都是未经测试的伪代码。我来自JavaScript，我的C++也超级生疏，所以我为我的无知道歉:)

我正在尝试使用 CUDA 回测许多不同的外汇策略。

使用 Thrust，我从 class（伪代码）中实例化了 1000 个对象：

#include <stdio.h>
#include <thrust/device_ptr.h>
#include <thrust/device_new.h>

#define N 1000

typedef struct dataPoint {
    ...
} dataPoint;

class Strategy {
    public:
        __device__ __host__ Strategy(...) {
            ...
        }

        __device__ __host__ void backtest(dataPoint data) {
            ...
        }
};

int main() {
    dataPoint data[100000];
    thrust::device_ptr<Strategy> strategies[1000];
    int i;

    // Instantiate 5000 strategies.
    for (i=0; i<1000; i++) {
        thrust::device_ptr<Strategy> strategies[i] = thrust::device_new<Strategy>(...);
    }

    // Iterate over all 100000 data points.
    for (i=0; i<100000; i++) {
        // Somehow run .backtest(data[j]) on each strategy here.
        // i.e. Run backtest() in parallel for all 1000
        // strategy objects here.
    }
}

现在假设我想运行 data 中每个项目的每个对象的 .backtest() 方法。在程序上我会做以下事情：

// Iterate over all 100000 data points.
for (j=0; j<100000; j++) {
    // Iterate over all 1000 strategies.
    for (i=0; i<1000; i++) {
        strategies[i].backtest(data[j]);
    }
}

我如何使用 CUDA 完成此操作，以便每次迭代 j 通过数据对所有策略并行 .backtest() 运行？

如果我必须完全重新构建所有内容，那就这样吧——我对任何必要的事情都持开放态度。如果 classes 无法做到这一点，那就这样吧。

Answer 1

典型的推力代码经常使用某些 C++ 习语（例如仿函数），因此如果您的 C++ 生锈了，您可能想阅读有关 C++ 仿函数的信息。您可能还想查看 thrust quick start guide 以讨论函子以及我们目前将使用的奇特迭代器。

总的来说，至少从表达式的角度来看，我认为推力很适合你的问题描述。考虑到针对这些类型问题的推力表达的灵活性，可能有很多方法可以解决问题。我将尝试向您的伪代码展示一些 "close" 的东西。但毫无疑问，有很多方法可以实现这一点。

首先，在推力上，我们通常会尽量避免 for 循环。这些将非常缓慢，因为它们通常在每次迭代时涉及主机代码和设备代码交互（例如，在每次迭代时调用 CUDA 内核，在引擎盖下）。如果可能的话，我们更喜欢使用推力算法，因为这些算法通常 "translate" 到一个或少量的 CUDA 内核，在引擎盖下。

thrust 中最基本的算法之一是 transform。它有多种风格，但基本上是获取输入数据并对其应用任意操作，逐个元素。

使用基本的推力变换操作，我们可以初始化您的数据和策略，而无需诉诸 for 循环。我们将为每种类型的对象（dataPoint、Strategy）构造一个适当长度的设备向量，然后我们将使用thrust::transform初始化每个向量。

这给我们留下了针对每个 Strategy 执行每个 dataPoint 的任务。理想情况下，我们也想并行执行它；不仅针对您建议的每个 for 循环迭代，而且针对 every Strategy 与 every dataPoint，所有 "at once"（即在单个算法调用中）。

实际上，我们可以考虑一个矩阵，一个轴由 dataPoint（在您的示例中为 100000 维）组成，另一个轴由 Strategy（在您的示例中为 1000 维）组成).对于这个矩阵中的每个点，我们设想它保存 Strategy 对 dataPoint.

的应用结果

总的来说，我们往往更喜欢将二维的概念理解为一维的。因此我们的结果 space 等于 dataPoint 的数量乘以 Strategy 的数量。我们将创建一个 result device_vector 这个大小（在您的示例中为 100000*1000）来保存结果。

为了演示，由于您几乎没有提供有关您想要执行的算术类型的指导，我们假设如下：

对dataPoint应用Strategy的结果是float。
dataPoint 是一个由 int（dtype - 本例忽略）和 float（dval）组成的结构。 dval 将包含 dataPoint(i)、1.0f + i*10.0f.
Strategy由一个multiplier和一个adder组成，用法如下：
```
Strategy(i) = multiplier(i) * dval + adder(i);
```
对 dataPoint 应用 Strategy 包括检索与 dataPoint 关联的 dval，并将其代入给出的等式上面的第 3 项。该等式在 class Strategy 的 backtest 方法中被捕获。 backtest 方法将类型为 dataPoint 的对象作为其参数，它将从中检索适当的 dval.

我们还需要介绍几个概念。二维结果矩阵的一维实现需要我们提供一个合适的索引手段，这样在二维矩阵中的每个点，给定它的线性维度，我们可以确定哪个Strategy和哪个dataPoint将用于计算该点的 result。总之，我们可以使用花式迭代器的组合来做到这一点。

简而言之，从 "inside out" 开始，我们将从一个变换迭代器开始，它采用索引映射仿函数和 thrust::counting_iterator 提供的线性序列，以便为每个索引（每个矩阵维度）。每个映射仿函数中的算法会将 result 的线性索引转换为矩阵的行和列的适当重复索引。给定此转换迭代器来创建重复的行或列索引，我们将该索引传递给置换迭代器，该置换迭代器为每个指示的 row/column 选择适当的 dataPoint 或 Strategy。然后将这两项 (dataPoint、Strategy) 压缩到 zip_iterator 中。然后将 zip_iterator 传递给 run_strat 仿函数，后者实际计算应用于给定 dataPoint.

的给定 Strategy

下面是概述上述概念的示例代码：

#include <iostream>
#include <thrust/device_vector.h>
#include <thrust/host_vector.h>
#include <thrust/transform.h>
#include <thrust/iterator/counting_iterator.h>
#include <thrust/iterator/permutation_iterator.h>
#include <thrust/iterator/zip_iterator.h>
#include <math.h>

#define TOL 0.00001f


// number of strategies
#define N 1000
// number of data points
#define DSIZE 100000

// could use int instead of size_t here, for these problem dimensions
typedef size_t idx_t;


struct dataPoint {

  int dtype;
  float dval;
};

class Strategy {

    float multiplier;
    float adder;
    idx_t id;

    public:

        __device__ __host__ Strategy(){
          id = 0;
          multiplier = 0.0f;
          adder = 0.0f;
          }
        __device__ __host__ Strategy(idx_t _id) {
          id = _id;
          multiplier = 1.0f + ((float)id)/(float)N;
          adder = (float)id;
        }

        __device__ __host__ float backtest(dataPoint data) {
          return multiplier*data.dval+adder;
        }
};

// functor to initialize dataPoint
struct data_init
{
  __host__ __device__
  dataPoint operator()(idx_t id){
    dataPoint temp;
    temp.dtype = id;
    temp.dval = 1.0f + id * 10.0f;
    return temp;
    }
};

// functor to initialize Strategy
struct strat_init
{
  __host__ __device__
  Strategy operator()(idx_t id){
    Strategy temp(id);
    return temp;
    }
};

// functor to "test" a Strategy against a dataPoint, using backtest method
struct run_strat
{
  template <typename T>
  __host__ __device__
  float operator()(idx_t id, T t){
    return (thrust::get<0>(t)).backtest(thrust::get<1>(t));
    }
};

// mapping functor to generate "row" (Strategy) index from linear index
struct strat_mapper : public thrust::unary_function<idx_t, idx_t>
{
  __host__ __device__
  idx_t operator()(idx_t id){
    return id/DSIZE;
    }
};

// mapping functor to generate "column" (dataPoint) index from linear index
struct data_mapper : public thrust::unary_function<idx_t, idx_t>
{
  __host__ __device__
  idx_t operator()(idx_t id){
    return id%DSIZE;
    }
};



int main() {
    // initialize data
    thrust::device_vector<dataPoint> data(DSIZE);
    thrust::transform(thrust::counting_iterator<idx_t>(0), thrust::counting_iterator<idx_t>(DSIZE), data.begin(), data_init());

    // initialize strategies
    thrust::device_vector<Strategy> strategies(N);
    thrust::transform(thrust::counting_iterator<idx_t>(0), thrust::counting_iterator<idx_t>(N), strategies.begin(), strat_init());

    // test each data point against each strategy

        // Somehow run .backtest(data[j]) on each strategy here.
        // i.e. Run backtest() in parallel for all 1000
        // strategy objects here.

    // allocate space for results for each datapoint against each strategy
    thrust::device_vector<float> result(DSIZE*N);
    thrust::transform(thrust::counting_iterator<idx_t>(0), thrust::counting_iterator<idx_t>(DSIZE*N), thrust::make_zip_iterator(thrust::make_tuple(thrust::make_permutation_iterator(strategies.begin(), thrust::make_transform_iterator(thrust::counting_iterator<idx_t>(0), strat_mapper())), thrust::make_permutation_iterator(data.begin(), thrust::make_transform_iterator(thrust::counting_iterator<idx_t>(0), data_mapper())))), result.begin(), run_strat());

    // validation
    // this would have to be changed if you change initialization of dataPoint
    // or Strategy
    thrust::host_vector<float> h_result = result;
    for (int j = 0; j < N; j++){
      float m =  1.0f + (float)j/(float)N;
      float a = j;
      for (int i = 0; i < DSIZE; i++){
        float d =  1.0f + i*10.0f;
        if (fabsf(h_result[j*DSIZE+i] - (m*d+a))/(m*d+a) > TOL) {std::cout << "mismatch at: " << i << "," << j << " was: " << h_result[j*DSIZE+i] << " should be: " << m*d+a << std::endl; return 1;}}}
    return 0;
}

备注：

如前所述，这是一种可能的实现方式。我认为它应该是 "reasonably" 高效的，但在 thrust 中可能会有更高效的实现。在尝试解决优化问题之前，对您的实际策略和回测方法进行更完整的分析可能是合适的。
最后的 transform 操作使用 counting_iterator 作为第一个参数（和第二个），但这实际上被忽略了 "dummy" 用法，只是为了大小适当的问题。它可以通过更简单的实现来消除，但在我看来，最简单的方法（不会进一步混淆代码）是使用 C++11 auto 来定义 zip_iterator，然后通过它本身加上它的一个偏移版本，thrust::transform，使用只接受一个输入向量而不是 2 个输入向量的版本。我认为这不会对性能产生太大影响，我觉得这有点更容易解析，但也许不是。

如何在 CUDA 中为多个 Thrust 对象成员函数调用内核函数？

How can I call a kernel function for multiple Thrust object member functions in CUDA?

c++

cuda

nvidia

thrust