# How to obtain large order perturbation series for cubic anharmonic oscillator?

Consider the potential

$$V(x)= \frac{x^2}{2} + gx^3.\tag{1}$$

Then the time-independent Schrödinger equation becomes

$$\left(-\frac{1}{2}\frac{d^2}{dx^2} + \frac{x^2}{2} + gx^3 \right)\psi = E(g) \psi.\tag{2}$$

Where $$E(g)$$ is the energy eigenvalues as a function of parameter $$g$$. One obtains the following perturbation series expansion

$$E(g) = \frac{1}{2} - \frac{11}{8} g^2 - \frac{465}{32} g^4 - \frac{39708}{128} g^6 -\ldots\tag{3}$$ for the ground state energy eigenvalue.

Now, how does one obtain such a series? I am aware of the standard Rayleigh–Schrödinger perturbation series which allows one to compute till the first few orders after which calculations become too messy. I am also aware of Feynman Path Integral way of computing ground state energy eigenvalues but no luck finding large order perturbation series terms. A brief outline of proceeding will be helpful.

• Hi @Suroj Dey: Series expansion taken from which reference? Commented May 29, 2020 at 21:00
• Table (3.1) of the book 'Path Integrals in Quantum Mechanics, Statistics, Polymer Physics, and Financial Markets' by Hagen Kleinert Commented Jun 14, 2022 at 19:20

1. Here we will for fun try to reproduce the first few terms in a perturbative series for the ground state energy $$E_0$$ of the 1D TISE \begin{align} \hat{H}\psi_0~=~&E_0\psi_0, \cr H~=~&\frac{p^2}{2}+\frac{\omega^2}{2}q^2+V_{\rm int}(q), \cr V_{\rm int}(q)~=~&gq^3, \cr g~=~&\frac{\lambda}{6},\end{align} \tag{A} using an Euclidean path integral in 0+1D \begin{align} e^{W_c[J]/\hbar}~=~&Z[J]\cr ~=~&\int_{q(0)=q(T)}^{\dot{q}(0)=\dot{q}(T)}\!{\cal D}q ~\exp\left\{ \frac{1}{\hbar}\int_{[0,T]}\!dt~(-L_E +J q) \right\} \cr ~=~&\exp\left\{-\frac{1}{\hbar}\int_{[0,T]}\!dt~V_{\rm int}\left(\hbar \frac{\delta}{\delta J}\right) \right\} Z_2[J],\end{align} \tag{B} cf. Refs. 1-3. The Euclidean Lagrangian is \begin{align} L_E~=~&\frac{1}{2}\dot{q}^2+\frac{\omega^2}{2}q^2+V_{\rm int}(q), \cr q(T)~=~&q(0),\qquad \dot{q}(0)~=~\dot{q}(T), \end{align} \tag{C} with periodic boundary conditions$$^1$$.

2. The free quadratic part is the harmonic oscillator (HO) \begin{align}Z_2[J]&~=~\cr Z_2[J\!=\!0]&\exp\left\{\frac{1}{2\hbar}\iint_{[0,T]^2}\!dt~dt^{\prime} J(t) \Delta(t,t^{\prime})J(t^{\prime}) \right\}.\end{align} \tag{D} The partition function for the HO can be calculated either via path integrals (e.g. as a functional determinant) or via its definition in statistical physics:
\begin{align} Z_2[J\!=\!0] ~=~&{\rm Tr}\left(e^{-\hat{H}(g=0)~T/\hbar}\right)\cr ~=~&\sum_{n\in\mathbb{N}_0}e^{-(n+1/2)\omega T}\cr ~=~&\left(2\sinh\frac{\omega T}{2}\right)^{-1}.\end{align} \tag{E} The free propagator/Greens function with periodic boundary conditions is \begin{align} \Delta(t,t^{\prime}) ~=~&\frac{1}{2\omega}\sum_{n\in\mathbb{Z}} e^{-\omega |t-t^{\prime}+nT|}\cr ~=~&\frac{1}{2\omega}e^{-\omega |t-t^{\prime}|}+\frac{1}{\omega}\frac{\cosh\omega(t-t^{\prime})}{e^{\omega T}-1}\cr ~=~&\frac{1}{2\omega}\left(\frac{e^{-\omega |t-t^{\prime}|}}{1-e^{-\omega T}} +\frac{e^{\omega |t-t^{\prime}|}}{e^{\omega T}-1}\right), \cr \left(-\frac{d^2}{dt^2}+\omega^2\right)\Delta(t,t^{\prime}) ~=~&\delta(t\!-\!t^{\prime}+T\mathbb{Z})~=:~{\rm III}_T(t\!-\!t^{\prime}),\end{align}\tag{F} where $${\rm III}_T$$ is the Dirac comb/shah function. Eq. (F) becomes in the $$T\to\infty$$ limit \begin{align} \Delta(t,t^{\prime})~=~&\frac{1}{2\omega}e^{-\omega |t-t^{\prime}|}, \cr \left(-\frac{d^2}{dt^2}+\omega^2\right)\Delta(t,t^{\prime}) ~=~&\delta(t\!-\!t^{\prime}).\end{align}\tag{G}

