Florent Renac
email: florent.renac@onera.fr

Entropy stable DGSEM for nonlinear hyperbolic systems in nonconservative form with application to two-phase flows

Keywords: 2000 MSC: 65M12, 65M70, 76T10 nonconservative hyperbolic systems, entropy stable schemes, discontinuous Galerkin method, summation-by-parts, two-phase flows

Entropy stable DGSEM for nonlinear hyperbolic systems in nonconservative form with application to two-phase flows

Florent Renac

Introduction

The discussion in this paper focuses on the high-order discretization of the Cauchy problem for nonlinear hyperbolic systems in nonconservative form:

∂ t u + A(u)∂ x u = 0, in R × (0, ∞), (1a) u(•, 0) = u 0 (•), in R, (1b)
where u(x, t) represents the vector of unknowns with values in the set of states Ω a ⊂ R m and A : Ω a u → A(u) ∈ R m×m is a smooth matrix-valued function. We assume that system (1a) is strictly hyperbolic over the set of states. If there exists a flux function f : Ω a → R m such that A(u) = f (u) for all u in Ω a , (1a) can be written in conservative form for which the concept of weak solutions in the sense of distributions is used to define admissible solutions.

In the general case where A is not the Jacobian of a flux function, the theory of distributions does not apply which makes difficult to give a meaning to the nonconservative product A(u)∂ x u at a point of discontinuity of the solution.

The work by Dal Maso, Lefloch, and Murat [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF] generalizes the notion of weak solutions from conservation laws to [START_REF] Abgrall | How to prevent pressure oscillations in multicomponent flow calculations: a quasi conservative approach[END_REF] and allows to define the nonconservative product for functions of bounded variations by extending the definition by Volpert [START_REF] Volpert | The space BV and quasilinear equations[END_REF]. The definition is based on a family of consistent and Lipschitz paths φ : [0, 1] × Ω a × Ω a → Ω a with bounded derivative ∂ s φ(s; •, •). We refer to [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF] for the complete theory and to [44, section 2] for an introduction. Across a discontinuity of speed σ, the nonconservative product A(u)∂ x u is then defined as the unique Borel measure defined by the so-called generalized Rankine-Hugoniot condition

σ[[u]] = 1 0 A φ(s; u -, u + ) ∂ s φ(s; u -, u + )ds, (2)
where [[u]] = u +u -, u -and u + are the left and right limits of u across the discontinuity. Note that the notion of weak solutions now depends on the family of paths under consideration [START_REF] Ph | Shock waves for nonlinear hyperbolic systems in nonconservative form[END_REF]Lemma 2.1].

Admissible weak solutions have to satisfy an entropy inequality

∂ t η(u) + ∂ x q(u) ≤ 0, (3)
for the smooth entropy-entropy flux pair (η, q) with η(•) a strictly convex function such that η (u) A(u) = q (u) for all u in Ω a . Many problems in science and engineering such as turbulence models, multiphase flows, shallow water flows, magnetohydrodynamics (MHD), etc. contain both conservative and nonconservative terms of the form

∂ t u + ∂ x f(u) + c(u)∂ x u = 0, (4)
which require different approaches for their discretizations. For smooth solutions we have A ≡ f + c and the entropy pair satisfies η (u) f (u) + c(u) = q (u) for u in Ω a . The objective of this work is to develop a general method to design arbitrary high-order schemes for (1) that satisfy the entropy inequality (3) at the semi-discrete level. We propose to use the discontinuous Galerkin spectral element method (DGSEM) based on the collocation between interpolation and quadrature points defined from Gauss-Lobatto quadrature rules [START_REF] Kopriva | On the quadrature and weak form choices in collocation type discontinuous Galerkin spectral element methods[END_REF]. Using diagonal norm summation-by-parts (SBP) operators and the entropy conservative numerical fluxes from Tadmor [START_REF] Tadmor | The numerical viscosity of entropy stable schemes for systems of conservation laws[END_REF], semi-discrete entropy conservative finite-difference and spectral collocation schemes have been derived in [START_REF] Fisher | High-order entropy stable finite difference schemes for nonlinear conservation laws: Finite domains[END_REF][START_REF] Carpenter | Entropy stable spectral collocation schemes for the Navier-Stokes equations: discontinuous interfaces[END_REF] for nonlinear conservation laws. An energy stable DGSEM for scalar nonlinear equations has been derived in [START_REF] Gassner | A Skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF] using the same framework and latter an entropy stable DGSEM for the compressible Euler equations on curvilinear hexahedral meshes has been proposed in [START_REF] Gassner | Split form nodal discontinuous Galerkin schemes with summation-by-parts property for the compressible Euler equations[END_REF]. The particular form of the SBP operators allows to take into account the numerical quadrature that approximates integrals in the numerical scheme compared to other techniques that require their exact evaluation to satisfy the entropy inequality [START_REF] Jiang | On a cell entropy inequality for discontinuous Galerkin methods[END_REF][START_REF] Hiltebrand | Entropy stable shock capturing spacetime discontinuous Galerkin schemes for systems of conservation laws[END_REF]. The work in [START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF] extends the DGSEM to triangular meshes and convection-diffusion equations. The DGSEM thus provides a general framework for the design of entropy conservative and entropy stable schemes for the discretization of nonlinear systems of conservation laws. Numerical experiments in [START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF] highlight the benefits on stability and robustness of the computations, though this not guarantees to preserve neither the entropy stability at the fully discrete level, nor positivity of the numerical solution which is necessary to define the entropy. Designs of fully discrete entropy stable and positive DGSEM have been proposed in [START_REF] Després | Entropy inequality for high order discontinuous Galerkin approximation of Euler equations[END_REF][START_REF] Després | Discontinuous Galerkin method for the numerical solution of Euler equations in axisymmetric geometry[END_REF][START_REF] Renac | A robust high-order Lagrange-projection like scheme with large time steps for the isentropic Euler equations[END_REF][START_REF] Renac | A robust high-order discontinuous Galerkin method with large time steps for the compressible Euler equations[END_REF]. A general framework for the design of entropy conservative and entropy stable schemes for steady-state conservation laws has been recently proposed in [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations[END_REF][START_REF] Abgrall | On the Connection between Residual Distribution Schemes and Flux Reconstruction[END_REF] that encompasses residual distribution schemes, discontinuous and continuous Galerkin methods with general quadrature formulas on simplex elements and flux reconstruction schemes on polygons. Some works rely on the discontinuous Galerkin (DG) approximation of nonconservative systems of shallow water flows [START_REF] Gassner | A well balanced and entropy conservative discontinuous Galerkin spectral element method for the shallow water equations[END_REF][START_REF] Dumbser | A staggered semi-implicit spectral discontinuous Galerkin scheme for the shallow water equations[END_REF][START_REF] Tassi | A discontinuous Galerkin finite element model for river bed evolution under shallow flows[END_REF][START_REF] Ranocha | Shallow water equations: split-form, entropy stable, well-balanced, and positivity preserving numerical methods[END_REF], MHD [START_REF] Liu | Entropy stable high order discontinuous Galerkin methods for ideal compressible MHD on structured meshes[END_REF][START_REF] Dumbser | A staggered semi-implicit spectral discontinuous Galerkin scheme for the shallow water equations[END_REF][START_REF] Bohm | An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations. Part I: Theory and Numerical Verification[END_REF], two-phase flows [START_REF] Van Zwieten | Efficient simulation of one-dimensional two-phase flow with a high-order h-adaptive space-time discontinuous Galerkin method[END_REF][START_REF] Franquet | Runge-Kutta discontinuous Galerkin method for the approximation of Baer and Nunziato type multiphase models[END_REF][START_REF] Franquet | Runge-Kutta discontinuous Galerkin method for reactive multiphase flows[END_REF][START_REF] Rhebergen | Discontinuous Galerkin finite element methods for hyperbolic nonconservative partial differential equations[END_REF][START_REF] Henry De Frahan | A new limiting procedure for discontinuous Galerkin methods applied to compressible multiphase flows with shocks and interfaces[END_REF][START_REF] Tokareva | HLLC-type Riemann solver for the Baer-Nunziato equations of compressible two-phase flow[END_REF][START_REF] Fraysse | Upwind methods for the Baer-Nunziato equations and higher-order reconstruction using artificial viscosity[END_REF][START_REF] Dumbser | A staggered semi-implicit spectral discontinuous Galerkin scheme for the shallow water equations[END_REF], etc. Note that the works in [START_REF] Gassner | A well balanced and entropy conservative discontinuous Galerkin spectral element method for the shallow water equations[END_REF] and [START_REF] Liu | Entropy stable high order discontinuous Galerkin methods for ideal compressible MHD on structured meshes[END_REF] use the DGSEM as discretization method and derive, respectively, high-order entropy conservative and well balanced discretization of the shallow water equations through skew-symmetric splitting techniques, and entropy stable schemes for the ideal compressible MHD equations by using two-point numerical fluxes from [START_REF] Chandrashekar | Entropy stable finite volume scheme for ideal compressible MHD on 2-D Cartesian meshes[END_REF] at element interfaces and treating the nonconservative product as source terms without particular treatment. Though not exhaustive, we also refer to the works in [START_REF] Abgrall | Numerical approximation of a compressible multiphase system[END_REF][START_REF] Dumbser | High-order unstructured Lagrangian one-step WENO finite volume schemes for non-conservative hyperbolic systems: applications to compressible multi-phase flows[END_REF][START_REF] Dumbser | FORCE schemes on unstructured meshes II: non-conservative hyperbolic systems[END_REF][START_REF] Dumbser | A simple robust and accurate a posteriori sub-cell finite volume limiter for the discontinuous Galerkin method on unstructured meshes[END_REF][START_REF] Tokareva | A flux splitting method for the Baer-Nunziato equations of compressible two-phase flow[END_REF]] and references therein as alternative techniques for high-order approximations of two-phase flows.

Here, we extend the works in [START_REF] Fisher | High-order entropy stable finite difference schemes for nonlinear conservation laws: Finite domains[END_REF][START_REF] Carpenter | Entropy stable spectral collocation schemes for the Navier-Stokes equations: discontinuous interfaces[END_REF][START_REF] Gassner | A Skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF][START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF] to nonconservative products by using the two-point entropy conservative numerical fluxes in fluctuation form introduced in [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF]. This extension is clarified through the direct link between fluctuation fluxes and conservative fluxes in the case of conservation laws. One difficulty in the design of an entropy stable DGSEM lies in the treatment of the integrals over discretization elements which contain space derivatives of test functions whose sign cannot be controlled. The use of entropy conservative numerical fluxes in those integrals allows however to remove their contribution to the global entropy production in the element. The properties of high-order accuracy and approximation of the cell averaged numerical solution are also difficult to derive due to the specific form of the fluctuation fluxes. Indeed, the consistency condition has less physical meaning for fluctuation fluxes compared to conservation fluxes which require homogeneity properties in closed form. Moreover, even in the case of path-conservative fluxes [START_REF] Parès | Numerical methods for non-conservative hyperbolic systems: a theoretical framework[END_REF] they require a priori knowledge of the underlying path. We thus introduce some assumptions on the form of the entropy conservative fluctuation fluxes and derive conditions on the scheme to keep high-order accuracy and a same semi-discrete scheme for the cell averaged approximate solution as in the original DGSEM. The method is fairly general and we provide examples of entropy conservative fluxes for nonconservative systems in various fields such as spray dynamics, gas dynamics, shallow water flows or two-phase flows. A deeper analysis is given for the discretization of the compressible Euler equations in nonconservative form [START_REF] Ch | A new comment on the computation of non-conservative products using Roe-type path conservative schemes[END_REF] and for two two-phase flow models in one space-dimension: a 2×2 system with a nonconservative product associated to a linearlydegenerate (LD) characteristic field, and the isentropic Baer-Nunziato model. We provide numerical examples where the original DGSEM fails to capture the entropy weak solution. The use of an entropy stable DGSEM scheme is here necessary to capture the correct solution and improve robustness of the computations. For the Baer-Nunziato model, we further analyze the properties of the discrete scheme and derive conditions on the time step to keep positivity of the partial densities and a maximum principle on the void fractions. These properties hold for the cell averaged numerical solution and motivate the use of a posteriori limiters [START_REF] Zhang | On maximum-principle-satisfying high order schemes for scalar conservation laws[END_REF][START_REF] Zhang | On positivity-preserving high order discontinuous Galerkin schemes for compressible Euler equations on rectangular meshes[END_REF] to extend them to nodal values within elements. Again, numerical experiments highlight stability and robustness improvement with the entropy stable scheme.

