Eigen JacobiSVD cuda編譯錯誤

Question

強調文本我有一個錯誤，關於在我的cuda函數中調用JacobiSVD。

這是導致錯誤的代碼的一部分。

Eigen::JacobiSVD<Eigen::Matrix3d> svd(cov_e, Eigen::ComputeThinU | Eigen::ComputeThinV); 

並且這是錯誤消息。

CUDA_voxel_building.cu(43): error: calling a __host__ function("Eigen::JacobiSVD , (int)2> ::JacobiSVD") from a __global__ function("kernel") is not allowed

我用下面的命令來編譯它。

nvcc -std=c++11 -D_MWAITXINTRIN_H_INCLUDED -D__STRICT_ANSI__ -ptx CUDA_voxel_building.cu 

我在ubuntu 16.04上使用代碼8.0和eigen3。 似乎像其他功能，如特徵值分解也給出了相同的錯誤。

任何人都知道解決方案嗎？我在下面附上我的代碼。

//nvcc -ptx CUDA_voxel_building.cu 
#include </usr/include/eigen3/Eigen/Core> 
#include </usr/include/eigen3/Eigen/SVD> 
/* 
#include </usr/include/eigen3/Eigen/Sparse> 

#include </usr/include/eigen3/Eigen/Dense> 
#include </usr/include/eigen3/Eigen/Eigenvalues> 

*/ 





__global__ void kernel(double *p, double *breaks,double *ind, double *mu, double *cov, double *e,double *v, int *n, char *isgood, int minpts, int maxgpu){ 
    bool debuginfo = false; 
    int idx = threadIdx.x + blockIdx.x * blockDim.x; 
    if(debuginfo)printf("Thread %d got pointer\n",idx); 
    if(idx < maxgpu){ 


     int s_ind = breaks[idx]; 
     int e_ind = breaks[idx+1]; 
     int diff = e_ind-s_ind; 

     if(diff >minpts){ 
      int cnt = 0; 
      Eigen::MatrixXd local_p(3,diff) ; 
      for(int k = s_ind;k<e_ind;k++){ 
       int temp_ind=ind[k]; 

       //Eigen::Matrix<double, 3, diff> local_p; 
       local_p(1,cnt) = p[temp_ind*3]; 
       local_p(2,cnt) = p[temp_ind*3+1]; 
       local_p(3,cnt) = p[temp_ind*3+2]; 
       cnt++; 
      } 

      Eigen::Matrix3d centered = local_p.rowwise() - local_p.colwise().mean(); 
      Eigen::Matrix3d cov_e = (centered.adjoint() * centered)/double(local_p.rows() - 1); 

      Eigen::JacobiSVD<Eigen::Matrix3d> svd(cov_e, Eigen::ComputeThinU | Eigen::ComputeThinV); 
    /*   Eigen::Matrix3d Cp = svd.matrixU() * svd.singularValues().asDiagonal() * svd.matrixV().transpose(); 


      mu[idx]=p[ind[s_ind]*3]; 
      mu[idx+1]=p[ind[s_ind+1]*3]; 
      mu[idx+2]=p[ind[s_ind+2]*3]; 

      e[idx]=svd.singularValues()(0); 
      e[idx+1]=svd.singularValues()(1); 
      e[idx+2]=svd.singularValues()(2); 

      n[idx] = diff; 
      isgood[idx] = 1; 

      for(int x = 0; x < 3; x++) 
      { 
       for(int y = 0; y < 3; y++) 
       { 
        v[x+ 3*y +idx*9]=svd.matrixV()(x, y); 
        cov[x+ 3*y +idx*9]=cov_e(x, y); 
        //if(debuginfo)printf("%f ",R[x+ 3*y +i*9]); 
        if(debuginfo)printf("%f ",Rm(x, y)); 
       } 
      } 
*/ 

     } else { 
      mu[idx]=0; 
      mu[idx+1]=0; 
      mu[idx+2]=0; 

      e[idx]=0; 
      e[idx+1]=0; 
      e[idx+2]=0; 

      n[idx] = 0; 
      isgood[idx] = 0; 

      for(int x = 0; x < 3; x++) 
      { 
       for(int y = 0; y < 3; y++) 
       { 
        v[x+ 3*y +idx*9]=0; 
        cov[x+ 3*y +idx*9]=0; 
       } 
      } 
     } 




    } 
} 

Answer 1

首先，Ubuntu 16.04提供了Eigen 3.3-beta1，這並不真正推薦使用。我會建議升級到更新的版本。此外，包括艾根，寫（例如）：

#include <Eigen/Eigenvalues> 

與-I /usr/include/eigen3編譯（如果使用操作系統提供的版本），或更好的-I /path/to/local/eigen-version。

然後，正如talonmies所說，你不能從內核調用主機函數，（我不確定此刻，爲什麼JacobiSVD沒有標記爲設備函數），但在你的情況下，它會做得更多無論如何，感覺要使用Eigen::SelfAdjointEigenSolver。既然你是分解矩陣是固定大小的3×3，你應該實際使用的優化computeDirect方法：

Eigen::SelfAdjointEigenSolver<Eigen::Matrix3d> eig; // default constructor 
eig.computeDirect(cov_e); // works for 2x2 and 3x3 matrices, does not require loops 

看來computeDirect甚至適用於被Ubuntu提供的測試版（我還是建議更新）。

一些不相關的注意事項：

1. 下列哪項是錯誤的，因爲你應該與指數0開始：

   local_p(1,cnt) = p[temp_ind*3]; 
   local_p(2,cnt) = p[temp_ind*3+1]; 
   local_p(3,cnt) = p[temp_ind*3+2]; 
   

   此外，您還可以在一行寫：

   local_p.col(cnt) = Eigen::Vector3d::Map(p+temp_ind*3); 
   

2. 此行不適合（除非diff==3）：

   Eigen::Matrix3d centered = local_p.rowwise() - local_p.colwise().mean(); 
   

   什麼你大概的意思是（local_p實際上是3XN不NX3）

   Eigen::Matrix<double, 3, Eigen::Dynamic> centered = local_p.colwise() - local_p.rowwise().mean(); 
   

   和計算cov_e當你需要.adjoint()第二個因素，不是第一次。

3. 你可以同時避免 '大' 矩陣local_p和centered，直接累加Eigen::Matrix3d sum2和Eigen::Vector3d sum與sum2 += v*v.adjoint()和sum +=v和計算

   Eigen::Vector3d mu = sum/diff; 
   Eigen::Matrix3d cov_e = (sum2 - mu*mu.adjoint()*diff)/(diff-1);