Nonreciprocal total cross section of quantum metasurfaces

Nikita Nefedkin; Michele Cotrufo; Andrea Alù

doi:10.1515/nanoph-2022-0596

Article Open Access

Nonreciprocal total cross section of quantum metasurfaces

Nikita Nefedkin , Michele Cotrufo and Andrea Alù

Published/Copyright: January 9, 2023

Published by

Become an author with De Gruyter Brill

Submit Manuscript Author Information Explore this Subject

From the journal Nanophotonics Volume 12 Issue 3

Abstract

Nonreciprocity originating from classical interactions among nonlinear scatterers has been attracting increasing attention in the quantum community, offering a promising tool to control excitation transfer for quantum information processing and quantum computing. In this work, we explore the possibility of realizing largely nonreciprocal total cross sections for a pair of quantum metasurfaces formed by two parallel periodic arrays of two-level atoms. We show that large nonreciprocal responses can be obtained in such nonlinear systems by controlling the position of the atoms and their transition frequencies, without requiring that the environment in which the atoms are placed is nonreciprocal. We demonstrate the connection of this effect with the asymmetric population of a slowly decaying dark state, which is critical to obtain large nonreciprocal responses.

Keywords: atomic arrays; nonreciprocity; quantum metasurface; subradiance

1 Introduction

Combining the classical concepts of metasurfaces and metamaterials with quantum effects has attracted increasing attention over the past few years [1–5]. Recently, ensembles of cold atoms in free space have emerged as a novel platform to control light–matter interactions at the few-photon level. In these systems, often termed quantum metamaterials, exotic phenomena emerge due to interplay between the field radiated and scattered by multiple atoms [1, 6], [7], [8], [9], [10]. Moreover, these systems offer exciting opportunities for applications in quantum information processing and metrology in free space [6, 11], [12], [13].

Among the many desirable features for quantum and classical computation, there is a strong need for devices that support nonreciprocal wave propagation, i.e., systems in which light is transmitted along opposite directions with different efficiencies. Several approaches for nonreciprocity have been recently suggested in both the classical [14–16] and quantum realms [17–19]. In particular, in recent works it has been shown that nonreciprocal propagation or scattering can be obtained in nonlinear systems if they couple asymmetrically to different input/output channels. In classical electromagnetism, this approach has been demonstrated with resonators (such as optical cavities or radio-frequency circuits) coupled asymmetrically to two ports and loaded with Kerr-like nonlinearities, such that their resonant frequencies depend on the stored energy [20–24]. As a consequence, the transmission can be markedly different along the two directions. In parallel, similar behaviors have been shown in quantum devices by suitably combining quantum nonlinearities with breaking of spatial symmetry [25–32]. Fratini et al. have shown [27] that a system of two atoms, optimally detuned and positioned in a waveguide, can lead to large asymmetries between the signals transmitted along the two directions for a large range of input powers. Importantly, both classical and quantum approaches have the important advantage of not requiring any form of external bias, which makes their practical implementation much easier and compatible, for instance, with conventional platforms for quantum computing and processing. However, so far nonlinearity-induced nonreciprocity in free-space quantum systems has not been investigated to the best of our knowledge.

In this work, we investigate whether it is possible to obtain free-space nonreciprocity in quantum metasurfaces, adding an important tool to this rapidly emerging area of research. We consider a system formed by two parallel 2D arrays, each including a large number of atoms, excited by classical plane waves. We demonstrate that, by tailoring the system geometry, it is possible to leverage the quantum nonlinearity of the atoms and obtain large nonreciprocal scattered fields for power levels at which the atomic saturation becomes non-negligible. In particular, we show that the total cross section of the system can be made strongly dependent on the direction of the impinging wave, a feature that cannot be obtained in reciprocal systems. Importantly, here we obtain nonreciprocity by solely engineering the position and detuning of the atoms, without requiring that the electromagnetic environment in which the atoms are placed (i.e., free space) supports nonreciprocal wave propagation [33].

2 System and model

We consider a system composed of two squared 2D arrays of two-level atoms (Figure 1). The arrays are parallel to each other, orthogonal to the z direction, and separated by a distance L. We assume that each array contains N_⊥ × N_⊥ identical two-level atoms (where N_⊥ is the number of atoms along each dimension). In each array (i = 1, 2) the atoms have the same transition frequency ω_i and they are separated by a lattice pitch a_i. Hence, the total number of atoms in the system is N ≡ 2 × N_⊥ × N_⊥. For simplicity, we assume that all atoms have the same transition dipole moment d = de_x, directed along the x axis. We denote by |0_j⟩ and |1_j⟩ the ground and excited states of the jth atom (j = 1, …, N). In the following, we will find it useful to normalize the decay rates of the atomic collective modes by the decay rate of a single atom in free space (with frequency ω_i and dipole moment magnitude d), which is given by γ0i=ωi3d23πϵ0ℏc3 [34], where c is the speed of light and ϵ₀ is the vacuum permittivity.

Figure 1:

Sketch of two periodic atomic arrays separated by the distance L. Both arrays are parallel to the xy plane, and are composed of N_⊥ × N_⊥ identical atoms with frequency ω_i, spatially separated by a lattice constant a_i.

We are interested in modeling the response of this system to impinging classical monochromatic plane waves of arbitrary intensities. For simplicity, we focus here on the case of waves propagating perpendicular to the plane of the arrays, but the formalism can be easily generalized to tilted excitation directions. The impinging electric field is Ein(r)e−iω0t=E0exexp(ik0r)e−iω0t, where the wave vector k0=±|k0|ẑ has magnitude |k₀| = ω₀/c = 2πλ₀/c and it points towards either the +z or −z direction. For brevity, we define as forward (f) the direction corresponding to propagation along the positive z direction (k0=|k0|ẑ), and as backward (b) the opposite direction (k0=−|k0|ẑ).

A full quantum description of the system response would require to explicitly describe the temporal dynamics of the infinite free-space modes that the atoms can interact with, which would make any analytical or numerical solution very challenging. Following standard approaches [35, 36], it is possible to drastically simplify the system description by applying the Born-Markov approximation. Assuming that the typical relaxation time of the atoms is much slower than the time required by light to propagate between atoms, the photonic degrees of freedom can be traced out from the full Hamiltonian of the system, resulting into an effective master equation which involves only the atomic degrees of freedom. In our scenario, the validity of the Born-Markov approximation requires that 1/γ_max ≫ a_max/c, where γ_max is the largest atomic decay rate and a_max is the largest inter-atomic distance. The effective master equation reads [35, 36]

(1)ρ̂̇=iℏρ̂,ĤS+∑i,jΓij22σ̂jρ̂σ̂i+−ρ̂σ̂i+σ̂j−σ̂i+σ̂jρ̂

where we have assumed that, in the frequency range of interest, thermal excitations in the environment are negligible, and thus dissipation can only occur from the atomic system into the environment. For a general arrangement of atoms in free space excited by a classical EM field, the effective Hamiltonian (in the frame rotating at the frequency ω₀ of the incident field) is [35, 37, 38]

(2)ĤS=ℏ∑k=1NΔkσ̂k+σ̂k−ΩRkσ̂k−ΩR*kσ̂k++ℏ∑i≠jΩijσ̂i+σ̂j,

where the indices k, i, j run over all the atoms. We defined the atom-field detuning, Δ_k ≡ ω_k − ω₀, and the interaction constant between the incident field and the kth atom, ΩRk≡d⋅Ein+(rk)/ℏ. For a fixed dipole moment and atom position, ΩRk is proportional to the impinging field amplitude. In the system considered in Figure 1, due to the planar symmetry of the atomic system and to the plane wave excitation, the coefficients in Eq. (2) simplify. In particular, atoms that belong to the same array will have the same detuning (denoted Δ⁽¹⁾ for the first array and Δ⁽²⁾ for the second array), and the same interaction constants, denoted ΩR(1) and ΩR(2). Moreover, ΩR(1) and ΩR(2) only differ by a phase factor eik0L, and we define |ΩR|≡|ΩR(1)|=|ΩR(2)|. The terms proportional to Γ_ij in Eq. (1) and Ω_ij in Eq. (2) account for the dipole-dipole dissipative and coherent interactions, respectively, and they read

