# Regular Bohr-Sommerfeld quantization rules for a *h*-pseudo-differential operator

## The method of positive commutators

Abdelwaheb Ifa* — Michel Rouleux**

* Université de Tunis El-Manar, Département des mathématiques
2092 Manar II
TUNISIA
abdelwaheb.ifa@fsm.rnu.tn

** Aix Marseille Univ, Univ Toulon, CNRS, CPT
FRANCE
rouleux@univ-tln.fr

**ABSTRACT.** We revisit in this Note the well known Bohr-Sommerfeld quantization rule (BS) for a 1-D Pseudo-differential self-adjoint Hamiltonian within the algebraic and microlocal framework of Helffer and Sjöstrand; BS holds precisely when the Gram matrix consisting of scalar products of some WKB solutions with respect to the "f ux norm" is not invertible.

**RÉSUMÉ.** Dans le cadre algébrique et microlocal élaboré par Helffer et Sjöstrand, on propose une ré-écriture de la règle de quantification de Bohr-Sommerfeld (BS) pour un opérateur auto-adjoint *h*-Pseudo-différentiel 1-D; elle s'exprime par la non-inversibilité de la matrice de Gram d'un couple de solutions WKB dans une base convenable, pour le produit scalaire associé à la "norme de f ux".

**KEYWORDS :** Semi-classical spectral asymptotics, quantization rules.

**MOTS-CLÉS :** Analyse spectrale semi-classique, règles de quantification.
# 1. Introduction

Let $p(x, \xi; h)$ be a smooth real classical Hamiltonian on $T^*\mathbb{R}$; we will assume that $p$ belongs to the space of symbols $S^0(m)$ for some order function $m$ (for example $m(x, \xi) = (1 + |\xi|^2)^M$) with

$$S^N(m) = \{p \in C^\infty(T^*\mathbb{R}) : \forall \alpha \in \mathbb{N}^2, \exists C_\alpha > 0, |\partial^\alpha p(x, \xi; h)| \le C_\alpha h^N m(x, \xi)\} \quad (1)$$

This allows to take Weyl quantization $P = p^w(x, h D_x; h)$ of $p$

$$P(x, hD_x; h)u(x; h) = (2\pi h)^{-1} \iint e^{\frac{i}{\hbar}(x-y)} \eta p(\frac{x+y}{2}, \eta; h) u(y) dy d\eta \quad (2)$$

so that $p^w(x, h D_x; h)$ is self-adjoint. We also assume that $p+i$ is elliptic (in the classical sense). We call as usual $p_0$ the principal symbol, and $p_1$ the sub-principal symbol; in case of Schrödinger operator $P(x, h D_x; h) = (h D_x)^2 + V(x)$, $p(x, \xi; h) = p_0(x, \xi) = \xi^2 + V(x)$. We make the hypothesis of [7], namely:

Fix some compact interval $I = [E_-, E_+]$, $E_- < E_+$ and assume that there exists a topological ring $\mathcal{A} \subset p_0^{-1}(I)$ such that $\partial \mathcal{A} = A_+ \cup A_-$ with $A_\pm$ a connected component of $p_0^{-1}(E_\pm)$. Assume also that $p_0$ has no critical point in $\mathcal{A}$, and $A_-$ is included in the disk bounded by $A_+$ (if this is not the case, we can always change $p$ to $-p$). We define the microlocal well $W$ as the disc bounded by par $A_+$.

For $E \in I$, let $\gamma_E \subset W$ be a periodic orbit in the energy surface $\{p_0(x, \xi) = E\}$ (so that $\gamma_E$ is an embedded Lagrangian manifold).

Let $K_h^N(E)$ be the microlocal kernel of $P - E$ of order $N$, i.e. the space of local solutions of $(P - E)u = \mathcal{O}(h^{N+1})$ in the distributional sense, microlocalized on $\gamma_E$. This is a smooth complex vector bundle over $\pi_x(\gamma_E)$. Here we address the problem of finding the set of $E = E(h)$ such that $K_h^N(E)$ contains a global section, i.e. of constructing a sequence of quasi-modes (QM) $(u_n(h), E_n(h))$ of a given order $N$. As usual we denote by $K_h(E)$ the microlocal kernel of $P - E \mod \mathcal{O}(h^\infty)$; since the distinction between $K_h^N(E)$ and $K_h(E)$ plays no important role here, we shall content to write $K_h(E)$.

Then if $E_+ < E_0 = \liminf_{|x, \xi| \to \infty} p_0(x, \xi)$, all eigenvalues of $P$ in $I$ are indeed given by Bohr-Sommerfeld quantization condition (BS) $S_h(E_n(h)) = 2\pi nh$, where the semi-classical action $S_h(E)$ has the asymptotics

$$S_h(E) \sim S_0(E) + hS_1(E) + h^2S_2(E) + \dots$$

We determine BS at any accuracy by computing quasi-modes. There are a lot of ways to derive BS: the method of matching of WKB solutions [3], known also as Liouville-Green method [15], which has received many improvements (see [21]), the method of the monodromy operator (see [12] and references therein), the method of quantization deformation based on Functional Calculus and Trace Formulas [14], [7], [6], [11], [1]. Note that the latter one already assumes BS, it only gives a very convenient way to derive it. In the real analytic case, BS rule, and also tunneling expansions, can be obtained using the so-called “exact WKB method” see e.g. [8], [9] when $P = -h^2\Delta + V(x)$ is Schrödinger operator.

Here we present another way to construct quasi-modes of order 2, based on [18], [13]. We stress that our method in the present scalar case, when carried to second order, is a bit
more intricated than [14], [7] and its refinements [11]; it is most useful for matrix valued operators with double characteristics such as Bogoliubov-de Gennes Hamiltonian ([10], [4], [5]), or Born-Oppenheimer type Hamiltonians ([2], [16]).

## 2. The microlocal Wronskian

The best algebraic and microlocal framework for computing 1-D quantization rules in the self-adjoint case, cast in the fundamental works [18], [13], is based on Fredholm theory, and the classical "positive commutator method" using conservation of some quantity called a "quantum flux".

Bohr-Sommerfeld quantization rules result in constructing quasi-modes by WKB approximation along a closed Lagrangian manifold $\Lambda_E \subset \{p_0 = E\}$, i.e. a periodic orbit of Hamilton vector field $H_p$ with energy $E$. This can be done locally according to the rank of the projection $\Lambda_E \to \mathbb{R}_x$.

Thus the set $K_h(E)$ of asymptotic solutions to $(P - E)u = 0$ along $\Lambda_E$ can be considered as a bundle over $\mathbb{R}$ with a compact base, corresponding to the "classically allowed region" at energy $E$. The sequence of eigenvalues $E = E_n(h)$ is determined by the condition that the resulting quasi-mode, gluing together asymptotic solutions from different coordinates patches along $\Lambda_E$, be single-valued, i.e. $K_h(E)$ have trivial holonomy.

