Skip to content

15. Resolvent Operator

PDF pages 653–684

15.2.1 Energy-Space Green Functions

Chapter 15 Resolvent Operator We will now develop a method of calculation based on the resolvent operator.1 This method is quite formal but also quite powerful, especially in scattering problems and problems involving transitions to continua. 15.1 Definition The resolvent operator is defined in the complex plane in terms of the time-independent Hamiltonian H by

\[ G(z) := \]

z −H . (15.1) (resolvent operator) Expanding into a basis |\alpha\rangle of eigenstates of H (i.e., multiplying by the identity in this basis on either side of the resolvent), we see that

\[ G(z) = \]

X \alpha

\[ |\alpha\rangle \langle \alpha| \]

z −E\alpha . (15.2) Thus, we see that G(z) has a simple pole at each eigenvalue E\alpha of H. In the case of a continuum of states, where the eigenvalues cluster together into a continuous interval along the real line, G(z) has instead a branch cut. For example, for a bound system with ionized states, G(z) has poles at each bound-state energy, but a branch cut beginning at the smallest ionization energy and extending to +\inftyalong the real axis. There are some subtleties in dealing with branch cuts in G(z) that we will return to later. However, away from the real axis, the resolvent operator is analytic (in the sense that all its matrix elements are analytic functions off the real axis). 15.2 Green Functions for the Schrödinger Equation 15.2.1 Energy-Space Green Functions Consider the following function, defined by the integral expression

\[ G+(E) := −lim \]
\[ \delta\rightarrow 0+ \]

i ¯h Z \infty

\[ d\tau ei(E−H)\tau/¯he−\tau\delta/¯h, \]

(15.3) 1The oft-cited canonical references are Marvin L. Goldberger and Kenneth M. Watson, Collision Theory (Wiley, 1964); and Albert Messiah, Quantum Mechanics (Wiley, 1958) (ISBN: 0486409244). However, a detailed description with many examples appears in Claude Cohen-Tannoudji, Jacques Dupont-Roc, and Gilbert Grynberg, Atom–Photon Interactions: Basic Processes and Applications (Wiley, 1992), Chapter III (ISBN: 0471625566). Also, in the context of spontaneous emission, see the section on the ‘‘Goldberger–Watson method’’ in G. S. Agarwal, ‘‘Quantum Statistical Theories of Spontaneous Emission and their Relation to Other Approaches,’’ Springer Tracts Modern Physics 70, 1 (1974) (doi: 10.1007/BFb0042382). Another nice, more general introduction is given by Doron Cohen, ‘‘Lecture Notes in Quantum Mechanics,’’ arXiv.org preprint (arXiv: quant- ph/0605180v3) (2008). For the resolvent method applied to the master equation, see P. Lambropoulos, ‘‘Spectral Line Shape in the Presence of Weak Collisions and Intense Fields,’’ Physical Review 164, 84 (1967) (doi: 10.1103/PhysRev.164.84).

15.2.2 Time-Dependent Green Functions and Propagators

Chapter 15. Resolvent Operator where the \delta factor is inserted to guarantee convergence of the integral. Carrying out this integral,

\[ G+(E) = −lim \]
\[ \delta\rightarrow 0+ \]

i ¯h Z \infty

\[ d\tau ei(E−H+i\delta)\tau/¯h \]

= −lim

\[ \delta\rightarrow 0+ \]
\[ ei(E−H+i\delta)\tau/¯h \]
\[ E −H + i\delta \]

\infty = lim

\[ \delta\rightarrow 0+ \]
\[ E −H + i\delta \]

=

\[ E −H + i0+ \]
\[ = G(E + i0+). \]

(15.4) Similarly, we can define the function

\[ G−(E) := lim \]
\[ \delta\rightarrow 0+ \]

i ¯h Z 0 −\infty

\[ d\tau ei(E−H)\tau/¯he+\delta\tau/¯h, \]

(15.5) which becomes

\[ G−(E) = lim \]
\[ \delta\rightarrow 0+ \]

i ¯h Z 0 −\infty

\[ d\tau ei(E−H−i\delta)\tau/¯h \]

= E −H −i0+

\[ = G(E −i0+). \]

(15.6) For reasons we will see, G+(E) is called the retarded Green function, in energy (frequency) space, while G−(E) is called the advanced Green function in energy space. We have thus shown that both Green functions are related to the resolvent via

\[ G\pm(E) = G(E \pm i0+) = \]

E −H \pm i0+ . (retarded and advanced Green functions) (15.7) That is, they are essentially the resolvent along the line displaced infinitesimally above and below the real axis. These differ due to the singular nature of the resolvent along the real axis. 15.2.2 Time-Dependent Green Functions and Propagators Now note that the definition (15.3) can be rewritten

\[ G+(E) = 1 \]

i¯h Z \infty −\infty

\[ d\tau ei(E+i0+)\tau/¯h U(\tau, 0) \Theta(\tau), \]

(15.8)

\[ where \Theta(\tau) is the Heaviside step function and U(\tau, 0) = e−iH\tau/¯h is the unitary time-evolution operator from \]

time 0 to \tau for evolution under the time-independent Hamiltonian H. That is, G+(E) is (within a spe- cific normalization convention) the Fourier transform of ‘‘half’’ of the time-evolution operator, U(\tau, 0)\Theta(\tau). Similarly, from the definition (15.5), we see that

\[ G−(E) = 1 \]

i¯h Z \infty −\infty

\[ d\tau ei(E−i0+)\tau/¯h [−U(\tau, 0) \Theta(−\tau)], \]

(15.9) so that up to a minus sign, G−(E) is the Fourier transform of the ‘‘other half’’ of the time-evolution operator, U(\tau, 0) \Theta(−\tau). Thus, defining the time-dependent Green functions

\[ G\pm(t, t0) := \pmU(t, t0) \Theta[\pm(t −t0)], \]

(retarded and advanced Green functions) (15.10)

15.2 Green Functions for the Schrödinger Equation these are related by the above Green functions by a Fourier transform:

\[ G\pm(E) = 1 \]

i¯h Z \infty −\infty

\[ d\tau ei[E+(sgn\tau)i0+]\tau/¯h G\pm(\tau, 0). \]

(Green-function Fourier transform) (15.11) Since the energy-space Green function G\pm(E) is the Fourier transform of the time-domain function G\pm(\tau, 0), the former is commonly written ˜G\pm(E). However, since we will usually be in the energy domain, here we will just use the argument of the Green function to determine in which space it lives. Note that since we are dealing with time-independent systems, G\pm(t, t0) only depends on t −t0. Why the terminology of advanced and retarded Green functions? First, recall [Eq. (4.36)] that the time-evolution operator satisfies the Schrödinger equation,

\[ i¯h\partial tU(t, t0) = HU(t, t0). \]

(15.12)

\[ We can use this relation and \partial \tau\Theta(\tau) = \delta(\tau) to differentiate the Green functions: \]
\[ i¯h\partial tG\pm(t, t0) = \pmi¯h\partial t \]

 U(t, t0)\Theta[\pm(t −t0)]

\[ = \pmHU(t, t0)\Theta[\pm(t −t0)] \pm i¯hU(t, t0) \]

\pm \delta(t −t0)

\[ = HG\pm(t, t0) + i¯h\delta(t −t0). \]

(15.13)

\[ In the last step, we used U(t0, t0) = 1. Thus, we have shown that \]
\[ (i¯h\partial t −H)G\pm(t, t0) = i¯h\delta(t −t0). \]

(Green functions for the Schrödinger equation) (15.14)

\[ Thus, G\pm(t, t0) is the solution to the Schrödinger equation, ‘‘driven’’ by a delta-function impulse at t = t0. In \]

particular, G+(t, t0) is the ‘‘retarded’’ Green function, because the ‘‘source’’ is in the past, and the response follows after the impulse, t > t0. Similarly, G−(t, t0) is the ‘‘advanced’’ Green function, because the source is in the future, and the response comes before the impulse, t < t0. Both Green functions obey the same equation, but correspond to different boundary conditions. Inverting the Fourier-transform relation (15.11) gives

\[ G\pm(\tau, 0) = −1 \]

2\pii Z \infty −\infty dE e−iE\tau/¯h G\pm(E). (Green-function inverse Fourier transform) (15.15) In particular, the case of G+(\tau, 0) is important, as it gives the time-evolution operator for evolving the system forward from t = 0 to \tau:

\[ U(\tau, 0) = −1 \]

2\pii Z \infty −\infty

\[ dE e−iE\tau/¯h G+(E) \]

(\tau > 0). (Green-function Fourier relation) (15.16) This relation is particularly useful in that it shows that matrix elements of the evolution operator,

\[ K(\beta, t; \alpha, t0) := \langle \beta|U(t, t0)|\alpha\rangle , \]

(15.17) (propagator) collectively called the propagator, can be computed from matrix elements of the retarded Green function. The propagator gives the transition amplitude, or the probability amplitude for the system to be in state \beta at time t, given that it was in state \alpha at the (earlier) time t0. Note that in the case of forward propagation t > t0, the propagator can also be regarded as the collection of matrix elements of of the retarded Green function G+(t, t0) in view of the definition (15.10).

15.2.3 Relation to Laplace Transform

Chapter 15. Resolvent Operator 15.2.3 Relation to Laplace Transform We have seen what amounts to this formalism before, when solving the optical Bloch equations via the Laplace-transform method. Recall (Section 5.5.2.1) that the Laplace transform of a time derivative is

\[ L [ ˙y(t)] = sL [y(t)] −y(0), \]

(15.18) where the Laplace transform is defined in general by

\[ L [y](s) := \]

Z \infty dt e−st y(t). (15.19) The Laplace transform of \partial tU(t, t0) is then

\[ L [\partial tU(t, t0)] = sL [U(t, t0)] −U(t0, t0) = sL [U(t, t0)] −1. \]

(15.20) Then the Laplace transform of Eq. (15.12) reads

\[ i¯hsL [U(t, t0)] −i¯h = HL [U(t, t0)], \]

(15.21) or

\[ (i¯hs −H)L [U(t, t0)] = i¯h, \]

(15.22) so that the Laplace transform of the evolution operator becomes

\[ i¯hL [U(t, t0)] = \]

i¯hs −H . (15.23)

\[ Comparing this to the definition (15.1) of the resolvent, we can identify z = i¯hs as the rescaled coordinate, \]

and G(z) is proportional to the Laplace transform of the evolution operator, but with a rescaled coordinate, rotated along the imaginary axis. The propagator relation (15.16) is essentially the inverse Laplace transform for \tau > 0, representing a convenient method to algebraically solve the initial-value problem for the propagator. 15.3 Transitions Between Discrete States Consider a quantum system described by Hamiltonian

\[ H = H0 + V, \]

(15.24) where H0 is the unperturbed Hamiltonian, with interaction V causing transitions among the eigenstates of H0. We then have two resolvents, one corresponding to the perturbed Hamiltonian and given by the original definition (15.1), and one corresponding to the unperturbed system,

\[ G0(z) := \]

z −H0 . (15.25) (unperturbed resolvent operator)

\[ It is convenient to relate these two operators. Starting with B = A + (B −A) for arbitrary A and B, we \]

multiply on the left by B−1 and on the right by A−1 to obtain the identity A = 1 B + 1 B (B −A) 1 A. (15.26) Letting A = z −H0 −V and B = z −H0, we find

\[ G(z) = G0(z) + G0(z)V G(z). \]

(15.27) (perturbed resolvent)

15.3.1 Example: Rabi Oscillations

15.3 Transitions Between Discrete States This relation may be used directly, or it may be iterated to obtain a perturbation series in terms of the unperturbed resolvent:

\[ G = G0 + G0V G0 + G0V G0V G0 + G0V G0V G0V G0 + \cdot \cdot \cdot . \]

(15.28) (perturbed resolvent) In the nonperturbative case, we may take matrix elements of the relation (15.27) in the basis of eigenstates of H0. The diagonal matrix elements become

\[ \langle \alpha|G(z)|\alpha\rangle = \langle \alpha|G0(z)|\alpha\rangle + \langle \alpha|G0(z)V G(z)|\alpha\rangle \]

= z −E\alpha + z −E\alpha X j

\[ \langle \alpha|V |j\rangle \langle j|G(z)|\alpha\rangle , \]

(15.29) while the off-diagonal elements become

\[ \langle \beta|G(z)|\alpha\rangle = \langle \beta|G0(z)|\alpha\rangle + \langle \beta|G0(z)V G(z)|\alpha\rangle \]

= z −E\beta X j

\[ \langle \beta|V |j\rangle \langle j|G(z)|\alpha\rangle , \]

(15.30) for \beta̸ = \alpha. Writing these more compactly,

\[ (z −E\alpha)G\alpha\alpha(z) = 1 + \]

X j

\[ V\alphajGj\alpha(z) \]
\[ (z −E\beta)G\beta\alpha(z) = \]

X j

\[ V\betajGj\alpha(z) \]
\[ (\beta̸ = \alpha). \]

(15.31) (resolvent matrix elements) The strategy then is to solve the algebraic equations, and then compute the inverse Fourier transform for the resolvent matrix elements to obtain the propagator. 15.3.1 Example: Rabi Oscillations As an example, we revisit the Rabi-flopping problem in the two-level atom from Section 5.2.2. In the rotating frame of the laser field, the free atomic Hamiltonian is

\[ H0 = −¯h∆|e\rangle \langle e|, \]

(15.32) while the laser-driven interaction is V = ¯hΩ 

\[ \sigma + \sigma\dagger \]

= ¯hΩ 

\[ |g\rangle \langle e| + |e\rangle \langle g| \]

 . (15.33) Then the relations Eqs. (15.31) for the resolvent matrix elements become

\[ (z + ¯h∆)Gee = 1 + ¯hΩ \]

2 Gge

\[ zGgg = 1 + ¯hΩ \]

2 Geg

\[ (z + ¯h∆)Geg = ¯hΩ \]

2 Ggg zGge = ¯hΩ 2 Gee. (15.34) Taking the third relation and using the second relation to eliminate Ggg decouples the relations and gives a relation for Geg alone:

\[ (z + ¯h∆)Geg = ¯hΩ \]

2 Ggg = ¯hΩ 2z  1 + ¯hΩ 2 Geg  . (15.35)

Chapter 15. Resolvent Operator Solving for Geg, we find

\[ Geg(z) = \]

¯hΩ/2

\[ z2 + ¯h∆z −¯h2Ω2/4. \]

(15.36) Factoring the quadratic denominator,

\[ Geg(z) = \]

¯hΩ/2

z + ¯h∆ 2 + ¯h˜Ω ! z + ¯h∆ 2 −¯h˜Ω !, (15.37) we see that this matrix element has poles at z = −¯h∆ 2 \pm ¯h˜Ω 2 , (15.38) where as before ˜Ω:= p Ω2 + ∆2 (15.39) is the generalized Rabi frequency. Now the corresponding matrix element of the retarded Green function is G+

\[ eg(E) = Geg(E + i0+) = \]

¯hΩ/2

E + ¯h∆ 2 + ¯h˜Ω

\[ 2 + i0+ \]

! E + ¯h∆ 2 −¯h˜Ω

\[ 2 + i0+ \]

!. (15.40) In this expression, the i0+ terms have the effect of shifting both poles infinitesimally below the real axis. Now according to Eq. (15.16), the relevant propagator is given by the inverse Fourier transform

\[ \langle e|U(\tau, 0)|g\rangle = −1 \]

2\pii Z \infty −\infty

\[ dE e−iE\tau/¯h G+ \]

eg(E) (\tau > 0). (15.41) Due to the form of the exponential factor, in the lower complex plane E = x −iy gives a damping factor of the form e−y\tau/¯h. Thus to evaluate this integral, we complete the contour around the lower half-plane. Re[z] Im[z] The contribution from the lower great half-circle vanishes due to the exponential factor and the asymptotic E−2 dependence of the remaining part of the integrand. However, the contour encloses both poles, and the Cauchy integral formula [Eq. (14.79)] states that the integral is 2\pii times the sum of the two residues; however, because the contour runs clockwise, we should also include another minus sign.. We thus find

\[ \langle e|U(\tau, 0)|g\rangle = (¯hΩ/2)ei(∆+˜Ω)\tau/2 \]

−¯h˜Ω

\[ + (¯hΩ/2)ei(∆−˜Ω)\tau/2 \]

¯h˜Ω = −iΩ ˜Ω

\[ ei∆\tau/2 sin \]

˜Ω\tau 2 . (15.42) This gives the probability amplitude of finding the atom in the excited state at time \tau, given that it was in the ground state at time 0. Thus, it agrees with our earlier solution in Section 5.2.2.1.

15.4 Level-Shift Operator 15.4 Level-Shift Operator Again, let’s start with a system with Hamiltonian H = H0 + V , where H0 is the unperturbed Hamiltonian. Working with the eigenstates of H0, we may be interested in particular in the matrix element G\alpha\alpha(z), which reflects the survival probability of state |\alpha\rangle due to the perturbation V . Now suppose we define projection operators for |\alpha\rangle and everything but |\alpha\rangle :2

\[ P\alpha := |\alpha\rangle \langle \alpha|, \]
\[ Q\alpha := 1 −P\alpha = \]

X

\[ j̸=\alpha \]
\[ |j\rangle \langle j|. \]

(15.43) Now recalling that the resolvent is defined by

\[ (z −H)G(z) = (z −H0 −V )G(z) = 1, \]

(15.44) we can insert a (P\alpha + Q\alpha) before the resolvent, and then operate on the right with P\alpha:

\[ (z −H0 −V )(P\alpha + Q\alpha)G(z)P\alpha = P\alpha. \]

(15.45) This then becomes

\[ (z −H0 −V )P\alpha[P\alphaG(z)P\alpha] + (z −H0 −V )Q\alpha[Q\alphaG(z)P\alpha] = P\alpha. \]

(15.46) where we used P 2

\[ \alpha = P\alpha and Q 2 \]

\alpha = Q\alpha. Operating on the left with P\alpha, and using P\alphaQ\alpha = Q\alphaP\alpha = 0 [in

\[ particular, P\alpha(z −H0)Q\alpha = 0], \]
\[ P\alpha(z −H)P\alpha[P\alphaG(z)P\alpha] −P\alphaV Q\alpha[Q\alphaG(z)P\alpha] = P\alpha. \]

(15.47) Operating on the left of Eq. (15.46) instead with Q\alpha,

\[ −Q\alphaV P\alpha[P\alphaG(z)P\alpha] + Q\alpha(z −H)Q\alpha[Q\alphaG(z)P\alpha] = 0. \]

(15.48) Solving for Q\alphaG(z)P\alpha,

\[ Q\alphaG(z)P\alpha = \]

Q\alpha

\[ Q\alpha(z −H)Q\alpha \]
\[ V P\alpha[P\alphaG(z)P\alpha], \]

(15.49) and putting this into Eq. (15.47) to decouple these two equations,

\[ P\alpha(z −H)P\alpha[P\alphaG(z)P\alpha] −P\alphaV Q\alpha \]

Q\alpha

\[ Q\alpha(z −H)Q\alpha \]
\[ V P\alpha[P\alphaG(z)P\alpha] = P\alpha. \]

(15.50) Simplifying and expanding out the Hamiltonian, we find P\alpha  z −H0 −V −V Q\alpha

\[ z −Q\alphaH0Q\alpha −Q\alphaV Q\alpha \]

V 

\[ P\alphaG(z)P\alpha = P\alpha. \]

(15.51) The part of the quantity in braces involving the interaction potential is important, and is called the level- shift operator:

\[ R(z) := V + V \]

Q\alpha

\[ z −Q\alphaH0Q\alpha −Q\alphaV Q\alpha \]

V. (15.52) (level-shift operator) Then Eq. (15.51) becomes

\[ P\alpha [z −H0 −R(z)] P\alphaG(z)P\alpha = P\alpha, \]

(15.53) or

\[ P\alphaG(z)P\alpha = \]

P\alpha

\[ z −P\alphaH0P\alpha −P\alphaR(z)P\alpha \]

. (15.54) (projection of the resolvent) 2We are essentially following the derivation of C. Cohen-Tannoudji et al., op. cit., Section III.B.2, p. 174. Essentially the same derivation is also given by P. Lambropoulos and P. Zoller, ‘‘Autoionizing states in strong laser fields,’’ Physical Review A 24, 379 (1981) (doi: 10.1103/PhysRevA.24.379).

15.4.1 Decomposition of the Level-Shift Operator

Chapter 15. Resolvent Operator In particular, by matching the coefficient of P\alpha, this means that the relevant matrix element of the resolvent is given by

\[ G\alpha\alpha(z) = \]
\[ z −E\alpha −R\alpha\alpha(z). \]

(resolvent matrix element in terms of level-shift operator) (15.55) Note that so far this is an exact expression. Furthermore, this formalism is easy to generalize to the case of survival among a subspace of multiple states: P\alpha is redefined as the sum of projectors over the relevant states, and we still have Q\alpha = 1 −P\alpha, with all the algebra so far carrying through except for the last expression above. 15.4.1 Decomposition of the Level-Shift Operator In analyzing the resolvent G(z), recall that we in general require the resolvent just next to the real axis, G(E \pm i0+). We will similarly require the level-shift operator at the same locations, R(E \pm i0+). Using the relation (Problem 14.6)

\[ x \pm i0+ = P 1 \]
\[ x ∓i\pi\delta(x), \]

(15.56) where the ‘‘P’’ denotes that an integral taken over that term is interpreted as a Cauchy principal value, we can rewrite Eq. (15.52) just off the real axis as

\[ R(E \pm i0+) = V + P V \]

Q\alpha

\[ E −Q\alphaHQ\alpha \]
\[ V ∓i\piV Q\alpha\delta(E −Q\alphaHQ\alpha)Q\alphaV. \]

(15.57) Now defining the Hermitian operators

\[ ∆(E) := P 1 \]

¯hV Q\alpha

\[ E −Q\alphaHQ\alpha \]

V (Hermitian part of level-shift operator) (15.58) and

\[ \Gamma(E) := 2\pi \]
\[ ¯h V Q\alpha\delta(E −Q\alphaHQ\alpha)Q\alphaV, \]

(anti-Hermitian part of level-shift operator) (15.59) we see that these are the Hermitian and anti-Hermitian parts of the level-shift operator, which we can now write as

\[ R(E \pm i0+) = V + ¯h∆(E) ∓i¯h\Gamma(E) \]

. (15.60) (level-shift operator) The operators ∆and \Gamma correspond to dispersive and dissipative effects (i.e., Stark shifts and decay) due to the perturbation V . To see this explicitly, we can write down the retarded Green function as G+

\[ \alpha\alpha(E) = G\alpha\alpha(E + i0+) = \]
\[ E −E\alpha −V\alpha\alpha −¯h∆\alpha\alpha(E) + i[¯h\Gamma\alpha\alpha(E)/2 + 0+]. \]

(15.61) In this form, the resolvent appears to have a pole where E satisfies

\[ E = E\alpha + V\alpha\alpha + ¯h∆\alpha\alpha(E) −i¯h\Gamma\alpha\alpha(E) \]

