Petsc FFTW

#include <math.h>
#include <petscvec.h>
#include <petscmat.h>

int main(int argc,char **args)
{
        int out = 0;
        int dimensions[] = {2000};
        int ndims = sizeof(dimensions)/sizeof(int);
        int N = 1;
        int i;
        PetscScalar v;
        Vec inputvec, backtransform;
        Mat fftwa;
        Vec fftwin,fftwout, fftwbackout;

    PetscInitialize(&argc,&args,(char *)0,"");

    for (i = 0; i < ndims; i++)
        N *= dimensions[i];

    VecCreateSeq(PETSC_COMM_SELF,N,&inputvec);
    VecCreateSeq(PETSC_COMM_SELF,N,&backtransform);
    MatCreateFFT(PETSC_COMM_WORLD,ndims,dimensions,MATFFTW,&fftwa);
    MatGetVecsFFTW(fftwa,&fftwin,&fftwout,&fftwbackout);

    for (i = 0; i < N; i++)
    { // fill the input vector with a sinewave
        v = sin(double(i)/60.0);
        VecSetValues(inputvec,1,&i,&v, INSERT_VALUES);
    }
    VecAssemblyBegin(inputvec);
    VecAssemblyEnd(inputvec);

    VecScatterPetscToFFTW(fftwa,inputvec,fftwin);
    MatMult(fftwa,fftwin,fftwout);
    MatMultTranspose(fftwa,fftwout,fftwbackout);
    VecScatterFFTWToPetsc(fftwa,fftwbackout,backtransform);

    VecView(inputvec, PETSC_VIEWER_DRAW_WORLD); PetscSleep(-1);
    VecView(fftwout, PETSC_VIEWER_DRAW_WORLD); PetscSleep(-1);
    VecView(backtransform, PETSC_VIEWER_DRAW_WORLD); PetscSleep(-1);

    VecDestroy(&fftwin); VecDestroy(&fftwout); VecDestroy(&fftwbackout);
    VecDestroy(&inputvec); VecDestroy(&backtransform); MatDestroy(&fftwa);

    PetscFinalize();
    return out;
}