3. The main idea is to use the fact that the ground state energy can be inferred from the connected vacuum bubbles$$^2$$
$$W_c[J\!=\!0]~\stackrel{(I)}{\sim}~ -E_0(g) T\quad\text{for}\quad T\to\infty.\tag{H}$$ Here we are using the linked-cluster theorem$$^3$$ \begin{align} e^{W_c[J=0]/\hbar} ~=~&Z[J\!=\!0] ~=~{\rm Tr}\left(e^{-\hat{H}T/\hbar}\right)\cr ~=~&\sum_{n\in\mathbb{N}_0}e^{-E_n(g)T/\hbar} \end{align} \tag{I} (We should stress that the free vacuum part $$e^{W_{c,2}[J=0]/\hbar} =Z_2[J\!=\!0]$$ in eq. (E) does not have an interpretation in terms of Feynman diagrams.) The Feynman rule for the cubic vertex is $$-\frac{\lambda}{\hbar}\int_{[0,T]}\!dt. \tag{J}$$

4. The unperturbed zero-mode energy for the HO is the well-known $$E_0(g\!=\!0)~=~E_0^{(0)}~=~\frac{\hbar\omega}{2},\tag{K}$$ cf. eq. (E). Set $$\hbar=1=\omega$$ to compare with OP's eq. (3).

5. There are two 2-loop vacuum bubbles: the dumbbell diagram $$O\!\!-\!\!O$$ and the sunset diagram $$\theta$$ with symmetry factor $$S=8$$ and $$S=12$$, respectively, cf. Fig. 1.

$$\uparrow$$ Fig. 1. (From Ref. 4.) The two 2-loop vacuum bubbles: the dumbbell diagram $$O\!\!-\!\!O$$ and the sunset diagram $$\theta$$.

They make up the 2nd-order contribution $$g^2E_0^{(2)}~=~-\frac{11\hbar^2 g^2}{8\omega^4} \tag{L}$$ to the ground state energy $$E_0(g)$$, cf. OP's eq. (3) & Ref. 5.

Proof of eq. (L): The dumbbell Feynman diagram $$O\!\!-\!\!O$$ is$$^4$$ \begin{align} \frac{(-\lambda/\hbar)^2}{8} \left(\frac{\hbar}{2\omega}\right)^{3}& \iint_{[0,T]^2}\!dt~dt^{\prime} \left(\frac{e^{-\omega |t-t^{\prime}|}}{1-e^{-\omega T}} +\frac{e^{\omega |t-t^{\prime}|}}{e^{\omega T}-1}\right) \coth^2\frac{\omega T}{2}\cr ~=~& \frac{\hbar\lambda^2}{32\omega^4}T\coth^2\frac{\omega T}{2}\cr ~=~& \frac{9\hbar g^2}{8\omega^4}T\coth^2\frac{\omega T}{2}.\end{align}\tag{M} The sunset Feynman diagram $$\theta$$ is$$^3$$ \begin{align} \frac{(-\lambda/\hbar)^2}{12} \left(\frac{\hbar}{2\omega}\right)^{3}& \iint_{[0,T]^2}\!dt~dt^{\prime} \left(\frac{e^{-\omega |t-t^{\prime}|}}{1-e^{-\omega T}} +\frac{e^{\omega |t-t^{\prime}|}}{e^{\omega T}-1}\right)^3\cr ~=~& \frac{\hbar\lambda^2}{144 \omega^4}T\left(1+\frac{3}{\sinh^2\frac{\omega T}{2}}\right)\cr ~=~& \frac{\hbar g^2}{4\omega^4}T\left(1+\frac{3}{\sinh^2\frac{\omega T}{2}}\right).\end{align}\tag{N} The sum of the dumbbell and sunset diagrams becomes $$\frac{\hbar g^2}{8\omega^4}T\left(11+\frac{15}{\sinh^2\frac{\omega T}{2}}\right) ~=~ \frac{11\hbar g^2}{8\omega^4}T + {\cal O}(T^0) .\tag{O}$$ To obtain eq. (L) compare eqs. (H) and (O). $$\Box$$

