Lattice Boltzmann method for warm fluid simulations of plasma wakefield acceleration

Simeoni, Daniele; Parise, Gianmarco; Guglietta, Fabio; Rossi, Andrea Renato; Rosenzweig, James; Cianchi, Alessandro; Sbragaglia, Mauro

doi:10.1063/5.0175910

A comprehensive characterization of lattice Boltzmann (LB) schemes to perform warm fluid numerical simulations of particle wakefield acceleration (PWFA) processes is discussed in this paper. The LB schemes we develop hinge on the moment matching procedure, allowing the fluid description of a warm relativistic plasma wake generated by a driver pulse propagating in a neutral plasma. We focus on fluid models equations resulting from two popular closure assumptions of the relativistic kinetic equations, i.e., the local equilibrium and the warm plasma closure assumptions. The developed LB schemes can, thus, be used to disclose insights on the quantitative differences between the two closure approaches in the dynamics of PWFA processes. Comparisons between the proposed schemes and available analytical results are extensively addressed.

I. INTRODUCTION

The process of particle acceleration plays a role of primary importance at the interface between fundamental and applied physics.¹ The growing costs and dimensions of conventional large-scale accelerators demand for new and more efficient technologies to push the energies reached by particle beams beyond the state of the art of modern day capabilities. In this context, plasma acceleration is a promising technique that would enable the construction of compact particle accelerators while retaining the same (or superior) energy gains obtained with conventional methods.^2–4 A ionized gas is perturbed via the injection of relativistic charged particles (particle wakefield acceleration, PWFA)⁵ or intense lasers (laser wakefield acceleration, LWFA),⁶ generally named driver: the interaction between the neutral plasma and the injected driver creates a wave like dynamics of positive and negative charges, and hence, strong accelerating fields (up to 100 GV/m) are developed; the interested reader might look into^7,8 to go into detail on the topic. The described process involves a large number of “actors”: the injected driver components, whether particles moving near the speed of light—typically electrons—or laser fields, and both the ions and electrons that make up the plasma. All of them interact with each other via electromagnetic forces, thus making it really difficult to predict and control the final behavior of a plasma acceleration experiment. Theoretical modeling and numerical simulations are, therefore, a powerful tool to help guide the design of new experiments.

On the side of numerical simulations, the most commonly used techniques in the field are represented by particle in cell (PIC) methods^9–13 that employ single particle dynamics to describe both the particles in the driver (in this paper, we will focus on PWFA) and the plasma components. These schemes are deeply fine-grained and, while they can capture the most refined phenomena that happen at the microscopic scale, they bring in low-statistics numerical noise.^9,14 An alternative numerical modeling can be proposed by recurring to continuum fluid descriptions of the system. Fluid models describe the plasma via macroscopic fields such as particle number density and fluid velocity, and they do so by solving the inviscid relativistic Euler equations that can be systematically derived from kinetic theory of charged gases via a coarse graining of the relativistic Maxwell–Vlasov system^15–17 (in their most general warm formulation, these equations are not closed and, hence, require additional constrains—more on this later on). Despite losing the ability to describe some kinetically pertinent features, numerical methods relying on the fluid description are still able to capture non trivial features of the PWFA system and are set up by construction not to show statistical noise. Some examples of fluid solvers used in the realm of PWFA are Architect,¹⁸ QFLUID,¹⁹ MARPLE,²⁰ FLASH,²¹ and the code used for hydrodynamic optically field-ionized plasma channels in Ref. 22.

From the theoretical point of view, fluid models have traditionally been developed by neglecting thermal effects, i.e., by neglecting pressure terms in the Euler equations. This choice has been motivated first and foremost by the fact that initial electron thermal energy in the plasma is expected to be of the order of $k_{b} T_{i} \sim 10$ eV,^23,24 which is a small value when compared with the electron rest energy $m_{e} c^{2} = 0.511$ MeV (the initial thermal energy normalized to the electron rest energy, $μ_{i} = k_{b} T_{i} / (m_{e} c^{2})$ ⁠, is the control parameter that is usually used to assess the importance of thermal effects), and second by the fact that these cold fluid models are easier to treat theoretically, providing even analytical results in some simplified cases.^25–27

Nevertheless, there is a series of important reasons that drive the development of warm fluid models for PWFA. First, the wave-like solutions to cold fluid equations become singular in the presence of wide and highly charged driver pulses (Wave Breaking^28,29) The presence of thermal effects may be one of the regularizing mechanisms that mitigate the singularity,^30–33 by altering the wakefield properties and allowing electron trapping in the wakefield.³⁴ A significant heating is also expected in the post-wavebreaking dynamic.³⁵ Second, studies targeted at the late stage dynamics of the process^36–40 point to the importance of electrons temperature (with particular emphasis on thermally driven ion-acoustic motion^39–43) for the restoration of the equilibrium conditions after a driver pulse has passed in the plasma channel. Understanding the restoration conditions is pivotal to enable the possibility of having high repetition rates of driver pulses in order to create sustained accelerating fields. Third, although the aforementioned initial temperatures would not lead to meaningful divergences in behavior from the cold case (at least in the first wave periods^34,35), a different situation is expected as many consecutive pulses are injected into the system. The energy deposited by every pulse would be partially transferred to the plasma,^36,38 leading to significant increases in temperature [some estimates provide O(1) keV increases in post wave-breaking situations^39,40]. Finally, we mention that temperature effects might also be relevant for the study of positron acceleration in quasi-hollow warm plasma channels.^41,44–47

If one wants to develop a warm fluid theory, often the only viable option when trying to derive analytical results, an immediate problem has to be tackled: additional fields are now present in the set of equations (namely, the pressure tensor fields), and therefore, a suitable closure has to be carefully selected. The closure problem has been studied in the community, and various models have been proposed^48–53 (see Ref. 54 for an outlook).

In this paper, we explore two closures. The first one⁴⁸ is based on the assumption that the underlying probability distribution solving for the Vlasov equation is a Maxwellian equilibrium, which is described in the relativistic framework by a Maxwell–Jüttner distribution.⁵⁵ Hereafter, we will refer to it as the local equilibrium closure (LEC). This choice leads to no entropy production and, therefore, grants the adoption of an isentropic equation of state to close for the pressure tensor field that in this framework can be described as a single scalar quantity. The LEC is in principle not well suited for the description of early stage dynamics, as the restoring mechanism that drives the plasma back to an eventual initial equilibrium state (particle collisions) happens on longer time scales than the ones which are typical of the early phases of PWFA. Nevertheless, this closure would still retain a physical relevance when studying late stage dynamics (when particle collisions start to become relevant⁴⁰). Furthermore, side-by-side comparisons against other closure schemes (or PIC solvers) might, indeed, reveal that the LEC is also helpful for qualitative assessments on early dynamics.

The second closure (hereafter named warm closure—WARMC) is based on the idea proposed in Refs. 49, 50, and 52 and later on reconsidered by Refs. 32, 56, and 57: the centralized moments equations obtained from the coarse-graining of the Vlasov equation are closed by neglecting the third order centralized moment, choice motivated by the assumption of weakly warm systems, and hence small momentum distribution variances. This leads to a closed set of equations that can be solved without having to make hypotheses on the underlying momentum probability distribution, except for the smallness of its variance. This gives also the possibility to evolve independently the various components of the pressure tensor, and hence to evaluate the expected momentum spread anisotropies: in fact, as the dynamics of the system is strongly focused in the acceleration direction, and no collisions are taken into account on these short timescales to regularize the process, momentum spread anisotropies are to be expected.^58,59

The main target of this paper is to develop novel numerical schemes for the simulation of warm fluid models in the context of PWFA, for both the LEC and the WARMC. The schemes rely on the lattice Boltzmann (LB) method^60,61 to solve the fluid equations. LB is a popular numerical scheme commonly used in computational fluid dynamics as an alternate scheme to direct hydrodynamical solvers. Its formulation is rooted in the kinetic theory of gases, and this provides a strong physical basis to the method. Furthermore, the space locality of the calculations involved in the method makes this solver prone to multi CPUs and/or multi GPUs parallelization.⁶⁰ In the past years, LB has been generalized to work in many fields other than classical fluid dynamics,⁶¹ and it is now widely accepted as a numerical solver for many physics problems that rely on a set of continuum equations. The LB formulation used in this paper is the so-called moment matching LB that solves systems described by advection diffusion equations.^60,61 This is the very same formulation used in Ref. 62 to solve for the cold fluid models in PWFA. Here, we extend such formulation to warm fluid models; hence, we are pushing further the usage of the LB method in the context of PWFA. The LB schemes for warm fluids are coupled to a finite difference time domain (FDTD) scheme that solves for the electromagnetic fields:^63,64 we refer to Ref. 62 for a more in depth explanation of the coupling. The development of the LB schemes for two different warm fluids closures will enable us to look for thermal effects that depend on the adopted closure scheme and, furthermore, to assess the importance of thermal spreads anisotropies showing side-by-side comparisons between the two closures.

This paper is organized as follows: in Sec. II, we provide a basic introduction to the LB method, with particular emphasis on the moment matching procedure; in Sec. III, basic equations for collisionless relativistic plasmas are reviewed, and the two fluid closures LEC and WARMC are discussed in Secs. IV and V, respectively; numerical results will be presented in Sec. VI, and an outlook and a discussion on future perspectives are given in Sec. VII.

II. LATTICE BOLTZMANN (LB) METHOD

In this section, we introduce the lattice Boltzmann (LB) method, which we use to solve the fluid equations both in the LEC and WARMC models. We first present the basics of the method, that are directly drawn from the kinetic theory of gases and in their original formulation tasked to the reproduction of classical (non relativistic) Navier–Stokes equations. We then illustrate how the whole procedure can be adapted to the simulation of generic advection equations (moment matching LB). We will see in the following sections how the relativistic warm fluid equations, in both the LEC (Sec. IV) and WARMC (Sec. V), can be recast into a set of advection equations and, hence, can be numerically solved via moment matching LB. The interested reader might look into^60,61 for more detailed discussions on both LB in its original formulation and its moment matching variant.

A. Basic LB

As the theoretical cornerstone of LB is kinetic theory, the natural starting point in its algorithmic derivation is the Boltzmann transport equation

(\partial_{t} + ξ \cdot \nabla_{x}) f (x, ξ, t) = Σ (x, ξ, t),

(1)

which describes the evolution of the phase space density

f (x, ξ, t)

referring to the number of fluid particles with velocity

ξ

at position x at time t. On the RHS of the equation,

Σ (x, ξ, t)

represents a production term that is usually expressed as the sum of two components

Σ (x, ξ, t) = S (x, ξ, t) + Ω (x, ξ, t),

(2)

S (x, ξ, t)

representing the action of volume body forces F on the particles, and

Ω (x, ξ, t)

the action of binary collisions

S (x, ξ, t) = - F \cdot \nabla_{ξ} f (x, ξ, t),

(3)

Ω (x, ξ, t) = - \frac{1}{τ} (f (x, ξ, t) - f^{e q} (x, ξ, t)),

(4)

where the latter has been here written by recurring to the customary Bhatnagar–Groos–Krook operator,⁶⁵ which expresses the tendency of the distribution function f to relax toward an equilibrium f ^eq in a typical relaxation time τ. The moments of the distribution function f that solves Eq.

M^{(n)} (x, t) = \int f (x, ξ, t) ξ^{n} d ξ,

(5)

are then used to retrieve the hydrodynamic fields. In this basic formulation, the zeroth- and first-order moments (i.e., n = 0, 1) provide the mass density and the momentum density of the fluid, respectively.⁶⁰

It can then be verified through the Chapman–Enskog expansion^60,61,66 that the moments obtained according to Eq. (5) verify the target continuum equations—again, the Navier–Stokes equations in this original formulation of the LB method—provided that f is sufficiently close to f ^eq. The most pivotal step in the development of the LB method is the realization that one needs only a truncated version of the distribution functions f and f ^eq to properly recover the desired moments $M^{(n)} (x, t)$ that solve the field equations.^67,68 For this reason, it proves expedient to expand both f and f ^eq into series of orthogonal Hermite polynomials in the variable $ξ$ and then to truncate this expansion up to the point where one recovers the macroscopic observables of interest.⁶⁹

The second most important ingredient in the development of the LB is the velocity discretization procedure. Because of the truncated Hermite polynomials expansion, one usually adopts Gauss-Hermite quadrature rules to prune the continuum velocity space and select a discrete set of

ξ_{i}

(i = 0, \dots, N_{pop} - 1)

velocities (stencil), each of them having a statistical weight w_i.^70,71 In Fig. 1, we report the discrete stencil adopted in this paper. This velocity discretization procedure splits Eq. into a set of N_pop equations, one for every discrete distribution function

f_{i} = f (x, ξ_{i}, t)

