为啥在这个例子中 PyCUDA 比 C CUDA 快

Posted 2023-03-23

【中文标题】为啥在这个例子中 PyCUDA 比 C CUDA 快

【英文标题】：Why PyCUDA is faster than C CUDA in this example为什么在这个例子中 PyCUDA 比 C CUDA 快

【发布时间】：2021-05-31 07:20:17

【问题描述】：

我正在探索从 OpenCL 迁移到 CUDA，并进行了一些测试来衡量 CUDA 在各种实现中的速度。令我惊讶的是，在下面的示例中，PyCUDA 实现比 C CUDA 示例快约 20%。

我阅读了很多关于 C CUDA 代码“发布构建”的帖子。我确实尝试在 makefile 中添加-Xptxas -O3，但这并没有什么不同。我还尝试调整执行内核的块大小。不幸的是，它也无助于提高速度。

我的问题是：

可能是什么原因导致 C CUDA 和 PYCUDA 之间的速度差异？ 如果 PYCUDA 中的“高级”（没有更好的词）编译是原因之一，我该如何优化我的 C CUDA 代码的编译？ 在这种情况下，还有其他方法可以提高 C CUDA 的速度吗？

虽然我很欣赏一般的 cmets，但我正在寻找可以在我的机器上验证的可行建议。谢谢！

import pycuda.autoinit
import pycuda.driver as drv
import numpy as np

from pycuda.compiler import SourceModule
import time


mod = SourceModule(
    """
__global__ void saxpy(int n, const float a, float *x, float *y)

    int i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n)
        y[i] = a * x[i] + y[i];
    

"""
)

saxpy = mod.get_function("saxpy")

N = 1 << 25
time_elapse = 0.0

for i in range(100):
    # print(i)
    # print(N)

    x = np.ones(N).astype(np.float32)
    y = 2 * np.ones(N).astype(np.float32)
    start = time.time()
    saxpy(
        np.int32(N),
        np.float32(2.0),
        drv.In(x),
        drv.InOut(y),
        block=(512, 1, 1),
        grid=(int(N / 512) + 1, 1),
    )
    time_elapse += (time.time() - start)


print(time_elapse )
print(y[-100:-1])
print(y.sum())
print(N * 4.0)

#include <stdio.h>
#include <time.h>
#define DIM 512



__global__ void saxpy(int n, float a, float *x, float *y)

    int i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n)
        y[i] = a * x[i] + y[i];


int main(int num_iterations)

    double start;
    double cputime;
    int N = 1 << 25;
    float *x, *y, *d_x, *d_y;
    int i, j;
    for (j = 0; j < num_iterations; j++)
    
        x = (float *)malloc(N * sizeof(float));
        y = (float *)malloc(N * sizeof(float));

        cudaMalloc(&d_x, N * sizeof(float));
        cudaMalloc(&d_y, N * sizeof(float));

        for (i = 0; i < N; i++)
        
            x[i] = 1.0f;
            y[i] = 2.0f;
        

        cudaMemcpy(d_x, x, N * sizeof(float), cudaMemcpyHostToDevice);
        cudaMemcpy(d_y, y, N * sizeof(float), cudaMemcpyHostToDevice);

        // Perform SAXPY on 1M elements
        start = clock();
        saxpy<<<(N + DIM) / DIM, DIM>>>(N, 2.0f, d_x, d_y);
        cputime += ((double)(clock() - start) / CLOCKS_PER_SEC);
        cudaMemcpy(y, d_y, N * sizeof(float), cudaMemcpyDeviceToHost);

        // float maxError = 0.0f;
        // for (int i = 0; i < N; i++)
        //     maxError = max(maxError, abs(y[i] - 4.0f));
        //     //printf("y[%d]: %f\n", i,y[i]);
        // 
        // printf("Max error: %f\n", maxError);

        cudaFree(d_x);
        cudaFree(d_y);
        free(x);
        free(y);
    

 
    printf("cpu time is %f\n", cputime);
    return 0;

我将上述文件保存为cuda_example.cu，并在makefile 中使用以下命令对其进行编译：

nvcc -arch=sm_61 -Xptxas -O3,-v -o main cuda_example.cu

【问题讨论】：

-Xptxas -O3 是nvcc 使用的默认值。您的计时框架对主机活动和设备活动的混合进行计时，这通常是一个坏主意。您可能希望单独比较各个内核的执行时间，然后比较各个主机代码的执行时间。我怀疑你会发现时间差异在后者。

在 pycuda 代码中，您永远不会在计算结果后将数据复制回来，您只需再次循环内核。在 C 版本中，分配、初始化、复制、计算、复制和释放都在同一个循环体中。这简直是​​低效的，与内核或编译器所做的优化无关。

另外，您是否尝试将 number_iterations 作为命令行输入参数？ int main(int num_iterations) 不是这样做的方法。

@geebert 哦，对不起，我实际上将函数重命名为 main 却没有注意到这一点。我会改变它。当我尝试计时saxpy<<<(N + DIM) / DIM, DIM>>>(N, 2.0f, d_x, d_y) 时，它给了我一个不切实际的小数字。我猜它只是排队的时间（不是真正的执行时间）？

cuda 中的内核启动是异步的，因此如果你只测量一个内核执行，你应该有一个运行时 api 调用，如 cudamemcpy，或者在内核调用之后有 cudadeviceSynchronize()，否则你是对的，它可能不是真正的执行时间。

【参考方案1】：

如果我按原样执行您的 CUDA-C 代码，并将 num_iterations 设置为 300，如下所示：

int num_iterations =300;

那么在 Geforce GTX 1650 上执行您的程序大约需要 60 秒。您的代码效率极低，因为您在每次迭代时都在 GPU 和设备之间来回复制数据。 因此，让我们将循环限制为内核执行：

#include <stdio.h>
#include <time.h>
#define DIM 512

__global__ void saxpy(int n, float a, float *x, float *y)

int i = blockIdx.x * blockDim.x + threadIdx.x;
if (i < n)
    y[i] = a * x[i] + y[i];


int main()

double start = clock();
int N = 1 << 25;
float *x, *y, *d_x, *d_y;
int i, j;

int num_iterations = 300;
x = (float *)malloc(N * sizeof(float));
y = (float *)malloc(N * sizeof(float));

cudaMalloc(&d_x, N * sizeof(float));
cudaMalloc(&d_y, N * sizeof(float));

for (i = 0; i < N; i++)

   x[i] = 1.0f;
   y[i] = 2.0f;

cudaMemcpy(d_x, x, N * sizeof(float), cudaMemcpyHostToDevice);
cudaMemcpy(d_y, y, N * sizeof(float), cudaMemcpyHostToDevice);

for (j = 0; j < num_iterations; j++)
    saxpy<<<(N + DIM) / DIM, DIM>>>(N, 2.0f, d_x, d_y);
    cudaDeviceSynchronize();

cudaMemcpy(y, d_y, N * sizeof(float), cudaMemcpyDeviceToHost);
cudaFree(d_x);
cudaFree(d_y);
free(x);
free(y);

double cputime = ((double)(clock() - start) / CLOCKS_PER_SEC);
printf("cpu time is %f\n", cputime);
return 0;

如果我这样做，那么执行时间将变为 1.36 秒。做类似于 PyCUDA 代码的事情，我得到了大约 19 秒的执行时间。

【讨论】：

运行 100 次迭代，C CUDA 执行内核耗时 0.001648 秒，PCUDA 耗时 4.6 秒。我的 GPU 是 NVIDIA P2000。

使用上述代码进行 100 次迭代，我得到大约 0.6 秒，这对我来说似乎是合理的。