R Abgrall

E Le Mélédo

P Öffner

On the Connection between Residual Distribution Schemes and Flux Reconstruction

In this short paper, we are considering the connection between the Residual Distribution Schemes (RD) and the Flux Reconstruction (FR) approach. We demonstrate that flux reconstruction can be recast into the RD framework and vice versa. Because of this close connection we are able to apply known results from RD schemes to FR methods. In this context we propose a first demonstration of entropy stability for the FR schemes under consideration and show how to construct entropy stable numerical schemes based on our FR methods. Simultaneously, we do not restrict the mesh to tensor structures or triangle elements, but rather allow polygons. The key of our analysis is a proper choice of the correction functions for which we present an approach here.

Introduction

We are interested in the approximation of mainly non-linear hyperbolic problems like Euler equations or the MHD equations. The construction of high order methods for these problems is widely studied in current research with the aim to find good ways to build preferable schemes. All of these methods have in common that they are based of either a finite difference (FD) or a finite element (FE) approach. In the last years, great efforts have been made to transform numerical schemes from one to another and to use techniques which are originally used in a different framework. Here, the summation-by-parts (SBP) operators are a good example to mention. SBP operators originate in the FD framework [START_REF] Kreiss | Finite element and finite difference methods for hyperbolic partial differential equations[END_REF] and lead to an ansatz to prove stability in a way similar to the continuous analysis (see [START_REF] Svärd | Review of summation-by-parts schemes for initial-boundary-value problems[END_REF][START_REF] Fernández | Review of summation-by-parts operators with simultaneous approximation terms for the numerical solution of partial differential equations[END_REF][START_REF] Hicken | Multidimensional summation-by-parts operators: General theory and application to simplex elements[END_REF] and the references therein). In [START_REF] Gassner | A skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF] the author transforms the technique to a Discontinuous Galerkin (DG) spectral element method (DGSEM) using the nodes of Lobatto-Legendre quadrature, and in [START_REF] Öffner | Error boundedness of correction procedure via reconstruction / flux reconstruction[END_REF][START_REF] Öffner | Stability of correction procedure via reconstruction with summation-by-parts operators for burgers' equation using a polynomial chaos approach[END_REF][START_REF] Ranocha | Generalised summation-by-parts operators and variable coefficients[END_REF][START_REF] Ranocha | Stability of artificial dissipation and modal filtering for flux reconstruction schemes using summation-by-parts operators[END_REF][START_REF] Ranocha | Summation-by-parts operators for correction procedure via reconstruction[END_REF][START_REF] Ranocha | Extended skew-symmetric form for summation-by-parts operators and varying Jacobians[END_REF] SBP operators are applied to the Correction Procedure via Reconstruction (CPR) or Flux Reconstruction (FR) method to extend stability proofs in a more general framework. In this paper, we also deal with the reinterpretation or transformation of two classes of numerical methods, meaningly the Residual Distribution (RD) schemes and the FR methods. Both lead to a general framework which contains several numerical schemes. The FR creates a unifying framework for several high-order methods such as DG, spectral difference (SD) and spectral volume (SV) methods. These connection are already pointed out in [START_REF] Huynh | High-order methods for computational fluid dynamics: A brief review of compact differential formulations on unstructured grids[END_REF] and in the references therein. However, since the early work of Roe, Deconinick and Struijs [START_REF] Roe | Characteristic-based schemes for the euler equations[END_REF][START_REF] Struijs | Fluctuation splitting schemes for the 2D Euler equations[END_REF], the RD schemes have been further developed in a series of papers, e.g. [START_REF] Abgrall | Toward the ultimate conservative scheme: following the quest[END_REF][START_REF] Abgrall | High order fluctuation schemes on triangular meshes[END_REF][START_REF] Deconinck | Status of multidimensional upwind residual distribution schemes and applications in aeronautics[END_REF][START_REF] Ricchiuto | Application of conservative residual distribution schemes to the solution of the shallow water equations on unstructured meshes[END_REF]. A connection between RD to DG is explained in [START_REF] Abgrall | Development of residual distribution schemes for the discontinuous galerkin method: The scalar case with linear elements[END_REF] and, because of the close relation between RD and FR to DG, it seems natural to study the link between RD and FR. Besides accuracy and robustness of the numerical scheme, another desirable property of numerical schemes is entropy stability. Recently, efforts have been made to construct numerical methods enjoying entropy stability. So far, linear stability remains mainly investigated in the context of FR, as e.g. [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF][START_REF] Castonguay | Energy stable flux reconstruction schemes for advection-diffusion problems[END_REF][START_REF] Jameson | On the non-linear stability of flux reconstruction schemes[END_REF][START_REF] Vincent | An extended range of stablesymmetric-conservative flux reconstruction correction functions[END_REF][START_REF] Wang | A review of flux reconstruction or correction procedure via reconstruction method for the Navier-Stokes equations[END_REF][START_REF] Vermeire | On the properties of energy stable flux reconstruction schemes for implicit large eddy simulation[END_REF][START_REF] Williams | Energy stable flux reconstruction schemes for advection-diffusion problems on tetrahedra[END_REF]. By embedding FR into the RD framework we are able to follow the steps of [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF] and construct FR schemes that are also entropy conservative. The key is a proper choice of the correction functions, which is discussed in this paper. Another beautiful consequence of our abstract approach it that we do not need to restrict our mesh in two dimensions on triangles or tensor structures. Our approach is valid for general polygons, extending the current results. The paper is organized as follows. In the next section we shortly repeat the main idea of RD schemes and FR schemes. We introduce the notations that will be used later in this work. After that we explain the flux reconstruction approach and formulate the schemes in the RD context. Therefore, the definition of the correction functions in FR are essential and two conditions guaranteeing conservation are derived. Furthermore, in the section 4, we transform theoretical convergence and stability results from RD to our FR schemes. Simultaneously, we make some preparations for discussing entropy stability. In the next section 5 we follow the steps from [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF] and construct an entropy conservative/stable numerical scheme based on our FR approach. Simultaneously, by bringing our investigations together we are able to derive straight conditions on our correction functions so that the resulting FR methods are naturally entropy conservative. We then summarize everything and conclude on the admissibility criteria for the correction functions. In the appendix we further extend the investigation on entropy stability and correction functions. Entropy stability is usually associated with the condition of Tadmor on the numerical flux. Following the study of [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF][START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF][START_REF] Abgrall | Some remarks about conservation for residual distribution schemes[END_REF] where similar conditions for RD schemes are proposed and making use of the link between RD and FR, we are able to derive conditions on FR schemes that guarantee entropy stability. It is shown that the correction functions have to fulfill an inequality which is derived from an entropy inequality. However, if this more a theoretical condition can be checked, it does not lead to a construction method for entropy conservative/stable schemes yet. Nevertheless, as it may be interesting for future research, we detail it in the appendix for the sake of completeness.