⁠, also called population. The careful selection of the velocity stencil is performed with the goal of exactly preserving the coarse-graining of the moments when passing from the continuous integrals of Eq. to discrete summations

M^{(n)} (x, t) = \int f (x, ξ, t) ξ^{n} d ξ = \sum_{i = 0}^{N_{pop} - 1} ξ_{i}^{n} f_{i} .

(6)

FIG. 1.

View large Download slide

Lattice Boltzmann discretization. Through a pruning procedure, the continuum velocity space is discretized in a minimal set of N_pop velocities. This discrete set (stencil) is selected in such a way as to guarantee the hydrodynamic equations for coarse grained fields (e.g., density, momentum). On the right, we report one of the most popular LB stencils, with N_pop = 19 velocities in 3D, which we have also used in this work. In the table, we report the adimensional velocity components of the stencil, together with corresponding statistical weights. The different shades of blue group velocities with the same magnitude.

The performed velocity discretization, when joined with an explicit time marching discretization of step

Δ t

⁠, finally delivers the Lattice Boltzmann equation

f_{i} (x + ξ_{i} Δ t, t + Δ t) = f_{i} (x, t) + Δ t Σ_{i} (x, t),

(7)

where space is discretized on a regular lattice of characteristic length

Δ x = ξ_{i} Δ t

⁠. The production term

Σ_{i} (x, ξ_{i}, t)

is a discretization of the continuous counterpart appearing in Eq.

Σ_{i} (x, t) = Ω (x, ξ_{i}, t) + S_{i} (x, t),

(8)

with

S_{i} (x, t)

being discretized according to one of the many forcing schemes available in the literature⁶⁰ to properly reproduce Eq. .

The algorithmic steps of the LB scheme are now clearly outlined. A set of N_pop versions of Eq. (7), one for every f_i, is evolved on a regular spatial grid: the f_i are updated at every node with the source term Eq. (8) (Source step) and then stream to neighboring nodes with their corresponding velocity $ξ_{i}$ (Streaming step). At every iteration, the hydrodynamic fields are obtained through the discrete summation appearing in Eq. (6).

B. Moment matching procedure for forced advection diffusion equation

In this section, we will make use of the framework established in Sec. II A to explain how the moment matching LB, tasked for the solution of advection diffusion equations (ADEs), works. This variant was initially developed as an alternative, more straightforward formulation of thermal LBs^72–74 and at the same time applied to model different physics phenomena, such as diffusive chemical reactions,⁷⁵ combustion problems,^76,77 dissolution in porous media^77,79 (more use cases can be found in Refs. 60 and 61). The target equation for a moment matching LB is a forced ADE for the scalar field

A (x, t)

\partial_{t} A + \nabla_{x} \cdot (A u) = κ \nabla_{x}^{2} A + F,

(9)

where the quantity A is advected with velocity

u (x, t)

and diffused via the diffusion parameter κ. The previously established framework (Sec. II A) can then be used by considering a discrete probability distribution function

f_{i} (x, t)

whose zeroth-order discrete moment is exactly the field

A (x, t)

A (x, t) = \sum_{i = 0}^{N_{pop} - 1} f_{i} (x, t) .

(10)

The velocity stencil described in Fig. 1 can then be used to evolve via Eq. the distribution functions f_i, and the Chapman-Enskog analysis grants that the coarse graining Eq. solves for Eq. . The last remaining ingredient to complete the algorithmic discussion of the method is the form of the Hermite truncated expansion

f_{i}^{eq}

⁠. It has been shown that a suitable form that recovers the correct equations for A and minimizes the computational error is⁸⁰^,

f_{i}^{eq} = w_{i} A (1 + \frac{u \cdot ξ_{i}}{ν^{2}} + \frac{{(u \cdot ξ_{i})}^{2}}{ν^{4}} - \frac{u \cdot u}{2 ν^{2}}),

(11)

where

ν = \frac{1}{\sqrt{3}} \frac{Δ x}{Δ t}

is a reference lattice velocity. This completes the description of the main algorithmic steps of the moment matching LB. The diffusion coefficient κ is shown via Chapman–Enskog procedure to be dependent on the relaxation time τ in the following way:

κ = ν^{2} (τ - \frac{Δ t}{2}) .

(12)

While in this work we keep τ as small as possible (taking into account the known inferior limit

τ > Δ t / 2

imposed by LB bulk stability conditions⁶⁰) in order to recover the physics of the problem, we note from Eq. that by adjusting the parameter τ in the simulations one might control the diffusive effect in the ADE, obtaining, therefore, a tunable regularizing effect that might be helpful when dealing with stiff advection equations. Finally, we report a couple of additional details for our implementation of the method. First, the source term S_i is chosen to be of the form

S_{i} (x, t) = w_{i} F (x, t),

(13)

in order to correctly reproduce the forcing appearing in Eq. .^60,61 Second, open boundary conditions are imposed at the extrema of the domain via the enforcing of the variables f_i

f_{i} (x_{b}, t) = f_{i} (x_{f}, t),

(14)

where

x_{b}

and

x_{f}

are, respectively, the boundary position and the nearest bulk fluid node. Finally, we remark that the whole procedure showed so far can be adapted to a 3D axisymmetric geometry: when doing so, after switching to cylindrical coordinates, we postulate that every quantity is independent from the azimuthal variable

A = A (r, z)

⁠. Furthermore, no angular motion is requested, hence the azimuthal component of the advection velocity u is zero.

Following^81,82 one can develop an axisymmetric moment matching LB by considering a 2D Cartesian geometry (r, z) where all the differential operators—the divergences and laplacians appearing in Eq. —are carefully adapted to the cylindrical geometry.⁸³ To this extent, we remark that differently from the cited sources we will consider A to be either a scalar or the coordinate component of some hydrodynamic tensor field. Hence, particular care will have to be taken when converting the differential operators to the cylindrical geometry. Finally, when adopting the cylindrical geometry, we have to modify the conditions applied at the r = 0 boundary node. If the LB advected quantity A is symmetric with respect to radial reflections, we keep employing Eq. . Instead, if A changes sign under radial reflections (the radial momentum component has this property), we employ the following condition:

f_{i} (x_{b}, t) = f_{\bar{i}} (x_{f}, t),

(15)

with

\bar{i}

being the mirrored direction to i with respect to the radial axis.

III. RELATIVISTIC KINETIC EQUATIONS FOR COLLISIONLESS PLASMA

In this section, we review the basic equations that can be used to build hydrodynamic models of warm plasmas starting from the kinetic theory of gases, all expressed within the framework of special relativity.^15–17 In Sec. IV and Sec. V, we will see how the obtained set of equations can be closed via either a local equilibrium assumption or a warm closure and then explain how to recast them into a set of advection equations.

In the following, we will work in a flat space-time with Minkowski metric signature $η^{α β} = diag (+, -)$ ⁠. When expressing formulas in a manifestly covariant form, we will adopt Einstein's summation convention, with Greek indexes running from 0 to $D$ and Latin indexes from 1 to $D$ ⁠. When not explicitly stated, we will use $\partial_{α} \equiv \frac{\partial}{\partial x^{α}}$ ⁠.

The starting point in a relativistic fluid theory of warm plasmas is the Relativistic Vlasov equation

ρ^{α} \frac{\partial f}{\partial x^{α}} + \frac{q}{c} F^{α β} ρ_{β} \frac{\partial f}{\partial ρ^{α}} = 0,

(16)

where

x^{α} = (c t, x)

is the space-time coordinate vector,

ρ^{α} = (ρ^{0}, ρ)

is the relativistic kinetic momentum vector,

F^{α β}

is Maxwell's electromagnetic field tensor, and

f = f (x^{α}, ρ^{α})

is the single-particle probability distribution function, describing the number of particles of mass m and charge q (m = m_e and

q = - e

in case of electrons) in the Lorentz invariant momentum space volume

\frac{d ρ}{ρ^{0}}

⁠.

The fields subject of hydrodynamic theories emerge from the coarse graining of the distribution function f in relativistic momentum space

M^{α_{1} \dots α_{n}} = c \int f ρ^{α_{1}} \dots ρ^{α_{n}} \frac{d ρ}{ρ^{0}},

(17)

with the first moments being assigned the following specific names: Invariant density h, Particle flow

N^{α}

⁠, and Energy-momentum tensor

T^{α β}

⁠. Equation implies the following conservation equations for the moments:^15,17,57

0 = \partial_{α} N^{α},

(18)

0 = \partial_{α} T^{α β} - \frac{q}{c} F^{β α} N_{α},

(19)

0 = \partial_{α} M^{α β γ} - \frac{q}{c} (F^{β α} T_{α}^{γ} + F^{γ α} T_{α}^{β}) .

(20)

Equations are the constitutive equations of a relativistic fluid, which must be coupled to Maxwell equations⁸⁴

0 = \partial_{α} F^{α β} - μ_{0} c q N^{β},

(21)

0 = \partial_{α} F_{β γ} + \partial_{β} F_{γ α} + \partial_{γ} F_{α β},

(22)

to provide a complete description of a plasma. Equations express conservation of mass, energy, and momentum once a proper decomposition for

N^{α}

and

T^{α β}

is provided: as an example, the cold fluid equations can be obtained from Eqs. and by assuming the following:

N^{α} = n_{0} U^{α}, T^{α β} = n_{0} m U^{α} U^{β},

(23)

U^{β} = γ (c, u)

being the fluid four velocity and n₀ the particle number density as measured in the fluid rest frame (density transforms under Lorentz boosts as

n = γ n_{0}

⁠). γ is the Lorentz factor associated with u. At variance with the cold fluid treatment, any other fluid model involving a non-zero temperature requires additional closing equations.

In Sec. IV and Sec. V, we will show two possible closures: the first one (LEC) postulates that the underlying distribution function describing the warm plasma [and hence solving Eq. ] is a local equilibrium distribution. This grants the use of an isentropic equation of state to close for the system; the second one (WARMC) rephrases Eqs. as centered moments equations and then assumes the third order centered moment to be zero on the excuse of small thermal spreads, i.e., by assuming small values of the initial thermal energy (⁠

i

subscripts stand for the initial values) to electron rest energy

μ_{i} = \frac{k_{b} T_{i}}{m_{e} c^{2}} .

(24)

This is the control parameter that is usually used to assess the importance of thermal effects in warm fluid models and will be extensively used in Sec. VI to label our results.

IV. LOCAL EQUILIBRIUM CLOSURE (LEC)

The key assumption of the LEC hypothesis is that the distribution f solving for Eq. is the relativistic version of the Maxwell–Boltzmann distribution, the Maxwell–Jüttner distribution.^55,85 This has a number of important consequences, first and foremost the fact that under this assumption the plasma can be considered to be an ideal fluid. It can be shown¹⁵ that, in this case, the particle flow and the energy momentum tensor assume the following forms:

N^{α} = n_{0} U^{α},

(25)

T^{α β} = (P_{0} + ϵ_{0}) \frac{U^{α} U^{β}}{c^{2}} - P_{0} η^{α β},

(26)

where ϵ₀ and P₀ are, respectively, the plasma energy density and pressure, as measured in the fluid's rest frame (all quantities written with the 0 subscript have to be intended as they are measured in this frame). Consequently, the conservation equations Eqs. and can be re-expressed as the relativistic counterpart of Euler's equations

\begin{matrix} \partial_{t} n + \nabla_{x} \cdot (n u) = 0, \\ \partial_{t} (\frac{h_{e} p}{m c^{2}}) + u \cdot \nabla_{x} (\frac{h_{e} p}{m c^{2}}) = - \frac{\nabla_{x} P_{0}}{n} + q (E + u \land B), \end{matrix}

(27)

where

h_{e} = (P_{0} + ϵ_{0}) / n_{0}

is the relativistic enthalpy per particle and

p = m γ u

is the relativistic fluid momentum. Equation is not yet closed, but another important feature of the LEC assumption is that there is no entropy production in the fluid. Therefore, to close the set of equations instead of a statement of energy conservation, we consider the entropy per particle

s_{0} = k_{b} log [\frac{1}{n_{0}} {(\frac{m k_{b} T}{2 π ℏ^{2}})}^{D / 2}] + k_{b} (1 + \frac{D}{2}),

(28)

here written in the small temperature limit,¹⁵ to be constant. This statement leads to the following scaling for temperature:

T = T_{i} {(\frac{n_{0}}{n_{i}})}^{2 / D} .

(29)

Furthermore, this grants the use of the Synge equation of state⁸⁶ that we consider here for small temperatures

P_{0} = n_{0} k_{b} T = n_{0} m c^{2} μ_{i} {(\frac{n_{0}}{n_{i}})}^{2 / D},

(30)

\begin{matrix} \int_{0} = n_{0} m c^{2} + \frac{D}{2} n_{0} k_{b} T = n_{0} m c^{2} [1 + \frac{D}{2} μ_{i} {(\frac{n_{0}}{n_{i}})}^{2 / D}] . \end{matrix}