(3)Ωij=−3πcωγ0iγ0jdi*ReG⃡(ri,rj,ω)dj,

(4)Γij=6πcωγ0iγ0jdi*ImG⃡(ri,rj,ω)dj,

where G⃡ is the free-space total Green’s tensor,

(5)G⃡(r,r′,k)=I⃡+1k2∇∇eik|r−r′|4π|r−r′|,

and I⃡ is a unit tensor.

The total electric field, given by the sum of the impinging field and the field scattered by the atoms, can be calculated using the input-output relations [37, 39],

(6)Êtot+(r)=Êin+(r)+Êsc+(r).

The incident field operator is defined as Êin+(r)=Ein+(r)Î, where Î is an identity operator acting in the Hilbert space of the atoms and Ein+(r) is the classical field defined above. The scattered field operator is given by

(7)Êsc+(r)=ω2ϵ0c2∑j=1NG⃡(r,rj,ω)djσ̂j.

Equations (6) and (7) allow to calculate the operators associated to the electric field. Following standard approaches, the expectation value of any operator Â in a given state of the system, described by the density matrix ρ̂, is obtained by ⟨Â⟩=Tr(Âρ̂). Equations (6) and (7) can be used to calculate either the field generated by the ensemble of atoms upon external excitation, or the field generated when the system is prepared in any arbitrary state and there is no external excitation (Êin+(r)=0). If the system is prepared in one of its eigenstates, denoted as |ψ_n⟩, the density operator is ρ̂n=|ψn〉〈ψn|.

3 Spectral analysis and fields generated by atomic collective states

In order to understand the scattering properties of this system, it is useful to first study its eigenstates. We initially focus on the single-excitation subspace of the full Hamiltonian. In order to study also the dissipation in the system, we recast the master equation Eq. (1) into a non-Hermitian Hamiltonian,

(8)Ĥeff=ĤS−iℏ∑j,kΓjk2σ̂j+σ̂k

where the imaginary part accounts for dissipation. By introducing the so-called “quantum jump” terms, σ̂j+σ̂k, this effective non-Hermitian Hamiltonian becomes formally equivalent to the master Eq. (1), as shown in Ref. [40].

3.1 Symmetric infinite arrays

For the general case in which the system does not obey parity symmetry (i.e. ω₁ ≠ ω₂ and/or a₁ ≠ a₂), it is challenging to diagonalize the Hamiltonian in Eq. (8) analytically. In order to gain insights on the system dynamics, we initially assume that the system is symmetric, i.e., ω₁ = ω₂ = ω, and a₁ = a₂ = a. Such symmetric double-array system was also investigated in ref. [11]. Assuming infinitely extended lattices, the eigenvalues and eigenstates of the Hamiltonian in Eq. (8) can be calculated analytically [11].

Following [11], we label each atom by a vector index j = (j_⊥, j_z), where j_⊥ = [j_x, j_y] are integer numbers identifying the in-plane position of each atom, R‖(jx,jy)=a(jxêx+jyêy), whereas j_z = 1, 2 specifies which array each atom belongs to (see Figure 1). The general single-excitation eigenstate of the Hamiltonian in Eq. (8), Ĥeff|ψn〉=En|ψn〉, can be written as

(9)|ψn〉=∑jψn,j|ej〉

where the sum runs over all the atoms in the system, |e_j⟩ is the excited state of the atom identified by j, and ψ_n,j is a complex coefficient denoting the contribution of the j-th atom to the nth eigenstate of the system. Due to parity symmetry and translational invariance, the eigenstates of Ĥeff have a plane-wave nature [11, 37]. In particular, each eigenstate can be associated to an in-plane vector q_n defined within the first Brillouin zone of the lattice, and to a parity index p_n = ±1. Specifically,

(10)ψn,j=eiaj⊥⋅qneiπpn−12(jz−1)/2N

(11)En=−i3πγ02(ka)2∑Qk2−|(qn−Q)⋅d|2kk2−|qn−Q|2×1+pneik2−|qn−Q|2L.

The sums in Eq. (11) run over all the vectors Q = (2π/a)(m_x + m_y) of the reciprocal lattice, where m_x and m_y are integers, and k ≡ ω/c. The decay rate of each eigenstate can be calculated from the imaginary part of the expression in Eq. (11), γn≡2Im[En]. It is easy to verify that γ_n is different from zero only when |q_n − Q| < k. In other words, only a finite number of terms in the infinite sum contributes to the decay rate γ_n. The decay rate of the nth eigenstate is therefore given by

(12)γn=3πγ0(ka)2∑Q:|qn−Q|<kk2−|(qn−Q)⋅d|2kk2−|qn−Q|2×1+pn⁡cosk2−|qn−Q|2L.

Equation (12) can be further simplified when the lattice constant is smaller than the atomic wavelength, a < λ, and for eigenstates with zero momentum q_n = 0. In this case the two eigenstates, corresponding to the two parity values p = ±1, have decay rates γ±1=3πγ0(ka)21±cos(kL). In the specific case where kL = mπ with m an integer, we obtain γ₋₁ = 0 and γ+1=6πγ0(ka)2=32πλa2γ0. Here, the system has a dark eigenstate, characterized by zero decay rate, and a bright eigenstate, with a decay rate which is finite and larger than γ₀ (because of the a < λ assumption). In the following, we will use the subscripts D and B to denote quantities pertaining to the dark and bright states, respectively. The dark and bright states are the result of collective effects in the atomic array, leading to constructive and destructive interference of the collective emission. As we will discuss later, the existence of these states is the key property paving the way to nonreciprocal effects in the scattered fields.

3.2 Symmetric finite arrays

In the previous section, we have considered quantum metasurfaces formed by infinitely extended periodic arrays, which allowed us to calculate analytically the eigenstates of the atomic ensemble and to show the existence of dark and bright states. In the remainder of the paper, we will focus on systems of finite size. This is due to two reasons: first, we want to be sure that the predicted effects are observable in experimentally feasible atomic systems, which are currently limited to squared arrays with N_⊥ × N_⊥ ≤ 14 × 14 [41]. Second, considering finite and small arrays will facilitate the analysis of their nonlinear dynamics. When multiple collective excitations are created in the atomic system, analytical results are very challenging to obtain and we will therefore resort to numerical calculations, which are only feasible for small systems.

As discussed in [11], in finite symmetric arrays the collective decay rates are generally different from the ones calculated in infinite arrays Eq. (12). One can still identify bright and dark states by looking at eigenstates whose decay rate are much larger or smaller, respectively, than the decay rate of a single atom in free space. However, the decay rate of dark states in finite systems remains always larger than zero. In Ref. [11], the authors addressed this issue and minimized the decay rate of a dark state by assuming that the arrays are not planar but instead possess a Gaussian-like curvature, calculated based on the phase profile of a Gaussian mode. This allowed them to find analytical expressions for decay rates and frequency shifts of the eigenstates, and to express the eigenmodes as Hermit–Gaussian modes with |q_n| close to 0.

In order to keep our system more realistic, we restrict ourselves to planar arrays, considering collective modes with arbitrarily large |q_n|, and we calculate numerically the single-excitation spectrum of the effective Hamiltonian in Eq. (8). As an example, in Figure 2 we consider a system of two identical arrays, each of them with N_⊥ × N_⊥ = 10 × 10 atoms and we plot the amplitudes of the dark state (defined here as the eigenstate with the smallest decay rate) and the bright state (defined here as the eigenstate with the largest decay rate) at each atomic location in the first array.

Figure 2:

Dark (left panels) and bright (right panel) eigenstates of a two-array system with N_⊥ × N_⊥ = 10 × 10 atoms in each array obeying mirror symmetry (a₁ = a₂ = λ₀/4, ω₁ = ω₂) for L = 0.6λ₀. Panel (a) and (c) show the amplitude |ψ_j| and phase arg(ψ_j) of the dark eigenstate at each atom in the first array. Panels (b) and (d) show the same quantities for the bright state. The amplitudes are the same in the second array, while all phases are shifted by π due to the parity parameter p_n in Eq. (10) (p_B,D = −1 in this case).

