Stochastic variational principles for the collisional Vlasov–Maxwell and Vlasov–Poisson equations

In this work, we recast the collisional Vlasov–Maxwell and Vlasov–Poisson equations as systems of coupled stochastic and partial differential equations, and we derive stochastic variational principles which underlie such reformulations. We also propose a stochastic particle method for the collisional Vlasov–Maxwell equations and provide a variational characterization of it, which can be used as a basis for a further development of stochastic structure-preserving particle-in-cell integrators.


Introduction
The collisional Vlasov equation  probability density function, and therefore we will use the normalization f (x, v, t)d 3 vd 3 x = 1 instead. A self-consistent model of plasma is obtained by coupling (1.1) with the Maxwell equations ∇ x · E = ρ, (1.2a) and denote the charge density and the electric current density, respectively, and the factor N tot is due to our normalization. The system (1.1)-(1.3) is usually referred to as the Vlasov-Maxwell equations. It will also be convenient to express the electric and magnetic fields in terms of the scalar ϕ(x, t) and vector A(x, t) potentials (1.4b) as is typical in electrodynamics. The Vlasov-Poisson equations are an approximation of the Vlasov-Maxwell equations in the non-relativistic zero-magnetic field limit (see §6). The main goal of this work is to provide a variational characterization of the Vlasov-Maxwell and Vlasov-Poisson equations via a stochastic Lagrange-d'Alembert type of a principle. Variational principles have proved extremely useful in the study of nonlinear evolution partial differential equations (PDEs). For instance, they often provide physical insights into the problem being considered; facilitate discovery of conserved quantities by relating them to symmetries via Noether's theorem; allow one to determine approximate solutions to PDEs by minimizing the action functional over a class of test functions (e.g. [1]); and provide a way to construct a class of numerical methods called variational integrators [2,3]. A variational principle for the collisionless Vlasov-Maxwell equations was first proposed in [4]. It has been used to derive various particle discretizations of the Vlasov-Maxwell and Vlasov-Poisson equations [5][6][7][8][9], including structure-preserving variational particle-in-cell (PIC) methods [10][11][12]. It has also been applied to gyrokinetic theory (e.g. [13,14]). For other formulations and extensions, see also [15].
A structure-preserving description of collisional effects is far less developed. A metriplectic framework for the Vlasov-Maxwell-Landau equations has been presented in [16,17]. More recently, a stochastic variational principle has been proposed in [18] to describe collisional effects for the Vlasov equation with a fixed external electric field. To the best of our knowledge, to date no variational principle has been derived for the collisional Vlasov-Maxwell and Vlasov-Poisson equations. In this work, we extend the notion of the stochastic Lagrange-d'Alembert principle presented in [18] to plasmas evolving in self-consistent electromagnetic fields. The main idea of our approach is to interpret the Vlasov equation (1.1) as a Fokker-Planck equation and consider the associated stochastic differential equations.
Main content. The main content of the remainder of this paper is, as follows.
In §2, we recast the collisional Vlasov-Maxwell equations as a system of coupled stochastic and partial differential equations. In §3, we discuss the relationship between particle methods and stochastic modelling. We formulate a stochastic particle discretization for the collisional Vlasov-Maxwell equations and cast it in a form that allows the derivation of a variational principle. In §4, we describe the variational structure underlying the stochastic particle discretization of the Vlasov-Maxwell system. The main result of this section is theorem 4.2, in which a stochastic Lagrange-d'Alembert principle for the particle discretization is proved. In §5, we generalize the ideas from §4 to the original undiscretized equations. The main result of this section is theorem 5.1, in which a stochastic Lagrange-d'Alembert principle is proved for a class of the collisional Vlasov-Maxwell equations. In §6, we prove a stochastic Lagrange-d'Alembert principle applicable to the Vlasov-Poisson equations. The main result of this section is theorem 6.1. Section 7 contains the summary of our work.