Assuming $\Lambda_E$ is smoothly embedded in $T^*\mathbb{R}^2$, it can be always be parametrized by a non degenerate phase function. Of particular interest are the critical points of the phase functions, or focal points which are responsible for the change in Maslov index. Recall that $a(E) = (x_E, \xi_E) \in \Lambda_E$ is a focal point if $\Lambda_E$ "turns vertical" at $a(E)$, i.e. $T_{a(E)}\Lambda_E$ is no longer transverse to the fibers $x = \text{Const.}$ in $T^*\mathbb{R}$.

In any case, however, $\Lambda_E$ can be parametrized locally either by a phase $S = S(x)$ (spatial representation) or a phase $S = \tilde{S}(\xi)$ (Fourier representation). Choose an orientation on $\Lambda_E$ and for $a \in \Lambda_E$ (not necessarily a focal point), denote by $\rho = \pm 1$ its oriented segments near $a$. Let $\chi^a \in C_0^\infty(\mathbb{R}^2)$ be a smooth cut-off equal to 1 near $a$, and $\omega_\rho^a$ a small neighborhood of supp[$P, \chi^a$] $\cap \Lambda_E$ near $\rho$. Here $\chi^a$ holds for $\chi^a(x, hD_x)$ as in (2), and we shall equally write $P(x, hD_x)$ (spatial representation) or $P(-hD_\xi, \xi)$ (Fourier representation).

**Definition 2.1.** Let $P$ be self-adjoint, and $u^a, v^a \in K_h(E)$ be supported on $\Lambda_E$. We call

$$W_{\rho}^{a}(u^{a}, \bar{v}^{a}) = \left(\frac{i}{\hbar} [P, \chi^{a}]_{\rho} u^{a} | v^{a}\right) \qquad (3)$$

the microlocal Wronskien of $(u^a, \bar{v}^a)$ in $\omega_\rho^a$. Here $\frac{i}{\hbar} [P, \chi^a]_\rho$ denotes the part of the commutateur supported microlocally on $\omega_\rho^a$.

To understand that terminology, let $P = -h^2 \Delta + V$, $x_E = 0$, and change $\chi$ to Heaviside unit step function $\chi(x)$, depending on $x$ alone. Then in distributional sense, we have $\frac{i}{\hbar} [P, \chi] = -i h \delta' + 2 \delta h D_x$, where $\delta$ denotes the Dirac measure at 0, and $\delta'$ its derivative, so that $(\frac{i}{\hbar} [P, \chi] u|u) = -i h (u'(0) \overline{u(0)} - u(0) \overline{u'(0)})$ is the usual Wronskien of $(u, \bar{u})$.
**Proposition 2.1.** Let $u^a, v^a \in K_h(E)$ as above, and denote by $\hat{u}$ the h-Fourier (unitary) transform of $u$. Then

$$ \left(\frac{i}{h}[P, \chi^a]u^a|v^a\right) = \left(\frac{i}{h}[P, \chi^a]\hat{u}^a|\hat{v}^a\right) = 0 \qquad (4) $$

$$ \left(\frac{i}{h}[P, \chi^a]_+ u^a|v^a\right) = -\left(\frac{i}{h}[P, \chi^a]_- u^a|v^a\right) \qquad (5) $$

(all equalities being understood mod $\mathcal{O}(h^\infty)$, resp $\mathcal{O}(h^{N+1})$ when considering $u^a, v^a \in K_h^N(E)$). Moreover, $W_\rho^a(u^a, \overline{v^a})$ doesn't depend mod $\mathcal{O}(h^\infty)$ (resp $\mathcal{O}(h^{N+1})$) of the choice of $\chi^a$ as above.

*Proof.* Since $u^a, v^a \in K_h(E)$ are distributions in $L^2$, the first equality (4) follows from Plancherel formula and the regularity of microlocal solutions in $L^2$, $p + i$ being elliptic. If $a$ is not a focal point, $u^a, v^a$ are smooth WKB solutions near $a$, so we can expand the commutator in $W = \left(\frac{i}{h}[P, \chi^a]u^a|v^a\right)$ and use that $P$ is self-adjoint to show that $W = \mathcal{O}(h^\infty)$. If $a$ is a focal point, $u^a, v^a$ are smooth WKB solutions in Fourier representation, so again $W = \mathcal{O}(h^\infty)$. Then (5) follows from Definition 2.1. $\square$

We can find a linear combination of $W_\pm^a$, (depending on $a$) which defines a sesquilinear form on $K_h(E)$, so that this Hermitean form makes of $K_h(E)$ a metric bundle, endowed with the gauge group $U(1)$. This linear combination is prescribed as the construction of Maslov index: namely we take $W^a(u^a, \overline{u^a}) = W_+^a(u^a, \overline{u^a}) - W_-^a(u^a, \overline{u^a}) > 0$ when the critical point $a$ of $\pi_{\Lambda_E}$ is traversed in the $-\xi$ direction to the right of the fiber (or equivalently $W^a(u^a, \overline{u^a}) = -W_+^a(u^a, \overline{u^a}) + W_-^a(u^a, \overline{u^a}) > 0$ when $a$ is traversed in the $+\xi$ direction to the left of the fiber). Otherwise, just exchange the signs. When $\gamma_E$ is a convex curve, there are only 2 focal points. In general there may be many focal points $a$, but each jump of Maslov index is compensated at the next focal point which is traversed to the other side of the fiber (co-cycle property). Maslov index is computed mod 4. Our method consists in constructing Gram matrix of a generating system of $K_h(E)$ in a suitable dual basis; its determinant vanishes precisely at the eigenvalues $E(h)$.

### 3. QM and BS in the case of a Schrödinger operator

