Spontaneous emission of a sodium Rydberg atom close to an optical nanofibre

Within the last two decades, the strong dipole–dipole interaction experienced by two neighbouring Rydberg-excited atoms [1] has become the main ingredient for many of the atomic quantum information protocol proposals (see [2] and references therein). In particular, this interaction can be so large as to even forbid the simultaneous resonant excitation of two atoms if their separation is less than a specific distance, called the blockade radius [3], which typically depends on the intensity of the laser excitation and the interaction between the Rydberg atoms [4]. The discovery of this 'Rydberg blockade' phenomenon [3, 5–9] paved the way for a new encoding scheme using atomic ensembles as collective quantum registers [5, 10–12] and repeaters [13–15]. In this novel framework, information is stored in collective spin-wave-like symmetric states, which contain fully delocalized atomic excitations. Qubits are more easily manipulated and more robust in this collective approach than in the usual single-particle paradigm.

Scalability is one of the crucial requirements for quantum devices [16] and interfacing atomic ensembles into a quantum network is a possible way to reach this goal. Photons naturally appear as ideal information carriers and the photon-based protocols considered so far include free-space [17], or guided propagation through optical fibres [13]. The former has the advantage of being relatively easy to implement, but presents the drawback of strong losses. The latter requires a cavity quantum electrodynamics setup, which is experimentally more involved. An alternative option would be to resort to optical nanofibres. Such fibres have recently received much attention [18, 19] because the coupling to the evanescent (resp. guided) modes of a nanofibre allows for easy-to-implement atom trapping [20, 21] (resp. detection [22]). This coupling increases in strength as the fibre diameter reduces and the atoms approach the fibre surface. It was also even shown that energy could be exchanged between two distant atoms via the guided modes of the fibre [23]. This strongly suggests that optical nanofibres could play the role of a communication channel between the nodes of an atomic quantum network consisting of Rydberg-excited atomic ensembles.

In this article, we make a first step towards this goal and investigate the emission rate of a highly-excited (Rydberg) sodium atom in the neighbourhood of an optical nanofibre made of silica. In the perspective of building a quantum network, we are particularly interested in quantifying how much spontaneously emitted light can be captured and guided along the fibre. Here, we study the influence of the atom to fibre distance, the radius of the fibre, and the symmetry of the Rydberg state, on the emission rates into the guided and radiative fibre modes. Extending the treatments of [24, 25], we find that up to ≈13%, of the spontaneously emitted light can be captured and guided along both directions of the fibre, which is comparable with the ratio of ≈30% obtained with a cesium atom initially in its lowest excited state $6{P}_{3/2}$ and located on the surface of a 200 nm diameter nanofibre [25]. Although the theoretical framework we use here is the same, numerical calculations are more complex than in [25] due to the larger number of transitions considered. Contrary to [25], we do not take into account the atomic hyperfine structure in the excited state, which is very small for Rydberg states [26].

The article is organized as follows. In section 2 we briefly present the system and introduce the expressions of the spontaneous emission rates. In section 3, we present the results of our numerical calculations and discuss the different behaviours observed when the atom is initially in an s or p Rydberg state. Finally, in section 4, we conclude and give perspectives of our work. Appendices A and B provide details about the guided and radiative electromagnetic modes, appendix C sketches the derivation of the spontaneous emission rates of the atom in the presence of the nanofibre and appendix D displays the atomic data we used in our calculations.

We consider a sodium atom, initially prepared in the highly-excited (Rydberg) level $n\leqslant 10$, in the vicinity of a silica nanofibre, whose radius is denoted by a and whose axis is conventionally taken as the z-axis, see figure 1. Our goal is to investigate how the presence of the fibre modifies the spontaneous emission rate of the atom: in particular, we want to study the influence of the radius of the fibre, the distance of the atom to the fibre as well as the symmetry of the Rydberg state $\left|{{nl}}_{j},{m}_{j}\right\rangle$ considered and the principal quantum number n on the spontaneous emission rate. Note that, though the configuration is the same as in [25], in this work, the atom is (relatively) highly excited and, in contrast to [25], several transition frequencies must therefore be considered which complicates the numerical work. The choice of the sodium atom and the maximal principal quantum number n_max = 10 is motivated by the fact that, for the relevant transitions $10\to n\,=\,3,\,\cdots ,\,9$, the fibre can be approximately considered as a non absorbing medium of respective refractive indices n₁ = (1.467, 1.450, 1.438, 1.399, 1.112, 1.615, 2.021) [27]. Such constraints may, however, be alleviated by resorting to the formalism of macroscopic quantum electrodynamics and the Green's function approach [28]. These techniques allow to take the absorption of the medium into account and therefore to deal with higher Rydberg states. This formalism and its application to the calculation of energy shifts will be investigated in a future work. Moreover, the choice of sodium, rather than rubidium or cesium which are more commonly used in nanofibre experiments, was made to allow us to neglect relativistic effects on the electronic wavefunctions and therefore simplify our treatment. The case of cesium will also be tackled in a future work.

Zoom In Zoom Out Reset image size

As recalled in [25], the free electromagnetic field in the presence of a cylindrical fibre can be decomposed into guided and radiative modes which respectively correspond to energy propagation along the fibre and radially to it (see appendices A and B).

Guided modes are characterized by their frequency ω > 0 and order m, which is a positive integer fixing the periodicity of the field with respect to φ. Due to the continuity conditions at the core-cladding interface of the fibre, the norm of the projection of the wavevector onto the z axis, denoted by β, can only take a discrete set of values which are the solutions of the so-called characteristic equation, equation (A.1) [29, 30]. The corresponding modes have different cutoff frequencies. In particular, if ω is sufficiently low, only the (so-called 'hybrid') mode HE₁₁, corresponding to m = 1, can propagate along the fibre. Since a given mode can propagate either in the positive or negative z-direction, an extra index f = ±1 is introduced, such that β × f is the (algebraic) projection of the wavevector onto the z-axis. To complete the description, one also allows for two different polarization directions labelled by p = ±1. For simplicity, we gather the characteristic numbers $\left(\omega ,m,f,p\right)$ into one symbol μ and replace the discrete/continuous sums ${\sum }_{{mfp}}{\int }_{0}^{\infty }d\omega$ by ${\sum }_{\mu }$. Finally the general form of the quantized guided field component is

$$\begin{eqnarray*}&&{\vec{E}}_{g}\left(\vec{r}\right)={\rm{i}}\displaystyle \sum _{\mu }\,\sqrt{\displaystyle \frac{{\hslash }\omega \beta ^{\prime} }{4\pi {\epsilon }_{0}}}{a}_{\mu }{\vec{e}}^{\left(\mu \right)}\left(r,\varphi \right){e}^{{\rm{i}}\left(f\beta z+p\varphi \right)}+{\rm{h.c}}.\end{eqnarray*}$$

