Non-Hermitian wave packet approximation of Bloch optical equations

0.49 MB·English

by  Eric Charron

Non-Hermitian wave packet approximation of Bloch optical equations
Eric Charron1 and Maxim Sukharev2
1Universit´e Paris-Sud, Institut des Sciences Mol´eculaires d’Orsay,
ISMO, CNRS, F-91405 Orsay, France
2Department of Applied Sciences and Mathematics,
Arizona State University, Mesa, Arizona 85212, USA
(Dated: September 20, 2012)
We introduce a new non-Hermitian approximation of Bloch equations. This approximation pro-
videsacompletedescriptionoftheexcitation,relaxationanddecoherencedynamicsofsingleaswell
as ensembles of coupled quantum emitters (atoms or molecules) in weak laser fields, taking into
account collective effects and dephasing. This is demonstrated by computing the numerical wave
packet solution of a time-dependent non-Hermitian Maxwell-Schro¨dinger equation describing the
interaction of electromagnetic radiation with a quantum (atomic or molecular) nano-structure and
by comparing the calculated transmission, reflection, and absorption spectra with those obtained
2 from the numerical solution of the Maxwell-Liouville-von Neumann equation. We provide the key
1 ingredientsforeasy-to-useimplementationoftheproposedschemeusingtime-dependentdecayrates
0 and identify the limits and error scaling of this approximation.
2
PACSnumbers: 42.50.Ct,78.67.-n,32.30.-r,33.20.-t,36.40.Vz,03.65.Yz
p
e
S
I. INTRODUCTION seriesofworksbyNeuhauseretal. [32–34]clearlydemon-

