One of the simplest model that is currently applied in plasma physics simulations is the Vlasov-Poisson system. This system describes the evolution of charged particles in the case where the only interaction considered is the mean-field force created through electrostatic effects. The system consists in Vlasov equations for phase space density of each particle species of charge and mass
coupled to its self-consistent electric field which satisfies the Poisson equation
where is the vacuum permittivity. On the one hand, for a smooth and nonnegative initial data , the solution to (1.1) remains smooth and nonnegative for all . On the other hand, for any function , we have
which leads to the conservation of mass, norms, for and kinetic entropy,
We also get the conservation of momentum
and total energy
Numerical approximation of the Vlasov-Poisson system has been addressed since the sixties. Particle methods (PIC), consisting in approximating the plasma by a finite number of macro particles, have been widely used . They allow to obtain satisfying results with a few number of particles, but a well-known drawback of this class of methods is their inherent numerical noise which only decreases in when the number of particles increases, preventing from getting an accurate description of the distribution function for some specific applications. To overcome this difficulty, Eulerian solvers, that is methods discretizing the Vlasov equation on a mesh of the phase space, can be considered. Many authors explored their design, and an overview of many different techniques with their pros and cons can be found in . Among them, we can mention finite volume methods 
which are a simple and inexpensive option, but in general low order. Fourier-Fourier transform schemes are based on a Fast Fourier Transform of the distribution function in phase space, but suffer from Gibbs phenomena if other than periodic conditions are considered. Standard finite element methods [31, 32] have also been applied, but may present numerical oscillations when approximating the Vlasov equation. Later, semi-Lagrangian schemes have also been proposed 
, consisting in computing the distribution function at each grid point by following the characteristic curves backward. Despite these schemes can achieve high order allowing also for large time steps, they require high order interpolation to compute the origin of the characteristics, destroying the local character of the reconstruction.
In the present article, using Hermite polynomials in the velocity variable, we write the Vlasov-Poisson system (1.1)–(1.2) as an hyperbolic system. This idea of using Galerkin methods with a small finite set of orthogonal polynomials rather than discretizing the distribution function in velocity space goes back to the 60’s [1, 18]. More recently, the merit to use rescaled orthogonal basis like the so-called scaled Hermite basis has been shown [8, 16, 28, 26, 30]. In , Holloway formalized two possible approaches. The first one, called symmetrically-weighted (SW), is based on standard Hermite functions as the basis in velocity and as test functions in the Galerkin method. It appears that this SW method cannot simultaneously conserve mass and momentum. It makes up for this deficiency by correctly conserving the norm of the distribution function, ensuring the stability of the method. In the second approach, called asymmetrically-weighted (AW), another set of test functions is used, leading to the simultaneous conservation of mass, momentum and total energy. However, the AW Hermite method does not conserve the norm of the distribution function and is then not numerically stable. In addition, the asymmetric Hermite method exactly solves the spatially uniform problem without any truncation error, a property not shared by either the traditional symmetric Hermite expansion or by finite difference methods. The aim of this work is to present a class of numerical schemes based on the AW Hermite methods and to provide a stability analysis.
In what follows, we consider two types of spatial discretizations for the Vlasov-Poisson system (1.1)-(1.2), written as an hyperbolic system using Hermite polynomials in the velocity variable: spectral methods and a discontinuous Galerkin (DG) method. Concerning spectral methods, the Fourier basis is the natural choice for the spatial discretization when considering periodic boundary conditions. Spectral Galerkin and spectral collocation methods for the AW Fourier-Hermite discretization have been proposed in [8, 21, 24]. In , authors study a time implicit AW Fourier-Hermite method allowing exact conservation of charge, momentum and energy, and highlight that for some test cases, this scheme can be significantly more accurate than the PIC method.
For the SW Fourier-Hermite method, a convergence theory was proposed in . In , authors study conservation and stability properties of a generalized Hermite-Fourier semi-discretization, including as special cases the SW and AW approaches. Concerning discontinuous Galerkin methods, they are similar to finite elements methods but use discontinuous polynomials and are particularly well-adapted to handling complicated boundaries which may arise in many realistic applications. Due to their local construction, this type of methods provides good local conservation properties without sacrificing the order of accuracy. They were already used for the Vlasov-Poisson system in [14, 6]
. Optimal error estimates and study of the conservation properties of a family of semi-discrete DG schemes for the Vlasov-Poisson system with periodic boundary conditions have been proved for the one and multi-dimensional  cases. In all these works, the DG method is employed using a phase space mesh.
Here, we adopt this approach only in physical space, as in , with a Hermite approximation in the velocity variable. In , such schemes with discontinuous Galerkin spatial discretization are designed in such a way that conservation of mass, momentum and total energy is rigorously provable.
As mentionned before, the main difficulty is to study the stability of approximations based on asymmetrically-weighted Hermite basis. Indeed, this choice fails to preserve norm of the approximate solution, and therefore to ensure long-time stability of the method. In this framework, the natural space to be considered is
and there is no estimate of the associated norm for the solution to the Vlasov-Poisson system (1.1)-(1.2). To overcome this difficulty, we introduce a new weighted space, with a time-dependent weight, allowing to prove global stability of the solution in this space. Actually, this idea has been already employed in [22, 23] to stabilize Hermite spectral methods for linear diffusion equations and nonlinear convection-diffusion equations in unbounded domains, yielding stability and spectral convergence of the considered methods. Here, we define the weight as
where is a nonincreasing positive function which will be designed in such a way that a global stability estimate can be established in the following weighted space:
with the corresponding norm, that is
Applying now Young inequality on the first term, we get for ,
We notice that if is such that
the first and third terms cancel each other, yielding
Finally, applying Grönwall’s Lemma gives
It now remains to define satisfying (1.4). We remark that
is a suitable choice, where is the initial value
of at and corresponds to the initial scale of the
distribution function, whereas is a free parameter. This function is positive, nonincreasing and the parameter will practically be chosen sufficiently small to ensure that does not decrease too fast towards as goes to infinity.
Summarizing, we have established the following result.
Notice that the weight depends on the solution to the Vlasov-Poisson, hence it is mandatory to control the norm of in order that this estimate makes sense. In the next section, we introduce the formulation of the Vlasov equation using the Hermite basis in velocity, and prove the analogous of Proposition 1.1 for the obtained system. Then in Section 3, we discuss conservation and stability properties for a class of spatial discretizations including Fourier spectral method and discontinuous Galerkin approximations. In Section 4, we introduce the time discretization that will be used to compute numerical results with the discontinuous Galerkin method. Finally in Section 5 we present numerical results for two stream instability, bump-on-tail problem and ion acoustic wave to highlight conservations and stability of our discretization.
2. Hermite spectral form of the Vlasov equation
For simplicity, we now set all the physical constants to one and consider the one dimensional Vlasov-Poisson system for a single species with periodic boundary conditions in space, where the Vlasov equation reads
with , position and velocity . Moreover, the self-consistent electric field is determined by the Poisson equation
where the density is given by
and is a constant ensuring the quasi-neutrality condition of the plasma
2.1. Hermite basis
where is the number of modes. The issue is then to determine a well-suited class of basis functions and to find the expansion coefficients . Here, our aim is to obtain a control on a norm of , in the spirit of that established in Proposition 1.1. We then choose the following basis of normalized scaled time dependent asymmetrically weighted Hermite functions:
where is the positive nonincreasing function given by
with an approximation of the electric field that will be defined below. Functions are the Hermite polynomials defined by , and for , has the following recursive relation
Let us also emphasize that for all , and that the Hermite basis (2.4) satisfies the following orthogonality property
where is the Kronecker delta function.
Now the objective is to obtain an evolution equation for each mode by substituting the expression (2.3) for into the Vlasov equation (2.1) and using the orthogonality property (2.6). Thanks to the definition (2.4) of and the properties of , we compute the different terms of (2.1). The time derivative term is given by
the transport term is
and finally the nonlinear term is
Then taking as test function and using orthogonality property (2.6), we obtain an evolution equation for each mode , :
with the understanding that for and . Meanwhile, we first observe that the density satisfies
and then the Poisson equation becomes
In what follows, we rather consider a generic weak formulation of (2.7)–(2.8). Indeed, this will allow us to straightforwardly apply the obtained results to spatial discretizations with spectral methods (see Section 3.1).
Let be a subspace of . We look for such that and for all ,
and for a couple such that , and for all and , we have
2.2. Conservation properties
We consider the first three equations of (2.9): for all ,
which will be related to the conservation of mass, momentum and energy. Thus, let us look at the conservation properties from (2.11).
First the total mass is defined as
hence the conservation of mass directly comes from (2.11) by taking as test function in the first equation.
Then the momentum is given by
which taking in the second equation of (2.11) yields
Now, choosing in the second equation of (2.10) gives
which gives the conservation of momentum.
Finally to prove the conservation of total energy , we compute the variation of the kinetic energy defined as
Finally, taking in the first equation of (2.10), we obtain
We now compute the time derivative of the potential energy defined by
Using as test function in the first equation of (2.10) and an integration by parts, we have
This concludes the proof since the total energy is the sum of the kinetic and potential energies.
2.3. Weighted stability of
We now establish the discrete counterpart of Proposition 1.1, namely the stability in an appropriately weighted norm. More precisely, with defined by (2.5), we consider the weight , and denote by the corresponding weighted norm. Using the definition (2.3) of , we have
which gives, using orthogonality property (2.6), the following expression for the weighted norm of :
We compute the time derivative of defined by (2.13):
Thanks to equation (2.9) with , we then have to estimate
First of all, the transport term vanishes. Indeed, reindexing the sum over and using that , we have
Then on the one hand, using again that for all and , we have
On the other hand, we notice that
Gathering these two identities in (2.14), we get
Applying Young inequality in the second sum, this provides:
Yet, using definition (2.5) of , we have that
Using Grönwall’s lemma, this gives the expected result. ∎
2.4. Control of
Using the definition (2.13) of , the stability result established in Proposition 2.2 gives a stability result in for the coefficients , provided that is bounded from below by a positive constant. Due to definition (2.5) of , this is achieved as soon as is controlled. This result is given in the following proposition.
Since we are in one space dimension, by Sobolev inequality, there exists a constant such that for all ,
Moreover, since , Poincaré-Wirtinger inequality applies and then there exists such that
Taking in the second equation of (2.10) and applying Cauchy-Schwarz inequality gives
and then one obtains
Hence, we need to control