Cooperative Decay of an Ensemble of Atoms in a One-Dimensional Chain with a Single Excitation

1. Introduction

Cooperative spontaneous emission by N excited two-level atoms has been extensively studied, since the seminal work by Dicke in 1954 [1] and Lehmberg in 1970 [2]. Whereas superradiance i.e. enhanced spontaneous emission due to positive interference between the emitters, has been well understood [3,4], subradiance, i.e. inhibited emission due to negative interference between the emitters, is more elusive and difficult to observe [5,6,7]. Nevertheless, the number of studies on subradiance has seen a large increase in the last few years, as it offers the opportunity of storing photons in emitter ensembles for times longer than the single emitter lifetime [8,9,10].

In disordered system, the cooperative decay must be studied numerically or by solving the dynamics of an initially excited ensemble or by calculating the eigenvalues of the system. Things are more simple when atoms form ordered arrays, easier to treat theoretically.

In this paper we consider a one-dimensional array of atoms in the single-excitation configuration, i.e. with only one atom excited among N [11,12]. Starting from the effective non-hermitial Hamiltonian, which include an imaginary part describing the spontaneous decay, we study the spectrum of the decay rates, obtaining the exact solution for a infinite chain and approximated expressions for a finite chain of N atoms. We compare the results for a infinite and finite chains obtained from the scalar light and the vectorial light models. We see that the scalar light model still provides a good approximation to the vectorial model for large lattice spacing.

2. Scalar Model

We consider N two-level atoms with the same atomic transition frequency ${\omega }_a$, linewidth $\mathsf{\Gamma }$ and dipole d. The atoms are prepared in a single-excitation state, by absorbing a photon with wave vector ${\mathbf{k}}_0$; $|g_j\rangle$ and $|e_j\rangle$ are the ground and excited states, respectively, of the j-th atom, $j=1,\dots ,N$, which is placed at position ${\mathbf{r}}_j$. We consider the single excitation effective Hamiltonian in the scalar approximation. If we assume that only one photon is present, when tracing over the radiation degrees of freedom, the dynamics of the atomic system can be described by the non-Hermitian Hamiltonian

$$\begin{array}{ccc}\hfill \widehat{H}& =& -i\frac{\hslash }{2}\sum _{j,m}{G}_{jm}{\widehat{\sigma }}_j^{†}{\widehat{\sigma }}_m,\hfill \end{array}$$

(1)

where ${\widehat{\sigma }}_j=|g_j\rangle \langle e_j|$ and ${\widehat{\sigma }}_j^{†}=|e_j\rangle \langle g_j|$ are the lowering and raising operators and

