Low Regularity Exponential-Type Integrators for Semilinear Schrödinger Equations

We introduce low regularity exponential-type integrators for nonlinear Schrödinger equations for which first-order convergence only requires the boundedness of one additional derivative of the solution. More precisely, we will prove first-order convergence in Hr\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H^r$$\end{document} for solutions in Hr+1\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H^{r+1}$$\end{document} (with r>d/2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$r > d/2$$\end{document}) of the derived schemes. This allows us lower regularity assumptions on the data than for instance required for classical splitting or exponential integration schemes. For one-dimensional quadratic Schrödinger equations, we can even prove first-order convergence without any loss of regularity. Numerical experiments underline the favorable error behavior of the newly introduced exponential-type integrators for low regularity solutions compared to classical splitting and exponential integration schemes.

Semilinear Schrödinger equations, in particular those of type (1) i∂ t u = −∆u + µ|u| 2p u, p ∈ N with µ ∈ R are nowadays extensively studied numerically.In this context, splitting methods (where the right-hand side is split into the kinetic and nonlinear part, respectively) as well as exponential integrators (including Lawson-type Runge-Kutta methods) contribute particularly attractive classes of integration schemes.For an extensive overview on splitting and exponential integration methods we refer to [15,16,18,26], and for their rigorous convergence analysis in the context of semilinear Schrödinger equations we refer to [3,4,6,8,9,11,12,25,31] and the references therein.However, within the construction of all these numerical methods the stiff part (i.e., the terms involving the differential operator ∆) is approximated in one way or another which generally requires the boundedness of two additional spatial derivatives of the exact solution.In particular, convergence of a certain order only holds under sufficient additional regularity assumptions on the solution.In the following we will illustrate the local error behavior of classical splitting and exponential integration methods and the thereby introduced smoothness requirements.
Splitting schemes.The Strang splitting for the cubic Schrödinger equation (i.e., p = 1) , where the right-hand side is split into the kinetic T (u) = i∆u and nonlinear part V (u) = −iµ|u| 2 u is rigorously analyzed in [25].In particular, its second-order convergence in H r for solutions in H r+4 and its first-order convergence in H r for solutions in H r+2 holds for all r ≥ 0, see [25] and [10], respectively.The result follows from the fact that the local error of the splitting scheme can be expressed through the double Lie commutator [T, [T, V ]](u) and the Lie commutator [T, V ](u) for second and first-order methods, respectively.The latter reads see [25,Section 4.2].Due to the appearance of ∆u in the local error, the boundedness of at least two additional derivatives of the exact solution is required.Based on [25], fractional error estimates for the Strang splitting were established in [10] which require the boundedness of 2 + 2γ additional derivatives for convergence of order 1 + γ, with 0 < γ < 1.Furthermore, in [19] first-order convergence of a filtered Lie splitting method for Schrödinger equations of type (1) with 1 ≤ 2p < 4/d was shown in L 2 (R d ) for solutions in H 2 (R d ).
Its construction is based on Duhamel's formula u(t n + τ ) = e iτ ∆ u(t n ) − iµ τ 0 e i(τ −s)∆ |u(t n + s)| 2p u(t n + s)ds by applying the approximation in the integral terms and solving the remaining integral over the free Schrödinger group e i(τ −s)∆ exactly.Due to the approximation (3) the above scheme introduces a local error of the form see [16].Hence, first-order convergence also requires the boundedness of at least two derivatives of the exact solution due to the appearance of ∆u and its complex conjugate counterpart.
The main novelty in this work lies in the development and analysis of a first-order exponentialtype integrator for Schrödinger equations of type (1) (4) with ϕ 1 (z) = e z − 1 z which only requires the boundedness of one additional derivative of the exact solution.The construction is based on Duhamel's formula looking at the twisted variable v(t) = e −it∆ u(t) and treating the dominant term triggered by u 2p in an exact way.This idea of twisting the variable is widely used in the analysis of partial differential equations in low regularity spaces (see for instance [1]) and also well known in the context of numerical analysis, see [24] for the introduction of Lawson-type Runge-Kutta methods.However, instead of approximating the appearing integrals with a Runge-Kutta method (see for instance [20]) we integrate the dominant second-order stiff parts exactly.
A similar approach has been successfully applied to the one dimensional KdV equation, see [17].However, due to the Burgers-type nonlinearity, additional smoothness assumptions on the exact solution are necessary.More precisely, first-order convergence in H 1 is only guaranteed for solutions in H 3 in the KdV setting.The introduced exponential integrators for nonlinear Schrödinger equations are in contrast first-order convergent in H r for solutions in H r+1 (with r > d/2).Furthermore, we show that our approach yields methods for quadratic nonlinearities of type u 2 and |u| 2 in one dimension which are convergent in H r for solutions in H r , i.e., no additionalsmoothness of the solutionis required.Such generous convergence results for nonlinear Schrödinger equations are up to our knowledge not yet know in the numerical analysis literature.
For practical implementation issues we impose periodic boundary conditions and refer to [1] for local wellposedness (LWP) results of nonlinear Schrödinger equations in H r (T d ).Note that in the case of a cubic nonlinearity we obtain LWP for r > 1/2 if d ≤ 3.In the case of a quadratic nonlinearity of type u 2 (or u 2 ) LWP can be pushed down to r > −1/2 on the one dimensional torus, see [21].In the following let r > d/2.We denote by • r the standard H r = H r (T d ) Sobolev norm.In particular, we exploit the well-known bilinear estimate (5) f g r ≤ c r,d f r g r which holds for some constant c r,d > 0. Furthermore, we denote by c a generic constant which may depend on r, d, p and γ.