Residual Distribution Schemes and Flux Reconstruction

Residual Distribution Schemes -Basic Formulation

We follow the ideas and notations for the RD methods of [START_REF] Abgrall | Residual distribution schemes: current status and future trends[END_REF][START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF][START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF][START_REF] Abgrall | Some remarks about conservation for residual distribution schemes[END_REF][START_REF] Abgrall | A high-order nonconservative approach for hyperbolic equations in fluid dynamics[END_REF]. As already described in [START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF], for example, the Discontinuous Galerkin Method (DG) can be interpreted as a RD scheme. Thus, by the close connection between DG and FR (see [START_REF] Huynh | A flux reconstruction approach to high-order schemes including discontinuous Galerkin methods[END_REF]), such an interpretation is also possible for FR schemes.

In this paper we are considering the steady state problem

div f (u) = d j=1 ∂f j ∂x j (u) = 0 for x ∈ Ω ⊂ R d (1)
together with the initial condition

(∇ u f (u) • n(x)) -(u -u b ) = 0 on ∂Ω, (2)
where n(x) is the outward normal vector at x ∈ ∂Ω and u b is a regular enough function. An extension of (1) to unsteady problems is straightforward by following the steps from [START_REF] Abgrall | A high-order nonconservative approach for hyperbolic equations in fluid dynamics[END_REF]. The flux function is given by

f j =    f 1,j . . . f p,j    and u =    u 1 . . . u p    ∈ D ⊂ R p
is the conserved variable.

Later on, we will focus on the entropy U , which fulfils the condition

∇ u U • ∇ u f j = ∇ u g j , ∀j = 1, . . . , d (3)
where

g = (g 1 , • • • , g d ), g j ∈ C 1 (Ω)
and is called the entropy flux function. If u is smooth, then the additional conservation relation

div g(u) = d j=1 ∂g j ∂x j (u) = 0 (4)
holds. If u is a weak entropy solution of ( 1)-( 2), then div g(u) ≤ 0 holds in the sense of distribution. We split Ω in a partition of elements K, and approximate our solution in each element by a polynomial of degree k. The numerical solution is denoted u h . Therefore, the numerical solution lies in the space

V h = K {u h ∈ L 2 (K) d , u h | K ∈ P k (K)}.
In addition, we set {φ σ } σ∈ K a set of basis functions for the space P k (K), where K is a set of degrees of freedom of linear forms acting on the set P k (K), which will be used in every element to express u h . Considering all the elements covering Ω, the set of degrees of freedom is denoted by S and a generic degree of freedom by σ. Furthermore, for any K we have that

∀x ∈ K, σ∈K φ σ (x) = 1.
The key of the RD schemes is to define residuals Φ K σ (u h ) on every element K, satisfying element-wise the following conservation relation

σ∈K Φ K σ (u h ) = ∂K f (u h , u h,-) d γ, (5)
where u h,-is the approximated solution on the other side of the local edge/face of K, f is a consistent numerical flux, i.e. f (u, u) = f (u) • n, and K is the boundary integral evaluated by a numerical quadrature rule. Simultaneously, we have to consider the residuals on the boundary elements Γ. For any degree of freedom belonging to the boundary Γ, we assume that Φ Γ σ (u h ) fulfils the conservation relation

σ∈Γ Φ Γ σ (u h ) = ∂Γ f (u h , u b ) -f (u h ) • n d γ. ( 6
)
The discretisation of (1)-( 2) is given by the following formula. For any σ ∈ S, it reads K⊂Ω, σ∈K

Φ K σ (u h ) + Γ⊂∂Ω, σ∈Γ Φ Γ σ (u h ) = 0. ( 7
)
We are able to embed the discretisation (7) into several numerical methods like finite element or DG depending on the solution space V h and the definition of the residuals, see [START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF] for details. Here, we repeat it shortly for the DG scheme before we extend it to the FR methods.

A weak formulation of DG reads: Find u h ∈ V h such that for any

v h ∈ V h , a(u h , v h ) := K⊂Ω - K ∇v h • f (u h ) d x + ∂K v h • f (u h , u h,-) d γ + Γ⊂∂Ω Γ v h • f (u h , u b ) -f (u h ) • n d γ = 0. (8)
Here we have defined for the boundary faces u h,-= u h and used the fact that the expression ∇v h • f (u h ) implies

∇v h • f (u h ) = d j=1 ∂ ∂x j v h • f j (u h ). ( 9
)
The strong version of DG is obtained by applying integration-by-parts (summation-by-parts in the discrete sense) another time. We obtain the corresponding RD scheme's residuals by comparing [START_REF] Abgrall | High order fluctuation schemes on triangular meshes[END_REF] and [START_REF] Abgrall | Development of residual distribution schemes for the discontinuous galerkin method: The scalar case with linear elements[END_REF]. For the inner elements, we get

Φ K,DG σ (u h ) = - K ∇φ σ • f (u h ) d x + ∂K φ σ • f (u h , u h,-) d γ. ( 10
)
The boundary residuals are given by

Φ Γ,DG σ (u h ) = Γ φ σ • f (u h , u b ) -f (u h ) • n d γ. (11)
Note that the expressions [START_REF] Castonguay | Energy stable flux reconstruction schemes for advection-diffusion problems[END_REF] and [START_REF] Deconinck | Status of multidimensional upwind residual distribution schemes and applications in aeronautics[END_REF] satisfy the conservation relations ( 5) and [START_REF] Abgrall | A high-order nonconservative approach for hyperbolic equations in fluid dynamics[END_REF].

In view of later use, we set an expression for the average of the left and right states of a on the boundary of K and a jump condition that for any function ω as follows.

{a} := 1 2 a K + a K,- [ω] := ω |K --ω |K (12)

Flux Reconstruction Approach on Triangles

In our research we focus on two dimensional problems and target the use general polygonal meshes. Therefore, we start by introducing the FR approach directly on triangles following the explanations of [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF]. For a detailed introduction into FR we strongly recommend the review article [START_REF] Huynh | High-order methods for computational fluid dynamics: A brief review of compact differential formulations on unstructured grids[END_REF] and the references therein. Instead of considering our steady state problem (1), we focus in this subsection on a two dimensional scalar conservation law

∂ t u + div(f ) = 0 ( 13
)
within an arbitrary domain Ω. The flux is now f = (f (u), g(u)) and div is divergence in the space variables x, y. The domain is splitted into N non-overlapping elements Ω k such that Ω = K k=1 Ω k . We use here a conforming triangulation of Ω. Note that a quadrangulation would also be possible. Both u and f are approximated by polynomials in every element Ω k and their total approximation in Ω is given by