$$G_{jm}=\left\{\begin{array}{cc}{\mathsf{\Gamma }}_{jm}-i{\mathsf{\Omega }}_{jm}\hfill & \mathrm{if}j\ne m,\hfill \\ 1\hfill & \mathrm{if}j=m,\hfill \end{array}\right.$$

(2)

where

$${\mathsf{\Gamma }}_{jm}=\mathsf{\Gamma }\frac{sin\left(k_0{r}_{jm}\right)}{k_0{r}_{jm}}\mathrm{and}{\mathsf{\Omega }}_{jm}=\mathsf{\Gamma }\frac{cos\left(k_0{r}_{jm}\right)}{k_0{r}_{jm}}.$$

(3)

$\widehat{H}$ contains both real and imaginary parts, which takes into account that the excitation is not conserved since it can leave the system by emission. We focus our attention on the decay term ${\mathsf{\Gamma }}_{jm}$. It can be obtained as the angular average of the radiation field propagating between the two atomic positions $r_j$ and $r_m$ with wave-vector $\mathbf{k}=k_0(sin\theta cos\varphi ,sin\theta sin\varphi ,cos\theta )$,

$${\mathsf{\Gamma }}_{jm}=\frac{\mathsf{\Gamma }}{2}{〈e^{-i\mathbf{k}·({\mathbf{r}}_j-{\mathbf{r}}_m)}+\mathrm{c}.\mathrm{c}.〉}_{\mathsf{\Omega }}$$

(4)

where the angular average is defined as

$${〈f(\theta ,\varphi )〉}_{\mathsf{\Omega }}=\frac{1}{4\pi }{\int }_0^{2\pi }d\varphi {\int }_0^{\pi }sin\theta f(\theta ,\varphi )d\theta .$$

Eq.(4) allows to factorize ${\mathsf{\Gamma }}_{jm}$ in the product of two terms, before doing the angular average.

Let now consider N atoms placed along a linear chain with period d, i.e. ${\mathbf{r}}_j=d(j-1){\widehat{\mathbf{e}}}_z$, with $j=1,\cdots ,N$. Then, we can write

$$e^{-i\mathbf{k}·{\mathbf{r}}_j}=E_j=e^{-i{k}_{0}d(j-1)cos\theta }$$

(5)

and

$${\mathsf{\Gamma }}_{jm}=\frac{\mathsf{\Gamma }}{2}{〈E_j{E}_m^{*}+\mathrm{c}.\mathrm{c}.〉}_{\mathsf{\Omega }}.$$

(6)

We calculate the Fourier transform of ${\mathsf{\Gamma }}_{jm}$ as

$${\mathsf{\Gamma }}_k=\frac{1}{N}\sum _{j=1}^N\sum _{m=1}^N{\mathsf{\Gamma }}_{jm}{e}^{ikd(j-m)}$$

(7)

with $k\in [0,2\pi /d]$. While (7) is the discrete Fourier transform of ${\mathsf{\Gamma }}_{jm}$ only in the case of the infinite chain, such that the system is periodic, nevertheless the wave vector k is still a good label for the modes when N is sufficiently large [11]. The true spectrum of the system can be obtained only by calculating numerically the eigenvalues of the matrix ${\mathsf{\Gamma }}_{jm}$, as done in all the published studied on the cooperative decay. Instead, (7) gives an approximated expression of the spectrum as a function of the continuous parameter k, from which an analytic expression may be obtained, providing a valuable and explicit information on the cooperative decay rate. By proceeding with this approach, we use (4) and write

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k& =& \frac{\mathsf{\Gamma }}{N}\langle |F_k\left(\theta \right){{|}^2\rangle }_{\mathsf{\Omega }}\hfill \end{array}$$

(8)

where

$$\begin{array}{ccc}\hfill |F_k{\left(\theta \right)|}^2& =& {\left|\sum _{j=1}^N{e}^{i[k-k_{0}cos\theta ]d(j-1)}\right|}^2=\frac{{sin}^2[(k-k_{0}cos\theta )dN/2]}{{sin}^2[(k-k_{0}cos\theta )d/2]}\hfill \end{array}$$

(9)

and

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k& =& \frac{\mathsf{\Gamma }}{4\pi N}{\int }_0^{2\pi }d\varphi {\int }_0^{\pi }sin\theta {\left|F_k\left(\theta \right)\right|}^{2}d\theta =\frac{\mathsf{\Gamma }}{k_{0}dN}{\int }_{(k-k_0)d/2}^{(k+k_0)d/2}\frac{{sin}^2\left(Nt\right)}{{sin}^{2}t}dt\hfill \end{array}$$

(10)

where we changed the integration variable from $\theta$ to $t=(k-k_{0}cos\theta )d/2$. For $k=k_0$ and in the limit $k_{0}d\ll \pi$ we have ${\mathsf{\Gamma }}_{k_0}\approx \mathsf{\Gamma }N$ i.e. the superradiant decay. For large N, we can approximate in the integral of Eq.(10)

$$\frac{{sin}^2\left(Nt\right)}{{sin}^{2}t}\approx N^2\sum _{m=0}^{\infty }{\mathrm{sinc}}^2\left[\left(t-m\pi \right)N\right]$$

(11)

where $\mathrm{sinc}\left(x\right)=sinx/x$, so that

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k& =& =\frac{\mathsf{\Gamma }N}{k_{0}d}\sum _{m=0}^{\infty }{\int }_{(k-k_0)d/2}^{(k+k_0)d/2}{\mathrm{sinc}}^2\left[\left(t-m\pi \right)N\right]dt\hfill \end{array}$$

(12)

In the limit $N\to \infty$, we approximate

$${\mathrm{sinc}}^2\left[\left(t-m\pi \right)N\right]\to \frac{\pi }{N}\delta (t-m\pi )$$

(13)

where $\delta \left(x\right)$ is the Dirac’s delta function. Thus,

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k& =& \frac{\mathsf{\Gamma }\pi }{k_{0}d}\sum _{m=0}^{\infty }\mathsf{\Pi }[2m\pi -k_{0}d<kd<2m\pi +k_{0}d]\hfill \end{array}$$

(14)

where $\mathsf{\Pi }(a<x<b)$ is the rectangular function, equal to 1 for $a<x<b$ and 0 elsewhere. Eq.(14) is the solution for an infinite chain, where k is the true index of the modes and the sum over m represents the decomposition in Brillouin zones. We can see that if $k_{0}d<\pi$ (i.e. $d/{\lambda }_0<0.5$) then ${\mathsf{\Gamma }}_k=0$ for $k_0<k<2\pi /d-k_0$ (i.e.subradiance). Instead, for $0<k<k_0$ and $2\pi /d-k_0<k<2\pi /d$, ${\mathsf{\Gamma }}_k=\mathsf{\Gamma }(\pi /k_{0}d)>\mathsf{\Gamma }$ i.e. enhanced radiance[11,12]) (see Figure 1a). For larger lattice spacing, for $\pi <k_{0}d<2\pi$ we have ${\mathsf{\Gamma }}_k=\mathsf{\Gamma }(\pi /k_{0}d)<\mathsf{\Gamma }$ for $0<k<2\pi /d-k_0$ and $k_0<k<2\pi /d$, whereas ${\mathsf{\Gamma }}_k=\mathsf{\Gamma }(2\pi /k_{0}d)>\mathsf{\Gamma }$ for $2\pi /d-k_0<k<k_0$ (see Figure 1b). More generally, for $m\pi <k_{0}d<(m+1)\pi$ with $m=2,3,\cdots$, ${\mathsf{\Gamma }}_k=\mathsf{\Gamma }(m\pi /k_{0}d)<\mathsf{\Gamma }$ for $0<k<2m\pi /d-k_0$ and ${\mathsf{\Gamma }}_k=\mathsf{\Gamma }[(m+1)\pi /k_{0}d]>\mathsf{\Gamma }$ for $2m\pi /d-k_0<k<2\pi$ (see (Figure 1c) and (Figure 1d)). Hence, for each subsequent interval of $k_{0}d$ of width equal to $\pi$, we have two regions of the spectrum, one where ${\mathsf{\Gamma }}_k$ is less than $\mathsf{\Gamma }$ and the other one where is larger than $\mathsf{\Gamma }$. The difference between these two values is $\Delta \mathsf{\Gamma }=\pi /k_{0}d$. The only case of pure subradiance, with zero decay rate, is for $k_{0}d<\pi$.

