并行编程与 GPU 编程核心笔记 + C++示例代码

#include <iostream>
#include <vector>
#include <chrono>
#include <cuda_runtime.h>

using namespace std;
using namespace chrono;

// 矩阵大小（需为块大小的整数倍，简化示例）
const int MATRIX_SIZE = 2048;
const int BLOCK_SIZE = 32;  // CUDA块大小，适配SM的并行能力

// CPU矩阵乘法（串行，用于对比）
void cpu_matrix_mult(const vector<float>& A, const vector<float>& B, vector<float>& C) {
    for (int i = 0; i < MATRIX_SIZE; ++i) {
        for (int j = 0; j < MATRIX_SIZE; ++j) {
            float sum = 0.0f;
            for (int k = 0; k < MATRIX_SIZE; ++k) {
                sum += A[i * MATRIX_SIZE + k] * B[k * MATRIX_SIZE + j];
            }
            C[i * MATRIX_SIZE + j] = sum;
        }
    }
}

// GPU核函数：矩阵乘法（使用共享内存优化）
__global__ void gpu_matrix_mult(const float* d_A, const float* d_B, float* d_C) {
    // 共享内存：存储A的子块和B的子块（快于全局内存）
    __shared__ float shared_A[BLOCK_SIZE][BLOCK_SIZE];
    __shared__ float shared_B[BLOCK_SIZE][BLOCK_SIZE];

    // 当前线程在矩阵中的全局坐标
    int row = blockIdx.y * blockDim.y + threadIdx.y;
    int col = blockIdx.x * blockDim.x + threadIdx.x;

    float sum = 0.0f;

    // 分块加载A和B到共享内存（流水线并行思想）
    for (int tile = 0; tile < MATRIX_SIZE / BLOCK_SIZE; ++tile) {
        // 加载A的子块到共享内存（线程分工）
        shared_A[threadIdx.y][threadIdx.x] = d_A[row * MATRIX_SIZE + tile * BLOCK_SIZE + threadIdx.x];
        // 加载B的子块到共享内存
        shared_B[threadIdx.y][threadIdx.x] = d_B[(tile * BLOCK_SIZE + threadIdx.y) * MATRIX_SIZE + col];

        __syncthreads();  // 等待块内所有线程加载完成（同步）

        // 子块内乘法累加（SIMD：同一指令作用于不同共享内存数据）
        for (int k = 0; k < BLOCK_SIZE; ++k) {
            sum += shared_A[threadIdx.y][k] * shared_B[k][threadIdx.x];
        }

        __syncthreads();  // 等待块内所有线程计算完成，避免覆盖共享内存
    }

    // 将结果写入全局内存
    d_C[row * MATRIX_SIZE + col] = sum;
}

int main() {
    // 1. 初始化CPU端矩阵（宿主内存）
    vector<float> h_A(MATRIX_SIZE * MATRIX_SIZE, 1.0f);  // 矩阵A全1
    vector<float> h_B(MATRIX_SIZE * MATRIX_SIZE, 1.0f);  // 矩阵B全1
    vector<float> h_C_cpu(MATRIX_SIZE * MATRIX_SIZE, 0.0f);
    vector<float> h_C_gpu(MATRIX_SIZE * MATRIX_SIZE, 0.0f);

    // 2. CPU矩阵乘法（串行）
    auto start = high_resolution_clock::now();
    cpu_matrix_mult(h_A, h_B, h_C_cpu);
    auto cpu_time = duration_cast<seconds>(high_resolution_clock::now() - start).count();
    cout << "CPU矩阵乘法耗时：" << cpu_time << "s" << endl;

    // 3. GPU矩阵乘法（并行）
    float *d_A, *d_B, *d_C;
    // 分配GPU设备内存（全局内存）
    cudaMalloc((void**)&d_A, MATRIX_SIZE * MATRIX_SIZE * sizeof(float));
    cudaMalloc((void**)&d_B, MATRIX_SIZE * MATRIX_SIZE * sizeof(float));
    cudaMalloc((void**)&d_C, MATRIX_SIZE * MATRIX_SIZE * sizeof(float));

    // 从CPU内存拷贝数据到GPU全局内存（文档：需手动管理内存传输）
    cudaMemcpy(d_A, h_A.data(), MATRIX_SIZE * MATRIX_SIZE * sizeof(float), cudaMemcpyHostToDevice);
    cudaMemcpy(d_B, h_B.data(), MATRIX_SIZE * MATRIX_SIZE * sizeof(float), cudaMemcpyHostToDevice);

    // 定义CUDA网格和块维度
    dim3 block(BLOCK_SIZE, BLOCK_SIZE);  // 块内线程数：32x32=1024（适配SM能力）
    dim3 grid(MATRIX_SIZE / BLOCK_SIZE, MATRIX_SIZE / BLOCK_SIZE);  // 网格大小

    // GPU核函数执行（并行计算）
    start = high_resolution_clock::now();
    gpu_matrix_mult<<<grid, block>>>(d_A, d_B, d_C);
    cudaDeviceSynchronize();  // 等待GPU执行完成
    auto gpu_time = duration_cast<milliseconds>(high_resolution_clock::now() - start).count();
    cout << "GPU矩阵乘法耗时：" << gpu_time << "ms" << endl;

    // 从GPU全局内存拷贝结果到CPU内存
    cudaMemcpy(h_C_gpu.data(), d_C, MATRIX_SIZE * MATRIX_SIZE * sizeof(float), cudaMemcpyDeviceToHost);

    // 4. 验证结果（简单对比前10个元素）
    bool result_valid = true;
    for (int i = 0; i < 10; ++i) {
        if (abs(h_C_cpu[i] - h_C_gpu[i]) > 1e-5) {
            result_valid = false;
            break;
        }
    }
    cout << "结果验证：" << (result_valid ? "正确" : "错误") << endl;
    cout << "GPU加速比：" << (double)cpu_time * 1000 / gpu_time << "x" << endl;

    // 释放GPU内存
    cudaFree(d_A);
    cudaFree(d_B);
    cudaFree(d_C);

    return 0;
}