# Asymptotic expansions for the linear PDEs with oscillatory input terms;
Analytical form and error analysis

Karolina Kropielnicka
Institute of Mathematics
Polish Academy of Sciences
A. Abrahama 18, 81-825 Sopot
Poland  Rafał Perczyński
Institute of Mathematics
University of Gdańsk
W. Stwosza 57, 80-308 Gdańsk
Poland

###### Abstract

Partial differential equations with highly oscillatory input terms are hardly
ever solvable analytically and their numerical treatment is difficult.
Modulated Fourier expansion used as an ansatz is a well known and extensively
investigated tool in asymptotic numerical approach for this kind of problems.
Although the efficiency of this approach has been recognised, its error
analysis has not been investigated rigorously for general forms of linear
PDEs. In this paper, we start such kind of investigations for a general form
of linear PDEs with an input term characterised by a single high frequency.
More precisely we derive an analytical form of such an expansion and provide a
formula for the error of its truncation. Theoretical investigations are
illustrated by computational simulations.

##  1 Introduction

In this paper we are concerned with the linear PDEs featured by highly
oscillatory input term $f(x,t)$:

$\displaystyle\partial_{t}u(x,t)=\mathcal{L}u(x,t)+f(x,t)u(x,t),\qquad
t\in[0,t^{\star}],\ x\in\Omega\subset\mathbb{R}^{m},$ (1.1) $\displaystyle
u(x,0)=u_{0}(x),$

where $\Omega$ is a bounded subset of $\mathbb{R}^{m}$ with smooth boundary
$\partial\Omega$ while $\mathcal{L}$ is a differential operator of degree $p$
acting in the following way

$\displaystyle\mathcal{L}v(x)=\sum_{|\textbf{p}|\leq
p}a_{\textbf{p}}(x)\frac{\partial^{|\mathbf{p}|}v(x)}{\partial
x_{1}^{p_{1}}\dots\partial
x_{m}^{p_{m}}},\quad\mathbf{p}=(p_{1},\dots,p_{m}),\ p>0,\ \quad x\in\Omega.$
(1.2)

The type of the boundary conditions does not need to be specified for the
purpose of this paper and depends on the applications. We require, however,
the whole initial boundary problem of considerations to be well-posed. The
main assumption is the oscillatory in time nature of the input term, more
precisely in this manuscript we assume that

$f(x,t)=\alpha(x){\mathrm{e}}^{{\mathrm{i}}\omega t},\quad\omega\gg 1,$ (1.3)

where parameter $\alpha(x)$ is independent of $\omega$. Let us observe that
the presented theory is also applicable to function
$f(x,t)=\beta(x)+\alpha(x){\mathrm{e}}^{{\mathrm{i}}\omega t}$, which can be
obtained by incorporating (also independent of $\omega$) $\beta(x)$ into the
linear operator $\mathcal{L}$.

The aim of this paper is to express the solution of (1.1) as an asymptotic
series

$u(x,t)\thicksim
p_{0,0}(x,t)+\sum_{r=1}^{\infty}\frac{1}{\omega^{r}}\sum_{n=0}^{\infty}p_{r,n}(x,t){\mathrm{e}}^{{\mathrm{i}}n\omega
t},$ (1.4)

which is also known as modulated Fourier expansion. In this paper we make the
first step towards the derivation of an asymptotic expansion of the solution
of (1.1) for a more general highly oscillatory function $f(x,t)$ of the
following form

$f(x,t)=\sum_{n=-N}^{N}\alpha_{n}(x,t){\mathrm{e}}^{{\mathrm{i}}\omega_{n}t},\quad|\omega_{n}|\gg
1,\quad N\in\mathbb{N},$ (1.5)

where $\alpha_{n}(x,t)$ are independent of $\omega_{m}$, $n,m\in\mathbb{N}$,
and where possibly occurring resonances will require additional attention.

The numerical-asymptotic approach based on (1.4) for problem of type (1.1)
with (1.5) has been widely applied to many differential problems. For example
paper [Condon et al.(2009)Condon, Deaño, and Iserles] engages this methodology
for problems in electronic engineering, [Condon et al.(2010)Condon, Deaño, and
Iserles] for systems of differential equations with extrinsic oscillation,
[Condon et al.(2012)Condon, Deaño, Iserles, and Kropielnicka] for ODEs with
time delay, while wave-like problems were threated in [Faou and
Schratz(2014)], [Condon et al.(2019)Condon, Iserles, Kropielnicka, and Singh]
and [Condon et al.(2021)Condon, Kropielnicka, Lademann, and Perczyński]. This
methodology bases on establishing the ansatz solution of type (1.4) for the
investigated problem. Incorporating the ansatz to the differential equation
allows for establishing recurrent or simple, non-oscillatory differential
relations for the coefficients $p_{r,n}(x,t)$. Final calculations of the
coefficients often require numerical treatment of simple, non-oscillatory
PDEs. For this reason methods based on the described approach are called
numerical-asymptotic solvers. Numerical-asymptotic methods stem from the
theory described in [Hairer et al.(2010)Hairer, Lubich, and Wanner], where the
theoretical fundaments has been presented for slightly simpler case of
expansions (where summation over powers of magnitude $\omega^{-r}$ is not
taken into account).

The novelty of the paper is twofold: We derive (1.4) purely analytically
(presenting the formulas for $p_{r,n}$) and we provide the formula for the
error of the approximation. Unlike in the cited publications we do not assume
that ansatz (1.4) solves (1.1) but derive it making use of Duhamel formula,
Neumann series, and of integration by parts.

To present the error of the approximation we observe that in practice instead
of an infinite series (1.4) we rather resort to partial sums:

$u(x,t)=p_{0,0}(x,t)+\sum_{r=1}^{R}\frac{1}{\omega^{r}}\sum_{n=0}^{N}p_{r,n}(x,t){\mathrm{e}}^{{\mathrm{i}}n\omega
t}+\mathcal{R}_{R,N},\qquad t\in[0,t^{\star}],\
x\in\Omega\subset\mathbb{R}^{m},$ (1.6)

and are interested in expression $\mathcal{R}_{R,N}$. According to our
knowledge no papers concerned with the numerical-asymptotic solvers based on
modulated Fourier expansions provide neither formula $\mathcal{R}_{R,N}$ nor
its error bounds. All that has been concluded by now is that the leading error
term (not the error!) of such truncation scales like $\omega^{-R-1}$. This is
however misleading, because the error constant depends on the differential
operator $\mathcal{L}$, on the initial condition $u_{0}(x)$ and on the
parameters $\alpha_{n}(x,t)$. Let us emphasize that the error constant of the
leading error term has never been investigated and the form of dependence on
the aforementioned parameters has not been known. In this manuscript we
provide analytical formulas not only for functions $p_{r,n}(x,t)$, but also
for the error term $\mathcal{R}_{R,N}$. In contrast to previous research we
investigate the whole error term of the truncation (not only the leading error
term). Moreover knowledge of the exact formula $\mathcal{R}_{R,N}$ allows for
the realistic estimate of the error. This is important, because it is widely
known that modulated Fourier expansions do not have to be convergent. This may
happen once error constant $C_{R}$ of the error term $\mathcal{R}_{R,N}$
exceeds the magnitude of $\omega^{-R-1}$. For this particular reason the
knowledge of the formula $\mathcal{R}_{R,N}$ is so important.

Modulated Fourier expansions of type (1.4) have been successfully applied to
the problems, where oscillatory nature of solutions has intrinsic origin
rather than extrinsic, see for example [Ernst and Christian(2001), Gauckler et
al.(2018)Gauckler, Hairer, and Lubich, Wang and Zhao(2021), Wang and
Zhao(2023)]. For obvious reasons our study does not cover these problems.

In Section 2 we present solution of the underlying problem with the Duhamel
formula and will present it as a series of multivariate integrals. Section 3
is devoted to the asymptotic expansion of these multivariate integrals
derived. The main theorem on the asymptotic expansion of the solution to (1.1)
is presented in Section 4. Numerical simulations are presented in Section 5.

##  2 The solution as a series of the multivariate oscillatory integrals

In this section we will show that function $u\in
C^{1}([0,t^{\star}),H^{p}(\Omega))$ solving problem (1.1) with $u_{0}\in
H^{p}(\Omega)$ can be represented as an infinite sum of multivariate highly
oscillatory integrals which, as will be mentioned in Remark 2.2 and later
claimed in Theorem 3.1, scale like $\mathcal{O}(\omega^{-d}),\ d=0,1,2,\dots$.

Let us recall that for the sake of generality of applications of the presented
methodology we do not impose any boundary conditions, but require well
possedness of the problem under considerations. The following assumption is
necessary for reasonings presented in this section.

###### Assumption 1

We assume that problem (1.1) equipped with suitable boundary conditions is
well posed and seek solutions $u[t]:=u(\cdot,t)\in
C^{1}([0,t^{\star}],H^{p}(\Omega))$ while $u_{0}\in H^{p}(\Omega)$. Moreover
we assume that the linear, differential operator
$\mathcal{L}:D(\mathcal{L})\supseteq H^{p}(\Omega)\rightarrow L^{2}(\Omega)$
defined with (1.2) is a generator of $C_{0}$ (strongly continuous) semigroup
$e^{t\mathcal{L}}:L^{2}(\Omega)\rightarrow L^{2}(\Omega)$, that
$f[t]:=f(\cdot,t)\in C^{1}([0,t^{\star}),L^{\infty}(\Omega))$ and that the
functions $a_{\textbf{p}}\in C_{b}^{\infty}(\Omega)$, $|\textbf{p}|\leq p$,
meaning that their all derivatives are bounded on $\Omega$

###### Remark 2.1