, (15.62) though of course the nature of the singularities of the resolvent depends on the exact forms of ∆\alpha\alpha(E) or \Gamma\alpha\alpha(E). In the pole approximation, we assume that the interaction V leads only to a very weak perturbation and thus we can replace E by E\alpha on the right-hand side of this expression:

\[ E \approx E\alpha + V\alpha\alpha + ¯h∆\alpha\alpha(E\alpha) −i¯h\Gamma\alpha\alpha(E\alpha) \]

. (15.63)

15.4.2 Perturbation Expansion

15.4 Level-Shift Operator Then the propagator from Eq. (15.16) becomes

\[ U\alpha\alpha(\tau, 0) = −1 \]

2\pii Z \infty −\infty

\[ dE e−iE\tau/¯h G+(E) \]

(\tau > 0)

\[ = e−i(E\alpha+V\alpha\alpha)\tau/¯h e−i∆\alpha\alpha(E\alpha)\tau e−\Gamma\alpha\alpha(E\alpha)\tau/2. \]

(15.64) Thus, we see that the energy of the state |\alpha\rangle has been shifted by the first-order shift V\alpha\alpha as well as the higher-order shift ∆\alpha\alpha(E\alpha). The population of the state |\alpha\rangle also decays at the rate \Gamma\alpha\alpha(E\alpha). 15.4.2 Perturbation Expansion Though exact, the expression we have is not yet so illuminating, particularly in the form of the level-shift operator (15.52). However, we can first note that Q\alpha z −H0

\[ Q\alpha = \]

Q\alpha

\[ z −Q\alphaH0Q\alpha \]

= Q\alpha z −H0 . (15.65) That is, since Q\alpha is a sum of projectors of eigenstates of H0, the operators Q\alpha and (z −H0)−1 commute, and (z −H0)−1 acts as a scalar when multiplied by each of the projectors in Q\alpha. The level-shift operator (15.60)

\[ has the form R = V + V Q\alphaGQ\alphaV ; then using the series (15.28) for the resolvent, we see the expanded \]

level-shift operator has the expansion

\[ R(z) = V + V \]

Q\alpha z −H0 V + V Q\alpha z −H0 V Q\alpha z −H0

\[ V + \cdot \cdot \cdot . \]

(series expansion of level-shift operator) (15.66) Then the matrix element we need for Eq. (15.55) is

\[ R\alpha\alpha(z) = V\alpha\alpha + \]

X

\[ \beta̸=\alpha \]
\[ V\alpha\betaV\beta\alpha \]

z −E\beta + X

\[ \beta̸=\alpha \]
\[ \gamma̸=\alpha \]
\[ V\alpha\betaV\beta\gammaV\gamma\alpha \]
\[ (z −E\beta)(z −E\gamma) + \cdot \cdot \cdot . \]

(series expansion of level-shift matrix element) (15.67) This series can be truncated at any order in V . Note that, at least up to the second-order term, setting z = E\alpha gives the energy shift of the unperturbed level \alpha in time-independent perturbation theory. To gain a bit more insight into this perturbation expansion, note that Eq. (15.55) can be expanded as

\[ G\alpha\alpha(z) = \]

z −E\alpha 

\[ 1 + R\alpha\alpha(z) \]

z −E\alpha + R 2

\[ \alpha\alpha(z) \]
\[ (z −E\alpha)2 + \cdot \cdot \cdot \]

 . (15.68) Thus, even when R\alpha\alpha(z) is truncated to some order in the perturbation, G\alpha\alpha(z) in the form (15.55) still contains contributions at all orders in V , and so corresponds to some nonperturbative expansion of the resolvent. We can compare this to the direct series expansion (15.28) of the resolvent, which gives

\[ G\alpha\alpha(z) = \]

z −E\alpha  1 +

\[ V\alpha\alpha \]

z −E\alpha + X \beta

\[ V\alpha\betaV\beta\alpha \]
\[ (z −E\alpha)(z −E\beta) + \]

X \beta\gamma

\[ V\alpha\betaV\beta\gammaV\gamma\alpha \]
\[ (z −E\alpha)(z −E\beta)(z −E\gamma) + \cdot \cdot \cdot \]

 . (15.69) [Note that the summation indices are not restricted as in Eq. (15.67).] The two series expansions (15.67) and (15.69) have exactly the same terms, but correspond to different orderings of the term. The direct expansion (15.69) is strictly in powers of the perturbation V . However, in the expression (15.67), we note that R\alpha\alpha never contains (z −E\alpha) in any denominator. Thus, the expansion (15.68) for G\alpha\alpha(z) corresponds to a power series (Laurent series) in (z −E\alpha)−1. This focus on the pole at E\alpha is sensible if we are interested in the effects of the perturbation on the state |\alpha\rangle , whose properties are determined by the resolvent in the vicinity of the perturbed energy of state |\alpha\rangle . Given a small perturbation, the pole should not have shifted too much,

15.5.1 Pole Approximation

Chapter 15. Resolvent Operator so we are still maintaining accuracy in the relevant region.3 (For example, the residue of the shifted pole gives the perturbed dynamics, and thus accuracy in this neighborhood is important.) It is also convenient to give the series expansion in terms of the Hermitian and anti-Hermitian parts of R(E \pm i0+). The Hermitian part ∆(E) from Eq. (15.58) expands in the same way as the resolvent in Eq. (15.67):

\[ ∆(E) = V + V \]

Q\alpha E −H0 V + V Q\alpha E −H0 V Q\alpha E −H0

\[ V + \cdot \cdot \cdot . \]

(expansion of level-shift operator, Hermitian part) (15.70) The diagonal matrix elements are then

\[ ¯h∆\alpha\alpha(E) = V\alpha\alpha + \]

X

\[ \beta̸=\alpha \]
\[ V\alpha\betaV\beta\alpha \]

E −E\beta + X

\[ \beta̸=\alpha \]
\[ \gamma̸=\alpha \]
\[ V\alpha\betaV\beta\gammaV\gamma\alpha \]
\[ (E −E\beta)(E −E\gamma) + \cdot \cdot \cdot . \]

(matrix-element expansion of level-shift operator, Hermitian part) (15.71) The same matrix element of the anti-Hermitian part, from Eq. (15.59), is

\[ \Gamma\alpha\alpha(E) = 2\pi \]

¯h X

\[ \beta̸=\alpha \]
\[ \gamma̸=\alpha \]
\[ V\alpha\betaV\gamma\alpha\delta(E −\delta\beta\gammaE\beta −V\beta\gamma), \]

(anti-Hermitian part of level-shift operator) (15.72) Notice the similarity to Fermi’s Golden Rule [Eq. (11.58)], except here we are summing over all possible decay paths. 15.5 Spontaneous Decay 15.5.1 Pole Approximation As in our previous analysis of spontaneous emission in Chapter 11, we take the state |e\rangle to be coupled to |g, 1k,\zeta\rangle , with the relevant matrix element of the interaction Hamiltonian reading [Eq. (11.63)]

\[ \langle e|HAF|g, 1k,\zeta\rangle = \]

r ¯h\omegak

\[ 2ϵ0V (ˆ\epsilonk,\zeta \cdot dge) eik\cdotr. \]

(15.73) Since we wish to examine the survival probability, we consider the diagonal matrix element of the resolvent, which from Eq. (15.55) we can write in terms of the level-shift operator as

\[ Gee(z) = \]

z −Ee −Ree(z). (15.74) Near the real axis, we can use Eq. (15.60) to write

\[ Gee(E \pm i0+) = \]

E −Ee −¯h∆ee(E) \pm i¯h\Gammaee(E)/2, (15.75)

\[ where we have used (HAF)ee = 0. We showed in Section 15.4.1 that in the pole approximation, \Gammaee(Ee) is \]

the decay rate of the excited state, and ∆ee(Ee) is the energy (Stark) shift of the excited state. We now must evaluate the dispersive matrix element from Eqs. (15.70),

\[ ¯h∆ee(E) = \]

X

\[ \beta̸=e \]

(HAF)e\beta(HAF)\betae E −E\beta + X

\[ \beta̸=e \]
\[ \gamma̸=e \]
\[ (HAF)e\beta(HAF)\beta\gamma(HAF)\gammae \]
\[ (E −E\beta)(E −E\gamma) \]
\[ + \cdot \cdot \cdot . \]

(15.76) 3This point is discussed in detail by Claude Cohen-Tannoudji, Jacques Dupont-Roc, and Gilbert Grynberg, op. cit., Section III.B.1, p. 172.

15.5.2 Line Shape of Spontaneous Decay

15.5 Spontaneous Decay

\[ Truncating this expression to lowest order and making the pole approximation (E \approx Ee), we have \]
\[ ¯h∆ee(E) \approx \]

X

\[ \beta̸=e \]

(HAF)e\beta(HAF)\betae Ee −E\beta , (15.77) which is just the Lamb shift of the excited state |e\rangle in second-order perturbation theory (evaluated in Section 13.12). Evidently, solving the full expansion (15.76) self-consistently for E leads to the exact Lamb- shifted energy

\[ ˜Ee = Ee + ¯h∆ee( ˜Ee). \]

(15.78) We also need the absorptive matrix element

\[ \Gammaee(E) = 2\pi \]

¯h X

\[ \beta̸=e \]
\[ \gamma̸=e \]
\[ (HAF)e\beta(HAF)\gammae\delta(E −\delta\beta\gammaE\beta −(HAF)\beta\gamma). \]

(15.79) Note that in the pole approximation, we can drop the off-diagonal terms where \beta̸ = \gamma, since the interaction HAF is assumed to be weak enough to justify perturbation theory, and thus the delta function can never reach resonance without the presence of the E\beta term in the argument. Thus,

\[ \Gammaee(Ee) \approx 2\pi \]

¯h X

\[ \beta̸=e \]
\[ |(HAF)e\beta|2 \delta(Ee −E\beta −(HAF)\beta\beta). \]

(15.80) But there is no first-order shift due to the dipole interaction HAF, so

\[ \Gammaee(Ee) \approx 2\pi \]

¯h X

\[ \beta̸=e \]
\[ |(HAF)e\beta|2 \delta(Ee −E\beta). \]

(15.81) This expression is equivalent to Fermi’s Golden Rule, as in Eq. (11.58), explicitly summed over the continuum of final states |g, 1k,\zeta\rangle , and we know from Section 11.6.1 that this leads to the correct spontaneous decay rate in free space. Again, not making any perturbative expansion would yield the same expression for the decay rate, but with Ee replaced by ˜Ee (i.e., the transition frequency \omega0 in the decay-rate formula (11.68) is the exact value, including the Lamb shift computed to all orders). Henceforth, we will simply absorb the Lamb shift into the bare-state energy E\alpha, since we will assume that when applying the results of any calculation, we will use the observed atomic energies, which already include the correct Lamb shift. 15.5.2 Line Shape of Spontaneous Decay To compute the line shape of spontaneous decay, we will use the matrix element \langle g, 1k,\zeta|G(z)|e\rangle of the

\[ resolvent operator, which will give the rate to create a photon in mode (k, \zeta), of frequency \omegak = ck. Starting \]

with the second identity in Eqs. (15.31), we have

\[ \langle g, 1k,\zeta|G(z)|e\rangle = \]

z −¯h\omegak

\[ \langle g, 1k,\zeta|HAF|e\rangle \langle e|G(z)|e\rangle . \]

(15.82) The excited-state matrix element of the resolvent is given by Eq. (15.61), with the shift and decay rate determined in the last section:

\[ \langle e|G+(E)|e\rangle = \]
\[ E −¯h\omega0 −¯h∆ee + i¯h\Gamma/2 + i0+ . \]