The Vlasov-Maxwell-Fokker-Planck equations (a) Stochastic reformulation
Various collision models and various forms of the collision operator C[f ] are considered in the plasma physics literature (e.g. [50,51]). A key step towards a stochastic variational principle is a probabilistic interpretation of the Vlasov equation (1.1). Therefore, in this work we will be interested only in those collision operators for which (1.1) takes the form of a linear or strongly nonlinear Fokker-Planck equation (e.g. [52][53][54]). Namely, we will assume that the collision operator can be expressed as for some symmetric positive semi-definite matrix D ij (x, v; f ) and vector K i (x, v; f ) functions, where the dependence of D ij and K i on f may in general be nonlinear, and may involve differential and integral forms of f . In that case (1.1) is an integro-differential equation, the so-called strongly nonlinear Fokker-Planck equation [52]. In case D ij and K i are independent of f , that is, for a vector function G(x, v; f ), and a family of vector functions g ν (x, v; f ) with ν = 1, . . . , M. Note that given a symmetric positive semi-definite matrix D ij , a decomposition (2.2) can always be found, but it may not be unique. For instance, one may take M = 3 and assume that g i ν = g ν i for i, ν = 1, 2, 3. Then the first equation in (2.2) implies that the family of functions g i ν can be determined by calculating the square root of the matrix D ij , and the second equation in (2.2) can be used to calculate the function G. If (1.1) has the form of a Fokker-Planck equation, then the particle density function f can be interpreted as the probability density function for a stochastic process (X(t), V(t)) ∈ R 3 × R 3 . This stochastic process then satisfies the Stratonovich stochastic differential equation [52][53][54][55] where W 1 (t), . . . , W M (t) denote the components of the standard M-dimensional Wiener process and • denotes Stratonovich integration. Note that the terms G and g ν can be interpreted as external forces, and that in their absence the equations (2.3) reduce to the equations of motion of a charged particle in an electromagnetic field. We will therefore refer to G and g ν as forcing terms. The electric and magnetic fields are coupled via the Maxwell equations (1.2). It should also be noted that unless (1.1) is linear, the right-hand side of (2. 3) depends on f . In order to obtain a self-consistent system, one can express f in terms of the stochastic processes X and V where E denotes the expected value, and δ is Dirac's delta. This can be further plugged into (1.3). Together, we get Equations (1.2), (2.3) and (2.4) form a self-consistent system of stochastic and partial differential equations whose solutions are the stochastic processes X(t), V(t), and the functions E(x, t), B(x, t).
Remark. Upon substituting (2.4a), the forcing terms G and g ν become functionals of the processes X and V, that is, G(x, v; f ) = G(x, v; X, V) and g ν (x, v; f ) = g ν (x, v; X, V). However, for convenience and simplicity, throughout this work we will stick to the notation G(x, v; f ) and g ν (x, v; f ), understanding that the probability density is given by (2.4a) (or by (3.2a) for particle discretizations; see §3).

(b) Examples
Below we list a few examples of collision operators that fit the description presented in §2a.

(i) Lenard-Bernstein operator
The Lenard-Bernstein collision operator where ν c > 0, μ > 0 and γ > 0 are parameters, models small-angle collisions, and was originally used to study longitudinal plasma oscillations [50,51,56]. It can be easily verified that an example decomposition (2.2) for M = 3 is given by the functions

(ii) Lorentz operator
The Lorentz collision operator models electron-ion interactions via pitch-angle scattering and is given by the formula where ν c (|v|) is the collisional frequency as a function of the absolute value of velocity, I is the 3 × 3 identity matrix and ⊗ denotes tensor product. The primary effect of this type of scattering is a change of the direction of the electron's velocity with negligible energy loss. More information about the Lorentz collision operator, including the exact form of the collision frequency, can be found in, e.g. [50,51,57,58]. It can be verified by a straightforward calculation that an example decomposition (2.2) for M = 3 is given by the functions Note that these functions do not explicitly depend on f , therefore also in this case (1.1) is a linear Fokker-Planck equation.

