Solitonic behaviour in coupled multi atom–cavity systems

In view of current research effort on semiconducting photonic-crystal defect microcavities, we consider the dynamics of an array of coupled optical cavities, each containing an ensemble of qubits. By concentrating on the strong coupling regime, we analytically prove that the nonlinearity inherent in the dynamics of each ensemble coupled to the respective cavity field allows the formation of solitons. We further show how the use of the Holstein–Primakoff transformation and the large-detuning limit with the cavity allows one to recover the Bose–Hubbard model.


2
The fields of nano-and micro-photonics have extensively taken advantage of the properties of versatility exhibited by photonic-crystal nano-or micro-cavities to design and fabricate compact semiconductor lasers and coupled-cavity waveguides. Recently, interest in the physics of photonic-crystal cavities has been extended to their coupling with light emitters for the construction of on-chip photonic sources and guides at the nano-scale. On the other hand, fundamental research has recently been conducted towards the use of coupled-cavity systems as potential candidates for the observation of quantum cooperative phenomena in strongly correlated many-body systems [1]. Superfluid-to-Mott-insulator quantum phase transitions and glassy phases have been theoretically predicted in arrays of mutually coupled resonators interacting with quantum particles. At the basis of these interesting effects is the nonlinear character of the so-called polariton, a combined state of a photon and a two-level quantum system. For easiness of language, we refer to such two-level systems as qubits until we address the physical setup we suggest. The flexibility of these systems, which makes them exploitable quantum simulators, can be used to realize effective multiple-spin dynamics useful for quantum information processing [2]. Intriguing possibilities are offered by arrays of cavities interacting with qubits: photonic Mott insulators can be achieved from qubit-mediated nonlinear photon-photon interactions [3], whereas the combination of intracavity qubit-photon interactions and intercavity tunnelling is powerful for nonlinear-optical operations. The general scheme of mutually coupled systems with embedded nonlinearities arising from interactions described at the quantum level will certainly be useful to describe quantum transport across nanostructures. In this respect, the consideration of these systems would set the ground, at a fundamental level, for future applicability in nanotechnology.
Here we explore nonlinearity originated from a cooperative effect in an ensemble of qubits, which can be embodied by nanoparticles such as quantum dots, and report for the first time the formation of solitary waves (or solitons) [4]-[6] in a periodic array of coupledcavity-nanoparticle systems [7]. Solitons, i.e. localized nonlinear waves characterized by striking stability against dynamical perturbations, are considered as one of the most remarkable effects in nonlinear optics [8]. The generation of solitons in classical systems of nonlinear coupled resonators was proposed in [6], extending the idea of linear waveguides realized through coupled cavities given in [9]. The situation we consider is different from a classical array of media with the Kerr nonlinearity [9] as we do not expand the qubit-medium polarization as a power series of the field or assume fast response of the qubit system relative to the field. Differently, we treat the dynamics of both cavity field and qubit ensembles on equal footing. Our source of nonlinearity is intrinsic to each ensemble of qubits. The identification of solitonic behaviour in a system of current theoretical and experimental interest is an important step forward in the grounding of coupled-cavity systems as exploitable quantum simulators. Our results are derived using a simple technique permitting the consideration of any desired/appropriate order of nonlinearity; nevertheless, catching the most salient features of the full dynamics at hand. The extension to bidimensional configurations and the treatment of the fully discrete case are possible. Semiconducting photonic-crystal microcavities (doped with multiple impurities) are a foreseeable setup for our proposal [10]. By using a few qubits per cavity, we address an interesting working point between the Bose-Hubbard model (found when each cavity is off-resonantly interacting with many qubits) and the fully microscopic polaritonic dynamics [1].
The remainder of this paper is organized as follows. In section 1, we introduce the physical model we consider and discuss in detail the steps needed in order to arrive at a cubic nonlinear Schrödinger equation. Section 2 is devoted to the discussion of a few important points related to the generality of our results and their possible extensions. In section 3, we briefly address the experimental setting where our proposal can be implemented, while section 4 summarizes our findings and points out possible directions for future investigations.

