cuda2 向量加法

向量加法

向量加法程序解读

#include<stdio.h>
#include<cuda.h>

typedef float FLOAT;
#define USE_UNIX 1  区别不同系统 

get thread id 1D block and 2D grid
#define get_tid() (block)

get block id， 2D grid

warm up 可选的，让gpu先运作起来，
函数前后来一个时间戳，获取程序时间长度
后面一个是纯粹在cpu上运行的代码，用于时间比较。

关于指针

在 CUDA 中，函数的参数需要明确指定其在设备端（称为 kernel 函数）还是在主机端执行。在 __global__ 修饰的 kernel 函数中，所有参数都会被显式地传递到设备端，并在该函数中进行操作。因此，在 kernel 函数参数声明中，需要使用指针类型 FLOAT * 来传递数据。

在您提供的函数定义中，FLOAT *x 表示一个指向存储输入数组 x 的内存空间的指针，FLOAT *y 表示一个指向存储输出数组 y 的内存空间的指针，int N 则表示数组的大小。这样定义的好处是在 kernel 函数中可以直接通过指针来访问数组元素，而不需要进行数据拷贝等额外操作，从而提高计算效率。

然而，由于 kernel 函数只能在设备端执行，因此在主机端调用 kernel 函数时，需要先将数据从主机端拷贝到设备端，再在设备端执行 kernel 函数。通常情况下，我们会使用类似于下面的代码来调用 kernel 函数：

FLOAT *d_x, *d_y; // 定义设备端指针
cudaMalloc(&d_x, N * sizeof(FLOAT)); // 在设备端分配空间
cudaMemcpy(d_x, x, N * sizeof(FLOAT), cudaMemcpyHostToDevice); // 将 x 数组拷贝到设备端
cudaMalloc(&d_y, N * sizeof(FLOAT)); // 在设备端分配空间
add<<<1, 1>>>(d_x, d_y, N); // 调用 kernel 函数
cudaMemcpy(y, d_y, N * sizeof(FLOAT), cudaMemcpyDeviceToHost); // 将 y 数组拷贝回主机端
cudaFree(d_x); // 释放设备端空间
cudaFree(d_y);

在这个例子中，x 和 y 是主机端指针，d_x 和 d_y 则是设备端指针。我们先在设备端使用 cudaMalloc() 分配空间，再使用 cudaMemcpy() 将数据从主机端拷贝到设备端，然后调用 kernel 函数，最后再将结果从设备端拷贝回主机端。通过这种方式，就可以在 CUDA 中实现高效的数据并行计算了。

关于N

在 CUDA 中，如果函数参数需要传递的是一个单一的值，那么就不需要使用指针类型。在您提供的 add() kernel 函数中，N 是表示数组大小的一个整数参数，不需要使用指针类型。

当您调用这个 kernel 函数时，可以直接使用下面的语句来传递参数：

add<<<num_blocks, num_threads>>>(x, y, N);

其中，N 是一个整数变量，不需要使用取地址符号 &。在设备端执行 kernel 函数时，N 可以直接被访问并使用，因为它是一个值参数。

需要注意的是，虽然 N 不需要使用指针类型，但它仍然需要在 kernel 函数参数列表中进行声明，因为 CUDA 在编译时需要确定每个 kernel 函数的参数列表。

代码

获取线程全局编号
如果编号小于向量长度，那么就进行相加，其实就是并行计算

后面是cpu的代码用于比较
这里是循环的串行计算

代码段2

这里就是区分unix 还是 windows下的时间

代码段3

热身代码调用了8次
开辟了一个块，含有256个线程

main

