Lie Symmetry Preservation by Finite Difference Schemes for the Burgers Equation

Invariant numerical schemes possess properties that may overcome the numerical properties of most of classical schemes. When they are constructed with moving frames, invariant schemes can present more stability and accuracy. The cornerstone is to select relevant moving frames. We present a new algorithmic process to do this. The construction of invariant schemes consists in parametrizing the scheme with constant coefficients. These coefficients are determined in order to satisfy a fixed order of accuracy and an equivariance condition. Numerical applications with the Burgers equation illustrate the high performances of the process.


Introduction
The partial differential equations used in fluid mechanics problems, and more generally in transfer phenomena, admit transformations conserving the whole set of solutions.The set of these transformations forms the symmetries of equations; it has the structure of a smooth local group, which is a Lie group.It is known as the symmetry group of the equation.This group plays an important role in the physical properties described by the equations.In fluid dynamics, some important flow properties such as the scaling laws, the Kolmogorov −5/3 law or the existence of special self-similar solutions, can be deduced from the symmetry group of the Navier-Stokes equations [1][2][3][4].In case of developed turbulence, Lindgren [5] shows that side laws computed from the symmetries of Navier-Stokes equations are more effective to predict high accuracy experimental measurements than the classical ones.Razafindralandy et al. has constructed a turbulence model for Large Eddy Simulation (LES) preserving the symmetry group of the Navier-Stokes equations [6][7][8].
Some models for anisothermal flows have also been studied using the Lie-symmetry group of the non-isothermal Navier-Stokes equations [9].
The success of this geometrical point of view naturally leads to the application of this approach to discrete space for numerical design.It is mainly crucial not to destroy the geometrical properties of geometric models when they are used in practical computation.More generally, one can expect numerical methods to reproduce geometrical properties of the continuous equation.Such numerical methods belong to the class of geometric integrators.Some of the most popular geometric integrators are symplectic integrators, adapted to differential hamiltonian systems [10].It is important to recall that they naturally reproduce the qualitative behaviour of the physical solutions (volume preservation in phase space, bounded energy conservation).Symplectic methods possess a numerical robustness that is significantly more important than those of classical methods over a long time integration [11,12].This performance is produced by conserving a geometric property called symplecticity of the differential system.Some processes allow to apply symplectic techniques on many partial differential equations (PDE), thanks to the fact that those PDE can be expressed in a formal way as Hamiltonian systems [13].There is also a generalization of symplectic schemes called multisymplectic schemes [14].But all of those extensions of numerical techniques of integration to PDE-some of them-are based on conservation laws [15,16], such as the energy conservation (for time independent systems) or the symplecticity conservation (for systems whose flows are symplectomorphisms).Those laws of conservation are tied to particular symmetries through Noether's theorem.But what if the PDE has other symmetries?How can we reproduce the gain of symplectic integrators for any PDEs admitting any symmetries?
Our aim is to present an approach which allows us to systematically construct numerical schemes preserving symmetries of the continuous equations.Such geometric integrators have been developed in many different ways.Some recent significant results have been obtained for the finite dimensional case [17,18].It has been shown that, e.g., in the case of first order differential equations, the numerical solutions of symmetry-preserving schemes coincide exactly with the solutions of the ODE [19].Particular attention has thus been paid to the symmetry analysis of differential-difference equations [20][21][22].For the PDE case, some pioneer analysis have been realized by Yanenko and Shokin.Their work on the subject yield to construct a semi-invariant numerical scheme by considering only finite order transformations [23,24].The implementation of this method has been realized by Hoarau et al. for the Burgers equation [25].Bakirova et al. propose to compute discrete difference invariants to get the expression of the numerical scheme [26].Budd et al. propose to modify the mesh grid in order to preserve symmetries [27,28].The symmetry-preserving method to which we propose a contribution is called the invariantization of numerical schemes.It is based on the concept of moving frame introduced by Élie Cartan, and recently reformulated by Fels and Olver [29][30][31].This method allows to construct invariant schemes, starting from classical scheme in an algorithmic way.The obtained set of invariant schemes could be infinite as the number of moving frames is.But all of the schemes are not competitive, nor at least as effective as the original classical scheme.The first applications of invariantization have been presented by Kim for ordinary differential equations (ODE) and some partial differential equations to convince that this process can lead to invariant schemes that outperform classical schemes [32,33].But the way to determine relevant moving frames is not well established yet [34].
This paper is devoted to present an efficient process to compute moving frames such that the invariant scheme has the guarantee to be at least as accurate as the original scheme.First of all we shall recall some basics and definitions and introduce invariant schemes.We will illustrate the performances of such schemes with the linear convection diffusion equation.We will then consider an algorithm consisting in the parametrization of the transformed scheme.The parametrization coefficients are computed in order to satisfy the equivariance condition and a fixed order of accuracy.The illustration of this process will be detailed with the non-linear Burgers equation.