6. Five 3-loop vacuum bubbles make up the 4th-order contribution, cf. Fig. 2.

$$\uparrow$$ Fig. 2. (From Ref. 4.) The five 3-loop vacuum bubbles.

7. It is in principle possible to calculate to any order by drawing Feynman diagrams. An $$n$$-loop-integral is analytically doable by breaking the integration region $$[0,T]^n$$ into $$n$$-simplexes.

References:

1. M. Marino, Lectures on non-perturbative effects in large $$N$$ gauge theories, matrix models and strings, arXiv:1206.6272; section 3.1.

2. R. Rattazzi, The Path Integral approach to Quantum Mechanics, Lecture Notes for Quantum Mechanics IV, 2009; subsection 2.3.6.

3. R. MacKenzie. Path Integral Methods and Applications, arXiv:quant-ph/0004090, section 6.

4. M. Srednicki, QFT, 2007; figures 9.1 + 9.2. A prepublication draft PDF file is available here.

5. I. Gahramanov & K. Tezgin, A resurgence analysis for cubic and quartic anharmonic potentials, arXiv:1608.08119, eqs. (3.1) + (3.2).

--

$$^1$$ We need 2 boundary conditions since the EOM is a 2nd-order ODE. It is natural to assume that the trajectories $$t\mapsto q(t)$$ are differentiable.

$$^2$$ On one hand, in the frequency space a connected Feynman diagram with $$n$$ external legs is of the form

$${\cal M}(k^0_1,\ldots, k^0_n)\delta(k^0_1\!+\!\ldots\!+\!k^0_n)\tag{P}$$

due to time translation invariance. Moreover $${\cal M}(k^0_1,\ldots, k^0_n)$$ is finite in QM, cf. e.g. my Phys.SE answer here. In particular each vacuum bubble is proportional to the time period $$\delta(0)\sim T$$. On the other hand, in the time domain the propagator is bounded

$$|\Delta(t,t^{\prime})|\leq \Delta(0,0)=\frac{1}{2\omega}\coth\frac{\omega T}{2},\tag{Q}$$

so Feynman diagrams in the time domain are finite if the time region $$[0,T]$$ is compact.

$$^3$$ It is in principle possible to derive all the energy-levels

$$E_n(g)~=~\underbrace{E_n^{(0)}}_{=(n+1/2)\hbar\omega}+ g\underbrace{E_n^{(1)}}_{=0} + g^2E_n^{(2)} +g^3\underbrace{E_n^{(3)}}_{=0} +{\cal O}(g^4)\tag{R}$$

from the identity (I)

\begin{align}\sum_{n\in\mathbb{N}_0}&e^{-E_n(g)T/\hbar} \cr ~\stackrel{(I)+{O}}{=}&\left(2\sinh\frac{\omega T}{2}\right)^{-1}\exp\left\{\frac{\hbar g^2}{8\omega^4}T\left(11+\frac{15}{\sinh^2\frac{\omega T}{2}}\right)+{\cal O}(g^4)\right\}\cr ~=~&e^{-\omega T/2}\left(\sum_{k\in\mathbb{N}_0}e^{-k\omega T}\right) \left\{1+\frac{\hbar g^2}{8\omega^4}T\left(11+ 60\sum_{\ell\in\mathbb{N}_0}\ell e^{-\ell\omega T}\right)+{\cal O}(g^4)\right\} \end{align}\tag{S}

