Chapter 3 Momentum Space CESE Method
3.1 The core scheme of the momentum space CESE method
The momentum space CESE method is developed to conquer the problems of wave propagation out of the boundary and ease the boundary treatment etc.. When the wave progresses out of the boundary, no matter how boundary condition we use the information of the wave has lost.
To treat the scattering problem the momentum space CESE method shows its useful capability.
The original CESE method is second-order accuracy method. It approximates the numerical solution by second-order Taylor expansion. The higher derivative term is truncated as an error.
It shows that the convection-diffusion equation Eq. (2.13) has loss its’ capability of accuracy when the viscosity coefficient is large. However, it is natural problem of the second order scheme.
Because of the higher derivative term is truncated, the accuracy is only to first derivative term.
The higher-order CESE method is to be developed for more diffusive cases in the future. In this section, we propose the momentum space CESE method. The higher derivative term is transformed to the source term. The second order scheme is enough to keep the accuracy if the source term is accurate enough. The momentum space CESE method uses the iterating process to keep the source term accuracy. For the explicit scheme, it has to keep the stability, i.e. the△t can not too large. It is time consuming for the explicit scheme. The iterating process added in each time step. Even if the iterations is few, it is still time consuming. But it is not a big problem for us, because of the space domain is smaller than the coordinate space method. For the quantum problem we describe in chapter 6, the domain for the momentum space CESE method is only [−5, 5]. For coordinate space method, the box size is must ≈ 3300 (see Ref.[13]). The total computing time of the time evolution is comparable for the coordinate space method, even it is shorter than it. However the momentum space CESE method is now a useful scheme to treat the
large diffusive problems. But for the future the higher-order CESE method will be developed. The higher-order CESE method can be easily extended to the momentum space CESE method. The basic momentum space CESE method is proposed to see its capability to keep accuracy even the equation is more diffusive.
3.1.1 Convection equation Eq. (2.1) can be transformed into
∂eu
∂t + iapeu = 0, (3.1)
by making the Fourier transformation, the system has the coordinate and the momentum representation alternatively:
For simplicity, we just use the same symbol u instead ofeu. The Eq. (3.1) becomes ∂u/∂t+iapu = 0. This is simply an ordinary differential equation. Its analytic solution is straightforward. With initial condition u(p, t = 0), the solution at any time is u(p, t) = u(p, t = 0)e−ipat. Obviously, the amplitude of the solution u(p, t) is stationary at any time in the momentum space. Though the equation and its solution in momentum space are rather simple, they serve as the calibration example for the momentum space CESE method. The Eq. (3.1) equivalently be written as
∇ · (0, u) = −iapu, (3.2)
where the operator∇ = (∂/∂p, ∂/∂t). Consider p and t as the coordinates of a two-dimensional Euclidean space E2. The conservation laws becomes
I
Let E2be divided into non-overlapping rectangular regions referred to as conservation elements CEs. The CEs with the mesh point (j, n) ∈ Ω are denoted by CE−(j, n) and CE+(j, n),
respectively. Let SE(j, n) be the rhombus shaped area DEGF depicted in Figure 3.1. The function value at the the center of CE±can be approximated by un−
1 2
j±12, (up)n−
1 2
j±12 and (ut)n−
1 2
j±12.
Definitions of CE and SE in E2.
For any (p, t) ∈ SE(j, n), u(p, t) and −→
h (p, t) are approximated by u(p, t; j, n) and
−
→h (p, t; j, n), respectively. We define
u(p, t; j, n) = unj + (up)nj(p− pj) + (ut)nj(t− tn), (3.4) where (pj, tn) is the coordinate of the mesh point (j, n) .
Note that unj, (up)nj and (ut)nj are constants in SE(j, n). We also have
−
→h (p, t; j, n) = (0, u(x, t; j, n)). (3.5)
Requiring that u = u(p, t; j, n) satisfies Eq. (3.1) within SE(j, n), one has
n −iap n
The approximation of the total flux leaving the boundary of CE±(j, n) is
where we use Taylor expansion to estimate the function value at the center of CE± and we also designate up = △p4 up. In this section, the value of (ut)nj is easy to solve. Using the transform in Eq. (3.6), we can easily obtain the value of (ut)nj. So the iterating process is not necessary. Here we just show that the iterating process also works in linear cases. For the following nonlinear cases the value of (ut)nj has convolution integral with unknown value. It is not easy to solve like the simple case here. In order to prepare for nonlinear situation (in chapter 4). With the aid of Eqs. (3.6) and (3.8), unj and (up)nj can be solved in terms of (u)n−
approximating the source term to the Taylor expansion of the center of¯¯BD¯¯and¯¯CD¯¯. Here the index ℓ is the number of time that Eq. (3.9) has been iterated, and unj solved by Eq. (3.8) can be denoted by unj,0. Using the Cauchy criterion, we define the convergence as
¯¯unj,ℓ− unj,ℓ−1¯¯< ϵ. (3.10)
We stop the iterations if the convergence criterion is reached for a plausible small ϵ. The criterion is usually arrived within iterations less than ten times.
3.1.2 Convection-diffusion equation
In the momentum space, Eq. (2.13) can be transformed into
∂eu
∂t + (iap + µp2)eu = 0, (3.11)
by making the Fourier transformation. For simplicity, we just use the same symbol u instead ofeu. Eq. (3.11) within SE(j, n), one has
(ut)nj =−(iapj+ µp2j)unj. (3.13) The approximation of the total flux leaving the boundary of CE±(j, n) is
F±(j, n) =
where the conditions are similar to the previous simple wave problem. With the aid of Eqs. (3.15) and (3.13), unj and (up)nj can be solved iteratively in terms of (u)n− time-marching scheme is derived similar to the previous simple wave case. The iterating process can be added in the source term, too. But for the equation, it is not necessary to use the iterating process. So we do not go detail to describe it. For more details, we will show in the nonlinear chapter.
3.2 Boundary treatment
The boundary treatment in the momentum space CESE method is quite easy. Because the momentum is directly related to the kinetic energy, extremely large energy for a system is usually unphysical. So, only a moderate momentum region will be sufficient for numerical modeling.
Also, the wave will simply vanish at the numerical boundary and cause no trouble like the methods in coordinate space. The non-reflecting boundary condition for momentum space CESE method is obtained in a simple way. That is when the numerical solution at the boundary, we add
the ghost cell. As depicted in Figure 2.5, we let un−
1 2
jb++12 = 0; un−
1 2
pjb++12 = 0, (3.16)
un−
1 2
jb−−12 = 0; un−
1 2
pjb−−12 = 0.
, then there are CE+(with the ghost cell) and CE−for unj
b+ and unpj
b+ (unj
b−and unpj
b−).