Background
First and foremost, it is necessary to recall the main definitions and notations used in this paper.Consider a (p + q) dimensional manifold M whose elements are noted in local coordinates z = (x 1 , ..., x p , u 1 , ..., u q ), where x i , i = 1, ..., p are the independent variables and u j , j = 1, ..., q are the dependent variables.The joint product M ⋄n of M is defined to be the cartesian product M ×n of M without the n-tuples of identical points: A numerical scheme for a PDE E = 0 can be formally seen as a couple of functions (N, ϕ) defined over the joint product that verifies some condition of coalescent limit over a set of exact points [35]: a set of n points (z 1 ; ...; z n ) ∈ M ⋄n is exact for E = 0 if an exact solution of the equation whose graph contains z 1 , ..., z n exists.Then a numerical scheme for an equation E = 0 can formally be defined as a couple of numerical applications (N, ϕ) which verifies: where S is the set of all exact points for E = 0, and The application N represents the discrete equation-the difference scheme, and ϕ determines the equations of the stencil.Now let us consider a partial differential equation E = 0, its symmetry group G and a numerical scheme (N, ϕ) for the PDE.The action of G can be naturally prolonged to the joint product as a product action, for g ∈ G and z ∈ M ⋄n .
We say that (N, ϕ) is a G-symmetric numerical scheme if, ∀g ∈ G : A stronger property [36] for numerical scheme is the G-invariance, which can be defined as ∀g ∈ G : In Section (3.) a development of a method of the construction of invariant schemes based upon moving frames is presented: Theorem 1.A moving frame exists if and only if the action of the group G is regular and free.
The proof of this statement is in [31].
Moving frames allow to construct invariant counterpart of any numerical application as numerical schemes defined below.A straightforward computation using the property of equivariance of ρ leads to the following proposition: Proposition 2. Let (N, ϕ) be a numerical scheme for E = 0, G a symmetry group, and ρ : M ⋄n → G a moving frame.Then the scheme ( N , φ) defined as is an invariant numerical scheme for E = 0.
This fundamental result says that it is enough to construct equivariant moving frame to get invariance.Up to the present, the construction of moving frames is made according to the Cartan's method of normalization [37].This algorithm reduces the determination of a moving frame to the choice of a cross-section to the group orbit.But as an infinity of sub-manifold transverse to the group orbit could exist, one could have such many choices of moving frames.Therefore, one could obtain a very large choice of invariant schemes.Their performances are not equal at all [33], some invariant scheme could even give very spurious results.One simple explanation is that the process of invariantization does not take into account numerical properties.The moving frame constructed by geometrical arguments on the cross-section could give rise to an invariant scheme with low order of accuracy for example.
Our aim is to propose a method of construction of an invariant numerical scheme using moving frames that assures the order of accuracy, so that the invariant scheme transports the fundamental geometric structure of the continuous equation and possesses good numerical properties.

Scheme Parametrization Algorithm
The algorithm we propose is not based on the normalization process.It consists in parametrizing the classical scheme following the symmetry transformation.Then we suppose an algebraic expression of the symmetry parameters depending on constant real coefficients.On one hand, we compute the equivariance equations in function of the constant coefficients in order to get the invariance of the transformed scheme.On the other hand, we compute the consistency order of the transformed scheme in function of the real constant coefficients to insure the accuracy order of the scheme.It can be formalized as follow: 1. Start with a numerical scheme: where U represents the stencil of the discrete equation 2. From the prolonged action of the r-dimensional symmetry group G on M ⋄n : where ε 1 , ..., ε r are the symmetry parameters, write down the expression of the transformed scheme: where the discrete points z are in the transformed stencil 3. Then for each symmetry parameters, suppose an algebraic form that depends on constant real coefficients, that is for i = 1, ..., r: The transformed discrete equation has an expression that depends on those real constant coefficients: for z ∈ Ū .
4. For each symmetry parameter, compute the equivariance relation: 5. Compute conditions over the constant coefficients from the order of accuracy of the scheme: This process gives a computational fashion for determining moving frames, among a family, that are relevant for numerical aims.The main difference with the process of normalization as used by P. Kim is that we do not have to fix a cross-section to the group orbit to get the expression of moving frames.We directly deal with moving frames.This can be summarized by the algorithm described as: ⇓ ←→ Equivariance condition and Order of accuracy conditions

