Simulating pump-probe photo-electron and absorption spectroscopy on the attosecond time-scale with time-dependent density-functional theory

1.9 MB·English

by  Umberto De Giovannini

Simulating pump-probe photo-electron and absorption spectroscopy on the
attosecond time-scale with time-dependent density-functional theory
Umberto De Giovannini,1,a) Gustavo Brunetto,1,2 Alberto Castro,3,b) Jessica Walkenhorst,1 and Angel
Rubio1,4,5,c)
1)Nano-Bio Spectroscopy Group and ETSF Scientific Development Center, University of the Basque Country
UPV/EHU, Avenida de Tolosa 72, 20018 San Sebastian, Spain
2)Instituto de F´ısica ”Gleb Wataghin”, Universidade Estadual de Campinas, 13083-970, Campinas, Sa˜o Paulo, Brazil
3)ARAID Foundation – Institute for Biocomputation and Physics of Complex Systems (University of Zaragoza),
Zaragoza, Spain
4)Centro de F´ısica de Materiales CSIC-UPV/EHU-MPC and DIPC, Av. Tolosa 72, E-20018 San Sebasti´an, Spain
5)Fritz-Haber-Institut der Max-Planck-Gesellschaft,Faradayweg 4-6, D-14195 Berlin, Germany
(Dated: 10 January 2013)
3
1 Molecular absorption and photo-electron spectra can be efficiently predicted with real-time time-dependent
0
density-functionaltheory(TDDFT).Weshowherehowthesetechniquescanbeeasilyextendedtostudytime-
2
resolved pump-probe experiments in which a system response (absorption or electron emission) to a probe

n pulse, is measured in an excited state. This simulation tool helps to interpret the fast evolving attosecond
a time-resolved spectroscopic experiments, where the electronic motion must be followed at its natural time-
J
scale. Weshowhowtheextradegreesoffreedom(pumppulseduration,intensity,frequency,andtime-delay),
9
which are absent in a conventional steady state experiment, provide additional information about electronic