as a power series in $$g$$ and $$e^{-\omega T}$$. It's a fun exercise to check that the identity (S) agrees with the 2nd-order correction from standard Rayleigh–Schrödinger (RS) perturbation theory \begin{align} g^2 E_n^{(2)} ~=~& \sum_{k\neq n} \frac{|\langle k^{(0)}| \hat{V}_{\rm int}| n^{(0)}\rangle |^2}{E_n^{(0)}-E_k^{(0)}}\cr ~=~& \frac{\hbar^2g^2}{8\omega^4} \left(\frac{n(n-1)(n-2)}{3} + 9n^3 -9(n+1)^3 -\frac{(n+3)(n+2)(n+1)}{3} \right)\cr ~=~& -\frac{\hbar^2g^2}{8\omega^4} \left(30n^2 +30n +11\right), \end{align} \tag{T} as it should.

$$^4$$ It is straightforward to check that \begin{align} \iint_{[0,T]^2}\!dt~dt^{\prime}e^{-\omega |t-t^{\prime}|} ~=~& 2\int_0^T \!dt\int_0^t\!dt^{\prime}e^{-\omega (t-t^{\prime})}\cr ~=~&\frac{2}{\omega}\left(T -\frac{1-e^{-\omega T}}{\omega} \right)\cr ~=~&\frac{2T}{\omega} + {\cal O}(T^0).\end{align} \tag{U}

• Notes for later: Strong Dirichlet BC $q(-T/2)=0=q(T/2)$: Greens function: $2\omega\Delta_D(t,t^{\prime})=e^{-\omega |t-t^{\prime}|}+\frac{e^{-\omega T}\cosh\omega(t-t^{\prime})-\cosh\omega(t+t^{\prime})}{\sinh\omega T}$. Commented Jun 29, 2020 at 15:39
• Notes for later: 2nd-order ODE. Periodic BC $q(0)=q(T)$ +no edge $\Rightarrow$ Differentiable periodic profile $q(t)=q(t+T) \Leftrightarrow q(0)=q(T) \wedge \dot{q}(0)=\dot{q}(T)$: Greens function: $2\omega\Delta_P(t,t^{\prime}) =\sum_{n\in\mathbb{Z}} e^{-\omega |t-t^{\prime}+nT|} =e^{-\omega |t-t^{\prime}|}+\frac{e^{-\omega T/2}\cosh\omega(t-t^{\prime})}{\sinh(\omega T/2)}$ $=e^{-\omega |t-t^{\prime}|}+\frac{2\cosh\omega(t-t^{\prime})}{e^{\omega T}-1}$. Commented Nov 17, 2022 at 14:18
• Notes for later: Difference: $\Delta_P(t,t^{\prime})-\Delta_D(t,t^{\prime}) =\frac{\cosh\omega(t+t^{\prime})+\cosh\omega(t-t^{\prime})}{2\omega\sinh\omega T} =\frac{\cosh\omega t\cosh\omega t^{\prime}}{\omega\sinh\omega T}$. Commented Nov 17, 2022 at 17:55
• physics.stackexchange.com/q/328439/2451 $\quad K^E(q_2,t_2;q_1,t_1)_J =\exp\left\{-\frac{1}{\hbar}\int_{[t_1,t_2]}\!dt~V_{\rm int}\left(\hbar \frac{\delta}{\delta J}\right) \right\} K^E_2(q_2,t_2;q_1,t_1)_J$; $\quad K^E_2(q_2,t_2;q_1,t_1)_J = K^E_2(0,t_2;0,t_1)_J\exp\left\{-\frac{1}{\hbar}S^E_{2,J}[q_{\rm cl}] \right\}$; $\quad K^E_2(0,t_2;0,t_1)_J =K^E_2(0,t_2;0,t_1)_{J=0}$ $\exp\left\{\frac{1}{2\hbar}\iint_{[t_1,t_2]^2}\!dt~dt^{\prime} J(t) \Delta_D(t,t^{\prime})J(t^{\prime}) \right\}$; $\quad K^E_2(0,t_2;0,t_1)_{J=0}=\sqrt{\frac{m\omega}{2\pi\hbar \sinh(\omega\Delta t) }}$; Commented Feb 19, 2023 at 12:21
• Notes for later: $\quad q_{\rm cl}(t) = \frac{q_2\sinh \omega (t-t_1)+q_1\sinh \omega (t_2-t)}{\sinh (\omega \Delta t)}$; $\quad S^E_{2,J}[q_{\rm cl}] = \frac{m\omega}{2}\left((q_2^2+q_1^2)\coth(\omega\Delta t)-\frac{2q_2q_1}{\sinh(\omega\Delta t)}\right)-\int_{[t_1,t_2]}\!dt~J(t)q_{\rm cl}(t)$; Commented Jun 23, 2023 at 6:32

