SYLLABUS Previous: 4.1 FFT with the
Up: 4 FOURIER TRANSFORM
Next: 4.3 Aliasing, filters and
4.2 Linear equations.
Slide : [
Scheme 
Code 
Run 
VIDEO
login]
Albeit slower, with
operations per timestep,
than the numerical schemes in
we have seen so far,
FFTs can be used to compute the Fourier representation of the
advectiondiffusion problem (1.3.2#eq.2) and describe the
evolution of each Fourier component
separately

(1) 
The equation can be integrated analytically without any restriction on
the timestep; starting directly from the initial condition, this yields
for each Fourier component

(2) 
After the initialization, where the initial condition is
discretized
by
sampling on a regular mesh and stored for subsequent transformation,
the scheme
implemented in JBONE
reads:
int N = mesh_.size(); //A power of 2
double boxLen = mesh_.size()*mesh_.interval(0); //Periodicity
double k = 2*Math.PI/boxLen; //Notation
Complex ik1 = new Complex( 0., k );
Complex ik2 = new Complex(k*k, 0. );
Complex advection = new Complex(ik1.scale(velocity)); //Variables
Complex diffusion = new Complex(ik2.scale(diffusCo));
Complex[] s0 = new Complex[f.length]; //FFT real to KSpace
s0=keepFFT.getFromKSpace(FFT.firstPart,boxLen); // only once
s[0] = s0[0];
for (int m=1; m<N/2+1; m++) { //Propagate directly
total= advection.scale((double)(m )); // from initial
total=total.add(diffusion.scale((double)(m*m))); // condition
exp=(total.scale(timeStep*(double)(jbone.step))).exp();
s[m ] = s0[m].mul(exp); // s0 contains the IC
s[Nm] = s[m].conj(); // f is real
}
FFT ffts = new FFT(s,FFT.inKSpace); //Initialize in Kspace
f = ffts.getFromXSpacePart(FFT.firstPart,boxLen); //FFT real back for plot
If you read carefully, you probably wonder now about the sign of the phase
factor, which is opposite from (4.2#eq.2); the reason is that
the phase chosen here for the spatial harmonics is exactly opposite to the
one that has been assumed in the FFT routine from
Numerical Recipes
[24].
This makes the scheme look like as if it evolves backwards in time  it
doesn't! Also note that once the initial condition has been transformed
to Kspace, subsequent transformations back to Xspace are in fact only
required for plotting.
The example below
illustrates the advection of a box calculated with the same time step
as previously.
JBONE applet: press Start/Stop
to simulate the advection of a box function, showing what appears to
be a perfect solution that alternates with a strong Gibbs phenomenon.
In fact, the oscillations are always present in the solution, but
disappear from the plot when straight lines are used to connect
adjacent mesh points.

The Gibbs phenomenon (an artifact from using
harmonic functions for the discretization of sharp edges discussed in
sect.1.4.5) is clearly visible except when
is an integer; but don't be fooled by the plot
routine, which (for simplicity, but still incorrectly) draws straight
lines between the grid points when it should in fact draw the complete
oscillations that are visible in figure 1.4.5#fig.1.
The power of the method is evident when taking large timesteps:
edit the parameter to TimeStep=64 , add some diffusion
Diffusion=0.3 and compare the solution obtained with a single step
with the result computed using your favorite FD or FEM scheme from the
previous sections^{20}.
This is very nice indeed, but remember that dealing with an inhomogeneous
medium
,
or complicated boundary conditions is problematic
in Fourier space.
SYLLABUS Previous: 4.1 FFT with the
Up: 4 FOURIER TRANSFORM
Next: 4.3 Aliasing, filters and
