Source-linked AI summary
Exactly Energy Conserving Semi-Implicit Particle in Cell Formulation
Giovanni Lapenta
TL;DR
Further improvements to PIC continue to be studied, with energy conservation a prominent goal. ECSIM modifies the semi-implicit particle mover and field equations to retain an explicit-like computational cycle without inner iteration. ECSIM conserves energy to round-off and remains stable with grid spacing up to 16 orders of magnitude larger than the Debye length in the reported tests.
Problem
Further improvements to PIC continue to be studied, with energy conservation a prominent goal.
Source: §1 · p. 1Method
ECSIM modifies the semi-implicit particle mover and field equations to retain an explicit-like computational cycle without inner iteration.
Source: §1 · p. 1, §6 · p. 10Results
ECSIM conserves energy to round-off and remains stable with grid spacing up to 16 orders of magnitude larger than the Debye length in the reported tests.
Source: §1 · p. 2, §6 · p. 10Takeaways & Limitations
The reported results indicate ECSIM can model macroscopically large systems with coarse resolution while conserving energy to round-off.
Source: §6 · p. 10ECSIM does not conserve momentum exactly; the reported error is about one part per thousand and does not accumulate in time.
Source: §4.1 · p. 7
Abstract
from arXiv · showhide
We report a new particle in cell (PIC) method based on the semi-implicit approach. The novelty of the new method is that unlike any of its semi-implicit predecessors at the same time retains the explicit computational cycle and conserves energy exactly. Recent research has presented fully implicit methods where energy conservation is obtained as part of a non linear iteration procedure. The new method (referred to as Energy Conserving Semi-Implicit Method, ECSIM), instead, does not require any non-linear iteration and its computational cycle is similar to that of explicit PIC. The properties of the new method are: i) it conserves energy exactly to round-off for any time step or grid spacing; ii) it is unconditionally stable in time, freeing the user from the need to resolve the electron plasma frequency and allowing the user to select any desired time step; iii) it eliminates the constraint of the finite grid instability, allowing the user to select any desired resolution without being forced to resolve the Debye length; iv) the particle mover has a computational complexity identical to that of the explicit PIC, only the field solver has an increased computational cost. The new ECSIM is tested in a number of benchmarks where accuracy and computational performance are tested.
1. Introduction
The paper asks whether PIC can conserve energy exactly while retaining explicit PIC’s step order and avoiding nonlinear iteration. It introduces ECSIM, a semi-implicit method, and reports exact energy conservation alongside stability and resolution benefits.
1. Introduction
- Explicit PIC can exhibit growing energy error, while fully implicit energy-conserving schemes require coupled particle-field solves with nonlinear iteration.Explicit PIC can achieve acceptable energy error with sufficiently small time steps.Source: §1 · p. 1
- ECSIM retains explicit PIC’s particle-then-field cycle without nonlinear iteration while conserving energy exactly.The authors present it as a new semi-implicit implementation.Source: §1 · p. 1
- Previous semi-implicit methods offered stability and multiple-time-scale handling but did not conserve energy; ECSIM changes the mover and field-equation derivation.Its mover needs no inner iteration and has complexity comparable to the explicit PIC mover, while its field solver costs more CPU effort.Source: §1 · p. 2
- ECSIM conserves energy to round-off with a direct solver, remains stable with grid spacing up to 16 orders of magnitude above the Debye length, and avoids anisotropic heating and cooling.The introduction also reports operation counts of the same order as previous semi-implicit PIC.Source: §1 · p. 2
- The paper derives ECSIM, proves its energy conservation, tests benchmarks, and examines finite-grid instability and anisotropic heating and cooling.Appendices cover stability analysis, test discretization, implementation accuracy, and a computational-cycle summary.Source: §1 · p. 2
2. Energy Conserving Moment Implicit Method
ECSIM modifies the semi-implicit mover and field formulation to retain an explicit-like computational cycle without nonlinear iteration while enabling exact energy conservation. Its linear current–field relation and mass matrices support the energy-conserving formulation.
2. Energy Conserving Moment Implicit Method
- The mover combines DIM D1 position advancement with θ-scheme velocity advancement, using a known midpoint position instead of an iterated predictor-corrector position.The two midpoint positions are described as conceptually the same; the proposed one is computed explicitly.Source: §2 · p. 2
- A rigorously derived current, without approximation or linearization, is linked linearly to the advanced electric field through symmetric, sparse mass matrices.For full three-component vectors, there are nine such matrices, each with 2d diagonals.Source: §2 · pp. 3–4
- The mass-matrix relation preserves symmetry between field interpolation to particles and current interpolation to the grid, making energy exchange identical in either calculation.Source: §2 · p. 4
- The linear particle response in Maxwell’s equations allows stepping over wave timescales and removes stability requirements.The passage notes this stability advantage is also present in previously known semi-implicit methods.Source: §2 · p. 4
- ECSIM’s mover and field formulation enable exact energy conservation without nonlinear iteration while retaining linear particle–field coupling.The approach preserves the semi-implicit computational cycle.Source: §2 · pp. 2, 4
- In a typical iPic3D run, the particle-mover step costs one-third as much, while evaluating mass matrices adds CPU time and memory use.The comparison includes particle motion and field interpolation to particles; the passage describes the extra matrix cost as manageable.Source: §2 · p. 4
3. Formal Proof of Energy Conservation
The energy proof equates particle and field energy exchange, under θ = 1/2 and a grid curl discretization that preserves the continuum curl identity. The resulting energy conservation is exact to round-off.
3. Formal Proof of Energy Conservation
- Energy conservation holds only for θ = 1/2, when particle and field energy exchange are evaluated consistently.Source: §3 · pp. 4–5
- Matching the particle and field J·E exchange terms yields energy conservation to round-off, confirmed by the reported results.Source: §3 · p. 5
- The proof requires current based on mean particle speeds and a grid curl operator that preserves the continuum identity; the examples use such a discretization.Source: §3 · p. 5
4. Benchmarks
The benchmarks test ECSIM on electrostatic and electromagnetic instabilities against explicit PIC, fully implicit ECPIC, and IMM implementations. They report energy conservation and agreement with linear theory, alongside momentum behavior and computational cost.
4. Benchmarks
- The MATLAB benchmarks use electrostatic 1D1V and electromagnetic 1D3V codes, follow electrons against fixed ions, and test two-stream beams in a periodic domain of length 2π.The electromagnetic case uses three velocity, current, and field components.Source: §4 · p. 5, §4.1 · p. 5
4.1. Electrostatic Example: Two-Stream Instability
- ECSIM and ECPIC conserve electrostatic energy, with ECSIM at machine precision; all four runs agree with linear theory, while explicit PIC has less than 1% energy error and IMM fares worse.The benchmark also compares the four methods’ electric-field-energy evolution.Source: Figure 1 · p. 6, §4.1 · p. 6
- The electrostatic timing comparison is indicative only: methods were unoptimized, and preconditioning can improve ECPIC efficiency.The authors note that future work should optimize mass-matrix implementation for modern architectures.Source: §4.1 · p. 6
- Momentum is not conserved exactly: its error is about one part per thousand and does not accumulate, though momentum-priority applications may favor explicit momentum-conserving methods.The numerical collision frequency oscillates without a secular trend and remains below one part per thousand.Source: §4.1 · p. 7, Figure 3 · p. 7
4.2. Electromagnetic Example: Transverse Beam (Weibel) Instability
- ECSIM costs less than IMM in the transverse benchmark despite about 30 seconds for nine mass matrices, because it avoids a particle-corrector mover.The timing comparisons are unoptimized; preconditioning can make ECPIC more efficient.Source: §4.1 · p. 6, §4.2 · pp. 7–8
- All transverse-instability runs reproduce the linear-theory growth rate; ECSIM conserves energy perfectly, ECPIC is limited by Newton tolerance, and IMM conserves energy well but not precisely.The comparison follows magnetic-field mode m = 3.Source: §4.2 · p. 7, Figure 4 · p. 8
5. Elimination of the Finite Grid Instability
ECSIM remains stable with exactly conserved total energy and constant electron temperature across simulations spanning grid spacings up to 16 orders of magnitude above the Debye length. The section also examines energy fluctuations and artificial anisotropic heating or cooling.
5. Elimination of the Finite Grid Instability
- Across the tested cases, ECSIM remains stable with exactly conserved total energy and constant electron temperature at grid spacings up to 16 orders of magnitude above the Debye length.The simulations used a stable, thermal-equilibrium plasma; the authors report runs spanning ∆x/λDe from 10 to 10^16.Source: §5 · p. 8
- Thermal and electric-field fluctuations nearly overlap and diminish as system size increases, while momentum conservation also improves toward noiseless fluid behavior.Thermal energy reaches a minimum near the 10^-15 direct-solver precision; beyond that, it is limited by total-energy error.Source: §5 · p. 8
- Larger CFL values shift the transition toward reduced fluctuations and macroscopic fluid behavior to smaller spatial scales, with otherwise similar qualitative results.The compared ensembles used CFL values of 0.0127 and 0.127.Source: §5 · pp. 8–9, Figure 6 · p. 9, Figure 7 · p. 9
5.1. Artificial Cooling and Artificial Anisotropy
- ECSIM and fully implicit methods eliminate the artificial anisotropy tested in a magnetized plasma setup.The comparison uses a typical Earth-space-environment condition with a uniform, magnetized initial electron population and immobile ions.Source: §5.1 · p. 9
- The IMM comparison is made within its stable range, whereas explicit PIC would exceed its stability limits at ωpe∆t = 10 and ∆x/λDe = 50.Figure 8 compares total energy and electron temperature by direction for standard IMM and ECSIM.Source: §5.1 · p. 9, Figure 8 · p. 10
- The standard IMM loses over 10% of its energy, entirely in the perpendicular plane, in a deliberately pathological case not representative of prior simulations.The authors note that in a production run they would reduce the time step; they report only a few examples and say complex production runs remain to be assessed.Source: §5.1 · pp. 9–10
6. Conclusions
The conclusions present ECSIM as a semi-implicit PIC method built around a new mover and mass-matrix representation of particle response. It conserves energy exactly and is unconditionally stable, while the authors identify limits in other conservation properties and further development needs.
6. Conclusions
- ECSIM combines a new particle mover with mass matrices for plasma response, producing a self-consistent linear field problem.The mover hybridizes the D1 and θ movers, advancing position in a leap-frog-like manner and using an implicit electric-field expression for velocity.Source: §6 · p. 10
- The method is exactly energy conserving and unconditionally stable, subject to the accuracy constraint that particles should not move too far in one time step.The stated mover accuracy condition is vthe∆t < ∆x; violating it can make the solution inaccurate without disrupting exact energy conservation.Source: §6 · p. 10
- ECSIM eliminates the Debye-length grid constraint in the tested simple Maxwellian plasma, with stable, accurate behavior reported at cells up to 16 orders of magnitude larger than the Debye length.The authors call for further testing of how this property holds in more complex setups that include drifts.Source: §6 · p. 10
- ECSIM lacks momentum and charge conservation, so users are advised to monitor those quantities in realistic applications.Source: §6 · p. 10
- Compared with standard IMM, ECSIM reduces particle-mover complexity by a factor of 3, while adding mass-matrix computation and linear Maxwell solves.The authors report that reduced mover cost outweighs matrix-computation cost; a full ECPIC performance comparison remains open because their ECPIC run was unpreconditioned.Source: §6 · pp. 10–11
- Further development includes improving momentum and local electric-field-divergence conservation with substep velocity updates and implementing ECSIM in a 3D production code.Source: §6 · p. 11
Appendix A. Von Neumann Stability Analysis of the New Mover
The Von Neumann analysis compares Langmuir-wave frequencies for ECSIM, the θ-scheme, and explicit leapfrog. ECSIM and the θ-scheme remain stable for all time steps, while ECSIM is slightly more accurate than the θ-scheme at intermediate steps.
Appendix A. Von Neumann Stability Analysis of the New Mover
- The analysis assumes harmonic time dependence and derives Langmuir-wave dispersion for the new mover with θ = 1/2.Source: Appendix A · p. 11
- ECSIM and the θ-scheme have real Langmuir-wave frequencies for all time steps, unlike leapfrog, which becomes unstable when ωpe∆t > 2.The stability comparison is shown for the three movers in Figure A.9.Source: Appendix A · pp. 11–12, Figure A.9 · p. 12
- At large ωpe∆t, both implicit schemes approach Nyquist frequency ω∆t = π, while explicit leapfrog requires ωpe∆t < 2.Source: Appendix A · p. 12
- At intermediate ωpe∆t, ECSIM stays closer to the physical Langmuir-wave frequency than the θ-scheme, with little expected practical consequence.Source: Appendix A · p. 12
Appendix B. Spatial Discretization and its Properties
The appendix specifies a staggered grid for the spatial discretization and shows how periodic boundaries make the discrete Poynting-flux balance exactly zero.
Appendix B. Spatial Discretization and its Properties
- The implementation places magnetic fields at cell centers and electric fields at vertices on a uniform, periodic 1D grid.The interval [0, L] is divided into N cells of width ∆x.Source: Appendix B · p. 12
- The grid notation ∇g × E and ∇g × B is given concrete meaning through discrete curl expressions for the staggered fields.Source: Appendix B · p. 12
- With periodic boundary conditions, interior terms cancel and matching endpoint terms cancel, leaving an exactly zero Poynting-flux balance.The appendix identifies this balance as physically meaningful.Source: Appendix B · p. 12
Appendix C. Rate of Convergence in Space and Time of the ECSIM
The mover’s discretization error decreases quadratically with both time step and grid spacing. The appendix uses separate module tests and notes that exact energy conservation and second-order accuracy can help diagnose coding errors.
Appendix C. Rate of Convergence in Space and Time of the ECSIM
- The convergence analysis tests the Maxwell solver with given sources and the mover with given electric fields separately, rather than analyzing their full coupling.The authors say full coupling also depends on particle number and shape, and magnetic-field coupling is first order in time.Source: Appendix C · p. 12
- Centered spatial and temporal discretizations are expected to yield second-order convergence as ∆t and ∆x decrease.This expectation follows from the stated second-order accuracy of the Maxwell and motion equations.Source: Appendix C · p. 12
- The mover’s L2 discretization error decreases with the second power of both ∆t and ∆x.The tests vary ∆t at fixed ∆x and then ∆x at fixed ∆t, comparing against a converged solution; Figure C.10 reports the results.Source: Appendix C · p. 12, Figure C.10 · p. 13
- Exact energy conservation together with second-order accuracy provides a debugging check, since even small boundary-condition errors can break energy conservation.Source: Appendix C · p. 13
Appendix D. Summary of the ECSIM Algorithm
ECSIM retains an explicit-like computational cycle without iteration between particles and fields. Its sequence constructs particle response matrices, solves linear Maxwell equations, and then advances particles using the new field.
Appendix D. Summary of the ECSIM Algorithm
- The cycle advances positions using known velocities, interpolates the known magnetic field, computes rotation operators, and forms the hatted current.Source: Appendix D · p. 13
- It constructs 3v mass matrices representing the plasma response and solves a linear set of Maxwell equations.The matrices may be stored or computed with matrix-free methods.Source: Appendix D · p. 13
- After the field solve, the algorithm interpolates the new field to particles and advances their velocity.Source: Appendix D · p. 13
- The computational cycle retains the nature of explicit or other semi-implicit PIC, with no iteration between particles and field.Source: Appendix D · p. 13