In Figure 2 we plot ${\mathsf{\Gamma }}_k$ vs $k_{0}d$ for $k=0$ (blue line) and $k=\pi /d$ (red line), obtained from Eq.(12) with $N=10$ (Figure 2a) and $N=100$ (Figure 2b). We see that for an infinite chain ($N\to \infty$) and $k_{0}d<2\pi$, ${\mathsf{\Gamma }}_{k=0}=\mathsf{\Gamma }\pi /k_{0}d$, whereas for $\pi <k_{0}d<2\pi$, ${\mathsf{\Gamma }}_{k=\pi /d}=2\pi /k_{0}d$. The dashed line in Figure 2a is the analytic solution (14) for infinite chain (i.e. $N\to \infty$).

The case a finite chain (i.e. with N finite) requires the evaluation of the integral in Eq.(12). An approximated analytic expression can be obtained substituting the function ${\mathrm{sinc}}^2\left(x\right)$ with the Lorentzian function $1/(1+x^2)$, which has the same peak value of unity and the same normalization value of $\pi$. Then, the integral in Eq.(12) yields

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k& =& \frac{\mathsf{\Gamma }}{k_{0}d}\sum _{m=0}^{\infty }\left[arctan\left(b_m\right)-arctan\left(a_n\right)\right]\hfill \end{array}$$