u σ = N k=1 u σ k , f σ = N k=1 f σ k ,
where σ represents the DOF in the FR context, i.e. the solution points to evaluate the polyomials. In place of doing every calculation in each element Ω k a reference element Ω s is chosen. Each element Ω k is mapped to the reference element Ω s and all calculation are done in Ω s . The initial equation ( 13) can be transformed to the following governing equation in the reference domain;

∂ t ũ + div( f ) = 0 ⇐⇒ ∂ t ũ + ∇ rs • f = 0, ( 14
)
where div is the divergence in the variables r, s in the computational space. To clarify the notion we use ∇ rs • until the end of this section. ũ, f are the transformed u and f . The quantities ũ and f can be directly calculated using the element mapping. We suppress this dependence in the following to simplify the notation, see [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF] for details. Let P p (Ω S ) denote the space of polynomials of degree less than p on Ω s and R p (Γ S ) be the polynomial space on the edges given by

R P (Γ s ) = v ∈ L 2 (Γ s ), v| Γ f ∈ P p (Γ f ), ∀Γ f
where Γ f stand for the edge f of the reference element Ω s . The approximation ũσ of the solution u within the reference element Ω s is done through a multi-dimensional polynomial of degree p, using the values of ũ at N p =1 2 (p + 1)(p + 2) solution points. The solution approximation then reads;

ũσ (r, s, t) = Np i=1 ũσ i l i (r, s)
where ũσ i is the value 1 of ũ at the solution point i and l i (r, s) is the multidimensional Lagrange polynomial associated with the solution points i in the reference element Ω s . We now detail the main idea of FR. A simple approach is to also approximate the flux function f by a polynomial f σ = ( f σ , gσ ). To build this approximation, a first polynomial decomposition of f σ is set as

f σ,D = Np i=1 f σ i l i (r, s), gσ,D = Np i=1 gσ i l i (r, s), (15)
where the coefficients of these polynomials, respectively denoted by f σ i and gσ i , are again evaluated at the solution points. Since f σ,D is always discontinuous at the boundary, it is called discontinuous flux. To overcome / reduce that problem a further term f σ,C is then added to f σ,D . f σ,C is set to work directly on the boundaries of each element and corrects f σ,D such that information of two neighboring elements interacts and properties like conservation still hold in the discretisation. We obtain for the approximation of the flux;

f σ = f σ,D + f σ,C .
This gives Flux Reconstruction its name. The selection/definition of these correction functions is essential. We now detail a possible construction for our special case. On each edge of the triangle, a set of N fp = p + 1 flux points are defined. These flux points are applied to couple the solution between neighboring elements. The correction function is then constructed as follows

f σ,C = 3 f =1 N fp j=1 f • n σK f,j -f σ,D • n f,j h f,j (r, s) (16)
The indices f, j correspond to a quantity at the flux point j of face f . Thus, in our case it is 1

≤ f ≤ 3 and 1 ≤ j ≤ N fp . The term f σ,D • n f,j
is the normal component of the transformed discontinuous flux at the flux point j, j, whereas f • n σK f,j is a normal transformed numerical flux computed at flux point f, j. We compute it by evaluating the multiply defined values of u σ at each flux point. More precisely, we first define by u σ,-the value of u σ computed in the current element and by u σ,+ its value computed using the information from the adjoint element that shares the same flux point j. This couples two neighboring elements and the information between them. We then evaluate u σ,+ and u σ,-at each flux point and compute f (u σ,-, u σ,+ ) • n σK f,j

. Finally, h f,j (r, s) has to be explained. This is a vector correction function associated with the flux points f, j and that lie in the Raviart-Thomas space RT p (Ω s ) of order p. Thus h fulfills the following two properties

∇ rs • h f,j ∈ P p (Ω s ) h f,j • n ∈ R p (Γ S ) (17)
and has also to satisfy

h f,j (r f2,j2 ) • n f2,j2 = 1 if f = f 2 and j = j 2 , 0 if f = f 2 or j = j 2 . ( 18
)
Because of ( 18) it follows that

f σ,C (r f,j ) • n f,j = f • n σK f,j -f σ,D • n f,j =: α f,j .
We also get f σ (r f,j )

• n f,j = f • n σK f,j
at each flux point f, j. Combining our results, the approximate solution values to the problem ( 13) can be updated at the solution points from

d u d t = -∇ rs • f σ ri = -∇ rs • f σ,D ri -∇ rs • f σ,C ri = - Np k=1 f σ k ∂l k (r, s) ∂r ri - Np k=1 gσ k ∂l k (r, s) ∂s ri - 3 f =1 N fp j=1 f • n σK f,j -f σ,D • n f,j ∇ rs h f,j (r) = - Np k=1 f σ k ∂l k (r, s) ∂r ri - Np k=1 gσ k ∂l k (r, s) ∂s ri - 3 f =1 N fp j=1 α f,j ∇ rs h f,j (r).
Defining our FR scheme reduces to select the distributions of flux points and solutions points, as well as the form of our correction functions. The choice leads to several numerical methods with different properties.

In [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF] special attention is paid on conservation and linear stability which restricts again the set of correction functions, but we do not go further into details here. Finally, we want to mention that some of the most famous schemes are embedded in this framework by a right choice of correction functions and point distributions. To give a concrete example, the nodal Discontinuous Galerkin Spectral Element Method of Gassner et al. [START_REF] Gassner | A comparison of the dispersion and dissipation errors of Gauss and Gauss-Lobatto discontinuous Galerkin spectral element methods[END_REF][START_REF] Gassner | A skew-symmetric discontinuous Galerkin spectral element discretization and its relation to SBP-SAT finite difference methods[END_REF] can be named.

div f h + f -f h • n ∇ψ = 0. ( 20
)
We derive conditions on ∇ψ so that this approach fits in the RD framework and that our methods have the desirable properties of conservation and stability. First, let us focus on FR schemes. The main idea of the FR schemes is that the numerical flux at the boundaries will be corrected by functions in such manner that information of two neighboring elements interacts and properties like conservation hold also in their discretisations. Let us consider our discretisation [START_REF] Mengaldo | Dealiasing techniques for high-order spectral element methods on regular and irregular grids[END_REF]. If we apply a Galerkin approach in every element K, then we obtain that for any v h ∈ V h the relation

K v h • div f h + f (u h , u h,-) -f h • n ∇ψ d x = 0 (21)
has to be fulfilled. Using the Gauss theorem in the above equation yields

