CUDA笔记3

最新推荐文章于 2024-07-25 14:12:59 发布
最新推荐文章于 2024-07-25 14:12:59 发布
阅读量1k 收藏 16
点赞数 18
分类专栏： CUDA 文章标签： 笔记 算法
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/qq_44576434/article/details/135730581
版权

CUDA编程-入门、加法、乘法

例子1：入门

#include <stdio.h>

__global__ void hello_from_gpu()
{
    const int bid = blockIdx.x;
    const int tid = threadIdx.x;
    printf("Hello world from block %d and thread %d!\n",bid,tid);
}
int main(void)
{
    hello_from_gpu<<<2,3>>>();
    cudaDeviceSynchronize();
    return 0;
}

Hello world from block 1 and thread 0!
Hello world from block 1 and thread 1!
Hello world from block 1 and thread 2!
Hello world from block 0 and thread 0!
Hello world from block 0 and thread 1!
Hello world from block 0 and thread 2!

每个线程执行一次“Hello world”

例子2：数组相加

#include <math.h>
#include <stdio.h>

void __global__ add(const double *x, const double *y, double *z, int count)
{
    const int n = blockDim.x * blockIdx.x + threadIdx.x;
	if( n < count)
	{
        printf("Block ID: %d, blockDim.x: %d,Thread ID: %d, n: %d\n",
               blockIdx.x, blockDim.x, threadIdx.x, n);
	    z[n] = x[n] + y[n];
        printf("z[%d]: %f\n", n, z[n]);
	}
    // 第二种xy输出
    // Add a barrier to ensure all threads have printed before exiting the kernel
    // __syncthreads();

    // // Only have one thread (threadIdx.x == 0) print the entire arrays
    // if (threadIdx.x == 0 && blockIdx.x == 0)
    // {
    //     for (int i = 0; i < count; ++i)
    //     {
    //         printf("Final Values - x[%d]: %f, y[%d]: %f, z[%d]: %f\n", i, x[i], i, y[i], i, z[i]);
    //     }
    // }

}
// void check(const double *z, const int N)
// {
//     bool error = false;
//     for (int n = 0; n < N; ++n)
//     {
//         if (fabs(z[n] - 3) > (1.0e-10))
//         {
//             error = true;
//         }
//     }
//     printf("%s\n", error ? "Errors" : "Pass");
// }
void check(const double *z, const double *h_x, const double *h_y, const int N)
{
    bool error = false;
    for (int n = 0; n < N; ++n)
    {
        // 修改此处的期望值计算方式
        // 不会执行，即error=false，返回Pass
        if (fabs(z[n] - (h_x[n] + h_y[n])) > (1.0e-10))
        {
            error = true;
        }
    }
    printf("%s\n", error ? "Errors" : "Pass");
}

int main(void)
{
    const int N = 6;
    const int M = sizeof(double) * N;
    double *h_x = (double*) malloc(M);
    double *h_y = (double*) malloc(M);
    double *h_z = (double*) malloc(M);

    for (int n = 0; n < N; ++n)
    {
        h_x[n] = 1;
        h_y[n] = n;
        //第一种xy输出
        printf("h_x[%d]: %f, h_y[%d]: %f\n", n, h_x[n], n, h_y[n]);
    }

    double *d_x, *d_y, *d_z;
    cudaMalloc((void **)&d_x, M);
    cudaMalloc((void **)&d_y, M);
    cudaMalloc((void **)&d_z, M);
    cudaMemcpy(d_x, h_x, M, cudaMemcpyHostToDevice);
    cudaMemcpy(d_y, h_y, M, cudaMemcpyHostToDevice);

    const int block_size = 2;
    /*N是数组大小
    如果N能够被block_size整除，那么grid_size就等于N/block_size,
    那么每个线程块处理一个数据。
    如果N不能被block_size整除，那么(N+block_size-1)/block_size
    的计算会确保生成的grid_size足够大，
    保所有的输入数据都能够被覆盖到，并且没有数据被漏掉
    如果 N = 10，block_size = 3，则 (N + block_size - 1) / block_size = 4。
    这确保了有足够的线程块来处理所有的10个数据元素。*/
    const int grid_size = (N + block_size - 1) / block_size;
    add<<<grid_size, block_size>>>(d_x, d_y, d_z,N);

    cudaMemcpy(h_z, d_z, M, cudaMemcpyDeviceToHost);
    // check(h_z, N);
    check(h_z, h_x, h_y, N);

    free(h_x);
    free(h_y);
    free(h_z);
    cudaFree(d_x);
    cudaFree(d_y);
    cudaFree(d_z);
    return 0;
}