(15.83) In the pole approximation, recall that \Gamma is the usual decay rate of |e\rangle , and ∆ee is the Lamb shift of the excited state, which we will absorb into the excited-state energy ¯h\omega0:

\[ \langle e|G+(E)|e\rangle = \]
\[ E −¯h\omega0 + i¯h\Gamma/2 + i0+ . \]

(15.84)

Chapter 15. Resolvent Operator Thus, from (15.82) and (15.84), we can write down the Green function

\[ \langle g, 1k,\zeta|G+(E)|e\rangle = \]
\[ \langle g, 1k,\zeta|HAF|e\rangle \]
\[ (E −¯h\omegak + i0+)(E −¯h\omega0 + i¯h\Gamma/2 + i0+). \]

(15.85) Transforming to find the propagator,

\[ \langle g, 1k,\zeta|U(\tau, 0)|e\rangle = −1 \]

2\pii Z \infty −\infty

\[ dE e−iE\tau/¯h\langle g, 1k,\zeta|G+(E)|e\rangle \]
\[ = −\langle g, 1k,\zeta|HAF|e\rangle \]

2\pii Z \infty −\infty dE

\[ e−iE\tau/¯h \]
\[ (E −¯h\omegak + i0+)(E −¯h\omega0 + i¯h\Gamma/2 + i0+). \]

(15.86) We carry out this integral as in Section 15.3.1 by completing the contour around the lower half-plane. The result is

\[ \langle g, 1k,\zeta|U(\tau, 0)|e\rangle = \langle g, 1k,\zeta|HAF|e\rangle \]



\[ e−i\omegak\tau \]
\[ ¯h(\omegak −\omega0) + i¯h\Gamma/2 + \]
\[ e−i\omega0\taue−\Gamma\tau/2 \]
\[ −¯h(\omegak −\omega0) −i¯h\Gamma/2 \]

 , (15.87) where the first term is the residue of the pole at E = ¯h\omegak and the second term is the residue of the pole at E = ¯h\omega0 −i¯h\Gamma/2. Simplifying this expression, we find

\[ \langle g, 1k,\zeta|U(\tau, 0)|e\rangle = \]
\[ \langle g, 1k,\zeta|HAF|e\rangle \]
\[ ¯h(\omegak −\omega0 + i\Gamma/2) \]

h

\[ e−i\omegak\tau −e−i\omega0\taue−\Gamma\tau/2i \]

. (15.88) Therefore, the probability to decay into mode (k, \zeta) is the squared modulus of this amplitude:

\[ P(k, \zeta, \tau) = |\langle g, 1k,\zeta|U(\tau, 0)|e\rangle |2 = \]
\[ |\langle g, 1k,\zeta|HAF|e\rangle |2 \]
\[ ¯h2[(\omegak −\omega0)2 + \Gamma2/4] \]

h

\[ 1 + e−\Gamma\tau −2e−\Gamma\tau/2 cos(\omegak −\omega0)\tau \]

i . (15.89) Let’s expand out the matrix element, using Eq. (15.73), which becomes

\[ |\langle g, 1k,\zeta|HAF|e\rangle |2 = ¯h\omegak \]

6ϵ0V |dge|2 (15.90) for an isotropic atom in quantization volume V . In passing to the continuum limit, recall that we make the replacement X k −\rightarrow V (2\pi)3 Z d3k. (15.91) Converting to spherical coordinates, in free space this is isotropic, so we may carry out the angular integral, and change variables \omega = ck, with the result X k −\rightarrow V 2\pi2c3 Z \infty

\[ d\omega \omega2. \]

(15.92) Since the summations are implicit in a later calculation of a probability, for now we can let V −1 −\rightarrow \omega2

\[ 2\pi2c3 d\omega \]

(15.93) so that

\[ |\langle g, 1k,\zeta|HAF|e\rangle |2 = \]

¯h\omega 3

\[ 12\pi2ϵ0c3 |dge|2d\omega = ¯h2\Gamma \]
\[ 4\pi d\omega, \]

(15.94) where

\[ \Gamma := \omega 3 \]

0 |dge|2 3\piϵ0¯hc3 (15.95)

15.5.3 Branches of the Resolvent

15.5 Spontaneous Decay

\[ is the usual free-space decay rate, and we have taken \omega \approx \omega0. The we can rewrite Eq. (15.89) as a continuous \]

probability density for emission after changing variables: at frequency \omega = \omegak and polarization \zeta at time \tau:

\[ P(\omega, \zeta, \tau) d\omega = \]

\Gamma

\[ 4\pi[(\omega −\omega0)2 + \Gamma2/4] \]

h

\[ 1 + e−\Gamma\tau −2e−\Gamma\tau/2 cos(\omega −\omega0)\tau \]

i d\omega. (15.96) If we don’t care about the polarization of the emitted light, we can sum over the two orthogonal polarizations to find

\[ P(\omega, \tau) d\omega = \]

\Gamma

\[ 2\pi[(\omega −\omega0)2 + \Gamma2/4] \]

h

\[ 1 + e−\Gamma\tau −2e−\Gamma\tau/2 cos(\omega −\omega0)\tau \]

i d\omega. (time-dependent spontaneous-emission line shape) (15.97) In the long-time limit, when a photon has certainly been emitted, this becomes

\[ P(\omega) d\omega = \]

\Gamma

\[ 2\pi[(\omega −\omega0)2 + \Gamma2/4] d\omega, \]

(15.98) (spontaneous-emission line shape) which is a properly normalized Lorentzian line shape of full width at half maximum of \Gamma. Alternately, integrating the transition probability density over all frequencies gives

\[ P(\tau) = \]

Z \infty

\[ d\omega P(\omega, \tau) \approx \]

Z \infty −\infty

\[ d\omega P(\omega, \tau) = 1 −e−\Gamma\tau, \]

(15.99) which is the expected exponential behavior. Of course, the above treatment predicts some nontrivial time dependence to the line shape, including an oscillatory component. The line shape is shown here for several interaction times. Gtoáo" Gto=o4 Gto=o3 Gto=o2 Gto=o1

\[ (w-w0)/G \]

-5 pPoo(w)/2 The departures from the long-time Lorentzian only occur if \Gamma\tau is not too large. However, to resolve any such difference would require a measurement time much longer than \Gamma−1, so these differences would be difficult to detect experimentally. Intuitively, the oscillations arise here when the exponential decay is ‘‘interrupted’’ at time \tau. The spectral response should thus be the convolution of the long-time Lorentzian with a sinc function that becomes increasingly narrow in time. Unless, of course, there is a way to ‘‘freeze’’ the interaction after some short time. This is not normally possible in spontaneous emission, but is possible in the spontaneous Raman problem we consider below. 15.5.3 Branches of the Resolvent In the previous section, we found the following situation: an atom in the excited state |e\rangle decays to the ground state |g\rangle an energy ¯h\omega0 lower, but the energy photon emitted has an uncertainty that doesn’t necessarily

15.5.4 Nonexponential Decay

Chapter 15. Resolvent Operator match, as it can be emitted in a range of width ¯h\Gamma around the transition energy ¯h\omega0. Clearly, energy is conserved on average, but we have a time-independent Hamiltonian for the coupled quantum atom–field system, so energy should be conserved in each process individually. So what gives? Well, first of all, we will see below that the coupled state |e\rangle no longer has a well-defined energy. We might expect this: since it decays, it is not an eigenstate of the full Hamiltonian, so it should not have a definite energy. Furthermore, the energy that goes into preparing the atom in |e\rangle must be yet more uncertain, since in view of the unstable nature of the state, the preparation must take place on a time scale much shorter than 1/\Gamma. Thus, the uncertainty in the energy of the emitted photon is easily accounted for by the uncertainty in the energy of setting up the problem. But now let’s explore the idea that the coupled state |e\rangle has no well-defined energy. The matrix element

\[ \langle e|G0(z)|e\rangle = \]

z −¯h\omega0 (15.100) of the unperturbed resolvent has a single pole at the excited-state energy ¯h\omega0. However, in Eq. (15.75) we wrote an expression for the same matrix element of the coupled resolvent

\[ \langle e|G(E \pm i0+)|e\rangle = \]

E −¯h[\omega0 + ∆ee(E)] \pm i¯h\Gammaee(E)/2, (15.101) in the vicinity of the real axis, where we see that the value of G(z) jumps as we cross the real axis. In coupling the atom to the continuum, the pole at the bare atomic energy appears to have changed into a branch cut, reflecting the unstable nature of the state, and the fact that it has no well-defined energy. Note that this was not the case in the Rabi-flopping example in Eq. (15.35)—rather, it is a consequence of coupling the excited state to a continuum. In particular, note from Eq. (15.75) that the retarded Green function G+

\[ ee(E) = Gee(E + i0+) appears to have a pole at E = ¯h[\omega0 + ∆ee(E)] −i¯h\Gammaee(E)/2, below the real \]

axis. However, G+ ee(E) is only defined above the real axis; just below the real axis, the resolvent is given by the different value G−

\[ ee(E) = Gee(E −i0+). So while the pole appears to have changed into a branch cut, \]

we can view it as having ‘‘disappeared’’ behind the branch cut that formed. That is, we may still find the pole if we analytically continue G+

\[ ee(E) = Gee(E + i0+) into the lower half plane. In this case, rather than \]

suffer the discontinuity in the resolvent in crossing the branch cut, one can think of crossing continuously through the real axis, but ending up in a different Riemann sheet, or the second Riemann sheet, since the function value in the lower half plane defined in this way differs from the function value given with the branch cut. Then we have the function G II

\[ ee(z) = \]
\[ z −¯h[\omega0 + ∆ee(z)] + i¯h\Gammaee(z)/2, \]

(15.102) which has the same functional form as G+

\[ ee(E) = Gee(E +i0+), but unlike Gee(z), it is defined in this way in \]

the lower half of the complex plane Im[z] < 0, and thus corresponds to the function extended to the second Riemann sheet. The usual example for extending through branch cuts is log(z), which, when defined as a function, has

\[ a branch cut along the negative real axis. Thus, log(−x\pmi\delta) = log(x)\pmi(\pi+\delta) if x > 0, which is a difference \]

of 2\pii + 2\delta that does not vanish as \delta −\rightarrow 0. This preserves the continuity (analyticity) and single-valued nature of log z, basically by excluding the discontinuity from the function’s domain to | arg z| < \pi (i.e., such that the discontinuity is never ‘‘detected’’ along a continuous path). Of course, adding any integer multiple of

\[ 2\pii to log z is still valid as a logarithm, since when inverted, exp(log z + 2\piin) = z. Each n thus corresponds \]

to an ‘‘extension’’ of the logarithm function to a different Riemann sheet. 15.5.4 Nonexponential Decay Thus, in making the pole approximation in Section 15.5.1 to arrive at the rate \Gamma of exponential decay of the excited state, we were somewhat sloppy, as it turns out. Implicitly, we solved an integral by closing a contour around the lower half-plane as in Section 15.3.1. But to do so, we needed to cross a branch cut, and we need to be more careful about this.