Assuming that the system is prepared in one of these eigenstates (and in absence of any impinging field), we can calculate the electric field intensity generated by the quantum metasurfaces. By using Eq. (7), the field intensity generated at any point is I(r)∝∑j,j′dj*G⃡*(r−rj,ω)G⃡r−rj′,ωdj′⟨σ̂j+σ̂j′⟩. The field distribution in the xy plane is expected to follow the pattern shown in Figure 2, whereas the distribution in the xz plane gives us new information about the character of emission of the system in the dark and bright states.

Figure 3 shows the field intensity generated by the bright and dark states in the plane of the array (xy) and in a plane perpendicular to it (xz). The in-plane field profiles (Figure 3a–b) bears similarity to a Hermit–Gaussian mode [11]. As expected, the dark mode (i.e., an eigenstate with almost-zero decay rate) is associated to an almost-zero field far from the arrays (Figure 3c). On the other hand, the bright mode generates a strong field along the two directions perpendicular to the arrays (Figure 3d).

Figure 3:

Field intensity in the dark and bright state: (a), (b) slices in xy plane; (c), (d) slices in xz plane.

While the decay rate of dark states in a finite atomic system will always be non-zero, it can be made arbitrarily small by increasing the number of atoms. In particular, as discussed in [37] the decay rate of the dark state is expected to scale as γ_D ∼ N⁻³ for collective modes at the edge of the Brillouin zone, where N is the total number of atoms in the system. In Figure 4, we show the dependence of γ_D on N, while all other parameters are left constant (see caption to Figure 2). We note that the decay rates of bright states fulfill the condition of applicability of the Born-Markov approximation discussed above, i.e., 1/γ_D ≫ a_max/c, for all the values of N considered here.

Figure 4:

Dependence of the decay rate of the dark state on the number of atoms in the system. The black line is the fit γ_D ∼ N^−α, α = 2.5.

For the largest arrays considered here (N_⊥ × N_⊥ = 20 × 20, N = 800), γ_D is about 10⁻⁷ smaller than the decay rate of a single atom in free space (γ₀). The black curve in Figure 4 shows the fit of the decay rates of the dark state, γ_D ∼ N^−α, with α ≈ 2.5. This value is smaller than the one found in the case of one array [37], which we tentatively attribute to the interactions between the arrays at the parameters of choice. In a different geometry, L = λ₀ and a = λ₀/5, we found that α = 3. As suggested in Refs. [11, 37], such extremely subradiant atomic systems can be used as a quantum memory in network and information processing applications, and protocols for writing and reading states in these systems have been proposed [6, 11].

4 Nonreciprocal total cross section in asymmetric systems

4.1 Arrays with few atoms – quantum analysis

After having elucidated the spectral properties of symmetric systems, we now investigate under what conditions a two-array system can lead to nonreciprocity. As discussed in the introduction, we expect, based on general considerations [23], that a system can become electromagnetically nonreciprocal when geometric asymmetry is combined with a nonlinear optical response. When restricted to the single-excitation subspace, the dynamics of the two-array system is fully linear, and it is analogous to a system of classical dipoles. The two-level atoms, however, feature a saturable absorption. A nonlinear response is therefore obtained at sufficiently high excitation powers. The required geometric asymmetry can be obtained by relaxing the assumption that a₁ = a₂ and/or ω₁ = ω₂. Note that, when considering impinging fields with arbitrarily large power, we cannot use low-power approximations [37, 38] which neglect nonlinear terms.

In order to quantify the degree of nonreciprocity for a given geometry and input power, we calculate the total extinction cross section of the system (for both excitation directions), which can be readily obtained via the optical theorem. Following the formalism described above, we assume that the system is excited by a plane wave with spatial field distribution E_in(r) = E₀e₀exp(ik₀r), where the impinging propagation direction is dictated by the wave vector k0=±|k0|ẑ. At large distances (r ≫ λ₀), the field scattered along the propagation direction, Esc(r=±rẑ), can be written as an asymptotic spherical wave,

(13)Esc(r=±rẑ)→eik0rrE0f(k=k0),

where f(k) is the scattering amplitude along the direction defined by k. By applying the optical theorem [42], the total extinction cross section of the system can be then calculated as

(14)σtot(k0)=4πk0Im(e0*⋅f(k=k0)).

Note that in the following text we use the term ‘‘total cross section’’ to refer to the total extinction cross section. In our model, the scattering amplitude can be found by solving numerically the master Eq. (1) and calculating the expectation values of the dipole operators in steady-state. After that, the expectation value of the scattered field Esc+ at a point very far from the quantum metasurfaces is found using Eq. (7), and then the scattering amplitude is extracted using Eq. (13).

We are interested in quantifying the nonreciprocal extinction of the system, i.e., how differently the system scatters and absorbs waves propagating along opposite directions. Following the forward (f)/backward (b) notation defined above, we define the forward total cross section, σtotf≡σtot(|k0|ẑ) and the backward total cross section, σtotb≡σtot(−|k0|ẑ). This allows us to the quantify the degree of nonreciprocity by the nonreciprocal extinction efficiency, defined as

(15)M≡maxσtotf,σtotbσtotf−σtotbσtotf+σtotb.

For a reciprocal system, M=0 always. In order to find an optimal geometry, we employed an optimization algorithm based on gradient descent to find the optimal system parameters which maximize the function M2 (we used the square of M in order to get a smooth objective function and its derivative). In the optimization, we fix N_⊥ to a given value, and we optimize the following parameters: the lattice constant of the array 1, a₁, the difference between the array lattice constants, δ ≡ a₂ − a₁, the distance between arrays L, the frequency detuning relative to the frequency of the incident field Δ = 2(ω₀ − ω₁) = 2(ω₂ − ω₀) = ω₂ − ω₁ (where we assume that the incident field frequency always lies in the middle between ω₁ and ω₂) and the amplitude of the incident field E₀. For numerical convenience, in the following plots, we normalize the total cross sections and the function M by the total cross section of a single atom with transition frequency ω₀, spontaneous emission rate γ₀ and dipole moment d, placed in free space, which is denoted σtot0 (see Appendix A for the derivation).

We start by optimizing the geometry and input power to achieve large nonreciprocal extinction in a simple system composed of two dimers. That is, each array is composed of only two atoms separated along the x axis by a distance a_i, and the two dimers are separated by a distance L along the z axis, see Figure 5c. We found that a maximum of M2 occurs for Δ = −γ₀, δ = a₂ − a₁ = λ₀/10, a₁ = λ₀/3, L = λ₀/10, |Ω_R| = γ₀/2. Figure 5a shows, for a system with these optimized parameters, the total cross sections for both excitation directions versus the impinging power. A clear nonreciprocal region, whereby the cross sections σtotf and σtotb are markedly different, is visible. As expected, at low impinging powers the total cross sections σtotf and σtotb are the same. This is due to the fact that, for low powers, the population inversion of the atoms is negligible ⟨σ̂jz⟩≈−1. In this scenario, the atoms behave as a system of classical linear dipoles, which is bound to be reciprocal. Moreover, the total cross sections are identical also at high powers. In this case the atoms of both arrays are completely saturated ⟨σ̂jz⟩≈0, and the whole system becomes fully transmissive along both directions. The black dashed curve in Figure 5a shows that the maximum of the function M is obtained for impinging powers such that |Ω_R| ≈ γ₀/2.

$Figure 5: Two-dimer system. (a) Total cross section, (b) population of the dark state and (d) von Neumann entropy of the two-dimer system (with total number of atoms N = 4) when exciting along the forward direction (red line) and the backward direction (blue line), versus the amplitude of the incident field. The two-dimer system is sketched in panel (c). The black dashed line in panel (a) shows the corresponding value of the nonreciprocal extinction efficiency M$\mathcal{M}$. The system parameters are Δ = −γ0, δ = a2 − a1 = λ0/10, a1 = λ0/3, L = λ0/10.$

Figure 5:

Two-dimer system. (a) Total cross section, (b) population of the dark state and (d) von Neumann entropy of the two-dimer system (with total number of atoms N = 4) when exciting along the forward direction (red line) and the backward direction (blue line), versus the amplitude of the incident field. The two-dimer system is sketched in panel (c). The black dashed line in panel (a) shows the corresponding value of the nonreciprocal extinction efficiency M. The system parameters are Δ = −γ₀, δ = a₂ − a₁ = λ₀/10, a₁ = λ₀/3, L = λ₀/10.

To gain more insight on the system dynamics, we look at the steady-state density matrices of the dimer system for the two opposite propagation directions, denoted ρ̂f and ρ̂b, and for a value of impinging powers corresponding to the largest nonreciprocal extinction. The elements of these matrices are shown in Figure 6a and b, with respect to the basis of the uncoupled atomic states. Specifically, the vector |s_i, s_ii, s_iii, s_iv⟩ denotes a state where each of the 4 atoms of the two-dimer system is either in the ground (s_j = 0) or the excited (s_j = 1) state. Figure 6a shows that, for forward excitation, the density matrix has a non-trivial structure and it is a superposition of several states of the uncoupled atomic basis. The density matrix for the backward excitation direction (Figure 6b), and for the same impinging power, has, in contrast, a simple structure in which the ground state |0, 0, 0, 0⟩ is the most populated one. This gives a first hint to explain the nonreciprocal total cross sections observed in Figure 5a. Similar to situations involving few atoms in waveguides [29, 32, 43], a larger ground state population leads to stronger reflection when the atomic system is excited on resonance. A larger reflection implies a larger total cross section, thus explaining the asymmetry observed in Figure 5a.

Figure 6:

The elements of the steady-state density matrices when the two-dimer system is excited along the forward (a) and backward (b) direction. In panels (a) and (b), the density matrix is plotted in the basis of the uncoupled atomic states (as described in the text). In panel (c), we plot the same density matrix as in panel (a), but with respect to the basis of the system eigenstates, |ψ_n⟩, enumerated from 1 to 16. The dark state corresponds to the eigenstate |ψ₅⟩. The system parameters are Δ = −γ₀, δ = a₂ − a₁ = λ₀/10, a₁ = λ₀/3, L = λ₀/10, |Ω_R| = γ₀/2.

Further insights can be gained by plotting the density matrix for forward excitation, ρ̂f, with respect to the basis of the eigenstates of the system, denoted |ψ_j⟩ (j = 1, …, 16), which are obtained by diagonalizing the corresponding Hamiltonian. The result, displayed in Figure 6c, clearly shows that for forward excitation the system populates almost exclusively the eigenstate |ψ₅⟩. By calculating the corresponding eigenvalues of (8) (not shown here), we verified that |ψ₅⟩ is a dark state, with decay rate γ_D = 0.6γ₀. Thus, nonreciprocity is intimately linked to an asymmetric state population. In a certain power range, the steady state of the system evolves towards a dark state when excited along the forward direction, while it evolves into the ground state when excited from the opposite direction. This mechanism is analogous to what described in Refs. [19, 29, 32] for the case of two atoms in a waveguide, whereby different population transfer paths are realized due to interference.

Such asymmetric state population is more clearly displayed in Figure 5b, which shows the population of the dark state, c_D, versus impinging power for the forward and backward directions of excitation. For forward excitation direction the population of the dark state reaches 0.6 within the range of powers where nonreciprocity is large, while it drops to low values at both low and high impinging powers. In the same range of powers, the dark state population is below 0.1 for backward excitation, and it remains low throughout the range of impinging powers.

The dynamics of the system can be further analyzed with the help of the von Neumann entropy which, for any arbitrary state ρ̂, is defined as

(16)SvN=−Tr(ρ̂ln(ρ̂)).

The entropy S_vN quantifies the degree of mixing of the state ρ̂. In particular, it is equal to zero when the system is in a pure state, and it increases as the state ρ̂ becomes a more and more complicated mixture of pure states. In our system the maximum value of entropy is S_vN = ln(2^N), where N is the total number of atoms, corresponding to a maximally mixed state. The von Neumann entropy versus impinging power for both excitation directions is shown in Figure 5d. For low powers, S_vN ≈ 0 for both directions, corresponding to the fact that system remains in the (pure) ground state. As the impinging power increases, the entropy increases as well, but in a different way for the two impinging directions. For forward excitation (red line in Figure 5d) the dark state gets partially populated at intermediate powers (c_D = 0.6). However, the steady-state of the system in this scenario is not pure, since additional states are populated. This increases the entropy. When exciting along the backward direction (blue line in Figure 5d), at the same intermediate power level, the system is instead almost entirely localized in the ground state, which results in a lower entropy. At large impinging powers, the entropy tends to its maximum value, S_vN = ln(2^N) ≈ 2.77 for both directions of excitation, highlighting that the system is in a maximally mixed state, see also Figure 15 in Appendix B.

The system considered in Figures 5 and 6 was obtained by optimizing the geometry and parameters of the atomic ensemble for a given input frequency ω₀. As discussed above, here nonreciprocity is due to the asymmetrical excitation of a dark state. Thus, the nonreciprocal behaviour is expected to depend strongly on the excitation frequency (for a given set of system parameters). This is confirmed by the calculations shown in Appendix B, where we calculate the nonreciprocal efficiency M for different values of the input frequency ω₀.

Having verified the emergence of nonreciprocal extinction in a small system composed of two dimers, we now consider a slightly bigger system, composed of two squared arrays with N_⊥ × N_⊥ = 2 × 2 atoms each. For this system we are unable to find the global maximum of the M2 function in the entire parameter space, due to the higher numerical complexity. However, we found that the same optimal parameters which optimize the two-dimer system (see the caption of Figure 5) lead to large nonreciprocal effects also in this larger system. In Figure 7 we show the total field intensity, the real part of scattered field (x-component) and the real part of total field (x-component) when this larger system is exited along the forward or backward direction, for the power level corresponding to |Ω_R| = γ₀/2 (a similar plot for the two-dimer system can be found in the Appendix B, Figure 14). A clear nonreciprocal behaviour can be seen in Figure 7: for a wave incident along the forward direction the magnitude of the scattered field is much smaller than the case in which the same wave propagates along the backward direction.

$Figure 7: Total field intensity, real part of Escx${E}_{\text{sc}}^{x}$ and real part of Etotx${E}_{\text{tot}}^{x}$ of 2-array system, with each array containing N⊥ × N⊥ = 2 × 2 atoms, for forward (panel a) and backward (panel b) excitation direction parameters: Δ = −γ0, δ = a2 − a1 = λ0/10, a1 = λ0/3, L = λ0/10, |Ω R | = γ0/2.$

Figure 7:

Total field intensity, real part of Escx and real part of Etotx of 2-array system, with each array containing N_⊥ × N_⊥ = 2 × 2 atoms, for forward (panel a) and backward (panel b) excitation direction parameters: Δ = −γ₀, δ = a₂ − a₁ = λ₀/10, a₁ = λ₀/3, L = λ₀/10, |Ω_R| = γ₀/2.

We emphasize that, in all the systems considered in this work, the difference between the total cross sections (i.e. the fact that σtotf>σtotb within a certain power range) is solely dictated by the geometry of the arrays and by the distance between them, which in turns introduces a phase shift in the transmitted field exp(ik₀L). In particular, the direction along which the total cross section is the largest can be switched by simply changing the value of L. For example, for the two-dimer system the total cross section along the backward excitation becomes larger than the one along the forward direction when, for example, L > λ₀/2.

In conclusion, nonreciprocity arises in the system due to the realization of different routes of population transfer depending on the excitation direction. When the system is excited along the forward direction the population is partially trapped in the dark state, whereas for the opposite direction the system remains in the ground state. This asymmetric population in the system’s steady state (Figure 6) leads to nonreciprocal extinction, which manifests itself in different total cross sections, Figure 5.

4.2 Large arrays – semiclassical approach

