# Efficient approximation of high-frequency Helmholtz solutions by Gaussian coherent states

Théophile Chaumont-Frelet, Victorita Dolean, Maxime Ingremeau

► To cite this version:

Théophile Chaumont-Frelet, Victorita Dolean, Maxime Ingremeau. Efficient approximation of high-frequency Helmholtz solutions by Gaussian coherent states. Numerische Mathematik, 2024, 10.1007/s00211-024-01411-0. hal-03747290

HAL Id: hal-03747290

https://inria.hal.science/hal-03747290v1

Submitted on 8 Aug 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
EFFICIENT APPROXIMATION OF HIGH-FREQUENCY HELMHOLTZ
SOLUTIONS BY GAUSSIAN COHERENT STATES

THÉOPHILE CHAUMONT-FRELET*, VICTORITA DOLEAN†,‡, AND MAXIME INGREMEAU†

ABSTRACT. We introduce new finite-dimensional spaces specifically designed to approximate the solutions to high-frequency Helmholtz problems with smooth variable coefficients in dimension d. These discretization spaces are spanned by Gaussian coherent states, that have the key property to be localised in phase space. We carefully select the Gaussian coherent states spanning the approximation space by exploiting the (known) micro-localisation properties of the solution. For a large class of source terms (including plane-wave scattering problems), this choice leads to discrete spaces that provide a uniform approximation error for all wavenumber k with a number of degrees of freedom scaling as $k^{d-1/2}$, which we rigorously establish. In comparison, for discretization spaces based on (piecewise) polynomials, the number of degrees of freedom has to scale at least as $k^d$ to achieve the same property. These theoretical results are illustrated by one-dimensional numerical examples, where the proposed discretization spaces are coupled with a least-squares variational formulation.

KEY WORDS. Gabor frames, Helmholtz equation, high-frequency problems

# 1. INTRODUCTION

Time-harmonic wave propagation is a mechanism at the center of a large number of physical and industrial applications. We may cite, among many, radar imaging [8], or seismic prospection [36]. In practice, numerical methods are required to approximately simulate the propagation of waves, and although several methods are available, it is still very challenging to compute accurate approximations in high-frequency regime.

Here, we consider the scalar Helmholtz equation, which is probably the simplest model for this kind of problems. Specifically, given a compactly supported right-hand side $f: \mathbb{R}^d \to \mathbb{C}$, our model problem is to find $u: \mathbb{R}^d \to \mathbb{C}$ such that

$$ (1.1a) \qquad -k^2\mu u - \nabla \cdot (A\nabla u) = f \quad \text{in } \mathbb{R}^d, $$

where $\mu$ and $A$ are (given) smooth coefficients that are respectively equal to 1 and $I$ outside a ball of radius $R > 0$, and $k > 0$ is the (given) wavenumber. This equation is supplemented with the Sommerfeld radiation condition at infinity. Namely, we require that

$$ (1.1b) \qquad \frac{\partial u}{\partial |\mathbf{x}|}(\mathbf{x}) - iku(\mathbf{x}) = o\left(|\mathbf{x}|^{-(d+1)/2}\right) \quad \text{as } |\mathbf{x}| \to +\infty. $$

A particularly important scenario covered by (1.1) is the scattering of a plane-wave, where the right-hand side takes the form

$$ (1.2) \qquad f := (k^2\mu + \nabla \cdot (A\nabla \cdot)) e^{ik\mathbf{d} \cdot \mathbf{x}}, $$

where $\mathbf{d} \in \mathbb{R}^d$, $|\mathbf{d}| = 1$ is the incidence direction.

As we propose a “volumic” method, we will actually replace the Sommerfeld radiation condition (1.1b) by a Perfectly Matched Layer (PML). This approach is entirely standard [2, 6, 14], and amounts to slightly modifying the coefficients $\mu$ and $A$ in (1.1). This process is detailed in Section 4.2.

<footnote>*University Côte d'Azur, Inria, CNRS, LJAD</footnote>

<footnote>†University Côte d'Azur, CNRS, LJAD</footnote>

<footnote>‡Department of Mathematics and Statistics, University of Strathclyde</footnote>
In this work, we investigate the use of discretization spaces based on Gaussian coherent states (GCS), that is, functions of the form

$$ \Phi_{k,x_0,\xi_0}(\mathbf{x}) := \left(\frac{k}{\pi}\right)^{d/4} e^{-\frac{k^2}{2}|\mathbf{x}-\mathbf{x}_0|^2} e^{-ik\xi_0 \cdot (\mathbf{x}-\mathbf{x}_0)}, $$

where $x_0, \xi_0 \in \mathbb{R}^d$ are user-selected parameters. The idea of decomposing a function as a discrete sum of Gaussian coherent states goes back to [13]. Here, following [4, 7], we focus on a lattice of phase-space points $[\mathbf{x}^{k,m}, \mathbf{\xi}^{k,n}] := \sqrt{(kR)^{-1}\pi}[\mathbf{m}, \mathbf{n}], [\mathbf{m}, \mathbf{n}] \in \mathbb{Z}^{2d}$, that are spaced by $\sim (kR)^{-1/2}$, and we write

$$ \Psi_{k,m,n} := \Phi_{k,x^{k,m},\xi^{k,n}}. $$

Then, our a discretization space is of the form

$$ W_{\Lambda} := \mathrm{Vect}\{\Psi_{k,m,n}; [\mathbf{m}, \mathbf{n}] \in \Lambda\}, $$

where $\Lambda \subset \mathbb{Z}^{2d}$ is a carefully chosen set of indices.

For the sake of simplicity, we will assume in the introduction that the domain is non-trapping. Our first result is that if $\Lambda$ is chosen as

$$ \Lambda := \{[\mathbf{m}, \mathbf{n}] \in \mathbb{Z}^{2d}; |\mathbf{x}^{k,m}|^2 + |\mathbf{\xi}^{k,n}|^2 \leq \rho k R\}, $$

for $\rho > 0$, we have

$$ (1.3) \qquad \dim W_{\Lambda} \approx \rho^d (kR)^d \quad \text{and} \quad \min_{w \in W_{\Lambda}} \|u-w\|_{\hat{H}^1(\mathbb{R}^d)} \leq C \rho^{-1} \|f\|_{L^2(\mathbb{R}^d)}, $$

for a general $f \in L^2(\mathbb{R}^d)$ with supp $f \subset B(0, R)$, where $\|\cdot\|_{\hat{H}^1(\mathbb{R}^d)}$ is a $H^1(\mathbb{R}^d)$-norm including a weight at infinity (see (4.1) below). As we describe in more length afterwards, this result is not very impressive on its own. Specifically, it is a standard approximation result similar to polynomial approximations: we need a fixed number of points per wavelength to achieve a constant accuracy. Our second result, which is key, deals with the case where $f$ takes the particular form (1.2). In this case, we select

$$ \Lambda := \{[\mathbf{m}, \mathbf{n}] \in \mathbb{Z}^{2d}; |p(\mathbf{x}^{k,m}, \mathbf{\xi}^{k,n})| \leq k^{-1/2+\epsilon}\} $$

where $\epsilon > 0$ can be selected arbitrarily small and

$$ p(\mathbf{x}, \mathbf{\xi}) := \mathbf{A}(\mathbf{x})\mathbf{\xi} \cdot \mathbf{\xi} - \mu(\mathbf{x}) \quad \forall \mathbf{x}, \mathbf{\xi} \in \mathbb{R}^d, $$

is the principal symbol associated with the differential operator in (1.1), and we have

$$ (1.4) \qquad \dim W_{\Lambda} \approx (kR)^{d-1/2+\epsilon} \quad \text{and} \quad \min_{w \in W_{\Lambda}} \|u-w\|_{\hat{H}^1(\mathbb{R}^d)} \leq C_{\epsilon,m} (kR)^{-m} \quad \forall m \in \mathbb{N}. $$

It means that for right-hand sides corresponding to scattering problems (and actually, a wider family of right-hand sides), Gaussian coherent states provide an accurate solution with $O((kR)^{d-1/2+\epsilon})$ DOFs. In fact, the convergence is even super-algebraic as the frequency increases.

To put (1.3) and (1.4) into perspective, we compare them with other standard methods. Actually, there are several options to numerically solve (1.1) (either with the Sommerfeld condition (1.1b) or with a PML approximation), that we review below.

The most versatile approach is probably the finite element method (FEM). The method hinges of a triangulation of the domain into elements of size $h$, and piecewise polynomial basis functions of degree $p$. It can be shown that if $p$ grows logarithmically with $k$, then the condition that $kh/p$ is constant provides (at least) a constant accuracy as $k$ increases [26, 29, 30]. As a result, high-order FEM essentially requires $O((kR)^d)$ degrees of freedom (DOFs) to achieve a constant accuracy. The resulting matrix is sparse.

Trefftz-like methods are similar to FEM in that they also rely on a mesh of the domain, but the polynomial shape functions are replaced by local solutions to the Helmholtz problem, such as plane-waves [23]. There are many ways to “glue” these local solution together, including partition of unity methods [28], least squares methods [31], the ultra weak variational method [3], the discontinuous enrichment method [12] or the variational theory of complex rays [34]. While these methods typically induce a large reduction of the number of DOFs as compared to FEM, they
<table><caption>TABLE 1.1. Comparison of commonly used discretization techniques</caption><thead><tr><th></th><th>FEM</th><th>GCS</th><th>BEM</th><th>asymptotic</th></tr></thead><tbody><tr><th>Cost</th><td>(kR)<sup>d</sup></td><td>(kR)<sup>d-1/2</sup></td><td>(kR)<sup>d-1</sup></td><td>1</td></tr><tr><th>Heterogeneous media</th><td>yes</td><td>yes</td><td>no</td><td>yes</td></tr><tr><th>Convergence at fixed k</th><td>yes</td><td>yes</td><td>yes</td><td>no</td></tr></tbody></table>

usually still need at least $O((kR)^d)$ DOFs, see, e.g., [5, 17, 22]. Similar to FEM, the resulting matrix is sparse.

The next family of methods we want to mention are boundary element methods (BEM) [35]. These methods rely on boundary integral equations which, strictly speaking, are not available for smoothly varying coefficients, since the expression of Green's function must be available. It is nevertheless interesting to include them in the comparison. These methods typically provide a constant accuracy with only $O((kR)^{d-1})$ DOFs [15]. However, the resulting matrix is dense and its entries are costly to compute. These issues can be mitigated using compression techniques, such as the fast multi-pole method [20] or hierarchical matrices [21].

Finally, asymptotic methods rely on the fact when the frequency is very large, it is sometimes possible to simplify the search for the solution of Helmholtz equation and the properties of the solution itself to computations involving only the underlying classical dynamics [10]. This is done using tools of semi-classical analysis, such as the WKB method and can lead to discrete problems with a number of DOFs independent of $k$. The main drawback of these approaches is that they are only asymptotically valid: they do not converge for fixed value of $k$. Besides, it is not always clear from which range of $k$ they are relevant.

As compared to FEM, the proposed GCS method thus gains “half a dimension” at high-frequencies, but it is still half a dimension higher than BEM. As compared to BEM however, our methodology has the advantage to apply in a generic framework where the Green’s function is not available. Another important comment is that in the (very) high-frequency regime, our method is more expensive than asymptotic methods. However, asymptotic methods cannot converge at fixed $k$, which our method does. This is summarized in Table 1.1.

In addition to the approximability results (1.3) and (1.4), we also present a least-squares method based on Gaussian coherent states for Problem (1.1). As we show, the convergence of the method is easily established. Besides, although the matrix is dense, we show that the entries decay algebraically away from the diagonal. As a result, the matrix is essentially banded, and we believe that efficient linear solvers can be devised. This will be analysed in more depth in future works.

We finally present a set of one-dimensional numerical experiments using the proposed least-squares method. Although the setting is rather simple, the results perfectly fit the theoretical analysis and readily shows that proposed approach allows for a drastic reduction of the number of DOFs in the high-frequency regime.

To the best of our knowledge, our micro-locally adapted spaces of Gaussian coherent states appear to be entirely original, but we would like to mention that similar basis functions have already been employed to discretize PDE problems. In particular, generalised coherent states like Hagerdorn wavepackets were used to describe the solution of time-dependent Schrödinger equation in [11, 18, 19, 27].

The remainder of our work is organised as follows. In Section 2, we precise the setting and state our key approximation result (1.4) in its more general form. Section 3 contains the proof of our findings. In Section 4, we apply the general theory of Sections 2 and 3 to our scattering model problem. Numerical examples are reported in Section 5. Finally, Appendix A collects technical results concerning Gaussian coherent states.

## 2. SETTING AND MAIN RESULTS

2.1. **Notations.** Throughout this work $\hbar \in \mathcal{H} \subset (0, 1]$ will denote a small parameter. When applying our general results to the Helmholtz equation, we will have $\hbar \sim (kR)^{-1}$, so that considering the set $(0, 1]$ amounts to ignoring low frequencies, and focusing on high frequencies when $\hbar \to 0$.
For the sake of generality, we restrict our analysis to a subset $\mathcal{H} \subset (0, 1]$ for reasons that will become apparent in Section 4. Notice that the case $\mathcal{H} = (0, 1]$ is not excluded.

2.1.1. Basic notation. The canonical basis of $\mathbb{R}^d$ or of $\mathbb{C}^d$ will be denoted by $(e_1, ..., e_d)$. If $\mathbf{x}, \mathbf{y} \in \mathbb{C}^d$, we write

$$ \mathbf{x} \cdot \mathbf{y} := \sum_{j=1}^{d} x_j y_j $$

without complex conjugation on the second argument, and $|\mathbf{x}|^2 = \mathbf{x} \cdot \bar{\mathbf{x}}$ is the usual Euclidean norm.

For a multi-index $\alpha \in \mathbb{N}^d$, $[\alpha] := \alpha_1 + \dots + \alpha_d$ denotes its usual $\ell_1$ norm. If $v : \mathbb{R}^d \to \mathbb{C}$, the notation

$$ \partial^{\alpha} v := \frac{\partial^{\alpha_1}}{\partial x_1} \cdots \frac{\partial^{\alpha_d}}{\partial x_d} v $$

is employed for the partial derivatives in the sense of distributions, whereas $\mathbf{x}^{\alpha} := x_1^{\alpha_1} \cdots x_d^{\alpha_d}$. Finally, if $\beta \in \mathbb{N}^d$ is another multi-index, we will sometimes need the notation

$$ \begin{pmatrix} \alpha \\ \beta \end{pmatrix} = \prod_{j=1}^{d} \begin{pmatrix} \alpha_j \\ \beta_j \end{pmatrix}, $$

and the notation $\alpha \le \beta$ means that $\alpha_j \le \beta_j$ for all $j \in \{1, \dots, d\}$.

If $\mathbf{n} \in \mathbb{Z}^d$, we employ the notation $|\mathbf{n}|^2 := n_1^2 + \dots + n_d^2$ for its $\ell_2$ norm. Finally, if $\Lambda \subset \mathbb{Z}^{2d}$, $\ell^2(\Lambda)$ has its usual definition, and we denote by $\| \cdot \|_{L^2(\mathbb{R}^d)}$ its usual norm.

2.1.2. Key functional spaces. In what follows, $L^2(\mathbb{R}^d)$ is the usual Lebesgue space of complex-valued square integrable functions over $\mathbb{R}^d$. The usual norm and inner products of $L^2(\mathbb{R}^d)$ are respectively denoted by $\| \cdot \|_{L^2(\mathbb{R}^d)}$ and $(\cdot, \cdot)$.

Following [4], since we are dealing with the (unbounded) $\mathbb{R}^d$ space, our analysis will require the weighted Sobolev spaces

$$ \hat{H}^p(\mathbb{R}^d) := \{ v \in L^2(\mathbb{R}^d) \mid \mathbf{x}^{\alpha} \partial^{\beta} v \in L^2(\mathbb{R}^d) \quad \forall \alpha, \beta \in \mathbb{N}^d; [\alpha], [\beta] \le p \}, $$

that we equip with the family of equivalent $\hbar$-weighted norms given by