(31)

This provides sufficient information for the closure of Eq. and paves the way to a successful numerical treatment of this set via moment matching LB.

A. Local equilibrium closure lattice Boltzmann (LEC-LB)

We can now express Eq. as a series of forced advection equations

\partial_{t} A + \nabla_{x} \cdot (uA) \equiv D_{u} (A) = F,

(32)

where

A = (\begin{matrix} n \\ [1 + (1 + \frac{D}{2}) μ_{i} {(\frac{n}{γ n_{i}})}^{2 / D}] n p \end{matrix}),

(33)

F = (\begin{matrix} 0 \\ - n_{i} m c^{2} μ_{i} \nabla_{x} {(\frac{n}{γ n_{i}})}^{1 + 2 / D} + q n (E + u \land B) \end{matrix}) .

(34)

This set can then be numerically solved via moment matching LB, using the techniques explained in Sec. II B. As already stated, we can adapt the method to work in a 3D axisymmetric environment, and particular care has to be taken when translating the equations to said geometry. As the final expressions are rather bulky, we report them in full in Appendix A.

The two-way coupling of the fluid with the electromagnetic fields needs no particular explanation: the plasma hydrodynamic quantities (particle density n and the transport velocity u) are obtained from the LB evolved quantities (the various components of A) and then fed to the FDTD Maxwell solver together with the driving bunch terms. The evolved electromagnetic fields are then plugged into Eq. (34) as source terms, and the next LB iteration can be performed.

The only detail worth of discussion is the determination of the transport velocity u from the LB-advected quantities A. In fact, due to the appearance of the γ Lorentz factor in Eq. (33), one cannot easily obtain the relationship between u and the second component of the vector A (the one containing momentum $p)$ ⁠, and therefore, there is the need for a specific iterative algorithm to determine it. Initially, set $p^{(0)}$ as the zeroth order moment of the distribution function coming out of the LB iteration, divided by n: Then,

Starting from k = 0, compute $u^{(k)}$ as

$u^{(k)} = \frac{p^{(k)}}{m \sqrt{1 + | p^{(k)} / (m c) |^{2}}} .$

(35)
Compute $γ^{(k)} = γ (u^{(k)})$ as

$γ^{(k)} = γ (u^{(k)}) = \frac{1}{\sqrt{1 - {(u^{(k)} / c)}^{2}}} .$

(36)
Set $p^{(k + 1)}$ as

$p^{(k + 1)} = \frac{p^{(0)}}{1 + (1 + \frac{D}{2}) μ_{i} {(\frac{n}{γ n_{i}})}^{2 / D}} .$

(37)
Repeat until convergence is reached.

V. WARM PLASMA CLOSURE (WARMC)

The derivation given in this section closely follows:⁵⁷ to explain how to derive the Warm Closure from the conservation equations Eqs. , one has to first define the thermal momentum

w^{μ}

as

w^{μ} = (\frac{n_{0}}{h}) U^{μ} \Rightarrow N^{μ} = n_{0} U^{μ} = h w^{μ} .

(38)

The second and third order

centered

moments of the distribution function f, respectively,

θ^{μ ν}

and

Q^{μ ν λ}

⁠, are then defined as

θ^{μ ν} = c \int (ρ^{μ} - w^{μ}) (ρ^{ν} - w^{ν}) f \frac{d ρ}{ρ^{0}},

(39)

Q^{μ ν λ} = c \int (ρ^{μ} - w^{μ}) (ρ^{ν} - w^{ν}) (ρ^{λ} - w^{λ}) f \frac{d ρ}{ρ^{0}} .

(40)

Therefore, one obtains

T^{μ ν} = θ^{μ ν} + h w^{μ} w^{ν},

(41)

M^{μ ν λ} = Q^{μ ν λ} + w^{μ} θ^{ν λ} + w^{ν} θ^{μ λ} + w^{λ} θ^{μ ν} + h w^{μ} w^{ν} w^{λ},

(42)

and due to the mass-shell condition

ρ^{μ} ρ_{μ} = m^{2} c^{2}

it follows that

T_{μ}^{μ} = h m^{2} c^{2} \Rightarrow θ_{μ}^{μ} = h (m^{2} c^{2} - w^{μ} w_{μ}),

(43)

M_{ν}^{μ ν} = N^{μ} m^{2} c^{2} \Rightarrow Q_{ν}^{μ ν} = 2 w_{ν} θ^{μ ν} .

(44)

The conservation equations Eqs. become, therefore,

0 = \partial_{α} (h w^{α}),

(45)

0 = \partial_{α} θ^{α β} + \partial_{α} (h w^{α} w^{β}) - \frac{h q}{c} F^{β α} w_{α},

(46)

\begin{matrix} 0 = \partial_{α} Q^{α β γ} + \partial_{α} (h w^{α} \frac{θ^{β γ}}{h}) + θ^{γ α} \partial_{α} w^{β} + θ^{α β} \partial_{α} w^{γ} - \frac{q}{c} (F^{β α} θ_{α}^{γ} + F^{γ α} θ_{α}^{β}) . \end{matrix}

(47)

The warm closure consists in taking

Q^{α β γ} = 0

based on the assumption of a small momentum spread

\partial_{α} (n_{0} U^{α}) = 0,

(48)

\partial_{α} (n_{0} U^{α} \frac{n_{0} U^{β}}{h}) = - \partial_{α} θ^{α β} + \frac{n_{0} q}{c} F^{β α} U_{α},

(49)

\begin{matrix} \partial_{α} (n_{0} U^{α} \frac{θ^{β γ}}{h}) = - θ^{γ α} \partial_{α} (\frac{n_{0} U^{β}}{h}) - θ^{α β} \partial_{α} (\frac{n_{0} U^{γ}}{h}) \\ + \frac{q}{c} (F^{β α} θ_{α}^{γ} + F^{γ α} θ_{α}^{β}), \end{matrix}

(50)

and this also provides an additional condition, obtained from Eq.

θ_{μ}^{μ} = h c^{2} [m^{2} - {(\frac{n_{0}}{h})}^{2}],

(51)

U_{ν} θ^{μ ν} = 0 .

(52)

Note that the previous expressions have been written using Eq. to express everything in terms of the more familiar quantities n₀ and

U^{α}

⁠. Equations constitute a closed set of fluid equations, which may be solved for the plasma evolution through moment matching LB.

A. Warm plasma closure lattice Boltzmann (WARMC-LB)

Equations can be rephrased as advection equations by realizing that, for any given quantity

G = {1, n_{0} U^{β} / h, θ^{β γ} / h}

⁠, one can write

\partial_{α} (n_{0} U^{α} G) = \partial_{t} (n G) + \nabla_{x} \cdot (n G u) = D_{u} (n G) .

(53)

Therefore, one has

\partial_{t} A + \nabla_{x} \cdot (uA) \equiv D_{u} (A) = F,

(54)

where

A = (\begin{matrix} n \\ n \frac{n_{0} U^{β}}{h} \\ n \frac{θ^{β γ}}{h} \end{matrix}),

(55)

F = (\begin{matrix} 0 \\ - \partial_{α} θ^{α β} + \frac{n_{0} q}{c} F^{β α} U_{α} \\ - θ^{γ α} \partial_{α} (\frac{n_{0} U^{β}}{h}) - θ^{α β} \partial_{α} (\frac{n_{0} U^{γ}}{h}) + + \frac{q}{c} (F^{β α} θ_{α}^{γ} + F^{γ α} θ_{α}^{β}) \end{matrix}) .

(56)

The system of equations can now be simulated via moment matching LB, interpreting all terms on the RHS of equations as external forcings. The Fluid-Maxwell coupling discussed in Sec. IV A applies also, in this case, with the exception of the iterative scheme for the transport velocity u since, in this framework, u can be naturally identified from the LB advected quantities.

Also in this case, we employ the axisymmetric description: the passages that are needed to adapt Eqs. (51), (52), (55), and (56) to this framework are rather lengthy but simple and we report the final full expressions in Appendix A.

As a final remark, it can be worth mentioning how the tensor

θ^{α β}

can be decomposed in its transversal and longitudinal pressure components (w.r.t. the direction of the driving bunch, here chosen to be the

z -

axis), respectively,

P_{⊥}

and

P_{∥} = P_{⊥} + Δ P

⁠, as to define a term of comparison against the isotropic pressure case provided by the LEC. We postulate^17,87 that in the fluid rest frame the anisotropic energy momentum tensor and the thermal momentum would read as

T_{0}^{μ ν} = diag (ϵ, P_{⊥}, P_{⊥}, P_{∥}),

(57)

w_{0}^{μ} = \frac{n_{0}}{h} (c, 0, 0, 0) .

(58)

It follows then from Eq. that

θ_{0}^{μ ν} = diag (ϵ_{0} - \frac{n_{0}^{2} c^{2}}{h}, P_{⊥}, P_{⊥}, P_{⊥}) + Δ P δ_{z}^{μ} δ_{z}^{ν} .

(59)

It is sufficient to apply a Lorentz boost

Λ_{ν}^{μ}

to derive this tensor in a generic reference frame

\begin{matrix} θ^{μ ν} = Λ_{α}^{μ} Λ_{β}^{ν} θ_{0}^{α β} = (P_{⊥} + ϵ_{0} - \frac{n_{0}^{2} c^{2}}{ρ}) \frac{U^{μ} U^{ν}}{c^{2}} - P_{⊥} η^{μ ν} + Δ P Λ_{z}^{μ} Λ_{z}^{ν} . \end{matrix}

(60)

Other than showing how to obtain the rest frame values for both pressure terms, this form of the tensor is actually useful to determine that some of its components are null in the axisymmetric framework. This effectively reduces the number of advection equations that have to be solved via moment matching LB.

VI. NUMERICAL RESULTS

In this section, we present the current capabilities of the method by reproducing known analytical results and we also show side-by-side comparisons between the two closures, highlighting the emergence of momentum spreads anisotropies. We divide the showcase of our results in three different parts: in the first, we start by considering temperature effects in a completely 1D scenario, where our equations of motions are made mono-dimensional by imposing translational symmetry along the radial directions (all derivatives w.r.t. transversal coordinates are, therefore, zero), imposing $D = 1$ in Eqs. (33) and (34) (this is equivalent to considering a 1D kinetic momentum space) and considering $(1 + 1) -$ dimensional tensors in the WARMC case. In this 1D1V set-up, we compare our numerical result with known analytical solutions.^30–32,56 Then, we move to the discussion of dispersion relations in a 1D3V setup:^57,88 there is still translational invariance along the transversal directions, but a 3D kinetic momentum space is considered [ $D = 3$ in Eqs. (33) and (34) and $θ^{μ ν}$ tensors are (3 + 1) dimensional]. Finally, we consider full spatially resolved simulations in a 3D axisymmetric environment (3D3V).

In this work, plasma ions are considered as an immobile background with constant plasma density

n_{i} = 10^{16}

{cm}^{- 3}

⁠, as they are several order of magnitude more massive than plasma electrons. Furthermore, the driving electron pulse is modeled by a rigid bi-Gaussian density n_b, with rms-sizes

σ_{z} = σ_{r} = 25

μ m

and moving at the speed of light c (moving from right to left)

n_{b} (r, z) = α exp (- \frac{{(z - c t - z_{0})}^{2}}{σ_{z}^{2}}) exp (- \frac{r^{2}}{σ_{r}^{2}}) .

(61)

The driving bunch is initially centered in z₀ and with a peak amplitude α such that the Normalized Charge Parameter

\tilde{Q}

\tilde{Q} = {(2 π)}^{D / 2} σ_{z} σ_{r}^{D - 1} α (\frac{k_{p}^{D}}{n_{i}}),

(62)

representing the strength of the perturbation, is kept at a fixed desired value. Here, k_p is the plasma wave number,

k_{p} = ω_{p} / c

⁠, with

ω_{p} = \sqrt{e^{2} n_{i} / (m ϵ_{0})}

the cold plasma frequency (m being the electron mass and ϵ₀ the vacuum permittivity). Unless explicitly stated, all the results presented from now on will consider a computational domain of size

L_{r} = 6 / k_{p}

and

L_{z} = 30 / k_{p}

⁠, with cell resolutions

Δ z = Δ r = 0.01 / k_{p}

and computational time step

Δ t = 0.001 / ω_{p}

⁠, for a simulated physical time of

30 / ω_{p}

⁠. The value of τ is tuned in every setup in order to gauge between the LB numerical stability conditions⁶⁰ and the smallness of the parameter κ introduced in Eq. , as to properly recover the physics of the problem: we have chosen a value of

τ = 0.52 Δ t

in the 1D setups, and

τ = 0.53 Δ t

in the 3D setups, and both these two values perfectly reproduce the available analytic solutions.