int main()
{
    int N = 20000000; 
    int nbytes = N * sizeof(FLOAT);  定义两千万的浮点数 的内存

    /* 1D block */
    int bs = 256;		块的大小 一维度的

    /* 2D grid */
    int s = ceil(sqrt((N + bs - 1.) / bs)); 不同的GPU在x y z 维度上有限制，不论是块还是网格 65536 z维度1024
    dim3 grid = dim3(s, s);			取平方根，表示二维的

    FLOAT *dx = NULL, *hx = NULL; d表示device h表示host 分别为两者申请内存
    FLOAT *dy = NULL, *hy = NULL;  三个向量 指针
    FLOAT *dz = NULL, *hz = NULL;

    int itr = 30;
    int i;
    double th, td;

    /* allocate GPU mem */ 申请GPU内存 三个向量  你这里是加法运算 用到3个向量
    cudaMalloc((void **)&dx, nbytes);   取地址 相当于给指针分配内存，并用二重指针接收
    cudaMalloc((void **)&dy, nbytes);
    cudaMalloc((void **)&dz, nbytes);

    if (dx == NULL || dy == NULL || dz == NULL) {     指针是否为空，是否正确分配到内存
        printf("couldn't allocate GPU memory\n");
        return -1;
    }

    printf("allocated %.2f MB on GPU\n", nbytes / (1024.f * 1024.f));

    /* alllocate CPU mem */
    hx = (FLOAT *) malloc(nbytes);
    hy = (FLOAT *) malloc(nbytes);
    hz = (FLOAT *) malloc(nbytes);

    if (hx == NULL || hy == NULL || hz == NULL) {
        printf("couldn't allocate CPU memory\n");
        return -2;
    }
    printf("allocated %.2f MB on CPU\n", nbytes / (1024.f * 1024.f));

    /* init */
    for (i = 0; i < N; i++) {   全部cpu上初始化
        hx[i] = 1;
        hy[i] = 1;
        hz[i] = 1;
    }

    /* copy data to GPU */
    cudaMemcpy(dx, hx, nbytes, cudaMemcpyHostToDevice);  dx目的地 hx从哪里拷贝， 拷贝多少字节 ，拷贝方向
    cudaMemcpy(dy, hy, nbytes, cudaMemcpyHostToDevice);
    cudaMemcpy(dz, hz, nbytes, cudaMemcpyHostToDevice);

    /* call GPU */
    cudaDeviceSynchronize(); 卡住，必须等Gpu所有跑完才跑 ，因为是异步的，后面的程序调用完，不等执行完，就继续往下执行了
    td = get_time();
    
    for (i = 0; i < itr; i++) vec_add<<<grid, bs>>>(dx, dy, dz, N);

    cudaDeviceSynchronize(); 等GPU跑完才开始往下执行，相当于计算gpu运行时间
    td = get_time() - td;

    /* CPU */
    th = get_time();  cpu是同步的，直接计算
    for (i = 0; i < itr; i++) vec_add_host(hx, hy, hz, N);
    th = get_time() - th;

    printf("GPU time: %e, CPU time: %e, speedup: %g\n", td, th, th / td);

    cudaFree(dx);
    cudaFree(dy);
    cudaFree(dz);

    free(hx);
    free(hy);
    free(hz);

    return 0;
}

这行代码使用了 C++ 标准库中的 ceil() 函数和数学库中的 sqrt() 函数，用于计算一个整数 N 被分为块大小为 bs 的若干块时，需要划分的块数。

具体地，表达式 (N + bs - 1.) / bs 表示将 N 加上一个块大小 bs 减一后再除以 bs，得到的结果是一个小数，但因为加上了一个浮点数 bs - 1.，所以结果会自动向上取整。例如，当 N=1000、bs=256 时，表达式的结果为 4.894，但向上取整后就变成了 5。

接着，函数 sqrt() 被用于对结果取平方根，得到的值表示将 N 划分为平方块数时每边的块数。最后，又调用了标准库中的 ceil() 函数，将该数向上取整并赋值给变量 s，即为划分后的块数。

换句话说，这行代码的作用是计算需要将一个长度为 N 的数据分割成块大小为 bs 的若干块时，一共需要划分多少块。这个划分结果会被保存在变量 s 中，供后续使用。

为什么要 +bs-1再除，直接除不好

由于加上 bs-1 的操作是为了确保整数除法可以向上取整，避免数据划分误差，而 (N + bs - 1) / bs 可以达到这个效果，因此使用 (N + bs - 1) / bs 进行整数除法可以得到正确的结果。当 N=1000、bs=256 时，使用 (N + bs - 1) / bs 进行整数除法得到的商是 (1000 + 256 - 1) / 256 = 4，即实际需要划分成 4 个块。再次感谢您的指正。

在进行整数除法时，为了避免舍入误差可能引起的结果错误，我们经常需要将被除数加上除数减一得到一个偏移值，使得偏移后的被除数能够确保被恰当地划分成为等大小的整数块。

在这个具体的例子中，我们需要将长度为 N 的数据分割成块大小为 bs 的若干块，并计算需要划分多少块。假设 N 不能被 bs 整除，则使用 (N + bs - 1) / bs 进行整数除法后的商会比真实的商小 1，即有可能少算一个块。例如，当 N=1000、bs=256 时，使用 (N ) / bs 进行整数除法得到的商是 3，而不是真实的商 4。

因此，我们需要对被除数进行调整，使得整数除法得到的商向上取整，即确保所有的数据都被正确地划分到某个块中。而加上 bs-1 的目的正是为了解决这个问题，确保整数除法得到的商是向上取整的结果，即每个块装满数据，并且不会有数据被遗漏。