$$ \|v\|_{\hat{H}_h^p(\mathbb{R}^d)}^2 := \sum_{[\alpha] \le p} \sum_{q \le p - [\alpha]} \hbar^{2[\alpha]} \|\mathbf{x}|^q \partial^\alpha v\|_{L^2(\mathbb{R}^d)}^2 $$

for all $p \in \mathbb{N}$.

$C^0(\mathbb{R}^d)$ is the set of complex-valued continuous functions defined over $\mathbb{R}^d$, and $C^\ell(\mathbb{R}^d)$ is the set of functions $v : \mathbb{R}^d \to \mathbb{C}$ such that $\partial^\alpha v \in C^0(\mathbb{R}^d)$ for all $\alpha \in \mathbb{N}^d$ with $[\alpha] \le \ell$. We introduce the notation

$$ \|v\|_{C^\ell(\mathbb{R}^d)} := \max_{\substack{\alpha \in \mathbb{N}^d \\ [\alpha] \le \ell}} \max_{\mathbf{x} \in \mathbb{R}^d} |(\partial^\alpha v)(\mathbf{x})|, \quad \forall v \in C^\ell(\mathbb{R}^d) $$

and $C_b^\ell(\mathbb{R}^d)$ is the subset of functions $v \in C^\ell(\mathbb{R}^d)$ such that $\|v\|_{C^\ell(\mathbb{R}^d)} < +\infty$. We also set

$$ C_b^\infty(\mathbb{R}^d) := \bigcap_{\ell \in \mathbb{N}} C_b^\ell(\mathbb{R}^d). $$

Finally, if $\Omega \subset \mathbb{R}^d$ is an open set, we denote by $C_c^\infty(\Omega)$ the set of smooth functions whose support is a compact subset of $\Omega$.
2.2. **The frame of Gaussian coherent states.** The goal of this work is to efficiently approximate the solution $u_\hbar$ to the equation $P_\hbar u_\hbar = f$ with a finite span of Gaussian coherent states. For $[m, n] \in \mathbb{Z}^{2d}$, we thus consider the Gaussian state

$$ \Psi_{\hbar, m, n}(\mathbf{x}) := (\pi\hbar)^{-d/2} e^{-\frac{1}{2\hbar}|\mathbf{x}-\mathbf{x}^{\hbar, m}|^2} e^{\frac{i}{\hbar} \boldsymbol{\xi}^{\hbar, n} \cdot (\mathbf{x}-\mathbf{x}^{\hbar, m})}, $$

where $\mathbf{x}^{\hbar, m} := \sqrt{\pi\hbar m}$ and $\boldsymbol{\xi}^{\hbar, n} := \sqrt{\pi\hbar n}$. The family of Gaussian coherent states $(\Psi_{\hbar, m, n})_{[m, n] \in \mathbb{Z}^{2d}}$ actually forms a frame over $L^2(\mathbb{R}^d)$, meaning there exists two constants $0 < \alpha < \beta < +\infty$ solely depending on $d$ such that

$$ \alpha \|v\|_{L^2(\mathbb{R}^d)}^2 \le \sum_{[m,n] \in \mathbb{Z}^{2d}} |(v, \Psi_{\hbar,m,n}^*)(\mathbf{x})|^2 \le \beta \|v\|_{L^2(\mathbb{R}^d)}^2. $$

This result was first proved in [7], but the idea of decomposing a function as a discrete sum of Gaussian states goes back to [13], where it was proved that the span of $(\Psi_{\hbar, m, n})_{[m, n] \in \mathbb{Z}^{2d}}$ is dense in $L^2(\mathbb{R}^d)$.

Actually, the frame property implies that there exists another family of functions $(\Psi_{\hbar, m, n}^*)([m, n] \in \mathbb{Z}^{2d})$ called the dual frame such that

$$ (2.1) \qquad v = \sum_{[m,n] \in \mathbb{Z}^{2d}} (v, \Psi_{\hbar,m,n}^*) \Psi_{\hbar,m,n} $$

for all $v \in L^2(\mathbb{R}^d)$. It is thus clear that any $v \in L^2(\mathbb{R}^d)$ may be well-approximated by (a large number of) Gaussian states. As we are going to develop hereafter, when considering the solution to a high-frequency PDE problem, a good approximation may be obtained with few Gaussian states, by carefully selecting the indices $[m, n]$ in (2.1).

**Remark 2.1** (General expansions in the Gaussian frame). The family $(\Psi_{\hbar, m, n})$ is not a Riesz basis, so that the expansion (2.1) of $v$ as a sum of $\Psi_{\hbar, m, n}$ is not unique. However, a crucial property of (2.1) is that this expansion is stable, in the sense that

$$ \sum_{[m,n] \in \mathbb{Z}^{2d}} |(v, \Psi_{\hbar,m,n}^*)|^2 \le \gamma \|v\|_{L^2(\mathbb{R}^d)}^2, $$

where $\gamma$ only depends on $d$, the dual frame being itself a frame. This is especially important at the numerical level in the presence of round-off errors [1].

2.3. **Settings and key assumptions.** Throughout this work, we consider a second order differential operator on $\mathbb{R}^d$ depending on $\hbar$, and taking the form

$$ (2.2) \qquad (P_\hbar v)(\mathbf{x}) = \hbar^2 \sum_{j,\ell=1}^{d} a_{j\ell}^\hbar(\mathbf{x}) \frac{\partial^2 v}{\partial x_j \partial x_\ell}(\mathbf{x}) + i\hbar \sum_{j=1}^{d} b_j^\hbar(\mathbf{x}) \frac{\partial v}{\partial x_j}(\mathbf{x}) + c^\hbar(\mathbf{x})v(\mathbf{x}), $$

where $a_{j\ell}^\hbar, b_j^\hbar, c^\hbar \in C_b^\infty(\mathbb{R}^d)$ for $1 \le j, \ell \le d$. For the sake of simplicity, we introduce

$$ C_{\text{coef},p} := \sup_{\hbar \in \mathcal{H}} \left( \sum_{j,\ell=1}^{d} \|a_{j,\ell}^{\hbar}\|_{C^p(\mathbb{R}^d)} + \sum_{j=1}^{d} \|b_j^{\hbar}\|_{C^p(\mathbb{R}^d)} + \|c^{\hbar}\|_{C^p(\mathbb{R}^d)} \right) \quad \forall p \in \mathbb{N}, $$

and assume that $C_{\text{coef},p} < +\infty$ for all $p \in \mathbb{N}$.

The principal symbol of $P_\hbar$ is the function $p_\hbar \in C^\infty(\mathbb{R}^{2d})$ defined by

$$ (2.3) \qquad p_\hbar(\mathbf{x}, \boldsymbol{\xi}) := \sum_{j,\ell=1}^{d} a_{j,\ell}^\hbar(\mathbf{x}) \xi_i \xi_j + \sum_{j=1}^{d} b_j^\hbar(\mathbf{x}) \xi_j + c^\hbar(\mathbf{x}). $$

For the sake of shortness, we will often write $p_\hbar(m, n) := p_\hbar(x^{\hbar, m}, \boldsymbol{\xi}^{\hbar, n})$ for $[m, n] \in \mathbb{Z}^{2d}$.

Along with the smoothness of the coefficients, we make two key assumptions. First, we assume that $P_\hbar : \hat{H}^p(\mathbb{R}^d) \to \hat{H}^p(\mathbb{R}^d)$ is invertible with the norm of $P_\hbar^{-1}$ being polynomially bounded in $\hbar$. Specifically, we assume that for all $f \in L^2(\mathbb{R}^d)$ there exists a unique $u_\hbar \in L^2(\mathbb{R}^d)$ such that $P_\hbar u_\hbar = f$. In addition, there exists $N \in \mathbb{N}$ such that for all $p \in \mathbb{N}$

$$ (2.4) \qquad \|u_\hbar\|_{\hat{H}_\hbar^p(\mathbb{R}^d)} \le C_{\text{sol},p} \hbar^{-N} \|f\|_{\hat{H}_\hbar^p(\mathbb{R}^d)} \quad \forall \hbar \in \mathcal{H} $$
for some constant $C_{\text{sol},p}$ independent of $\hbar$. Our second assumption is that there exists a value $\delta_0 > 0$ such that

$$ (2.5) \qquad \sup_{\mathbf{h} \in \mathcal{H}} \operatorname{diam}\{(\mathbf{x}, \boldsymbol{\xi}) \in \mathbb{R}^{2d}; |p_{\mathbf{h}}(\mathbf{x}, \boldsymbol{\xi})| < \delta_0\} \le D_0 < +\infty $$

for some $D_0 \in \mathbb{R}$.

In the remainder of this work, we allow generic constants $C$ to depend on $\{C_{\text{coef},p}\}_{p \in \mathbb{N}}$, $\{C_{\text{sol},p}\}_{p \in \mathbb{N}}$, $N$ and $D_0$. We also employ the notation $C_{\alpha,\beta,...}$ if the constant $C$ is additionally allowed to depend on other previously introduced quantity $\alpha, \beta, ...$

**2.4. Statement of the approximability result.** Our main result is that, if $f$ is micro-localised near the set $\{(\mathbf{x}, \boldsymbol{\xi}) \in \mathbb{R}^{2d}; p_{\mathbf{h}}(\mathbf{x}, \boldsymbol{\xi}) = 0\}$, then so is the solution $u_{\mathbf{h}}$ to $P_{\mathbf{h}}u_{\mathbf{h}} = f$. This is a standard result when micro-localisation is understood in terms of pseudo-differential operators (see for instance [37, Theorem 6.4]), but here, by micro-localisation properties, we mean that $f$ can be approached by a linear combination of $\Psi_{h,m,n}$ with $[m, n]$ in a certain region of $\mathbb{Z}^{2d}$. Hence, our results may not be easily recovered from [37, Theorem 6.4].

**Theorem 2.2** (Approximability for Gaussian state right-hand sides). Let $\varepsilon > 0$ and $0 \le \alpha < \delta_0/2$. For $\mathbf{h} \in \mathcal{H}$, consider the right-hand side

$$ f_{\mathbf{h}} := \sum_{[m,n] \in \Lambda_{\mathbf{h},\text{rhs}}} F_{m,n}^{\mathbf{h}} \Psi_{\mathbf{h},m,n}. $$

where $F^{\mathbf{h}} \in \Lambda_{\mathbf{h},\text{rhs}}$ with

$$ (2.6) \qquad \Lambda_{\mathbf{h},\text{rhs}} := \{[m, n] \in \mathbb{Z}^{2d} \mid |p_{\mathbf{h}}(\mathbf{x}^{\mathbf{h},m}, \boldsymbol{\xi}^{\mathbf{h},n})| \le \alpha + 2\mathbf{h}^{1/2-\varepsilon}\}. $$

Then, if $u_{\mathbf{h}}$ is the solution to $P_{\mathbf{h}}u_{\mathbf{h}} = f_{\mathbf{h}}$, we have

$$ \left\| u_{\mathbf{h}} - \sum_{[m,n] \in \Lambda_{\mathbf{h},\text{sol}}} (u_{\mathbf{h}}, \Psi_{\mathbf{h},m,n}^{*}) \Psi_{\mathbf{h},m,n} \right\|_{\hat{H}_{\mathbf{h}}^{p}(\mathbb{R}^{d})} \le C_{\varepsilon,p,m} \mathbf{h}^{-m} \|F_{\mathbf{h}}\|_{\ell^{2}(\mathbb{Z}^{2d})} \quad \forall m \in \mathbb{N}. $$

for all $p \in \mathbb{N}$ with

$$ \Lambda_{\mathbf{h},\text{sol}} := \{[m, n] \in \mathbb{Z}^{2d} \mid |p_{\mathbf{h}}(\mathbf{x}^{\mathbf{h},m}, \boldsymbol{\xi}^{\mathbf{h},n})| \le \alpha + 4\mathbf{h}^{1/2-\varepsilon}\}. $$

In practice, the right-hand side of the problem is not a finite linear combination of Gaussian coherent states. However, many right-hand sides of interest become well-approximated by such combination in the high-frequency regime. This is in particular the case when considering scattering by a plane-wave (see Lemma 4.4 below).

**Corollary 2.3** (Approximability for micro-localised right-hand sides). Let $p \ge 0$. Consider a set of right-hand sides $(f_{\mathbf{h}})_{\mathbf{h} \in \mathcal{H}} \subset \hat{H}^p(\mathbb{R}^d)$ and assume that there exists a set of sequences $(F^{\mathbf{h}})_{\mathbf{h} \in \mathcal{H}} \subset \mathbb{Z}^{2d}$ such that

$$ \begin{gathered} \|F^{\mathbf{h}}\|_{\ell^2(\mathbb{Z}^{2d})} \le C \\ \|f_{\mathbf{h}} - \sum_{[m,n] \in \Lambda_{\mathbf{h},\text{rhs}}} F_{m,n}^{\mathbf{h}} \Psi_{\mathbf{h},m,n}\|_{\hat{H}_{\mathbf{h}}^{p}(\mathbb{R}^{d})} \le C_m \mathbf{h}^{-m} \quad \forall m \in \mathbb{N} \end{gathered} $$

for all $\mathbf{h} \in \mathcal{H}$. Then, we have

$$ \left\| u_{\mathbf{h}} - \sum_{[m,n] \in \Lambda_{\mathbf{h},\text{sol}}} (u_{\mathbf{h}}, \Psi_{\mathbf{h},m,n}^{*}) \Psi_{k,m,n} \right\|_{\hat{H}_{\mathbf{h}}^{p+2}(\mathbb{R}^{d})} \le C_{\varepsilon,p,m} \mathbf{h}^{-m} \quad \forall m \in \mathbb{N}. $$

### 3. PROOF OF THEOREM 2.2

This section is devoted to the detailed proof of Theorem 2.2 and Corollary 2.3.
**3.1. Preliminary results on Gaussian states.** We start by stating key preliminary results on the Gaussian and dual frames. We first point out that the converse to (2.1) is true, namely, that

$$ (3.1) \qquad v = \sum_{[m,n] \in \mathbb{Z}^{2d}} (v, \Psi_{h,m,n}^* \Psi_{h,m,n}) \quad \forall v \in L^2(\mathbb{R}^d), $$

see, e.g. [1, 4, 7]. We also record that the following bound

$$ (3.2) \qquad \|\Psi_{h,m,n}\|_{\hat{H}_h^s(\mathbb{R}^d)} \le C_s(1 + |[m, n]|^s) $$

holds true for all $p \in \mathbb{N}$, see [4, Lemma C.1]. Finally, we will need the following expansion result.

**Proposition 3.1** (Tight expansion). For all $[m, n] \in \mathbb{Z}^{2d}$, there exists a sequence of coefficients $U^{m,n} \subset \ell^2(\mathbb{Z}^{2d})$ such that

$$ \|U^{m,n}\|_{\ell^p(\mathbb{Z}^{2d})} \le C_p \quad \forall p \in [1, +\infty], $$

and for all $\varepsilon > 0$ and $s \in \mathbb{N}$, we have

$$ (3.3) \qquad \left\| \Psi_{h,m,n}^* - \sum_{\substack{[m',n'] \in \mathbb{Z}^{2d} \\ |[m,n]-[m',n']| \le \hbar^{-\varepsilon}}} U_{m',n'}^{m,n} \Psi_{h,m',n'} \right\|_{\hat{H}_h^s(\mathbb{R}^d)} \le C_{\varepsilon,s,m} \hbar^m \quad \forall m \in \mathbb{N}. $$

*Proof.* We start by applying (3.1) to $v = \Psi_{h,m,n}^*$, leading to

$$ \Psi_{h,m,n}^* = \sum_{[m',n'] \in \mathbb{Z}^{2d}} (\Psi_{h,m,n}^*, \Psi_{h,m',n'}^*) \Psi_{k,m,n}. $$

Next, we recall from [4, Proposition 4.2] that

$$ |(\Psi_{h,m,n}^*, \Psi_{h,m',n'}^*)| \le C e^{-|[m,n]-[m',n']|^{1/2}}. $$

As a result, we define $U_{m',n'}^{m,n} := (\Psi_{h,m,n}^*, \Psi_{h,m',n'}^*)$, so that $U^{m,n}$ indeed belongs to $\ell^p(\mathbb{Z}^{2d})$ for $1 \le p \le +\infty$, and