As a warm-up, we derive the well known BS quantization rule using microlocal Wronskians in case of a potential well, i.e. $\gamma_E$ is convex. Consider the spectrum of Schrödinger operator $P(x, hD_x) = (hD_x)^2 + V(x)$ near the energy level $E_0 < \liminf_{|x|\to\infty} V(x)$, when $\{V \le E\} = [x'_E, x_E]$ and $x'_E, x_E$ are simple turning points, $V(x'_E) = V(x_E) = E$, $V'(x'_E) < 0$, $V'(x_E) > 0$. It is convenient to start the construction from the focal points $a$ or $a'$. We set $a' = x'_E$, $a = x_E$, identifying the focal point $a = a_E = (x_E, 0)$ with its projection $x_E$. We know that microlocal solutions $u$ of $(P - E)u = 0$ near $a$ are of the form

$$ u^a(x; h) = \frac{C}{\sqrt{2}} \left( e^{i \frac{\pi}{4}} (E - V(x))^{-\frac{1}{4}} e^{\frac{i}{\hbar} S(a,x)} + e^{-i \frac{\pi}{4}} (E - V(x))^{-\frac{1}{4}} e^{-\frac{i}{\hbar} S(a,x)} + \mathcal{O}(h) \right) \qquad (6) $$
where $C \in \mathbb{C}$, $S(y, x) = \int_y^x \xi_+(t) dt$, and $\xi_+(t)$ is the positive root of $\xi^2 + V(t) = E$. In the same way, the microlocal solutions of $(P - E)u = 0$ near $a'$ have the form, with $C' \in \mathbb{C}$

$$ u^{a'}(x; h) = \frac{C'}{\sqrt{2}} \left( e^{-i \frac{\pi}{4}} (E-V(x))^{-\frac{1}{4}} e^{\frac{i}{h} S(a', x)} + e^{i \frac{\pi}{4}} (E-V(x))^{-\frac{1}{4}} e^{-\frac{i}{h} S(a', x)} + \mathcal{O}(h) \right) \quad (7) $$

These expressions result in computing by the method of stationary phase the oscillatory integral that gives the solution of $(P(-hD_\xi, \xi) - E)\hat{u} = 0$ in Fourier representation. The change of phase factor $e^{\pm i \frac{\pi}{4}}$ accounts for Maslov index. For the sake of simplicity, we omit henceforth $\mathcal{O}(h)$ terms, but the computations below extend to all order in $h$ (practically, at least for $N = 2$), thus giving the asymptotics of BS. This will be elaborated in the next Section.

The semi-classical distributions $u^a$, $u^{a'}$ span the microlocal kernel $K_h$ of $P - E$ in $(x, \xi) \in]a', a[ \times \mathbb{R}$; they are normalized using microlocal Wronskians as follows.

Let $\chi^a \in C_0^\infty(\mathbb{R}^2)$ as in the Introduction be a smooth cut-off equal to 1 near $a$. Without loss of generality, we can take $\chi^a(x, \xi) = \chi_1^a(x) \chi_2^a(\xi)$, so that $\chi_2 \equiv 1$ on small neighborhoods $\omega_\pm^a$, of $\text{supp}(\frac{i}{h}[P, \chi^a]) \cap \{\xi^2 + V = E\}$ in $\pm\xi > 0$. We define $\chi^{a'}$ similarly. By (6) and (7) we have, mod $\mathcal{O}(h)$:

$$ \frac{i}{h}[P, \chi^a]u^a = \sqrt{2} C (\chi_1^a)'(x) \left( e^{i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{\frac{i}{h} S(a,x)} - e^{-i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{-\frac{i}{h} S(a,x)} \right) $$

$$ \frac{i}{h}[P, \chi^a]u^{a'} = \sqrt{2} C' (\chi_1^{a'})'(x) \left( e^{-i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{\frac{i}{h} S(a',x)} - e^{i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{-\frac{i}{h} S(a',x)} \right) $$

Let

$$ F_{\pm}^{a}(x; h) = \frac{i}{h} [P, \chi^{a}]_{\pm} u^{a}(x, h) = \pm \sqrt{2} C (\chi_{1}^{a})'(x) e^{\pm i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{\pm \frac{i}{h} S(a,x)} $$

so that:

$$ \begin{aligned} (u^a | F_+^a - F_-^a) &= |C|^2 \left( e^{i \frac{\pi}{4}} (E-V)^{-\frac{1}{4}} e^{\frac{i}{h} S(a,x)} \left| (\chi_1^a)'(x) e^{i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{\frac{i}{h} S(a,x)} \right. \right. \\ & \qquad \left. + |C|^2 \left( e^{-i \frac{\pi}{4}} (E-V)^{-\frac{1}{4}} e^{-\frac{i}{h} S(a,x)} \left| (\chi_1^a)'(x) e^{-i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{-\frac{i}{h} S(a,x)} \right. \right. \\ & \qquad \left. \left. = |C|^2 \left( \int (\chi_1^a)'(x) dx + \int (\chi_1^a)'(x) dx \right) + \mathcal{O}(h) \right) = 2|C|^2 + \mathcal{O}(h) \end{aligned} $$

(the mixed terms such as $(e^{i \frac{\pi}{4}} (E-V)^{-\frac{1}{4}} e^{\frac{i}{h} S(a,x)})(\chi_1^a)'(x) e^{-i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{-\frac{i}{h} S(a,x)}$) are $\mathcal{O}(h^\infty)$ because the phase is non stationary), thus $u^a$ is normalized mod $\mathcal{O}(h)$ if we choose $C = 2^{-1/2}$. In the same way, with

$$ F_{\pm}^{a'}(x; h) = \frac{i}{h} [P, \chi^{a'}]_{\pm} u^{a'}(x, h) = \pm \sqrt{2} C' (\chi_1^{a'})'(x) e^{\mp i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{\pm \frac{i}{h} S(a', x)} $$

we get

$$ (u^{a'} | F_{+}^{a'} - F_{-}^{a'}) = |C'|^2 \left( \int (\chi_{1}^{a'})'(x) dx + \int (\chi_{1}^{a'})'(x) dx \right) + \mathcal{O}(h) = -2|C'|^2 + \mathcal{O}(h) $$

and we choose again $C' = C$ which normalizes $u^{a'} \mod \mathcal{O}(h)$. Normalization carries to higher order, as is shown in the next Section.
So there is a natural duality product between $K_h(E)$ and the span of functions $F_+^a - F_-^a$ and $F_+^{a'} - F_-^{a'}$ in $L^2$. As in [18], [13] we can show that this space is microlocally transverse to $\text{Im}(P - E)$ on $(x, \xi) \in ]a', a[ \times \mathbb{R}$, and thus identifies with the microlocal co-kernel $K_h^*(E)$ of $P - E$; in general $\dim K_h(E) = \dim K_h^*(E) = 2$, unless $E$ is an eigenvalue, in which case $\dim K_h = \dim K_h^* = 1$ (showing that $P - E$ is of index 0 when Fredholm.)

Microlocal solutions $u^a$ and $u^{a'}$ extend as smooth solutions on the whole interval $]a', a[$; we denote them by $u_1$ and $u_2$. Since there are no other focal points between $a$ and $a'$, they are expressed by the same formulae (which makes the analysis particularly simple) and satisfy:

$$ (u_1|F_+^a - F_-^a) = 1, \quad (u_2|F_+^{a'} - F_-^{a'}) = -1 $$

Next we compute (still modulo $\mathcal{O}(h)$)

$$ \begin{aligned} (u_1|F_+^{a'} - F_-^{a'}) &= \frac{1}{2} (e^{i \frac{\pi}{4}} (E-V)^{-\frac{1}{4}} e^{\frac{i}{h} S(a,x)} |(\chi_1^{a'})'(x) e^{-i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{\frac{i}{h} S(a',x)}) \\ &\quad + \frac{1}{2} (e^{-i \frac{\pi}{4}} (E-V)^{-\frac{1}{4}} e^{-\frac{i}{h} S(a,x)} |(\chi_1^{a'})'(x) e^{i \frac{\pi}{4}} (E-V)^{\frac{1}{4}} e^{-\frac{i}{h} S(a',x)}) \\ &= \frac{i}{2} e^{-\frac{i}{h} S(a',a)} \int (\chi_1^{a'})'(x) dx - \frac{i}{2} e^{\frac{i}{h} S(a',a)} \int (\chi_1^{a'})'(x) dx \\ &= -\sin(S(a',a)/h) \end{aligned} $$

(taking again into account that the mixed terms are $\mathcal{O}(h^\infty)$). Similarly $(u_2|F_+^a - F_-^a) = \sin(S(a',a)/h)$. Now we define Gram matrix

$$ G^{(a,a')}(E) := \begin{pmatrix} (u_1|F_+^a - F_-^a) & (u_2|F_+^a - F_-^a) \\ (u_1|F_+^{a'} - F_-^{a'}) & (u_2|F_+^{a'} - F_-^{a'}) \end{pmatrix} \quad (8) $$

whose determinant $-1 + \sin^2(S(a',a)/h) = -\cos^2(S(a',a)/h)$ vanishes precisely on eigenvalues of $P$ in $I$, so we recover the well known BS quantization condition

$$ \oint \xi(x) dx = 2 \int_{a'}^a (E - V(x))^{\frac{1}{2}} dx = 2 \pi h \left(k + \frac{1}{2}\right) + \mathcal{O}(h); \quad k \in \mathbb{Z} \quad (9) $$

and $\det G^{(a,a')}(E)$ is nothing but Jost function which is computed e.g. in [DeDi], [9] by another method.

# 4. The general case

By the discussion after Proposition 2.1, it clearly suffices to consider the case when $\gamma_E$ contains only 2 focal points which contribute to Maslov index.

## 4.1. Well normalized quasi-modes mod $\mathcal{O}(h^2)$

Let $a = (x_E, \xi_E)$ be such a focal point. Following a well known procedure we can trace back to [17], we first seek for asymptotic solutions in Fourier representation near $a$ of the form $\hat{u}(\xi) = e^{i\psi(\xi)/h}b(\xi; h)$. Here the phase $\psi = \psi_E$ solves Hamilton-Jacobi
equation $p_0(-\psi'(\xi), \xi) = E$, and can be normalized by $\psi(\xi_E) = 0$; the amplitude $b(\xi; h) = b_0(\xi) + hb_1(\xi) + \dots$ has to be found recursively together with $a(x, \xi; h) = a_0(x, \xi) + ha_1(x, \xi) + \dots$, such that

$$hD_{\xi}(e^{i(x\xi+\psi(\xi))/h}a(x, \xi; h)) = e^{i(x\xi+\psi(\xi))/h}b(\xi; h) \\ \times (p_0(x, \xi) - E + h\tilde{p}_1(x, \xi) + h^2\tilde{p}_2(x, \xi) + \mathcal{O}(h^3))$$

$p_0$ being the principal symbol of $P$, $\tilde{p}_1$ its sub-principal symbol for the standard (Feynman) quantization, etc... Define $\lambda(x, \xi)$ by $p_0(x, \xi) - E = \lambda(x, \xi)(x + \psi'(\xi))$, we get first

$$\lambda(-\psi'(\xi), \xi) = \alpha(\xi) = \partial_x p_0(-\psi'(\xi), \xi)$$

This yields $a_0(x, \xi) = \lambda(x, \xi)b_0(\xi)$ and solving a first order ODE $L(\xi, D_\xi)b_0 = 0$, with $L(\xi, D_\xi) = \alpha(\xi)D_\xi + (\frac{1}{2i}\alpha'(\xi) - p_1(\psi'(\xi), \xi))$ we get

$$b_0(\xi) = C_0|\alpha(\xi)|^{-1/2}e^{i\int_{\alpha}^{p_1}}$$

with an arbitrary constant $C_0$, we take independent of $E$. This gives in turn $a_1(x, \xi) = \lambda(x, \xi)b_1(\xi) + \lambda_0(x, \xi)$, with

$$\lambda_0(x, \xi) = \frac{b_0(\xi)\tilde{p}_1 + i\frac{\partial a_0}{\partial \xi}}{x + \partial_\xi \psi}$$

and $b_1(\xi)$ solution of $L(\xi, D_\xi)b_1 = \tilde{p}_2b_0 + i\partial_\xi\lambda_0|_{x=-\psi'(\xi)}$. We eventually get

$$b_0(\xi) + hb_1(\xi) = (C_0 + hC_1 + hD_1(\xi))|\alpha(\xi)|^{-1/2}e^{i\int_{\alpha}^{p_1}}$$

where we have set

$$D_1(\xi) = \alpha(\xi_E) \int_{\xi_E}^{\xi} (i\tilde{p}_2b_0 - \partial_\xi\lambda_0|_{x=-\psi'(\xi')}) |\alpha(\xi')|^{-1/2} e^{-i\int_{\alpha}^{p_1}} d\xi' \quad (10)$$

The integration constants $C_0, C_1, \dots$ will be determined by normalizing the microlocal Wronskians as follows.

We compute $W_{\rho}^a(u^a, \overline{u^a})$ in Fourier representation, with $\hat{u}(\xi; h) = e^{i\psi(\xi)/h}b(\xi; h)$. Recall $\chi^a \in C_0^\infty(\mathbb{R}^2)$, $\chi^a \equiv 1$ near $a_E$; without loss of generality, we can take $\chi^a(x, \xi) = \chi_1(x)\chi_2(\xi)$, so that $\chi_2 \equiv 1$ on small neighborhoods $\omega_\pm^a$, of $\text{supp}\frac{i}{h}[P, \chi^a] \cap \gamma_E$ in $\pm(\xi - \xi_E) > 0$. Thus we need only consider the variations of $\chi_1$. Weyl symbol of $\frac{i}{h}[P, \chi^a]$ is given by $c(x, \xi; h) = (\partial_\xi p_0(x, \xi) + h\partial_\xi p_1(x, \xi))\chi'_1(x) + \mathcal{O}(h^2)$, so

$$\frac{i}{h}[P, \chi^a]\hat{u}(\xi) = (2\pi h)^{-1} \iint e^{i(-(\xi-\eta)y+\psi(\eta))/h} c(y, \frac{\xi+\eta}{2}; h)(b_0 + hb_1)(\eta) dy d\eta$$

Evaluating by stationary phase, we find $\frac{i}{h}[P, \chi^a]\hat{u}(\xi) = e^{i\psi(\xi)/h}d(\xi; h)$, where $d(\xi; h) = d_0 + hd_1 + \mathcal{O}(h^2)$ with $d_0 = c_0b_0$ and $d_1$ a function of $c_0, c_1, b_0, b_1$ we have determined so far. It follows

$$(\frac{i}{h}[P, \chi^a] + \hat{u}|\hat{u}) = \int_{\xi_E}^{\infty} d_0(\xi)\overline{b_0(\xi)} d\xi + \mathcal{O}(h)$$
The leading term is simply

$$ \int_{\xi_E}^{\infty} |b_0|^2 c_0 d\xi = |C_0|^2 \operatorname{sgn}(\alpha(\xi_E)) \int_{\xi_E}^{\infty} \chi_1'(-\psi'(\xi))\psi''(\xi) d\xi = |C_0|^2 \operatorname{sgn}(\alpha(\xi_E)) \quad (11) $$

so $u^a$ is normalized mod $\mathcal{O}(h)$, provided $\alpha(\xi_E) > 0$, when $2|C_0|^2 = 1$; we take $C_0 = 1/\sqrt{2}$ as in Schrödinger case. Next step in normalization involves the term $D_1(\xi)$ defined in (10); integrating by parts, we can remove $\chi_1(\xi)$ and its second derivative, so to end up with a simple integral like (11). The computation being rather lengthy, we only state the final result:

**Lemma 4.1.** With the hypotheses above, the microlocal Wronskian near a focal point $a_E$ is given by

$$ \begin{aligned} W^a(u^a, \overline{u^a}) &= \left( \frac{i}{h} ([P, \chi^a]_+ - [P, \chi^a]_-) \hat{u}^a \middle| \hat{u}^a \right) \\ &= 2 \operatorname{sgn}(\alpha(\xi_E)) (|C_0|^2 + h (2 \operatorname{Re}(\overline{C_0}C_1) + |C_0|^2 \partial_x \left( \frac{p_1}{\partial_x p_0} \right)(\xi_E)) + \mathcal{O}(h^2)) \end{aligned} $$

The condition that $u^a$ be normalized mod $\mathcal{O}(h^2)$ (once we have chosen $C_0$ to be real), is then

$$ C_1(E) = -\frac{1}{2} C_0 \partial_x \left( \frac{p_1}{\partial_x p_0} \right) (\xi_E) $$

so that now $W^a(u^a, \overline{u^a}) = 2\operatorname{sgn}(\alpha(\xi_E))C_0^2(1 + \mathcal{O}(h^2))$. This procedure carries to any order, we say then that $u^a$ is well-normalized. It can be formalized by considering $\{a_E\}$ as a Poincaré section, and Poisson operator the operator that assigns to the initial condition $C_0$ on $\{a_E\}$, in a unique way, the well-normalized (forward) solution $u^a$ to $(P-E)u^a = 0$; see next Section.

The next task consists in extending the solutions away from $a_E$ in the spatial representation. First we expand $u^a(x) = (2\pi h)^{-1/2} \int e^{i\psi(\xi)/h} b(\xi; h) d\xi$ near $a$ by stationary phase, selecting the 2 critical points $\xi_\pm(x) = \xi_\pm(x)$, that correspond to the phase functions $\varphi_\rho(x) = x\xi_\rho(x) + \psi(\xi_\rho(x))$. So we have

$$ \begin{aligned} u^a(x; h) &= u_+^a(x; h) + u_-^a(x; h) \\ &= \frac{1}{\sqrt{2}} \sum_{\pm} \left( \frac{\partial_\xi p_0(x, \xi_\pm(x))}{i} \right)^{-\frac{1}{2}} e^{\frac{i}{h} S_\pm(x_E, x; h)} (1 + \mathcal{O}(h)) \end{aligned} \quad (12) $$

where

$$ S_{\pm}(x_E, x; h) = \varphi_{\pm}(x_E, x) + h \int_{\xi_E}^{\xi_{\pm}(x)} \frac{p_1(-\psi'(\zeta), \zeta)}{\alpha(\zeta)} d\zeta + \sqrt{2} h^2 \operatorname{Im} \left( D_1(\xi_{\pm}(x)) \right) \quad (13) $$

and $\varphi_\pm(x_E, x) = x_E \xi_E + \int_{x_E}^x \xi_\pm(y) dy$. Then we use standard WKB theory with Ansatz $u_\rho^a(x) = a_\rho(x; h)e^{i\varphi_\rho(x)/h}$. Let $\beta_0(x) = \partial_\xi p_0(x, \varphi'(x)) = -\frac{\alpha(\xi(x))}{\xi'(x)}$, $\beta_1(x) = \partial_\xi p_1(x, \varphi'(x))$. Omitting the index $\rho$, we find the usual half-density

$$ a_0(x) = \tilde{C}_0 |\beta_0(x)|^{-\frac{1}{2}} \exp \left( -i \int \frac{p_1(x, \varphi'(x))}{\beta_0(x)} dx \right) $$
with a new constant $\tilde{C}_0 \in \mathbf{R}$; the next term is

$$a_1(x) = (\tilde{C}_1 + \tilde{D}_1(x)) |\beta_0(x)|^{-\frac{1}{2}} \exp \left( -i \int \frac{p_1(x, \varphi'(x))}{\beta_0(x)} dx \right)$$

and $\tilde{D}_1(x)$ a complex function with

$$\mathrm{Re}(\tilde{D}_1(x)) = -\frac{\tilde{C}_0}{2} \frac{\beta_1(x)}{\beta_0(x)} + \mathrm{Const.}$$

and

$$\mathrm{Im}(\tilde{D}_1(x)) = \tilde{C}_0 \left( \int \frac{\beta_1(x)}{\beta_0^2(x)} p_1(x, \varphi'(x)) dx - \int \frac{p_2(x, \varphi'(x))}{\beta_0(x)} dx \right) \quad (14)$$

## 4.2. The homology class of the generalized action

Here we identify the various terms in (10) and (14). First on $\gamma_E$ we have $\psi(\xi) = \int -x d\xi + \mathrm{Const.}$, and $\varphi(x) = \int \xi dx + \mathrm{Const.}$. By Hamilton equations

$$\dot{\xi}(t) = -\partial_x p_0(x(t), \xi(t)), \quad \dot{x}(t) = \partial_\xi p_0(x(t), \xi(t))$$

so on $\gamma_E$ we have $\int \frac{p_1}{\alpha} d\xi = -\int \frac{p_1}{\partial_\xi p_0} dx = -\int_{\gamma_E} p_1 dt$. The form $p_1 dt$ is called the subprincipal 1-form. Next we consider $D_1(\xi)$ as the integral over $\gamma_E$ of the 1-form, defined near $a$ in Fourier representation as

$$\Omega_1 = T_1 d\xi = \mathrm{sgn}(\alpha(\xi)) (i \tilde{p}_2 b_0 - \partial_\xi \lambda_0) |\alpha|^{-1/2} e^{-i \int \frac{p_1}{\alpha} d\xi}$$

Using WKB constructions, $\Omega_1$ can also be extended in the spatial representation. Since $\gamma_E$ is Lagrangian, $\Omega_1$ is a closed form that we are going to compute modulo exact forms. Using integration by parts, the integral of $\Omega_1(\xi)$ in Fourier representation simplifies to

$$\sqrt{2} \mathrm{Re}(D_1(\xi)) = -\frac{1}{2} \left[ \partial_x \left( \frac{p_1}{\partial_x p_0} \right) (-\psi'(\zeta), \zeta) \right]_{\xi_E}^{\xi} \quad (15)$$

$$\sqrt{2} \mathrm{Im}(D_1(\xi)) = \int_{\xi_E}^{\xi} T_1(\zeta) d\zeta + \left[ \frac{\psi''}{6\alpha} \frac{\partial^3 p_0}{\partial x^3} + \frac{\alpha'}{4\alpha^2} \frac{\partial^2 p_0}{\partial x^2} \right]_{\xi_E}^{\xi} \quad (16)$$

$$T_1(\zeta) = \frac{1}{\alpha} \left( p_2 - \frac{1}{8} \frac{\partial^4 p_0}{\partial x^2 \partial \zeta^2} + \frac{\psi''}{12} \frac{\partial^4 p_0}{\partial x^3 \partial \zeta} + \frac{(\psi'')^2}{24} \frac{\partial^4 p_0}{\partial x^4} \right) + \frac{1}{8} \frac{(\alpha')^2}{\alpha^3} \frac{\partial^2 p_0}{\partial x^2} \\ + \frac{1}{6} \frac{\psi''}{\alpha^2} \frac{\partial^3 p_0}{\partial x^3} - \frac{p_1}{\alpha^2} \left( \partial_x p_1 - \frac{p_1}{2\alpha} \frac{\partial^2 p_0}{\partial x^2} \right) \quad (17)$$

Eq. (15) already shows that $\mathrm{Re}(\Omega_1)$ is exact. We can carry the integration in $x$-variable between the focal points $a_E$ and $a'_E$, and in $\xi$-variable again near $a'_E$. Now let $\Omega(x, \xi) = f(x, \xi) dx + g(x, \xi) d\xi$, where $f(x, \xi), g(x, \xi)$ are any smooth functions on $\mathcal{A}$. By Stokes formula

$$\int_{\gamma_E} \Omega(x, \xi) = \int \int_{\{p_0 \le E\}} (\partial_x g - \partial_\xi f) dx \wedge d\xi$$
where, following [CdV], we have extend $p_0$ inside the disk bounded by $A_-$ so that it coincides with a harmonic oscillator in a neighborhood of the origin ($p_0(0) = 0$, say). Making the symplectic change of coordinates $(x, \xi) \mapsto (t, E)$:

$$ \iint_{\{p_0 \le E\}} (\partial_x g - \partial_\xi f) \, dx \wedge d\xi = \int_0^E \int_0^{T(E')} (\partial_x g - \partial_\xi f) \, dt \wedge dE' \quad (18) $$

where $T(E')$ is the period of the flow of Hamilton vector field $H_{p_0}$ at energy $E'$ ($T(E')$ being a constant near 0). Using these expressions, we recover the well known action integrals (see e.g. [7]):

**Lemma 4.2.** Let $\Gamma dt$ be the restriction to $\gamma_E$ of the 1-form

$$ \omega_0(x, \xi) = \left( \frac{\partial^2 p_0}{\partial x^2} \frac{\partial p_0}{\partial \xi} - \frac{\partial^2 p_0}{\partial x \partial \xi} \frac{\partial p_0}{\partial x} \right) dx + \left( \frac{\partial^2 p_0}{\partial x \partial \xi} \frac{\partial p_0}{\partial \xi} - \frac{\partial^2 p_0}{\partial \xi^2} \frac{\partial p_0}{\partial x} \right) d\xi $$

We have $\mathrm{Re} \oint_{\gamma_E} \Omega_1 = 0$, whereas

$$ \mathrm{Im} \oint_{\gamma_E} \Omega_1 = \frac{1}{48} \left(\frac{d}{dE}\right)^2 \int_{\gamma_E} \Gamma dt - \int_{\gamma_E} p_2 dt - \frac{1}{2} \frac{d}{dE} \int_{\gamma_E} p_1^2 \quad (19) $$

### 4.3. BS quantization rule

Recall from (12) the asymptotic solution $u^a(x; h)$ near $a = a_E$. In the last term of (13) we can substitute (16) with $T_1$ as in (17), that is

$$ \sqrt{2} \mathrm{Im}(D_1(\xi_{\pm}(x))) = \int_{x_E}^x T_1(\xi_{\pm}(y)) \xi'_{\pm}(y) \, dy $$

Similarly, the asymptotic solution near $a' = a'_E$ is given by

$$ \begin{aligned} u^{a'}(x; h) &= u_{+}^{a'}(x; h) + u_{-}^{a'}(x; h) \\ &= \frac{1}{\sqrt{2}} \sum_{\pm} \left( \frac{\partial_{\xi} p_0(x, \xi_{\pm}(x))}{-i} \right)^{-\frac{1}{2}} e^{\frac{i}{\hbar} S_{\pm}(x', x; h)} (1 + \mathcal{O}(h)) \end{aligned} \quad (20) $$

where as in (12), using (16) and (17)

$$ S_{\pm}(x', x; h) = \varphi_{\pm}(x', x) - h \int_{x'_E}^x \frac{p_1(y, \xi_{\pm}(y))}{\partial_{\xi} p_0(y, \xi_{\pm}(y))} dy + h^2 \int_{x'_E}^x T_1(\xi_{\pm}(y)) \xi'_{\pm}(y) dy \quad (21) $$

and similarly for $S_{\pm}(x_E, x; h)$. The semi-classical distributions $u^a, u^{a'}$ are well normalized as in Lemma 4.1. We compute $F_{\pm}^a = \frac{i}{\hbar}[P, \chi^a]_{\pm} u_{\pm}^a$. Still mod $\mathcal{O}(h)$

$$ F_{\pm}^{a}(x; h) = \pm \frac{1}{\sqrt{2}} e^{\pm i \frac{\pi}{4}} \left( \pm \partial_{\xi} p_{0}(x, \xi_{\pm}(x)) \right)^{\frac{1}{2}} e^{\frac{i}{\hbar} S_{\pm}(x_{E}, x; h)} \frac{d\chi_{1}^{a}(x)}{dx}(x) \quad (22) $$

and using that the mixed terms $(u_{\pm}^a | F_{\pm}^a)$ are $\mathcal{O}(h^{\infty})$, we find $(u^a | F_{+}^a - F_{-}^a) \equiv 1 \pmod{\mathcal{O}(h)}$. In the same way, near $a'$ we have $(u^a | F_{+}^a - F_{-}^a) \equiv -1$. The normalized microlocal solutions $u^a$ and $u^{a'}$, uniquely extended along $\gamma_E$, are now called $u_1$ and $u_2$. They verify

$$ (u_1 | F_{+}^{a'} - F_{-}^{a'}) = \frac{i}{2} \left( e^{\frac{i}{\hbar} A_{-}(x_E, x'_E; h)} - e^{\frac{i}{\hbar} A_{+}(x_E, x'_E; h)} \right) \quad (23) $$
$$ (u_2|F_+^a - F_-^a) = \frac{i}{2} \left( e^{-\frac{i}{h} A_-(x_E, x'_E; h)} - e^{-\frac{i}{h} A_+(x_E, x'_E; h)} \right) \quad (24) $$

where the generalized actions are given by

$$ A_{\pm}(x_E, x'_E; h) = S_{\pm}(x_E, x; h) - S_{\pm}(x'_E, x; h) = (x_E - x'_E) \xi_E + \int_{x_E}^{x'_E} \xi_{\pm}(y) dy - h \int_{x_E}^{x'_E} \frac{p_1(y, \xi_{\pm}(y))}{\partial_{\xi} p_0(y, \xi_{\pm}(y))} dy + h^2 \int_{x_E}^{x'_E} T_1(\xi_{\pm}(y)) \xi'_{\pm}(y) dy $$

We have

$$ \int_{x'_E}^{x_E} (\xi_+(y) - \xi_-(y)) dy = \oint_{\gamma_E} \xi(y) dy \\ \int_{x'_E}^{x_E} \left( \frac{p_1(y, \xi_+(y))}{\partial_{\xi} p_0(y, \xi_+(y))} - \frac{p_1(y, \xi_-(y))}{\partial_{\xi} p_0(y, \xi_-(y))} \right) dy = \int_{\gamma_E} p_1 dt \\ \int_{x'_E}^{x_E} \left( T_1(\xi_+(y)) \xi'_+(y) - T_1(\xi_-(y)) \xi'_-(y) \right) dy = \mathrm{Im} \oint_{\gamma_E} \Omega_1(\xi(y)) $$

On the other hand, Gram matrix as in (8) has determinant

$$ - \cos^2 \left( \frac{A_-(x_E, x'_E; h) - A_+(x_E, x'_E; h)}{2h} \right) $$

which vanishes precisely when BS holds. Summing up, we eventually obtain:

**Proposition 4.1.** With the notations and hypotheses stated in the Introduction, BS is given in the interval I by $S_h(E) = 2\pi n h$, $n \in \mathbb{Z}$, where the semi-classical action $S_h(E) \sim S_0(E) + h S_1(E) + h^2 S_2(E) + \cdots$ consists of:

(i) the classical action

$$ S_0(E) = \oint_{\gamma_E} \xi(x) dx = \int \int_{\{p_0 \le E\} \cap W} d\xi \wedge dx $$

(ii) Maslov index and the integral of the sub-principal 1-form $p_1 dt$

$$ S_1(E) = -\pi - \int_{\gamma_E} p_1(x(t), \xi(t)) dt $$

(iii) the second order term

$$ S_2(E) = \frac{1}{24} \frac{d}{dE} \int_{\gamma_E} \Delta dt - \int_{\gamma_E} p_2 dt - \frac{1}{2} \frac{d}{dE} \int_{\gamma_E} p_1^2 dt $$

where

$$ \Delta(x, \xi) = \frac{\partial^2 p_0}{\partial x^2} \frac{\partial^2 p_0}{\partial \xi^2} - \left( \frac{\partial^2 p_0}{\partial x \partial \xi} \right)^2 $$

We recall that $S_3(E) = 0$. Note that the signs in front of the first and third term of our formula for $S_2(E)$ differ from those in [7].
5. The discrete spectrum of P in I

Here we recover the fact that BS determines asymptotically all eigenvalues of P in I. We adapt the argument of [19]. It is convenient to think of {$a_E$} and {$a'_E$} as zero-dimensional "Poincaré sections" of $\gamma_E$. Let $\mathcal{K}^a(E)$ be the operator (Poisson operator) that assigns to its "initial value" $C_0 \in L^2(\{a_E\}) \approx \mathbb{R}$ the well normalized solution $u(x; h) = \int e^{\frac{i}{h}(x\xi + \psi(\xi))} b(\xi; h) d\xi$ to $(P - E)u = 0$ near $\{a_E\}$. By construction, we have:

$$ \pm \mathcal{K}^a(E)^* \frac{i}{h} [P, \chi^a]_{\pm} \mathcal{K}^a(E) = \mathrm{Id}_{a_E} = 1 \quad (25) $$

We define objects "connecting" $a$ to $a'$ along $\gamma_E$ as follows: let $T = T(E) > 0$ such that $\exp TH_{p_0}(a) = a'$. Choose $\chi_f^a$ ($f$ for "forward") be a cut-off function supported microlocally near $\gamma_E$, equal to 0 along $\exp tH_{p_0}(a)$ for $t \le \epsilon$, equal to 1 along $\gamma_E$ for $t \in [2\epsilon, T + \epsilon]$, and back to 0 next to $a'$, e.g. for $t \ge T + 2\epsilon$. Let similarly $\chi_b^a$ ($b$ for "backward") be a cut-off function supported microlocally near $\gamma_E$, equal to 1 along $\exp tH_{p_0}(a)$ for $t \in [-\epsilon, T - 2\epsilon]$, and equal to 0 next to $a'$, e.g. for $t \ge T - \epsilon$. By (25) we have

$$ \mathcal{K}^a(E)^* \frac{i}{h} [P, \chi^a]_+ \mathcal{K}^a(E) = \mathcal{K}^a(E)^* \frac{i}{h} [P, \chi_f^a] \mathcal{K}^a(E) = 1 \quad (26) $$

$$ -\mathcal{K}^a(E)^* \frac{i}{h} [P, \chi^a]_- \mathcal{K}^a(E) = -\mathcal{K}^a(E)^* \frac{i}{h} [P, \chi_b^a] \mathcal{K}^a(E) = 1 \quad (27) $$

which define a left inverse $R_+^a(E) = \mathcal{K}^a(E)^* \frac{i}{h} [P, \chi_f^a]$ to $\mathcal{K}^a(E)$ and a right inverse $R_-^a(E) = -\frac{i}{h} [P, \chi_b^a] \mathcal{K}^a(E)$ to $\mathcal{K}^a(E)^*$. We define similar objects connecting $a'$ to $a$, $T' = T'(E) > 0$ such that $\exp T'H_{p_0}(a) = a'$ ($T = T'$ if $p_0$ is invariant by time reversal), in particular a left inverse $R_+^{a'}(E) = \mathcal{K}^{a'}(E)^* \frac{i}{h} [P, \chi_f^{a'}]_+$ to $\mathcal{K}^{a'}(E)$ and a right inverse $R_-^{a'}(E) = -\frac{i}{h} [P, \chi_b^{a'}] \mathcal{K}^{a'}(E)$ to $\mathcal{K}^{a'}(E)^*$, with the additional requirement

$$ \chi_b^a + \chi_b^{a'} = 1 \quad (28) $$

near $\gamma_E$. Define now the pair $R_+(E)u = (R_+^a(E)u, R_+^{a'}(E)u)$, $u \in L^2(\mathbb{R})$ and $R_-(E)$ by $R_-(E)u_- = R_-^a(E)u_- + R_-^{a'}(E)u_-$, $u_- = (u_-^a, u_-^{a'}) \in \mathbb{C}^2$, we call Grushin operator $\mathcal{P}(z)$ the operator defined by the linear system

$$ \frac{i}{h} (P - z)u + R_-(z)u_- = v, \quad R_+(z)u = v_+ \quad (29) $$

From [19], we know that the problem (29) is well posed, and

$$ \mathcal{P}(z)^{-1} = \begin{pmatrix} E(z) & E_+(z) \\ E_-(z) & E_-(z) \end{pmatrix} $$

with $(P - z)^{-1} = E(z) - E_+(z) E_-(z)^{-1} E_-(z)$. Actually one can show that the effective Hamiltonian $E_+(z)$ is Gram matrix (8). There follows that the spectrum of P in I is precisely the set of E we have determined by BS quantization rule
6. References

[1] P. ARGYRES, “The Bohr-Sommerfeld quantization rule and Weyl correspondence”, *Physic*, vol. 2, 1965.

[2] H. BAKLOUTI, “Asymptotique des largeurs de résonances pour un modèle d’effet tunnel microlocal”, *Ann.IHP (Phys.Th)*, vol. 68, num. 2, 1998.

[3] C. BENDER, S.ORZSAG, “Advanced Mathematical Methods for Scientists and Engineers”, Springer, 1979.

[4] A. BENSOUISSI, A. Ifa, M. ROULEUX, “Andreev reflection and the semi-classical Bogoliubov-de Gennes Hamiltonian”, *Proceedings "Days of Diffraction 2009"*, Saint-Petersburg p.37-42, IEEE 2009.

[5] A.BENSOUISSI, N. M'HADBI, M. ROULEUX, “Andreev reflection and the semi-classical Bogoliubov- de Gennes Hamiltonian: resonant states”, *Proceedings Days of Diffraction 2011*, Saint-Petersburg. p.39-44. IEEE 101109/DD.2011.6094362, 2011.

[6] M. CARGO, A. GRACIA-SAZ, R. LITTLEJOHN, M. REINSCH, P. DE RIOS "Moyal star product approach to the Bohr-Sommerfeld approximation", J.Phys.A:Math and Gen.38, 2005.

[7] Y. COLIN DE VERDIÈRE, “Bohr-sommerfeld Rules to All orders”, Ann. H.Poincaré, vol. 6, 2005.

[8] E. DELABAERE, F. PHAM "Resurgence methods in semi-classical asymptotics", Ann. Inst. H.Poincaré, vol. 71, num. 1, 1999.

[9] E. DELABAERE, H. DILLINGER, F. PHAM "Exact semi-classical expansions for 1-D quantum oscillators", J.Math.Phys, vol. 38, num. 12, 1997.

[10] K.P. DUNCAN, B.L. GYÖRFFY "Semiclassical theory of quasiparticles in superconducting state", *Annals of Phys*, vol. 298, 2002.

[11] A. GRACIA-SAZ "The symbol of a function of a pseudo-differential operator", Ann. Inst. Fourier, vol. 55, num. 7, 2005.

[12] B. HELFFER, D. ROBERT "Puits de potentiel généralisés et asymptotique semi-classique", Annales Inst. H.Poincaré, vol. 41, num. 3, 1984.

[13] B. HELFFER, J. SJÖSTRAND "Semi-classical analysis for Harper's equation III", Mémoire No 39, Soc. Math. de France, vol. 117, num. 4, 1988.

[14] R. LITTLEJOHN "Lie Algebraic Approach to Higher-Order Terms III", Preprint 17p, 2013.

[15] F.W.J. OLVER "Introduction to Asymptotics and Special Functions", Academic Press, New York, 1974.

[16] M. ROULEUX "Tunneling effects for h-Pseudodifferential Operators, Feshbach resonances and the Born-Oppenheimer approximation in: Evolution Equations, Feshbach Resonances, Singular Hodge Theory", Adv. Part. Diff. Eq Wiley-VCH, 1999.

[17] J. SJÖSTRAND "Analytic singularities of solutions of boundary value problems", Proc. NATO ASI on Singularities in boundary value problems, D.Reidel, 1980.

[18] J. SJÖSTRAND "Density of states oscillations for magnetic Schrödinger operators, in: Bennewitz (ed.)", Diff. Eq. Math. Phys, Univ. Alabama, Birmingham, 1990.

[19] J. SJÖSTRAND, M. ZWORSKI "Quantum monodromy and semi-classical trace formula", J. Math. Pures Appl, vol. 81, num. 9, 2002.

[20] A. VOROS "Asymptotic h-expansions of stationary quantum states", Ann. Inst. H. Poincaré Sect. A(N.S), vol. 26, 1977.

[21] D.R. YAFAEV "The semiclassical limit of eigenfunctions of the Schrödinger equation and the Bohr-Sommerfeld quantization condition, revisited. ", Algebra i Analiz, vol. 2, num. 6, 2010.