The paper is organized as follows. Section 2 presents the DGSEM for the space discretization of nonconservative systems [START_REF] Abgrall | How to prevent pressure oscillations in multicomponent flow calculations: a quasi conservative approach[END_REF] and its entropy stable version through the use of entropy conservative and entropy stable numerical fluxes. In section 3, we derive the semi-discrete entropy inequality and give conditions on the numerical fluxes to keep highorder accuracy and the semi-discrete scheme for the cell averaged numerical solution. Four examples of entropy conservative fluxes are given in section 4 for different nonconservative systems and some other examples are given in Appendix A. We further investigate the stability and robustness properties of an entropy stable DGSEM for the isentropic Baer-Nunziato model in section 5. Numerical experiments for the four examples of section 4 are given in section 6. Finally, concluding remarks about this work are given in section 7.

DGSEM formulation

The DG method consists in defining a semi-discrete weak formulation of problem [START_REF] Abgrall | How to prevent pressure oscillations in multicomponent flow calculations: a quasi conservative approach[END_REF]. The domain is discretized with a grid Ω h = ∪ j∈Z κ j with cells κ j = [x j-1 2 , x j+ 1 2 ], x j+ 1 2 = jh and h > 0 the space step (see Figure 1) that we assume to be uniform without loss of generality.

Numerical solution

We look for approximate solutions in the function space of discontinuous polynomials

V p h = {v h ∈ L 2 (Ω h ) : v h | κ j ∈ P p (κ j ), κ j ∈ Ω h }
, where P p (κ j ) denotes the space of polynomials of degree at most p in the element κ j . The approximate solution to (1) is sought under the form

u h (x, t) = p l=0 ϕ l j (x)U l j (t), ∀x ∈ κ j , κ j ∈ Ω h , t ≥ 0, (5)
where U 0≤l≤p j are the degrees of freedom (DOFs) in the element κ j . The subset (ϕ 0 j , . . . , ϕ p j ) constitutes a basis of V p h restricted onto a given element. In this work we will use the Lagrange interpolation polynomials 0≤k≤p associated to the Gauss-Lobatto nodes over the segment [-1, 1]:

s 0 = -1 < s 1 < • • • < s p = 1: k (s l ) = δ k,l , 0 ≤ k, l ≤ p, (6)
with δ k,l the Kronecker symbol. The basis functions with support in a given element κ j thus write ϕ k j (x) = k (σ j (x)) where σ j (x) = 2(xx j )/h and x j = (x j+ 1 2 + x j-1 2 )/2 denotes the center of the element. 3 The DOFs thus correspond to the point values of the solution: given 0 ≤ k ≤ p, j in Z, and t ≥ 0, we have u h (x k j , t) = U k j (t) for x k j = x j + s k h/2. The left and right traces of the numerical solution at interfaces x j± 1 2 of a given element hence read (see Figure 1):

u - j+ 1 2 (t) := u h (x - j+ 1 2 , t) = U p j (t), ∀t ≥ 0, (7a)
u + j-1 2 (t) := u h (x + j-1 2 , t) = U 0 j (t), ∀t ≥ 0. ( 7b
)
It is convenient to introduce the discrete derivative matrix with entries

D kl = l (s k ) = h 2 d x ϕ l j (x k j ), 0 ≤ k, l ≤ p. (8)
In the DGSEM, the integrals over elements are approximated by using a Gauss-Lobatto quadrature rule with nodes collocated with the interpolation points of the numerical solution

κ j f (x)dx h 2 p l=0 ω l f (x l j ), (9)
with ω l > 0, x l j = x j + s l h/2 the weights and nodes of the quadrature rule which impose p l=0 ω l = 1 -1 ds = 2, and s l the node positions in [-1, 1]. This leads to the definition of the discrete inner product in the element κ j f, g p j :=

h 2 p l=0 ω l f (x l j )g(x l j ).
As noticed in [START_REF] Gassner | A Skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF], the DGSEM satisfies the summation-by-parts property:

ω k D kl + ω l D lk = δ kl (δ kp -δ k0 ), 0 ≤ k, l ≤ p. ( 10
)
Note also that the property p l=0 l ≡ 1 implies

p l=0 D kl = 0, 0 ≤ k ≤ p. (11)

Space discretization

The semi-discrete form of the DG discretization in space of problem (1) reads [START_REF] Franquet | Runge-Kutta discontinuous Galerkin method for the approximation of Baer and Nunziato type multiphase models[END_REF][START_REF] Rhebergen | Discontinuous Galerkin finite element methods for hyperbolic nonconservative partial differential equations[END_REF]: find u h in (V p h ) m such that

Ω h v h ∂ t u h dx + κ j ∈Ω h κ j v h A(u h )∂ x u h dx + j∈Z v - j+ 1 2 D -u - j+ 1 2 (t), u + j+ 1 2 (t) + j∈Z v + j-1 2 D + u - j-1 2 (t), u + j-1 2 (t) = 0, ∀v h ∈ V p h , t > 0, (12)
where the numerical fluxes D ± (•, •) in fluctuation form will be defined below.

The projection of the initial condition (1b) onto (V p h ) m reads

Ω h v h (x)u h (x, 0)dx = Ω h v h (x)u 0 (x)dx, ∀v h ∈ V p h .
Substituting v h for the Lagrange interpolation polynomials (6) and using the Gauss-Lobatto quadrature (9) to approximate the volume integrals, [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF] becomes

ω k h 2 dU k j dt + ω k A(U k j ) p l=0 U l j D kl + δ kp D -(U p j , U 0 j+1 ) + δ k0 D + (U p j-1 , U 0 j ) = 0, ∀ j ∈ Z, 0 ≤ k ≤ p, t > 0. ( 13
)
In section 2.3, we propose to modify the volume integrals in (13) so as to satisfy an entropy balance. Note that (13) satisfies some conservation property:

h d u j dt + A(u h ), d x u h p j + D -(U p j , U 0 j+1 ) + D + (U p j-1 , U 0 j ) = 0, (14)
for the cell averaged solution

u j (t) := 1 h κ j u h (x, t)dx = 1 2 p k=0 ω k U k j (t).
The numerical fluxes in fluctuation form satisfy the following consistency property

D ± (u, u) = 0, ∀u ∈ Ω a , (15)
and may also satisfy the path-conservative property [START_REF] Parès | Numerical methods for non-conservative hyperbolic systems: a theoretical framework[END_REF] D

-(u -, u + ) + D + (u -, u + ) = 1 0 A φ(s; u -, u + ) ∂ s φ(s; u -, u + )ds, (16)
for a given path φ.

Entropy stable numerical fluxes

In the following, we use the usual terminology from [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF] and denote by entropy conservative for the entropyentropy flux pair (η, q) in (3), the numerical fluxes D ± ec satisfying:

η (u -) D - ec (u -, u + ) + η (u + ) D + ec (u -, u + ) = q(u + ) -q(u -), ∀u ± ∈ Ω a . (17)
Furthermore, we will assume that the numerical fluxes at interfaces in [START_REF] Castro | Why many theories of shock waves are necessary: Convergence error in formally path-consistent schemes[END_REF] are entropy stable in the following sense:

η (u -) D -(u -, u + ) + η (u + ) D + (u -, u + ) ≥ q(u + ) -q(u -), ∀u ± ∈ Ω a . ( 18
)
As done in [START_REF] Gassner | A Skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF][START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF] for hyperbolic conservation laws, we modify the volume integrals in [START_REF] Castro | Why many theories of shock waves are necessary: Convergence error in formally path-consistent schemes[END_REF] to satisfy the entropy inequality at the semi-discrete level. The semi-discrete scheme now reads

ω k h 2 dU k j dt + R k j (u h ) = 0, ∀ j ∈ Z, 0 ≤ k ≤ p, t > 0, (19)
with

R k j (u h ) = ω k p l=0 D(U k j , U l j )D kl + δ kp D -(U p j , U 0 j+1 ) + δ k0 D + (U p j-1 , U 0 j ), (20)
and

D(u -, u + ) := D - ec (u -, u + ) -D + ec (u + , u -), ∀u ± ∈ Ω a , (21)
where D ± ec (•, •) are some entropy conservative fluctuation fluxes [START_REF] Cockburn | Runge-Kutta discontinuous Galerkin methods for convection-dominated problems[END_REF]. The reason of using D(U k j , U l j ) instead of A(U k j )U l j in (20) will be made clear below in equation ( 25) and constitutes the main ingredient for entropy stability.

3. Properties of the semi-discrete scheme 3.1. Entropy stable scheme Theorem 3.1 proves a semi-discrete entropy inequality for the scheme [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] together with entropy stable fluxes at interfaces, while Theorem 3.2 establishes high-order accuracy and the preservation of equation ( 14) for the cell averaged solution.

Theorem 3.1 (entropy stable DGSEM). Let D(•, •) defined in [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF] with D ± ec (•, •) consistent ( 15) and entropy conservative [START_REF] Cockburn | Runge-Kutta discontinuous Galerkin methods for convection-dominated problems[END_REF] fluctuation fluxes, and let D ± (•, •) be consistent [START_REF] Chandrashekar | Entropy stable finite volume scheme for ideal compressible MHD on 2-D Cartesian meshes[END_REF] and entropy stable (18) fluctuation fluxes. Then, the semi-discrete DGSEM [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] satisfies the following entropy inequality for the pair (η, q) in (3)

h d η j dt + Q(U p j , U 0 j+1 ) -Q(U p j-1 , U 0 j ) ≤ 0, ( 22
)
with η j = p k=0 ω k 2 η(U k j ) and either Q(U p j , U 0 j+1 ) = q(U p j ) + η (U p j ) D -(U p j , U 0 j+1 ), (23) or
Q(U p j , U 0 j+1 ) = q(U 0 j+1 ) -η (U 0 j+1 ) D + (U p j , U 0 j+1 ). (24)
Proof. Left multiplying [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] with η (U k j ) and adding up over 0 ≤ k ≤ p, we obtain

h d η j dt + k,l ω k η (U k j ) D(U k j , U l j )D kl + η (U p j ) D -(U p j , U 0 j+1 ) + η (U 0 j ) D + (U p j-1 , U 0 j ) = 0,
where the second term may be transformed into

k,l ω k η (U k j ) D(U k j , U l j )D kl (21) = k,l ω k η (U k j ) D - ec (U k j , U l j ) -D + ec (U l j , U k j ) D kl (10) = k,l ω k η (U k j ) D - ec (U k j , U l j )D kl + ω l η (U k j ) D + ec (U l j , U k j )D lk -δ kl (δ kp -δ k0 )η (U k j ) D + ec (U l j , U k j ) (15) = k↔l k,l ω k η (U k j ) D - ec (U k j , U l j ) + η (U l j ) D + ec (U k j , U l j ) D kl (17) = k,l ω k q(U l j ) -q(U k j ) D kl (11) = k,l ω k q(U l j )D kl (10) = (11) q(U p j ) -q(U 0 j ), (25)
where k ↔ l indicates an inversion of indices k and l in some of the terms. We thus obtain

h d η j dt + q(U p j ) -q(U 0 j ) + η (U p j ) D -(U p j , U 0 j+1 ) + η (U 0 j ) D + (U p j-1 , U 0 j ) = 0,
and using ( 23) we get [START_REF] Després | Discontinuous Galerkin method for the numerical solution of Euler equations in axisymmetric geometry[END_REF] to the right-hand-side, we finally deduce

h d η j dt + Q(U p j , U 0 j+1 ) -q(U 0 j ) + η (U 0 j ) D + (U p j-1 , U 0 j ) = 0, and adding Q(U p j-1 , U 0 j ) -q(U p j-1 ) -η (U p j-1 ) D -(U p j-1 , U 0 j ) = 0 from
h d η j dt + Q(U p j , U 0 j+1 ) -Q(U p j-1 , U 0 j ) = q(U 0 j ) -η (U 0 j ) D + (U p j-1 , U 0 j ) -q(U p j-1 ) -η (U p j-1 ) D -(U p j-1 , U 0 j ) (18)
≤ 0.