Model and results
The salient features of the mechanism under consideration are elucidated by addressing the case of a linear array of coupled microcavities. For simplicity, we describe the cavities as M equally spaced and independent resonators. Each cavity contains N identical qubits whose ground (excited) state we denote as |g (|e ). The cavities are sufficiently close to each other to allow for their mutual coupling through evanescent-photon hopping (see figure 1). The cavity coupling drops exponentially with the distance so that we take only nearest-neighbour interactions. Each qubit interacts with the respective cavity via the electric-dipole coupling. The volume of each cavity can be made very small (we call λ c the wavelength of the field) so that the coupling within each cavity can be large [11]. Moreover, as discussed later, we consider fixed positions of the qubits within a cavity, as in the case of quantum dots grown in semiconducting matrices. The free energy of the system isĤ 0 = M j=1 (ω câ † where ω c is the frequency of the cavity fields, ω eg the energy spacing of the qubits,â † j the creation operator of mode j andσ z q j the z-Pauli matrix for the qth qubit interacting with the jth cavity of the array (we seth = 1). At resonance ω c ω eg and, in the interaction picture with respect toĤ 0 , the coupling Hamiltonian readsĤ where h.c. stands for the hermitian conjugate, J j is the photon-hopping strength for the jth pair of consecutive cavities and j is the coupling rate between the jth cavity field and the 4 corresponding qubit ensemble [12]. In analogy with the common approach to cooperative effects such as optical bistability [7], phase factors accounting for the spatial distribution (and/or finite dimensions) of the elements of a qubit system are incorporated using the collective operatorŜ ± j = qσ ± q j e −ik c x q with k c being the wavevector of the cavity field, x q the position of the qth element inside a cavity with respect to a chosen reference frame andσ − q j = (σ + q j ) † = |g q j e|. From here on, we assume linear dimensions of each qubit system small with respect to the cavity wavelength λ c . Moreover, in order to fix the ideas, the parameters are taken as real and homogeneous so that J j = J and j = , for any j, although the generalization is straightforward. On the other hand, the seminal works on the behaviour of multi-atom systems in cavities conducted by Bonifacio et al and by Haake and Glauber [12] have shown that atomic cooperative effects are granted as far as the cooperative length (i.e. 'the maximum distance over which atoms in the radiating medium can cooperate in the generation of a superradiant pulse' [12]) is larger than the wavelength of the radiation interacting with the atomic system. This condition is largely satisfied in the range of frequencies we are interested in. The statistical properties of the system, under this condition, can thus be gathered without considering the effects inŜ ± j induced by the relative positions of the qubits within a given cavity. Our assumption of small qubit sample, having long cooperative length, thus allows us to safely neglect any position dependence in the collective operators, which becomeŜ − j = N q=1σ − q j . One could relax the assumption of equal cavity frequencies, keeping in mind that, in the interaction picture with respect to the correspondingĤ 0 , both the inter-cavity and the qubit-light coupling rates acquire time-dependent phase factors. The condition of quasiresonance within each qubit-cavity systems allows us to keep the latter as time-independent. This would also be the case for the inter-cavity couplings in a situation of quasi-homogeneous cavity frequencies. However, an explicit, non-negligible time dependence of J i 's would result in a more complicated physical picture requiring a numerical assessment. As this is not the focus of the present discussion and the model encompassed byĤ is a good approximation for experimentally achievable physical situations (due to the progress achieved in nanolithography, which allows to set very precisely the geometrical factors determining the cavity frequencies in an array [10]), we restrict our attention to cavities having all the same frequency.
The first term (and its hermitian conjugate) inĤ of equation (1) describes the inter-cavity coupling while the second term accounts for in situ exchange of excitations between a qubit system and the respective field (for ω c , we use the rotating wave approximation). The extended Hilbert space of the whole system makes the goal of tackling the dynamics ruled bŷ H formidable. However, we can exploit the collective behaviour of each ensemble by using the Holstein-Primakoff (HP) transformation, which maps a physical collective spin into an effective boson [13]. Assuming N 1, we introduce the HP operatorsb † guarantees that the collective operatorsŜ ±,z j satisfy the necessary algebraic structure. In terms of effective HP bosons, equation (1) takes the form The enhancement of the intracavity coupling due to the collective nature of the coupling is clearly seen. The nonlinear interaction between the two bosons entering the dynamics is encompassed byÂ j 's, which depends on the number of excitations in the state of the HP boson. The strength of any nonlinear effect results from a trade-off between b † jb j and N . Ressayre and Tallet [14] treated collective atomic effects using the HP transformation and showed the emergence of superradiance effects under the condition that the size of the sample is smaller than the cooperative length.
For our discussion, it is useful to expandÂ j in power series with respect to 1/N . Although any order in this expansion could be taken, for clarity of presentation we consider the situation of a 'mesoscopic' number of qubits per ensemble, i.e. N is such thatÂ j 1 −b † jb j /2N . Physically, this implies that the number of qubits per cavity is large enough for the HP transformation to be valid but sufficiently small not to blur any nonlinearity. We thus find the form of the nonlinear coupling, up to this order of approximation The structure of this Hamiltonian is interesting. The term (b † jâ j + h.c.) shows linear interaction, exchanging a quantum between the field and the ensemble of qubits. On the other hand, the extra dependence onb † jb j in the termsb † jâ j (b † jb j ) + h.c. makes the Hamiltonian manifest nonlinear behaviour which, noticeably, is not of the cross-Kerr form [7]. It is worth stressing that, while in the usual approach to coupled-cavity problems, the third term in equation (3), i.e. the photonhopping term, is diagonalized by means of a canonical transformation introducing collective cavity modes [1], here we keep its structure which is instrumental to the following discussion. Interestingly, by assuming the conditions for adiabatic elimination of the qubit degrees of freedom, eachb j would effectively become proportional toâ j and (3) would reduce to the Bose-Hubbard model [1]. We explicitly consider the case of cavities having small dissipation rates with respect to the effective rate of nonlinear dynamics so that the Heisenberg equations for the two bosonic species involved in the dynamics are The coupled nature and the discreteness of these nonlinear equations make their direct solution extremely demanding. In order to tackle the challenge, we put equations (4) into a manageable form using plausible physical assumptions. First, we eliminate the difficulties related to the non-commutativity of the operators in equations (4). We neglect the quantum fluctuations of the operatorsb j andb † j calculated over each qubit's state and concentrate on their mean values, thus providing a mean-field analysis, which is nevertheless sufficient to gather the crucial points of our investigation. Moreover, we take each cavity as prepared in a coherent state (by coupling the array to an external source). This allows the operators in the dynamical equations to be replaced with complex scalar functions as 4â l → α l andb l → β l . Equations (4) 6 now become As anticipated, our approach is designed so as to allow the generalization of these equations to the case of arbitrary-order expansion in the operatorÂ j . So far, we retained the discreteness of the equations of motion. For a clear understanding of our results (preserving the important aspects of the dynamics), we take the continuous limit where the function γ depends on a position-variable x along the array of cavities (here γ = α, β and their hermitian conjugates). This simply means that the inter-distance d between two consecutive qubit-cavity systems is much smaller than the scale represented by the wavelength λ of the wave-like excitation propagating across the array. That is, by taking k = 2π/λ, it must be kd 2π so that γ j → γ (x, t). The photon-hopping contribution to the first of equations (3) is thus modified reminiscent of a Schrödinger equation with a cubic nonlinearity (SCN) [4], which admits solitonic solutions. However, the equation at hand differs from an SCN as the nonlinear term couples the dynamics of real and effective bosons. Instead of using numerical methods for the solution of the problem (such as Newtonian relaxation techniques [16]), we would like to get a clear picture of the physical process. We thus decide to use an approach that can capture the pivotal features of the physics behind this problem. In [17], a multiple-scale technique has been used to investigate gap solitons in nonlinear periodic structures. The flexibility of this technique and its successful applications so far suggest its adaptation to the present situation. The method is based on the expansion of both time and space derivatives in terms of mutually independent time-and length scales v p = µ p v (v = t, x) according to ∂ v = p 0 µ p ∂ v p . Similarly γ = p 1 µ p γ ( p) , where we have introduced the small parameter µ. In its essence, the multiple-scale method utilizes a separation of spatial and temporal scales analogous to other standard techniques used in quantum optics, such as the slowly varying envelope approximation [7]. Following [17], we stop the expansion at order p = 3 in γ , which is enough, in the conditions of small nonlinearity at hand, to encompass the important aspects of the system's dynamics. We now solve the problem corresponding to a scale of order p − 1 and use it in the one of order p. By replacing the scale expansion in the equations of motion (we drop any explicit dependence on x and t) and collecting terms corresponding to the same power of µ, we get the universal structures α,β depending on higher-order derivatives of γ ( p−k) (k ∈ Z). We proceed stepwise: at order µ we have (1) α,β = 0 so that α (1) and β (1) satisfy linear equations in the slow variables x 0 and t 0 which are combined to give In order to solve this equation, we take the plane-wave ansatz α (1) = Ee i(kx 0 −ωt 0 ) with E an envelope function depending on space and time scales faster than x 0 and t 0 and ω the frequency of such a carrier excitation. The associated solvability condition leads to the dispersion relation ω ± (k) = −J (1 − d 2 k 2 /2) ± J 2 (1 − d 2 k 2 /2) 2 + 2 N , thus defining 7 effective 'optical' and 'acoustic' branches (corresponding to ω + (k) and ω − (k), respectively), in analogy with a diatomic crystal. The difference with our case is that the two species in the crystal (i.e. true and effective bosons) share the same site location. We remark that ω ± (k) are the long-wavelength approximations of the branches coming from the exact dispersion relation for the discrete equations (4) in the linear limit. These readω ± (k) = −J cos(kd) ± J 2 cos 2 (kd) + 2 N , which becomes ω ± (k) for kd 1 and exhibits a band gap at the edge of the first Brillouin zone [19]. Here, for definiteness, we concentrate on the optical branch. The previous level of solution leaves E unknown. Its form is determined by going to higher order in p and imposing the appropriate solvability conditions (the nullity of any secular terms, as required in order for the multiple-scale approach to hold [17]). At order µ 2 this results in the equation (∂ t 1 + v g + ∂ x 1 )E = 0, so that the envelope function must depend on the variable ξ = x 1 − v g + t 1 with the group velocity v g + = ∂ k ω + = (2k J d 2 ω 2 + )/(ω 2 + + 2 N ). The iteration of this approach leads to the general SCN where , E = ε/µ, χ = ξ/µ and t 2 = µ 2 t, as implied by the multiple-scale definition. Equation (7) can be exactly solved by means of well-known inverse-scattering methods (ISM) and is known to have the solitonic solution [18] , as in our case. We have thus found a solitary-wave behaviour in the array being studied, which was our central aim. Here, η, σ, ν and φ 0 are integration constants to be determined from the boundary conditions associated with a given physical problem.

A few remarks
Some remarks are in order: the use of complex scalar quantities in going from equations (4) to (5) corresponds to neglecting quantum correlations between qubits and fields. This, however, does not imply the classical nature of the predicted solitonic behaviour. In fact, the non-classical features of these excitations depend only on the quantum nonlinear interaction between the ensemble of qubits and cavity fields. The signature of such nonlinearity is in the quantum Heisenberg equations (4), whose form is dictated byĤ hp and the commutation rules of the involved operators. The use of complex numbers (the mean values of the corresponding operators) amounts to neglecting any operator fluctuations and this does not affect the form of the dynamical equations. The next step would be the inclusion of linear fluctuations in equations (5), which can be done following [21]. In particular, Nagasako et al in [21] have shown that the influence of fluctuations on the solitonic solution gathered for the mean expectation values can be rather small (even negligible). The use of similar quantitative analyses in our specific case is the object of ongoing investigation. On a different level, the form of the nonlinearity considered here requires some discussions. First, the inclusion of higher order terms in the expansion ofÂ j 's results in quite small modifications. For instance, for N 10 and including terms of order N −2 , the rate in front of a nonlinear term in equation (3) is modified by less than 3%. Moreover, a very small extra nonlinearity is added (which is ∼2.5% of the rate obtained by retaining only order N −1 ), giving rise to a cubic-quintic nonlinear Schrödinger equation. The solution of the latter has solitonic character as well and is very close to what we have found (due to the negligible rate of quintic nonlinearity). Multi-soliton solutions can also be found through ISM [18]. The transition from nonlinear to linear regime occurs for c 1 c 2 . In the long-wavelength approximation and in the realistic situation of J , this corresponds to √ N /J , which bounds the number of qubits per cavity, in agreement with our analysis. The continuous limit can be bypassed using the discrete dynamical equations so as to find a solution valid within the whole Brillouin zone (the discrete solitons predicted in [20]). This can be done by taking γ l = p µ p γ ( p) (ξ l , µt, φ l ) with ξ l = µ(ld −λt), φ l being a phase difference between each system andλ which is found via solvability conditions. The relaxation of the long-wavelength assumption would give insight into the dynamics within the whole Brillouin zone, thus allowing the study of gap solitons. However, this will require a detailed numerical analysis which is the focus of current investigation. Finally, we mention that the required tradeoff between the number of atoms per cavity and the degree of nonlinearity is a result of the contingent situation we consider and the HP transformation. Obviously, this is not the only option available and it would be interesting to study the effects of different forms of nonlinearity arising from coherent matter-light interactions. In particular, the possibility of having nonlinear effects which grow with the number of atoms per cavity, as in [23], would be a particularly intriguing case to explore. Finally, we briefly discuss the behaviour of the collective-qubit excitation β(x, t). Once the solitonic nature of the envelope function E is ascertained, one can determine the analogous character of β(t) from the condition β(t) = ( √ N /ω)Ee ikx 0 −ωt 0 , which reveals the solitonic nature of the collective qubit-like excitations with a modified amplitude.

The physical setup
We identify in an array of photonic-crystal microcavities a potential candidate for the observation of our predictions [10] (see figures 1(b) and 2). The linear arrangement addressed here can be generalized to a bidimensional configuration by replacing the cavity adjacency matrix, which determines the hopping term in equations (4), with a proper tensor. Each cavity is a single defect (created, for instance, every two or three sites) in the crystal's lattice and is doped with N ∼10 qubits (substitutional Si donor impurities or InAs quantum dots in a GaAs photonic crystal [10], with linear dimensions in the range of tens of nm, the density being determined by adjusting the bulk doping density) [11]. Importantly, the possibility of cooperative behaviour (in the form of superradiance) of quantum dots in a system close to the building block of the one addressed here has been experimentally demonstrated [22].
The inter-cavity spacing can be as small as tenths of the wavelength of the radiation confined in the crystal [10]. At cavity frequencies of hundreds of THz, ∼40 GHz is realistic 9 (for both Si impurities and InAs dots) because of the small mode volume of microcavity fields (linear dimension in the range of 100 nm, which in turn gives an indication of the required energy density). A proper pattern of the array allows /J ∼10 or smaller. For cavity quality factors around 10 5 -10 6 (see [11]), which are within the experimental reality, the inter-cavity tunnelling time is shorter than the photon lifetime in our system, therefore securing a coherent transient period where equation (7) is rigorously valid.
For the spontaneous emission from the excited state of the qubits, we can take advantage of the relatively long coherence time of donor impurities in a bulk semiconductor matrix (1 ns for Si at the frequency considered here and up to 2 ns for InAS quantum dots in GaAs structures [11]). This is important: a solitonic behaviour is known to be possible in conditions of negligible dissipation, which is the usual assumption in most of the analytical approaches. Spontaneous emission can also be quenched by detuning the two-level systems from the respective cavity field by . It is possible to see that equation (7) still holds with c 1,2 and v g+ modified by taking ω + → ω + + , where is the modified dispersion relation of the optical branch. The solitonic nature of the corresponding solution is preserved for / ∼0.1.

Conclusions and perspectives
We have addressed the nonlinear dynamics of an array of coupled cavities from a quantummechanical perspective. Built-in nonlinearities arise from the coupling of the cavities with ensembles of qubits and compete with photon tunnelling so as to create solitons. We have used a simple method to demonstrate this effect, from a quantum standpoint, in a system of coupled photonic-crystal microcavities. Our study adds to the affirmation of coupledcavity systems as flexible quantum simulators and is a potentially useful approach to future studies of quantum solitonic transport in nanotechnology. The mean-field description of our coupled cavity-quantum two-level systems yields a behaviour which takes us smoothly from the polaritonic to the solitonic regime. Future studies will concentrate on the behaviour of quantum fluctuations around the mean-field solutions.