15.5 Spontaneous Decay To locate the branch cut in the resolvent, recall Eq. (15.81) for \Gammaee(E). There we can see that the branch cut exists anywhere that \Gammaee(E) is nonzero, since it represents the discontinuity across the branch cut. This function is in fact nonzero for any energy where there exists a possible decay energy E\beta. But since the decay is to state |g, 1k,\zeta\rangle , where the atomic ground state has zero energy, the possible decay energies are ¯h\omegak, that is to say, any positive real number. So the branch cut extends from the origin, and along the entire positive real axis. When we close the contour, including the branch cut, we can do so as shown here. Re[z] Im[z] Eo=o0 The contour starts into the lower half of the complex plane, but in crossing the branch cut, we enter the second Riemann sheet (where the contour is shown by a dotted line). To compensate for this, we return along the negative real axis, turn around the branch point at the origin, and then continue along the great semicircle. This contour still encloses the pole, which is again in the second Riemann sheet (and thus shown as a grey dot), and the residue here results in the exponential decay that we have already derived. The extra contribution is due to the two portions along the negative imaginary axis, which comes from extending Eq. (15.16) to these portions, along with the Eqs. (15.101) and (15.102) for the appropriate Green-function expressions, with the result

\[ \langle e|U ∥(\tau, 0)|e\rangle = −1 \]

2\pii Z 0 −\infty d(iy) ey\tau/¯h [G II ee(iy) −Gee(iy)] = −1 2\pi Z 0 −\infty

\[ dy ey\tau/¯h \]



\[ iy −¯h\omega0 + i¯h\Gamma/2 − \]
\[ iy −¯h\omega0 −i¯h\Gamma/2 \]

 = 2\pii Z \infty

\[ dy e−y\tau/¯h \]



\[ y −i¯h\omega0 −¯h\Gamma/2 − \]
\[ y −i¯h\omega0 + ¯h\Gamma/2 \]

 , (15.103) where we have absorbed the Lamb shift into \omega0. To compute the integral here, first note that the exponential integral (Problem 13.2) E1(z) is defined by

\[ E1(z) := \]

Z \infty z dy e−y y . (15.104) Letting y −\rightarrow \alphay (with Re[\alpha] > 0 to guarantee convergence),

\[ E1(z) = \]

Z \infty

\[ z/\alpha \]

dy e−\alphay y , (15.105)

\[ and letting y −\rightarrow y + \beta (\beta not on the negative real axis, and \beta̸ = 0, as we will see from the expansion for \]

E1(z) later on),

\[ E1(z) = e−\alpha\beta \]

Z \infty

\[ z/\alpha−\beta \]

dy e−\alphay

\[ y + \beta . \]

(15.106) Now taking z = \alpha\beta, we find the integral formula Z \infty e−\alphay

\[ y + \beta dy = e\alpha\betaE1(\alpha\beta) \]
\[ (Re[\alpha] > 0, \beta /\in R−, \beta̸ = 0). \]

(15.107) Then Eq. (15.103) becomes

\[ \langle e|U ∥(\tau, 0)|e\rangle = e−i\omega0\tau \]

2\pii h

\[ e−\Gamma\tau/2E1(−i\omega0\tau −\Gamma\tau/2) −e\Gamma\tau/2E1(−i\omega0\tau + \Gamma\tau/2) \]

i (15.108)

15.5.5 Frequency-Dependent Decay Rate

Chapter 15. Resolvent Operator Using the asymptotic expansion (Problem 15.3)

\[ E1(z) = e−z \]

z  1 −1 z + 2 z2 −3!

\[ z3 + \cdot \cdot \cdot + \]

n!

\[ (−z)n + \cdot \cdot \cdot \]

 (15.109) to lowest order,

\[ \langle e|U ∥(\tau, 0)|e\rangle = \]

2\pii 

\[ −i\omega0\tau −\Gamma\tau/2 − \]
\[ −i\omega0\tau + \Gamma\tau/2 \]



\[ + O(\tau −2) \]

= i\Gamma

\[ 2\pi(\omega 2 \]
\[ 0 + \Gamma2/4)\tau + O(\tau −2). \]

(15.110) Thus, we see that at long times, this calculation predicts a slow power-law decay of the excited-state pop- ulation as t−2, which dominates the exponential decay at long times. Unfortunately, this calculation is not quite right, as we have ignored the dependence of \Gammaee(E) on E, replacing it by its pole-approximation value of \Gamma. 15.5.5 Frequency-Dependent Decay Rate To derive an improved expression for \Gammaee(E), we return to Eq. (15.81), but without making the pole approx- imation:

\[ \Gammaee(E) \approx 2\pi \]

¯h X

\[ \beta̸=e \]
\[ |(HAF)e\beta|2 \delta(E −E\beta). \]

(15.111) In our treatment of spontaneous emission using Fermi’s Golden rule (Section 11.6), we showed that the sum over the continuum modes and the delta function are replaced by the density of states \rho(E), so that

\[ \Gammaee(E) = 2\pi \]
\[ ¯h |(HAF)e\beta|2 \rho(E). \]

(15.112) The density of states is [Eq. (11.66)]

\[ \rho(E) = E 2V \]

\pi2¯h3c3 , (15.113) and we have the usual matrix element [Eq. (11.63)]

\[ |(HAF)e\beta|2 = |\langle e|HAF|g, 1k,\zeta\rangle |2 = ¯h\omegakd 2 \]

ge 6ϵ0V = Ed 2 ge 6ϵ0V , (15.114)

\[ upon identifying the initial and final-state energies as the same (due to the delta function), so that E = ¯h\omegak. \]

Combining Eqs. (15.111)-(15.114), we find

\[ \Gammaee(E) = \]

E 3d 2 ge

\[ 3\piϵ0¯h4c3 = \Gamma \]

E 3 (¯h\omega0)3 , (15.115) after using the usual expression \Gamma = \omega 3 0 d 2 ge 3\piϵ0¯hc3 . (15.116) This is, in fact, not the answer we want. This result comes from using the electric-dipole Hamiltonian

\[ HAF = −d \cdot E. \]

On the other hand, we can use the alternative Coulomb-gauge interaction Hamiltonian

\[ HAF = (e/me)pe \cdot A from Section 9.3. As we discussed before in Section 9.3.2, the ratio of the dipole to the \]

Coulomb Hamiltonians is \omega/\omega0, or in this case E/¯h\omega0. Since the square of the matrix element enters our calculation, Eq. (15.115) instead becomes

\[ \Gammaee(E) = \]

\Gamma ¯h\omega0 E (15.117) (decay-function matrix element)

15.5 Spontaneous Decay in the Coulomb gauge. Of course, in the pole approximation none of this matters, since we identify \omega \approx \omega0 anyway. Why prefer the Coulomb gauge to the dipole gauge now? Recall that they should give the same result, but provided that all atomic levels are included in the interaction, which we are not doing. But in a practical sense, putting the cubic expression (15.115) into Eq. (15.103) in place of \Gamma means we have to factor a cubic polynomial to find the poles and thus carry out the contour integration, and frankly, who wants to do that? The other reason comes from the origin of \Gamma(E) as the singular part of ∆(E), as in the factorization starting with Eq. (15.56). To see this, we can start with the second-order truncation of Eq. (15.66)

\[ Ree(z) = Vee + \]

X

\[ \beta̸=e \]
\[ Ve\betaV\betae \]

z −E\beta = X

\[ \beta̸=e \]
\[ Ve\betaV\betae \]

z −E\beta , (15.118) where we have dropped the first-order term, which we will comment on in a bit. The usual mode sum becomes X

\[ \beta̸=e \]

|Ve\beta|2 −\rightarrow ¯h2\Gamma 2\pi\omega 3 Z \infty

\[ d\omega \omega3, \]

(15.119) with E\beta −\rightarrow ¯h\omega, so that

\[ Ree(z) = ¯h2\Gamma \]

2\pi\omega 3 Z \infty d\omega \omega3 z −¯h\omega . (15.120) This is the (divergent) expression for the Lamb shift in the rotating-wave approximation (before separating out the decay). Then applying Eq. (15.56),

\[ Ree(E \pm i0+) = ¯h2\Gamma \]

2\pi\omega 3 – Z \infty d\omega \omega3 E −¯h\omega ∓i ¯h2\Gamma 2\omega 3 Z \infty

\[ d\omega \omega3\delta(E −¯h\omega) \]

= ¯h2\Gamma 2\pi\omega 3 – Z \infty d\omega \omega3 E −¯h\omega ∓i ¯h\GammaE3 2(¯h\omega0)3 , (15.121) where we identify the last term as (∓i¯h/2)\Gammaee(E), gives us the result (15.115). But recall that the first (Lamb- shift) term here, which diverges as \omega2, must be renormalized by adding the dipole self-energy contribution [Eq. (13.153), Section 13.12.2.1] HP⊥= 2ϵ0 Z d3r P ⊥2(r), (15.122) which comes in at first order here, and then the result must be mass-renormalized. The result then agrees with the Coulomb-gauge calculation, which is only linearly divergent:

\[ Ree(E \pm i0+) = ¯h2\Gamma \]

2\pi\omega0 – Z \infty d\omega \omega E −¯h\omega ∓i ¯h\GammaE 2¯h\omega0 . (15.123) Since the renormalization is necessary to get a physical answer anyway, here we will prefer the Coulomb gauge (particularly in that the renormalization is not straightforward when z is not fixed to \omega0).

15.5.7 Pole Contribution

Chapter 15. Resolvent Operator 15.5.6 Branch Contribution Now we return to calculating corrections to exponential decay, using our improved expression (15.117) for \Gammaee(E). Retracing the derivation above, starting with Eq. (15.103), we find

\[ \langle e|U ∥(\tau, 0)|e\rangle = −1 \]

2\pii Z 0 −\infty d(iy) ey\tau/¯h [G II ee(iy) −Gee(iy)] = −1 2\pi Z 0 −\infty

\[ dy ey\tau/¯h \]



\[ iy −¯h\omega0 + i¯h\Gamma(iy)/2 − \]

iy −¯h\omega0 −i¯h\Gamma(iy)/2  = 2\pii Z \infty

\[ dy e−y\tau/¯h \]



\[ y −i¯h\omega0 −iy\Gamma/2\omega0 \]

\[ y −i¯h\omega0 + iy\Gamma/2\omega0 \]



\[ = \beta− \]

2\pii Z \infty dy

\[ e−y\tau/¯h \]
\[ y −i¯h\omega0\beta− \]
\[ −\beta+ \]

2\pii Z \infty dy

\[ e−y\tau/¯h \]
\[ y −i¯h\omega0\beta+ \]

= 2\pii h

\[ \beta−ei\beta−\omega0\tauE1(i\beta−\omega0\tau) −\beta+ei\beta+\omega0\tauE1(i\beta+\omega0\tau) \]

i , (15.124) where we have introduced the notation

\[ \beta\pm := \]

1 \pm i \Gamma 2\omega0

\[ = \omega0(\omega0 ∓i\Gamma/2) \]

\omega 2

\[ 0 + \Gamma2/4 \]

, (15.125) where \beta\pm \approx 1 for the typical case of \omega0 ≫\Gamma. The asymptotic expansion to lowest order vanishes,

\[ \langle e|U ∥(\tau, 0)|e\rangle = \]

2\pii  \beta−

\[ i\beta−\omega0\tau − \]

\beta+

\[ i\beta+\omega0\tau \]



\[ + O(\tau −2) = O(\tau −2), \]

(15.126) so keeping the second-order term in the expansion (15.109), we find

\[ \langle e|U ∥(\tau, 0)|e\rangle = \]

2\pii  − \beta−

\[ (i\beta−\omega0\tau)2 + \]

\beta+

\[ (i\beta+\omega0\tau)2 \]



\[ + O(\tau −3) \]

=

\[ 2\pii\omega 2 \]

0 \tau 2  1 \beta− −1 \beta+ 

\[ + O(\tau −3) \]

= − \Gamma 2\pi\omega 3

\[ 0 \tau 2 + O(\tau −3). \]