A similar result holds with [START_REF] Dumbser | High-order unstructured Lagrangian one-step WENO finite volume schemes for non-conservative hyperbolic systems: applications to compressible multi-phase flows[END_REF].

Remark 3.1. Entropy fluxes ( 23) and ( 24) are different if ( 18) is satisfied with a strict inequality in contrast to entropy conservative DGSEM, see Corollary 3.1. We stress that both hold for a given choice of numerical fluxes D ± (•, •), but lead to different expressions of the entropy dissipation rate evaluated at the left face of the cell with (23) (see equation ( 26)) or at the right face of the cell with [START_REF] Dumbser | High-order unstructured Lagrangian one-step WENO finite volume schemes for non-conservative hyperbolic systems: applications to compressible multi-phase flows[END_REF]. We do not know a priori if either [START_REF] Després | Discontinuous Galerkin method for the numerical solution of Euler equations in axisymmetric geometry[END_REF] or ( 24) leads to a stronger entropy dissipation and this should be application dependent. Equation ( 26) is also used in section 4.3 for discussing the entropy dissipation rate at the discrete level for the Euler equations [START_REF] Jiang | On a cell entropy inequality for discontinuous Galerkin methods[END_REF].

Entropy conservation then results as an immediate consequence.

Corollary 3.1 (entropy conservative fluxes). Under the assumptions of Theorem 3.1, the semi-discrete DGSEM [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] is entropy conservative iff. the numerical fluxes at interfaces are entropy conservative [START_REF] Cockburn | Runge-Kutta discontinuous Galerkin methods for convection-dominated problems[END_REF]. The numerical entropy flux reads

Q(U p j , U 0 j+1 ) = q(U p j ) + η (U p j ) D - ec (U p j , U 0 j+1 ), = q(U 0 j+1 ) -η (U 0 j+1 ) D + ec (U p j , U 0 j+1 ).
High-order accuracy and the conservation-like property ( 14) require further assumptions on the form of the entropy conservative fluxes [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF] which are summarized in Theorem 3.2 below. We stress that this form of fluctuation fluxes is fairly general and includes for instance skew-symmetric splittings (see Corollary 3.2). Theorem 3.2. Under the assumptions of Theorem 3.1 and further assuming that the entropy conservative fluctuation fluxes [START_REF] Cockburn | Runge-Kutta discontinuous Galerkin methods for convection-dominated problems[END_REF] have the following form

D ± ec (u -, u + ) = A ± (u -, u + )[[u]], (27a) A(u -, u + ) := A -(u -, u + ) + A + (u -, u + ), (27b) A(u -, u + ) + A(u + , u -) = A(u -) + A(u + ), (27c) A(u, u) = A(u), (27d)
for all u ± and u in Ω a , where

[[u]] = u + -u -.
Then, the semi-discrete DGSEM ( 19) is a high-order approximation in space of smooth solutions of the nonconservative equation (1a) which satisfies [START_REF] Ch | A new comment on the computation of non-conservative products using Roe-type path conservative schemes[END_REF].

Proof. First, to prove accuracy, it is sufficient to prove that the volume integral in ( 19) is a high-order approximation of

A(u)∂ x u at points x k j , 0 ≤ k ≤ p, for smooth enough solutions u. Let π p h : L 2 (Ω h ) ∩ C(Ω h ) u → π p h (u) ∈ V p h be the Lagrange projection onto V p
h associated to nodes [START_REF] Ambroso | The drift-flux asymptotic limit of barotropic two-phase two-pressure models[END_REF]. Since the Lagrange interpolation error is of order O(h p+1 ), we have for u and v in C p+1 (Ω h ):

d x π p h (uv)(x) = u(x)d x v(x) + v(x)d x u(x) + O(h p ), ∀x ∈ Ω h . (28)
Let t > 0, introducing the interpolation polynomial

a k h (x) := p l=0 A -(U k j , U l j )ϕ l j (x), we have a k h (x k j ) = A -(U k j , U k j ) and d x a k h (x k j ) = p l=0 A -(U k j , U l j )d x ϕ l j (x k j )
. Using [START_REF] Fisher | High-order entropy stable finite difference schemes for nonlinear conservation laws: Finite domains[END_REF] for the product a k h u h at point x k j , we obtain

2 h p l=0 A -(U k j , U l j )U l j D kl = A -(U k j , U k j )d x u h (x k j ) + 2 h p l=0 A -(U k j , U l j )U k j D kl + O(h p ).
Applying the same rule for a k h (x) := p l=0 A + (U l j , U k j )ϕ l j (x), we finally obtain

2 h p l=0 A -(U k j , U l j )(U l j -U k j )D kl = A -(U k j , U k j )d x u h (x k j ) + O(h p ), (29a) 2 h
p l=0 A + (U l j , U k j )(U l j -U k j )D kl = A + (U k j , U k j )d x u h (x k j ) + O(h p ). ( 29b
)
We thus have

2 h p l=0 D(U k j , U l j )D kl (27a) = (21) 2 h p l=0 A -(U k j , U l j ) + A + (U l j , U k j ) (U l j -U k j )D kl (29) = A -(U k j , U k j ) + A + (U k j , U k j ) d x u h (x k j ) + O(h p ) (27b) = A(U k j , U k j )d x u h (x k j ) + O(h p ) (27d) = A(U k j )d x u h (x k j ) + O(h p ).
Then, to obtain ( 14), we add up (19) over 0 ≤ k ≤ p and obtain

h d u j dt + p k=0 ω k p l=0 D(U k j , U l j )D kl + D -(U p j , U 0 j+1 ) + D + (U p j-1 , U 0 j ) = 0,
where the second term may be transformed into

k,l ω k D(U k j , U l j )D kl (21) = k,l ω k D - ec (U k j , U l j ) -D + ec (U l j , U k j ) D kl (10) = k,l ω k D - ec (U k j , U l j )D kl + ω l D + ec (U l j , U k j )D lk - δ kl (δ kp -δ k0 )D + ec (U l j , U k j ) (15) = k↔l k,l ω k D - ec (U k j , U l j ) + D + ec (U k j , U l j ) D kl (27a,b) = k,l ω k A(U k j , U l j )(U l j -U k j )D kl (10) = (27d) k,l ω k A(U k j , U l j )U l j D kl + ω l A(U k j , U l j )U k j D lk - A(U p j )U p j + A(U 0 j )U 0 j (27c) = k↔l k,l ω k A(U k j ) + A(U l j ) U l j D kl -A(U p j )U p j + A(U 0 j )U 0 j (10) = k,l ω k A(U k j )U l j D kl , = A(u h ), d x u h p j ,
which completes the proof. Now, we consider splittings of the nonconservative product for smooth solutions of the form

A∂ x u = αA∂ x u + (1 -α) ∂ x (Au) -(∂ x A)u , 0 ≤ α ≤ 1. ( 30
)
Entropy stable schemes based on the above decomposition fall into the assumptions of Theorem 3.2 as stated below.

Corollary 3.2 (skew-symmetric splitting). Entropy conservative fluxes in [START_REF] Cockburn | Runge-Kutta discontinuous Galerkin methods for convection-dominated problems[END_REF] for the splitting (30) read

D ± ec (u -, u + ) = A ± (u -, u + )[[u]], A ± (u -, u + ) = 1 2 αA(u ± ) + (1 -α)A(u ∓ ) , (31)
and constitute particular cases of the high-order entropy conservative fluxes (27) of Theorem 3.2.

Proof. First, using ( 31) and ( 21) to evaluate the volume integral in [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF], we obtain

D(u -, u + ) = αA(u -) + (1 -α)A(u + ) [[u]] = αA(u -)(u + -u -) + (1 -α)A(u + )(u + -u -) = αA(u -)u + + (1 -α) A(u + )u + -A(u + )u -,
since from (11) the term A(u -)u -has no contribution to the volume integral in [START_REF] Coquel | Adjoint approximation of nonlinear hyperbolic systems with non-conservative products[END_REF], i.e., ω k p l=0 A(U k j )U k j D kl = 0. The above relation implies that ( 21) is a volume discretization of the right-hand-side of [START_REF] Franquet | Runge-Kutta discontinuous Galerkin method for the approximation of Baer and Nunziato type multiphase models[END_REF].

Then, from (27b) we have

A(u -, u + ) = A -(u -, u + ) + A + (u -, u + ) = 1 2 A(u -) + A(u + ) , which indeed satisfies (27c,d).

Entropy conservative fluxes for conservation laws

In the particular case where (1) reduces to a conservation law, i.e., A(u) = f (u), it has been shown in [START_REF] Fisher | High-order entropy stable finite difference schemes for nonlinear conservation laws: Finite domains[END_REF][START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF] that it is possible to satisfy the entropy inequality ( 22) by using the entropy conservative fluxes h ec (u -, u + ) from Tadmor [START_REF] Tadmor | The numerical viscosity of entropy stable schemes for systems of conservation laws[END_REF] which satisfy

[[η ]] h ec (u -, u + ) = [[η f -q]], ∀u ± ∈ Ω a , h ec (u, u) = f(u), ∀u ∈ Ω a .
The link between the fluctuation fluxes D ± ec (u -, u + ) and the conservative flux h ec (u -, u + ) reads

h ec (u -, u + ) = f(u -) + D - ec (u -, u + ) = f(u + ) -D + ec (u -, u + ),
from which we deduce that

D - ec (u -, u + ) -D + ec (u + , u -) = h ec (u -, u + ) + h ec (u + , u -) -2f(u -), (32)
and using [START_REF] Carpenter | Entropy stable spectral collocation schemes for the Navier-Stokes equations: discontinuous interfaces[END_REF] the volume integral in [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] becomes

ω k p l=0 D(U k j , U l j )D kl = ω k p l=0 h ec (U k j , U l j ) + h ec (U l j , U k j ) D kl . (33)
In [START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF], a slightly different choice has been made: D(u -, u + ) := 2h ec (u -, u + ), where h ec (•, •) is assumed to be symmetric. In fact, it may be easily verified that the properties of Theorem 3.3 in [START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF] also hold with [START_REF] Gassner | A Skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF] which may be seen as a generalization of the framework of entropy stable DGSEM to nonsymmetric entropy conservative fluxes by using the symmetrizer (h ec (u -, u + ) + h ec (u + , u -))/2.

Examples

In this section we consider four different nonconservative scalar equations and systems in one space dimension and provide examples of entropy conservative and entropy stable numerical fluxes that fall into the category considered in Theorems 3.1 and 3.2 in section 2. [START_REF] Abgrall | On the Connection between Residual Distribution Schemes and Flux Reconstruction[END_REF]. These examples will be used in the numerical experiments of section 6 and some other examples are given in Appendix A for gas dynamics, spray dynamics, shallow water flows, etc. In the following, it is convenient to introduce the average operator u := u -+u + 2 .

Burgers' equation

The inviscid Burgers' equation in nonconservative form reads

∂ t u + u∂ x u = 0, (34)
with entropy η(u) = u 2 2 and entropy flux q(u) = u 3 3 . Entropy conservative fluctuation fluxes of the form ( 27) read

D - ec (u -, u + ) = 2u -+ u + 6 [[u]], D + ec (u -, u + ) = u -+ 2u + 6 [[u]]. (35)
Using [START_REF] Fraysse | Upwind methods for the Baer-Nunziato equations and higher-order reconstruction using artificial viscosity[END_REF], with f (u) = u 2 2 , and looking for an equivalent symmetric entropy conservative flux for conservative equations, we obtain

h ec (u -, u + ) = D - ec (u -, u + ) -D + ec (u -, u + ) + 2 f (u -) 2 = (u -) 2 + u -u + + (u + ) 2 6 ,
which corresponds to the entropy conservative skew-symmetric splitting of the Burgers' equation [START_REF] Tadmor | Skew-selfadjoint form for systems of conservation law[END_REF][START_REF] Gassner | A Skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF]. Note that this skew-symmetric splitting already consists in a nonconservative formulation of the space derivative

∂ x u 2 2 = 2 3 ∂ x u 2 2 + 1 3 u∂ x u [33].

Nonconservative product associated to a LD field

Let us introduce the following nonlinear hyperbolic system representative of two-phase flow problems [START_REF] Coquel | Adjoint approximation of nonlinear hyperbolic systems with non-conservative products[END_REF] where the LD characteristic field plays the role of interface velocity [START_REF] Coquel | Closure laws for a two-fluid two-pressure model[END_REF]:

∂ t u + g(u)∂ x u = 0, (36a) ∂ t v + ∂ x f (u) = 0, (36b) with g(u) = u + v and f (u) = v 2 -u 2 2 .
The eigenvalues are g(u) associated to the LD field and v associated to a genuinely nonlinear field so the system is strictly hyperbolic over the set of states Ω a = {(u, v) ∈ R 2 : u > 0}. It satisfies an entropy inequality for the pair η(u) = (u+v) 2 2 and q(u) = (u+v) 3 3 . Entropy conservative fluctuation fluxes are

D - ec (u -, u + ) = 1 6 2g(u -) + g(u + ) [[u]] (2v -+ v + )[[v]] -(2u -+ u + )[[u]] , (37a)
D + ec (u -, u + ) = 1 6 g(u -) + 2g(u + ) [[u]] (v -+ 2v + )[[v]] -(u -+ 2u + )[[u]] . (37b)
Note that the regularized system

∂ t u + g(u)∂ x u = ∂ 2 xx u, ∂ t v + ∂ x f (u) = ∂ 2 xx v, with > 0 gives ∂ t η(u) + ∂ x q(u) -∂ 2 xx η(u) = -(∂ x u) 2 + (∂ x v) 2 ≤
0, so the associated viscous profiles will give the physically admissible solutions in the limit = 0 + . Using this result for numerical purposes, we design the following entropy stable flux

D ± (u -, u + ) = 2g(u ± )+g(u ∓ ) 6 [[u]] ± f (u ± ) -ĥ(u -, u + ) ± v [[u]], ĥ(u -, u + ) = f (u -)+ f (u + ) 2 -β s 2 [[v]], (38)
with numerical parameters v ≥ 0 and β s ≥ 0. Setting v = 0, it may be checked that the fluctuations fluxes in [START_REF] Ismail | Affordable, entropy-consistent Euler flux functions II: Entropy production at shocks[END_REF] are entropy conservative providing that

β s = ([[v]] -[[u]] 2 /[[v]])/6. In practice, we set β s = max |v ± |, |g(u ± )|, ([[v]] - [[u]] 2 /[[v]]
)/6, 0 and v > 0 to get an entropy stable flux.

Euler equations in nonconservative form

he compressible Euler equations may be written in nonconservative form as

∂ t ρ + ∂ x (ρu) = 0, (39a) ∂ t (ρu) + ∂ x (ρu 2 + p) = 0, (39b) ∂ t (ρe) + ∂ x (ρeu) + p∂ x u = 0, (39c)
with ρ the density, u the velocity, e the specific internal energy. The equations are supplemented with an equation of states for ideal gas, p = (γ -1)ρe, with γ = C p C v the ratio of specific heats, and admissible solutions satisfy an entropy inequality for the pair η(u) = -ρs, q(u) = -ρus, with s = C v ln(Tρ 1-γ ) and T the temperature. System (39) satisfies the additional conservation law for the total specific energy E = e + u 2 2 :

∂ t (ρE) + ∂ x (ρEu + pu) = 0. ( 40
)
We now consider a class of fluctuations fluxes of the form

D -(u -, u + ) =           h ρ -ρ -u - h ρu + p -ρ -u -2 -p - h ρe -ρ -u -e -+ β -[[u]]           , D + (u -, u + ) =           -h ρ + ρ + u + -h ρu -p + ρ + u + 2 + p + -h ρe + ρ + u + e + + β + [[u]]           , (41)
where the unknown functions h ρ , h ρu , h ρe , p , and β ± depend on u ± . We look for entropy stable fluxes that further conserve the total energy, so we apply condition [START_REF] Cockburn | Runge-Kutta discontinuous Galerkin methods for convection-dominated problems[END_REF], with η(u) = ρE and q(u) = (ρE + p)u, and condition [START_REF] Coquel | Closure laws for a two-fluid two-pressure model[END_REF], with η(u) = -ρs and q(u) = -ρus, thus getting the two following relations

uh ρ -h ρu -p + β + + β -[[u]] = 0, (42a)
-C v h ρ [[ln θ]] + (γ -1)[[ln ρ]] + (pθ -β -θ --β + θ + )[[u]] + u[[pθ]] + h ρe [[θ]] ≥ 0, (42b)
where θ = T -1 . Looking for fluxes that satisfy the above relations for all jump values, it can be easily checked that the following choice is total energy preserving and entropy stable

h ρ = u ρ -k ρ [[ρ]], h ρu = uh ρ , h ρe = C v θ h ρ -k e [[e]], p = θp θ -k p [[u]], β ± = p 2 , ( 43
)
where â = [[a]] [[ln a]] denotes the logarithmic mean [START_REF] Ismail | Affordable, entropy-consistent Euler flux functions II: Entropy production at shocks[END_REF], providing that the parameters k ρ , k e , and k p are all non negative. Moreover, the fluxes are entropy conservative when k ρ = k e = k p = 0. In the numerical experiments, we will use the following choice of parameters k ρ = v ρc 2ρ and k e = k p = v ρc 2 with v ≥ 0 and c the speed of sound leading to the following dissipation rate of entropy

d η j dt + 1 h Q(U p j , U 0 j+1 ) -Q(U p j-1 , U 0 j ) = - v ρc j-1 2 2h C p -C v ρ [[ln ρ]][[ρ]] + θ[[u]] 2 -C v [[T]][[θ]] x j-1 2 ,
where Q(•, •) is defined from [START_REF] Després | Discontinuous Galerkin method for the numerical solution of Euler equations in axisymmetric geometry[END_REF], and may be compared to the dissipation rate for the formal equations

∂ t η(u) + ∂ x q(u) = - 4µ 3 θ(∂ x u) 2 - 3γC v 4Pr ∂ x T∂ x θ ,
where µ > 0 is the dynamic viscosity coefficient and Pr denotes the Prandtl number. Thus, v ρch plays the role of viscosity and setting k ρ = 0 in the numerical fluxes [START_REF] Munkejord | Comparison of Roe-type methods for solving the two-fluid model with and without pressure relaxation[END_REF], one recovers a dissipation rate which scales with the physical one for an equivalent Prandtl number 3γ/4 close to the physical value Pr = 0.72. We recall that the jump conditions depend on the underlying choice of path in (2) which should be consistent with the viscous profile leading to the physically relevant solution [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF][START_REF] Ph | Shock waves for nonlinear hyperbolic systems in nonconservative form[END_REF]. The design of the dissipation introduced in the numerical scheme is thus of great importance as many schemes fail to capture the relevant solution [START_REF] Abgrall | A comment on the computation of non-conservative products[END_REF][START_REF] Ch | A new comment on the computation of non-conservative products using Roe-type path conservative schemes[END_REF].

Isentropic Baer-Nunziato model

We finally consider the two-pressure two-velocity isentropic model [START_REF] Baer | A two-phase mixture theory for the deflagration-to-detonation transition (DDT) in reactive granular materials[END_REF][START_REF] Ambroso | The drift-flux asymptotic limit of barotropic two-phase two-pressure models[END_REF] with void fractions α i , densities ρ i , velocities u i , and a general equation of states p i = p i (ρ i ) with p i (ρ i ) > 0 and p i (ρ i ) < 0 for phases i = 1, 2. It is useful to introduce the specific internal energy e i and enthalpy h i of both phases defined by ρ 2 i e i (ρ i ) = p i (ρ i ) and

ρ i h i (ρ i ) = ρ i e i (ρ i ) + p i (ρ i ). Likewise, we introduce the speeds of sound c 2 i (ρ i ) = p i (ρ i ).

Two-phase flow model

Neglecting source terms modeling relaxation mechanisms, the governing equations have the form (4) with

u =                    α 1 α 1 ρ 1 α 1 ρ 1 u 1 α 2 ρ 2 α 2 ρ 2 u 2                    , f(u) =                    0 α 1 ρ 1 u 1 α 1 (ρ 1 u 2 1 + p 1 ) α 2 ρ 2 u 2 α 2 (ρ 2 u 2 2 + p 2 )                    , c(u)∂ x u =                    u 2 0 -p 1 0 p 1                    ∂ x α 1 , (44)
where u 2 and p 1 have been chosen as closure laws for the interface velocity and pressure, respectively. Both phases are assumed to satisfy the saturation condition

α 1 + α 2 = 1. ( 45
)
The set of states is Ω a = {u ∈ R 5 : ρ i > 0, α i > 0, i = 1, 2} and the system satisfies an entropy inequality (3) for the pair

η(u) = 2 i=1 α i ρ i u 2 i 2 + e i (ρ i ) , q(u) = 2 i=1 α i ρ i u 2 i 2 + h i (ρ i ) u i . (46)
We stress that the Baer-Nunziato system is only weakly hyperbolic and the assumptions in the introduction exclude resonance effects [START_REF] Berthon | Why many theories of shock waves are necessary: kinetic relations for non-conservative systems[END_REF], though the numerical experiments in section 6 will consider solutions close to resonance.

Note that given δ in R, the nonnegative functions ψ + δ (α 2 ) = max(α 2 -δ, 0) and ψ - δ (α 2 ) =min(α 2 -δ, 0) constitute convex entropies for the Baer-Nunziato system and satisfy the following inequalities in the sense of distributions

∂ t ψ ± δ (α 2 )ρ 2 + ∂ x ψ ± δ (α 2 )ρ 2 u 2 ≤ 0, (47)
which are easily obtained by combining the equations for the void fraction α 2 and the partial density α 2 ρ 2 .

Following the lines of Tadmor's proof of a minimum entropy principle for the equations of gas dynamics [55, section 3], together with the work by Serre on invariant domains for conservation laws [START_REF] Serre | Domaines invariants pour les systèmes hyperboliques de lois de conservation[END_REF], a maximum principle holds for the void fractions. This is summarized in the following lemma. Lemma 4.1 (maximum principle). The following estimates hold for solutions of the isentropic Baer-Nunziato model ( 4), [START_REF] Parès | Numerical methods for non-conservative hyperbolic systems: a theoretical framework[END_REF]:

ess inf |x|≤X+tu max 2 α 0 i (x) ≤ α i (x, t) ≤ ess sup |x|≤X+tu max 2 α 0 i (x), for almost all |x| ≤ X, t > 0, (48)
for i = 1, 2, where u max

2 = max C |u 2 | over C = {(x, τ) : |x| ≤ X + (t -τ)u max 2 , 0 ≤ τ ≤ t} and α 0 i (•) = α i (•, 0).
Proof. Integrating (47) over C, we get ∂C ψ ± δ (α 2 )ρ 2 (n t + u 2 n x )ds ≤ 0, where (n t , n x ) = 1 + (u max 2 ) 2 -1/2 (u max 2 , 1) denotes the unit normal on ∂C pointing outward C, so

n t + u 2 n x = 1 + (u max 2 ) 2 -1/2 (u max 2 + u 2 ) ≥ 0 by definition of u max 2 . We thus obtain 0 ≤ |x|≤X ψ ± δ (α 2 )ρ 2 dx ≤ |x|≤X+tu max 2 ψ ± δ (α 0 2 )ρ 0 2 dx,
where the left inequality follows from positivity of the integrand. Note that the integrations in the above inequalities are justified in the case of weak solutions, see [51, Sec. I.4].

Using the above result with the nonnegative functions ψ - δ (α) with δ = ess inf |x|≤X+tu max 2 α 0 2 , then ψ + δ (α) with δ = ess sup |x|≤X+tu max 2 α 0 2 , satisfying ψ ± δ (α 0 2 ) = 0, we get ψ ± δ (α 2 ) ≡ 0 almost everywhere in {(x, t) : |x| ≤ X}. We thus obtain (48) for α 2 . The same result also holds for α 1 through the saturation condition (45).

Entropy conservative numerical fluxes

Let us consider the following fluxes

D - ec (u -, u + ) = h(u -, u + ) -f(u -) + d -(u -, u + ), (49a) D + ec (u -, u + ) = f(u + ) -h(u -, u + ) + d + (u -, u + ), (49b)
with