// xy值
h_x[0]: 1.000000, h_y[0]: 0.000000
h_x[1]: 1.000000, h_y[1]: 1.000000
h_x[2]: 1.000000, h_y[2]: 2.000000
h_x[3]: 1.000000, h_y[3]: 3.000000
h_x[4]: 1.000000, h_y[4]: 4.000000
h_x[5]: 1.000000, h_y[5]: 5.000000
//线程块，线程，n的值
Block ID: 1, blockDim.x: 2, Thread ID: 0, n: 2
Block ID: 1, blockDim.x: 2, Thread ID: 1, n: 3
Block ID: 2, blockDim.x: 2, Thread ID: 0, n: 4
Block ID: 2, blockDim.x: 2, Thread ID: 1, n: 5
Block ID: 0, blockDim.x: 2, Thread ID: 0, n: 0
Block ID: 0, blockDim.x: 2, Thread ID: 1, n: 1
//输出的z值
z[2]: 3.000000
z[3]: 4.000000
z[4]: 5.000000
z[5]: 6.000000
z[0]: 1.000000
z[1]: 2.000000

Block ID: 0, blockDim.x: 2, Thread ID: 0, n: 0，即第1线程块，第1个线程，根据以下公式：

  const int n = blockDim.x * blockIdx.x + threadIdx.x; 
  z[n] = x[n] + y[n];

算出n=0，进一步，z[0]=x[0]+y[0]=1+0=1
即第一步算出索引，第二步带入索引值求出z，其余同理。
请注意，由于线程和块的索引计算，每个线程处理不同的数据元素，从而实现了并行计算。在实际的大规模并行计算中，这种并行性能够带来显著的加速。
如Cuda笔记1所讲的，所有线程执行相同的核函数，且并行执行。
每个线程负责处理一个计算任务。在GPU中，大量的线程可以同时执行相同的指令，每个线程独立地处理不同的数据或任务，这被称为SIMD（Single Instruction, Multiple Data）架构。（分治思想，将大的任务分成多个一样的小任务，每个thread处理一个计算任务）
例如，假设有一个大数组需要排序，可以将数组分割成若干段，每个线程负责对其中一段进行排序。这样，每个线程独立执行排序任务，最后再将各段有序的数组合并成整体有序的数组。

通过规约(Reduction)操作理解GPU内存体系（源于知乎-紫气东来）

reduce 算法可以描述为
在这里插入图片描述
其中 ×○可表示为求 sum，min，max，avg 等操作，最后获得的输出相比于输入一般维度上会递减。
在GPU中，reduce采用了一种树形的计算方式，例如下面的求和问题：
在这里插入图片描述
由于GPU没有针对global数据的同步操作，只能针对block的数据进行同步。所以，一般而言将reduce分为两个阶段，其示意图如下：
在这里插入图片描述
假设给定一个长度为N的数组，需要计算该数组的所有元素之和。首先需要将数组分为m个小份。而后，在第一阶段中，开启m个block计算出m个小份的reduce值。最后，在第二阶段中，使用一个block将m个小份再次进行reduce，得到最终的结果。
chatgpt的回答：
在GPU编程中，Reduce（归约）是一种并行计算模式，用于将一个大型数据集的元素合并为一个单一的值。Reduce操作通常用于计算总和、平均值、最大值、最小值等。在GPU中，由于存在大量的并行处理单元（CUDA核心），可以通过归约操作充分利用并行性，加速数据处理。

典型的GPU中的Reduce流程如下：

数据分割： 将大型数据集划分成小块，每个块由一个线程块（Thread Block）处理。

并行计算： 每个线程块中的线程并行地对其分配的数据块执行归约操作。这可以是对数据进行相加、取最大值、取最小值等。

层级归约： 由于线程块的数量通常大于硬件所支持的线程数，需要进行多层级的归约。每一层级的归约将当前层级的中间结果传递到下一层级，直到最终得到一个单一的结果。

最终结果： 最后一层的归约结果即为整个数据集的最终结果。

CUDA编程中，可以使用**共享内存（shared memory）**来加速归约操作，因为共享内存是线程块内的线程共享的。这减少了对全局内存的访问次数，提高了性能。

归约操作的实现可能因问题的不同而有所差异，但其核心思想是通过并行计算和层级归约来充分利用GPU的并行计算能力。

例子3：矩阵乘法（全局内存vs共享内存）

1）全局内存

#include <stdio.h>
#include <math.h>
#include "error.cuh"

#define BLOCK_SIZE 16