Some boundary conditions may enforce additional assumptions on the domain of
the operator $\mathcal{L}$. For example zero boundary conditions would require
$H^{p}(\Omega)\cap H^{1}_{0}(\Omega)$ instead of $H^{p}(\Omega)$ in Assumption
1. Careful changes of the domains of operators and/or coefficients, however,
would not affect the suggested methodology.

Having Assumption 1 satisfied, application of the Duhamel principle (also
known as variation-of-constants) to equation (1.1) results in its equivalent
integral form

$\displaystyle
u[t]={\mathrm{e}}^{t\mathcal{L}}u_{0}+\int_{0}^{t}{\mathrm{e}}^{(t-\tau)\mathcal{L}}f[\tau]u[\tau]\,\mathrm{d}\tau.$
(2.7)

Let us define space $Y:=C([0,t^{\star}];L^{2}(\Omega))$ with standard norm
$\|v\|_{Y}=\max_{0\leq t\leq t^{\star}}\|v[t]\|_{L^{2}(\Omega)}$. Then for
function $f$ and operator $\mathcal{L}$ satisfying Assumption 1 we can define
operator $K(u)[t]:Y\rightarrow Y$ as an integral

$K(u)[t]=\int_{0}^{t}{\mathrm{e}}^{(t-\tau)\mathcal{L}}f[\tau]u[\tau]\,\mathrm{d}\tau,\qquad
0\leq t\leq t^{\star},$ (2.8)

and analyse its multiple compositions

$\displaystyle
K^{d}(u)[t]=\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}f[\tau_{1}]\int_{0}^{\tau_{1}}{\mathrm{e}}^{(\tau_{1}-\tau_{2})\mathcal{L}}f[\tau_{2}]\dots\int_{0}^{\tau_{d-1}}{\mathrm{e}}^{(\tau_{d-1}-\tau_{d})\mathcal{L}}f[\tau_{d}]u[\tau_{d}]\,\mathrm{d}\tau_{d}\dots\,\mathrm{d}\tau_{1}.$
(2.9)

Highly oscillatory nature of function $f$, see definition (1.3), indicates
that these are $d$ times nested highly oscillatory integrals. Obviously
equation (2.7) can be expressed as

$u[t]={\mathrm{e}}^{t\mathcal{L}}u_{0}+Ku[t],\quad t\in[0,t^{\star}],$ (2.10)

and, as often happens, solved in means of the iterative approximations

$u^{[n]}[t]={\mathrm{e}}^{t\mathcal{L}}u_{0}+Ku^{[n-1]}[t],\ {\rm where}\
u^{[0]}={\mathrm{e}}^{t\mathcal{L}}u_{0}.$ (2.11)

We can easily observe that formula (2.11) leads to the well recognised and
widely applied truncated sums of Neumann series, also known as Dyson series or
Time Ordering  operator

$u^{[n]}[t]=\sum_{d=0}^{n}K^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0},\ {\rm where}\
K^{0}={\rm Id}.$ (2.12)

The following theorem is crucial for further investigations as it not only
justifies forthcoming expansions but also clearly establishes error bounds for
approximation in means of Neumann series for the particular form of (1.1),
thus of (2.7) as they are equivalent in our setting.

###### Theorem 2.1

Let us make Assumption 1. Then the elements of an infinite sequence
$\left\\{u^{[n]}\right\\}_{n=0}^{\infty}\in Y$ defined as in (2.12) with the
operator $K$ given with (2.8) converge in norm $\|\,\cdot\,\|_{Y}$ to the
unique solution $u^{\star}$ of (1.1) as $n\rightarrow\infty$. Moreover, the
following estimate holds

$\qquad\|u^{\star}-u^{[n]}\|_{Y}\leq
C_{t^{\star}}\|u_{0}\|_{L^{2}(\Omega)}\exp\Big{(}C_{t^{\star}}\|f\|_{Y}{t^{\star}}\Big{)}\left\|K^{n+1}\right\|_{Y},\\\
$ (2.13)

where $C_{t^{\star}}=\max_{0\leq
t\leq{t^{\star}}}\|e^{t\mathcal{L}}\|_{L^{2}(\Omega)}$.

Proof Let us start with the observation that there exists $C_{t^{\star}}$ such
that $\|{\mathrm{e}}^{t\mathcal{L}}\|_{L^{2}(\Omega)}\leq C_{t^{\star}}$,
$t\in[0,{t^{\star}}]$, and that

$\displaystyle\|K^{n}\|_{Y}$ $\displaystyle\leq$
$\displaystyle\left(C_{t^{\star}}\|f\|_{Y}\right)^{n}\int_{0}^{t^{\star}}\int_{0}^{\tau_{1}}\dots\int_{0}^{\tau_{n-1}}\,\mathrm{d}\tau_{n}\dots\,\mathrm{d}\tau_{1}\leq\frac{(C_{t^{\star}}\|f\|_{Y}{t^{\star}})^{n}}{n!}.$
(2.14)

This means that
$\sum_{d=0}^{\infty}\left\|K^{d}\right\|_{Y}\leq\exp\Big{(}C_{t^{\star}}\|f\|_{Y}{t^{\star}}\Big{)}$
and that $\left\\{u^{[n]}\right\\}_{n=0}^{\infty}\in Y$ forms a Cauchy
sequence with its limit defined as

$\displaystyle
u^{\star}[t]:=\lim_{n\to\infty}u^{[n]}[t]=\sum_{d=0}^{\infty}K^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}.$
(2.15)

Now we can notice that $u^{\star}$ solves (1.1). Indeed, since
$\lim_{n\rightarrow\infty}\|K^{n}\|_{Y}=0$ we obtain

$(I-K)u^{\star}[t]=(I-K)\sum_{d=0}^{\infty}K^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}=\left(I-\lim_{n\rightarrow\infty}K^{n}\right){\mathrm{e}}^{t\mathcal{L}}u_{0}={\mathrm{e}}^{t\mathcal{L}}u_{0},$

which coincides with the fixed point of (2.10) thus solves (2.7) and (1.1).
Observation that

$u^{\star}[t]-u^{[n]}[t]=K^{n+1}\sum_{d=0}^{\infty}K^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$
(2.16)

together with (2.14) results in the estimate (2.13) which completes the poof.
$\Box$

###### Remark 2.2

Estimates (2.13) and (2.14) lead to

$\left\|u^{\star}-u^{[n]}\right\|_{Y}\leq
C_{t^{\star}}\|u_{0}\|_{L^{2}(\Omega)}\exp\Big{(}C_{t^{\star}}\|f\|_{Y}{t^{\star}}\Big{)}\frac{\left(C_{t^{\star}}\|f\|_{Y}{t^{\star}}\right)^{n+1}}{(n+1)!}$
(2.17)

which obviously specifies certain rate of convergence. Given, however that $f$
is highly oscillatory in time, i.e. $f(x,t)={\mathrm{e}}^{{\mathrm{i}}\omega
t}\alpha(x),\ \omega\gg 1$ this estimate is suboptimal. It is well known that
integrals of highly oscillatory functions are of small magnitude. To
illustrate it let us consider the integral of a continuously differentiable
function $g(s)=c(s){\mathrm{e}}^{{\mathrm{i}}\omega s}$, $s\in[0,t]$. One can
easily observe that

$\left|\int_{0}^{t}c(s){\mathrm{e}}^{{\mathrm{i}}\omega
s}\,\mathrm{d}s\right|\leq\bar{C}\min\\{t,\omega^{-1}\\},$

where $\bar{C}$ depends on the coefficient $c(s)$ and on its derivative.

Following this argumentation and observing that estimate (2.13) utilises
multivariate integral of highly oscillatory functions, see (2.9), we can
deduce that under stronger assumptions on our problem the estimate of the rate
of convergence in (2.17) may be improved to

$\left\|u^{\star}-u^{[n]}\right\|\leq C^{n}\cdot\omega^{-n-1}.$

One needs to keep in mind, however that the constant is independent of
$\omega$ but might get dangerously large (even larger then $\omega^{-n-1}$ !)
as it depends on other few factors, namely
$C^{n}(t^{\star},\mathcal{L},u_{0},\alpha)$. For this reason we are interested
in the formula describing the error of truncation.

##  3 Multivariate oscillatory integral as an asymptotic expansion

In this section we are concerned with the operator (2.8) where
$u={\mathrm{e}}^{t\mathcal{L}}u_{0}$ and
$f(x,t)=\alpha(x){\mathrm{e}}^{i\omega t}$, $\omega\gg 1$. To distinguish this
form we will use letter $T$ instead of $K$:

$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}=\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\int_{0}^{\tau_{1}}{\mathrm{e}}^{(\tau_{1}-\tau_{2})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{2}}\dots\int_{0}^{\tau_{d-1}}{\mathrm{e}}^{(\tau_{d-1}-\tau_{d})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega
t_{d}}{\mathrm{e}}^{\tau_{d}\mathcal{L}}u_{0}\,\mathrm{d}\tau_{d}\dots\,\mathrm{d}\tau_{1}$
(3.18)

We are interested in an asymptotic expansion of integrals
$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0},d=1,\dots n$. More specifically for
given $d$ we will present $T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$ as a finite
series and an error term related to the asymptotic expansion. This will lead
us towards the expansions of $n$-th approximations of (2.7). Indeed, it is
enough to observe that $u^{[n]}$ is defined, see (2.12), as an infinite sum of
integrals $T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$. Let start with necessary
assumption and auxiliary lemmas to which we will refer to in the further part
of this manuscript.

###### Assumption 2

We assume that problem (1.1) equipped with suitable boundary conditions is
well posed and seek solutions $u[t]:=u(\cdot,t)\in
C^{n}([0,t^{\star}],H^{np}(\Omega))$ while $u_{0}\in H^{np}(\Omega)$, where
$n\in\mathbb{N}$ is a fixed number. Moreover we assume that the linear,
differential operator $\mathcal{L}:D(\mathcal{L})\supseteq
H^{np}(\Omega)\rightarrow L^{2}(\Omega)$ defined with (1.2) is a generator of
a differential (in strong operator topology) semigroup
$e^{t\mathcal{L}}:L^{2}(\Omega)\rightarrow L^{2}(\Omega)$, that $\alpha\in
W^{np,\infty}(\Omega)$