h(u -, u + ) =                      0 α 1 u 1 ĥ1 (ρ - 1 , ρ + 1 ) α 1 u 1 2 ĥ1 (ρ - 1 , ρ + 1 ) + p 1 α 2 u 2 ĥ2 (ρ - 2 , ρ + 2 ) α 2 u 2 2 ĥ2 (ρ - 2 , ρ + 2 ) + p 2                      , d ± (u -, u + ) = [[α 1 ]] 2                     u ± 2 ± β s ±β s ĥ1 (ρ - 1 , ρ + 1 ) -p ± 1 ± β s u 1 ĥ1 (ρ - 1 , ρ + 1 ) ∓β s ĥ2 (ρ - 2 , ρ + 2 ) p ± 1 ∓ β s u 2 ĥ2 (ρ - 2 , ρ + 2 )                     , (50)
where β s > 0 is a measure of the spectral radius of A(u h ) and will be evaluated in Lemma 5.1.

Lemma 4.2. The fluctuation fluxes ( 49), together with the following choice for the numerical fluxes for the partial densities in (50):

ĥi (ρ - i , ρ + i ) = [[p i (ρ i )]] [[h i (ρ i )]] if ρ - i ρ + i , ρ i if ρ - i = ρ + i = ρ i , i = 1, 2, (51)
are entropy conservative.

Proof. Inserting ( 46) into ( 17) and using the Leibniz identities [[

u 2 i 2 ]] = u i [[u i ]], [[α i p i u i ]] = α i (p i [[u i ]] + u i [[p i ]]) + p i u i [[α i ]], i = 1, 2, (52)
we obtain

(η -) D - ec + (η + ) D + ec -[[q]] (49)
=

[[η (f -h)]] + (η -) d -+ (η + ) d + -[[q]] (52)
=

-α

1 u 1 ĥ1 [[h 1 - ¡ ¡ u 2 1 2 ]] -α 1 ( ¨ü 1 2 ĥ1 + p 1 )[[u 1 ]] + [[ $ $ $ $ $ $ $ α 1 ρ 1 u 1 (h 1 - u 2 1 2 )]] + [[α 1 ( & & & ρ 1 u 2 1 + p 1 )u 1 ]] -α 2 u 2 ĥ2 [[h 2 - ¡ ¡ u 2 2 2 ]] -α 2 ( ¨ü 2 2 ĥ2 + p 2 )[[u 2 ]] + [[ $ $ $ $ $ $ $ α 2 ρ 2 u 2 (h 2 - u 2 2 2 )]] + [[α 2 ( & & & ρ 2 u 2 2 + p 2 )u 2 ]] -[[ $ $ $ $ $ $ $ α 1 ρ 1 ( u 2 1 2 + h 1 )u 1 + $ $ $ $ $ $ $ α 2 ρ 2 ( u 2 2 2 + h 2 )u 2 ]] + [[α 1 ]] u 2 (p 2 -p 1 ) -u 1 p 1 + u 2 p 1 + β s 2 [[p 2 -p 1 ]] + ĥ1 [[h 1 - u 2 1 2 ]] + u 1 [[u 1 ]]) -ĥ2 [[h 2 - u 2 2 2 ]] + u 2 [[u 2 ]] (51) = -α 1 (u 1 [[p 1 ]] + p 1 [[u 1 ]]) + [[α 1 u 1 p 1 ]] -u 1 p 1 [[α 1 ]] -α 2 (u 2 [[p 2 ]] + p 2 [[u 2 ]]) + [[α 2 u 2 p 2 ]] -u 2 p 2 [[α 2 ]] (52)
= 0, which completes the proof. Some remarks are in order. The numerical conservation flux h(•, •) in ( 49) is symmetric, consistent and differentiable, while the fluctuation fluxes have the form (27a) with A(u -, u + ) = (u 2 , 0, -p 1 , 0, p 1 ) and therefore satisfy (27c,d) and are path-conservative ( 16) for a linear path in u 2 and p 1 . Due to the presence of the nonlinear fluxes ĥi , the d ± are examples of fluctuation fluxes in non-splitting form. Finally, the DGSEM with the fluxes ( 49) is by construction conservative for the mixture density and momentum.

High-order DGSEM for the isentropic Baer-Nunziato model

Entropy stable fluxes

We now focus on the design of a positive and entropy stable DG scheme for the two-pressure two-velocity isentropic model ( 4) with [START_REF] Parès | Numerical methods for non-conservative hyperbolic systems: a theoretical framework[END_REF]. For that purpose, we introduce the fully discrete scheme for a one-step first-order explicit time discretization and analyze its properties. High-order time integration will be done by using strong-stability preserving explicit Runge-Kutta methods [START_REF] Shu | Efficient implementation of essentially non-oscillatory shock-capturing schemes[END_REF] that keep the properties of the first-order in time scheme. Let t (n) = n∆t, with ∆t > 0 the time step, set λ = ∆t h , and use the notations n) ). The DGSEM scheme for solving the isentropic Baer-Nunziato equations reads

