Albeit slower, with operations per time-step, than the numerical schemes in we have seen so far, FFTs can be used to compute the Fourier representation of the advection-diffusion problem (1.3.2#eq.2) and describe the evolution of each Fourier component separately
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[N-m] = 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 K-space, subsequent transformations back to X-space are in fact only required for plotting.
The example below illustrates the advection of a box calculated with the same time step as previously.
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 sections20.
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