时间戳

首页 > 解决方案 > 减少cuda内核运行时间：内核中矩阵的动态内存分配

c++ - 减少cuda内核运行时间：内核中矩阵的动态内存分配

问题描述

我想通过在 GPU 上并行运行矩阵运算来为大量较小的矩阵执行 OLS 拟合。我编写了似乎可以正常运行的代码，但是它比预期的要慢。目前，尽管在 GPU 上进行并行计算，但在 CPU 上的单个线程上运行它需要更短的时间。Nvidia Visual Profiler 似乎表明内存分配占用了大量时间。我怀疑罪魁祸首是内核内部不同大小矩阵的动态内存分配。我需要有关加快内核运行时的建议和帮助。

我尝试对循环中创建的每个矩阵使用 new 和 delete 。

这是内核：

__global__
void comb_ols(double *y, double *X, double *R2 ,const unsigned int M, const unsigned int N, int* sub_col, int *sub_size, int* cumulative_size, const unsigned int numberOfCalculations){

    int size;   
    int start_index;

    int index = blockIdx.x*blockDim.x+threadIdx.x;
    int stride = blockDim.x*gridDim.x;  
    for(int i = index; i < numberOfCalculations; i+=stride){    

        size = sub_size[i];
        start_index = cumulative_size[i];             

        double *sub_matrix = new double[M*(1+size)];


            for(int j = 0; j < size; j++){
            for(int k  = 0; k<M; k++){
                sub_matrix[k] = 1;
                sub_matrix[k + M * (1 +  j)] = X[k + M * (sub_col[start_index+j]+1)];                                           
                                            }       
            }
        }

        R2[i] = getR2(y,sub_matrix,M,size+1);


        delete [] sub_matrix;
    }
}

在设备函数 getR2 中，我们有以下内容：

__device__
double getR2(double *y, double *X ,const unsigned int M, const unsigned int N) {

    // Initilize values
    double R2, numerator;
    double* A = new double[N*N];
    double* IA = new double[N*N];
    double* yX = new double[N];  
    // Generate all components
    XtX(X, A, M, N);
    LUPDecompose(A, N);
    LUPInvert(A, N, IA);
    yTX(y, X, yX, M, N);
    // Calc R2
    numerator = olsR2numerator(yX, IA, N);
    R2 = numerator / yTy(y, M);
    //R2 = yTy(y,M);

    delete[] A;
    delete[] IA;
    delete[] yX;

    return R2;
}

实际的内核调用是这样的：

com_ols<<<numBlocks, blockSize >>>(Y,X,R2,M,N,sub_columns, sub_size, cumulative_size, numberOfCalculations);

目前，内核运行时间大约为 1.4 秒，而在单线程 CPU 上为 0.7 秒。我预计内核运行时间会快得多，因为它只是循环矩阵运算的许多迭代，这应该适用于 gpu。不同大小的矩阵的内存分配方式有些低效。你们对在内核中动态存储各种大小的矩阵有什么看法？这应该如何以最有效的方式完成？

对给定代码的任何其他反馈表示赞赏。

标签： c++cudagpudynamic-memory-allocation

解决方案

在我看来，这里适用三个非常简单的经验法则：

  1. 无论您在什么平台上编程，动态内存分配总是很昂贵。
  2. 除非绝对必要，否则高性能代码从不使用动态内存分配。
  3. 如果动态内存分配是绝对必要的，请预先分配内存并尽可能地重用它

如果您查看您的代码，它违反了所有这三个概念。

sub_size您清楚地知道（或可以简单地计算）内核启动之前的最大值是多少。使用先验知识对您有利 - 为计算预先分配堆内存，该计算足够大以处理数据集中最大的问题并在线程的生命周期内重新使用它。你的内核很容易看起来像这样：

__global__
void comb_ols(double *y, double *X, double *R2 ,const unsigned int M, 
             const unsigned int N, int* sub_col, int *sub_size, int* cumulative_size, 
             const unsigned int numberOfCalculations, const int max_size){

    int size;   
    int start_index;

    int index = blockIdx.x*blockDim.x+threadIdx.x;
    int stride = blockDim.x*gridDim.x;

    double *sub_matrix = new double[M*(1+max_size)];
    R2scratch temp(1+max_size);

    for(int i = index; i < numberOfCalculations; i+=stride){    

        size = sub_size[i];
        start_index = cumulative_size[i];             
        for(int j = 0; j < size; j++){
            for(int k  = 0; k<M; k++){
                sub_matrix[k] = 1;
                sub_matrix[k + M * (1 +  j)] = X[k + M * (sub_col[start_index+j]+1)];                                           
                                            }       
            }
        }
        R2[i] = getR2(y,sub_matrix,M,size+1,temp);
    }
    delete [] sub_matrix;
}

设备功能如下：

struct R2scratch
{
    double* A;
    double* IA;
    double* yX;  

    __device__
    R2scratch(int N) {
        A = new double[N*N];
        IA = new double[N*N];
        yX = new double[N];  
    };

    __device__
    ~R2scratch() {
        delete[] A;
        delete[] IA;
        delete[] yX;
    };
};

__device__
double getR2(double *y, double *X ,const unsigned int M, const unsigned int N, 
             R2scratch &scratch) {

    // Initilize values
    double R2, numerator;
    double* A = scratch.A;
    double* IA = scratch.IA;
    double* yX = scratch.yX;

    // Generate all components
    XtX(X, A, M, N);
    LUPDecompose(A, N);
    LUPInvert(A, N, IA);
    yTX(y, X, yX, M, N);
    // Calc R2
    numerator = olsR2numerator(yX, IA, N);
    R2 = numerator / yTy(y, M);
    //R2 = yTy(y,M);

    return R2;
}

[代码显然是在浏览器中编写的，从未编译和测试，使用风险自负]。

通过这样做，您可以在多次计算中分摊一次性内存分配的成本，这应该比您当前的方法更有效。

推荐阅读