$$ \mathcal{E} := \Psi_{h,m,n}^* - \sum_{\substack{[m',n'] \in \mathbb{Z}^{2d} \\ |[m,n]-[m',n']| \le \hbar^{-\varepsilon}}} U_{m',n'}^{m,n} \Psi_{h,m',n'} = \sum_{\substack{[m',n'] \in \mathbb{Z}^{2d} \\ |[m,n]-[m',n']| > \hbar^{-\varepsilon}}} U_{m',n'}^{m,n} \Psi_{h,m',n'}. $$

We then observe that

$$ \begin{aligned} \|\mathcal{E}\|_{\hat{H}_h^s(\mathbb{R}^d)} &\le \sum_{\substack{[m',n'] \in \mathbb{Z}^{2d} \\ |[m,n]-[m',n']| > \hbar^{-\varepsilon}}} |U_{m',n'}^{m,n}| \cdot \|\Psi_{h,m',n'}\|_{\hat{H}_h^s(\mathbb{R}^d)} \\ &\le C \sum_{\substack{[m',n'] \in \mathbb{Z}^{2d} \\ |[m,n]-[m',n']| > \hbar^{-\varepsilon}}} (1 + |[m,n]|)^s e^{-|[m,n]-[m',n']|^{1/2}}, \end{aligned} $$

due to (3.2), and (3.3) follows since

$$ \sum_{\substack{[m',n'] \in \mathbb{Z}^{2d} \\ |[m,n]-[m',n']| > \hbar^{-\varepsilon}}} (1 + |[m,n]|)^s e^{-|[m,n]-[m',n']|^{1/2}} \le C_{\varepsilon,s,m} \hbar^{-m}. $$

□

We close this section with two technical results. As we believe they are of independent interest, and because their proof require tedious computations, they are presented later in Appendix A.

**Proposition 3.2** (Quasi orthogonality). Consider two sets of indices $\Lambda, \Lambda' \subset \mathbb{Z}^{2d}$ with

