0

我正在尝试复制 matlab fft 功能，它对矩阵进行逐行（或逐列）fft。每行将是袖带计划中的批次之一。

我可以使用 cufftExecC2C（下面代码中注释掉的部分有效）使其工作，但不能使用 cufftExecR2C。我的代码使用 cufftPlan1d，但理想情况下我想使用 cufftPlanMany 来实现它。

我想知道我做错了什么，以及是否有更好的方法来做到这一点。谢谢你。

// linker -> input -> additional dependencies -> add 'cufft.lib'
// VC++ Directories -> include directories - > add 'C:\ProgramData\NVIDIA Corporation\CUDA Samples\v6.0\common\inc'

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

#include <iostream>

#define NX 6
#define NY 5

void printArray(float *my_array);
void printComplexArray(float2 *my_array);

int main(){

/************************************************************ C2C ************************************************************/
/*  
    float2 *initial_array = (float2 *)malloc(sizeof(float2) * NX * NY);
    for (int h = 0; h < NX; h++){
        for (int w = 0; w < NY; w++){
            initial_array[NY * h + w].x = 0;
            initial_array[NY * h + w].y = 0;
        }
    }
    initial_array[NY*3 + 0].x = 1;
    initial_array[NY*5 + 0].x = 1;

    printComplexArray(initial_array);

    float2 *transformed_array= (float2 *)malloc(sizeof(float2) * NX * NY);

    cufftComplex *gpu_initial_array;
    cufftComplex *gpu_transformed_array;

    cudaMalloc((void **)&gpu_initial_array, NX*NY*sizeof(cufftComplex));
    cudaMalloc((void **)&gpu_transformed_array, NX*NY*sizeof(cufftComplex));

    cudaMemcpy(gpu_initial_array, initial_array, NX*NY*sizeof(float2), cudaMemcpyHostToDevice);

    cufftHandle plan;
    cufftPlan1d(&plan, NY, CUFFT_C2C, NX);

    cufftExecC2C(plan, gpu_initial_array, gpu_transformed_array, CUFFT_FORWARD);

    cudaMemcpy(transformed_array, gpu_transformed_array, NX*NY*sizeof(cufftComplex), cudaMemcpyDeviceToHost);

    printComplexArray(transformed_array);
*/
/************************************************************ C2C ************************************************************/

/************************************************************ R2C ************************************************************/

    float *initial_array = (float *)malloc(sizeof(float) * NX * NY);
    for (int h = 0; h < NX; h++){
        for (int w = 0; w < NY; w++)
            initial_array[NY * h + w] = 0;
    }

    initial_array[NY*3 + 0] = 1;

    printArray(initial_array);

    float2 *transformed_array= (float2 *)malloc(sizeof(float2) * (NY/2+1) * NX);

    cufftReal *gpu_initial_array;
    cufftComplex *gpu_transformed_array;

    cudaMalloc((void **)&gpu_initial_array, NX*NY*sizeof(cufftReal));
    cudaMalloc((void **)&gpu_transformed_array, (NY/2+1)*NX*sizeof(cufftComplex));

    cudaMemcpy(gpu_initial_array, initial_array, NX*NY*sizeof(float), cudaMemcpyHostToDevice);

    cufftHandle plan;
    cufftPlan1d(&plan, NY, CUFFT_R2C, NX);

    //                       ***** cufftPlanMany *****
    //int n[2] = {NX, NY};
    //cufftPlanMany(&plan,1,n,NULL,1,0,NULL,1,0,CUFFT_R2C,NX);

    cufftExecR2C(plan, gpu_initial_array, gpu_transformed_array);

    cudaMemcpy(transformed_array, gpu_transformed_array, NX*(NY/2+1)*sizeof(cufftComplex), cudaMemcpyDeviceToHost);

    printComplexArray(transformed_array);

/************************************************************ R2C ************************************************************/

    cufftDestroy(plan);
    free(initial_array);
    free(transformed_array);
    cudaFree(gpu_initial_array);
    cudaFree(gpu_transformed_array);

    std::system("pause");
    return 0;
}

void printArray(float *my_array){
    for (int h = 0; h < NX; h++){
        for (int w = 0; w < NY; w++)
            std::cout << my_array[NY * h + w] << " | ";
        std::cout << std::endl; 
    }
    std::cout << std::endl;     
}

void printComplexArray(float2 *my_array){
    for (int h = 0; h < NX; h++){
        for (int w = 0; w < NY; w++)
            std::cout << my_array[NY * h + w].x << " + " << my_array[NY * h + w].y << " | ";
        std::cout << std::endl;
    }
    std::cout << std::endl; 
}
4

1 回答 1

3

看来您的问题在于您打印结果的方式。您不能使用相同的例程来打印 和 的两种CUFFT_R2C情况CUFFT_C2C。在前一种情况下，你有一个(NY/2+1)*NX大小的输出，而后一种情况你有一个NY*NX大小的输出。下面的固定代码应该可以工作。

此外，添加适当的CUDA 错误检查和CUFFT 错误检查也很好，我也将其添加到下面的代码中。

