Numerical stability of 2nd order Runge-Kutta integration algorithms for use in particle-in-cell codes

Abstract

An essential ingredient of particle-in-cell (PIC) codes is a numerically accurate and stable integration scheme for the particle equations of motion. Such a scheme is the well known time-centered leapfrog (LF) method [1] accurate to 2 nd order with respect to the timestep Δt. However, this scheme can only be used for forces independent of velocity unless a simple enough implicit implementation is possible. The LF scheme is therefore inapplicable in Monte-Carlo treatments of particle collisions [2] and/or interactions with radio-frequency fields [3]. We examine here the suitability of the 2 nd order Runge-Kutta (RK) method. We find that the basic RK scheme is numerically unstable, but that conditional stability can be attained by an implementation which preserves phase space area. Examples are presented to illustrate the performance of the RK schemes. We compare analytic and computed electron orbits in a traveling nonlinear wave and also show selfconsistent PIC simulations describing plasma flow in the vicinity of a lower hybrid antenna.PACS : 52.65.Rr, 52.65.Pp Numerical stability of 2 nd order Runge-Kutta integration schemesIn (PIC) simulations with many particles and evolving on long time scales it is essential to use a sufficiently simple, accurate and stable integration scheme for the electron and ion equations of motion dv e,i /dt ≡v e,i = qE z /m + a e,i ; dz/dt ≡ż = v ,where a e,i are the external acting forces per unit mass, q/m is the particle charge to mass ratio and E z is the self-consistent electric field determined from the Poisson equation. Accuracy and stability of an integration scheme do not go hand-in-hand and need to be discussed separately [4]. Accuracy has to do with the order and magnitude of the lowest order truncation terms, while stability refers to the evolution of numerical perturbations. Here we concentrate on stability. Accuracy is evident from the scheme lowest-order truncation term.As discussed in detail in [5], E z is only a function of time and particle position. If, in addition, the external forces a e,i acting on the particles are also independent C100