$$ \mu := \mathrm{diam}(\Lambda) < +\infty \qquad \rho := \mathrm{dist}(\Lambda, \Lambda') > 0. $$
Consider $L \in \mathbb{N}$, smooth coefficients $(A_a)_{a \in \mathbb{N}^d} \subset C_b^\infty(\mathbb{R}^d)$ and the differential operator

$$ P_{h,L,A} := \sum_{\substack{\alpha \in \mathbb{N}^d \\ [\alpha] \le L}} \hbar^{[\alpha]} \Lambda_\alpha \partial^\alpha. $$

Then, for all $q > 0$ and $\alpha \in \mathbb{N}^d$, we have

$$ (3.4) \quad \sum_{[m,n] \in \Lambda} \sum_{[m',n'] \in \Lambda'} |(x^\alpha P_{h,L,A} \Psi_{h,m,n}, \Psi_{h,m',n'})|^q \le C_{L,A,q,m} (1 + (\hbar\mu)^{[\alpha]+L})^q (1 + \rho)^{-m} $$

for all $m \in \mathbb{N}$.

**Proposition 3.3** (Control of $(P_h - p)^L$). For all $[m, n] \in \mathbb{Z}^{2d}$, we have

$$ (3.5a) \quad \left\| \left( P_h - p_h(x^{k,m}, \xi^{k,n}) \right)^L \Psi_{h,m,n} \right\|_{L^2(\mathbb{R}^d)} \le C_L (1 + |n|^2)^L \hbar^{L/2}, $$

and

$$ (3.5b) \quad \left\| \left( P_h^* - \bar{p}_h(x^{k,m}, \xi^{k,n}) \right)^L \Psi_{h,m,n} \right\|_{L^2(\mathbb{R}^d)} \le C_L (1 + |n|^2)^L \hbar^{L/2}. $$

**3.2. Main proof.** We focus here on the proof of Theorem 2.2. Hence, we fix an $\varepsilon > 0$, and consider a right-hand side micro-localised near $\{p_h = 0\}$. Specifically, we will assume that

$$ f_h := \sum_{[m,n] \in \Lambda_{h,rhs}} F_{m,n}^{\hbar} \Psi_{h,m,n}, $$

where $F^\hbar \in \ell^2(\Lambda_{rhs})$ with $\Lambda_{h,rhs} := \Lambda(0, \alpha + 2\hbar^{1/2-\varepsilon})$. Our goal is to show that the associated solution $u_h$ is essentially micro-localised near $\{p_h = 0\}$ as well. Specifically, setting $\Lambda_{h,near} := \Lambda(0, \alpha + 4\hbar^{1/2-\varepsilon})$, our goal will be to show that

$$ u_h^{\text{near}} := \sum_{[m,n] \in \Lambda_{h,near}} (u_h, \Psi_{h,m,n}^*) \Psi_{h,m,n} $$

is “close” to $u_h$.

The key idea is to separate the set of indices of $u_h$ into $\Lambda_{h,near}$, $\Lambda_{h mid} := \Lambda(\alpha + 4\hbar^{1/2-\varepsilon}, \alpha + 6\hbar^{1/2-\varepsilon})$ and $\Lambda_{h,far} := \Lambda(\alpha + 6\hbar^{1/2-\varepsilon}, +\infty)$. We shall also need the “enlarged” sets

$$ \Lambda_{h,near}^* := \Lambda(0, \alpha + 5\hbar^{1/2-\varepsilon}), \quad \Lambda_{h mid}^* := \Lambda(\alpha + 3\hbar^{1/2-\varepsilon}, \alpha + 7\hbar^{1/2-\varepsilon}), \quad \Lambda_{h,far}^* := \Lambda(\alpha + 5\hbar^{1/2-\varepsilon}, +\infty), $$
for the test functions.

We first state some elementary properties of these sets of indices. We do not report the (straight-forward) proofs for the sake of shortness.

**Lemma 3.4** (Index sets). Assume that $\alpha + 7\hbar^{1/2-\varepsilon} \le \delta_0$, then we have

$$ (3.6) \quad \operatorname{dist}(\Lambda_{h,far}, \Lambda_{h,near}^*) \ge C\hbar^{-\varepsilon} \quad \operatorname{dist}(\Lambda_{h,rhs}, \Lambda_{h mid}^*) \ge C\hbar^{-\varepsilon} \quad \operatorname{dist}(\Lambda_{h,near}, \Lambda_{h,far}^*) \ge C\hbar^{-\varepsilon} $$

and

$$ (3.7) \quad \operatorname{diam}(\Lambda_{h,rhs}^*) \le C\hbar^{-1/2} \quad \operatorname{diam}(\Lambda_{h,near}^*) \le C\hbar^{-1/2} \quad \operatorname{diam}(\Lambda_{h mid}^*) \le C\hbar^{-1/2} $$

In addition, if $\hbar$ is small enough, we have

$$ (3.8) \quad \left\{ [m', n'] \in \mathbb{Z}^{2d} \mid \exists [m, n] \in \Lambda_{h, mid}; \quad |[m, n] - [m', n']| \le \hbar^{-\varepsilon/2} \right\} \subset \Lambda_{h, mid}^*. $$

**Lemma 3.5** (Quasi orthogonality away from RHS micro-support). For $F^\hbar \in \ell^2(\Lambda_{h,rhs})$, consider the right-hand side

$$ f_h = \sum_{[m,n] \in \Lambda_{h,rhs}} F_{m,n}^{\hbar} \Psi_{h,m,n} $$

and the associated solution $u_h$. Then, if $[m, n] \in \Lambda_{h mid}^*$, we have

$$ |(u_h, \Psi_{h,m,n})| \le C_{\varepsilon,m} \hbar^m \|F^\hbar\|_{\ell^2(\mathbb{Z}^{2d})} $$

for all $m \in \mathbb{N}$.
Proof. Throughout the proof, we fix a pair of indices $[m, n] \in \Lambda_{h, \text{mid}}^*$. By definition of $\Lambda_{h, \text{mid}}^*$, the assumption that $\alpha + 7\hbar^{1/2-\epsilon} \le \delta_0$ and (2.5), we have

$$ (3.9) \qquad c\hbar^{1/2-\epsilon} \le |p_h(m, n)| \le C. $$

In particular, we can write

$$ u_h = \frac{1}{p_h(m, n)} f_h + \frac{1}{p_h(m, n)} (P_h - p_h(m, n)) u_h. $$

Since $f_h$ is smooth, so is $u_h$, and we can iterate this relation, leading to

$$ u_h = \sum_{\ell=1}^{r} \frac{1}{p_h(m, n)^\ell} (P_h - p_h(m, n))^{\ell-1} f_h + \frac{1}{p_h(m, n)^r} (P_h - p_h(m, n))^r u_h, $$

and

$$ (3.10) \qquad (u_h, \Psi_{h,m,n}) = \sum_{\ell=1}^{r} \frac{1}{p_h(m, n)^\ell} (f_h, (P_h^* - \bar{p}_h(m, n))^{\ell-1} \Psi_{h,m,n}) \\ \qquad \qquad + \frac{1}{p_h(m, n)^r} (u_h, (P_h^* - \bar{p}_h(m, n))^r \Psi_{h,m,n}) $$

for all $r \in \mathbb{N}$.

Then, if $[m', n'] \in \Lambda_{h, \text{rhs}}$, using the upper-bound in (3.9), we have

$$ \begin{aligned} |(\Psi_{k,m',n'}, (P_h^* - \bar{p}_h(m,n))^{\ell-1} \Psi_{h,m,n})| &= \left| \sum_{j=0}^{\ell-1} \binom{\ell-1}{j} \bar{p}_h(m,n)^{\ell-1-j} (\Psi_{k,m',n'}, (P_h^*)^j \Psi_{h,m,n}) \right| \\ &\le C_\ell \sum_{j=0}^{\ell-1} |(\Psi_{h,m',n'}, (P_h^*)^j \Psi_{h,m,n})|. \end{aligned} $$

Recalling, (3.7), we have $\hbar^{\ell-1} \text{diam}(\Lambda_{h, \text{mid}}^*) \le C\hbar^{\ell-3/2} \le C$ for all $\ell \ge 2$. As a result, applying (3.4) gives

$$ |(\Psi_{h,m',n'}, (P_h^* - \bar{p}_h(m,n))^{\ell-1} \Psi_{h,m,n})| \le C_{\ell,n} (1 + |[m,n] - [m'n']|)^{-n} \le C_{\ell,n} \hbar^{\epsilon n}, $$

for all $n \in \mathbb{N}$, since $|[m, n] - [m', n']| \ge \hbar^{-\epsilon}$ due to (3.6). The case $\ell = 1$ also easily follows by Proposition 3.2. We then write that

$$ \begin{aligned} |(f_h, (P_h^* - \bar{p}_h(m, n))^{\ell-1} \Psi_{h,m,n})| &\le \sum_{[m', n'] \in \Lambda_{h, \text{rhs}}} |F_{m,n}^h| |(\Psi_{h,m',n'}, (P_h^* - \bar{p}_h(m, n))^{\ell-1} \Psi_{h,m,n})| \\ &\le C_{\ell,n} \hbar^{\epsilon n} \sum_{[m', n'] \in \Lambda_{h, \text{rhs}}} |F_{m,n}^h| \\ &\le C_{\ell,n} \hbar^{\epsilon n-d} \|F^h\|_{\ell^2(\mathbb{Z}^{2d})}, \end{aligned} $$

where we used that $|\Lambda_{h, \text{rhs}}| \le C \text{diam}(\Lambda_{h, \text{rhs}})^{2d} \le C\hbar^d$ due to (3.7). As a result, using the lower-bound in (3.9), we have

$$ \frac{1}{|p(m, n)|^\ell} |(f_h, (P_h^* - \bar{p}_h(m, n))^{\ell-1} \Psi_{h,m,n})| \le C_{\ell,n} \hbar^{\epsilon n - d - \ell/2} \|F^h\|_{\ell^2(\mathbb{Z}^{2d})}, $$

and

$$ \sum_{\ell=1}^{r} \frac{1}{|p(m, n)|^\ell} |(f_h, (P_h^* - \bar{p}_h(m, n))^{\ell-1} \Psi_{h,m,n})| \le C_{m,n} \hbar^{\epsilon n - d - r/2} \|F^h\|_{\ell^2(\mathbb{Z}^{2d})}, $$

for all $n \in \mathbb{N}$. Thus, for any $m \in \mathbb{N}$, we can select $n = n(m, d, r, \epsilon)$ such that $\epsilon n - d - r/2 \ge m$, leading to

$$ \sum_{\ell=1}^{r} \frac{1}{|p(m, n)|^\ell} |(f_h, (P_h^* - \bar{p}_h(m, n))^{\ell-1} \Psi_{h,m,n})| \le C_{\epsilon,m,r} \hbar^m \|F^h\|_{\ell^2(\mathbb{Z}^{2d})}. $$
On the other hand, using again (3.7) and applying (3.5b), we have

$$ \|(P_h^* - \bar{p}(m, n))^T \Psi_{h,m,n} \|_{L^2(\mathbb{R}^d)} \le C_r \hbar^{r/2}, $$

and the lower bound in (3.9) shows that

$$ \frac{1}{|p(m, n)|^r} \|(P_h^* - \bar{p}(m, n))^T \Psi_{h,m,n} \|_{L^2(\mathbb{R}^d)} \le C_r \hbar^{\varepsilon r}. $$

We then write that

$$ \frac{1}{|p(m, n)|^m} |(u_h, (P_h^* - \bar{p}(m, n))^T \Psi_{h,m,n})| \le C_r \hbar^{\varepsilon r} \|u_h\|_{L^2(\mathbb{R}^d)} \le C_r \hbar^{\varepsilon r - N} \|f_h\|_{L^2(\mathbb{R}^d)} \le C_r \hbar^{\varepsilon r - N} \|F^h\|_{\ell^2(\mathbb{Z}^{2d})}, \le C_{\varepsilon, m} \hbar^m \|F^h\|_{\ell^2(\mathbb{Z}^{2d})}, $$

up to picking $r$ such that $\varepsilon r - N \ge m$. $\square$

*Proof of Theorem 2.2.* We expand $u_h$ in the frame $(\Psi_{h,m,n})_{[m,n]\in\mathbb{Z}^{2d}}$ as

$$ u_h = u_h^{\text{near}} + u_h^{\text{mid}} + u_h^{\text{far}} $$

where

$$ u_h^{\text{near}} := \sum_{[m,n] \in \Lambda_{h,\text{near}}} (u_h, \Psi_{h,m,n}^*) \Psi_{h,m,n} $$

$$ u_h^{\text{mid}} := \sum_{[m,n] \in \Lambda_{h,\text{mid}}} (u_h, \Psi_{h,m,n}^*) \Psi_{h,m,n} $$

$$ u_h^{\text{far}} := \sum_{[m,n] \in \Lambda_{h,\text{far}}} (u_h, \Psi_{h,m,n}^*) \Psi_{h,m,n}. $$

The proof then consists in showing that $u_h^{\text{mid}}$ and $u_h^{\text{far}}$ are small.

**Step 1.** We first treat the $u_h^{\text{mid}}$ term. To do so, we start by introducing the approximation

$$ (3.11) \qquad \tilde{u}_h^{\text{mid}} := \sum_{[m,n] \in \Lambda_{h,\text{mid}}} \sum_{[m',n'] \in \mathbb{Z}^{2d} \atop |[m,n]-[m',n']| \le \hbar^{-\varepsilon/2}} U_{m',n'}^{m,n} (u_h, \Psi_{h,m',n'}) \Psi_{h,m,n}. $$

Recalling (3.8), all the $[m', n']$ indices in the sum belong to the enlarged set $\Lambda_{h,\text{mid}}^*$, so that

$$ (3.12) \qquad |(u_h, \Psi_{h,m',n'})| \le C_{\varepsilon,m} \hbar^m \|F^h\|_{\ell^2(\mathbb{Z}^{2d})} $$

by Lemma 3.5. Recalling (3.7), $|[m, n]| \le C\hbar^{-1/2}$ for all $[m, n] \in \Lambda_{h,\text{mid}}$, and we have from (3.2)

$$ (3.13) \qquad \| \Psi_{h,m,n} \|_{\hat{H}_h^p(\mathbb{R}^d)} \le C_p. $$

Thus, plugging (3.12) and (3.13) into (3.11), we have

$$ \|\tilde{u}_h^{\text{mid}}\|_{L^2(\mathbb{R}^d)} \le C_{p,\varepsilon,m} \hbar^m \|F^h\|_{\ell^2(\mathbb{Z}^{2d})} \sum_{[m,n] \in \Lambda_{h,\text{mid}}} \sum_{[m',n'] \in \mathbb{Z}^{2d} \atop |[m,n]-[m',n']| \le \hbar^{-\varepsilon/2}} 1 $$

$$ (3.14) \qquad \le C_{p,\varepsilon,m} \hbar^{m-\varepsilon d^2/2} \|F^h\|_{\ell^2(\mathbb{Z}^{2d})}. $$

We now estimate the difference between $u_h^{\text{mid}}$ and $\tilde{u}_h^{\text{mid}}$

$$ u_h^{\text{mid}} - \tilde{u}_h^{\text{mid}} = \sum_{[m,n] \in \Lambda_{h,\text{mid}}} \left( u_h, \Psi_{h,m,n}^* - \sum_{[m',n'] \in \mathbb{Z}^{2d} \atop |[m,n]-[m',n']| \le \hbar^{-\varepsilon/2}} U_{m',n'}^{m,n} \Psi_{h,m',n'} \right) \Psi_{h,m,n}, $$
so that

$$ \|u_h^{\text{mid}} - \tilde{u}_h^{\text{mid}}\|_{\hat{H}_h^p(\mathbb{R}^d)} \le C_p \|u_h\|_{L^2(\mathbb{R}^d)} \sum_{[m,n] \in \Lambda_{h,\text{mid}}} \left\| \Psi_{h,m,n}^* - \sum_{[m',n'] \in \mathbb{Z}^{2d} \setminus [m,n]-[m',n'] \le \hbar^{-\epsilon/2}} \Psi_{h,m',n'} \right\|_{L^2(\mathbb{R}^d)} $$

and it follows from Proposition 3.1 that

$$ (3.15) \qquad \|u_h^{\text{mid}} - \tilde{u}_h^{\text{mid}}\|_{\hat{H}_h^p(\mathbb{R}^d)} \le C_{p,\epsilon,m} \hbar^{m-d} \|u_h\|_{\hat{H}_h^p(\mathbb{R}^d)} \le C_{p,\epsilon,m} \hbar^{m-d-N} \|F^\hbar\|_{\ell^2(\mathbb{Z}^{2d})}. $$

Then, it follows from (3.14) and (3.15) that

$$ \|u_h^{\text{mid}}\|_{\hat{H}_h^p(\mathbb{R}^d)} \le C_{p,\epsilon,m} \hbar^m \|F^\hbar\|_{\ell^2(\mathbb{Z}^{2d})} \quad \forall m \in \mathbb{N}, $$

up to redefining $m$.

**Step 2.** We then turn our attention to $u_h^{\text{far}}$. On the one hand, we can apply Proposition 3.2 with $P_{\hbar,L,A} := \hbar^{[\beta]}\partial^\beta \circ P_\hbar$. This gives

$$ \sum_{[m,n] \in \Lambda_{h,\text{near}}^*} |(x^\alpha \partial^\beta (P_\hbar u_\hbar^{\text{far}}), \Psi_{\hbar,m,n})| = \hbar^{-[\beta]} \left| \sum_{[m,n] \in \Lambda_{h,\text{near}}^*} \sum_{[m',n'] \in \Lambda_{h,\text{far}}} (u_\hbar, \Psi_{\hbar,m',n'}^*(x^\alpha P_{h,L,A} \Psi_{\hbar,m',n'}, \Psi_{\hbar,m,n})) \right| \\ \le C_{\alpha,\beta,r} \hbar^{r-[\beta]} \|u_\hbar\|_{L^2(\mathbb{R}^d)} \quad \forall r \in \mathbb{N}. $$

As a result, since $|\Lambda_{h,\text{near}}^*| \le C(\text{diam}(\Lambda_{h,\text{near}}^*))^{2d} \le C\hbar^{-d}$ due to (3.7), we have

$$ \sum_{[m,n] \in \Lambda_{h,\text{near}}^*} |(x^\alpha \partial^\beta (P_h u_h^{\text{far}}), \Psi_{h,m,n})|^2 \le C\hbar^{-d} \left( \sum_{[m,n] \in \Lambda_{h,\text{near}}^*} |(x^\alpha \partial^\beta (P_h u_h^{\text{far}}), \Psi_{h,m,n})| \right)^2 \\ \le C_{\alpha,\beta,r} \hbar^{2r-2[\beta]-d} \|u_h\|_{L^2(\mathbb{R}^d)}^2 $$

which we rewrite as

$$ \sum_{[m,n] \in \Lambda_{h,\text{near}}^*} |(x^\alpha \partial^\beta (P_h u_h^{\text{far}}), \Psi_{h,m,n})|^2 \le C_{\alpha,\beta,m} \hbar^{2m} \|u_h\|_{L^2(\mathbb{R}^d)}^2 \quad \forall m \in \mathbb{N} $$

after changing variables.

On the other hand, if $[m,n] \in \Lambda_{h,\text{far}}^*$, we write that

$$ P_h u_h^{\text{far}} = f_h - P_h u_h^{\text{near}} + P_h u_h^{\text{mid}}, $$

so that

$$ |(x^\alpha \partial^\beta (P_h u_h^{\text{far}}), \Psi_{h,m,n})|^2 \le C \left( |(x^\alpha \partial^\beta f_h, \Psi_{h,m,n})|^2 + |(x^\alpha \partial^\beta (P_h u_h^{\text{near}}), \Psi_{h,m,n})|^2 + |(x^\alpha \partial^\beta (P_h u_h^{\text{mid}}), \Psi_{h,m,n})|^2 \right). $$

We then have

$$ \sum_{[m,n] \in \Lambda_{h,\text{far}}^*} |(x^\alpha \partial^\beta (P_h u_h^{\text{mid}}), \Psi_{h,m,n})|^2 \le C \|x^\alpha \partial^\beta (P_h u_h^{\text{mid}})\|_{L^2(\mathbb{R}^d)}^2 \le C \|u_h^{\text{mid}}\|_{\hat{H}_h^{[\beta]+2}(\mathbb{R}^d)}^2 \le C_{\beta,m} \hbar^m \|F^\hbar\|_{\ell^2(\mathbb{Z}^{2d})}^2 $$

due to **Step 1**. Next, thanks to Proposition 3.2,

$$ \sum_{[m,n] \in \Lambda_{h,\text{far}}^*} |(x^\alpha \partial^\beta f_h, \Psi_{h,m,n})|^2 \le \sum_{[m,n] \in \Lambda_{h,\text{far}}^*} \sum_{[m',n'] \in \Lambda_{h,\text{rhs}}^*} |F_{m,n}^\hbar|^2 |(x^\alpha \partial^\beta \Psi_{h,m',n'}, \Psi_{h,m,n})|^2 \\ \le C_{\alpha,\beta,m}^2 \hbar^{2m} \|F^\hbar\|_{\ell^2(\mathbb{Z}^{2d})}^2. $$

Finally, still by Proposition 3.2, it holds that

$$ \sum_{[m,n] \in \Lambda_{h,\text{far}}^*} |(x^\alpha \partial^\beta (P_h u_h^{\text{near}}), \Psi_{h,m,n})|^2 \le \sum_{[m,n] \in \Lambda_{h,\text{far}}^*} \sum_{[m',n'] \in \Lambda_{h,\text{near}}^*} |(u, \Psi_{h,m',n'}^*)|(x^\alpha \partial^\beta (P_h \Psi_{h,m',n'}, \Psi_{h,m,n})|^2 \\ \le C_{\alpha,\beta,r} \hbar^{2r} \|u_h\|_{L^2(\mathbb{R}^d)}^2 \le C_{\alpha,\beta,m} \hbar^{2m} \|F^\hbar\|_{\ell^2(\mathbb{Z}^{2d})}^2. $$
Since $\Lambda_{\hbar, \text{near}}^*$ and $\Lambda_{\hbar, \text{far}}^*$ form a partition of $\mathbb{Z}^{2d}$, we have thus shown that

$$ (3.16) \qquad \|\mathbf{x}^{\alpha} \partial^{\beta} (P_{\hbar} u_{\hbar}^{\text{far}})\|_{L^2(\mathbb{R}^d)} \le C_{\alpha, \beta, m} \hbar^m \|F^{\hbar}\|_{\ell^2(\mathbb{Z}^{2d})}. $$

Letting $f_{\hbar}^{\text{far}} = P_{\hbar} u_{\hbar}^{\text{far}}$, we see from (3.16) that $u_{\hbar}^{\text{far}}$ solves $P_{\hbar} u_{\hbar}^{\text{far}} = f_{\hbar}^{\text{far}}$ with a right-hand side $f_{\hbar}^{\text{far}} \in \hat{H}_{\hbar}^{p}(\mathbb{R}^d)$ such that

$$ \|f_{\hbar}^{\text{far}}\|_{\hat{H}_{\hbar}^{p}(\mathbb{R}^d)} \le C_{p,m} \hbar^m \|F^{\hbar}\|_{\ell^2(\mathbb{Z}^{2d})} \quad \forall m \in \mathbb{N}. $$

Then, we conclude the proof with (2.4). □

### 4. APPLICATION TO THE HELMHOLTZ EQUATION

We now turn our attention to the model problem of the Helmholtz equation (1.1), with a particular focus on the case of plane-wave scattering.

4.1. **Notation.** In this section $k$ will denote the wavenumber in the Helmholtz problem, For the sake of simplicity, we assume that $kR \ge 1$. We will apply the results of Section 2 with $\hbar \sim (kR)^{-1}$. As a result, the norms

$$ (4.1) \qquad \|v\|_{\hat{H}_k^p(\mathbb{R}^d)} := \sum_{\substack{\alpha \in \mathbb{N}^d \\ [\alpha] \le p}} \sum_{q \le p - [\alpha]} k^{-2[\alpha]} \left\| \left|\frac{\mathbf{x}}{R}\right|^q \partial^{\alpha} v \right\|_{L^2(\mathbb{R}^d)}^2 $$

will be convenient. Notice that if

$$ \mathcal{F}(v)(\xi) := (2\pi)^{-d/2} \int_{\mathbb{R}^d} v(\mathbf{x}) e^{-i\mathbf{x} \cdot \xi} d\mathbf{x}, \quad \text{a.e. } \xi \in \mathbb{R}^d. $$

is the Fourier transform of $v \in L^2(\mathbb{R}^d)$, and if we define the “reverse” norm by

$$ (4.2) \qquad \|v\|_{\hat{H}_k^p(\mathbb{R}^d)}^2 := \sum_{\substack{\alpha \in \mathbb{N}^d \\ [\alpha] \le p}} \sum_{q \le p - [\alpha]} k^{-2q} \left\| \left|\frac{\mathbf{x}}{R}\right|^q \partial^{\alpha} v \right\|_{L^2(\mathbb{R}^d)}^2, $$

then we have

$$ (4.3) \qquad c(R) \|\mathcal{F}v\|_{\hat{H}_k^p(\mathbb{R}^d)} \le \|v\|_{\hat{H}_k^p(\mathbb{R}^d)} \le C(R) \|\mathcal{F}v\|_{\hat{H}_k^p(\mathbb{R}^d)}. $$

We will also use the following “standard” norm

$$ \|v\|_{\hat{H}_k^p(\mathbb{R}^d)}^2 := \sum_{\substack{\alpha \in \mathbb{N}^d \\ [\alpha] \le p}} k^{-2[\alpha]} \|\partial^{\alpha} v\|_{L^2(\mathbb{R}^d)}^2. $$

4.2. **Model problem.** We consider smooth coefficients $\mu$, $A \in C_b^\infty(\mathbb{R}^d)$ that are respectively equal to 1 and $I$ outside $B(0, R)$. Given $f: \Omega \to \mathbb{C}$ our model problem is to find $u: \Omega \to \mathbb{C}$ such that

$$ (4.4a) \qquad -k^2 \mu u - \nabla \cdot (A \nabla u) = f \quad \text{in } \mathbb{R}^d $$

and

$$ (4.4b) \qquad \frac{\partial u}{\partial |\mathbf{x}|}(\mathbf{x}) - iku(\mathbf{x}) = o\left(|\mathbf{x}|^{(-d-1)/2}\right) \quad \text{as } |\mathbf{x}| \to +\infty. $$

Problem (4.4) is well-posed in the sense that for all $f \in L_{\text{comp}}^2(\mathbb{R}^d)$, there exists a unique $R_k f := u \in L_{\text{loc}}^2(\mathbb{R}^d)$ such that (4.4a) and (4.4b) hold true. Here, we will further assume that the problem has polynomially bounded resolvent, meaning that there exists $C, N$ such that

$$ (4.5) \qquad \|\chi R_k \chi\|_{L^2(\mathbb{R}^d) \to L^2(\mathbb{R}^d)} \le C(kR)^N $$

where $\chi$ is a smooth cut-off function that takes the value 1 on $B(0, R)$ and 0 outside $B(0, 2R)$.

**Remark 4.1** (When does the polynomial bound actually hold?). The bound (4.5) is known to hold in several situations:
* When the dynamics induced by the Hamiltonian p has no trapped trajectory, i.e., when every trajectory leaves any compact set in finite time, the assumption holds with $K = [1, +\infty)$. See for instance [16]. This situation is often referred to as "non-trapping".

* When the dynamics induced by p has a trapped set, and the dynamics is hyperbolic, close to this trapped set, it has been conjectured in [38] that (2.4) always holds with $K = [1, +\infty)$. Actually, this is already known when the trapped set is "filamentary enough", see [32, 33].

* Without any assumption on the dynamics, (2.4) holds taking $K$ to be $[1, +\infty)$ from which we exclude a set of frequencies $k$ whose intersection with $\{k \ge k_0\}$ has a length going to zero as $k_0 \to +\infty$. We refer the reader to [25] for more details.

### 4.3. Perfectly matched layers.
As advertised in the introduction, the formulation (4.4) is not suited for immediate discretization by "volume" methods, as the radiation condition is hard to take into account. We will thus rely on an equivalent formulation that uses perfectly matched layers (PML).

Specifically, given $f: \Omega \to \mathbb{C}$, we consider the problem to find $u: \Omega \to \mathbb{C}$ such that $P_k u = f$ where

$$ (4.6) \qquad P_k u := -\frac{1}{k^2} ((I + iM)^{-1}\nabla) \cdot ((I + iM)^{-1}A\nabla u) - \mu u. $$

In (4.6), the (SPD) matrix function $M$ is given by

$$ M(\mathbf{x}) := \frac{g(|\mathbf{x}|)}{|\mathbf{x}|^3} (|\mathbf{x}|^2 \mathrm{Id} - \mathbf{x} \otimes \mathbf{x}) + \frac{g'(|\mathbf{x}|)}{|\mathbf{x}|^2} \mathbf{x} \otimes \mathbf{x}, $$

where $g: \mathbb{R} \to \mathbb{R}$ is a user-defined function such that $g(r) = 0$ if $r \le R$, $g(r) = r$ if $r \ge R_0 > R$, and $g'(r) \ge 0$. In what follows, we will assume that $g$ a smooth function to satisfy the assumptions of Section 2, but many results about PML still hold with less regular $g$ (see, e.g., [14]).

Notice that, if $|\mathbf{x}| \le R$, $M = 0$, so that the original operator is not modified on the support of $\mu - 1$ and $A - I$. On the other hand, $M = I$ if $|\mathbf{x}| \ge R_0$, so that dissipation is introduced away from the origin, where the operator simply reads

$$ P_k v = -\frac{2i}{k^2} \Delta v - \mu v \quad \text{whenever} \quad \mathrm{supp} v \cap B(0, 2R) = \emptyset. $$

This transformation can be naturally interpreted as a complex deformation of coordinates. It is also often called the complex scaling technique. Crucially, the PML is designed in such way that

$$ (4.7) \qquad ((P_k)^{-1} f)|_{B(0,R)} = (R_k f)|_{B(0,R)} $$

whenever $\mathrm{supp} f \subset B(0, R)$. We refer the reader to [9, §4.5] or [14] for more information.

### 4.4. Abstract setting.
We now verify that the Helmholtz problem formulated with PML indeed fits the abstract setting of Section 2.3. The only non-trivial facts to establish are the polynomial resolvent estimates in $\hat{H}^q(\mathbb{R}^d)$ in (2.4) and the boundedness of the energy layer (2.5).

**Lemma 4.2 (Resolvent estimates).** Let $q \in \mathbb{N}$. For all $f \in \hat{H}^q(\mathbb{R}^d)$, there exists a unique $u \in L^2(\mathbb{R}^d)$ such that $P_k u = f$. In addition we have $u \in \hat{H}^{q+2}(\mathbb{R}^d)$ with

$$ \|u\|_{\hat{H}_k^{q}(\mathbb{R}^d)} \le C(kR)^N \|f\|_{\hat{H}_k^{q}(\mathbb{R}^d)}. $$

Furthermore, if $\mathrm{supp} f \subset B(0, R)$, then $u \in \hat{H}^{q+2}(\mathbb{R}^d)$ with

$$ (4.8) \qquad \|u\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)} \le C(kR)^N \|f\|_{\hat{H}_k^q(\mathbb{R}^d)}. $$

*Proof.* We first invoke Theorem 1.6 of [14] which states that

$$ \|P_k^{-1}\|_{L^2(\mathbb{R}^d) \to H_k^2(\mathbb{R}^d)} \le C \|\chi R_k \chi\|_{L^2(\mathbb{R}^d) \to L^2(\mathbb{R}^d)}. $$

Then, using a usual bootstrap technique, we easily show that

$$ (4.9) \qquad \|P_k^{-1}\|_{H_k^q(\mathbb{R}^d) \to H_k^{q+2}(\mathbb{R}^d)} \le C \|\chi R_k \chi\|_{L^2(\mathbb{R}^d) \to L^2(\mathbb{R}^d)} $$

for all $q \in \mathbb{N}$.
We then need to take care of the weights in the $\hat{H}^q(\mathbb{R}^d)$ norms. To do so, we observe that if $u = P_k f$, then we may write

$$2ik^2u + \Delta u = g,$$

with

$$g := -2k^2P_k u + (2k^2P_k u + \Delta u + 2ik^2)u = -2k^2(f + Q_k u),$$

where $Q_k$ is differential operator of degree 2 with smooth coefficients supported in $B(0, R_0)$. Let $\alpha, \beta \in \mathbb{N}^d$ with $|\alpha|, |\beta| \le 2$. Since $Q_k$ is compactly supported, we have

$$(4.10) \qquad \left\| \left(\frac{\mathbf{x}}{R}\right)^{\alpha} g \right\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)} \le C_{\alpha} k^2 \left( \|u\|_{H_k^{q+2}(\mathbb{R}^d)} + \|f\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)} \right).$$