所以，在当前的情况下，将 N 加上 bs-1 再除以 bs，可以确保划分的块数为向上取整后的结果，避免了直接进行整数除法所可能带来的数据划分误差。

ceil

是的，ceil 这个函数是向上取整的。它的全称是 “ceiling”，表示将一个实数向上取整到最近的整数。具体来说，如果给定一个实数 x，那么 ceil(x) 的值就是大于或等于 x 的最小整数。

例如，ceil(2.3) 的结果是 3，ceil(-1.5) 的结果是 -1，ceil(4) 的结果是 4。在计算机科学中，ceil 函数常用于需要对数据进行分块或分页的操作中，确保每一块或每一页都有足够的容量，避免数据划分不均匀或发生越界问题。

3

这段代码使用了 CUDA 中的动态内存分配函数 cudaMalloc，其作用是在设备端（GPU）上分配一块指定大小的内存，并返回指向该内存的指针。具体来说，这段代码中的 dx 是一个指向 FLOAT 类型数据的指针，nbytes 是需要分配的内存大小，单位是字节。

在这段代码中，cudaMalloc((void **)&dx, nbytes) 表示将大小为 nbytes 的内存块分配给 dx 指针，并且内存的起始地址被存储在 dx 指向的内存单元中。因为 cudaMalloc 函数返回值类型是 void，所以需要将其返回的指针强制转换为 FLOAT * 类型，即 cudaMalloc((void **)&dx, nbytes) 中的 (FLOAT **) 部分，表示将 &dx 的类型从 FLOAT ** 转换为 void **，以便与 cudaMalloc 函数的参数类型匹配。

值得注意的是，在使用完动态分配的内存后，程序应该调用 cudaFree 函数释放内存，以避免内存泄露。例如，可以使用 cudaFree(dx) 来释放 dx 所指向的内存块。

为什么要 x y z

在 CUDA 编程中，通常需要在主机端（CPU）和设备端（GPU）之间传输数据。为了方便这些数据的管理和操作，CUDA 引入了统一内存（Unified Memory）的概念，即将主机端和设备端的内存作为同一块内存进行管理。在使用统一内存时，我们可以使用 cudaMallocManaged 函数来分配一块统一内存，该函数会在主机端和设备端上都分配一部分内存，然后将这两部分内存映射到同一块物理内存上。

因此，在使用统一内存时，我们只需要分配一个指向 FLOAT 类型数据的指针，然后将其在主机端和设备端之间传递即可。具体来说，我们可以通过下面的代码实现：

FLOAT *x = NULL;
cudaMallocManaged(&x, nbytes);

其中，nbytes 表示需要分配的内存大小，单位是字节。这段代码会在主机端和设备端同时分配一块大小为 nbytes 的统一内存，然后将 x 指向这块内存的起始地址。在程序运行过程中，我们可以通过在主机端对 x 进行操作来修改设备端的内存数据，或者通过在设备端对 x 进行操作来修改主机端的内存数据。

至于为什么需要申请三个指针 dx、dy、dz，我猜测可能是在某些情况下需要同时传递多个数组到 GPU 中进行并行计算，这时可以分别使用三个指针分别指向这些数组的起始地址。具体而言，例如在进行向量加法的操作时，我们可以将两个向量分别存储在 dx 和 dy 中，然后将结果保存在 dz 中。在 GPU 中进行并行计算时，可以通过下面的代码实现：

__global__ void add_vectors(FLOAT *dx, FLOAT *dy, FLOAT *dz, int n) {
    int i = threadIdx.x + blockIdx.x * blockDim.x;
    if (i < n) {
        dz[i] = dx[i] + dy[i];
    }
}

其中，add_vectors 是一个 GPU 内核函数，用于实现向量加法的并行计算。dx、dy 和 dz 分别表示两个输入向量和一个输出向量的起始地址，n 表示向量的长度。每个线程计算一部分数据，通过线程索引 threadIdx.x 和块索引 blockIdx.x 来计算需要处理的数据的起始位置 i，然后通过 dz[i] = dx[i] + dy[i] 计算结果并保存到输出向量中。

malloc

malloc() 是一个在 C 和 C++ 中用于动态分配内存的函数，它的原型如下：

void* malloc(size_t size);

它接受一个参数 size，表示需要分配内存块的大小（以字节为单位），返回值是一个指向分配内存块的 void 指针。如果分配失败，则返回空指针（NULL）。

例如，下面的代码使用了 malloc() 函数动态分配了一个长度为 10 的整数数组，并将其初始化为 0：

#include <stdio.h>
#include <stdlib.h>

int main() {
    int* a = (int*) malloc(10 * sizeof(int));
    if (a == NULL) {
        printf("Failed to allocate memory.\n");
        exit(EXIT_FAILURE);
    }

    for (int i = 0; i < 10; i++) {
        a[i] = 0;
    }

    // 访问并处理 a 数组

    free(a); // 释放内存块
    return 0;
}