__global__ void gpu_matrix_mult(int *a,int *b, int *c, int m, int n, int k)
{ 
    //不是thread的一维索引，而是复合线程索引值，更符合矩阵运算等
    //索引就是变值，不要想复杂了，就相当于一维数组，二维数组中的索引一回事
    //一个索引对应一个线程（一个线程的位置，线程是最小单位），然后用这个线程计算，GPU只是很多线程同时计算而已。
    //计算当前线程在矩阵中的行索引。
    int row = blockIdx.y * blockDim.y + threadIdx.y; 
    int col = blockIdx.x * blockDim.x + threadIdx.x;
    int sum = 0;
    //确保当前线程的索引在矩阵的有效范围内，以避免越界访问。
    if( col < k && row < m) 
    {
        //对于每一列（或每一行，具体取决于矩阵乘法的顺序），
        //将对应位置的元素相乘并累加到 sum 中。
         //读取m的一行，n的一列
        for(int i = 0; i < n; i++) 
        {
            /*在实际的 GPU 编程中，选择线程块和网格的配置通常取决于算法和数据的特性，
            以及对硬件架构的优化。
            这儿矩阵相乘，因此选择二维的线程块和网格，线程块的 y 索引和 x 索引
            可以自然地映射到结果矩阵中的行和列。
            */ 
            sum += a[row * n + i] * b[i * k + col];
            //即每个(row,col)位置的线程处理一个计算任务
            /*
            与CPU的对比，
            void cpu_matrix_mult(int *h_a, int *h_b, int *h_result, int m, int n, int k) {
              for (int i = 0; i < m; ++i) 
              {
                 for (int j = 0; j < k; ++j) 
                 {
                    int tmp = 0.0;
                    for (int h = 0; h < n; ++h) 
                    {
                        tmp += h_a[i * n + h] * h_b[h * k + j];
                    }
                    h_result[i * k + j] = tmp;
                }
            }
        }
            事实上，外面的两个for循环用线程的动态索引取代了，而线程不是按顺序执行的
            而是多个线程同时执行的
            cpu中的i和j分别相当于GPU中的row，col 
            */
        }
        //row，col确定当前线程在输出矩阵 C 中的二维索引
        //用的一维表示二维a[i*M+j] i是行下标，M是一行的元素个数，j是列下标
        c[row * k + col] = sum;
    }
} 

__global__ void gpu_matrix_mult_shared(int *d_a, int *d_b, int *d_result, int n) 
{
    __shared__ int tile_a[BLOCK_SIZE][BLOCK_SIZE];
    __shared__ int tile_b[BLOCK_SIZE][BLOCK_SIZE];

    int row = blockIdx.y * BLOCK_SIZE + threadIdx.y;
    int col = blockIdx.x * BLOCK_SIZE + threadIdx.x;
    int tmp = 0;
    int idx;

    for (int sub = 0; sub < gridDim.x; ++sub) 
    {
        idx = row * n + sub * BLOCK_SIZE + threadIdx.x;
        tile_a[threadIdx.y][threadIdx.x] = row<n && (sub * BLOCK_SIZE + threadIdx.x)<n? d_a[idx]:0;
        idx = (sub * BLOCK_SIZE + threadIdx.y) * n + col;
        tile_b[threadIdx.y][threadIdx.x] = col<n && (sub * BLOCK_SIZE + threadIdx.y)<n? d_b[idx]:0;

        __syncthreads();
        for (int k = 0; k < BLOCK_SIZE; ++k) 
        {
            tmp += tile_a[threadIdx.y][k] * tile_b[k][threadIdx.x];
        }
        __syncthreads();
    }
    if(row < n && col < n)
    {
        d_result[row * n + col] = tmp;
    }
}

void cpu_matrix_mult(int *h_a, int *h_b, int *h_result, int m, int n, int k) {
    for (int i = 0; i < m; ++i) 
    {
        for (int j = 0; j < k; ++j) 
        {
            int tmp = 0.0;
            for (int h = 0; h < n; ++h) 
            {
                tmp += h_a[i * n + h] * h_b[h * k + j];
            }
            h_result[i * k + j] = tmp;
        }
    }
}

