Simon Baumstark

Erwan Faou

Katharina Schratz

UNIFORMLY ACCURATE EXPONENTIAL-TYPE INTEGRATORS FOR KLEIN-GORDON EQUATIONS WITH ASYMPTOTIC CONVERGENCE TO THE CLASSICAL NLS SPLITTING by

We introduce efficient and robust exponential-type integrators for Klein-Gordon equations which resolve the solution in the relativistic regime as well as in the highly-oscillatory non-relativistic regime without any step-size restriction under the same regularity assumptions on the initial data required for the integration of the corresponding nonlinear Schrödinger limit system. In contrast to previous works we do not employ any asymptotic/multiscale expansion of the solution. This allows us to derive uniform convergent schemes under far weaker regularity assumptions on the exact solution. In addition, the newly derived first-and second-order exponential-type integrators converge to the classical Lie, respectively, Strang splitting in the nonlinear Schrödinger limit. *

with ϕ 1 given in Definition 6. Note that the definition of the initial value u 0 * follows from (9). For complex-valued functions z (i.e., for u ≡ v) we similarly derive the exponential-type integration scheme (28)

, where the scheme in v n+1 * is obtained by replacing u n * ↔ v n * on the right-hand side of (28) with initial value

which after a Fourier pseudo-spectral space discretization only requires the usage of two Fast Fourier transforms (and its corresponding inverse counter parts) instead of three.

) , where the constant K r,M depends on r and M , but can be chosen independently of c.

Global error bound: Plugging the stability bound (32) as well as the local error bound (31) into (30) yields by a bootstrap argument that

) we have that

Together with the bound in (33) this completes the proof.

) is always satisfied for initial values

as then thanks to (18) we have

uniformly in c. For (u n * , v n * ) defined in (28) we set

Then, there exists a T r > 0 and τ 0 > 0 such that for all τ ≤ τ 0 and t n ≤ T r we have for all c > 0 that

, where the constants K 1,r,M 2γ , K 2,r,M and K * r,M,M 2γ ,tn can be chosen independently of c. Proof. -The proof follows the line of argumentation to the proof of Theorem 11 using "fractional estimates" of the operator A c .

Next we point out an interesting observation: For sufficiently smooth solutions the exponential-type integration scheme (28) converges in the limit c → ∞ to the classical Lie splitting of the corresponding nonlinear Schrödinger limit (2).

Remark 14 (Approximation in the non relativistic limit c → ∞)

The exponential-type integration scheme (28) corresponds for sufficiently smooth solutions in the limit (

), essentially to the Lie Splitting ([17, 8])

) for some constant k r that depends on M 2δ+ε and T , but is independent of c.

Proof. -In the following fix r > d/2, 0 < δ ≤ 2 and ε > 0:

1. Initial value approximation: Thanks to (37) we have by the definition of the initial value u 0 * in (28), respectively,

Regularity of the numerical solutions (u n

* , v n * ): Thanks to the regularity assumption (36) we have by Theorem 13 that there exists a T > 0 and τ 0 > 0 such that for all τ ≤ τ 0 we have *

and a similar bound on v n * -v n * ,∞ . Solving the recursion yields the assertion.

).

Remark 16. -Note that for l = ±2, -4 we have (see Definition 6)

Furthermore, (42) yields that

Thus, for sufficiently large values of c, more precisely if

and under the same regularity assumption (34) we may take instead of (28) the scheme

as a first-order numerical approximation to (u * (t n+1 ), v * (t n+1 )) in (12).

However, note that in the strongly non-relativistic limit regime (such that in particular cτ > 1) we may immediately take the Lie splitting scheme proposed in [9] as a suitable first-order approximation to (12) thanks to the following observation:

) for the nonlinear Schrödinger limit equation (2), namely,

). This assertion follows from [9] thanks to the approximation

Introduction

Cubic Klein-Gordon equations [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] c -2 ∂ tt z -∆z + c 2 z = |z| 2 z, z(0, x) = z 0 (x), ∂ t z(0, x) = c 2 z ′ 0 (x) are extensively studied numerically in the relativistic regime c = 1, see [START_REF] Gauckler | Error analysis of trigonometric integrators for seimilinear wave equations[END_REF][START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF] and the references therein. In contrast, the so-called "non-relativistic regime" c ≫ 1 is numerically much more involved due to the highly-oscillatory behavior of the solution. We refer to [START_REF] Engquist | Highly Oscillatory Problems[END_REF][START_REF] Hairer | Geometric Numerical Integration[END_REF] and the references therein for an introduction and overview on highly-oscillatory problems.

Analytically, the non-relativistic limit regime c → ∞ is well understood nowadays: The exact solution z of (1) allows (for sufficiently smooth initial data) the expansion z(t, x) = 1 2 e ic 2 t u * ,∞ (t, x) + e -ic 2 t v * ,∞ (t, x) + O(c -2 )

2000 Mathematics Subject Classification. -35C20 & 65M12 & 35L05.
on a time-interval uniform in c, where (u * ,∞ , v * ,∞ ) satisfy the cubic Schrödinger limit system

(2)

i∂ t u * ,∞ = 1 2 ∆u * ,∞ + 1 8 |u * ,∞ | 2 + 2 |v * ,∞ | 2 u * ,∞ u * ,∞ (0) = ϕ -iγ i∂ t v * ,∞ = 1 2 ∆v * ,∞ + 1 8 |v * ,∞ | 2 + 2 |u * ,∞ | 2 v * ,∞ , v * ,∞ (0) = ϕ -iγ with initial values z(0, x) c→∞ -→ γ(x) and c -1 c 2 -∆ -1/2 ∂ t z(0, x) c→∞ -→ ϕ(x),
see [START_REF] Masmoudi | From nonlinear Klein-Gordon equation to a system of coupled nonlinear Schrödinger equations[END_REF]Formula (1.3)] and for the periodic setting [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF]Formula (37)].

Also numerically, the non-relativistic limit regime c ≫ 1 has recently gained a lot of attention: Gautschi-type methods (see [START_REF] Hochbruck | A Gautschi-type method for oscillatory second-order differential equations[END_REF]) are analyzed in [START_REF] Bao | Analysis and comparison of numerical methods for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF]. However, due to the difficult structure of the problem they suffer from a severe time-step restriction as they introduce a global error of order c 4 τ 2 which requires the CFL-type condition c 2 τ < 1. To overcome this difficulty socalled limit integrators which reduce the highly-oscillatory problem to the corresponding nonoscillatory limit system (i.e., c → ∞ in (1)) as well as uniformly accurate schemes based on multiscale expansions were introduced in [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] and [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF][START_REF] Chartier | Uniformly accurate numerical schemes for highly oscillatory Klein-Gordon and nonlinear Schrödingier equations[END_REF]. In the following we give a comparison of these methods focusing on their convergence rates and regularity assumptions:

Limit integrators: Based on the modulated Fourier expansion of the exact solution (see [START_REF] Cohen | Modulated Fourier expansions of highly oscillatory differential equations[END_REF][START_REF] Hairer | Geometric Numerical Integration[END_REF]) numerical schemes for the Klein-Gordon equation in the strongly non-relativistic limit regime c ≫ 1 were introduced in [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF]. The benefit of this ansatz is that it allows us to reduce the highlyoscillatory problem [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] to the integration of the corresponding non-oscillatory limit Schrödinger equation [START_REF] Bao | Uniformly accurate multiscale time integrators for highly oscillatory second order differential equations[END_REF]. The latter can be carried out very efficiently without imposing any c-dependent step-size restriction. However, as this approach is based on the asymptotic expansion of the solution with respect to c -2 , it only allows error bounds of order O(c -2 + τ 2 ) when integrating the limit system with a second-order method. Henceforth, the limit integration method only yields an accurate approximation of the exact solution for sufficiently large values of c.

Uniformly accurate schemes based on multiscale expansions: Uniformly accurate schemes, i.e., schemes that work well for small as well as for large values of c were recently introduced for Klein-Gordon equations in [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF][START_REF] Chartier | Uniformly accurate numerical schemes for highly oscillatory Klein-Gordon and nonlinear Schrödingier equations[END_REF]. The idea is thereby based on a multiscale expansion of the exact solution. We also refer to [START_REF] Bao | Uniformly accurate multiscale time integrators for highly oscillatory second order differential equations[END_REF] for the construction and analysis in the case of highlyoscillatory second-order ordinary differential equations. The multiscale time integrator (MTI) pseudospectral method derived in [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] allows two independent error bounds at order O(τ 2 + c -2 ) and O(τ 2 c 2 ) for sufficiently smooth solutions. These error bounds immediately imply that the MTI method converges uniformly in time with linear convergence rate at O(τ ) for all c ≥ 1 thanks to the observation that min(c -2 , τ 2 c 2 ) ≤ τ . However, the optimal quadratic convergence rate at O(τ 2 ) is only achieved in the regimes when either 0 < c = O(1) (i.e., the relativistic regime) or 1 τ ≤ c (i.e., the strongly non-relativistic regime). In the context of ordinary differential equations similar error estimates were established for MTI methods in [START_REF] Bao | Uniformly accurate multiscale time integrators for highly oscillatory second order differential equations[END_REF]. The first-order uniform convergence of the MTI-FP method [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] holds for sufficiently smooth solutions: First-order convergence in time holds in H 2 uniformly in c for solutions in H 7 with sup 0≤t≤T z(t)