In this expression, β' stands for the derivative $\left(\tfrac{d\beta }{d\omega }\right)$, ${\vec{e}}_{\mu }$ is the electric-field profile function of the mode $\left(\mu \right)$ whose expression is given in appendix A, while a_μ is the annihilation operator of the mode, satisfying the bosonic commutation rules $\left[{a}_{\mu },{a}_{\mu ^{\prime} }^{\dagger }\right]=\delta \left(\omega -\omega ^{\prime} \right){\delta }_{{mm}^{\prime} }{\delta }_{{ff}^{\prime} }{\delta }_{{pp}^{\prime} }$.

Radiative modes are characterized by their frequency ω > 0, their (positive integer) order m and the projection of the wavevector on the nanofibre axis β which can now vary continuously between −ω n₂/c and ω n₂/c. Here, the negative or positive sign indicates the direction of the propagation of the radiation mode along the z-axis. A last number is needed to fully determine a radiative mode, i.e. the polarization number p = ±1. The two values of p correspond to two modes of orthogonal polarizations, see appendix B. For simplicity, we gather the characteristic numbers $\left(\omega ,\beta ,m,p\right)$ into one symbol ν and replace the discrete/continuous sums ${\sum }_{{mp}}{\int }_{0}^{\infty }d\omega {\int }_{-{{kn}}_{2}}^{{{kn}}_{2}}d\beta$ by ${\sum }_{\nu }$. The general form of the quantized radiative field component is

$$\begin{eqnarray*}&&{\vec{E}}_{r}\left(\vec{r}\right)={\rm{i}}\displaystyle \sum _{\nu }\,\sqrt{\displaystyle \frac{{\hslash }\omega }{4\pi {\epsilon }_{0}}}{a}_{\nu }{\vec{e}}^{\left(\nu \right)}\left(r,\varphi \right){e}^{{\rm{i}}\left(\beta z+m\varphi \right)}+{\rm{h.c}}.\end{eqnarray*}$$

In this expression ${\vec{e}}_{\nu }$ is the electric-field profile function of the mode $\left(\nu \right)$ whose expression is given in appendix B, while a_ν is the annihilation operator of the mode, satisfying the bosonic commutation rules $\left[{a}_{\nu },{a}_{\nu ^{\prime} }^{\dagger }\right]=\delta \left(\omega -\omega ^{\prime} \right){\delta }_{{mm}^{\prime} }{\delta }_{{pp}^{\prime} }$.

In the presence of the nanofibre, the spontaneous emission rate Γ_M of an atom from a state $\left|M\right\rangle$ is the sum of the rates γ_MN from $\left|M\right\rangle$ to all lower states $\left|N\right\rangle$, i.e. ${{\rm{\Gamma }}}_{M}\equiv {\sum }_{N\lt M}{\gamma }_{{MN}}$ with

$$\begin{eqnarray}&&{\gamma }_{{MN}}\equiv 2\pi \displaystyle \sum _{\lambda }{\left|{G}_{\lambda {MN}}\right|}^{2}\delta \left({\omega }_{\lambda }-{\omega }_{{MN}}\right).\end{eqnarray} \tag{ 1 }$$

In the expression above, the sum is performed over all electromagnetic modes denoted by λ, whether they be guided $\left(\lambda =\mu \right)$ or radiative $\left(\lambda =\nu \right)$; we moreover introduce the quantities

$$\begin{eqnarray*}\begin{array}{rcl}{G}_{\mu {MN}} & \equiv & -\sqrt{\displaystyle \frac{\omega \beta ^{\prime} }{4\pi {\epsilon }_{0}{\hslash }}}\left({\vec{d}}_{{MN}}\cdot {\vec{e}}^{\left(\mu \right)}\right){e}^{{\rm{i}}\left(f\beta z+p\varphi \right)},\\ {G}_{\nu {MN}} & \equiv & -\sqrt{\displaystyle \frac{\omega }{4\pi {\epsilon }_{0}{\hslash }}}\left({\vec{d}}_{{MN}}\cdot {\vec{e}}^{\left(\nu \right)}\right){e}^{{\rm{i}}\left(\beta z+m\varphi \right)},\end{array}\end{eqnarray*}$$

characterizing the coupling of the different electromagnetic modes to the atomic transition $\left|M\right\rangle \to \left|N\right\rangle$ of frequency ${\omega }_{{MN}}\equiv \left({E}_{M}-{E}_{N}\right)/{\hslash }$ and dipole matrix element ${\vec{d}}_{{MN}}$. Finally the decoherence rate between states $\left|M\right\rangle$ and $\left|N\right\rangle$ is given by

$$\begin{eqnarray}&&{{\rm{\Gamma }}}_{{MN}}\equiv \displaystyle \frac{1}{2}\left({{\rm{\Gamma }}}_{M}+{{\rm{\Gamma }}}_{N}\right).\end{eqnarray} \tag{ 2 }$$

For a detailed derivation of equations (1), (2), see appendix C.

In this section, we present the numerical results we obtained for the spontaneous emission rate of a sodium atom $\left(Z=11\right)$ initially prepared either in $\left|{{ns}}_{1/2},{m}_{j}\right\rangle$ or $\left|{{np}}_{j},{m}_{j}\right\rangle$ states with $n\leqslant 10$ and j = 1/2 or 3/2. We study the influence of the principal quantum number, n, and the distance from the atom to the fibre surface on the emission rate. We also show how the fibre's radius modifies the relative weights of the different transitions' contributions to the total rate. For simplicity, we consider the contributions of the guided and radiative modes separately. The atomic data we used can be found in appendix D.

3.1. Guided modes