and that the functions $a_{\textbf{p}}\in C_{b}^{\infty}(\Omega)$,
$|\textbf{p}|\leq p$, meaning that their all derivatives are bounded on
$\Omega$

###### Lemma 3.1

Let us make Assumption 2. The $k$-th derivative of expression
${\mathrm{e}}^{(t-\tau)\mathcal{L}}\alpha{\mathrm{e}}^{\tau\mathcal{L}}$ can
be presented in the following way

$\displaystyle\partial_{\tau}^{k}\left({\mathrm{e}}^{(t-\tau)\mathcal{L}}\alpha{\mathrm{e}}^{\tau\mathcal{L}}\right)={\mathrm{e}}^{(t-\tau)\mathcal{L}}\underbrace{\Big{[}\dots\big{[}[\alpha,\mathcal{L}],\mathcal{L}\big{]},\dots\Big{]}}_{\mathcal{L}\
appears\ k\
\text{times}}{\mathrm{e}}^{\tau\mathcal{L}}=(-1)^{k}{\mathrm{e}}^{(t-\tau)\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{k}\big{(}\alpha\big{)}{\mathrm{e}}^{\tau\mathcal{L}},\quad
k\leq n,$

where
$\mathrm{ad}\,_{\mathcal{L}}^{0}\big{(}\alpha\big{)}=\alpha,\quad\mathrm{ad}\,_{\mathcal{L}}^{k}\big{(}\alpha\big{)}=\big{[}\mathcal{L},\mathrm{ad}\,_{\mathcal{L}}^{k-1}\big{(}\alpha\big{)}\big{]}$
and $[X,Y]\equiv XY-YX$ is the commutator of $X$ and $Y$.

###### Lemma 3.2

Let us make Assumption 2. The following integral can be expressed as a sum of
$n+1$ components scaling like $\omega^{-m},\quad m=1,\ldots n.$

$\displaystyle\int_{0}^{t}{\mathrm{e}}^{{\mathrm{i}}\omega\tau}{\mathrm{e}}^{(t-\tau)\mathcal{L}}\alpha{\mathrm{e}}^{\tau\mathcal{L}}u_{0}\,\mathrm{d}\tau=\sum_{m=0}^{n-1}\frac{1}{({\mathrm{i}}\omega)^{m+1}}\Big{[}{\mathrm{e}}^{{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{m}\big{(}\alpha\big{)}{\mathrm{e}}^{t\mathcal{L}}-{\mathrm{e}}^{t\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{m}\big{(}\alpha\big{)}\Big{]}u_{0}$
(3.19) $\displaystyle\hskip
108.12054pt+\frac{1}{({\mathrm{i}}\omega)^{n}}\int_{0}^{t}{\mathrm{e}}^{{\mathrm{i}}\omega\tau}{\mathrm{e}}^{(t-\tau)\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{n}\big{(}\alpha\big{)}{\mathrm{e}}^{\tau\mathcal{L}}u_{0}\,\mathrm{d}\tau$

Lemma 3.1 can be easily proved by induction while Lemma 3.2 by $n$-timeas
integration by parts. Let us observe that in case $u_{0}\in
H^{(n+1)p}(\Omega)$ and $\alpha\in W^{(n+1)p,\infty}(\Omega)$, the last
component of expansion (3.19) scales like $\omega^{-(n+1)}$. Indeed, the last
term can be integrated by parts again.

To derive the formula for the asymptotic series of the multivariate highly
oscillatory integral $T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$, we will heavily
utilise Lemma 3.2. We will start with the asymptotic expansion of the
innermost integral, then we will expand asymptotically each successive
integrals up to the outermost. Integration by parts of each integral will
result in appearance of the two terms
$\sum_{m=0}^{n-1}\frac{1}{({\mathrm{i}}\omega)^{m+1}}{\mathrm{e}}^{{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{m}\big{(}\alpha\big{)}{\mathrm{e}}^{t\mathcal{L}}u_{0}$
and
$\sum_{m=0}^{n-1}\frac{-1}{({\mathrm{i}}\omega)^{m+1}}{\mathrm{e}}^{t\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{m}\big{(}\alpha\big{)}u_{0}$.
As a consequence, the number of components in the series of expansion of
$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$ will grow exponentially in $d$,
therefore to facilitate the notation we introduce the following definitions.

###### Definition 3.1

Let us define multi-index ${\bf k^{s}}=(k_{1},\ldots,k_{s})\in\mathbb{N}^{s}$,
and understand that $|{\bf k^{s}}|=|{\bf k^{s}}|_{1}=k_{1}+\ldots+k_{s}$.

###### Definition 3.2

Let us make Assumption 2 and define operator
$\mathcal{N}_{k_{l}}^{k_{r}}[\mathcal{L},\alpha]$ as follows

$\mathcal{N}_{k_{l}}^{k_{r}}[\mathcal{L},\alpha]=\left\\{\begin{array}[]{ll}\mathrm{ad}\,_{\mathcal{L}}^{k_{r}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{r-1}}(\dots\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{l}}(\alpha))\big{)},&l<r\\\
\\\ \mathrm{ad}\,_{\mathcal{L}}^{k_{l}}(\alpha),&l=r\\\ \\\ {\rm Id},&l>r\\\
\end{array}\right.,$

where $k_{l},k_{l+1},\dots,k_{r-1},k_{r}$, $l\leq r$ are inputs of multi-index
${\bf k^{s}}$, and Id is an identity function on space $W$. For clarity of
exposition we will alternate notation
$\mathcal{N}_{k_{l}}^{k_{r}}[\mathcal{L},\alpha]$ with
$\mathcal{N}_{k_{l}}^{k_{r}}$ in more complicated expressions.

###### Remark 3.1

$\mathcal{N}_{k_{l}}^{k_{r}}$ can be also defined by the recursive formula

$\mathcal{N}_{k_{l}}^{k_{l}}=\mathrm{ad}\,_{\mathcal{L}}^{k_{l}}\big{(}\alpha\big{)},\quad\mathcal{N}_{k_{l}}^{k_{r}}=\mathrm{ad}\,_{\mathcal{L}}^{k_{r}}\big{(}\alpha\mathcal{N}_{k_{l}}^{k_{r-1}}\big{)},\quad\text{for}\quad
r>l\quad\text{and}\quad\mathcal{N}_{k_{l}}^{k_{r}}={\rm Id},\quad l>r.$

###### Definition 3.3

Let us make Assumption 2. Having operator
$\mathcal{N}_{k_{l}}^{k_{r}}[\mathcal{L},\alpha]$ introduced in Definition 3.2
we define function $\mathcal{F}_{\ell}[\mathcal{L},\alpha]$, independent of
$t$, which satisfy the following recursive relationship

$\displaystyle\mathcal{F}_{0}[\mathcal{L},\alpha]={\rm Id},$
$\displaystyle\mathcal{F}_{\ell}[\mathcal{L},\alpha]=-\sum_{m=0}^{\ell-1}\Big{(}a_{\ell-m}^{\ell}\mathcal{N}_{k_{m+1}}^{k_{\ell}}[\mathcal{L},\alpha]\mathcal{F}_{m}[\mathcal{L},\alpha]\Big{)},$
(3.20)

where the coefficients $a_{j}^{\ell}\in\mathbb{Q}$, $j=0,1,2,\ldots,\ell$ are
numbers satisfy:

$a_{j}^{\ell}=\frac{1}{\prod_{r=0}^{j-1}(j-r)^{k_{\ell-r}+1}},\quad
a_{0}^{\ell}=1.$ (3.21)

As it was in the case of $\mathcal{N}_{k_{l}}^{k_{r}}[\mathcal{L},\alpha]$ we
will write $\mathcal{F}_{\ell}$ instead of
$\mathcal{F}_{\ell}[\mathcal{L},\alpha]$ if only it simplifies notation and
does not lead to misunderstandings.

###### Example 3.1

To illustrate the Definition 3.3 we present the exact formulas for
$\mathcal{F}_{j},j=1,2,3$.

$\displaystyle\mathcal{F}_{1}[\mathcal{L},\alpha]$ $\displaystyle=$
$\displaystyle-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}\big{(}\alpha\big{)}$
$\displaystyle\mathcal{F}_{2}[\mathcal{L},\alpha]$ $\displaystyle=$
$\displaystyle-\frac{1}{2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}+\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)$
$\displaystyle\mathcal{F}_{3}[\mathcal{L},\alpha]$ $\displaystyle=$
$\displaystyle\frac{-1}{3^{k_{3}+1}\cdot
2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha))\big{)}+\frac{1}{2^{k_{3}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\big{)}\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)$
$\displaystyle+\frac{1}{2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}-\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)$

Subsequent corollary, which can be easily proved by induction, will be
referred to in the proof of the main theorem while calculating nested
integrals with leading frequencies $\frac{1}{({\mathrm{i}}\omega)^{k_{1}+1}}$
and $\frac{1}{({\mathrm{i}}\omega)^{k_{2}+1}}$, where $k_{1},k_{2}$ will be
inputs of multi-index ${\textbf{k}^{2}}$

###### Corollary 3.1

The following identity holds

$\displaystyle\sum_{k_{1}=0}^{n}\frac{1}{({\mathrm{i}}\omega)^{k_{1}+1}}a_{k_{1}}\sum_{k_{2}=0}^{n-k_{1}-1}\frac{1}{({\mathrm{i}}\omega)^{k_{2}+1}}b_{k_{2}}=\sum_{k_{1}+k_{2}=0}^{n-1}\frac{1}{({\mathrm{i}}\omega)^{k_{1}+k_{2}+2}}a_{k_{1}}b_{k_{2}}$
(3.22)

with assumption that $\sum_{k=0}^{-1}a_{k}=0$.

###### Theorem 3.1

Let us make Assumption 2. The highly oscillatory integral
$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$ defined with (3.18) can be presented
as a sum of finite series $\mathcal{S}^{d}_{n}[\mathcal{L},\alpha](t)$ and the
remaining term $E^{d}_{n}[\mathcal{L},\alpha](t)$

$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}=\mathcal{S}^{d}_{n}[\mathcal{L},\alpha](t)+E^{d}_{n}[\mathcal{L},\alpha](t),\quad\
n\geq d,\quad{\rm where}$ (3.23)

$\displaystyle\mathcal{S}^{d}_{n}[\mathcal{L},\alpha](t)$
$\displaystyle=\sum_{|{\bf
k^{d}}|=0}^{n-d}\frac{1}{({\mathrm{i}}\omega)^{d+|{\mathbf{k}^{d}}|}}\sum_{\ell=0}^{d}a_{d-\ell}^{d}{\mathrm{e}}^{(d-\ell){\mathrm{i}}\omega
t}\mathcal{N}_{k_{\ell+1}}^{k_{d}}{\mathrm{e}}^{t\mathcal{L}}\mathcal{F}_{\ell}u_{0},$
(3.24) $\displaystyle E^{d}_{n}[\mathcal{L},\alpha](t)$
$\displaystyle=\frac{1}{({\mathrm{i}}\omega)^{n}}\sum_{|{\mathbf{k}^{d}}|=n-d+1}\
\
\sum_{\ell=0}^{d-1}b_{d-\ell}^{d}\int_{0}^{t}{\mathrm{e}}^{{(t-\tau_{1})\mathcal{L}}}\mathcal{N}_{k_{\ell+1}}^{k_{d}}{\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{(d-\ell){\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}\mathcal{F}_{\ell}u_{0}$
(3.25)
$\displaystyle+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}E^{d-1}_{n}[\mathcal{L},\alpha](\tau_{1})u_{0}\,\mathrm{d}\tau_{1},$

and $E_{n}^{0}[\mathcal{L},\alpha](t)\equiv 0$ while
$E_{n}^{d-1}[\mathcal{L},\alpha](t)$ is the remaining term of the analogous
integral $T^{d-1}{\mathrm{e}}^{t\mathcal{L}}u_{0}$ and coefficients
$b_{j}^{\ell}$, $j=0,1,2,\ldots,\ell$ satisfy:

$b_{j}^{\ell}=\frac{1}{j^{k_{\ell}}\cdot\prod_{r=1}^{j-1}(j-r)^{k_{\ell-r}+1}},\quad
b_{0}^{\ell}=1,$ (3.26)

with $\prod_{r=1}^{0}=1$. Moreover

$\|T^{d}\|_{Y}\leq\omega^{-d}C_{T}^{d}(t^{\star},\mathcal{L},u_{0},\alpha),$
(3.27)

were $C_{T}^{d}(t^{\star},\mathcal{L},u_{0},\alpha$) is a constant depending
on $d$, $t^{\star}$, operator $\mathcal{L}$, initial condition $u_{0}$ and
function $\alpha$.

###### Remark 3.2

To shorten notation will write $S_{n}^{d}(t)$ and $E_{n}^{d}(t)$ instead of
$S_{n}^{d}[\mathcal{L},\alpha](t)$ and $E_{n}^{d}[\mathcal{L},\alpha](t)$
respectively whenever it does not lead to misunderstanding. We will write also
$T^{d}{\mathrm{e}}^{t\mathcal{L}}$ instead of
$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$.

Proof of Theorem 3.1

We will prove (3.23) by induction on $d$. Thus taking $d=1$ we consider the
single integral
$T^{1}{\mathrm{e}}^{t\mathcal{L}}=\int_{0}^{t}{\mathrm{e}}^{{\mathrm{i}}\omega\tau}{\mathrm{e}}^{(t-\tau)\mathcal{L}}\alpha{\mathrm{e}}^{\tau\mathcal{L}}\,\mathrm{d}\tau$.
In this simplest case multi-index ${\mathbf{k}^{1}}\in\mathbb{N}^{1}$ consists
of one input $k_{1}$ only, and $a_{1}^{1}=b_{1}^{1}=1$ while
$\mathcal{F}_{1}=-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}\big{(}\alpha\big{)}$.
With these inputs the considered integral is given with the following formula

$\displaystyle
T^{1}{\mathrm{e}}^{t\mathcal{L}}=\sum_{k_{1}=0}^{n-1}\frac{1}{({\mathrm{i}}\omega)^{1+k_{1}}}\left({\mathrm{e}}^{{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha){\mathrm{e}}^{t\mathcal{L}}+{\mathrm{e}}^{t\mathcal{L}}\left(-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\right)\right)+\frac{1}{({\mathrm{i}}\omega)^{n}}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{n}(\alpha){\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1},$

which fully coincides with (3.19) from Lemma 3.2 and proves the first step of
induction. Suppose now that the formula is valid for
$T^{d}{\mathrm{e}}^{t\mathcal{L}}$. Then

$\displaystyle T^{d+1}{\mathrm{e}}^{t\mathcal{L}}=$
$\displaystyle\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}T^{d}{\mathrm{e}}^{\tau_{1}\mathcal{L}}\,\mathrm{d}\tau_{1}\stackrel{{\scriptstyle(1)}}{{=}}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\Big{(}\mathcal{S}^{d}_{n}(\tau_{1})+E^{d}_{n}(\tau_{1})\Big{)}\,\mathrm{d}\tau_{1}$
$\displaystyle\stackrel{{\scriptstyle(2)}}{{=}}$ $\displaystyle\sum_{|{\bf
k^{d}}|=0}^{n-d}\frac{1}{({\mathrm{i}}\omega)^{d+|{\bf
k^{d}}|}}\sum_{\ell=0}^{d}a_{d-\ell}^{d}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{(d-\ell+1){\mathrm{i}}\omega\tau_{1}}\mathcal{N}_{k_{\ell+1}}^{k_{d}}{\mathrm{e}}^{\tau_{1}\mathcal{L}}\mathcal{F}_{\ell}\,\mathrm{d}\tau_{1}+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}E^{d}_{n}(\tau_{1})\,\mathrm{d}\tau_{1}$
$\displaystyle\stackrel{{\scriptstyle(3)}}{{=}}$ $\displaystyle\sum_{|{\bf
k^{d}}|=0}^{n-d}\frac{1}{({\mathrm{i}}\omega)^{d+|{\bf
k^{d}}|}}\sum_{k_{d+1}=0}^{n-d-|{\bf
k^{d}}|-1}\frac{1}{({\mathrm{i}}\omega)^{k_{d+1}+1}}\sum_{\ell=0}^{d}a_{d-\ell}^{d}\frac{1}{(d-\ell+1)^{k_{d+1}+1}}\cdot$
$\displaystyle\hskip
113.81102pt\cdot\Bigg{(}{\mathrm{e}}^{(d-\ell+1){\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{d+1}}\big{(}\alpha\mathcal{N}_{k_{\ell+1}}^{k_{d}}\big{)}{\mathrm{e}}^{t\mathcal{L}}\mathcal{F}_{\ell}-{\mathrm{e}}^{t\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{k_{d+1}}\big{(}\alpha\mathcal{N}_{k_{\ell+1}}^{k_{d}}\big{)}\mathcal{F}_{\ell}\Bigg{)}$
$\displaystyle+$ $\displaystyle\sum_{|{\bf
k^{d}}|=0}^{n-d}\frac{1}{({\mathrm{i}}\omega)^{d+|{\bf
k^{d}}|}}\frac{1}{({\mathrm{i}}\omega)^{n-d-|{\bf
k^{d}}|}}\sum_{\ell=0}^{d}a_{d-\ell}^{d}\frac{1}{(d-\ell+1)^{n-d-|{\bf
k^{d}}|}}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{n-d-|{\bf
k^{d}}|}\big{(}\alpha\mathcal{N}_{k_{\ell+1}}^{k_{d}}\big{)}{\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{(d-\ell+1){\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}\mathcal{F}_{\ell}$
$\displaystyle\quad+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}E^{d}_{n}(\tau_{1})\,\mathrm{d}\tau_{1}$
$\displaystyle\stackrel{{\scriptstyle(4)}}{{=}}$ $\displaystyle\sum_{|{\bf
k^{d+1}}|=0}^{n-d-1}\frac{1}{({\mathrm{i}}\omega)^{d+1+|{\bf
k^{d+1}}|}}\sum_{\ell=0}^{d}a_{d+1-\ell}^{d+1}\Bigg{(}{\mathrm{e}}^{(d-\ell+1){\mathrm{i}}\omega
t}\mathcal{N}_{k_{\ell+1}}^{k_{d+1}}{\mathrm{e}}^{t\mathcal{L}}\mathcal{F}_{\ell}-{\mathrm{e}}^{t\mathcal{L}}\mathcal{N}_{k_{\ell+1}}^{k_{d+1}}\mathcal{F}_{\ell}\Bigg{)}$
$\displaystyle\quad+\frac{1}{({\mathrm{i}}\omega)^{n}}\sum_{|{\bf
k^{d+1}}|=n-d}\sum_{\ell=0}^{d}b_{d+1-\ell}^{d+1}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathcal{N}_{k_{\ell+1}}^{k_{d+1}}{\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{(d-\ell+1){\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}\mathcal{F}_{\ell}+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}E^{d}_{n}(\tau_{1})\,\mathrm{d}\tau_{1}$
$\displaystyle\stackrel{{\scriptstyle(5)}}{{=}}$ $\displaystyle\sum_{|{\bf
k^{d+1}}|=0}^{n-d-1}\frac{1}{({\mathrm{i}}\omega)^{d+1+|{\bf
k^{d+1}}|}}\sum_{\ell=0}^{d+1}a_{d+1-\ell}^{d+1}{\mathrm{e}}^{(d-\ell+1){\mathrm{i}}\omega
t}\mathcal{N}_{k_{\ell+1}}^{k_{d+1}}{\mathrm{e}}^{t\mathcal{L}}\mathcal{F}_{\ell}$
$\displaystyle\quad+\frac{1}{({\mathrm{i}}\omega)^{n}}\sum_{|{\bf
k^{d+1}}|=n-d}\sum_{\ell=0}^{d}b_{d+1-\ell}^{d+1}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathcal{N}_{k_{\ell+1}}^{k_{d+1}}{\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{(d-\ell+1){\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}\mathcal{F}_{\ell}+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha
E^{d}_{n}(\tau_{1})\,\mathrm{d}\tau_{1}$
$\displaystyle\stackrel{{\scriptstyle(6)}}{{=}}$
$\displaystyle\mathcal{S}^{d+1}_{n}(t)+E^{d+1}_{n}(t)$

In (1) we exploit inductive assumption, in (2) formula for
$\mathcal{S}^{d}_{n}(t)$, in (3) we refer to Lemma 3.2 and perform integrating
by parts $n-d-|{\bf k_{d}}|$ times. Equivalence (4) is valid due to Corollary
3.1, Definition 3.2 and formulas (3.21), (3.26) where coefficients
$a_{j}^{\ell},\ b_{j}^{\ell}$ are defined. Definition 3.3 is utilised in (5),
more specifically we observe that
$\mathcal{F}_{d+1}=-\sum_{\ell=0}^{d}a_{d+1-\ell}^{d+1}\mathcal{N}_{\ell+1}^{d+1}\mathcal{F}_{\ell}$.
Expressions preceding equivalence (6) coincide with (3.23) which proves the
first thesis of the theorem.

To prove (3.27) we start with the observation that there exists such a
constant $\bar{C}_{n}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)$ that
$\|\mathcal{S}^{d}_{n}[\mathcal{L},\alpha]\|_{Y}\leq\omega^{-d}\bar{C}_{n}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)$.
By induction on $n$ we can easily prove that there exists also such a constant
$\tilde{C}_{n}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)$ that
$\|E^{d}_{n}[\mathcal{L},\alpha]\|_{Y}\leq\omega^{-n}\tilde{C}_{n}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)$.
Since $n\geq d$ we obtain

$\|T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}\|_{Y}=\|\mathcal{S}^{d}_{n}[\mathcal{L},\alpha](t)+E^{d}_{n}[\mathcal{L},\alpha](t)\|_{Y}\leq\omega^{-d}\Big{(}\bar{C}_{n}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)+\tilde{C}_{n}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)\Big{)},$