(15)

where $a_m=[(k-k_0)d/2-m\pi ]N$ and $b_m=[(k+k_0)d/2-m\pi ]N$. This expression allows to evaluate the subradiant decay rate in the limit of large N, using the identity $arctan\left(z\right)=\pm \pi /2-arctan(1/z)$ (where the positive sign is for $z>0$ and the negative sign for $z<0$). Then, for $k_{0}d<\pi$ and the band interval $k_0<k<2\pi /d-k_0$, in the first two terms of the sum in Eq.(15) $a_0,b_0>0$, and $a_1,b_1<0$, so that, for large N and far from the band edges,

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k& =& \frac{4\mathsf{\Gamma }}{N}\left[\frac{1}{{\left(kd\right)}^2-{\left(k_{0}d\right)}^2}+\frac{1}{{(kd-2\pi )}^2-{\left(k_{0}d\right)}^2}\right].\hfill \end{array}$$

(16)

We observe the dependence on $1/N$, typical of subradiance [6]. In particular, at the center of the band, $k=\pi /d$,

$${\mathsf{\Gamma }}_{k=\pi /d}=\frac{8\mathsf{\Gamma }}{N[{\pi }^2-{\left(k_{0}d\right)}^2]},$$

(17)

so for $k_{0}d\ll \pi$ the subradiant decay is independent on the lattice constant d. Conversely, superradiance is obtained for $k=0$,

$${\mathsf{\Gamma }}_{k=0}=\frac{2\mathsf{\Gamma }}{k_{0}d}arctan\left(\frac{k_{0}dN}{2}\right).$$

(18)

In the limit $k_{0}d\ll 2\pi /N$, ${\mathsf{\Gamma }}_{k=0}\approx \mathsf{\Gamma }N$, as it can be observed in Figure 2a.

Finally, it is interesting to see how ${\mathsf{\Gamma }}_k$ is approximating the true spectrum of a finite chain. Figure shows ${\mathsf{\Gamma }}_k$, as calculated from Eq.(12), and the N eigenvalues ${\lambda }_i$ of the $N\times N$ matrix ${\mathsf{\Gamma }}_{jm}$, ordered from the largest to the smallest, as function of $kd=\pi (i-1/2)/N$, for $i=1,\cdots ,N$ for $k_{0}d=\pi /2$, $N=10$ (a) and $N=50$ (b). We see that ${\mathsf{\Gamma }}_k$ reproduces rather satisfactorily the features of the spectrum.

3. Vectorial Model

We now extend the previous analysis to the vectorial model, taking into account the polarization of the electromagnetic field. The non-hermitian Hamiltonian is now

$$\begin{array}{ccc}\hfill \widehat{H}& =& -i\frac{\hslash }{2}\sum _{\alpha ,\beta }\sum _{j,m}{G}_{\alpha ,\beta }({\mathbf{r}}_j-{\mathbf{r}}_m){\widehat{\sigma }}_{j,\alpha }^{†}{\widehat{\sigma }}_{m,\beta }.\hfill \end{array}$$

(19)

where $\alpha ,\beta =(x,y,z)$, ${\widehat{\sigma }}_{j,\alpha }=|g_j\rangle \langle e_{j,\alpha }|$ and ${\widehat{\sigma }}_{j,\alpha }^{†}=|e_{j,\alpha }\rangle \langle g_j|$ are the lowering and raising operators between the ground and the three excited states $|e_{j,\alpha }\rangle$ and [2,13]