Figure 2 displays the spontaneous emission rates, ${{\rm{\Gamma }}}_{g}^{10s}$ and ${{\rm{\Gamma }}}_{g}^{10p}$, of an atom initially prepared in the states $\left|10{s}_{1/2},{m}_{j}\right\rangle$ and $\left|10{p}_{j},{m}_{j}\right\rangle$ with j = 1/2 or 3/2, respectively, into the guided modes as a function of the distance r to the fibre axis, see figure 1. Note that the rates are presented relative to the spontaneous emission rates ${{\rm{\Gamma }}}_{0}^{10s}$, ${{\rm{\Gamma }}}_{0}^{10p}$ in vacuum and r is expressed in units of the fibre radius with a = 100 nm. As expected, in both cases, the influence of the guided modes vanishes as r increases, and therefore ${{\rm{\Gamma }}}_{g}^{10s},{{\rm{\Gamma }}}_{g}^{10p}\to 0$ when $r\to +\infty$. The maximal value is obtained for r = a, i.e. when the atom is on the fibre surface. More precisely, we have ${{\rm{\Gamma }}}_{g}^{10s}\approx 0.18{{\rm{\Gamma }}}_{0}^{10s}$ for an atom initially prepared in $| 10{s}_{1/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$ and ${{\rm{\Gamma }}}_{g}^{10p}\approx (0.027,0.035,0.044)$ $\times {{\rm{\Gamma }}}_{0}^{10p}$ for an atom initially prepared in $(| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{1/2}$, ${m}_{j}\,=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{3}{2}\rangle )$. In these calculations we assumed that the electronic wave-function of the Rydberg atom is not affected by the nanofibre, which deserves further study. As a more realistic configuration, we shall consider that the Rydberg atom is located at a distance from the fibre surface which is much larger than its radius ${r}_{{Na}}\approx 5\ \mathrm{nm}=\tfrac{a}{20}$. For $r\,=a+10{r}_{{Na}}\approx 1.5\,\times \,a$, we obtain the spontaneous rate ${{\rm{\Gamma }}}_{g}^{10s}\approx 0.066\times {{\rm{\Gamma }}}_{0}^{10s}$ for an atom initially prepared in $| 10{s}_{1/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$ and ${{\rm{\Gamma }}}_{g}^{10p}\approx \left(0.006,0.011,0.015\right)\,\times {{\rm{\Gamma }}}_{0}^{10p}$ for an atom initially prepared in $(| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{1/2}$, ${m}_{j}\,=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{3}{2}\rangle )$. Moreover, we note that in general, ${{\rm{\Gamma }}}_{g}^{10p}\ll {{\rm{\Gamma }}}_{g}^{10s}$, and ${{\rm{\Gamma }}}_{g}^{10{p}_{3/2},{m}_{j}=\pm \tfrac{1}{2}}\lt {{\rm{\Gamma }}}_{g}^{10{p}_{1/2}}\,\lt {{\rm{\Gamma }}}_{g}^{10{p}_{3/2},{m}_{j}=\pm \tfrac{3}{2}}$. The latter relation can be qualitatively understood by geometric arguments on the coupling of guided modes with the atomic orbitals. The more a state is polarized along z, the less it couples to the guided modes which are essentially polarized orthogonally to the fibre axis. This is consistent with what we observe, since the states $\left|10{p}_{3/2},{m}_{j}=\pm \tfrac{1}{2}\right\rangle$ are better aligned along z than the states $\left|10{p}_{1/2},{m}_{j}=\pm \tfrac{1}{2}\right\rangle$ which themselves are more aligned along z than $\left|10{p}_{3/2},{m}_{j}=\pm \tfrac{3}{2}\right\rangle$. This can be seen on their relation with the decoupled basis states⁷ .

Zoom In Zoom Out Reset image size

Figure 3 shows the influence of the principal quantum number n on the spontaneous emission rate ${{\rm{\Gamma }}}_{g}^{{ns}}$ into the guided modes for an atom initially prepared in the state $\left|{{ns}}_{1/2},{m}_{j}=\pm \tfrac{1}{2}\right\rangle$ for n = 5–10. The higher the value of n, the more ${{\rm{\Gamma }}}_{g}^{{ns}}$ is peaked as a function of r/a around 1. Moreover, the plots get closer and closer as n increases: the curves n = 9, 10 cannot be distinguished and for clarity, the curve n = 8 has not been plotted.

Zoom In Zoom Out Reset image size

Finally, figures 4 and 5 illustrate the influence of the fibre radius, a, on the spontaneous emission rate from the state $\left|{{ns}}_{1/2},{m}_{j}=\pm \tfrac{1}{2}\right\rangle$ into the guided modes. More precisely, figure 4 displays the partial spontaneous emission rates along the specific transition $10s\to 3p$ (Note that n = 3 corresponds to the ground state of the sodium atom) into different guided modes HE_mn, EH_mn, TE_mn and TM_mn. Two cases are considered: (i) the atom is located on the fibre surface, i.e. at a distance r = a from the z-axis, and (ii) the atom is placed at a fixed distance of 150 nm from the fibre surface, i.e. at a distance r = a + 150 nm from the z-axis. As expected, case (ii) gives rise to much weaker relative rates than case (i), since the atom is further away from the fibre and therefore the guided modes are strongly attenuated. Moreover, as a increases, the cutoff frequencies of higher modes become smaller: when the cutoff frequency of one mode passes below the frequency of the transition $10s\to 3p$, this mode starts to contribute to the spontaneous emission rate. The peaked structure observed on the different plots results from the peaked shape of the mode intensity profile itself with respect to a.

Zoom In Zoom Out Reset image size

Figure 5 displays the partial spontaneous emission rates ${\gamma }_{g}^{10s\to {np}}$ into the guided modes along the respective transitions $10s\to {np}$ as well as the total rate ${{\rm{\Gamma }}}_{g}^{10s}={\sum }_{3\leqslant n\leqslant 10}{\gamma }_{g}^{10s\to {np}}$ as functions of the fibre radius, a, in the same two cases (i, ii) as above. One observes that, due to the range chosen for a, only the transitions $10s\to {np}$ for $n\,=\,3,4,5$ give relevant contributions to the total rate. It also appears that only the transition $10s\to 3p$ substantially couples to higher-order guided modes, while the other transitions couple only to the fundamental guided mode HE₁₁. On the range chosen for a, the peak structure observed for the total emission rate is therefore mainly due to the partial rate ${\gamma }_{g}^{10s\to 3p}$, while the other transitions smoothly modify the value of ${{\rm{\Gamma }}}_{g}^{10s}$. Note that the intensity profiles of the guided modes relative to the different transition frequencies are expected to coincide up to a rescaling of the a-axis: this scaling factor is given by the ratio of the frequencies. The positions of the peaks of the different partial rates ${\gamma }_{g}^{10s\to {np}}$ should therefore also coincide up to a simple scaling. The heights of the peaks, however, are expected to be different since, for instance, the dipole matrix element is not the same for the different transitions.

Zoom In Zoom Out Reset image size

3.2. Radiative modes

We now turn to the contribution of the radiative modes to the total spontaneous emission rates. Figure 6 displays the spontaneous emission rates ${{\rm{\Gamma }}}_{r}^{10s}$ and ${{\rm{\Gamma }}}_{r}^{10p}$ of an atom initially prepared in the states $\left|10{s}_{1/2},{m}_{j}\right\rangle$, $\left|10{p}_{1/\mathrm{2,3}/2},{m}_{j}\right\rangle$, respectively into the radiative modes as a function of the distance r to the fibre axis, see figure 1. Note that the rates are renormalized by the spontaneous emission rates in vacuum ${{\rm{\Gamma }}}_{0}^{10s}$, resp. ${{\rm{\Gamma }}}_{0}^{10p}$, and r is expressed in units of the fibre radius a = 100 nm. As expected, in both cases, the influence of the fibre vanishes as r increases, i.e. ${{\rm{\Gamma }}}_{r}^{10s},{{\rm{\Gamma }}}_{r}^{10p}\to 1$ for $r\to +\infty$. The maximal value is observed for r = a, i.e. when the atom is on the fibre surface. More precisely, we have ${{\rm{\Gamma }}}_{r}^{10s}\,\approx 1.24\,\times {{\rm{\Gamma }}}_{0}^{10s}$ for an atom initially prepared in $\left|10{s}_{1/2},{m}_{j}=\pm \tfrac{1}{2}\right\rangle$ and ${{\rm{\Gamma }}}_{r}^{10p}\approx \left(1.19,1.23,1.29\right)\times {{\rm{\Gamma }}}_{0}^{10p}$ for an atom initially prepared in $(| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{1/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{3}{2}\rangle )$. For an atom at r ≈ 1.5 × a, i.e. at a distance from the fibre surface, we obtain the spontaneous rate ${{\rm{\Gamma }}}_{r}^{10s}\approx 1.041\times {{\rm{\Gamma }}}_{0}^{10s}$ for an atom initially prepared in $| 10{s}_{1/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$ and ${{\rm{\Gamma }}}_{r}^{10p}\approx \left(1.028,1.044,1.062\right)\,\times {{\rm{\Gamma }}}_{0}^{10p}$ for an atom initially prepared in $(| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{1/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$, $| 10{p}_{3/2}$, ${m}_{j}=\pm \tfrac{3}{2}\rangle )$. This allows us to compute the proportion of light which is emitted into the guided and radiative modes. For instance, for an atom initially prepared in the state $| 10{s}_{1/2}$, ${m}_{j}=\pm \tfrac{1}{2}\rangle$, ${{\rm{\Gamma }}}_{g}/\left({{\rm{\Gamma }}}_{g}+{{\rm{\Gamma }}}_{r}\right)\approx 13 \%$ when the atom is located on the fibre surface $\left(r=a\right)$, and ${{\rm{\Gamma }}}_{g}/\left({{\rm{\Gamma }}}_{g}+{{\rm{\Gamma }}}_{r}\right)\,\approx 6 \%$ when the atom is located at 50 nm from the fibre surface $\left(r=1.5\times a\right)$. Since light is mostly spontaneously emitted into the radiative modes, it seems quite challenging to efficiently interface a Rydberg atom with a guided mode of the nanofibre and, thence, to build a valuable quantum network. The use of atomic ensembles might alleviate this concern, since, as already demonstrated in free-space, their spontaneous emission could be made highly directional and their coupling strength is enhanced [17]. These issues and the perspectives they offer will be addressed in a future work.

Zoom In Zoom Out Reset image size

Finally, in figure 6, one observes a damped semi-oscillatory behaviour for ${{\rm{\Gamma }}}_{r}^{10s}$ and ${{\rm{\Gamma }}}_{r}^{10p}$ as functions of r, and for ${{\rm{\Gamma }}}_{r}^{10p}$ the oscillations of the different contributions j = 1/2, 3/2 are not in phase. These features result from the behaviour of the different transition components ${\gamma }_{{nl}\to n^{\prime} l^{\prime} }$ shown in figure 7 for ${nl}={{np}}_{3/2},{m}_{j}=\pm \tfrac{1}{2}$, which is itself due to the oscillatory behaviour of the radiative field. For a transition of frequency ω, the frequency of oscillation with r is approximately given by 2ω/c.

Zoom In Zoom Out Reset image size

We have investigated the influence of an optical nanofibre on the spontaneous emission rate of a sodium atom prepared in a Rydberg state. The respective contributions of the guided and radiative modes to the total rate were numerically determined, for different principal quantum numbers and different symmetries, and their remarkable features were physically discussed.

Though the radiative modes' contribution is dominant, a small fraction of the spontaneously emitted light is transferred into the guided mode of the nanofibre. This effect might be enhanced by resorting to atomic ensembles which could offer stronger and more directional collective coupling. Using thicker fibres, with more than one guided mode, may also yield for a higher ratio of spontaneous emission into the guided modes. This potentially paves the way towards the implementation of a quantum network based on Rydberg atomic ensembles linked by nanofibres, which will be further addressed in a future work.

This research was supported by the Centre National de la Recherche Scientifique (CNRS) via the grant 'PICS QuaNet'. SNC acknowledges support from OIST Graduate University. The authors thank Antoine Browaeys, Tridib Ray and Fam Le Kien for fruitful discussions.

A guided mode is characterized by a set $\mu \,\equiv \left(\omega ,\beta ,m,f=\pm ,\,p=\pm \right)$. β is the projection of the wavevector onto the axis of the nanofibre whose value is determined by the eigenvalue equation

$$\begin{eqnarray}&&\begin{array}{r}\begin{array}{l}\left(\displaystyle \frac{{n}_{1}^{2}{J}_{m}^{{\rm{{\prime} }}}(\kappa a)}{a\kappa {J}_{m}(\kappa a)}+\displaystyle \frac{{n}_{2}^{2}{K}_{m}^{{\rm{{\prime} }}}(\gamma a)}{a\gamma {K}_{m}(\kappa a)}\right)\left(\displaystyle \frac{{J}_{m}^{{\rm{{\prime} }}}(\kappa a)}{a\kappa {J}_{m}(\kappa a)}+\displaystyle \frac{{K}_{m}^{{\rm{{\prime} }}}(\gamma a)}{a\gamma {K}_{m}(\gamma a)}\right)\\ \,=\,{\left(\displaystyle \frac{mc\beta }{\omega }\right)}^{2}{\left(\displaystyle \frac{1}{{\left(\gamma a\right)}^{2}}+\displaystyle \frac{1}{{\left(\kappa a\right)}^{2}}\right)}^{2}.\end{array}\end{array}\end{eqnarray} \tag{ A.1 }$$

Here we introduced $\kappa \equiv \sqrt{{n}_{1}^{2}{k}^{2}-{\beta }^{2}}$, $\gamma \equiv \sqrt{{\beta }^{2}-{n}_{2}^{2}{k}^{2}}$ and $k\equiv \tfrac{\omega }{c}$. a is the radius of the fibre, n₁ is the core index, n₂ ≈ 1 is the index of the surrounding vacuum. J_m and K_m denote the Bessel functions of the first kind and the modified Bessel functions of the second kind, respectively. Note that, when the monomode conditions are fulfilled, only the hybrid modes HE₁₁ with m = 1 exist, and are fully characterized by $\mu \,\equiv \left(\omega ,f,p\right)$.

The polarization vectors of the guided mode $\left(\mu \right)$ for r < a are given by

$$\begin{eqnarray*}\begin{array}{rcl}{e}_{r}^{\left(\mu \right)} & = & \displaystyle \frac{\beta C}{2{\rm{i}}\kappa }\displaystyle \frac{{K}_{m}(\gamma a)}{{J}_{m}(\kappa a)}\left({J}_{m-1}(\kappa r)(1-{ms})-{J}_{m+1}(\kappa r)(1+{ms})\right)\\ {e}_{\varphi }^{\left(\mu \right)} & = & \displaystyle \frac{p\beta C}{2\kappa }\displaystyle \frac{{K}_{m}(\gamma a)}{{J}_{m}(\kappa a)}\left({J}_{m-1}(\kappa r)(1-{ms})+{J}_{m+1}(\kappa r)(1+{ms})\right),\\ {e}_{z}^{\left(\mu \right)} & = & C\displaystyle \frac{{K}_{m}(\gamma a)}{{J}_{m}(\kappa a)}{J}_{m}(\kappa r),\end{array}\end{eqnarray*}$$

while, for r > a, they are

$$\begin{eqnarray*}\begin{array}{rcl}{e}_{r}^{\left(\mu \right)} & = & \displaystyle \frac{\beta C}{2{\rm{i}}\gamma }\left({K}_{m-1}(\gamma r)(1-{ms})+{K}_{m+1}(\gamma r)(1+{ms})\right),\\ {e}_{\varphi }^{\left(\mu \right)} & = & \displaystyle \frac{\beta {pC}}{2\gamma }\left({K}_{m-1}(\gamma r)(1-{ms})-{K}_{m+1}(\gamma r)(1+{ms})\right),\\ {e}_{z}^{\left(\mu \right)} & = & {{CK}}_{m}(\gamma r),\end{array}\end{eqnarray*}$$

where

$$\begin{eqnarray*}&&s=\displaystyle \frac{\tfrac{1}{{\gamma }^{2}{a}^{2}}+\tfrac{1}{{\kappa }^{2}{a}^{2}}}{\tfrac{{J}_{m}^{\prime} (\kappa a)}{\kappa {{aJ}}_{m}(\kappa a)}+\tfrac{{K}_{m}^{\prime} (\gamma a)}{\gamma {{aK}}_{m}(\gamma a)}}.\end{eqnarray*}$$

Using the normalization condition

$$\begin{eqnarray*}&&{\int }_{0}^{2\pi }{d}\varphi {\int }_{0}^{+\infty }n{\left(r\right)}^{2}{\left|{{\bf{e}}}^{\left(\mu \right)}\right|}^{2}\,r{d}r=1,\end{eqnarray*}$$

we deduce that

$$\begin{eqnarray*}\begin{array}{rcl}| C| & = & \displaystyle \frac{2}{a\beta {K}_{m}\left(\gamma a\right)\sqrt{2\pi ({n}_{1}^{2}{A}_{1}+{n}_{2}^{2}{A}_{2})}}\end{array}\end{eqnarray*}$$

with the abbreviations

$$\begin{eqnarray*}\begin{array}{rcl}{A}_{1} & = & {\left(\displaystyle \frac{1}{\kappa {J}_{m}\left(\kappa a\right)}\right)}^{2}\times \left({\left(1-{ms}\right)}^{2}\left({J}_{m-1}^{2}\left(\gamma a\right)-{J}_{m}\left(\gamma a\right){J}_{m-2}\left(\gamma a\right)\right)\right.\\ & & +{(1+{ms})}^{2}\left({J}_{m+1}^{2}\left(\gamma a\right)-{J}_{m}\left(\gamma a\right){J}_{m+2}\left(\gamma a\right)\right)\\ & & \left.+2\displaystyle \frac{{\kappa }^{2}}{{\beta }^{2}}\left({J}_{m}^{2}\left(\gamma a\right)-{J}_{m-1}\left(\gamma a\right){J}_{m+1}\left(\gamma a\right)\right)\right),\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl}{A}_{2} & = & {\left(\displaystyle \frac{1}{\gamma {K}_{m}\left(\gamma a\right)}\right)}^{2}\times \left({\left(1-{ms}\right)}^{2}\left(-{K}_{m-1}^{2}\left(\gamma a\right)+{K}_{m}\left(\gamma a\right){K}_{m-2}\left(\gamma a\right)\right)\right.\\ & & +{\left(1+{ms}\right)}^{2}\left(-{K}_{m+1}^{2}\left(\gamma a\right)+{K}_{m}\left(\gamma a\right){K}_{m+2}\left(\gamma a\right)\right)\\ & & \left.+2{\left(\displaystyle \frac{\gamma }{\beta }\right)}^{2}\left(-{K}_{m}^{2}\left(\gamma a\right)+{K}_{m-1}\left(\gamma a\right){K}_{m+1}\left(\gamma a\right)\right)\right).\end{array}\end{eqnarray*}$$

A radiative mode is characterized by a set $\nu \,\equiv \left(\omega ,\beta ,m,p=\pm \right)$, where m is the order of the mode, and the meaning of p will be explained below.

Defining the quantities $\kappa \equiv \sqrt{{n}_{1}^{2}{k}^{2}-{\beta }^{2}}$, $\sigma \,\equiv \sqrt{{n}_{2}^{2}{k}^{2}-{\beta }^{2}}$ and k ≡ ω/c, one can write the polarization vectors of the radiative mode $\left(\nu \right)$ for r < a:

$$\begin{eqnarray*}\begin{array}{rcl}{e}_{r}^{\left(\nu \right)} & = & \displaystyle \frac{1}{{\rm{i}}\kappa }\left(\beta {{AJ}}_{m}^{\prime} (\kappa r)+{\rm{i}}B\displaystyle \frac{\omega m}{r\kappa }{J}_{m}(\kappa r)\right),\\ {e}_{\varphi }^{\left(\nu \right)} & = & \displaystyle \frac{1}{{\rm{i}}\kappa }\left({\rm{i}}A\displaystyle \frac{m\beta }{\kappa r}{J}_{m}(\kappa r)-\omega {{BJ}}_{m}^{\prime} (\kappa r)\right),\\ {e}_{z}^{\left(\nu \right)} & = & {{AJ}}_{m}(\kappa r),\end{array}\end{eqnarray*}$$

while for r > a:

$$\begin{eqnarray*}\begin{array}{rcl}{e}_{r}^{\left(\nu \right)} & = & \displaystyle \frac{1}{{\rm{i}}\sigma }\left[\beta \left({{CJ}}_{m}^{\prime} (\sigma r)+{{EY}}_{m}^{\prime} (\sigma r)\right)+\displaystyle \frac{{\rm{i}}m\omega }{\sigma r}\left({{DJ}}_{m}\left(\sigma r\right)\right.\right.\\ & & \left.\left.+{{FY}}_{m}\left(\sigma r\right)\right)\right],\\ {e}_{\varphi }^{\left(\nu \right)} & = & \displaystyle \frac{1}{{\rm{i}}\sigma }\left[{\rm{i}}\displaystyle \frac{\nu \beta }{\sigma r}\left({{CJ}}_{m}\left(\sigma r\right)+{{EY}}_{m}\left(\sigma r\right)\right)-\omega \left({{DJ}}_{m}^{\prime} (\sigma r)\right.\right.\\ & & \left.\left.+{{FY}}_{m}^{\prime} (\sigma r)\right)\right],\\ {e}_{z}^{\left(\nu \right)} & = & {{CJ}}_{m}\left(\sigma r\right)+{{EY}}_{m}\left(\sigma r\right),\end{array}\end{eqnarray*}$$

where Y_m denote Bessel functions of the second kind. The coefficients C, D, E and F are related to A and B as follows

$$\begin{eqnarray*}\begin{array}{rcl}C & = & -\displaystyle \frac{\pi a{\sigma }^{2}}{2{\epsilon }_{0}{n}_{2}^{2}}\left({{AL}}_{2}+{\rm{i}}{{BV}}_{2}\right),\\ D & = & {\rm{i}}\displaystyle \frac{\pi a{\sigma }^{2}}{2}\left(A{\mu }_{0}{V}_{2}+{\rm{i}}{{BM}}_{2}\right),\\ E & = & \displaystyle \frac{\pi a{\sigma }^{2}}{2{\epsilon }_{0}{n}_{2}^{2}}\left({{AL}}_{1}+{\rm{i}}{{BV}}_{1}\right),\\ F & = & -{\rm{i}}\displaystyle \frac{\pi a{\sigma }^{2}}{2}\left(A{\mu }_{0}{V}_{1}+{\rm{i}}{{BM}}_{1}\right),\end{array}\end{eqnarray*}$$

with

$$\begin{eqnarray*}\begin{array}{rcl}{V}_{1} & = & \displaystyle \frac{m\beta }{a\omega {\mu }_{0}{\kappa }^{2}{\sigma }^{2}}\left({n}_{2}^{2}-{n}_{1}^{2}\right){J}_{m}\left(\sigma a\right){J}_{m}\left(\kappa a\right),\\ {V}_{2} & = & \displaystyle \frac{m\beta }{a\omega {\mu }_{0}{\kappa }^{2}{\sigma }^{2}}\left({n}_{2}^{2}-{n}_{1}^{2}\right){Y}_{m}\left(\sigma a\right){J}_{m}\left(\kappa a\right),\\ {M}_{1} & = & \displaystyle \frac{1}{\kappa }{J}_{m}\left(\sigma a\right){J}_{m}^{\prime} \left(\kappa a\right)-\displaystyle \frac{1}{\sigma }{J}_{m}^{\prime} \left(\sigma a\right){J}_{m}\left(\kappa a\right),\\ {M}_{2} & = & \displaystyle \frac{1}{\kappa }{Y}_{m}\left(\sigma a\right){J}_{m}^{\prime} \left(\kappa a\right)-\displaystyle \frac{1}{\sigma }{Y}_{m}^{\prime} \left(\sigma a\right){J}_{m}\left(\kappa a\right),\\ {L}_{1} & = & \displaystyle \frac{{\epsilon }_{0}{n}_{1}^{2}}{\kappa }{J}_{m}\left(\sigma a\right){J}_{m}^{\prime} \left(\kappa a\right)-\displaystyle \frac{{\epsilon }_{0}{n}_{2}^{2}}{\sigma }{J}_{m}^{\prime} \left(\sigma a\right){J}_{m}\left(\kappa a\right),\\ {L}_{2} & = & \displaystyle \frac{{\epsilon }_{0}{n}_{1}^{2}}{\kappa }{Y}_{m}\left(\sigma a\right){J}_{m}^{\prime} \left(\kappa a\right)-\displaystyle \frac{{\epsilon }_{0}{n}_{2}^{2}}{\sigma }{Y}_{m}^{\prime} \left(\sigma a\right){J}_{m}\left(\kappa a\right).\end{array}\end{eqnarray*}$$

In the single-mode approximation, a guided mode is completely specified by the frequency ω, the direction of propagation f = ±1 and the polarization p = ±1. By contrast, at first glance, this is not the case for the radiative modes any longer. Once β, ω and m are fixed, we are left with two constants A and B, and a normalization condition will only determine one constant. We must therefore separate these into two modes. For instance, we can just set A = 0 for one mode and B = 0 for the other one. We want, however, the two modes to be orthogonal to each other. An alternative method consists in setting B = piη A with the parameter p = ±1, then imposing an orthogonality condition between ${{\bf{e}}}^{\left(p=+1\right)}$ and ${{\bf{e}}}^{\left(p=-1\right)}$. Explicitly, this condition is written:

$$\begin{eqnarray*}&&\begin{array}{l}{\displaystyle \int }_{0}^{2\pi }d\varphi {\displaystyle \int }_{0}^{\infty }n{\left(r\right)}^{2}\left[{{\bf{e}}}^{\left(\nu \right)}\cdot {{\bf{e}}}^{\left(\nu ^{\prime} \right)}\right]{rdr}\\ \quad =\,{\delta }_{{pp}^{\prime} }{\delta }_{{mm}^{\prime} }\delta \left(\omega -\omega ^{\prime} \right)\delta \left(\beta -\beta ^{\prime} \right).\end{array}\end{eqnarray*}$$

If we consider the vacuum surrounding with the index n₂ = 1, this leads to:

$$\begin{eqnarray*}\begin{array}{rcl}\eta & = & \sqrt{\displaystyle \frac{{L}_{1}^{2}+{L}_{2}^{2}+{\epsilon }_{0}{\mu }_{0}\left({V}_{1}^{2}+{V}_{2}^{2}\right)}{{V}_{1}^{2}+{V}_{2}^{2}+\tfrac{{\epsilon }_{0}}{{\mu }_{0}}\left({M}_{1}^{2}+{M}_{2}^{2}\right)}},\\ 1 & = & \displaystyle \frac{2\pi \omega }{{\sigma }^{2}}\left[\left({\left|C\right|}^{2}+{\left|E\right|}^{2}\right)+{c}^{2}\left({\left|D\right|}^{2}+{\left|F\right|}^{2}\right)\right].\end{array}\end{eqnarray*}$$

The second normalization equation allows us to calculate the form of $\left|A\right|$

$$\begin{eqnarray*}&&\left|A\right|=\displaystyle \frac{1}{\sigma a}{\left[\displaystyle \frac{2}{\omega {\pi }^{3}}\displaystyle \sum _{j=1,2}\left(\displaystyle \frac{1}{{\epsilon }_{0}}{\left|{L}_{j}-p\eta {V}_{j}\right|}^{2}+{\left|{\mu }_{0}{V}_{j}-p\eta {M}_{j}\right|}^{2}\right)\right]}^{-\tfrac{1}{2}}.\end{eqnarray*}$$

This shows that $\nu =\left(\omega ,\beta ,m,p\right)$ completely determines a radiative mode.

With the definitions ${\sigma }_{{MN}}\equiv \left|M\right\rangle \left\langle N\right|$, ${\omega }_{{MN}}\equiv \left({E}_{M}-{E}_{N}\right)/{\hslash }$, k ≡ ω/c, the Hamiltonian of the full system consisting of the atom and the electric field takes the form H = H_at + H_f +H_int with

$$\begin{eqnarray*}\begin{array}{rcl}{H}_{f} & = & \displaystyle \sum _{\mu }{\hslash }\omega {a}_{\mu }^{\dagger }{a}_{\mu }+\displaystyle \sum _{\nu }{\hslash }\omega {a}_{\nu }^{\dagger }{a}_{\nu },\\ {H}_{{at}} & = & \displaystyle \sum _{m}{\hslash }{\omega }_{M}{\sigma }_{{MM}},\\ {H}_{{int}} & = & -\vec{D}\cdot \vec{E},\end{array}\end{eqnarray*}$$

where $\vec{D}={\sum }_{M,N}{\vec{d}}_{{MN}}{\sigma }_{{MN}}$ and $\vec{E}\left(\vec{r}\right)={\vec{E}}_{g}\left(\vec{r}\right)+{\vec{E}}_{r}\left(\vec{r}\right)$ are the atomic dipole operator and the total electric field operator, respectively. Switching to the interaction picture relative to ${H}_{0}\equiv {H}_{{at}}+{H}_{f}$, and resorting to the rotating wave approximation we get the interaction Hamiltonian (note that the states $\left|M\right\rangle$ are ordered by increasing energies)

$$\begin{eqnarray*}\begin{array}{rcl}{\tilde{H}}_{{int}}\left(t\right) & \approx & {\rm{i}}{\hslash }\displaystyle \sum _{N\lt M}\displaystyle \sum _{\mu }\,{G}_{\mu {MN}}{\sigma }_{{MN}}{a}_{\mu }{e}^{-{\rm{i}}\left(\omega -{\omega }_{{MN}}\right)t}+{\rm{h}}.{\rm{c}}.\\ & & +{\rm{i}}{\hslash }\displaystyle \sum _{N\lt M}\displaystyle \sum _{\nu }\,{G}_{\nu {MN}}{\sigma }_{{MN}}{a}_{\nu }{e}^{-{\rm{i}}\left(\omega -{\omega }_{{MN}}\right)t}+{\rm{h}}.{\rm{c}}.\end{array}\end{eqnarray*}$$

where we introduced

$$\begin{eqnarray*}\begin{array}{rcl}{G}_{\mu {MN}} & \equiv & -\sqrt{\displaystyle \frac{\omega \beta ^{\prime} }{4\pi {\epsilon }_{0}{\hslash }}}\left({\vec{d}}_{{MN}}\cdot {\vec{e}}^{\left(\mu \right)}\right){e}^{{\rm{i}}\left(f\beta z+p\varphi \right)},\\ {G}_{\nu {MN}} & \equiv & -\sqrt{\displaystyle \frac{\omega }{4\pi {\epsilon }_{0}{\hslash }}}\left({\vec{d}}_{{MN}}\cdot {\vec{e}}^{\left(\nu \right)}\right){e}^{{\rm{i}}\left(\beta z+m\varphi \right)}.\end{array}\end{eqnarray*}$$

For simplicity, from now on we shall use λ to denote either guided modes, i.e. $\lambda =\left(\omega ,f=\pm ,p=\pm \right)$, or radiative modes, i.e. $\lambda =\left(\omega ,\beta ,m,p\right)$, and use ${\sum }_{\lambda }$ to represent the sum, either discrete or continuous, of these modes, whence

$$\begin{eqnarray}\begin{array}{rcl}{\tilde{H}}_{{int}}\left(t\right) & \approx & {\rm{i}}{\hslash }\displaystyle \sum _{N\lt M}\displaystyle \sum _{\lambda }\,{G}_{\lambda {MN}}{\sigma }_{{MN}}{a}_{\lambda }{e}^{-{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MN}}\right)t}\\ & & -{\rm{i}}{\hslash }\displaystyle \sum _{N\lt M}\displaystyle \sum _{\lambda }\,{G}_{\lambda {MN}}^{* }{a}_{\lambda }^{\dagger }{\sigma }_{{NM}}{e}^{+{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MN}}\right)t}.\end{array}\end{eqnarray} \tag{ C.1 }$$