which proves (3.27) and completes the proof. $\Box$

###### Remark 3.3

In a more general settings, but without specifying the terms of asymptotic
series and error related to the truncation, multivariate highly oscillatory
integrals were investigated for example in [Iserles and Norsett(2006)],
[Olver(2008)]. Authors of these manuscripts shown that

$\int_{\Omega}f(\textbf{y}){\mathrm{e}}^{{\mathrm{i}}\omega
g(\textbf{y})}\,\mathrm{d}V\sim\mathcal{O}(\omega^{-d}),$

for $\Omega\subset\mathbb{R}^{d}$ having a piecewise smooth boundary where
$\nabla g\notin\bar{\Omega}$ and $\nabla g\not\perp\partial{\Omega}$. Also the
non-resonance condition is assumed.

###### Example 3.2

To illustrate Theorem 3.1 we will present the first three integrals, namely
$T^{1}{\mathrm{e}}^{t\mathcal{L}}$, $T^{2}{\mathrm{e}}^{t\mathcal{L}}$ and
$T^{3}{\mathrm{e}}^{t\mathcal{L}}$.

$\displaystyle T^{1}{\mathrm{e}}^{t\mathcal{L}}$ $\displaystyle=$
$\displaystyle\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}=\mathcal{S}^{1}_{n}(t)+E^{1}_{n}(t)$
$\displaystyle\mathcal{S}^{1}_{n}(t)$ $\displaystyle=$
$\displaystyle\sum_{k_{1}=0}^{n-1}\frac{1}{({\mathrm{i}}\omega)^{k_{1}+1}}\Big{[}{\mathrm{e}}^{{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}\big{(}\alpha\big{)}{\mathrm{e}}^{t\mathcal{L}}+{\mathrm{e}}^{t\mathcal{L}}\big{(}\underbrace{-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}\big{(}\alpha\big{)}}_{\mathcal{F}_{1}}\big{)}\Big{]}$
$\displaystyle E^{1}_{n}(t)$ $\displaystyle=$
$\displaystyle\frac{1}{({\mathrm{i}}\omega)^{n}}\int_{0}^{t}{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{n}\big{(}\alpha\big{)}{\mathrm{e}}^{\tau_{1}\mathcal{L}}\,\mathrm{d}t\tau_{1}$

$\displaystyle T^{2}{\mathrm{e}}^{t\mathcal{L}}$ $\displaystyle=$
$\displaystyle\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\int_{0}^{\tau_{1}}{\mathrm{e}}^{(\tau_{1}-\tau_{2})\mathcal{L}}\alpha{\mathrm{e}}^{\tau_{2}\mathcal{L}}{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{2}}\,\mathrm{d}\tau_{2}\,\mathrm{d}\tau_{1}=\mathcal{S}^{2}_{n}(t)+E^{2}_{n}(t)$
$\displaystyle\mathcal{S}^{2}_{n}(t)$ $\displaystyle=$
$\displaystyle\sum_{|{\bf k^{2}}|=0}^{n-2}\frac{1}{({\mathrm{i}}\omega)^{|{\bf
k^{2}}|+2}}\Bigg{[}\frac{1}{2^{k_{2}+1}}{\mathrm{e}}^{2{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}{\mathrm{e}}^{t\mathcal{L}}+{\mathrm{e}}^{{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha){\mathrm{e}}^{t\mathcal{L}}\big{(}\underbrace{-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)}_{\mathcal{F}_{1}}\big{)}$
$\displaystyle\hskip
68.28644pt+{\mathrm{e}}^{t\mathcal{L}}\Big{(}\underbrace{-\frac{1}{2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}+\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)}_{\mathcal{F}_{2}}\Big{)}\Bigg{]}$
$\displaystyle E^{2}_{n}(t)$ $\displaystyle=$
$\displaystyle\frac{1}{({\mathrm{i}}\omega)^{n}}\sum_{|{\bf
k^{2}}|=n-1}\Bigg{[}\frac{1}{2^{k_{2}}}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}{\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{2{\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}$
$\displaystyle\hskip
68.28644pt+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha){\mathrm{e}}^{\tau_{1}\mathcal{L}}\Big{(}-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\Big{)}{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}\Bigg{]}$
$\displaystyle+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}E^{1}_{n}(\tau_{1})\,\mathrm{d}\tau_{1}$

