The Green's function of the Lax-Wendroff and Beam-Warming schemes

We prove a sharp uniform generalized Gaussian bound for the Green's function of the Lax-Wendroff and Beam-Warming schemes. Our bound highlights the spatial region that leads to the well-known (rather weak) instability of these schemes in the maximum norm. We also recover uniform bounds in the maximum norm when these schemes are applied to initial data of bounded variation.

For 1 ≤ q < +∞, we let q (Z; C) denote the Banach space of complex valued sequences indexed by Z and such that the norm: , is finite.We also let ∞ (Z; C) denote the Banach space of bounded complex valued sequences indexed by Z and equipped with the norm: We let N * denote the set {1, 2, 3, . . .} of positive integers.The letter C, resp.c, denotes some large, resp.small, positive constant that may vary from one line to the other and possibly within the same line (for instance, we use the conventions C + C = C, C/c = C and so on).The dependence of the constants on the various involved parameters is made precise in the statement of the results and in the proofs.

Introduction and main result 1.Introduction
The Lax-Wendroff and Beam-Warming schemes are two second order finite difference approximations of the transport equation: For simplicity, the velocity a is a constant here which we choose to be positive (otherwise the Beam-Warming scheme defined in (3) below should be upwinded the other way round).The solution to (1) corresponds to a transport along the characteristics which implies that any L p norm of the initial condition u 0 is preserved: A desirable feature of any finite difference approximation of (1) is to satisfy a similar stability property in the discrete setting.Finite difference approximations of the solution to (1) amount to replacing u by a piecewise constant function with respect to both space and time.We thus introduce time and space steps ∆t > 0 and ∆x > 0 which we choose such that the ratio ∆t/∆x is a fixed positive constant.From now on, we let λ > 0 denote the fixed constant a ∆t/∆x.The solution u to (1) is then approximated by a sequence of functions defined as: where it remains to define inductively the sequence (u n j ) (n,j)∈N×Z .For the linear transport equation (1), the so-called Lax-Wendroff scheme reads: and the Beam-Warming scheme reads1 : In both (2) and (3), the sequence (f j ) j∈Z is meant to provide with an approximation of the initial condition u 0 for (1).A typical choice is: We recall that in (2) and (3), the fixed constant λ stands for a ∆t/∆x.Both schemes (2) and (3) are known to lead (at least formally) to second order approximations of the solution to (1) and to exhibit a dispersive behavior.This phenomenon is evidenced in Figure 1 where we show the propagation of the step function: at time t = 2 by the numerical schemes (2) and (3) (we choose a = 1 in the computations).Oscillating wavetrains are generated either behind or ahead of the discontinuities in the solution.Our goal in this article is to give an accurate description of the so-called Green's function (a.k.a. the fundamental solution) for the recurrence relations (2) and (3).Our result gives a quantitative description of the oscillating wavetrains displayed in Figure 1.Moreover it also gives a detailed analysis of the (rather weak) instability of the schemes (2) and (3) in ∞ (Z; C).We now briefly recall some bibliographic references on this instability phenomenon and state our main result.