H 7 + c -2 ∂ t z(t) H 6 ≤ 1 (see [1, Theorem 4.1]
). First-order uniform convergence also holds in H 1 under weaker regularity assumptions, namely for solutions in H 6 satisfying sup 0≤t≤T z(t)

H 6 + c -2 ∂ t z(t) H 5 ≤ 1 if an additional CFL-type condition is imposed in space dimensions d = 2, 3 (see [1, Theorem 4.9]).
A second-order uniformly accurate scheme based on the Chapman-Enskog expansion was derived in [START_REF] Chartier | Uniformly accurate numerical schemes for highly oscillatory Klein-Gordon and nonlinear Schrödingier equations[END_REF] for the Klein-Gordon equation. Thereby, to control the remainders in the expansion, second-order uniform convergence in H r (r > d/2) requires sufficiently smooth solutions with in particular z(0) ∈ H r+10 . Also, due to the expansion, the problem needs to be considered in d + 1 dimensions.

We establish exponential-type integrators which converge with second-order accuracy in time uniformly in all c > 0. In comparison, the multiscale time integrators (MTI) derived in [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF][START_REF] Bao | Uniformly accurate multiscale time integrators for highly oscillatory second order differential equations[END_REF] only converge with first-order accuracy uniformly in all c ≥ 1. This is due to the fact that the MTI methods are based on the multiscale decomposition

z(t, x) = e itc 2 z n + (t, x) + e -itc 2 z n -(t, x) + r n (t, x)
which leads to a coupled second-order system in time in the c 2 -frequency waves z n ± and the rest frequency waves r n (cf. [1, System (2.4)]) and only allows numerical approximations at order O(τ 2 + c -2 ) and O(τ 2 c 2 ).

In contrast to [START_REF] Bao | A uniformly accurate multiscale time integrator pseudospectral method for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF][START_REF] Chartier | Uniformly accurate numerical schemes for highly oscillatory Klein-Gordon and nonlinear Schrödingier equations[END_REF][START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] we do not employ any asymptotic/multiscale expansion of the solution, but construct exponential-type integrators based on the following strategy:

1. In a first step we reformulate the Klein-Gordon equation (1) as a coupled first-order system in time via the transformations

u = z -i c -∆ + c 2 -1 ∂ t z, v = z -i c -∆ + c 2 -1 ∂ t z.
2. In a second step we rescale the coupled first-order system in time by looking at the so-called "twisted variables"

u * (t) = e ic 2 t u(t), v * (t) = e -ic 2 t v(t).
This essential step will later on allow us to treat the highly-oscillatory phases e ±ic 2 t and their interaction explicitly. 3. Finally, we iterate Duhamel's formula in (u * (t), v * (t)) and integrate the interactions of the highly-oscillatory phases exactly by approximating only the slowly varying parts.

This strategy in particular allows us to construct uniformly accurate exponential-type integrators up to order two which in addition asymptotically converge to the classical splitting approximation of the corresponding nonlinear Schrödinger limit system (2) given in [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF]. More precisely, the second-order exponential-type integrator converges for c → ∞ to the classical Strang splitting scheme

(3) u n+1 * ,∞ = e -i τ 2 ∆ 2 e -iτ 3 8 |e -i τ 2 ∆ 2 v n * ,∞ | 2 e -i τ 2 ∆ 2 u n * ,∞ , u 0
* ,∞ = ϕ -iγ associated to the nonlinear Schrödinger limit system (2) (see also Remark 31) where for simplicity we assumed that z is real-valued such that u * = v * . A similar result holds for the asymptotic convergence of the first-order exponential-type integration scheme towards the classical Lie splitting approximation (see also Remark 15).

In [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] the Strang splitting (3) is precisely proposed for the numerical approximation of non-relativistic Klein-Gordon solutions. However, in contrast to the uniformly accurate exponential-type integrators derived here, the scheme in [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] only yields second-order convergence in the strongly non-relativistic regime c > 1 τ due to its error bound at order O(τ 2 + c -2 ).

The main novelty in this work thus lies in the development and analysis of efficient and robust exponential-type integrators for the cubic Klein-Gordon equation (1) which

• allow second-order convergence uniformly in all c > 0 without adding an extra dimension to the problem. • resolve the solution z in the relativistic regime c = 1 as well as in the non-relativistic regime c → ∞ without any c-dependent step-size restriction under the same regularity assumptions as needed for the integration of the corresponding limit system. • in addition to converging uniformly in c, converge asymptotically to the classical Lie, respectively, Strang splitting for the corresponding nonlinear Schrödinger limit system (2) in the non-relativistic limit c → ∞.

Our strategy also applies to general polynomial nonlinearities f (z) = |z| 2p z with p ∈ N. However, for notational simplicity, we will focus only on the cubic case p = 1. Furthermore, for practical implementation issues we impose periodic boundary conditions, i.e., x ∈ T d .

We commence in Section 2 with rescaling the Klein-Gordon equation (1) which then allows us to construct first-, and second-order schemes that converge uniformly in c, see Section 3 and 4, respectively.

Scaling for uniformly accurate schemes

In a first step we reformulate the Klein-Gordon equation (1) as a first-order system in time which allows us to resolve the limit-behavior of the solution, i.e., its behavior for c → ∞ (see also [START_REF] Masmoudi | From nonlinear Klein-Gordon equation to a system of coupled nonlinear Schrödinger equations[END_REF][START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF]).

For a given c > 0, we define the operator

∇ c = -∆ + c 2 . ( 4
)
With this notation, equation (1) can be written as ( 5)

∂ tt z + c 2 ∇ 2 c z = c 2 f (z)
with the nonlinearity f (z) = |z| 2 z. In order to rewrite the above equation as a first-order system in time, we set

(6) u = z -ic -1 ∇ -1 c ∂ t z, v = z -ic -1 ∇ -1 c ∂ t z such that in particular (7) z = 1 2 (u + v). Remark 1. -If z is real, then u ≡ v.
A short calculation shows that in terms of the variables u and v equation ( 5) reads ( 8)

i∂ t u = -c ∇ c u + c ∇ -1 c f ( 1 2 (u + v)), i∂ t v = -c ∇ c v + c ∇ -1 c f ( 1 2 (u + v)) with the initial conditions (see (1)) (9) u(0) = z(0) -ic -1 ∇ -1 c z ′ (0), and v(0) = z(0) -ic -1 ∇ -1 c z ′ (0)
. Formally, the definition of ∇ c in (4) implies that c ∇ c = c 2 + "lower order terms in c". [START_REF] Gauckler | Error analysis of trigonometric integrators for seimilinear wave equations[END_REF] This observation motivates us to look at the so-called "twisted variables" by filtering out the highly-oscillatory parts explicitly: More precisely, we set u * (t) = e -ic 2 t u(t), v * (t) = e -ic 2 t v(t). [START_REF] Hairer | Geometric Numerical Integration[END_REF] This idea of "twisting" the variable is well known in numerical analysis, for instance in the context of the modulated Fourier expansion [START_REF] Cohen | Modulated Fourier expansions of highly oscillatory differential equations[END_REF][START_REF] Hairer | Geometric Numerical Integration[END_REF], adiabatic integrators [START_REF] Lorenz | Adiabatic integrators for highly-oscillatory second-order linear differential equations with time-varying eigendecomposition[END_REF][START_REF] Hairer | Geometric Numerical Integration[END_REF] as well as Lawson-type Runge-Kutta methods [START_REF] Lawson | Generalized Runge-Kutta processes for stable systems with large Lipschitz constants[END_REF]. In the case of "multiple high frequencies" it is also widely used in the analysis of partial differential equations in low regularity spaces (see for instance [START_REF] Bourgain | Fourier transform restriction phenomena for certain lattice subsets and applications to nonlinear evolution equations[END_REF]) and has been recently successfully employed numerically for the construction of low-regularity exponential-type integrators for the KdV and Schrödinger equation, see [START_REF] Hofmanova | An exponential-type integrator for the KdV equation[END_REF][START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF].

In terms of (u * , v * ) system (8) reads (cf. [18, Formula (2.1)]) [START_REF] Hochbruck | A Gautschi-type method for oscillatory second-order differential equations[END_REF] i∂