u (n) h (•) = u h (•, t (n) ) and U k,n j = U k j (t (
ω k h 2 U k,n+1 j -U k,n j ∆t + R k j (u (n) h ) = 0, ( 53
)
with R k j (•) defined in [START_REF] Coquel | Adjoint approximation of nonlinear hyperbolic systems with non-conservative products[END_REF] and where the entropy conservative fluxes [START_REF] Rhebergen | Discontinuous Galerkin finite element methods for hyperbolic nonconservative partial differential equations[END_REF] are used in the definition of [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF]. We follow the strategy in [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF] to design entropy stable fluxes at interfaces:

D ± (u -, u + ) = D ± ec (u -, u + ) ± v D v (u -, u + )[[η (u)]], (54)
with v ≥ 0 and the positive diagonal matrix

D v (u -, u + ) = max(c 1 , c 2 ) diag 0, α 1 ρ 1 u 1 2 +c 1 2 , α 1 ρ 1 , α 1 ρ 1 u 1 2 +c 1 2 , α 2 ρ 2 . ( 55
)

Properties of the discrete scheme

We have the following results that guarantee positivity of the solution and the maximum principle [START_REF] Renac | A robust high-order discontinuous Galerkin method with large time steps for the compressible Euler equations[END_REF] for general equation of states for the fully discrete solution of the DGSEM.

Theorem 5.1. Assume that ρ 0≤k≤p,n i, j∈Z > 0 and α 0≤k≤p,n i, j∈Z > 0 for i = 1, 2, then under the CFL conditions

λ max j∈Z max 0≤k≤p 1 ω k u (n) 2,h , d x ϕ k j p j + δ k,p β s -u p,n 2, j 2 + δ k,0 β s +u 0,n 2, j 2 < 1 2 , (56a)
λ max j∈Z max i=1,2 (β s +ũ n i, j ) ĥi (U 0,n j ,U p,n j )-ρ 0,n i, j ci (U 0,n j ,U p,n j ) ρ 0,n i, j + (β s -u n i, j-1 2 ) ĥi (U p,n j-1 ,U 0,n j )+ρ 0,n i, j ci (U p,n j-1 ,U 0,n j ) ρ 0,n i, j , (β s +u n i, j+ 1 2
) ĥi (U p,n j ,U 0,n j+1 )-ρ p,n i, j ci (U p,n j ,U 0,n j+1 ) ρ p,n i, j

+ (β s -ũ n i, j ) ĥi (U 0,n j ,U p,n j )+ρ p,n i, j ci (U 0,n j ,U p,n j ) ρ p,n i, j < 2 p(p + 1) , ( 56b
)
we have for the cell averages at time t (n+1)

α i ρ i (n+1) j > 0, α i (n+1) j > 0, i = 1, 2, j ∈ Z, (57)
and

α 1 (n+1) j = p k=0 ω k 2 -λ u (n) 2,h , d x ϕ k j p j + δ k,p β s -u p,n 2, j 2 + δ k,0 β s +u 0,n 2, j 2 α k,n 1, j + λ β s -u p,n 2, j 2 α 0,n 1, j+1 + λ β s +u 0,n 2, j 2 α p,n 1, j-1 , j ∈ Z, ( 58
)
is a convex combination of DOFs at time t (n) where

β s = max j∈Z max i=1,2 -u 0,n 2, j , u p,n 2, j , ρ p,n i, j-1 ci (U p,n j-1 ,U 0,n j ) ĥi (U p,n j-1 ,U 0,n j ) -u n i, j-1 2 , ũn i, j - ρ p,n i, j ci (U 0,n j ,U p,n j ) ĥi (U 0,n j ,U p,n j ) , ρ 0,n i, j ci (U 0,n j ,U p,n j ) ĥi (U 0,n j ,U p,n j ) -ũn i, j , u n i, j+ 1 2 - ρ 0,n i, j+1 ci (U p,n j ,U 0,n j+1 ) ĥi (U p,n j ,U 0,n j+1 ) , (59)
ũn i, j = u 0,n i, j +u p,n i, j 2 , u n i, j+ 1 2 = u p,n i, j +u 0,n i, j+1 2
, and ci (u

-, u + ) = v max(c 1 , c 2 ) [[h i - u 2 i 2 ]] u 2 i +c 2 i .
Proof. We know by Theorem 3.2 that the cell averaged semi-discrete DGSEM [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] holds. Thus, summing the first component in ( 53) over 0 ≤ k ≤ p gives

α 1 (n+1) j := p k=0 ω k 2 α k,n+1 1, j (53)
= (54) α 1 (n) j -λ u (n) 2,h , ∂ x α (n) 1,h p j -λ d - 1 (U p,n j , U 0,n j+1 ) + d + 1 (U p,n j-1 , U 0,n j ) (5) = (50) p k=0 ω k 2 α k,n 1, j -λ p k=0 u (n) 2,h , d x ϕ k j p j α k,n 1, j -λ u p,n 2, j -β s 2 (α 0,n 1, j+1 -α p,n 1, j ) + u 0,n 2, j +β s 2 (α 0,n 1, j -α p,n 1, j-1 ) = p-1 k=1 ω k 2 -λ u (n) 2,h , d x ϕ k j p j α k,n 1, j + ω 0 2 -λ u (n) 2,h , d x ϕ 0 j p j -λ β s +u 0,n 2, j 2 α 0,n 1, j + λ β s +u 0,n 2, j 2 α p,n 1, j-1 + ω p 2 -λ u (n) 2,h , d x ϕ p j p j -λ β s -u p,n 2, j 2 α p,n 1, j + λ β s -u p,n 2, j 2 α 0,n 1, j+1 ,
which is indeed a convex combination under the conditions (56a) and β s ≥ max(-u 0,n 2, j , u p,n 2, j ) satisfied by [START_REF] Tokareva | A flux splitting method for the Baer-Nunziato equations of compressible two-phase flow[END_REF]. The positivity of the cell averaged partial densities relies on the techniques introduced in [START_REF] Perthame | On positivity preserving finite volume schemes for Euler equations[END_REF][START_REF] Zhang | On positivity-preserving high order discontinuous Galerkin schemes for compressible Euler equations on rectangular meshes[END_REF] to rewrite a conservative high-order scheme for the mean value as a convex combination of positive first-order schemes. Since the equations for partial densities of phases i = 1, 2 are conservative, we get by summing the 2i-th component in ( 53) over 0 ≤ k ≤ p:

α i ρ i (n+1) j = α i ρ i (n) j -λ H 2i (U p,n j , U 0,n j+1 ) -H 2i (U p,n j-1 , U 0,n j ) = p-1 k=1 ω k 2 α k,n i, j ρ k,n i, j + ω p 2 α p,n i, j ρ p,n i, j -2λ ω p H 2i (U p,n j , U 0,n j+1 ) -H 2i (U 0,n j , U p,n j ) + ω 0 2 α 0,n i, j ρ 0,n i, j -2λ ω 0 H 2i (U 0,n j , U p,n j ) -H 2i (U p,n j-1 , U 0,n j ) (60)
which is also a convex combination with coefficients ω k 2 of DOFs and updates of a three-point scheme with numerical flux 54) and [START_REF] Sainsaulieu | Ondes progressives solutions de systèmes convectifs-diffusifs et systèmes hyperboliques non conservatifs[END_REF]. It remains to show that this scheme is positive. Consider the scheme

H 2i (u -, u + ) = α i u i -β s 2 [[α i ]] ĥi -ci (u -, u + )α 1 ρ 1 from (
αρ n+1 j = αρ n j -λ FV H(U n j , U n j+1 ) -H(U n j-1 , U n j ) , (61)
where we have removed indices on phases for the sake of clarity. This scheme may be easily recast into the form

αρ n+1 j = 1 -λ FV (C - j+ 1 2 + C + j-1 2 ) αρ n j + λ FV C + j+ 1 2 αρ n j+1 + λ FV C - j-1 2 αρ n j-1 with coefficients C ± = 1 2ρ ± (β s ∓ u) ĥ ± ρ ± c for which positivity follows from conditions C ± j± 1 2 ≥ 0 and λ FV (C - j+ 1 2 + C + j-1 2 ) ≤ 1. Identifying with (60) with λ = ω 0 2 λ FV = ω p 2 λ FV , ω 0 = ω p = 2 p(p+1)
, we obtain conditions (56b) and (59).

Remark 5.1. In practice, we replace ( 56) and ( 59) by the more affordable conditions

λ max j∈Z max 0≤k≤p 1 ω k u (n) 2,h , d x ϕ k j p j + δ k,p β s -u p,n 2, j 2 + δ k,0 β s +u 0,n 2, j 2 < 1 2 , (62)
and

β s = max j∈Z max i=1,2 |u p,n i, j-1 | + c p,n i, j-1 , |u 0,n i, j+1 | + c 0,n i, j+1 . (63)
The rationale for this choice is as follows. First, [START_REF] Tokareva | HLLC-type Riemann solver for the Baer-Nunziato equations of compressible two-phase flow[END_REF] still is a convex combination under ( 62) and [START_REF] Zhang | On positivity-preserving high order discontinuous Galerkin schemes for compressible Euler equations on rectangular meshes[END_REF]. Then, let go back to the three-point scheme [START_REF] Volpert | The space BV and quasilinear equations[END_REF] and consider polytropic ideal gas: p = (γ -1)ρe, h = γe, and c 2 = (γ -1)h with γ > 1. The positivity condition reads β s ≥ |u| + ρ ± c ĥ with ρ ± c ĥ

= v ρ ± c [[h]] [[p]] [[h-u 2 /2]] u 2 +c 2 . Observe that [[p]] ρ ± [[h]] = (γ-1)ρ γρ ± 1 + e[[ρ]] ρ[[e]] ≥ γ-1 2γ , [[h- u 2 2 ]] u 2 +c 2 ≤ [[h]] (γ-1)h + [[u 2 ]] 2u 2 ≤ γ+1 γ-1 ,
from positivity of the thermodynamic variables and the fact that e (ρ) > 0. We thus get ρ ± c ĥ ≤ 2γ(γ+1) (γ-1) 2 v c, so ( 63) implies (59) for v ≤ (γ-1) 2 2γ(γ+1) . Though this condition seems restrictive, we used v 1 2p+1 during the numerical experiments without violating positivity (57) of Theorem 5.1.

The following result is useful to prevent spurious oscillations in the numerical solution. Indeed, the present schemes satisfies the Abgrall's criterion [START_REF] Abgrall | How to prevent pressure oscillations in multicomponent flow calculations: a quasi conservative approach[END_REF] that states that uniform velocity and pressure must remain uniform at all time.

Theorem 5.2 (Abgrall's criterion). Assume that the velocity and pressure are uniform and equal at time t (n) :

u k,n i, j = u, p k,n i, j = p, i = 1, 2, ∀ j ∈ Z, 0 ≤ k ≤ p, (64)
then they remain uniform and equal at time t (n+1) .

Proof.

The assumption (64) on pressures require uniform densities ρ k,n i, j = ρ i , i =

]] = 0. Then, the entropy conservative fluxes [START_REF] Rhebergen | Discontinuous Galerkin finite element methods for hyperbolic nonconservative partial differential equations[END_REF] and entropy stable fluxes (54) reduce to

D ± (u -, u + ) = D ± ec (u -, u + ) = [[α 1 ]] 2                    u ± β s ρ 1 (u ± β s ) ρ 1 u(u ± β s ) -ρ 2 (u ± β s ) -ρ 2 u(u ± β s )                    , so the explicit residuals in (53) become R k j (u h ) = R k 1, j (1, ρ 1 , ρ 1 u, -ρ 2 , -ρ 2 u) with R k 1, j = ω k u p l=0 α k 1, j D kl + δ kp (u-β s ) 2 [[α 1 ]] j+ 1 2 + δ k0 (u+β s ) 2 [[α 1 ]] j-1 2 , (65a) (45)
=

-

ω k u p l=0 α k 2, j D kl - δ kp (u-β s ) 2 [[α 2 ]] j+ 1 2 -δ k0 (u+β s ) 2 [[α 2 ]] j-1 2 . (65b)
We thus rewrite [START_REF] Stewart | Two-phase flow: models and methods[END_REF] as

ω k h 2 α k,n+1 1, j -α k,n 1, j ∆t + R k,n 1, j = 0, (66a)
ω k h 2 α k,n+1 1, j (ρ k,n+1 1, j -ρ 1 )+ρ 1 (α k,n+1 1, j -α k,n 1, j ) ∆t + ρ 1 R k,n 1, j = 0, ( 66b
)
ω k h 2 (α 1 ρ 1 ) k,n+1 j (u k,n+1 1, j -u)+u (α 1 ρ 1 ) k,n+1 -ρ 1 α k,n 1, j ∆t + ρ 1 uR k,n 1, j = 0, (66c)
ω k h 2 α k,n+1 2, j (ρ k,n+1 2, j -ρ 2 )+ρ 2 (α k,n+1 2, j -α k,n 2, j ) ∆t -ρ 2 R k,n 1, j = 0, (66d)
ω k h 2 (α 2 ρ 2 ) k,n+1 j (u k,n+1 2, j -u)+u (α 2 ρ 2 ) k,n+1 -ρ 2 α k,n 2, j ∆t -ρ 2 uR k,n 1, j = 0. (66e)
Then, (66b)-ρ 1 ×(66a) implies ρ k,n+1

1, j = ρ 1 , while (66c)-ρ 1 u×(66a) gives u k,n+1

1, j = u. Then, using (65) and (66d,e) we obtain ρ k,n+1 2, j = ρ 2 and u k,n+1 2, j = u.

Limiting strategy

The properties in Theorem 5.1 hold only for the cell averaged value of the numerical solution at time t (n+1) , which is not sufficient for robustness and stability of numerical computations. However, these results motivate the use of a posteriori limiters introduced in [START_REF] Zhang | On maximum-principle-satisfying high order schemes for scalar conservation laws[END_REF][START_REF] Zhang | On positivity-preserving high order discontinuous Galerkin schemes for compressible Euler equations on rectangular meshes[END_REF]. These limiters aim at extending preservation of invariant domains [START_REF] Zhang | On positivity-preserving high order discontinuous Galerkin schemes for compressible Euler equations on rectangular meshes[END_REF] or maximum principle [START_REF] Zhang | On maximum-principle-satisfying high order schemes for scalar conservation laws[END_REF] from mean values to nodal values within elements.

We enforce positivity of nodal values of partial densities and the maximum principle (48) by using the linear limiter Ũk,n+1

j = θ j (U k,n+1 j -u (n+1) j ) + u (n+1) j , 0 ≤ k ≤ p, j ∈ Z, (67)
with 0 ≤ θ j ≤ 1 defined by θ j := min(θ ρ i j , θ α i j : i = 1, 2) where

θ ρ i j = min ρ i (n+1) j - ρ i (n+1) j -ρ min i, j , 1 , ρ min i, j = min 0≤k≤p ρ k,n+1 i, j , θ α i j = min α i (n+1) j -m α i j α i (n+1) j -α min i, j , α i (n+1) j -M α i j α i (n+1) j -α max i, j , 1 , α min i, j = min 0≤k≤p α k,n+1 i, j , α max i, j = max 0≤k≤p α k,n+1 i, j , 0 <
1 is a parameter (we set = 10 -8 in our numerical tests), and

m α i j = min α p,n i, j-1 , α 0≤k≤p,n i, j , α 0,n i, j+1 , M α i j = max α p,n i, j-1 , α 0≤k≤p,n i, j , α 0,n i, j+1 .
The limiter (67) guarantees that ρ0≤k≤p,n+1 j ≥ together with a discrete maximum principle on the void fractions Likewise, phase densities and velocities remain unchanged by the limiter (67) so uniform velocity and pressure profiles are conserved.

m α i j ≤ α0≤k≤p,n+1 i, j ≤ M α i j . It

Numerical experiments

In the following, we consider Riemann problems for the nonconservative systems introduced in section 4 associated to initial conditions

u 0 (x) = u L , x < 0, u R , x > 0.
The set of initial conditions u 0 is given in Table 1. For the time integration, we use the 3-stage third-order strongstability preserving Runge-Kutta time integration scheme of Shu and Osher [START_REF] Shu | Efficient implementation of essentially non-oscillatory shock-capturing schemes[END_REF].

Entropy conservative schemes

We first validate the entropy conservation property of the modified numerical scheme ( 19)-( 21) stated in Corollary 3.1 which constitutes the basis of the theoretical results of section 3. We follow the experimental setup introduced in [START_REF] Bohm | An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations. Part I: Theory and Numerical Verification[END_REF]. We consider some examples of section 4 and use the corresponding entropy conservative fluxes in [START_REF] Coquel | Adjoint approximation of nonlinear hyperbolic systems with non-conservative products[END_REF] and [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF] and choose an initial condition resulting in the development of an isolated shock wave with a moderate strength (see Table 1) on a domain of unit length with periodic boundary conditions. As a result of entropy conservation of the space discretization, only the time integration scheme should modify the global entropy budget at the discrete level. We thus evaluate the difference Table 1: Initial conditions and physical parameters of Riemann problems with U = u for the Burgers' equation [START_REF] Gassner | A well balanced and entropy conservative discontinuous Galerkin spectral element method for the shallow water equations[END_REF], U = (u, v) for the 2 × 2 system (36), U = (ρ, u, p) for the Euler system [START_REF] Jiang | On a cell entropy inequality for discontinuous Galerkin methods[END_REF], and U = (α 1 , ρ 1 , u 1 , ρ 2 , u 2 ) for the isentropic Baer-Nunziato system (4), [START_REF] Parès | Numerical methods for non-conservative hyperbolic systems: a theoretical framework[END_REF]

E Ω (t) := j η(u h ) j (t) -η(u 0 ) j , (68)
which quantifies the difference between the discrete entropy at final time and the initial entropy over the domain Ω h . Results are shown in Table 2, for examples 4.1 to 4.4. We use a space grid of N = 100 elements and fourthorder accuracy in space p = 3. For the Euler equations, we use the initial condition from [START_REF] Bohm | An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations. Part I: Theory and Numerical Verification[END_REF] and test the entropy conservation with numerical fluxes [START_REF] Ph | Shock waves for nonlinear hyperbolic systems in nonconservative form[END_REF] with numerical parameter v = 0 for η(u) = -ρs, and v = 1 for η(u) = ρE. We observe that the errors (68) decrease to machine accuracy when refining the time step with an experimental order of convergence of order three corresponding to the theoretical approximation order of the time integration scheme, thus validating the entropy conservative property of the numerical schemes.

Riemann problems

We now consider more challenging test cases and test the effect of the present method on the robustness and stability of numerical computations of discontinuous solutions. Figures 3 to 7 compare the numerical solution in symbols with the exact solution in lines.

We evaluate the time step with a safety factor of ∆t = 0.9 × λh, where λ is evaluated from

λ max j∈Z max 0≤k≤p ρ s (A U k,n k ) , 2 v h 2p+1 ≤ 1 2p+1 ,
where ρ s (•) denotes the spectral radius, for examples 4.1 to 4.3 and from [START_REF] Zhang | On maximum-principle-satisfying high order schemes for scalar conservation laws[END_REF] for the isentropic Baer-Nunziato model.

Entropy stable schemes at element interfaces are obtained by adding viscosity operators that mimic, at the discrete level, a physical parabolic regularization in the same way as done in [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF]. These operators are introduced in sections 4.2, 4.3 and 5. Let us stress that we here consider systems [START_REF] Hiltebrand | Entropy stable shock capturing spacetime discontinuous Galerkin schemes for systems of conservation laws[END_REF] and ( 4), [START_REF] Parès | Numerical methods for non-conservative hyperbolic systems: a theoretical framework[END_REF] having nonconservative products associated with LD characteristic fields for which finite difference schemes have been shown to converge to the physically relevant solution [START_REF] Castro | Why many theories of shock waves are necessary: Convergence error in formally path-consistent schemes[END_REF]. However the present strategy may fail for strong shocks where the agreement between regularizations at discrete and continuous levels may not be satisfied [START_REF] Castro | Why many theories of shock waves are necessary: Convergence error in formally path-consistent schemes[END_REF].

Burgers' equation

We use the entropy conservative fluxes [START_REF] Gassner | Split form nodal discontinuous Galerkin schemes with summation-by-parts property for the compressible Euler equations[END_REF] in the volume integral [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF], while we use the following fluctuation fluxes at interfaces

D ± (u -, u + ) = ∓h(u -, u + ) ± u ± , h(u -, u + ) = min{ v 2 2 : v ∈ [u -, u + ]}, if u -≤ u + , max{ v 2 2 : v ∈ [u + , u -]}, if u -> u + . (69)
where h(•, •) is the Godunov numerical flux and is entropy stable. Figure 2 displays solutions at different times obtained with the nonconservative Burgers' equation [START_REF] Gassner | A well balanced and entropy conservative discontinuous Galerkin spectral element method for the shallow water equations[END_REF] and a smooth initial condition u 0 (x) = sin(2πx). Results highlight the good resolution of the shock wave and the decrease of oscillations within cells while p is increased with the entropy stable scheme. In contrast, using the same numerical flux at interfaces (69) with the original form of the volume integral in (13) leads to large oscillations which amplify in time when the solution becomes discontinuous as observed in Fig. 2d.

Nonconservative product associated to a LD field

Figure 3 shows results for a 1-shock, 2-contact problem (RP0 in Table 1) for system [START_REF] Hiltebrand | Entropy stable shock capturing spacetime discontinuous Galerkin schemes for systems of conservation laws[END_REF]. We compare solutions obtained with the entropy stable scheme [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF], with (21) evaluated from the entropy conservative fluxes [START_REF] Henry De Frahan | A new limiting procedure for discontinuous Galerkin methods applied to compressible multiphase flows with shocks and interfaces[END_REF], or with the original DGSEM [START_REF] Castro | Why many theories of shock waves are necessary: Convergence error in formally path-consistent schemes[END_REF]. In both cases, we use the same entropy stable numerical fluxes [START_REF] Ismail | Affordable, entropy-consistent Euler flux functions II: Entropy production at shocks[END_REF] at interfaces. The results highlight the importance of the modification of the volume integral in [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] to satisfy the entropy inequality. The second order solution without this modification does not tend to the exact weak solution and contains non-physical waves even when the mesh is refined. We note that higher-order computations for p ≥ 2 without the correction [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF] were seen to blow up due to a change of sign in the u component of the solution which induces a loss of strict hyperbolicity of system [START_REF] Hiltebrand | Entropy stable shock capturing spacetime discontinuous Galerkin schemes for systems of conservation laws[END_REF]. The correction [START_REF] Maso | Definition and weak stability of nonconservative products[END_REF] successfully stabilizes the computations and the numerical solution now tends to the exact entropy solution.

Nonconservative Euler equations

Again, we first analyze the effect of the formulation of the volume integral. The results are shown in Fig. 4 for a third-order discretization of an isolated shock wave [START_REF] Ch | A new comment on the computation of non-conservative products using Roe-type path conservative schemes[END_REF] (test "shock" in Table 1) where we compare two different 19),( 41), (b) scheme ( 13), [START_REF] Ph | Shock waves for nonlinear hyperbolic systems in nonconservative form[END_REF].