In the previous section we have focused on systems composed of few atoms. Keeping the number of atoms small (N ≤ 10) allowed us to numerically solve the full master Eq. (1) without any approximation. However, as the number of atoms increases, the exponential increase of the dimension of the Hilbert space prevents us from using the full quantum formalism to calculate the dynamics of larger arrays. This is particularly challenging for the scenario considered in this work: as we are considering large excitation powers, we cannot truncate the Hilbert space to the single-excitation subspace. A common way to address this challenge, and to be able to describe systems composed of hundreds of atoms, is to employ the so-called semiclassical approximation. The core idea of this approximation is to neglect quantum correlations between atoms. In particular, one begins by writing down the Langevin–Heisenberg equations for the operators σ̂j, σ̂j+, σ̂jz, and then take the average values of both sides of the each equation. Solving the time-dependent differential equations for ⟨σ̂j⟩, ⟨σ̂j+⟩, ⟨σ̂jz⟩, however, requires knowing the time-dependent values of two-operator expectation values, such as ⟨σ̂i+σ̂j⟩ and others, which, in turn, depend on three-operator expectation values and so on. This results in a cascaded hierarchy of a large number of differential equations. Within the semiclassical approximation, one assumes that the expectation value of a product of operators is always factorizable, i.e., ⟨ÂB̂⟩≈⟨Â⟩⟨B̂⟩ for any pair of operators Â and B̂ [34, 44]. This approximation allows to describe a system of N atoms by solving a system of 3N coupled differential equations instead of 2^2N. In particular, the equations describing the evolution of the system are [35]:

(17)(17)σ˙k=−iΔk−Γkk2σk+∑j≠kiΩjk+Γjk2σkzσj−idkEin+(rk)ℏσkz

(18)(18)σ˙k+=iΔk−Γkk2σk++∑j≠k−iΩjk+Γjk2σj+σkz+idk*Ein−(rk)ℏσkz

(19)(19)σ˙kz=−Γkkσkz+1−2∑j≠kiΩkj+Γkj2σk+σj +−iΩjk+Γjk2σj+σk+2idkEin+(rk)ℏσk+−dk*Ein−(rk)ℏσk,

where the first and second equations are the complex conjugated of each other. The semiclassical approximation allows us to use the optimization procedure described above (based on maximizing M2, Eq. (15)) for very large arrays. In Figure 8, we show an example of an optimized geometry (parameters in caption) with N_⊥ = 10 (total number of atoms N = 200).

$Figure 8: Total field intensity, real part of Escx${E}_{\text{sc}}^{x}$ and real part of Etotx${E}_{\text{tot}}^{x}$ upon forward (panel a) and backward (panel b) excitation of a 2-array system, each of them composed by N⊥ × N⊥ = 10 × 10 atoms. Parameters: δ = a2 − a1 = 0, a1 = λ0/4, Δ = 2γ0, L = 0.6λ0, |Ω R | = 2γ0.$

Figure 8:

Total field intensity, real part of Escx and real part of Etotx upon forward (panel a) and backward (panel b) excitation of a 2-array system, each of them composed by N_⊥ × N_⊥ = 10 × 10 atoms. Parameters: δ = a₂ − a₁ = 0, a₁ = λ₀/4, Δ = 2γ₀, L = 0.6λ₀, |Ω_R| = 2γ₀.

From Figure 8, we see a similar behavior of the scattered field as for the small arrays (compare with Figure 7) – the significant difference being the scattered field for the forward and backward directions of the excitation. When excited along the forward direction, the scattered field (Figure 8a, center panel) has a larger magnitude than when excited from the opposite direction. The perturbation to the total field is also more noticeable in the case of the forward direction of excitation (Figure 8a, left panel). As mentioned earlier, the direction along which the scattering is larger is set by the phase shift between the arrays and the specific in-plane geometry. For the particular set of optimal parameters considered here, differently from the small arrays considered earlier, for intermediate powers the system is trapped in the dark state when excited along the backward direction.

We can understand better the details of the scattering process by looking at 1D slices of the field distributions shown in Figure 8. Let us consider first the case of forward excitation direction (Figure 9a): the scattered field (red solid lines) has a high amplitude, comparable to the one of the incident field (E₀, black dashed lines). For z < 0 the real parts of the incident and scattered fields are in phase, while the imaginary parts cancel each other. This leads to the emergence of a standing-wave pattern of the total field for z < 0 and to a large reflection level. For z > 0 the incident and scattered fields interfere destructively, and the total field is suppressed for λ₀ < z < 3λ₀. Beyond this region (i.e. for z > 3λ₀), the total field increases and its magnitude becomes comparable to the incident field. This is due to the fact that, since the incident field is assumed to be an infinitely extended plane wave, diffraction of the impinging field from the edges of the array dominates the total field at large distances, where the scattered field is instead almost zero. For the backward excitation direction (Figure 9b) the incident and scattered fields interfere destructively for z < 0. However, due to the large difference between the amplitudes of these two fields, the total field is not fully canceled. For z > 0 the incident and scattered fields have a mutual phase shift of ≈π/2, which does not lead to either destructive or constructive interference.

Figure 9:

Real part, imaginary part and absolute value of the x component of incident, total and scattered fields at x = 0, y = 0, for the same system as in Figure 8.

The nonreciprocal extinction in this large system is also confirmed by the total cross sections. In Figure 10, we plot the total cross sections for forward (red line) and backward (blue line) impinging directions versus the amplitude of the incident field. Similarly to what found in Figure 5a for the case of small arrays, we found that σtotf and σtotb are markedly different for intermediate powers. Interestingly, we note that increasing the size of the system does not increase the nonreciprocal extinction efficiency M. The maximum ratio between the two cross sections is σtotf/σtotb≈1.5, which is slightly smaller than the one obtained in Figure 5 for the two-dimer system. Moreover, we note that the maximum value of M occurs at higher values of input amplitude, Ω_R ≈ 2γ₀.

$Figure 10: Total cross section associated with the forward and backward direction of excitation versus the absolute interaction constant between the incident field and the atoms, |Ω R |/γ0, normalized by the total cross section of a single atom, σtot0${\sigma }_{\text{tot}}^{0}$, N⊥ × N⊥ = 100. Parameters are the same as in Figure 8.$

Figure 10:

Total cross section associated with the forward and backward direction of excitation versus the absolute interaction constant between the incident field and the atoms, |Ω_R|/γ₀, normalized by the total cross section of a single atom, σtot0, N_⊥ × N_⊥ = 100. Parameters are the same as in Figure 8.

The connection between the nonreciprocal behavior and the excitation of a dark state can be further proved by analyzing the distribution of the amplitudes and phases of the dipole moments of each atom (given by arg(σ_j)) in the steady state and at the power level where nonreciprocity is maximized. Figure 11 shows the phases of the dipoles in polar plots, while Figure 12 shows the amplitude and phase of each atom in cartesian color-coded scatter plots.

Figure 11:

Phases of the atomic dipole moments in the steady state at optimal parameters for forward (panel a) and backward (panel b) direction of excitation. Red stars and blue crosses denote the average phases in the array 1 and array 2, respectively, N_⊥ × N_⊥ = 100. Parameters are the same as in Figure 8.

Figure 12:

Absolute values and phases of the atomic dipole moments in the steady state at optimal parameters for forward (panel a) and backward (panel b) direction of excitation. N_⊥ × N_⊥ = 100, other parameters are the same as in Figure 9.

For forward impinging direction (Figure 11a), the phases of the dipoles are distributed over a large interval. On the contrary, for backward impinging direction (Figure 11b), the phases are tightly localized in a small interval around 3π/2. To quantify this phenomena, we define the average phase arg(σ)kd̄=∑jkargσjkd/N⊥2 and corresponding variance Vararg(σ)kd, where j_k denotes the jth atom of the kth array (k = 1, 2), the subscript d = f, b denotes the excitation direction, and the sum runs over all atoms in the array specified by k. For forward excitation direction, the phases of the dipoles in array 1 are grouped around arg(σ)1f̄≈π/10 with a relatively small variance, Vararg(σ)1f=0.025, whereas for array 2 we find arg(σ)2f̄≈4π/5, and the variance Vararg(σ)2f=0.25 is larger by an order of magnitude. The difference between the average dipole phases in the two arrays is arg(σ)1f̄−arg(σ)2f̄|≈7π/10, which is equals ≈17π/10 when taking into account the phase shift between arrays (≈π).

