CUDA中嵌套复数向量GPU复制与FFT正确实现及异常排查

Windows平台C++项目，使用CUDA 10.1的cuFFT库执行设备端FFT计算。现有嵌套复数向量std::vector<std::vector<std::complex<float>>> unitedVector，需将其复制到cudaMalloc分配的GPU指针std::complex<float>* gpu_data中执行FFT。当前代码前3-4次计算结果与Matlab一致，但之后会触发异常，怀疑数据复制或FFT实现存在问题。

原初始化代码：

cufftHandle plan;
int data_size = unitedVector[0].size();  // 每个子向量的数据长度
int num_channels = unitedVector.size();  // 通道数（子向量数量）
std::complex<float>* gpu_data;

cufftPlan1d(&plan, data_size, CUFFT_C2C, num_channels); // 创建FFT计划
cudaMalloc((void**)&gpu_data, num_channels * data_size * sizeof(std::complex<float>)); // GPU内存分配

原数据复制代码：

cudaMemcpy(gpu_data, unitedVector[0].data(), num_channels * data_size * sizeof(std::complex<float>), cudaMemcpyHostToDevice);

原FFT执行与fftshift代码：

cufftResult result = cufftExecC2C(plan, reinterpret_cast<cufftComplex*>(gpu_data), reinterpret_cast<cufftComplex*>(gpu_data), CUFFT_FORWARD);

// 频谱旋转（fftshift）
int half_size = data_size / 2;
for (int i = 0; i < num_channels; ++i) {
    cudaMemcpy(unitedVector[i].data(), gpu_data + i * data_size + half_size, half_size * sizeof(std::complex<float>), cudaMemcpyDeviceToHost);
    cudaMemcpy(unitedVector[i].data() + half_size, gpu_data + i * data_size, half_size * sizeof(std::complex<float>), cudaMemcpyDeviceToHost);
}

1. 核心数据复制错误 #

原代码直接使用unitedVector[0].data()作为源地址，一次性复制所有通道的数据，但std::vector<std::vector<T>>的子向量在内存中是非连续存储的，这会导致除第一个通道外，其余通道的数据都是无效的垃圾值，多次运行后因内存访问越界触发异常。

修复方式：逐个复制每个通道的数据到GPU的对应内存区域：

// 逐个复制每个通道的宿主数据到GPU
for (int i = 0; i < num_channels; ++i) {
    cudaError_t err = cudaMemcpy(
        gpu_data + i * data_size,
        unitedVector[i].data(),
        data_size * sizeof(std::complex<float>),
        cudaMemcpyHostToDevice
    );
    if (err != cudaSuccess) {
        fprintf(stderr, "cudaMemcpy HostToDevice failed for channel %d: %s\n", i, cudaGetErrorString(err));
        cudaFree(gpu_data);
        cufftDestroy(plan);
        return;
    }
}

2. 缺失错误检查 #

原代码未对任何CUDA/cuFFT操作的返回值做检查，无法定位异常触发的具体环节（比如内存分配失败、计划创建失败、FFT执行错误等），多次运行后可能因资源未正确释放导致内存耗尽。

修复方式：为所有CUDA/cuFFT操作添加错误检查：

// 检查GPU内存分配
cudaError_t cuda_err = cudaMalloc((void**)&gpu_data, num_channels * data_size * sizeof(std::complex<float>));
if (cuda_err != cudaSuccess) {
    fprintf(stderr, "cudaMalloc failed: %s\n", cudaGetErrorString(cuda_err));
    return;
}

// 检查cuFFT计划创建
cufftResult plan_err = cufftPlan1d(&plan, data_size, CUFFT_C2C, num_channels);
if (plan_err != CUFFT_SUCCESS) {
    fprintf(stderr, "cufftPlan1d failed: %d\n", plan_err);
    cudaFree(gpu_data);
    return;
}

// 检查FFT执行结果
cufftResult exec_err = cufftExecC2C(plan, reinterpret_cast<cufftComplex*>(gpu_data), reinterpret_cast<cufftComplex*>(gpu_data), CUFFT_FORWARD);
if (exec_err != CUFFT_SUCCESS) {
    fprintf(stderr, "cufftExecC2C failed: %d\n", exec_err);
    cudaFree(gpu_data);
    cufftDestroy(plan);
    return;
}

3. FFTShift的奇数长度处理缺陷 #

原代码中half_size = data_size / 2仅适用于偶数长度的FFT，若data_size为奇数，会导致数据复制不完整，残留未覆盖的内存区域，引发未定义行为。

修复方式：适配奇数长度的fftshift：

int half_size = (data_size + 1) / 2;
int remaining_size = data_size - half_size;