Numerical Illustration: The Burgers Equation
Consider the Burgers equation: The symmetry group of the Equation 19 is generated by 5 one-parameters subgroups given by the following infinitesimal generators [38]: Symmetries are given by the one-parameter transformations: • Spatial translation: • Time translation: • Projection: • Galilean boost: Most classical numerical schemes are invariant under the spatial and the time translations, as the scale transformation.But projection and Galilean boost are generally not respected.We will focus on the numerical properties induced by conserving them.

Construction of Invariant Numerical Scheme
We deal with the following transformation, depending on parameters ε 1 , ε 2 , ε 3 and ε 5 : Remark.This transformation corresponds to all Burger's symmetries except the scale transformation.
It would have been possible to take explicitly into account the scaling transformations, but for simplicity they are being omitted.Indeed, scale transformations would have added a multiplicative factor to the transformations ( 25), ( 26), (27).This factor vanishes in the computation of the invariant scheme.The scale transformations stay a symmetry during our process of invariantization.
Remark.As noted above, our interest is focused on the preservation of the projection and the Galilean transformations.The order they act on the discrete variables is not important because these two transformations commute.
Let us construct a numerical scheme invariant from the explicit forward time and centered space (FTCS) numerical scheme, in its non-conservative form: The transformed scheme is then: with x n+1 j − x n j = 0 because the grid is time-space orthogonal.We now have to specify the expression of ε 3 and ε 5 to make invariant the transformed scheme.
The fundamental result in Prop. 2 of the construction of invariant numerical scheme using moving frames guarantees that if ε 3 and ε 5 are moving frames for the actions ( 25), ( 26) and ( 27) then the transformed scheme ( 29) is an invariant scheme.

Construction of the moving frames
The condition of equivariance for ρ is: that is, for the left member of (30), and applied to points of the stencil used by the FTCS scheme: ) and for the right member: So that the equivariance condition (30) is equivalent to: In order to determine a family of moving frames allowing us to make the invariant numerical scheme, we suppose an algebraic expression of ε 3 .
As [ε 3 ] = s −1 , one chooses a combination of discrete variables which is dimensionally relevant.To keep the explicit form of the numerical scheme, there must not be any term in u n+1 .On the other hand, to keep the order of degree of the convective term, the parameter of symmetrization must be at most degree one.Finally, the construction of a moving frame from the application ε 3 requires that there is no term in ∆x and ∆t alone.Suppose then the algebraic form: Similar arguments for ε 5 (of dimension m.s −1 ) allow us to suppose a general form like: Conditions over a, b, c and d, e, f such that (32) and ( 33) become moving frames, and so verify (31), are: The determination of the expressions of the family of moving frames associated to a given group action is, of course, independent of the considered numerical scheme.But in order to get a couple of moving frames that can make an invariant numerical scheme conserving effective numerical properties, one has to impose some numerical arguments.

Order of Accuracy
A fundamental result assures that the invariantization of a numerical scheme using moving frames preserves the consistency [32].But the order of consistency may vary.
We consider the temporal term T t , the convective term T c and the diffusive term T d : The condition of consistency of the invariant scheme is: where Tt , Tc and Td are the invariant counterpart of the time, convection and diffusion terms: with coefficients a, b, c, d, e and f satisfying the equivariance condition (34).
Computing Taylor expansion gives: In order to keep the same degree of order, it is necessary to have: But we do not have any constrain over coefficients (d, e, f ) of ε 5 ; in fact ε 5 only appears in the convective term.Or the Taylor expansion of Tc indicates that: It means that the symmetrized convective term is identically null, therefore the convective phenomena are produced by the temporal term.
The moving frame associated to the projection is then determined in one way: The expression ε 5 is arbitrary.
Remark.It would have been possible to obtain this expression for ε 3 with the normalization process by taking the cross-section: Remark.The parameter of symmetrization associated to the Galilean boost does not occur when the numerical scheme is expressed in the regular and orthogonal original mesh.It is not the case anymore if the mesh grid is not orthogonal, nor if the numerical solution is expressed in the transformed frame of reference.

Numerical Applications
As mentioned before, most classical schemes do not respect the projection transformation and the Galilean boost.We expect from invariant schemes to have specially better behaviour to translate solutions related to those symmetries.