A reminder on the instability of dispersive schemes in the maximum norm
Let us recall a few facts about Laurent (or convolution) operators.If a ∈ 1 (Z; C), we let L a denote the so-called Laurent operator associated with the sequence a [TE05, Nik02], which is defined by: whenever the defining formula (4) for the sequence L a u makes sense.For instance, the scheme (2) corresponds to the (finitely supported) sequence: and the scheme (3) corresponds to the (finitely supported) sequence: We always consider a sequence a ∈ 1 (Z; C) in what follows, and we shall even assume for our main result that the sequence a is finitely supported.Young's inequality shows that the operator L a acts boundedly on q (Z; C) for any q ∈ [1, +∞].The spectrum of L a is well-understood since the socalled Lévy-Wiener Theorem [New75] characterizes the invertible elements of 1 (Z; C) for the convolution product (and we have the morphism property L a • L b = L a b , where denotes the convolution product on 1 (Z; C)).Namely, the spectrum of L a as an operator acting on q (Z; C) does not depend on q and is nothing but the image of the Fourier transform of the sequence a, see [TE05]: Since a belongs to 1 (Z; C), its Fourier transform is continuous on R. It actually belongs to the so-called Wiener algebra, see [Kah70].
Let us now take a closer look at the norm of L a when acting on either 2 (Z; C) or ∞ (Z; C).For future use, we let F a denote the Fourier transform of the sequence a: (5) By Fourier analysis, see [GKO95, SW97, TE05], we have: Since the maps (a → L a ) and (a → F a ) are morphisms, this immediately gives: which gives the well-known von Neumann necessary and sufficient condition for 2 -stability (see In the context of finite difference schemes, F a is usually referred to as the amplification factor.Under the von Neumann condition, the operator L a is a contraction on 2 (Z; C) as is any of its powers.For the Lax-Wendroff (2) and Beam-Warming (3) schemes, we obtain (with rather self-explanatory notation): We then compute : and we thus find that the Lax-Wendroff scheme is2 -stable for λ ∈ (0, 1] and that the Beam-Warming scheme is 2 -stable for λ ∈ (0, 2].In what follows, we avoid the trivial cases λ = 1 in (2) and λ ∈ {1, 2} in (3) for which the schemes reduce to simple shift operators.We thus consider from now on : • the Lax-Wendroff scheme (2) with λ ∈ (0, 1), • or the Beam-Warming scheme (3) with λ ∈ (0, 1) ∪ (1, 2).
In particular, there holds the so-called dissipation condition: as well as the conservativity condition The stability of the Laurent operator L a on ∞ (Z; C) is less easy to study since the norm 2 : is usually larger than 1.The only favorable case is that when the coefficients a are nonnegative real numbers, which corresponds to monotone schemes or equivalently to the case where the a 's are the jump probabilities of a random walk on Z.There is no obvious reason why, even under a dissipation condition as in (7) for the Fourier transform F a , the operator L a should be power bounded when acting on ∞ (Z; C).This stability property in the maximum norm has been thoroughly investigated in the fundamental contribution [Tho65], see also [Des08,DSC14] for some recent developments.The main result in [Tho65] shows that stability in the maximum norm is related to the Taylor expansion of the Fourier transform F a at any point where the modulus of F a attains the value 1.For the Lax-Wendroff and Beam-Warming schemes, the modulus of the Fourier transform equals 1 only at the frequency 0, see (7), and we compute the Taylor expansions: as θ tends to 0. Since a purely imaginary term of the form i θ 3 arises in (8a) and (8b) (with a power 3 that is less than the first even power with a coefficient of negative real part), [Tho65,Theorem 3] shows that the associated convolution operators are not power bounded on ∞ (Z; C).In other words, the schemes (2) and (3) are not stable in ∞ (Z; C) uniformly in time.Namely, there holds (with, again, quite self-explanatory notation): Actually, [Tho65, Theorem 3] gives the lower bounds: for a suitable constant c > 0 (that only depends on λ).The sharp growth rate: as n tends to infinity, was proved in [Hed65] (see also [Hed66] for a general classification).We shall recover this growth rate later on by making it even more precise (see Section 4 for more details).

Main result
From now on, we consider a convolution operator L a where a ∈ 1 (Z; C) is a given sequence whose associated Fourier transform is denoted F a , see (5).We are interested in giving an accurate description of the powers L n a , as the integer n becomes large, in a case that exhibits a dispersive behavior as reported in Figure 1.We thus make the following assumption.
Assumption 1.The complex valued sequence a is finitely supported and sums to 1: Its Fourier transform F a (that is a trigonometric polynomial) satisfies the dissipation condition: Moreover, there exist a real number α, a nonzero real number c 3 , a positive real number c 4 and a holomorphic function ϕ defined on a neighborhood of 0 such that, as θ ∈ C tends to zero, there holds: Let us observe that since F a is a trigonometric polynomial, the definition (5) shows that it extends to a holomorphic function on the whole complex plane.Hence the holomorphy of the remainder ϕ in (10) is automatic.
In order to state our main result, we let from now on δ denote the discrete Dirac mass defined by: For any n ∈ N and j ∈ Z, we then use the notation G n j for the so-called Green's function associated with the operator L a , that is: Since L a is a convolution operator, the Green's function G n is nothing but the sequence a convolved with itself n − 1 times: G 0 = δ, G 1 = a, G 2 = a a and so on.Equivalently, we have: Our main result is the following.
Theorem 1. Assume that the constant c 3 in Assumption 1 is positive (which holds for the Lax-Wendroff scheme with λ ∈ (0, 1) and for the Beam-Warming scheme with λ ∈ (1, 2)).Then there exist two constants C > 0 and c > 0 such that the Green's function (G n j ) (n,j)∈N×Z satisfies the uniform bounds: and: where G n j is defined for n ∈ N * and j ∈ Z as: If c 3 is negative (which holds for the Beam-Warming scheme with λ ∈ (0, 1)), the bounds depending on the sign of j − α n should be switched and c 3 should be replaced by In what follows, only the values of G n j for j − α n < 0 will matter, but there is no problem defining G n j for any (n, j) ∈ N * × Z.An immediate consequence of Theorem 1 is the following corollary.
Corollary 1. Assume that the constant c 3 in Assumption 1 is positive.Then there exists a constant C > 0 such that for any n ∈ N * , the Green's function (G n j ) j∈Z satisfies: and where the approximate Green's function G n j is defined in (14).
In other words, the reason for the nonuniform integrability of the Green's function (G n j ) j∈Z is located in the region {j − α n < 0} (assuming c 3 > 0), and the resulting instability is entirely described by means of the approximate Green's function (G n j ) j∈Z .This observation will be used systematically in Section 4 below.
Let us clarify the position of our main result with some recent advances on local limit theorems for complex valued sequences.More precisely, in probability theory (the case where the a 's are nonnegative real numbers), local limit theorems describe the asymptotic behavior of G n j as n becomes large, see for instance [Pet75, Chapter VII].There have been recent developments of this theory in the case of complex valued sequences, of which a culminating point is [RSC15] (see [RSC17] for results on multidimensional problems).Going more into the details, when we apply [RSC15, Theorem 1.2] to sequences a satisfying Assumption 1, we obtain3 : uniformly with respect to j ∈ Z, and the function H i c 3 3 is defined as the oscillatory integral4 : Theorem 1 above does not give such an accurate (universal) description of G n j , but it gives global uniform bounds that cannot be obtained (at least, not in a straightforward way) from [RSC15].This difference between Theorem 1 and [RSC15] should be expected since the analysis in [RSC15] only retains the first two terms in the Taylor expansion (10) while we make here a very strong assumption on the third term (the O(θ 4 ) term in (10) that produces most of the dissipation mechanism).Last, we observe that the uniform bound in Theorem 1 for the case c 3 (j − α n) ≥ 0 is compatible with the well-known fast decaying behavior of the Airy function on R + .Theorem 1 and [RSC15] should therefore be seen as complementary.
We provide in Figure 2 with a representation of the Green's function (G n j ) j∈Z for the Lax-Wendroff scheme (2) for λ = 3/4 at various time iterations.Entirely similar pictures may be generated for the Beam-Warming scheme (3).The right picture in Figure 2 compares the exact Green's function (G n j ) j∈Z , that is depicted in blue, with yet another approximate Green's function from the one given in Theorem 1. Namely, we plot in red the approximate Green's function (H n j ) j∈Z defined by: where Ai denotes the Airy function.The values of c 3 and c 4 for the Lax-Wendroff scheme (2) are: The Gaussian like factor in the definition of H n j fits very well with the observed damping in the oscillations of G n j .With no Gaussian factor, the very slow decay of the Airy function on R − would not fit with the observed fast decaying behavior of G n j .However, we have not been able so far to obtain accurate bounds for the difference G n j − H n j .This is postponed to a future work.
We now give the proof of Theorem 1 in the case c 3 > 0 and leave the case c 3 < 0 to the interested reader.The rest of this article is organized as follows.In Section 2, we prove the sharp bound (12).We then give the proof of the bound (13) in Section 3. Corollary 1 as well as two other consequences of Theorem 1 are discussed in Section 4 in connection with [Hed65,Hed66] and [ELR90]. 2 Proof of the main result.I. Uniform bound and the fast decaying tail

Preliminary facts and notation
First of all, let us fix a constant δ 0 > 0 such that the function ϕ in Assumption 1 is holomorphic on the open square {z ∈ C / max ( |Re z|, |Im z| ) < 2 δ 0 }.We then define the constant C 0 > 0 as: In many arguments below, we shall use some contours that are located within the closed square {z ∈ C / max ( |Re z|, |Im z| ) ≤ δ 0 } in order to be able to bound from above the modulus of ϕ by the constant C 0 .The terms that involve ϕ will always be dealt with as remainders.The constants δ 0 and C 0 are fixed once and for all.We recall the (Fourier based) formula for the Green's function G n j .This is the same starting point as in [Tho65,Des08,DSC14,RSC15]. By the standard properties of the Fourier transform, we have: and, more generally, since G n is the convolution of a with itself n − 1 times, we have: We thus obtain the expression: where the final change of variables has been performed in order to stick as much as possible to the notation in [Hed65].From now, we use the notation: where the real number α corresponds to the first coefficient in the Taylor expansion (10).With this definition, we can rewrite (19) as: Since the sequence a is finitely supported, the sequence G n is also finitely supported for any n ∈ N.More precisely, there exists a positive integer M ∈ N * such that G n j = 0 for any j ∈ Z satisfying | j | > M n.The integer M measures the size of the support of a. From now on, we thus only consider the case | j | ≤ M n.Up to choosing M even larger, we can thus assume that the parameter ω in (20 The estimate of the Green's function G n j is split in several steps, depending on the size of ω.In some regimes, ω plays the role of a small parameter in (21), while n plays the role of a large parameter.In other regimes, ω is basically treated as a positive or negative constant and the integer n is the only (large) parameter in (21).Of course, our main attention below is to obtain estimates where all constants are independent of ω and n ∈ N * , which is similar to being independent of j ∈ Z and n ∈ N * .

The uniform bound
In this paragraph, we recall the argument in [RSC15] that yields a uniform O(n − 1/3 ) bound for the Green's function.The result follows from two lemmas on oscillatory integrals.
The second result is due to van der Corput.
Lemma 2 (van der Corput).There exists a numerical constant C > 0 such that for any real numbers a ≤ b, for any real valued function Following [RSC15], the combination of Lemma 1 and Lemma 2 yields a uniform bound for the Green's function, as we now recall.
Proposition 1.Under Assumption 1, there exists a constant C > 0 such that the Green's function in (11) satisfies: Proof.With the radius δ 0 > 0 fixed above and the constant C 0 in (18), we choose δ > 0 such that δ ≤ δ 0 and δ ≤ π , and where c 4 > 0 is the constant associated with the θ 4 term in the Taylor expansion (10).We then use the expression (21) and split the integral as: with: The estimate of ε n j follows from the dissipation assumption (9), from which we get: for some suitable constant C > 0 that is independent of n ∈ N * .It remains to focus on the oscillatory integral on the interval [−δ, δ].As in [RSC15], we use the Taylor expansion (10) (hence the restriction δ ≤ δ 0 ) and write: By applying the van der Corput Lemma (Lemma 2), there exists a constant C > 0 that is independent of ω and n such that: Furthermore, with our choice for the parameter δ, we have: and, differentiating the expression for g n (θ), we also get the bound: for some uniform constant C. We thus obtain that the quantity: is finite.Applying Lemma 1, we get the final estimate: with a constant C that does not depend on ω ∈ R nor on n ∈ N * .Going back to the decomposition (22) of G n j , the claim of Proposition 1 follows.

The fast decaying tail
In view of Proposition 1, we now consider n ∈ N * and j ∈ Z such that j − α n ≥ n 1/3 or, in other words, ω ≥ n − 2/3 .We are going to use a contour deformation argument in order to prove the generalized Gaussian bound stated in Theorem 1.The contour is guessed by following the so-called saddle point method, see [dB81].To determine the location of the saddle point as well as the path direction through the saddle point, we use the (truncated) phase function: For ω > 0, the two saddle points are ±i ω/(3 c 3 ) and the one with "lowest altitude" is i ω/(3 c 3 ), hence the choice made below in the proof of Proposition 2. The negative direction through this saddle point corresponds to the real axis (see Figure 3).We verify below that this choice of contour deformation, which is associated with the truncated phase, handles well the complete phase in (21) that includes the θ 4 term as well as the O(θ5 ) remainder.Our first result for the regime ω > 0 is the following.
Proposition 2. Under Assumption 1, there exists ω 0 > 0 and there exist two constants C > 0 and c > 0 such that the Green's function in (11) satisfies: , as long as n ∈ N * and the parameter ω defined in (20 Proof.In the regime considered in Proposition 2, the parameter ω in ( 21) is positive and small (but cannot go arbitrarily close to 0) and the integer n is thought as being large (at least large enough so that n − 2/3 ≤ ω 0 with ω 0 > 0 fixed as above).We start from the formula (21) and split again G n j as: with δ ∈ (0, π) to be fixed and: For reasons that will be made clear in the following lines, we choose the parameters δ and ω 0 such that the following inequalities hold 5 : With this choice of δ (that is fixed once and for all), we use Assumption 1 to write G n j as: with for suitable constants C > 0 and c > 0. This is the same first step as in the proof of Proposition 1.For later use, we define: and now focus on this term, which is the second term on the right hand side in the decomposition (24).By our choice of δ and ω 0 in (23), and the restriction n −2/3 ≤ ω ≤ ω 0 , we can use the contour deformation depicted in Figure 3.This contour remains within the closed square [−δ 0 , δ 0 ] × [−δ 0 , δ 0 ] on which ϕ is a holomorphic function and we can bound its modulus by C 0 .Applying Cauchy's formula [Rud87], we thus get: , where ε n j (1) corresponds to the integral on the left vertical segment, ε n j (2) corresponds to the integral on the right vertical segment, and the leading contribution H n j corresponds to the integral on the horizontal segment (these contributions are depicted in red in Figure 3).We obtain the expressions: and Let us start with ε n j (1).From now on, when we write exp( i • • • ), the dots always stand for a real number whose expression is useless since in the end a modulus of this expression will be simply estimated by 1. Going back to the defining expression for ε n j (1), we expand the various expressions within the integral and compute: We take the modulus of each side of the equality and apply the triangle inequality to get: The integration contour in the case n −2/3 ≤ ω ≤ ω 0 .The (approximate) saddle point is represented with a red bullet.The black bullets represent the end points of the three segments along which we compute the integrals ε n j (1), ε n j (2) and H n j .
where we used the fact that the integration contour is located within the region where the modulus of ϕ is less than C 0 .We now apply the Hölder inequality in C 2 to get (y is nonnegative here): which gives: imply that the terms with δ 5 and y 5 can be absorbed by half the ones with the power 4, namely: Let us now note that on the interval [0, ω/(3 c 3 )], we have: and we also use the last inequality in (23) to get: for suitable constants C > 0 and c > 0.
The estimate of the integral ε n j (2) along the right vertical segment is entirely similar.At this stage, we can collect (25) and ( 27) to show that for n −2/3 ≤ ω ≤ ω 0 , the Green's function G n j satisfies: where the expression of the (presumably leading) contribution H n j is given in (26).Let us therefore turn to the study of H n j .
We expand the expressions within the integral on the right hand side of (26) and obtain: We take the modulus on each side of the equality and apply the triangle inequality to get (the same Hölder inequality as above is used to deal with the remainder term on the second line): Again, our restrictions on ω 0 and δ in (23) imply that the final remainder terms can be absorbed by half of some already arising with a "good" sign, and we get the estimate: At last, we use the bound from above (see (23)): to get: and we therefore end up with our final estimate: We now combine (28) and (29) to get: , as long as j and n satisfy n −2/3 ≤ ω ≤ ω 0 .For such integers, it can be easily seen that the leading contribution on the right hand side of this last inequality is the second one.Namely, given C, c and ω 0 positive, we can always find other positive constants C and c such that for 0 < j − α n ≤ ω 0 n, there holds: .
We have thus obtained the estimate of the Green's function G n j as claimed in Proposition 2.
The final case to deal with in this section is when the parameter ω belongs to the interval [ω 0 , 2 M ].
Proposition 3.Under Assumption 1, with the same ω 0 > 0 as in Proposition 2, there exist two constants C > 0 and c > 0 such that the Green's function in (11) satisfies: , as long as n ∈ N * and the parameter ω defined in (20 Proof.The proof follows similar lines as that of Proposition 2. The difference is that we can no longer choose the same contour as in Figure 3 since for too large values of ω, there is no reason why the contour would remain within the holomorphy region of ϕ.Nevertheless, we still fix the parameters δ and ω 0 as in (23) and decompose G n j as in (24), where the remainder ε n j is uniformly exponentially small, see (25).We now use the contour depicted in Figure 4 where the "height" has been kept fixed equal to ω 0 /(3 c 3 ), independently of ω, in order to remain within the region where ϕ is holomorphic and bounded by C 0 .
The integration contour in the case ω 0 ≤ ω ≤ 2 M .The (approximate) saddle point is represented with a red bullet.The black bullets represent the end points of the three segments along which we compute the integrals Keeping the same notation as in the proof of Proposition 2 (see Figure 4), we have thus decomposed the Green's function G n j as: , where the first remainder ε n j is estimated as in (25), and the other terms are given by 6 : Let us start with the remainder term ε n j (1).As in the proof of Proposition 2, we expand the quantities within the integral and then apply the triangle inequality to get: e 16 n C 0 y 5 dy .
By using the restrictions (23) on δ and ω 0 , we can still absorb the δ 5 and y 5 terms to get: By using again (23) and ω ≥ ω 0 ≥ ω 0 /3, we end up with: The estimate of the other remainder term ε n j (2) is similar and we still get: where H n j is now given by (30).
6 For ε n j (1) and ε n j (2), the difference with the proof of Proposition 2 is in the upper bound of the interval over which we integrate.
We expand the various terms in (30) and obtain the expression: We take the modulus on each side of the inequality and apply the already used Hölder inequality to absorb the final remainder, which yields: We now use ω ≥ ω 0 as well as the restriction (23) on ω 0 to get: , where we recall that ω 0 > 0 has been fixed.Combining with (31), we have thus obtained the uniform exponential bound: for ω ≥ ω 0 , and we can convert this bound into: for ω ∈ [ω 0 , 2 M ] where ω 0 and M have already been fixed.Going back to the definition of ω, this gives the result of Proposition 3.

Conclusion
In this short paragraph, we explain why the above preliminary results imply the validity of (12).If n ≥ 1 and n 1/3 ≤ j − α n ≤ 2 M n, we use Propositions 2 and 3 to obtain the existence of positive constants C and c (independent of j and n) such that: This proves the validity of (12) for n 1/3 ≤ j − α n ≤ 2 M n.The constants C and c are now fixed.
For 2 M n < j − α n, the validity of ( 12) is even more clear since G n j is zero.It therefore remains to treat the case 0 ≤ j − α n ≤ n 1/3 .We use Proposition 1 to obtain: for some other constant C (possibly larger than the above constant C ), and (12) follows by using the inequalities: , for 0 ≤ j − α n ≤ n 1/3 .The bound (12) follows by choosing c := c and C := max(C , C e c ).We now turn to the case j − α n < 0 which is where the oscillations in the Green's function will arise.
3 Proof of the main result.II.The oscillations

The oscillations
We are now interested in the regime ω < 0 and start with the case where ω is small.This is the most difficult region where the Green's function exhibits oscillations.This is also the unique region that is the cause for the ∞ instability phenomenon which we have recalled in the introduction.As in Section 2, the regime −n − 2/3 ≤ ω ≤ 0 will de dealt with by Proposition 1, so we consider from now on ω ≤ −n − 2/3 , that is j − α n ≤ − n 1/3 .Proposition 4.Under Assumption 1, there exist ω 0 > 0 and there exist two constants C > 0 and c > 0 such that for any j ∈ Z and n ∈ N * , with G n j defined as in (14), there holds: as long as n ∈ N * and the parameter ω defined in (20) satisfy −ω 0 ≤ ω ≤ −n −2/3 (hence j − α n < 0).
Proof.Following the proof of Proposition 2, we introduce a parameter δ ∈ (0, π) to be fixed later on and split the quantity G n j as: with: For reasons that will be made clear in the following lines, we choose some parameters δ and ω such that the following inequalities hold: The parameter δ in the above decomposition of G n j is fixed once and for all, and we consider ω ∈ [− ω , − n − 2/3 ].We shall further need to restrict the possible values of ω later on but the restrictions (33) are a starting point for several terms that arise below.
With the choice (33) for δ, we use Assumption 1 to write G n j as (ω is negative here): with for suitable constants C > 0 and c > 0. For the integral on the right hand side of (34), we use Cauchy's formula and choose the contour depicted in Figure 5 which consists in: • A segment (with slope π/4) from the point −i Figure 5: The integration contour in the case −ω 0 ≤ ω ≤ −n −2/3 .The two red bullets represent the approximate saddle points and the black bullets represent the end points of the six segments along which we compute the contributions ε n j (1), ε n j (3), H j,n , H j,n , ε n j (4) and ε n j (2).
Thanks to the restrictions (33), this contour is included in the closed square [−δ 0 , δ 0 ] × [−δ 0 , δ 0 ] on which ϕ is a holomorphic function and we can bound its modulus by C 0 .According to our choice of contour, we decompose G n j in (34) as: where the six contributions ε n j (1), ε n j (3), H j,n , H j,n , ε n j (4), ε n j (2) correspond to the integral of the holomorphic function: along each of the six segments that make the above defined contour (we refer to Firgure 5 for an illustration).
Let us start with the contribution ε n j (1) that corresponds to the integral along the first vertical segment: Expanding the terms within the integral, we get: e − n (−δ+i y) 5 ϕ(δ−i y) dy .
We take the modulus on each side of the equality and apply the triangle inequality as well as the (already used) Hölder inequality for the remainder term to obtain (the arguments here are the same as in the proof of Proposition 2): Let us now note that on the interval [0, | ω |/(3 c 3 )], we have: and we also use the inequality (see (33) and use to get: for suitable constants C > 0 and c > 0 that do not depend on ω and n.The estimate of the integral ε n j (2) along the final vertical segment is entirely similar so we can collect ( 35) and (37) to obtain: We now turn to the contribution ε n j (3) which corresponds to the integral along the horizontal segment from −δ where the upper bound Ξ(ω) in the integral is defined as: We collect the real and purely imaginary contributions in the exponential functions of the expression for ε n j (3), which yields: The situation is less favorable than what occurred in the proof of Proposition 2 because we have a large factor in front of the integral (since ω is small, the | ω | 3/2 term is dominant with respect to the second term which scales as ω 2 ).Applying the triangle inequality as well as Hölder's inequality, we use (33) and get: Using (33) again, we have: and we thus obtain: where we used: We now use the following inequality which is valid for any couple (a, X) of positive numbers: and obtain: The estimate of the contribution ε n j (4) is entirely similar and we have thus obtained, for suitable constants C and c: Comparing with the estimate in (38), we see that the upper bound e − c n in (38) is always smaller than the right hand side in (39) (up to choosing appropriate constants C and c).We can thus add (38) and (39) and obtain, for new constants C and c: for −ω ≤ ω ≤ −n −2/3 and ω satisfying (33).
We now focus on the integrals H j,n and H j,n which correspond to the leading contributions arising from the segments that pass through the (approximate) saddle points.These contributions will give rise to the term G n j defined in (14).In what follows, we shall further need to restrict the possible values of ω since the condition |ω| ≤ ω will not be sufficient to absorb several terms that arise during the analysis.
Let us start with H j,n .We parametrize the corresponding segment by: and we therefore obtain the expression: where Θ(t) is given in (41) and the interval bounds t min , t max correspond to the end points of the interval in (41): We use the expression of Θ(t) in (41) to compute: where the complex valued functions p 0 , . . ., p 4 depend on ω only, are continuous on R and satisfy the following expansions as the real variable ω tends to zero: Re At this stage, we have written H j,n under the form: In what follows, we use the expansions (43) to simplify the expression of H j,n by getting rid of the terms in the exponential that contribute for remainders only.It is important to understand that here, a remainder is meant to be a quantity υ j,n that is uniformly (in n) summable with respect to j.Up to now, we have assumed that ω is negative and satisfies n −2/3 ≤ | ω | ≤ ω where ω > 0 has been fixed in order to meet the restrictions (33).The further restrictions on ω will help us absorb several terms.Namely, with ω fixed as in (33), there exists a constant C > 0 such that, for any ω ∈ [−ω , 0), there holds (see (43)): Using the inequality: as well as the estimate (see ( 41) and ( 42)): for a suitable constant C > 0 that does not depend on ω ∈ [−ω , 0), we obtain the estimate (for another suitable constant C independent of ω and n, at least for | ω | ≤ ω ): The crucial observation is now the following.In the integral on the right hand side of (45), either t is negative and t 3 is also negative, so we have: Or t is nonnegative and we have: In both cases, we can use a bound of the form: for a suitable constant c > 0. Going back to (45), we can thus choose ω 0 > 0 such that, for any ω ∈ [− ω 0 , 0), there holds: Using (43a), we find that the real part of p 0 (ω) can absorb the O(| ω | 5/2 ) remainder term, and we are eventually led to the estimate: It now remains to compute the integral on the right hand side and we obtain our first main simplification: We now simplify the term p 0 (ω) in the left hand side of (46) by using (43a) and (43b).Namely, by using the triangle inequality and (46), we get: for suitable constants C > 0 and c > 0 and | ω | ≤ ω 0 for a sufficiently small constant ω 0 > 0. The final integral is dealt with as above, meaning that we can absorb the t 3 term within the O(t 2 ) on the interval [t min , t max ] and we end up with the estimate: for suitable constants C > 0 and c > 0 that are independent of ω and n, and for | ω | ≤ ω 0 for a sufficiently small constant ω 0 > 0.
There are still two steps to further simplify the contribution H j,n .The first step consists in restricting to a symmetric interval for the integral.We recall once again the definition (42) of the interval bounds, and use again (47) with the triangle inequality to obtain: The length of the interval [t min , −t max ] is O(| ω |), see (42), and the function that is integrated is increasing with respect to t on the considered interval so we get: Once we have restricted to a symmetric interval in t, the final step consists in getting rid of the O(t 3 ) term in the integrated function.To achieve this, we use the inequality: and the fact that the function: is odd to get: For | t | ≤ t max , we use the bound: and we can thus use (49) to get: we end up with our final estimate: (50) where (50) holds for any ω ∈ [− ω 0 , −n − 2/3 ] and ω 0 > 0 is a sufficiently small constant.
We can now use exactly the same arguments to deal with the last contribution H j,n .Leaving the details to the interested reader, we end up with the estimate: (51) that is entirely similar to (50).Eventually, adding (50) and (51) and using the triangle inequality again, we end up with: where we recall the definition (14) of G n j .Adding (52) with (40), we see that the largest term on the right hand side is that of (40), which gives: and the proof of Proposition 4 is thus complete.

The tail
It remains to consider the regime where the parameter ω in (21) is negative and not small.Our result is summarized in the following Proposition.
Proof.As in the proof of Proposition 2, we decompose the Green's function into: with δ ∈ (0, π) to be fixed and: The parameter δ ∈ (0, π) is fixed here in order to satisfy7 : From Assumption 1, we have the uniform bound: so we focus on the second term in the above decomposition of G n j , which we denote H n j as in the proof of Proposition 2. Thanks to our choice for δ, we can use Cauchy's formula and use the segments [− δ, −i δ] ∪ [−i δ, δ] rather than the interval [− δ, δ] as an integration contour (see Figure 6).Those two segments are included in the closed square {z ∈ C / max ( |Re z|, |Im z| ) ≤ δ 0 } on which ϕ is holomorphic and its modulus is bounded by C 0 .We thus obtain a decomposition: , with (ω is negative here): We deal with the estimate of the first integral H n j (1) and leave the similar estimate of H n j (2) to the interested reader.Expanding the terms within the exponential function, we compute: We take the modulus on each side of the equality sign and use the uniform bound on ϕ to get: Using | ω | ≥ ω 0 and (53), we have: We then get: With the restriction (53), we can see that the function: achieves its minimum at t = 0, and we thus get: where we used (53) one last time.In the end, we have proved that each term in the decomposition of G n j is exponentially small in the considered regime −2 M ≤ ω ≤ −ω 0 .

Conclusion
We now discuss the proof of the bound (13) in Theorem 1.For | j | ≤ M n, we use the definition (14) of G n j to get: In particular, with ω 0 > 0 given as in Proposition 4, there exist two constants C > 0 and c > 0 such that the approximate Green's function in (14) satisfies: as long as j and n satisfy −2 M ≤ ω ≤ −ω 0 (recall that the parameter ω defined in (20)).Combining with Proposition 5, we obtain: Hence the bound (32) of Proposition 4 holds not only for −ω 0 ≤ ω ≤ −n − 2/3 but for the wider interval −2 M ≤ ω ≤ −n − 2/3 (up to choosing suitable constants C and c).This proves the validity of (13) for −2 M n ≤ j − α n ≤ −n 1/3 .Let us now discuss the case −n 1/3 ≤ j − α n ≤ 0. We go back to the definition (14) of G n j and use the following estimate from above: since we have |j − α n| ≤ n 1/3 .Adding with the uniform estimate of Proposition (1), we get: and this proves the validity of (13) for −n 1/3 ≤ j − α n ≤ 0. We now discuss the remaining values j < − 2 M n for which the Green's function G n j vanishes.For j < − 2 M n, we are in the case j − α n < − M n and G n j = 0, hence the bound (13) reduces to proving: for some suitable constants C and c.For j − α n < − M n, the definition (14) of G n j gives: Hence (54) is valid for j < − 2 M n since the decay of G n j is actually even faster than the right hand side of (54).This completes the proof of Theorem 1.

Consequences
This section is devoted to the proof of several consequences of Theorem 1.We first prove Corollary 1.We then state and prove two other consequences, the first of which makes the divergence of (G n j ) j∈Z in the 1 norm precise.Our last result deals with the action of the Laurent operator L a on sequences with bounded variations and gives an alternative proof to one of the main results in [ELR90].

On the instability in the maximum norm. I
We first prove Corollary 1.In many occurrences below, we shall compare Riemann sums with integrals of nonincreasing functions.Namely, for a ∈ R and f ∈ L 1 ((a, +∞); R + ) a nonincreasing function, for any nondecreasing sequence (x k ) k∈N with values in [a, +∞), there holds: Proof of Corollary 1.For n ∈ N * , we define the integer: and our first goal is to show a uniform upper bound for the sum: . The uniform bound for the term G n Jn−1 follows from Proposition 1 (actually, this term is not only bounded but tends to zero).We now use (12) to obtain (for suitable positive constants C and c): and we observe that the condition j ∈ J ε,n is equivalent to requiring: Applying the triangle inequality, we thus obtain: Choosing n ≥ N ε for some appropriate N ε ∈ N * , we thus obtain: and the claim (62) follows.
The conclusion of Lemma 3 on the equivalence between (58) and (63) follows from the triangle inequality.
In view of Lemma 3, we now wish to prove that (63) holds, where the sequence (H n j ) j≤Jn is defined in Lemma 3. A final simplification amounts to substituting a discrete sum by an integral, which will allow us to perform change of variables and other algebraic operations more easily.A precise statement is the following.
Lemma 4. With the sampling points y j,n given in (59) and the sequence (H n j ) j≤Jn defined in Lemma 3, there holds: Consequently, (63) holds if and only if: e − β 0 y 2 y 1/4 cos β 1 n 1/4 y 3/2 − π 4 dy = . (65) Let us assume for a moment that the conclusion of Lemma 4 holds and let us complete the proof of Corollary 2. For ease of reading, we define the sequence: ∀ n ∈ N * , I n := 3 β 1 2 π +∞ 0 e − β 0 y 2 y 1/4 cos β 1 n 1/4 y 3/2 − π 4 dy , of which we aim at computing the limit.We perform a first change of variable x = β 1 y 3/2 in the integral and obtain the relation: Using the definition (60) of β 0 and β 1 and using the functional relation x Γ(x) = Γ(x + 1), we end up proving (recall the definition (56)): which means that (65) holds.Going back all the way up, we have thus proved that (55) holds.At this stage, it just remains to prove Lemma 4.
Proof of Lemma 4. The multiplicative factor (3 β 1 )/(2 π) in the definition of H n j and in front of the integral in (64) is harmless so we omit it from now on.We let ε > 0 and consider two positive numbers A ε , B ε that satisfy: e − β 0 y 2 j,n ≤ ε 3 .
In particular, we get: We therefore focus on the difference: since the remaining parts on the left hand side of (64) are estimated (in absolute value) by 2 ε/3.We recall that throughout this proof we omit the multiplicative factor (3 β 1 )/(2 π) in the definition of H n j .
Then f n is Lipshitzean on the interval [A ε , B ε + 2] and, more precisely, there exists a constant C ε > 0 that depends on ε but not on n ∈ N * , and such that: We thus get the bound: f n (y) dy = as long as y j,n satisfies: Of course, the constant C ε in (66) is independent of n ∈ N * .Summing (66) over j, we get: Let us see what Corollary 2 gives in the case of the Lax-Wendroff scheme (2).Recalling (8a), we have: which means that for any λ ∈ (0, 1), the Green's function of the Lax-Wendroff scheme (2) satisfies: (1 − λ 2 ) 1/8 λ 1/4 .

Figure 1 :
Figure 1: Propagation of a step function by the Lax-Wendroff scheme (2) (left) and Beam-Warming scheme (3) (right) schemes with λ = 3/4.The initial condition is depicted in blue and the numerical solution at time t = 2 is depicted in red.

Figure 2 :
Figure 2: The Green's function of the Lax-Wendroff scheme (2) at various time iterations (left) and a zoom on the oscillations for n = 2400 (right).The chosen parameter is λ = 3/4.On the right: the Green's function G n j in blue and the approximate Green's function H n j defined in (17) in red dots.

e
By comparing the Riemann sums below with integrals, we therefore also get: j≤Jn / y j,n ≤Aε − β 0 y 2 j,n + j≤Jn / y j,n ≥Bε+n − 1/2 j≤Jn / Aε<y j,n <Bε+n − 1/2H n j − y J min −1,n y Jmax,n f n (y) dy ≤ C ε n − 1/4 ,where the interval of indices [J min , J max ] ∩ Z corresponds to the set of integers j such that (67) is satisfied.It is then easy to verify that the quantity:Bε Aε f n (y) dy − y J min −1,n y Jmax,nf n (y) dy tends to zero as n tends to infinity so overall, we can fix an integer N ε ∈ N * such that for any n ≥ N ε , there holds: j≤Jn / Aε<y j,n <Bε+n − 1/2 y) dy ≤ ε , for n ≥ N ε so Lemma 4 is proved.