- K ∇v h • f h + α∇ψ d x + ∂K v h • f h • n + f -f h • n ∇ψ • n d γ = 0. ( 22
)
To guarantee conservation we demand that the flux over the element boundaries should be expressed only by the numerical flux of elements sharing this boundary. Therefore, we require that

f h • n + f (u h , u h,-) -f h • n ∇ψ • n = f (u h , u h,-), which implies ∇ψ • n ≡ 1 (23)
on the boundary. The relation [START_REF] Ranocha | Generalised summation-by-parts operators and variable coefficients[END_REF] yields us a first property on our correction function ∇ψ.

Remark 3.1. The condition (23) can be further weaken. Since we are using quadrature rules to evaluate the integrals (21) or ( 22), ( 23) has to be be only fulfilled at the quadrature points. We can guarantee the property [START_REF] Ranocha | Generalised summation-by-parts operators and variable coefficients[END_REF] in two dimension (d = 2) when using functions ∇ψ lying in the lowest order Raviart-Thomas space [START_REF] Raviart | A mixed finite element method for 2-nd order elliptic problems[END_REF], up to some scaling. The relation [START_REF] Ranocha | Generalised summation-by-parts operators and variable coefficients[END_REF] is then automatically fulfilled as a basic property of this function space. In [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF] the authors already considered the Raviart-Thomas elements focusing on triangles. We are considering the more general case of polygons.

To demonstrate the connection between FR and RD and to build an numerical scheme we have to apply quadrature formulas to evaluate the continious integrals.

From Flux Reconstruction to Residual Distribution Schemes

In this part of the paper we show the connection between Flux Reconstruction and Residual Distribution Schemes. The key is a proper definition of the residuals. If one has again a look on the formulation of the residuals ( 8)-( 11) from section 2.1 and compare them now with the formulations of ( 21)-( 22), it can be noticed that the equations share similar structures. By passing from integrals to quadrature formulas and splitting v h along {φ σ } σ∈S , we can define the residuals in the following manner.

Φ K,F R σ (u h ) := K φ σ • div f h + f (u h , u h,-) -f h • n ∇ψ d x (24)
An other approach is to use Gauss formula (integration-by-parts/summation-by-parts), leading to

Φ K,F R σ (u h ) := - K ∇φ σ • f h d x + ∂K φ σ α∇ψ • n + f h • n d γ - K ∇φ σ • α∇ψ d x (23) = - K ∇φ σ • f h d x + ∂K φ σ f (u h , u h,-) d γ - K ∇φ σ • α∇ψ d x :=rσ . ( 25
)
Recalling property [START_REF] Ranocha | Generalised summation-by-parts operators and variable coefficients[END_REF], the boundary residuals reads

Φ K,F R σ (u h ) = Γ φ σ f (u h , u b ) -f h • n d γ. ( 26
)
Comparing the residuals ( 24) and ( 25) with the residuals (10) of the DG scheme2 , we can write

Φ K,F R σ (u h ) = Φ K,DG σ (u h ) + r σ . (27)
Furthermore, the conservation relation ( 5) directly provides a second property on ∇ψ,

explicitly σ∈K r σ = - σ∈K K ∇φ σ • α∇ψ d x = 0. ( 28
)
If we apply the residuals ( 24)-( 26) on our underlying steady state problem ( 1)-( 2), we are able to write our model problem in the shape of [START_REF] Abgrall | High order fluctuation schemes on triangular meshes[END_REF]. For any σ ∈ S, it reads

K⊂Ω, σ∈K Φ K,F R σ (u h ) + Γ⊂∂Ω, σ∈Γ Φ Γ,F R σ (u h ) = 0. ( 29
)
With the definitions of the residuals ( 24)-( 27) and the discretisation (29), the Flux Reconstruction is embedded within the RD framework. By ensuring that conditions ( 23) and ( 28) hold, the conservation relation (5) for our residuals is guaranteed and we are now able to use the theoretical results of RD [START_REF] Abgrall | Residual distribution schemes: current status and future trends[END_REF][START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF][START_REF] Abgrall | Some remarks about conservation for residual distribution schemes[END_REF][START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF][START_REF] Abgrall | A high-order nonconservative approach for hyperbolic equations in fluid dynamics[END_REF] for the FR schemes under consideration. Naturally, the conservation properties of Flux Reconstruction schemes also hold and also stability results will transfer from the RD framework to FR.

Remark 3.2. If we are considering a two dimensional problem (1), our approach does not restrict the splitting of the domain Ω to a specific geometric structure like triangles or rectangles. The results are valid more generally for all polygons. This approach then extends the results of [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF][START_REF] Huynh | High-order methods for computational fluid dynamics: A brief review of compact differential formulations on unstructured grids[END_REF][START_REF] Mengaldo | Dealiasing techniques for high-order spectral element methods on regular and irregular grids[END_REF][START_REF] Vincent | A new class of high-order energy stable flux reconstruction schemes[END_REF] on FR to general grids.

Before we focus on entropy stability of our FR methods we shortly repeat some well-known results of RD schemes from [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF][START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF][START_REF] Abgrall | Some remarks about conservation for residual distribution schemes[END_REF] and references therein.

Transformation Results to Flux Reconstruction

As it is described inter alia in [START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF] for the RD schemes, a generalization of the classical Lax-Wendroff theorem is valued. It transfers naturally to our FR formulation in RD ( 24)- [START_REF] Raviart | A mixed finite element method for 2-nd order elliptic problems[END_REF].

Theorem 4.1 (Theorem 2.2 of [START_REF] Abgrall | High order fluctuation schemes on triangular meshes[END_REF]). Assume the family of meshes T = (T h ) is shape regular. We assume that the residuals {Φ K σ } σ∈K for K an element or a boundary element of T h satisfy: • For any M ∈ R + , there exists a constant C which depends only on the family of meshes T h and M such that for any

u h ∈ V h with ||u h || ∞ ≤ M , then ||Φ K σ (u h | K )|| ≤ C σ,σ ∈K ||u h σ -u h σ ||.
• The conservation relations (5) and (6) hold.

If there exists a Constant C max such that the solutions of the scheme (7) (or (29)) satisfy ||u h || ∞ ≤ C max and a function v ∈ L 2 (Ω) d such that (u h ) h or at least a sub-sequence converges to v in L 2 (Ω) d , then v is a weak solution of (1).

In view of the proof, the following relation is essential and can be derived from the conservation relations ( 5) and [START_REF] Abgrall | A high-order nonconservative approach for hyperbolic equations in fluid dynamics[END_REF]. For any v h ∈ V h written as