schemes with entropy stable fluxes (41) at interfaces: (ES) with the modified volume integral [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF] built with entropy conservative fluxes from (41); (no ES) with the original volume integral [START_REF] Castro | Why many theories of shock waves are necessary: Convergence error in formally path-consistent schemes[END_REF]. As a consequence, the latter scheme is not entropy stable. The entropy stable DGSEM scheme is seen to converge to the physically relevant shock solution, while a spurious discontinuity is observed with the original volume integral for which more artificial viscosity in the numerical fluxes, v = 2 in the parameters of [START_REF] Munkejord | Comparison of Roe-type methods for solving the two-fluid model with and without pressure relaxation[END_REF], was necessary to stabilize the computation. The rather simplicity of the numerical fluxes (41) may be compared to the inability of the path-conservative Roe-type fluxes [START_REF] Parès | Numerical methods for non-conservative hyperbolic systems: a theoretical framework[END_REF] to capture the entropy solution as pointed out in [START_REF] Abgrall | A comment on the computation of non-conservative products[END_REF][START_REF] Ch | A new comment on the computation of non-conservative products using Roe-type path conservative schemes[END_REF] and to the relative complexity of the correction proposed by [START_REF] Ch | A new comment on the computation of non-conservative products using Roe-type path conservative schemes[END_REF] relying on Glimm's random sampling strategy in the context of first-order finite volume methods. Figure 5 displays the solutions obtained for the three Riemann problems of Table 1 with a fourth-order entropy stable scheme with fluxes (41) and an artificial viscosity parameter v = 1. We observe a sharp capture of shock and contact waves as well as a very good resolution of the rarefaction waves. However, high-frequency oscillations develop in the neighborhood of contacts and cover a large region for the Lax problem. These oscillations are not damped by the scheme even with a large dissipation coefficient v and their amplitude increases with p (results not shown here). They may be related to the fact that the dissipation is introduced in the numerical scheme only through interfaces, so internal DOFs U 0<k<p k may suffer from a lack of stabilization mechanisms. As a consequence, entropy stability holds for the cell-averaged scheme and only large scale components of the error are damped efficiently.

Isentropic Baer-Nunziato model

For the numerical experiments on the isentropic Baer-Nunziato model ( 4),( 44), we consider polytropic ideal gas with equations of state of the form p i (ρ i ) = κρ γ i i with κ > 0 and γ i > 1, i = 1, 2. Computations are done with the entropy stable numerical scheme [START_REF] Stewart | Two-phase flow: models and methods[END_REF] and fourth-order accuracy, p = 3. The limiter (67) is applied at the end of each stage unless stated otherwise. Problems RP1 and RP2 come from [START_REF] Coquel | A robust entropy-satisfying finite volume scheme for the isentropic Baer-Nunziato model[END_REF], while RP3 is adapted from [START_REF] Munkejord | Comparison of Roe-type methods for solving the two-fluid model with and without pressure relaxation[END_REF].

We first consider the advection of a discontinuity in the void fraction with uniform velocities, u 1,0 = u 2,0 = 1 and pressures, p 1,0 = p 2,0 = 1.25, so the mass and momentum equations in ( 4),( 44) are trivially satisfied. The pressure law parameters are κ = 1, γ 1 = 1.4, and γ 2 = 1.2. Figure 6 presents the solution obtained at time t = 0.1 with and without limiter. In both cases, the velocities and pressures remain uniform as expected from Theorem 5.2, but the limiter is seen to introduce numerical dissipation that smears the contact discontinuity. The design of a sharp limiter would help to improve the solution but is beyond the scope of the present study where we rather focus on stability and robustness issues.

Figure 7 presents the solution of Riemann problems associated to the initial conditions of Table 1. For RP1 and RP2, we use κ = 1, γ 1 = 3 and γ 2 = 1.5. RP2 considers solutions close to resonance with a vanishing phase 2 and the contact discontinuity separates a mixture region where the two phases coexist from a single phase region. The shock and rarefaction waves are well captured, while the contact wave is slightly diffused as an effect of the limiter as observed in the precedent experiment. RP3 is adapted from the experiment with large relative velocity for one pressure models in [START_REF] Munkejord | Comparison of Roe-type methods for solving the two-fluid model with and without pressure relaxation[END_REF] and we set κ = 10 5 and γ 1 = γ 2 = 1.4. Spurious oscillations of low amplitude are observed in the neighborhood of the strong shocks as a result of the limiter (67) which only guarantees positivity by reducing the amplitude of these oscillations, but the results are in good quantitative agreement with the exact solution. We stress that our experiments show that the correction (21) of the volume integral is strongly needed for stabilizing the computations which would blow up otherwise.

Concluding remarks

In this work, we introduce a general framework for the design of entropy stable DGSEM for the discretization of nonlinear hyperbolic systems in nonconservative form. The framework relies on the use of SBP operators and twopoint entropy conservative fluctuation fluxes [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF] to evaluate the integral over discretization elements, thus removing its contribution to the global entropy production within the element, together with entropy stable fluxes at element interfaces. This contribution may be seen as a generalization of the framework on entropy stable spectral collocation schemes for conservation laws introduced in [START_REF] Fisher | High-order entropy stable finite difference schemes for nonlinear conservation laws: Finite domains[END_REF][START_REF] Carpenter | Entropy stable spectral collocation schemes for the Navier-Stokes equations: discontinuous interfaces[END_REF][START_REF] Gassner | A Skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF][START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF]. In particular, the generalizations to multiple space dimensions with quadrangles, hexahedra, or simplex elements; the use of bound-preserving or TVD limiters; and the discretization of viscous terms will keep the entropy inequality for the DGSEM as shown in [START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF]. Applications show that the method proves to be robust, stable and entropy satisfying for the high-order discretization of the nonconservative Burgers and Euler equations, and two-phase flow models: a 2 × 2 system with a nonconservative product associated to a LD field and the isentropic Baer-Nunziato model. Numerical experiments show that the proposed modification of the integral over discretization elements has a very beneficial effect on the stability and robustness of the computations. However, spurious oscillations still remain in the numerical solutions in the neighborhood of discontinuities which may be attributed to some of the drawbacks of the present approach: semidiscrete entropy stability; entropy stability for the cell-averaged DGSEM only; stabilization introduced only through numerical fluxes that does not stabilize internal DOFs; positivity of the solution not guaranteed with the Lagrange interpolation polynomials even if nodal values are positive; integration error associated to mass lumping; etc. Future work will concern the improvement of the limiter to preserve contact discontinuities, the analysis of the well-balanced property, and the extension of the method to the Baer-Nunziato model with general equations of states including the transport equations for partial energies [START_REF] Baer | A two-phase mixture theory for the deflagration-to-detonation transition (DDT) in reactive granular materials[END_REF].

D - ec (u -, u + ) =                 ĥg -α -ρ - g u - g ĥg u g + α p -α -(ρ - g (u - g ) 2 + p -) -p - 2 [[α]] ĥl -(1 -α -)ρ l u -

Figure 1 :

 1 Figure 1: Mesh with definition of left and right traces at interfaces x j± 1 2

  also keeps the entropy inequality (3) at the discrete level in the sense that[START_REF] Chen | Entropy stable high order discontinuous Galerkin methods with suitable quadrature rules for hyperbolic conservation laws[END_REF] Lemma 3.1] for η convex we have η

Figure 2 :

 2 Figure 2: Burgers: solution obtained at different times t = 0.1, 0.2, 0.3, 0.4, and 0.5 for a polynomial degree p, N = 20 cells and entropy stable (ES) modification (21) or not (no ES).

Figure 3 : 2 ×

 32 Figure 3: 2 × 2 system: RP0 discretized with polynomial degree p, N cells and entropy stable (ES) modification (21) or not (no ES).

Figure 4 :

 4 Figure 4: Euler equations: density profiles for the shock problem discretized with a polynomial degree p = 2, N = 100 or N = 400 cells and entropy stable numerical fluxes (41) at interfaces; (a) scheme (19),(41), (b) scheme (13),[START_REF] Ph | Shock waves for nonlinear hyperbolic systems in nonconservative form[END_REF].

Figure 5 :Figure 6 :

 56 Figure 5: Euler equations: Riemann problems discretized with a polynomial degree p = 3, N = 100 cells and entropy stable scheme.

Figure 7 :

 7 Figure 7: Isentropic Baer-Nunziato model: Riemann problems discretized with a polynomial degree p = 3, N = 100 cells and entropy stable scheme.

  .

	test	model	left state U L	right state U R	t
	EC1	(34)	1.1	1	0.1
	EC2	(36)	0.5, 1.1	0.5, 1	0.15
	RP0	(36)	3, 0.5	0.75, 1	0.15
	EC3	(39)	1, 0, 1	1.125, 0, 1.1	0.5
	shock	(39)	0.4765625, 2.3046638387921279, 1	0.125, 0, 0.1	0.1
	Sod	(39)	1, 0, 1	0.125, 0, 0.1	0.2
	Lax	(39)	0.445, 0.698, 3.528	0.5, 0, 0.571	0.13