int main(int argc, char const *argv[])
{
    int m=1000;
    int n=1000;
    int k=1000;

    int *h_a, *h_b, *h_c, *h_cc, *h_cs;
    //void **) &h_a将分配的内存地址赋值给指针h_a 
    //sizeof(int)*m*n)分配一个 int 类型的数组，其大小为 m*n 个元素
    //即主机端分配一个 int 类型的二维数组，大小为 m 行、n 列。
    CHECK(cudaMallocHost((void **) &h_a, sizeof(int)*m*n));
    CHECK(cudaMallocHost((void **) &h_b, sizeof(int)*n*k));
    CHECK(cudaMallocHost((void **) &h_c, sizeof(int)*m*k));
    CHECK(cudaMallocHost((void **) &h_cc, sizeof(int)*m*k));
    CHECK(cudaMallocHost((void **) &h_cs, sizeof(int)*m*k));
    

    /*CUDA事件
    cudaEvent_t: 这是CUDA库中表示事件的数据类型。
    CUDA事件通常用于测量 GPU 执行时间，用于确定 GPU 执行某个任务所花费的时间。
    cudaEventCreate表示用于创建 CUDA 事件
    原型“cudaError_t cudaEventCreate(cudaEvent_t* event);”
    
    */
    cudaEvent_t start, stop,stop_share;
    CHECK(cudaEventCreate(&start));
    CHECK(cudaEventCreate(&stop));
    CHECK(cudaEventCreate(&stop_share));

    //将整个二维矩阵h_a赋值为1
    //数组是在内存按顺序存储，因此二维其实也是一维
    for (int i = 0; i < m; ++i) {
        for (int j = 0; j < n; ++j) {
            h_a[i * n + j] = 1;
        }
    }

    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < k; ++j) {
            h_b[i * k + j] = 0;
        }
    }

    int *d_a, *d_b, *d_c, *d_c_share;
    CHECK(cudaMalloc((void **) &d_a, sizeof(int)*m*n));
    CHECK(cudaMalloc((void **) &d_b, sizeof(int)*n*k));
    CHECK(cudaMalloc((void **) &d_c, sizeof(int)*m*k));
    CHECK(cudaMalloc((void **) &d_c_share, sizeof(int)*m*k));
    /*  cudaDeviceSynchronize()目的是确保设备上的所有操作都已完成。
    这段程序片段的目的是确保CUDA分配的设备内存已经成功完成。*/
    cudaDeviceSynchronize(); 
    // cudaEventRecord是CUDA中记录事件的函数，用于将当前时间戳记录到由start表示的事件中。
    CHECK(cudaEventRecord(start));
    // copy matrix A and B from host to device memory(cudaMemcpy全局内存)
    CHECK(cudaMemcpy(d_a, h_a, sizeof(int)*m*n, cudaMemcpyHostToDevice));
    CHECK(cudaMemcpy(d_b, h_b, sizeof(int)*n*k, cudaMemcpyHostToDevice));

    unsigned int grid_rows = (m + BLOCK_SIZE - 1) / BLOCK_SIZE;
    unsigned int grid_cols = (k + BLOCK_SIZE - 1) / BLOCK_SIZE;
    /*
    dim3三维网格和线程块，z=1，blockidx.z=0,threadidx.z=0
    创建了一个二维的线程块与线程网格
    */
    dim3 dimGrid(grid_cols, grid_rows);
    dim3 dimBlock(BLOCK_SIZE, BLOCK_SIZE);
    

    gpu_matrix_mult<<<dimGrid, dimBlock>>>(d_a, d_b, d_c, m,n,k);    

    CHECK(cudaMemcpy(h_c, d_c, (sizeof(int)*m*k), cudaMemcpyDeviceToHost));
    //cudaThreadSynchronize();
    CHECK(cudaEventRecord(stop));
    CHECK(cudaEventSynchronize(stop));
    
    gpu_matrix_mult_shared<<<dimGrid, dimBlock>>>(d_a, d_b, d_c_share, n);
    CHECK(cudaMemcpy(h_cs, d_c_share, (sizeof(int)*m*k), cudaMemcpyDeviceToHost));
    
    CHECK(cudaEventRecord(stop_share));
    CHECK(cudaEventSynchronize(stop_share));
    
    float elapsed_time, elapsed_time_share;
    CHECK(cudaEventElapsedTime(&elapsed_time, start, stop));
    CHECK(cudaEventElapsedTime(&elapsed_time_share, stop, stop_share));
    printf("Time_global = %g ms.\n", elapsed_time);
    printf("Time_share = %g ms.\n", elapsed_time_share);

    CHECK(cudaEventDestroy(start));
    CHECK(cudaEventDestroy(stop));    

    //cpu_matrix_mult(h_a, h_b, h_c, m, n, k);

    int ok = 1;
    /*for (int i = 0; i < m; ++i)
    { 
        for (int j = 0; j < k; ++j)
        {
            if(fabs(h_cs[i*k + j] - 0)>(1.0e-10))
            {
                printf("hcs: %d hc: %d  ",h_cs[i*k + j], h_c[i*k + j]);
                ok = 0;
            }
        }
    }

    if(ok)
    {
        printf("Pass!!!\n");
    }
    else
    {
        printf("Error!!!\n");
    }
    */
    // free memory
    CHECK(cudaFree(d_a));
    CHECK(cudaFree(d_b));
    CHECK(cudaFree(d_c));
    CHECK(cudaFreeHost(h_a));
    CHECK(cudaFreeHost(h_b));
    CHECK(cudaFreeHost(h_c));
    return 0;
}

重点：
①在实际的 GPU 编程中，选择线程块和网格的配置通常取决于算法和数据的特性，以及对硬件架构的优化。这儿矩阵相乘，因此选择二维的线程块和网格，线程块的 y 索引和 x 索引可以自然地映射到结果矩阵中的行和列。
②即每个(row,col)位置的线程处理一个计算任务。事实上，CPU外面的两个for循环用GPU线程的动态索引取代了，而线程不是按顺序执行的。而是多个线程同时执行的。CPU中的i和j分别相当于GPU中的row，col
③定义MK个thread，刚好动态映射_a，d_b的mn、nk，以及输出结果d_c中的mk个值
④copy matrix A and B from host to device memory(cudaMemcpy全局内存)
⑤逻辑上是二维的，实际内存物理上是一维的，因此要用一维索引形式表示

2）共享内存

理论：
在这里插入图片描述
整体思路：
申请很多线程，每个线程求取pd中的每个值，就需要从M和N读取一行和一列的值。那么要求出pdsut中每个值，就要读取m和n中by=1行，bx=1中的所有值。即每一行需要多很多次，放到共享内存。但是共享内存小，因此分块。做子矩阵（蓝色）线程。
1、在共享内存申请空间
2、pd相当于一个block，这里面的每个thread只读取一个元素，放入share memory中，所有线程同时行动
3、计算时，每个线程读取share memory读取行的元素（蓝色块中的黄色行）,结果存在temp中。然后黄色矩阵的行，结果放在temp中。计算完后，放到pd中。