structure and dynamics that improve a system characterization. As an extension of this approach, time-
]
s dependent 2D spectroscopies can also be simulated, in principle, for large-scale structures and extended
u
systems.
l
c
-
m
I. INTRODUCTION visiblelonger(fewfemtoseconds)pulseusedforitsgener-
t ation. Combining two XUV attosecond pulses is in prin-
a
ciplepossible(andhasbeentheoreticallyanalysed7),but
. Pump-probe experiments are the preferred tech-
s
unfortunately the low outputs of current XUV attosec-
c nique to study the dynamical behaviour of atoms and
i molecules: the dynamics triggered by the pump pulse ondpulsesleadtomuchtooweaksignals. Anotherchoice
s
to make is the final observable, i.e. what kind of system
y can be monitored by the time-dependent reaction of the
h system to the probe pulse, a reaction that can be mea- reaction is to be measured as a function of the time de-
p suredintermsof,forexample,theabsorptionofthepulse lay. In this work we focus on two common choices. One,
[ intensity, or of the emission of electrons.1 The time reso- observingthe emissionof electrons(their energies, angu-
lar distribution, or total yield) from the pumped system
1 lution of these experiments is mainly limited by the du-
due to the probe pulse. This can be called time-resolved
v ration of the pulses – although it is also limited by the
8 abilityoftheexperimentertoascertaintheirrelativetime photo-electron spectroscopy (TRPES). Two, observing
5 the optical absorption of the probe signal, which can be
delay and shape. In order to precisely fix this delay, the
9 calledtimeresolvedabsorptionspectroscopy,ortransient
two pulses are coherently synchronized – in fact, they
1 absorption spectroscopy (TAS).
have the same origin, or one of them is used to generate
.
1 theother–,sothatthedelayisgaugedbyanopticalpath Bothtechniquescanofcoursebeusedtolookatlonger
0 difference.2 The electron dynamics has a natural time- timeresolutions–andthereisalreadyasubstantialbody
3
scale in the range of attoseconds, and therefore could of literature describing such experiments. If we look at
1
notbestudiedbypump-probespectroscopyuntilthead- molecular reaction on the scale of tens or hundreds of
:
v vent of attosecond-pulse laser sources one decade ago.3,4 femtoseconds, the atomic structure will have time to re-
i
X Nowadays time-resolved spectroscopy can be utilized to arrange. These techniques are thus mainly employed to
monitor electron dynamics in real time, giving birth to observe modification, creation, or destruction of bonds,
r
a the field of attosecond physics.5,6 a field now named femtochemistry.8
A wealth of possibilities exists, depending on the fre- TAS, for example, has been succesfully employed to
quencies, durations and intensities of the two pulses. A watch the first photo-synthetic events in cholorophylls
common set-up in attosecond physics employs a XUV and carotenoids,9 that transform the energy gained by
attosecond pulse and the relatively more intense NIR or lightabsorptionintomolecularrearrangements. Areview
describing the essentials of this technique can be found
in Ref. 10. Note, however, that in addition to following
chemical reactions, femtosecond-long pulses may also be
a)Electronicmail: [email protected] used for example for characterizing the final electronic
b)Electronicmail: acastro@bifi.es quantum state of ionized atoms.11
c)Electronicmail: [email protected] In TRPES, the probe pulse generates free electrons
2
throughphoto-ionization, andonemeasurestheirenergy systems. Finally, the above-mentioned experiment of
orangulardistributionasafunctionoftime;ifthistimeis Goulielmalkis et al22 was analized with the model de-
inthefemtosecondscaleonecanfollowmoleculardynam- scribedinRef.32,whichtreatedthepumpIRpulsenon-
ics in the gas phase, as demonstrated already in the mid perturbatively.
1990s,12,13 althoughthistechniquehadalreadybeenem- Indeed,itwouldbedesirabletoanalyzetheseprocesses
ployed to follow electronic dynamics on surfaces.14 This with a non-perturbative theory (since at least one of the
methodology is well documented in recent articles.15–19 pulses is usually very intense), which at the same time
If the goal –as in this work– is to study the electronic is capable of going beyond the SAE and accounting for
dynamics only, disentangling them from the vibronic de- many-electron interaction effects. This last fact is rel-
grees of freedom, then one must move down these spec- evant since the attosecond time resolution obtained in
troscopic methods to the attosecond regime.20 In this this type of experiments is able to unveil the fast dy-
regime,bothTASandTRPEShaverecentlybeendemon- namical electron-electron interaction effects. The SAE,
strated. Regarding TAS, we may cite as a prototypical whichessentiallyassumesthatonlyoneelectronactively
example the recent experiment of Holler et al21 where, respondstothelaserpulse, hasbeensuccessfullyusedto
thetransientabsorptionofanattosecondpulsetrain(cre- interpretmanystrong-fieldprocesses. However,itsrange
ated by high harmonic generation) by a Helium gas tar- ofvalidityislimited, androughlyspeakingitisexpected
get, was studied in the presence of an intense IR pulse. to fail whenever the energies of multielectron excitations
The absorption was observed to oscillate as a function become comparable to the laser frequencies or the single
of the time-delay of pump and probe. Another example electron excitations.33
is the real-time observation of valence electron motion Time-dependentdensityfunctionaltheory(TDDFT)34
reported by Goulielmakis et al.22 in principle meets all requirements: it may be used non-
Several cases of use of TRPES with attosecond pulses perturbatively, includes the electron-electron interaction
have also been recently reported. For example, Uiber- andcanhandleout-of-equilibriumsituations. Ithasbeen
acker et al.23 could observe in real time the light in- routinely used in the past decades to study the electron
ducedelectrontunnelingprovokedbyastrongNIRpulse, dynamics in condensed matter in equilibrium. By this
demonstrating how this electron tunneling can be used we mean that, usually, one computes the linear or non-
to probe short lived electronic states. Smirnova et al24 linear response properties of systems in the ground state
also studied the ionization of an atom by an attosecond (oratthermalequlibrium). Inpump-probeexperiments,
XUV pulse in the presence of an intense laser pulse, as however,onemustcomputetheresponseofasystemthat
a function of the time delay between both. Jonnsson et is being driven out of equilibrium by an initial pulse. In
al25 employed attosecond pulse trains and a Helium tar- this work, we will explore the usability of TDDFT for
get, and attosecond photoelectron spectroscopy was also this purpose, and show how, at least for the two cases of
demonstrated to yield useful information for condensed TAS and TRPES, the extension is straightforward.
matter systems.26
All these advances demand an appropriate theoretical
modeling. Theuseofmorethanonepulseoflightintrin- II. THEORY
sically requires to go beyond any “linear spectroscopy”
technique–althoughifthepulsesareweakaperturbative Density-functional theory35 establishes a one-to-one
treatment may still be in order. This non-linear behav- correspondence between the ground-state density and
ior, provides much more information about the system the external potential of a many-electron system. This
at the cost of a increasingly difficult analysis. The use of implies that any system property is, in principle, a
two (or more) coherent pulses of light, with fine control ground-state density functional. The computation of
over their shape (sometimes called a “multidimensional the ground-state density usually follows the Kohn-Sham
analysis”), permits a deeper characterization. This fact (KS)36 scheme, in which one utilizes a fictitious system
wasalreadyacknowledgedinthefieldofnuclearmagnetic of non-interacting electrons that has the same ground
resonance, or later in femtochemistry – see for example state density. For excited states, properties, however,
Refs. 27 and 28 for theoretical treatments of these cases. or in order to simulate the behavior of the system in
A recent theoretical analysis of attosecond TAS based time-dependent external fields, one must use its time-
on perturbation theory was given by Baggesen et al.7 dependent version, TDDFT.34,37,38
Gaarde et al29 presented a study for relatively weak In the case of TDDFT, a one-to-one correspondence
pumping IR pulses in combination with XUV ultrafast also exists between the time-dependent densities and
probes,forHeliumtargetsandbasedonthesingleactive potentials. One also uses an auxiliary fictitious sys-
electron approximation (SAE). Very recently, the exper- tem of non-interacting electrons that produces the same
imented reported by Ott et al.30, in which the ultrafast time-dependent density. This substitution is the source
TAS of Helium displayed characteristic beyond-SAE fea- of the great computational simplification, since a non-
tures, was theoretically analyzed in Ref. 31, utilizing an interacting system of electrons can in general be repre-
exact solution of the time-dependent Schr¨odinger equa- sented by a single Slater determinant formed by a set
tion, that cannot however be easily extended to larger of “Kohn-Sham” orbitals, ϕ (i = 1,...,N/2). We will
i
3
assume a spin-compensated system of N electrons dou- be that of the ground state. If we want to analyze a
bly occupying N/2 spatial orbitals. If the real system TAS experiment, the objective is to obtain the response
is irradiated with an external field characterized by a of the excited states that are visited by the system as it
scalar potential v((cid:126)r,t) (the extension to vector poten- is driven by the pump pulse (i.e. the response function
tials is also possible), the “time-dependent Kohn-Sham” of a system out of equilibrium). This extension will be
(TDKS) equations that characterize the evolution of the treated in Section IIA.
fictitioussystemare(atomicunitswillbeusedhereafter): Likewise, TDDFT can be used to compute strong
field non-linear photo-electron spectra of atoms and
∂ 1
i ϕ ((cid:126)r,t)=− ∇2ϕ ((cid:126)r,t)+v [n]((cid:126)r,t)ϕ ((cid:126)r,t), (1) molecules, for example with the method recently devel-
∂t i 2 i KS i oped by some of us.40 These spectra, however, are also
N/2 characteristic of the ground state, although, as it will
(cid:88)
n((cid:126)r,t)=2 |ϕi((cid:126)r,t)|2. (2) be shown in Section IIB, the methodology can be easily
i=1 extended to tackle the pump-probe case (time-resolved
photoelectron spectroscopy19).
The time-dependent density n((cid:126)r,t) is the central object,
and is identical for the real and for the KS systems. The
KS potential v is a functional of this density, and is
KS
A. Attosecond transient absorption spectroscopy
defined as:
v [n]((cid:126)r,t)=v ((cid:126)r)+v((cid:126)r,t)+v [n]((cid:126)r,t)+v [n]((cid:126)r,t), When an electromagnetic pulse passes through a gas
KS 0 H xc
(3) sample, the molecules polarize, and this polarization
wheretheHartreepotentialv correspondstoaclassical modifies the otherwise free propagation of light – one
H
electrostatic term of the consequences being the partial absorption of it.
In a dilute gas, assuming the electric dipole approxima-
(cid:90) n((cid:126)r,t)
v [n]((cid:126)r,t)= d3r(cid:48) , (4) tion and a sufficiently weak pulse, the dipole-dipole lin-
H |(cid:126)r(cid:48)−(cid:126)r|
ear dynamical polarizability entirely determines the po-
larization of the medium, and therefore the amount of
v ((cid:126)r) is the static external potential that characterizes
0 absorption. This is usually understood at equilibrium:
the system in its ground state (in a molecule, originated
the gas is formed by molecules at thermal equilibrium
by a set of nuclei), and the “exchange and correlation”
(perhaps at sufficiently low temperature so that they all
potential is v [n]. The exchange correlation potential is
xc can be considered to be at their ground state), and the
also a functional of the density and accounts for all the
onlylightpulsepresentisthatwhoseabsorptionwewant
intricatemany-electroneffects. Itisinpracticeunknown
to measure.
and must be approximated.34,39
In the pump-probe situation discussed here, however,
TheTDKSequationscanbeutilized,eitherdirectlyor
one wishes to compute the absorption of a probe pulse
inappropriatelytransformedmanner,tocomputethere-
by a set of molecules that is also irradiated by a pump,
sponseofamany-electronsystemtoaperturbation,weak
either simultaneously or with a given delay. The task is
or strong. In the perturbative regime, ideally one wishes
therefore to compute the response of the electric dipole
toobtaintheresponsefunctions[(hyper)-polarizabilities,
with and without the probe pulse – the difference being
optical and magnetic susceptibilities, ...], since (i) these
the excess of polarization, responsible for the absorption
objects then permit to predict any reaction in the ap-
oftheprobe. Wewillassume,asitisoftenthecase,that
propriate order, and (ii) experiments typically provide
thepumppulseisintense,whereastheprobeisweakand
spectra that are directly related to the response func-
can be treated in first order perturbation theory.
tions – e.g. the optical absorption cross section of a gas
Thissituationisamenabletoageneralizeddefinitionof
isproportionaltotheimaginarypartofthedipole-dipole
responsefunctions,suchastheonegivenintheappendix
molecular polarizability. In contrast, in the strong-field
ofRef.41andalsodiscussedindetailinRef.42. Wewill
regime, where perturbative treatments become cumber-
reviewherethisdefinition,adaptingittothepump-probe
some, one normally computes the particular response of
situation. Let us depart from a Hamiltonian in the form
the system to the perturbation of interest by directly
(we only treat the electric part neglecting the magnetic
propagating the TDKS equations in real time.
term of the electromagnetic field):
The vast majority of TDDFT applications have ad-
dressed the first-order response of the ground-state sys- Hˆ [E](t)=Hˆ+E(t)Vˆ (5)
0
tem to weak electric fields – which can provide the ab-
sorption spectrum, the optically-allowed excitation ener- where Hˆ is the static Hamiltonian that defines the sys-
gies and oscillator strengths, etc. This can be performed temitself,andE(t)Vˆ isthecouplingtothe“pump”laser
by linearizing the TDKS equations in the frequency do- pulse. This is the “unperturbed” Hamiltonian, that con-
main and casting the result into matrix-eigenvalue form, tains only the pump pulse; the full Hamiltonian results
or by propagating the same equations in real time ap- of the addition of the probe pulse f(t)Vˆ:
plying a sufficiently weak dipole perturbation. In any
case,theresponsefunctioncomputedinthismannerwill Hˆ(t)=Hˆ+E(t)Vˆ +f(t)Vˆ . (6)
4
The evolution of the system is given by: Hˆ) and of the pump shape E, as opposed to the con-
ventional response functions, which are only system de-
i∂ ρˆ(t)=(cid:104)Hˆ(t),ρˆ(t)(cid:105) , (7) pendent. Note also its dependence on two times t and
∂t t(cid:48), which cannot be reduced to only one by making use
of time-translational invariance, as it is customary when
andinitially(t=t ,sometimebeforethearrivalofboth
0 working at equilibrium.
pump or probe), the system is at equilibrium:
The response itself, δA(t), will be a functional of both
(cid:104) (cid:105) pump and probe pulses. If we take its Fourier transform
Hˆ,ρˆeq =0. (8) we may write it as:
(cid:90) ∞
ForafixedpumpE,wemayassumethesystemevolution δA[E,f](ω)= dt(cid:48)f(t(cid:48))χ [E](ω,t(cid:48)). (16)
Aˆ,Vˆ
to be a functional of the probe shape: ρˆ= ρˆ[f], and we −∞
may expand ρˆin a Taylor series (in the functional sense)
In order to compute the response function, one can use
around f =0:
as a probe a delta perturbation, i.e. f(t(cid:48)) = λδ(t(cid:48) −τ),
∞ which permits to identify:
(cid:88)
ρˆ[f]= ρˆ [f], (9)
n 1
n=0 χAˆ,Vˆ[E](ω,τ)= λδA[E,λδτ](ω). (17)
where ρˆ is the unperturbed system evolution (i.e. the
0 The action of such a delta-perturbation applied at the
evolution of the system in the presence of the pump
instant τ on the system is given by:
only: f = 0), and ρˆ is n-th order in the perturbation:
n
ρˆ[λf] = λnρˆ[f]. The system response to this perturba- |Φ(t→τ+)(cid:105)=e−iλVˆ|Φ(τ)(cid:105). (18)
tion is measured in terms of the expectation value of an
observable Aˆ, which can likewise be expanded: Fromnowonwewillrestrictthediscussiontopurestates,
since ensembles are not needed for the results that will
∞
A(t)=Tr(cid:104)ρˆ(t)Aˆ(cid:105)= (cid:88)A (t), (10) be shown below. However it can easily be extended to
n general ensembles. We also restrict the discussion to a
n=0 specific response function: the dipole-dipole polarizabil-
ity α[E](t,t(cid:48)) = χ [E](t,t(cid:48)), where both Aˆ and Vˆ are
where Dˆ,Dˆ
the atomic or molecular dipole operator Dˆ– taking into
(cid:104) (cid:105)
A (t)=Tr ρˆ (t)Aˆ . (11) accountthat,forthefrequenciesthatwearedealingwith,
n n
the dipole of interest is that of the electrons, and the
Forsufficientlyweakprobes,weareonlyinterestedinthe clamped nuclei approximation can be used. Moreover,
first term: wechoosetoworkwithlightpolarizedinthexdirection,
so that:
(cid:104) (cid:105)
δA(t)=A(t)−A (t)≈A (t)=Tr ρˆ (t)Aˆ , (12)
0 1 1 N
Dˆ =−(cid:88)xˆ , (19)
i
whichislinearlyrelatedtof throughapump-dependent
i=1
response function:
where N is the number of electrons. The expectation
(cid:90) ∞ value of this electronic dipole is an explicit functional of
A1(t)= dt(cid:48) f(t(cid:48))χAˆ,Vˆ[E](t,t(cid:48)). (13) the time-dependent density, and so is its variation:
−∞
(cid:90)
The response function is given by:41 δD[E,f](ω)=− d3r δn((cid:126)r,ω)x, (20)
(cid:110) (cid:104) (cid:105)(cid:111)
χ [E](t,t(cid:48))=iθ(t−t(cid:48))Tr ρˆ Aˆ [E](t),Vˆ [E](t(cid:48)) where δn((cid:126)r,ω) is the Fourier-transformed difference be-
Aˆ,Vˆ eq H H
tweentheelectronicdensitiesobtainedwithandwithout
(14)
the probe pulse. This straightforward formula in terms
Inside the commutator, the operators appear in the
of the density is what makes TDDFT specially suited
Heisenberg representation:
for these computations: we may safely utilize the Kohn-
Oˆ [E](t)=Uˆ[E](t ,t)Oˆ Uˆ[E](t,t ), (15) Shamsystemofnon-interactingelectrons. Thedeltaper-
H 0 0
turbation, Eq. (18), must be applied to each one of the
where Uˆ[E] is the time propagation operator in the pres- Kohn-Sham orbitals, and takes the following form:
enceofthepumponly–hencethefunctionaldependence ϕ((cid:126)r,t→τ+)=eiλxϕ ((cid:126)r,τ). (21)
i
on E. We keep this functional dependence explicit in the
notation for χAˆ,Vˆ[E](t,t(cid:48)), to stress that it is a property The absorption of a particular probe pulse f, is deter-
of both the system (defined by the static Hamiltonian minedbytheinducedpolarization,givenbyδD[E,f](ω).
5
Wewillcomputethedynamicalpolarizabilityα[E](ω,τ), (a) 1 M(r)
whichisthepolarizationinducedbyadeltaperturbation,
i.e.:
1
α[E](ω,τ)= δD[E,λδ ](ω), (22)
λ τ
C
since it allows to compute any particular response
through the integration of Eq. (16), as long as the probe
(b)
isweak. Inparticular,wewilllookattheimaginarypart
ofα[E](ω,τ),whichisthepartresponsibleforabsorption.
Note, finally, that in a 3D situation the polarizability is
not a scalar but a tensor, since there are three possible A
light polarization directions, and three components for
the system dipole moment. In most cases, one is inter-
ested in the trace of this tensor, an averaged quantity
that corresponds to the absorption of a randomly ori-
ented sample of molecules.
B
B. Time-resolved photoelectron spectroscopy
The photoelectron spectra presented in this work are
FIG. 1. Schematic description of the space partition imple-
produced within TDDFT using the recently introduced
mented by the mask method. A mask function (a) is used
Mask Method.40 This method is based on a geometrical to implement the spatial partitions (b). In region A (inter-
partitioningandamixedreal-andmomentum-spacetime action region) the TDKS equations are numerically solved in
evolution scheme.43 In the following, we summarize the real space while in B (free propagation region) electrons are
maintraitsofthetechnique(wereferthereadertoRef.40 evolvedanalyticallyasfreeparticlesinmomentumspace. Re-
for a complete description), and demonstrate how it can gion C is where ϕA and ϕB overlap.
be straightforwardly applied to the non-equilibrium sit-
uation required by pump-probe experiments.
(b). The space is divided into two regions, A and B;
In photoemission processes a light source focused on
the inner region A, containing the system with enough
a sample transfers energy to the system. Depending on
empty space around, is where electrons are allowed to
the light intensity electrons can absorb one or more pho-
interact with each other and with the system, and re-
tonsandescapefromthesampleduetothephotoelectric
gion B, defined as the complement of A, is where elec-
effect. In experiments, electrons are detected and their
trons are non-interacting and freely propagating. Every
momentum is measured. By repeating measurements on
KS orbital ϕ (r) can be decomposed accordingly with
similarly prepared samples it is possible to estimate the i
ϕ (r) = ϕA(r)+ϕB(r), so that ϕA(r) resides mainly in
probability to measure an electron with a given momen- i i i i
region A and ϕB(r) mainly in region B.
tum. From a computational point of view, the descrip- i
Thegeometricalpartitionisimplementedbyasmooth
tion of such processes for complex systems is a challeng-
maskfunctionM(r)definedtobeonedeepintheinterior
ingproblem. Themaindifficultyarisesfromthenecessity
of A and zero outside (see Fig. 1 (a)):
of describing properly electrons in the continuum.
In typical experimental setups, detectors are situated ϕA(r)=M(r)ϕ (r), (23)
i i
far away from the sample and electrons overcoming the ϕB(r)=(1−M(r))ϕ (r). (24)
i i
ionizationbarriertravelalongwaybeforebeingdetected.
The distances, that electrons travel are usually orders Themaskfunctiontakescareoftheboundaryconditions
of magnitude larger than the typical interaction length in A by forcing every function to be zero at the border.
scales in the sample. During their journey towards the Inordertogiveagooddescriptionoffunctionsextending
detector,andfarawayfromtheparentsystem,theyprac- over the whole space it is convenient to represent ϕBi (r)
ticallyevolveasfreeparticlesdrivenbyanexternalfield. orbitals in momentum space ϕ˜Bi (p).
The solution of the Schr¨odinger equation for free elec- A mixed real and momentum-space time evolution
tronsinatimedependentexternalfieldisknownanalyt- schemecanthenbeeasilyderivedfollowingthegeometri-
ically in terms of plane waves as Volkov states. It seems cal splitting. Given a set of orbitals at time t their value
therefore a waste of resources to solve the Schr¨odinger at a successive time t+∆t is provided by
equationnumericallyinthewholespaceifaconsiderable (cid:40)
ϕA(r,t+∆t)=M(r)e−iHˆ∆tϕA(r,t)
part of the wave function can be described analytically. i i
Inordertotakeadvantageofthepreviousobservations ϕ˜Bi (p,t+∆t)=e−i(p−A2(t))2∆tϕBi (p,t)+ϕ˜Ai (p,t+∆t)
we partition the space according to the scheme in Fig. 1 (25)
6
withHˆ beingtheeffectivesingle-particleTDDFTHamil- (PAD), can easily be evaluated by expressing P(p) in
tonian, A(t) the total external time dependent vector polar coordinates with respect to a given azimuth axis.
potential (the coupling with the external field is conve- It is noteworthy that during the evolution defined in
niently expressed in the velocity gauge), and Eq. 25 the part of the density contained in A transferred
1 (cid:90) toB isnotallowedtoreturn. Clearly, incaseswherethe
ϕ˜A(p,t+∆t)= dr(1−M(r))e−iHˆ∆tϕA(r,t)eip·r, externalfieldisstrongenoughtoproduceelectronorbits
i (2π)3/2 i
crossing the boundary of A and backscattering to the
(26)
core the mask method provides a poor approximation.
constituting the portion of electrons leaving the system
InthesecasesabiggerregionAoramorerefinedscheme
attimet+∆t. Ateachiterationintheevolutiontheout-
must be employed40. The laser fields employed in this
goingcomponentsofϕA(r)aresuppressedintheinterac-
i work are weak enough that we can safely assume region
tion region by the multiplication with M(r) while being
Atoalwaysbesufficientlylargetocontainalltherelevant
collected as plane waves in ϕ˜B(p) via ϕ˜A(p). The re-
i i electron trajectories.
sultingmomentumspacewavefunctionsarethenevolved
analytically simply by a phase multiplication.
Theadvantageofusingsuchanapproachresidesinthe
fact that we can conveniently store the wavefunctions III. RESULTS
on a spatial grid inside A while treating wavefunctions
in B (and therefore the tails extending to infinity) as The theory described in the previous section has been
free-electrons in momentum space. Moreover the mask implemented in the octopus code44; we refer the reader
function introduces a region C, where the wavefunctions toRefs.45and46fortheessentialpointsofthenumerical
inAandB overlap(seeFig.1)andthatactsasmatching methodology.
layer. Inspiteofthefactthat,fromatheoreticalpointof In order to simplify the illustration of the results a
view,thematchingbetweeninnerandouterregioncould clamped ion approximation has been used in the cal-
be performed on a single surface, from a numerical point culations for the molecular case. Further inclusion of
of view, having a whole region to perform the matching theionicmotioncouldbedoneatthesemi-classicallevel
is more stable and less influenced by different choices of with Ehrenfest dynamics45,47 – already implemented in
spatial grids. the code – and without any essential modification to the
FromthemomentumcomponentsoftheorbitalsinBit theory presented.
ispossibletoevaluatethemomentum-resolvedphotoelec-
tron probability distribution as a sum over the occupied
orbitals
A. One-dimensional model Helium
occ
(cid:88)
P(p)≈ lim |ϕ˜B(p,t)|2, (27)
t→∞ i As first example we study the absorption spectrum of
i=1
an excited one-dimensional soft-Coulomb Helium atom.
thelimitt→∞ensuringthatalltheionizedcomponents This is an exactly solvable model that provides a useful
are collected. This scheme is entirely non-perturbative; benchmarktotestdifferentapproximations. Wewillfirst
in a pump-probe setup, it does not assume linearity in discuss the exact solution, and later apply TDDFT. A
either pump or probe. Therefore, it can be applied in more realistic 3D model will be presented in the next
the same manner when two pulses are present than with section.
one pulse only, as it was shown in Ref. 40. Like in The 1D model of the Helium atom is defined by the
Sec.IIAwecangeneralizethepreviousderivationtoad- following Hamiltonian:
dress transient photoelectron spectroscopy (spin-, angle-
and energy-resolved) in practice by employing a pump- H(t)=T +V (t)+V , (29)
ext ee
probeschemeandperformingnumericalsimulationswith
twotimedelayedexternalpulses. ATRPESmapisthen where
generatedbyperformingacomputationforeachdifferent
time delay. 2 2
V =− − +E(t)(x +x ) (30)
From P(p) several relevant quantities can be calcu- ext (cid:112)1+x2 (cid:112)1+x2 1 2
1 2
lated. The energy-resolved photoelectron probability
P(E), usually referred to as photoelectron spectrum is the external potential: the electron Coulomb interac-
√
(PES), can be obtained by integrating P(p) over solid tion 1/|x| is softened to 1/ 1+x2. The coupling with
angles the external time-dependent field E(t) is expressed in
(cid:90) length gauge, and electrons are confined to move along
P(E =p2/2)= dΩ P(p). (28) the x direction only. Finally,
p
4π
The angular- and energy-resolved photoelectron proba- 1(cid:18) ∂2 ∂2 (cid:19)
T =− + (31)
bility P(θ,φ,E), or photoelectron angular distribution 2 ∂x2 ∂x2
1 2
7
30
20
a.u. 16 20 LDA
on/ 12 a.u.
cti n/
Cross-se 48 0.59 ss-sectio 10 EXX
o
0 00..5577 Cr

00..5555 Exact
0 0
00..11
00..22 00..33 00..5533
00..44
00..55
00.66 0.7 0.51 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
FIG. 2. Out of equilibrium absorption spectrum as function
ofthepumplaserfrequencyforone-dimensionalHelium. The FIG. 3. Comparison of absorption spectra calculated in dif-
system is driven out of equilibrium by 45 cycle sin2 envelope ferent approximations for a one-dimensional Helium model.
laser pulses of intensity I = 5.26×1011 W/cm2, at differ- Thefilledcurvesarethespectrafortheunperturbedsystems
ent carrier frequencies and then probed right after. Maximal whilethesolidlinesarethespectraofthesystemexcitedbya
response is observed for frequencies close to the first optical laserasinFig.2resonantwiththefirstallowedopticaltransi-
transition ωP =0.533 a.u.. tion: exacttime-dependentSchro¨dingerequationωP =0.533
a.u. (in blue), EXX ω = 0.549 a.u. (in red), and LDA
P
ω = 0.475 a.u. (in green). The dashed blue line is the
P
isthekineticenergy,andtheelectron-electroninteraction absorption of the system perturbed by a 180 cycle laser and
is probed at t = 30.62 fs, where the population on the excited
state is maximal. The lines have been shifted by a vertical
1
V = . (32) constant to facilitate the comparison between results.
ee (cid:112)
1+(x −x )2
1 2
This model is numerically solvable given the exact map-
pingdiscussedinRefs.48and49,whereitisprovedthat In orderto analyze thispoint further, a cutat the res-
the many-body problem of N electrons in one dimension onant frequency is displayed in the lower (blue) curves
isequivalenttothatofoneelectroninN dimensions. The of Fig. 3. The filled curve represents the spectrum ob-
wave functions and other necessary functions are repre- tainedfromthesysteminitsgroundstatewhilethesolid
sented on a real space regular grid; a squared (linear for line corresponds to the spectra of the system excited by
one-dimensional TDDFT) box of size L = 200 a.u. and a laser pulse with a frequency resonant with the first
spacing ∆x = 0.2 a.u. has been employed in all the cal- optical transition and probed after the perturbation. By
culations. directcomparisonofthetwospectraitiseasytodiscrim-
In order to illustrate how an external field can modify inatethepeaksassociatedwiththeground-stateabsorp-
theopticalpropertiesofasysteminFig.2weshowascan tion from the ones characterizing the absorption from
ofthenon-equilibriumabsorptionspectrumgeneratedby the excited state. In particular, the peaks related to
a 45 cycle sin2 envelope pulse with intensity I = 5.26× the ground-state absorption are located at the energies
1011 W/cm2 atdifferentangularfrequencies0.51[a.u.]≤ ω0→1 =(cid:15)1−(cid:15)0 =0.533 a.u. corresponding to the transi-
ωP ≤ 0.59 [a.u.]. The plot displays Im[α[E,τ](ω,ωP)], tion from the ground ((cid:15)0 = −2.238 a.u.) to the first ex-
choosing τ right at the end of the pump pulse. citedstate((cid:15)1 =−1.705)andω0→3 =(cid:15)3−(cid:15)0 =0.672a.u.
It can be seen how the absorption around the first ex- corresponding to the third excited state ((cid:15)3 = −1.566
citationfrequency0.533a.u. isstronglydiminishedwhen a.u.) – the direct excitation of the second excited state
thefrequencyofthepumpresonateswiththatfrequency. is forbidden by symmetry. The solid curves show finger-
Inthatsituation,anabsorptionpeakalsoappearsaround printsofthepopulationofthefirstexcitedstate,namely
the excitation frequency corresponding to the transition thepeakscorrespondingtotransitionsfromthatfirstex-
from the first to the second excited state, in our case at citedstatetoothers: inparticular,thepeakappearingat
0.076 a.u. This behavior is a direct consequence of the the low energy ω1→2 =(cid:15)2−(cid:15)1 =0.076 a.u. is associated
fact that the laser is pumping the system to the first ex- with the transition from the first excited state (cid:15)1 to the
cited state and this process is more efficient for a field second one (cid:15)2 =−1.629 a.u.
tuned with the excitation energy. The absorption spec- These spectra contain information that is not con-
trum is therefore a mixture of the one corresponding to tained in the equilibrium ones. For example, let us con-
theunperturbedgroundstateandthatofthefirstexcited sider the spectra that would be produced by each sin-
state. gle eigenstate, given by the state-dependent dynamical
8
polarizabilities, which may be written in the sum-over- 1
states form as:
(cid:34)
α(i)(ω)=(cid:88) |(cid:104)Ψi|Dˆ|Ψj(cid:105)|2 (33)
ω−((cid:15) −(cid:15) )+i0+
j i
j(cid:54)=i n
o
− |(cid:104)Ψi|Dˆ|Ψj(cid:105)|2 (cid:35) . pulati 45 cycles 180 cycles
ω+((cid:15) −(cid:15) )+i0+ o
j i P
The poles of this function provide us with the eigenvalue
differences(cid:15) −(cid:15) ;ifthisvalueispositive,thecorrespond-
j i
ing term is associated with a photon-absorption process;
if it is negative, with a stimulated emission term. The
0
weightassociatedwitheachoneofthesepolesprovideus 0 20 40 60
t/fs
with the dipole coupling matrix elements (cid:104)Ψ |Dˆ|Ψ (cid:105).
i j
During the time evolution, the wavefunction can be FIG. 4. Exact population on the ground |η (t)|2 =
0
expanded on the basis of eigenstates of the unperturbed |(cid:104)Ψ |Ψ(t)(cid:105)|2 (solid lines) and first excited |η (t)|2 =
(cid:80) 0 1
systemΨ(t)= iηi(t)Ψi. Whenthesystemisprobedat |(cid:104)Ψ1|Ψ(t)(cid:105)|2 (dashed lines) states as a function of time for
acertaintimet,theresultingspectrumcanbethoughtas different laser pulses. In red a 45 cycles pulse with parame-
alinearcombinationofthespectraproducedbyeachsin- ters as in Fig. 3, and in blue a longer 180 cycles pulse with
gleeigenstate. Ananalysisofthetransientspectrummay the same parameters.
therefore provide information about the mixing weights
η , and about excitation energies and dipole couplings
i
between excited states – information that is absent in A Rabi oscillation is a fluctuation behavior of states
equilibrium ground-state linear response. occupationoccurringduetotheinteractionofanoscilla-
In our case, we find by direct projection of the tory optical field in resonance with a two-level system.50
time-dependent wave function onto the eigenstates, that The occupation probability alternates with the Rabi fre-
the system after the pulse is composed mainly of the quency Ω(t)=f(t)µ0→1, where µ0→1 is the dipole tran-
ground and the first excited state with weights, |η |2 = sition matrix element between the states and f(t) is the
0
|(cid:104)Ψ |Ψ(t)(cid:105)|2 =0.7120 and |η |2 =|(cid:104)Ψ |Ψ(t)(cid:105)|2 =0.2876. electric field envelope. Extremal points of the popula-
0 1 1
The same information can be recovered by comparing tions should be located at times where the pulse area
(cid:82)t
the perturbed and the unperturbed spectrum at E0→1. Θ(t) = −∞dτf(t)µ0→1 is an integer multiple of π,
At this energy we only have the contribution coming Θ(t) = nπ. With the numerically calculated matrix el-
from Ψ → Ψ and its inverse Ψ → Ψ . The peak ement µ = 1.11 a.u. the first maximal population
0 1 1 0 0→1
height of the perturbed spectrum after the laser pulse of the excited state is expected at t = 30.65 fs, in good
h is therefore a combination of the heights associated agreement with what is observed. The absorption spec-
t
with the ground h and the excited h states: h = trum at this time, as shown in Fig 3 (dashed blue line),
0 1 t
|η |2h +|η |2h . At this energy h = −h due to the displays a considerable enhancement at ω and a neg-
0 0 1 1 0 1 1→2
nature of the transition 1 → 0 the ratio α = h /h = ative emission peak at ω as expected from a pure
t 0 0→1
|η |2 − |η |2 = 0.4258 thus gives direct information excited state.
0 1
on the difference of the mixing weights. Complement- ItisinterestingtostudythesamemodelwithTDDFT
ing this information with a two-level system assumption insteadofwithanexacttreatmentinordertoaddressthe
|η |2 +|η |2 = 1 we obtain |η |2 = (1+α)/2 = 0.7129 performance of available (mainly static) xc-functionals.
0 1 0
and |η |2 = (1−α)/2 = 0.2871 in good agreement with In Fig. 3 we display results obtained with TDDFT, em-
1
the results calculated by direct projection of the wave- ployingtwodifferentexchange-correlation(xc)functional
function. approximations: exact exchange51 (EXX) in red, and
In Fig. 4 we display the population weights for two one-dimensional local density approximation52 (LDA) in
different laser pulses. The red lines correspond to the green. The calculations were performed in the adiabatic
same laser pulse as in Fig. 3 while the blue lines pertain approximation using the same parameters as in the ex-
to a four times longer laser with the same parameters act case. The laser frequency was tuned to match the
(intensity,envelopeshapeandcarrierfrequency)and180 first optical transition appearing at: ω =0.549 a.u. for
P
opticalcycles. Forbothlasersthepopulationsofboththe EXX, and ω =0.475 a.u. for LDA.
P
ground and first excited state at each time almost sum The unperturbed spectrum (solid curve) provided by
toone,indicatinganessentialtwo-leveldynamics. Inthe EXX is in good agreement with the exact calculation,
caseofthelongpulseweobserveamaximum(minimum) and the perturbed one qualitatively reproduces the ex-
of the population over the excited (ground) state at t = act result. In particular the new peak appearing at low
30.62 fs. This behavior can be understood in terms of energyassociatedwiththetransition1→2iswellrepre-
Rabi oscillations. sented. In contrast LDA is only capable of reproducing
9
2 2.5
a.u. 0.8
/
n
Cross-sectio 1 00..46 s-section/a.u.
s
o
r
C
0.2
0
0
0 0.2 0.4 0.6 0.8 0
0 2 4 6 8 10
FIG.5. Comparisonoftheabsorptionspectraofunperturbed FIG. 6. Helium transient absorption spectrum scan for dif-
(filled curve) and perturbed He atom probed at τ = 5.39 ferenttimedelays. Pumplaserpicturedintheupperpanelis
(solidline)andaftertheendofthepulseτ =8.68fs(dashed the same as in Fig. 5.
line). The spectrum range is below the ionization threshold.
The atom is excited by a 45 cycle sin2 envelope laser pulse
polarized along the x-axis with carrier ω = 0.79 a.u. reso-
P
nantwiththefirstopticaltransition,intensityI =2.6×1012 trumfortheunperturbedatom(filledcurve)andtheper-
W/cm2. turbed one probed with a delta perturbation right after
the pump pulse at τ = 8.68 fs (dashed line). The com-
parisonpresentsmanytraitssimilartotheonesdiscussed
one peak both for the perturbed and unperturbed cases. in Sec. IIIA for the one-dimensional Helium model. In
Thisisduetotheknownproblemofasymptoticexponen- particular, fingerprints of the population of the first ex-
tial decay of the functional that in this one-dimensional cited state can be observed in the appearance of a peak
example supports only a single bound excited state. in the gap at ω2p→3s = 0.079 a.u. associated with the
A common feature of both approximations is consti- transition 1s2p → 1s3s. The second peak, associated
tutedbythepresenceofnegativevaluesintheperturbed with the transition 1s2 → 1s2p, presents height changes
spectra. Thiscanbetrackeddowntothelackofmemory correlatedwiththeformerone. Wealsoobtainthesmall
intheadiabaticxc-functionalapproximation53. Thelack artifacts,suchastheenergyshiftsandthenegativevalues
orwrongmemorydependenceinthefunctionalresultsin attributedtothexc-kernelmemorydependencediscussed
slightly displaced absorption and emission peaks associ- previously.
atedwiththesametransition. Thisfact, analyzedinthe Additional detailson the excitationprocess can beac-
lightofEq.33,results,atthetransitionenergy,inasum quired by expanding the time dimension of the absorp-
of two Lorentzian curves with different sign and slightly tion spectrum. In Fig. 6 the time resolved absorption
different centers. This explains why we get two inverted spectrum (TAS) map is displayed. The map was pro-
peakswhereweshouldhaveonlyasingleonegoingfrom duced by probing the system at different time delays.
positive to negative strength as we populate the excited As the delay is increased we observe the build-up of the
state – as shown by the exact (blue) curves in Fig. 3. peakassociatedwiththestatebeingpumpedbythelaser
pulse at ω . This changes are reflected in the oscil-
2p→3s
lations of the ground state first optical peak at ω .
1s→2p
B. Helium atom in 3D In TDDFT the knowledge of the wavefunction is lost in
favorofthedensity,whichdoesnotallowustodoapop-
In this section we study the real Helium atom. We ulationanalysisbasedonsimplewavefunctionprojection
employed the EXX functional and discretized TDDFT Thetransientabsorptionspectrum,ontheotherhand,is
equations on a spherical box of radius R=14 a.u., spac- an explicit density functional, and its computation with
ing ∆r =0.4 a.u. and absorbing boundaries 2 a.u. wide. TDDFT may help us to understand the evolution of the
We begin by investigating the changes in absorption state populations.
of He under the influence of an external UV laser field The peak appearing in the gap presents a maximum
drivingthe system withthe frequencyof thefirst dipole- at τ = 5.39 fs that emerges before the end of the pump
allowed excitation. To this end we used a 45 cycle sin2 pulse (τ = 8.68 fs). This peak is associated only with
laser pulse in velocity gauge with carrier ω = 0.79 a.u. the transition from the 1s2p → 1s3s and therefore its
P
resonant with the 1s2 → 1s2p transition, of intensity heightisproportionaltothe2pexcitedstatepopulation.
I = 2.6 × 1012 W/cm2 polarized along the x-axis. In The oscillation can then be interpreted in terms of Rabi
Fig. 5 we show a comparison of the absorption spec- physics as discussed in Sec. IIIA.
10
180 1
(a) (b)
1.8 -2
1.4
)
u. E) 0 0
E/a. P(( (c) (d) 1
1 log -3
E))
P(
g(
0.6 lo
-9
-1 1 3 5 7 9
-7 0
0.6 1.7 0 360
E/a.u.
FIG.7. Heliumtransientphotoelectronspectruminlogarith-
mic scale. The pump laser (upper panel) is the same as in
FIG. 8. Energy- and angular- resolved photoelectron spec-
Fig. 5 and the probe is a 40 cycles trapezoidal laser pulse
tra for Helium at fixed delay τ = 8.99 fs. Panel (c) dis-
with 8 cycles ramp, ω = 1.8 a.u., I = 5.4×109 W/cm2
p plays a logarithmic scale PES P(E) comparison at fixed de-
aligned with the pump pulse.
lays τ = −1.69 fs (red) and τ = 8.99 fs (green). The other
panels depict normalized PADs P(θ,φ,E) with polar coordi-
natesreferredtoaxisz atfixeddelayτ =8.99fsandenergy:
Further insight can be achieved investigating the pho- (a)E =0.66a.u.,(b)E =0.88a.u.,and(d)E =1.67a.u..
1 2 3
toemission properties of the system. In Fig. 7 we show Whitecrossesmarktheintersectionbetweentheprobepolar-
the TRPES map, as calculated in a pump-probe set up. ization axis and the cutting sphere.
Photoelectrons are calculated with the technique out-
lined in Sec. IIB. The pump pulse is the same as the
one employed for TAS. The probe is a 40 cycles trape- begins to emerge. This peak corresponds to emission
zoidal laser pulse (8 cycles ramp) with carrier frequency from the pump-excited 2p state E = ω + ω − I .
3 P p P
ω = 1.8 a.u., intensity I = 5.4 × 109 W/cm2, polar- It is a process where the atom, initially in the ground
p
ized along the x-axis and is weak enough to discard non- state, absorbs a photon from the pump and gets excited
linear effects. We performed a scan for different time to the 2p bound state. The subsequent absorption of a
delays, measuring each delay as the difference from the probe photon frees the electron into the continuum. The
probe center to the beginning of the pump. Negative peakatE isunderstoodintermsofpumpphotonsonly:
1
delays correspond to the situation where the probe pre- E = 2ω −I . The ionization mechanism shares the
3 P P
cedes the pump. Moreover, in order to include all the first step with the E process, namely a 2s→2p excita-
3
relevant trajectories a spherical box of R = 30 a.u. was tion produced by the absorption of a ω photon. In the
P
employed, and photoelectrons were recorded only during second step the electron is directly excited to a contin-
the up-time of the probe pulse. uum state by the absorption of a second ω photon. In
P
The TRPES map in Fig 7 shows three main features thelinearregime,thedirectphotoionizationcross-section
at E = 0.66 .a.u, E = 0.88 a.u. and E = 1.67 a.u.. decays exponentially with energy54. For this reason and
1 2 3
In our case the probe pulse is weak, and photoelectrons duetothedisparityinintensitybetweenpumpandprobe
escaping the system undergo photoelectric-effect energy this ionization channel is by far the most favorable one.
conservation. A bound electron can absorb a single pho- In direct photoemission processes, the photoelectron
ton and escape from the atom with a maximum kinetic angular distribution (PAD) contains information about
energy E = ω −I , where ω is the probe carrier fre- the electronic configuration of the ionized state.55 In or-
p P p
quency and I is the field-free ionization energy. The der to support the energetic arguments PADS P(θ,φ,E)
P
ionizationpotentialcanbeevaluatedinDFTastheneg- atτ =8.99fsarepresentedinFig8(a),(b),(d)together
ativeenergyofthehighestoccupiedKSorbital(HOMO) withcutsonTRPESmapatτ =0fsandτ =8.99fs(b).
I = (cid:15) = 0.92 a.u. Thus the peak appearing at E For each energy marked in Fig. 8 (c) we perform spheri-
P 2s 2
is energetically compatible with photoelectrons emitted calcutsofthephotoemissionprobabilityonenergyshells
from the 2s level: E2 =ωp−IP. Consistently, this peak at E =E1, E2, E3. Each cut is then plotted in polar co-
is the only one appearing at negative delays where the ordinates with θ being the angle from the z-axis and φ
pulses do not overlap. Moreover, the peak strength is the angle in the xy-plane measured from the x-axis. In-
weakly varying with the delay while slightly shifting to- tersection of the lasers polarization axis with the sphere
wards lower values around 3 fs in accordance with TAS are marked with a white cross.
findings. At about the same delay time the peak at E Fig 8 shows clearly that photoelectrons at E (a) and
3 1