When considering the backward direction, instead, we observe a strong phasing of the dipole moments [45, 46] both within each array and between array 1 and array 2. The average phases are equal for both arrays, arg(σ)1b̄≈arg(σ)2b̄≈−π/2, and the variances are Vararg(σ)1b=2⋅10−3 and Vararg(σ)2b=2⋅10−2 for array 1 and array 2, respectively. Combined with the fact that the arrays are separated by the distance L ≈ λ₀/2, this shows that the dipoles in the two arrays have almost the opposite phase. Opposite phases lead to suppression of the scattered field (see Figure 10b), which is a manifestation of trapping the system’s population in the dark state. Importantly, despite the similar phases and their distributions across the arrays, the amplitudes of the atomic dipole moments are distributed differently on array 1 and array 2, see Figure 12b.

Thus, also within a semiclassical approximation approach, we found a clear connection between nonreciprocity and excitation of a dark state. This trapping manifests itself as phasing of the dipoles in both arrays when the system is excited from one direction. In addition, the dipoles in array 1 and array 2 have opposite phases, which also highlights the connection to a subradiant state.

5 Conclusions

In this work, we demonstrated large nonreciprocal extinction from asymmetric pairs of quantum metasurfaces composed of two-level atoms. Nonreciprocity is obtained by properly combining geometrical asymmetry and nonlinear response, hereby provided by the atoms saturable absorption. For a given number of atoms in the system, we optimized several parameters in the metasurface design to achieve the largest asymmetry between the total cross sections calculated when the system is excited from opposite directions. We have investigated systems composed of few atoms, which can be described by a fully quantum master equation approach, and also systems composed of hundreds of atoms, whereby the semiclassical approximation is required to keep the numerical calculations feasible. In both cases, we showed that the occurrence of nonreciprocal extinction is intimately related to the asymmetric population of a dark state, and it results in nontrivial features in the spatial distribution of the scattered and total fields. Moreover, we show that the nonreciprocal extinction, and the excitation of a dark state, also manifests peculiar phasing phenomena, whereby the phases of atomic dipoles condense near some common average phase within array 1 and array 2, and the phase shift between arrays is ≈π.

Our results show that the dark state plays a key role in the emergence of nonreciprocity in complex free-space atomic systems, where atoms interact with each other through dipole-dipole interactions. We expect these results to stimulate experimental studies in atomic lattices of trapped cold atoms and quantum metasurfaces. Dark-state-induced nonreciprocity in atomic arrays may contribute to the development of quantum technologies that requiring efficient and tunable state transfer and state management at microscopic scale.

Corresponding authors: Nikita Nefedkin, Photonics Initiative, Advanced Science Research Center, City University of New York, New York, NY 10031, USA, E-mail: nnefedkin@gc.cuny.edu; and Andrea Alù, Photonics Initiative, Advanced Science Research Center, City University of New York, New York, NY 10031, USA; and Physics Program, Graduate Center, City University of New York, New York, NY 10016, USA, E-mail: aalu@gc.cuny.edu

Funding source: Simons Foundation

Award Identifier / Grant number: Unassigned

Funding source: Air Force Office of Scientific Research

Award Identifier / Grant number: Unassigned

Author contributions: All the authors have accepted responsibility for the entire content of this submitted manuscript and approved submission.
Research funding: This work was supported by the Simons Foundation and the Air Force Office of Scientific Research.
Conflict of interest statement: The authors declare no conflicts of interest regarding this article.

Appendix A: Derivation of the total cross section for one atom

We consider a system consisting of a two-level atom, with transition frequency ω_a and dipole moment d, located at r₀ and excited by a classical monochromatic EM field with an arbitrary profile and frequency ω₀. The Hamiltonian of this system reads

(A1)Ĥ=ℏωaσ̂+σ̂−dEin+(r0)e−iω0tσ̂−d*Ein−(r0)eiω0tσ̂+.

In order to get rid of the explicit time dependence, we rewrite the Hamiltonian in a frame rotating at the frequency of the incident field ω₀. The can be obtained via a unitary transformation, described by the operator

(A2)Û=exp−iℏω0σ̂+σ̂

(A3)H̄̂=Û+ĤÛ−iÛ+∂Û∂t,

and the Hamiltonian in the rotating frame is

(A4)H̄̂=ℏΔσ̂+σ̂−ℏΩRσ̂+ΩR*σ̂+,

where Δ = ω_a − ω₀, and ΩR=dEin+(r0)/ℏ. The atom interacts with the EM modes of free space and, therefore, spontaneously emits at the rate

(A5)γ0=ω03d23πϵ0ℏc3.

After applying Born and Markov approximations, we obtain the master equation in Lindblad form [34],

ρ̂̇=−iℏH̄̂,ρ̂+γ021+n(ωa,T)(2σ̂ρσ̂+−ρ̂σ̂+σ̂−σ̂+σ̂ρ̂)+γ02n(ωa,T)(2σ̂+ρ̂σ̂−ρ̂σ̂σ̂+−σ̂σ̂+ρ̂),

where n(ωa,T)=1/eℏωa/kBT−1 is the Bose–Einstein distribution. By looking for the steady-state of the master equation, we find the steady-state density matrix,

(A6)ρ11s.s.=4|ΩR|2γ022nωa,T+12+8|ΩR|2+4Δ2+nωa,T2nωa,T+1,

(A7)ρ12s.s.=2iΩR*γ02nωa,T+1+2iΔ2nωa,T+1γ022nωa,T+12+8|ΩR|2+4Δ2,

(A8)ρ21s.s.=−2ΩR2Δ+iγ02nωa,T+12nωa,T+1γ022nωa,T+12+8|ΩR|2+4Δ2,

(A9)ρ22s.s.=1−4|ΩR|2γ022nωa,T+12+8|ΩR|2+4Δ2+nωa,T2nωa,T+1.

Following standard procedures, we can calculate the average value of the atom polarization as ⟨p̂⟩=Tr(ρ̂s.s.p̂) assuming that n(ω_a, T) ≈ 0, resulting into

(A10)〈p̂〉=d〈σ̂〉=d2i(γ0−2iΔ)ΩR*γ02+8|ΩR|2+4Δ2.

We then substitute the average polarization ⟨p̂⟩ into Eq. (7) and using Eq. (14) we find the total cross section of 1 atom, σtota.

Appendix B: Further analysis of the dimer system

B.1 Spectral properties of the system

In the main text, the frequency of the incident field was fixed to a certain value ω₀ and the system geometry and parameters where then optimized to maximize nonreciprocity. In this section we discuss the dependence of the total cross section, the von Neumann entropy of the system and the dark state population on the input frequency ω₀, for the same parameters considered in the main text (see Section 4.1. and the caption of Figure 6). In the main text, the atoms were symmetrically detuned with respect to the incident EM field frequency, i.e., the input frequency was fixed to ω₀ = ω_1,2 ≡ (ω₁+ω₂)/2.

Figure 13 shows the dependence of the total cross sections, dark state populations and von Neumann entropy for the forward and backward directions of excitation on the frequency of the incident field (normalized by the average atomic frequency, ω_1,2).

$Figure 13: Two-dimer system: spectral properties. (a) Total cross section and nonreciprocal extinction efficiency, (b) population of the dark state, and (c) von Neumann entropy of the two-dimer system (with total number of atoms N = 4) when exciting along the forward direction (red line) and the backward direction (blue line), versus the frequency of the incident field. The black dashed line in panel (a) shows the corresponding value of the nonreciprocal extinction efficiency M$\mathcal{M}$. The system parameters are Δ = −γ0, δ = a2−a1 = λ0/10, a1 = λ0/3, L = λ0/10.$

Figure 13:

Two-dimer system: spectral properties. (a) Total cross section and nonreciprocal extinction efficiency, (b) population of the dark state, and (c) von Neumann entropy of the two-dimer system (with total number of atoms N = 4) when exciting along the forward direction (red line) and the backward direction (blue line), versus the frequency of the incident field. The black dashed line in panel (a) shows the corresponding value of the nonreciprocal extinction efficiency M. The system parameters are Δ = −γ₀, δ = a₂−a₁ = λ₀/10, a₁ = λ₀/3, L = λ₀/10.

