High-order exponential solver method for particle-in-cell simulations

Abstract

Outstanding advances in solid-state laser technology, employing the optical parametric chirped-pulse-amplification (OPCPA) technique, have led physicists to focus laser pulses to highly-relativistic intensities which led to novel schemes for charged-particle acceleration and radiation generation in laser-driven plasmas. Microscopic understanding of these highly nonlinear processes is possible via accurate modeling of the laser-plasma interaction using particle-in-cell (PIC) simulations. Numerous codes are available and they rely on finite difference time domain methods on Yee-grids or on the analytical solution of the Maxwell-equations in spectral space. In this work, we present a solution bridging these two methods, which we call finite difference exponential time domain solution. This method could provide a very high accuracy even in 3D, but with improved locality, similar to the pseudospectral analytical methods without relying on transformation to special basis functions. We verified the accuracy and the convergence of the method in various benchmarks, including laser propagation in vacuum and in underdense plasma. We also simulated electron injection in a non-linear laser-plasma wakefield acceleration and surface high-harmonic generation in the overdense regime. The results are then compared with those obtained from standard PIC codes.

Figures (20)

• Figure 1: The $\omega(\mathdutchcal{k}')$ dispersion curves of the spectral advection equation Eq. (\ref{['eq:exponential_advect']}) using 4th $(a)$, 6th $(b)$ and 8th $(c)$ order Taylor expansion of the exponential with time steps $\Delta t = 0.8$ (blue and orange curves) and $\Delta t = 0.4$ (purple and red curves). For the derivative, the ideal spectral representation $\mathdutchcal{k}' = \mathdutchcal{k}$ was used. We plot two sets of curves on each figure where one set is corresponding to the real part (in blue and purple) and to imaginary part (in orange and red) of $\omega(\mathdutchcal{k})$, the latter is up shifted by the value of 0.5.
• Figure 2: The $\mathdutchcal{k}' = \sqrt{g^2(\mathdutchcal{k})}$ spectral curves Eq. (\ref{['eq:spatial_spectrum2']}) for centered Eq. (\ref{['eq:spatial_derivativeC']}) and staggered Eq. (\ref{['eq:spatial_derivativeSt']}) finite difference operators using 4th (red), 10th (blue), 20th (green) and 30th (purple) order accurate formulas. We also show the attenuation function of the binomial filter with dashed gray lines for reference. Here $\mathdutchcal{k}_{\max}=\pi/\Delta x$ (i.e. $\lambda_{\max}= 2 \Delta x$) corresponds to the largest frequency representable on the grid. The laser pulse is typically discretized below $0.2\mathdutchcal{k}_{\max}$ (which equals $\Delta x = \lambda_0/10$). We show the spectral range $> 0.5\mathdutchcal{k}_{\max}$ on the left ($b$) which could contain various numerical defects depending on the finite-difference formula. Temporally, the exponential propagation inherits the real part of $\omega(\mathdutchcal{k})$ curve from the $\mathdutchcal{k}'(\mathdutchcal{k})$ spectral functions of the differences shown here.
• Figure 3: Spectral attenuation functions of the current filters from Table \ref{['tab:spatial_filter']}: binomial and squared binomial filter (gray and brown), our lowpass filters (red, blue, green), divergence correction filter (orange). The worst case spectral functions of the 2nd and 3rd order PIC particle shapes are also shown (purple dotted lines). For reference, we also show the spectral functions of the 20th order staggering (dashed gray) and difference operator (gray).
• Figure 4: Illustration of Yee-staggering near the grid point $(x_k,y_j,z_i)$. We show subcell corresponding to staggered coordinates with dashed lines. Our $2\times$ supersampling reconstructs all field values at all 8 corners of the subcell (marked with dots of various color).
• Figure 5: Sketch of our 2D checkered parallel deposition algorithm. At the beginning of the simulation we calculate table of parallel zones which are $4\times4$ cell blocks, and we assign each zone to $2\times2$ density groups (shown here by blue, red, white and green). The particles corresponding to each cell are accessed through the particle map which provides the starting point of a linked list. In the example we show here that the cell $(i,k)$ stores the particle index $\mathdutchcal{p}_0$, and the particle at $\mathdutchcal{p}_0$ points to $\mathdutchcal{p}_1$, and the latter points to $\mathdutchcal{p}_n$, which terminates the list with an invalid particle index. Parallel density deposition can be done independently in each zone by depositing into one of the 4 components in the density buffer with matching group color and summing the density afterwards.