t u * = -A c u * + c ∇ -1 c e -ic 2 t f 1 2 (e ic 2 t u * + e -ic 2 t v * ) i∂ t v * = -A c v * + c ∇ -1 c e -ic 2 t f 1 2 (e ic 2 t v * + e -ic 2 t u * ) with the leading operator (13) A c := c ∇ c -c 2 .
Remark 2. -The advantage of looking numerically at (u * , v * ) instead of (u, v) lies in the fact that the leading operator -c ∇ c in system (8) is of order c 2 (see [START_REF] Gauckler | Error analysis of trigonometric integrators for seimilinear wave equations[END_REF]) whereas its counterpart -A c in system [START_REF] Hochbruck | A Gautschi-type method for oscillatory second-order differential equations[END_REF] is "of order one in c" (see Lemma 3 below).

In the following we construct integration schemes for [START_REF] Hochbruck | A Gautschi-type method for oscillatory second-order differential equations[END_REF] based on Duhamel's formula ( 14)

u * (t n + τ ) = e iτ Ac u * (t n ) -ic ∇ -1 c τ 0 e i(τ -s)Ac e -ic 2 (tn+s) f 1 2 (e ic 2 (tn+s) u * (t n + s) + e -ic 2 (tn+s) v * (t n + s)) ds, v * (t n + τ ) = e iτ Ac v * (t n ) -ic ∇ -1 c τ 0 e i(τ -s)Ac e -ic 2 (tn+s) f 1 2 (e ic 2 (tn+s) v * (t n + s) + e -ic 2 (tn+s) u * (t n + s)) ds.
Thereby, to guarantee uniform convergence with respect to c we make the following important observations. We define the Sobolev norm on T d by the formula

u 2 r = k∈Z d (1 + |k| 2 ) r |û k | 2 , where ûk = 1 (2π) d T d u(x)e ik•x dx, where for k = (k 1 , . . . , k d ) ∈ Z d , we set k • x = k 1 x 1 + • • • k d x d and |k| 2 = k • k.
Moreover, for a given linear bounded operator L we denote by L r its corresponding induced norm.

Lemma 3 (Uniform bound on the operator

A c ). -For all c ∈ R we have that A c u r ≤ 1 2 u r+2 . (15)
Proof. -The operator A c acts a the Fourier multiplier

(A c ) k = c 2 -c c 2 + |k| 2 , k ∈ Z d .

Thus, the assertion follows thanks to the bound

A c u 2 r = k∈Z d (1 + |k| 2 ) r c c 2 + |k| 2 -c 2 2 |û k | 2 ≤ k∈Z d (1 + |k| 2 ) r |k| 2 2 2 |û k | 2 ,
where we have used that

√ 1 + x 2 ≤ 1 + 1 2 x 2 for all x ∈ R. Lemma 4.
-For all t ∈ R we have that

(16) e itAc r = 1 and e -itAc -1 u r ≤ 1 2 |t| u r+2 .
Proof. -The first assertion is obvious and the second follows thanks to the estimate |(e ix -1)| ≤ |x| which holds for all x ∈ R together with the essential bound on the operator A c given in [START_REF] Lawson | Generalized Runge-Kutta processes for stable systems with large Lipschitz constants[END_REF].

In particular, the time derivatives (u ′ * (t), v ′ * (t)) can be bounded uniformly in c.

Lemma 5 (Uniform bounds on the derivatives

(u ′ * (t), v ′ * (t))) Fix r > d/2. Solutions of (12) satisfy (17) u * (t n + s) -u * (t n ) r ≤ 1 2 |s| u * (t n ) r+2 + 1 8 |s| sup 0≤ξ≤s u * (t n + ξ) r + v * (t n + ξ) r 3 , v * (t n + s) -v * (t n ) r ≤ 1 2 |s| v * (t n ) r+2 + 1 8 |s| sup 0≤ξ≤s u * (t n + ξ) r + v * (t n + ξ) r 3 .
Proof.

-The assertion follows thanks to Lemma 4 together with the bound

(18) c ∇ -1 c r ≤ 1
which implies by Duhamel's perturbation formula ( 14) that

u * (t n + s) -u * (t n ) r ≤ |s| A c u * (t n ) r + 1 8 |s| c ∇ -1 c r sup 0≤ξ≤s u * (t n + ξ) r + v * (t n + ξ) r 3 ≤ 1 2 |s| u * (t n ) r+2 + 1 8 |s| sup 0≤ξ≤s u * (t n + ξ) r + v * (t n + ξ) r 3 .
Similarly we can establish the bound on the derivative v ′ * (t).

We will also employ the so-called "ϕ j functions" given in the following Definition.

Definition 6 (ϕ j functions [START_REF] Hochbruck | Exponential integrators[END_REF]). -Set

ϕ 0 (z) := e z and ϕ k (z) := 1 0 e (1-θ)z θ k-1 (k -1)! dθ, k ≥ 1 such that in particular ϕ 0 (z) = e z , ϕ 1 (z) = e z -1 z , ϕ 2 (z) = ϕ 0 (z) -ϕ 1 (z) z .
In the following we assume local-wellposedness (LWP) of ( 12) in H r .

Assumption 7. -Fix r > d/2 and assume that there exists a T r = T > 0 such that the solutions (u * (t), v * (t)) of ( 12) satisfy

sup 0≤t≤T u * (t) r + v * (t) r ≤ M uniformly in c. Remark 8.
-The previous assumption holds under the following condition on the initial data

z(0) r + c -1 ∇ -1 c z ′ (0) r ≤ M 0
where M 0 does not depend on c as can be easily proved from the formulation (14).

A first-order uniformly accurate scheme

In this section we derive a first-order exponential-type integration scheme for the solutions (u * , v * ) of ( 12) which allows first-order uniform time-convergence with respect to c. The construction is thereby based on Duhamel's formula [START_REF] Hofmanova | An exponential-type integrator for the KdV equation[END_REF] and the essential estimates in Lemma 3, 4 and 5. For the derivation we will for simplicity assume that z is real, which (by Remark 1) implies that u = v such that system (12) reduces to [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF] i∂

t u * = -A c u * + 1 8 c ∇ -1 c e -ic 2 t e ic 2 t u * + e -ic 2 t u * 3 with mild-solutions (20) u * (t n + τ ) = e iτ Ac u * (t n ) - i 8 c ∇ -1 c τ 0 e i(τ -s)Ac e -ic 2 (tn+s) e ic 2 (tn+s) u * (t n + s) + e -ic 2 (tn+s) u * (t n + s) 3 ds.
3.1. Construction. -In order to derive a first-order scheme, we need to impose additional regularity assumptions on the exact solution u * (t) of ( 19).

u * (t n + τ ) = e iτ Ac u * (t n ) - i 8 c ∇ -1 c e iτ Ac τ 0 e -ic 2 (tn+s) e ic 2 (tn+s) u * (t n ) + e -ic 2 (tn+s) u * (t n ) 3 ds + R(τ, t n , u * ),
where the remainder R(τ, t n , u * ) satisfies thanks to the bounds ( 16), ( 17) and ( 18) that

(22) R(τ, t n , u * ) r ≤ τ 2 k r (M 2 ),
for some constant k r (M 2 ) which depends on M 2 (see Assumption 9) and r, but is independent of c. Solving the integral in (21) (in particular, integrating the highly-oscillatory phases exp(±ilc 2 s) exactly) furthermore yields by adding and subtracting the term

τ 3i 8 e iτ Ac |u * (t n )| 2 u * (t n ) (see Remark 17 below for the purpose of this manipulation) that (23) u * (t n + τ ) = e iτ Ac 1 -τ 3i 8 |u * (t n )| 2 u * (t n ) -τ 3i 8 c ∇ -1 c -1 e iτ Ac |u * (t n )| 2 u * (t n ) -τ i 8 c ∇ -1 c e iτ Ac e 2ic 2 tn ϕ 1 (2ic 2 τ )u 3 * (t n ) + e -2ic 2 tn ϕ 1 (-2ic 2 τ )3|u * (t n )| 2 u * (t n ) + e -4ic 2 tn ϕ 1 (-4ic 2 τ )u * 3 (t n ) + R(τ, t n , u * )
with ϕ 1 given in Definition 6.

As the operator e itAc is a linear isometry in H r and by Taylor series expansion

|1 -x -e -x | = O(x 2 ) we obtain for r > d/2 that (24) e iτ Ac 1 -τ 3i 8 |u * (t n )| 2 u * (t n ) -e iτ Ac e -τ 3i 8 |u * (tn)| 2 u * (t n ) r ≤ k r 3τ 2 u * (t n ) 3 r
for some constant k r independent of c.

The bound in (24) allows us to express (23) as follows (25)