(15.127) Thus, this dominates the exponential decay at long times, so the probability amplitude decreases asymptot- ically as t−4. That is, the survival probability P(t) eventually becomes4 P(t) ∼  \Gamma 2\pi\omega 3 2 1 t4 , (15.128) but for typical transitions where \omega0 ≫\Gamma, the crossover from exponential to this power-law behavior happens only at extremely long times, where the survival probability is essentially undetectably small. 15.5.7 Pole Contribution Returning to Eq. (15.75) for the retarded Green function, G+

\[ ee(E) = \]

E −Ee + i¯h\Gammaee(E)/2, (15.129) 4J. Mostowski and K. Wódkiewicz, ‘‘On the Decay Law of Unstable States,’’ Bulletin L’Académie Polonaise des Science, Série des Sciences Mathématiques, Astronomiques et Physiques 21, 1027 (1973); P. L. Knight and P. W. Milonni, ‘‘Long- Time Deviations from Exponential Decay in Atomic Spontaneous Emission Theory,’’ Physics Letters 56A 275 (1976) (doi:

\[ 10.1016/0375-9601(76)90306-6). \]

15.5.8 Short Times

15.5 Spontaneous Decay we can use the expression (15.117) for the function \Gammaee(E) to go beyond the pole approximation for the contribution of the pole, which previously just gave exponential decay at rate \Gamma. We now have G+

\[ ee(E) = \]
\[ E −¯h\omega0 + i\GammaE/2\omega0 \]

= \beta+

\[ E −\beta+¯h\omega0 \]

, (15.130) and thus, with Eq. (15.16), the pole contribution to the propagator becomes

\[ \langle e|U •(\tau, 0)|e\rangle = −1 \]

2\pii Z \infty −\infty

\[ dE e−iE\tau/¯h G+ \]

ee(E)

\[ = \beta+e−i\beta+\omega0\tau \]
\[ = \omega0(\omega0 −i\Gamma/2) \]

\omega 2

\[ 0 + \Gamma2/4 \]
\[ e−i˜\omega0\taue−˜\Gamma\tau/2 = (˜\omega0 −i˜\Gamma/2) \]

\omega0

\[ e−i˜\omega0\taue−˜\Gamma\tau/2, \]

(15.131) where the contour completed around the lower half-plane encloses the single pole at \beta+, and we have defined the shifted resonance frequency

\[ ˜\omega0 := \]

\omega 3 \omega 2

\[ 0 + \Gamma2/4 = \]

\omega0

\[ 1 + (\Gamma/2\omega0)2 \]

(15.132) (shifted resonance frequency) and the shifted decay rate ˜\Gamma := \Gamma

\[ 1 + (\Gamma/2\omega0)2 , \]

(15.133) (shifted decay rate) both of which have very small corrections of order (\Gamma/\omega0)2 (typically ∼10−16 for alkali dipole transitions) as a result of the more precise treatment of the pole, accounting for the frequency dependence of the decay rate (15.117). 15.5.8 Short Times

\[ To focus on short times, we will need the series expansion of the exponential integral around z = 0 (see \]

Problem 13.2),

\[ E1(x) = −\gamma −log x − \]

\infty X j=1 (−1)jxj jj! , (15.134)

\[ where \gamma \approx 0.577 215 664 901 532 860 607 is the Euler–Mascheroni constant. From Eq. (15.124), we had the \]

branch contribution

\[ \langle e|U ∥(\tau, 0)|e\rangle = \]

2\pii h

\[ \beta−ei\beta−\omega0\tauE1(i\beta−\omega0\tau) −\beta+ei\beta+\omega0\tauE1(i\beta+\omega0\tau) \]

i (15.135) to the propagator, where

\[ \beta\pm := \]

1 \pm i \Gamma 2\omega0

\[ = \omega0(\omega0 ∓i\Gamma/2) \]

\omega 2

\[ 0 + \Gamma2/4 \]
\[ = ˜\omega0 ∓i˜\Gamma/2 \]

\omega0 . (15.136) Sadly, this expression diverges at short times. Namely, the logarithmic term in Eq. (15.134) leads to a short-time scaling of

\[ \langle e|U ∥(\tau, 0)|e\rangle = −1 \]

2\pii h

\[ \beta−log(i\beta−\omega0\tau) −\beta+ log(i\beta+\omega0\tau) \]

i . (15.137) Since the terms do not exactly cancel, there is a logarithmic divergence at short times. Evidently, there is a problem at short times with taking the integral over all frequencies.

Chapter 15. Resolvent Operator 15.5.8.1 Hard Cutoff To handle this, we can introduce a high-energy cutoff \Lambda for the integrals. This echoes the strategy in the Lamb shift (Section 13.12), where the argument was that the logarithmically divergent integral should be cut off at the large energy \Lambda = mec2, where relativistic effects should take over and naturally cut off the integral.5 To evaluate the integral with the frequency cutoff, Z \Lambda e−\alphay

\[ y + \beta dy = \]

Z \infty e−\alphay

\[ y + \beta dy − \]

Z \infty \Lambda e−\alphay

\[ y + \beta dy \]

= Z \infty e−\alphay

\[ y + \beta dy −e−\alpha\Lambda \]

Z \infty e−\alphay

\[ y + \beta + \Lambda dy \]

(15.138) where we have let y −\rightarrow y + \Lambda in the second step. Then using Eq. (15.138), we have the integral formula Z \Lambda e−\alphay

\[ y + \beta dy = e\alpha\betaE1(\alpha\beta) −e\alpha(\beta+\Lambda)E1[\alpha(\beta + \Lambda)] \]
\[ (Re[\alpha] > 0, \beta /\in R−, \beta + \Lambda /\in R−, \beta̸ = 0, \beta + \Lambda̸ = 0). \]

(15.139) Then retracing the derivation of Eq. (15.124), but cutting off the integral at energy \Lambda,

\[ \langle e|U ∥(\tau, 0)|e\rangle = \beta− \]

2\pii Z \Lambda dy

\[ e−y\tau/¯h \]
\[ y −i¯h\omega0\beta− \]
\[ −\beta+ \]

2\pii Z \Lambda dy

\[ e−y\tau/¯h \]
\[ y −i¯h\omega0\beta+ \]

= 2\pii h

\[ \beta−ei\beta−\omega0\tauE1(i\beta−\omega0\tau) −\beta+ei\beta+\omega0\tauE1(i\beta+\omega0\tau) \]

i − 2\pii h

\[ \beta−e(i\beta−\omega0−\Lambda/¯h)\tauE1[(i\beta−\omega0 −\Lambda/¯h)\tau] −\beta+e(i\beta+\omega0−\Lambda/¯h)\tauE1[(i\beta+\omega0 −\Lambda/¯h)\tau] \]

i . (15.140) at short times, using the expansion (15.134), this becomes

\[ \langle e|U ∥(\tau, 0)|e\rangle = \]

2\pii h

\[ \beta−log (1 −\Lambda/i\beta−¯h\omega0) −\beta+ log (1 −\Lambda/i\beta+¯h\omega0) \]

i

\[ + O(\tau), \]

(15.141) which is finite. Notice that even for large \Lambda, the logarithms are still comparatively of order unity and almost the same because the \beta\pm are both close to unity. The two terms then nearly cancel, with the difference of order

\[ \beta−−\beta+ = \]

\omega0\Gamma \omega 2

\[ 0 + \Gamma2/4 \approx \Gamma \]

\omega0 , (15.142) which is much smaller than unity. Thus, the contribution at \tau = 0 is negligible compared to the pole

\[ contribution. In principle, the decay rate should vanish at \tau = 0 (Section 11.7.1), but this does not appear \]

to be the case here with this cutoff or the soft cutoff below. The problem with this cutoff procedure is that it modifies the long-time scaling. The asymptotic calculation to O(\tau −1) from Eq. (15.126) should now have the extra contribution

\[ \langle e|U ∥\Lambda(\tau, 0)|e\rangle = \]

2\pii  \beta−

\[ i\beta−\omega0\tau −\Lambda/¯h − \]

\beta+

\[ i\beta+\omega0\tau −\Lambda/¯h \]



\[ + O(\tau −2) \]

(15.143) from the new cutoff terms. However, due to the presence of \Lambda, these terms no longer cancel, but

\[ \langle e|U ∥\Lambda(\tau, 0)|e\rangle = \]

\chi

\[ \pi[\chi2 + (1 −i\xi2)](\Lambda/¯h)\tau + O(\tau −2) \]

\approx \chi

\[ \pi(\Lambda/¯h)\tau + O(\tau −2), \]

(15.144) 5For issues regarding this hard cutoff and the dipole approximation in this calculation, see J. Seke and W. N. Herfort, ‘‘Deviations from exponential decay in the case of spontaneous emission from a two-level atom,’’ Physical Review A 38, 833 (1988) (doi: 10.1103/PhysRevA.38.833); J. Seke and W. Herfort, ‘‘Finite-time deviations from exponential decay in the case of spontaneous emission from a two-level hydrogenic atom,’’ Physical Review A 40, 1926 (1989) (doi: 10.1103/PhysRevA.40.1926). Note that they obtain a slightly different, cutoff-dependent asymptotic scaling.

15.5.9 Intermediate Times

15.5 Spontaneous Decay where \chi := \Gamma/2\omega0 ≪1 and \xi := ¯h\omega0/\Lambda ≪1. Thus, the long-time scaling behavior is P(t) ∼ \chi2

\[ \pi2(\Lambda/¯h)2\tau 2 = \]

\Gamma2

\[ 4\pi2\omega 2 \]
\[ 0 (\Lambda/¯h)2\tau 2 , \]

(15.145) which is not the \tau −4 behavior we expect from Eq. (15.128). Note that there was no divergence problem at long times before, so this new scaling is a symptom that indicates this cutoff procedure is not quite right: the long-time scaling behavior should be cutoff-independent. Since it scales as \Lambda−2, we might expect the numerical coefficient to be small, but it will still eventually dominate. 15.5.8.2 Soft Cutoff A different scenario for cutting off the integral is to smoothly bring the integral to zero at large frequencies. Physically, this represents the fact that the effects of short wavelengths should be attenuated by the finite size of the atom, since the atom ‘‘smooths’’ out the wave on this length scale. Thus, it is appropriate to take a cutoff energy of \Lambda ∼2\pic/a, where a is the atomic radius (Bohr radius). We explicitly miss this effect in the dipole approximation, which treats the atom as a point. A simple functional form for the cutoff is an exponential of the form e−y/\Lambda to cut off large energies y. This corresponds to assuming a Lorentzian shape for the atom, with the cutoff modeling the convolution of the atomic profile with the field modes of different frequencies. Thus, Eq. (15.124) becomes

\[ \langle e|U ∥(\tau, 0)|e\rangle = \]

2\pii Z \infty

\[ dy e−y\tau/¯h \]



\[ y −i¯h\omega0 −iy\Gamma/2\omega0 \]

\[ y −i¯h\omega0 + iy\Gamma/2\omega0 \]

 e−y/\Lambda = 2\pii Z \infty

\[ dy e−y(\tau+¯h/\Lambda)/¯h \]



\[ y −i¯h\omega0 −iy\Gamma/2\omega0 \]

\[ y −i¯h\omega0 + iy\Gamma/2\omega0 \]

 = 2\pii n

\[ \beta−ei\beta−\omega0(\tau+¯h/\Lambda)E1[i\beta−\omega0(\tau + ¯h/\Lambda)] −\beta+ei\beta+\omega0(\tau+¯h/\Lambda)E1[i\beta+\omega0(\tau + ¯h/\Lambda)] \]

o , (15.146) which is exactly the same as the result without any cutoff, but with the time displaced forward \tau −\rightarrow \tau +¯h/\Lambda. This avoids the singularity at \tau = 0 since there