$\displaystyle T^{3}{\mathrm{e}}^{t\mathcal{L}}$ $\displaystyle=$
$\displaystyle\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\int_{0}^{\tau_{1}}{\mathrm{e}}^{(\tau_{1}-\tau_{2})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{2}}\int_{0}^{\tau_{2}}{\mathrm{e}}^{(\tau_{2}-\tau_{3})\mathcal{L}}\alpha{\mathrm{e}}^{\tau_{3}\mathcal{L}}{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{3}}\,\mathrm{d}\tau_{3}\,\mathrm{d}\tau_{2}\,\mathrm{d}\tau_{1}=\mathcal{S}^{3}_{n}(t)+E^{3}_{n}(t)$
$\displaystyle\mathcal{S}^{3}_{n}(t)$ $\displaystyle=$
$\displaystyle\sum_{|{\bf k^{3}}|=0}^{n-3}\frac{1}{({\mathrm{i}}\omega)^{|{\bf
k^{3}}|+3}}\Bigg{(}\frac{1}{3^{k_{3}+1}\cdot
2^{k_{2}+1}}{\mathrm{e}}^{3{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha))\big{)}{\mathrm{e}}^{t\mathcal{L}}+$
$\displaystyle\hskip
28.45274pt+\frac{1}{2^{k_{3}+1}}{\mathrm{e}}^{2{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\big{)}{\mathrm{e}}^{t\mathcal{L}}\Big{(}\underbrace{-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)}_{\mathcal{F}_{1}}\Big{)}$
$\displaystyle\hskip 28.45274pt+{\mathrm{e}}^{{\mathrm{i}}\omega
t}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}(\alpha){\mathrm{e}}^{t\mathcal{L}}\Big{(}\underbrace{-\frac{1}{2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}+\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha(x))\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)}_{\mathcal{F}_{2}}\Big{)}$
$\displaystyle\hskip
28.45274pt+{\mathrm{e}}^{t\mathcal{L}}\Big{(}\underbrace{\frac{-1}{3^{k_{3}+1}\cdot
2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha))\big{)}+\frac{1}{2^{k_{3}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\big{)}\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)}$
$\displaystyle\hskip
28.45274pt+\underbrace{\frac{1}{2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}-\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)}_{\mathcal{F}_{3}}\Big{)}\Bigg{)}+$
$\displaystyle E^{3}_{n}(t)$ $\displaystyle=$
$\displaystyle\frac{1}{({\mathrm{i}}\omega)^{n}}\sum_{|{\bf
k^{3}}|=n-2}\Bigg{(}\frac{1}{3^{k_{3}}}\frac{1}{2^{k_{2}+1}}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha))\big{)}{\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{3{\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}$
$\displaystyle+\frac{1}{2^{k_{3}}}\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\big{)}{\mathrm{e}}^{\tau_{1}\mathcal{L}}\Big{(}-\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\Big{)}{\mathrm{e}}^{2{\mathrm{i}}\omega\tau_{1}}\,\mathrm{d}\tau_{1}$
$\displaystyle+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{k_{3}}(\alpha){\mathrm{e}}^{\tau_{1}\mathcal{L}}{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}\Big{(}-\frac{1}{2^{k_{2}+1}}\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}\big{(}\alpha\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\big{)}+\mathrm{ad}\,_{\mathcal{L}}^{k_{2}}(\alpha)\mathrm{ad}\,_{\mathcal{L}}^{k_{1}}(\alpha)\Big{)}\,\mathrm{d}\tau_{1}\Bigg{)}$
$\displaystyle+\int_{0}^{t}{\mathrm{e}}^{(t-\tau_{1})\mathcal{L}}\alpha{\mathrm{e}}^{{\mathrm{i}}\omega\tau_{1}}E^{2}_{n}(\tau_{1})\,\mathrm{d}\tau_{1}$

