使用可重定位设备代码构建 R 包

Build R package with relocatable device code

我正在编写一个 R 包,它使用 Thrust 来处理内存分配并避免编写我自己的 CUDA 内核。

在某些情况下,我从设备代码而不是主机代码调用 cuBLAS 例程。这改变了编译要求。虽然代码使用下面的 nvcc 命令编译,但可能需要显式调用主机 linker (g++)。我怎样才能修改当前的构建过程来完成这个?

我使用的步骤是:

  1. 使用-dc开关

  2. 编译包含设备可重定位代码的输出文件(max.o
  3. 创建一个库(libmax.a)到link
  4. 使用-c开关

  5. 编译包含包装函数(somePackage.o)的输出文件
  6. 使用 -shared 开关

  7. 创建 link 到 libmax.a 的共享库 (somePackage.so)

工作示例如下所示:

iterator.h:定义了一些类型,包括strideAccessor.

max.h: max.cu

中的函数声明

max.cu:定义一个函数,该函数在每个 n 维度的串联数组中查找最大元素的索引 d.

somePackage.cu:处理 R/C++ 接口的包装器

$ cat iterator.h
#ifndef ITER_H
#define ITER_H

#include <thrust/host_vector.h>
#include <thrust/device_vector.h>
#include <thrust/iterator/transform_iterator.h>
#include <thrust/iterator/permutation_iterator.h>
#include <thrust/tuple.h>
#include <thrust/iterator/zip_iterator.h>

typedef thrust::device_vector<int> ivec_d;
typedef thrust::device_vector<double> fvec_d;
typedef thrust::device_vector<int>::iterator intIter;
typedef thrust::device_vector<double>::iterator realIter;
typedef thrust::host_vector<int> ivec_h;
typedef thrust::host_vector<double> fvec_h;

typedef thrust::counting_iterator<int> countIter;

//Used for generating rep( (1:len)*incr, times=infinity)
struct stride: public thrust::unary_function<int, int>{

  int incr;

  __host__ __device__ stride(int incr=1): incr(incr){}

  __host__ __device__ int operator()(int x){

    return x*incr;
  }
};

typedef thrust::transform_iterator<stride, countIter> strideIter;
typedef thrust::permutation_iterator<realIter, strideIter> strideAccessor;


#endif

$ cat max.h
#include "iterator.h"

void cublas_max(fvec_d &x, ivec_d &result, int n, int d);

$ cat max.cu
#include "iterator.h"
#include <thrust/functional.h>
#include <thrust/transform.h>
#include <cublas_v2.h>
#include <iostream>

struct whichMax : thrust::unary_function<double, int>{
  int dim;

  __host__ __device__ whichMax(int dim): dim(dim){}

  __host__ __device__ int operator()(double &vec){

    cublasHandle_t handle;
    cublasCreate_v2(&handle);
    int incx=1, n = dim, result =0;
    double *vec_ptr = thrust::raw_pointer_cast(&vec);

    //find the first index of a maximal element
    cublasIdamax(handle, n, vec_ptr, incx, &result);
    cublasDestroy_v2(handle);
    return result;
  }
};

void cublas_max(fvec_d &x, ivec_d &result, int n, int d){

  stride f(d);
  strideIter siter = thrust::transform_iterator<stride, countIter>(thrust::make_counting_iterator<int>(0), f);
  strideAccessor stridex = thrust::permutation_iterator<realIter, strideIter>(x.begin(), siter);

  whichMax g(d);

  //find the index of maximum for each of n subvectors
  thrust::copy(result.begin(), result.end(), std::ostream_iterator<int>(std::cout, " "));
  std::cout << std::endl;
  thrust::transform(stridex, stridex + n, result.begin(),  g);
  thrust::copy(result.begin(), result.end(), std::ostream_iterator<int>(std::cout, " "));
  std::cout << std::endl;
}

$ cat somePackage.cu
#include "iterator.h"
#include "max.h"
#include <thrust/host_vector.h>
#include <thrust/device_vector.h>
#include <R.h>
#include <Rinternals.h>
#include <Rmath.h>
#include <iostream>

extern "C" SEXP Rcublas_max(SEXP x, SEXP n, SEXP dim){

  double *xptr = REAL(x);
  int N = INTEGER(n)[0], D = INTEGER(n)[0];

  fvec_d dx(xptr, xptr+N*D);
  ivec_d dresult(N);

  cublas_max(dx, dresult, N, D);

  ivec_h hresult(N);
  thrust::copy(dresult.begin(), dresult.end(), hresult.begin());

  SEXP indices = PROTECT(allocVector(INTSXP, N));

  for(int i=0; i<N; ++i)
    INTEGER(indices)[i] = hresult[i];

  UNPROTECT(1);
  return indices;
}

$ make
nvcc -dc -arch=sm_35 -Xcompiler -fPIC -lcublas_device -lcublas_device max.cu -o max.o
nvcc -lib -arch=sm_35 -Xcompiler -fPIC -lcublas_device -lcublas_device max.o -o libmax.a
nvcc -c -arch=sm_35 -Xcompiler -fPIC -lcublas_device somePackage.cu -lmax -I/home/emittman/src/R-3.3.1/builddir/include -I. -o somePackage.o
nvcc -shared -arch=sm_35 -Xcompiler -fPIC -lcublas_device somePackage.o -I/home/emittman/src/R-3.3.1/builddir/include -I. -L. -lcublas_device -lmax -o somePackage.so
ptxas info    : 'device-function-maxrregcount' is a BETA feature

我用 Rcpp 创建了一个 R 包,从 C++ 共享库调用一些外部函数,然后调用 CUDA 内核来执行所需的计算。

你在这里试图做的是将你的 CUDA 代码编译成静态库,然后 link 它到你的 R 包(它本身将被编译成一个共享库)。我的方法和你的不同,我提供我的方法的描述只是为了给你一个不同的想法。

这是一个简化的例子。

kernels.cu 包含 CUDA 代码的共享库:

__global__
void my_cuda_kernel( ... ) {
    // ......
}

main.cu 包含 CUDA 代码的共享库:

extern "C" {
    void do_cuda_work( ... ) {
        thrust::copy( ... );
        my_cuda_kernel <<< ... >>> ( ... );
    }
}
R包中的

package.cpp:

extern void do_cuda_work( ... );

// [[Rcpp::export]]
void call_cuda_code( ... ) {
    do_cuda_work( ... );
}

要将CUDA代码编译成共享库,需要使用:

nvcc -arch=sm_35 -dc ... kernels.cu -o kernels.o
nvcc -arch=sm_35 -dc ... main.cu -o main.o
nvcc --shared -arch=sm_35 ... kernels.o main.o ... libMyCUDALibrary.so

请注意,要使单独的编译工作,您需要为编译器指定 -arch=sm_35,为编译器指定 linker 和 -dc。一旦您成功创建了共享库,link将您的 R 包连接到它就相当简单了。 (但是,您可能需要在 R 包的 src 文件夹下创建一个 Makevars 文件以指定包含和库路径,也许还有 RPATH):

CXX_STD= CXX11
PKG_CPPFLAGS= -I../../../CPP/include
PKG_LIBS= -L../../../CPP/bin/Release -lMyCUDALibrary -Wl,-rpath=$$HOME/MyCUDALibrary/CPP/bin/Release `$(R_HOME)/bin/Rscript -e "Rcpp::LdFlags()"`