In your example, you can do this analytically as the unperturbed potential in the harmonic oscillator, for which there are analytical solutions for the eigenenergies and eigenfunctions.

## General formulae

Use perturbation theory and a sensible choice of your unperturbed basis.

The energy $$E_n$$ will be written as: $$E_n = E_n^{(0)} + E_n^{(1)} + E_n^{(2)} + \dots ,$$ where the LHS is the true value and the LHS terms are $$n^{\mathrm{th}}$$ corrections.

Same thing for the wavefunction $$\psi_n$$: $$\psi_n = \psi_n^{(0)} + \psi_n^{(1)} + \psi_n^{(2)} + \dots .$$

The first order correction for the energy $$E_n^{(1)}$$ is given by: $$E_n^{(1)} = \int \psi_n^{(0)*} \hat{H}^{(1)} \psi_n^{(0)},$$ the second order is: $$E_n^{(2)} = \int \psi_n^{(0)*} \hat{H}^{(1)} \psi_n^{(1)},$$ and so on.

The first order correction for the wavefunction is: $$\psi_n^{(1)} = \sum_{i\neq n} \psi_n^{(0)} \frac{\int \psi_i^{(0)*} \hat{H}^{(1)} \psi_n^{(0)}}{E_n^{(0)} - E_i^{(0)}} .$$

$$H^{(1)}$$ is the perturbation to the Hamiltonian.

You can combine the previous two formulae to re-write the second-order energy correction as: $$E_n^{(2)} = \sum_{i\neq n} \frac{| \int \psi_i^{(0)*} \hat{H}^{(1)} \psi_n^{(0)}|^2}{E_n^{(0)} - E_i^{(0)}} .$$

I am going to use $$\hat{H}_1 = gx^3$$ for the perturbed potential. The unperturbed potential is the harmonic potential so will use the analytical eigenergies and eigenfunctions for $$E_n^0$$ and $$\psi^0_n$$.

In your example, let's take $$x^2/2$$ to the the "basis" potential, which is none other than the harmonic oscillator which known and analytical wavefunctions (so that you can compute the integrals easily). The perturbed contribution $$H^{(1)}$$ is $$gx^3$$.

So the total energy $$E$$ will be the energy from the harmonic oscillator plus the correction due to the perturbed potential.

For the ground state then ($$n=0$$):

• Zero-th order: $$E^{(0)} = \frac{1}{2},$$ which is the harmonic (unperturbed potential) contribution.

• First order: $$E^{(1)} = 0,$$ because $$gx^3$$ is odd.

• Second order: I did this with Mathematica using the last equation from the previous section. It converges after 3 terms in the sum and it gives: $$E^{(2)} = -1.375 g^2 = -\frac{11}{8}g^2.$$

Then you go on.

# --- Addition: ---

## Higher-order

I then copied all the 3rd, 4th and higher terms from wikipedia into Mathematica and got:

• Third order: $$E^{(3)} = 0.$$

• Fourth order (converged in 8 terms): $$E^{(4)} = -14.5313 g^4 = -\frac{465}{32}g^4.$$

Have you read the basic paper: "Anharmonic Oscillator. II. A Study of Perturbation Theory in Large Order," by Carl M. Bender and Tai Tsun Wu, Phys. Rev. D 7, 1620 (1973)?