##  4 The solution as a modulated Fourier expansion

In this Section we will formulate the main result on the form of the
asymptotic expansion for problem (1.1) and on the accuracy of its truncation.
Assumption 2 is sufficient for the derivation of the coefficients of magnitude
$\omega^{-r}$, $r=0,\ldots,n$ and for the proof that the error $R_{n}(t)$
scales like $\omega^{-n}$. To obtain error $R_{n}(t)$ scaling like
$\omega^{-(n+1)}$ a stronger regularity of the problem is required, because
one more integration by parts will be necessary, so let us start with the
formulation of the assumption.

###### Assumption 3

We assume that problem (1.1) equipped with suitable boundary conditions is
well posed and seek solutions $u[t]:=u(\cdot,t)\in
C^{n+1}([0,T],H^{(n+1)p}(\Omega))$ while $u_{0}\in H^{(n+1)p}(\Omega)$, where
$n\in\mathbb{N}$ is a fixed number. Moreover we assume that the linear,
differential operator $\mathcal{L}:D(\mathcal{L})\supseteq
H^{(n+1)p}(\Omega)\rightarrow L^{2}(\Omega)$ defined with (1.2) is a generator
of a differential (in strong operator topology) semigroup
$e^{t\mathcal{L}}:L^{2}(\Omega)\rightarrow L^{2}(\Omega)$, that $\alpha\in
W^{(n+1)p,\infty}(\Omega)$ and that the functions $a_{\textbf{p}}\in
C_{b}^{\infty}(\Omega)$, $|\textbf{p}|\leq p$, meaning that their all
derivatives are bounded on $\Omega$

###### Theorem 4.1

Let us make Assumption 3. The solution of (1.1) can be presented as a
modulated Fourier series of form (1.6) where

$p_{0,0}(t)={\mathrm{e}}^{t\mathcal{L}}u_{0},\quad
p_{r,s}(t)=\sum_{d=s}^{r}\sum_{|{\bf
k^{d}}|=r-d}\frac{a_{s}^{d}}{({\mathrm{i}})^{r}}\mathcal{N}_{k_{d-s+1}}^{k_{d}}{\mathrm{e}}^{t\mathcal{L}}\mathcal{F}_{d-s}u_{0},\quad
r\neq 0,\quad{\rm with}\quad\sum_{\mathbf{|k^{0}|}}=0,$ (4.28)

and

$R_{n}(t)=\sum_{d=1}^{n}E^{d}_{n}[\mathcal{L},\alpha](t)+\sum_{d=n+1}^{\infty}\Big{(}\mathcal{S}^{d}_{d}[\mathcal{L},\alpha](t)+E^{d}_{d}[\mathcal{L},\alpha](t)\Big{)},$
(4.29)

with $S_{n}^{d}[\mathcal{L},\alpha](t)$ and
$E_{n}^{d}[\mathcal{L},\alpha](t)$, $n\geq d$, defined as (3.24) and (3.25),
respectively. Moreover $R_{n}(t)$ satisfies the following relation

$\|R_{n}\|_{Y}\leq\omega^{-(n+1)}C_{R}^{n}(t^{\star},\mathcal{L},u_{0},\alpha),$
(4.30)

were $C_{R}^{n}(t^{\star},\mathcal{L},u_{0},\alpha$) is a constant depending
on $n$, $t^{\star}$, operator $\mathcal{L}$, initial condition $u_{0}$ and
function $\alpha$.

Proof Based on the Theorem 2.1 and its proof, see equations (2.15) and (2.16),
we conclude that for $f(x,t)=\alpha(x){\mathrm{e}}^{i\omega t}$ solution
$u^{\star}$ of (1.1) satisfies

$u^{\star}[t]=\sum_{d=0}^{n}T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}+T^{n+1}\sum_{d=0}^{\infty}T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0},$
(4.31)

where $T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}$ is defined with (3.18). For the
completeness of exposition let us rewrite the thesis of Theorem 3.1

$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}=\mathcal{S}_{n}^{d}[\mathcal{L},\alpha](t)+E_{n}^{d}[\mathcal{L},\alpha](t),$
(4.32)

and sort the summands of $\mathcal{S}_{n}^{d}[\mathcal{L},\alpha](t)$
according to their magnitude $\frac{1}{\omega^{r}}$ and frequency
${\mathrm{e}}^{s{\mathrm{i}}\omega t}$ as follows

$\mathcal{S}_{n}^{d}[\mathcal{L},\alpha](t)=\frac{1}{\omega^{d}}\sum_{s=0}^{d}S_{d,s}^{d}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}u_{0}+\frac{1}{\omega^{d+1}}\sum_{s=0}^{d}S_{d+1,s}^{d}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}u_{0}+\ldots+\frac{1}{\omega^{n}}\sum_{s=0}^{d}S_{n,s}^{d}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}u_{0},\quad d=1,\dots n.$ (4.33)

Based on Theorem 3.1 we can conclude that expression $S_{r,s}^{d}$ is of the
following form

$S_{r,s}^{d}(t)=\sum_{\mathbf{|k^{d}|}=r-d}\frac{a_{s}^{d}}{({\mathrm{i}})^{r}}\mathcal{N}_{k_{d-s+1}}^{k_{d}}{\mathrm{e}}^{t\mathcal{L}}\mathcal{F}_{d-s}.$
(4.34)

Let us write the first $n$ oscillatory elements of expansion (4.31) in terms
of (4.32) and (4.34) in a form of a table

$\begin{array}[]{ccccccc}T^{1}{\mathrm{e}}^{t\mathcal{L}}=&\frac{1}{\omega}\sum_{s=0}^{1}S_{1,s}^{1}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}&+\frac{1}{\omega^{2}}\sum_{s=0}^{1}S_{2,s}^{1}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}&+\frac{1}{\omega^{3}}\sum_{s=0}^{1}S_{3,s}^{1}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}+&\ldots&+\frac{1}{\omega^{n}}\sum_{s=0}^{1}S_{n,s}^{1}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}&+E^{1}_{n}(t)\\\
T^{2}{\mathrm{e}}^{t\mathcal{L}}=&&\quad\frac{1}{\omega^{2}}\sum_{s=0}^{2}S_{2,s}^{2}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}&+\frac{1}{\omega^{3}}\sum_{s=0}^{2}S_{3,s}^{2}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}+&\ldots&+\frac{1}{\omega^{n}}\sum_{s=0}^{2}S_{n,s}^{2}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}&+E^{2}_{n}(t)\\\
T^{3}{\mathrm{e}}^{t\mathcal{L}}=&&&\quad\frac{1}{\omega^{3}}\sum_{s=0}^{3}S_{3,s}^{3}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t+}&\ldots&+\frac{1}{\omega^{n}}\sum_{s=0}^{3}S_{n,s}^{3}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}&+E^{3}_{n}(t)\\\ \vdots&&&&\ddots&\vdots\\\
T^{n}{\mathrm{e}}^{t\mathcal{L}}=&&&&&\quad\frac{1}{\omega^{n}}\sum_{s=0}^{n}S_{n,s}^{n}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}&+E^{n}_{n}(t)\\\ \end{array}$

Summing up $r$-th column of the above table we obtain

$\frac{1}{\omega^{r}}\sum_{d=0}^{r}\sum_{s=0}^{d}S_{r,s}^{d}{\mathrm{e}}^{s{\mathrm{i}}\omega
t}=\frac{1}{\omega^{r}}\sum_{s=0}^{r}{\mathrm{e}}^{s{\mathrm{i}}\omega
t}\sum_{d=s}^{r}S_{r,s}^{d},\quad{\rm where}\quad S_{r,s}^{0}(t)=0.$

Hence, based also on (4.32) and (4.33), equation (4.31) is equivalent to

$\begin{split}&u^{\star}[t]\\\
&={\mathrm{e}}^{t\mathcal{L}}u_{0}+\sum_{r=0}^{n}\frac{1}{\omega^{r}}\sum_{s=0}^{r}{\mathrm{e}}^{s{\mathrm{i}}\omega
t}\sum_{d=s}^{r}S_{r,s}^{d}u_{0}+\sum_{d=1}^{n}E^{d}_{n}u_{0}+\sum_{d=n+1}^{\infty}T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}\\\
&={\mathrm{e}}^{t\mathcal{L}}u_{0}+\sum_{r=1}^{n}\frac{1}{\omega^{r}}\sum_{s=0}^{r}{\mathrm{e}}^{s{\mathrm{i}}\omega
t}\sum_{d=s}^{r}\sum_{\mathbf{|k^{d}|}=r-d}\frac{a_{s}^{d}}{({\mathrm{i}})^{r}}\mathcal{N}_{k_{d-s+1}}^{k_{d}}{\mathrm{e}}^{t\mathcal{L}}\mathcal{F}_{d-s}u_{0}+\sum_{d=1}^{n}E^{d}_{n}u_{0}+\sum_{d=n+1}^{\infty}\Big{(}\mathcal{S}^{d}_{d}[\mathcal{L},\alpha](t)+E^{d}_{d}[\mathcal{L},\alpha](t)\Big{)}\\\
&={\mathrm{e}}^{t\mathcal{L}}u_{0}+\sum_{r=0}^{n}\frac{1}{\omega^{r}}\sum_{s=0}^{r}{\mathrm{e}}^{s{\mathrm{i}}\omega
t}p_{r,s}(t)u_{0}+R_{n}(t)\\\ \end{split}$ (4.35)

