计算科学 - 多次运行的 OpenCL C 矩阵乘法 - 吾爱随笔录

我正在尝试通过 OpenCL C 将 dgemm / MPI 矩阵乘法器转换到 GPU 上。我的问题是，下面的代码给出了 9x6 矩阵乘以 6x450 矩阵的正确输出，直到 6x450 矩阵略有变化并且内核输出为移动了 315 次。包括我试图模仿以供参考的 dgemm。

DGEMM：

    umAxBtrans(A, 450, 6, B, 9, 6, &fQ1);
    umAxBtrans(double **A, int Arows, int Acols, double *B, int Brows, int Bcols, double **C){
    op(A) = 'T';
    op(B) = 'F';
    dgemm(&opA, &opB, &rowsB, &rowsA, &colsA, &one, *B, &colsB, *A, &colsA, &zero, *C, &rowsB);

dgemm 手册页： http: //www.math.utah.edu/software/lapack/lapack-blas/dgemm.html

现在，如果该 dgemm 实现不是有点令人困惑，那么我不知道是什么。我已经在这个矩阵乘法代码之外处理了 B 矩阵的转置，因为 B 是一个静态矩阵。然而，矩阵 A 在函数调用的每个循环中都会更新，这意味着每次使用下面的函数和内核代码将新的 A 按值传递给 openCL。

开放式CL：

    clAxBtrans(A, 6, 450, Bcl, 9, 6, &C, context, program, queue, device_counter);
    clAxBtrans(double* A, int Arows, int Acols, double* Bcl, int Brows, int Bcols, double** C, cl_context context, cl_program program, cl_command_queue* queue, cl_uint device_counter){

在这个函数中，我将 A 的值与前一个 dgemm 调用中使用的尺寸一起传递（我已经尝试通过引用和值传递 A，没有区别），以及已经转置的静态 Bcl。然后我设置内核并定义：

    globalWorkSize[0]=450;
    globalWorkSize[1]=9;

然后像这样编写内核：

    __kernel void clAxBtrans(__global double* A,
                     __global double* B,
                     int rowsA,
                     int colsA,
                     int rowsB,
                     int colsB,
                     __global double* C)
    {
    int globalx = get_global_id(0);
    int globaly = get_global_id(1);
    double tmp = 0;
    double tmp0 = 0;
    double tmp1 = 0;
    double tmp2 = 0;
    double tmp3 = 0;
    double tmp4 = 0;
    double tmp5 = 0;
    tmp0 = B[globaly * colsB + 0] * A[0 * colsA + globalx];
    tmp1 = B[globaly * colsB + 1] * A[1 * colsA + globalx];
    tmp2 = B[globaly * colsB + 2] * A[2 * colsA + globalx];
    tmp3 = B[globaly * colsB + 3] * A[3 * colsA + globalx];
    tmp4 = B[globaly * colsB + 4] * A[4 * colsA + globalx];
    tmp5 = B[globaly * colsB + 5] * A[5 * colsA + globalx];
    tmp = tmp0 + tmp1 + tmp2 + tmp3 + tmp4 + tmp5;
    barrier(CLK_GLOBAL_MEM_FENCE);
    C[globaly * colsA + globalx] = tmp;
    }

主机代码到处都有阻塞调用，因为这是在 MPI 并行函数中，并且 localWorkSize 设置为 null 以允许 gpu 选择最“最佳”的大小。

输出是一个 6x450 C 矩阵，它与之前的实现相匹配，直到 A 更新。我的意思是，对于前 5 次运行，A 是一个静态矩阵，但之后每次运行都会发生非常小的变化，具体取决于生成的 C 矩阵。这让我相信我在记忆中错误地解释了矩阵 A。我正在使用下面的 clBuffer 调用将其写入内核：

    Acl = clCreateBuffer(context, CL_MEM_READ_ONLY | CL_MEM_COPY_HOST_PTR, sizeof(double)*Arows*Acols, A, &err);

对于参考矩阵 B 是下面的 9x6，A 是一个完整的 6x450，填充了 .111111 双精度数，结果矩阵 C 是一个 9x450，前 5 次运行填充了 .111111 双精度数。

矩阵 B：

    [1, 0 ,0 ,0 ,0 ,0]
    [0, 1 ,0 ,0 ,0 ,0]
    [0, 0 ,0 ,1 ,0 ,0]
    [0, 1 ,0 ,0 ,0 ,0]
    [0, 0 ,1 ,0 ,0 ,0]
    [0, 0 ,0 ,0 ,1 ,0]
    [1, 0 ,0 ,0 ,0 ,0]
    [0, 0 ,1 ,0 ,0 ,0]
    [0, 0 ,0 ,0 ,0 ,1]

因为我正在复制主机 ptr，所以我相信我的问题要么在于我认为 dgemm 如何使用我的实际内核代码（最有可能）解释矩阵上的维度翻转，要么在于我如何将其编写为缓冲区。我已经检查过，opencl 和以前工作的代码上的第一个新 A 是相同的，但是生成的 C 矩阵是不同的（指向错误的内核代码或内存管理）；如果需要，我可以包含生成的 C 比较，但是第一个 C 结果错误地移动了 315 个 indeces，这导致不同的解决方案真正失控。

我非常感谢你能给我的任何帮助，即使它指出我做了一些非常愚蠢的事情！谢谢！