From equation (C.1), we get the Heisenberg equations for the field and atomic operators, a_λ and σ_PQ

$$\begin{eqnarray}&&{\partial }_{t}{a}_{\lambda }=-\displaystyle \sum _{N\lt M}\,{G}_{\lambda {MN}}^{* }{\sigma }_{{NM}}{e}^{+{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MN}}\right)t},\end{eqnarray} \tag{ C.2 }$$

$$\begin{eqnarray}\begin{array}{rcl}{\partial }_{t}{\sigma }_{{PQ}} & = & \displaystyle \sum _{M\lt Q}\displaystyle \sum _{\lambda }\,{G}_{\lambda {QM}}{\sigma }_{{PM}}\left(t\right){a}_{\lambda }\left(t\right){e}^{-{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{QM}}\right)t}\\ & & -\displaystyle \sum _{P\lt M}\displaystyle \sum _{\lambda }\,{G}_{\lambda {MP}}{\sigma }_{{MQ}}\left(t\right){a}_{\lambda }\left(t\right){e}^{-{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MP}}\right)t}\\ & & -\displaystyle \sum _{Q\lt M}\displaystyle \sum _{\lambda }\,{G}_{\lambda {MQ}}^{* }{a}_{\lambda }^{\dagger }\left(t\right){\sigma }_{{PM}}\left(t\right){e}^{+{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MQ}}\right)t}\\ & & +\displaystyle \sum _{M\lt P}\displaystyle \sum _{\lambda }\,{G}_{\lambda {PM}}^{* }{a}_{\lambda }^{\dagger }\left(t\right){\sigma }_{{MQ}}\left(t\right){e}^{+{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{PM}}\right)t}.\end{array}\end{eqnarray} \tag{ C.3 }$$

