RCPP和CULA：分段錯誤

Question

我提取的相關位從gputoolsř -package到使用RCPP通過動態加載鏈接到culatools一個共享庫我的GPU上運行的QR分解。在我的Mac上，所有內容都可以在終端和R.app中順利運行。該結果與- [R的QR（）功能，但問題是，分割故障發生在離開R.app（使用終端時，不會發生錯誤）：

*** caught segfault *** 
address 0x10911b050, cause 'memory not mapped' 

我想我縮小了問題的指針「a」和「頭」在鏈接到culatools .c文件：

#include<cula.h> 

void gpuQR(const int *m, const int *n, float *a, const int *lda, float *tau) 
{ 
    culaInitialize(); 
    culaSgeqrf(m[0], n[0], a, lda[0], tau); 
    culaShutdown(); 
} 

我編譯使用我的Mac .c文件：

/usr/local/cuda/bin/nvcc -gencode arch=compute_10,code=sm_10 -gencode arch=compute_11,code=sm_11 -gencode arch=compute_12,code=sm_12 -gencode arch=compute_13,code=sm_13 -gencode arch=compute_20,code=sm_20 -c -I. -I/usr/local/cula/include -m64 -Xcompiler -fPIC gpuQR.c -o gpuQR.o 
/usr/local/cuda/bin/nvcc -gencode arch=compute_10,code=sm_10 -gencode arch=compute_11,code=sm_11 -gencode arch=compute_12,code=sm_12 -gencode arch=compute_13,code=sm_13 -gencode arch=compute_20,code=sm_20 -shared -m64 -Xlinker -rpath,/usr/local/cula/lib64 -L/usr/local/cula/lib64 -lcula_core -lcula_lapack -lcublas -o gpuQR.so gpuQR.o 

我寫使用RCPP和動態加載共享庫gpuQR.so .cpp文件：

#include <Rcpp.h> 
#include <dlfcn.h> 

using namespace Rcpp; 
using namespace std; 

typedef void (*func)(int*, int*, float*, int*, float*); 

RcppExport SEXP gpuQR_Rcpp(SEXP x_, SEXP n_rows_, SEXP n_cols_) 
{  
    vector<float> x = as<vector<float> >(x_); 
    int n_rows = as<int>(n_rows_); 
    int n_cols = as<int>(n_cols_); 
    vector<float> scale(n_cols); 

    void* lib_handle = dlopen("path/gpuQR.so", RTLD_LAZY); 
    if (!lib_handle) 
    { 
     Rcout << dlerror() << endl; 
    } else { 
     func gpuQR = (func) dlsym(lib_handle, "gpuQR"); 
     gpuQR(&n_rows, &n_cols, &(x[0]), &n_rows, &(scale[0])); 
    } 

    dlclose(lib_handle); 

    for(int ii = 1; ii < n_rows; ii++) 
    { 
     for(int jj = 0; jj < n_cols; jj++) 
     { 
      if(ii > jj) { y[ii + jj * n_rows] *= scale[jj]; } 
     } 
    } 

    return wrap(x); 
} 

我使用編譯cpp文件中- [R：

library(Rcpp) 
PKG_LIBS <- sprintf('%s $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)', Rcpp:::RcppLdFlags()) 
PKG_CPPFLAGS <- sprintf('%s', Rcpp:::RcppCxxFlags()) 
Sys.setenv(PKG_LIBS = PKG_LIBS , PKG_CPPFLAGS = PKG_CPPFLAGS) 
R <- file.path(R.home(component = 'bin'), 'R') 
file <- 'path/gpuQR_Rcpp.cpp' 
cmd <- sprintf('%s CMD SHLIB %s', R, paste(file, collapse = ' ')) 
system(cmd) 

並舉了一個例子：

dyn.load('path/gpuQR_Rcpp.so') 

set.seed(100) 
A <- matrix(rnorm(9), 3, 3) 
n_row <- nrow(A) 
n_col <- ncol(A) 

res <- .Call('gpuQR_Rcpp', c(A), n_row, n_col) 
matrix(res, n_row, n_col) 

      [,1]  [,2]  [,3] 
[1,] 0.5250958 -0.8666927 0.8594266 
[2,] -0.2504899 -0.3878644 -0.1277837 
[3,] 0.1502908 0.4742033 -0.8804248 

qr(A)$qr 

      [,1]  [,2]  [,3] 
[1,] 0.5250957 -0.8666925 0.8594266 
[2,] -0.2504899 -0.3878643 -0.1277838 
[3,] 0.1502909 0.4742033 -0.8804247 

有沒有人有一個想法如何解決分段故障？

Answer 1

問題是由從.cpp文件刪除

dlclose(lib_handle); 

解決。這產生以下：

#include <Rcpp.h> 
#include <dlfcn.h> 

using namespace Rcpp; 
using namespace std; 

typedef void (*func)(int*, int*, float*, int*, float*); 

RcppExport SEXP gpuQR_Rcpp(SEXP x_, SEXP n_rows_, SEXP n_cols_) 
{  
    vector<float> x = as<vector<float> >(x_); 
    int n_rows = as<int>(n_rows_); 
    int n_cols = as<int>(n_cols_); 
    vector<float> scale(n_cols); 

    void* lib_handle = dlopen("path/gpuQR.so", RTLD_LAZY); 
    if (!lib_handle) 
    { 
     Rcout << dlerror() << endl; 
    } else { 
     func gpuQR = (func) dlsym(lib_handle, "gpuQR"); 
     gpuQR(&n_rows, &n_cols, &(x[0]), &n_rows, &(scale[0])); 
    } 

    for(int ii = 1; ii < n_rows; ii++) 
    { 
     for(int jj = 0; jj < n_cols; jj++) 
     { 
      if(ii > jj) { x[ii + jj * n_rows] *= scale[jj]; } 
     } 
    } 

    return wrap(x); 
} 

的cpp文件可以在使用編譯ř：

library(Rcpp) 
PKG_LIBS <- sprintf('%s $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)', Rcpp:::RcppLdFlags()) 
PKG_CPPFLAGS <- sprintf('%s', Rcpp:::RcppCxxFlags()) 
Sys.setenv(PKG_LIBS = PKG_LIBS , PKG_CPPFLAGS = PKG_CPPFLAGS) 
R <- file.path(R.home(component = 'bin'), 'R') 
file <- 'path/gpuQR_Rcpp.cpp' 
cmd <- sprintf('%s CMD SHLIB %s', R, paste(file, collapse = ' ')) 
system(cmd) 

實際。C文件鏈接到culatools是：

#include<cula.h> 

void gpuQR(const int *m, const int *n, float *a, const int *lda, float *tau) 
{ 
    culaInitialize(); 
    culaSgeqrf(m[0], n[0], a, lda[0], tau); 
    culaShutdown(); 
} 

它可以使用被編譯：

dyn.load('path/gpuQR_Rcpp.so') 

set.seed(100) 

n_row <- 3 
n_col <- 3 
A <- matrix(rnorm(n_row * n_col), n_row, n_col) 

res <- .Call('gpuQR_Rcpp', c(A), n_row, n_col) 
matrix(res, n_row, n_col) 

      [,1]  [,2]  [,3] 
[1,] 0.5250958 -0.8666927 0.8594266 
[2,] -0.2504899 -0.3878644 -0.1277837 
[3,] 0.1502908 0.4742033 -0.8804248 

qr(A)$qr 

      [,1]  [,2]  [,3] 
[1,] 0.5250957 -0.8666925 0.8594266 
[2,] -0.2504899 -0.3878643 -0.1277838 
[3,] 0.1502909 0.4742033 -0.8804247 

這裏：

gcc -c -I/usr/local/cula/include gpuQR.c 
gcc -shared -Wl,-rpath,/usr/local/cula/lib64 -L/usr/local/cula/lib64 -lcula_lapack -o gpuQR.so gpuQR.o 

QR分解然後可以在- [R使用進行是使用具有16個CUDA核心的NVIDIA GeForce 9400M GPU進行基準測試的結果：

n_row <- 1000; n_col <- 1000 
A <- matrix(rnorm(n_row * n_col), n_row, n_col) 
B <- A; dim(B) <- NULL 

res <- benchmark(.Call('gpuQR_Rcpp', B, n_row, n_col), qr(A), columns = c('test', 'replications', 'elapsed', 'relative'), order = 'relative') 

            test replications elapsed relative 
1 .Call("gpuQR_Rcpp", B, n_row, n_col)   100 38.037 1.000 
2        qr(A)   100 152.575 4.011 

Answer 2

實際上並不需要動態加載鏈接到culatools的共享庫。我最初想到了這一點，但我沒有得到使用Rcpp編譯的.cpp文件。總之，新的.cpp文件是：

#include<Rcpp.h> 
#include<cula.h> 

using namespace Rcpp; 
using namespace std; 

RcppExport SEXP gpuQR_Rcpp(SEXP x_, SEXP n_rows_, SEXP n_cols_) 
{  
    vector<float> x = as<vector<float> >(x_); 
    int n_rows = as<int>(n_rows_); 
    int n_cols = as<int>(n_cols_); 
    vector<float> scale(n_cols); 

    culaInitialize(); 
    culaSgeqrf(n_rows, n_cols, &(x[0]), n_rows, &(scale[0])); 
    culaShutdown(); 

    for(int ii = 1; ii < n_rows; ii++) 
    { 
     for(int jj = 0; jj < n_cols; jj++) 
     { 
      if(ii > jj) { x[ii + jj * n_rows] *= scale[jj]; } 
     } 
    } 

    return wrap(x); 
} 

.cpp文件使用編譯：

library(Rcpp) 
PKG_LIBS <- sprintf('-Wl,-rpath,/usr/local/cula/lib64 -L/usr/local/cula/lib64 -lcula_lapack %s $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)', Rcpp:::RcppLdFlags()) 
PKG_CPPFLAGS <- sprintf('-I/usr/local/cula/include %s', Rcpp:::RcppCxxFlags()) 
Sys.setenv(PKG_LIBS = PKG_LIBS , PKG_CPPFLAGS = PKG_CPPFLAGS) 
R <- file.path(R.home(component = 'bin'), 'R') 
file <- 'path/gpuQR_inc.cpp' 
cmd <- sprintf('%s CMD SHLIB %s', R, paste(file, collapse = ' ')) 
system(cmd) 

在哪裏設置相應的路徑culatools。整個過程不會跑的更快，但是不需要再編譯連接到culatools的共享庫並動態加載它。

我認爲，這是一個很好的替代gputoolsř -package延長- [R與C++和在GPU上執行線性代數運算。