9 strated that even a single molecule in a close proximity
1 of a plasmonic material can significantly alter scattering
Optics of nanoscale materials has attracted consider-
spectra.
] able attention in the past several years [1–4] due to vari-
ll ous important applications [5]. Exploring electrodynam- In many applications (such as optics of molecular lay-
a ics of near-fields associated with subwavelength systems, ers coupled to plasmonic materials [20–23, 29], for in-
h
researchers are now truly dwelling into nanoscale [5–9]. stance) self-consistent modeling relies heavily on the nu-
-
s Owing to both new materials processing techniques [10] merical integration of the corresponding Maxwell-Bloch
e
andcontinuousprogressinlaserphysics[11]theresearch equations assuming that static emitter-emitter interac-
m
in nano-optics is currently transitioning from linear sys- tions can be neglected (which is true for systems at rel-
t. tems, where materials and their relative arrangement atively low densities). Such an approximation results in
a control optical properties [12], to the nonlinear regime expressingthelocalpolarizationintermsofaproductof
m
[13]. The latter expands optical control capabilities far thelocaldensityofquantumemittersandthelocalaver-
- beyond conventional linear optics as in the case of active agedsingleemitter’sdipolemoment[35]. Oneofthefirst
d
n plasmonic materials [14], for instance, combining highly efficient numerical schemes for simulations of nonlinear
o localized electromagnetic (EM) radiation driven by sur- optical phenomena of quantum media driven by exter-
c face plasmon-polaritons (SPP) with non-linear materi- nal classical EM radiation was proposed by Ziolkowski
[ als [15]. et al. [36]. Using a one-dimensional example of ensem-

blesoftwo-levelatomsitwasshownthatthecorrespond-
2 Yet another promising research direction, namely
v optics of highly coupled exciton-polariton systems, is ing Maxwell-Bloch equations can be successfully inte-
7 emerging [16, 17]. It basically reincarnates a part grated using an iterative scheme based on the predictor-
2 corrector method (strongly coupled method). Later on
of research in semiconductors [18, 19], bringing it to
3 this approach has been extended to two- [37] and three-
nanoscale via deposition of ensembles of quantum emit-
3
dimensional systems [38]. Although such a scheme accu-
. ters (molecules [20–23], quantum dots [24–27]) directly
9 ratelycapturesthesystem’sdynamics,itcanbecomeex-
on to plasmonic materials. Even in the linear regime,
0 tensively slow for multidimensional systems [30]. More-
when the external EM radiation is not significantly ex-
2
over this method is limited to two-level systems only. A
citing the quantum sub-system, SPP near-fields can be
1
moreefficienttechniquebasedonthedecouplingofBloch
: strongly coupled to quantum emitters. This manifests
v equationsfrom the Ampere law was proposed in2003 by
itself as a Rabi splitting widely observed in transmission
i Bid´egaray [39]. This latter method, usually referred to
X experiments [28]. Moreover, a new phenomenon, namely
as a weakly coupled method, noticeably improves the ef-
collectivemolecularmodesdrivenbySPPnear-fields,has
r
a beenobserved[21]andrecentlyexplained[29]. Itwasalso ficiency of the numerical integration of Maxwell-Bloch
equations, allowing to consider multilevel quantum me-
shownthatnanoscaleclusterscomprisedofopticallycou-
dia [30].
pled quantum emitters exhibit collective scattering and
absorption [30] similar to Dicke superradiance [31]. It The approach proposed in the present paper is based
is hence important to be able to account for collective on this weakly coupled method [30]. It further improves
effects in a self-consistent manner. We also note that a the numerical efficiency of the Maxwell-Bloch integrator
2
for ensembles of multilevel quantum emitters. By incor- The approach we propose here is different since it is
porating a new non-Hermitian wave packet propagation basedonperturbationtheory. Italsosignificantlyspeeds
technique into the weakly coupled method, we demon- up calculations since it only requires the propagation of
strate that our approach can be successfully applied to a single wave function under the action of an easy-to-
ensembles of multilevel atoms and diatomic molecules. implementtime-dependenteffectiveHamiltonianinorder
The paper is organized as follows. We first introduce to reproduce accurately the full dynamics. We demon-
ourmodelinsectionIIA,using,asanexample,anensem- strate the efficiency and accuracy of our method using,
ble of two-level atoms. The non-Hermitian wave packet first, ensembles of two-level atoms, and then, interacting
approximation is then described in section IIB. Applica- diatomic molecules including both the vibrational and
tions of the proposed method to the case of a nano-layer rotational degrees of freedom. We also discuss the con-
comprisedoftwo-levelatomsarediscussedinsectionIIC. ditions under which our method is no longer valid. All
We then generalize our method to the case of interact- results are compared with data obtained via direct inte-
ingmulti-levelemittersinsectionIID.Finally,usingthis gration of Maxwell-Bloch equations using a weakly cou-
generalizedapproach,wecalculateinsectionIIEtheop- pled method [30].
ticalpropertiesofanano-layerofcoupledmolecules,tak-
ing into account both the vibrational and rotational de-
grees of freedom, and we reveal several interesting fea- A. Atomic two-level system
tures in the absorption spectra. Last section III summa-
rizes our work. Let us first consider a two-level quantum system in-
teracting with EM radiation. We label the two levels as
|0(cid:105) and |1(cid:105), with associated energy eigenvalues (cid:126)ω and
0
(cid:126)ω ,respectively. Thecorrespondingdensitymatrixρˆ(t)
II. THEORETICAL MODEL AND 1
APPLICATIONS satisfies the dissipative Liouville-von Neumann equation
[50]
Aspointedoutintheprevioussection,theimportance ∂ρˆ
of collective effects manifested in recent experiments i(cid:126) =[Hˆ,ρˆ]−i(cid:126)Γˆρˆ (1)
∂t
and theoretical papers calls for the development of self-
consistent models capable of taking into account mutual whereHˆ =Hˆ +Vˆ(t)isthetotalHamiltonianandΓˆ isa
0 i
EMinteractioninensemblesofquantumemitters. While superoperator, taken in the Lindblad form [51], describ-
direct numerical integration of corresponding Maxwell- ingrelaxationanddephasingprocessesunderMarkovap-
Bloch equations can be based on either a strong coupled proximation. The field free Hamiltonian reads
method [36] or a weakly coupled one [39], for multi-level
systems in two- and especially three-dimensions such a Hˆ0 =(cid:126)ω0|0(cid:105)(cid:104)0|+(cid:126)ω1|1(cid:105)(cid:104)1|. (2)
brute-force approach becomes numerically very expen-
The interaction of the two-level system with EM radia-
sive, both in terms of computation (CPU) times and
tion is taken in the form
memory requirements. Indeed, the main disadvantage of
such approaches is that both the CPU times and mem- Vˆ(t)=(cid:126)Ω(t)(cid:0)|1(cid:105)(cid:104)0|+|0(cid:105)(cid:104)1|(cid:1). (3)
ory requirements scale generally at least as N2, where i
h
Nh is the dimension of the Hilbert space of the system. Ω(t) denotes here the instantaneous Rabi frequency as-
This unfavorable scaling is directly related to the size of sociated with the coupling between the quantum system
the reduced density matrix used to describe the system and an external field. In Eq.(1), the non-diagonal ele-
and to the associated number of time-dependent equa- ments of the operator Γˆ include a pure dephasing rate
tions one has to solve to follow the system’s dynamics. γ∗, and the diagonal elements of this operator consist of
By contrast, one has only to solve a reduced set of N
h theradiationlessdecayrateΓoftheexcitedstate. Equa-
time-dependent equations when the system can be de-
tions (1)-(3) lead to the well-known Bloch optical equa-
scribed by a single wave function. It would therefore
tions [35] describing the quantum dynamics of a coupled
be extremely useful to be able to derive a “Schr¨odinger-
two-level system
type” approximation which could describe on an equal
footing the field-induced coherent dynamics of a multi- ρ˙ = iΩ(t)(ρ −ρ )+Γρ (4a)
00 01 10 11
level quantum system and the associated relaxation and (cid:20) 2γ∗+Γ(cid:21)
decoherence processes. This goal has been achieved in ρ˙ = iΩ(t)(ρ −ρ )+ iω − ρ (4b)
01 00 11 B 2 01
the past using different approaches, the three main con-
tributions being the stochastic Schr¨odinger method used (cid:20) 2γ∗+Γ(cid:21)
ρ˙ = iΩ(t)(ρ −ρ )− iω + ρ (4c)
inconjunctionwithaMonteCarlointegrator[40–45],the 10 11 00 B 2 10
Gadzuk jumping wave packet scheme [46, 47] and the ρ˙ = iΩ(t)(ρ −ρ )−Γρ (4d)
11 10 01 11
variational wave packet method [48, 49]. These methods
require the propagation of N different wave functions whereω =ω −ω istheBohrtransitionfrequencyand
f B 1 0
and are therefore mainly attractive when N (cid:28)N . where the dot denotes the time derivative.
f h
3
We assume in the following that the system is initially sity matrix ρs(t)
inthegroundstate|0(cid:105),andwewillshowthat,undersome
assumptions, the subsequent induced excitation and re- ρ˙s00 = iΩ(t)(ρs01−ρs10)+γ1ρs11 (11a)
laxationdynamicscanbeaccuratelydescribedbyanon- (cid:20) γ −γ (cid:21)
ρ˙s = iΩ(t)(ρs −ρs )+ iω − 1 0 ρs (11b)
Hermitian wave packet approximation. 01 00 11 B 2 01
(cid:20) (cid:21)
γ −γ
ρ˙s = iΩ(t)(ρs −ρs )− iω + 1 0 ρs (11c)
10 11 00 B 2 10
B. Non-Hermitian two-level wave packet
ρ˙s = iΩ(t)(ρs −ρs )−γ ρs (11d)
approximation 11 10 01 1 11
to be compared with the exact equations(4a)-(4d).
Within the aforementioned simplified model, the sys- Two obvious choices can then be made for the empiri-
tem’s wave packet |Ψ(t)(cid:105) can be expanded as cal decay parameters γ and γ :
0 1
|Ψ(t)(cid:105)=c0(t)|0(cid:105)+c1(t)|1(cid:105). (5) • The first choice γ1 = Γ allows to reproduce cor-
rectly the equations(4a) and(4d) describing the
populations at the cost of degrading the descrip-
Let us now assume that the ground and excited states
tion of the coherences ρ (t) and ρ (t).
energies include an imaginary part that we denote as 01 10
+(cid:126)γ /2 and −(cid:126)γ /2, respectively. Inserting expan-
0 1 • The second choice, γ −γ = 2γ∗ +Γ, allows to
1 0
sion(5) into the time-dependent Schr¨odinger equation
reproduce correctly the equations(4b) and(4c) de-
scribing the coherences, at the cost of degrading
∂
i(cid:126) |Ψ(t)(cid:105)=Hˆ |Ψ(t)(cid:105) (6) the description of the populations.
∂t
The optimal choice for applications in linear optics of
and projecting it onto states |0(cid:105) and |1(cid:105) yields the fol- nano-materials should be as follows: in weak fields, the
lowing set of coupled equations for the coefficients cn(t) variation of the populations, as described by perturba-
tion theory, is a second order term with respect to the
coupling amplitude while the variation of the coherences
(cid:16) γ (cid:17)
ic˙ = ω +i 0 c +Ω(t)c (7a) is a first order term. For a correct description of the
0 0 2 0 1
quantum dynamics in weak fields, it is important to de-
(cid:16) γ (cid:17)
ic˙ = Ω(t)c + ω −i 1 c . (7b) scribefirstordertermsaccurately. Thereforeweproceed
1 0 1 2 1
with the second choice
We will now derive the differential equations describ-
γ −γ =2γ∗+Γ. (12)
ing the temporal dynamics of the products ρs = c c∗, 1 0
ij i j
where the subscript s corresponds to this simplified non-
FromEqs.(10)and(12)weobtainasetoftwoempirical
Hermitian Schr¨odinger-type model. Our goal is to ob-
time dependent decay rates γ (t) and γ (t) that we can
0 1
tain the decay rates γ and γ which will allow for an
0 1 insert into the Schr¨odinger-type approximation
approximate description of the system’s excitation and
relaxation dynamics. For the evolution of the popula- (2γ∗+Γ)|c (t)|2
γ (t) = 1 (13a)
tions ρs00(t) and ρs11(t), one gets 0 |c0(t)|2−|c1(t)|2
(2γ∗+Γ)|c (t)|2
ρ˙s = iΩ(t)(ρs −ρs )+γ ρs (8a) γ (t) = 0 (13b)
00 01 10 0 00 1 |c (t)|2−|c (t)|2
ρ˙s = iΩ(t)(ρs −ρs )−γ ρs . (8b) 0 1
11 10 01 1 11
Wefinallyobtainthefollowingsetoftime-dependentcou-
The conservation of the total norm pled equations
ρ˙s00(t)+ρ˙s11(t)=0 (9) ic˙ = (cid:18)ω +i(γ∗+Γ/2)|c1|2(cid:19)c +Ω(t)c (14a)
0 0 |c |2−|c |2 0 1
0 1
thus results in (cid:18) (γ∗+Γ/2)|c |2(cid:19)
ic˙ = Ω(t)c + ω −i 0 c (14b)
1 0 1 |c |2−|c |2 1
γ ρs (t)=γ ρs (t) (10) 0 1
0 00 1 11
This system is solved numerically using a fourth order
It therefore appears that, since the populations ρs (t) Runge-Kutta algorithm [52]. Eqs. (14a), (14b) include
00
and ρs (t) are generally time dependent, at least one of twonon-linearnon-Hermitiantermswhichcan(aswewill
11
the two rates γ0 or γ1, which are not yet fully defined, illustratebelow)accuratelyreproducethedissipativedy-
must be taken as a time-dependent function. namics in weak fields, i.e. for |c |2 (cid:28) |c |2 ≈ 1. Indeed,
1 0
Taking Eq.(10) into account, one finally obtains the asonecaneasilyshowthat,inthislimit,Eqs.(11a)-(11d)
followingsetofBlochequationsfortheapproximateden- are strictly equivalent to the exact Bloch Eqs.(4a)-(4d).
4
C. Application to a uniform nano-layer of atoms absorption A(E) spectra of an atomic layer of thickness
∆z =400nmasafunctionoftheincidentphotonenergy
Whileourapproximationisequallyapplicabletothree- E for an atomic transition energy of EB = (cid:126)ωB = 2eV.
dimensional systems, we have chosen, for the sake of The results are shown in Fig.1.
simplicity, to test it on a simplified one-dimensional sys- ThetransmissionT(E)andreflectionR(E)spectraare
tem consisting of a uniform infinite layer of atoms whose obtained from the normalized Poynting vector
thickness ∆z lies in the range of a few hundred nanome- (cid:12) (cid:12)
ters. An incident radiation field propagating in the pos- (cid:12)(cid:12)E(cid:101)xH(cid:101)y(cid:12)(cid:12)
itive z-direction is represented by a transverse-electric S = (cid:12) (cid:12) (18)
mode with respect to the propagation axis. It is char- (cid:12)(cid:12)E(cid:101)x,incH(cid:101)y,inc(cid:12)(cid:12)
acterized by a single in-plane electric field component
at a specific location under and above the atomic layer,
E (z,t) and a single in-plane magnetic field component
x
Hy(z,t). To account for the symmetry of the atomic po- respectivelly, where E(cid:101)x, H(cid:101)y and E(cid:101)x,inc, H(cid:101)y,inc are the
larization response, the atoms in the layer are described Fourier components of the total and incident EM fields.
as two-level systems with the following states: an s-type The absorption spectrum is then simply obtained as
ground state and a p -type excited state. This model
x A(E)=1−T(E)−R(E). (19)
is a one-dimensional simplification of a more general ap-
proach used in Ref. [30] and we refer the reader to the
In Fig. 1 the solution of Maxwell-Liouville-von Neu-
body of this paper for details.
mann equations is shown as a blue solid line while the
The time-domain Maxwell’s equations for the dynam-
open red squares are from the solution of our approxi-
ics of the electromagnetic fields
mate non-Hermitian Schr¨odinger model. One can notice
a perfect agreement of the two methods irrespective of
µ ∂ H = −∂ E (15a)
0 t y z x the atomic density. These results clearly demonstrate
(cid:15)0∂tEx = −∂zHy−∂tPx (15b) that, in the linear regime where the atomic excitation
probabilityremainssmallandvarieslinearlywiththein-
are solved using a generalized finite-difference time-
cidentfieldintensity,asimplewavepacketpropagationis
domain technique where both the electric and magnetic
sufficient to mimic the excitation and relaxation dynam-
fields are propagated in discretized time and space [53].
ics of an ensemble of quantum emitters. Therefore in
In these equations, µ and (cid:15) denote the magnetic per-
0 0 theweakfieldregime, thepropagationofthefulldensity
meability and dielectric permittivity of free space. The
matrix is unnecessary.
macroscopic polarization of the atomic system
At low density, most of the incident radiation is sim-
ply transmitted and a small part of the incident energy
P (z,t)=n(cid:104)µ (cid:105) (16)
x x is absorbed by the atomic ensemble at energies close
to the atomic transition energy of 2 eV. The absorp-
istakenastheexpectationvalueoftheatomictransition
tion spectrumshows theconventional Lorentzian profile.
dipole moment µ , where n is the atomic density.
x
At higher densities, the absorption spectrum is strongly
Aself-consistentmodelisbasedonthenumericalinte-
modified due to the appearance of collective excitation
gration of Maxwell’s equations(15a) and(15b), coupled
modes [30], and a large part of the incident energy is ei-
viaEq.(16)tothequantumdynamics. Inthemean-field
ther absorbed or reflected from the atomic nano-layer at
approximationemployedhereitisassumedthattheden-
photon energies close to the transition energy.
sity matrix of the ensemble is expressed as a product of
The blue lines with squares in Fig. 2 show, as a
density matrices of individual quantum emitters driven
function of the incident field intensity, the relative error
byalocalEMfield. Inordertoaccountfordipole-dipole
|A (E )−A (E )|/A (E )obtainedinthecalculation
interactions within a single grid cell we follow [54] and S B L B L B
of the absorption spectrum at the transition energy E
introduce Lorentz-Lorenz correction for a local electric B
usingtheSchr¨odingerapproximationA (E )whencom-
field term into quantum dynamics according to S B
pared to the solution of the full Liouville-von Neumann
Ex,local =Ex+ 3Pεx , (17) eaqtoumatiicondeAnsLi(tyEBn).=T2h.5e×so1l0id25bmlu−e3lainneditshfeordatshheedlowbleuset
0
line is for the highest density n=2.5×1027m−3.
where E is the solution of Maxwell’s equations (15a), One can see that irrespective of the atomic density
x
(15b) and macroscopic polarization is evaluated accord- the relative error of the non-Hermitian wave packet ap-
ing to Eq. (16). proximation scales linearly with the field intensity. This
The quantum dynamics is evaluated by computing is not surprising. Indeed, solving coupled Schr¨odinger
the atomic dipole moment either using the single- equations(14a)-(14b)isstrictlyequivalenttosolvingthe
atom density matrix(4a)-(4d) or the single-atom wave approximateBlochequations(11a)-(11d). Thelatterdif-
packet(14a)-(14b). To compare two approaches in the fer from the exact Bloch equations (4a)-(4d) by a term
linear regime (i.e. for ρ (cid:28) ρ ≈ 1, |c |2 (cid:28) |c |2 ≈ 1), proportional to the excited state population ρ . The
11 00 1 0 11
we calculate the transmission T(E), reflection R(E) and maximum excited state population ρmax is also shown
11
5
n = 2.5 1025 m-3 n = 2.5 1026 m-3 n = 2.5 1027 m-3 100
1.00 1.0 1.0
0.8 0.8
0.95 10-1 Maximum excited state population
E) 0.6 0.6
T( 0.90 0.4 0.4
0.85 (a) 0.2 (d) 0.2 (g) 10-2 Relative error
0.0 0.0
10-3 10-1 0.4
10-4 10-2 0.3 10-3
R(E) 1100--65 1100--43 00..21 10-4
(b) (e) (h)
10-7 10-5 0.0
0.15 1.0 1.0 10-5
0.8 0.8
) 0.10 0.6 0.6 n = 2.25 1025 m-3
(AE 0.05 0.4 0.4 10-6 n = 2.25 1027 m-3
0.2 0.2
(c) (f) (i)
0.00 1.9 E (2e.0V) 2.1 0.0 1.9 E (2e.0V) 2.1 0.0 1.9 E (2e.0V) 2.1 10-170-6 10-5 10-4 10-3 10-2 10-1
I (a.u.)
0
FIG.1: (Coloronline)TransmissionT(E)(panelsa,d,g),re-
flectionR(E)(panelsb,e,h)andabsorptionA(E)(panelsc,f,i) FIG.2: (Coloronline)Log-Logplotofthemaximumexcited
spectraofanatomiclayerofthickness∆z=400nmasafunc- state population (red solid line with circles) and the relative
tion of the incident photon energy E. The atomic density is error (blue lines with squares) in the calculation of the ab-
n=2.5×1025m−3 inthefirstcolumn,n=2.5×1026m−3 in sorption spectrum A(E) at the transition energy E = (cid:126)ω
B B
thesecondcolumn,andn=2.5×1027m−3inthelastcolumn. using the Schro¨dinger approximation when compared to the
ThedecayrateandpuredephasingrateareΓ=1012s−1 and solutionofthefullLiouville-vonNeumannequationasafunc-
γ∗ = 1015s−1, respectively. The atomic transition energy is tion of the incident field intensity in atomic units. The blue
(cid:126)ω =2eV and the transition dipole moment is 2D. The so- solid line is for the lowest atomic density n=2.5×1025m−3
B
lutionofMaxwell-Liouville-vonNeumannequationsisshown and the blue dashed line is for the highest atomic density
asabluesolidlinewhiletheopenredsquaresarefromtheso- n=2.5×1027m−3. All other parameters are as in Fig. 1.
lutionofourapproximatenon-HermitianSchro¨dingermodel.
The electron coordinates are denoted by the vector r ,
in Fig. 2 as a function of the field intensity. It can be e
and the vector R ≡ (R,Rˆ) represents the internuclear
seen that it also varies linearly in the present weak field
vector.
regime. We can conclude (from Fig. 2 and numerous
calculations we have performed) that as long as the ex- We now separate the global electronic coordinate r
e
citedstatepopulationremainssmallerthan1%ourwave of all electrons into the coordinate r of the core elec-
c
packet approximation can be used safely. trons and the coordinate r of the active electron [55].
The ground X(cid:0)1Σ+(cid:1) electronic state is considered as a
g
2sσ state, and the electronic wave function Φ (r |R) is
g e
D. Generalization to multi-level systems expressedinthemolecularframe(Hund’scase(b)repre-
sentation) as the product
In order to generalize our Schr¨odinger-type approxi-
mation of the excitation and dissipation dynamics of an Φ (r |R)=φ (r |R)R (r)Y (rˆ) (21)
g e g c X 00
ensemble of quantum emitters to a multi-level system
we now consider the case of neutral diatomic molecules.
whereR (r)andY (rˆ)aretheradialandangularparts
More specifically, we have chosen the particular case of X 00
of the electronic wave function associated with the ac-
the ground and first excited electronic states of the Li
2 tiveelectron. Similarly,the2pσ excitedstateofA(cid:0)1Σ+(cid:1)
molecule and we follow the electronic dynamics and the u
symmetry is expressed in the molecular frame as
nuclear motion by expanding the total molecular wave
functionΨ(r ,R,t)usingtheBorn-Oppenheimerexpan-
e
sion Φ (r |R)=φ (r ,r|R)R (r)Y (rˆ). (22)
e e e c A 10
Ψ(r ,R,t)=χ (R,t)Φ (r |R)+χ (R,t)Φ (r |R) (20)
e g g e e e e
Due to the Σ symmetry of both electronic states, the
fwuhnecrteioΦngs(aress|Roc)iaatneddΦwei(trhe|tRhe)dgeronuontedtXhe(cid:0)e1lΣec+tr(cid:1)oannicdwfiarvset aron-dvibχra(tRio,nt)alcatinmbe-edeepxpenadnednetd wonavae lfiumnictteidonssetχogf(Rno,tr)-
g e
excited A(cid:0)1Σ+(cid:1) electronic states of Li , respectively. malized Wigner rotation matrices in order to take into
u 2
6
account the rotational degree of freedom, following [56] where c (t) and c (t) are time-dependent complex coef-
g v
χ (R,t) = (cid:88) χg (R,t)DN∗(Rˆ) (23a) ficients. ϕgv/,Ne(R) denote here the bound ro-vibrational
g N,M M,0 eigenstates of the g/e electronic potentials [58]. It is
N,M
not necessary here to take into account the dissociative
χ (R,t) = (cid:88) χe (R,t)DN∗(Rˆ), (23b) nuclear eigenstates associated with the excited potential
e N,M M,0
N,M since their coupling with the ground vibrational level of
the ground electronic state is negligible.
whereN denotesthemolecularrotationalquantumnum-
We thus arrive at a multi-level system which is very
ber while M denotes its projection on the electric field
similar to the atomic case described in sections IIA and
polarization axis x of the laboratory frame.
IIB, except that the initial ground state is now coupled
Introducing these expansions in the time-dependent
withalargesetofexcitedlevels. Forconvenienceandfor
Schr¨odinger equation (6) describing the molecule-field
an easy comparison with the atomic case, we will label
interaction and projecting onto the electronic and rota-
thegroundstateasstatenumber0andtheexcitedstates
tionalbasisfunctionsyields,inthedipoleapproximation,
as states number j (cid:62) 1. Our reference calculations will
the following set of coupled differential equations for the
nuclear wave packets χg (R,t) and χe (R,t) bebasedonthenumericalsolutionsofthecorresponding
N,M N,M Liouville-von Neumann equations
∂
i(cid:126)∂tχgN,M = HˆNg (R)χgN,M −Ex(t)µAX(R) ρ˙00 = (cid:88)(cid:2)iΩj(t)(ρ0j −ρj0)+Γρjj(cid:3) (28a)
× (cid:88) MN(cid:48),M(cid:48)χe (24a) j(cid:62)1
N,M N(cid:48),M(cid:48)
ρ˙ = iΩ (t)(ρ −ρ )
N(cid:48),M(cid:48) 0j j 00 jj
∂ (cid:20) 2γ∗+Γ(cid:21)
i(cid:126)∂tχeN,M = HˆNe (R)χeN,M −Ex(t)µAX(R) + i(ωj −ω0)− 2 ρ0j (28b)
× (cid:88) MNN,(cid:48)M,M(cid:48)∗χgN(cid:48),M(cid:48) (24b) ρ˙j0 = iΩj(t)(ρjj −ρ00)
N(cid:48),M(cid:48) (cid:20) 2γ∗+Γ(cid:21)
− i(ω −ω )+ ρ (28c)
where µ (R)=(cid:104)R |r|R (cid:105) is the electronic transition j 0 2 0j
AX A X r
dipole. Thero-vibrationalnuclearHamiltoniansHˆg/e(R) ρ˙jj = iΩj(t)(ρj0−ρ0j)−Γρjj (28d)
N
are defined as
where (cid:126)ω is the total energy of the excited state j and
(cid:126)2 (cid:20) ∂2 N(N +1)(cid:21) j
Hˆg/e(R)=− − +V (R) (25) where Ωj(t) denotes the instantaneous Rabi frequency
N 2µ ∂R2 R2 g/e associated with the molecule-field interaction as defined
in Eqs.(24a) and(24b). γ∗ and Γ denote the pure de-
where µ denotes the molecular reduced mass. The po-
phasing rate and the relaxation rate associated with the
tential energy curves V (R) and V (R) associated with
g e
excited states, respectively.
the ground and first excited electronic states of the
In comparison with this “exact” model, our Schr¨odin-
molecule are taken from Ref. [57] and the matrix ele-
ments MN(cid:48),M(cid:48) which couple the nuclear wave packets ger-type approximation will be based on the numerical
N,M solutionsofthecoupledequationsforthetime-dependent
evolvingontheseelectronicpotentialcurvescanbewrit-
expansion coefficients
ten using 3j-symbols as
MNN,(cid:48)M,M(cid:48) = (−1)M(2N√4+π1)N(cid:88)(cid:48),M(cid:48)(cid:18)N0(cid:48) 10 N0 (cid:19) ic˙0 = (cid:20)ω0+iγ02(t)(cid:21)c0+(cid:88)j(cid:62)1Ωj(t)cj (29a)
×(cid:18)MN(cid:48)(cid:48) 01 −NM (cid:19) (26) ic˙j = Ωj(t)c0+(cid:20)ωj −iγj2(t)(cid:21)cj, (29b)
We assume that the molecules are prepared at time wherethetime-dependentdecayratesγ (t)andγ (t)are
0 j
t = 0 in the ro-vibrational level v = 0 and N = 0 of the now defined as
ground electronic state. In weak linearly polarized fields
and except for a phase factor this ground state compo- γ (t) = (2γ∗+Γ)(cid:80)j(cid:62)1|cj(t)|2 (30)
nent of the molecular wave function remains unaffected 0 |c (t)|2−(cid:80) |c (t)|2
0 j(cid:62)1 j
and the excited state component is limited to N = 1
(2γ∗+Γ)|c (t)|2
and M =0. The ground and excited nuclear wave pack- γ (t) = 0 (31)
ets are thus finally expanded in terms of ro-vibrational j |c0(t)|2−(cid:80)j(cid:62)1|cj(t)|2
eigenstates as
One can show that with such a definition of γ (t)
χg (R,t) = c (t)ϕg (R) (27a) 0
0,0 g 0,0 and γj(t), Eqs.(29a) and(29b) are strictly equivalent to
χe (R,t) = (cid:88)c (t)ϕe (R) (27b) Eqs.(28a)-(28d) in the limit of weak couplings. In the
1,0 v v,1
next section we demonstrate that our Schr¨odinger-type
v
7
n = 2.5 1025 m-3 n = 2.5 1027 m-3 As in the case of two-level atoms, one can notice a
0.9 0.6 perfect agreement of the two methods irrespective of the
23 (a) (c)
4 molecular density. This shows that, in the linear regime
1
5 wherethemolecularexcitationprobabilityremainssmall
0.6 0.4
) 6 andvarieslinearlywiththeincidentfieldintensity,asim-
E 0 Zoom
A( 7 ple wave packet propagation is sufficient to mimic the
0.3 0.2
8 excitation, dissipation and decoherence dynamics of an
9 ensemble of multi-level quantum emitters. The propaga-
00..00 00..00 tionofthefulldensitymatrixisthen,again,unnecessary.
It is important to note that at low density, we observe
0.3 0.2 a series of overlapping vibrational resonances which re-
E) flects the vibrational structure of the excited molecular
A( potentialandwhichfollowstheFranck-Condonprinciple
0.6 0.4
[59, 60]. The green labels seen in panel (a) of Fig. 3 in-
(b) (d) dicate the excited state vibrational level responsible for
0.91.6 1.8 2.0 2.2 0.6 1.0 2.0 3.0 the observed resonance.
E (eV) E (eV) At higher density, the absorption spectrum is strongly
distorted and one observes, just like in the atomic case
FIG. 3: (Color online) Absorption spectra A(E) of a Li2 [30], the appearance of collective excitation modes, the
molecular layer of thickness ∆z = 400nm as a function of difference being that a vibrational structure may still be
the incident photon energy E. The molecular density is n= present in some of these collective molecular excitation
2.5×1025m−3 in the left column (panels (a) and (b)) and
modes. A detailed analysis of the physics underlying the
n=2.5×1027m−3 in the right column (panels (c) and (d)).
appearanceoftheseintriguingmolecularcollectivemodes
The solutions of Maxwell-Liouville-von Neumann equations
will be presented in another paper. For high densities,
are shown as blue solid lines in the first raw (panels (a) and
wecouldalsoobservethattheradiationfieldisnottrans-
(c))whiletheredsolidlines(invertedspectra,panels(b)and
(d))arefromthesolutionsofourapproximatenon-Hermitian mitted anymore through the molecular layer in the ab-
Schr¨odinger model. All other parameters are as in Fig. 1. sorption window 1.2eV (cid:54) E (cid:54) 2.5eV of the molecule.
Within this window, the field is either absorbed or re-
flected. The small inset seen in panel (c) of Fig. 3 shows
amagnificationoftheabsorptionspectrumintheenergy
model accurately reproduces the excitation and dissipa- range 1.7eV (cid:54) E (cid:54) 2.1eV corresponding to the “nor-
tiondynamicsofmulti-levelquantumsystemsinthelimit
mal” low-density aborption spectrum. One can see in
of weak couplings, i.e. for (cid:80) |c |2 (cid:28)|c |2 ≈1.
j 0 this inset, and by comparing with panel (a) of Fig. 3,
j(cid:62)1
that the absorption spectrum is not strongly modified at
high molecular density in this energy region.
Figure(4)finallyshowstheCPUtimenecessaryforthe
E. Application to a uniform nano-layer of
molecules calculation of these absorption spectra as a function of
the number of quantum states introduced in this multi-
levelmodelforboththeMaxwell-Liouville-vonNeumann
We consider a simplified one-dimensional system sim-
(bluelinewithcircles)andMaxwell-Schr¨odinger(redline
ilar to the one discussed in case of two-level atoms. A
with squares) approaches. One can observe a substan-
uniform infinite layer with a thickness of ∆z = 400 nm
tial difference in CPU times which can prove of crucial
comprised of Li molecules is exposed to incident lin-
2 importance when one has to deal with realistic three-
early polarized radiation. The incident field propagates
dimensional systems. This origin of the observed gain in
in the positive z-direction. We calculate the absorption
CPU time relies on the necessity of propagating a single
spectrum A(E) of this molecular layer just as we did in
wave function instead of a full density matrix.
section IIC for atoms.
To ascertain the validity of the proposed Schr¨odinger-
type approximation, absorption spectra are represented
in Fig. 3 as a function of the incident photon energy III. SUMMARY
E, in the linear regime, for two different molecular den-
sities: n = 2.5 × 1025m−3 in the left column (panels We proposed a new non-Hermitian approximation of
(a) and (b)) and n=2.5×1027m−3 in the right column Bloch optical equations. Our method provides an accu-
(panels(c)and(d)). Thesolutionsobtainedviaintegrat- rate, complete description of the excitation, relaxation
ingMaxwell-Liouville-vonNeumannequationsareshown and decoherence dynamics of single as well as ensem-
in the upper panels (a) and (c) as blue solid lines while bles of coupled quantum emitters (atoms or molecules)
the ones obtained from our approximate non-Hermitian in weak EM fields, taking into account collective effects
Schr¨odinger model are shown in red in the lower panels anddephasing. Wedemonstratedtheapplicabilityofthe
(b) and (d) (inverted spectra). method by computing optical properties of thin layers
8
600 comprised of two-level atoms and diatomic molecules. It
was shown that, in the limit of weak incident fields, the
Liouville - Maxwell
dynamics of interacting quantum emitters can be suc-
500 Schrödinger - Maxwell
cessfullydescribedbyoursetofapproximatedequations,
which result in a substantial gain both in CPU time and
400 computer memory requirements. These calculations also
)
(s reveal some intriguing new collective molecular excita-
e
m tion modes which will be presented in detail in another
-ti 300 publication. The proposed approach was demonstrated
U
P to provide a substantial increase in numerical efficiency
C
200 for self-consistent simulations.
100
0
2 3 4 5 6 7 8 9 10
Number of states IV. ACKNOWLEDGEMENTS
FIG. 4: (Color online) CPU process time necessary on a
Intel Xeon E5-1650 processor for the calculation of the ab- E.C. would like to acknowledge useful and stimulating
sorption spectrum of a Li2 molecular nano-layer of thickness discussions with O. Atabek and A. Keller from Univer-
∆z =400nm as a function of the number of quantum levels
sit´e Paris-Sud (Orsay). M.S. is grateful to the Univer-
included in the calculation. The blue line with circles is for
sit´e Paris-Sud (Orsay) for financial supports through an
the solution obtained from Maxwell-Liouville-von Neumann
invited Professor position in 2011. E.C. acknowledges
equations while the red line with squares is for our proposed
supportsfromANR(contractAttowaveANR-09-BLAN-
Schr¨odinger-typeapproximation. Thespatialgridhasatotal
size of 2560nm with a spatial step of 1nm. The time prop- 0031-01), and from the EU (Project ITN-2010-264951,
agation is performed on a temporal grid of total size 1.7ps CORINF). We also acknowledge the use of the comput-
with a time step of 1.7as. ing facility cluster GMPCS of the LUMAT federation
(FR LUMAT 2764).
[1] S.A.MaierandH.A.Atwater,J.Appl.Phys.98,011101 N. I. Zheludev, Appl. Phys. Lett. 89, 031118 (2006).
(2005). [16] J. Bellessa, C. Bonnand, J. C. Plenet, and J. Mugnier,
[2] S. Lal, S. Link, and N. J. Halas, Nat. Photon. 1, 641 Phys. Rev. Lett. 93, 036404 (2004).
(2007). [17] D. E. Chang, A. S. Sorensen, P. R. Hemmer, and M. D.
[3] E. Ozbay, Science 311, 189 (2006). Lukin, Phys. Rev. Lett. 97, 053002 (2006).
[4] P. Berini, Adv. Opt. Photon. 1, 484 (2009). [18] V.M.AgranovichandG.C.LaRocca,SolidStateCom-
[5] M. I. Stockman, Opt. Express 19, 22029 (2011). mun. 135, 544 (2005).
[6] S. A. Maier, M. L. Brongersma, P. G. Kik, S. Meltzer, [19] G. Khitrova, H. M. Gibbs, M. Kira, S. W. Koch, and
A. A. G. Requicha, and H. A. Atwater, Adv. Mater. 13, A. Scherer, Nat. Phys. 2, 81 (2006).
1501 (2001). [20] J. Dintinger, S. Klein, F. Bustos, W. L. Barnes, and
[7] F.J.GarciadeAbajo,Rev.Mod.Phys.79,1267(2007). T. W. Ebbesen, Phys. Rev. B 71, 035424 (2005).
[8] T. W. Ebbesen, C. Genet, and S. I. Bozhevolnyi, Phys. [21] T. W. Ebbesen, A. Salomon, and C. Genet, Angew.
Today 61, 44 (2008). Chem. Int. Edit. 48, 8748 (2009).
[9] J.A.Schuller,E.S.Barnard,W.S.Cai,Y.C.Jun,J.S. [22] D. D. Lekeufack, A. Brioude, A. W. Coleman, P. Miele,
White, and M. L. Brongersma, Nature Materials 9, 193 J. Bellessa, L. D. Zeng, and P. Stadelmann, Appl. Phys.
(2010). Lett. 96, 253107 (2010).
[10] E. Hutter and J. H. Fendler, Adv. Mater.ials 16, 1685 [23] A.Berrier,R.Cools,C.Arnold,P.Offermans,M.Crego-
(2004). Calama, S. H. Brongersma, and J. Gomez-Rivas, ACS
[11] R. E. Slusher, Rev. Mod. Phys. 71, S471 (1999). Nano 5, 6226 (2011).
[12] W. L. Barnes and W. A. Murray, Adv. Mater.ials 19, [24] H. Park, A. V. Akimov, A. Mukherjee, C. L. Yu, D. E.
3771 (2007). Chang, A. S. Zibrov, P. R. Hemmer, and M. D. Lukin,
[13] N. Kroo, S. Varro, G. Farkas, P. Dombi, D. Oszetzky, Nature 450, 402 (2007).
A. Nagy, and A. Czitrovszky, J. Mod. Optic. 55, 3203 [25] V. K. Komarala, A. L. Bradley, Y. P. Rakovich, S. J.
(2008). Byrne, Y. K. Gun’ko, and A. L. Rogach, Appl. Phys.
[14] S.Wuestner,A.Pusch,K.L.Tsakmakidis,J.M.Hamm, Lett. 93, 123102 (2008).
and O. Hess, Philos. T. R. Soc. A 369, 3525 (2011). [26] M. L. Andersen, S. Stobbe, A. S. Sorensen, and P. Lo-
[15] A.V.Krasavin,K.F.MacDonald,A.S.Schwanecke,and dahl, Nat. Phys. 7, 215 (2011).
9
[27] N.Livneh,A.Strauss,I.Schwarz,I.Rosenberg,A.Zim- Physics 111, 10126 (1999).
ran, S. Yochelis, G. Chen, U. Banin, Y. Paltiel, and [46] J. Gadzuk, L. Richter, S. Buntin, D. King, and R. Ca-
R. Rapaport, Nano Lett. 11, 1630 (2011). vanagh, Surface Science 235, 317 (1990).
[28] D. E. Gomez, K. C. Vernon, P. Mulvaney, and T. J. [47] K. Finger and P. Saalfrank, Chemical Physics Letters
Davis, Nano Lett. 10, 274 (2010). 268, 291 (1997).
[29] A. Salomon, R. J. Gordon, Y. Prior, T. Seideman, and [48] T. Gerdts and U. Manthe, The Journal of Chemical
M. Sukharev, Phys. Rev. Lett. 109, 073002 (2012). Physics 106, 3017 (1997).
[30] M. Sukharev and A. Nitzan, Phys. Rev. A 84 (2011). [49] L.Pesce,T.Gerdts,U.Manthe,andP.Saalfrank,Chem-
[31] V. I. E. A. V. Andreev and Y. A. Il’inski˘ı, Sov. Phys. ical Physics Letters 288, 383 (1998).
Usp. 23, 493 (1980). [50] K. Blum, Density matrix theory and applications
[32] K. Lopata and D. Neuhauser, J. Chem. Phys. 127, (Plenum Press, New York, 1996), 2nd ed.
154715 (2007). [51] H.-P. Breuer and F. Petruccione, The Theory of Open
[33] K.LopataandD.Neuhauser,J.Chem.Phys.130(2009). Quantum Systems (Oxford University Press, New York,
[34] K. Lopata and D. Neuhauser, J. Chem. Phys. 131, 2002).
014701 (2009). [52] W. H. Press, B. P. Flannery, S. A. Teukolsky, and
[35] L.AllenandJ.H.Eberly,Opticalresonanceandtwo-level W. T. Vetterling, Numerical Recipes: The Art of Sci-
atoms (Wiley, New York,, 1975). entific Computing (Cambridge University Press, Cam-
[36] R.W.Ziolkowski,J.M.Arnold,andD.M.Gogny,Phys. bridge, 2002).
Rev. A 52, 3082 (1995). [53] A. Taflove and S. C. Hagness, Computational elec-
[37] G. Slavcheva, J. M. Arnold, I. Wallace, and R. W. Zi- trodynamics: the finite-difference time-domain method
olkowski, Phys. Rev. A 66, 063418 (2002). (Artech House, Boston, 2005), 3rd ed.
[38] A. Fratalocchi, C. Conti, and G. Ruocco, Phys. Rev. A [54] C.M.BowdenandJ.P.Dowling,Phys.Rev.A47,1247
78, 013806 (2008). (1993).
[39] B. Bid´egaray, Numer. Meth. Part. D. E. 19, 284 (2003). [55] E. Charron and M. Raoult, Phys. Rev. A 74, 033407
[40] J.Dalibard,Y.Castin,andK.Mølmer,Phys.Rev.Lett. (2006).
68, 580 (1992). [56] E.Charron,A.Giusti-Suzor,andF.H.Mies,Phys.Rev.
[41] R.Dum,P.Zoller,andH.Ritsch,Phys.Rev.A45,4879 A 49, R641 (1994).
(1992). [57] I.Schmidt-Mink,W.Mu¨ller,andW.Meyer,Chem.Phys.
[42] C.W.Gardiner,A.S.Parkins,andP.Zoller,Phys.Rev. 92, 263 (1985).
A 46, 4363 (1992). [58] E. Charron and A. Suzor-Weiner, J. Chem. Phys. 108,
[43] R. Dum, A. S. Parkins, P. Zoller, and C. W. Gardiner, 3922 (1998).
Phys. Rev. A 46, 4382 (1992). [59] J. Franck and E. G. Dymond, Trans. Faraday Soc. 21,
[44] N.GisinandI.C.Percival,JournalofPhysicsA:Math- 536 (1926).
ematical and General 25, 5677 (1992). [60] E. Condon, Phys. Rev. 28, 1182 (1926).
[45] D. E. Makarov and H. Metiu, The Journal of Chemical