where coefficients $p_{r,s}(t)$ and error $R_{n}(t)$ coincide with (4.28) and
(4.29), respectively, which proves the first part of the thesis. To prove
estimate (4.30) let us observe that

$R_{n}(t)=\sum_{d=1}^{n}E^{d}_{n}u_{0}+T^{n+1}\sum_{d=0}^{\infty}T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}.$

As it was pointed out in the proof of Theorem 3.1
$\|E^{d}_{n}[\mathcal{L},\alpha]\|_{Y}\leq\omega^{-n}\tilde{C}_{n}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)$
which does not suffice for the desired estimate of $R_{n}(t)$. Assumption 3
allows for a longer expansion of $T^{d}$

$T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}=\mathcal{S}_{n+1}^{d}[\mathcal{L},\alpha](t)+E_{n+1}^{d}[\mathcal{L},\alpha](t),$

where $\mathcal{S}_{n+1}^{d}[\mathcal{L},\alpha](t)$ and
$E_{n+1}^{d}[\mathcal{L},\alpha](t)$ are defined as previously which means
that $E^{r}_{n}(t)$ can be replaced with a sum
$\frac{1}{\omega^{n+1}}\sum_{s=0}^{1}S_{n+1,s}^{d}(t){\mathrm{e}}^{s{\mathrm{i}}\omega
t}+E^{d}_{n+1}(t)$ for $d=1,\ldots n$ which is bounded in $Y$ norm by
$\omega^{-(n+1)}C_{E}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)$. Based on
Theorem 2.1 we claim that

$\|T^{n+1}\sum_{d=0}^{\infty}T^{d}{\mathrm{e}}^{t\mathcal{L}}u_{0}\|_{Y}\leq
C_{t^{\star}}\|u_{0}\|_{L^{2}(\Omega)}\exp\Big{(}C_{t^{\star}}\|\alpha\|_{L^{2}(\Omega)}{t^{\star}}\Big{)}\left\|T^{n+1}\right\|_{Y},$

where $C_{t^{\star}}=\max_{0\leq
t\leq{t^{\star}}}\|e^{t\mathcal{L}}\|_{L^{2}(\Omega)}$, while from Theorem 3.1
we know that

$\|T^{n+1}\|_{Y}\leq\omega^{-(n+1)}C_{T}^{n+1}(t^{\star},\mathcal{L},u_{0},\alpha).$

Altogether we get that

