# The Klein-Gordon Equation
## Background

The focusing/defocusing nonlinear Klein-Gordon equation describes the evolution of a possible complex scalar field $u$ according to,
    ${\frac {\partial ^{2}u}{\partial t^{2}}}-\Delta u+u=\pm |u|^{2}u$,         ( 1)
where $+$ is the focusing case and $-$ the defocusing case in a similar manner to the nonlinear Schrödinger equation. Blow up of three dimensional radially symmetric real solutions to this equation have recently been numerically studied by Donninger and Schlag. Two dimensional simulations of the Klein-Gordon equation can be found in Yang. The linear Klein-Gordon equation occurs as a modification of the linear Schrödinger equation that is consistent with special relativity, see for example Landau or Grenier. At the present time, there have been no numerical studies of blow up of solutions to this equation without the assumption of radial symmetry. This equation has generated a large mathematical literature and is still poorly understood. Most of this mathematical literature has concentrated on analyzing the equation on an infinite three dimensional space with initial data that either decays exponentially as one tends to infinity or is nonzero on a finite set of the domain. Here, we will simulate this equation in a periodic setting. Since this equation is a wave equation, it has a finite speed of propagation of information, much as a sound wave in air takes time to move from one point to another. Consequently for short time simulations, a simulation of a solution that is only nonzero on a finite part of the domain is similar to a simulation on an infinite domain. However, over long times, the solution can spread out and interact with itself on a periodic domain, whereas on an infinite domain, the interaction over long times is significantly reduced and the solution primarily spreads out. Understanding the interactions in a periodic setting is an interesting mathematical problem. The Klein-Gordon equation has a conserved energy given by
$\int {\frac {1}{2}}\left({\frac {\partial u}{\partial t}}\right)^{2}+{\frac {u^{2}}{2}}+{\frac {1}{2}}\left|\nabla u\right|^{2}\mp {\frac {\left|u\right|^{4}}{4}}\mathrm {d} \mathbf {x}$.
The equation is also time reversible. For long time simulations, one wants to construct numerical methods that approximately conserve this energy and are also time reversible. When using Fourier spectral methods, we primarily need to ensure that the time discretization preserves these properties, since the spectral spatial discretization will typically automatically satisfy these properties. Following Donninger and Schlag, we use two schemes. First, an implicit-explicit time stepping scheme which is time reversible but only conserves the energy approximately and is given by
    ${\frac {u^{n+1}-2u^{n}+u^{n-1}}{(\delta t)^{2}}}-\Delta {\frac {u^{n+1}+2u^{n}+u^{n-1}}{4}}+{\frac {u^{n+1}+2u^{n}+u^{n-1}}{4}}=\pm \left|u^{n}\right|^{2}u^{n}$         ( 2)
and second, a fully implicit time stepping scheme with fixed point iteration
    ${\frac {u^{n+1,k+1}-2u^{n}+u^{n-1}}{(\delta t)^{2}}}-\Delta {\frac {u^{n+1,k+1}+2u^{n}+u^{n-1}}{4}}+{\frac {u^{n+1,k+1}+2u^{n}+u^{n-1}}{4}}$ $=\pm {\frac {\left|u^{n+1,k}\right|^{4}-\left|u^{n-1}\right|^{4}}{4(u^{n+1,k}-u^{n-1})}}$         ( 3)
which conserves a discrete energy exactly
    $\int {\frac {1}{2}}\left({\frac {u^{n+1}-u^{n}}{\delta t}}\right)^{2}+{\frac {1}{2}}\left({\frac {u^{n+1}+u^{n}}{2}}\right)^{2}+{\frac {1}{2}}\left|\nabla {\frac {u^{n+1}+u^{n}}{2}}\right|^{2}\mp {\frac {\left|{u}^{n+1}\right|^{4}+\left|{u}^{n}\right|^{4}}{8}}$.         ( 4)
As before, the superscript $n$ denotes the time step and $k$ denotes the iterate in the fixed point iteration scheme. Iterations are stopped once the difference between two successive iterates falls below a certain tolerance.