That paper shows you how to obtain the series for the quartic case. We were able to find several hundred terms for other QM perurbation series by their method.

See the appendix in: J.Reeve, M. Stone "Late terms in the asymptotic expansion for the energy levels of a periodic potential" Phys.Rev.D 18 (1978) 4746

I have not done it for the cubic, but it is probably not too hard.

• Commented May 29, 2020 at 20:55
1. In this answer, we will sketch the derivation of a large-$$n$$ expression$$^1$$ for the $$a_n$$-coefficients in the perturbative series $$E_0(g)~=~\sum_{n\in\mathbb{N_0}} a_n g^n \tag{A}$$ for the ground state energy by considering a single instanton bounce.

2. The Euclidean Lagrangian is \begin{align}L_E~=~&\frac{1}{2}\dot{q}^2+V(q), \cr q(T)~=~&q(0),\end{align} \tag{B} with periodic boundary condition. Let us be a bit more general and assume that the potential is of the form \begin{align} V(q)~=~&V_2(q)+V_{\rm int}(q),\cr V_2(q)~=~&\frac{1}{2}\omega^2q^2, \cr V_{\rm int}(q)~=~&g~{\rm sgn}(q) {\cal V}(|q|)\cr ~=~&-V_{\rm int}(-q),\cr 0~\leq~& {\cal V}(|q|)~=~{\cal O}(|q|^3).\end{align} \tag{C} This means the potential is unstable for $$q\to\pm\infty$$, and therefore has an instanton bounce for $$q>0$$ (or $$q<0$$) if $$g<0$$ (or $$g>0$$), respectively, cf. e.g. Refs. 1-3 and this related Phys.SE post.

3. The Euclidean partition function $$Z$$ is defined by analytical continuation of the position variable $$q$$ in the complex plane. Remarkably we only need to assume that this is possible; the details are largely irrelevant. Since initially all the parameters in the Euclidean partition function $$Z$$ are manifestly real, we expect the analytical continuation to satisfy the Schwarz reflection principle $$2i{\rm Im} Z(g)~=~Z(g+i0^+)-Z(g-i0^+).\tag{D}$$

4. Let us focus on the case $$g<0$$ so that the instanton $$q_{\rm cl}\geq 0$$ is in the positive $$q$$-direction. The classical energy is conserved: \begin{align}\frac{1}{2}\dot{q}_{\rm cl}^2-V(q_{\rm cl})~=~&E_{\rm cl} \cr\cr \Updownarrow ~& \cr\cr \pm\dot{q}_{\rm cl}~=~&v(q_{\rm cl})\cr ~=~&\sqrt{2(V(q_{\rm cl})+E_{\rm cl})}.\end{align}\tag{E} The classical action contribution for the instanton bounce is \begin{align} \frac{S_{\rm cl}}{2} ~\stackrel{(B)}{=}~& \int_0^{T/2}\! dt\left(\frac{1}{2}\dot{q}_{\rm cl}^2+V(q_{\rm cl})\right)\cr ~\stackrel{(E)}{=}~& \int_{q_-}^{q_+} \! \frac{dq}{v(q)}\left(2V(q)+E_{\rm cl}\right) \cr ~\stackrel{(E)}{=}~& \int_{q_-}^{q_+} \! dq\frac{2V(q)+E_{\rm cl}}{\sqrt{2(V(q)+E_{\rm cl})}}\cr ~\stackrel{E_{\rm cl}\lesssim 0}{\simeq}&~\int_0^{q_+} \! dq\sqrt{2V(q)}.\end{align}\tag{F}

5. A simple scaling argument shows that the coupling constant $$g$$ effectively plays the role of $$\hbar$$ in the WKB/stationary phase approximation. The leading asymptotic WKB approximation for $$g\to 0$$ of the Euclidean partition function is \begin{align} {\rm Re}Z(T) ~\stackrel{\rm WKB}{\sim}&~ \frac{1}{2\sinh\frac{\omega T}{2}}\cr~\approx~& e^{-\omega T/2}\end{align}\tag{G} and \begin{align} {\rm Im}Z(T) ~\stackrel{\rm WKB}{\sim}&~\cr\cr-\frac{T}{2}\sqrt{\frac{1}{2\pi \hbar}\frac{\partial E_{\rm cl}}{\partial T}}& \exp\left(-\frac{S_{\rm cl}-TE_{\rm cl}}{\hbar}\right)~\ll ~1.\end{align}\tag{39.95} Eq. (39.95) from Ref. 1 follows essentially via the Gelfand-Yaglom formula for functional determinants, cf. e.g. this related Phys.SE post.