for (int i = 0; i < num_channels; ++i) {
    // 复制后半段到前半段
    cudaError_t err1 = cudaMemcpy(
        unitedVector[i].data(),
        gpu_data + i * data_size + half_size,
        remaining_size * sizeof(std::complex<float>),
        cudaMemcpyDeviceToHost
    );
    // 复制前半段到后半段
    cudaError_t err2 = cudaMemcpy(
        unitedVector[i].data() + remaining_size,
        gpu_data + i * data_size,
        half_size * sizeof(std::complex<float>),
        cudaMemcpyDeviceToHost
    );
    if (err1 != cudaSuccess || err2 != cudaSuccess) {
        fprintf(stderr, "cudaMemcpy DeviceToHost failed for channel %d\n", i);
        cudaFree(gpu_data);
        cufftDestroy(plan);
        return;
    }
}

4. 资源泄漏问题 #

原代码未在使用完毕后释放GPU内存和cuFFT计划，多次运行会导致资源耗尽，触发异常。

修复方式：在计算完成后释放资源：

// 释放GPU内存
cudaFree(gpu_data);
// 销毁cuFFT计划
cufftDestroy(plan);

#include <cufft.h>
#include <vector>
#include <complex>
#include <cstdio>

void performFFT(std::vector<std::vector<std::complex<float>>>& unitedVector) {
    if (unitedVector.empty() || unitedVector[0].empty()) {
        fprintf(stderr, "Input vector is empty\n");
        return;
    }

    int data_size = unitedVector[0].size();
    int num_channels = unitedVector.size();
    std::complex<float>* gpu_data = nullptr;
    cufftHandle plan;

    // 分配GPU内存并检查错误
    cudaError_t cuda_err = cudaMalloc((void**)&gpu_data, num_channels * data_size * sizeof(std::complex<float>));
    if (cuda_err != cudaSuccess) {
        fprintf(stderr, "cudaMalloc failed: %s\n", cudaGetErrorString(cuda_err));
        return;
    }

    // 创建cuFFT计划并检查错误
    cufftResult plan_err = cufftPlan1d(&plan, data_size, CUFFT_C2C, num_channels);
    if (plan_err != CUFFT_SUCCESS) {
        fprintf(stderr, "cufftPlan1d failed: %d\n", plan_err);
        cudaFree(gpu_data);
        return;
    }

    // 逐个复制通道数据到GPU
    for (int i = 0; i < num_channels; ++i) {
        cuda_err = cudaMemcpy(gpu_data + i * data_size, unitedVector[i].data(), data_size * sizeof(std::complex<float>), cudaMemcpyHostToDevice);
        if (cuda_err != cudaSuccess) {
            fprintf(stderr, "cudaMemcpy HostToDevice failed for channel %d: %s\n", i, cudaGetErrorString(cuda_err));
            cudaFree(gpu_data);
            cufftDestroy(plan);
            return;
        }
    }

    // 执行FFT并检查错误
    cufftResult exec_err = cufftExecC2C(plan, reinterpret_cast<cufftComplex*>(gpu_data), reinterpret_cast<cufftComplex*>(gpu_data), CUFFT_FORWARD);
    if (exec_err != CUFFT_SUCCESS) {
        fprintf(stderr, "cufftExecC2C failed: %d\n", exec_err);
        cudaFree(gpu_data);
        cufftDestroy(plan);
        return;
    }

    // 执行FFTShift，适配奇偶长度
    int half_size = (data_size + 1) / 2;
    int remaining_size = data_size - half_size;
    for (int i = 0; i < num_channels; ++i) {
        cuda_err = cudaMemcpy(unitedVector[i].data(), gpu_data + i * data_size + half_size, remaining_size * sizeof(std::complex<float>), cudaMemcpyDeviceToHost);
        if (cuda_err != cudaSuccess) {
            fprintf(stderr, "cudaMemcpy DeviceToHost failed for channel %d (first part): %s\n", i, cudaGetErrorString(cuda_err));
            cudaFree(gpu_data);
            cufftDestroy(plan);
            return;
        }
        cuda_err = cudaMemcpy(unitedVector[i].data() + remaining_size, gpu_data + i * data_size, half_size * sizeof(std::complex<float>), cudaMemcpyDeviceToHost);
        if (cuda_err != cudaSuccess) {
            fprintf(stderr, "cudaMemcpy DeviceToHost failed for channel %d (second part): %s\n", i, cudaGetErrorString(cuda_err));
            cudaFree(gpu_data);
            cufftDestroy(plan);
            return;
        }
    }

    // 释放资源
    cudaFree(gpu_data);
    cufftDestroy(plan);
}

内容的提问来源于stack exchange，提问作者joe_dowson