The total cross sections have a spectral bandwidth determined by the dark state decay rate, γ_D and, therefore, the nonreciprocity effect is observable within a comparable bandwidth Δω ≈ γ_D,. At the same time, due to the optimization procedure described in the main text, the maximum of M is reached when the incident field frequency is approximately equal to the average atomic frequency, ω₀ = ω_1,2. Moreover, both the population of the dark state, c_D, and the von Neumann entropy have the maximum contrast for forward and backward directions of excitation at ω₀ = ω_1,2 as well. The doublet of peaks visible in some of the panels of Figure 13 is due to the detuning between atoms in array 1 and array 2.

B.2 Total and scattered fields

In the main text we studied nonreciprocity emerging in the simple dimer system (whose geometry is sketched in Figure 5c). Using the approach developed in the main text, we found the total field and the scattered field for the dimer system at the optimal parameters which maximize the nonreciprocal extinction efficiency (15), see the caption of Figure 6. In Figure 14, we show additional results for this system, namely the total intensity, scattered field and total field.

$Figure 14: Total field intensity, real part of Escx${E}_{\text{sc}}^{x}$ and real part of Etotx${E}_{\text{tot}}^{x}$ upon forward (panel a) and backward (panel b) excitation of the two-dimer system, each of them composed by 2 atoms. Parameters: Δ = −γ0, δ = a2 − a1 = λ0/10, a1 = λ0/3, L = λ0/10.$

Figure 14:

Total field intensity, real part of Escx and real part of Etotx upon forward (panel a) and backward (panel b) excitation of the two-dimer system, each of them composed by 2 atoms. Parameters: Δ = −γ₀, δ = a₂ − a₁ = λ₀/10, a₁ = λ₀/3, L = λ₀/10.

From Figure 14, one can see that all the characteristic features of the fields behavior matches the case of a larger arrays shown in Figure 8.

B.3 Steady-state density matrix at the high power

Further let us have a look at the steady state of the system at high powers of excitation, i.e., when |Ω_R| ≫ γ₀. The corresponding density matrix is shown in Figure 15.

Figure 15:

The elements of the steady-state density matrices when the two-dimer system is excited along either the forward and backward direction. The density matrix is plotted in the basis of the uncoupled atomic states. The system parameters are Δ = −γ₀, δ = a₂ − a₁ = λ₀/10, a₁ = λ₀/3, L = λ₀/10, |Ω_R| = 20γ₀.

Here the density matrix is close to a diagonal form, i.e. all the nondiagonal elements are negligible. Moreover, the diagonal elements which are responsible for the populations of the correspondent states are all equal. This fact is reflected in the behavior of the von Neumann entropy, which reaches its maximum at the high power of the incident field, see Figure 5d.

References

[1] R. Bekenstein, I. Pikovski, H. Pichler, E. Shahmoon, S. F. Yelin, and M. D. Lukin, “Quantum metasurfaces with atom arrays,” Nat. Phys., vol. 16, p. 676, 2020. https://doi.org/10.1038/s41567-020-0845-5.Search in Google Scholar

[2] A. S. Solntsev, G. S. Agarwal, and Y. S. Kivshar, “Metasurfaces for quantum photonics,” Nat. Photonics, vol. 15, p. 327, 2021. https://doi.org/10.1038/s41566-021-00793-z.Search in Google Scholar

[3] K. Wang, J. G. Titchener, S. S. Kruk, et al.., “Quantum metasurface for multiphoton interference and state reconstruction,” Science, vol. 361, p. 1104, 2018. https://doi.org/10.1126/science.aat8196.Search in Google Scholar PubMed

[4] M. Sakhdari, M. Hajizadegan, M. Farhat, and P. Y. Chen, “Efficient, broadband and wide-angle hot-electron transduction using metal-semiconductor hyperbolic metamaterials,” Nano Energy, vol. 26, p. 371, 2016. https://doi.org/10.1016/j.nanoen.2016.05.037.Search in Google Scholar

[5] P. Y. Chen and M. Farhat, “Modulatable optical radiators and metasurfaces based on quantum nanoantennas,” Phys. Rev. B, vol. 91, p. 035426, 2015. https://doi.org/10.1103/physrevb.91.159904.Search in Google Scholar

[6] M. Reitz, C. Sommer, and C. Genes, “Cooperative quantum phenomena in light-matter platforms,” PRX Quantum, vol. 3, p. 010201, 2022. https://doi.org/10.1103/prxquantum.3.010201.Search in Google Scholar

[7] K. E. Ballantine and J. Ruostekoski, “Quantum single-photon control, storage, and entanglement generation with planar atomic arrays,” PRX Quantum, vol. 2, p. 040362, 2021. https://doi.org/10.1103/prxquantum.2.040362.Search in Google Scholar

[8] S. J. Masson, I. Ferrier-Barbut, L. A. Orozco, A. Browaeys, and A. Asenjo-Garcia, “Many-body signatures of collective decay in atomic chains,” Phys. Rev. Lett., vol. 125, p. 263601, 2020. https://doi.org/10.1103/physrevlett.125.263601.Search in Google Scholar PubMed

[9] K. E. Ballantine and J. Ruostekoski, “Unidirectional absorption, storage, and emission of single photons in a collectively responding bilayer atomic array,” Phys. Rev. Res., vol. 4, p. 033200, 2022. https://doi.org/10.1103/physrevresearch.4.033242.Search in Google Scholar

[10] D. Fernández-Fernández and A. González-Tudela, “Tunable directional emission and collective dissipation with quantum metasurfaces,” Phys. Rev. Lett., vol. 128, p. 113601, 2022. https://doi.org/10.1103/physrevlett.128.113601.Search in Google Scholar

[11] P. O. Guimond, A. Grankin, D. V. Vasilyev, B. Vermersch, and P. Zoller, “Subradiant Bell states in distant atomic arrays,” Phys. Rev. Lett., vol. 122, p. 093601, 2019. https://doi.org/10.1103/physrevlett.122.093601.Search in Google Scholar

[12] D. E. Chang, J. S. Douglas, A. González-Tudela, C. L. Hung, and H. J. Kimble, “Colloquium: quantum matter built from nanoscopic lattices of atoms and photons,” Rev. Mod. Phys., vol. 90, p. 031002, 2018. https://doi.org/10.1103/revmodphys.90.031002.Search in Google Scholar

[13] M. Manzoni, M. Moreno-Cardoner, A. Asenjo-Garcia, J. V. Porto, A. V. Gorshkov, and D. Chang, “Optimization of photon storage fidelity in ordered atomic arrays,” New J. Phys., vol. 20, p. 083048, 2018. https://doi.org/10.1088/1367-2630/aadb74.Search in Google Scholar PubMed PubMed Central

[14] C. Caloz, A. Alu, S. Tretyakov, D. Sounas, K. Achouri, and Z. L. Deck-Léger, “Electromagnetic nonreciprocity,” Phys. Rev. Appl., vol. 10, p. 047001, 2018. https://doi.org/10.1103/physrevapplied.10.047001.Search in Google Scholar

[15] H. Nassar, B. Yousefzadeh, R. Fleury, et al.., “Nonreciprocity in acoustic and elastic materials,” Nat. Rev. Mater., vol. 5, p. 667, 2020. https://doi.org/10.1038/s41578-020-0206-0.Search in Google Scholar

[16] S. Manipatruni, J. T. Robinson, and M. Lipson, “Optical nonreciprocity in optomechanical structures,” Phys. Rev. Lett., vol. 102, p. 213903, 2009. https://doi.org/10.1103/physrevlett.102.213903.Search in Google Scholar PubMed

[17] A. Mahoney, J. Colless, S. Pauka, et al.., “On-chip microwave quantum hall circulator,” Phys. Rev. X, vol. 7, p. 011007, 2017. https://doi.org/10.1103/physrevx.7.011007.Search in Google Scholar

[18] S. Zhang, Y. Hu, G. Lin, et al.., “Thermal-motion-induced non-reciprocal quantum optical system,” Nat. Photonics, vol. 12, p. 744, 2018. https://doi.org/10.1038/s41566-018-0269-2.Search in Google Scholar