#include <stdio.h>
#include <math.h>
#include "error.cuh"

#define BLOCK_SIZE 16


__global__ void gpu_matrix_mult(int *a,int *b, int *c, int m, int n, int k)
{ 
    //不是thread的一维索引，而是复合线程索引值，更符合矩阵运算等
    //索引就是变值，不要想复杂了，就相当于一维数组，二维数组中的索引一回事
    //一个索引对应一个线程（一个线程的位置，线程是最小单位），然后用这个线程计算，GPU只是很多线程同时计算而已。
    //计算当前线程在矩阵中的行索引。
    int row = blockIdx.y * blockDim.y + threadIdx.y; 
    int col = blockIdx.x * blockDim.x + threadIdx.x;
    int sum = 0;
    //确保当前线程的索引在矩阵的有效范围内，以避免越界访问。
    if( col < k && row < m) 
    {
        //对于每一列（或每一行，具体取决于矩阵乘法的顺序），
        //将对应位置的元素相乘并累加到 sum 中。
         //读取m的一行，n的一列
        for(int i = 0; i < n; i++) 
        {
            /*在实际的 GPU 编程中，选择线程块和网格的配置通常取决于算法和数据的特性，
            以及对硬件架构的优化。
            这儿矩阵相乘，因此选择二维的线程块和网格，线程块的 y 索引和 x 索引
            可以自然地映射到结果矩阵中的行和列。
            */ 
            sum += a[row * n + i] * b[i * k + col];
            //即每个(row,col)位置的线程处理一个计算任务
            /*
            与CPU的对比，
            void cpu_matrix_mult(int *h_a, int *h_b, int *h_result, int m, int n, int k) {
              for (int i = 0; i < m; ++i) 
              {
                 for (int j = 0; j < k; ++j) 
                 {
                    int tmp = 0.0;
                    for (int h = 0; h < n; ++h) 
                    {
                        tmp += h_a[i * n + h] * h_b[h * k + j];
                    }
                    h_result[i * k + j] = tmp;
                }
            }
        }
            事实上，外面的两个for循环用线程的动态索引取代了，而线程不是按顺序执行的
            而是多个线程同时执行的
            cpu中的i和j分别相当于GPU中的row，col 
            */
        }
        //row，col确定当前线程在输出矩阵 C 中的二维索引
        //用的一维表示二维a[i*M+j] i是行下标，M是一行的元素个数，j是列下标
        c[row * k + col] = sum;
    }
} 

__global__ void gpu_matrix_mult_shared(int *d_a, int *d_b, int *d_result, int M, int N, int K) 
{
    /*
    __shared__: 这是一个CUDA关键字，
    表示声明的变量是共享内存中的变量。共享内存是一块对于线程块中的所有线程可见的内存。
    int tile_a[BLOCK_SIZE][BLOCK_SIZE];这是一个二维整数数组，
    表示在共享内存中分配的二维数组。tile_a 被用于存储一个矩阵的子矩阵。
    怎么理解被用于存储一个矩阵的子矩阵？
    矩阵乘法通常涉及到将大矩阵分解为子矩阵，
    然后每个线程块负责计算一个或多个子矩阵的乘法。
    在矩阵乘法的上下文中，tile_a 存储了一个矩阵的子矩阵。当每个线程块负责计算某个子矩阵的乘法时，
    将这个子矩阵加载到 tile_a 中，以便线程块中的每个线程可以快速访问和共享这部分数据。
    
    */
    __shared__ int tile_a[BLOCK_SIZE][BLOCK_SIZE];
    __shared__ int tile_b[BLOCK_SIZE][BLOCK_SIZE];

    int row = blockIdx.y * BLOCK_SIZE + threadIdx.y;
    int col = blockIdx.x * BLOCK_SIZE + threadIdx.x;
    int tmp = 0;
    int idx;
    //读取每个矩阵小块
    // for (int sub = 0; sub < gridDim.x; ++sub) 
    for (int sub = 0; sub <= N/BLOCK_SIZE; ++sub) 
    {
        // idx = row * n + sub * BLOCK_SIZE + threadIdx.x;
        // tile_a[threadIdx.y][threadIdx.x] = row<n && (sub * BLOCK_SIZE + threadIdx.x)<n? d_a[idx]:0;
        // idx = (sub * BLOCK_SIZE + threadIdx.y) * n + col;
        // tile_b[threadIdx.y][threadIdx.x] = col<n && (sub * BLOCK_SIZE + threadIdx.y)<n? d_b[idx]:0;
        //详细版
        int r = row;
        int c = sub * BLOCK_SIZE + threadIdx.x;
        idx = r * N + c;

        if (r >= M || c >= N)
        {
            tile_a[threadIdx.y][threadIdx.x]=0;
        }
        else
        {
            tile_a[threadIdx.y][threadIdx.x]=d_a[idx];
        }

        r = sub * BLOCK_SIZE + threadIdx.y;
        c = col;
        idx = r * K + c;
        if (c >= K || c >= N)
        {
            tile_b[threadIdx.y][threadIdx.x]=0;
        }
        else
        {
            tile_b[threadIdx.y][threadIdx.x]=d_b[idx];
        }

        __syncthreads();

        //相乘累加到一个temp中
        for (int k = 0; k < BLOCK_SIZE; ++k) 
        {
            tmp += tile_a[threadIdx.y][k] * tile_b[k][threadIdx.x];
        }
        __syncthreads();
    }

    if(row < M && col < K)
    {
        d_result[row * K + col] = tmp;
    }
}