### Matlab and Python Programs

Listings A , B , C and D demonstrate Matlab implementations of these time stepping schemes. In one dimension, the Klein-Gordon equation has easily computable exact solutions, (see for example Nakanishi and Schlag) which can be used to test the accuracy of the numerical schemes. These equations seem to display three possibilities for the behavior of solutions which are dependent on the initial conditions:
* the solutions could disperse or thermalize, that is a given localized initial condition spreads out over the entire space
* the solutions blow up or become infinite
* a portion of the solution travels around as a localized particle while the rest of the solution disperses.
Since the equations are reversible, there is also the possibility that a solution which is initially distributed over the spatial domain localizes itself.

### A Two-Dimensional OpenMP Fortran Program

The programs that we have developed in Fortran have become rather long. Here we add subroutines to make the programs shorter and easier to maintain. Listing E is the main Fortran program which uses OpenMP to solve the 2D Klein-Gordon equation. Notice that by using subroutines, we have made the main program significantly shorter and easier to read. It is still not as simple to read as the Matlab program, but is significantly better than some of the previous Fortran programs. It is also much easier to maintain, and once the subroutines have been written and debugged, they may be reused in other programs. The only drawback in using too many subroutines is that one may encounter a slight decrease in performance due to the overhead of calling a subroutine and passing data to it. The subroutines are in listings F , G , H , I , J , K and an example makefile is in listing L . Finally listing M contains a Matlab program which produces pictures from the binary files that have been computed. One can then use another program to take the images and create a video.

### A Three-Dimensional MPI Fortran Program using 2DECOMP&FFT

For this section, we will use the library 2DECOMP&FFT available from http://www.2decomp.org/index.html. The website includes some examples which indicate how this library should be used, in particular the sample code at http://www.2decomp.org/case_study1.html is a very helpful indication of how one converts a code that uses FFTW to one that uses MPI and the aforementioned library.
We now give a program for the three-dimensional nonlinear Klein-Gordon equation. The program uses the same subroutine structure as the two-dimensional code. To make the program easy to reuse, the subroutine listed in listing U has been created to read an INPUTFILE which specifies the parameters to use for the program and so the program does not need to be recompiled every time it is run. To enable the program to scale better, the arrays which hold the Fourier frequencies and grid points have also been decomposed so that only the portions of the arrays used on each processor are created and stored on the processor. The program is listed in listing W , to use this program simply compile it using MPI wrappers and gfortran (other compilers such as Intel fortran, PGI fortran or IBM XL fortran should also work), and then run it in the directory from which the INPUTFILE and data are stored. An example submission script is in listing Y, this requires an environment variable named "inputfile" to be set with the location of the INPUTFILE that specifies problem parameters, such as grid size and number of timesteps.
A further addition is a short postprocessing program to create header files to allow one to use the bov (brick of values) format that allows one to use the parallel visualization software VisIt. The program VisIt can be downloaded from https://wci.llnl.gov/codes/visit/home.html. This program also run on laptops, desktops as well as parallel computer clusters. Documentation on using VisIt is available here https://wci.llnl.gov/codes/visit/manuals.html and here http://www.visitusers.org/index.php?title=Main_Page. A short video tutorial on how to use VisIt remotely is available at http://cac.engin.umich.edu/resources/software/visit.html. An alternative program to create images is Paraview which can be downloaded from http://www.paraview.org/. To use Paraview, xdmf format header files can be used to get Paraview to read in the binary data output files. A short video tutorial on rendering with Paraview can be found at http://www.youtube.com/watch?v=rqCXFr73sfk&list=PL4FFYwZhtT0BlDBtjsSkQ4LcM1v9K24No. A Fortran program to create these header files is in listing X.

### Exercises