We can eliminate the field degree of freedom by inserting the formal solution of equation (C.2)

$$\begin{eqnarray*}&&{a}_{\lambda }\left(t\right)={a}_{\lambda }\left({t}_{0}\right)-\displaystyle \sum _{N\lt M}\,{G}_{\lambda {MN}}^{* }{\int }_{{t}_{0}}^{t}{\sigma }_{{NM}}\left(s\right){e}^{+{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MN}}\right)s}\end{eqnarray*}$$

into equation (C.3). Then performing Markov approximation [31] and using

$$\begin{eqnarray*}&&{\int }_{0}^{\infty }{{dse}}^{\pm {\rm{i}}{\rm{\Delta }}s}=\pi \delta \left({\rm{\Delta }}\right)\pm {\rm{i}}{ \mathcal P }\left(\displaystyle \frac{1}{{\rm{\Delta }}}\right)\end{eqnarray*}$$

we get

$$\begin{eqnarray*}\begin{array}{rcl}{\partial }_{t}{\sigma }_{{PQ}} & \approx & -\left\{{{\rm{\Gamma }}}_{{PQ}}+{\rm{i}}{{\rm{\Delta }}}_{{PQ}}\right\}{\sigma }_{{PQ}}\left(t\right)\\ & & +\,{\delta }_{{PQ}}\displaystyle \sum _{P\lt M}{\gamma }_{{MP}}{\sigma }_{{MM}}\left(t\right)+{\xi }_{{PQ}}\left(t\right),\end{array}\end{eqnarray*}$$

