After hours of browsing the DSP posts and resources online, I still struggle to understand why my code diverges when I activate zero-padding dealiasing. When I deactivate it, everything 'works well', and I pinned down the issue to the points I raise below.

I'm doing a loop where I

  1. retrieve from the solution of a linear system some Fourier coefficients in frequency space - say $\hat{u}_k$ and $\hat{v}_k$ - for each $k \in [-K,K]$
  2. apply backward DFT to go to physical space.
  3. form product $uv$
  4. apply forward DFT to go back to frequency space and broadcast the de-aliased product to update the r.h.s. of the linear system for the next iteration

The Orszag "3/2" zero-padding rule ensures dealiasing with such quadratic linearities $uv$. I extended accordingly my DFT input/output arrays at steps 2. and 4 (from $N$ to $N_\text{pad}$) and made sure the zeros are inserted at the right location.

  • Is full dealiasing achieved at step 3? In other words, should I only apply zero-padding during my backward DFT at step 2? Should my forward DFT call be of size $N$ only?

  • It seems that there (p18) and there (p4) they apply a factor 3/2 = $N_\text{pad}/N$. Can someone explain mathematically or in physical terms the reason for such a scaling? I grasp that the 'numerical' FFT is of size $N_\text{pad}$ whereas you're interested in a convolution of a size $N$ only but I'm still puzzled. In this open-source library there's no such scaling*.

*I'm aware FFTW returns an unnormalized DFT in contrast to MATLAB or numpy where ifft is normalized.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.