As already mentioned, just like other common PWFA solvers,¹⁸ the electromagnetic fields are solved via an FDTD scheme^63,64 through numerical integration of the curl Maxwell equations [Faraday's and Ampère's laws in Eqs. (21) and (22)]. It can be shown in fact⁸⁹ that when these equations are considered together with the continuity equation, the divergence Maxwell counterparts [Gauss laws in Eqs. (21) and (22)] are automatically satisfied, provided that the fields are correctly initialized. For this reason, we properly initialize the electromagnetic fields by solving analytically the Maxwell system with just a rigid Gaussian density (the driving bunch) in its rest-frame and then boost the quantities to the lab-frame.^90,91

A. 1D1V results

We adapt our scheme to work in a 1D environment by employing a value of σ_r which is way bigger than the computational domain along the radial direction L_r so that the driving bunch of Eq. (61) reduces effectively to a single Gaussian profile along the z coordinate, and the system assumes a translational invariance along the r coordinate. Equation (62) is, therefore, employed by imposing $D = 1$ ⁠. Furthermore, the 1V condition is reached for the LEC by setting $D = 1$ in Eqs. (33) and (34), while in the WARMC case is sufficient to initialize to zero the transversal components of the $θ^{μ ν}$ tensor.

We present in Fig. 2 the results for a numerical benchmark of the method obtained by comparing both the two fluid solvers against the analytic solutions that can be obtained by considering Eq. (27) and Eqs. (48)–(52) in a 1D1V environment. The corresponding equations are well studied in the literature: for the LEC, see Sec. VI in Ref. 31 or equivalently;³⁰ for the WARMC, see Ref. 32 and replace the laser driven case with a particle driven case. The theoretical solutions are obtained via finite differences integration of Eq. (6.4) in Ref. 31 (LEC) and Eq. (8) in Ref. 32 (WARMC). Although coming from different set of equations, these solutions are equivalent in the limit of small temperature. Being this the case, they are represented by a single solid black curve in Fig. 2. We observe that both the two methods are perfectly able to reproduce the theoretical results. In this simplified 1D scenario, it is also possible to start appreciating some effects of temperature on the dynamics (the choice of μ_i is done in order to evidence thermal effects in the first periods of the wave): by comparing against the analytic solutions of the cold fluid model (light purple in Fig. 2), it is possible to see that temperature leads to a decrease in the thermal wavelength of the plasma wave. As it will be seen in Sec. VI B, this effect is in stark contrast to what can be observed in a 3V environment, where the wavelength increases or decreases depending on the selected fluid closure.

FIG. 2.

View large Download slide

Comparison between the numerical results obtained from the moment matching LB with the LEC (LEC-LB), the moment matching LB with the WARMC (WARMC-LB) and their respective warm (⁠ $μ_{i} = 0.04 \to k_{b} T_{i} = 20$ $keV$ ⁠) 1D theory^31,32 (which coincide for this choice of the parameters and hence is represented here by a single solid black line). We also show for reference the analytic solution that can be derived from the cold fluid model (light purple). From top to bottom, we show results for the electron plasma density n, the electric field E (normalized w.r.t. $E_{p} = \frac{m c ω_{p}}{e}$ ⁠), and the plasma velocity u. All the curves are plotted with the co-moving variable $ζ = z - c t$ on the $x -$ axis. We show two significant regimes: $\tilde{Q} = 0.01$ on the left column and $\tilde{Q} = 0.5$ on the right column. The Gaussian driving bunch (green line) appearing in the top panels is vertically shifted by a unit factor for visualization purposes.

B. 1D3V results

An increase in complexity with respect to the 1D1V case is represented by the 1D3V case, that is targeted specifically at the recovery of dispersion relations for the plasma wave, while still retaining a 1D simulation framework. In this context, it is actually possible to extract for both closures a dispersion relation linking the frequency ω of the plasma wave to its wavenumber k. In the literature, such result has been obtained with various levels of approximations.^88,92 In our case, a theoretical study can be conducted by considering Eqs. and , linearly perturbing them and then obtaining the dispersion relation via Fourier analysis.^57,93,94 The following formulas are obtained at order

O (μ_{i})

⁠:

{(\frac{ω}{ω_{p}})}^{2} = 1 - \frac{5}{2} μ_{i} + \frac{5}{3} μ_{i} {(\frac{k}{k_{p}})}^{2} LEC,

(63)

{(\frac{ω}{ω_{p}})}^{2} = 1 - \frac{5}{2} μ_{i} + 3 μ_{i} {(\frac{k}{k_{p}})}^{2} WARMC .

(64)

Note that the WARMC result in the previous equation coincides with findings in Ref. 88 and represents the relativistic Bohm-Gross relationship for warm Langmuir waves, where the 5/2 correction term stands for thermal inertia, i.e., hydrodynamic transport correction due to temperature.

The 1D3V framework is such that only longitudinal waves are present in the system, and therefore, one has phase velocity equal to the bunch velocity

v_{p h} = \frac{ω}{k} = c

⁠. Therefore, in this context, one obtains at order

O (μ_{i})

{(\frac{ω}{ω_{0}})}^{2} = 1 - \frac{5}{6} μ_{i} LEC,

(65)

{(\frac{ω}{ω_{0}})}^{2} = 1 + \frac{1}{2} μ_{i} WARMC .

(66)

At variance with what was observed in Sec. VI A, Eq. foresees a closure dependent behavior of the wave: with respect to the cold case (where

ω = ω_{p}

⁠), LEC plasmas exhibit wavelengths decreasing with temperature, whereas WARMC plasmas wavelengths increase with temperature. This can be observed in Fig. 3, where we report a comparison between the theoretical predictions and the numerical results. In both the two closure schemes, we observe a very good agreement.

FIG. 3.

View large Download slide

Dispersion relation dependency on temperature. Comparison of the numerics obtained from the moment matching LB simulations vs the theory predictions Eqs. (65) and (66). In both the two closures, our numerical codes well reproduce the theoretical results. We remark that theory predictions are operated in the linear regime (small $\tilde{Q}$ perturbation) and in the assumption of small initial temperatures (small μ_i values).

C. 3D3V results

We now move our analysis to fully spatially resolved simulations. In Fig. 4, we provide an initial comparison of the two fluid closures in a 3D3V axisymmetric environment. The temperature is set to an initial value of $k_{b} T_{i} = 20$ $keV$ ⁠, and the chosen regime of the driving pulse is quasi-linear, $\tilde{Q} = 0.5$ ⁠. In the figure, we show color plots for four significant quantities: the plasma number density n, the longitudinal fluid velocity u_z, the longitudinal accelerating field E_z, and the transversal wakefield $E_{r} - c B_{ϕ}$ ⁠. It is possible to appreciate that at variance with the usual dynamics of the cold case,⁶² the typical wave-like pattern of peaks and valleys is perturbed by the emergence of an acoustic behavior that promotes electron motion out of the radial axis. The presence of pressure waves can be appreciated by the appearance of Mach cone structures⁹⁵ that can be easily spotted by looking at anyone of the panels of Fig. 4: as the driving bunch travels at the speed of light c, it perturbs the plasma medium. The perturbation propagates at a specific velocity c_s (that depends on the initial temperature μ_i) smaller than c, and the wave-fronts of the perturbation, that in a sub-sonic flow would radially propagate all over the space, are instead constrained in a conical structure around the perturbation: a supersonic Mach cone is observed.

FIG. 4.

View large Download slide

Qualitative comparison of 3D3V LEC-LB simulation against the WARMC-LB counterpart. We show four significant quantities, all properly adimensionalized w.r.t. previously defined values: the particle density n, the longitudinal fluid velocity u_z, and both the two components of the wakefields, E_z and $E_{r} - c B_{ϕ}$ ⁠. The chosen initial temperature value is $μ_{i} = 0.04 \to k_{b} T_{i} = 20$ $keV$ (again, selected for enhancing thermal effects in the first wave periods) while $\tilde{Q} = 0.5$ ⁠.

From a qualitative point of view, at this stage of the analysis, there is no clear difference in behavior between the two fluid closures. A more quantitative characterization can be performed by analyzing these acoustic structures in the linear regime (i.e., for small values of the parameter

\tilde{Q}

⁠), where it is possible to elaborate an analytical theory for both the two fluid models: in a nutshell, by using perturbation theory on, respectively, Eqs. and , it is possible to derive a forced Klein–Gordon equation for the electron plasma density n, where the only difference between the fluid closures appear in the temperature dependent coefficients of the equation (see the Appendix B for details). By looking at this equation, the acoustic nature of the dynamics becomes clear, and it is possible to derive some analytical results to compare with the numerical simulations. In Fig. 5, we show the result of numerical simulations (top half of the panels) vs the analytical solutions (bottom half of the panel). There is good agreement between the two. By inspecting the two cited panels, it is evident that the two fluid models, although exhibiting similar qualitative features, differ in the inclination Θ (plotted for reference in the figure) of the conical structure. The theory of supersonic flows⁹⁵ tells us that this inclination is linked to the intensity of the velocity c_s via

c_{s} = c sin (Θ),

(67)

which means that the wavefront of the perturbation propagates in the two models with different velocities. The warm linear theory provides, in fact, a dependency (look at Appendix B for more analytical details) of c_s from the initial temperature

{(\frac{c_{s}}{c})}^{2} = \frac{5}{3} μ_{i} LEC,

(68)

{(\frac{c_{s}}{c})}^{2} = 3 μ_{i} WARMC .

(69)

FIG. 5.

View large Download slide

Comparison of simulations (top half of both panels) against analytic solutions of the warm linear theory (bottom half of both panels), for both the two closures (top panel LEC, while bottom panel WARMC). The quantity shown is the particle density n, normalized w.r.t. the initial plasma density n_i. We also plot the Mach cone structure, together with its inclination Θ for reference.

One can then run a set of simulations with different initial temperatures μ_i and study the inclination of the conical envelope (from which we extract c_s) as a function of the temperature. In Fig. 6, we show the results of this analysis that again show good match between the simulations and expectations from theory.

FIG. 6.

View large Download slide

Acoustic velocity c_s dependency on the initial plasma temperature μ_i. Comparison of the numerics obtained from the moment matching LB simulations vs the warm linear theory predictions Eqs. (68) and (69). In both the two closures, our numerical codes well reproduce the theoretical results. We remark that theory predictions are operated in the linear regime (small $\tilde{Q}$ perturbation) and in the assumption of small initial temperatures (small μ_i values).

As a final element of discussion, we present the results of Fig. 7, where we consider the rest frame pressures ratio

P_{∥} / P_{⊥}

for the WARMC model. This quantity is expected to be equal to one in the case of an isotropic pressure tensor, as it is the case in the LEC model; therefore, in Fig. 7, we consider the following quantity:

\frac{P_{∥}}{P_{⊥}} - 1 .

(70)

FIG. 7.

View large Download slide

Study of the thermal spread anisotropies in the WARMC model. The pressure ratio $P_{∥} / P_{⊥}$ ⁠, whose value is obtained in the numerics via inversion of Eq. (60), is compared against its isotropic value 1. Differences of order 10% are observed in proximity of the bunch. The chosen initial temperature value is $μ_{i} = 0.04 \to k_{b} T_{i} = 20$ $keV$ (again, selected for enhancing thermal effects in the first wave periods) while $\tilde{Q} = 0.5$ ⁠.

The fluid rest frame pressures $P_{∥}$ and $P_{⊥}$ are obtained from the numerics by inverting Eq. (60) (the lab-frame components of the $θ^{μ ν}$ tensor are obtained from the simulations). We note thermal spread anisotropies are mostly relevant in the proximity of the bunch (located at $ζ k_{p} = 3.0$ ⁠) where they become of the order of 10% and become less and less important as one proceeds further away from the perturbation. To the best of our knowledge, this is the first time that such analysis is conducted on the WARMC model for a spatially resolved plasma (a study on a laser excited, 1D restricted plasma can be found in Ref. 96) and this preliminary results indicate that the LEC model, that is missing this pressure anisotropy feature by construction, could still be used to characterize plasma behavior in late-stage dynamics studies. We reserve further investigation, in particular studies on the dependency of this anisotropic feature on driving bunch parameters and initial temperatures, for later research.

VII. OUTLOOK AND PERSPECTIVES

