降低CUDA内核运行时:内核中矩阵的动态内存分配

Decrease cuda kernel runtime: dynamic memory allocation of matrices in kernel

本文关键字:内核 动态 内存 分配 CUDA 运行时 降低      更新时间:2023-10-16

我想通过在GPU上并行运行矩阵操作来适合大量较小的矩阵。我写了似乎正在运行的代码,但是它比预期的要慢。当前,尽管GPU上有平行的计算,但是在CPU上的单个线程上运行它需要较短的时间。Nvidia Visual Profiler似乎表明内存分配需要大量时间。我怀疑这是罪魁祸首是内核内不同尺寸矩阵的动态内存分配。我需要建议并帮助加快内核运行时。

我尝试使用循环中创建的每个矩阵使用新的和删除。

这是内核:

__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。如何分配不同大小的矩阵的记忆效率低下。你们如何在内核中动态存储各种尺寸的矩阵?这应该如何以最有效的方式完成?

对给定代码上的任何其他反馈表示感谢。

在我看来,这是适用的三个非常简单的经验规则:

  1. 动态内存分配是始终昂贵的,无论您编程如何。
  2. pertarant代码从不使用动态内存分配,除非绝对必要。
  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;
}

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

通过执行此操作,您可以在许多计算上摊销一次内存分配的成本,这应该比您当前的方法更有效。