6. The (imaginary part of the) ground state energy is \begin{align}-{\rm Im}E_0 ~=~&\lim_{T\to \infty}\frac{{\rm Im} W_c(T)}{T}\cr ~=~&\lim_{T\to \infty}\frac{\hbar~{\rm Im} \ln Z(T)}{T}\cr ~\approx~&\lim_{T\to \infty}\frac{\hbar~{\rm Im}Z(T)}{T~{\rm Re}Z(T)}.\end{align}\tag{H} The sought-for $$a_n$$-coefficients can then in principle be extracted from \begin{align} a_n ~=~&\frac{1}{2\pi i}\oint_{\mathbb{R}} \! \mathrm{d}g \frac{E_0(g)}{g^{n+1}}\cr ~=~&\frac{1}{2\pi i}\int_{\mathbb{R}} \! \mathrm{d}g \frac{E_0(g-i0^+)-E_0(g+i0^+)}{g^{n+1}}\cr ~=~& -\frac{1}{\pi}\int_{\mathbb{R}} \! \mathrm{d}g \frac{{\rm Im} E_0(g)}{g^{n+1}}.\end{align}\tag{I} This is the main answer to OP's question.

7. Let us provide some more details. We are interested in the limit $$T\to \infty$$ where $$q_-\to 0$$ and \begin{align} 0 ~<~& -E_{\rm cl}\cr ~=~& V(q_+)\cr ~=~& V(q_-) \cr ~\simeq~& V_2(q_-)\cr ~=~& \frac{1}{2}\omega^2 q_-^2\to 0.\end{align}\tag{J} We calculate \begin{align} \frac{T}{2} ~=~& \int_{q_-}^{q_+} \! \frac{dq}{v(q)}\cr ~\stackrel{(E)}{=}~& \int_{q_-}^{q_+} \! \frac{dq}{\sqrt{2(V(q)+E_{\rm cl})}}\cr ~\stackrel{(L)+(M)}{\simeq}&~ \frac{I+I_2}{\omega},\end{align}\tag{K} where \begin{align} I ~=~& \int_0^{q_+}\!dq~\left(\frac{\omega}{\sqrt{2V(q)}} -\frac{\omega}{\sqrt{2V_2(q)}}\right)\cr ~=~& \int_0^{q_+}\!dq~\left(\frac{\omega}{\sqrt{2V(q)}} -\frac{1}{q}\right) \end{align}\tag{L} and where \begin{align} I_2 ~=~& \omega\int_{q_-}^{q_+} \! \frac{dq}{\sqrt{2(V_2(q)+E_{\rm cl})}}\cr ~=~& \int_{q_-}^{q_+} \! \frac{dq}{\sqrt{q^2 -q_-^2}} \cr ~=~& \left[ \ln\left(\sqrt{q^2 -q_-^2}+q \right)\right]_{q_-}^{q_+}\cr ~\simeq~&\ln\frac{2q_+}{q_-}.\end{align}\tag{M} This implies \begin{align} q_-~\stackrel{(M)}{\simeq}~&2q_+e^{-I_2} \cr\cr \Downarrow ~& \cr\cr -E_{\rm cl} ~\simeq~& \frac{1}{2}\omega^2 q_-^2\cr ~\stackrel{(K)}{\simeq}~& 2\omega^2q_+^2e^{2I-\omega T} \cr\cr \Downarrow ~& \cr\cr \frac{\partial E_{\rm cl}}{\partial T} ~\simeq~&-\omega E_{\rm cl}.\end{align}\tag{N}
This leads to the main formula \begin{align} {\rm Im}E_0~\stackrel{(39.95)+(N)}{\sim}&~\cr\cr -\frac{q_+\omega^{3/2}}{2}\sqrt{\frac{\hbar}{\pi}} & e^I\exp\left(-\frac{1}{\hbar} S_{\rm cl}\right).\end{align} \tag{39.103}