In the development of a fluid model for the simulations of PWFA processes, a multitude of physical ingredients has to be taken into account to provide a realistic description of the process. In an earlier paper,⁶² some of the authors started the exploration of lattice Boltzmann (LB) schemes for the construction of fluid models for PWFA and considered the simplified case of cold fluid models. This paper represents a step forward, in that we have explored LB fluid schemes accounting for thermal effects.^{30–32,34,39,41,44–46} The inclusion of thermal effects is a rather non trivial task due to the theoretical complication represented by the choice of a proper closure to the set of equations.⁵⁴ We have handled this problem by selecting two of the most popular closure models that have been discussed in the literature: the first one relies on the assumption of a local equilibrium (LEC),⁴⁸ while the second one involves truncating the hierarchy of centered equations at an arbitrary order (WARMC).⁵⁷ We have then shown how to successfully adapt the LB schemes to both closures. If from one side the LEC is nominally not appropriate for a collisionless warm plasma, from the other side, the WARMC is obtained under the assumption of an asymptotically small temperature. Any finite temperature, however small, can raise the question on what is the right closure scheme to obtain the correct fluid model for collisionless warm plasma dynamics. The preliminary comparisons shown in this work, although not presenting strong qualitative differences in the dynamics, give a clear indication that the selection of a closure scheme is pivotal for the quantitative assessment of PWFA experiments. To this aim, a one-to-one comparison between the predictions of fluid models and the PIC simulations (or a numerical solution of the Vlasov equations) could be helpful to shed lights on the matter. Work is in progress along this direction.

The next logical step in the development of the method is the inclusion of ions dynamics. Since plasma ions are way more massive then electrons, their dynamics happens on longer time scales than the ones examined in PWFA studies that target the early stage evolution of system, unless a strongly dense driving bunch is considered.⁹⁷ There is though a growing research interest for studies that inquire the late time evolution of the system where this kind of dynamics cannot be ignored anymore:^36,38,40,98 unlocking the movement of ions brings in a new wealth of physical processes such as soliton dynamic⁴⁰ and ion channels formation.³⁹ Furthermore, most of these studies invoke thermal effects to explain the acoustic motion of the ions,^38,40,41,43 and it is then clear how the development of a method that would be able of handling both temperature and ions dynamics is of the uttermost importance. The inclusion of ion motion in the numerical scheme presented in this paper is theoretically trivial and only brings in a bigger computational effort (more equations need to be integrated): the aforementioned characterization is, therefore, completely within the reach of the LB method.

We would also like to mention a couple of advancements that would make the method presented in this paper a more complete numerical tool for the simulation of wakefield acceleration processes: the first is represented by the possibility of including a laser source field as the plasma perturbating mechanism, as the LWFA has always been an alternative (w.r.t. PWFA) route to wakefield acceleration;⁸ the second is the possibility of handling non-rigid particle bunches, as the tracking of the driving bunch properties is both an important element of diagnostic and also a feature to be kept under control in modern day PWFA experiments (see, for example, Ref. 99) Again, work in this direction is in progress. On the computational side, we would like to mention that the code developed for this work is amenable to a number of numerical optimizations (such for example porting on GPUs) and it is a firm intention of the authors to realize these improvements in future works.

As the most important design choice presented in this paper is the selection of the LB method as a fluid solver, a few comments are in place. In this paper, we have shown that the LB method is well capable of recovering analytic results, and we would like to mention further advantages that could make LB a convenient and viable option in the context of PWFA simulations. Although PIC methods are able to model the dynamics at the particle level, and hence are able to model kinetic effects, they constantly have to keep at bay their inherent numerical noise.^9,100 Therefore, in some situations, it would be preferable to dispose of alternative tools for prototyping and to reserve more complete PIC simulations for later stages. Solvers based on a strict discretization of the Vlasov equations would not suffer of the noise problem and still be able to capture mesoscopic effects, but for any dimensionality more than 1D, they would be too computationally demanding and hungry for memory resources. Fluid solvers (and the LB method among them) present then a viable option for the realization of quick PWFA simulations. They are coarse-grained and hence do not suffer of numerical noise, and still capture a good amount of the physical effects that are encountered in PWFA. Furthermore, the LB method lends itself exceptionally well to parallelization on both CPUs and GPUs.⁶⁰ In our implementation, we achieved parallelization on multi CPUs using the message passing interface (MPI). This approach involves dividing the computational domain into multiple rectangular sub-domains, corresponding to the number of processors. Leveraging the local nature of computations [Eq. (7), r.h.s.], communications are solely required to exchange populations between neighboring processors during the “streaming” process [Eq. (7), l.h.s.]. By effectively disentangling “compute” and “communicate,” this strategy greatly enhances the parallelization process.^101–104 Simulations in this study were conducted on an Intel Xeon E5-2695@2.40 GHz processor. A representative simulation (like those shown in Figs. 4 and 5) run on 96 processors requires approximately 182 min for the LEC-LB model and 368 min for the WARMC-LB model for $3 \cdot 10^{4}$ time steps. Memory requirements for such simulations are ∼1.6 GB for the LEC-LB model and ∼4.3 GB for the WARMC-LB model although these numbers could be decreased by further optimizations of the code. As a final note, we would like to remark that LB takes roughly $50 %$ of the compute time in the simulation, whereas the Maxwell solver for the electromagnetic fields takes 1%. This is to be expected as the amount of computation appearing in LB is significantly higher. Figure 8 presents the execution time per iteration and speedup data for varying numbers of processors (strong scaling).

FIG. 8.

View large Download slide

Top panel: execution time per iteration (in [s]) as a function of the number of processors for both LEC-LB (blue dotted) and WARMC-LB (orange solid). The scale is logarithmic. Bottom panel: corresponding speedups (data from the two closures are mostly overlapping). The simulation setup is the same as for results in Figs. 4 and 5.

Furthermore, an advantage of the LB method over Vlasov solvers is that the first adopts a smart pruning of the velocity space,^60,61 thus improving the computational efficiency. It is possible to increase the number of discrete velocities N_pop (see Sec. II A). However, this leads to a proportional increase in both the computational cost and memory requirements, making it a crucial factor to consider when choosing the LB stencil. Common practice is, indeed, to choose the minimum value N_pop that ensures the recovery of the hydrodynamic properties.⁶⁰ We hasten to remark, however, that at variance with usual hydrodynamic solvers, LB's theoretical formulation is strongly grounded in kinetic theory, and its fluid behavior is obtained via the smart discretization of the velocity space cited above. This is a remarkable feature that makes the method a strong candidate for the inclusion of kinetic effects into PWFA fluid solvers. In fact, recent studies^105–107 in other research fields show that LB is, indeed, capable of capturing behaviors beyond hydrodynamics just by increasing the number of discrete kinetic velocities. This could open novel perspectives for developing more refined numerical schemes for simulations of PWFA processes, with the obvious need of a precise comparison/benchmark against some reference PIC/Vlasov simulations. Further investigation on this is in progress.

We finally would like to mention that an alternative route to LB wakefield simulation could be represented by the Relativistic Lattice Boltzmann method.¹⁰⁸ This is an extension of LB hydrodynamic schemes (like the ones exposed in Sec. II A) to the theory of special relativity, originally created for simulation in astrophysical¹⁰⁹ and condensed matter^110,111 contexts, and, therefore, constitutes a promising tool for the simulation of warm plasmas within the LEC assumption. Adapting it to a PWFA framework is though more technical and less immediate w.r.t. the moment matching LB used in this paper, therefore, the authors reserve further development on this line for future works.

ACKNOWLEDGMENTS

The authors gratefully acknowledge Fabio Bonaccorso for his technical support. This work was supported by the Italian Ministry of University and Research (MUR) under the FARE program (No. R2045J8XAW), project “Smart-HEART.” MS gratefully acknowledges the support of the National Center for HPC, Big Data and Quantum Computing, Project CN_00000013 - CUP E83C22003230001, Mission 4 Component 2 Investment 1.4, funded by the European Union - NextGenerationEU.

AUTHOR DECLARATIONS

Conflict of Interest

The authors have no conflicts to disclose.

Author Contributions

Daniele Simeoni: Conceptualization (lead); Formal analysis (lead); Investigation (lead); Methodology (equal); Software (equal); Visualization (lead); Writing – original draft (lead). Gianmarco Parise: Conceptualization (supporting); Formal analysis (supporting); Software (equal); Writing – original draft (supporting). Fabio Guglietta: Conceptualization (supporting); Software (supporting); Writing – original draft (supporting). A. R. Rossi: Conceptualization (supporting); Writing – original draft (supporting). James Rosenzweig: Conceptualization (supporting); Writing – original draft (supporting). Alessandro Cianchi: Conceptualization (supporting); Writing – original draft (supporting). Mauro Sbragaglia: Conceptualization (supporting); Methodology (equal); Writing – original draft (supporting).

DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

APPENDIX A: FULL EXPRESSIONS FOR THE FLUID ADVECTION EQUATIONS

We report here the full form of the advection equations solved via the moment matching LB method in both the LEC and WARMC model. The equations are expressed in a 3D3V axisymmetric environment. In this context, we slightly modify the advection operator

D_{u}

defined in the main text to adopt its cylindrical counterpart

{\tilde{D}}_{u} (A) \equiv \partial_{t} A + \partial_{r} (A u_{r}) + \partial_{z} (A u_{z}) .

(A1)

All cylindrical correction terms deriving from the divergence transformation are moved to the RHS of the equations and are treated as source terms. The axisymmetry condition provides some useful simplifications:

The azimuthal fluid velocity is zero

$0 = u_{ϕ} .$

(A2)
Axisymmetric Maxwell equations lead to some electromagnetic field components to be zero

$0 = E_{ϕ} = B_{r} = B_{z} .$

(A3)
In the WARMC model, thanks to the strategy explained in the main text (end of Sec. V A) one can derive rest frame quantities

$0 = θ^{0 ϕ} = θ^{r ϕ} = θ^{z ϕ},$

(A4)

$P_{⊥} = θ^{ϕ ϕ},$

(A5)

$P_{∥} = - \frac{c^{2} θ_{μ}^{μ} + θ^{00} u^{2} + θ^{ϕ ϕ} u_{r}^{2}}{c^{2} + u_{z}^{2}} .$

(A6)

1. LEC equations

{\tilde{D}}_{u} (n) = - \frac{u_{r} n}{r},

(A7)

{\tilde{D}}_{u} (n A_{r}) = - \frac{u_{r} A_{r}}{r} - n_{i} m_{e} c^{2} μ_{i} \partial_{r} {(\frac{n}{γ n_{i}})}^{5 / 3} - e n (E_{r} - u_{z} B_{ϕ}),

(A8)

{\tilde{D}}_{u} (n A_{z}) = - \frac{u_{r} A_{z}}{r} - n_{i} m_{e} c^{2} μ_{i} \partial_{z} {(\frac{n}{γ n_{i}})}^{5 / 3} - e n (E_{z} + u_{r} B_{ϕ}),

(A9)

with

(\begin{matrix} A_{r} \\ A_{z} \end{matrix}) = [1 + \frac{5}{2} μ_{i} {(\frac{n}{γ n_{i}})}^{2 / 3}] n (\begin{matrix} p_{r} \\ p_{z} \end{matrix}) .

(A10)

2. WARMC equations

Define

λ = n / h

{\tilde{D}}_{u} (n) = - \frac{n u_{r}}{r},

(A11)

{\tilde{D}}_{u} (n λ c) = - \frac{n λ c u_{r}}{r} - \frac{1}{c} \partial_{t} θ^{00} - \partial_{r} θ^{0 r} - \partial_{z} θ^{0 z} - \frac{θ^{0 r}}{r} - \frac{e n}{c} (E_{r} u_{r} + E_{z} u_{z}),

(A12)

{\tilde{D}}_{u} (n λ u_{r}) = - \frac{n λ u_{r}^{2}}{r} - \frac{1}{c} \partial_{t} θ^{0 r} - \partial_{r} θ^{r r} - \partial_{z} θ^{r z} - \frac{θ^{r r} - θ^{ϕ ϕ}}{r} - e n (E_{r} - u_{z} B_{ϕ}),

(A13)

{\tilde{D}}_{u} (n λ u_{z}) = - \frac{n λ u_{z} u_{r}}{r} - \frac{1}{c} \partial_{t} θ^{0 z} - \partial_{r} θ^{r z} - \partial_{z} θ^{z z} - \frac{θ^{r z}}{r} - e n (E_{z} + u_{r} B_{ϕ}),

(A14)

{\tilde{D}}_{u} (λ θ^{r r}) = - \frac{λ θ^{r r} u_{r}}{r} - 2 \frac{θ^{0 r}}{c} \partial_{t} (λ u_{r}) - 2 θ^{r r} \partial_{r} (λ u_{r}) - 2 θ^{r z} \partial_{z} (λ u_{r}) - 2 \frac{e}{c} (E_{r} θ^{0 r} - c B_{ϕ} θ^{r z}),

(A15)

{\tilde{D}}_{u} (λ θ^{ϕ ϕ})) = - 3 \frac{λ θ^{ϕ ϕ} u_{r}}{r},

(A16)

{\tilde{D}}_{u} (λ θ^{z z}) = - \frac{λ θ^{z z} u_{r}}{r} - 2 \frac{θ^{0 z}}{c} \partial_{t} (λ u_{z}) - 2 θ^{r z} \partial_{r} (λ u_{z}) - 2 θ^{z z} \partial_{z} (λ u_{z}) - 2 \frac{e}{c} (E_{z} θ^{0 z} + c B_{ϕ} θ^{r z}),