void cpu_matrix_mult(int *h_a, int *h_b, int *h_result, int m, int n, int k) {
    for (int i = 0; i < m; ++i) 
    {
        for (int j = 0; j < k; ++j) 
        {
            int tmp = 0.0;
            for (int h = 0; h < n; ++h) 
            {
                tmp += h_a[i * n + h] * h_b[h * k + j];
            }
            h_result[i * k + j] = tmp;
        }
    }
}

int main(int argc, char const *argv[])
{
    int m=1000;
    int n=1000;
    int k=1000;

    int *h_a, *h_b, *h_c, *h_cc, *h_cs;
    //void **) &h_a将分配的内存地址赋值给指针h_a 
    //sizeof(int)*m*n)分配一个 int 类型的数组，其大小为 m*n 个元素
    //即主机端分配一个 int 类型的二维数组，大小为 m 行、n 列。
    CHECK(cudaMallocHost((void **) &h_a, sizeof(int)*m*n));
    CHECK(cudaMallocHost((void **) &h_b, sizeof(int)*n*k));
    CHECK(cudaMallocHost((void **) &h_c, sizeof(int)*m*k));
    CHECK(cudaMallocHost((void **) &h_cc, sizeof(int)*m*k));
    CHECK(cudaMallocHost((void **) &h_cs, sizeof(int)*m*k));
    

    /*CUDA事件
    cudaEvent_t: 这是CUDA库中表示事件的数据类型。
    CUDA事件通常用于测量 GPU 执行时间，用于确定 GPU 执行某个任务所花费的时间。
    cudaEventCreate表示用于创建 CUDA 事件
    原型“cudaError_t cudaEventCreate(cudaEvent_t* event);”
    
    */
    cudaEvent_t start, stop,stop_share;
    CHECK(cudaEventCreate(&start));
    CHECK(cudaEventCreate(&stop));
    CHECK(cudaEventCreate(&stop_share));

    //将整个二维矩阵h_a赋值为1
    //数组是在内存按顺序存储，因此二维其实也是一维
    for (int i = 0; i < m; ++i) {
        for (int j = 0; j < n; ++j) {
            h_a[i * n + j] = 1;
        }
    }

    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < k; ++j) {
            h_b[i * k + j] = 0;
        }
    }

    int *d_a, *d_b, *d_c, *d_c_share;
    CHECK(cudaMalloc((void **) &d_a, sizeof(int)*m*n));
    CHECK(cudaMalloc((void **) &d_b, sizeof(int)*n*k));
    CHECK(cudaMalloc((void **) &d_c, sizeof(int)*m*k));
    CHECK(cudaMalloc((void **) &d_c_share, sizeof(int)*m*k));
    /*  cudaDeviceSynchronize()目的是确保设备上的所有操作都已完成。
    这段程序片段的目的是确保CUDA分配的设备内存已经成功完成。*/
    cudaDeviceSynchronize(); 
    // cudaEventRecord是CUDA中记录事件的函数，用于将当前时间戳记录到由start表示的事件中。
    CHECK(cudaEventRecord(start));
    // copy matrix A and B from host to device memory(cudaMemcpy全局内存)
    CHECK(cudaMemcpy(d_a, h_a, sizeof(int)*m*n, cudaMemcpyHostToDevice));
    CHECK(cudaMemcpy(d_b, h_b, sizeof(int)*n*k, cudaMemcpyHostToDevice));

    unsigned int grid_rows = (m + BLOCK_SIZE - 1) / BLOCK_SIZE;
    unsigned int grid_cols = (k + BLOCK_SIZE - 1) / BLOCK_SIZE;
    /*
    dim3三维网格和线程块，z=1，blockidx.z=0,threadidx.z=0
    创建了一个二维的线程块与线程网格
    */
    dim3 dimGrid(grid_cols, grid_rows);
    dim3 dimBlock(BLOCK_SIZE, BLOCK_SIZE);
    

    gpu_matrix_mult<<<dimGrid, dimBlock>>>(d_a, d_b, d_c, m,n,k);    

    CHECK(cudaMemcpy(h_c, d_c, (sizeof(int)*m*k), cudaMemcpyDeviceToHost));
    /*使用这两个函数可以计算CUDA代码的执行时间，
    方法是在cudaEventRecord(start)和cudaEventRecord(stop)之间执行CUDA代码
    cudaEventSynchronize(stop)：等待stop事件完成，确保GPU执行到此处时
    已完成前面的所有操作。
    这是一种同步操作，保证我们能够准确地测量CUDA代码的执行时间。
    */
    CHECK(cudaEventRecord(stop));
    CHECK(cudaEventSynchronize(stop));
    
   
    gpu_matrix_mult_shared<<<dimGrid, dimBlock>>> (d_a, d_b, d_c_share, m,n,k);
    CHECK(cudaMemcpy(h_cs, d_c_share, (sizeof(int)*m*k), cudaMemcpyDeviceToHost));
    
    CHECK(cudaEventRecord(stop_share));
    CHECK(cudaEventSynchronize(stop_share));
    
    //通过cudaEventElapsedTime计算两个事件之间的时间差。
    float elapsed_time, elapsed_time_share;
    CHECK(cudaEventElapsedTime(&elapsed_time, start, stop));
    CHECK(cudaEventElapsedTime(&elapsed_time_share, stop, stop_share));
    printf("Time_global = %g ms.\n", elapsed_time);
    printf("Time_share = %g ms.\n", elapsed_time_share);
    
    /*
cudaEventDestroy函数用于销毁CUDA事件。在CUDA程序中，为了避免内存泄漏，
应该在不再需要使用CUDA事件的地方调用cudaEventDestroy来释放相应的资源。*/
    CHECK(cudaEventDestroy(start));
    CHECK(cudaEventDestroy(stop));    

    //cpu_matrix_mult(h_a, h_b, h_c, m, n, k);

    // int ok = 1;
    /*for (int i = 0; i < m; ++i)
    { 
        for (int j = 0; j < k; ++j)
        {
            if(fabs(h_cs[i*k + j] - 0)>(1.0e-10))
            {
                printf("hcs: %d hc: %d  ",h_cs[i*k + j], h_c[i*k + j]);
                ok = 0;
            }
        }
    }

    if(ok)
    {
        printf("Pass!!!\n");
    }
    else
    {
        printf("Error!!!\n");
    }
    */
    // free memory
    CHECK(cudaFree(d_a));
    CHECK(cudaFree(d_b));
    CHECK(cudaFree(d_c));
    CHECK(cudaFreeHost(h_a));
    CHECK(cudaFreeHost(h_b));
    CHECK(cudaFreeHost(h_c));
    return 0;
}