u * (t n + τ ) = e iτ Ac e -τ 3i 8 |u * (tn)| 2 u * (t n ) -τ 3i 8 c ∇ -1 c -1 e iτ Ac |u * (t n )| 2 u * (t n ) -τ i 8 c ∇ -1 c e iτ Ac e 2ic 2 tn ϕ 1 (2ic 2 τ )u 3 * (t n ) + e -2ic 2 tn ϕ 1 (-2ic 2 τ )3|u * (t n )| 2 u * (t n ) + e -4ic 2 tn ϕ 1 (-4ic 2 τ )u * 3 (t n ) + R(τ, t n , u * ),
where the remainder R(τ, t n , u * ) satisfies thanks to ( 22) and ( 24) that

(26) R(τ, t n , u * ) r ≤ τ 2 k r (M 2 ),
for some constant k r (M 2 ) which depends on M 2 (see Assumption 9) and r, but is independent of c.

The expansion (25) of the exact solution u * (t) builds the basis of our numerical scheme: As a numerical approximation to the exact solution u * (t) at time t n+1 = t n + τ we choose the exponential-type integration scheme (27)

u n+1 * = e iτ Ac e -τ 3i 8 |u n * | 2 u n * -τ 3i 8 c ∇ -1 c -1 e iτ Ac |u n * | 2 u n
Remark 10 (Practical implementation). -To reduce the computational effort we may express the first-order scheme (28) in its equivalent form

u n+1 * = e iτ Ac e -τ i 8 |u n * | 2 +2|v n
In Section 3.2 below we prove that the exponential-type integration scheme (28) is first-order convergent uniformly in c for sufficiently smooth solutions. Furthermore, we give a fractional convergence result under weaker regularity assumptions and analyze its behavior in the nonrelativistic limit regime c → ∞. In Section 3.3 we give some simplifications in the latter regime.

Convergence analysis. -

The exponential-type integration scheme (28) converges (by construction) with first-order in time uniformly with respect to c, see Theorem 11. Furthermore, a fractional convergence bound holds true for less regular solutions, see Theorem 13. In particular, in the limit c → ∞ the scheme converges to the classical Lie splitting applied to the nonlinear Schrödinger limit system, see Lemma 15.

Theorem 11 (Convergence bound for the first-order scheme)

Fix r > d/2 and assume that

z(0) r+2 + c -1 ∇ -1 c z ′ (0) r+2 ≤ M 2 (29)
uniformly in c. For (u n * , v n * ) defined in (28) we set

z n := 1 2 e ic 2 tn u n * + e -ic 2 tn v n * .
Then, there exists a T r > 0 and τ 0 > 0 such that for all τ ≤ τ 0 and t n ≤ T r we have for all c > 0 that

z(t n ) -z n r ≤ τ K 1,r,M 2 e tnK 2,r,M ≤ τ K * r,M,M 2 ,tn
, where the constants K 1,r,M 2 , K 2,r,M and K * r,M,M 2 ,tn can be chosen independently of c.

Proof. -Fix r > d/2. First note that the regularity assumption on the initial data in (29) implies the regularity Assumption 9 on (u * , v * ), i.e., there exists a T r > 0 such that

sup 0≤t≤Tr u * (t) r+2 + v * (t) r+2 ≤ k(M 2 )
for some constant k that depends on M 2 and T r , but can be chosen independently of c.

In the following let (φ t u * , φ t v * ) denote the exact flow of ( 12) and let (Φ τ u * , Φ τ v * ) denote the numerical flow defined in (28), i.e.,

u * (t n+1 ) = φ τ u * (u * (t n ), v * (t n )), u n+1 * = Φ τ u * (u n * , v n *
) and a similar formula for the functions v * (t n ) and v n * . This allows us to split the global error as follows (30)

u * (t n+1 ) -u n+1 * = φ τ u * (u * (t n ), v * (t n )) -Φ τ u * (u n * , v n * ) = Φ τ u * (u * (t n ), v * (t n )) -Φ τ u * (u n * , v n * ) + φ τ u * (u * (t n ), v * (t n )) -Φ τ u * (u * (t n ), v * (t n )).
Local error bound: With the aid of (26) we have by the expansion of the exact solution in (25) and the definition of the numerical scheme (28) that (31)

φ τ u * (u * (t n ), v * (t n )) -Φ τ u * (u * (t n ), v * (t n )) r = R(τ, t n , u * , v * ) r ≤ τ 2 k r (M 2
) for some constant k r which depends on M 2 and r, but can be chosen independently of c.

Stability bound: Note that for all l ∈ Z we have that

ϕ 1 (iτ c 2 l) r ≤ 2.
Thus, as e itAc is a linear isometry for all t ∈ R we obtain together with the bound (18) that as long as