$\|R_{n}\|_{Y}\leq\omega^{-(n+1)}\Big{(}\sum_{d=1}^{n}C_{E}^{d}(t^{\star},\mathcal{L},u_{0},\alpha)+C_{t^{\star}}\|u_{0}\|_{L^{2}(\Omega)}\exp\Big{(}C_{t^{\star}}\|\alpha\|_{L^{2}(\Omega)}{t^{\star}}C_{T}^{n+1}(t^{\star},\mathcal{L},u_{0},\alpha)\Big{)},$

which completes the proof. $\Box$

###### Example 4.1

The first few coefficients of the derived example are given with the following
formulas

$\displaystyle p_{0,0}(t)={\mathrm{e}}^{t\mathcal{L}}u_{0},$ $\displaystyle
p_{1,0}(t)=-\frac{1}{{\mathrm{i}}}{\mathrm{e}}^{t\mathcal{L}}\alpha u_{0},$
$\displaystyle
p_{1,1}(t)=\frac{1}{{\mathrm{i}}}\alpha{\mathrm{e}}^{t\mathcal{L}}u_{0},$
$\displaystyle
p_{2,0}(t)=-\frac{1}{2}{\mathrm{e}}^{t\mathcal{L}}\alpha^{2}u_{0}+{\mathrm{e}}^{t\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{1}\left(\alpha\right)u_{0},$
$\displaystyle
p_{2,1}(t)=-\mathrm{ad}\,_{\mathcal{L}}^{1}\left(\alpha\right){\mathrm{e}}^{t\mathcal{L}}u_{0}+\alpha{\mathrm{e}}^{t\mathcal{L}}\alpha
u_{0},$ $\displaystyle
p_{2,2}(t)=-\frac{1}{2}\alpha^{2}{\mathrm{e}}^{t\mathcal{L}}u_{0};$
$\displaystyle
p_{3,0}(t)={\mathrm{i}}{\mathrm{e}}^{t\mathcal{L}}\left(-\mathrm{ad}\,_{\mathcal{L}}^{2}(\alpha)-\frac{1}{4}\mathrm{ad}\,_{\mathcal{L}}^{1}(\alpha^{2})+\mathrm{ad}\,_{\mathcal{L}}^{1}(\alpha)\alpha+\frac{1}{2}\alpha\mathrm{ad}\,_{\mathcal{L}}^{1}(\alpha)-\frac{1}{6}\alpha^{3}\right)u_{0},$
$\displaystyle
p_{3,1}(t)={\mathrm{i}}\left(\mathrm{ad}\,_{\mathcal{L}}^{2}(\alpha){\mathrm{e}}^{t\mathcal{L}}-\mathrm{ad}\,_{\mathcal{L}}^{1}(\alpha){\mathrm{e}}^{t\mathcal{L}}\alpha-\alpha{\mathrm{e}}^{t\mathcal{L}}\mathrm{ad}\,_{\mathcal{L}}^{1}(\alpha)+\frac{1}{2}\alpha{\mathrm{e}}^{t\mathcal{L}}\alpha^{2}\right)u_{0},$
$\displaystyle
p_{3,2}(t)={\mathrm{i}}\left(\frac{1}{4}\mathrm{ad}\,_{\mathcal{L}}^{1}(\alpha^{2}){\mathrm{e}}^{t\mathcal{L}}+\frac{1}{2}\alpha\mathrm{ad}\,_{\mathcal{L}}^{1}(\alpha){\mathrm{e}}^{t\mathcal{L}}-\frac{1}{2}\alpha^{2}{\mathrm{e}}^{t\mathcal{L}}\alpha\right)u_{0},$
$\displaystyle
p_{3,3}(t)={\mathrm{i}}\frac{1}{6}\alpha^{3}{\mathrm{e}}^{t\mathcal{L}}u_{0}.$

##  5 Numerical examples

To illustrate the proposed approach we analyse two problems chosen in such a
way, that their exact solutions are known. This way we present the exact
errors between the analytical solutions $u^{\star}[t]$ and their first four
approximations

$U^{[n]}:=p_{0,0}(t)+\sum_{r=1}^{n}\frac{1}{\omega^{r}}\sum_{s=0}^{r}{\mathrm{e}}^{s{\mathrm{i}}\omega
t}p_{r,s}(t),\quad n=0,\ldots,3,\quad t\in[0,1],$ (5.36)

Example 1.

Let us start with damped heat equation equipped with initial–boundary
conditions

$\displaystyle\partial_{t}u(x,t)=\Delta
u(x,t)+{\mathrm{e}}^{{\mathrm{i}}\omega t}\cdot u(x,t),\quad x\in[0,\pi],\quad
t\in[0,1],$ $\displaystyle
u(x,0)=u_{0}(x)=\sin(x)\cdot\exp\left(-\frac{{\mathrm{i}}}{\omega}\right),$
(5.37) $\displaystyle u(0,t)=u(\pi,t)=0,$

with known solution

$u^{\star}_{1}[t](x)={\mathrm{e}}^{-t-{\mathrm{i}}{\mathrm{e}}^{{\mathrm{i}}\omega
t}/\omega}\sin(x).$

(a) $\omega=10$

(b) $\omega=100$

Figure 5.1: The base-10 logarithmic errors
$R^{[n]}_{1}=\|u^{\star}_{1}-U^{[n]}_{1}\|_{Y}$, $n=0,\ldots,3$ obtained for
the heat equation (5).

Let us observe that the initial condition $u_{0}(x)$ and function $\alpha(x)$
do not affect error constant. Indeed, $u_{0}$ is non-oscillatory while
$\alpha$ is constant - their derivatives are small and vanish, respectively.
Thus, according to the expectations the error decays as the number of summands
in (5.36) grows. Moreover we can see that the method improves with the grow of
the oscillatory parameter $\omega$. This phenomena is presented in Figure
5.1(a) \- 5.1(b), where errors
$R^{[n]}_{1}:=\|u^{\star}_{1}-U^{[n]}_{1}\|_{Y}$, $n=0,\ldots,3$ are displayed
in logarithmic scale.

Example 2.

Let us consider now transport equation

$\displaystyle\partial_{t}u(x,t)=\nabla
u(x,t)+\cos(c(x^{2}+x))\cdot{\mathrm{e}}^{{\mathrm{i}}\omega t}\cdot
u(x,t),\quad c={\rm const},\quad x\in[-1,1],\quad t\in[0,1],$ (5.38)
$\displaystyle u(x,0)=u_{0}(x)=x,$

with known solution

$\displaystyle u^{\star}_{2}[t](x)=(t+x)\cdot$
$\displaystyle\cdot\exp\left(\frac{(-1)^{3/4}\sqrt{\pi}\exp(f_{1}(x,t))\left(e^{\frac{i\left(c^{2}+\omega^{2}\right)}{2c}}\left(\text{erf}\left(g_{1}(x)\right)-\text{erf}\left(f_{2}(x,t)\right)\right)-\text{erfi}\left(f_{3}(x,t)\right)+\text{erfi}\left(g_{2}(x)\right)\right)}{4\sqrt{c}}\right)$

where

$\displaystyle f_{1}(x,t)$ $\displaystyle=$
$\displaystyle-\frac{i\left(c^{2}-2c\omega(2t+2x+1)+\omega^{2}\right)}{4c},$
$\displaystyle g_{1}(x)$ $\displaystyle=$
$\displaystyle\frac{\sqrt[4]{-1}(2cx+c+\omega)}{2\sqrt{c}},$ $\displaystyle
f_{2}(x,t)$ $\displaystyle=$
$\displaystyle\frac{\sqrt[4]{-1}(2ct+2cx+c+\omega)}{2\sqrt{c}},$
$\displaystyle f_{3}(x,t)$ $\displaystyle=$
$\displaystyle\frac{\sqrt[4]{-1}(2c(t+x)+c-\omega)}{2\sqrt{c}},$
$\displaystyle g_{2}(x)$ $\displaystyle=$
$\displaystyle\frac{\sqrt[4]{-1}(2cx+c-\omega)}{2\sqrt{c}},$

and $\text{erf}(x)$, $\text{erfi}(x)$ are the error function and the imaginary
error function respectively, i.e
$\text{erf}(x)=\frac{2}{\sqrt{\pi}}\int_{0}^{x}{\mathrm{e}}^{-x^{2}}\,\mathrm{d}x$,
$\text{erfi}(x)=-{\mathrm{i}}\text{erf}({\mathrm{i}}x)$. In this example we
wish to show how function $\alpha(x)=\cos(c(x^{2}+x))$ may affect the error
constant. In the non-oscillatory case of $\alpha$, where $c=1$, the error
decays with number $n$ in (5.36), see Figure 5.2(a). An interesting situation
takes place when $\alpha$ becomes oscillatory, for example $c=31$. Derivatives
of $\alpha$ grow fast resulting in $R^{[2]}_{2}>\omega^{-2}$ and
$R^{[3]}_{2}>\omega^{-3}$, what can be observed in Figure 5.2(b).

(a) $c=1,\quad\omega=100$

(b) $c=31,\quad\omega=100$

Figure 5.2: The base-10 logarithmic errors
$R^{[n]}_{2}=\|u^{\star}_{2}-U^{[n]}_{2}\|_{Y}$, $n=0,\ldots,3$ obtained for
the heat equation (5.38)

In the table below we present first few functions $p_{r,n}$ of both examples.

$p_{r,s}$ |  $\mathcal{L}=\nabla$, $\alpha(x)=\cos(c(x^{2}+x)),\ u_{0}(x)=x$ |  $\mathcal{L}=\Delta$, $\alpha(x)\equiv 1$, $u_{0}(x)=\sin(x){\mathrm{e}}^{-{\mathrm{i}}/\omega}$
---|---|---
$p_{0,0}$ | $u_{0}(x+t)$ | ${\mathrm{e}}^{-t}\sin(x){\mathrm{e}}^{-{\mathrm{i}}/\omega}$
$p_{1,0}$ | $-\frac{1}{{\mathrm{i}}}\alpha(x+t)u_{0}(x+t)$ | $-\frac{1}{{\mathrm{i}}}{\mathrm{e}}^{-t}\sin(x){\mathrm{e}}^{-{\mathrm{i}}/\omega}$
$p_{1,1}$ | $\frac{1}{{\mathrm{i}}}\alpha(x)u_{0}(x+t)$ | $\frac{1}{{\mathrm{i}}}{\mathrm{e}}^{-t}\sin(x){\mathrm{e}}^{-{\mathrm{i}}/\omega}$
$p_{2,0}$ | $-\frac{1}{2}\alpha^{2}(x+t)u_{0}(x+t)+\alpha^{\prime}(x+t)u_{0}(x+t)$ | $-\frac{1}{2}{\mathrm{e}}^{-t}\sin(x){\mathrm{e}}^{-{\mathrm{i}}/\omega}$
$p_{2,1}$ | $-u_{0}(x+t)\alpha^{\prime}(x)+\alpha(x)\alpha(x+t)u_{0}(x+t)$ | ${\mathrm{e}}^{-t}\sin(x){\mathrm{e}}^{-{\mathrm{i}}/\omega}$
$p_{2,2}$ | $-\frac{1}{2}\alpha^{2}(x)u_{0}(x+t)$ | $-\frac{1}{2}{\mathrm{e}}^{-t}\sin(x){\mathrm{e}}^{-{\mathrm{i}}/\omega}$
Table 1: Coefficients $p_{r,n}$ of asymptotic expansion of Example 1. and
Example 2.

##  6 Acknowledgement

We would like to thank Arieh Iserles from the University of Cambridge and
Marcus Webb from the University of Manchester for numerous stimulating
discussions which were very helpful and inspirational.

We wish also to thank Anna Kamont from Institute of Mathematics Polish Academy
of Sciences for her support, deep observations and rigorous approach to the
contents of this paper, Jakub Skrzeczkowski from University of Warsaw for
fruitful discussions about semi-group theory and Karolina Lademann from the
University of Gdańsk for contribution to the performance of numerical
examples.

The work of KK in this project was supported by The National Center for
Science (NCN), based on Grant No. 2019/34/E/ST1/00390. This work was also
partially supported by the Simons Foundation Award No 663281 granted to the
Institute of Mathematics of the Polish Academy of Sciences for the years
2021-2023. Numerical simulations were carried out at the Centre of Informatics
Tricity Academic Supercomputer and networK (CI TASK) in Gdańsk.

## References

  * [Condon et al.(2009)Condon, Deaño, and Iserles] M. Condon, A. Deaño, A. Iserles,  On highly oscillatory problems arising in electronic engineering,  M2AN Math. Model. Numer. Anal. 43 (2009) 785–804. URL: https://doi.org/10.1051/m2an/2009024. doi:10.1051/m2an/2009024.
  * [Condon et al.(2010)Condon, Deaño, and Iserles] M. Condon, A. Deaño, A. Iserles,  On systems of differential equations with extrinsic oscillation,  Discrete Contin. Dyn. Syst. 28 (2010) 1345–1367. URL: https://doi.org/10.3934/dcds.2010.28.1345. doi:10.3934/dcds.2010.28.1345.
  * [Condon et al.(2012)Condon, Deaño, Iserles, and Kropielnicka] M. Condon, A. Deaño, A. Iserles, K. Kropielnicka,  Efficient computation of delay differential equations with highly oscillatory terms,  ESAIM Math. Model. Numer. Anal. 46 (2012) 1407–1420. URL: https://doi.org/10.1051/m2an/2012004. doi:10.1051/m2an/2012004.
  * [Faou and Schratz(2014)] E. Faou, K. Schratz,  Asymptotic preserving schemes for the Klein-Gordon equation in the non-relativistic limit regime,  Numer. Math. 126 (2014) 441–469. URL: https://doi.org/10.1007/s00211-013-0567-z. doi:10.1007/s00211-013-0567-z.
  * [Condon et al.(2019)Condon, Iserles, Kropielnicka, and Singh] M. Condon, A. Iserles, K. Kropielnicka, P. Singh,  Solving the wave equation with multifrequency oscillations,  J. Comput. Dyn. 6 (2019) 239–249. URL: https://doi.org/10.3934/jcd.2019012. doi:10.3934/jcd.2019012.
  * [Condon et al.(2021)Condon, Kropielnicka, Lademann, and Perczyński] M. Condon, K. Kropielnicka, K. Lademann, R. Perczyński,  Asymptotic numerical solver for the linear Klein-Gordon equation with space- and time-dependent mass,  Appl. Math. Lett. 115 (2021) Paper No. 106935, 7. URL: https://doi.org/10.1016/j.aml.2020.106935. doi:10.1016/j.aml.2020.106935.
  * [Hairer et al.(2010)Hairer, Lubich, and Wanner] E. Hairer, C. Lubich, G. Wanner, Geometric numerical integration, volume 31 of Springer Series in Computational Mathematics, Springer, Heidelberg, 2010. Structure-preserving algorithms for ordinary differential equations, Reprint of the second (2006) edition.
  * [Ernst and Christian(2001)] H. Ernst, L. Christian,  Long-time energy conservation of numeri- cal methods for oscillatory differential equations,  SIAM J. Numer. Anal. 38 (2001) 414–441. URL: https://www.jstor.org/stable/3062063.
  * [Gauckler et al.(2018)Gauckler, Hairer, and Lubich] L. Gauckler, E. Hairer, C. Lubich,  Dynamics, numerical analysis, and some geometry,  in: Proceedings of the International Congress of Mathematicians—Rio de Janeiro 2018. Vol. I. Plenary lectures, World Sci. Publ., Hackensack, NJ, 2018, pp. 453–485.
  * [Wang and Zhao(2021)] B. Wang, X. Zhao,  Error estimates of some splitting schemes for charged-particle dynamics under strong magnetic field,  SIAM J. Numer. Anal. 59 (2021) 2075–2105. URL: https://doi.org/10.1137/20M1340101. doi:10.1137/20M1340101.
  * [Wang and Zhao(2023)] B. Wang, X. Zhao,  Geometric two-scale integrators for highly oscillatory system: uniform accuracy and near conservations,  SIAM J. Numer. Anal. 61 (2023) 1246–1277. URL: https://doi.org/10.1137/21M1462908. doi:10.1137/21M1462908.
  * [Iserles and Norsett(2006)] A. Iserles, S. Norsett,  Quadrature methods for multivariate highly oscillatory integrals using derivatives,  Math. Comput. 75 (2006) 1233–1258. doi:10.1090/S0025-5718-06-01854-0.
  * [Olver(2008)] S. Olver, Numerical approximation of highly oscillatory integrals, Ph.D. thesis, University of Cambridge Cambridge, 2008\.