Low regularity exponential integrators
We consider nonlinear Schrödinger equations (NLS) of type with µ ∈ R and where we have set ∆ = d k=1 ∂ 2 x k .Furthermore, we employ the so-called ϕ 1 function defined as ϕ 1 (z) = e z − 1 z .
In the following we derive exponential-type integrators for Schrödinger equations of type (6) based on iterating Duhamel's formula in the twisted variable v(t) = e −it∆ u(t).Note that the twisted variable satisfies with mild solution given by ( 8) In order to derive our numerical scheme we proceed as follows.As e it∆ is a linear isometry on H r for all t ∈ R and the bilinear estimate (5) holds we obtain that ( 9) In this sense we have for |s| ≤ τ for a small time step τ .
Remark 1.1.-Note that the twisted variable allows the formal expansion such that the approximation (10) holds without additional regularity assumptions (see (9) for the rigorous estimate).The classical exponential integrator (2), in contrast, is based on the approximation (3) for the original solution u.This approximation error is small only under additional smoothness assumptions, since Plugging (10) into (8) yields the approximation which is the basis of our numerical scheme.Hence, we are left with deriving a numerical approximation to the integral e −i(tn+s)∆ |e i(tn+s)∆ w| 2p e i(tn+s)∆ w ds.
To illustrate the idea we will first consider the cubic case p = 1 in Section 1.1 below.In Section 1.2 we will deal with general nonlinearities p ≥ 1.The special case of quadratic nonlinearities of type u 2 and |u| 2 , respectively, will be treated in Section 3.
1.1.Cubic nonlinearities p = 1.-For notational simplicity we first illustrate the idea in one dimension.In Section 1.1.2we will give the generalization to arbitrary dimensions d ≥ 1.
. Then we will denote its Fourier expansion by f (x) = k∈Z fk e ikx .Furthermore, we define a regularization of ∂ −1 x through its action in Fourier space by and, by continuity, In the case of a cubic nonlinearity (i.e., p = 1 in ( 6)) in one spatial dimension the integral in (11) can be expressed in terms of the Fourier expansion as follows ( 13) Naturally, the next step would be the exact integration of the appearing integral ( 14) However, the obtained convolution can not be converted into the physical space straightforwardly and therefore only yields a practical scheme in dimension one as the iteration needs to be carried out fully in Fourier space.In order to obtain an efficient and practical low regularity implementation we therefore only treat the dominant quadratic term 2k 2 1 exactly.For this purpose, we write the integrand in ( 14) as 15) into ( 13) and integrating the first term yields that ( 16) where Since |β| −γ (e iβ − 1) is uniformly bounded for β ∈ R and 0 ≤ γ ≤ 1, the integral in the above formula is of the order τ 1+γ .This shows that ( 17) for some constant c > 0 and 0 ≤ γ ≤ 1.
Next we define the auxiliary function g(x) = k∈Z ĝk e ikx through its Fourier coefficients which allows us to express the bound on the remainder given in (17) as follows With the aid of the bilinear estimate ( 5) we can thus conclude that for r > 1/2 , where we have used the definition of the Fourier coefficients in (18).Hence, we can conclude the following bound on the remainder for some constant c > 0.
Inserting the dominant term of ( 16) into (11) and twisting the solution back again, i.e., setting yields the following exponential-type integration scheme for the cubic NLS ( 20) Next we generalize the above scheme to arbitrary dimensions d ≥ 1.
-In the following we use the notation In the case of a cubic nonlinearity p = 1 in dimension d the integral in ( 11) can be expressed in terms of the Fourier expansion as follows Note that we have In order to obtain an efficient and practical low regularity implementation we again treat the dominant quadratic term 2j 2 exactly.This yields that (21) where similarly to above (see (19)) we have that for some constant c > 0.
Inserting the approximation ( 21) into (11) and twisting the solution back again, i.e., setting yields the generalization of the exponential-type integrator (20) to dimensions d ≥ 1 Next we consider the general case of (6) with p ∈ N.
1.2.Nonlinearities with integer p. -In the case of a general nonlinearity p ∈ N in (6) the integral in (11) can be expressed in terms of the Fourier expansion as follows e −i(tn+s)∆ |e i(tn+s)∆ w| 2p e i(tn+s)∆ w ds.
Similar observations to above (see also [28]) yield the following numerical scheme for some constant c > 0.
In the original variable u the exponential-type integration scheme for the nonlinear Schrödinger equation ( 6) reads ( 25) In the following section we give an error analysis for the above scheme.