v h = σ∈S v σ φ σ : 0 = - Ω ∇v h • f (u h ) d x + ∂Ω v h f (u h , u b ) -f (u h ) • n d γ + e∈E h e {v h } f (u h , u h,-) d γ + K∈Ω 1 #K   σ,σ ∈K (v σ -v σ ) Φ K σ (u h ) -Φ K,DG σ (u h )   + Γ∈∂Ω 1 #Γ   σ,σ ∈Γ (v σ -v σ ) Φ Γ σ (u h ) -Φ Γ,DG σ (u h )   , (30)
with Φ •,DG σ as in [START_REF] Castonguay | Energy stable flux reconstruction schemes for advection-diffusion problems[END_REF] and [START_REF] Deconinck | Status of multidimensional upwind residual distribution schemes and applications in aeronautics[END_REF]. A consequence of ( 30) is the following entropy inequality: Proposition 4.2 (Proposition 3.2 from [START_REF] Abgrall | Some remarks about conservation and entropy stability for residual distribution schemes[END_REF]). Let (U, g) be an couple entropy-flux for (1) and ĝ a numerical entropy flux consistent with g • n. Assume that the residuals satisfy:

σ∈K ∇ u U (u σ ), Φ K σ ≥ ∂K ĝ(u h , u h,-) d γ for any element K, σ∈e ∇ u U (u σ ), Φ e σ ≥ ∂e ĝ(u h , u b ) -g(u h ) • n) d γ for any boundary edge e. (31)
Under the assumption of the theorem 4.1, the limit weak solution satisfies the following entropy inequality:

For any τ ∈ C 1 (Ω), τ ≥ 0, - Ω ∇τ • g(u) d x + ∂Ω - τ g(u b ) • n d γ ≤ 0.
The theorem 4.1 and the proposition 4.2 ensure entropy stability. Therefore, assumption (31) is essential. In [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF], the property (31) is further analysed and compared with the theory of Tadmor about entropy conservative/stable numerical flux functions [START_REF] Tadmor | The numerical viscosity of entropy stable schemes for systems of conservation laws[END_REF][START_REF] Tadmor | Entropy stability theory for difference approximations of nonlinear conservation laws and related time-dependent problems[END_REF]. We can transfer this investigation, yielding a further condition on the correction functions. One can show that they have to satisfy an inequality in some sense to guarantee entropy stability. However, this is more a theoretical condition and until now we do not see how this helps us to construct entropy conservative/stable FR schemes. For reasons of completeness and for future perspectives we develop this in the appendix. Another approach will be followed in the next section.

Entropy Conservative/Stable Flux Reconstruction Schemes

We show in this section how to construct an entropy conservative scheme staring from our FR schemes Φ K,F R σ as the defined in ( 25) and [START_REF] Raviart | A mixed finite element method for 2-nd order elliptic problems[END_REF]. From now on v represents the entropy variable ∇ u U (u). Note that since the entropy is strictly convex, the mapping u → v(u) is one-to-one. Here, we concentrate only on the the inner elements. For a detailed analysis and the study of the boundary we strongly recommend [START_REF] Abgrall | Some remarks about conservation for residual distribution schemes[END_REF]. We give further conditions on our correction functions in [START_REF] Raviart | A mixed finite element method for 2-nd order elliptic problems[END_REF] to get an entropy conservative/stable numerical FR scheme and presented a way to construct entropy stable FR schemes. This is done for the first time whereas [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF][START_REF] Vincent | A new class of high-order energy stable flux reconstruction schemes[END_REF][START_REF] Vincent | An extended range of stablesymmetric-conservative flux reconstruction correction functions[END_REF] derive only linear stability. Now, let ĝ be the associated numerical entropy flux for the entropy variable v consistent with g • n. As it is described in [START_REF] Abgrall | Some remarks about conservation for residual distribution schemes[END_REF] the entropy conditions can be formulated in the RD setting for every inner element by the following formula.

σ∈K v σ , Φ CS σ = ∂K ĝ(v h , v h,-) d γ (32)
To construct a new scheme Φ CS σ which fulfils the condition (32), we consider

Φ CS σ = Φ K,F R σ + τ σ (27) = Φ K,DG σ + r σ + τ σ , (33)
with some τ σ that has to be built. In order to guarantee the conservation property, we ask

σ∈K (Φ K,F R σ + τ σ ) = ∂K f (v h , v h,-) d γ, (34)
which implies

σ∈K τ σ = 0. ( 35
)
The question is now to construct τ σ under this constraint. When (32) holds, we have

σ∈K v σ , τ σ = ∂K ĝ(v h , v h,-) d γ - σ∈K v σ , Φ K,F R σ := E. ( 36
)
Thus, a solution to ( 35)-( 36) is given by:

τ σ = α σ (v σ -v), α σ = E σ∈K (v σ -v) 2 , v = 1 #K σ∈K v σ . ( 37
)
This can be seen by the following short calculation

σ∈K τ σ = σ∈K α σ (v σ -v) (37) = σ∈K α σ v σ - σ∈K α σ v σ = 0, E = σ∈K v σ , τ σ = σ∈K v σ , α σ (v σ -v) = σ∈K α σ ( v σ , v σ -v σ , v σ ) = σ∈K α σ ( v σ , v σ -2 v σ , v σ + v σ , v σ ) = σ∈K α σ (v σ -v) 2
The scheme [START_REF] Tadmor | Entropy stability theory for difference approximations of nonlinear conservation laws and related time-dependent problems[END_REF] is entropy stable by construction, but one can wonder about its accuracy. We are using the second approach presented in [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF], since we can preserve the accuracy. The entropy flux g(u) ⊂ R d and the normal flux f fulfill the following relation.

3 g = v, f -θ with ∇ v θ = f (38)
The crucial point is the error E which has to be approximated as accurately as possible. For simplicity reason we are considering

Ẽ := -E = - ∂K ĝ(v h , v h,-) d γ + σ∈K v σ , Φ K,F R σ = - ∂K ĝ(v h , v h,-) d γ + σ∈K v σ , Φ K,DG σ + r σ (39)
The entropy error Ẽ for DG was already investigated in [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF]. We follow the steps which are analogous except for the additional terms r σ . The numerical entropy ĝ will be defined later on. We hide the dependence of f h on v h in the following to simplify the notation. We further need the following two relations. Using [START_REF] Castonguay | A new class of high-order energy stable flux reconstruction schemes for triangular elements[END_REF] we have:

∇v h • f h = div v h , f h -v h • ∇f h , (40)
which can be also written in components. It would read

d j=1 p i=1 ∂ ∂x j v h i f i,j = d j=1 ∂ ∂x j p i=1 (v h i f i,j ) - d j=1 p i=1 v h i ∂ ∂x j f i,j .
We also get from ( 3) and ( 4)