u n * r ≤ 2M and u(t n ) r ≤ M we have that (32) Φ τ u * (u * (t n ), v * (t n ))-Φ τ u * (u n * , v n
Under weaker regularity assumptions on the exact solution we obtain uniform fractional convergence of the formally first-order scheme (28).

Theorem 13 (Fractional convergence bound for the first-order scheme)

Fix r > d/2 and assume that for some 0 < γ ≤ 1

z(0) r+2γ + c -1 ∇ -1 c z ′ (0) r+2γ ≤ M 2γ (34) z(0) c→∞ -→ γ and c -1 ∇ -1 c z ′ (0) c→∞ -→ ϕ.
More precisely, the following Lemma holds.

Lemma 15. -Fix r > d/2 and let 0 < δ ≤ 2. Assume that (36) z(0) r+2δ+ε + c -1 ∇ -1 c z ′ (0) r+2δ+ε ≤ M 2δ+ε
for some ε > 0 uniformly in c and let the initial value approximation (there exist functions ϕ, γ such that)

z(0) -γ r + c -1 ∇ -1 c z ′ (0) -ϕ r ≤ k r c -δ (37)
hold for some constant k r independent of c.

Then, there exists a T > 0 and τ 0 > 0 such that for all τ ≤ τ 0 the difference of the first-order scheme (28) for system [START_REF] Hochbruck | A Gautschi-type method for oscillatory second-order differential equations[END_REF] and the Lie splitting (35) for the limit Schrödinger equation (2) satisfies for t n ≤ T and all c > 0 with

(38) τ c 2-δ ≥ 1 that u n * -u n (39) u n * r+2δ + v n * r+2δ
≤ m 2δ as long as t n ≤ T for some constant m 2δ depending on M 2δ+ε and T , but not on c.

3. Regularity of the numerical solutions (u n * ,∞ , v n * ,∞ ): Thanks to the regularity assumption (36) we have by (37) and the global first-order convergence result of the Lie splitting for semilinear Schrödinger equations (see for instance [START_REF] Faou | Geometric numerical integration and Schrödinger equations[END_REF][START_REF] Lubich | On splitting methods for Schrödinger-Poisson and cubic nonlinear Schrödinger equations[END_REF]) that there exists a T > 0 and τ 0 > 0 such that for all τ ≤ τ 0 we have (40)

u n * ,∞ r + v n * ,∞ r ≤ m 0 as long as t n ≤ T for some constant m 0 depending on M r and T , but not on c.

4. Approximations: Using the following bounds, γ > 1

(41) 1 + x 2 -1 - 1 2 x 2 ≤ x 2γ and 1 √ 1 + x 2 -1 ≤ x 2γ-2 ,
together with the Definition of ϕ 1 (see Definition 6) we have for every f ∈ H r+2+2δ ,

A c + ∆ 2 f r + c ∇ -1 c -1 f r+2 + ϕ 1 (ilc 2 τ )f r+2+δ ≤ k r c -δ f r+2+2δ (42)
for l = ±2, -4 and for some constant k r independent of c, where we used (38) for the last estimate.

5. Difference of the numerical solutions: Thanks to the a priori regularity of the numerical solutions (39) and (40) we obtain with the aid of (42) under assumption (38) for the difference

u n * -u n * ,∞ that (43) u n+1
4. A second-order uniformly accurate scheme

In this section we derive a second-order exponential-type integration scheme for the solutions (u * , v * ) of ( 12) which allows second-order uniform time-convergence with respect to c. For notational simplicity we again assume that z is real, which reduces the coupled system [START_REF] Hochbruck | A Gautschi-type method for oscillatory second-order differential equations[END_REF] to equation [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF] with mild-solutions [START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF] (see also Remark 1).

The construction of the second-order scheme is again based on Duhamel's formula [START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF] and the essential estimates in Lemma 3, 4 and 5. However, the second-order approximation is much more involved due to the fact that

u ′ * (t) = O(1), but u ′′ * (t) = O(c 2 ).
The latter observation prevents us from simply applying the higher-order Taylor series expansion

u * (t n + s) = u * (t n ) + su ′ * (t n ) + O s 2 u ′′ * (t n + ξ)
in Duhamel's formula [START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF] as this would lead to the "classical" c-dependent error at order O(τ 2 c 2 ). Therefore we need to carry out a much more careful frequency analysis by iterating Duhamel's formula [START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF] twice and controlling the appearing highly-oscillatory terms e ±ic 2 t and their interactions e ilc 2 t (l ∈ Z) precisely.

4.1. Construction of a second-order uniformly accurate scheme. -In Section 4.1.1 we state the necessary regularity assumptions on the solution u * and derive two useful expansions. In Section 4.1.2 we collect some useful lemmata on highly-oscillatory integrals and their approximations. These approximations will then allow us to construct a uniformly accurate second-order scheme in Section 4.1.3. The rigorous convergence analysis is given in Section 4.2.

Regularity and expansion of the exact solution.

-In order to derive a second-order scheme, we need to impose additional regularity on the exact solution u * (t) of [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF]. In Lemma 20 below we derive two useful expansions of the exact solution u * of [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF]. For this purpose we introduce the following definition.

Definition 19. -For some function v and t n , t ∈ R we set (45) Ψ c 2 (t n , t, v) := te 2ic 2 tn ϕ 1 2ic 2 t v 3 + 3te -2ic 2 tn ϕ 1 -2ic 2 t |v| 2 v + te -4ic 2 tn ϕ 1 -4ic 2 t v 3 .
The above defintion allows us the following expansions of the exact solution u * .

Lemma 20. -Fix r > d/2. Then the exact solution of (19) satisfies the expansions

u * (t n + s) = e isAc u * (t n ) - 3i 8 c ∇ -1 c s 0 e i(s-ξ)Ac e iξAc u * (t n ) 2 e iξAc u * (t n ) dξ - i 8 c ∇ -1 c Ψ c 2 (t n , s, u * (t n )) + R 1 (t n , s, u * )
and

u * (t n + s) = e isAc u * (t n ) - i 8 c ∇ -1 c 3s |u * (t n )| 2 u * (t n ) + Ψ c 2 (t n , s, u * (t n )) + R 2 (t n , s, u * )
with Ψ c 2 defined in (45) and where the remainders satisfy

(46) R 1 (t n , s, u * ) r + R 2 (t n , s, u * ) r ≤ s 2 k r (M 2 )
for some constant k r (M 2 ) which depends on M 2 , but is independent of c.

Proof. -Note that by Duhamel's perturbation formula [START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF] we have that (47)

u * (t n + s) = e isAc u * (t n ) - i 8 c ∇ -1 c s 0 e i(s-ξ)Ac 3 |u * (t n + ξ)| 2 u * (t n + ξ) + e 2ic 2 (tn+ξ) u * (t n + ξ) 3 +3e -2ic 2 (tn+ξ) |u * (t n + ξ)| 2 u * (t n + ξ) + e -4ic 2 (tn+ξ) u * (t n + ξ) 3 dξ.
Therefore, the bound on c ∇ -1 c given in [START_REF] Masmoudi | From nonlinear Klein-Gordon equation to a system of coupled nonlinear Schrödinger equations[END_REF] in particular implies that for ξ ∈ R u * (t n + ξ) -e iξAc u * (t n ) r ≤ ξk r (1 + M 0 ) 3 for some constant k r which is independent of c. Together with Lemma 4 and 5 the assertion then follows by integrating the highly-oscillatory phases exp ±ilc 2 ξ exactly.

In the next section we collect some important definitions and useful lemmata on highlyoscillatory integrals.

Preliminary lemmata on highly-oscillatory integrals.

-The construction of a secondorder approximation to u * based on the iteration of Duhamel's formula [START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF] that holds uniformly in all c > 0 leads to interactions of the highly-oscillatory phases e ic 2 t . More precisely, we need to handle highly-oscillatory integrals of type (48)

τ 0 e is(δc 2 -Ac) e isAc v l e -isAc v m ds, δ ∈ {-4, -2, 2}.
In order to control these integrals we first need to distinguish the non-resonant case δ ∈ {-4, -2}

where ∀c > 0, k ∈ N : (δc 2 -A c ) k = δc 2 -c c 2 + k 2 + c 2 = 0
from the resonant case δ = 2 in which the operator δc 2 -A c may become singular.

In Lemma 21 we outline how to control the non-resonant case δ ∈ {-4, -2}. Lemma 23 treats the resonant case δ = 2.

Lemma 21. -Fix r > d/2. Then we have for δ 1 = -2 and δ 2 = -4 that for j = 1, 2 and l, m ∈ N * , (49)

τ 0 e is(δ j c 2 -Ac) e isAc v l e -isAc v m ds = τ ϕ 1 iτ (δ j c 2 -A c ) v l v m + iτ 2 ϕ 2 iτ (δ j c 2 -A c ) lv l-1 v m A c v -mv l v m-1 A c v + R(t n , s, v),
where the remainder satisfies

(50) R(t n , s, v) r ≤ k r τ 3 v r+4 v l+m-1 r
for some constant k r which is independent of c.

Proof. -By Taylor series expansion of e isAc and noting [START_REF] Lawson | Generalized Runge-Kutta processes for stable systems with large Lipschitz constants[END_REF] we obtain that (51)

τ 0 e -isAc e iδ j c 2 s e isAc v l e -isAc v m ds = τ 0 e is(δ j c 2 -Ac) v l v m + is lv l-1 v m A c v -mv l v m-1 A c v ds + R(t n , s, v),
where thanks to [START_REF] Lawson | Generalized Runge-Kutta processes for stable systems with large Lipschitz constants[END_REF] we have for r > d/2 that (50) holds for the remainder. The assertion then follows by the definition of the ϕ j functions given in Definition 6.

As our numerical scheme will be built on the approximation in (49) we need to guarantee that the constructed term

τ 2 ϕ 2 iτ (δ j c 2 -A c ) lv l-1 v m A c v -mv l v m-1 A c v
is uniformly bounded with respect to c in H r for all functions v ∈ H r . This stability analysis is carried out in Remark 22 below, where we in particular exploit the bilinear estimate

(52) vw r ≤ k v r 1 w r 2 for all r ≤ r 1 + r 2 -d 2 with r 1 , r 2 , -r = d 2 and r 1 + r 2 ≥ 0.
Remark 22 (Stability in Lemma 21). -Note that for δ 1 = -2, respectively, δ 2 = -4 we have that

(53) 0 = δ j c 2 -A c = δ j c 2 -c ∇ c + c 2 = -(c 2 + c ∇ c ) if j = 1 -(3c 2 + c ∇ c ) if j = 2 .
Thanks to (53) which in particular implies that

( ∇ c ) k = c 2 + |k| 2 ≤ √ c 2 + |k| 2 = c + |k| and 1 c 2 + c ( ∇ c ) k ≤ min |c| -2 , |c c 2 + k 2 | -1 ≤ min |c| -2 , (c|k|) -1
we obtain together with the bilinear estimate (52) that for δ j = -2, -4 (54)

τ 2 ϕ 2 iτ (δ j c 2 -A c ) (vA c w) r = τ ϕ 0 (iτ (δ j c 2 -A c )) -ϕ 1 (iτ (δ j c 2 -A c )) (δ j c 2 -A c ) (vA c w) r ≤ 2τ 1 (c 2 + c ∇ c ) (vA c w) r ≤ 2τ 1 (c 2 + c ∇ c ) v2c 2 w r + 2τ 1 (c 2 + c ∇ c ) (vc ∇ 0 w) r ≤ 4k r τ v r w r
for all r > d/2 and all functions v and w and some constant k r > 0. The estimate (54) guarantees stability of our numerical scheme built on the approximation in (49).

A simple manipulation allows us to treat the resonant case, i.e., δ = 2 in (48), similarly to Lemma 21.

Lemma 23. -Fix r > d/2 and let c = 0. Then we have that

(55) τ 0 e is(2c 2 -Ac) e isAc v l e -isAc v m ds = τ ϕ 1 iτ (2c 2 -1 2 ∆) v l v m + iτ 2 ϕ 2 iτ (2c 2 -1 2 ∆) ( 1 2 ∆ -A c ) v l v m + lv l-1 v m A c v -mv l v m-1 A c v + R(t n , s, v),
where the remainder satisfies

(56) R(t n , s, v) r ≤ k r τ 3 v r+4 v l+m-1 r
for some constant k r which is independent of c.

Proof. -Note that as 2c

2 -A c = 2c 2 -1 2 ∆ + 1 2 ∆ -A c we obtain (57) τ 0 e is(2c 2 -Ac) e isAc v l e -isAc v m ds = τ 0 e is(2c 2 -1 2 ∆) e is( 1 2 ∆-Ac) e isAc v l e -isAc v m ds = τ 0 e is(2c 2 -1 2 ∆) 1 + is( 1 2 ∆ -A c ) v l v m + is lv l-1 v m A c v -mv l v m-1 A c v ds + R(t n , s, v),
where thanks to [START_REF] Lawson | Generalized Runge-Kutta processes for stable systems with large Lipschitz constants[END_REF] we have for r > d/2 that (56) holds for the remainder. The assertion then follows by the definition of the ϕ j functions given in Definition 6.

Again we need to verify that the constructed term 55) can be bounded uniformly with respect to c in H r for all functions v ∈ H r . This is done in the following remark.