$$G_{\alpha ,\beta }\left(\mathbf{r}\right)=\frac{3}{2}\frac{e^{i{k}_{0}r}}{i{k}_{0}r}\left\{{\delta }_{\alpha ,\beta }-{\widehat{n}}_{\alpha }{\widehat{n}}_{\beta }+[{\delta }_{\alpha ,\beta }-3{\widehat{n}}_{\alpha }{\widehat{n}}_{\beta }][i/\left(k_{0}r\right)-1/{\left(k_{0}r\right)}^2]\right\}$$

(20)

with $r=\left|\mathbf{r}\right|$ and ${\widehat{n}}_{\alpha }$ being the components of the unit vector $\widehat{\mathbf{n}}=\mathbf{r}/r$. We consider the linear chain with period d, i.e. ${\mathbf{r}}_j=d(j-1){\widehat{\mathbf{e}}}_z$, with $j=1,\cdots ,N$, and all the dipoles aligned and with an angle $\delta$ with respect to the chain’s axis, so that ${\widehat{n}}_{\alpha }={\widehat{n}}_{\beta }=cos\delta$ and

$$G_{\alpha ,\alpha }\left(\mathbf{r}\right)=\frac{3\mathsf{\Gamma }}{2}\frac{e^{i{k}_{0}r}}{i{k}_{0}r}\left\{{sin}^2\delta +[1-3{cos}^2\delta ][i/\left(k_{0}r\right)-1/{\left(k_{0}r\right)}^2]\right\}.$$

(21)

Notice that if the dipoles are randomly oriented in the 3D space, ${\langle {cos}^2\delta \rangle }_{\mathsf{\Omega }}=1/3$ and $G_{\alpha ,\alpha }\left(\mathbf{r}\right)=e^{i{k}_{0}r}/i{k}_{0}r$, i.e. the scalar model. The decay rate for the vectorial model is

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}^{\left(\delta \right)}\left(r_{jm}\right)& =& \frac{3\mathsf{\Gamma }}{2}\left\{{sin}^2\delta j_0\left(k_0{r}_{jm}\right)+(3{cos}^2\delta -1)\frac{j_1\left(k_0{r}_{jm}\right)}{k_0{r}_{jm}}\right\}\hfill \end{array}$$

(22)

where $j_0\left(x\right)=sinx/x$ and $j_1\left(x\right)=sinx/x^2-cosx/x$ are the spherical Bessel functions of order $n=0$ and $n=1$. As before, it is possible to write ${\mathsf{\Gamma }}^{\left(\delta \right)}\left(r_{jm}\right)$ as angular average of the radiation field between the two atoms. By using the identities

$$\begin{array}{ccc}\hfill {\langle e^{ixcos\theta }\rangle }_{\mathsf{\Omega }}& =& \frac{1}{4\pi }{\int }_0^{2\pi }d\varphi {\int }_0^{\pi }sin\theta e^{ixcos\theta }d\theta =j_0\left(x\right)\hfill \end{array}$$

(23)

$$\begin{array}{ccc}\hfill {\langle {cos}^2\theta e^{ixcos\theta }\rangle }_{\mathsf{\Omega }}& =& \frac{1}{4\pi }{\int }_0^{2\pi }d\varphi {\int }_0^{\pi }sin\theta {cos}^2\theta e^{ixcos\theta }d\theta =j_0\left(x\right)-2\frac{j_1\left(x\right)}{x}.\hfill \end{array}$$

(24)

we can write

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}^{\left(\delta \right)}\left(r_{jm}\right)& =& \frac{3\mathsf{\Gamma }}{4}\left\{\frac{1+{cos}^2\delta }{2}{\langle E_j{E}_m^{*}\rangle }_{\mathsf{\Omega }}+\frac{1-3{cos}^2\delta }{2}{\langle {cos}^2\theta E_j{E}_m^{*}\rangle }_{\mathsf{\Omega }}+\mathrm{c}.\mathrm{c}.\right\}\hfill \end{array}$$

(25)

where $E_j$ has been defined in the previous section. The Fourier transform of ${\mathsf{\Gamma }}^{\left(\delta \right)}\left(r_{jm}\right)$ is