where we introduced the different decay rates and energy shifts due to spontaneous emission into the modes of the fibre

$$\begin{eqnarray*}\begin{array}{rcl}{{\rm{\Gamma }}}_{{PQ}} & \equiv & \displaystyle \frac{1}{2}\left({{\rm{\Gamma }}}_{P}+{{\rm{\Gamma }}}_{Q}\right)\\ {{\rm{\Delta }}}_{{PQ}} & \equiv & {{\rm{\Delta }}}_{P}-{{\rm{\Delta }}}_{Q}\\ {{\rm{\Gamma }}}_{P} & \equiv & \displaystyle \sum _{M\lt P}{\gamma }_{{PM}}\\ {\gamma }_{{PQ}} & \equiv & 2\pi \displaystyle \sum _{\lambda }{\left|{G}_{\lambda {PQ}}\right|}^{2}\delta \left({\omega }_{\lambda }-{\omega }_{{PQ}}\right)\\ {{\rm{\Delta }}}_{P} & \equiv & { \mathcal P }\left(\displaystyle \sum _{M\lt P}\displaystyle \sum _{\lambda }\displaystyle \frac{{\left|{G}_{\lambda {PM}}\right|}^{2}}{{\omega }_{\lambda }-{\omega }_{{PM}}}-\displaystyle \frac{{\left|{G}_{\lambda {PM}}^{\left(0\right)}\right|}^{2}}{{\omega }_{\lambda }-{\omega }_{{PM}}}\right)+{{\rm{\Delta }}}_{P}^{\left(0\right)}\end{array}\end{eqnarray*}$$

