An explicit method can suffer from numerical instability in case when Courant number is sufficiently large. The Courant number is inversely proportional to square of grid spacing. This case shows the occurrence of numerical instability that starts at spots where the grid spacing is exceptionally small. Finally, the code gives "floating point overflow" numerical error.
Equation to solve
i
∂Ψ
∂t
+ 10∆Ψ = 0.
Domain
Cylinder. The radius of the cylinder is R=1 and the length L=4.
Mesh
A mesh for a cylinder constructed with use of blockMesh utility.
The cylinder is divided in four prism-like blocks and a central block, which is a rectangular parallelepiped with a 16×16×32 grading. Four side blocks are curvilinear hexahedrons with a similar grading.
Initial condition
The function at zero time is Gaussian distribution,
Ψ(x,y,z) = exp(−16 [(x−0.4)^{2}−y^{2}−z^{2}]).
Blue color corresponds to Ψ = 0, and red to Ψ = 1. The cylinder is sliced in order to show the interior.
Boundary conditions
Neumann (zero gradient normal to the boundary).
Calculations
The time step is ∆t=1.8·10^{−5}. Calculation time is around 1 minute.
Results
Numerical divergence caused by large Courant number. The sub-figures show calculated absolute value of the wave function. Numerical instability starts to develop at t=0.002 in a region where the grid spacing is smallest. Then, it spreads to an entire region.