8. Example: odd quartic interaction potential: \begin{align} {\cal V}(|q|)~=~& \frac{1}{4}|q|^4 \cr\cr \Downarrow ~& \cr\cr q_+~=~&\omega\sqrt{\frac{2}{|g|}},\end{align}\tag{O} \begin{align} S_{\rm cl}~\stackrel{(F)}{=}~& 2\int_0^{q_+}\!dq~\sqrt{2V(q)} \cr ~=~& \int_0^{q^2_+}\!d(q^2)~\sqrt{\omega^2 - \frac{|g|}{2}q^2} \cr ~=~& \frac{2\omega^2}{|g|}\int_0^1\!du~\sqrt{1-u}\cr ~=~& \frac{4\omega^2}{3|g|},\end{align}\tag{P} \begin{align} I ~\stackrel{(L)}{=}~& \int_0^{q_+}\!dq~\left(\frac{\omega}{\sqrt{2V(q)}} -\frac{1}{q}\right)\cr ~=~ & \int_0^{q_+}\!\frac{dq}{q}\left(\frac{1}{\sqrt{1 -\frac{|g|}{2\omega^2}q^2}} -1\right)\cr ~=~& \left[-\ln\left(\sqrt{1 -\frac{|g|}{2\omega^2}q^2} +1\right)\right]_0^{q_+} \cr ~=~& \ln 2,\end{align}\tag{Q} \begin{align} {\rm Im}E_0~\stackrel{(39.103)}{\sim}&~\cr\cr -\sqrt{\frac{2\hbar}{\pi |g|}}\omega^{5/2}&\exp\left(-\frac{4\omega^2}{3|g|\hbar}\right).\end{align}\tag{R}

9. Example: odd cubic interaction potential: \begin{align} {\cal V}(|q|)~=~& |q|^3 \cr\cr \Downarrow ~& \cr\cr q_+~=~&\frac{\omega^2}{2|g|}, \end{align}\tag{S} \begin{align} S_{\rm cl}~\stackrel{(F)}{=}~& 2\int_0^{q_+}\!dq~\sqrt{2V(q)} \cr ~=~&2 \int_0^{q_+}\!dq~q\sqrt{\omega^2 -2|g|q} \cr ~=~& \frac{1}{2g^2} \int_0^{\omega^2}\!du~(\omega^2 - u)\sqrt{u}\cr ~=~& \frac{2\omega^5}{15g^2}, \end{align}\tag{T} \begin{align} I ~\stackrel{(L)}{=}~& \int_0^{q_+}\!dq~\left(\frac{\omega}{\sqrt{2V(q)}} -\frac{1}{q}\right)\cr ~=~& \int_0^{q_+}\!\frac{dq}{q}\left(\frac{1}{\sqrt{1 -\frac{2|g|}{\omega^2}q}} -1\right) \cr ~=~& \left[-2\ln\left(\sqrt{1 -\frac{2|g|}{\omega^2}q} +1\right)\right]_0^{q_+}\cr ~=~& \ln 4,\end{align}\tag{U} \begin{align} {\rm Im}E_0~\stackrel{(39.103)}{\sim}&~\cr\cr -\frac{\omega^{7/2}}{|g|}\sqrt{\frac{\hbar}{\pi}}&\exp\left(-\frac{2\omega^2}{15|g|^2\hbar}\right).\end{align}\tag{V}

References:

1. J. Zinn-Justin, QFT & Critical Phenomena, 2002; chapter 39.

2. M. Marino, Instantons and large $$N$$; chapter 2.

3. M. Marino, Lectures on non-perturbative effects in large $$N$$ gauge theories, matrix models and strings, arXiv:1206.6272; section 3.1.

4. I. Gahramanov & K. Tezgin, A resurgence analysis for cubic and quartic anharmonic potentials, arXiv:1608.08119; section 3.

--

$$^1$$ It is in principle possible to derive exact expressions for finite $$n$$ by considering multi-instantons and their interactions, cf. e.g. Ref. 4.