需要注意的是，使用 malloc() 函数分配的内存块需要手动调用 free() 函数进行释放，否则会导致内存泄漏。

cudaMalloc((void **)&dx, nbytes);

在 CUDA 中，cudaMalloc() 函数用于在设备端分配内存。其中，void ** 是一个指向指针的指针，它的作用是使函数能够修改传入的指针的值。该函数原型如下：

cudaError_t cudaMalloc(void** devPtr, size_t size);

第一个参数 devPtr 是一个指向指针的指针，用于存储函数分配的设备端内存的地址。需要注意的是指针本身也需要进行分配，因此需要使用 void ** 来接收该指针的地址。

因此，cudaMalloc((void **)&dx, nbytes) 中的 &dx 表示 dx 变量的地址，&dx 实际上是一个指向 dx 的指针，因此需要将其转换为 void ** 类型，并且在前面加上 & 取其地址，以便能够存储函数分配的设备端内存的地址。

例如，下面的代码使用了 cudaMalloc() 函数分配了一个大小为 nbytes 字节的设备内存空间，并将其地址存储到 dx 指针中：

#include <stdio.h>
#include <cuda_runtime.h>

int main() {
    int *dx;
    size_t nbytes = 1024;

    cudaError_t err = cudaMalloc((void **)&dx, nbytes);
    if (err != cudaSuccess) {
        printf("Failed to allocate device memory: %s\n", cudaGetErrorString(err));
        exit(EXIT_FAILURE);
    }

    // 在设备中使用 dx 指针分配的内存块

    cudaFree(dx); // 释放设备内存
    return 0;
}

在使用完 cudaMalloc() 分配的内存空间后，需要使用 cudaFree() 函数将其释放，否则会导致设备端内存泄漏。

cudaDeviceSynchronize();

在CUDA程序中，cudaDeviceSynchronize()函数用于线程同步，它可以使当前主机线程等待所有设备上的CUDA核函数执行完成，并且在执行完成后返回是否有错误发生。该函数没有任何参数。

当CUDA主机代码调用一个CUDA核函数时，该调用是异步的，即主机代码不会等待CUDA核函数执行完成而是会立即返回并继续执行下一条语句。因此，如果主机代码需要使用CUDA核函数的输出结果，那么就需要使用 cudaDeviceSynchronize() 函数来等待CUDA核函数执行完成。

例如，下面的代码展示了 cudaDeviceSynchronize() 函数的使用场景：

#include <stdio.h>
#include <cuda_runtime.h>

__global__ void kernel(int *a, int *b, int *c) {
    int i = threadIdx.x;
    c[i] = a[i] + b[i];
}

int main() {
    int n = 256;
    int *a, *b, *c;
    int *d_a, *d_b, *d_c;

    // 分配内存
    a = (int*) malloc(n * sizeof(int));
    b = (int*) malloc(n * sizeof(int));
    c = (int*) malloc(n * sizeof(int));
    cudaMalloc((void**)&d_a, n * sizeof(int));
    cudaMalloc((void**)&d_b, n * sizeof(int));
    cudaMalloc((void**)&d_c, n * sizeof(int));

    // 初始化输入数据
    for (int i = 0; i < n; i++) {
        a[i] = i;
        b[i] = 2 * i;
    }

    // 将输入数据从主机内存拷贝到设备内存中
    cudaMemcpy(d_a, a, n * sizeof(int), cudaMemcpyHostToDevice);
    cudaMemcpy(d_b, b, n * sizeof(int), cudaMemcpyHostToDevice);

    // 调用 CUDA 核函数
    kernel<<<1, n>>>(d_a, d_b, d_c);

    // 等待 CUDA 核函数执行完成
    cudaDeviceSynchronize();

    // 将计算结果从设备内存拷贝回主机内存
    cudaMemcpy(c, d_c, n * sizeof(int), cudaMemcpyDeviceToHost);

    // 输出计算结果
    for (int i = 0; i < n; i++) {
        printf("%d + %d = %d\n", a[i], b[i], c[i]);
    }

    // 释放内存
    free(a); free(b); free(c);
    cudaFree(d_a); cudaFree(d_b); cudaFree(d_c);
    return 0;
}

在这个例子中，cudaDeviceSynchronize()函数的作用是等待 kernel() 函数执行完成，并且等待所有设备上的CUDA核函数执行完成。如果不使用该函数调用 kernel() 后立即将结果从设备拷贝到主机可能会导致程序输出以及计算结果出现错误。

总结

10倍以上速度