τ 2 ϕ 2 iτ (2c 2 -1 2 ∆) ( 1 2 ∆ -A c ) v l v m + lv l-1 v m A c v -mv l v m-1 A c v in (
Remark 24 (Stability in Lemma 23). -Note that the operator 2c 2 -1 2 ∆ satisfies the bounds c|k|

2c 2 -1 2 ∆ k = c|k| 2c 2 + 1 2 |k| 2 ≤ 2, c 2 2c 2 -1 2 ∆ k = c 2 2c 2 + 1 2 |k| 2 ≤ 1 2
and furthermore

(A c ) k = c c 2 + |k| 2 -c 2 ≤ 2c 2 + c|k|.
The above estimates together with the bilinear estimate (52) imply that for r > d/2 (58)

τ 2 ϕ 2 iτ (2c 2 -1 2 ∆) (vA c w) 2 r ≤ τ k (1 + |k| 2 ) r (2c 2 + 1 2 |k| 2 ) 2 k=k 1 +k 2 v k 1 (A c ) k 2 w k 2 2 ≤ τ m r k (1 + |k| 2 ) r c 4 (2c 2 + 1 2 |k| 2 ) 2 k=k 1 +k 2 |v k 1 ||w k 2 | 2 + τ m r k (1 + |k| 2 ) r c 2 (2c 2 + 1 2 |k| 2 ) 2 k=k 1 +k 2 |v k 1 ||k 2 ||w k 2 | 2 ≤ τ m r k (1 + |k| 2 ) r k=k 1 +k 2 |v k 1 ||w k 2 | 2 + τ m r k (1 + |k| 2 ) r-1 k=k 1 +k 2 |v k 1 ||k 2 ||w k 2 | 2 ≤ τ m r v 2 r w 2 r + τ k r v 2 r ∂ x w 2 r-1 ≤ τ km r v 2 r w 2 r
for some constant m r > 0 which guarantees stability of the numerical method built on the approximation in Lemma 23.

Next we need to analyze integrals involving the highly-oscillatory function Ψ c 2 defined in [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF]. The following lemma yields a uniform approximation.

Lemma 25. -Fix r > d/2. Then for any polynomial p(v) in v and v we have that

(59) τ 0 e i(τ -s)Ac p e isAc v c ∇ -1 c Ψ c 2 (t n , s, v)ds = τ 2 p(v)c ∇ -1 c ϑ c 2 (t n , τ, v) + R(t n , τ, v) with (60)
ϑ c 2 (t n , τ, v) : = e 2ic 2 tn ϕ 1 2ic 2 τ -1 2iτ c 2 v 3 + 3e -2ic 2 tn ϕ 1 -2ic 2 τ -1 -2iτ c 2 |v| 2 v + e -4ic 2 tn ϕ 1 -4ic 2 τ -1 -4iτ c 2 v 3
and where the remainder satisfies

(61) R(t n , τ, v) r ≤ k r τ 3 (1 + v r+2 ) 5
for some constant k r independent of c.

Proof. -Thanks to the approximation ( 16) and the fact that Ψ c 2 (t n , s, u * (t n )) is of order one in s uniformly in c we have that

τ 0 e i(τ -s)Ac p e isAc v c ∇ -1 c Ψ c 2 (t n , s, v)ds = p (v) c ∇ -1 c τ 0 Ψ c 2 (t n , s, v)ds + R(t n , τ, v),
where the remainder satisfies for r > d/2 the bound (61).

Finally, we need to handle the interaction of highly-oscillatory phases e ilc 2 t with the highlyoscillatory function Ψ c 2 defined in [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF].

Lemma 26. -Let c = 0. Then, we have for l ∈ N that (62) Ω c 2 ,l (t n , τ, v) := 1 τ 2 τ 0 e ilc 2 s Ψ c 2 (t n , s, v)ds = e 2ic 2 tn ϕ 1 (l + 2)ic 2 τ -ϕ 1 lic 2 τ 2iτ c 2 v 3 + 3e -2ic 2 tn ϕ 1 (l -2)ic 2 τ -ϕ 1 lic 2 τ -2iτ c 2 |v| 2 v + e -4ic 2 tn ϕ 1 (l -4)ic 2 τ -ϕ 1 lic 2 τ -4iτ c 2 v 3
as well as that τ 0 e ilc 2 s sds = τ 2 ϕ 2 (ilc 2 τ ).

Proof. -Note that by Definition 19 we have that Ψ c 2 (t n , s, v) = e 2ic 2 tn e (l+2)ic 2 s -e lic 2 s 2ic 2 v 3 + 3e -2ic 2 tn e -2ic 2 s -e lic 2 s -2ic 2 |v| 2 v + e -4ic 2 tn e (l-4)ic 2 s -e lic 2 s -4ic 2 v 3 which implies the assertion by Definition 6 of ϕ 1 and ϕ 2 .

With the above lemmata at hand we can commence the construction of the second-order uniformly accurate scheme. 4.1.3. Uniform second-order discretization of Duhamel's formula. -Our starting point is again Duhamel's perturbation formula (see [START_REF] Strauss | Numerical solution of a nonlinear Klein-Gordon equation[END_REF])

u * (t n + τ ) = e iτ Ac u * (t n ) - i 8 c ∇ -1 c τ 0 e i(τ -s)Ac e -ic
I c 2 (τ, t n , u * ) := τ 0 e i(τ -s)Ac e 2ic 2 (tn+s) u 3 * (t n + s) + 3e -2ic 2 (tn+s) |u * (t n + s)| 2 u * (t n + s) + e -4ic 2 (tn+s) u * 3 (t n + s) ds.
In order to obtain a second-order uniformly accurate scheme based on the decomposition (63) we need to carefully analyze the highly-oscillatory phases in I * (τ, t n , u * ) and I c 2 (τ, t n , u * ). We commence with the analysis of I * (τ, t n , u * ).

I * (τ, t n , u * ) = e iτ Ac u * (t n ) - 3i 8 c ∇ -1 c τ 0 e i(τ -s)Ac |u * (t n + s)| 2 u * (t n + s)ds = e iτ Ac u * (t n ) - 3i 8 c ∇ -1 c I 1 * (τ, t n , u * ) + 3i 8 c ∇ -1 c i 8 τ 0 e i(τ -s)Ac 2 e isAc u * (t n ) 2 c ∇ -1 c Ψ c 2 (t n , s, u * (t n )) -e isAc u * (t n ) 2 c ∇ -1 c Ψ c 2 (t n , s, u * (t n )) ds + R(τ, t n , u * ),
where we have set

I 1 * (τ, t n , u * ) := τ 0 e i(τ -s)Ac e isAc u * (t n ) 2 e isAc u * (t n ) - 3i 4 e isAc u * (t n ) 2 c ∇ -1 c s 0 e i(s-ξ)Ac |e iξAc u * (t n )| 2 e iξAc u * (t n )dξ + 3i 8 e isAc u * (t n ) 2 c ∇ -1 c s 0 e -i(s-ξ)Ac |e iξAc u * (t n )| 2 e -iξAc u * (t n )dξ ds
and the remainder satisfies

(68) R(τ, t n , u * ) r ≤ τ 3 k r (M 4 )
for some constant k r (M 4 ) which depends on M 4 , but is independent of c. Lemma 25 allows us to handle the highly-oscillatory integrals involving the function Ψ c 2 in (67). Thus, in order to obtain a uniform second-order approximation of I * (τ, t n , u * ) it remains to derive a suitable second-order approximation to I 1 * (τ, t n , u * ). 1.1.) Approximation of I 1 * (τ, t n , u * ): The midpoint rule yields the following approximation (69)

I 1 * (τ, t n , u * ) = τ e i τ 2 Ac e i τ 2 Ac u * (t n ) 2 e i τ 2 Ac u * (t n ) - 3i 4 e i τ 2 Ac u * (t n ) 2 c ∇ -1 c τ /2 0 e i( τ 2 -ξ)Ac |e iξAc u * (t n )| 2 e iξAc u * (t n )dξ + 3i 8 e i τ 2 Ac u * (t n ) 2 c ∇ -1 c τ /2 0 e -i( τ 2 -ξ)Ac |e iξAc u * (t n )| 2 e -iξAc u * (t n )dξ + R(τ, t n , u * (t n )),
where the remainder satisfies thanks to [START_REF] Lawson | Generalized Runge-Kutta processes for stable systems with large Lipschitz constants[END_REF] and ( 18) that

(70) R(τ, t n , u * (t n )) r ≤ τ 3 k r (M 4 )
with k r independent of c.

Next we approximate the two remaining integrals in (69) with the right rectangular rule, i.e., (71)

τ /2 0 e i( τ 2 -ξ)Ac |e iξAc u * (t n )| 2 e iξAc u * (t n )dξ = τ 2 |e i τ 2 Ac u * (t n )| 2 e i τ 2 Ac u * (t n ) + R(τ, t n , u * (t n )),
where the remainder satisfies again thanks to (15) that