$${\mathsf{\Gamma }}_k^{\left(\delta \right)}=\frac{1}{N}\sum _{j=1}^N\sum _{m=1}^N{\mathsf{\Gamma }}^{\left(\delta \right)}\left(r_{jm}\right)e^{ikd(j-m)}$$

(26)

with $k\in [0,2\pi /d]$. Then,

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k^{\left(\delta \right)}& =& \frac{3\mathsf{\Gamma }}{4N}\left[(1+{cos}^2\delta )\langle |F_k\left(\theta \right){|}^2{\rangle }_{\mathsf{\Omega }}+(1-3{cos}^2\delta )\langle {cos}^2\theta |F_k\left(\theta \right){{|}^2\rangle }_{\mathsf{\Omega }}\right]\hfill \end{array}$$

(27)

where $|F_k{\left(\theta \right)|}^2$ is defined in Eq.(9). By evaluating the angular average and using the approximation (11),

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k^{\left(\delta \right)}& =& \frac{3\mathsf{\Gamma }N}{2k_{0}d}\sum _{m=0}^{\infty }{\int }_{(k-k_0)d/2}^{(k+k_0)d/2}\left\{{sin}^2\delta +\frac{1}{2}(1-3{cos}^2\delta )\frac{{(kd-2t)}^2-{\left(k_{0}d\right)}^2}{{\left(k_{0}d\right)}^2}\right\}\hfill \\ & \times & {\mathrm{sinc}}^2\left[\left(t-m\pi \right)N\right]dt.\hfill \end{array}$$

(28)

In the limit $N\to \infty$, using (13),

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k^{\left(\delta \right)}& =& \frac{3\mathsf{\Gamma }\pi }{2k_{0}d}\sum _{m=0}^{\infty }\left\{{sin}^2\delta +\frac{1}{2}(1-3{cos}^2\delta )\frac{{(kd-2\pi m)}^2-{\left(k_{0}d\right)}^2}{{\left(k_{0}d\right)}^2}\right\}\hfill \\ & \times & \mathsf{\Pi }[2m\pi -k_{0}d<kd<2m\pi +k_{0}d]\hfill \end{array}$$

(29)

We observe that we still have subradiance, i.e. ${\mathsf{\Gamma }}_k^{\left(\delta \right)}=0$ for $k_{0}d<\pi$ and in the interval $k_0<k<2\pi /d-k_0$. However, the expression of ${\mathsf{\Gamma }}_k^{\left(\delta \right)}$ changes with respect to the scalar model solution, i.e. it is no more constant in k. In particular, for $k_{0}d<\pi$ and in the intervals $0<k<k_0$ and $2\pi /d-k_0<k<2\pi /d$,

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k^{\left(\delta \right)}& =& \frac{3\mathsf{\Gamma }\pi }{2k_{0}d}\left\{{sin}^2\delta +\frac{1}{2}(1-3{cos}^2\delta )\left(\frac{k^2}{k_0^2}-1\right)\right\}\hfill \end{array}$$

(30)

At the edge $k=k_0$ we have ${\mathsf{\Gamma }}_{k_0}^{\left(\delta \right)}=(3\mathsf{\Gamma }\pi /2k_{0}d){sin}^2\delta$. Eq.(30) is in agreement with the results of ref.[11]. In Figure 4 we compare the results of the vectorial and the scalar models for an infinite chain, for the same values of $k_{0}d$ as in Figure 1 and for $\delta =\pi /2$. We observe that increasing $k_{0}d$ the differences between the scalar and the vectorial models become less pronounced, and the spectrum becomes more flat. Also the behavior of ${\mathsf{\Gamma }}_k^{\left(\delta \right)}$ vs $k_{0}d$ for $k=0$ and $k=\pi /d$, shown in Figure 5, is similar to that obtained from the scalar model (see Figure 2).