and the associated Langevin forces

$$\begin{eqnarray}\begin{array}{rcl}{\xi }_{{PQ}}\left(t\right) & \equiv & \displaystyle \sum _{M\lt Q}\displaystyle \sum _{\lambda }\,{G}_{\lambda {QM}}{\sigma }_{{PM}}\left(t\right){a}_{\lambda }\left({t}_{0}\right){e}^{-{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{QM}}\right)t}\\ & & -\displaystyle \sum _{P\lt M}\displaystyle \sum _{\lambda }\,{G}_{\lambda {MP}}{\sigma }_{{MQ}}\left(t\right){a}_{\lambda }\left({t}_{0}\right){e}^{-{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MP}}\right)t}\\ & & -\displaystyle \sum _{Q\lt M}\displaystyle \sum _{\lambda }\,{G}_{\lambda {MQ}}^{* }{a}_{\lambda }^{\dagger }\left({t}_{0}\right){\sigma }_{{PM}}\left(t\right){e}^{+{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{MQ}}\right)t}\\ & & +\displaystyle \sum _{M\lt P}\displaystyle \sum _{\lambda }\,{G}_{\lambda {PM}}^{* }{a}_{\lambda }^{\dagger }\left({t}_{0}\right){\sigma }_{{MQ}}\left(t\right){e}^{+{\rm{i}}\left({\omega }_{\lambda }-{\omega }_{{PM}}\right)t}.\end{array}\end{eqnarray} \tag{ C.4 }$$