Self-Similar Solutions
Self-similar solutions have a particular importance in the fact that they allow to describe long-time behaviours or singularity.They can be computed by firstly reducing the PDE in some differential equation, and then solve the reduced equation.As standard numerical schemes (like FTCS) do not preserve every single symmetry, we expect they do not acquire the physical solutions as good as the invariant numerical schemes, unless having very thin grid.
Let us consider global invariants of the projection: The partial derivatives of u in terms of new coordinates are: Injecting those transformations into the burgers equation gives the reduced differential equation associated to the projection.It corresponds to the steady Burgers equation: The general solution of this ODE is: where D 0 and D 1 are arbitrary real constants.
Substituting the reduction transformation, one obtains general self similar solution of Burgers' equation for the projection on the upper half plan {t > 0}: )) In numerical example, boundary and initial conditions are given by a smooth solution (K 1 = 0 and D 0 = (2ν) −1/2 ) corresponding to a viscous shock: This viscous shock tends to a shock when the viscosity ν becomes closer to 0.
Let us consider the problem over the domain Ω =] − 1, 1[ whose boundary and initial conditions are given by the shock solution: Figure 1a shows the numerical behaviour of both classical and invariant scheme.The classical scheme presents some oscillations around the shock, although the invariant scheme reproduces the exact solution accurately.When the viscosity becomes smaller (but with constant CFL), the classical solution blows up, which is not the case for the invariant numerical scheme as shown in Figure 1b.The process of invariance has made the scheme more stable without any loss of accuracy.

Galilean Invariance
A fundamental principle of physics is the Galilean invariance.The systems are equivalent if they differ by a Galilean translation.¿From the numerical point of view, it corresponds to apply a Galilean boost of velocity λ (Figure 2): The frame moves at constant velocity while the solution is translated up to λ.
In order to test how a numerical scheme can reproduce this property, consider the problem over the domain Ω =] − 2; 10[, whose boundary conditions are given by the exact solution: This solution represents a damped sinusoidal weave.In Figure 3a, the numerical solution associated to the classical FTCS scheme shows an extra convection velocity growing up to the Galilean boost.The Galilean shift has produced error in evaluation of the velocity of the travelling weave.When comparing with the invariant scheme under Galilean transformation in Figure 3b, we observe that the behaviour is insensitive up to the Galilean shift as expected.The accuracy is compared in function of the value of the Galilean boost λ for the classical and the invariant scheme in Table 1, with the discrete errors: The error of the invariant scheme is the same at 10 −6 , whatever the value of λ is.When looking at the L 2 error the accuracy of the invariant scheme is unchanged until 10 −3 .The exact solution Equation 45 is enough smooth for low viscosity.The loss of Galilean invariance produces error in the convective velocity for classical scheme.But when the solution allows shock as the self-similar solution associated to the projection Equation 44, some numerical oscillations may appear even if the numerical conditions are effective enough in the original referential.Figure 4a represents the classical solution for different Galilean boost λ.Unphysical oscillations appear around the shock as λ grows.The numerical solution stays advected by an additional velocity, so that the velocity of the shock depends on the Galilean boost too.Figure 4b shows the solution associated to the invariant scheme.Even if some dissipation is produced around the shock, the solution stays the same whatever the value of λ is.For this example, the invariance under the Galilean transformation and the projection seems to be required to represent the physical solution numerically.

Conclusions
We have presented a new method of construction of invariant numerical schemes using moving frames.The process is based on geometric considerations and include as a parameter of construction the order of accuracy of the invariant scheme.The main motivation is to get rational arguments to choose moving frames leading to accurate schemes.This problem is one of the major limitation of the normalization process.It also points out that attributing geometrical features to schemes can break fundamental numerical properties.In this paper, we propose to fix the order of accuracy of invariant scheme to be the same as the original classical scheme.One could include other fundamental numerical conditions as stability conditions, dispersive relations, etc.Through computed examples over the Burgers equation, we have illustrated how the conservation of symmetry by numerical scheme has its importance for the capture of self-similar solutions and the respect of Galilean invariance.One of the main trump of such schemes is in succeeding in naturally capturing global features of equations described in the symmetry group.Comparable results could be obtained with classical schemes by dealing with a very refined grid at greater computational cost.

Figure 2 .Figure 3 .
Figure 2. Evolution of the frame.Straight and uniform displacement of the frame for a fixed boost velocity λ.

Table 1 .
Errors in function of Galilean boost.λ