\[ \langle e|U ∥(\tau = 0, 0)|e\rangle = −1 \]

2\pii h

\[ \beta−log[i\beta−¯h\omega0/\Lambda] −\beta+ log[i\beta+¯h\omega0/\Lambda] \]

i , (15.147) which is again finite and negligible compared to unity. At long times, \tau + ¯h/\Lambda \approx \tau, so we obtain the correct \tau −4 scaling at long times. 15.5.9 Intermediate Times The total survival probability is then given by combining the pole amplitude from Eq. (15.131) and the branch amplitude from Eq. (15.146):

\[ P(\tau) = |\langle e|U(\tau, 0)|e\rangle |2 = \]
\[ \langle e|U •(\tau, 0)|e\rangle + \langle e|U ∥(\tau, 0)|e\rangle \]

=

\[ \beta+e−i\beta+\omega0\tau + \]

2\pii n

\[ \beta−ei\beta−\omega0(\tau+¯h/\Lambda)E1[i\beta−\omega0(\tau + ¯h/\Lambda)] −\beta+ei\beta+\omega0(\tau+¯h/\Lambda)E1[i\beta+\omega0(\tau + ¯h/\Lambda)] \]

o . (15.148) We have seen that for short times, the pole contribution dominates, and thus the decay is exponential. For long times, the branch contribution dominates, and the decay crosses over to a power law. At intermediate times, when both contributions are important, the behavior is somewhat more complicated. The pole contri- bution always oscillates at optical frequencies, but we have seen that asymptotically, the branch contribution does not. Thus, there can be optical-frequency beating between the two contributions. Unfortunately, it is difficult to visualize these high-frequency beats on the long decay time scales, except for the unrealistic case where \Gamma is not too different from \omega0. This plot shows the case \Gamma/\omega0 = 10−1, with a cutoff \Lambda/¯h\Gamma = 102, along with the exponential pole decay alone and the asymptotic \tau −4 decay from Eq. (15.128). The oscillations at the crossover are clear here.

Chapter 15. Resolvent Operator asymptotic pole combined Gt survival probability Poo(t) 10-20 10-4 10-8 10-12 10-16 This next plot shows the case \Gamma/\omega0 = 10−2, with a cutoff \Lambda/¯h\Gamma = 103, again, along with the exponential pole decay alone and the asymptotic \tau −4 decay from Eq. (15.128). Here, the oscillations are already too rapid to meaningfully plot, so we instead plot the envelope of the oscillations (calculated from the sum and difference of the absolute values of the two contributions). The increases in both parameters shift the crossover to a later time and smaller survival probability. asymptotic pole combined Gt survival probability Po(t) 10-25 10-5 10-10 10-15 10-20 For the realistic case of 87Rb on the 780 nm D2 transition, we have \Gamma/2\pi = 6.07 MHz and \omega0/2\pi = 384.23 THz. We can also take \Lambda = hc/a, where the atomic radius a = 2.99 a0 comes from the dipole moment of 2.99 ea0 for the |F = 2, mF = 2\rangle −\rightarrow |F ′ = 3, m′ F = 3\rangle stretched-state hyperfine transition, and a0 \approx 0.529 Å is the Bohr radius. Thus, the parameters we need are \Gamma/\omega0 = 1.58\times10−8 and \Lambda/¯h\Gamma = 5\times1010. Note that to obtain the correct asymptotic behavior numerically, a cancellation between the different terms is necessary to get a smaller number, so arbitrary-precision arithmetic is required in this regime (standard double-precision, floating-point arithmetic gives an error-dominated asymptotic scaling as \tau −2.) Also, note that had we instead used the hard cutoff, the asymptotic scaling (15.145) makes a substantial difference in the long-time region of this plot, even with the relativistic cutoff.

15.6 Spontaneous Raman Scattering asymptotic pole combined Gt survival probability Poo(t) 10-70 10-10 10-20 10-30 10-40 10-50 10-60 The crossover occurs after some 130 lifetimes, with a survival probability well below 10−50. Long-time nonexponential decay of atomic spontaneous emission is unlikely to ever be seen in an experiment. This is the way the numbers work out in this problem but keep in mind that long-time nonexponential decay is a generic phenomenon, and the method here is a good way to get the full time dependence of the decay. 15.5.10 Interpretation of Nonexponential Decay What is the meaning of this long, nonexponential tail of the decay curve? First of all recall that exponential decay follows from having a constant rate of decay,

\[ \partial tP = −\GammaP(t) \]

−\rightarrow

\[ P(t) = P(0) e−\Gammat. \]

(15.149) This fundamentally means that the system decays in exactly the same way at each instant in time, inde- pendent of its past history. This solution is unique, so any deviation from exponential decay points to a ‘‘memory’’ in the system, or a breakdown of the Markov approximation (the Born–Markov master equation of Section 4.5, or equivalently the Lindblad master equation of Section 19.1 assume the Markov approxima- tion and thus cannot predict this kind of nonexponential decay). The ‘‘memory’’ of the atom of the emitted photon is somewhat counterintuitive, however: evidently the photon-emission amplitude, even though it propagates rapidly away from the atom to infinity, has some long tail that interacts with the atom and interferes with the remaining decay amplitude. 15.6 Spontaneous Raman Scattering As another example of the resolvent method, consider spontaneous Raman scattering in a three-level \Lambda atom, where the transition |g\rangle −\rightarrow |e\rangle is coupled by a laser field with detuning ∆, and spontaneous decay occurs on the |e\rangle −\rightarrow |f\rangle transition at rate \Gamma. D W G |og‚ |oe‚ |ofo‚

Chapter 15. Resolvent Operator Technically, |e\rangle must also decay to |g\rangle if the transition can be coupled by the laser, but we assume that this decay route is much slower than the decay to |f\rangle (see Problem 15.5). This model also describes quenching of a metastable state by coupling to a quickly decaying state, or influence on the metastability of a state by coupling to another decaying level.6 The free Hamiltonian is

\[ H0 = ¯h∆|g\rangle \langle g| −¯h\omegaef|f\rangle \langle f| \]

(15.150) in the rotating frame of the laser field (Section 5.1.5), taking the energy Ee of |e\rangle to be zero, and where

\[ \omegaef = (Ee −Ef)/¯h. The atom–field coupling is given in the rotating-wave approximation by \]

V = ¯hΩ

\[ \sigma + \sigma\dagger \]
  • X k,\zeta ¯h h
\[ gk,\zeta(r)\sigma\dagger \]
\[ f ak,\zeta + H.c. \]

i , (15.151) where \sigma := |g\rangle \langle e|, \sigmaf := |f\rangle \langle e|, Ωis the usual Rabi frequency for the laser field, and gk,\zeta are the free-space coupling coefficients for the vacuum field [Eq. (11.7)]. Now we can focus on the coupling between |g\rangle and |e\rangle . Defining the projector P := |g\rangle \langle g| + |e\rangle \langle e| and

\[ the orthogonal projector Q := 1 −P, we can use the result (15.54) in terms of the level-shift operator, \]
\[ PG(z)P = \]

P z −PH0P −PR(z)P . (15.152) The resolvent in the subspace of |g\rangle and |e\rangle can then be written in matrix form as  Gee(z) Geg(z) Gge(z) Ggg(z)  =  z −Ee −Ree(z) −Reg(z) −Rge(z) z −Eg −Rgg(z) −1 . (15.153) Since we want to analyze the survival probability of |g\rangle , we can use the inversion formula  a b c d −1 = ad −bc  d −b −c a  (15.154) to write

\[ Ggg(z) = \]

z −Ee −Ree(z) [z −Ee −Ree(z)][z −Eg −Rgg(z)] −Rge(z)Reg(z), (15.155) which we will now evaluate. Using Eq. (15.60) for the level-shift operator, we can compute the matrix element

\[ Ree(E + i0+) = Vee + ¯h∆ee(E) −i¯h\Gammaee(E) \]

. (15.156)

\[ In the pole approximation, we take E = Ee, and then ¯h∆ee(Ee) is the Lamb shift—which we absorb into \]
\[ the excited-state energy—of |e\rangle due to the coupling to the vacuum continuum, and \Gamma = \Gammaee(Ee) represents \]

the spontaneous decay of |e\rangle −\rightarrow |f\rangle due to the vacuum coupling, and thus

\[ Ree(E + i0+) = −i¯h\Gamma \]

2 . (15.157) Similarly,

\[ Rgg(E + i0+) = 0, \]

(15.158) since |g\rangle is not coupled to the vacuum continuum. To get the off-diagonal matrix elements, we can use the perturbative expansion (15.66) up to second order,

\[ R(z) = V + V \]

Q z −H0 V, (15.159) 6Claude Cohen-Tannoudji, Jacques Dupont-Roc, and Gilbert Grynberg, Atom–Photon Interactions: Basic Processes and Applications (Wiley, 1992), Section III.C.3.

15.6.1 Weak Pumping

15.6 Spontaneous Raman Scattering so that to second order

\[ Reg(z) = Veg = ¯hΩ \]

2 , (15.160) with the same result for Rge(z), since

\[ \langle e|V QV |g\rangle = \]

X k,\zeta

\[ \langle e|V |f, 1k,\zeta\rangle \langle f, 1k,\zeta|V |g\rangle = 0, \]

(15.161) again since |g\rangle is not coupled to the vacuum. Now that we have the level-shift operator in the subspace of |g\rangle and |e\rangle , we note that a nice interpre- tation of Eq. (15.152) is that PG(z)P is the resolvent operator of the effective subspace Hamiltonian

\[ P[H0 −R(z)]P = \]

 Ee −i¯h\Gamma/2 ¯hΩ/2 ¯hΩ/2 Eg  , (15.162) which is now no longer Hermitian due to the decay. Returning now to the resolvent matrix element (15.155), which now becomes G+

\[ gg(E) = Ggg(E + i0+) = \]
\[ E −Ee + i¯h\Gamma/2 \]
\[ (E −Ee + i¯h\Gamma/2)(E −Eg) −(¯hΩ/2)2 , \]

(15.163) which has poles E\pm = 1  Ee + Eg −i¯h\Gamma \pm s Ee −Eg −i¯h\Gamma 2 + (¯hΩ)2   (15.164) (shifted energies) corresponding to the eigenvalues of the effective Hamiltonian (15.162). Thus, the propagator from the inversion formula (15.16) gives the survival amplitude

\[ \langle g|U(\tau, 0)|g\rangle = −1 \]

2\pii Z \infty −\infty

\[ dE e−iE\tau/¯h G+ \]

gg(E). (15.165) We can do this integral via a contour around the lower half-plane, which encloses both poles, since the square root of Eq. (15.164) always has an imaginary part smaller in magnitude than i¯h\Gamma/2 (this is apparent when visualizing the squaring and square root operations as respectively doubling and halving the complex angle). Then with G+

\[ gg(E) = \]
\[ E −Ee + i¯h\Gamma/2 \]

(E −E+)(E −E−), (15.166) the propagator becomes

\[ \langle g|U(\tau, 0)|g\rangle = \]

E+ −E− 

\[ E+ −Ee + i¯h\Gamma \]



\[ e−iE+\tau/¯h − \]

 E−−Ee + i¯h\Gamma 

\[ e−iE−\tau/¯h \]

 , (survival amplitude) (15.167) which is a fairly complicated expression, which we can analyze more intuitively in the the limits of weak and strong pumping. 15.6.1 Weak Pumping For weak pumping, Ωis small, and thus we can expand the square root in Eq. (15.164) to lowest order in Ω: E\pm \approx 1  Ee + Eg −i¯h\Gamma \pm  Ee −Eg −i¯h\Gamma      1 + (¯hΩ)2  Ee −Eg −i¯h\Gamma 2       , (15.168)