(72) R(τ, t n , u * (t n )) r ≤ τ 2 k r (M 4 )
with k r independent of c. Plugging (71) into (69) yields, with the notation

U * (t n ) = e i τ 2 Ac u * (t n ) (73) that I 1 * (τ, t n , u * ) = e i τ 2 Ac τ |U * (t n )| 2 U * (t n ) - τ 2 2 3i 4 |U * (t n )| 2 c ∇ -1 c |U * (t n )| 2 U * (t n ) + τ 2 2 3i 8 U * (t n ) 2 c ∇ -1 c |U * (t n )| 2 U * (t n ) + R(τ, t n , u * (t n )),
where thanks to (68), ( 70) and (72) the remainder satisfies the bound R(τ, t n , u * (t n )) r ≤ τ 2 k r (M 4 ) with k r independent of c.

In order to obtain asymptotic convergence to the classical Strang splitting scheme (3) associated to the nonlinear Schrödinger limit (2) we add and subtract the term

e i τ 2 Ac τ 2 2 3i 8 |U * (t n )| 4 U * (t n )
in the above approximation of I 1 * (τ, t n , u * ). This yields that (74)

I 1 * (τ, t n , u * ) = e i τ 2 Ac τ |U * (t n )| 2 U * (t n ) - τ 2 2 3i 8 |U * (t n )| 4 U * (t n ) - τ 2 2 3i 4 |U * (t n )| 2 c ∇ -1 c -1 |U * (t n )| 2 U * (t n ) + τ 2 2 3i 8 U * (t n ) 2 c ∇ -1 c -1 |U * (t n )| 2 U * (t n ) + R(τ, t n , u * (t n )).
The above decomposition allows us a second-order approximation of I * (τ, t n , u * ) which holds uniformly in all c:

1.2.) Final approximation of I * (τ, t n , u * ): Plugging (74) into (67) yields with the aid of Lemma 25 that

I * (τ, t n , u * ) = e i τ 2 Ac U * (t n ) - 3i 8 τ |U * (t n )| 2 U * (t n ) + - 3i 8 2 τ 2 2 |U * (t n )| 4 U * (t n ) -τ 3i 8 c ∇ -1 c -1 e i τ 2 Ac |U * (t n )| 2 U * (t n ) + τ 2 θ c ∇ c-1 (t n , τ, U * (t n )) -τ 2 3 32 c ∇ -1 c |u * (t n )| 2 c ∇ -1 c ϑ c 2 (t n , τ, u * (t n )) + τ 2 3 64 c ∇ -1 c (u * (t n )) 2 c ∇ -1 c ϑ c 2 (t n , τ, u * (t n )) + R(τ, t n , u * )
with a remainder R of order O(τ 3 ) uniformly in c. The Taylor series expansion 1 + x + x 2

2e x = O(x 3 ) furthermore allows us the following final representation of I * :

(75)

I * (τ, t n , u * ) = e i τ 2 Ac exp - 3i 8 τ |U * (t n )| 2 U * (t n ) -τ 3i 8 c ∇ -1 c -1 e i τ 2 Ac |U * (t n )| 2 U * (t n ) + τ 2 θ c ∇ c-1 (t n , τ, U * (t n )) -τ 2 3 32 c ∇ -1 c |u * (t n )| 2 c ∇ -1 c ϑ c 2 (t n , τ, u * (t n )) + τ 2 3 64 c ∇ -1 c (u * (t n )) 2 c ∇ -1 c ϑ c 2 (t n , τ, u * (t n )) + R(τ, t n , u * ) with (76) θ c ∇ c-1 (t n , τ, v) := - 1 2 9 64 e i τ 2 Ac c ∇ -1 c -1 |v| 4 v - 1 2 9 32 c ∇ -1 c e i τ 2 Ac |v| 2 c ∇ -1 c -1 |v| 2 v + 1 2 9 64 c ∇ -1 c e i τ 2 Ac v 2 c ∇ -1 c -1 |v| 2 v
and where ϑ c 2 is defined in (60) and the remainder R(τ, t n , u * ) satisfies

(77) R(τ, t n , u * (t n )) r ≤ τ 3 k r (M 4 )
with k r independent of c. The approximation of I * (τ, t n , u * ) given in (75) yields the first terms in our numerical scheme. In order to obtain a full approximation to u * (t n + τ ) in (63) we next derive a second-order approximation to I c 2 (τ, t n , u * ). where Ω c 2 ,l is defined in Lemma 26 and the remainder satisfies (80) R(t n , τ, u * ) r ≤ τ 3 k r (M 4 ), with k r independent of c.

2.) Second term I

3.) Final approximation of u * (t n + τ ): Plugging (75) as well as (78) into (63) builds the basis of our second-order scheme: As a numerical approximation to the exact solution u * at time t n+1 we take the second-order uniform accurate exponential-type integrator: U n * = e i τ 2 Ac u n * and

(81)

u n+1 * = e i τ
2 Ac e -iτ 3 8 |U n

Numerical experiments

In this section we numerically confirm first-, respectively, second-order convergence uniformly in c of the exponential-type integration schemes (28) and (81). In the numerical experiments we use a standard Fourier pseudospectral method for the space discretization with the largest Fourier mode K = 2 10 (i.e., the spatial mesh size ∆x = 0.0061) and integrate up to T = 0.1. In Figure 1 we plot (double logarithmic) the time-step size versus the error measured in a discrete H 1 norm of the first-order scheme (28) and the second-order scheme ( 81 c=1.00 c=5.00 c=10.00 c=50.00 c=100.00 c=500.00 c=1000.00 c=5000.00 c=10000.00

Figure 1. Error of the first-, respectively, second-order exponential-type integration scheme (28) and (81). The slope of the dotted and dashed line is one and two, respectively.

Assumption 9 .

 9 -Fix r > d/2 and assume that u * ∈ C([0, T ]; H r+2 (T d )) and in particular sup 0≤t≤T u * (t) r+2 ≤ M 2 uniformly in c. Applying Lemma 4 and Lemma 5 in (20) allows us the following expansion (21)

Assumption 18 .

 18 -Fix r > d/2 and assume that u * ∈ C([0, T ]; H r+4 (T d )) and in particular sup 0≤t≤T u * (t) r+4 ≤ M 4 uniformly in c.

  c 2 (τ, t n , u * ): Applying the second approximation in Lemma 20 yields together with Lemma 4 and by the definition of I c 2 (τ, t n , u * ) in (65) thatI c 2 (τ, t n , u * ) = τ 0 e i(τ -s)Ac e 2ic 2 (tn+s) e isAc u * (t n ) 3 + 3e -2ic

2 -

 2 cos(x) for different values c = 1, 5, 10, 50, 100, 500, 1000, 5000, 10000.

  2 (tn+s) e ic 2 (tn+s) u * (t n + s) + e -ic 2 (tn+s) u * (t n + s) Ac |u * (t n + s)| 2 u * (t n + s)ds and the terms involving u 3 * , u * 3 and |u * | 2 u * (65)

	which we split into two parts by separating the linear plus classical cubic part |u * | 2 u * from the terms involving u 3 * , u * 3 and |u * | 2 u * . More precisely, we set
	(63)	u * (t n + τ ) = I * (τ, t n , u * ) -	i 8	c ∇ -1 c I c 2 (τ, t n , u * )
	with the linear as well as classical cubic part |u * | 2 u *
	(64)	I * (τ, t n , u * ) := e iτ Ac u * (t n ) -	3i 8	c ∇ -1 c	0	τ	e i(τ -s)
								3	ds

  Ψ c 2 (t n , s, u * (t n )) + R 1 (t n , s, u * )with Ψ c 2 defined in (45) and where the remainder R 1 is of order O(s 2 ) uniformly in c. Plugging the approximation (66) into I * (τ, t n , u * ) defined in (64) yields that

	1.) First term I * (τ, t n , u * ): By Lemma 20 we have that
	u * (t n + s) = e isAc u * (t n ) --i 8 c ∇ -1 c (67) (66)	3i 8	c ∇ -1 c	0	s	e i(s-ξ)Ac e iξAc u * (t n )	2	e iξAc u * (t n ) dξ

  (t n )| 2 u * (t n ) + Ψ c 2 (t n , s, u * (t n )) (t n )| 2 u * (t n ) + Ψ c 2 (t n , s, u * (t n )) (t n )| 2 u * (t n ) + Ψ c 2 (t n , s, u * (t n )) (tn+s) (u * (t n )) 2 c ∇ -1 c 3s|u * (t n )| 2 u * (t n ) + Ψ c 2 (t n , s, u * (t n )) ds + R(t n , τ, u * )with Ψ c 2 defined in (45) and where thanks to Lemma 4, 20 and the fact that Ψ c 2 is of order one in s uniformly in c the remainder satisfies R(τ, t n , u * (t n )) r ≤ τ 3 k r (M 4 ) with k r independent of c.Lemma 21, 23 together with Lemma 26 thus allow us the following expansion of I c 2 : We have (78)I c 2 (τ, t n , u * ) = I 1 c 2 (τ, t n , u * ) + R(t n , τ, u * ) with the highly-oscillatory term (79) I 1 c 2 (τ, t n , u * ) := τ e 2ic 2 tn e iτ Ac ϕ 1 iτ (2c 2 -1 2 ∆) u 3 * (t n ) + iτ 2 e 2ic 2 tn e iτ Ac ϕ 2 iτ (2c 2 -1 2 ∆) ( 1 2 ∆ -A c )u 3 * (t n ) + 3u 2 * (t n )A c u * (t n ) + 3τ e -2ic 2 tn e iτ Ac ϕ 1 (iτ (-2c 2 -A c )) |u * (t n )| 2 u * (t n ) + 3iτ 2 e -2ic 2 tn e iτ Ac ϕ 2 (iτ (-2c 2 -A c )) u * 2 (t n )A c u * (t n ) -2|u * (t n )| 2 A c u * (t n ) + τ e -4ic 2 tn e iτ Ac ϕ 1 (iτ (-4c 2 -A c ))u * 3 (t n ) -iτ 2 e -4ic2 tn e iτ Ac ϕ 2 (iτ (-4c 2 -A c ))3u * 2 (t n )A c u * (t n ) -τ 2 3i 8 e 2ic 2 tn (u * (t n )) 2 c ∇ -1

	3 ds
	+ e 2ic 2 (tn+s) (u + 3e -2ic 2 (tn+s) -τ 0 -3i 8 i 8 (u + 2i 8 |u * (t n )| 2 c ∇ -1 c 3s|u + 3i 8 e -4ic 2

2

(tn+s) e isAc u * (t n ) 2 e -isAc u * (t n )

+ e -4ic 2 (tn+s) e -isAc u * (t n ) * (t n )) 2 c ∇ -1 c 3s|u * * (t n )) 2 c ∇ -1 c 3s|u * * c 3ϕ 2 (2ic 2 τ )|u * (t n )| 2 u * (t n ) + Ω c 2 ,2, (t n , τ, u * (t n )) -τ 2 3i 8 e -2ic 2 tn (u * (t n )) 2 c ∇ -1 c 3ϕ 2 (-2ic 2 τ )|u * (t n )| 2 u * (t n ) + Ω c 2 ,-2 (t n , τ, u * (t n )) + τ 2 6i 8 e -2ic 2 tn |u * (t n )| 2 c ∇ -1 c 3ϕ 2 (-2ic 2 τ )|u * (t n )| 2 u * (t n ) + Ω c 2 ,-2 (t n , τ, u * (t n )) + τ 2 3i 8 e -4ic 2 tn (u * (t n )) 2 c ∇ -1 c 3ϕ 2 (-4ic 2 τ )|u * (t n )| 2 u * (t n ) + Ω c 2 ,-4 (t n , τ, u * (t n ))