(iii) Coulomb/Landau operator
The more general Coulomb collision operator has the form (2.1) with where N tot appears due to our normalization of f , δ ij is Kronecker's delta, and Γ = (4π q 4 /m 2 ) ln Λ, with ln Λ denoting the so-called Coulomb logarithm. The Coulomb operator describes collisions in which the fundamental two-body force obeys an inverse square law, and makes the assumption that small-angle collisions are more important than collisions resulting in large momentum changes [50,51,59]. A decomposition (2.2) can be found, for example, via the procedure outlined in §2a. However, the expressions for G and g ν are complicated, therefore we are not stating them here explicitly. Note that D ij and K i explicitly depend on f . Therefore, for the Coulomb operator the Vlasov equation (1.1) is a strongly nonlinear Fokker-Planck equation. Note also that D ij and K i can be explicitly written as functionals of the stochastic processes X and V as and The collision operator (2.1) with D ij and K i as in (2.9) can also be expressed in an equivalent, although more symmetric form, known as the Landau form of the Coulomb operator, or simply the Landau collision operator (e.g. [50]).

Stochastic particle discretization of the Vlasov-Maxwell equations
Particle modelling is one of the most popular numerical techniques for solving the Vlasov equation (e.g. [60,61]). In this section, we discuss the connections between particle methods and stochastic modelling. The standard particle method for the collisionless Vlasov equation for the particle density function, and deriving the corresponding ordinary differential equations satisfied by the 'particle' positions X a (t) and velocities V a (t), which turn out to be the characteristic equations. Note that we did a qualitatively similar thing in §2a, where we turned the original collisional Vlasov equation into the system of stochastic differential equations (2.3), which in the absence of the forcing terms G and g ν have the same form as the characteristic equations, and in fact the 'particles' X a (t) and V a (t) can be interpreted as realizations of the stochastic processes X(t) and V(t) for different elementary events ω ∈ Ω.
When the right-hand side of (2.3) does not depend on f , then (2.3) can in principle be solved numerically with the help of any standard stochastic numerical method (e.g. [55]), and each realization of the stochastic processes can be simulated independently of others. When the righthand side of (2.3) depends on f , then all realizations of the stochastic processes have to be solved for simultaneously, so that at each time step the probability density function f can be numerically approximated (e.g. [52]). Such an approach, however, does not quite lend itself to a geometric formulation. Therefore, in order to be able to introduce a variational principle in §4, let us consider 2N stochastic processes and Wiener processes. Note that the systems (3.1) are decoupled from each other for different values of a, and each system is driven by an independent Wiener process W a . Therefore, the pairs (X a , V a ) for a = 1, . . . , N are independent identically distributed stochastic processes, each with the probability density function f that satisfies the original Fokker-Planck equation (1.1). In that sense (3.1) is equivalent to (2.3). The advantage is that instead of considering N realizations of the six-dimensional stochastic process (X, V) in (2.3), one can consider one realization of the 6Ndimensional process (X 1 , V 1 , . . . , X N , V N ) in (3.1). Such a reformulation will allow us to identify an underlying stochastic variational principle in §4. The last step leading to the stochastic particle discretization is approximating the probability density function f in (3.1). This can be done with the help of the law of large numbers, namely, one can approximate (2.4) for large N as and It is easy to see that (3.2a) coincides with the standard Ansatz used in particle modelling (with the weights w a = 1/N). Therefore, the system of stochastic differential equations ( (3.2), and with the electromagnetic field coupled via the Maxwell equations (1.2), can be considered as a stochastic particle discretization of the collisional Vlasov-Maxwell equations.
Remark. Upon substituting (3.2a), the forcing terms G and g ν become functionals of the processes X 1 , . . . , X N and V 1 , . . . , V N . Similar to the discussion in §2a, for convenience and simplicity, throughout this work we will stick to the notation G(x, v; f ) and g ν (x, v; f ), understanding that the probability density is given by (3.2a) for particle discretizations.

Variational principle for the particle discretization
In this section, we propose an action functional which can be understood as a stochastic version of the Low action functional [4], and we prove a variational principle underlying the particle discretization introduced in §3, akin to the stochastic Lagrange-d'Alembert principle first introduced in [18].

(a) Function spaces
Before we introduce the action functional, we need to identify suitable function spaces on which it will be defined. For simplicity, let our spatial domain be the whole three-dimensional space R 3 , and let us consider the time interval The stochastic processes X a (t) and V a (t) satisfy (3.1), so they are in particular F t -adapted semimartingales, and have almost surely continuous paths [62]. We also notice that there is no diffusion term in (3.1), therefore we even have that the processes X a (t) are almost surely of class C 1 . We introduce the notation Note that this set is a vector space [62]. The potentials ϕ and A satisfy the Maxwell equations (1.2) and (1.4), therefore we require them to be of class C 2 . However, since our spatial domain is unbounded, we further need to assume that the vector fields E and B are square integrable. We introduce the notation and where X 0 (R n ) is simply the space of compactly supported elements of X(R n ).

(b) Action functional
Let us consider the action functional whereẊ a denotes the time derivative of X a , and the electric and magnetic fields E and B are expressed in terms of the partial derivatives of the potentials ϕ and A as in (1.4). Following the standard convention in stochastic analysis, we will omit writing elementary events ω ∈ Ω as arguments of stochastic processes unless otherwise needed, i.e. X a (t) ≡ X a (ω, t). The action functional (4.4) resembles the Low action functional introduced in [4]. In fact, it can be viewed as a particle discretization of the Low action functional, written in terms of stochastic processes [5,6,8,[10][11][12]. The term P a · (Ẋ a − V a ) is the so-called Hamilton-Pontryagin kinematic constraint (e.g. [63,64]) that enforces thatẊ a = V a using the Lagrange multiplier P a , which turns out to be the conjugate momentum. In principle, this constraint is not necessary in our context-we could omit it and replace V a withẊ a in (4.4). We will, however, keep it in order to make a clear connection with the theory developed in [35]. It also makes the notation in the proof of the stochastic Lagrange-d'Alembert principle in §4c more convenient and elegant. Note that the action functional S is itself a random variable, as ω ∈ Ω is one of its arguments. The variations of S with respect to its arguments are given by (see appendix A for the details of the derivations) and formulate a stochastic variational principle, we need the following lemma, whose proof is given in appendix B.
Remark. Equation (4.8) means that P(t), X(t) and V(t) satisfy a stochastic differential equation, which can be written in the differential form as (4.9) We are now in a position to formulate and prove a stochastic variational principle that generalizes the deterministic Lagrange-d'Alembert principle for forced Lagrangian and Hamiltonian systems, akin to the stochastic variational principle introduced in [18].
Proof. Let us first consider the variations with respect to A in (4.11). Given the boundary conditions for δA, from the standard calculus of variations we have that δ A S = 0 (see equation δϕ if and only if (1.2a) holds. Further, for variations with respect to V a we have that δ V a S = 0 (see equation (4.5b)) for all δV a if and only if (4.10b) is satisfied almost surely, which follows from the standard theorem of the calculus of variations, since the integral in (4.5b) is a standard Lebesgue integral, and the integrands are almost surely continuous. Similarly, δ P a S = 0 (see equation (4.5c)) for all δP a if and only if (4.10a) is satisfied almost surely. Finally, for variations with respect to X a , equations (4.5a) and (4.11) give Remark. Equation (4.10) is expressed in terms of the Lagrange multipliers P a , which, as can be seen in (4.10b), turn out to be the conjugate momenta. The conjugate momenta can be eliminated, and equation (4.10) can be recast as equation (3.1b), which is shown in the following theorem. Proof. By calculating the stochastic differential on both sides of (4.10b) and substituting (4.10a), we obtain Remark. Theorems 4.2 and 4.3 provide a variational formulation of the stochastic particle method from §3. One can further perform a variational discretization of the electromagnetic fields A and ϕ, for instance along the lines of [10,65] or [66], thus obtaining a stochastic PIC discretization of the collisional Vlasov-Maxwell equations. The resulting structure-preserving numerical methods will be investigated in a follow-up work.

Variational principle for the Vlasov-Maxwell equations
The form of the action functional (4.4) and of the Lagrange-d'Alembert principle (4.11) suggests that it should be possible to formulate a similar variational principle for the stochastic reformulation of the Vlasov-Maxwell system discussed in §2a. In this section, we provide such a variational principle for a class of collision operators. denotes the expected value of the random variable Y. Note that unlike S in (4.4), the action functionalS is not a random variable, as the dependence on ω ∈ Ω is integrated out with respect to the probability measure by calculating the expected value. In fact, S could be regarded as a Monte Carlo approximation ofS when the processes X 1 , . . . , X N are independent and identically distributed as X, and similarly for V and P. An important issue to consider is the domain of this action functional. In a similar manner to (4.3), one may want to take as the domain the set on which the formula (5.1) is well defined. This domain, however, turns out to be too big, in the sense that, as will be discussed below, due to the presence of the expected value the variations ofS do not uniquely determine the set of stochastic evolution equations (2.3). It is therefore necessary to restrict (5.2) to a smaller subspace or submanifold which is compatible with the considered collision operator. Below we will demonstrate how this can be done for a class of collision operators (2.1) for which D ij (x, v; f ) = const, that is, we have This class encompasses, for instance, the Lenard-Bernstein operator (2.6), or the more general nonlinear energy and momentum preserving Dougherty collision operator and its modifications [67][68][69][70][71][72][73][74]. For a given collision operator of the form (5.3), we define a compatible subset of C 0 Ω,T , namely, Note that for any P 1 , P 2 ∈ C col we have that d(P 1 − P 2 ) = (Z 1 − Z 2 ) dt, that is, P 1 − P 2 ∈ C 1 Ω,T . Therefore, the pair (C col , C 1 Ω,T ) is an affine subspace of C 0 Ω,T . The action functionalS can now be defined asS Similar to the calculations in §4b, the variations ofS with respect to V and P are given by, respectively, except that here δP ∈ C 1 Ω,T , so that P + δP ∈ C col . For the variation ofS with respect to X we have Since P ∈ C col , we have that dP = Z dt + m M ν=1 χ χ χ ν dW ν (t). Furthermore, the variations δX are almost surely of class C 1 , and therefore have sample paths of almost surely finite variation. Consequently, the quadratic covariation [χ χ χ ν · δX, W ν ] T 0 = 0 almost surely [62]. Since the expected By plugging this in (5.8), we finally obtain δ XS = N tot · E P(T) · δX(T) − P(0) · δX(0) The variations with respect to A and ϕ are the same as in (4.5d,e), respectively, only with the charge and electric current densities given by (2.4) rather than (3.2). The total variation ofS with respect to the variations of all arguments is given by
Proof. Similar to the proof of theorem 4.2, the equations δ ϕS = 0 and δ AS = 0 are equivalent to (1.2a) and (1.2d), respectively. Note that C 0 Ω,T is a subspace of L 2 (Ω × [0, T], R 3 ), and Y 1 , is an inner product on that space. Therefore, by substituting equations (5.6), (5.7) and (5.10) in equation (5.13), and using the fact that the variations are arbitrary, we establish equivalence with equations (5.12a,b), as well as with the equation for i = 1, 2, 3, which in turn is equivalent to equation (5.12c), given the assumption P ∈ C col . Proof. Similar to the proof of theorem 4.3, by calculating the stochastic differential on both sides of equation (5.12b) and comparing with equation (5.12c), one eliminates P and obtains equation (2.3b).
Remark. Note that the forcing terms g ν do not explicitly appear in the variational equation (5.13). By comparing theorem 4.2 and theorem 5.1, one could intuitively expect that the relevant variational principle should read However, due to the presence of the expected value in this equation, part or all of the information about the Stratonovich integral term is lost, as we saw in (5.9) for instance. Therefore, if the domain (5.2) is chosen forS, then the variational equations (5.13) or (5.15) do not determine a unique set of stochastic differential equations that need to be satisfied by the considered stochastic processes. Consequently, it is necessary to encode the missing information about the forcing terms g ν in the definition of the action functionalS by restricting its domain to a subset compatible with the considered collision operator. For the class of collision operators (5.3) a suitable choice of the domain is proposed in (5.5). For other collision operators appropriate domains will be nonlinear subspaces of (5.2), and they will be investigated in a follow-up work.

Variational principle for the Vlasov-Poisson equations
In the full Vlasov-Maxwell system, the scalar ϕ and vector A potentials are independent dynamic variables, and as such have to appear explicitly in the action functional alongside the stochastic processes X, V and P. In order to ensure the correct coupling between the stochastic processes and the electromagnetic field, an expected value was necessary in the definition of the action functional (5.1). This created a difficulty in deriving a variational principle, as pointed out in the remark following theorem 5.2. This difficulty can be circumvented for the Vlasov-Poisson equations because in this case the electrostatic potential ϕ is uniquely determined by the stochastic process X, as will be demonstrated below.

(a) The collisional Vlasov-Poisson equations
The collisional Vlasov-Poisson equations where E = −∇ x ϕ (6.2a) and x ϕ = −ρ, (6.2b) and the charge density ρ is given by (1.3), are an approximation of the Vlasov-Maxwell equations in the non-relativistic zero-magnetic field limit. The associated stochastic differential equations take the form dX = V dt (6.3a) and dV = q m E(X, t) + G(X, V; f ) dt + principle can be derived just like in § §3 and 4, respectively. Also, a variational principle analogous to the Lagrange-d'Alembert principle presented in §5 can be derived in a similar fashion. However, by doing so, one encounters the same difficulty with including the Stratonovich integral. In the case of the Vlasov-Poisson equations a different variational principle can be obtained by observing that the electrostatic potential ϕ can be expressed as a functional of the stochastic process X, ϕ : by solving Poisson's equation (6.2b). Given the charge density function (2.4b) and specific boundary conditions, the solution of Poisson's equation can be written using an appropriate Green's function for the Laplacian. Assuming the spatial domain is unbounded, the standard Green's function yields From (6.2a) we have the electric field (6.6)

(b) Action functional
Let us consider the action functional defined by the formulâ where the electrostatic potential ϕ is given by (6.5). Note that similar to S in (4.4), the functional S is itself random, and can be viewed as the action functional of particles represented by the process X which are moving in the electric field generated by particles represented by the process Y. Similar to the calculations in §4b, the variations ofŜ with respect to X, V and P are given by, respectively, where the electric field E is given by (6.6). Note that we are not considering variations with respect to Y. Let us for convenience define the joint variation ofŜ with respect to X, V and P as δ (X,V,P)Ŝ = δ XŜ + δ VŜ + δ PŜ . (6.10)

(c) The stochastic Lagrange-d'Alembert principle
In the following theorem, we formulate a variational principle for the system of equations (2.4b), (6.2) and (6.3). Note that E(X(t), t, X) is the electric field generated by a distribution of charged particles represented by the process X at time t, and evaluated at the random point x = X(t) in space. Furthermore, the notation δ XŜ [X, X, V, P] means that the variation ofŜ is evaluated for the arguments X, Y, V, P with Y = X. and V, P ∈ C 0 Ω,T be stochastic processes, and let ϕ(·, ·, X) ∈ X(R) be given by (6.5). Assume that G(·, ·; f ) and g ν (·, ·; f ) for ν = 1, . . . , M are C 1 functions of their arguments, where f is given by (2.4a). Then X, V and P satisfy the system of stochastic differential equationṡ X(t) = V(t), (6.11a) and dP(t) = qE X(t), t, X + m G X(t), V(t); f dt for arbitrary variations δX ∈ C 1 Ω,T , and δV, δP ∈ C 0 Ω,T , with δX(0) = δX(T) = 0 almost surely, where the action functionalŜ is given by (6.8).
Proof. Analogous to the proof of theorem 4.2.
Remark. It is straightforward to see that equations (6.11), together with (6.5) and (6.6), are equivalent to the system of equations (2.4b), (6.2) and (6.3). The Lagrange-d'Alembert principle (6.12) is unusual in that the variations of the action functionalŜ with respect to the argument Y are omitted. Thanks to such a form, however, the action functional does not require an expected value, and the collisional effects can be correctly included. A similar idea to solve Poisson's equation and plug the solution into the action functional was presented in [15], where the authors proposed a variational principle for the collisionless Vlasov-Poisson equations. In that approach the energy of the electric field was also included in the variational principle, and the variations were taken with respect to all arguments of the action functional. This approach could be adapted to the stochastic reformulation of the Vlasov-Poisson equations, but the corresponding action functional would have a form similar to (5.1), that is, it would need to contain an expected value, and therefore we would face a similar difficulty as for the Vlasov-Maxwell equations in §5b.

Summary and future work
In this work, we have considered novel stochastic formulations of the collisional Vlasov-Maxwell and Vlasov-Poisson equations, and we have identified new stochastic variational principles underlying these formulations. We have also proposed a stochastic particle method for the Vlasov-Maxwell equations and proved the corresponding stochastic variational principle.
Our work can be extended in several ways. The stochastic variational principle introduced in §4 can be used to construct stochastic variational PIC numerical algorithms for the collisional Vlasov-Maxwell and Vlasov-Poisson equations. Variational integrators are an important class of geometric integrators. This type of numerical scheme is based on discrete variational principles and provides a natural framework for the discretization of Lagrangian systems, including forced, dissipative or constrained ones. These methods have the advantage that they are symplectic when applied to systems without forcing, and in the presence of a symmetry, they satisfy a discrete version of Noether's theorem. For this reason, they demonstrate superior performance in long-time simulations; see [3,[75][76][77][78][79][80][81][82][83][84]. Variational integrators were introduced in the context of finite-dimensional mechanical systems, but were later generalized to Lagrangian field theories [ and applied in many computations, for example in elasticity, electrodynamics, fluid dynamics, or plasma physics; see [10][11][12]65,72,[85][86][87]. Stochastic variational integrators were first introduced in [35] and further studied in [18,34,37,39,40].
In §5, we have proposed a general action functional for the collisional Vlasov-Maxwell equations. However, we have also determined that in order to prove a relevant variational principle, the domain of this action functional has to be restricted in a way compatible with the collision operator of interest. We have shown that for a class of collision operators with constant diffusion terms, a suitable subdomain is an affine subspace (i.e. a submanifold). A natural continuation of our work would be to investigate submanifolds of (5.2) which are suitable for other collision operators.
Another aspect worth a more detailed investigation is the issue of existence and uniqueness of the solutions of the stochastic reformulations presented in this work, which are non-trivial systems of coupled stochastic and partial differential equations. This question is closely connected to the issue of existence and uniqueness of the solutions of the original collisional Vlasov-Maxwell and the Vlasov-Poisson equations. General results are available in the collisionless case (e.g. [88][89][90]), but the theory for the collisional equations is less developed (see [29,[91][92][93][94][95][96] and references therein).
Finally, as is typical for particle methods in general, the stochastic particle discretization proposed in §3 will require a large number of particles for accurate numerical simulations, which is computationally expensive. Structure-preserving model reduction methods [105,106] have been recently successfully applied to particle discretizations of the collisionless Vlasov equation [107]. It would be of great practical interest to combine our results with model reduction techniques in order to develop new efficient structure-preserving data-driven numerical methods for the collisional Vlasov-Maxwell equations. V j ∂A j ∂x i (X a , t)δX i a + P a · δẊ a dt.
Since δX a is almost surely differentiable, we have that its stochastic differential is simply dδX a = δẊ a dt. Furthermore, both δX a and P a are almost surely continuous semimartingales, therefore using the integration by parts formula for semimartingales [62] we can write T 0 P a · δẊ a dt = T 0 P a • dδX a = P a (t) · δX a (t) where the Stratonovich integrals are understood in the sense that δX a • dP a = i δX i a • dP i a . By substituting (A 3) in (A 2), we obtain (4.5a). Variations with respect to δV a , δP a ∈ C 0 Ω,T are defined analogously to (A 1). Similar computations (note that integration by parts is not necessary) yield (4.5b) and (4.5c), respectively.
The variation of S with respect to the variation δA ∈ X 0 (R 3 ) of the vector potential A is defined as Switching the order of differentiation and integration, integrating by parts, and using the fact that δA is compactly supported, one arrives at (4.5d), where in the derivations we have used (3.2c) and and the remaining calculations are standard, and can be found in, e.g. [108,109]. The variation of S with respect to the variation δϕ ∈ X 0 (R) of the scalar potential ϕ is defined in a similar fashion, and after similar calculations one obtains (4.5e).