Now, we note that

$$(4.11) \qquad \left(\frac{\mathbf{x}^{\alpha}}{R}\right) \partial^{\beta} u = i^{\alpha-\beta} \mathcal{F}^{-1} \left( \left(\frac{1}{R}\partial\right)^{\alpha} \left( \frac{\xi^{\beta}}{-|\xi|^2 + 2ik^2} \mathcal{F}(g) \right) \right),$$

Let us write

$$\frac{\xi^{\beta}}{-|\xi|^2 + 2ik^2} = k^{[\beta]-2} \frac{(\xi/k)^{\beta}}{-|\xi/k|^2 + 2i},$$

so that the map

$$\xi \mapsto \frac{\xi^{\beta}}{-|\xi|^2 + 2ik^2}$$

has $C^\ell$ norm bounded by $C_\ell k^{[\beta]-2-\ell} \le C_\ell k^{[\beta]-2}$, since $k \ge 1$.

We deduce from (4.3) and (4.10) that

$$\left\| \left(\frac{1}{R}\partial\right)^{\alpha} \left( \frac{\xi^{\beta}}{-|\xi|^2 + 2ik^2} \mathcal{F}(g) \right) \right\|_{\hat{H}_k^q(\mathbb{R}^d)} \le C k^{[\beta]} \left( \|u\|_{H_k^{q+2}(\mathbb{R}^d)} + \|f\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)} \right)$$

and hence, thanks to (4.11),

$$\|u\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)} \le C \sum_{[\alpha]+[\beta] \le 2} k^{-[\beta]} \left\| \left(\frac{\mathbf{x}}{R}\right)^{\alpha} \partial^{\beta} u \right\|_{\hat{H}_k^p(\mathbb{R}^d)} \le C \sum_{[\alpha]+[\beta] \le 2} k^{-[\beta]} \left\| \mathcal{F} \left( \left(\frac{\mathbf{x}}{R}\right)^{\alpha} \partial^{\beta} u \right) \right\|_{\hat{H}_k^p(\mathbb{R}^d)} \\ \le C \sum_{[\alpha]+[\beta] \le 2} \left( \|u\|_{H_k^{q+2}(\mathbb{R}^d)} + \|f\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)} \right) \le C \|\chi R_k \chi\|_{L^2(\mathbb{R}^d) \to L^2(\mathbb{R}^d)} \|f\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)},$$

thanks to (4.9). The first part of the result follows from (4.5).

Now, if we further assume that $f$ is supported in $B(0, R)$, $g$ is compactly supported and equation (4.10) may be replaced by

$$\left\| \left(\frac{\mathbf{x}}{R}\right)^{\alpha} g \right\|_{\hat{H}_k^q(\mathbb{R}^d)} \le C_{\alpha} k^2 \left( \|f\|_{H_k^q(\mathbb{R}^d)} + \|u\|_{H_k^{q+2}(\mathbb{R}^d)} \right),$$

and the same reasoning as above leads to

$$\|u\|_{\hat{H}_k^{q+2}(\mathbb{R}^d)} \le C \left( \|u\|_{H_k^{q+2}(\mathbb{R}^d)} + \|f\|_{H_k^q(\mathbb{R}^d)} \right) \le C \|\chi R_k \chi\|_{L^2(\mathbb{R}^d) \to L^2(\mathbb{R}^d)} \|f\|_{\hat{H}_k^q(\mathbb{R}^d)},$$

as announced.

□

**Lemma 4.3** (Boundedness of the energy layer). For all $\delta \in (0, 1)$, the set

$$\{[\mathbf{x}, \xi] \in \mathbb{R}^{2d} \mid |p(\mathbf{x}, \xi)| \le \delta\}$$

is bounded.

*Proof.* Fix $0 < \delta < 1$, and let $U := \{[\mathbf{x}, \xi] \in \mathbb{R}^{2d}; |p(\mathbf{x}, \xi)| \le \delta\}$. Let $\mathbf{x} \in \mathbb{R}^d$. We first assume that $|\mathbf{x}| \le R_0$. Since we know that

$$C |\xi|^2 - C' \le |p(\mathbf{x}, \xi)| \le \delta,$$

we see that $\{\xi \in \mathbb{R}^d \mid |p(\mathbf{x}, \xi)| \le \delta\}$ is bounded, and thus $U \cap (B_\mathbf{x}(0, R_0) \times \mathbb{R}^d)$ is bounded. On the other hand, if $|\mathbf{x}| \ge R_0$, we have

$$p(\mathbf{x}, \xi) = (1+i)^2 |\xi|^2 - 1 = 2i|\xi|^2 - 1$$
so that

$$|p(\mathbf{x}, \boldsymbol{\xi})|^2 = 4|\boldsymbol{\xi}|^4 + 1 \ge 1 > \delta,$$

which implies that $U \setminus (B_\mathbf{x}(0, R_0) \times \mathbb{R}^d) = \emptyset$.

<figure><img src="image_2.png" /></figure>

We finally show that the right-hand side associated with plane-wave scattering are indeed well-approximated by Gaussian coherent states in order to apply Corollary 2.3 later on.

**Lemma 4.4** (Approximability of plane-wave right-hand sides). For $k > 0$, consider the right-hand side $f_k := \tilde{\chi}e^{ikd \cdot n}$ where $\tilde{\chi} \in C_c^\infty(B_R)$ and $d \in \mathbb{R}^d$ with $|d| = 1$. Then, there exists $F^k \in \ell^2(\mathbb{Z}^{2d})$ such that

$$(4.12) \qquad \left\| f_k - \sum_{[m,n] \in \Lambda_{k,\text{rhs}}} F_{m,n}^k \Psi_{k,m,n} \right\|_{\hat{H}_k^p(\mathbb{R}^d)} \le C_{\varepsilon,p,m} (kR)^{-m},$$

where

$$\Lambda'_{k, \text{rhs}} := \{ [m, n] \in \mathbb{Z}^{2d} \mid \text{dist}(x^{k,m}, \text{supp}(\chi)) \le (kR)^{-1/2+\varepsilon} \text{ and } |\xi^{k,n}| - 1 \le (kR)^{-1/2+\varepsilon} \}.$$

In particular, if $p$ is the symbol of the function appearing in section 4.3 and if $\chi$ is supported in the region where $p(\mathbf{x}, \boldsymbol{\xi}) = |\boldsymbol{\xi}|^2 - 1$, then $\Lambda'_{k,\text{rhs}}$ is of the form (2.6), with $\alpha = 0$.

*Proof.* From now on, we fix $p, m$ and $\varepsilon$. We will first show that, if $[m, n] \notin \Lambda'_{k,\text{rhs}}$, then we have

$$(4.13) \qquad |(f_k, \Psi_{k,m,n})| \le C_m(kR)^{-m} \quad \forall m \in \mathbb{N}.$$

The quantity $(f_k, \Psi_{k,m,n})$ is of the form $\int_{\mathbb{R}^d} \chi(\mathbf{x})e^{ik\varphi_{m,n}(\mathbf{x})}d\mathbf{x}$, with

$$\varphi_{m,n}(\mathbf{x}) = \mathbf{x} \cdot (\boldsymbol{\xi} - \boldsymbol{\xi}^{k,n}) + \frac{i}{2}|\mathbf{x} - \mathbf{x}^{k,m}|^2,$$

so that

$$\nabla \varphi_{m,n}(\mathbf{x}) = \boldsymbol{\xi} - \boldsymbol{\xi}^{k,n} + i(\mathbf{x} - \mathbf{x}^{m,k}).$$

In particular, we have

$$|\nabla \varphi_{m,n}(\mathbf{x})| \ge |\boldsymbol{\xi}^{k,n}| - 1 + \text{dist}(\mathbf{x}^{k,m}, \text{supp}(\chi)) \ge \frac{1}{2}(kR)^{-\frac{1}{2}+\varepsilon}.$$

we may use the method of non-stationary phase (i.e., integrate by parts several times, as in [37, Lemma 3.14]) to deduce (4.13).

Now, combining (4.13) with Proposition 3.1, we see that there exists $C$ such that, for any $[m, n] \in \Lambda_{k,\text{rhs}}$, we have

$$(4.14) \qquad |(f_k, \Psi_{k,m,n}^*)| \le C_m(kR)^{-m} \quad \forall m \in \mathbb{N}.$$

On the other hand, it follows from [4, Theorem 3.1] that we have

$$(4.15) \qquad \left\| f_k - \sum_{\substack{[m,n] \in \mathbb{Z}^{2d} \\ |m,n| \le k}} (f_k, \Psi_{k,m,n}^*) \Psi_{k,m,n} \right\|_{\hat{H}_k^p(\mathbb{R}^d)} \le C_{p,m} (kR)^{-m}.$$

Writing

$$f_k = \sum_{[m,n] \in \Lambda_{k,\text{rhs}}} (f_k, \Psi_{k,m,n}^*) \Psi_{k,m,n} + \sum_{\substack{[m,n] \in \mathbb{Z}^{2d} \setminus \Lambda_{k,\text{rhs}} \\ |m,n| \le k}} (f_k, \Psi_{k,m,n}^*) \Psi_{k,m,n} + \left( f_k - \sum_{\substack{[m,n] \in \mathbb{Z}^{2d} \\ |m,n| \le k}} (f_k, \Psi_{k,m,n}^*) \Psi_{k,m,n} \right),$$

we deduce from (4.14) and (4.15) that the last two terms have a $\hat{H}_k^p$ norm bounded by $C_{\varepsilon,p,m}(kR)^{-m}$, and the result follows. □

<figure><img src="image_3.png" /></figure>
4.5. **Approximability estimates.** We are now ready to present our approximability estimates for the Helmholtz problem. We start with an approximation result for general right-hand sides that does not hinge on Section 2, but rather on the results established in [4]. Recall that $N$ was introduced in (4.5).

**Theorem 4.5** (Approximability at a fixed frequency). Assume $f \in \hat{H}^p(\mathbb{R}^d)$ with $\text{supp } f \subset B(0, R)$. If

$$ \Lambda := \{[\mathbf{m}, \mathbf{n}] \in \mathbb{Z}^{2d} \mid \|[\mathbf{m}, \mathbf{n}]\| \le \sqrt{\rho(kR)^N}\}, $$

then, we have

$$ \left\| u - \sum_{[\mathbf{m}, \mathbf{n}] \in \Lambda} (u, \Psi_{\mathbf{k}, \mathbf{m}, \mathbf{n}}^*) \Psi_{\mathbf{k}, \mathbf{m}, \mathbf{n}} \right\|_{\hat{H}_k^q(\mathbb{R}^d)} \le C_p \rho^{-(2+p-q)/2} \|f\|_{\hat{H}_k^{p}(\mathbb{R}^d)} $$

for all $p, q \in \mathbb{N}$ with $q \le p$.

*Proof.* Let us set $D := (kR/\pi)^{-1/2}(\rho(kR)^N)^{1/2}$. We start with [4, Theorem 3.1], showing that

$$ \left\| u - \sum_{[\mathbf{m}, \mathbf{n}] \in \Lambda} (u, \Psi_{\mathbf{k}, \mathbf{m}, \mathbf{n}}^*) \Psi_{\mathbf{k}, \mathbf{m}, \mathbf{n}} \right\|_{\hat{H}_k^q(\mathbb{R}^d)} \le C_p D^{q-p-2} \|u\|_{\hat{H}_k^{p+2}(\mathbb{R}^d)} $$

The result then follows using (4.8), since