v(v h ) • ∇f (v h ) = d j=1 v(v h ) • ∂ ∂x j f j (v h ) = d j=1 ∂ ∂x j g j (v h ) = div g(v h ). ( 41
)
Here, v(v h ) is used to emphasis that condition (3) is only valued for the entropy variable v and we can not assume directly that it holds for its interpolation. Therefore, we have to use the flux itself in this context and not the interpolated one. Since we use the approximated flux function f h in our FR schemes we get a slight different investigation as in [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF]. However it does not change the major result. With (40), ( 41) and Gauss's theorem 4 we are able to rewrite (39) as

Ẽ (40) = - K div v h , f h d x + K v h • ∇f h d x + ∂K v h , f (v h , v h,-) -ĝ(v h , v h,-) d γ + σ∈K v σ , r σ = - K div v h , f h d x + K div v h , f h d x :=SU R1 - K div v h , f h d x + K v h -v(v h ) • ∇f (v h ) d x :=SU R2 + K v h • ∇ f h -f (v h ) d x :=SU R3 + K v(v h ) • ∇f (v h ) d x - K ∇v h • (αψ) d x :=CO + ∂K v h , f (v h , v h,-) -ĝ(v h , v h,-) d γ Gauss&(41) = SU R 1 + SU R 2 + SU R 3 + CO - ∂K v h , f (v h ) • n d γ + ∂K v h , f (v h ) • n d γ :=BO - ∂K v h , f (v h ) • n d γ + ∂K g(v h ) • n d γ + ∂K v h , f (v h , v h,-) -ĝ(v h , v h,-) d γ = SU R 1 + SU R 2 + CO + BO + ∂K -v h , f (v h ) • n + v h , f (v h , v h,-) + g(v h ) • n -ĝ(v h , v h,-) d γ
Assuming a smooth exact solution, we can use a quadrature formula of order k for the three first volume terms (SU R 1 , SU R 2 , SU R 3 ) and obtain

SU R 1 = O(h k+d+1 ), SU R 2 = O(h k+d+1 ), SU R 3 = O(h k+d+1 ). ( 42
)
For the boundary term BO, we can get BO = O(h k+d+1 ) using a quadrature formula of order k + 1. The last term has to be investigated. We use here for the numerical entropy flux

ĝ(v h , v h,-) = {v h }, f (v h , v h,-) -θ({v h }). (43)
Applying [START_REF] Williams | Energy stable flux reconstruction schemes for advection-diffusion problems on tetrahedra[END_REF] and (43) in the last term yields after some calculations to

-v h , f (v h ) • n + v h , f (v h , v h,-) + g(v h ) • n -ĝ(v h , v h,-) = O(h 2(k+1) )
and in total

∂K -v h , f (v h ) • n + v h , f (v h , v h,-) + g(v h ) • n -ĝ(v h , v h,-) d γ = O(h d+2k+1 )
with a suitable quadrature formula. We now consider the last element term CO, which has to approximate zero at least as O(h k+d+1 ). Indeed, if we assume that a quadrature of order k for the volume integrals leads to CO = O(h d+k+1 ), then we can merge the errors and retrieve

Ẽ = σ∈K v σ , Φ σ - ∂K ĝ(v h , v h,-) d γ = O(h d+k+1 ), (44)
provided that a quadrature formula of order k for the volume integrals and order k + 1 for the boundary integrals are applied. This last equality furnishes a new condition for admissible correction functions and/or their quadrature rule. This yields us to an extension of Proposition 3.3.

Remark 5.1. To obtain an entropy stable scheme, i.e. to have

σ∈K v σ , Φ ST σ ≥ ∂K ĝ(v h , v h,-) d γ, (45)
we can combine the above approach and add an additional term to Φ CS σ . More explicitly we can set:

Φ ST σ = Φ CS + Ψ σ ,
where Ψ σ satisfy σ∈K Ψ σ = 0 and σ∈K v σ , Ψ σ ≥ 0. In [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF] two expressions for Ψ can be found. One contains jumps and the other one streamlines.

However, we are interested here in entropy conservative/stable FR schemes. Therefore, the inequality (45) should already hold for FR schemes of type [START_REF] Ranocha | Summation-by-parts operators for correction procedure via reconstruction[END_REF]. Let E DG be the entropy error [START_REF] Vincent | An extended range of stablesymmetric-conservative flux reconstruction correction functions[END_REF] of DG. Up to assuming

σ∈K v σ , r σ ≥ E DG , (46)
we get automatically

σ∈K v σ , Φ K,F R σ - ∂K ĝ(v h , v h,-) d γ = σ∈K v σ , r σ - ∂K ĝ(v h , v h,-) d γ - σ∈K v σ , Φ K,DG σ E DG ≥ 0.
Thus, we are now able to present a way to build functions such that (32) is naturally fulfilled. To simplify our notation we include the scaling term α in our correction functions ∇ψ and use from now on ∇ ψ := α∇ψ. Therefore, conditions ( 23) and ( 28) transfer to

∇ ψ • n = f (u h , u h,-) -f h • n (23*) σ∈K r σ = - σ∈K K ∇φ σ • ∇ ψ d x = 0. ( 28*)
For the construction of an entropy conservative scheme, we first add the term τ to our FR schemes. A solution of (37) is then obtained. A way to determine a conservation preserving τ is to read its formulation directly from the the ansatz (37) used together with the above mentioned conditions, especially (28*). On the one side we have

r σ = K ∇φ σ • ∇ ψ d x,
but on the other side we also ask the entropy condition [START_REF] Tadmor | The numerical viscosity of entropy stable schemes for systems of conservation laws[END_REF]. In order to fulfill this last condition, we may write r σ as

r σ = α σ (E) (v -v) ,
where α σ (E) is defined as in [START_REF] Wang | A review of flux reconstruction or correction procedure via reconstruction method for the Navier-Stokes equations[END_REF] and where we emphasized the dependence on the entropy error E. Bringing the two relations together we define our correction function by solving the following discrete Neumann problem.

K ∇φ σ • ∇ ψ d x = α σ (E) (v -v) , ∇ ψ • n = f (u h , u h,-) -f h • n (47)
The equation ( 32) is then fulfilled and we naturally get entropy conservation for our FR schemes. We can add jump or streamline terms to the schemes, as described in remark 5.1 or more explicitly in [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF], to get entropy stability.

Summary

In this short paper we demonstrated the connection between Residual Distribution Schemes and the Flux Reconstruction approach. We saw that the FR schemes can be written as RD schemes. This link enables us to transform the well-known results about RD schemes into the FR framework. The crucial point is to derive suitable correction functions in FR. From our previous analysis we can formulate the following main result.