解释部分代码：
在这里插入图片描述

参考链接：
https://zhuanlan.zhihu.com/p/654027980
https://www.bilibili.com/video/BV17T4y117vK?p=5

tensor.shape
关注
  • 18
    点赞
  • 踩
  • 16
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
专栏目录
cuda 排序算法笔记
jacke121的专栏
03-30 2691
Thrust是cuda自带的c++库，cuda安装好之后，这个库也默认安装好了。 这个库基本是采用类似STL的接口方式，因此对于开发者非常友好，开发者不再需要关注内存与显存相关的问题了。 排序入门例子 #include <thrust/host_vector.h> #include <thrust/device_vector.h> #include <thrust/generate.h> #include <thrust/sort.h> #incl..
CUDA 笔记本电脑 安装详细说明
03-16
这是我在做毕业设计中用笔记本安装CUDA的详细说明，本人笔记本式索尼VAIO系列的电脑，NVIDIA公司的CUDA笔记本驱动不支持该系列电脑，故本人做了一个安装说明。谢谢。本人还会上传一些附件的安装软件。
Cuda笔记1
qq_44576434的博客
01-21 425
培训001 1 1+…+100，CPU是串行执行，GPU是分成几部分同时计算，如1+2+3,4+5+6… 培训002 一来一回 每种定义有对应的调用位置，和执行位置，不对会报错。 下图是用NVPROF时间分析 下图是资源分析 CUDA线程层次 一个SM执行多个block 所有线程执行同一个kernel，每个线程处理的数据不一样 线程在cuda core里面，block以块为单位分配到sm中，grid在device中 如图描述的是：9个block分给3个sm 这张ppt的流程与下面这个
cudaMemcpy学习笔记
jacke121的专栏
07-31 1667
cudaMemcpy学习笔记
CUDA学习笔记
weixin_47885125的博客
12-31 1056
期末抱佛脚
CUDA笔记-卷积计算
黑山老妖的博客
04-26 4362
CUDA 3D convolution - ijpq - 博客园overview https://www.cnblogs.com/ijpq/p/15405106.html
我的CUDA学习笔记
08-29
"CUDA学习笔记" 本资源摘要信息中，我们将对CUDA学习笔记进行详细的知识点梳理，涵盖了CUDA程序架构、Thread管理、存储器体系、混合编程、CUDA程序接口描述等方面的知识点。 1. CUDA 程序架构 CUDA 程序架构中...
CUDA-:CUDA C学习笔记
03-16
这篇学习笔记将深入探讨CUDA C的基础知识和关键概念。 1. **CUDA架构** - **CUDA核心**: GPU内部包含大量的流处理器（Streaming Multiprocessors, SM），每个SM由多个CUDA核心组成，它们执行并行计算任务。 - **...
CUDA学习笔记，GPU编程，C++编程
06-08
本篇笔记将深入探讨CUDA编程的基础概念，包括GPU架构、并行计算原理以及CUDA编程的关键要素。 首先，我们需要理解为什么需要GPU进行并行计算。传统的CPU（中央处理器）擅长于串行处理任务，但对于大规模并行计算，...
【yolov4】学习笔记一
12-21
本篇学习笔记主要围绕YOLOv4在PyTorch环境中的实现进行探讨，同时也涉及CUDA的使用。 首先，YOLOv4的核心改进包括引入了多种技术，如 Mish 激活函数、SPP-Block（Spatial Pyramid Pooling）、CBAM（Channel ...
deepseek-vl 论文阅读笔记
samoyan的博客,记录技术成长~
07-22 1223
我们的语言模型基于DeepSeek LLM（DeepSeek-AI，2024），其微设计大体遵循LLaMA（Touvron等，2023a,b）的设计，采用带有RMSNorm（Zhang和Sennrich，2019）函数的Pre-Norm结构，并使用SwiGLU（Shazeer，2020）作为前馈网络（FFN）的激活函数，中间层维度为8/3模型维度。此外，我们引入了一种新的“模态预热”策略。为了促进创新并支持广泛的应用需求，我们公开了两个版本的模型，分别为1.3B和7B，以满足不同计算能力的需求。
vite5+vue3开发阅读APP实战笔记20240725
Python私教
07-25 386
修改src/page/home/index.vue，使用组件和API。注意，这里无论是组件还是API，都不需要我们手动引入，我们配置了自动导入以后，vite会帮我们自动导入。创建 src/page/home/index.vue。修改 vite.config.js。修改vite.config.js。修改 src/main.js。修改src/App.vue。
Golang学习笔记20240725，Go语言基础语法
最新发布
Python私教
07-25 531
如果我们在遍历的时候，不想要key或者value，也可以用下划线替代，下划线叫做匿名变量。go语言里面的switch语句每个case天生就是独立的，不需要加break。指针可以直接对内存地址进行操作。使用*表示指针，使用&取地址。切片在go语言里面非常常用，因为其有动态扩展的特性。运行方式3：goland右键运行。使用加号可以对字符串进行拼接。
学习笔记-系统框图传递函数公式推导
weixin_45875994的博客
07-24 789
输出C(s)是前向路径（通过G(s)）的输出与反馈路径（通过H(s)）的输出之差（在负反馈系统中）。因此，我们可以写出以下关系式：这里，G(s)R(s)是前向路径的输出，而G(s)H(s)C(s)是反馈路径的输出（注意，反馈路径的输出先经过H(s)，再与前向路径的输出相减）。
谷粒商城实战笔记-50-51-商品分类的删除
epitomizelu的专栏
07-23 558
这两节的内容是开发商品分类的删除功能。
Typora笔记上传到CSDN
m0_48362854的博客
07-25 524
我花了一个小时弄这个 找了很多csdn的笔记 都弄完感觉最方便的还是直接在csdn上写（因为有的之前不是用的阿里云oss 而是gitee 免费仓库 但是好像后来不能使用了 就会导致你上传的笔记图片不显示 很麻烦 使用阿里云或者腾讯云的话 需要花钱 虽然两块钱用一年 但是感觉也是很麻烦 不如直接在网页上写写算了 当然如果笔记图片很多的 还是建议配置一下 ）参考这个配置完阿里云oss 没什么问题 但是PicGo也是参考他的 因为版本是旧的 可以点击我下边的官方链接下载。图片的链接是本地的 当然没法显示。
从零开始学Java（超详细韩顺平老师笔记梳理）08——面向对象编程中级（上）IDEA常用快捷键、包、封装、继承
apple_68589597的博客
07-23 994
本篇文章是韩顺平老师课程P264-P296（第8章：面向对象编程中级）上半部分重点内容的详细总结，包括IDEA常用快捷键与模板、包的使用与规范、访问修饰符、封装（重点）、继承（重点）
SAP PP学习笔记31 - 计划运行的步骤2 - Scheduling(日程计算)，BOM Explosion(BOM展开）
shi_ly的专栏
07-24 764
本章主要讲Scheduling(日程计算),其中包括 外部调达中的调度，生产订单中的调度，修改需求日期(其实就是在物料中增加安全时间和BOM中的提前期偏置量)。另外还讲了BOM展开的内容，详细讲了BOM选择的各种方法，除了按数量选BOM之外，还有其他比如 按日期，按生产版本等，而其中按日期又分为按订单开始日期，订单终了日期等，越趋复杂。
InternLM学习笔记
haveanybody的博客
07-24 222
【代码】InternLM学习笔记。
笔记本电脑安装cuda
09-09
3. 安装过程中，按照安装向导的指示进行操作。可能需要选择自定义安装选项来选择需要安装的组件。确保选择正确的组件和版本。 4. 完成安装后，重启笔记本电脑以使更改生效。 请注意，这只是一个概述的步骤，具体的...
评论 1
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值