我想通过在 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。如何分配不同大小的矩阵的内存效率低下。你们如何看待在内核中动态存储各种大小的矩阵?这应该如何以最有效的方式完成?
对给定代码的任何其他反馈表示赞赏。
最佳答案
在我看来,三个非常简单的经验法则适用于此:
- 无论您在什么平台上编程,动态内存分配总是是昂贵的。
- 除非绝对必要,否则高性能代码从不使用动态内存分配。
- 如果动态内存分配是绝对必要,预分配内存并尽可能重复使用
如果您查看您的代码,它违反了所有这三个概念。
在内核启动之前,您清楚地知道(或者可以简单地计算)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;
}
[代码显然是用浏览器编写的,从未编译和测试,使用风险自负]。
通过这样做,您可以在多次计算中分摊一次性内存分配的成本,这应该比您当前的方法更有效。
关于c++ - 减少 cuda 内核运行时 : dynamic memory allocation of matrices in kernel,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/55915999/