(A17)

\begin{matrix} {\tilde{D}}_{u} (λ θ^{r z}) = - \frac{λ θ^{r z} u_{r}}{r} - \frac{θ^{0 z}}{c} \partial_{t} (λ u_{r}) - \frac{θ^{0 r}}{c} \partial_{t} (λ u_{z}) - θ^{r z} \partial_{r} (λ u_{r}) - θ^{r r} \partial_{r} (λ u_{z}) - θ^{z z} \partial_{z} (λ u_{r}) - θ^{r z} \partial_{z} (λ u_{z}) \\ - \frac{e}{c} [E_{r} θ^{0 z} + E_{z} θ^{0 r} + c B_{ϕ} (θ^{r r} - θ^{z z})] . \end{matrix}

(A18)

The mass-shell conditions Eqs. and provide the remaining non-zero components of the tensor

θ^{μ ν}

θ^{0 r} = \frac{u_{r}}{c} θ^{r r} + \frac{u_{z}}{c} θ^{r z},

(A19)

θ^{0 z} = \frac{u_{r}}{c} θ^{r z} + \frac{u_{z}}{c} θ^{z z},

(A20)

θ^{00} = \frac{u_{r}}{c} θ^{0 r} + \frac{u_{z}}{c} θ^{0 z} .

(A21)

APPENDIX B: WARM LINEAR THEORY

The fluid equations, respectively, Eqs. and for the two closures, can be linearly perturbed w.r.t. the initial rest state when coupled with Maxwell equations Eqs. and . One then obtains a forced Klein–Gordon equation for the density perturbation

n_{1} = n - n_{i}

\partial_{t}^{2} n_{1} - c_{s}^{2} \nabla^{2} n_{1} = - a_{e}^{2} n_{1} - a_{b}^{2} n_{b},

(B1)

where the coefficients c_s, a_e, and a_b depend on the initial temperature in a way that's dictated by the selected closure scheme. At order

O (μ_{i})

⁠, one has

