c - 在 FFTW3 中进行二阶导数

Question

我已经使用前向 FFT 和 IFFT(对结果进行归一化)测试了我的代码的一些实际功能，这工作正常。

不过，我想对实函数求二阶导数。为了简单起见，我将 sin(2*pi*t) 作为测试用例。这是我使用的相关代码(库中的 FFT 函数):

int main(void)
{
    int i;
    int nyh = (N/2) + 1;

    double result_array[nyh][2];

    double x_k[nyh][2];
    double x_r[N];

    FILE* psit;
    psit=fopen("psitest.txt","w");  

    init();

    fft(x, result_array);  //function in a library, this has been tested
    psi(result_array, x_k);  
    ifft(x_k, x_r);        //function in a library, this has been tested

    for(i=0;i<N;i++)
    {
        fprintf(psit, "%f\n", x_r[i]);
    }


    fclose(psit);

    return 0;
}

void psi(double array[nyh][2], double out[nyh][2])
{
    int i;

    for ( i = 0; i < N/2; i++ )
    {
        out[i][0] = -4.0*pi*pi*i*i*array[i][0];
        out[i][1] = -4.0*pi*pi*i*i*array[i][1];
    }   
    out[N/2][0]=0.0;
    out[N/2][1]=0.0;
}

void init()
{
    int i;

    for(i=0;i<N;i++) 
    {
        x[i] = sin(2.0*pi*i/N);
    }
}

现在问题来了:该算法对 sin( 2*pi*t*K) 形式的任何函数都能完美运行，其中 K 是一个整数，但如果我将其作为测试function sin(3*pi*t)，算法失败。我看不到编码中的错误。

请注意，因为函数是实数，所以我只需要取k值的一半。这不是问题。

Answer 1

我的猜测是 sin(3*pi*t) 引入了不连续性，因为它不会在您的样本间隔中给出整数个周期。对于大多数与 FFT 相关的应用程序，您将应用一个窗口函数来处理此类不连续性，但显然这会在您的导数中引入一个误差项，我不确定您是否能够纠正这一点。