As done for the scalar model, we can obtain an approximated expression for ${\mathsf{\Gamma }}_k^{\left(\delta \right)}$ for finite lattice, substituting the function ${\mathrm{sinc}}^2\left(x\right)$ with the Lorentzian function $1/(1+x^2)$ in Eq.(28) and solving the integral. The result is the following:

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k^{\left(\delta \right)}& =& \frac{3\mathsf{\Gamma }}{2k_{0}d}\sum _{m=0}^{\infty }\left\{\left[{sin}^2\delta +\frac{1}{2}(1-3{cos}^2\delta )\left(\frac{{(kd-2\pi m)}^2-{\left(k_{0}d\right)}^2}{{\left(k_{0}d\right)}^2}-\frac{4}{N^2{\left(k_{0}d\right)}^2}\right)\right]\right.\hfill \\ & \times & \left.(arctan{b}_m-arctan{a}_m)+\frac{2}{N{k}_{0}d}(1-3{cos}^2\delta )\left[1-\frac{kd-2\pi m}{k_{0}d}ln\sqrt{\frac{1+b_m^2}{1+a_m^2}}\right]\right\}\hfill \\ \hfill \end{array}$$

(31)

In particular, the subradiant decay, for $k_0<k<2\pi /d-k_0$ and $k_{0}d<\pi$, is, for $N\gg 1$,

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_k^{\left(\delta \right)}& =& \frac{3\mathsf{\Gamma }}{N}\left\{{sin}^2\delta \left[\frac{2}{{\left(kd\right)}^2-{\left(k_{0}d\right)}^2}+\frac{2}{{(kd-2\pi )}^2-{\left(k_{0}d\right)}^2}\right]+(1-3{cos}^2\delta )\right.\hfill \\ & \times & \left.\left[\frac{4}{{\left(k_{0}d\right)}^2}-\frac{1}{{\left(k_{0}d\right)}^3}\left(kdln\frac{kd+k_{0}d}{kd-k_{0}d}+(kd-2\pi )ln\frac{kd+k_{0}d-2\pi }{kd-k_{0}d-2\pi }\right)\right]\right\}\hfill \end{array}$$

(32)

At the center of the band, at $k=\pi /d$,

$$\begin{array}{ccc}\hfill {\mathsf{\Gamma }}_{k=\pi /d}^{\left(\delta \right)}& =& \frac{3\mathsf{\Gamma }}{N}\left[\frac{4{sin}^2\delta }{{\pi }^2-{\left(k_{0}d\right)}^2}+\frac{1-3{cos}^2\delta }{{\left(k_{0}d\right)}^2}\left(4+\frac{2\pi }{k_{0}d}log\frac{\pi -k_{0}d}{\pi +k_{0}d}\right)\right]\hfill \end{array}$$

(33)

For $k_{0}d\ll \pi$ the subradiant decay is ${\mathsf{\Gamma }}_{k=\pi /d}^{\left(\delta \right)}\approx (8\mathsf{\Gamma }/{\pi }^{2}N)(3{sin}^2\delta -1)$.

4. Conclusions

In this work, we have obtained analytic expressions for the cooperative spontaneous decay rate from a infinite and finite one-dimensional chain of N atoms in the single-excitation configuration. We approximate the true spectrum of the cooperative decay matrix by a function depending on a continuous index k which mimics the quasi-momentum in the first Brillouin zone of a periodic system. The calculations have been done both for the scalar and vectorial light models. The results for infinite lattice are in agreement with those of ref.[11]. However, our approach allows us to obtain an approximated expression for a finite lattice and for arbitrary lattice spacing. In particular, the results show that the ’average’ subradiant decay scales as $1/N$, up to now never demonstrated to our knowledge. This is also in odd with the result of ref.[11], where a dependence as $1/N^3$ is found for the most subradiant eigenvalues. However, this scaling refer to the discrete eigenvalues, whose behavior depends on the microscopic details, such as the polarization of the atoms. On the contrary, our results describe the ’continuous limit’ of the spectrum, neglecting the discreteness of the system. It will be interesting to extend the present ’continuous-limit’ approach of the spectrum from one-dimensional chains to two- and three- dimensional finite arrays, a field not yet completely understood up to now.