+ R(t n , τ, u * ),

Acknowledgement

The authors gratefully acknowledge financial support by the Deutsche Forschungsgemeinschaft (DFG) through CRC 1173. This work was also partly supported by the ERC Starting Grant Project GEOPARDI No 279389.

where I 1 c 2 (τ, t n , u n * ) is defined in (79) and with ϕ 1 , ϕ 2 given in Definition 6, θ c ∇ c-1 given in (76), ϑ c 2 in (60) and Ω c 2 ,l in (62).

4.2. Convergence analysis. -The exponential-type integration scheme (81) converges (by construction) with second-order in time uniformly with respect to c.

Theorem 27 (Convergence bound for the second-order scheme)

Fix r > d/2 and assume that (82) z(0) r+4 + c -1 ∇ -1 c z ′ (0) r+4 ≤ M 4 uniformly in c. For u n * defined in (81) we set

Then, there exists a T r > 0 and τ 0 > 0 such that for all τ ≤ τ 0 and t n ≤ T r we have for all c > 0 that

, where the constants K 1,r,M 2 , K 2,r,M and K * r,M,M 4 ,tn can be chosen independently of c. Proof. -First note that the regularity assumption on the initial data in (82) implies the regularity Assumption 18 on u * (t), i.e., there exists a T r > 0 such that

for some constant k that depends on M 4 and T r , but can be chosen independently of c.

In the following let φ t denote the exact flow of [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF], i.e., u * (t n+1 ) = φ τ (u * (t n )) and let Φ τ denote the numerical flow defined in (81), i.e.,

Taking the difference of ( 20) and (81) yields that

Local error bound: With the aid of the expansion (75) and (78) we obtain by the representation of the exact solution in (63) together with the error bounds (80) and ( 77) that

for some constant k r which depends on M 4 and r, but can be chosen independently of c.

Stability bound: Note that by the definition of ϕ 2 in Definition 6,

for some constant k r independent of c. Together with the bound [START_REF] Masmoudi | From nonlinear Klein-Gordon equation to a system of coupled nonlinear Schrödinger equations[END_REF], the definition of ϕ 1 in Definition 6 and the stability estimates (54) and (58) we thus obtain as long as u

, where the constant K r,M depends on r and M , but can be chosen independently of c.

Global error bound: Plugging the stability bound (86) as well as the local error bound (84) into (83) yields by a bootstrap argument that

where the constants are uniform in c. Note that as u = v we have by ( 7) and ( 11) that

Together with the bound in (87) this completes the proof.

Remark 28 (Fractional convergence and convergence in L 2 ) A fractional convergence result as Theorem 13 for the first-order scheme also holds for the second-order exponential-type integrator (81): Fix r > d/2 and let 0 ≤ γ ≤ 1. Assume that

Then, the scheme (81) is convergent of order τ 1+γ in H r uniformly with respect to c. Furthermore, for initial values satisfying

the exponential-type integration scheme (81) is second-order convergent in L 2 uniformly with respect to c by the strategy presented in [START_REF] Lubich | On splitting methods for Schrödinger-Poisson and cubic nonlinear Schrödinger equations[END_REF].

In analogy to Remark 15 we make the following observation: For sufficiently smooth solutions the exponential-type integration scheme (81) converges in the limit c → ∞ to the classical Strang splitting of the corresponding nonlinear Schrödinger limit equation [START_REF] Bao | Uniformly accurate multiscale time integrators for highly oscillatory second order differential equations[END_REF].

Remark 29 (Approximation in the non relativistic limit c → ∞)

The exponential-type integration scheme (81) corresponds for sufficiently smooth solutions in the limit u n * c→∞ -→ u n * ,∞ , essentially to the Strang Splitting ( [START_REF] Lubich | On splitting methods for Schrödinger-Poisson and cubic nonlinear Schrödinger equations[END_REF][START_REF] Faou | Geometric numerical integration and Schrödinger equations[END_REF])

for the cubic nonlinear Schrödinger limit system (2). More precisely, the following Lemma holds.

for some ε > 0 uniformly in c and let the initial value approximation (there exist functions ϕ, γ such that)

Then, there exists a T > 0 and τ 0 > 0 such that for all τ ≤ τ 0 the difference of the second-order scheme (81) for system [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF] and the Strang splitting (88) for the limit Schrödinger equation (2) satisfies for t n ≤ T and all c > 0 with τ c ≥ 1

for some constant k r that depends on M 3 and T , but is independent of c.

Proof. -The proof follows the line of argumentation to the proof of Lemma 15 by noting that for l = -4, -2 and n = -4, -2, 2

for some constant k r independent of c.

4.3.

Simplifications in the "weakly to strongly non-relativistic limit regime". -In the "weakly to strongly non-relativistic limit regime", i.e., for large values of c, we may again (substantially) simplify the second-order scheme (81) and nevertheless obtain a well suited, second-order approximation to u * (t n ) in [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF].

Remark 31 (Limit scheme [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF]). -For sufficiently large values of c and sufficiently smooth solutions, more precisely, if z(0) r+4 + c -1 ∇ -1 c z ′ (0) r+4 ≤ M 4 and τ c > 1 we may take instead of (81) the classical Strang splitting (see [START_REF] Lubich | On splitting methods for Schrödinger-Poisson and cubic nonlinear Schrödinger equations[END_REF][START_REF] Faou | Geometric numerical integration and Schrödinger equations[END_REF]) for the nonlinear Schrödinger limit equation ( 2 as a second-order numerical approximation to u * (t n ) in [START_REF] Ostermann | Low regularity exponential-type integrators for semilinear Schrödinger equations in the energy space[END_REF]. The assertion follows from [START_REF] Faou | Asymptotic preserving schemes for the Klein-Gordon equation in the nonrelativistic limit regime[END_REF] thanks to the approximation