2
$\begingroup$

While trying to solve a stochastic Gross-Pitaevskii equation I have found a problem that can be tracked down to something buggy occurring in the simplest Schrodinger equation possible:

$$\partial_t \psi = i \partial^2_x \psi$$

This has the very simple solution of a plane wave:

$$\psi(t) = A e^{i k x - i\omega t}, \text{with } \omega = k^2$$

and my problem is that in some cases I cannot recover this very simple solution numerically. This is due to the fact that my code uses a very standard Fast Fourier Transform method to go to k-space, where the evolution is trivial, and then does the inverse FFT to go back to real space. This method naturally enforces periodic boundary conditions in your simulation grid since you are expanding your solution as a sum of periodic functions. What occurs is that, when the initial condition is a plane wave with a periodicity that does not match the domain size, the algorithm is not able to provide the correct solution.

I show here an example. The problem is clear in $|\psi(t)|$, which should remain constant:

enter image description here

The effect is less dramatic in the phase but still present. One can see some ripples that come from the bottom edges.

enter image description here

I'm not surprised by the fact that, if I enforce periodic boundary conditions, computing the evolution of a wave which is not periodic in the simulation domain yields problems. However, the method is very standard and used extensively, with packages like XMDS (xmds.org) employing it by default. Therefore, I am surprised by the fact that I did not find any mention of this method failing to solve such an extremely simple example.

My question is, is there something that I'm missing? Should I just get over it and assume that this is not a good method if I expect a solution of this kind? Does someone know a reference where this is documented?

Note: Cross-posted in https://math.stackexchange.com/questions/1612508/solving-a-simple-schrodinger-equation-with-fast-fourier-transforms, and also in https://scicomp.stackexchange.com/questions/21821/solving-a-simple-schroedinger-equation-with-fast-fourier-transforms.

$\endgroup$
  • 1
    $\begingroup$ if you take periodic boundary conditions in the $x$-direction, the wave number $k$ can take only discrete values (in steps of $2\pi/L_x$); you just need to make sure that your initial condition satisfies this, there is not much more to say. $\endgroup$ – Carlo Beenakker Jan 14 '16 at 20:41
  • 1
    $\begingroup$ @CarloBeenakker That goes on the lines of what I thought. My concern then is the following: What if I use a non-problematic initial condition (like, k = 0) but add a term of white noise to make the equation stochastic? This noise could in principle excite all possible k's. Could this noise excite a wave like the one I used as initial condition and propagate pathologically as in this simple example? To put it short, does this result mean that I shouldn't use FFT to solve problems with delocalized waves? $\endgroup$ – Carlos_San Jan 15 '16 at 13:47
  • 1
    $\begingroup$ just make sure your white noise has only components at the discrete wave numbers $k_n = 2\pi n/L_x$; as you make $L_x$ longer and longer, this discretization will become less and less significant. $\endgroup$ – Carlo Beenakker Jan 15 '16 at 13:49
2
$\begingroup$

If I understand correctly, what you're doing amounts to:

  1. Starting with initial data $\psi(x,0)$ that is represented as a truncated Fourier series.
  2. Computing the exact (up to roundoff errors) time evolution of this initial condition on a periodic domain, using discrete Fourier transforms.

This is often referred to as a Fourier spectral method. Note that no matter what initial data you wish impose, by representing it using a finite, equispaced set of point values and then solving with Fourier transforms, you are effectively representing it as a truncated Fourier series. Therefore, "when the initial condition is a plane wave with a periodicity that does not match the domain size", by sampling that initial condition on your grid and imposing periodicity, you'll obtain an initial condition with a range of frequencies in it. Your PDE is dispersive, so you'll see those different frequencies separate, as shown in your first plot.

You claim that "$|\psi(t)|$ should remain constant"; I think you mean to refer to $|\psi(x,t)|$. But with a dispersive equation and an initial condition that is not monochromatic, this claim is not justified.

By the way, I believe that scicomp.SE is a more useful place to post questions like this.

$\endgroup$

Your Answer

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

Not the answer you're looking for? Browse other questions tagged or ask your own question.