[19] A. R. Hamann, C. Müller, M. Jerger, et al.., “Nonreciprocity realized with quantum nonlinearity,” Phys. Rev. Lett., vol. 121, p. 123601, 2018. https://doi.org/10.1103/physrevlett.121.123601.Search in Google Scholar

[20] D. L. Sounas and A. Alù, “Fundamental bounds on the operation of Fano nonlinear isolators,” Phys. Rev. B, vol. 97, p. 115431, 2018. https://doi.org/10.1103/physrevb.97.115431.Search in Google Scholar

[21] D. L. Sounas, J. Soric, and A. Alu, “Broadband passive isolators based on coupled nonlinear resonances,” Nat. Electron., vol. 1, p. 113, 2018. https://doi.org/10.1038/s41928-018-0025-0.Search in Google Scholar

[22] K. Y. Yang, J. Skarda, M. Cotrufo, et al.., “Inverse-designed non-reciprocal pulse router for chip-based LiDAR,” Nat. Photonics, vol. 14, p. 369, 2020. https://doi.org/10.1038/s41566-020-0606-0.Search in Google Scholar

[23] M. Cotrufo, S. A. Mann, H. Moussa, and A. Alù, “Nonlinearity-induced nonreciprocity-Part I,” IEEE Trans. Microwave Theory Tech., vol. 69, p. 3569, 2021. https://doi.org/10.1109/tmtt.2021.3079250.Search in Google Scholar

[24] M. Cotrufo, S. A. Mann, H. Moussa, and A. Alù, “Nonlinearity-induced nonreciprocity-Part II,” IEEE Trans. Microwave Theory Tech., vol. 69, p. 3584, 2021. https://doi.org/10.1109/tmtt.2021.3082192.Search in Google Scholar

[25] D. Roy, “Few-photon optical diode,” Phys. Rev. B, vol. 81, p. 155117, 2010. https://doi.org/10.1103/physrevb.81.155117.Search in Google Scholar

[26] D. Roy, “Cascaded two-photon nonlinearity in a one-dimensional waveguide with multiple two-level emitters,” Sci. Rep., vol. 3, p. 2337, 2013. https://doi.org/10.1038/srep02337.Search in Google Scholar PubMed PubMed Central

[27] F. Fratini, E. Mascarenhas, L. Safari, et al.., “Fabry-perot interferometer with quantum mirrors: nonlinear light transport and rectification,” Phys. Rev. Lett., vol. 113, p. 243601, 2014. https://doi.org/10.1103/physrevlett.113.243601.Search in Google Scholar

[28] F. Fratini and R. Ghobadi, “Full quantum treatment of a light diode,” Phys. Rev. A, vol. 93, p. 023818, 2016. https://doi.org/10.1103/physreva.93.023818.Search in Google Scholar

[29] C. Müller, J. Combes, A. R. Hamann, A. Fedorov, and T. M. Stace, “Nonreciprocal atomic scattering: a saturable, quantum Yagi-Uda antenna,” Phys. Rev. A, vol. 96, p. 053817, 2017. https://doi.org/10.1103/physreva.96.053817.Search in Google Scholar

[30] E. Mascarenhas, M. F. Santos, A. Auffèves, and D. Gerace, “Quantum rectifier in a one-dimensional photonic channel,” Phys. Rev. A, vol. 93, p. 043821, 2016. https://doi.org/10.1103/physreva.93.043821.Search in Google Scholar

[31] Y. L. L. Fang and H. U. Baranger, “Multiple emitters in a waveguide: nonreciprocity and correlated photons at perfect elastic transmission,” Phys. Rev. A, vol. 96, p. 013842, 2017. https://doi.org/10.1103/physreva.96.013842.Search in Google Scholar

[32] N. Nefedkin, M. Cotrufo, A. Krasnok, and A. Alù, “Dark-state induced quantum nonreciprocity,” Adv. Quantum Technol., vol. 5, p. 2100112, 2022. https://doi.org/10.1002/qute.202100112.Search in Google Scholar

[33] S. A. Hassani Gangaraj, L. Ying, F. Monticone, and Z. Yu, “Enhancement of quantum excitation transport by photonic nonreciprocity,” Phys. Rev. A, vol. 106, p. 033501, 2022. https://doi.org/10.1103/physreva.106.033501.Search in Google Scholar

[34] H. J. Carmichael, Statistical Methods in Quantum Optics 1: Master Equations and Fokker-Planck Equations, vol. 1, New York, Springer Science & Business Media, 1999.10.1063/1.883009Search in Google Scholar

[35] R. J. Bettles, M. D. Lee, S. A. Gardiner, and J. Ruostekoski, “Quantum and nonlinear effects in light transmitted through planar atomic arrays,” Commun. Phys., vol. 3, p. 141, 2020. https://doi.org/10.1038/s42005-020-00404-3.Search in Google Scholar

[36] D. S. Wild, “Algorithms and platforms for quantum science and technology,” Ph.D. thesis, Harvard University, 2020.Search in Google Scholar

[37] A. Asenjo-Garcia, M. Moreno-Cardoner, A. Albrecht, H. J. Kimble, and D. E. Chang, “Exponential improvement in photon storage fidelities using subradiance and “selective radiance” in atomic arrays,” Phys. Rev. X, vol. 7, p. 031024, 2017. https://doi.org/10.1103/physrevx.7.031024.Search in Google Scholar

[38] E. Shahmoon, D. S. Wild, M. D. Lukin, and S. F. Yelin, “Cooperative resonances in light scattering from two-dimensional atomic arrays,” Phys. Rev. Lett., vol. 118, p. 113601, 2017. https://doi.org/10.1103/physrevlett.118.113601.Search in Google Scholar PubMed

[39] H. T. Dung, L. Knöll, and D. G. Welsch, “Resonant dipole-dipole interaction in the presence of dispersing and absorbing surroundings,” Phys. Rev. A, vol. 66, p. 063810, 2002. https://doi.org/10.1103/physreva.66.063810.Search in Google Scholar

[40] F. Minganti, A. Miranowicz, R. W. Chhajlany, and F. Nori, “Quantum exceptional points of non-Hermitian Hamiltonians and Liouvillians: the effects of quantum jumps,” Phys. Rev. A, vol. 100, p. 062131, 2019. https://doi.org/10.1103/physreva.100.062131.Search in Google Scholar

[41] J. Rui, D. Wei, A. Rubio-Abadal, et al.., “A subradiant optical mirror formed by a single structured atomic layer,” Nature, vol. 583, p. 369, 2020. https://doi.org/10.1038/s41586-020-2463-x.Search in Google Scholar PubMed

[42] J. D. Jackson, Classical Electrodynamics, New York, American Association of Physics Teachers, 1999.Search in Google Scholar

[43] D. E. Chang, A. S. Sørensen, E. A. Demler, and M. D. Lukin, “A single-photon transistor using nanoscale surface plasmons,” Nat. Phys., vol. 3, p. 807, 2007. https://doi.org/10.1038/nphys708.Search in Google Scholar

[44] M. O. Scully and M. S. Zubairy, Quantum Optics, New York, American Association of Physics Teachers, 1999.10.1119/1.19344Search in Google Scholar

[45] N. Nefedkin, E. Andrianov, A. Zyablovsky, et al.., “Superradiance of a subwavelength array of classical nonlinear emitters,” Opt. Express, vol. 24, p. 3464, 2016. https://doi.org/10.1364/oe.24.003464.Search in Google Scholar

[46] N. Nefedkin, E. Andrianov, A. Pukhov, and A. Vinogradov, “Superradiance enhancement by bad-cavity resonator,” Laser Phys., vol. 27, p. 065201, 2017. https://doi.org/10.1088/1555-6611/aa6f6d.Search in Google Scholar

Received: 2022-09-29

Accepted: 2022-12-02

Published Online: 2023-01-09

This work is licensed under the Creative Commons Attribution 4.0 International License.

Articles in the same Issue

https://doi.org/10.1515/nanoph-2022-0596

Keywords for this article

atomic arrays; nonreciprocity; quantum metasurface; subradiance

Creative Commons

BY 4.0