Note that due to the normal operator ordering in equation (C.4), $\left\langle {\xi }_{{PQ}}\right\rangle =0$.

From the relation ${\rho }_{{QP}}\left(t\right)=\mathrm{Tr}\left[\rho \left({t}_{0}\right){\sigma }_{{PQ}}\left(t\right)\right]$ one immediately deduces the evolution equation for the density matrix

$$\begin{eqnarray*}&&{\partial }_{t}{\rho }_{{QP}}\approx -\left({{\rm{\Gamma }}}_{{PQ}}+{\rm{i}}{{\rm{\Delta }}}_{{PQ}}\right){\rho }_{{QP}}\left(t\right)+{\delta }_{{PQ}}\displaystyle \sum _{P\lt M}{\gamma }_{{MP}}{\rho }_{{MM}}\left(t\right).\end{eqnarray*}$$

In particular, for coherences $\left(P\ne Q\right)$, we obtain

$$\begin{eqnarray*}&&{\rho }_{{QP}}\left(t\right)\approx {e}^{-\left({{\rm{\Gamma }}}_{{PQ}}+{\rm{i}}{{\rm{\Delta }}}_{{PQ}}\right)\left(t-{t}_{0}\right)}{\rho }_{{QP}}\left({t}_{0}\right).\end{eqnarray*}$$

In order to calculate the rates of spontaneous emission for levels $10{s}_{1/2}$ and $10{p}_{1/\mathrm{2,3}/2}$, we need energies and transition dipole moments involving s, p and d lower levels. Regarding energies, we take experimental values from the NIST database [32]. Transition dipole moments are calculated using the Cowan codes [33].

The vector associated with the dipole operator is expressed as irreducible tensors ${\hat{d}}_{q}$ (q = 0, ±1), such that ${\hat{d}}_{0}={\hat{d}}_{z}$ and ${\hat{d}}_{\pm 1}=\mp ({\hat{d}}_{x}\pm {\rm{i}}{\hat{d}}_{y})\sqrt{2}$. Their matrix elements in the coupled atomic basis $\{| n{{\ell }}_{j}{m}_{j}\rangle \}$ read

$$\begin{eqnarray}&&\begin{array}{l}\left\langle n^{\prime} {\ell }{{\prime} }_{j^{\prime} }{m}_{j}^{\prime} \right|{\hat{d}}_{q}\left|n{{\ell }}_{j}{m}_{j}\right\rangle =e{(-1)}^{j+{\ell }^{\prime} +s}\sqrt{(2j+1)}\\ \quad \times \,\left\{\begin{array}{ccc}{\ell } & s & j\\ j^{\prime} & 1 & {\ell }^{\prime} \end{array}\right\}\left\langle n^{\prime} {\ell }^{\prime} \parallel \hat{{\bf{r}}}\parallel n{\ell }\right\rangle {C}_{{{jm}}_{j}1q}^{j^{\prime} {m}_{j}^{\prime} },\end{array}\end{eqnarray} \tag{ D.1 }$$

where e is the absolute value of the electron charge, $\left\{\begin{array}{ccc}a & b & c\\ d & e & f\end{array}\right\}$ is a Wigner 6j symbol, and ${C}_{a\alpha b\beta }^{c\gamma }$ a Clebsch–Gordan coefficient [34]. The quantity $\langle n^{\prime} {\ell }^{\prime} \parallel \hat{{\bf{r}}}\parallel n{\ell }\rangle$ is the reduced matrix element of the position operator of the outermost electron. In our calculations, it is supposed to be independent from j and $j^{\prime}$.

Table D1 contains the quantities $\langle n^{\prime} {\ell }^{\prime} \parallel \hat{{\bf{r}}}\parallel 10s\rangle$ and $\langle n^{\prime} {\ell }^{\prime} \parallel \hat{{\bf{r}}}\parallel 10p\rangle$ relevant for our calculation. They give radiative lifetimes of 0.855, 8.58 and 8.56 μs for $10s$, $10{p}_{1/2}$ and $10{p}_{3/2}$ respectively, which are in correct agreement with the values reported in [35].

Table D1.  Reduced matrix elements $\langle n^{\prime} {\ell }^{\prime} \parallel \hat{{\bf{r}}}\parallel n{\ell }\rangle$ in units of ea₀ (e is the absolute value of the electron charge and a₀ is the Bohr radius), for $n{\ell }=10s$ and $10p$.

$n^{\prime}$	$\langle n^{\prime} p\parallel \hat{{\bf{r}}}\parallel 10s\rangle$	$\langle n^{\prime} s\parallel \hat{{\bf{r}}}\parallel 10p\rangle$	$\langle n^{\prime} d\parallel \hat{{\bf{r}}}\parallel 10p\rangle$
3	0.145 809 5	0.020 545 8	−0.015 973 4
4	0.357 453 2	0.091 422 7	−0.085 674 7
5	0.739 883 4	0.240 707 9	−0.263 861 7
6	1.488 610 2	0.538 127 5	−0.694 403 3
7	3.236 754 6	1.174 545 4	−1.859 137 1
8	9.104 099 1	2.803 803 7	−6.129 887 2
9	71.148 579 0	8.958 641 1	160.001 178 7
10		93.200 142 5