15.6.2 Strong Pumping

Chapter 15. Resolvent Operator or E+ \approx Ee −i¯h\Gamma + (¯hΩ)2  Ee −Eg −i¯h\Gamma  E−\approx Eg − (¯hΩ)2  Ee −Eg −i¯h\Gamma . (15.169) Note that the eigenvalues here are only small corrections to the original eigenvalues. Recalling that Ee = 0 and Eg = ¯h∆, E+ ¯h \approx −i\Gamma 2 − Ω2  ∆+ i\Gamma  = −i\Gamma 2 −Ω2 (∆−i\Gamma/2)

\[ ∆2 + \Gamma2/4 \]

 = −i\Gamma 2 −˜∆+ i˜\Gamma E− ¯h \approx ∆+ Ω2  ∆+ i\Gamma

\[  = ∆+ Ω2 (∆−i\Gamma/2) \]
\[ ∆2 + \Gamma2/4 \]
\[  = ∆+ ˜∆−i˜\Gamma \]

2 , (15.170) where we have defined ˜∆:= " Ω2

\[ ∆2 + \Gamma2/4 \]

 # ∆ (15.171) (shift of |g\rangle ) and ˜\Gamma := " Ω2

\[ ∆2 + \Gamma2/4 \]

 # \Gamma. (15.172) (decay rate of |g\rangle ) Thus, the survival amplitude (15.167) becomes

\[ \langle g|U(\tau, 0)|g\rangle = \]

˜∆−i˜\Gamma !

\[ ei ˜∆\taue−(\Gamma−˜\Gamma)\tau/2 + \]
\[ ∆+ ˜∆+ i(\Gamma −˜\Gamma) \]

!

\[ e−i(∆+ ˜∆)\taue−˜\Gamma\tau/2 \]
\[ (∆+ i\Gamma/2) + 2( ˜∆−i˜\Gamma/2) \]

, (15.173) or noting that ˜\Gamma ≪\Gamma and ˜∆≪∆,

\[ \langle g|U(\tau, 0)|g\rangle \approx e−i(∆+ ˜∆)\taue−˜\Gamma\tau/2, \]

(15.174) (weak-pumping survival amplitude) This expression shows that the survival amplitude rotates at the natural (unperturbed) frequency of ∆, plus an ac Stark shift ˜∆due to the pumping laser. There is also the slow decay of |g\rangle at rate ˜\Gamma. Note that in the full expression (15.173) there is also a fast-decaying term, decaying at rate \Gamma −˜\Gamma, and shifted by −˜∆from zero energy. This is because the weak field mixes the ground and excited states slightly, so the part of |e\rangle mixed into |g\rangle decays essentially at the decay rate for |e\rangle , and has the opposite Stark shift as expected for a two-level system. One curious effect is that ˜\Gamma −\rightarrow 0 as \Gamma −\rightarrow \infty. Since a decay from |e\rangle to |f\rangle is a measurement of whether or not the atom is in |e\rangle (indicated by the detection of an emitted photon), \Gamma is essentially the rate at which the measurement is taking place. If this measurement is strong enough, the atom can never be promoted from |g\rangle to |e\rangle in the first place—an example of the quantum Zeno effect. 15.6.2 Strong Pumping

\[ In the limit of strong pumping (Ω≫\Gamma), the eigenvalues/poles from Eq. (15.164) become \]

E\pm \approx 1  Ee + Eg −i¯h\Gamma \pm ¯h˜Ω  = ¯h  ∆−i\Gamma 2 \pm ˜Ω  , (15.175)

15.6.3 General Case

15.6 Spontaneous Raman Scattering where ˜Ω:= p Ω2 + ∆2 (15.176) is the usual generalized Rabi frequency. Then the survival amplitude (15.167) becomes

\[ \langle g|U(\tau, 0)|g\rangle = 1 \]

˜Ω h ∆+ ˜Ω 

\[ e−i∆\tau/2e−\Gamma\tau/4e−i˜Ω\tau − \]

 ∆−˜Ω 

\[ e−i∆\tau/2e−\Gamma\tau/4ei˜Ω\taui \]

, (15.177) or

\[ \langle g|U(\tau, 0)|g\rangle = e−i∆\tau/2e−\Gamma\tau/4 \]

" cos ˜Ω\tau 2 −i∆ ˜Ω sin ˜Ω\tau # . (strong-pumping survival amplitude) (15.178) These are the usual generalized Rabi oscillations [cf. Eq. (5.59)], noting the sign difference in ∆], but now damped at rate \Gamma/2. Here the field mixes |g\rangle and |e\rangle together in equal parts, so |g\rangle decays at half the decay rate of |e\rangle . 15.6.3 General Case The general case interpolates between simple exponential decay and damped Rabi oscillations in a reasonable obvious way, as shown here for the on-resonance case ∆= 0. Do=o0

\[ W/Go=o0.1 \]
\[ W/Go=o0.2 \]
\[ W/Go=o0.5 \]
\[ W/Go=o1 \]
\[ W/Go=o2 \]
\[ W/Go=o5 \]

Gt survival probability Pgoo(t) Even for relatively weak pumping Ω/\Gamma = 0.1, when the decay is essentially exponential, the one obvious feature is the nonexponential decay at short times, since the whole process must start via a part of a Rabi oscillation from |g\rangle to |e\rangle . Of course, we already know that the decay must be nonexponential at short times in any case (Section 11.7.1).

Chapter 15. Resolvent Operator Do=o0

\[ W/Go=o0.05 \]
\[ W/Go=o0.1 \]
\[ W/Go=o0.2 \]
\[ W/Go=o0.5 \]
\[ W/Go=o1 \]

Gt survival probability Pgoo(t)

\[ For the off-resonance case (with Ω= 1), the Rabi oscillations are incomplete, and become more rapid, since \]

the oscillations occur around the generalized Rabi frequency. Obviously, the decay becomes slower for larger detunings, but also note that the fast oscillations damp out before a smooth decay takes over. Wo=o1

\[ D/Go=o0.1 \]
\[ D/Go=o0.5 \]
\[ D/Go=o1 \]
\[ D/Go=o5 \]
\[ D/Go=o20 \]

Gt survival probability Pgoo(t)

15.7 Exercises 15.7 Exercises Problem 15.1 Show that the resolvent operator

\[ G(z) := \]

z −H (15.179) for the Hamiltonian H is analytic off the real axis, in the sense that every matrix element \langle \psi|G(z)|\psi′\rangle for arbitrary states |\psi\rangle , |\psi′\rangle is an analytic function anywhere away from the real axis. State explicitly your criteria for analyticity. Problem 15.2 The inhomogeneous Helmholtz equation

\[ \nabla 2 + k2 \]
\[ \psi(r) = −f(r), \]

(15.180) where f(r) is an arbitrary source function, has a Green function (resolvent) defined by [see Eq. (14.57), noting that we are ditching the ϵ0 but keeping the minus sign] −\nabla 2 −k2

\[ G(r, r′; k2) = \deltad(r −r′) \]

(15.181) in d spatial dimensions. (a) If we assume the Helmholtz equation to be defined on a compact domain, show that the retarded ‘‘energy-space’’ Green function G(r, r′; k2) can be written in the form7

\[ G+(r, r′; k2) = \]

X n

\[ \psin(r)\psi∗ \]

n(r′) k 2 n −k2 −i0+ , (15.182) where \psin(x) are the eigenfunctions of the homogeneous version of Eq. (15.180) with (discrete) eigen- values k = kn. Be careful with the sign of the imaginary deformation here! (b) Show that in the continuum limit where kn −\rightarrow p, the (retarded) Green function may be written as

\[ G+(r, r′; k2) = \]

(2\pi)d Z

\[ ddp \psip(r)\psi∗ \]

p(r′) p2 −k2 −i0+ (15.183) in d spatial dimensions. Problem 15.3 Derive the asymptotic expansion

\[ E1(z) = e−z \]

z  1 −1 z + 2 z2 −3!

\[ z3 + \cdot \cdot \cdot + \]

n!

\[ (−z)n + \cdot \cdot \cdot \]

 . (15.184) Problem 15.4 Work out a formula for the inverse Laplace transform, using the integral formula for the propagator in terms of the resolvent operator. State any restrictions on the validity of your formula. Problem 15.5 In analyzing the spontaneous Raman problem, we ignored any decay back to the initial (ground) state |g\rangle . Suppose we modify the setup to explicitly include a decay rate of \Gamma′ from |e\rangle −\rightarrow |g\rangle . 7see, e.g., Marco Schäfer, Idrish Huet, and Holger Gies, ‘‘Energy-momentum tensors with worldline numerics,’’ International Journal of Modern Physics Conference Series 14, 511 (2012) (doi: 10.1142/S2010194512007647), arXiv.org preprint (arXiv: quant-ph/0605180v3).

Chapter 15. Resolvent Operator D W G Goo' |og‚ |oe‚ |ofo‚ (a) Why is the resolvent method not a natural approach to handle this new problem? (b) Derive a corrected formula for the decay rate of |g\rangle in the weak pumping limit, accounting for the new decay path. Hint: set up and solve Einstein-type rate equations for this system, generalizing the results from the resolvent approach as appropriate (e.g., introducing an auxiliary decay path). You need not retrace the derivation using the resolvent method if you can just indicate the appropriate changes. Problem 15.6 Consider an atom at a fixed location in an optical cavity. The optical cavity is initially in the vacuum state, and its resonance frequency \omega does not necessarily coincide with the atomic resonance frequency \omega0. The atom starts in the excited state. (a) Compute the decay rate for the atom, assuming the ‘‘bad-cavity’’ limit of large \kappa. Ignore decay into non-cavity modes. Hint: what is the level structure of this problem? (b) The enhancement of the atomic spontaneous emission rate by a cavity is called the Purcell effect. What is now known as the Purcell factor was given by Purcell8 as

\[ \etaP = 3Q\lambda3 \]

4\pi2V , (15.185) where Q is the quality factor of the cavity, \lambda is the emission wavelength, and V is the cavity volume. Purcell’s result was that multiplying the atomic decay rate by this factor gives the cavity-modified decay rate. Show that your result is consistent with Purcell’s for a cavity whose resonance matches that of the atom, under the assumption that the atomic dipole is aligned with the cavity-mode polarization

\[ (ˆ\epsilon \cdot dge = dge, without the factor of \]

\sqrt 3). Problem 15.7 Suppose the intensity of an optical cavity of resonant frequency \omega decays exponentially at rate \kappa. The cavity spectrum is bounded from below, and thus should decay nonexponentially at long times. For example, given that the cavity begins with exactly one photon, the photon’s survival probability should become nonexponential at long times. (a) Treating the cavity decay rate as approximately independent of frequency, give an expression for the asymptotic survival probability for long times. (b) Estimate the scaled time \kappat of crossover to nonexponential decay for a linear, two-mirror cavity of length 10 cm, assuming identical mirrors with 99% intensity reflection coefficients and a resonance wavelength of 532 nm. Also, estimate the survival probability at this crossover time. Problem 15.8 Consider the spontaneous-Raman problem, for which we derived the survival probability of |g\rangle in Section 15.6. 8E. M. Purcell, ‘‘Spontaneous Emission Probabilities at Radio Frequencies,’’ Physical Review 69, 681 (1946) (doi: 10.1103/PhysRev.69.674.2).

15.7 Exercises D W G |og‚ |oe‚ |ofo‚ Under the condition of weak excitation (small Ωor large |∆|), derive an expression for the spectral lineshape of the emitted light (assuming the long-time limit). Interpret your solution.