Table 2 :

 2 Entropy conservation error (68) and corresponding order of convergence O (see Table1for a description of the test cases).

	test case	EC1		EC2		EC3 (η(u) = -ρs)	EC3 (η(u) = ρE)	EC4
		E Ω (t)	O	E Ω (t)	O	E Ω (t)	O	E Ω (t)	O	E Ω (t)	O
	∆t	6.6310e-06	-	5.3463e-07	-	9.5312e-07	-	2.3712e-06	-	1.2134e-07	-
	∆t/2	1.9086e-06 1.80 7.0781e-08 2.92 2.4818e-07 1.94 6.1800e-07 1.94 1.5786e-08 2.94
	∆t/4	2.8396e-07 2.75 8.9489e-09 2.98 3.5602e-08 2.80 8.8700e-08 2.80 2.0335e-09 2.96
	∆t/8	3.6485e-08 2.96 1.1206e-09 3.00 4.5374e-09 2.97 1.1303e-08 2.97 2.5627e-10 2.99
	∆t/16	4.5811e-09 3.00 1.4004e-10 3.00 5.6837e-10 3.00 1.4137e-09 3.00 3.1668e-11 3.02
	∆t/32	5.7291e-10 3.00 1.7320e-11 3.02 7.0299e-11 3.02 1.7074e-10 3.05 3.1352e-12 3.34
	∆t/64	7.1240e-11 3.01 1.7950e-12 3.27 7.2221e-12 3.28 9.2604e-12 4.20 1.2636e-12 1.31

D ± roe (u -, u + ) = A(u -, u + ) ± [[u]] = A(ṽ)[[u]], with ṽ = (τ, u, p) .

Appendix A. Entropy conservative fluctuation fluxes

In the following, we provide further examples of entropy conservative numerical fluxes of the form [START_REF] Dumbser | A simple robust and accurate a posteriori sub-cell finite volume limiter for the discontinuous Galerkin method on unstructured meshes[END_REF] for different nonconservative systems. We recall the definitions of the jump and average operators: [[u]] := u +u -and u := u -+u + 2 .

Appendix A.1. coupled Burgers' equation The following nonconservative system was first proposed in [START_REF] Berthon | Nonlinear scheme for approximating a non-conservative hyperbolic system[END_REF]:

where we recover the Burgers' equation for the sum u + v. Entropy and entropy flux are therefore η(u) = (u+v) 2 2 and q(u) = (u+v) 3 3 . Entropy conservative fluctuation fluxes of the form (27) may also be derived:

which correspond to the path-conservative and entropy conservative fluxes derived in [START_REF] Castro | Entropy conservative and entropy stable schemes for nonconservative hyperbolic systems[END_REF].

Appendix A.2. Euler equations in Lagrangian coordinates

The Euler equations in Lagrangian coordinates may be written in nonconservative form:

with τ the specific volume, u the velocity, e the specific internal energy. The equations are supplemented with a general equation of states for the pressure p = p(τ, e) and admissible solutions satisfy the entropy inequality

with Tds = de + pdτ, and T the temperature. Entropy conservative fluctuation fluxes are

. Shallow water equations

The Shallow water equations with topography read

where h, u, b are the water height, depth-averaged velocity, bottom topography, and g denotes the gravitational acceleration. The entropy and entropy flux are η(u) = hu 2 +gh 2 Appendix A.4. One-pressure model of spray dynamics

We now consider the one-pressure two-velocity four equations system for modeling the dynamics of a spray of liquid droplets in a gas at thermodynamic equilibrium [START_REF] Sainsaulieu | Ondes progressives solutions de systèmes convectifs-diffusifs et systèmes hyperboliques non conservatifs[END_REF][START_REF] Stewart | Two-phase flow: models and methods[END_REF]. Let ρ g be the gas density, ρ l > 0 the constant and uniform liquid density, α the void fraction of the gas, and u g and u l the velocities of the gas and liquid phases. The variables obey the following hyperbolic system

over the set of states Ω a = {u ∈ R 4 : ρ g > 0, 0 < α < 1}. The gas pressure p = p(ρ g ) satisfies p (ρ g ) > 0, and θ(α) = θ 0 (1 -α) δ , with 1 < δ < 2, where θ 0 denotes the total pressure of the gas on a droplet. The system satisfies an entropy inequality (3) for the pair

where ρ 2 g e (ρ g ) = p(ρ g ) and h(ρ g ) = e(ρ g ) + p(ρ g )/ρ g . It can be checked that the following fluxes are entropy conservative: 

g or α -α + , ρ l (1 -α)u l if ρ - g = ρ + g = ρ g and α -= α + = α.


Fatimata Diaby

Antonio Clemente
email: antonio.clemente@cea.fr

Laurent Dussopt

Ronan Sauleau
email: ronan.sauleau@univ-rennes1.fr

Kien Pham

Erwan Fourn

Design of a 3-Facet Linearly-Polarized Transmitarray Antenna at Ka-band

Keywords: transmitarray, faceted arrays, Ka-band

This contribution presents the design of a 3-facet transmitarray operating at Ka-band. The proposed design methodology and initial numerical analysis are presented. The study is validated by electromagnetic (EM) simulation of a 400element linearly-polarized 3-facet transmitarray antenna.

I. INTRODUCTION

A transmitarray (TA) is typically composed of a planar, faceted or conformal array realized in PCB (Printed Circuit Board) technology and illuminated by a focal source or a focal array (Fig. 1). The phase shift associated to each array element (unit-cell) is tuned to collimate the radiated beam in a specific direction. Various types of TAs have been introduced at millimeter waves, e.g. with high-gain fixed beam [START_REF] Jouanlanne | Wideband linearly-polarized transmitarray antenna for 60 GHz backhauling[END_REF], [START_REF] Miao | 140 GHz high-gain LTCC-integrated transmit-array antenna using a wideband SIW aperture-coupling phase delay structure[END_REF], dualband / dual-polarization [START_REF] Pham | Dual-band transmitarrays with dual-linear polarization at Ka-band[END_REF] or beam-steering capabilities [START_REF] Palma | Circularly-polarized reconfigurable transmitarray in Kaband with beam scanning and polarization switching capabilities[END_REF]. Thanks to their spatial feeding architecture and beamforming capability, TAs are an excellent candidate for the development of innovative low-profile antennas for SATCOM-on-the-move (SOTM) terrestrial terminals at Ka-band. From the array theory, it is well known that element radiation patterns and array geometry have an impact on the scanning capability and performance. As in the case of phased arrays and reflectarrays, the scan range of a TA could be increased using a conformal or faceted radiating aperture. Furthermore, a faceted architecture can enhance the aerodynamics when the antenna is mounted on the fuselage of an aircraft or reduce the visual impact of the antenna. Previous works on reflectarrays concluded that the use of a convex conformal architecture could increase the antenna bandwidth [START_REF] Nayeri | Radiation analysis and characteristics of conformal reflectarray antennas[END_REF], [START_REF] Beccaria | On the bandwidth of convex conformal reflectarray antennas[END_REF].

The objective of this study is to investigate the feasibility of faceted TA and validate the proposed CAD tool through fullwave EM simulations.

II. FACETED TRANSMITARRAY: DESIGN METHODOLOGY

To illustrate the design methodology, let us consider a 3facet TA illuminated by a feed horn, as illustrated in Fig. 1. As well known, the relative phase of the electromagnetic wave radiated by unit-cell m can be written as , [START_REF] Jouanlanne | Wideband linearly-polarized transmitarray antenna for 60 GHz backhauling[END_REF] where is the phase radiated by the focal source towards the direction of unit-cell m, k 0 is the wave number in vacuum at the center frequency f 0 , is the distance between the phase centers of the focal source and the unit-cell m, F is the focal distance, are the coordinates of unit-cell m, is the phase of the radiation pattern of unit-cell m in the incident direction from the focal source, and is the phase shift introduced by unit-cell m. On the other hand, the main beam of a planar phased array can be steered in direction using a simple linear phase distribution across the array aperture, as given by ( 2) 3) is used to define the phase distribution over the TA aperture and scan the radiated beam in 2D over a large angular window. Ideally, a continuous 360° phase tuning range would be desired. Practically, the desired phase law is quantized to make a trade-off between unit-cell complexity (number of dielectrics and metal layers), fabrication tolerances, bandwidth, and aperture efficiency. In this paper, a 3-bit phase quantization is considered.

Table 1 reports the maximum gain of a 3-facet TA and compares its performance with the one of the planar array of same number of unit-cells. The TA is composed of 400 elements arranged on a three-panel 20×20 lattice (the central panel and lateral panels contain 8×20 and 6×20 elements, respectively). Both lateral panels are inclined by an angle α with respect to the central panel. Here 3 values are considered for α: 10°, 30° and 70°. All numerical results have been obtained assuming ideal unit-cells and focal source (cos n θ, n=4) with 100% efficiency. The bandwidth of the different configurations is plotted in Fig. 2. For a strongly faceted array (α = 70° here), the beam scanning range and the bandwidth of the TA are improved, but a gain loss is observed compared to planar configuration.

III. 3-FACET TRANSMITARRAY

A. Linearly-Polarized 3-Bit Unit-Cell

The architecture of the unit-cells used in the study is based on three-metal layers (top patch, ground plane, and bottom patch) printed on two dielectric substrates. Four different architectures (Types I -IV, Fig. 1) have been used to implement the required phase states. As presented in [START_REF] Jouanlanne | Wideband linearly-polarized transmitarray antenna for 60 GHz backhauling[END_REF] for each architecture, a 1-bit unit-cell is designed by physically rotating one patch by 180° around the metallized via used to connect the patch antennas printed on the top and bottom layer of the dielectric structure. The unit-cells exhibit transmission loss lower than 2 dB between 27 GHz and 32 GHz; the transmission phase between two consecutive states is around with a maximum error of 19°.

B. Transmitarray Simulation

Full-wave EM simulations have been performed to validate the theoretical methodology used to model faceted TAs. The 3faceted (α = 10°) TA with 400 elements is based on the linearly-polarized unit-cells described in Section III.A. Comparison between the radiation patterns obtained with the numerical tool (Matlab) and full-wave simulations are presented in Fig. 3, and show excellent agreement. Additional results for different geometries and scanning angles will be presented at the conference.

IV. CONCLUSIONS

In this paper, the principle of faceted transmitarrays has been described in detail. The numerical implementation of the 3-facet transmitarray was validated through 3-D electromagnetic simulations. Good agreement has been obtained between the theoretical and simulated results. We showed that for a certain angle of the faceting, the bandwidth and the scanning performance of the TA are improved at the cost of gain reduction. The objective of this paper, is to present and validate the proposed design methodology. Additional studies considering different array sizes will be presented during the conference.

Fig. 1 .

 1 Fig. 1. Schematic view of the 3-facet TA and corresponding unit-cells.By combining (1) and (2), the required phase shift for each unit-cell can be derived as

Fig. 2 .

 2 Fig. 2. Normalized gain versus frequency of the planar and 3-facet TAs.

Fig. 3 .

 3 Fig. 3. Simulated normalized radiation pattern of the 3-facet TA computed on the E-plane at 29 GHz.

TABLE I .

 I MAXIMUM GAIN OF THE 3-FACET TA AS A FUNCTION OF THE DESIRED SCAN ANGLE

			Gain dBi (actual scan angle)
	Desired scan angle (°)	Planar	3 facets α = 10°	3 facets α = 30°	3 facets α = 70°
	0	29.1 (0°)	29.2 (0°)	29.1 (0°)	27.4 (0°)
	10	29.0 (10°)	29.0 (10°)	28.9 (10°)	27.2 (10°)
	20	28.8 (20°)	28.8(20°)	28.7 (20°)	26.4 (19°)
	30	28.5 (30°)	28.5 (30°)	28.3 (30°)	26.1 (30°)
	40	27.9 (40°)	28.0 (40°)	27.7 (40°)	25.9 (40°)
	50	27.1 (49°)	27.2 (50°)	26.8 (50°)	25.5 (50°)
	60	26.1 (59°)	26.1 (59°)	25.2 (58°)	25.0 (60°)
	70	24.7 (68°)	24.6 (68°)	23.9 (69°)	24.1 (70°)
	80	22.9 (73°)	22.7 (74°)	22.0 (77°)	22.7 (78°)

This works is supported by the French National Research Agency (ANR) in the project named TRANSMIL (TRANSmitarray antenna at MILlimeterwave) under Grant ANR-14-CE28-0023.