#include <stdio.h>
#include <stdlib.h>
#include <cufft.h>
#include <cuda_runtime.h>
#include <assert.h>

#include <iostream>

#define NX 6
#define NY 5

void printArray(float *my_array);
void printComplexSymmetricArray(float2 *my_array);

/********************/
/* CUDA ERROR CHECK */
/********************/
#define gpuErrchk(ans) { gpuAssert((ans), __FILE__, __LINE__); }
inline void gpuAssert(cudaError_t code, char *file, int line, bool abort=true)
{
    if (code != cudaSuccess) 
    {
        fprintf(stderr,"GPUassert: %s %s %d\n", cudaGetErrorString(code), file, line);
        if (abort) exit(code);
    }
}

/*********************/
/* CUFFT ERROR CHECK */
/*********************/
static const char *_cudaGetErrorEnum(cufftResult error)
{
    switch (error)
    {
        case CUFFT_SUCCESS:
            return "CUFFT_SUCCESS";

        case CUFFT_INVALID_PLAN:
            return "CUFFT_INVALID_PLAN";

        case CUFFT_ALLOC_FAILED:
            return "CUFFT_ALLOC_FAILED";

        case CUFFT_INVALID_TYPE:
            return "CUFFT_INVALID_TYPE";

        case CUFFT_INVALID_VALUE:
            return "CUFFT_INVALID_VALUE";

        case CUFFT_INTERNAL_ERROR:
            return "CUFFT_INTERNAL_ERROR";

        case CUFFT_EXEC_FAILED:
            return "CUFFT_EXEC_FAILED";

        case CUFFT_SETUP_FAILED:
            return "CUFFT_SETUP_FAILED";

        case CUFFT_INVALID_SIZE:
            return "CUFFT_INVALID_SIZE";

        case CUFFT_UNALIGNED_DATA:
            return "CUFFT_UNALIGNED_DATA";
    }

    return "<unknown>";
}

#define cufftSafeCall(err)      __cufftSafeCall(err, __FILE__, __LINE__)

inline void __cufftSafeCall(cufftResult err, const char *file, const int line)
{
    if( CUFFT_SUCCESS != err) {
        fprintf(stderr, "CUFFT error in file '%s', line %d\n %s\nerror %d: %s\nterminating!\n",__FILE__, __LINE__,err, \
                            _cudaGetErrorEnum(err)); \
        cudaDeviceReset(); assert(0); \
    }
}

/********/
/* MAIN */
/********/
int main(){

    float *initial_array = (float *)malloc(sizeof(float) * NX * NY);
    for (int h = 0; h < NX; h++){
        for (int w = 0; w < NY; w++)
            initial_array[NY * h + w] = 0;
        }

    initial_array[NY*3 + 0] = 1;

    printArray(initial_array);

    float2 *transformed_array= (float2 *)malloc(sizeof(float2) * (NY/2+1) * NX);

    cufftReal *gpu_initial_array;
    cufftComplex *gpu_transformed_array;

    gpuErrchk(cudaMalloc((void **)&gpu_initial_array, NX*NY*sizeof(cufftReal)));
    gpuErrchk(cudaMalloc((void **)&gpu_transformed_array, (NY/2+1)*NX*sizeof(cufftComplex)));

    gpuErrchk(cudaMemcpy(gpu_initial_array, initial_array, NX*NY*sizeof(float), cudaMemcpyHostToDevice));

    cufftHandle plan;
    cufftSafeCall(cufftPlan1d(&plan, NY, CUFFT_R2C, NX));

    cufftSafeCall(cufftExecR2C(plan, (cufftReal*)gpu_initial_array, (cufftComplex*)gpu_transformed_array));

    gpuErrchk(cudaMemcpy(transformed_array, gpu_transformed_array, NX*(NY/2+1)*sizeof(cufftComplex), cudaMemcpyDeviceToHost));

    printComplexSymmetricArray(transformed_array);

    cufftSafeCall(cufftDestroy(plan));
    free(initial_array);
    free(transformed_array);
    gpuErrchk(cudaFree(gpu_initial_array));
    gpuErrchk(cudaFree(gpu_transformed_array));

    std::system("pause");
    return 0;
}

/***********************/
/* PRINTOUT REAL ARRAY */
/***********************/
void printArray(float *my_array){
    for (int h = 0; h < NX; h++){
        for (int w = 0; w < NY; w++)
            std::cout << my_array[NY * h + w] << " | ";
            std::cout << std::endl; 
        }
    std::cout << std::endl;     
}

/************************************/
/* PRINTOUT COMPLEX SYMMETRIC ARRAY */
/************************************/
void printComplexSymmetricArray(float2 *my_array){
    for (int h = 0; h < NX; h++){
        for (int w = 0; w < NY/2+1; w++)
            std::cout << my_array[(NY/2+1) * h + w].x << " + " << my_array[(NY/2+1) * h + w].y << " | ";
            std::cout << std::endl;
    }
    std::cout << std::endl; 
}
于 2014-07-09T06:23:58.780 回答