Recall from our previous work that we wish to solve the general energy eigenvalue problem with \(\hat{H} |n\rangle = E_n |n\rangle\), with the Hamiltonian able to be proken into an unperturbed and perturbed part via \(\hat{H} = \hat{H}_0 + \hat{V}\), and the unperturbed energy eigenvalue problem is given by \(\hat{H}_0 |n\rangle_0 = E_n^0 |n\rangle_0\). This implies that \(\hat{H} |n\rangle = E_n |n\rangle\) can be written as \[(E_n -\hat{H}_0) |n\rangle = \hat{V}|n\rangle.\] Therefore, let’s define, as before the two following projection operators (onto and perpendicular to the unperturbed ground state): \[\hat{P}_n |n\rangle_0 \, {}_0\langle n| ~~\text{and}~~ \hat{Q}_n = 1- \hat{P}_n.\] Now recall that we previously found that \[\hat{Q}_n \lvert n \rangle = \frac{\hat{Q}_n}{E_n - \hat{H}_0} \hat{V} \lvert n \rangle.\] However, beware that there is a caveat to this—namely that \(\hat{Q}_n\) projects \(\frac{1}{E_n - \hat{H}_0}\) perpendicular to the unperturbed ground state only when \(\hat{V} = 0\), because the denominator has an \(E_n\) in it instead of an \(E_N^{(0)}\). Now note that: \[\lvert n \rangle = (\hat{P}_n + \hat{Q}_n) \lvert n \rangle = \hat{P}_n \lvert n \rangle + \frac{\hat{Q}_n }{E_n - \hat{H}_0} \hat{V}\lvert n \rangle.\] So we write: \[\bigg(1 - \frac{\hat{Q}_n}{E_n - \hat{H}_0} \hat{V}\bigg) \lvert n \rangle = \hat{P}_n \lvert n \rangle = \lvert n \rangle_0 \, _0\langle n | n \rangle = |n \rangle_0,\] using the same normalization as we did before. Then we have that \[\lvert n \rangle = \left[1 - \frac{\hat{Q}_n \hat{V}}{E_n - \hat{H}_0}\right]^{-1} \lvert n \rangle_0= \lvert n \rangle_0 + \frac{\hat{Q}_n \hat{V} }{E_n - \hat{H}_0} \lvert n \rangle_0 + \frac{\hat{Q}_n \hat{V} \hat{Q}_n \hat{V} }{(E_n - \hat{H}_0)^2} \lvert n \rangle_0+ \cdots\]
To find \(E_n\), multiply \((E_n - \hat{H}_0) \lvert n \rangle = \hat{V} \lvert n \rangle\) by \(_0\langle n |\) from the left to find that \[E_n - E_n^{(0)} = {}_0\langle n \lvert \hat{V} \lvert n \rangle.\] So we find that \[\begin{aligned} E_n &= E_n^{(0)} + {}_0\langle n \lvert \hat{V} \lvert n \rangle _0 + {}_0\langle n \lvert \hat{V} \frac{\hat{Q}_n}{E_n - \hat{H}_0} \hat{V} \lvert n \rangle_0 + \cdots \\ &= E_n^{(0)} + V_{nn} + \sum_{m \neq n} \frac{V_{nm} V_{mn}}{E_n - E_m^{(0)}} + \cdots \end{aligned}\] Note that \(E_n\) appears on both the RHS and LHS. This generates a new equation for \(E_n\). But the series is much simpler than for the non-degenerate case, and we did not need to assume that the system was non-degenerate. On the downside, it is often less accurate than Rayleigh-Schrödinger perturbation theory. The full series for \(\lvert n \rangle\) becomes: \[\lvert n \rangle = \sum_{m=0}^\infty \left( \frac{\hat{Q}_n \hat{V}}{E_n - \hat{H}_0} \right)^m \lvert n\rangle_0.\] For \(E_n\), we have: \[\begin{aligned} E_n = &E_n^{(0)} + V_{nn} + \sum_{m \neq n} \frac{V_{nm} V_{mn}}{E_n - E_m^{(0)}} + \sum_{m \neq n} \sum_{m' \neq n} \frac{V_{nm} V_{mm'} V_{m'n}}{(E_n - E_m^{(0)})(E_n - E_{m'}^{(0)})} + \cdots \\ &+ \sum_{m_1 \neq n} \sum_{m_2 \neq n} \cdots \sum_{m_\ell \neq n} \frac{V_{nm_1} V_{m_1m_2} V_{m_2m_3} \cdots V_{m_{\ell-1}m_\ell} V_{m_\ell n}} {(E_n - E_{m_1}^{(0)})(E_n - E_{m_2}^{(0)}) \cdots (E_n - E_{m_\ell}^{(0)})}. \end{aligned}\] This leads to a high-order polynomial equation for \(E_n\). Usually, only one root is the physical root in this polynomial equation.
Let’s consider an example to see this process in action. Consider the following one-dimensional Hamiltonian: \[\hat{H}_0 = \frac{\hat{p}^2}{2m} + \frac{1}{2} k \hat{x}^2\] where \(\hat{V} = c \hat{x}\). This corresponds to a linear shift of the harmonic-oscillator Hamiltonian. The full Hamiltonian is then: \[\hat{H} = \frac{\hat{p}^2}{2m} + \frac{1}{2} k \hat{x}^2 + c \hat{x} = \frac{\hat{p}^2}{2m} + \frac{1}{2} k \left(\hat{x} + \frac{c}{k}\right)^2 - \frac{c^2}{2k}.\] Now, define \(\hat{x}' = \hat{x} + \frac{c}{k}\) to find that \[[\hat{x}', \hat{p}] = i\hbar\] and \[\hat{H} = \frac{\hat{p}^2}{2m} + \frac{1}{2} k \hat{x}'^2 - \frac{c^2}{2k}\] and we get that \[E_n = \hbar \omega \left(n + \frac{1}{2}\right) - \frac{c^2}{2k}.\] This means the energy is shifted to second order only!
Now, we calculate energies with perturbation theory. First we determine all of the quantities we will need. We rewrite \(\hat{V}\) as \(\hat{V}=c \hat{x} = c \sqrt{\frac{\hbar}{2m\omega}} \left(\hat{a} + \hat{a}^\dagger\right)\). Then, we find that \(V_{nn}=\langle n \lvert c \hat{x} \lvert n \rangle = 0\), \(V_{mn}=0\) unless \(m=n\pm 1\). Then, we have that \[V_{n,n+1} = c \sqrt{\frac{\hbar}{2m\omega}} \sqrt{n+1} ~\text{and}~ V_{n,n-1} = c \sqrt{\frac{\hbar}{2m\omega}} \sqrt{n}.\] We also have that \(E_n^{(0)} - E_{n+1}^{(0)} = -\hbar\omega\) and \(E_n^{(0)} - E_{n-1}^{(0)} = \hbar\omega\).
We now compute the Rayleigh-Schrödinger perturbation theory: \[\begin{aligned} E_n &= \hbar \omega \left(n + \frac{1}{2}\right) + 0 + \frac{\lvert V_{n,n+1} \rvert^2}{E_n^{(0)} - E_{n+1}^{(0)}} + \frac{\lvert V_{n,n-1} \rvert^2}{E_n^{(0)} - E_{n-1}^{(0)}} \\ &= \hbar \omega \left(n + \frac{1}{2}\right) + \frac{c^2}{2m\omega} \frac{\hbar}{\hbar\omega} \left(-(n+1) + n\right) \\ &= \hbar \omega \left(n + \frac{1}{2}\right) - \frac{c^2}{2m\omega^2}. \end{aligned}\] Hence, \(E_n = \hbar \omega (n +\tfrac{1}{2}) - \frac{c^2}{2k}\). This agrees with the exact answer.
As a further check, let’s look at the third-order correction: \[\Delta E_n^{(3)} = \sum_m' \sum_{m'}' \frac{V_{nm} V_{mm'} V_{m'n}}{(E_n^{(0)} - E_m^{(0)})(E_n^{(0)} - E_{m'}^{(0)})} - V_{nn} \sum_m' \frac{\lvert V_{nm} \rvert^2}{(E_n^{(0)} - E_m^{(0)})^2}.\] But, \(V_{nn}=0\) and \(V_{nm} V_{mm'} V_{m'n} =0\) since \(m=n\pm1\) and \(m' = n\pm 1\) in all cases implies that we must also have that \(V_{mm'}=0\). Thus, \(\Delta E_n^{(3)} =0\) and similarly for all \(m\geq 3\), we have \(\Delta E_n^{(m)} =0\).
Let’s now turn to Wigner-Brillouin perturbation theory. Here we have: \[E_n = \hbar \omega \left(n + \frac{1}{2}\right) + 0 + \frac{c^2 \hbar}{2m\omega} \left[\frac{n+1}{E_n - \hbar \omega (n+\frac{3}{2})} + \frac{n}{E_n - \hbar \omega (n-\frac{1}{2})}\right]\] Multiply by \((E_n - \hbar \omega (n+\frac{3}{2}))(E_n - \hbar \omega (n-\frac{1}{2}))\) to get: \[\begin{aligned} E_n (E_n - \hbar \omega (n+\tfrac{3}{2}))&(E_n - \hbar \omega (n-\tfrac{1}{2})) \\ &= \hbar \omega (n+\tfrac{1}{2}) (E_n - \hbar \omega (n+\tfrac{3}{2}))(E_n - \hbar \omega (n-\tfrac{1}{2})) \\ &+ \frac{c^2 \hbar}{2m\omega} \bigg[\left(n+1\right)(E_n - \hbar \omega (n-\tfrac{1}{2})) + n(E_n - \hbar \omega (n+\tfrac{3}{2}))\bigg] \end{aligned}\] Expanding, we find that \[\begin{aligned} E_n^3 &+ E_n^2[-\hbar \omega (2n+1) - \hbar\omega(n+\tfrac{1}{2})] \\ &- E_n \left[(\hbar \omega)^2 (n^2 + n - \tfrac{3}{4}+2n^2 +2n+\tfrac{1}{2}) - \frac{c^2 \hbar}{2m\omega}(2n+1)\right] \\ &-(\hbar \omega)^3 \left[n^3 + \tfrac{3}{2}n^2 - \tfrac{1}{4}n - \tfrac{3}{8}\right] + \frac{c^2 \hbar}{2m\omega} \hbar \omega \left[n^2 + \tfrac{1}{2}n - \tfrac{1}{2} + n^2 + \tfrac{3}{2}n\right] = 0. \end{aligned}\] Simplifying, we find that \[\begin{aligned} E_n^3 &+ E_n^2 \left[-\hbar \omega (3n + \tfrac{3}{2})\right] + E_n \left[(\hbar \omega)^2 \left(3n^2 + 3n - \tfrac{1}{4}\right) - \frac{c^2 \hbar}{2m\omega} (2n+1)\right] \\ &- (\hbar \omega)^3 \left[n^3 + \tfrac{3}{2}n^2 - \tfrac{1}{4}n - \tfrac{3}{8}\right] + \frac{c^2 \hbar^2}{2m\omega} (2n^2 + 2n - \frac{1}{2}) = 0. \end{aligned}\] This should factorize if it gives the exact answer. But, we have that \[\begin{aligned} \biggr(E_n &-\left. \hbar \omega \left(n + \frac{1}{2}\right) + \frac{c^2}{2m\omega^2}\right) \\ & \times \left(E_n^2 + E_n \left[-\hbar \omega (2n+1) - \frac{c^2}{2m\omega^2}\right] \right. \\ &\left. ~~~~~~~~~~~~~+ (\hbar \omega)^2 \left(n^2 + n - \frac{3}{4}\right) - \frac{c^2 \hbar}{2m\omega} \left(n + \frac{1}{2}\right) + \frac{c^4}{(2m\omega^2)^2}\right) \end{aligned}\] is off by an extra term \(\frac{c^6}{(2m\omega^2)^3}\), which means that the result will have an error of order \(\mathcal{O} \left(\frac{c^6}{(2m\omega^2)^3}\right)\). Also note that we have 3 roots, not one. This means that some roots ofthe Wigner-Brillouin perturbation theory are unphysical. In general, one often finds that Wigner-Brillouin perturbation theory is less accurate than Rayleigh-Schrödinger perturbation theory, as shown here, which further implies that higher order terms must contribute to cancel extra terms and ultimately give us the exact answer.
Lastly, let’s check the energy for \(n=0\), where we can solve the polynomial equation exactly: \[E_0 = \hbar \omega \frac{1}{2} + \frac{c^2 \hbar}{2m\omega} \frac{1}{E_0 - \frac{3}{2} \hbar \omega},\] This then becomes \[\begin{aligned} E_0^2 - E_0 \cdot \frac{3}{2} \hbar \omega &= \hbar \omega \left(E_0 - \frac{3}{2} \hbar \omega\right) + \frac{c^2 \hbar}{2m\omega}, \\ 0&=E_0^2 - E_0 \cdot 2\hbar \omega + \frac{3}{4}(\hbar \omega)^2 - \frac{c^2 \hbar}{2m\omega} , \end{aligned}\] The solutions are \[\begin{aligned} E_0 &= \hbar \omega + \frac{1}{2} \sqrt{4 (\hbar \omega)^2 - 3 (\hbar \omega)^2 + \frac{2c^2 \hbar}{m\omega}}, \\ &= \hbar \omega \pm \frac{1}{2} \sqrt{(\hbar \omega)^2 + \frac{2c^2 \hbar}{m\omega}}. \end{aligned}\] Now take the root: \[\begin{aligned} E_0 &= \hbar \omega - \frac{1}{2} \hbar \omega \sqrt{1 + \frac{2c^2}{k \hbar \omega}} \\ &= \hbar \omega - \frac{1}{2} \hbar \omega \left(1 + \frac{c^2}{k \hbar \omega} - \frac{4}{8} \left(\frac{c^2}{k \hbar \omega}\right)^2 + \cdots\right) \\ &= -\frac{1}{2} \hbar \omega - \frac{c^2}{2k} + \frac{1}{4} \frac{c^4}{k^2 \hbar \omega} + \cdots \end{aligned}\] Note the error that we have at order \(V^4\)! This will be canceled by higher-order terms.