2: /*
  3:     Provides an interface to the FFT packages.
  4: */

  6: #include <../src/mat/impls/fft/fft.h>   /*I "petscmat.h" I*/

 10: PetscErrorCode MatDestroy_FFT(Mat A)
 11: {
 13:   Mat_FFT        *fft = (Mat_FFT*)A->data;

 16:   if (fft->matdestroy){
 17:     (fft->matdestroy)(A);
 18:   }
 19:   PetscFree(fft->dim);
 20:   PetscFree(A->data);
 21:   PetscObjectChangeTypeName((PetscObject)A,0);
 22:   return(0);
 23: }

 27: /*@
 28:       MatCreateFFT - Creates a matrix object that provides FFT via an external package

 30:    Collective on MPI_Comm

 32:    Input Parameter:
 33: +   comm - MPI communicator
 34: .   ndim - the ndim-dimensional transform
 35: .   dim - array of size ndim, dim[i] contains the vector length in the i-dimension
 36: -   type - package type, e.g., FFTW or FFTCU

 38:    Output Parameter:
 39: .   A  - the matrix

 41:   Options Database Keys:
 42: + -mat_fft_type - set FFT type

 44:    Level: intermediate
 45:    
 46: @*/
 47: PetscErrorCode MatCreateFFT(MPI_Comm comm,PetscInt ndim,const PetscInt dim[],const MatType mattype,Mat* A)
 48: {
 50:   PetscMPIInt    size;
 51:   Mat            FFT;
 52:   PetscInt       N,i;
 53:   Mat_FFT        *fft;

 56:   if (ndim < 1) SETERRQ1(comm,PETSC_ERR_USER,"ndim %d must be > 0",ndim);
 57:   MPI_Comm_size(comm, &size);

 59:   MatCreate(comm,&FFT);
 60:   PetscNewLog(FFT,Mat_FFT,&fft);
 61:   FFT->data = (void*)fft;
 62:   N = 1;
 63:   for (i=0; i<ndim; i++){
 64:     if (dim[i] < 1) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"dim[%d]=%d must be > 0",i,dim[i]);
 65:     N *= dim[i];
 66:   }

 68:   PetscMalloc(ndim*sizeof(PetscInt),&fft->dim);
 69:   PetscMemcpy(fft->dim,dim,ndim*sizeof(PetscInt));

 71:   fft->ndim = ndim;
 72:   fft->n    = PETSC_DECIDE;
 73:   fft->N    = N;
 74:   fft->data = PETSC_NULL;
 75: 
 76:   MatSetType(FFT,mattype);
 77:   FFT->ops->destroy = MatDestroy_FFT;
 78: 
 79:   /* get runtime options */
 80:   PetscOptionsBegin(((PetscObject)FFT)->comm,((PetscObject)FFT)->prefix,"FFT Options","Mat");
 81:   PetscOptionsEnd();

 83:   *A = FFT;
 84:   return(0);
 85: }