Error analysis
In this section we carry out the error analysis of the exponential-type integrator (22).In the following, we set Then, for all t ∈ R we have where L depends on f r and g r .
Proof.-Note that for all t ∈ R and w ∈ H r it holds by (12) that Thus, as e it∆ is a linear isometry on H r and H r is an algebra for r > d/2 we obtain that This proves the assertion.
, where φ t denotes the exact flow of (7).Then, Proof.-Note that by (23) we have that Furthermore, (8) implies that Thus as e it∆ is a linear isometry on H r and H r is an algebra for r > d/2 we obtain that Together with the estimate on the difference v(t k + s) − v(t k ) given in (9) and the estimate on the remainder R τ p (v, t) given in (24) this yields for r > d/2 and 0 ≤ γ ≤ that r+γ which concludes the stated local error bound.
These two lemmata allow us to prove the following convergence bound.
Theorem 2.3.-Let r > d/2 and 0 < γ ≤ 1. Assume that the exact solution of (7) satisfies v(t) ∈ H r+γ for 0 ≤ t ≤ T .Then, there exists a constant τ 0 > 0 such that for all step sizes τ ≤ τ 0 and t n ≤ T we have that the global error of ( 22) is bounded by Proof.-The triangle inequality yields that . Thus, we obtain with the aid of Lemmata 2.2 and 2.1 that as long as v k ∈ H r we have that where c depends on sup t k ≤t≤t k+1 v(t) r+γ and L depends on v(t k ) r as well as on v k r .As long as v(t k ) ∈ H r+γ and v k ∈ H r for 0 ≤ k ≤ n we therefore obtain by iterating the estimate (26 The assertion then follows by induction, respectively, a Lady Windermere's fan argument (see, for example [11,14,25]).
The above theorem immediately yields a convergence result for the exponential-type integration scheme (22) approximating the solution of the nonlinear Schrödinger equation ( 6).
Corollary 2.4.-Let r > d/2 and 0 < γ ≤ 1. Assume that the exact solution of (6) satisfies u(t) ∈ H r+γ for 0 ≤ t ≤ T .Then, there exists a constant τ 0 > 0 such that for all step sizes τ ≤ τ 0 and t n ≤ T we have that the global error of (25) is bounded by u(t n ) − u n r ≤ cτ γ , where c depends on sup 0≤t≤T u(t) r+γ .
Proof.-As e it∆ is a linear isometry on H r for all t ∈ R we have that Thus, the assertion follows from Theorem 2.3.
Here, the key-relations (29) ( allow us to derive first-order exponential-type integrators for Schrödinger equations of type ( 27) and ( 28) without imposing any additional regularity assumptions on the solution.
3.1.Quadratic Schrödinger equations of type (27).-Set v(t) = e −it∂ 2 x u(t), where u is the solution of the Schrödinger equation (27).Then the twisted variable v(t) satisfies (31) i∂ which by Duhamel's formula yields that To construct our numerical scheme we proceed as above.We start from the following approximation and rewrite the above integral in Fourier space as follows The key-relation ( 29) now shows that Together with the approximation in (33) this yields the following integration scheme: (34) In order to obtain an approximation to the original solution u(t) of the quadratic Schrödinger equation ( 27) at time t n = nτ we twist the variable back again, i.e., we set u n := e itn∂ 2 x v n .This yields the following exponential-type integrator 3.2.Quadratic Schrödinger equations of type (28).-For the quadratic Schrödinger equation (28) the equation in the twisted variable v(t) = e −it∂ 2 x u(t) reads To construct our numerical scheme we proceed as above.We start from the following approximation and rewrite the above integral in Fourier space as follows The key-relation (30) now shows that Together with the approximation in (36) this yields (by twisting the variable back again) the following integration scheme in the original solution u(t) = e itn∂ 2 x v(t) (37) 3.3.Error analysis.-In this section we carry out the error analysis of the exponential-type integrators ( 35) and (37).In the following let r > 1/2.We commence with the quadratic Schrödinger equation of type (27).Let φ t denote the exact flow of (31), i.e., v(t) = φ t (v(0)).The following lemmata provide a local error bound of the scheme (34).Lemma 3.1.-Let r > 1/2 and assume that v(t k + t) = φ t (v(t k )) ∈ H r for 0 ≤ t ≤ τ , where φ t denotes the exact flow of (27).Then x is a linear isometry on H r for all t ∈ R we obtain by taking the difference of (32) with the approximation (33) and using the bilinear estimate (5) that where c depends on sup 0≤s≤τ φ s (v(t k )) r .This yields the assertion.
Next we state the stability estimate.
Lemma 3.2.-Let r > 1/2 and f, g ∈ H r .Then, for all t ∈ R we have Proof.-The assertion follows from the representation of the numerical flow given in (33) together with the bilinear estimate (5).
The above lemmata allow us to prove the following global convergence result.
Theorem 3.3.-Let r > 1/2.Assume that the exact solution of (31) satisfies v(t) ∈ H r for 0 ≤ t ≤ T .Then, there exists a constant τ 0 > 0 such that for all step sizes τ ≤ τ 0 and t n ≤ T we have that the global error of (34) is bounded by Proof.-The proof follows the line of argumentation to the proof of Theorem 2.3 thanks to the local error bound in Lemma 3.1 and the stability estimate of Lemma 3.2.
The above theorem immediately yields first-order convergence of the exponential-type integrator (35) (in H r ) without any additional smoothness assumptions on the exact solution of (27).
Corollary 3.4.-Let r > 1/2.Assume that the exact solution of (27) satisfies v(t) ∈ H r for 0 ≤ t ≤ T .Then, there exists a constant τ 0 > 0 such that for all step sizes τ ≤ τ 0 and t n ≤ T we have that the global error of (35) is bounded by u(t n ) − u n r ≤ cτ, where c depends on sup 0≤t≤T u(t) r .
Proof.-The bound follows directly from Theorem 3.3 as e it∂ 2 x is a linear isometry on H r for all t ∈ R.
Similarly to above we obtain the following first-order convergence result for the exponentialtype integrator (37) approximating the solution of the quadratic Schrödinger equation (28).
Corollary 3.5.-Let r > 1/2.Assume that the exact solution of (28) satisfies v(t) ∈ H r for 0 ≤ t ≤ T .Then, there exists a constant τ 0 > 0 such that for all step sizes τ ≤ τ 0 and t n ≤ T we have that the global error of (37) is bounded by u(t n ) − u n r ≤ cτ, where c depends on sup 0≤t≤T u(t) r .
Proof.-The proof is similar to that of Corollary 3.4 and will therefore be omitted here.

Numerical experiments
In this section we compare the newly developed exponential-type integrator (25) applied to Schrödinger equations with low regularity solutions in the energy space with the classical Lie splitting, the classical Strang splitting and the classical first-order exponential integrator.In the numerical experiments we use a standard Fourier pseudospectral method for the space discretization where we choose the largest Fourier mode K = 2 10 (i.e., the spatial mesh size ∆x = 0.0061).The numerical experiments show that the new integrator ( 25) is preferable over the classical splitting and exponential integration methods for low regularity solutions.In particular, the experiments indicate that the proposed integrator (25) is convergent of order one even for rough solutions whereas the classical methods suffer from order reduction.
Although the new integrator ( 25) was derived for nonlinearities with integer p ≥ 1, it can be applied to problem (6) with non-integer p as well.Preliminary numerical results for small non-integer values of p are presented in Section 4.3.(6).-In this section we consider nonlinear Schrödinger equations of type (6).The associated Lie splitting (of classical order one) reads (38)

Numerical experiments for Schrödinger equations of type
, and the associated Strang splitting (of classical order two) is given by (39) Furthermore, the classical first-order exponential integrator (of classical order one) is defined through In our numerical experiments we choose µ = 1 (unless otherwise noted) and set (41) where x j = (−1 + j K )π and rand(2K, 1) returns 2K uniformly distributed random numbers between 0 and 1.We compare the above schemes with the newly derived exponential-type integrator (25) for initial values (42) for different values of ϑ normalized in L 2 .(For typical initial values, see Figure 1.)All experiments are carried out with constant step sizes τ = j 512 for j = 1, . . ., 512.(6).The error at time T = 1 measured in a discrete H 1 norm of the exponential-type integrator (25), the Lie splitting scheme (38), the Strang splitting scheme (39) as well as the classical exponential integrator (40) for the initial value (42) normalized in L 2 is illustrated in Figure 2 for ϑ = 3/2, respectively, ϑ = 2 and in Figure 3 for ϑ = 3, respectively, ϑ = 5.As a reference solution we take the schemes themselves with a very small time step size if ϑ = 5.For ϑ = 5 we take the Strang splitting scheme (39) with a very small time step size as the reference solution for all schemes.The experiments show that Strang splitting is second-order convergent for ϑ ≥ 5, and that Lie splitting is first-order convergent for ϑ ≥ 3.These observations are in line with the convergence proofs from the literature [10,25], see also our discussion in the introduction.For smaller values of ϑ, both Lie and Strang splitting show a zig-zag behaviour, depending on whether the local errors accumulate or happy error cancellation occurs.The classical exponential integrator of order one is first-order convergent for ϑ ≥ 3 as expected (see the discussion in the introduction).For smaller values of ϑ the convergence behaviour of the exponential integrator gets less regular and order reduction occurs.On the other hand, the new exponential-type integrator introduced in this paper is first-order convergent for all ϑ ≥ 2. This behaviour is in line with Theorem 2.3.The integrator even shows first-order convergence and small errors for ϑ = 3 2 .4.1.2.Quintic Schrödinger equations.-In this section we consider the quintic Schrödinger equation, i.e., p = 2 in (6), which appears for instance as the mean field limit of a Boson gas with three-body interactions, see [7].The error at time T = 1 measured in a discrete H 1 norm of the exponential-type integrator (25), the Lie splitting scheme (38), the Strang splitting scheme (39) as well as the classical exponential integrator (40) for the initial value (42) normalized in L 2 is illustrated in Figure 4 for ϑ = 3/2, respectively, ϑ = 2 and in Figure 5 for ϑ = 3, respectively, ϑ = 5.As a reference solution we take the schemes themselves with a very small time step size if ϑ = 5.For ϑ = 5 we take the Strang splitting scheme (39) with a very small time step size as the reference solution for all schemes.The outcome of these numerical experiments is almost the same as for the cubic nonlinear Schrödinger equation.The only notable difference is   (27).In order to derive the splitting methods, we split the right-hand side of ( 27) into the linear and nonlinear part, respectively.Note that (for sufficiently small t) the exact solution of the subproblem iu ′ (t) = µu(t) 2 is given by .
Thus, the associated Lie splitting scheme (of classical order one) reads (43) and the associated Strang splitting scheme (of classical order two) is given by (44) In Examples 4.1 and 4.2 we compare the above splitting methods as well as the classical firstorder exponential integration scheme (45) with the newly derived exponential-type integrator (35) for smooth and non-smooth solutions, respectively, where we set µ = 1 and integrate up to T = 1.The numerical experiments in particular indicate that the exponential-type integrator (35) is preferable over the commonly used splitting methods for low regularity solutions.This is due to the fact that within the exponential integration scheme all the stiff parts (i.e., the terms involving the differential operator ∂ 2 x ) are solved exactly.
In Example 4.3 we consider a quadratic Schrödinger equation of type (27) with a small nonlinearity on a longer time interval T = 10.The numerical experiment indicates the favorable behavior of the exponential integrator (35) even in the smooth case.This is due to the fact that its error constant is triggered by the nonlinearity (and not by the differential operator ∂ 2 x ).-We choose the non-smooth initial value (41) normalized in L 2 .The error at time T = 1 measured in a discrete L 2 norm of the exponential-type integrator (35), the Lie splitting scheme (43), the Strang splitting scheme (44) and the classical exponential integrator (45) is illustrated in Figure 7.As a reference solution we take the schemes themselves with a very small time step size.For this example, the classical exponential integrator is not convergent.Lie and Strang splitting suffer from a strong order reduction (down to order one half).Here, the new exponential-type integrator is by far the best method.It is first-order convergent, as predicted by Theorem 3.3.
Example 4.3 (Small nonlinearity).-We consider the quadratic Schrödinger equation of type (27) with small nonlinearity, i.e., (47) i∂ t u = −∂ 2 x u + µu 2 , µ = 0.01, T = 10.The error at time T = 10 measured in a discrete L 2 norm of the exponential-type integrator (35), the Lie splitting scheme (43), the Strang splitting scheme (44) and the classical exponential integrator (45) for smooth and non-smooth solutions is illustrated in Figure 8 left and right, respectively.As a reference solution we take the Strang splitting scheme for smooth and the schemes themselves for non-smooth solutions with a very small time step size.For smooth initial data, all methods show their classical orders of convergence.The new exponential-type integrator, however, is by far the most accurate scheme (for the considered range of step sizes).For non-smooth initial data, the classical exponential integrator does not converge; Lie and Strang splitting suffer again from a strong order reduction down to order one half.Only the new exponential-type integrator is very accurate and first-order convergent, as expected.(magenta, star), classical exponential integrator (45) (green, dashed) and exponentialtype integration scheme (35) (red, dash dotted) for a small nonlinearity (47).Left picture: Smooth initial value (46).The slope of the continuous lines is one and two, respectively.Right picture: Non-smooth initial value (41).The slope of the continuous lines is one half, and one, respectively.

4.3.
Numerical experiments for semilinear Schrödinger equations with non-integer exponents p.-In this section we numerically investigate the convergence rates in the case of semilinear Schrödinger equations of type (6) with non-integer exponents p.In Figure 9 we plot the error at time T = 1 measured in a discrete H 1 norm of the exponential-type integrator (25), the Lie splitting scheme (38), the Strang splitting scheme (39) as well as the classical exponential integrator (40) for different values of p > 0. In all the simulations we take the smooth initial 4 (lower row, right).Error in the energy space H 1 of the Lie splitting (38) (blue, dotted circle), Strang splitting (39) (magenta, star), classical exponential integrator (40) (green, dashed) and exponential-type integration scheme (25) (red, dash dotted) with smooth initial value (48).The slope of the continuous and dotted line is one and two, respectively.

4. 2 .
Numerical experiments for quadratic Schrödinger equations.-In section we consider quadratic Schrödinger equations of type

Example 4 . 1 (
Smooth solutions).-We choose the smooth initial value(46) u(0, x) = sin x cos x ∈ C ∞ (T)normalized in L 2 .The error at time T = 1 measured in a discrete L 2 norm of the exponentialtype integrator (35), the Lie splitting scheme (43), the Strang splitting scheme (44), and the classical exponential integrator (45) is illustrated in Figure6.As a reference solution we take the Strang splitting scheme (44) with a very small time step size.All the methods show their classical orders of convergence, as expected.