I've been running into trouble updating the solution in Fourier space during the diffusion sub-step. Basically for each timestep, I do the following:
- Reaction substep using the exact solution in the assignment description
- Forward FFT on the output of step 1
- Update the solution in Fourier space with u^(k,t+1) = u^(k,t)exp(-k*k*dt)
- Inverse Fourier transform (scaled by N) on the output of step 3
If I run my code without doing Step 3, I just get the solution to the reaction equation, since step 4 just cancels out step 2 in this case (i.e., I get the vector u scaled by u.size() if I apply forward and then inverse FFT without doing the update in step 3). But when I include the update in Step 3, my solution looks like nonsense for all t>0 (basically just oscillating rapidly around zero with very small amplitude).
My approach for that step was to loop over all elements in u^, corresponding to each k-value, like this:
for(int i = 0; i < uhat.size(); i++){
uhat[i] *= exp(-1.0*k[i]*k[i]*dt);
}
where k is a vector containing 2*pi*q/L for -N/2 <= q < N/2. I assumed that uhat[i] corresponds to the k-value stored at k[i] since they are both the same size, but maybe this is not correct and could be the cause of my problems.