# Fourier pseudo-spectral method and numerical dissipation

Performing a direct numerical simulation of isotropic turbulence with Fourier pseudo-spectral method (Orzag & Patterson, PRL, 1972) using FFT. For a background of the method, which is widely used in turbulence community, you can see this course: http://www.math.ualberta.ca/~bowman/m655/lab3d.pdf

Using the so called $2/3$ rules for dealiasing which consists of doing $$\hat{u}(\mathbf{k},t)=0~~~~~~~~~~~~~~if~~~k > \frac{2}{3} k_{MAX}$$ where $t$ is the time, $\mathbf{k}$ is the wave number, $k_{MAX}$ is the maximum wave number and $\hat{u}$ is the spectral amplitude of the velocity.

Does the dealiasing act as a numerical dissipation ? In other words, is there a leak of energy due to the dealiasing ?

-
It's not clear what your question is asking. Could you clarify? –  Geoff Oxberry Oct 15 '12 at 19:51
what about now ? –  aberration Oct 15 '12 at 20:06
It's much better, thanks! –  Geoff Oxberry Oct 15 '12 at 20:07

The dealiasing of the convolution doesn't act as numerical dissipation. In fact, energy is conserved only if you kick out the aliased terms.

The idea behind dealiasing FFT-based convolutions is to get rid of extra terms that are added by the FFT. A convolution is just a sum, and you can compute it by just calculating the sum. However, this is really slow, so it's better to Fourier-transform the inputs and multiply the result and then invert the Fourier transform, which, by the convolution theorem, is the same thing as a convolution.

But the convolution theorem only works when the inputs are of infinite length; for finite length inputs, extra terms (aliases) show up which aren't physical. The point of dealiasing is to get back the original equation that you were trying to calculate, while still allowing you to use FFTs to speed up the computation.

-
Malcolm didn't mention the FFTW++ library that he co-authored, but it helps you perform dealiased convolutions quickly by handling the padding implicitly (among other things) . See fftwpp.sourceforge.net. –  Matthew Emmett Oct 16 '12 at 0:20
Oh yeah, thanks Matt! The convolution routines in FFTW++ really reduce the cost of computing dealiased convolutions. There is a reference on the sourceforge page and more on my web page. I worked on these with John Bowman, whose 3D turbulence lab is linked in the OP. –  Malcolm Oct 16 '12 at 6:49

Understanding what the dealiasing does requires understanding why you're taking FFTs in the first place. For me, the method of weighted residuals approach provides the simplest framework:

1. You're expanding a solution as a sum of Fourier modes and dropping it into the Navier--Stokes equations. Your test functions are also a Fourier basis.
2. When you multiply the trial and test functions together and take the inner product (i.e. integrate with an appropriate conjugation), you'll see that you have one integral you cannot simply evaluate because it is nonlinear.
3. You can approximate that integral using a uniform quadrature scheme where you simply sum the value at uniform collocation points.
4. Notice your quadrature is exact for a quadratic nonlinearity if you use 1.5 times as many quadrature points as you had Fourier modes. This is the three-halves factor.
5. Finally, notice the quadrature can be efficiently computed as a FFT.

All the FFT magic can be thought of as an efficient, approximate quadrature scheme which can be made exact for quadratic nonlinearities. Other dealiasing variants play other tricks to make the quadrature exact for the same class of nonlinearities.

So, to answer your question, dealiasing adds no dissipation. This can be seen because every step in the procedure was exact.

-

@Malcolm gave an expert explanation of the dealiasing of the FFT based convolution sums. My intuition tells me however that this is still not the complete answer @abberation needed, so I'll give few more sentences that reflect my own struggle for understanding of numerical methods.

Dealiasing doesn't act as a numerical dissipation, but it does the same thing as numerical dissipation - it prevents your simulation to blow-up.

There is an unphysical transfer of energy from high to low wavenumbers which causes, first the instability (often followed by a simptome - 2h waves, h-grid node distance, the smallest waves that can be represented on a grid).

The reason behind this is aliasing of high wave numbers to those low wavenumbers that can be represented by the grid. This causes "spectral blocking" - the unphysical increase of the amplitudes (read energy) of the high wavenumbers, which gets worse as you approach the upper limit of the represented wavenumbers.

First remedy (Philips, 1959): Filter out, or make equal to zero, upper half of the represented wavenumber spectrum.

Second remedy (Orszag, 1971): The remedy known as the Orszag's Two-Thirds rule. Steven said "filtering out half of the wavenumbers is a waste, we should keep 2/3N and filter out only the upper 1/3N."

The conclusion.

Our simulations sometimes tend to blow-up. If you live in second order, finite volume universe (like the author of these lines did), and know nothing about Spectral Methods, they tell you "Oh your simulation is blowing up? Try using first order upwind instead of central scheme, or blend your central scheme with 10% of upwind in deffered correction approach. That will add dissipation but will prevent your simulation from blowing-up! Oh you're doing LES? Then you shouldn't use upwind, unless you are performing ILES, use just central scheme and Smagorinsky will model the physical dissipation that happens on lowest scales." If you live in the Spectral Methods universe, you learn how unphysical transfer of energy takes place, and you find the remedy in dealiasing by Two-Thirds Rule or "padding".

My inspiration for this post is a wonderfull book by J. Boyd "Chebyshev and Fourier Spectral Methods", section 11.