$$ \begin{aligned} D^{q-p-2} \|u\|_{\hat{H}_k^{p+2}(\mathbb{R}^d)} &\le C_p ((\rho(kR)^{N-1})^{(q-p-2)/2} \|u\|_{\hat{H}_k^{p+2}(\mathbb{R}^d)} \\ &\le C_p ((\rho(kR)^{N-1})^{(q-p-2)/2} (kR)^N \|f\|_{\hat{H}_k^p(\mathbb{R}^d)} = C_p \rho^{(q-p-2)/2} \|f\|_{\hat{H}_k^p(\mathbb{R}^d)}. \end{aligned} $$

□

Our second approximability estimate applies specifically to high-frequency scattering problems. It is a direct consequence of Theorem 2.2 and Lemma 4.4.

**Theorem 4.6** (Approximability in the high-frequency regime). Fix $\varepsilon > 0$ and consider the index set

$$ \Lambda := \{[\mathbf{m}, \mathbf{n}] \in \mathbb{Z}^{2d} \mid |p(\mathbf{m}, \mathbf{n})| \le (kR)^{-1/2+\varepsilon}\}. $$

Then, if the right-hand side is of the form $f_k := \tilde{\chi} e^{ik \cdot \mathbf{x}}$ where $\mathbf{d} \in \mathbb{R}^d$ with $|\mathbf{d}| = 1$ and $\tilde{\chi} \in C_c^\infty(B_R)$ is supported in the region where $\mu = 1$ and $\mathbf{A} = \mathbf{I}$, for all $q \in \mathbb{N}$, we have

$$ \left\| u_k - \sum_{[\mathbf{m}, \mathbf{n}] \in \Lambda} (u, \Psi_{\mathbf{k}, \mathbf{m}, \mathbf{n}}^*) \Psi_{\mathbf{k}, \mathbf{m}, \mathbf{n}} \right\|_{\hat{H}_k^q(\mathbb{R}^d)} \le C_{\varepsilon, q, m} (kR)^{-m} \quad \forall \mathbf{m} \in \mathbb{N}. $$

4.6. **A least-squares method.** In this section, we introduce a least-squares method based on Gaussian coherent states. For a finite set $\Lambda \subset \mathbb{Z}^{2d}$, we consider the discretization space

$$ W_{\Lambda} := \text{Vect}\{\Psi_{\mathbf{k}, \mathbf{m}, \mathbf{n}}; [\mathbf{m}, \mathbf{n}] \in \Lambda\}. $$

Then, the least squares method consists in finding $u_\Lambda \in W_\Lambda$ such that

$$ (4.16) \qquad (P_k u_\Lambda, P_k w_\Lambda) = (f, P_k w_\Lambda) $$

for all $w_\Lambda \in W_\Lambda$. Classically, there exists a unique solution $u_\Lambda$, and we have

$$ (4.17) \qquad \|P_k(u - u_\Lambda)\|_{L^2(\mathbb{R}^d)} = \min_{w_\Lambda \in W_\Lambda} \|P_k(u - w_\Lambda)\|_{L^2(\mathbb{R}^d)}. $$

Using Theorem 4.5 we can easily show that the method converges for any fixed frequency, if sufficiently many DOFs per wavelength are employed, as stated in the following corollary:
**Corollary 4.7 (Convergence at a fixed frequency).** Consider the set

$$ \Lambda_{\rho} := \{[\mathbf{m}, \mathbf{n}] \in \mathbb{Z}^{2d} \mid |\mathbf{m}|^2 + |\mathbf{n}|^2 \leq \rho(kR)^N\}, $$

then the following error estimates holds true:

$$ \|u - u_{\Lambda}\|_{H^2(\mathbb{R}^d)} \leq C_p \rho^{-p/2} \|f\|_{\hat{H}^p(\mathbb{R}^d)}. \quad \forall p \in \mathbb{N}. $$

Relying on Theorem 4.6 instead, we can show a refined error estimate in the high-frequency regime.

**Corollary 4.8 (Convergence in the high-frequency regime).** Fix $\varepsilon > 0$,

$$ \Lambda := \{[\mathbf{m}, \mathbf{n}] \in \mathbb{Z}^{2d} \mid |p(\mathbf{x}^{k,m}, \xi^{k,n})| \leq (kR)^{-1/2+\varepsilon}\}, $$

and let $f_k := \tilde{\chi}e^{ik\mathbf{d}\cdot\mathbf{x}}$ with $\tilde{\chi} \in C_c^{\infty}(B_R)$ supported in the region where $\mu = 1$ and $\mathbf{A} = \mathbf{I}$. Then, if $k$ is sufficiently large, we have

$$ \|u - u_{\Lambda}\|_{\hat{H}^2(\mathbb{R}^d)} \leq C_{\varepsilon,m}(kR)^{-m} \quad \forall m \in \mathbb{N}. $$

## 5. NUMERICAL RESULTS

In this section, we provide numerical illustrations of the above theory in the one-dimensional case. The purpose of these examples is simply to illustrate and support our theoretical findings. Extension to higher dimensions and discussion about efficient linear system assembly and solve will be reported elsewhere.

**5.1. Setting.** We consider the one dimensional case, where the differential operator in (4.6) simplifies to

$$ (5.1) \qquad P_k v = -\mu \nu u - k^{-2}(\alpha \nu^{-1} v')' $$

where $\mu, \alpha > 0$ are smooth “physical” coefficients, and $\nu$ is the PML scaling. For the sake of simplicity we take $R := 1$, meaning that $\mu = \alpha = 1$ outside of $(-1, 1)$ and that our right-hand sides $f$ will be supported in $(-1, 1)$. Here, we select $\nu := 1 + i\sigma$ where $\sigma$ is a stretching function defined as

$$ \sigma(x) := a(\mathbf{1}_{(1,+\infty)}(x-1)^r + \mathbf{1}_{(-\infty,-1)}(-1-x)^r) $$

with $a := 1/10$ and $r := 4$. Notice that this choice slightly departs from our theoretical framework as the coefficients are only $C^3$ and not $C^{\infty}$.

In all our experiments, the grid in phase space is chosen as

$$ x^{k,m} = \sqrt{k^{-1}\pi m}, \quad \xi^{k,n} = \sqrt{k^{-1}\pi n}, \quad n, m \in \mathbb{Z} $$

and for a given values of $k$ and $\delta$, the set of indices is taken to be

$$ \Lambda := \{(m, n) \in \mathbb{Z}^2 \mid |p(x^{k,m}, \xi^{k,n})| < \delta\}. $$

**5.2. Homogeneous medium with analytical solution.** Our first example concerns the case where $\alpha = \mu = 1$. The right-hand side is given by

$$ f_k := P_k(\phi e^{ikx}) = (\phi'' + 2ik\phi')e^{ikx}, $$

where $\phi$ is the only even function in $C^3(\mathbb{R}^d)$ such that $\phi = 0$ on $(-\infty, -3/4)$, $\phi = 1$ on $(-1/2, 0)$, and $\phi$ is a polynomial of degree 7 in $(-3/4, -1/2)$, and the associated solution is $u_k(x) := \phi(x)e^{ikx}$. Results are reported in Table 5.1. Figure 5.1 represents the points $(x^{k,m}, \xi^{k,n})$ included in the space $\Lambda$ for different values of $\delta$ for the case $k = 400$.

On Figure 5.2, we present the convergence history of the method as $\delta$ is increased for different values of $k$. These curves illustrate the fact that the method converges for any fixed $k$ when increasing the number of coherent Gaussian states, as predicted by our theory. Besides, on Figure 5.3 we represent the values of $\delta$ and corresponding number of freedom $N_{\text{dofs}}$ required to achieve a value of about $2 \times 10^{-5}$ (second column of Table 5.1) for different frequencies. Interestingly, the expected rates, namely, $\delta \sim k^{-1/2}$ and $N_{\text{dofs}} \sim k^{1/2}$ are observed.
<figure><img src="image_4.png" /><figcaption>FIGURE 5.1. Representation of $\Lambda$ for $k = 400$ in the homogeneous example. Top-left: $\delta = 0.1$. Top-right: $\delta = 0.4$. Bottom-left: $\delta = 0.6$. The horizontal and vertical distances between dots are $\sqrt{\pi/k}$.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIGURE 5.2. Convergence histories in the homogeneous example</figcaption></figure>
<table><caption>TABLE 5.1. Numerical results for different k while varying δ</caption><thead><tr><th>k</th><th colspan="6">Relative H<sup>1</sup><sub>k</sub> error (δ) [N<sub>dofs</sub>]</th></tr></thead><tbody><tr><td>20</td><td>3.3539e-01 (0.5) [60]</td><td>2.2356e-05 (2.0) [177]</td><td>6.6722e-06 (4.0) [249]</td><td></td><td></td><td></td></tr><tr><td>50</td><td>2.3180e-02 (0.4) [38]</td><td>1.7831e-05 (1.0) [229]</td><td>4.4104e-06 (3.0) [419]</td><td></td><td></td><td></td></tr><tr><td>100</td><td>4.4322e-01 (0.3) [96]</td><td>3.9496e-05 (0.8) [278]</td><td>2.5610e-06 (2.0) [645]</td><td></td><td></td><td></td></tr><tr><td>200</td><td>3.9807e-02 (0.2) [58]</td><td>1.8555e-05 (0.6) [334]</td><td>4.4586e-06 (1.0) [825]</td><td></td><td></td><td></td></tr><tr><td>400</td><td>3.6922e-01 (0.1) [74]</td><td>3.6289e-05 (0.4) [532]</td><td>2.3452e-05 (0.6) [646]</td><td></td><td></td><td></td></tr></tbody></table>

<figure><img src="image_6.png" /><figcaption>FIGURE 5.3. Computational cost for a fixed accuracy in the homogeneous example</figcaption></figure>

**5.3. Scattering in an heterogeneous medium.** We now focus on the case where $\alpha = 1$ and $\mu$ is the only even $C^3$ functions that equals 2 on $[0, 0.7]$, equals 1 on $[0.8, +\infty]$ and is polynomial of degree 7 on $[0.7, 0.8]$. We select the right-hand side $f_k(x) := k^2(\mu - 1)e^{ikx}$. Here, the analytical solution is not available, and instead, we rely on a reference solution computed by a Lagrange finite element method of order 4 on a grid with $h = 0.02 \cdot k^{-9/8}$ in order to avoid the pollution effect (see, e.g., [24]). The results are listed on Table 5.2, and Figure 5.6 represents the phase-space points included in the space $\Lambda$ for the case $k = 200$.

<table><caption>TABLE 5.2. Numerical results for $k = 20, 50, 100, 200$ while varying $\delta$</caption><thead><tr><th>k</th><th colspan="6">Relative H<sup>1</sup><sub>k</sub> error (δ) [N<sub>dofs</sub>]</th></tr></thead><tbody><tr><td>20</td><td>1.7861e-01 (2.0) [183]</td><td>4.1982e-02 (4.0) [255]</td><td>1.6388e-03 (8.0) [399]</td><td>1.3463e-04 (12.0) [569]</td><td></td><td></td></tr><tr><td>50</td><td>3.0122e-01 (2.0) [337]</td><td>2.5263e-02 (4.0) [511]</td><td>1.0504e-03 (6.0) [695]</td><td>4.5593e-04 (8.0) [841]</td><td></td><td></td></tr><tr><td>100</td><td>2.0614e-01 (1.0) [376]</td><td>7.7805e-03 (2.0) [665]</td><td>1.7647e-04 (4.0) [995]</td><td>1.3300e-04 (6.0) [1309]</td><td></td><td></td></tr><tr><td>200</td><td>3.8799e-01 (0.5) [250]</td><td>8.5455e-03 (1.0) [748]</td><td>6.9635e-04 (2.0) [1233]</td><td>1.2300e-04 (4.0) [1879]</td><td></td><td></td></tr></tbody></table>

We provide the same figures than in the previous experiment, and arrive at similar observation. First, Figure 5.4 illustrates the convergence of the method for fixed values of $k$ as the number of degrees of freedom is increased. Second, we compare the values of $\delta$ and $N_{dofs}$ to achieve an accuracy about $10^{-4}$ (last column of Table 5.2) for different frequencies. In this example too, the expected rates are numerically observed.

## 6. CONCLUSION

We propose a new family of finite-dimensional spaces to approximate the solutions to high-frequency Helmholtz problems. These discretization spaces are spanned by Gaussian coherent states, which have the key property to be micro-localised in phase space. This unique feature allows to carefully select which Gaussian coherent states are included in the discretization, leading to a frequency-aware discretization space specifically tailored to approximate the solutions to scattering problems efficiently.
<figure><img src="image_7.png" /><figcaption>FIGURE 5.4. Convergence histories in the heterogeneous example</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>FIGURE 5.5. Computational cost for a fixed accuracy in the heterogeneous example</figcaption></figure>

Our key findings correspond to two types of approximability results. First, assuming for simplicity that the problem is non-trapping, for general $L^2$ right-hand sides, we show that the Gaussian state approximation converges and provides a uniform error for all frequency if the number of degrees of freedom grows as $(kR)^d$. This result is similar to approximation results available for finite element discretizations. Our second result applies when the right-hand side corresponds to a plane-wave scattering problem. In this case, we show that it is sufficient for the number of degrees of freedom to grow only as $(kR)^{d-1/2}$ to achieve a constant accuracy for increasing frequencies. To the best of the authors knowledge, this last estimate suggests that the proposed discretization space requires substantially less degrees of freedom than any available method in the literature for high-frequency scattering problems in general smooth heterogeneous media.
<figure><img src="image_9.png" /><figcaption>FIGURE 5.6. Representation of $\Lambda$ for $k = 200$ in the heterogeneous example. Top-left: $\delta = 0.5$. Top-right: $\delta = 1$. Bottom-left: $\delta = 2$. The horizontal and vertical distances between dots are $\sqrt{\pi/k}$.</figcaption></figure>

We also present a set of numerical examples where our Gaussian state spaces are coupled with a least-squares variational formulation. Although the setting is elementary, these examples successfully illustrate the key features of our abstract analysis.

While we believe that the proposed results are very encouraging for a further development of the proposed method, there still remain several challenges that we would like to address in future works. First, we have chosen to focus on a least-squares method, because it is simpler to analyse than a Galerkin formulation. However, least-squares methods are typically poorly-conditioned as compared to their Galerkin counterparts. While there is no reason to believe that the proposed discrete space would not work with a Galerkin variational formulation, the analysis appears to be substantially more complex. Second, especially for three-dimensional problems, computing the entries of the matrices associated with Gaussian coherent states is not a trivial task, and specific algorithms should be employed. Third, although the entries of the discrete matrices decrease super-algebraically away from their diagonal, they are still dense in principle. While we are convinced that efficient truncation or compression methods can be employed, this is still to be analysed. Finally, the resulting linear systems (after possible truncation and/or compression) will likely have different properties than usual volume methods. As a result, the design of specific preconditioners will probably be required.

## APPENDIX A. INTERPLAY BETWEEN GAUSSIAN STATES AND DIFFERENTIAL OPERATORS

**Lemma A.1 (Action of a differential operator).** Given $\beta \in \mathbb{N}^d$ and smooth coefficients $(A_\alpha)_{\alpha \in \mathbb{N}^d}$, consider the differential operator

$$ (A.1) \qquad \mathcal{P}_{\hbar, \beta, A} := \sum_{\substack{\alpha \in \mathbb{N}^d \\ \alpha \le \beta}} \hbar^{[\alpha]} A_\alpha \partial^\alpha. $$
For all $[x_0, \xi_0] \in \mathbb{Z}^{2d}$ there exists a smooth function $g_{h,\beta,A,x_0,\xi_0}$ such that

$$
(A.2) \qquad (P_{h,\beta,A}\Phi_{h,x_0,\xi_0})(x) = g_{h,\beta,A,x_0,\xi_0}(x)\Phi_{h,x_0,\xi_0}(x) \quad \forall x \in \mathbb{R}^d.
$$

In addition, we have

$$
(A.3) \quad \left| (\partial^\gamma g_{h,\beta,A,x_0,\xi_0})(\mathbf{x}) \right| \le C_{\beta,A,\gamma} \left( 1 + |\mathbf{x} - \mathbf{x}_0|^{[\beta]} + |\xi_0|^{[\beta]} \right) \quad \forall \gamma \in \mathbb{N}^d.
$$

Proof. Recalling [4, Lemma A.2], for $\alpha \in \mathbb{N}$, we have

$$
(\partial^{\alpha} \Phi_{h, x_0, \xi_0})(\mathbf{x}) = \hbar^{-[\alpha]/2} q_{\alpha}(\hbar^{-1/2}(\mathbf{x} - \mathbf{x}_0 + i\xi_0)) \Phi_{h, x_0, \xi_0}(\mathbf{x}),
$$

where $q_\alpha$ is a polynomial of degree less than or equal to $[\alpha]$. Hence, we readily see that (A.2)
holds true with

$$
g_{h,\beta,A,x_0,\xi_0}(\mathbf{x}) := \sum_{\substack{\alpha \in \mathbb{N}^d \\ \alpha \le \beta}} \hbar^{[\alpha]/2} A_{\alpha}(\mathbf{x}) q_{\alpha}(\hbar^{-1/2}(\mathbf{x} - \mathbf{x}_0 + i\xi_0)).
$$

To establish (A.3), we start with Leibniz’ rule

$$
(\partial^\gamma g_{h,\beta,A,x_0,\xi_0})(\mathbf{x}) := \sum_{\substack{\alpha \in \mathbb{N}^d \\ \alpha \le \beta}} \hbar^{[\alpha]/2} \sum_{\substack{\gamma \in \mathbb{N}^d \\ \delta \le \gamma}} \begin{pmatrix} \delta \\ \gamma \end{pmatrix} \hbar^{-[\delta]/2} (\partial^\gamma - \delta A_\alpha)(\mathbf{y}) (\partial^\delta q_\alpha)(\hbar^{-1/2}(\mathbf{x} - \mathbf{x}_0 + i\xi_0)),
$$

so that

$$
| (\partial^\gamma g_{h,\beta,A,x_0,\xi_0})(\mathbf{y}) | \le C_{\beta,A,\gamma} \sum_{\substack{\alpha \in \mathbb{N}^d \\ \alpha \le \beta}} \hbar^{[\alpha]/2} \sum_{\substack{\delta \in \mathbb{N}^d \\ \delta \le \gamma}} \hbar^{-[\delta]/2} | (\partial^\delta q_\alpha)(\hbar^{-1/2}(\mathbf{x} - \mathbf{x}_0 + i\xi_0)) |.
$$

Next, we employ the fact that $\partial^\delta q_\alpha = 0$ whenever $\delta > \alpha$, and the estimate

$$
| (\partial^\delta q_\alpha)(\hbar^{-1/2}(\mathbf{x} - \mathbf{x}_0 + i\xi_0)) | \le C_{\alpha,\delta} \left( 1 + \hbar^{-([\alpha]-[\delta])/2} |\mathbf{x} - \mathbf{x}_0 + i\xi_0|^{[\alpha]-[\delta]} \right),
$$

leading to

$$
\begin{align*}
| (\partial^\gamma g_{h,\beta,A,x_0,\xi_0})(\mathbf{x}) | &\le C_{\beta,A,\gamma} \sum_{\substack{\alpha \in \mathbb{N}^d \\ \alpha \le \beta}} \hbar^{[\alpha]/2} \sum_{\substack{\gamma \in \mathbb{N}^d \\ \delta \le \gamma}} (1 + (\hbar^{-1/2} |\mathbf{x} - \mathbf{x}_0 + i\xi_0|)^{[\alpha]-[\delta]}) \\
&\le C_{\beta,A,\gamma} \sum_{\substack{\alpha \in \mathbb{N}^d \\ \alpha \le \beta}} \hbar^{[\alpha]/2} (1 + (\hbar^{-1/2} |\mathbf{x} - \mathbf{x}_0 + i\xi_0|)^{[\alpha]}) \\
&\le C_{\beta,A,\gamma} \sum_{\substack{\alpha \in \mathbb{N}^d \\ \alpha \le \beta}} (1 + |\mathbf{x} - \mathbf{x}_0 + i\xi_0|^{[\alpha]}) \\
&\le C_{\beta,A,\gamma} (1 + |\mathbf{x} - \mathbf{x}_0 + i\xi_0|^{[\beta]})
\end{align*}
$$

Then, (A.3) follows since

$$
|\mathbf{x} - \mathbf{x}_0 + i\xi_0|^{[\beta]} \le C_{\beta} \left( |\mathbf{x} - \mathbf{x}_0|^{[\beta]} + |\xi_0|^{[\beta]} \right).
$$

□

For the reader’s convenience, we recall the following basic fact.

**Proposition A.2 (Moments of a Gaussian).** *Consider the Gaussian function*

$$
(A.4) \qquad G_{\hbar}(\mathbf{y}) := \hbar^{-d/2} e^{-\frac{1}{\hbar} |\mathbf{y}|^2}.
$$

Then, for $\beta, \gamma \in \mathbb{N}^d$, we have

$$
(A.5) \quad \int_{\mathbb{R}^d} |\mathbf{x}^{\beta} (\partial^{\gamma} G_{\hbar}) (\mathbf{x} - \mathbf{x}_0)| d\mathbf{x} \le C_{\beta, \gamma} \hbar^{-[\gamma]/2} (\hbar^{[\beta]/2} + |\mathbf{x}_0|^{[\beta]}) \quad \forall \mathbf{x}_0 \in \mathbb{R}^d.
$$

We now show that Gaussian states are localised in phase-space.
**Lemma A.3 (Quasi orthogonality).** For all $[x_0, \xi_0] \in \mathbb{R}^{2d}$ and $[x'_0, \xi'_0] \in \mathbb{R}^{2d}$, we have

$$|(\mathcal{P}_{\hbar,\beta,A} \Phi_{\hbar,x_0,\xi_0}, \Phi_{\hbar,x'_0,\xi'_0})| \le C_{\beta,A,m} \hbar^{m/2} \frac{1 + |\xi_0|^{[\beta]} + |x_0 - x'_0|^{[\beta]}}{1 + |\xi_0 - \xi'_0|^m} e^{-\frac{1}{4\hbar}|x_0 - x'_0|^2} \quad \forall m \in \mathbb{N}.$$

In particular,

$$|(\mathcal{P}_{\hbar,\beta,A} \Phi_{\hbar,x_0,\xi_0}, \Phi_{\hbar,x'_0,\xi'_0})| \le C_{\beta,A,m} \hbar^{m/2} \frac{1 + |\xi_0|^{[\beta]}}{1 + |[x_0, \xi_0] - [x'_0, \xi'_0]|^m} \quad \forall m \in \mathbb{N}.$$

*Proof.* We have

$$\Phi_{\hbar,x_0,\xi_0}(\boldsymbol{x}) \overline{\Phi_{\hbar,x'_0,\xi'_0}(\boldsymbol{x})} = \theta \hbar^{-d/2} e^{-\frac{1}{4\hbar}|\boldsymbol{x}_0 - \boldsymbol{x}'_0|^2} e^{-\frac{1}{\hbar}|\boldsymbol{x} - \frac{1}{2}(\boldsymbol{x}_0 + \boldsymbol{x}'_0)|^2} e^{\frac{i}{\hbar}(\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0) \cdot \boldsymbol{x}}$$

where $\theta \in \mathbb{C}$ with $|\theta| \le \pi^{-d/2}$. Recalling (A.2) and the notation (A.4) from Proposition A.2, we have

$$(\mathcal{P}_{\hbar,\beta,A} \Phi_{\hbar,x_0,\xi_0})(\boldsymbol{x}) \overline{\Phi_{\hbar,x'_0,\xi'_0}(\boldsymbol{x})} = \theta g_{\hbar,\beta,A,x_0,\xi_0}(\boldsymbol{x}) e^{-\frac{1}{4\hbar}|\boldsymbol{x}_0 - \boldsymbol{x}'_0|^2} G_{\hbar} \left( \boldsymbol{x} - \frac{1}{2}(\boldsymbol{x}_0 + \boldsymbol{x}'_0) \right) e^{\frac{i}{\hbar}(\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0) \cdot \boldsymbol{x}}$$

so that

$$|(\mathcal{P}_{\hbar,\beta,A} \Phi_{\hbar,x_0,\xi_0}, \Phi_{\hbar,x'_0,\xi'_0})| \le C e^{-\frac{1}{4\hbar}|\boldsymbol{x}_0 - \boldsymbol{x}'_0|^2} |I|,$$

with

$$I := \hbar^{-d/2} \int_{\mathbb{R}^d} g_{\hbar,\beta,A,x_0,\xi_0}(\boldsymbol{x}) G_{\hbar,x_0,x'_0}(\boldsymbol{x}) e^{\frac{i}{\hbar}(\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0) \cdot \boldsymbol{x}} d\boldsymbol{x},$$

where $G_{\hbar,x_0,x'_0}(\boldsymbol{x}) := G_{\hbar}(\boldsymbol{x} - \frac{1}{2}(\boldsymbol{x}_0 + \boldsymbol{x}'_0))$.

If $\xi_0 = \xi'_0$, this concludes the proof since

$$\begin{aligned} |I| &\le \hbar^{-d/2} \int_{\mathbb{R}^d} |g_{\hbar,\beta,A,x_0,\xi_0}(\boldsymbol{x})| G_{\hbar,x_0,x'_0}(\boldsymbol{x}) d\boldsymbol{x} \\ &\le C_{\beta} \int_{\mathbb{R}} (1 + |\boldsymbol{x} - \boldsymbol{x}_0|^{[\beta]}) G_{\hbar,x_0,x'_0}(\boldsymbol{x}) d\boldsymbol{x} \\ &= C_{\beta} \int_{\mathbb{R}^d} (1 + |\boldsymbol{x}|^{[\beta]}) G_{\hbar} \left( \boldsymbol{x} - \frac{1}{2}(\boldsymbol{x}_0 - \boldsymbol{x}'_0) \right) d\boldsymbol{x} \\ &\le C_{\beta} (1 + |\boldsymbol{x} - \boldsymbol{x}_0|^{[\beta]}). \end{aligned}$$

If $\xi_0 \ne \xi'_0$, then there exists a $j \in \{1, \dots, d\}$ such that $|(\xi_0 - \xi'_0)_j| \ge 1/\sqrt{d}|\xi_0 - \xi'_0|$. We integrate $m$ times by part with respect to $x_j$, leading to

$$\begin{aligned} I &= \hbar^{-d/2} \left( \frac{i}{\hbar} (\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0)_j \right)^{-m} \int_{\mathbb{R}^d} \partial_j^m (g_{\hbar,\beta,A,x_0,\xi_0} G_{\hbar,x_0,x'_0}(\boldsymbol{x})) e^{\frac{i}{\hbar}(\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0) \cdot \boldsymbol{x}} d\boldsymbol{x} \\ &\le C \hbar^m |\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0|^{-m} I'. \end{aligned}$$

where

$$\begin{aligned} |I'| &= \left| \hbar^{-d/2} \int_{\mathbb{R}^d} \partial_j^m (g_{\hbar,\beta,A,x_0,\xi_0} G_{\hbar,x_0,x'_0}(\boldsymbol{x})) e^{\frac{i}{\hbar}(\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0) \cdot \boldsymbol{x}} d\boldsymbol{x} \right| \\ &= \left| \hbar^{-d/2} \sum_{\ell=0}^{m} \binom{m}{\ell} \int_{\mathbb{R}^d} (\partial_j^{m-\ell} g_{\hbar,\beta,A,x_0,\xi_0}(\boldsymbol{x})) (\partial_j^{\ell} G_{\hbar,x_0,x'_0}(\boldsymbol{x})) e^{\frac{i}{\hbar}(\boldsymbol{\xi}_0 - \boldsymbol{\xi}'_0) \cdot \boldsymbol{x}} d\boldsymbol{x} \right| \\ &\le C_m \sum_{\ell=1}^{m} I_{\ell}'' . \end{aligned}$$

with

$$I_{\ell}'' := \hbar^{-d/2} \int_{\mathbb{R}^d} |\partial_j^{m-\ell} g_{\hbar,\beta,A,x_0,\xi_0}(\boldsymbol{x})| |\partial_j^{\ell} G_{\hbar,x_0,x'_0}(\boldsymbol{x})| d\boldsymbol{x}.$$
We then employ (A.3) and (A.5), showing that

$$
\begin{align*}
|I_{\ell}''| & \le C_{\beta,A,m} \hbar^{-d/2} \int_{\mathbb{R}^d} (1 + |\mathbf{x} - \mathbf{x}_0|^{[\beta]} + |\boldsymbol{\xi}_0|^{[\beta]}) |\partial_{\ell}^j G_{\hbar, \mathbf{x}_0, \mathbf{x}'_0}(\mathbf{x})| d\mathbf{x} \\
& = C_{\beta,A,m} \hbar^{-d/2} \int_{\mathbb{R}^d} (1 + |\mathbf{x}|^{[\beta]} + |\boldsymbol{\xi}_0|^{[\beta]}) |\partial_{\ell}^j G_{\hbar, \mathbf{x}_0, \mathbf{x}'_0}(\mathbf{x} + \mathbf{x}_0)| d\mathbf{x} \\
& = C_{\beta,A,m} \int_{\mathbb{R}^d} (1 + |\mathbf{x}|^{[\beta]} + |\boldsymbol{\xi}_0|^{[\beta]}) \left| \partial_{\ell}^j (G_{\hbar}) \left( \mathbf{x} - \frac{1}{2}(\mathbf{x}_0 - \mathbf{x}'_0) \right) \right| d\mathbf{x} \\
& \le C_{\beta,A,m} \hbar^{-\ell/2} \left( (1 + |\boldsymbol{\xi}_0|^{[\beta]}) + \left| \frac{1}{2}(\mathbf{x}_0 - \mathbf{x}'_0) \right|^{[\beta]} \right) \\
& \le C_{\beta,A,m} \hbar^{-\ell/2} (1 + |\boldsymbol{\xi}_0|^{[\beta]} + |\mathbf{x}_0 - \mathbf{x}'_0|^{[\beta]}).
\end{align*}
$$

□

**Lemma A.4** (Action of $(P_h - p)$ on Gaussian states). *Consider the second-order differential operator*

$$
P_h := \hbar^2 \sum_{k, \ell=1}^{d} (A_h)_{k\ell} \partial_{k\ell} + i\hbar b_h \cdot \nabla + c_h,
$$

where $A_h, b_h$ and $c_h$ are smooth functions, together with its symbol

$$
p(\mathbf{x}, \boldsymbol{\xi}) = A_h(\mathbf{x})\boldsymbol{\xi} \cdot \boldsymbol{\xi} + b_h(\mathbf{x}) \cdot \boldsymbol{\xi} + c_h(\mathbf{x}).
$$

For $\mathbf{x}_0, \boldsymbol{\xi}_0 \in \mathbb{R}^d$, we have

$$
(A.6) \qquad (P_h - p(\mathbf{x}_0, \boldsymbol{\xi}_0)) \Phi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0} = r_h, \mathbf{x}_0, \boldsymbol{\xi}_0 \Phi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0}
$$

with

$$
\begin{align*}
r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0}(\mathbf{x}) &:= \hbar \operatorname{tr} A_h(\mathbf{x}) \\
&\quad + i[(A_h + A_h^T)(\mathbf{x})\boldsymbol{\xi}_0 + b_h(\mathbf{x})] \cdot (\mathbf{x} - \mathbf{x}_0) \\
&\quad + (A_h(\mathbf{x}) - A_h(\mathbf{x}_0))\boldsymbol{\xi}_0 \cdot \boldsymbol{\xi}_0 + i(b_h(\mathbf{x}) - b_h(\mathbf{x}_0)) \cdot \boldsymbol{\xi}_0 + c_h(\mathbf{x}) - c_h(\mathbf{x}_0) \\
&\quad + A_h(\mathbf{x})(\mathbf{x} - \mathbf{x}_0) \cdot (\mathbf{x} - \mathbf{x}_0).
\end{align*}
$$

In particular, there exists smooth functions $\alpha_{h,x_0,\xi_0,\beta}$ such that

$$
(A.7) \qquad r_{h,x_0,\xi_0}(\mathbf{x}) = \sum_{\substack{\beta \in \mathbb{N}^d \\ |\beta| \le 2}} \hbar^{(1-|\beta|)/2} \alpha_{h,x_0,\xi_0,\beta}(\mathbf{x}) (\mathbf{x} - \mathbf{x}_0)^{\beta}
$$

and

$$
(A.8) \qquad \| \partial^\gamma \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, \beta} \|_{L^\infty(\mathbb{R}^d)} \le C_{\beta, \gamma} (1 + |\boldsymbol{\xi}_0|^2).
$$

*Proof.* We obtain (A.6) by straightforward computations. We easily identify that

$$
\alpha_{h,x_0,\xi_0,0} = \hbar^{1/2} \operatorname{tr} A_h \qquad \alpha_{h,x_0,\xi_0,e_j+e_\ell} = \hbar^{1/2} (A_h)_{j\ell}.
$$

For the terms with $[\beta] = 1$ we further write that for $\varphi \in \{(A_h)_{j\ell}, (b_h)_\ell, c_h\}$, there exist smooth functions $\varphi^{\mathbf{x}_0, k}$ such that

$$
\varphi(\mathbf{x}) - \varphi(\mathbf{x}_0) = \sum_{k=1}^{d} \varphi^{\mathbf{x}_0, k}(\mathbf{x}) \cdot (\mathbf{x} - \mathbf{x}_0),
$$

so that

$$
\alpha_{\mathbf{x}_0, \boldsymbol{\xi}_0, e_k} = A_{\hbar}^{\mathbf{x}_0, k} \boldsymbol{\xi}_0 \cdot \boldsymbol{\xi}_0 + i \mathbf{b}_{\hbar}^{\mathbf{x}_0, k} \cdot \boldsymbol{\xi}_0 + c_{\hbar}^{\mathbf{x}_0, k} + i \left( \sum_{j=1}^{d} ((A_{\hbar})_{kj} + (A_{\hbar})_{jk}) \boldsymbol{\xi}_{0,j} + (\mathbf{b}_{\hbar})_k \right)
$$

for $k \in \{1, \dots, d\}$.

□

We will now explore further the action of $P_h$ on $\psi_{h,\mathbf{x}_0,\boldsymbol{\xi}_0}$, and show that it is close to $p(\mathbf{x}_0, \boldsymbol{\xi}_0)\psi_{h,\mathbf{x}_0,\boldsymbol{\xi}_0}$
**Lemma A.5 (Action of $(P_h - p)^L$ on Gaussian states).** For $L \in \mathbb{N}$, we have

$$ (P_h - p(\mathbf{x}_0, \boldsymbol{\xi}_0))^L \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0} = r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L} \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0}, $$

where $r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}$ can be written in the form

$$ r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}(\mathbf{x}) = \sum_{\substack{\boldsymbol{\beta} \in \mathbb{Z}^d \\ |\boldsymbol{\beta}| \le 2L}} \hbar^{(L-[\boldsymbol{\beta}])/2} \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}}(\mathbf{x}) (\mathbf{x} - \mathbf{x}_0)^{\boldsymbol{\beta}}, $$

where the functions $\alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, \boldsymbol{\beta}}$ are smooth and satisfy

$$ (A.9) \qquad \|\partial^\gamma \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}}\|_{L^\infty(\mathbb{R}^d)} \le C_{L, \boldsymbol{\beta}, \gamma} (1 + |\boldsymbol{\xi}_0|^2)^L. $$