Theorem 6.1. We are considering the steady state problem (1) together with the initial condition (2). We are using a FR approach div f h + α∇ψ = 0. ( 20)

with α = f -f h • n .
Let us further assume that our correction functions α∇ψ fulfill the following two conditions

• ∇ψ • n ≡ 1 ( 23
)
• σ∈K r σ = - σ∈K K ∇φ σ • α∇ψ d x = 0. ( 28
)
Then, our FR schemes are conservative and we are able to recast them into the RD framework. All the results of section 4 apply automatically to our FR schemes.

If we further choose our correction functions so that

σ∈K v σ , r σ ≥ E DG , ( 46
)
where E DG is the entropy error (36) of DG, then our FR scheme is additionally entropy stable. Finally, the correction functions are determined by solving the discrete Neuman problem (47), leading to an entropy conservative FR scheme.

A naturally arising question is how to select this correction functions so that the conditions are fulfilled. As explained in remark 3.1, we search our correction functions ψ in the Raviart-Thomas space. Solving then the discrete Neumann problem (47) selects entropy conservative FR schemes.

In [START_REF] Vincent | An extended range of stablesymmetric-conservative flux reconstruction correction functions[END_REF] Vincent et al. already described and developed FR schemes on triangular grids. Their main problem was to describe the correction functions. They used the Raviart-Thomas space and could prove several conditions on their schemes focusing on linear stability. They used in their analysis/construction a direct representation with flux points and solutions points on triangles whereas we apply an abstract approach.

Our advantage is that we do not use the geometrical structure of the grid. Therefore, our results are valid for general polygons (in two space dimensions), and include the schemes of [START_REF] Vincent | An extended range of stablesymmetric-conservative flux reconstruction correction functions[END_REF]. However, this paper deals with the interpretation/transformation of FR schemes into the RD framework, allowing the use of theoretical results specific to RD in the context of FR. We are considering especially entropy stability and derive conditions which are linked to the selection of correction functions. We further give an idea to build entropy conservative FR schemes. In a forthcoming paper we construct these FR schemes for general polygons and test them numerically against some benchmarks.

Appendix

7.1 Entropy Stability-An Approach in the Sense of Tadmor

In this appendix we follow the steps from [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF] and have a closer look on the property (31) for our FR residuals ( 24)- [START_REF] Ranocha | Extended skew-symmetric form for summation-by-parts operators and varying Jacobians[END_REF]. For simplicity reasons we assume in the following that K is a fixed triangle. The results are nevertheless extendible to general polytopes with degrees of freedom on the boundary of K. We may consider a triangulation T K of K elements whose vertices are exactly the elements of S. We denote by fσ,σ the flux between two DOFs σ and σ and by n σ,σn the normal vector on the direct edge between σ and σ n . In [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF] it is shown that we can split the residual in the following way;

Φ K σ = edges[σ,σ ] fσ,σ + f b σ , (48)
with f b σ = ∂K φ σ fσ (u h , u h,-) d γ. Further properties and a detailed analysis can be found in [START_REF] Abgrall | A general framework to construct schemes satisfying additional conservation relations. application to entropy conservation and entropy dissipative schemes[END_REF]. As an additional example we derive the flux fσ,σ for FR.

Example 7.1 (Flux Reconstruction schemes in the P 1 case). The residuals are simply

Φ K,F R σ (u h ) = - K ∇φ σ f h + α∇ψ d x + ∂K φ σ f (u h , u h,-) d γ
The flux between two DOFs σ and σ is given by

fσ,σ (u h , u h,-) = ∂K (φ σ -φ σ ) f (u h , u h,-) d γ - K ∇(φ σ -φ σ ) • f h + α∇ψ d x,
where the equality ∇(φ σφ σ ) = n σσ |K| comes from the geometry (See figures 1 and 2 ). We obtain

fσ,σ (u

h , u h,-) = ∂K (φ σ -φ σ ) f (u h , u h,-) d γ - K f h + α∇ψ d x • n σσ |K| . ( 49
)
Remark 7.2. If the same quadrature formula is used on each edge, then ∂K (φ σφ σ ) d γ = 0 and by introducing the cell average u h of u h and the average f h of the flux f h , we can rewrite the flux as

fσ,σ (u h , u h,-) = ∂K (φ σ -φ σ ) f (u h , u h,-) -f h • n d γ - K f h + α∇ψ d x • n σσ |K| ,
and the first term can be interpreted as a dissipation.

In Tadmors's work [START_REF] Tadmor | The numerical viscosity of entropy stable schemes for systems of conservation laws[END_REF][START_REF] Tadmor | Entropy stability theory for difference approximations of nonlinear conservation laws and related time-dependent problems[END_REF] about entropy stability, sufficient conditions are derived for the schemes to be entropy stable/conservative. In particular, it is shown that the numerical entropy flux needs to fulfill some dissipation inequality. We deduce here equivalent conditions for our setting, paying attention to the correction function that plays an important role. Using equation (48) we can write our FR residuals as

Φ K,F R σ = edges[σ,σ ] f F R σ,σ + f b,F R σ , ( 50
) with f b,F R σ = ∂K φ σ f (u h , u h,-) d x.
From proposition 4.2 we know that the condition ( 31) is sufficient to guarantee entropy stability. Let us now insert (50) into it. We obtain

5 σ∈K v σ , Φ K,F R σ = σ∈K v σ , [σ,σ ] f F R σ,σ + f b,F R σ = σ∈K v σ , [σ,σ ] f F R σ,σ + σ∈K v σ , ∂K φ σ f (v h , v h,-) d x = 1 2 σ>σ v σ -v σ , f F R σ,σ + ∂K v h , f (v h , v h,-) d γ (51)
We recall that we want to tune the correction function such that the obtained scheme is stable. It reduces here to ask

σ∈K v σ , Φ K,F R σ ≥ ∂K ĝ(v h , v h,-) d γ. (52)
We now derive a condition that guarantees (52) to hold. Introducing the potential Θ h K in K by

Θ h K := σ∈K θ σ φ σ with θ σ = v σ , f (v σ ) -g(v).
where we used the fact that the entropy flux g is related to the flux f by g = v, f -Θ, and defining the numerical entropy flux ĝ in this spirit yields to the following expression.

ĝ(v h , v h,-) := {v h }, f (v h , v h,-) -{Θ h K } • n
Together with [START_REF] Fernández | Review of summation-by-parts operators with simultaneous approximation terms for the numerical solution of partial differential equations[END_REF] we may rewrite (51) as:

1 2 σ>σ v σ -v σ , f F R σ,σ + ∂K v h , f (v h , v h,-) d γ ≥ ∂K {v h }, f (v h , v h,-) -{Θ h K } • n d γ ⇐⇒ 1 2 σ>σ v σ -v σ , f F R σ,σ + ∂K v h , f (v h , v h,-) d γ + ∂K Θ h K • n d γ + 1 2 ∂K [Θ h K ] • n d γ ≥ ∂K {v h }, f (v h , v h,-) d γ ⇐⇒ 1 2 σ>σ v σ -v σ , f F R σ,σ + ∂K Θ h K • n d γ + 1 2 ∂K [Θ h K ] • n d γ - 1 2 ∂K [v h ], f (v h , v h,-) d γ ≥ 0
Rearranging the last inequality we get:

1 2 σ>σ v σ -v σ , f F R σ,σ + ∂K Θ h K • n d γ :=C K - 1 2 ∂K [v h ], f (v h , v h,-) -[Θ h K ] d γ :=B ∂K ≥ 0 ( 53
)
Taking now a closer look on our inequality (53), we may split this condition for entropy stability in two parts. One part is actually working on the boundary of K, denoted by B ∂K , and one in K, denoted by C K .

5 Here, we considered directly the entropy variable v similar to [START_REF] Tadmor | The numerical viscosity of entropy stable schemes for systems of conservation laws[END_REF] instead of u. Furthermore, we deal with an oriented graph. Given two vertices of this graph σ and σ , we write σ > σ a direct edge and we shorten the notation by Let's pull a connection to the work of [START_REF] Tadmor | The numerical viscosity of entropy stable schemes for systems of conservation laws[END_REF]. We recall that a numerical flux f is entropy stable in the sense of Tadmor if

[v h ], f (v h , v h,-) -[Θ h K ] ≤ 0 (54)
holds. Thus combining (54) and the use of an entropy stable numerical flux in the sense of Tadmor, we can guarantee that B ∂K ≤ 0. We have only left to consider C K . Let us study its first term. We have:

1 2 σ>σ v σ -v σ , f F R σ,σ = 1 2 σ>σ v σ -v σ , f F R σ,σ -(θ σ -θ σ ) • n σ,σ + 1 2 σ>σ (θ σ -θ σ ) • n σ,σ .
(θ σ -θ σ ) • n σ,σ = σ∈K θ σ • N σ = - σ∈K θ σ ∂K φ σ • n d γ = - ∂K Θ h K • n d γ. (56)
By (56) we get for C K :

1 2 σ>σ v σ -v σ , f F R σ,σ -(θ σ -θ σ ) • n σ,σ + ∂K Θ h K • n d γ - ∂K Θ h K • n d γ ≥ 0
Thus we finally get:

1 2 σ>σ v σ -v σ , f F R σ,σ -(θ σ -θ σ ) • n σ,σ ≥ 0. ( 57
)
This condition is analogous to the one of Tadmor about entropy stability. Therefore, for an entropy stable flux function f F R σ,σ , the entropy stability is guaranteed. However, note that f F R σ,σ depends on the correction functions. We show it on our example 7.1.

Example 7.3. Putting (49) into (57) yields to the condition

1 2 σ>σ v σ -v σ , f F R σ,σ -(θ σ -θ σ ) • n σ,σ = 1 2 σ>σ v σ -v σ , ∂K (φ σ -φ σ ) f (u h , u h,-) d γ - K f h + α∇ψ d x • n σσ |K| -(θ σ -θ σ ) • n σ,σ = 1 2 σ>σ v σ -v σ , ∂K (φ σ -φ σ ) f (u h , u h,-) d γ - K f h d x • n σσ |K| -(θ σ -θ σ ) • n σ,σ - 1 2 σ>σ v σ -v σ , K α∇ψ d x • n σσ |K| ≥0
The first line of the last expression is again the same as for DG whereas the last line gives us special conditions on the correction functions which can be tuned such that this inequality holds. Thus, the first step is the selection of α∇ψ such that

1 2 σ>σ v σ -v σ , K α∇ψ d x • n σσ |K| ≤ 0.
In the above example we have seen that the correction function has a direct influence on entropy stability.

We can now further restrict our correction functions so that the inequality (57) holds. A detailed analysis as well as more examples of those restrictions be found in [START_REF] Abgrall | Some remarks about conservation for residual distribution schemes[END_REF].

Remark 7.4 (Extension of Theorem 6.1). Our investigation yields us an additional condition, so that we can extend our Theorem 6.1 by the the following: If we further choose our correction functions so that f F R σ,σ is entropy stable in the sense of Tadmor (57) then our FR scheme is additionally entropy stable.

Figure 1 :

 1 Figure 1: Notations for the schemes. On the left: Definition of the control volume for the degree of freedom σ. The vertex σ plays the role of the vertex 1 on the left picture for the triangle K. The control volume C σ associated to σ = 1 is green on the right and corresponds to 1IGK on the left. The vectors n ij are normal to the internal edges scaled by the corresponding edge length.

1 Figure 2 :

 12 Figure 2: Representation of the control volume associated to DOF 1.

Furthermore, we have 1 2 σ>σ(

 2 θ σθ σ ) • n σ,σ = σ∈K θ σ • N σ ,(55)where the j-th component of N writesN j σ = edge[σ,σ ]σ,σ n j σ,σ and where σ,σ is defined by σ,σ := σ and σ are not on the same edge,1 if [σ, σ ] is an edge and σ → σ is direct, -1 if [σ, σ ] is an edge and σ → σ is indirect. Additionally, N fulfills f b σ = f (u h ) • N σ andwe can write N σ := -∂K φ σ • n d γ. Putting this in (55) we get 1 2 σ∈K σ>σ

We neglect the dependence on the mapping here again.

Connection between Flux Reconstruction and Residual DistributionInstead of using a variational or integral form like in DG, FR schemes are applied in their discretisation of the differential form (1) as it is described in subsection 2.2. The flux function is approximated by a polynomial of degree k + 1 denoted by f h . The discretisation of our underlying problem (1) readsdiv(f h + α∇ψ) = 0, (19)where α∇ψ is our correction function with the scaling term α = ff h • n. We change here the notation from subsection 2.2 on purpose to clarify that we are dealing now with the general case. We can translate it into triangles by setting ∇ψ := h. We get

Here, we neglect the fact that in our description of DG we did not approximate the flux function by a polynomial.

Here, the expression v, f := v T • f , where the flux is interpreted as a p × d matrix and the product is a d-row vector.

Here, we assume that the quadrature has sufficient order or accuracy such that also the discrete version of the theorem is fulfilled.

Acknowledgements

The second and third authors have been funded in by the the SNF project (Number 175784) " Solving advection dominated problems with high order schemes with polygonal meshes: application to compressible and incompressible flow problems".