LEC {\begin{matrix} {(\frac{c_{s}}{c})}^{2} = \frac{5}{3} μ_{i}, \\ {(\frac{a_{e}}{ω_{p}})}^{2} = 1 - \frac{5}{2} μ_{i}, \\ a_{b}^{2} = a_{e}^{2}, \end{matrix}

(B2)

WARMC {\begin{matrix} {(\frac{c_{s}}{c})}^{2} = 3 μ_{i}, \\ {(\frac{a_{e}}{ω_{p}})}^{2} = 1 - \frac{5}{2} μ_{i}, \\ {(\frac{a_{b}}{ω_{p}})}^{2} = 1 - \frac{1}{2} μ_{i} . \end{matrix}

(B3)

This result, aside from giving an explicit dependency of the acoustic velocity c_s w.r.t. temperature [result Eqs. and ], tells us that the behavior of the two fluid schemes is qualitatively the same. One can further inspect Eq. to determine a dispersion relation for the Langmuir wave. Assuming to be far away from the driving bunch (so that n_b can be ignored), one derives the following by recurring to Fourier transformation

(\partial_{t}, \nabla) \to (- i ω, i k)

⁠:

ω^{2} = a_{e}^{2} + c_{s}^{2} k^{2},

(B4)

which reduces to the result Eqs. and when considered at order

O (μ_{i})

⁠. Last thing to discuss is the strategy to the solution of Eq. . After performing the traveling wave ansatz, which states that every field is dependent on z and t only through the co-moving coordinate

ζ = z - c t

⁠, we obtain the following equation:

(c^{2} - c_{s}^{2}) \partial_{ζ}^{2} n_{1} - c_{s}^{2} \nabla_{⊥}^{2} n_{1} = - a_{e}^{2} n_{1} - a_{b}^{2} n_{b} .

(B5)

At this point, we introduce the Hankel Tranformation of order 0¹¹²

H_{0} [f (r)] = \hat{f} (w) = \int_{0}^{+ \infty} f (r) r J_{0} (w r) d r,

(B6)

where

J_{0} (w r)

is the Bessel Function of first kind and order 0. This transform behaves nicely in the presence of transversal laplacians of radial functions

(\nabla_{⊥}) \to (i w)

⁠. Therefore, if one Hankel transforms Eq. , a Forced Harmonic oscillator in the variable ζ is obtained

\begin{matrix} \partial_{ζ}^{2} {\hat{n}}_{1} (w, ζ) = - (\frac{a_{e}^{2} + w^{2} c_{s}^{2}}{c^{2} - c_{s}^{2}}) {\hat{n}}_{1} (w, ζ) - (\frac{a_{b}^{2}}{c^{2} - c_{s}^{2}}) {\hat{n}}_{b} (w, ζ), \end{matrix}

(B7)

which can be solved via the Green’s function method. The final anti-transformed solution, to be numerically integrated, reads, therefore, as

n_{1} (r, ζ) = - \frac{a_{b}^{2}}{\sqrt{c^{2} - c_{s}^{2}}} \int_{- \infty}^{ζ} \int_{0}^{+ \infty} \frac{{\hat{n}}_{b} (w, ζ') J_{0} (w r)}{\sqrt{a_{e}^{2} + w^{2} c_{s}^{2}}} sin [\sqrt{\frac{a_{e}^{2} + w^{2} c_{s}^{2}}{c^{2} - c_{s}^{2}}} (ζ - ζ')] dwd ζ' .

(B8)

REFERENCES

1.

S.

Gourlay

,

T.

Raubenheimer

, and

V.

Shiltsev

,

Front. Phys.

10

,

920520

(

2022

).

https://doi.org/10.3389/fphy.2022.920520

Google Scholar

Crossref

2.

J.

Faure

,

Y.

Glinec

,

A.

Pukhov

,

S.

Kiselev

,

S.

Gordienko

,

E.

Lefebvre

,

J.-P.

Rousseau

,

F.

Burgy

, and

V.

Malka

,

Nature

431

,

541

(

2004

).

https://doi.org/10.1038/nature02963

Google Scholar

Crossref

PubMed

3.

I.

Blumenfeld

,

C. E.

Clayton

,

F.-J.

Decker

,

M. J.

Hogan

,

C.

Huang

,

R.

Ischebeck

,

R.

Iverson

,

C.

Joshi

,

T.

Katsouleas

,

N.

Kirby

,

W.

Lu

,

K. A.

Marsh

,

W. B.

Mori

,

P.

Muggli

,

E.

Oz

,

R. H.

Siemann

,

D.

Walz

, and

M.

Zhou

,

Nature

445

,

741

(

2007

).

https://doi.org/10.1038/nature05538

Google Scholar

Crossref

PubMed

4.

M.

Litos

,

E.

Adli

,

W.

An

,

C. I.

Clarke

,

C. E.

Clayton

,

S.

Corde

,

J. P.

Delahaye

,

R. J.

England

,

A. S.

Fisher

,

J.

Frederico

,

S.

Gessner

,

S. Z.

Green

,

M. J.

Hogan

,

C.

Joshi

,

W.

Lu

,

K. A.

Marsh

,

W. B.

Mori

,

P.

Muggli

,

N.

Vafaei-Najafabadi

,

D.

Walz

,

G.

White

,

Z.

Wu

,

V.

Yakimenko

, and

G.

Yocky

,

Nature

515

,

92

(

2014

).

https://doi.org/10.1038/nature13882

Google Scholar

Crossref

PubMed

5.

P.

Chen

,

J. M.

Dawson

,

R. W.

Huff

, and

T.

Katsouleas

,

Phys. Rev. Lett.

54

,

693

(

1985

).

https://doi.org/10.1103/PhysRevLett.54.693

Google Scholar

Crossref

PubMed

6.

T.

Tajima

and

J. M.

Dawson

,

Phys. Rev. Lett.

43

,

267

(

1979

).

https://doi.org/10.1103/PhysRevLett.43.267

Google Scholar

Crossref

7.

R.

Bingham

,

J. T.

Mendonça

, and

P. K.

Shukla

,

Plasma Phys. Controlled Fusion

46

,

R1

(

2003

).

https://doi.org/10.1088/0741-3335/46/1/R01

Google Scholar

Crossref

8.

E.

Esarey

,

C. B.

Schroeder

, and

W. P.

Leemans

,

Rev. Mod. Phys.

81

,

1229

(

2009

).

https://doi.org/10.1103/RevModPhys.81.1229

Google Scholar

Crossref

9.

C.

Birdsall

and

A.

Langdon

,

Plasma Physics via Computer Simulation

, Series in Plasma Physics and Fluid Dynamics (

CRC Press

,

2018

).

Google Scholar

Crossref

10.

R. A.

Fonseca

,

L. O.

Silva

,

F. S.

Tsung

,

V. K.

Decyk

,

W.

Lu

,

C.

Ren

,

W. B.

Mori

,

S.

Deng

,

S.

Lee

,

T.

Katsouleas

, and

J. C.

Adam

, in

Computational Science—ICCS 2002

, edited by

P. M. A.

Sloot

,

A. G.

Hoekstra

,

C. J. K.

Tan

, and

J. J.

Dongarra

(

Springer

,

Berlin, Heidelberg

,

2002

), pp.

342

–

351

.

Google Scholar

Crossref

11.

R.

Lehe

,

M.

Kirchen

,

I. A.

Andriyash

,

B. B.

Godfrey

, and

J.-L.

Vay

,

Comput. Phys. Commun.

203

,

66

(

2016

).

https://doi.org/10.1016/j.cpc.2016.02.007

Google Scholar

Crossref

12.

H.

Burau

,

R.

Widera

,

W.

Hönig

,

G.

Juckeland

,

A.

Debus

,

T.

Kluge

,

U.

Schramm

,

T. E.

Cowan

,

R.

Sauerbrey

, and

M.

Bussmann

,

IEEE Trans. Plasma Sci.

38

,

2831

(

2010

).

https://doi.org/10.1109/TPS.2010.2064310

Google Scholar

Crossref

13.

C.

Benedetti

,

A.

Sgattoni

,

G.

Turchetti

, and

P.

Londrillo

,

IEEE Trans. Plasma Sci.

36

,

1790

(

2008

).

https://doi.org/10.1109/TPS.2008.927143

Google Scholar

Crossref

14.

R.

Hockney

and

J.

Eastwood

,

Computer Simulation Using Particles

(

CRC Press

,

2021

).

Google Scholar

Crossref

15.

C.

Cercignani

and

G. M.

Kremer

,

The Relativistic Boltzmann Equation: Theory and Applications

(

Birkhäuser Basel

,

2002

).

Google Scholar

Crossref

16.

S. R.

de Groot

,

W. A.

van Leeuwen

, and

C. G.

van Weert

,

Relativistic Kinetic Theory: Principles and Applications

(

North-Holland Pub. Co

.,

Amsterdam

,

1980

).

Google Scholar

17.

L.

Rezzolla

and

O.

Zanotti

,

Relativistic Hydrodynamics

(

Oxford University Press

,

2013

).

Google Scholar

Crossref

18.

A.

Marocchino

,

F.

Massimo

,

A.

Rossi

,

E.

Chiadroni

, and

M.

Ferrario

,

Nucl. Instrum. Methods Phys. Res., Sect. A

829

,

386

(

2016

).

https://doi.org/10.1016/j.nima.2016.03.005

Google Scholar

Crossref

19.

P.

Tomassini

and

A.

Rossi

,

Plasma Phys. Controlled Fusion

58

,

034001

(

2015

).

https://doi.org/10.1088/0741-3335/58/3/034001

Google Scholar

Crossref

20.

A.

Boldarev

,

V.

Gasilov

,

O.

Olkhovskaya

,

S.

Dyachenko

,

G.

Bagdasarov

,

S.

Boldyrev

,

I.

Gasilova

, and

E.

Dorofeeva

, “

Object-oriented code MARPLE3D: Simulations of radiative hydrodynamic/MHD effects at high-performance computer systems

,” ECCOMAS 2012–European Congress on Computational Methods in Applied Sciences and Engineering (

2012

).

Google Scholar

21.

N.

Cook

,

O.

Tresca

,

N. P.

Dover

,

C.

Maharjan

,

M. N.

Polyanskiy

,

Z.

Najmudin

,

P.

Shkolnikov

, and

I.

Pogorelsky

,

AIP Conf. Proc.

1777

,

090002

(

2016

).

https://doi.org/10.1063/1.4965677

Google Scholar

Crossref

22.

S. M.

Mewes

,

G. J.

Boyle

,

A. F.

Pousa

,

R. J.

Shalloo

,

J.

Osterhoff

,

C.

Arran

,

L.

Corner

,

R.

Walczak

,

S. M.

Hooker

, and

M.

Thévenet

,

Phys. Rev. Res.

5

,

033112

(

2023

).

https://doi.org/10.1103/PhysRevResearch.5.033112

Google Scholar

Crossref

23.

M.

Anania

,

E.

Chiadroni

,

A.

Cianchi

,

D.

Di Giovenale

,

M.

Ferrario

,

F.

Flora

,

G.

Gallerano

,

A.

Ghigo

,

A.

Marocchino

,

F.

Massimo

,

A.

Mostacci

,

L.

Mezi

,

P.

Musumeci

, and

M.

Serio

,

Nucl. Instrum. Methods Phys. Res., Sect. A

740

,

193

(

2014

).

https://doi.org/10.1016/j.nima.2013.10.053

Google Scholar

Crossref

24.

A. J.

Gonsalves

,

K.

Nakamura

,

J.

Daniels

,

C.

Benedetti

,

C.

Pieronek

,

T. C. H.

de Raadt

,

S.

Steinke

,

J. H.

Bin

,

S. S.

Bulanov

,

J.

van Tilborg

,

C. G. R.

Geddes

,

C. B.

Schroeder

,

C.

Tóth

,

E.

Esarey

,

K.

Swanson

,

L.

Fan-Chiang

,

G.

Bagdasarov

,

N.

Bobrova

,

V.

Gasilov

,

G.

Korn

,

P.

Sasorov

, and

W. P.

Leemans

,

Phys. Rev. Lett.

122

,

084801

(

2019

).

https://doi.org/10.1103/PhysRevLett.122.084801

Google Scholar

Crossref

PubMed

25.

R. D.

Ruth

,

A. W.

Chao

,

P. L.

Morton

, and

P. B.

Wilson

,

Part. Accel.

17

,

171

(

1985

).

Google Scholar

26.

P.

Chen

,

J.

Su

,

T.

Katsouleas

,

S.

Wilks

, and

J.

Dawson

,

IEEE Trans. Plasma Sci.

15

,

218

(

1987

).

https://doi.org/10.1109/TPS.1987.4316688

Google Scholar

Crossref

27.

W.

Lu

,

C.

Huang

,

M. M.

Zhou

,

W. B.

Mori

, and

T.

Katsouleas

,

Phys. Plasmas

12

,

063101

(

2005

).

https://doi.org/10.1063/1.1905587

Google Scholar

Crossref

28.

A. I.

Akhiezer

and

R. V.

Polovin

,

Sov. Phys. JETP

3

,

696

–

705

(

1956

) [Zh. Eksp. Teor. Fiz. 30(5), 915 (1956)].

Google Scholar

29.

J. M.

Dawson

,

Phys. Rev.

113

,

383

(

1959

).

https://doi.org/10.1103/PhysRev.113.383

Google Scholar

Crossref

30.

T.

Katsouleas

and

W. B.

Mori

,

Phys. Rev. Lett.

61

,

90

(

1988

).

https://doi.org/10.1103/PhysRevLett.61.90

Google Scholar

Crossref

PubMed

31.

J. B.

Rosenzweig

,

Phys. Rev. A

38

,

3634

(

1988

).

https://doi.org/10.1103/PhysRevA.38.3634

Google Scholar

Crossref

32.

C. B.

Schroeder

,

E.

Esarey

, and

B. A.

Shadwick

,

Phys. Rev. E

72

,

055401

(

2005

).

https://doi.org/10.1103/PhysRevE.72.055401

Google Scholar

Crossref

33.

N.

Jain

,

J.

Palastro

,

J.

Antonsen

,

T. M. W. B.

Mori

, and

W.

An

,

Phys. Plasmas

22

,

023103

(

2015

).

https://doi.org/10.1063/1.4907159

Google Scholar

Crossref

34.

E.

Esarey

,

C. B.

Schroeder

,

E.

Cormier-Michel

,

B. A.

Shadwick

,

C. G. R.

Geddes

, and

W. P.

Leemans

,

Phys. Plasmas

14

,

056707

(

2007

).

https://doi.org/10.1063/1.2714022

Google Scholar

Crossref

35.

K. V.

Lotov

,

Phys. Rev. ST Accel. Beams

6

,

061301

(

2003

).

https://doi.org/10.1103/PhysRevSTAB.6.061301

Google Scholar

Crossref

36.

R.

Gholizadeh

,

T.

Katsouleas

,

C.

Huang

,

W. B.

Mori

, and

P.

Muggli

,

Phys. Rev. ST Accel. Beams

14

,

021303

(

2011

).

https://doi.org/10.1103/PhysRevSTAB.14.021303

Google Scholar

Crossref

37.

M. F.

Gilljohann

,

H.

Ding

,

A.

Döpp

,

J.

Götzfried

,

S.

Schindler

,

G.

Schilling

,

S.

Corde

,

A.

Debus

,

T.

Heinemann

,

B.

Hidding

,

S. M.

Hooker

,

A.

Irman

,

O.

Kononenko

,

T.

Kurz

,

A.

Martinez de la Ossa

,

U.

Schramm

, and

S.

Karsch

,

Phys. Rev. X

9

,

011046

(

2019

).

https://doi.org/10.1103/PhysRevX.9.011046

Google Scholar

Crossref

38.

R.

D'Arcy

,

J.

Chappell

,

J.

Beinortaite

,

S.

Diederichs

,

G.

Boyle

,

B.

Foster

,

M. J.

Garland

,

P. G.

Caminal

,

C. A.

Lindstrøm

,

G.

Loisch

,

S.

Schreiber

,

S.

Schröder

,

R. J.

Shalloo

,

M.

Thévenet

,

S.

Wesch

,

M.

Wing

, and

J.

Osterhoff

,

Nature

603

,

58

(

2022

).

https://doi.org/10.1038/s41586-021-04348-8

Google Scholar

Crossref

PubMed

39.

R.

Zgadzaj

,

T.

Silva

,

V. K.

Khudyakov

,

A.

Sosedkin

,

J.

Allen

,

S.

Gessner

,

Z.

Li

,

M.

Litos

,

J.

Vieira

,

K. V.

Lotov

,

M. J.

Hogan

,

V.

Yakimenko

, and

M. C.

Downer

,

Nat. Commun.

11

,

4753

(

2020

).

https://doi.org/10.1038/s41467-020-18490-w

Google Scholar

Crossref

PubMed

40.

V. K.

Khudiakov

,

K. V.

Lotov

, and

M. C.

Downer

,

Plasma Phys. Controlled Fusion

64

,

045003

(

2022

).

https://doi.org/10.1088/1361-6587/ac4523

Google Scholar

Crossref

41.

T.

Silva

,

L. D.

Amorim

,

M. C.

Downer

,

M. J.

Hogan

,

V.

Yakimenko

,

R.

Zgadzaj

, and

J.

Vieira

,

Phys. Rev. Lett.

127

,

104801

(

2021

).

https://doi.org/10.1103/PhysRevLett.127.104801

Google Scholar

Crossref

PubMed

42.

A. A.

Sahai

,

T. C.

Katsouleas

,

F. S.

Tsung

, and

W. B.

Mori

, “

Long term evolution of plasma wakefields

,” in Proceedings of the North American Particle Accelerator Conference 2013 PAC2013, Pasadena, CA, Sep.-Oct. 2013 (JACoW Publishing, Geneva, Switzerland, 2014), pp. 90–92, see https://jacow.org/PAC2013/papers/MOPAC10.pdf.

43.

A. A.

Sahai

,

Phys. Rev. Accel. Beams

20

,

081004

(

2017

).

https://doi.org/10.1103/PhysRevAccelBeams.20.081004

Google Scholar

Crossref

44.

T.

Wang

,

V.

Khudik

, and

G.

Shvets

,

arXiv:2110.10290

(

2021

).

45.

S.

Diederichs

,

C.

Benedetti

,

E.

Esarey

,

M.

Thévenet

,

A.

Sinn

,

J.

Osterhoff

, and

C. B.

Schroeder

,

Phys. Plasmas

30

,

073104

(

2023

).

https://doi.org/10.1063/5.0155489

Google Scholar

Crossref

46.

S.

Diederichs

, “

Positron acceleration in a plasma column

,” Ph.D. thesis (

University of Hamburg

,

2023

).

Google Scholar

47.

G. J.

Cao

,

C. A.

Lindstrøm

,

E.

Adli

,

S.

Corde

, and

S.

Gessner

, “

Positron acceleration in plasma wakefields

,” arXiv:2309.10495 (2023).

48.

A. J.

Toepfer

,

Phys. Rev. A

3

,

1444

(

1971

).

https://doi.org/10.1103/PhysRevA.3.1444

Google Scholar

Crossref

49.

W. A.

Newcomb

,

Phys. Fluids

25

,

846

(

1982

).

https://doi.org/10.1063/1.863814

Google Scholar

Crossref

50.

P.

Amendt

and

H.

Weitzner

,

Phys. Fluids

28

,

949

(

1985

).

https://doi.org/10.1063/1.865066

Google Scholar

Crossref

51.

W. A.

Newcomb

,

Phys. Fluids

29

,

881

(

1986

).

https://doi.org/10.1063/1.865891

Google Scholar

Crossref

52.

J. G.

Siambis

,

Phys. Fluids

30

,

896

(

1987

).

https://doi.org/10.1063/1.866343

Google Scholar

Crossref

53.

S.

Pennisi

and

A. M.

Anile

,

Phys. Fluids B

3

,

1091

(

1991

).

https://doi.org/10.1063/1.859837

Google Scholar

Crossref

54.

O.

Muscato

,

Phys. Fluids B

5

,

3036

(

1993

).

https://doi.org/10.1063/1.860690

Google Scholar

Crossref

55.

F.

Jüttner

,

Annalen der Phys.

339

,

856

(

1911

).

https://doi.org/10.1002/andp.19113390503

Google Scholar

Crossref

56.

C. B.

Schroeder

,

J. Phys.: Conf. Ser.

169

,

012007

(

2009

).

https://doi.org/10.1088/1742-6596/169/1/012007

Google Scholar

Crossref

57.

C. B.

Schroeder

and

E.

Esarey

,

Phys. Rev. E

81

,

056403

(

2010

).

https://doi.org/10.1103/PhysRevE.81.056403

Google Scholar

Crossref

58.

W. A.

Newcomb

,

Phys. Fluids

29

,

1854

(

1986

).

https://doi.org/10.1063/1.865999

Google Scholar

Crossref

59.

P.

Amendt

,

Phys. Fluids

29

,

1458

(

1986

).

https://doi.org/10.1063/1.865664

Google Scholar

Crossref

60.

T.

Krüger

,

H.

Kusumaatmaja

,

A.

Kuzmin

,

O.

Shardt

,

G.

Silva

, and

E. M.

Viggen

,

The Lattice Boltzmann Method

(

Springer International Publishing

,

2017

).

Google Scholar

Crossref

61.

S.

Succi

,

The Lattice Boltzmann Equation: For Complex States of Flowing Matter

(

Oxford University Press

,

2018

).

Google Scholar

Crossref

62.

G.

Parise

,

A.

Cianchi

,

A.

Del Dotto

,

F.

Guglietta

,

A. R.

Rossi

, and

M.

Sbragaglia

,

Phys. Plasmas

29

,

043903

(

2022

).

https://doi.org/10.1063/5.0085192

Google Scholar

Crossref

63.

K.

Yee

,

IEEE Trans. Antennas Propag.

14

,

302

(

1966

).

https://doi.org/10.1109/TAP.1966.1138693

Google Scholar

Crossref

64.

J. F.

Olakangil

,

Finite-Difference Time-Domain Model of Plasma

(

Utah State University

,

2000

).

Google Scholar

65.

P. L.

Bhatnagar

,

E. P.

Gross

, and

M.

Krook

,

Phys. Rev.

94

,

511

(

1954

).

https://doi.org/10.1103/PhysRev.94.511

Google Scholar

Crossref

66.

S.

Chapman

and

T.

Cowling

,

The Mathematical Theory of Non-Uniform Gases

, 3rd ed., Cambridge Mathematical Library (

Cambridge University Press

,

1970

).

Google Scholar

67.

H.

Grad

,

Commun. Pure Appl. Math.

2

,

325

(

1949

).

https://doi.org/10.1002/cpa.3160020402

Google Scholar

Crossref

68.

H.

Grad

,

Commun. Pure Appl. Math.

2

,

331

(

1949

).

https://doi.org/10.1002/cpa.3160020403

Google Scholar

Crossref

69.

X.

Shan

,

X. F.

Yuan

, and

H.

Chen

,

J. Fluid Mech.

550

,

413

–

441

(

2006

).

https://doi.org/10.1017/S0022112005008153

Google Scholar

Crossref

70.

X.

Shan

,

Phys. Rev. E

81

,

036702

(

2010

).

https://doi.org/10.1103/PhysRevE.81.036702

Google Scholar

Crossref

71.

X.

Shan

,

J. Comput. Sci.

17

,

475

(

2016

).

https://doi.org/10.1016/j.jocs.2016.03.002

Google Scholar

Crossref

72.

X.

Shan

,

Phys. Rev. E

55

,

2780

(

1997

).

https://doi.org/10.1103/PhysRevE.55.2780

Google Scholar

Crossref

73.

X.

He

,

S.

Chen

, and

G. D.

Doolen

,

J. Comput. Phys.

146

,

282

(

1998

).

https://doi.org/10.1006/jcph.1998.6057

Google Scholar

Crossref

74.

I. V.

Karlin

,

D.

Sichau

, and

S. S.

Chikatamarla

,

Phys. Rev. E

88

,

063310

(

2013

).

https://doi.org/10.1103/PhysRevE.88.063310

Google Scholar

Crossref

75.

S.

Ponce Dawson

,

S.

Chen

, and

G. D.

Doolen

,

J. Chem. Phys.

98

,

1514

(

1993

).

https://doi.org/10.1063/1.464316

Google Scholar

Crossref

76.

K.

Yamamoto

,

X.

He

, and

G. D.

Doolen

,

J. Stat. Phys.

107

,

367

(

2002

).

https://doi.org/10.1023/A:1014583226083

Google Scholar

Crossref

77.

Z.

Wang

,

T.

Lei

, and

K. H.

Luo

,

Proc. Combust. Inst.

39

,

5365

(

2023

).

https://doi.org/10.1016/j.proci.2022.11.011

Google Scholar

Crossref

78.

F.

Verhaeghe

,

S.

Arnout

,

B.

Blanpain

, and

P.

Wollants

,

Phys. Rev. E

73

,

036316

(

2006

).

https://doi.org/10.1103/PhysRevE.73.036316

Google Scholar

Crossref

79.

Q.

Kang

,

D.

Zhang

, and

S.

Chen

,

J. Geophys. Res.: Solid Earth

108

,

ECV9

, https://doi.org/10.1029/2003JB002504 (

2003

).

Google Scholar

80.

B.

Chopard

,

J. L.

Falcone

, and

J.

Latt

,

Eur. Phys. J. Spec. Top.

171

,

245

(

2009

).

https://doi.org/10.1140/epjst/e2009-01035-5

Google Scholar

Crossref

81.

K. N.

Premnath

and

J.

Abraham

,

Phys. Rev. E

71

,

056706

(

2005

).

https://doi.org/10.1103/PhysRevE.71.056706

Google Scholar

Crossref

82.

S.

Srivastava

,

P.

Perlekar

,

J. H. M.

ten Thije Boonkkamp

,

N.

Verma

, and

F.

Toschi

,

Phys. Rev. E

88

,

013309

(

2013

).

https://doi.org/10.1103/PhysRevE.88.013309

Google Scholar

Crossref

83.

G. B.

Arfken

,

H. J.

Weber

, and

F. E.

Harris

,

Elsevier Mathematical Methods for Physicists

(

Elsevier India

,

2012

).

Google Scholar

84.

J. D.

Jackson

,

Classical Electrodynamics

, 3rd ed. (

John Wiley and Sons

,

Nashville, TN

,

1998

).

Google Scholar

85.

G.

Chacón-Acosta

,

L.

Dagdug

, and

H. A.

Morales-Técotl

,

Phys. Rev. E

81

,

021126

(

2010

).

https://doi.org/10.1103/PhysRevE.81.021126

Google Scholar

Crossref

86.

J.

Synge

,

The Relativistic Gas

, Series in Physics (

North-Holland Publishing Company

,

1957

).

Google Scholar

87.

P. S.

Letelier

,

Phys. Rev. D

22

,

807

(

1980

).

https://doi.org/10.1103/PhysRevD.22.807

Google Scholar

Crossref

88.

P. C.

Clemmow

,

A. J.

Willson

, and

J. A.

Ratcliffe

,

Proc. R. Soc. London, Ser. A

237

,

117

(

1956

).

https://doi.org/10.1098/rspa.1956.0165

Google Scholar

Crossref

89.

F.

Massimo

, “

Modeling acceleration of a high brightness electron beam by plasma wakefield

,” Ph.D. thesis (

University of Rome La Sapienza

,

2015

).

Google Scholar

90.

P.

Londrillo

,

C.

Gatti

, and

M.

Ferrario

,

Nucl. Instrum. Methods Phys. Res. Sect. A

740

,

236

(

2014

).

https://doi.org/10.1016/j.nima.2013.10.028

Google Scholar

Crossref

91.

F.

Massimo

,

A.

Marocchino

, and

A.

Rossi

,

Nucl. Instrum. Methods Phys. Res. Sect. A

829

,

378

(

2016

).

https://doi.org/10.1016/j.nima.2016.02.043

Google Scholar

Crossref

92.

D.

Bohm

and

E. P.

Gross

,

Phys. Rev.

75

,

1851

(

1949

).

https://doi.org/10.1103/PhysRev.75.1851

Google Scholar

Crossref

93.

I.

Hutchinson

and

J.

Freidberg

,

Intyroduction to Plasma Physics

(

Massachusetts Institute of Technology: MIT OpenCouseWares

,

2003

).

Google Scholar

94.

M. U.

Lee

,

G. S.

Yun

, and

J.-Y.

Ji

,

J. Korean Phys. Soc.

82

,

740

(

2023

).

https://doi.org/10.1007/s40042-023-00736-7

Google Scholar

Crossref

95.

L.

Landau

and

E.

Lifshitz

,

Fluid Mechanics, Second Edition: Volume 6 (Course of Theoretical Physics)

(

Butterworth-Heinemann

,

1987

).

Google Scholar

96.

E.

Khalilzadeh

,

J.

Yazdanpanah

,

J.

Jahanpanah

,

A.

Chakhmachi

, and

E.

Yazdani

,

Phys. Plasmas

22

,

113115

(

2015

).

https://doi.org/10.1063/1.4936276

Google Scholar

Crossref

97.

J. B.

Rosenzweig

,

A. M.

Cook

,

A.

Scott

,

M. C.

Thompson

, and

R. B.

Yoder

,

Phys. Rev. Lett.

95

,

195002

(

2005

).

https://doi.org/10.1103/PhysRevLett.95.195002

Google Scholar

Crossref

PubMed

98.

L. M.

Gorbunov

,

P.

Mora

, and

A. A.

Solodov

,

Phys. Plasmas

10

,

1124

(

2003

).

https://doi.org/10.1063/1.1559011

Google Scholar

Crossref

99.

R.

Pompili

,

D.

Alesini

,

M. P.

Anania

,

M.

Behtouei

,

M.

Bellaveglia

,

A.

Biagioni

,

F. G.

Bisesto

,

M.

Cesarini

,

E.

Chiadroni

,

A.

Cianchi

,

G.

Costa

,

M.

Croia

,

A.

Del Dotto

,

D.

Di Giovenale

,

M.

Diomede

,

F.

Dipace

,

M.

Ferrario

,

A.

Giribono

,

V.

Lollo

,

L.

Magnisi

,

M.

Marongiu

,

A.

Mostacci

,

L.

Piersanti

,

G.

Di Pirro

,

S.

Romeo

,

A. R.

Rossi

,

J.

Scifo

,

V.

Shpakov

,

C.

Vaccarezza

,

F.

Villa

, and

A.

Zigler

,

Nat. Phys.

17

,

499

(

2021

).

https://doi.org/10.1038/s41567-020-01116-9

Google Scholar

Crossref

100.

A. B.

Langdon

,

Phys. Fluids

22

,

163

(

1979

).

https://doi.org/10.1063/1.862452

Google Scholar

Crossref

101.

M.

Mazzeo

and

P.

Coveney

,

Comput. Phys. Commun.

178

,

894

(

2008

).

https://doi.org/10.1016/j.cpc.2008.02.013

Google Scholar

Crossref

102.

M.

Bernaschi

,

S.

Melchionna

,

S.

Succi

,

M.

Fyta

,

E.

Kaxiras

, and

J.

Sircar

,

Comput. Phys. Commun.

180

,

1495

(

2009

).

https://doi.org/10.1016/j.cpc.2009.04.001

Google Scholar

Crossref

103.

M.

Bernaschi

,

M.

Fatica

,

S.

Melchionna

,

S.

Succi

, and

E.

Kaxiras

,

Concurrency Comput.: Pract. Exper.

22

,

1

(

2010

).

https://doi.org/10.1002/cpe.1466

Google Scholar

Crossref

104.

S.

Succi

,

G.

Amati

,

M.

Bernaschi

,

G.

Falcucci

,

M.

Lauricella

, and

A.

Montessori

,

Comput. Fluids

181

,

107

(

2019

).

https://doi.org/10.1016/j.compfluid.2019.01.005

Google Scholar

Crossref

105.

V. E.

Ambruş

and

R.

Blaga

,

Phys. Rev. C

98

,

035201

(

2018

).

https://doi.org/10.1103/PhysRevC.98.035201

Google Scholar

Crossref

106.

L.

Bazzanini

,

A.

Gabbana

,

D.

Simeoni

,

S.

Succi

, and

R.

Tripiccione

,

J. Comput. Sci.

51

,

101320

(

2021

).

https://doi.org/10.1016/j.jocs.2021.101320

Google Scholar

Crossref

107.

V. E.

Ambruş

,

L.

Bazzanini

,

A.

Gabbana

,

D.

Simeoni

,

S.

Succi

, and

R.

Tripiccione

,

Nat. Comput. Sci.

2

,

641

(

2022

).

https://doi.org/10.1038/s43588-022-00333-x

Google Scholar

Crossref

PubMed

108.

A.

Gabbana

,

D.

Simeoni

,

S.

Succi

, and

R.

Tripiccione

, “

Relativistic lattice Boltzmann methods: Theory and applications

,”

Phys. Rep.

863

,

1

–

63

(

2020

).

https://doi.org/10.1016/j.physrep.2020.03.004

Google Scholar

Crossref

109.

L. R.

Weih

,

A.

Gabbana

,

D.

Simeoni

,

L.

Rezzolla

,

S.

Succi

, and

R.

Tripiccione

,

Mon. Not. Roy. Astron. Soc.

498

,

3374

(

2020

).

https://doi.org/10.1093/mnras/staa2575

Google Scholar

Crossref

110.

A.

Gabbana

,

M.

Mendoza

,

S.

Succi

, and

R.

Tripiccione

,

Comput. Fluids

172

,

644

(

2018

).

https://doi.org/10.1016/j.compfluid.2018.02.020

Google Scholar

Crossref

111.

D.

Simeoni

,

A.

Gabbana

, and

S.

Succi

,

Commun. Comput. Phys.

33

,

174

(

2023

).

https://doi.org/10.4208/cicp.OA-2022-0051

Google Scholar

Crossref

112.

A. C.

Offord

,

Proc. London Math. Soc.

s2-39

,

49

(

1935

).

https://doi.org/10.1112/plms/s2-39.1.49

Google Scholar

Crossref

Lattice Boltzmann method for warm fluid simulations of plasma wakefield acceleration

I. INTRODUCTION

II. LATTICE BOLTZMANN (LB) METHOD

A. Basic LB

B. Moment matching procedure for forced advection diffusion equation

III. RELATIVISTIC KINETIC EQUATIONS FOR COLLISIONLESS PLASMA

IV. LOCAL EQUILIBRIUM CLOSURE (LEC)

A. Local equilibrium closure lattice Boltzmann (LEC-LB)

V. WARM PLASMA CLOSURE (WARMC)

A. Warm plasma closure lattice Boltzmann (WARMC-LB)

VI. NUMERICAL RESULTS

A. 1D1V results

B. 1D3V results

C. 3D3V results

VII. OUTLOOK AND PERSPECTIVES

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

APPENDIX A: FULL EXPRESSIONS FOR THE FLUID ADVECTION EQUATIONS

1. LEC equations

2. WARMC equations

APPENDIX B: WARM LINEAR THEORY

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

Lattice Boltzmann method for warm fluid simulations of plasma wakefield acceleration

I. INTRODUCTION

II. LATTICE BOLTZMANN (LB) METHOD

A. Basic LB

B. Moment matching procedure for forced advection diffusion equation

III. RELATIVISTIC KINETIC EQUATIONS FOR COLLISIONLESS PLASMA

IV. LOCAL EQUILIBRIUM CLOSURE (LEC)

A. Local equilibrium closure lattice Boltzmann (LEC-LB)

V. WARM PLASMA CLOSURE (WARMC)

A. Warm plasma closure lattice Boltzmann (WARMC-LB)

VI. NUMERICAL RESULTS

A. 1D1V results

B. 1D3V results

C. 3D3V results

VII. OUTLOOK AND PERSPECTIVES

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

APPENDIX A: FULL EXPRESSIONS FOR THE FLUID ADVECTION EQUATIONS

1. LEC equations

2. WARMC equations

APPENDIX B: WARM LINEAR THEORY

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Related Content

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only