*Proof.* We shall prove this result by induction. The case $L = 0$ trivially holds with $\alpha_{\mathbf{x}_0, \boldsymbol{\xi}_0, 0, 0} = 1$. The case $L = 1$ is treated in Lemma A.4. Hence, let us assume that $L \in \mathbb{N}$ is such that (A.9) holds for all $L' \le L$. We have

$$ \begin{aligned} (P_h - p(\mathbf{x}_0, \boldsymbol{\xi}_0))^{L+1} \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0} &= (P_h - p(\mathbf{x}_0, \boldsymbol{\xi}_0)) (r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L} \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0}) \\ &= [P_h - p(\mathbf{x}_0, \boldsymbol{\xi}_0), r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}] \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0} + r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L} r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, 1} \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0}, \end{aligned} $$

Writing the operator $P_h$ in the form (2.2), we deduce that

$$ \begin{aligned} [P_h - p(\mathbf{x}_0, \boldsymbol{\xi}_0), r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}] \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0} &= \hbar^2 \sum_{j, \ell=1}^{d} a_{j\ell} \left( \frac{\partial^2 r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j \partial x_\ell} \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0} + \frac{\partial r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j} \frac{\partial \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0}}{\partial x_\ell} \right) \\ &\quad + i\hbar \sum_{j=1}^{d} b_j \frac{\partial r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j} \psi_{h, \mathbf{x}_0, \boldsymbol{\xi}_0} \\ &= \tilde{r}_{k, \mathbf{x}_0, \boldsymbol{\xi}_0, L+1} \psi_{k, \mathbf{x}_0, \boldsymbol{\xi}_0}, \end{aligned} $$