* Compare the accuracy of the implicit and semi-implicit time stepping schemes in eqs. 2 and 3 . Which scheme produces the most accurate results in the least amount of real time?
* Write serial Fortran programs to solve the two- and three-dimensional Klein-Gordon equations using the fully implicit time stepping scheme in eq. 3 .
* Write OpenMP parallel Fortran programs to solve the two- and three-dimensional Klein-Gordon equations using the fully implicit time stepping scheme in eq. 3 .
* The MPI command MPI_BCAST is used in the subroutine readinputfile, listed in list N . Look up this command (possibly using one of the references listed in the introduction to programming section) and explain what it does.
* Write an MPI parallel Fortran program to solve the two- and three-dimensional Klein-Gordon equations using the fully implicit time stepping scheme in eq. 3 .
* Compare the results of fully three-dimensional simulations with periodic boundary conditions ($\mathbb {T} ^{3}$) with analytical predictions for blow up on the entire real space ($\mathbb {R} ^{3}$) summarized in Donninger and Schlag.
* Grenier explains that the linear Klein-Gordon equation can be written as two coupled Schrödinger equations. One can extend this formulation to the nonlinear Klein-Gordon equation. If we let
    $u=\phi +\xi \quad {\text{and}}\quad {\frac {\partial u}{\partial t}}=\phi -\xi$         ( 5)
then the two coupled equations
    $i{\frac {\partial }{\partial t}}{\begin{bmatrix}\phi \\\xi \end{bmatrix}}={\begin{bmatrix}-\Delta -1&-\Delta \\\Delta &\Delta +1\end{bmatrix}}{\begin{bmatrix}\phi \\\xi \end{bmatrix}}\pm {\begin{bmatrix}1\\-1\end{bmatrix}}{\frac {|\phi +\xi |^{2}(\phi +\xi )}{2}}$         ( 6)
are equivalent to the nonlinear Klein-Gordon equation
    ${\frac {\partial ^{2}u}{\partial t^{2}}}-\Delta u+u=\pm u^{3}$.         ( 7)
** Fill in the details to explain why eqs. 5 and 6 are equivalent to eq. 7 . In particular show that by adding and subtracting the two equations in eqs. 5 and 7 , we get $i{\frac {\partial }{\partial t}}\left(\phi +\xi \right)=-\left(\phi -\xi \right)$ and $i{\frac {\partial }{\partial t}}\left(\phi -\xi \right)=-\Delta \left(\phi +\xi \right)-\left(\phi +\xi \right)\pm \left|\phi +\xi \right|^{2}\left(\phi +\xi \right)$. Differentiating the first of these equations and substituting it into the second, then recalling that we defined $u=\phi +\xi$ in eq. 5 gives us the Klein-Gordon equation in eq. 7 .
** Solve these two equations using either the implicit midpoint rule or the Crank Nicolson method.

## Notes

* ↑ An incomplete but easily accessible mathematical introduction to this equation can be found at http://wiki.math.toronto.edu/DispersiveWiki/index.php/Semilinear_NLW.
* ↑ Donninger and Schlag (2011)
* ↑ Yang (2006)
* ↑ Landau (1996)
* ↑ Grenier (1994)
* ↑ Donninger and Schlag (2011)
* ↑ Nakanishi and Schlag (2011)
* ↑ At the present time, Matlab’s video commands cannot reliably produce a single video from a very long simulation, so it is better to use Matlab to create still images.
* ↑ Donninger and Schlag (2011)
* ↑ Grenier (1994)

## References

Donninger, R.; Schlag, W. (2011). "Numerical study of the blowup/global existence dichotomy for the focusing cubic nonlinear Klein-Gordon equation". Nonlinearity. 24: 2547–2562. doi:10.1088/0951-7715/24/9/009.
Grenier, W. (1994). Relativistic Quantum Mechanics. Springer.
Landau, R.H. (1996). Quantum Mechanics II. Wiley.
Schlag, W. (2011). Invariant Manifolds and Dispersive Hamiltonian Evolution Equations. European Mathematical Society.
Yang, L. (2006), Numerical studies of the Klein-Gordon-Schrodinger equations, National University of Singapore {{citation}}: Cite has empty unknown parameter: |month= (help)