where

$$ \tilde{r}_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L+1} := \sum_{j, \ell=1}^{d} \left( \hbar^2 a_{j\ell} \frac{\partial^2 r_{k, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j \partial x_\ell} + \hbar \frac{\partial r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j} (i\boldsymbol{\xi}_0, \ell - (x_\ell - x_{0, \ell})) \right) + i\hbar \sum_{j=1}^{d} b_j \frac{\partial r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j}. $$

Using the induction hypothesis, we have

$$ \begin{aligned} \frac{\partial r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j}(\mathbf{x}) &= \sum_{\substack{\boldsymbol{\beta} \in \mathbb{Z}^d \\ |\boldsymbol{\beta}| \le 2L}} \hbar^{\max(L-[\boldsymbol{\beta}], 0)/2} \frac{\partial \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}}}{\partial x_j}(\mathbf{x}) (\mathbf{x} - \mathbf{x}_0)^{\boldsymbol{\beta}} \\ &\quad + \sum_{\substack{\boldsymbol{\beta} \in \mathbb{Z}^d \\ |\boldsymbol{\beta}| \le 2L}} \hbar^{\max(L-[\boldsymbol{\beta}], 0)/2} \beta_j \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}}(\mathbf{x}) (\mathbf{x} - \mathbf{x}_0)^{\boldsymbol{\beta}-e_j} \\ &= \sum_{\substack{\boldsymbol{\beta} \in \mathbb{Z}^d \\ |\boldsymbol{\beta}| \le 2L}} \hbar^{\max(L-[\boldsymbol{\beta}], 0)/2} \frac{\partial \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}}}{\partial x_j}(\mathbf{x}) (\mathbf{x} - \mathbf{x}_0)^{\boldsymbol{\beta}} \\ &\quad + \sum_{\substack{\boldsymbol{b} \in \mathbb{Z}^d \\ |\boldsymbol{b}+e_j] \le 2L}} \hbar^{\max(L-1-[\boldsymbol{b}], 0)/2} (\beta_j + 1) \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}+e_j}(\mathbf{x}) (\mathbf{x} - \mathbf{x}_0)^{\boldsymbol{\beta}}. \end{aligned} $$

Setting $\alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, \boldsymbol{\beta}} = 0$ if $[\boldsymbol{\beta}] > 2L$, this gives us

$$ \hbar \frac{\partial r_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L}}{\partial x_j}(\mathbf{x}) = \sum_{\substack{\boldsymbol{\beta} \in \mathbb{Z}^d \\ |\boldsymbol{\beta}| \le 2L}} \hbar^{\max(L+1-[\boldsymbol{\beta}], 0)/2} \left( \hbar \frac{\partial \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}}}{\partial x_j}(\mathbf{x}) + (\beta_j + 1) \alpha_{h, \mathbf{x}_0, \boldsymbol{\xi}_0, L, \boldsymbol{\beta}+e_j}(\mathbf{x}) \right) (\mathbf{x} - \mathbf{x}_0)^{\boldsymbol{\beta}}, $$

these terms are multiplied either by a smooth function, a term of the order $|\boldsymbol{\xi}_0|$ or a power of $(\mathbf{x} - \mathbf{x}_0)$, which always enters the induction for $L + 1$. The term with the second derivative is treated similarly.
For the remaining term, we simply write that

$$
\begin{aligned}
r_{\hbar, x_0, \xi_0, L} r_{\hbar, x_0, \xi_0, 1} &= \left( \sum_{\substack{\beta' \in \mathbb{Z}^{2d} \\ [\beta'] \le 2}} \hbar^{\max(1-[\beta'], 0)/2} \alpha_{\hbar, x_0, \xi_0, 1, \beta'}(x-x_0)^{\beta'} \right) \left( \sum_{\substack{\beta \in \mathbb{Z}^{2d} \\ [\beta] \le 2L}} \hbar^{\max(L-[\beta], 0)/2} \alpha_{\hbar, x_0, \xi_0, L, \beta}(x-x_0)^{\beta} \right) \\
&= \sum_{\substack{\beta \in \mathbb{Z}^{2d} \\ [\beta'] \le 2}} \sum_{\substack{\beta' \in \mathbb{Z}^{2d} \\ [\beta'] \le 2L}} \hbar^{\max(L-[\beta], 0)/2 + \max(1-[\beta'], 0)/2} \alpha_{\hbar, x_0, \xi_0, L, \beta} \alpha_{\hbar, x_0, \xi_0, 1, \beta'}(x-x_0)^{\beta+\beta'},
\end{aligned}
$$

and the result follows since

$$
\max(L - [\beta], 0) + \max(1 - [\beta'], 0) \geq \max(L + 1 - [\beta + \beta'], 0).
$$

<figure><img src="image_10.png" /></figure>

**Proposition A.6 (Control of $(P_h - p)^L$).** We have

$$
\| (P_h - p(\mathbf{x}_0, \boldsymbol{\xi}_0))^L \Phi_{\hbar, \mathbf{x}_0, \boldsymbol{\xi}_0} \|^2 \leq C_{\mathbf{A}, \mathbf{b}, c} (1 + |\boldsymbol{\xi}_0|^2) \hbar^{L/2}.
$$

for all $[\mathbf{x}_0, \mathbf{\xi}_0] \in \mathbb{R}^{2d}$.

*Proof.* We deduce from the previous lemma that

$$
\begin{aligned}
\|(P_\hbar - p(\mathbf{x}_0, \boldsymbol{\xi}_0))^L \psi_{\mathbf{x}_0, \boldsymbol{\xi}_0, k}\|^2 &\le C(1 + |\boldsymbol{\xi}_0|^2)^L (\pi\hbar)^{-d/2} \int_{\mathbb{R}^d} (|\mathbf{x} - \mathbf{x}_0| + \hbar)^L e^{-\frac{1}{\hbar}|\mathbf{x} - \mathbf{x}_0|^2} d\mathbf{x} \\
&= C(1 + |\boldsymbol{\xi}_0|^2)(\pi\hbar)^{-d/2} \int_{\mathbb{R}^d} (\hbar^{1/2}|\mathbf{y}| + \hbar)^L e^{-|\mathbf{y}|^2} d\mathbf{y} \\
&\le C(1 + |\boldsymbol{\xi}_0|^2)\hbar^{L/2}.
\end{aligned}
$$

<figure><img src="image_11.png" /></figure>

REFERENCES

1. B. Adcock and D. Huybrechs, *Frames and numerical approximation*, SIAM Review **61** (2019), no. 3, 443–473.

2. J. P. Bérenger, *A perfectly matched layer for the absorption of electromagnetics waves*, J. Comput. Phys. **114** (1994), 185–200.

3. O. Cessenat and B. Després, *Using plane waves as base functions for solving time harmonic equations with the ultra weak variational formulation*, J. Comput. Acous. **11** (2003), no. 02, 227–238.

4. T. Chaumont-Frelet and M. Ingremeau, *Decay of coefficients and approximation rates in Gabor Gaussian frames*, hal-03746979, 2022.

5. T. Chaumont-Frelet and F. Valentin, *A multiscale hybrid-mixed method for the Helmholtz equation in heterogeneous domains*, SIAM J. Numer. Anal. **58** (2020), no. 2, 1029–1067.

6. F. Collino and P. Monk, *The perfectly matched layer in curvilinear coordinates*, SIAM J. Sci. Comp. **19** (1998), no. 6, 2061–2090.

7. I. Daubechies, A. Grossman, and Y. Meyer, *Painless nonorthogonal expansions*, J. Math. Phys. **27** (1986), no. 5, 1271–1283.

8. R. C. Dorf, *Electronics, power electronics, optoelectronics, microwaves, electromagnetics and radar*, Taylor & Francis, 2006.

9. S. Dyatlov and M. Zworski, *Mathematical theory of scattering resonances*, vol. 200, American Mathematical Soc., 2019.

10. B. Engquist and O. Runborg, *Computational high frequency wave propagation*, Acta Numerica **2003** (2003), 181–266.

11. E. Faou, V. Gradinaru, and C. Lubich, *Computing semiclassical quantum dynamics with hagedorn wavepackets*, SIAM J. Sci. Comput. **31** (2009), no. 4, 3027–3041.

12. C. Farhat, I. Harari, and L. P. Franca, *The discontinuous enrichment method*, Comput. Methods Appl. Mech. Engrg. **190** (2001), no. 48, 6455–6479.

13. D. Gabor, *Theory of communication. part 1: The analysis of information*, Journal of the Institution of Electrical Engineers-Part III: Radio and Communication Engineering **93** (1946), no. 26, 429–441.

14. J. Galkowski, D. Lafontaine, and E. A. Spence, *Perfectly-matched-layer truncation is exponentially accurate at high frequency*, arXiv preprint arXiv:2105.07737 (2021).

15. J. Galkowski and E.A. Spence, *The Helmholtz boundary element method does not suffer from the pollution effect*, arXiv:2201.09721, 2022.
16. J. Galkowski, E.A. Spence, and J. Wunsch, *Optimal constants in nontrapping resolvent estimates and applications in numerical analysis*, Pure and Applied Analysis **2** (2019), no. 1, 157–202.

17. C.J. Gittelson, R. Hiptmair, and I. Perugia, *Plane wave discontinuous Galerkin methods: analysis of the h-version*, ESAIM Math. Model. Numer. Anal. **43** (2009), 297–331.

18. V. Gradinaru and G.A. Hagedorn, *Convergence of a semiclassical wavepacket based time-splitting for the Schrödinger equation*, Numer. Math. **126** (2014), no. 1, 53–73.

19. V. Gradinaru and O. Rietmann, *Hagedorn wavepackets and Schrödinger equation with time-dependent, homogeneous magnetic field*, J. Comput. Phys. **445** (2021), 110581.

20. L. Greengard and V. Rokhlin, *A fast algorithm for particle simulations*, J. Comput. Phys. **73** (1987), no. 2, 325–348.

21. W. Hackbusch, *Hierarchical matrices: algorithms and analysis*, vol. 49, Springer, 2015.

22. R. Hiptmair, A. Moiola, and I. Perugia, *Plane wave discontinuous Galerkin methods for the 2d Helmholtz equation: analysis of the p-version*, SIAM J. Numer. Anal. **49** (2011), no. 1, 264–284.

23. __________, *A survey of Trefftz methods for the Helmholtz equation*, pp. 237–279, Springer International Publishing, 2016.

24. F. Ihlenburg and I. Babuška, *Finite element solution of the Helmholtz equation with high wave number. Part II: The h-p version of the FEM*, SIAM J. Numer. Anal. **34** (1997), no. 1, 315–358.

25. D. Lafontaine, E.A. Spence, and J. Wunsch, *For most frequencies, strong tapping has a weak effect in frequency-domain scattering*, preprint arXiv:1903.12172, 2019.

26. __________, *Wavenumber-explicit convergence of the hp-fem for the full-space heterogeneous Helmholtz equation with smooth coefficients*, Comput. Math. Appl. **113** (2022), 59–69.

27. C. Lasser and C. Lubich, *Computing quantum dynamics in the semiclassical regime*, Acta Numerica **29** (2020), 229–401.

28. J. M. Melenk and I. Babuška, *The partition of unity finite element method: basic theory and applications*, Comput. Methods Appl. Mech. Engrg. **139** (1996), no. 1-4, 289–314.

29. J.M. Melenk and S. Sauter, *Convergence analysis for finite element discretizations of the Helmholtz equation with Dirichlet-to-Neumann boundary conditions*, Math. Comp. **79** (2010), no. 272, 1871–1914.

30. __________, *Wavenumber explicit convergence analysis for Galerkin discretizations of the Helmholtz equation*, SIAM J. Numer. Anal. **49** (2011), no. 3, 1210–1243.

31. P. Monk and D.-Q. Wang, *A least-squares method for the Helmholtz equation*, Comput. Methods Appl. Mech. Engrg. **175** (1999), no. 1-2, 121–136.

32. S. Nonnenmacher, *Spectral problems in open quantum chaos*, Nonlinearity **24** (2011), no. 12, R123.

33. S. Nonnenmacher and M. Zworski, *Quantum decay rates in chaotic scattering*, Acta mathematica **203** (2009), no. 2, 149–233.

34. H. Riou, P. Ladeveze, and B. Sourcis, *The multiscale VTCR approach applied to acoustics problems*, J. Comput. Acous. **16** (2008), no. 04, 487–505.

35. S.A. Sauter and C. Schwab, *Boundary element methods*, Springer, 2010.

36. A. Tarantola, *Inversion of seismic reflection data in the acoustic approximation*, Geophysics **49** (1984), no. 8, 1259–1266.

37. M. Zworski, *Semiclassical analysis*, vol. 138, American Mathematical Soc., 2012.

38. __________, *Mathematical study of scattering resonances*, Bulletin of Mathematical Sciences **7** (2017), no. 1, 1–85.