Theoretical Analysis and Numerical Simulation of
Attosecond Time Delays in Photoionization
by
Jing Su
B.S., University of Science and Technology of China, 2008
M.S., University of Colorado at Boulder, 2013

A thesis submitted to the
Faculty of the Graduate School of the
University of Colorado in partial fulfillmentof the requirements for the degree ofDoctor of Philosophy
Department of Physics
2014This thesis entitled:
Theoretical Analysis and Numerical Simulation of Attosecond Time Delays in Photoionization
written by Jing Su
has been approved for the Department of Physics

Prof. Andreas Becker

Dr. Agnieszka Jaroń-Becker

Date

The final copy of this thesis has been examined by the signatories, and we find that both thecontent and the form meet acceptable presentation standards of scholarly work in the above
mentioned discipline.

Su, Jing (Ph.D., Physics)
Theoretical Analysis and Numerical Simulation of Attosecond Time Delays in Photoionization
Thesis directed by Prof. Andreas Becker

Recent developments in laser technology, in particular the advances in high-harmonic generation, enable the generation of ultrashort extreme ultraviolet (XUV) pulses with attosecond (1as = 10−18 s) duration. Such tools open the opportunity to study electron dynamics in atoms
and molecules on its intrinsic time scale. As an example, the attosecond streaking technique wasrecently applied to time resolve the photoionization process in atomic and solid systems. In thistechnique, an isolated attosecond XUV pulse, that ionizes the electron in the target system, is superimposed with a few-cycle streaking pulse (usually of near-infrared wavelengths). The streaking
pulse modulates the final momentum (or energy) of the photoelectron. The measured streaking
trace, i.e., the final momentum (or energy) as a function of the relative delay between these twopulses, contains time information of the photoionization process. By comparing two streaking
traces measured for photoionization from the 2s and 2p orbitals in a neon atom, Schultze et al.
[Science 328, 1658 (2010)] found a temporal offset of 21 ± 5 as between them and interpreted thisvalue as the time delay between photoionization from the 2s and 2p orbitals. This experiment has
initiated a debate among theoreticians, in particular about the origin of the measured time delay. A
correct interpretation of the delay is extremely important for our understanding of the attosecond
streaking technique and an exact analysis of time resolved measurements of this and other ultrafast
processes.
In this thesis we systematically study the attosecond time delays in photoionization using
numerical simulations. We first propose a new method, based on the fundamental definition of atime delay, to theoretically study the photoionization process induced by an XUV pulse from a timedependent perspective. We then turn to analyze the time delays measured in streaking experiments.
Our results show that for single-photon ionization the observed streaking time delay arises fromiv
the finite-range propagation of the photoelectron in the coupled field of the ionic potential and thestreaking pulse. Consequently, we conclude that the photon absorption occurs instantaneously at
the center of the XUV pulse, i.e., with no time delay. Our analysis further reveals that the streaking
time delay can be interpreted as a finite sum of piecewise field-free time delays weighted by the
relative instantaneous streaking field strength and provides itself as a useful tool for imaging thepresence of an additional potential located at a distance from the ionic core. We finally extend ourtime delay studies to the two-photon ionization process and show that the absorption time delay issignificantly different for nonresonant and resonant two-photon ionization. Our results imply thatthe absorption of two photons in the nonresonant case occurs instantaneously, without time delay,at the center of the XUV pulse. However, in the resonant scenario we find a substantial absorption
time delay that changes linearly with the duration of the XUV pulse. Our further theoreticalanalysis shows that this absorption time delay can be related to the phase acquired by the electron
during its transition from the initial ground state to the continuum.

Dedication

To my love and family

Acknowledgements

First and foremost, I would like to thank my advisor Prof. Andreas Becker for his continuous
support and help in the last five years. Without his guidance and encouragement I would have nochance to finish my Ph.D. work. A lot of times when I met challenges and problems in our projectand did not have confidence about it, he taught and encouraged me with his patience, motivation,enthusiasm, and immense knowledge. He is a very nice and easygoing person, however with very
rigorous attitude towards science. I have really learned a lot from him in these years, in particular,
the ability of critical thinking and the way of presenting results with logic, which will help me
throughout my life.
Next, I would like to thank Dr. Agnieszka Jaroń-Becker, who is like the co-advisor for mein the group. She has a very broad knowledge and also much experience of researches in both
Chemistry and Physics, which inspire and bring us many interesting ideas in our project. She isalso very warmhearted and always be ready to help others in both research as well as everyday life.
I also owe many of my thanks to some previous members in our group: Dr. Shaohao Chen,
Dr. Antonio Picón, and Dr. Norio Takemoto. When I was a freshman in the group they taught mea lot about numerical simulations, parallel computations, mathematical libraries, and many otherthings. Talking with them was really helpful and always a lot of fun. I still remember how manytimes I went to disturb Shaohao in one afternoon when I was learning the numerical library of ourgroup. They are also willing to share their experience of research and life with me and taught me
how to make a balance between them.
Working in Becker’s ultrafast group is really an amazing experience. The atmosphere in theviigroup is so good that every one has a smiling face almost every day. We discussed not only Physics
but also everything interesting. Thanks to these great and humorous people: Hongcheng Ni, Dr.
Carlos Hernández-Garcı́a, Daniel Weflen, Yuqing Xia, Andrew Spott, Michelle Miller, and Cory
Goldsmith.
I also would like to take this opportunity to thank the other members of my thesis committee:
Prof. Henry Kapteyn, Prof. Steven Cundiff, and Prof. David Jonas. It is my great honor to havethem sit on my committee. As experimentalists, they provided me many interesting ideas anduseful suggestions from a different perspective during the discussions in my Comps III exam.
The Janus supercomputer and the clusters at JILA have given me continuous computational
support during these years. In particular, the Janus supercomputer has enabled me to do a lot of
calculations simultaneously, which has greatly reduced the time of my Ph.D. career. So I am very
grateful to all the people working for the Janus supercomputer as well as the staff members of the
JILA computing group.
Finally, and most importantly, I would like to thank my family and all of my friends forsupporting and encouraging me in so many years. My parents and sister have given me theirendless support and unconditional love no matter where I am and whether I am successful or not.
A special thank goes to my girlfriend Ye Gu, who has been in a long-distance relationship with me
for almost ten years. Without her and my family I would never be able to overcome the difficulties
and finally accomplish my Ph.D. study.

Contents

Chapter
1 Introduction

1

1.1

Electron dynamics on the attosecond time scale . . . . . . . . . . . . . . . . . . . . .

1

1.2

Attosecond pulses generated through high harmonic generation . . . . . . . . . . . .

3

1.3

Attosecond streaking technique and its application to photoionization . . . . . . . .

6

1.4Time delay in photoionization and previous theoretical interpretations . . . . . . . .

9

1.5

Organization of this thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14

2 Numerical Model and Method
2.1

2.2

Single-active-electron model for atoms and molecules . . . . . . . . . . . . . . . . . . 17
2.1.1

Potential average method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18

2.1.2

Optimized effective potential method . . . . . . . . . . . . . . . . . . . . . . . 22

Numerical methods of solving time-dependent Schrödinger equation . . . . . . . . . 25
2.2.1

Propagation method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

2.2.2

Methods of obtaining initial states . . . . . . . . . . . . . . . . . . . . . . . . 28

2.2.3

Absorbing boundaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

3 Time Delays in Single-Photon Ionization
3.1

17

35

Calculation of Wigner-Smith-like time delays from a time-dependent perspective . . 36
3.1.1

Wigner-Smith time delay . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36

3.1.2

Theoretical method and back-propagation technique . . . . . . . . . . . . . . 40

ix

3.2

3.1.3

Application to single-photon ionization by an XUV pulse . . . . . . . . . . . 46

3.1.4

Comparison with results of other theoretical methods . . . . . . . . . . . . . 51

3.1.5

Dependence of time delay on XUV pulse parameters . . . . . . . . . . . . . . 53

3.1.6

Application to laser-assisted XUV photoionization . . . . . . . . . . . . . . . 55

Interpretation of time delays measured in streaking experiments . . . . . . . . . . . . 60
3.2.1

Extracting streaking time delays in numerical simulations . . . . . . . . . . . 62

3.2.2

Calculating propagation time delays using classical approaches . . . . . . . . 64

3.2.3

Physical interpretation of streaking time delay in single-photon ionization . . 71

3.2.4

Finite-range time delay and its application to imaging . . . . . . . . . . . . . 78

3.2.5

Influences of pedestal and additional static field . . . . . . . . . . . . . . . . . 85

4 Time Delays in Two-Photon Ionization
4.1

4.2

4.3

89

Time delays extracted in two-photon ionization using streak camera . . . . . . . . . 90
4.1.1

Streaking time delay . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90

4.1.2

Propagation time delay . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94

4.1.3

Absorption time delay . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96

Alternative studies of absorption time delays in two-photon ionization . . . . . . . . 98
4.2.1

Trajectory analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

4.2.2

Phase derivative analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104

Control of two-photon processes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107
4.3.1

Short introduction to coherent control . . . . . . . . . . . . . . . . . . . . . . 108

4.3.2

Theoretical models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109

4.3.3

Control of excitation to dissociative state . . . . . . . . . . . . . . . . . . . . 111

4.3.4

Control of excitation to bound states . . . . . . . . . . . . . . . . . . . . . . . 113

5 Conclusions

116

Bibliography

118

x

TablesTable
2.1

Comparison of eigenenergies (in atomic units) from different methods for neon atom.
The experimental values are from Ref. [60], the HF results are from Ref. [61], andall the other results are from our calculations. . . . . . . . . . . . . . . . . . . . . . . 19

3.1

Results of numerical calculations for the times tΨg ,R , tΨ(0) ,R and the time delayg∆tΨg ,R for different spatial steps δx and a fixed time step of δt = 0.002. Results
are obtained for ionization from the ground state of the 1D Yukawa potential andR = [0, ±460]. The parameters of the XUV pulse were: peak intensity I = 1 × 1015
W/cm2 , frequency ω = 100 eV, pulse duration T = 400 as, and carrier-envelope
phase φ = −π/2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
3.2

Results of numerical calculations for the times tΨg ,R and tΨ(0) ,R and the time delayg∆tΨg ,R for different time steps δt and a fixed spatial step of δx = 0.02. All the otherparameters were the same as in Table 3.1. . . . . . . . . . . . . . . . . . . . . . . . . 48
3.3

Streaking time delays for different time intervals [t1 , t2 ) of the streaking trace. An
8-cycle streaking pulse and two 1D potentials [VC (z) and VCG (z) (z0 = 650)], havebeen used for the TDSE simulations. All values in this table are given in a.u.. . . . . 80

4.1Time delays for different parameters of the streaking field. We consider resonant TPI
in the 1D Coulomb potential streaked by a 3-cycle laser pulse. XUV laser parameters
are: ωXUV = 81.81 eV, NXUV = 45, IXUV = 1.0 × 1013 W/cm2 , and φXUV = −π/2. . 96

Figures

Figure1.1Time scales of atomic and electronic motions in atoms, molecules, nanostructures,and solids. Figure captured from Ref. [1]. . . . . . . . . . . . . . . . . . . . . . . . .

1.2

2

Semi-classical description of the HHG process: (a) electron is ionized by the driving
laser through tunneling, (b) and (c) electron propagates and accelerates in the laser
field, and (d) electron returns and combines with the parent ion by emitting a highenergy photon. Figure captured from Ref. [8]. . . . . . . . . . . . . . . . . . . . . . .

1.3

3

Generation of a single- or double-burst attosecond pulse: (a) and (c) The spectrogram(window Fourier transform) of the atomic HHG driven by a 5 fs, 750 nm pulse, (b)and (d) spectra of the emitted high-energy photons. (a) and (b) are calculated fora cosine-shaped driving pulse, while (c) and (d) are for a sine-shaped driving pulse.
Filtering out the high-energy part of the HHG spectrum yields a single (a) or double(c) attosecond XUV pulse. Figure captured from Ref. [8]. . . . . . . . . . . . . . . .

4

xii
1.4

Conventional (left, captured from Ref. [1]) and attosecond (right, captured from
Ref. [31]) streak cameras. The conventional streak camera uses a transverse voltage
ramp to deflect photoelectrons ionized by the incident light pulse. Electrons ionizedat different instants have different ionization probabilities and also different horizontal deflections on the detector screen. One can therefore reconstruct the temporal
structure of the incident pulse from the electron distribution on the screen. The attosecond streak camera uses an attosecond XUV pulse to initiate an ultrafast process(e.g., ionization) and a few-cycle IR pulse to streak the momentum of the electron.
By scanning the XUV pulse at different positions of the streaking pulse, one canmeasure the momentum as a function of the relative delay between these two pulses,
i.e., a streaking trace. The time information of the ultrafast process is stored in thestreaking trace. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1.5

7

Physical principle of the attosecond streak camera. Electrons released by the attosecond XUV pulse at different positions of the streaking pulse have different final
momenta, which are approximately given by Eq. (1.1). Figure captured from Ref. [1].

1.6

7Streaking traces from experimental measurements in neon atoms (left, captured from
Ref. [3]) and numerical simulations in hydrogen and helium ion atoms (right, captured from Ref. [37]). In the left figure the upper and lower traces correspond to the2s and 2p electrons respectively. In the right figure the central white line of eachtrace is the expectation value of the momentum. . . . . . . . . . . . . . . . . . . . . 10

1.7(a) Streaking traces obtained from numerical simulations for photoemission in the√
1D Coulomb potential VC (z) = −3.0/ z 2 + 2.0. A 8-cycle streaking pulse with
Is = 1 × 1012 W/cm2 , λs = 800 nm, and φs = −π/2 is used in the calculations. We
compare the streaking trace obtained from the TDSE (blue solid line) with that fromthe original streaking formula, Eq. (1.1), (green dashed line). (b) An enlargement of(a) to show the momentum shift ∆k, the streaking time delay ∆ts , and the relation
between them. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

xiii
2.1

Neon photoionization cross section as a function of photon energy. The blue solid linewith circles are numerical TDSE results by simulating the photoionization process
in the SAE potential produced using the PA method, and the green solid line withasterisks are data from Ref. [62]. The relative error between them are shown as thered dashed line with squares (right axis). . . . . . . . . . . . . . . . . . . . . . . . . . 19

2.2

Comparison of eigenstates calculated from the ITP [(a)-(d)] and spectral [(e)-(h)]methods. We have considered a hydrogen atom in (ρ, z) coordinates. The field-free
propagations of the 1s state are shown for four different propagation steps (N = 0,
3000, 6000, and 9000) with a time step of ∆t = 0.02. When using the spectral
method, the eigenenergies are obtained from the ITP calculations and then the
eigenfunctions are evaluated from Eq. (2.59). . . . . . . . . . . . . . . . . . . . . . . 29

3.1

Schematic diagram for the WS time delay in a scattering scenario. The WS timedelay ∆tWS is defined as the difference of two times: one is the time tp (R) a particle(or wave packet) spends in region R within the potential (upper row), the otherone is the time t0 (R) that a free particle spends in the same region R without
the presence of the potential (lower row), as the region R extends to infinity, i.e.,∆tWS = tp (R) − t0 (R)|R→∞ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37

3.2

Schematic diagram for understanding Eq. (3.8). We consider a normalized delta wave
packet that enters the region R at instant ti and leaves R at tf . In the lower panelR
we plot the probability to find the wave packet in the region R, i.e., R |Ψ(r, t)|2 dr,
as a function of time. Eq. (3.8) calculates the area of the red rectangular region, i.e.,
Sred = (tf − ti ) × 1 = tf − ti , which is essentially the time the wave packet spends inthe region R. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41

xiv
3.3

Schematic diagram for the back-propagation technique. In such a procedure wefirst solve the TDSE describing the interaction of the model potential with the laser
pulse forward in time and obtain the ionizing wave packet by separating it from thewhole wave function long after the laser pulse ceases (upper panel). We then backpropagate one side of the normalized ionizing wave packet (e.g., right side) once
within the potential and once as a free particle and absorb the corresponding wave
packet at the center of the grid (x = 0) by using ECS boundary at the other side(e.g., left side) of the grid (lower panel). The difference of these two propagation
times gives a time delay that relates to the propagation on one side. Repeating thisprocedure again for the other side of the ionizing wave packet and adding the two timedelay contributions together, we finally obtain a time delay for the corresponding
photoionization process accumulated in region R. . . . . . . . . . . . . . . . . . . . . 45

3.4Time delays ∆tΨi ,R and time difference ∆T (Ψi , Ψj , R) as a function of the outer
integration boundary xouter for two potentials: (a) short-range Yukawa potentialand (b) long-range Coulomb potential. Time delays obtained for the ground andfirst-excited states are represented by blue dashed lines and green dash-dotted lines,
respectively; while the red dotted lines show the results for the time difference between the delays. In (b) the black solid and blue dashed lines correspond to two
different forward propagation distances: hxforward i = 2000 and 3000, respectively,
for the ionization from the ground state. In all calculations we have used an XUVpulse with peak intensity I = 1 × 1015 W/cm2 , central frequency ω = 100 eV, pulse
duration T = 400 as, and carrier-envelope phase φ = −π/2 for the ionization. . . . . 49

3.5Time delays and difference between time delays as a function of inner integration
boundary xinner . Symbols and laser parameters are the same as in Fig. 3.4. We alsoplotted the WS time delays as black dots in this figure.

. . . . . . . . . . . . . . . . 49

xv
3.6

Extrapolation method to calculate the WS time delay. The solid blue lines show theouter integration boundary xouter as a function of time for (a) the Yukawa potentialand (b) the Coulomb potential, calculated from Eq. (3.12) by solving TDSE. The
green dashed and red dash-dotted lines are obtained by linearly fitting the solid bluelines in two regions: [150, 250) and [400, 500) respectively. The two insets show the
behaviors near t = 0 for each line. ∆t in panel (a) corresponds to the WS time delay
for a short-range potential. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52

3.7Time delays for ionization from the ground state of the Yukawa potential as functionsof (a) the XUV photon frequency (T = 400 as) and (b) the pulse duration of the
XUV pulse (ω = 100 eV). We have used three different methods to calculate the timedelay: back propagation method (blue diamonds), phase derivative method (green
dashed lines with asterisks), and trajectory extrapolation method (red open circles).
Other laser parameters are: I = 1 × 1015 W/cm2 and φ = −π/2. . . . . . . . . . . . 54

3.8Time delays (upper row) and time delay differences (lower row) as a function of
the outer boundary xouter of R for Yukawa potential (left column) and Coulomb
potential (right column). For each potential we have centered the XUV pulse at two
different positions, which correspond to the maximum (blue dash-dotted line) andzero (green solid line) of the vector potential of the streaking field, respectively. The
XUV parameters are: IXUV = 1 × 1015 W/cm2 , ωXUV = 100 eV, TXUV = 400 as, andφXUV = −π/2. The streaking parameters are: Is = 1 × 1012 W/cm2 , λs = 800 nm,
Ns = 3 cycle, and φs = −π/2. The small box in (d) shows the long-range behaviorof the two curves. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58

xvi
3.9

Time delays with (blue lines with squares) and without (blue lines with asterisks)
streaking field as well as relative differences between the results (green lines withtriangles) as a function of XUV central frequency for (a) Yukawa potential and (b)Coulomb potential. The XUV pulse is centered at the central zero point (maximum
vector potential point) of the streaking field. Other laser parameters are the sameas in Fig. 3.8. For the Coulomb case, the time delays are calculated at xouter = 800.

61

3.10 Relative differences of time delays as a function of the streaking intensity for (a)Yukawa potential and (b) Coulomb potential. The XUV pulse is centered at thecentral zero point of the streaking field. Laser parameters are the same as in Fig. 3.8
expect IIR changes. For the Coulomb case, the time delays are taken at xouter = 400,
800 and 1200. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
3.11 Comparison of streaking time delays from quantum streaking simulations (black
circles) with classical results from the perturbative approach (blue dashed lines) andthe full numerical solution (red solid lines). In our analysis we considered threepotentials: (a) 1D Coulomb potential [VC (z)], (b) the combination of 1D Coulomb
and Gaussian potentials [VCG (z)], (c) 3D Coulomb potential [V (r)]. . . . . . . . . . 67
3.12 Classical predictions for the time delay as a function of the delay ti between XUVionizing and IR streaking pulses for different fitting parameter α by using Eq. (3.37).
Calculations are performed for the 1D potential VCG (z) in Eq. (3.38) with Z = 3.0,
a = 2.0, V0 = −0.5, σ = 2.0, and z0 = 20. . . . . . . . . . . . . . . . . . . . . . . . . 69
3.13 Relative difference between momentum shift ∆k calculated classically from the perturbative approach and the full numerical solution for VC (z) (blue solid line withcircles) and VCG (z) (z0 = 200, green solid line with stars). The relative difference isdefined as |(∆k (pert) − ∆k (numerical) )/∆k (numerical) |. We have considered photoemission of an electron with final asymptotic momentum of 2.0, which is streaked by a
3-cycle 800 nm laser pulse. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69

xvii
3.14 Comparison of momenta calculated using different methods: full numerical solutionof the Newton’s equation (red solid line) and Eq. (3.42) (green dashed line). We haveconsidered an electron propagating in the combined potential of the 1D Coulomb
potential VC (z) and a 3-cycle 800 nm streaking pulse with Is = 1 × 1012 W/cm2 . Asa reference the field-free momentum k (0) (t) of an electron propagating in VC (z) isalso shown as blue solid line. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
3.15 Results of numerical simulations for the streaking time delays ∆ts (solid line withcircles and diamonds) as a function of (a) the range xp of a short-range potential,(b) the frequency of the ionizing XUV field, and (c) the wavelength of the streaking
pulse are compared with those for the Wigner-Smith time delay [dashed lines withsquares in (a) and (b) and solid circles in (c)]. Laser parameters are: IXUV = 1×1015
W/cm2 , TXUV = 600 as, ωXUV = 100 eV [(a) and (c)], φXUV = −π/2, Is = 1 × 1012
W/cm2 , Ns = 3 cycle [(a) and (b)], Ts = 32.02 fs [(c)], λs = 800 nm [(a) and (b)],and φs = −π/2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
3.16 Classical estimates for the ’field-free’ delay distance of an electron, released at thepeak of a three-cycle streaking pulse (wavelength of 800 nm) in the 1D Coulomb
potential, as a function of the XUV photon energy. This distance was determinedby the position of the electron at the time instant when the streaking field changed
to 99% (solid line) and 95% (dashed line) of the peak field strength. . . . . . . . . . 76

xviii3.17 (a) Streaking time delay ∆ts as a function of the position z0 of the Gaussian
potential. We compare results of TDSE simulations (stars, squares, and circles
for V0 = −0.5, −2.0, and −4.0, respectively) with those of classical calculations,Eq. (3.37) (red solid lines). Also shown is the TDSE result without Gaussian potential (diamond on the right end). In the classical calculations the parameter αvaries between 0.985 and 1.084. (b) Comparison of ∆ts from TDSE calculations
for different cycle numbers of the streaking pulse with (green solid line with circles)
and without (blue dashed line with squares) the Gaussian potential. The Gaussian
potential is located at z0 = 650. Other parameters are given in the text. . . . . . . . 79
3.18 2D model potentials, as defined in Eq. (3.47) [(a)] and Eq. (3.48) [(b)], plotted on alogarithmic scale as log[−V (x, y)]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
3.19 Time delay as a function of streaking pulse cycle number for the 2D potentials definedin Eq. (3.47) [(a)] and Eq. (3.48) [(b) and (c)] and polarizations of the streaking (andco-aligned XUV ionizing) field in x-direction (green line with stars), at 45◦ (cyan linewith crosses), and in y-direction (red line with squares). The results are comparedwith those for the pure 2D Coulomb potential without the additional potential (blueline with circles) and a streaking field polarized in x-direction. In (c), the original
streaking time delays for the 2D potential, Eq. (3.48) as shown in (b), have beenshifted to match the result for the shortest streaking pulse for the pure 2D Coulomb
potential. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84
3.20 Streaking field with a pedestal (red solid line): E0 = 5.34 × 10−3 (i.e., Is = 1 × 1012
W/cm2 ), Ts = 331, βp = 0.2, and Tp = 750. As a comparison a basic 3-cycle 800 nm
streaking pulse (blue dashed line) is also shown. . . . . . . . . . . . . . . . . . . . . . 86

xix
3.21 Streaking time delay as a function of (a) and (b) the pedestal length Tp (βp = 0.2),
and (c) and (d) the pedestal strength βp (Tp = 750). The left column shows theresults for Z = 1.0, while the right column is for Z = 3.0. The influence of anadditional static field on the streaking time delay is also present (green dashed lines
with asterisks and crosses) in (a) and (b). . . . . . . . . . . . . . . . . . . . . . . . . 87
3.22 Streaking time delays as a function of the strength of the additional static field fora basic 3-cycle 800 nm streaking field for two potentials: blue solid line with circles
for Z = 1.0 and green solid line with squares for Z = 3.0. . . . . . . . . . . . . . . . 88
4.1

Schematic diagram for obtaining time delays in two-photon ionization (TPI) using
attosecond streaking technique. An electron initially bound in the ground state |gi
is ionized by absorbing two photons from the XUV pulse either through a resonant
state |ei (resonant TPI) or not (nonresonant TPI), and then streaked by the longwavelength field in the continuum. By changing the relative delay between the XUVand streaking pulses, one obtains a streaking trace, from which a temporal shift (orstreaking time delay ∆ts ) can be extracted. After accounting for the propagation
time delay ∆tprop that the electron accumulates in the continuum, we are able toaccess the absorption time delay ∆tabs of the TPI process. . . . . . . . . . . . . . . . 92

4.2

Streaking time delay as a function of XUV pulse duration for (a) 1D Coulomb and(b) 3D He potentials. Three ionization processes are considered: nonresonant TPI(blue solid lines with circles), resonant TPI (green solid lines with squares), andone-photon ionization (red solid lines with asterisks). Laser parameters are given in
the text. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92

xx
4.3

Absorption time delay (red solid lines with asterisks) as a function of XUV pulse
duration for nonresonant TPI in 1D Coulomb [(a)] and 3D He [(b)] potentials, andresonant TPI in 1D Coulomb [(c)] and 3D He [(d)] potentials. Also shown are thestreaking time delays (blue solid lines with circles) and the propagation time delays(green solid lines with squares). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95

4.4

Absorption time delay (red dash-dotted lines) as a function of XUV pulse duration forresonant TPI in 1D Coulomb potential for different streaking parameters (diamonds:
λs = 2400 nm and Is = 1.0×1012 W/cm2 , asterisks: λs = 3200 nm and Is = 1.0×1012
W/cm2 , squares: λs = 4800 nm and Is = 1.0 × 1012 W/cm2 , and crosses: λs = 4800
nm and Is = 1.0 × 1011 W/cm2 ). Also shown are the streaking time delays (bluesolid lines) and the propagation time delays (green dashed lines). . . . . . . . . . . . 97

4.5

Excitation probability as a function of the duration of the XUV pulse for the nonresonant TPI in the 1D Coulomb potential. We have used two methods to calculate
the excitation probability in the first excited state: the numerical solution of TDSE
(black open circles) and the perturbation theory (red solid line). . . . . . . . . . . . 97

xxi
4.6

Trajectory analysis of time delays in TPI. We have shown the trajectories of thephotoelectron for three processes [(a) and (b), blue lines for TPI, green lines fororiginal SPI, and red lines for adjusted SPI] and the time delay as a function of
the position of the photoelectron [(c) and (d), green lines for original SPI, red lines
for adjusted SPI]. The trajectory is calculated using the expectation value hxi of theelectron position as a function of time by solving TDSE. The time delay is obtained as
the time difference for the same hxi between the TPI and SPI trajectories. We haveconsidered both resonant (left column) and nonresonant (right column) TPI. Herethe original SPI means the SPI process initiated by an XUV pulse with ωSPI = 2ωTPI
while the adjusted SPI stands for using an XUV pulse with ωSPI = 2ωTPI +∆ω, whichmakes the final energy of the SPI wave packet exactly the same as that of the TPI
wave packet. The inset in (a) shows how we obtain the time delay between the TPI
and SPI trajectories. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100

4.7Time delay obtained from trajectory analysis as a function of XUV pulse duration.
We have used two ways of calculating the time delays from the trajectories: one
relies on fitting the original time delay curve [e.g., green lines in Fig. 4.6(c) and (d)]to Eq. (4.5) (green lines with squares), the other one relates the time delay to theconstant delay line [e.g., red lines in Fig. 4.6(c) and (d)] (red lines with asterisks).
Results for both nonresonant [(a)] and resonant [(b)] TPI are shown. Also shown
are the absorption time delays (blue lines with open circles) of TPI extracted fromstreaking scenarios. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103

4.8

Comparison of absorption time delays from different methods: numerical streaking
experiment (blue dashed lines with circles) and phase analysis (red solid lines). As
before, both nonresonant [(a)] and resonant [(b)] TPI processes are studied. . . . . . 106

xxii
4.9

Potential energy curves of the ground and lowest lying excited states of our 2D
molecular model systems: (a) model with a dissociative excited state [Eq. (4.11)] and(b) model with bound excited states [Eq. (4.12)]. We are interested in controlling
excitations between states of the blue and red curves.

. . . . . . . . . . . . . . . . . 109

4.10 Two-photon coherent control of excitation to a dissociative molecular state. The
upper panels show that the electric field distribution as a function of time [i.e., theFourier transform of Eq. (4.14)] for (a) φ = 0 and (b) φ = π/2. The insets show
the central field cycle and the carrier-envelope phase (CEP) of each subpulse. In thelower panels the time evolution of the probabilities in the first excited state (bluesolid line) is presented for (c) φ = 0 and (d) φ = π/2. We also show the dissociation
probabilities (red dashed lines). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112
4.11 Same as Fig. 4.10, but for a molecular model system with fixed internuclear distanceR = R0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112
4.12 Excitation probabilities to the ground vibrational (ν2nd = 0, circles and dashedline) and the first excited vibrational states (ν2nd = 1, asterisks and dashed line) as
functions of φ. A spectral phase modulated field defined in Eq. (4.14) is used withparameters: α = 1.2024, β = 8Trevival = 85.6682 fs, ω/2 = 5.83 eV, ∆ω = 0.1 eV,
and I0 = 1 × 1011 W/cm2 . The symbols (circles and asterisks) are numerical results
obtained by solving the TDSE, while the curves (dashed lines) are obtained using
second-order perturbation theory. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115

xxiii4.13 Two-photon coherent control of excitations to a superposition of two vibrational
states. The upper panels show that the electric field distribution as a function oftime [i.e., the Fourier transform of Eq. (4.14)] for (a) φ = 0.96π (dark pulse) and (b)φ = 0.46π (bright pulse). The insets show the central field cycle and the CEP of each
subpulse. In the lower panels the time evolution of the probabilities in the ground(blue solid line) and first excited (red dashed line) vibrational state is presented
for the (c) dark and (d) bright pulses, respectively. Also shown is the ionization
probability (black dashed-dotted line). . . . . . . . . . . . . . . . . . . . . . . . . . . 115

Chapter 1

Introduction

1.1

Electron dynamics on the attosecond time scale
Understanding of quantum electron dynamics has played an extremely important role in thedevelopment of science and technology in the last 100 years. From chemical reactions between
different species to physical interactions of atoms and molecules with external fields, electrons
participate in nearly every aspect and their dynamics has been shown to be fundamental andessential for all processes. Based on the energy-time uncertainty principle, ∆E∆t ∼ ~, one expects
a time scale of tens to hundreds of femtosecond (1 fs = 10−15 s) for vibrational motions in moleculesand solids (Fig. 1.1), which defines the characteristic time scale for the motion of atoms in thesematerials. The time scale of electron dynamics is however even smaller since the energy spacing ofelectronic energy levels is usually of the order of electron volt (eV). An intuitive picture of this timescale is given by the simple Bohr model for a hydrogen atom. Within this model the kinetic energyof the electron in the ground state is Ek = v 2 /2 = 1/2 and therefore the period of the classicalorbital is T = 2π/v = 2π in Hartree atomic units (a.u.)1 , which is on the order of 100 attoseconds(1 as = 10−18 s and 1 a.u. = 24.19 as).
Many interesting physical processes occur in the attosecond time regime, including XUVphotoionization [2–4], strong field tunneling ionization [5, 6], molecular orbital rearrangement after
ionization [7], and many other examples. In order to time resolve these attosecond processes, it isnecessary to use measurement tools with an attosecond resolution (e.g., [1, 8, 9]). For example, if
1

Hartree atomic units, e = m = ~ = 1, are used throughout this thesis, unless otherwise stated.

2

Figure 1.1: Time scales of atomic and electronic motions in atoms, molecules, nanostructures, andsolids. Figure captured from Ref. [1].

3Figure 1.2: Semi-classical description of the HHG process: (a) electron is ionized by the driving
laser through tunneling, (b) and (c) electron propagates and accelerates in the laser field, and(d) electron returns and combines with the parent ion by emitting a high-energy photon. Figure
captured from Ref. [8].

one wants to study the real-time evolution of electronic orbitals after ionization in a molecule bytaking snapshorts of these orbitals, a “camera” with a shutter time on the attosecond time scale isrequired. Without such a fine resolution, the picture taken will have too much information and thefast motions will be smeared out. The developments of laser technology, in particular the advancesin high harmonic generation (HHG) in gases [10–20], have enabled us to generate isolated pulses
with attosecond duration, which therefore provide an extremely useful tool for observing attosecond
electron dynamics.

1.2

Attosecond pulses generated through high harmonic generation
Strong-field lasers with peak intensities on the order of 1013 to 1015 W/cm2 have fieldstrengths comparable to the internal electric field strength in atoms and molecules. Such a strong
field hence enables many nonlinear processes to occur, such as the HHG process, which can happen
in an atomic or molecular gas medium by focusing an intense laser beam [usually near-infrared(IR)] onto it. The HHG process can be well explained using a semi-classical model [12–14] within
three steps (Fig. 1.2): First, an electron initially bound in an atom or molecule is ionized by thestrong laser field through the tunneling mechanism. Then the ionized electron propagates and4

Figure 1.3: Generation of a single- or double-burst attosecond pulse: (a) and (c) The spectrogram(window Fourier transform) of the atomic HHG driven by a 5 fs, 750 nm pulse, (b) and (d) spectraof the emitted high-energy photons. (a) and (b) are calculated for a cosine-shaped driving pulse,
while (c) and (d) are for a sine-shaped driving pulse. Filtering out the high-energy part of the HHG
spectrum yields a single (a) or double (c) attosecond XUV pulse. Figure captured from Ref. [8].

5accelerates under the influence of the strong laser field in the continuum. In the last step, with acertain probability, the electron can be driven back and finally recaptured by the parent atomic ormolecular ion. During the recombination of the electron and the parent ion, the electron has toemit an energetic photon in order to conserve the total energy. Classical calculations show that
only electrons ionized in the second and fourth quarters of the optical cycle of the laser field canbe driven back to the parent ion and then may recombine. This temporal confinement of HHG
therefore usually results in the generation of a train of attosecond pulses spaced by half of the
optical period of the driving laser pulse.
Classical calculations also reveal that electrons ionized at different phases of the driving laser
pulse have different trajectories, returning energies and hence lead to the emission of photons withdifferent energies. A maximum energy of the returning electron occurs for liberation of the electron
at 17 ◦ with respect to the maximum of each half cycle of the driving pulse with the emission of aphoton of energy Ω = 3.17Up + Ip . Here, Ip is the ionization potential of the target system, andUp = I/4ω 2 is the ponderomotive energy that refers to the cycle-averaged energy of the electron in
the driving laser field with I and ω being its intensity and frequency, respectively. The maximum
energy of the HHG photon, given by Ω = 3.17Up + Ip , corresponds to the cutoff energy of the HHG
spectrum [e.g., Fig. 1.3(b)]. These conclusions can also be derived within quantum mechanicalanalysis [15]. As exemplified in Fig. 1.3(a) and (c), a time-frequency analysis [21] of the HHG
spectrum clearly shows when different photon frequencies are emitted.
For many applications one prefers to use an isolated attosecond pulse, of which the generation
is nowadays possible in a few different ways [19, 22]. Among them, one method is the so-called
amplitude gating [18, 23, 24] that makes use of a cosine-shaped laser field comprising merely
a few oscillation cycles (Fig. 1.3). In such a case only the electron wave packet ionized in thesecond quarter of the central half cycle (i.e., the half cycle with largest intensity) has significant
contributions to the HHG spectrum, in particular to the high energy part of the spectrum that isaround the cutoff energy [Fig. 1.3(a)]. As a result, after filtering out the high energy part of thespectrum, a single, isolated attosecond pulse is generated as illustrated in Fig. 1.3(a) and (b). In

6contrast, if one instead uses a sine-shaped laser field, the generated field has two identical bursts
after filtering out the same energy region [Fig. 1.3(c) and (d)]. As the first demonstration of thismethod, experimentalists generated a 650 as pulse with a central frequency of 90 eV by using a 7
fs driving laser at wavelength of 750 nm [23]. Currently, the limit of this method is the generationof a 80 as pulse centered at a photon energy of 80 eV using a 3.5 fs, 720 nm driving laser [24].
Very recently, another very promising method of generating isolated attosecond pulses basedon phase matching has been proposed and demonstrated by M.-C. Chen et al. [25]. It is well known
that electromagnetic fields emitted from different atoms (or molecules) that are phase matchedinterference constructively and add up coherently [26, 27]. The phase matching condition depends
on a few parameters of the experiment, such as the pressure of the gas, the frequency, intensity,and duration of the driving laser pulse [28–30]. In particular, the temporal window over which
phase matching occurs shrinks rapidly with increasing driving laser wavelength [25, 30]. In theexperimental demonstration a 14-cycle driving laser pulse at a central frequency of 2 µm has been
used, which is more stable and easier to generate than few-cycle 800 nm pulses. When the intensityof the driving laser is increased to 1.6 × 1014 W/cm2 (with a proper pressure of argon gas), thetemporal phase matching window is confined to the half peak cycle of the driving pulse and thus
a single, isolated soft X-ray pulse is generated. A pulse of duration of 300 as at a central photon
energy of 140 eV has been observed. Due to its simple and robust properties, this novel scheme
not only provides a new way of generating isolated attosecond pulse but also may make attosecond
science and technology accessible to a broader community.

1.3

Attosecond streaking technique and its application to photoionization
The concept of the streak camera has been around for a few decades with the original intention

to characterize the temporal structure of a short light pulse [32, 33]. On the left side of Fig. 1.4 aschematic view of a conventional streak camera is shown. In such an apparatus photoelectrons areknocked off a photocathode by a short light pulse. The time-dependent ionization probability is
proportional to the temporal intensity of the pulse. During the propagation to the detector screen,

7

Figure 1.4: Conventional (left, captured from Ref. [1]) and attosecond (right, captured from
Ref. [31]) streak cameras. The conventional streak camera uses a transverse voltage ramp to deflect
photoelectrons ionized by the incident light pulse. Electrons ionized at different instants have different ionization probabilities and also different horizontal deflections on the detector screen. One
can therefore reconstruct the temporal structure of the incident pulse from the electron distribution
on the screen. The attosecond streak camera uses an attosecond XUV pulse to initiate an ultrafast
process (e.g., ionization) and a few-cycle IR pulse to streak the momentum of the electron. By
scanning the XUV pulse at different positions of the streaking pulse, one can measure the momentum as a function of the relative delay between these two pulses, i.e., a streaking trace. The time
information of the ultrafast process is stored in the streaking trace.

Figure 1.5: Physical principle of the attosecond streak camera. Electrons released by the attosecond
XUV pulse at different positions of the streaking pulse have different final momenta, which areapproximately given by Eq. (1.1). Figure captured from Ref. [1].

8the electrons are deflected by a transverse microwave voltage ramp V (t) and the degree of deflection
is directly related to the release instants of the electrons. Therefore the spatial distribution of thephotoelectrons on the screen can be converted to the temporal structure of the incident pulse. Sucha conventional streak camera setup can reach a subpicosecond (1 ps = 10−12 s) resolution for a
voltage ramp that is within a fraction of a nanosecond (1 ns = 10−9 s).
The attosecond streak camera shares the concept with the conventional one and was originallyintroduced to measure the temporal structure of an attosecond XUV pulse [23, 34–36]. It consists of
two laser pulses: an attosecond XUV pulse that initiates an ultrafast process (e.g., ionization of theelectron in atoms and molecules), and a few-cycle streaking pulse (usually an IR pulse) that is used
to modulate the final momentum of the photoelectron according to its release instant (see the rightside of Fig. 1.4). Photoelectrons ionized by the XUV pulse at different phases of the streaking field
have different final asymptotic momenta, which can be easily estimated using classical mechanics
by neglecting the Coulomb potential of the parent ion. Assuming that the electron is released at
instant ti and then propagates in the streaking field Es (t), the final momentum is given within thisapproximation by
Z T
kf ' k0 −Es (t)dt = k0 − As (ti ),

(1.1)

ti

where k0 =p2(ω − Ip) is the field-free asymptotic momentum with ω being the XUV photon

energy and Ip being the ionization potential of the system and As (t) is the vector potential of thestreaking pulse with duration T and As (T ) = 0. As illustrated in Fig. 1.5, the final momentum
clearly depends on the ionization instant and therefore temporal information of the ionization (aswell as the XUV pulse) is mapped onto the momentum space. By scanning the application of the
XUV pulse as a function of the phase of the streaking pulse, one can obtain a so-called streaking
trace (Fig. 1.6), which shows the final energy (or momentum) of the photoelectron as a function of
the relative delay between the two pulses. As we will discuss in this thesis such a streaking tracerecords useful information about the ultrafast process initiated by the XUV pulse.

9

1.4

Time delay in photoionization and previous theoretical interpretations
As an important demonstration of the possible use of the attosecond streaking technique inattosecond science, Schultze et al. [3] applied it to study temporal dynamics of the XUV photoionization process in atoms. They intended to address a few very fundamental questions like: Is there
any time delay between ionization processes from different orbitals, for example, the 2s and 2p
orbitals in a neon atom? And when does the XUV photoionization occur in an atom? In theirexperiment two streaking traces are recorded simultaneously for electron emission from the 2s and2p orbitals in a neon atom2 by scanning a 200 as attosecond pulse with energy of 100 eV (energetic
enough to ionize either electron) with respect to a 4 fs IR streaking pulse (left side of Fig. 1.6).
The answer to the first question appears to be easy to find by directly comparing the two
streaking traces of the 2s and 2p photoelectrons. A time delay of 21 ± 5 as for the ionization of
the 2p electron with respect to the ionization of the 2s electron was found in the comparison of the
traces using several different data processing algorithms [3]. The authors hence claimed they haveobserved that the 2p electron takes more time to ionize than the 2s electron. In order to clearly
see that there is indeed a small temporal shift (or time delay) between the traces, which is not so
obvious in the original experimental data (left side of Fig. 1.6), we show also an example of timedelay streaking results from Ref. [37] on the right side of Fig. 1.6, which are based on numericalsimulations of streaking experiments in two atoms. By comparing the two traces from ground stateionization of the hydrogen atom and the helium ion atom, a time delay of about 30 as can be clearlyseen from the inset of this figure.
In order to address the second question, one common way is to use the streaking trace
predicted by the original streaking formula, Eq. (1.1), as a zero-time reference. Since the relative
delay between the XUV and streaking pulses is defined to be the time difference of the center ofthe XUV pulse with respect to that of the streaking pulse with a positive sign representing that the
XUV pulse is applied at the trailing edge of the streaking pulse, one actually uses the center of the2
The same idea has been also applied to a solid system to study the photoionization time delay between ionizationfrom a valence and a conduction band [2].

10Figure 1.6: Streaking traces from experimental measurements in neon atoms (left, captured from
Ref. [3]) and numerical simulations in hydrogen and helium ion atoms (right, captured from
Ref. [37]). In the left figure the upper and lower traces correspond to the 2s and 2p electrons
respectively. In the right figure the central white line of each trace is the expectation value of the
momentum.

1.98
2.05

TDSE

momentum (a.u.)

momentum (a.u.)

original streaking
2

1.95

TDSE
original streaking

1.975
1.97
∆k

1.965

∆ ts

1.96
1.955

1.9
1.95
1.85

(a)
−50
0
50
delay between XUV and streaking pulses (a.u.)

(b)
1.945
−30

−29
−28
−27
−26
delay between XUV and streaking pulses (a.u.)

Figure 1.7: (a) Streaking traces obtained
from numerical simulations for photoemission in the 1D√
Coulomb potential VC (z) = −3.0/ z 2 + 2.0. A 8-cycle streaking pulse with Is = 1 × 1012 W/cm2 ,λs = 800 nm, and φs = −π/2 is used in the calculations. We compare the streaking trace obtained
from the TDSE (blue solid line) with that from the original streaking formula, Eq. (1.1), (green
dashed line). (b) An enlargement of (a) to show the momentum shift ∆k, the streaking time delay∆ts , and the relation between them.

11
XUV pulse as the zero time reference when comparing the experimentally observed or theoreticallycalculated streaking traces to the predicted one. To exemplify this concept, in Fig. 1.7 we show thestreaking trace from the numerical solution of the time-dependent Schrödinger equation (TDSE)
corresponding to photoemission of an electron in a 1D model potential (blue solid line) as wellas the trace predicted by the original streaking formula, Eq. (1.1), (green dashed line). Here we
only show the expectation value of the numerically obtained final momentum in order to visualize
the relative shift between the two traces. A clear temporal offset (or time delay) can be seen in
Fig. 1.7(b) when we zoom into a certain part of the traces. As used by others (e.g., [37–39]), we
denote such a temporal offset as the streaking time delay ∆ts . Theoretically, the trace from theoriginal streaking formula provides a reference to study the origin of time delays and we will takeadvantage of this method in all of our studies within this thesis.
The photoionization time delay experiment by Schultze et al. [3] has initiated a lot of theoretical interest, in particular in view of the interpretation of the delay. In the original report [3],Schultze et al. gave a qualitative estimation and explanation of the delay by relating it to theso-called Wigner-Smith (WS) time delay [40, 41]. The WS time delay is a concept introducedin scattering studies and is defined to be the time difference between the time a particle spends
within the scattering potential and the time that a free particle with an equal final energy spends
in the same region under the condition that the region extends from −∞ to +∞. For a short-range
potential that decays faster than 1/r, the WS time delay can be evaluated by∆tWS =

dϕ
,
dE

(1.2)

with ϕ being the scattering phase shift. Since the WS delay diverges for a long-range Coulombinteraction [41], Schultze et al. [3] employed Eq. (1.2) to avoid the divergence problem, which later
raised concerns from some theorists including us (e.g., [42, 43]). This WS delay concept was later
revisited by Kheifets et al. in Ref. [44], in which an extension of this method was proposed totake account of the extra phase shift induced by the multielectron correlation effect3 . However,3

Although its relation to the streaking time delay remains somehow unclear, the WS time delay concept in the12neither the original WS calculation nor the extension could reproduce the time delay measured in
the experiment.
In order to seek further interpretation, it was proposed to compare the streaking time delay with the WS time delay for some single-electron potentials such as the hydrogen atom andthe helium ion atom, in which direct numerical simulations of streaking experiments without any
approximations are possible. Nagele et al. [37] and Zhang et al. [52] independently found that, generally, the streaking time delay is not given by the WS time delay calculated using Eq. (1.2) evenfor these simple single-electron systems, except for very short-ranged potentials. This discrepancy
however implies there should be some additional contribution to the streaking time delay. Nagele
et al. [37] pointed out that the (classical) propagation of the photoelectron in the coupled field of
the atomic potential and the streaking pulse can induce an extra momentum shift to the original
formula, Eq. (1.1), and therefore should give rise to an additional time delay. This coupling effectcan be viewed as a correction to the strong field approximation (that neglects the Coulomb potential) used in the derivation of Eq. (1.1). They hence proposed to use a classical-trajectory Monte
Carlo simulation, which treats the Coulomb potential on equal footing with the streaking field, tocalculate this coupling related time delay and showed that the obtained coupling delay ∆tCLC 4 isequal to the streaking time delay. However, the authors did not further elaborate on this result.
Instead, inspired by this result, Ivanov et al. [38] derived an analytical formula for the coupling time delay ∆tCLC via a quantum approach based on eikonal approximation. Based on thisapproximation, their analysis implies that for a hydrogen atom the streaking time delay can beapproximated as ∆ts = ∆tWS + ∆tCLC with ∆tWS given by Eq. (1.2). By separating the neon
potential into a −1/r potential and a remaining short-range part, they applied this relation to bothparts and obtained a total time delay of 18 as, which is within the uncertainty of the measured
streaking delay of 21 ± 5 as and seems to explain the experimental result. However, later Nagele etmean time has been widely used in theoretical analysis of time delays in atomic and molecular photoionizations(e.g., [45–51]).
4
Here, “CLC” refers to “Coulomb laser coupling”. Since the coupling delay ∆tCLC is induced by the classicalpropagation of the electron in the coupled field the ionic potential and the streaking pulse, most of the time we callit propagation time delay ∆tprop in this thesis.

13
al. [53] pointed out that the sum rule ∆ts = ∆tWS + ∆tCLC cannot be generalized to short-range
potentials since for short-range potentials the streaking time delay is found to be equal to the WS
time delay and hence there should be no coupling delay contribution. As a result they interpreted
that the coupling delay ∆tCLC solely arises from the coupling effect of the long-range tail of the
Coulomb potential with the streaking field. In other words, they suggested to separate the wholepotential into two parts, where for the short-range part ∆ts = ∆tWS while for the long-range part∆ts = ∆tCLC . Since the long-range Coulomb tail is the same for single-electron systems with thesame parent ion charge, they further proposed that one can calculate the coupling delay by employing ∆tCLC = ∆ts − ∆tWS for the hydrogen atom and then make use of it for all other more complex
systems. Similar to the work of Ivanov et al. [38], Zhang et al. [39] also obtained the coupling timedelay ∆tCLC using eikonal approximation and compared it with the streaking time delay, but they
did not further relate their results with the proposed sum rule. Therefore, up to now, the debate
about the interpretation of the streaking time delay even in single-electron systems is still goingon. In particular, the following questions remain not fully answered: Is the streaking time delay
related to the WS time delay? If yes, how are they related? For a short-range potential, is therea coupling time delay and how can it be calculated? Is the sum rule ∆ts = ∆tWS + ∆tCLC true?
And is it necessary to separate a Coulomb potential into two parts (short-range and long-range) inorder to apply the proposed sum rule accordingly?
Meanwhile some other effects including the polarization effect and the multielectron effectwere also investigated. For example, Baggesen et al. [54], Zhang et al. [39], and Nagele et al. [37]independently studied the influence of the polarization of the initial state due to the streaking pulseon the streaking time delay. Although they reach different conclusions for which initial state thiseffect is present and how the time delay related to this effect can be calculated, their works indicatethat this effect may influence the streaking time delays measured in streaking experiments. Anotherimportant effect in streaking experiments is the multielectron correlation in photoionization. Dueto the properties of noble gases, experimentalists prefer to use them in studying ultrafast processes,
which however may increase the theoretical difficulty in analyzing them. The neon atom used by

14
Schultze et al. in the time delay experiment is such an example. There are ten electrons in a neon
atom and at least eight of them (2s and 2p electrons) will rearrange during the photoionization
process in the delay experiment even though only one of them will be ionized. The correlation effecthas been shown (e.g., [44, 53, 55, 56]) to be important and may have a quite significant contribution
to the measured delay. However, due to the complexity of this multielectron system, solving the
problem from first principles without any approximations is still not feasible.
The time delay concept can also been implemented using the so-called reconstruction of
attosecond beating by interference of two-photon transition (RABITT) technique, in which an
attosecond pulse train is used instead of an isolated pulse (see Refs. [4, 57] for details). Theoretically,
it was shown (e.g., [58, 59]) that the RABITT technique gives rise to the same time delay as theattosecond streaking technique if applied to the same system.

1.5

Organization of this thesis
As discussed in the last section, the debate on the interpretation of the streaking time delay is

still ongoing and many related questions need to be addressed, even for single-electron systems likethe hydrogen atom. In order to contribute to this debate and provide an independent interpretationof the streaking time delay, we study the time delay concept as well as the streaking experimentin detail and systematically in this thesis. Inspired by the concerns about the sum rule ∆ts =∆tWS + ∆tCLC , we investigate the problem from two aspects: On the one hand, we start our studies
by revisiting the WS time delay concept including the proposal of a new time-dependent numericalmethod to calculate ∆tWS and its application in the photoionization process. On the other hand, we
employ a classical model to describe the electron dynamics in the coupled field of the ionic potential
and the streaking pulse and calculate the coupling (or propagation) time delay. Based on the resultsfrom these two aspects, we are able to conclude that: (a) The WS time delay indeed diverges fora long-range potential and therefore is not suitable to explain (even partially) the streaking timedelay. (b) For single-photon ionization (SPI) in a single-electron system (without polarization
effect), the coupling (or propagation) time delay is equal to the streaking time delay and hence it

15
is the only contribution to the streaking delay. Moreover, the coupling time delay, or equivalently
the streaking time delay, can be interpreted as the sum of piecewise field-free time delays weighted
by the ratio of the instantaneous streaking field strength relative to the field strength at the instantof ionization5 . With our understanding of the attosecond streaking technique, in particular of thestreaking time delay, we are able to use this tool to study other ultrafast processes such as thetwo-photon ionization (TPI) process, in which a few new and interesting phenomena are observedin our numerical simulations. The structure and organization of this thesis is given as follows.
In Chapter 2 we review the numerical model that we use to describe the atomic or molecularsystem of interest as well as the method that is utilized to numerically solve the correspondingtime-dependent Schrödinger equation (TDSE). We first introduce two different methods to obtaina single-active-electron (SAE) model for a multielectron system, which can significantly reduce ournumerical cost in simulating laser-matter interactions. Then in the second section we review thesolution of the TDSE in a grid representation using the finite difference method. The propagation
scheme, the method of obtaining the initial eigenstates, and the boundary conditions are describedin detail respectively in this section as well.
Chapter 3 is the central part of this thesis, in which we study the time delay in SPI from a
few different perspectives. In the first section, we introduce a new numerical method to calculate
the WS time delay for the photoionization process based on its fundamental definition by solving
the TDSE. This new method gives a time-dependent view on the WS time delay and thus clearlyshows its evolution during the propagation of the electron in the potential, which is not available
from the original definition of the WS time delay. The results reconfirm that the WS time delay
is an ill-defined concept for photoionization in a long-range potential and hence its application to
the streaking experiment is questionable. In the second section, we focus on the analysis of thestreaking experiment. We propose to consider the electron dynamics after the photon absorption
in the combined field of the ionic potential and the streaking field by a classical model. The
5

We define the instant of ionization as the instant of the transition of the electron into the continuum by absorbing
photons.

16
results for the propagation time delay given by this model indeed agree very well with the streaking
time delay obtained in TDSE simulation, which indicates that the propagation time delay is the
only contribution to the streaking time delay, i.e., the streaking time delay arises from the finite
propagation of the electron from the instant of ionization to the end of the streaking pulse in the
combined potential. We further show that one can interpret the propagation time delay as a finite
sum (or integral) of the piecewise field-free time delays weighted by a relative field strength. We
also note that the finite-range property of the streaking time delay opens the possibility to use thestreaking technique as a new tool of imaging. To relate our theoretical analysis to the experiment,
we also study the influence of some experimental parameters such as the presence of the additional
static field and the long tail of the streaking pulse on the streaking time delay.
Our studies of the streaking time delay in SPI imply that the photon absorption in SPI
occurs most likely (in a quantum mechanical sense) at the center of the attosecond XUV pulse. In
Chapter 4, we extend our studies to two-photon absorption processes. After properly accounting
for the propagation (or measurement-induced) time delay, we show that for the nonresonant TPI,
the two-photon absorption process occurs simultaneously with SPI, i.e., most likely at the centerof the XUV pulse. However time delays obtained for the resonant TPI are significantly different:
The absorption of the second photon occurs with a positive time delay with respect to the centerof the XUV pulse and this delay changes linearly with the XUV duration. In the second section,
we confirm these findings from the streaking simulations with a trajectory analysis based on TDSE
simulations and an analysis of the phase acquired in the two-photon absorption process using
second-order perturbation theory. These findings in TPI, in particular the linear dependence of
the absorption time delay on the XUV duration, suggest a new way of controlling the initiation ofultrafast processes, which for example may be useful in controlling chemical reactions. In the third
section of this chapter, we review our early work of two-photon coherent control of excitation inmolecular model potentials, which does not directly relate to the time delay work but may inspiresome future applications of the TPI delay. We end our thesis with a summary of our work in these
five years.

Chapter 2

Numerical Model and Method

It is well-known that ultrafast quantum processes, including the interaction of atoms andmolecules with laser fields, can be described by the time-dependent Schrödinger equation (TDSE),
that however usually does not have an analytical solution. Nowadays with the development of supercomputers, numerically solving the TDSE is becoming more and more popular. Therefore, before
we discuss certain quantum processes in detail, it is necessary to introduce and clarify the modeland method we use in our numerical simulations. In section 2.1, we describe the single-activeelectron (SAE) model, in which only one electron in a multielectron system is treated as active
while all the others are considered to be frozen as a background. In section 2.2, we present ournumerical approaches to solving TDSE, including the methods of obtaining initial field-free states,
propagating the wave function in the time domain, and applying different absorbing boundary
conditions.

2.1

Single-active-electron model for atoms and molecules
Consider a system consisting of N nuclei and n electrons interacting with an external fieldA(r, t). The Hamiltonian of such a system is given by
X 1
X 1
H(R1 , · · · , RN ; r1 , · · · , rn ; t) =[pα − Zα A(Rα , t)]2 +[pi + A(ri , t)]2
2m
2m
α
i
α
i

(2.1)

+ V (R1 , · · · , RN ; r1 , · · · , rn ),where pα and pi are the canonical momentum operator for the nucleus α with position Rα , mass
mα , and charge Zα and electron i with position ri , mass mi , and charge Zi = 1, respectively and18the last potential term consists of three contributions from the charge-charge Coulomb interactions,i.e.,
V (R1 , · · · , RN ; r1 , · · · , rn ) =

XX
α β>α

XX
XX
Zα Zβ
Zα
1
−
+
.
|Rα − Rβ ||Rα − ri ||ri − rj |
α
i

j

(2.2)

i>j

To study the electron and nuclei dynamics in such a system from first principles, one needs to solve
the corresponding TDSE, which is very challenging even for a two-particle (six-dimensional) system
via numerical approaches. Therefore, approximately describing a complex multielectron systemwith simple models is always of interest for theoreticians. In the case of investigating ultrafast
electron dynamics in atoms and molecules, the dynamics of the complex system can be often wellapproximated by the Born-Oppenheimer approximation, in which the nuclei barely change theirpositions and hence can be assumed as fixed in the time window of the electron dynamics. Thisapproximation reduces the system from an (N +n)-particle problem to an n-electron problem, which
effectively decreases the numerical complexity. Furthermore, when analyzing electron dynamics,
almost in all cases only a few electrons participate in the interactions with internal or external fields
and in many cases (e.g., the HHG process) only one electron is effectively involved. This leads to thedevelopment of single-active-electron (SAE) models, in which only one electron is treated as active
and the analysis is focused on the dynamics of this electron in the quantum process of interest.
Meanwhile, all the other electrons are considered to be frozen and their influence on the active
electron is modeled by an effective potential. There are different approaches of generating a SAE
potential for the same system. Here in this section, we introduce two of them, namely, the potential
average (PA) method, and the optimized effective potential (OEP) method.

2.1.1

Potential average method

The idea of the PA method is very intuitive. For an n-electron system with all nuclei fixed,
the wave function and the potential can be written as Ψ(r1 , · · · , rn ) and V (r1 , · · · , rn ), respectively,
with r1 to rn being the electron’s positions. Then the SAE potential can be constructed asR
dr2 · · · drn Ψ∗ (r1 , · · · , rn )V (r1 , · · · , rn )Ψ(r1 , · · · , rn )
R
,
(2.3)
V (r1 ) =dr2 · · · drn Ψ∗ (r1 , · · · , rn )Ψ(r1 , · · · , rn )

19

Table 2.1: Comparison of eigenenergies (in atomic units) from different methods for neon atom.
The experimental values are from Ref. [60], the HF results are from Ref. [61], and all the otherresults are from our calculations.
orbital1s2s2pexperimental−31.979−1.782−0.797

HF−32.772−1.930−0.850PA−32.520−1.244−0.457LSDA−30.306−1.271−0.442KLI+SIC+LSDA−31.129−1.662−0.805

0.3

0.06

0.2

0.04

0.1

0.02

0
900

1000

1100
1200
1300
photon energy (eV)1400

relative error

cross section (Mb)

SAE+TDSE
reference
relative error

01500

Figure 2.1: Neon photoionization cross section as a function of photon energy. The blue solidline with circles are numerical TDSE results by simulating the photoionization process in the SAE
potential produced using the PA method, and the green solid line with asterisks are data from
Ref. [62]. The relative error between them are shown as the red dashed line with squares (rightaxis).

20
where we assume the electron with coordinate r1 is the active electron and integrate out all the othercoordinates r2 to rn . For most of the atoms, the multielectron wave function can be respresented
in terms of Slater type functions, which can be found in some references. And for a complex
system like a molecule, one can use quantum chemistry softwares (e.g., GAUSSIAN and GAMESS)
to calculate its wave function. Usually the wave function is in a Hartree-Fock (HF) form, which
consists of one Slater determinant or a sum of several Slater determinants.
To simplify the model, we further assume that our wave function can be written in the Hartree
approximation, i.e.,Ψ(r1 , · · · , rn ) = φ1 (r1 ) · · · φn (rn ),(2.4)

where φi (ri ) is the normalized wave function of the ith electronic orbital, and for different orbitals,
φ(r)’s are orthogonal. As mentioned above, the potential of a complex system with electrons(labeled by i and j) and fixed nuclei (labeled by α and β) is given by
V (r1 , · · · , rn ) =

XX
α β>α

XX
XX
Zα Zβ
Zα
1
−
+
.
|Rα − Rβ ||Rα − ri ||ri − rj |
α
i

j

(2.5)

i>j

Since the first term here does not depend on r, it does not change when substituted to Eq. (2.3).
For the second term, if i = 1, we haveXαhX φ∗ (r1 )(Zα /|Rα − r1 |)φ1 (r1 )Zα1i =|Rα − r1 |φ∗1 (r1 )φ1 (r1 )
α
X
Zα
=
,
|Rα − r1 |
α

(2.6)
(2.7)

which also remains unchanged under Eq. (2.3). When i 6= 1, we have
XX
α

i6=1

XX
Zα
h
i =
|Rα − ri |
α
i6=1

R

φ∗i (ri )(1/|Rα − ri |)φi (ri )driR ∗φi (ri )φi (ri )dri

X X Z |φi (ri )|2=dri .
|Rα − ri |
α

(2.8)
(2.9)

i6=1

In our current model, all the nuclei are assumed to be fixed, so the double sum integral above isa constant, which is equivalent to a shift of all energy levels. Now we consider the third term, if

21
j = 1, then
X
i6=1

X
1
i =
h
|r1 − ri |

R

i6=1

φ∗i (ri )(1/|r1 − ri |)φi (ri )driR ∗φi (ri )φi (ri )dri

X Z |φi (ri )|2dri .
=|r1 − ri |(2.10)(2.11)

i6=1

For j 6= 1, we haveXX
i>j j6=1

XX
1
h
i =
|ri − rj |i>j j6=1RR

φ∗i (ri )φ∗j (rj )(1/|ri − rj |)φi (ri )φj (rj )dri drj
RR ∗φi (ri )φ∗j (rj )φi (ri )φj (rj )dri drj

X X Z Z |φi (ri )|2 |φj (rj )|2
dri drj .
=|ri − rj |(2.12)(2.13)

i>j j6=1

This double sum is also a constant, which does not influence the shape of the effective potential.
Thus, by neglecting the constant terms, we can write down the SAE potential as
V (r1 ) = −

X
α

X
Zα
+
|Rα − r1 |i6=1Z

|φi (ri )|2dri .
|r1 − ri |(2.14)The physical meaning of each term in Eq. (2.14) is clear: the first term is the attraction energy
between the active electron (associated with coordinate r1 ) and the nuclei, and the second termis the repulsion energy between the active electron and all the other frozen electrons, which now
are treated as an electron distribution. For an ultrafast process in which the Born-Oppenheimer
approximation applies, we can calculate an effective potential for each nuclei configuration {Rα },in case we want to include the nuclear motions.
To test this model, we calculate eigenenergies of several lowest energy levels for small atoms(e.g., He, Li, and Ne) and compare them with the accurate experimental results. Interestingly, thissimple model provides good results for the lowest inner shell orbital (for atoms, the 1s orbital),
where we find that the calculated energy is very close to the real one (e.g., see Table 2.1). This
indicates that the exchange-correlation effect, which is neglected by this model, is not important
for the most inner-shell electrons. We have also calculated the photoionization cross section as afunction of photon energy for several atoms. For the neon atom (Fig. 2.1) there is a good agreement
between results of our TDSE calculations based on the present SAE model and those from Ref. [62].

22
2.1.2

Optimized effective potential method

The PA method we proposed in the last subsection is easy to implement and works well for
inner-shell electrons. However, it ignores the important exchange-correlation effect and is therefore
not ideal for outer-shell electrons (see Table 2.1), which however are of importance in many lasermatter interaction processes (e.g., HHG). Fortunately, some more accurate models for outer-shellelectrons have been established and well-tested based on density functional theory (DFT). In thissection, we introduce one modified optimized effective potential (OEP) scheme using the selfinteraction-free DFT.
In the Kohn-Sham (KS) DFT formulation [63], the one-electron effective potential is givenasVσ (r) = −

X
α

Zα
+
|Rα − r|Z

ρ(r 0 )
dr 0 + Vxcσ (r),|r − r 0 |(2.15)

where Rα and r are defined as before, σ is the spin polarization, ρ(r) is the total electron density
ρ(r) =Nσ
XXσ|φiσ (r)|2 ,(2.16)

i=1

and Vxcσ (r) is the exchange-correlation potential given by
Vxcσ (r) =δExc [ρ↑ , ρ↓ ].
δρσ (r)(2.17)

Once the specific form of the exchange-correlation functional Exc [ρ↑ , ρ↓ ] is known, we can solve theKS equation
1 2− ∇ + Vσ (r) φiσ (r) = εiσ φiσ (r),2

(2.18)

in a self-consistent way.
However, in the traditional DFT, most forms of Exc [ρ↑ , ρ↓ ], such as the local spin density
approximation (LSDA) [63], contains spurious self-interaction contributions. This self-interaction
contribution can be seen from the total energy
E[ρ↑ , ρ↓ ] = T [ρ] −

X
α

Z
Zα

ρ(r)
1
dr +
|Rα − r|2

Z Z

ρ(r)ρ(r 0 )drdr 0 + Exc [ρ↑ , ρ↓ ],|r − r 0 |(2.19)

23
where T [ρ] is the total noninteracting kinetic energy. For a one-electron system the last two
terms should in principle cancel each other. However, in real calculations, one has to take someapproximations of Exc [ρ↑ , ρ↓ ], which usually do not satisfy this requirement. The self-interaction
contribution may result in an inaccuracy of the energies of the outer-shell orbitals, because it
causes an incorrect long-range behavior of the effective potential in Eq. (2.15). Thus, to construct
a good SAE model for outer-shell electrons, one has to find some efficient way to eliminate the selfinteraction contribution. Many methods have been suggested to solve the self-interaction problem.
Here we introduce two of them: the self-interaction correction (SIC) method [64] and the KriegerLi-Iafrate (KLI) method [65–67] based on the optimized effective potential (OEP) scheme [68, 69].
In the SIC approach [64], the following form of the exchange-correlation functional is used:
SIC
Exc[ρ↑ , ρ↓ ] = Exc [ρ↑ , ρ↓ ] −

Z Z
Nσ
XX
1
ρiσ (r)ρiσ (r 0 ){drdr 0 + Exc [ρiσ , 0]},
0|
2
|r
−
r
σ

(2.20)

i=1

SIC [ρ , ρ ] , which can effecwhere ρiσ is the electron density of the ith spin σ-orbital. This Exc↑ ↓tively remove the self-interaction contribution, however, leads to different potentials for different
orbitals [64]. So the wave functions for different orbitals are not orthogonal and orthogonalizationis required.
Another promising approach is the OEP scheme [68, 69], in which one uses the HF exchangeonly energy functional
X Z Z φ∗iσ (r)φiσ (r 0 )φ∗jσ (r 0 )φjσ (r)drdr 0 .
2|r − r 0 |

1
HF
Exc
=−

(2.21)ijσ

In the one-electron limit, this functional can cancel out with the third term in Eq. (2.19), andtherefore is self-interaction-free. To make use of this orbital-dependent (not density-dependent)
functional in a DFT approach, it was proposed [69] to choose a one-electron effective potential
VσOEP (r) such that it can minimize the total-energy functional, i.e.,
δE OEP [{φi↑ , φj↓ }]= 0,
δVσOEP (r)(2.22)

which can be converted to an integral equation and then solved together with Eq. (2.18). While
the physical idea of the OEP method is simple, its computational implementation is quite difficult

24
due to the integral equation derived from Eq. (2.22). Krieger, Li, and Iafrate [65–67] suggested an
approximative procedure to circumvent this difficulty, in which the OEP VσOEP (r) is determinedby solving linear equations instead of the complicated integral equation.
In the KLI approach, the exchange-correlation potential is shown to be
KLI
Vxcσ
(r) =

X ρiσ (r)i

ρσ (r)

vxciσ (r) +X ρiσ (r)i6=mρσ (r)KLI(V xciσ − v xciσ ),(2.23)

where
1 X
δExc [φiσ ]=−vxciσ (r) =φiσ δφ∗iσ
φiσ (r)jZ

φiσ (r 0 )φ∗jσ (r 0 )φjσ (r) 0dr ,|r − r 0 |KLI(2.24)KLI(r)|φiσ i,
V xciσ = hφiσ |Vxcσ(2.25)

v xciσ = hφiσ |vxciσ (r)|φiσ i,(2.26)and the second summation in Eq. (2.23) does not include the highest-occupied orbital mσ. The
KLI

term of V xciσ − v xciσ can be obtained by solving the following linear equations
NX
σ −1

KLI

S

(δji − (Mσ )ji )(V xciσ − v xciσ ) = V xcjσ − v xcjσ ,(2.27)

i=1

where
Z
(Mσ )ji =ρjσ (r)ρiσ (r)dr,ρσ (r)(2.28)andS

V xcjσ = hφjσ |NσX
ρiσ (r)vxciσ (r)i=1ρσ (r)|φjσ i.

(2.29)At the time of the development of the KLI method, the HF nonlocal exchange-only functional
in Eq. (2.21) was used. This functional is accurate for the exchange part of Exc , but is computationally more expensive than the DFT Exc ’s, such as LSDA. Thus, in some reference papers [67, 70]
it has been suggested to use a combined version of the KLI method and the SIC method. It completely follows the KLI procedure, while the only change is that, at the step of Eq. (2.24), one uses
the self-interaction-free Exc in Eq. (2.20) rather than that in Eq. (2.21). For the widely-used simple
LSDA functional [63], Eq. (2.24) then yields
 1/3Z
6
ρiσ (r 0 ) 0
vxciσ (r) ={−[ρσ (r)]1/3 + [ρiσ (r)]1/3 } −dr .
π|r − r 0 |(2.30)

25
This new KLI-SIC scheme is self-interaction-free, orbital-independent, and computationally not
expensive. The improvement on computational efficiency becomes particularly important when atime-dependent process is considered, in which one needs to reconstruct the SAE at each time step.
As a summary of the two methods to generate a SAE potential, we compare eigenenergies
of different orbitals obtained using different methods for a neon atom in Table 2.1, in which theexperimental eigenenergies from Ref. [60] serve as reference. We see that, for outer-shell orbitals,
the KLI-SIC method with the LSDA functional gives best results, while the PA method performs
best for the inner-shell orbital. Furthermore, it appears that all methods can give a reasonable
result for the inner-shell orbital. Therefore, it is usually a good choice to use the KLI-SIC method
if electrons from many orbitals participate in the process of interest. For some atoms, the SAE
potential based on KLI-SIC is fitted to an analytical formula with several fitting parameters andthese parameters are given in Ref. [71].

2.2

Numerical methods of solving time-dependent Schrödinger equation
To numerically simulate laser-matter interactions, we need to solve the TDSEii
∂Ψ(r, t) h
= Ĥ0 + ĤI (t) Ψ(r, t),∂t(2.31)

where r is the coordinate of the particle, Ĥ0 is the field-free Hamiltonian of the system, and ĤI (t)describes the interaction of the system with the laser field. In the SAE approximation and using
the length-gauge, this equation yields
 2
p∂Ψ(r, t)=+ V (r) + E(t) · r Ψ(r, t),i∂t
2(2.32)

where p = −i∇ is the momentum operator, V (r) is the SAE potential, and E(t) is the laser
field. Using the standard Crank-Nicolson method [72], this second-order differential equation canbe discretized on a grid in space and time and then numerically solved.

26
2.2.1

Propagation method

For a small time step ∆t, the solution of Eq. (2.31) can be approximated as
Ψ(t + ∆t) ' e−iĤ(t)∆t Ψ(t) '1 − i ∆t
2 Ĥ(t)1 + i ∆t
2 Ĥ(t)Ψ(t),(2.33)

where the Cayley unitary transformation [73] is used in the second step in order to preserve thenorm of the wave function. This equation can be further written as





∆t
∆t
1 + i Ĥ(t) Ψ(t + ∆t) = 1 − i Ĥ(t) Ψ(t).
22(2.34)

Moreover, the Hamiltonian Ĥ(t) can be decomposed into different coordinate directions, where
each component is a one-dimensional (1D) Hamiltonian given byĤx (t) = A(x, y, t)∂∂2+ B(x, y, t)+ V (x, y, t),∂x2∂x(2.35)

where x is the coordinate of interest and y stands for all the other coordinates. Discretizing this1D Hamiltonian using the second-order central scheme, Eq. (2.34) yields
Ψn+1 (t + ∆t) − 2Ψn (t + ∆t) + Ψn−1 (t + ∆t)∆t
An (t)
Ψn (t + ∆t) + i
2
(∆x)2

Ψn+1 (t + ∆t) − Ψn−1 (t + ∆t)+Bn (t)+ Vn (t)Ψn (t + ∆t)2∆xΨn+1 (t) − 2Ψn (t) + Ψn−1 (t)Ψn+1 (t) − Ψn−1 (t)∆t+ Bn (t)An (t)+ Vn (t)Ψn (t) ,= Ψn (t) − i
2
(∆x)2
2∆x
(2.36)

27
where n represents the grid point (xn , y) in the discretized space grid. It is not difficult to showthat this equation leads to the following system of linear equations:  2 X1 1 X2



 0

 .
 ..



 0


0

···

3 X1

0

0

0

0

2 X2

3 X2

0

···

0

0

1 X3

2 X3

3 X3

···

..
.

..
.

..
.

..

.

0
..
.

0
..
.

0

0

0

···

1 XN −1

2 XN −1

0

0

0

···

0

1 XN

0




0



0



..

.



3 XN −1 


2 XN



 Ψ1 (t + ∆t)  Ψ2 (t + ∆t) 






 Ψ3 (t + ∆t) 




.
..








 ΨN −1 (t + ∆t) ΨN (t + ∆t)(2.37)=

 2 Y1 1 Y2



 0

 .
 ..



 0


0


···

3 Y1

0

0

0

0

2 Y2

3 Y2

0

···

0

0

1 Y3

2 Y3

3 Y3

···

..
.

..
.

..
.

..

.

0
..
.

0
..
.

0

0

0

···

1 YN −1

2 YN −1

0

0

0

···

0

1 YN

0




0



0



..

.



3 YN −1 


2 YN





 Ψ1 (t)  Ψ2 (t) 






 Ψ3 (t) 
,



.
..








 ΨN −1 (t) ΨN (t)

where N is the total number of grid points in x direction and
1 Xn

=

2 Xn

=

3 Xn

=

1 Yn

=

2 Yn

=

3 Yn

=

i∆t
i∆t
An (t) −Bn (t),
2
2(∆x)
4∆x
i∆t
i∆t
1−
An (t) +Vn (t),
2
(∆x)
2
i∆t
i∆t
An (t) +Bn (t),
2(∆x)2
4∆x
i∆t
i∆t
−
An (t) +Bn (t),
2
2(∆x)
4∆x
i∆t
i∆t
1+
An (t) −Vn (t),
2
(∆x)
2
i∆t
i∆t
−
An (t) −Bn (t).
2
2(∆x)
4∆x

(2.38)(2.39)(2.40)(2.41)(2.42)(2.43)

28
The matrix on the left side of Eq. (2.37) is tridiagonal and the linear equations can hence be solvedby the typical LU decomposition method with a time complexity of O(N log N ).
For a 2D Hamiltonian in (x, y) coordinates, it can be decomposed asĤ(t) = Ĥx (t) + Ĥy (t),(2.44)

in which the mixing (x, y) terms are separated symmetrically. Then the wave function at timet + ∆t is given byΨ(t + ∆t) ' e−i[Ĥx (t)+Ĥy (t)]∆t Ψ(t)(2.45)−iĤy (t)∆t/2 −iĤx (t)∆t −iĤy (t)∆t/2

'e

e

e

Ψ(t),

where in the second approximation the Baker-Campbell-Hausdorff formula [74] is used. This approximation can be shown to have an accuracy of O(∆t3 ). So, for obtaining the wave function
Ψ(t + ∆t), we first propagate the wave function Ψ(t) in y direction by a time step of ∆t/2, thenpropagate the updated wave function Ψ0 (t) in x direction by ∆t, and finally propagate again thefurther-updated wave function Ψ00 (t) in y direction by ∆t/2. In each x- (or y-) direction propagation
step, Eq. (2.37) is solved for all y (or x) values. The propagation in 3D case is similar
Ψ(t + ∆t) ' e−i[Ĥx (t)+Ĥy (t)+Ĥz (t)]∆t Ψ(t)(2.46)−iĤz (t)∆t/2 −iĤy (t)∆t/2 −iĤx (t)∆t −iĤy (t)∆t/2 −iĤz (t)∆t/2

'e

e

e

e

e

Ψ(t),

where the order of (x, y, z) does not matter as long as the symmetry of the propagator is preserved.
The generalization to higher dimensions is straightforward and not limited to Cartesian coordinates.

2.2.2

Methods of obtaining initial states

We have reviewed the approach of propagating the wave function in time in the last subsection. However, before the start of the propagation, one needs to know the initial condition forthe differential equation, which is usually one of the eigenstates of the system of interest. Theseeigenstates can be obtained by solving the boundary condition problem of the time-independent
Schrödinger equationĤ0 Ψ(r) = EΨ(r),(2.47)

29

Figure 2.2: Comparison of eigenstates calculated from the ITP [(a)-(d)] and spectral [(e)-(h)]methods. We have considered a hydrogen atom in (ρ, z) coordinates. The field-free propagationsof the 1s state are shown for four different propagation steps (N = 0, 3000, 6000, and 9000) with atime step of ∆t = 0.02. When using the spectral method, the eigenenergies are obtained from the
ITP calculations and then the eigenfunctions are evaluated from Eq. (2.59).

30
with certain boundary conditions [e.g., Ψ(r = 0) = Ψ(r = ∞) = 0 for 3D Coulomb potential]. For
1D systems, Eq. (2.47) can be solved easily with the typical diagonalization or shooting method,which may however become extremely slow and difficult for systems with higher dimensions. Fortunately, it has been shown that the imaginary time propagation (ITP) method and the spectral
method are two very powerful and efficient tools to resolve this problem.
The idea of the ITP method is based on the propagation of a randomly-assumed wave function(usually a Gaussian wave packet for atoms) in the field-free Hamiltonian, i.e.,i∂Ψ(r, t)= Ĥ0 Ψ(r, t).
∂t(2.48)

Expanding the wave function in the Hilbert space of Ĥ0 , one obtains the time-dependent wavefunction Ψ(r, t) asΨ(r, t) =X

cj φj (r)e−iEj t ,(2.49)jwhere φj (r) is the jth eigenstate with energy Ej , and cj is the projection coefficient of Ψ(r, t = 0)
on φj (r). If we do an imaginary time propagation instead of the real time evolution, i.e., t → −it0 ,the wave function after N steps of propagations becomesΨ(r, t = N ∆t) =X

cj φj (r)e−Ej N ∆t .

(2.50)

j

When N is large enough, all the other terms are negligible small as compared to the first term,
which can be seen fromcj e−Ej N ∆t= (cj /c0 )e(E0 −Ej )N ∆t → 0.
c0 e−E0 N ∆t

(2.51)

Therefore, the ground state is given by
φ0 (r) = Ψ0 (r, t = N ∆t),(2.52)

where Ψ0 (r, t = N ∆t) here denotes the normalized wave function. The excited states can beobtained with the same idea, while the only modification is that those states that are lower inenergy than the calculated state need to be projected out at every propagation step. For degeneratestates, one has to use the symmetry of the wave function to distinguish them, which can be done

31
by using only half of the grid in the imaginary time propagation. As an example, for the 1s state of
an atom, it has a symmetry of Ψ(ρ, z) = Ψ(ρ, −z) in the (ρ, z) coordinates. We hence only use thehalf grid with ρ > 0 and z > 0, and set the first- and second-order derivatives at z = 0 according
to the symmetry, i.e., Ψ0 (z = 0) = 0 and Ψ00 (z = 0) = 2[Ψ(∆z) − Ψ(z = 0)]/∆z 2 .
The ITP method is easy to implement and the resulting eigenfunctions work well for mostof our simulations. However, for some calculations, especially when excited states are involved, its
accuracy is not ideal. This is due to the error accumulation in the calculation since the higher states
are calculated based on the lower states. Moreover, when two states have close eigenenergies but
are not degenerate, this method converges very slowly and its accuracy is also questionable. Even
for a ground state, we clearly see some probability flux emerging during the field-free propagationof the initial states calculated using ITP [see Fig. 2.2(b)-(d)], when the background probability is
chosen to be small (e.g., 10−16 in Fig. 2.2). We therefore like to introduce another method called
spectral method to improve the eigenfunctions for some high-precision calculations.
The spectral method [75] is also based on the propagation of an initially-guessed wave function. Let us define a correlation function C(t) as
Z
C(t) = hΨ(r, t = 0)|Ψ(r, t)i =Ψ∗ (r, t = 0)Ψ(r, t)dr,(2.53)

which can be further reduced toC(t) =

X|cj |2 e−iEj t ,(2.54)

jusing Eq. (2.49). The Fourier transform of this equation is
e
C(E)
=

X

|cj |2 δ(E − Ej ),(2.55)

jfrom which the eigenenergies Ej can be obtained from the peaks of this spectrum. With a known
Ej , the eigenfunction can be reconstructed fromZ ∞φj (r) = aΨ(r, t)eiEj t dt,(2.56)−∞

where a is determined by normalization. In practice, integrating to infinity is impossible, so the32following window function

w(t) =
1 − cos(2πt/T ), if 0 6 t 6 T,



0,

(2.57)

if t > T,is suggested to keep the integration within the region (0, T ), with T being the total propagation
time. We thus have
1
e
C(E)
=
T

Z T
C(t)w(t)dt,(2.58)

Ψ(r, t)w(t)eiEj t dt,(2.59)0andaφj (r) =T

Z T
0which both can be numerically evaluated easily.
To efficiently make use of Eq. (2.59), one needs very accurate eigenenergies (e.g., at least withan accuracy of 10−6 a.u.), which however may be not easy to extract from the correlation spectrum
calculated from Eq. (2.58) in simulations. Our numerical experience shows that the eigenenergies
calculated from the ITP method are usually very accurate and well beyond our requirement ofusing Eq. (2.59). We therefore propose a combined scheme of the ITP and the spectral methods,
in which we use the ITP method to calculate the eigenenergies and then obtain the wave functions
from the spectral method using Eq. (2.59). Results in Fig. 2.2 [(e)-(h)] show that the quality ofthe eigenstate is greatly improved using this combined scheme.

2.2.3

Absorbing boundaries

During laser-matter interactions, ionization is the most general consequence and also contains
a lot of useful and significant information about the system. To simulate ionization processes, wealways need to propagate the ionizing wave packet, i.e., the outgoing ionization probability, on our
numerical grid. Usually such a wave packet has a momentum in the outgoing direction, finally
it will hence hit and be reflected by the boundaries of our grid, if no boundary conditions areused. Such a reflection will cause some artificial physical processes, such as the interference withan outgoing wave packet that approaches the boundary at a later time, and therefore needs to33
be eliminated. While there exists many different carefully-designed boundary conditions that cansolve this problem, we introduce two most-widely used ones: the masking function method and theexterior complex scaling (ECS) method.
In the masking function method one simply multiplies the wave function by a cosine-based
function that smoothly decays from 1 to 0 in the absorbing region. Considering the 1D case, thismethod gives
Ψ(x) = Ψ(x)fmask (x) = Ψ(x) cosα (π x − x0),
2 x1 − x0(2.60)

where x0 and x1 are the starting and ending points of the absorbing region (x0 < x < x1 ),respectively, and α is a coefficient selected depending on the specific problem. The masking function
boundary is easy to code and works well for most cases, especially for those problems where a little
reflection is allowed to occur.
To suppress the reflection to an even lower level, one can use the ECS method ([76] andreferences therein), in which an exponentially decaying term in the propagator can absorb theoutgoing flux, if a complex coordinate axis is artificially assumed in the absorbing region asx = x0 + (x − x0 )eiη .

(2.61)

Here, η is the scaling angle with 0 < η < π/2, and we have again simplified our model to 1D
and assumed that x0 is the starting point of the absorbing boundary. Under these assumptions
and together with the length-gauge, the solution of the reduced 1D TDSE from Eq. (2.32) in theabsorbing region yields
2

Ψ(t + ∆t) = e−i[p̂ /2+V0 (x)+E(t)x]∆t Ψ(t)(2.62)−i cos(2η)p̂2 ∆t/2 − sin(2η)p̂2 ∆t/2 −iV0 (x)∆t −iE(t)[x0 +(x−x0 )eiη ]∆t

=e

e

e

e

Ψ(t),

with the potential term V0 (x) untransformed in this region andiηe−iE(t)[x0 +(x−x0 )e ]∆t = e−iE(t)[x0 +(x−x0 ) cos η]∆t eE(t)[x0 +(x−x0 ) sin η]∆t .

(2.63)

It is now clear that the second exponential term in the second row of Eq. (2.62) decays exponentiallyand hence can absorb the outgoing probability that enters the boundary region. However, one may

34note that the second exponential term on the right-hand side of Eq. (2.63) may become an undesired
source when the sign of the laser field E(t) is positive. In some extreme conditions, this term maycause the wave function to explode and therefore induce numerical problems. He et al. [76] hence
suggested to drop the eiη term in Eq. (2.63), i.e., to do the propagation using22

Ψ(t + ∆t) = e−i cos(2η)p̂ ∆t/2 e− sin(2η)p̂ ∆t/2 e−iV0 (x)∆t e−iE(t)x∆t Ψ(t),(2.64)

which means the untransformed standard field coupling is used over the entire grid including theabsorbing region. This simple modification of the traditional ECS method has been shown tostrongly absorb the outgoing flux without any undesired explosion and also effectively suppress
reflections as compared to the masking function method [76].

Chapter 3Time Delays in Single-Photon Ionization

Due to the great advances in the development in laser technologies, nowadays we are able
to use attosecond XUV pulses to observe, study, and even control electron dynamics in atoms andmolecules on their intrinsic time scale [1]. As an important application, the time information of
XUV-induced single-photon ionization has been extracted using the attosecond streak camera [2, 3]
as well as the RABITT technique [4]. A substantial time delay (a few tens of attoseconds) between
ionization from different orbitals has been reported in both neon [3] and argon [4] atoms, andtherefore initiated a lot of theoretical interest in interpreting the physical origin of such a time delay.
In this chapter, we want to contribute to this debate by exploring this problem from several different
perspectives and finally giving our understanding of this photoionization time delay. In section
3.1, we introduce a new numerical method of calculating the Wigner-Smith (WS) time delay inthe photoionization process using its fundamental definition from a time-dependent perspective. In
section 3.2, we turn to the experimental measurement of the streaking time delay and systematicallystudy problems related to the streaking experiment, including the physical interpretation of thestreaking time delay, the finite-range property of this delay and its application to imaging, as wellas the influence of pedestal and additional static field on this delay.

36

3.1

Calculation of Wigner-Smith-like time delays from a time-dependent
perspective1
In this section, we propose a time-dependent theoretical approach to calculate time delays inphotoionization, which not only is more general for different cases (e.g., short- and long-range potentials) but also provides more information about the process than the traditional time-independent
methods. This new numerical method is based on the fundamental definition of a time delay andrealized using the back-propagation technique. The time delays calculated using this new method
are then compared with those from the traditional methods, namely, the phase derivative method
and the trajectory method. We also apply this method to study the influence of the parametersof the XUV pulse in the XUV photoionization process as well as the parameters of an additionallaser pulse in the laser-assisted photoionization on the time delay.

3.1.1

Wigner-Smith time delay

The Wigner-Smith (WS) time delay concept was originally introduced and developed for ascattering process [40, 41]. It is defined to be the difference between the time tp (R) a particle
spends in a certain region R of a potential and the time t0 (R) that a free particle spends in thesame region without the presence of the potential (Fig. 3.1). Since a scattering process usually
starts and ends far away from the scattering center, one requires that the region R extends toinfinity in the definition of the WS time delay and therefore mathematically the WS time delay isgiven by∆tWS = tp (R) − t0 (R)|R→∞ .

(3.1)

For scattering in a short-range potential, the WS time delay is related to the scattering phase shift
ϕ by∆tWS =
1

dϕ
.
dE

(3.2)

The results of this section are presented in J. Su, H. Ni, A. Becker, and A. Jaroń-Becker, Phys. Rev. A 87,
033420 (2013) and J. Su, H. Ni, A. Becker, and A. Jaroń-Becker, J. Mod. Opt. 60, 1484 (2013).

37

Figure 3.1: Schematic diagram for the WS time delay in a scattering scenario. The WS time delay∆tWS is defined as the difference of two times: one is the time tp (R) a particle (or wave packet)
spends in region R within the potential (upper row), the other one is the time t0 (R) that a freeparticle spends in the same region R without the presence of the potential (lower row), as the region
R extends to infinity, i.e., ∆tWS = tp (R) − t0 (R)|R→∞ .

38
Wigner [40] has given a beautiful derivation of Eq. (3.2) by using a superposition of twomonoenergetic beams of slightly different energies as a substitute for a wave packet. Here, alternatively, we derive this formula by employing a wave packet analysis. To make our derivation simple
but without losing the central physical idea, we restrict the motion of the particle along the x
direction, i.e., considering a 1D scattering process. As initial condition, we assume that a particle
is propagating towards a short-range potential V (x) from x = −∞, given by the following incident
wave functionZΨ(x, t) =A(k)ei(kx−Et) dk,(3.3)

where k is the momentum and E = k 2 /2. To obtain the classical trajectory of the incident particle,
one uses the stationary phase approximation. According to this approximation, the rapid phasevariation in Eq. (3.3) will generally lead to destructive interference so that the integral is determinedby the point where the phase is stationaryxd(kx − Et) = − t = 0,
dE
k

(3.4)

which further reduces to a trajectory as
t=

x
.
k

(3.5)

After the particle passes through the potential, the wave function of the particle at x = ∞ can be
expressed byZ
Ψ(x, t) =A(k)ei[kx−Et+ϕ(E)] dk,(3.6)

where ϕ(E) is the scattering phase shift. Using the stationary phase analysis again, one obtains
the trajectory of the outgoing particle after the scattering ast=

x dϕ(E)
+
.
k
dE

(3.7)

According to its original definition, the WS time delay is the delay the particle accumulates during
its propagation in the scattering potential as compared to a freely-propagating particle. It is obviousthat Eqs. (3.7) and (3.5) describe the motions of the scattering and the corresponding free particles

39
respectively, and therefore the WS time delay associated with this scattering process is dϕ(E)dE bycomparing these two trajectories. While the above derivation for the WS time delay is based on a
1D assumption, one can generalize this procedure and show that Eq. (3.2) is also true for a moregeneral 3D case, for which a detailed discussion can be found e.g. in Ref. [77]. Due to the factthat a photoionization process is a half scattering process, the application of the WS time delay
formula for a photoionization process is straightforward: its WS time delay is 12 dϕ(E)dE , where again
ϕ(E) is the phase shift for the complete scattering process.
One may note that the above derivation for the WS time delay does not hold for a longrange potential, because the outgoing wave function of the particle can not simply be expressed by
Eq. (3.6), which is derived by assuming the potential vanishes quickly enough when x → ∞. For
long-range potentials, usually there are more terms in the phase of the scattered wave function.
As an example, for the Coulomb potential, there exists a logarithmic phase term of the form(1/k) ln(2kr). Consequently, the time delay is not only position dependent but also diverges whenr → ∞, which we will later discuss in detail.
The above method of calculating the WS time delay, i.e., using Eq. (3.2), is simple and notdifficult to evaluate for a photoionization process. However, such a method has a few drawbacks:
On the one hand, the WS time delay calculated via Eq. (3.2) for a short-range potential only
gives one time delay value for an ultrafast process, which is time-independent and does not allow
a time-dependent analysis of the process. On the other hand, such a method cannot be applied tolong-range potentials such as the Coulomb potential, which exists in all atoms and molecules. Asan alternative, it has been proposed (e.g., [37, 38, 44]) to still employ Eq. (3.2) to calculate theWS time delay by using the Coulomb scattering phase shift. Such an idea however does not matchthe original definition of the WS time delay and only accounts for the short-range part of the total
potential. Moreover, for a given problem it is not always obvious how one should separate the fullpotential into a short-range and a long-range part.
We therefore seek for an alternative time-dependent theoretical approach to calculate timedelays in photoionization, which addresses some of the concerns regarding the WS time delay and40
its determination. We further attempt to apply such an approach in time-dependent numerical
grid simulations, which are known to be a powerful tool in calculating and analyzing processes onan ultrashort time scale. The present theoretical analysis of a time delay is intended to be generaland not focused, in particular, on the recent streaking experiments. Once formulated, tested andestablished this may turn out in future as a useful step towards understanding the physics of timedelays in streaking experiments and other precise measurements of ultrashort time scales.
Our proposal is based on the quantum mechanical expression for the time a particle spends
inside a certain region R of a potential. By comparing this time to the corresponding time fora free particle, a time delay is given, which is well-defined for any finite region. To the best of
our knowledge, this fundamental definition of a time delay has not been applied in the analysisof time-dependent processes initiated or driven by ultrashort laser pulses. It however offers a few
interesting features: First, for any finite region R the time delay is well-defined for any physicalrelevant potential and independent whether or not the limit for an extension of the region towards
infinity exists. This enables a theoretical analysis in particular for long-range potentials without
any restriction of the potential. Second, in the limit to infinity, if well-defined, the time delay
should converge to the WS delay. Third, the time delay can be determined as a function of timeafter the emission of the photoelectron even during its interaction with some other external fields
such as the IR field present in the streaking experiment. This expands the options for a theoreticalanalysis of ultrashort time-dependent processes. Fourth, there is no a-priori separation of shortand long-range parts in the potential necessary and the influence of both contributions can bestudied.

3.1.2

Theoretical method and back-propagation technique

For a particle in a given normalized state Ψ(r, t), the time it spends inside a region R of apotential V (r) can be expressed as [78]Z ∞tΨ,R =

Z
dt

−∞

R

dr|Ψ(r, t)|2 .

(3.8)

41

Figure 3.2: Schematic diagram for understanding Eq. (3.8). We consider a normalized delta wave
packet that enters the region R at instant ti and leaves R Rat tf . In the lower panel we plot theprobability to find the wave packet in the region R, i.e., R |Ψ(r, t)|2 dr, as a function of time.
Eq. (3.8) calculates the area of the red rectangular region, i.e., Sred = (tf − ti ) × 1 = tf − ti , whichis essentially the time the wave packet spends in the region R.

42
In order to understand this interpretation, we consider a normalized delta wave function that enters
the left side of the region R (i.e., [xi , xf ]) at instant ti and then leaves the right side at tf (seeR
Fig. 3.2). One may note that the inner integration in Eq. (3.8), i.e., R |Ψ(r, t)|2 dr, is equal to theprobability to find the wave packet in the region R. Since the wave packet is in the region R only
during the time interval [ti , tf ], the result of the inner integration is always zero at other times.
Therefore the double integration in Eq. (3.8) can be determined by the area of the red rectangularregion that istΨ,R = Sred = (tf − ti ) × 1 = tf − ti .

(3.9)

This is just the time the wave packet spends in the region R. One major advantage of this time
formula is that it does not require time information about when the wave packet enters and leaves
the region R, as long as the wave packet is initially outside this region. Such a property may
become useful for calculating times when the start and end points for a process of interest are notavailable. One may also notice that for such an application the probability needs to be 1 when the
wave packet is in the region R. Thus one needs to normalize the wave packet when making use ofthis formula.
While tΨ,R is, in general, finite for finite regions and any Ψ(r, t), it is useful to compare tΨ,Rto the time spent by a free particle in R (or another reference time):
Z ∞ ZtΨ(0) ,R =dt
dr|Ψ(0) (r, t)|2 .
−∞(3.10)

R

Here, Ψ(0) (r, t) is the free particle state corresponding to Ψ(r, t). The difference between tΨ,R andtΨ(0) ,R defines the time delay associated with Ψ(r, t), the region R, and the potential V (r)∆tΨ,R = tΨ,R − tΨ(0) ,R .

(3.11)

Similar to the WS time delay, the quantity ∆tΨ,R should have a finite limit as the radius of R grows
to infinity if the interaction vanishes quickly enough. Thus, ∆tΨ,R→∞ and the associated quantum
mechanical operator are well-defined for short-range potentials V (r) only.
The above definition provides a useful concept to calculate time delays in time-dependent
processes, in particular on an ultrashort time scale. While it is known as the basis for deriving a43time delay in scattering scenarios, it has not been applied for the theoretical analysis of processes
initiated or driven by ultrashort intense laser pulses. We thus intend to apply this method to
obtain time delays for a photoionization process in our numerical simulations. Physically, we areinterested in the time that an initially-bounded electron needs to leave a certain region (centeredabout the position of the nucleus) following ionization due to the absorption of an XUV photon.
To this end, we note that the expressions above can be readily applied to a photoelectron described(ion)by the continuum wave function Ψi(r, t) in our adoption of Eq. (3.8)

1tΨi ,R =
Pion
(ion)

where Ψi

Z ∞

Z
dt

−∞

R

(ion)

dr|Ψi

(r, t)|2 ,(3.12)(r, t) is normalized by the ionization probability
Z ∞Pion =
−∞

(ion)

dr|Ψi

(r, t → ∞)|2 .

(3.13)

We can then define the time delay associated with the ionization from a specific initial state analogous to Eq. (3.11) as∆tΨi ,R = tΨi ,R − tΨ(0) ,R ,(3.14)i(0)

where Ψi (r, t) is the free-particle state corresponding to the continuum wave function after transition from the initial state Ψi (r, t = 0). According to this definition, we expect negative values
for the time delays, since a free wave packet should spend more time in a given region R than thecorresponding wave packet that has the same asymptotic momentum propagating in an attractive
potential. We also consider the difference in the time delays for the ionization from two different
initial states Ψi (r, t = 0) and Ψj (r, t = 0) as∆T (Ψi , Ψj ; R) = ∆tΨi ,R − ∆tΨj ,R .

(3.15)In order to use the above definitions in a numerical simulation of a photoionization process,(ion)

we need to identify the ionizing part of the wave function Ψi(r, t) as well as the corresponding(0)

free particle state Ψi (r, t). Moreover, to make use of Eq. (3.14) correctly, one needs to know thetime of ionization. Otherwise, for a wave function that is initially bound in region R, integrating the44time from −∞ to ∞ using Eq. (3.14) results with an infinite value. However, in traditional TDSE
simulations, it is not straightforward to obtain either the time of ionization or the form of the wave
packet after the transition into the continuum. Moreover, the information of the corresponding freeparticle state is not available in the simulations either. It hence appears to be difficult to makeuse of Eq. (3.14) in our calculations directly. We circumvent this obstacle by using the following
back-propagation technique.
We first solve the TDSE of the system, initially in the state Ψi (r, t = 0), under the interaction
with the external XUV field on a space-time grid:
 2

p
∂
i Ψ(r, t) =+ V (r) + Vlight (r, t) Ψ(r, t),∂t
2(3.16)

where p is the momentum operator, V (r) is the field-free potential of the system, and Vlight (t)represents the interaction with the ionizing light field. After the end of the interaction with thelight field we separate the ionizing part of the wavefunction from the remaining bound parts, either
via projection onto analytically or numerically known states or via spatial separation of the ionizing
part at large distances on the grid. After removal of the bound parts we propagate the remaining
ionizing part of the wave function backwards in time without taking account of the interaction withthe light field using two different Hamiltonians, once including the potential V (r) as 2
∂ (ion)p(ion)
i Ψi (r, t) =+ V (r) Ψi (r, t),∂t
2(3.17)

p2 (0)∂ (0)Ψi (r, t) =Ψ (r, t).
∂t
2 i

(3.18)and once as a free particlei

In order to calculate the time delay ∆tΨi ,R for a given region R, the wavepacket has to be locatedoutside of R at the start of the back-propagation and the propagation needs to be terminated as
the wave packet reaches the center of R, i.e., the location of the residual target ion. The latter
point will be further discussed in the application of the method below.

45

Figure 3.3: Schematic diagram for the back-propagation technique. In such a procedure we firstsolve the TDSE describing the interaction of the model potential with the laser pulse forward intime and obtain the ionizing wave packet by separating it from the whole wave function long after
the laser pulse ceases (upper panel). We then back propagate one side of the normalized ionizing
wave packet (e.g., right side) once within the potential and once as a free particle and absorb thecorresponding wave packet at the center of the grid (x = 0) by using ECS boundary at the otherside (e.g., left side) of the grid (lower panel). The difference of these two propagation times gives
a time delay that relates to the propagation on one side. Repeating this procedure again for theother side of the ionizing wave packet and adding the two time delay contributions together, wefinally obtain a time delay for the corresponding photoionization process accumulated in region R.

46
3.1.3

Application to single-photon ionization by an XUV pulse

The theoretical method outlined above is, in general, applicable to ionization of an atom ormolecule in any light field. Here we present results for the application to photoionization of anelectron initially bounded in two different model potentials. First, we use a short-range Yukawa
potential in 1D
VY (x) = − √Z|x|x2 + ae− b ,(3.19)

where Z is the effective nuclear charge, a is the soft-core parameter, and b is a parameter that
determines the effective range of this 1D potential. For our simulations we choose Z = 3.0, a = 2.0,
and b = 30.0, which relate to energies of −1.6742 and −1.0124 of the ground and first excitedstates. As a long-range interaction we make use of the Coulomb potential in 1D
VC (x) = − √Zx2 + a.

(3.20)For Z = 3.0 and a = 2.0, the energies of the lowest two states are −1.7117 and −1.0807, which areclose to the energies of the Yukawa potential.
For the interaction with the XUV laser pulse we use the length gauge, i.e.,
Vlight (t) = EXUV (t)x,(3.21)

where EXUV (t) represents a linearly polarized XUV pulse with a sin2 envelope, i.e.,
EXUV (t) = E0 sin2 (πt/T ) cos(ωt + φ),(3.22)

where E0 is the peak amplitude, T is the pulse duration, ω is the central frequency, and φ is thecarrier-envelope phase (CEP).
To solve the corresponding TDSE, we use the previously introduced Crank-Nicolson method
in a grid representation. We use a spatial step of δx = 0.02 and a time step of δt = 0.002, andthe grid extends from −4000 to 4000, which is large enough to hold the full wave function in our
simulations on the grid. The initial ground and first excited states are obtained by the ITP method.
We continue the propagation of the wave function after the interaction with the XUV pulse until

47
Table 3.1: Results of numerical calculations for the times tΨg ,R , tΨ(0) ,R and the time delay ∆tΨg ,Rgfor different spatial steps δx and a fixed time step of δt = 0.002. Results are obtained for ionizationfrom the ground state of the 1D Yukawa potential and R = [0, ±460]. The parameters of the XUVpulse were: peak intensity I = 1 × 1015 W/cm2 , frequency ω = 100 eV, pulse duration T = 400 as,
and carrier-envelope phase φ = −π/2.
spatial steptΨg ,RtΨ(0) ,R∆tΨg ,R

0.5
0.2
0.1
0.05
0.02
0.01

275.3389
240.6028
236.8736
235.9751
235.7258
235.6903

274.4569
241.5927
237.9043
237.0109236.7628
236.7274

0.8820−0.9899−1.0307−1.0358−1.0370−1.0371g

the ionizing parts of the wave packet reach a distance beyond |x| > 500 and hence can be well
separated from the remaining bound parts. This allows us to remove the latter parts from the gridand remain the ionizing parts of the wave function only. We then propagate the ionizing parts at
negative and positive x backwards in time independently, either under the influence of the potential,
VY or VC , or as a free particle. We determine the corresponding times tΨi ,R and tΨ(0) ,R for bothiparts of the ionizing wave function and added the two contributions. In the 1D calculations wedefined the region as R = [±xinner , ±xouter ], where xinner and xouter ≤ 500 are the inner and outerboundaries, respectively, and the ±-signs apply to back-propagation of the two parts of the ionizing
wave packet along the positive/negative x-axis, respectively. We absorb the wave function beyond
the inner boundary xinner using the ECS method [76].
Based on the results of recent experimental observations and theoretical calculations, we
expect that the time delay ∆tΨi ,R as well as the difference in the time delays for the ionizationfrom different initial states ∆T (Ψi , Ψj , R) are of the orders of a few tens of attoseconds. Resolutionof such small times requires an analysis of the time and spatial steps in the numerical simulations inorder to establish appropriate limits for grid parameters towards a convergence of the results in ourpresent studies. In Tables 3.1 and 3.2 we present a set of numerical results obtained for different
δx and δt in the case of the Yukawa potential. We see that a convergence of the time delay ∆tΨi ,R

48
Table 3.2: Results of numerical calculations for the times tΨg ,R and tΨ(0) ,R and the time delayg∆tΨg ,R for different time steps δt and a fixed spatial step of δx = 0.02. All the other parameters
were the same as in Table 3.1.
time steptΨg ,RtΨ(0) ,R∆tΨg ,R

0.1
0.05
0.02
0.01
0.005
0.002
0.001

237.0072
236.0391
235.7752
235.7378
235.7284
235.7258
235.7254

238.0431
237.0757
236.8121236.7747
236.7654
236.7628
236.7624−1.0359−1.0366−1.0369−1.0369−1.0370−1.0370−1.0370gwithin less than 0.001 a.u. (i.e., < 0.025 as) is reached for a time step of δt = 0.002 and a spatial
step of δx = 0.02. Similar conclusions hold for our studies with the Coulomb potential as well.
As mentioned above, the time delay ∆tΨi ,R depends on the size of the region R, and shouldbe negative and converge to a finite limit for short-range potentials only. To test these expectations,
we perform a set of simulations for the time delays for photoionization from the ground and excitedstates of both potentials as a function of xouter by fixing xinner = 0, i.e., for the region R = [0, xouter ].
As expected, the values for the time delays are negative and decrease as the size of R (i.e., xouter )increases for each of the results presented in Fig. 3.4. For the Yukawa potential (panel a) convergenceis found for outer boundaries xouter > 150. Consequently, for large values of the outer boundary
we obtain a well-defined value for the time difference ∆T (red dashed line) of the time delays forionization from the ground and the excited states.
In contrast, our results do not show a convergence for the time delays as a function of xouterin the case of the long-range Coulomb potential [see Fig. 3.4(b)]. This reflects the well-known
logarithmic divergence of the time delay for the Coulomb potential in scattering processes and,hence, for ionization from any bound state within the potential. Of course, in these cases a WS
time delay as the derivative of the phase shift [see Eq. (3.2)] cannot be defined as well, since
its derivation requires a finite limit of ∆tΨi ,R→∞ . It is interesting to point out that the results in
Fig. 3.4(b) further show that the logarithmic divergence is, in general, still present for the difference

49

0

∆ tf

∆tg−∆tf−0.6−0.8∆ tg−∆ tf−0.5
time delay (a.u.)

−0.4

−1

−1.5

−2

−1
(a)
0

∆ tg

∆tf

−0.2
time delay (a.u.)

0

∆tg

100

200
300
x_outer (a.u.)

400

500

−2.5
0

(b)
100

200

300
x_outer

400

500

Figure 3.4: Time delays ∆tΨi ,R and time difference ∆T (Ψi , Ψj , R) as a function of the outer
integration boundary xouter for two potentials: (a) short-range Yukawa potential and (b) long-range
Coulomb potential. Time delays obtained for the ground and first-excited states are represented
by blue dashed lines and green dash-dotted lines, respectively; while the red dotted lines show
the results for the time difference between the delays. In (b) the black solid and blue dashed lines
correspond to two different forward propagation distances: hxforward i = 2000 and 3000, respectively,
for the ionization from the ground state. In all calculations we have used an XUV pulse with peakintensity I = 1 × 1015 W/cm2 , central frequency ω = 100 eV, pulse duration T = 400 as, andcarrier-envelope phase φ = −π/2 for the ionization.

0

time delay (a.u.)

−0.2
−0.4
−0.6
−0.8
−1

∆tg

−1.2

∆tf

−1.4
0

∆tg − ∆tf1

2
3
x_inner (a.u.)

4

5

Figure 3.5: Time delays and difference between time delays as a function of inner integration
boundary xinner . Symbols and laser parameters are the same as in Fig. 3.4. We also plotted the
WS time delays as black dots in this figure.

50
between a pair of time delays obtained for the ionization from two different initial states. Thus, sucha time delay difference does not simply depend on the short-range character of the potential but also
contains information about the long-range part of the potential and is therefore not well-defined
as well. The present results agree well with the conclusions of early works on time delays [41]. We
may however reemphasize that any time delay obtained for a finite region via the present method isfinite and therefore well-defined, even in the case of the long-range Coulomb potential. As pointed
out above, this allows us to study certain aspects with respect to the parameters of the XUV pulse
and the effects of an additional field even in the Coulomb case.
It is necessary to point out that for the Coulomb potential the time delay also depends onthe distance that the ionizing wave packet is propagated in forward direction. This is due to thelong-range character of the Coulomb potential, since the central momentum of the ionizing wave
packet decreases with an increase of the forward propagation distance. Thus, the velocity of thefree particle during the back-propagation decreases as well. In Fig. 3.4(b) we show this effect byshowing results for the time delay from the ground state for two forward propagation distances:hxforward i = 2000 (black solid line) and 3000 (blue dashed line). As expected, the time delays for
hxforward i = 3000 are slightly smaller than those for hxforward i = 2000. This shows the need to use
rather large grids for the numerical simulations in the case of a Coulomb potential. However, thissmall dependence on the forward-propagation distance does not change our conclusions regarding
the convergence of the results towards infinite regions.
We also note from the results in Fig. 3.4 that the time delay increases most strongly inthe region close to the center of the potential, where the potential changes most strongly. This
indicates that the results should depend on the choice of the inner boundary xinner of the region
R. To study this feature, we fix the outer boundary of R at xouter = 500, which is large enoughto obtain converged results in the case of the Yukawa potential, and then vary the inner boundary
xinner . The results in Fig. 3.5 show the expected dependence on the choice of xinner : the absolutevalues of the time delays decrease by half as xinner increases from 0 to 5. In the remainder of thepresent studies we have chosen the xinner = 0 as the inner boundary, since this value corresponds

51
to the expectation value of x for all the bound states investigated here.

3.1.4

Comparison with results of other theoretical methods

Up to now, we have introduced the theoretical principle of our new method to calculate timedelays for photoionization processes and also established a numerical technique to demonstrate its
application. In order to further support our theory, we want to compare the time delays calculatedusing our new method with results from those of previously existing and well-tested methods.
To this end, for a short-range potential, we like to make use of the WS time delay formula, i.e.,
Eq. (3.2), and further review another method based on the analysis of the photoelectron trajectoryin the continuum.
To make use of the WS time delay formula for the 1D Yukawa potential of Eq. (3.19), we
calculate the scattering phase shift from a time-independent scattering approach. We assume√
that an electron with energy E (i.e., momentum k = 2E) is incident from x = −∞ towards
the Yukawa potential. We then solve the corresponding time-independent Schrödinger equation
numerically using the fourth order Runge-Kutta method up to |x| = 500, project the numericalsolution onto the appropriate plane-wave solutions for x → ±∞, and obtain the corresponding
phase shift ϕ(E). By repeating the same procedure for an electron with energy E + ∆E, we know
the phase shift ϕ(E + ∆E) and then evaluate the WS time delay approximately as(scat)∆tWS (E) =ϕ(E + ∆E) − ϕ(E).
∆E

(3.23)In order to take account of the energy spread of the ionizing wave packet in a specific photoionization(scat)

process, we averaged ∆tWS (E) over the energy spectrum of the wave packet, as obtained in ourtime-dependent numerical simulations. Finally, we consider the photoionization as a half-scattering
process and divide the result of the average by two. The resulting WS time delays are shown as
black dots in Fig. 3.5 and are in good agreement with our numerical results, obtained from thetime-dependent calculations, for xinner = 0 and xouter = 500. This is in support of the applicabilityof our approach to obtain time delays from the time-dependent numerical simulations.

52

500
TDSE
fitting in [150 250)
fitting in [400 500)

x_outer (a.u.)

400

(a)

300
10
200
∆t

5
100
0
−2
0

0

50

0

2

100
150
time (a.u.)

4

200

250

500
TDSE
fitting in [150 250)
fitting in [400 500)

x_outer (a.u.)

400

(b)

300
10
200
5
100
0
0

0

50

−2

100
150
time (a.u.)

0

2
200

4
250

Figure 3.6: Extrapolation method to calculate the WS time delay. The solid blue lines show
the outer integration boundary xouter as a function of time for (a) the Yukawa potential and (b)the Coulomb potential, calculated from Eq. (3.12) by solving TDSE. The green dashed and red
dash-dotted lines are obtained by linearly fitting the solid blue lines in two regions: [150, 250) and[400, 500) respectively. The two insets show the behaviors near t = 0 for each line. ∆t in panel (a)corresponds to the WS time delay for a short-range potential.

53
As discussed before, the derivation of the WS time delay is based on the comparison of two
trajectories, one of the free particle and the other of the scattering particle. This actually provides
another useful way of calculating the WS time delays since both trajectories are easily accessable inour numerical simulations of a photoionization process. In practice, the WS time delay is given as
the time difference between the time obtained by extrapolating the linear part of the photoelectron
trajectory back to the center of the potential (x = 0) and the time zero (chosen as the center ofthe XUV pulse). There is no need to numerically calculate the free particle trajectory since it isjust a line that goes through the origin (i.e., x = 0 and t = 0) with the same slope as that ofthe photoelectron trajectory. Usually, the photoelectron trajectory is calculated as the expectation
value of the ionizing wave packet as a function of time [44, 46]. Here, as shown in Fig. 3.6, we
instead obtain the trajectory by plotting the outer boundary xouter of the region R as a function oftΨ,R (blue solid lines), which requires no more calculations once the back-propagation simulations
are done. As an example, we show how the WS time delay is obtained using this method for the
short-range Yukawa potential in Fig. 3.6(a), where the fit of trajectories in two different regions
gives the same time delay. However, based on our above discussions, it is expected that this method
does not work for a long-range potential, such as the Coulomb potential considered in Eq. (3.20),
since the trajectory of the photoelectron ionized from the Coulomb potential has a logarithmic
term and cannot be fitted to a linear line. This is illustrated in Fig. 3.6(b), from which one sees
the fitted lines have different interceptions with the time axis for different fitting regions.

3.1.5

Dependence of time delay on XUV pulse parameters

In order to systematically test the applicability of our time-dependent method of calculatingtime delays and also study the dependence of the time delay introduced here on the parameters of
the XUV ionizing pulse, in Fig. 3.7 we present our results as functions of (a) the XUV frequency at
a fixed pulse duration of T = 400 as and (b) the duration of the XUV pulse at a fixed frequency ofω = 100 eV. The peak intensity is I = 1 × 1015 W/cm2 and the CEP is φ = −π/2 in each of thesesimulations. Other than the time-dependent back-propagation method (blue diamonds), we also

54

0

(a)

time delay (a.u.)

−2
−4
−6
−8
back propagation
phase derivative
trajectory extrapolation

−10
60

80
100
120
XUV central frequency (eV)

140

−0.9
(b)

time delay (a.u.)

−1
−1.1
−1.2
−1.3
back propagation
phase derivative
trajectory extrapolation

−1.4
200

400

600
800 1000 1200
XUV pulse duration (as)

1400

1600

Figure 3.7: Time delays for ionization from the ground state of the Yukawa potential as functionsof (a) the XUV photon frequency (T = 400 as) and (b) the pulse duration of the XUV pulse(ω = 100 eV). We have used three different methods to calculate the time delay: back propagation
method (blue diamonds), phase derivative method (green dashed lines with asterisks), and trajectory extrapolation method (red open circles). Other laser parameters are: I = 1 × 1015 W/cm2
and φ = −π/2.

55show time delays from the other two methods: the energy derivative with respect to the phase shift(green dashed lines with asterisks) and the trajectory extrapolation method (red open circles). The
results from different methods agree very well with each other and this again shows the applicationof the time-dependent back-propagation method on time delays is solid and convincing.
The trends of the time delays also agree well with qualitative expectations. The absolutevalue of the time delay decreases towards zero as the frequency of the ionizing XUV pulse and,hence, the final kinetic energy of the emitted electron increases [Fig. 3.7(a)]. This is due to thefact that the effect of the potential on the motion of the photoelectron becomes negligible in
the limit of infinitely large kinetic energy of the electron (i.e., infinite large XUV frequency) andtherefore the time spent in the potential approaches that of the free particle in this limit. We
further find that the absolute value of the time delay decreases with an increase of the XUV pulse
duration [Fig. 3.7(b)]. This dependence is closely related to that presented in Fig. 3.7(a) and canbe qualitatively understood as follows. It is found that the expectation value of the kinetic energyof the ionizing wave packet increases as the XUV pulse duration increases, which causes the timedelay to increase according to its dependence on photoelectron energy in Fig. 3.7(a). Furthermore,
due to the finite pulse duration the ionizing wave packet has a certain bandwidth about a central
kinetic energy. Consequently, the time delay obtained for the wave packet can be considered as anaverage over contributions at particular electron energies within the bandwidth (weighted by theionization probability at a given energy). As indicated by the results in Fig. 3.7(a) the time delay
does not change linearly with the kinetic energy. Therefore, the time delay obtained for a wave
packet will be smaller than its contribution at the central kinetic energy or the expectation value ofthe kinetic energy. This difference decreases and, thus, the time delay for the wave packet increases
as the energy bandwidth of the wave packet decreases, i.e., as the pulse duration increases.

3.1.6

Application to laser-assisted XUV photoionization

As introduced in Chapter 1, the HHG process is a popular and promising way of generating
attosecond XUV pulses. Usually a few-cycle intense IR laser field is used to drive the HHG process.

56
So in many applications using the attosecond XUV pulses, the strong laser field is present andtherefore influences the interaction of the XUV pulse and the systems of interest. In particular,
the interplay of the strong laser field and the ionic potential of the system is very importantfor many attosecond measurements, including the attosecond streak camera [34], laser-induced
tomography [79], diffraction imaging of molecules [80], and many other applications. In this section,as another application of our proposed method of calculating time delays in photoionization, we liketo extend this idea to the situation, in which the XUV photoionization is assisted by a streaking
field and the interplay of the streaking field and the ionic potential can be treated together as atime-dependent combined potential.
The attosecond streaking experiment resembles such a set up and we like to study the effectof a streaking field on WS-like time delays. To this end, we will first discuss how the streaking
field can be included in our numerical simulations of time delays and then study the impact of a
streaking field on the time delays for the short-range Yukawa as well as the long-range Coulomb
potentials by varying the parameters of the streaking field. We may stress that it is found (seeSection 3.2) that the time delay studied here does not correspond to those observed in streaking
experiments. The difference between these two delays are apparent and significant: the present one
is evaluated directly from the definition of a time delay from the time domain, while the streaking
time delay is extracted from a streaking trace in the momentum space.
The streaking field is considered as part of the potential V (r) in Eq. (3.16) and can be
represented byVs (r, t) = Es (t) · r,

(3.24)

where the common length-gauge is used. Thus, in the present calculations the streaking field istreated on equal footing with the atomic potential. After the forward propagation of the wave
function from its initial state and the separation of the bound and ionizing part of the wave
function, we then propagate the ionizing part of the wave function backwards in time once within
the combination of the atomic potential and the streaking field and once as a free particle.

57
As a result, we obtain the time delay associated with the ionizing part of the wave functionin the combined potential of the short- or long-range interaction and the streaking field as
(IR)

(IR)

∆tΨi ,R = tΨi ,R − tΨ(0) ,R(3.25)(IR)

where tΨi ,R is the time the ionizing wave packet spends in region R in the presence of both the(short- or long-range) atomic potential and the streaking field, and tΨ(0) ,R is the time for thecorresponding free particle.
Here, the forward propagation of the wave function has to be continued as long as thestreaking field is present. As pointed out above, for the long-range Coulomb potential, the timedelay depends on the distance that the wave packet is propagated in the forward direction. In
order to keep the corresponding error small in our current analysis we use a large grid of −13000
to 13000 and terminate the forward propagation when the expectation value of the ionizing wave
packet reaches 8000. We increase the spatial step to δx = 0.1 and the time step to δt = 0.02 ascompared to the previous calculations to save computation time. Test calculations show that therelative error of the present results is about 1%. In order to be consistent but also without loosing
generality, we make use of the Yukawa and the Coulomb potentials with the same parameters as
before, but choose the first excited state as the initial state. We have checked that the ionization
probability induced by the streaking field is negligible up to an intensity of 1 × 1013 W/cm2 .
We apply the ionizing XUV pulse at two positions of the streaking field: one centers at
the maximum point (zero of the vector potential) and the other centers at the central zero point(maximum of the vector potential). In the upper row of Fig. 3.8, the results for the time delays in
the presence of the streaking field (green solid lines: XUV centered at zero of the vector potential;
blue dashed lines: XUV centered at maximum of vector potential) are shown as a function of theouter boundary xouter of the region R for the Yukawa (left) and the Coulomb (right) potentials.
Similar as in the results without streaking field, we see that there is a well-defined limit of timedelays for the short-range Yukawa potential as the region R increases, while there is no convergencefound for the Coulomb potential.

58

time delay (a.u.)

−0.2

(a)

−0.6

0.1

−0.5 (b)

−1

−0.4

−1.5

−0.8
0
time delay difference (a.u.)

XUV at vector potential max
XUV at vector potential zero

400

800

0
0.1

(c)

0.05

(d)

0.05

400

800
−3

10

0

x 10

8
6

0

1200

2000

4000

0

400
800
x_outer (a.u.)

0

400
800
x_outer (a.u.)

1200

Figure 3.8: Time delays (upper row) and time delay differences (lower row) as a function of the outerboundary xouter of R for Yukawa potential (left column) and Coulomb potential (right column).
For each potential we have centered the XUV pulse at two different positions, which correspondto the maximum (blue dash-dotted line) and zero (green solid line) of the vector potential of thestreaking field, respectively. The XUV parameters are: IXUV = 1 × 1015 W/cm2 , ωXUV = 100 eV,
TXUV = 400 as, and φXUV = −π/2. The streaking parameters are: Is = 1 × 1012 W/cm2 , λs = 800
nm, Ns = 3 cycle, and φs = −π/2. The small box in (d) shows the long-range behavior of the two
curves.

59
To illustrate the effect of the streaking field, we present in the lower row of Fig. 3.8 thedifference between the time delays in the presence of the streaking field to those without streaking
field as a function of xouter , i.e.,(IR)∆T = ∆tΨi ,R − ∆tΨi ,R ,(3.26)

Although there is no well-defined limit of the time delays for infinite regions in the Coulomb case,neither with nor without streaking field, for any finite region the time delays introduced here are
well-defined and the effect of the streaking field can thus be analyzed. The same argument applies
to the weak dependence of the Coulomb results on the distance of forward propagation in our
simulations.
For both potentials, we see that the time delay difference ∆T oscillates for xouter < 400. Thisoscillation is due to the presence of the streaking field, since the ionizing wave packet propagatesup to about x ' 400 before the streaking field ceases in present simulations. We note that thedifferences ∆T are small, less than 3% for the Yukawa potential and less than the numerical errorof 1% for the Coulomb potential, as compared to the time delays induced by the atomic potentials
themselves. We therefore do not find a significant effect of the streaking field, neither for a shortrange nor for a long-range potential, on the time delays in our present simulations. As pointed out
before, this result however does not allow to draw conclusions for the streaking experiment, since
time delays introduced here are not equivalent to time delays measured in streaking experiments.
Before we continue to further study the influence of the streaking field, we note a subtle pointin the results obtained for the Coulomb potential, which are presented in Fig. 3.8. While neither
the time delays with and without streaking field converge as a function of the outer boundary
xouter , we find a converged result (within the numerical error) for the time delay difference ∆T , if
the XUV pulse is centered at the zero position of the vector potential of the streaking field [seegreen solid line in Fig. 3.8(d)]. This occurs since in this case the momentum distribution of theionizing wave packet at the end of the forward propagation is the same as that without streaking
field. In contrast, if the XUV pulse is applied at the maximum position of the vector potential, the60final momentum distribution is shifted and thus no convergence of the time delay difference within
the range of present boundaries is found [see blue dashed line in the inset of Fig. 3.8(d)].
The conclusion that the streaking field does not influence significantly the time delay introduced here holds over a large range of XUV frequencies as well as for streaking intensities up to
about 1013 W/cm2 . In Fig. 3.9 we present the results of time delays obtained for (a) the Yukawa
potential and (b) the Coulomb potential with (blue lines with squares) and without (blue lines
with asterisks) streaking field as a function of the XUV frequency. Since the results are in close
agreement, we also show the relative difference between them (green lines with triangles), which
does not exceed 5% and 2% in the Yukawa and Coulomb cases, respectively.
As one would expect, the relative difference between the results of time delays obtained withand without streaking field does increase with an increase of the streaking field intensity. This can
be clearly seen from the results shown in Fig. 3.10. It appears that for streaking intensities up
to 1012 W/cm2 the relative difference between the results is small enough such that there is no
significant effect on the time delay. While the relative difference quickly increases beyond 10% in
the case of the Yukawa potential with a further increase of the streaking intensity, the 10%-limit is
not reached for a streaking intensity of 1013 W/cm2 in the case of the Coulomb potential.

3.2

Interpretation of time delays measured in streaking experiments2
The photoemission time delay experiment in neon atom using the attosecond streaking tech-nique has initiated a lot of theoretical discussions about understanding this experiment, especially
how one should interpret the measured time delay (e.g., [37–39, 44, 81]). As mentioned in the
introduction, it was proposed that this measured time delay consists of two contributions [37, 38]:
one is the WS time delay, which is calculated as the energy derivative with respect to the scattering
phase shift induced by the electrostatic Coulomb potential, and the other one is the propagation
time delay that the photoelectron accumulates during its propagation in the combined potential of2

The results of this section are presented in J. Su, H. Ni, A. Becker, and A. Jaroń-Becker, Phys. Rev. A 88,
023413 (2013), J. Su, H. Ni, A. Becker, and A. Jaroń-Becker, Phys. Rev. A 89, 013404 (2014), and J. Su, H. Ni, A.
Becker, and A. Jaroń-Becker, Chin. J. Phys. 52, 404 (2014).

−0.5

0.04

−1

0.03

−1.5

0.02
XUV only

−2

−2.5
60

0.01

XUV+IR

0

0.02

−2

0.015

−4

0.01

XUV only−6

relative difference

(a)

80
100
120
XUV central frequency (eV)0.005

XUV+IR

0
140

−8
60

relative difference of time delay

0.05

time delay (a.u.)0

relative difference of time delay

time delay (a.u.)61

relative difference

(b)

80
100
120
XUV central frequency (eV)0140

Figure 3.9: Time delays with (blue lines with squares) and without (blue lines with asterisks)
streaking field as well as relative differences between the results (green lines with triangles) as afunction of XUV central frequency for (a) Yukawa potential and (b) Coulomb potential. The XUVpulse is centered at the central zero point (maximum vector potential point) of the streaking field.
Other laser parameters are the same as in Fig. 3.8. For the Coulomb case, the time delays arecalculated at xouter = 800.

0.08
relative difference of time delay

relative difference of time delay

0.3
0.25
0.2
0.15
0.1
0.05
(a)
0 11
10

12

10
IR intensity (W/cm2)

0.07

x_outer=400
x_outer=800
x_outer=1200

0.06
0.05
0.04
0.03
0.02
0.01
(b)

13

10

0 11
10

12

10
IR intensity (W/cm2)

13

10

Figure 3.10: Relative differences of time delays as a function of the streaking intensity for (a)Yukawa potential and (b) Coulomb potential. The XUV pulse is centered at the central zero pointof the streaking field. Laser parameters are the same as in Fig. 3.8 expect IIR changes. For the
Coulomb case, the time delays are taken at xouter = 400, 800 and 1200.

62
the streaking and the Coulomb fields, which is sometimes also called as the Coulomb-laser couplingtime delay. However, the sum rule ∆ts = ∆tWS + ∆tprop seems to be questionable for the following
reasons: First, for a long-range Coulomb interaction, the WS time delay is not a well-defined concept, as originally pointed out by Smith [41] and recently numerically reconfirmed by us [42, 43].
Evaluating it using the energy derivative with respect to the scattering phase shift, as done for ashort-range WS delay, only accounts for the short-range part of the Coulomb potential and does notobey the original definition for the WS time delay. Furthermore, while the sum rule appears to be agood empirical formula (e.g., [53]), a rigorous proof was never given. For instance, in Ref. [38] thisrule has been derived however based on the eikonal approximation, which is quite questionable forthe weak intensities used in streaking experiments. Furthermore, in Ref. [37] the propagation timedelay was calculated using the so-called classical-trajectory Monte Carlo method and the results do
not support the sum rule. Finally, even within the framework of the sum rule, applications differ
in some details: Ivanov et al. [38] apply the rule for both short- and long-range potentials, while
Nagele et al. [37] claim that the sum rule needs to be used for the long-range interaction only. In
the latter work, it has been assumed that the propagation time delay is zero. In order to contributeto the ongoing debate, we like to systematically study the streaking experiment in this section.

3.2.1

Extracting streaking time delays in numerical simulations

As introduced in Chapter 1, the streaking time delay is the temporal shift of the real streaking trace with respect to the trace predicted by the original streaking formula, kf = k0 − As (ti ).
To obtain such a delay in numerical simulations, we need to calculate the streaking trace accurately without making any approximations. To this end, we simulate a streaking experiment by
numerically solving the following TDSE 2
∂Ψ(r, t)pi=
+ V (r) + (EXUV (t) + Es (t)) · r Ψ(r, t),∂t
2(3.27)

where both the XUV and the streaking fields are assumed to have a sin2 envelope shape andbe linearly polarized in z-direction. We solve the TDSE on a grid in space and time by using

63
the common Crank-Nicolson method with spatial steps determined by the specific problem. In
each simulation we propagate the wave function on the grid for a sufficiently long time until bothlaser pulses cease and the ionizing wave packet is far from the ionic core. We confirm that in all
simulations the outgoing wave packets stay on the grid and do not reach its boundaries.
To obtain the momentum distribution we use one of the two following methods: one is toproject the total wave function onto the correct continuum eigenstate, the other is to spatiallyseparate the ionizing wave packet from the total wave function, which is possible due to the longpropagation times used in the simulations, and then perform a Fourier transform. For the latterchoice, we have to propagate the ionizing wave packet sufficiently far away from the nucleus so thatthe small error due to the projection onto the plane wave is negligible.
By varying the relative delay τ between the XUV and the streaking pulses, we obtain thestreaking trace, in which the momentum of the photoelectron kf is given as a function of τ . We
usually present the streaking trace in the form of the expectation value of the momentum as afunction of the relative delay in order to compare with other calculations. For model systems
with more than one dimensions, we calculate the expectation value of the momentum in the laser
polarization direction (z axis) with an opening angle of 5 to 10 degrees. An example for a streaking
trace, over about the central cycle of a 8-cycle streaking pulse at 800 nm and 1 × 1012 W/cm2 , isshown in Fig. 1.7(a). As already pointed out, the streaking time delay relates to a momentum shift∆k and can be understood as the temporal shift of the TDSE trace (blue solid line) with respect tothe trace predicted by the original streaking formula (green dashed line). In practice we however
do not compare these two traces directly to obtain the streaking time delay. Instead we extractsuch a delay by fitting the TDSE streaking trace tokf,z (τ ) = k0,z − αAs (τ + ∆ts )(3.28)

with k0,z , α, and ∆ts as three fitting parameters, which can be retrieved using the least-square
method.

64
3.2.2

Calculating propagation time delays using classical approaches

In view of the preexisting interpretation of the streaking time delay, it appears to be necessary to reinvestigate the original streaking formula, kf ' k0 − A(ti ), derived based on strong
field approximation (SFA). As motioned in the introduction, the original streaking formula canbe derived classically by assuming the photoelectron is set free in the continuum instantaneously
at ionization time ti with momentum k0 and then propagates only in the streaking laser field,
i.e., the ionic potential is neglected during the propagation of the photoelectron. While such an
approximation has been shown to work well in certain cases (e.g., [34]), however, for a streaking
time delay of the order of a few tens of attoseconds, which relates to a momentum shift of theorder of 10−3 a.u. (for a typical 800 nm streaking pulse with intensity of 1012 W/cm2 ), a more
careful method that also takes account of the ionic potential as well as the interplay between these
two fields is of course necessary. In this subsection, we calculate the propagation (or coupling)time delay in a streaking experiment using classical approaches that take account of both fields
and then compare this calculated propagation time delay with the streaking time delay obtained in
TDSE simulations. The comparison will lead us a new interpretation of the streaking time delay
for single-photon ionization.
It is known that in strong-field physics the propagation of an electron in the continuum canbe often well described by classical analysis (e.g., [82]). We therefore make use of this method
to analyze the time delay obtained in the streaking experiment. After the transition into thecontinuum due to XUV photon absorption the dynamics of the electron in the general 3D case isgiven bydk
= −Es (t) − ∇V (r),dt

(3.29)

where k is the momentum of the photoelectron, Es (t) is the streaking field, and V (r) is the atomic
potential. For a linearly polarized field the simultaneous interaction of the electron with the field
and the potential takes effect along the direction of the polarization, which we choose as the z-axis:dVdkz= −Es (t) −.
dtdz(3.30)

65
One approach to solve Eq. (3.30) is to use perturbation theory by assuming |Es |  | dV
dz |. To
first order, the final asymptotic momentum is then given by [assuming A(t) = 0 and V (z) = 0 fort → ∞]
(1)

(0)

(1)

kf,z (ti ) = kf,z (ti ) + ∆kf,z (ti )Z ∞Es (t)V z (0) (t)V (zi )= ki,z − As (ti ) +
+
dt
(0)
ki,z
ti
kz (t)2
Z ∞Es (t)V z (0) (t)dt,' k0,z − As (ti ) +(0)tikz (t)2
where k0,z =(3.31)q
2 + 2V (z ) ' kki,zii,z + V (zi )/ki,z is the streaking-field-free asymptotic momentum,
Z t(0)z (t) = zi + [ki − As (ti )] (t − ti ) +As (t0 )dt0 ,(3.32)tiandkz(0) (t) = ki,z − As (ti ) + As (t),(3.33)

are the zeroth order solutions for the position and momentum of the electron. ti , zi , and ki,z are
the initial time, position, and momentum of the electron after its transition into the continuum.
We note that Eq. (3.31) is equivalent to Eq. (25) in Ref. [83], in which the authors used
a quantum approach to calculate the momentum shift induced by the so-called Coulomb-laser
coupling effect in laser-assisted photoionization based on the eikonal approximation. This result
was later adopted in Ref. [38] to study the influence of the streaking field on the measured timedelay in a streaking experiment. However, one may notice that, similar to what we have done
here, the formula derived in Ref. [83] was based on the assumption that the laser field is rather
strong (e.g., 1.5 × 1014 W/cm2 in one of the numerical examples in Ref. [83]), which is usually notapplicable for a streaking experiment (streaking intensities are typically in the range of 1 × 1010
W/cm2 to 1 × 1012 W/cm2 ). Below we therefore compare results of the perturbative approach with
the numerical solution for the classical electron dynamics in order to test the applicability of theperturbative approach.
The use of perturbation theory can be avoided by integrating Eq. (3.30) directly. By multiplying dz to both sides of Eq. (3.30) and then integrating it, the solution of Eq. (3.30) for the

66
asymptotic momentum of the electron at z → ∞ can be written assZ∞

2 + 2V (z ) − 2ki,zikf,z (ti ) =

s
2 −2
k0,z

=

Es (t)kz (t)dtti

Z ∞Es (t)kz (t)dt.

(3.34)

ti

Making use of this equation requires knowing kz (t), which is available by numerically solving
Eq. (3.30). As we will see later, writing the final momentum in this form is very helpful forinterpretation of the propagation time delay.
To relate the final momentum of the electron to the propagation time delay, we set the final
momentum from Eq. (3.31) or (3.34) equal tokf,z (ti ) = k0,z − αAs (ti + ∆tprop )(3.35)' k0,z − αAs (ti ) + αEs (ti )∆tprop
as used in the fitting of the streaking results. Then the propagation time delay that the electron
accumulates in the combined field of the streaking field and the atomic potential is finally given by[for Es (ti ) 6= 0]∆tprop,pert 'R T Es (t)V (z (0) (t))
dt(α − 1)As (ti ) + ti s
(0)2[kz (t)]αEs (ti )for the perturbative approach, andq
R
2 − 2 Ts E (t)k (t)dt − kk0,zsz0,z + αAs (ti )ti∆tprop,exact 'αEs (ti )(3.36)(3.37)

for the exact solution. One may note that in both formulas the upper limit of the integrals has
been replaced by Ts , i.e., the pulse duration of the streaking field, which is valid since Es (t) = 0
for t > Ts .
In order to show the accuracy of the predictions of the two classical approaches and also
the relation between the propagation time delay and the streaking time delay, we compare in
Fig. 3.11 the results of quantum streaking simulations (black circles) for the streaking time delay∆ts with the classical predictions, Eq. (3.36) (blue dashed lines) and Eq. (3.37) (red solid lines), forphotoionization of an electron in (a) the 1D Coulomb potential VC (z), Eq. (3.20), with Z = 3.0 and

67

0

−18

(a)

(b)

1D TDSE
perturbative approachnumerical solution−20−22
time delay (as)time delay (as)

−20

−40

−60

−80

20

40

60
80
100
final electron energy (eV)

−26
−28
−30
−32

1D TDSE
perturbative approachnumerical solution

−100

−24

−34
−36
0

120

200

400
600
800
1000
Gaussian position z0 (a.u.)

1200

0
(c)
−10

H

time delay (as)

−20

He+

−30
−40
−50
−60
3D TDSE
perturbative approachnumerical solution

−70
−80
20

40

60
80
100
final electron energy (eV)120

Figure 3.11: Comparison of streaking time delays from quantum streaking simulations (black circles) with classical results from the perturbative approach (blue dashed lines) and the full numericalsolution (red solid lines). In our analysis we considered three potentials: (a) 1D Coulomb potential[VC (z)], (b) the combination of 1D Coulomb and Gaussian potentials [VCG (z)], (c) 3D Coulomb
potential [V (r)].

68
a = 2.0, (b) the combination of 1D Coulomb and Gaussian potentials (the electron was initiallybound in the ground state of the Coulomb potential)“|z|−z0Z−σ
V (z) = VCG (z) = − √+ V0 ez2 + a”2,(3.38)

with Z = 3.0 and a = 2.0 for the Coulomb potential, and V0 = −0.5, σ = 2.0 for the Gaussian
potential centered at z0 , and (c) a 3D Coulomb potential of the formZ
V (r) = − ,r(3.39)

with Z = 1.0 and 2.0 for H atom and He ion, respectively.
The streaking time delays ∆ts in the 1D quantum simulations are obtained by scanning thecentral cycle of the streaking field and then fitting the trace using Eq. (3.28). To ionize the electron
initially bounded in the 1D Coulomb potential, we have used an XUV pulse with IXUV = 1 × 1015
W/cm2 , ωXUV = 100 eV, TXUV = 600 as, and φXUV = −π/2. The parameters of the streaking pulse
were Is = 1 × 1012 W/cm2 , λs = 800 nm, Ns = 3 [Fig. 3.11(a)] or 8 [Fig. 3.11(b)], and φs = −π/2.
The 3D quantum results in Fig. 3.11(c) were extracted from Ref. [37]. In the classical calculations
we assumed a transition to the continuum at the center of the XUV pulse (i.e., ti = τ = TXUV /2).
It has been previously shown (e.g., [37, 38]) that the results of classical streaking simulations
depend on the choice of the initial position zi . We have chosen zi to be the most probable positionof the electron in the initial state [e.g., xi = 0 for VCG (z)]. Alternatively, one can sample theinitial conditions in Monte-Carlo calculations (see, e.g., [37]). To make use of the classical analysis
we further note that the propagation time delay ∆tprop in Eqs. (3.36) and (3.37) depends on the
choice of α. We determine α such that ∆tprop remains approximately constant while varying ti over
one field cycle (solid line with squares in Fig. 3.12). Please note that, independent of the choiceof α, our classical prediction for ∆tprop diverges for Es = 0. In practice, we can actually extract
the propagation time delay ∆tprop directly from the classical streaking trace by using the fitting
method. We have tested that this alternative method leads to the same time delays as those using
the time delay formula Eq. (3.36) or (3.37).

690

α = 1.025
α = 1.035
α = 1.045

time delay (a.u.)

−0.5

−1

−1.5

−2
−60

−40
−20
0
20
40
60
delay ti between XUV and streaking pulses (a.u.)

Figure 3.12: Classical predictions for the time delay as a function of the delay ti between XUVionizing and IR streaking pulses for different fitting parameter α by using Eq. (3.37). Calculations
are performed for the 1D potential VCG (z) in Eq. (3.38) with Z = 3.0, a = 2.0, V0 = −0.5, σ = 2.0,
and z0 = 20.

0.3

relative difference

0.25

0.2

0.15

0.1
Coulomb
Coulomb−Gaussian
0.05 11
10

12

13

14

10
10
10
streaking pulse intensity (W/cm2)

15

10

Figure 3.13: Relative difference between momentum shift ∆k calculated classically from theperturbative approach and the full numerical solution for VC (z) (blue solid line with circles)
and VCG (z) (z0 = 200, green solid line with stars). The relative difference is defined as|(∆k (pert) − ∆k (numerical) )/∆k (numerical) |. We have considered photoemission of an electron withfinal asymptotic momentum of 2.0, which is streaked by a 3-cycle 800 nm laser pulse.

70
The comparison in Fig. 3.11 reveals that the results of the quantum streaking simulations(black circles) and the full classical calculations (red solid lines) are in excellent agreement with eachother. In the present set of calculations the difference between the quantum and full classical results
does not exceed 2 as, which shows that the numerical streaking simulations can be well analyzed andinterpreted using a classical approach. Such a conclusion has also been reached based on the resultsof classical trajectory Monte Carlo calculations in Ref. [37]. This indicates that the transition of
the electron from the bound state into the continuum occurs instantaneously upon application of
the XUV pulse (i.e., most likely at the center of the XUV pulse) and the observed temporal offset∆ts , or time delay, arises due to the propagation of the photoelectron in the continuum.
On the other hand the comparison further shows that the results of the perturbative classicalapproach (blue dashed lines) agree neither with those of the quantum simulations (black circles)
nor with the full classical calculations (red solid lines). We find discrepancies of more than 30 as ascompared to the quantum results and the full classical results in some of our present calculations[e.g., the low energy part in Fig. 3.11(a)]. This degree of deviation between the two classicalapproaches persists for intensities of the streaking field up to 1014 W/cm2 . We exemplify thisby showing in Fig. 3.13 the relative difference between the results for the momentum shift ∆k,obtained from the two classical approaches for both VC (z) (blue solid line with circles) and VCG (z)(green solid line with stars). We note that, in general, results of the perturbative approach and thefull numerical solution only tend to agree for high streaking field intensities, which would induce
ionization from the target by the streaking field itself.
We can therefore further conclude that the discrepancy between the results of the perturbativeclassical approach and the quantum streaking simulations does arise from the first-order approximation of the Coulomb potential, but not from the use of classical theory itself. The transitionof the electron from the bound state into the continuum in photoionization occurs at distances atwhich the perturbation condition |Es |  | dV
dz | is not fulfilled for moderate streaking field strengths.
This conclusion is further supported by the results for the combination of Coulomb and Gaussian
potentials [VCG (z)]. In this case the final time delay has two contributions, one resulting from the71
Coulomb potential and the other from the rather weak Gaussian potential. In test calculations
we studied both contributions independently and found that in the present results the deviation
between the two classical results arises from the Coulomb potential alone. This can be also seenfrom the results in Fig. 3.11(b) as the difference between the two classical results is independent of
the position of the Gaussian potential.
Before proceeding, we may briefly discuss the option to extend the perturbative classicalresult for the analysis of streaking calculations. The deviations found for the first-order approximation can be significant. It remains to be studied if higher-order corrections can improve theresults sufficiently or a non-perturbative treatment becomes necessary. Instead, it might be moreinteresting to partition the space into an inner region close to the nucleus and an outer region, in
which |Es |  | dV
dz | is fulfilled for the outer one. Then, an application of the perturbative result
for the time delay related to the outer region appears to be satisfied. We may point out that the(perturbative) eikonal approximation was indeed initially introduced in the context of strong-field
ionization [84]. For strong fields, in contrast to weak-field photoionization, the perturbative condition is usually well satisfied for distances beyond the tunnel exit, at which the electron enters thecontinuum. Results of previous theoretical calculations [58] indicate that it might be then useful
to approximate the time delay in the inner region by the WS time delay for the short-range partof the potential. Of course, the quality of such an approximation should depend on the potential
of interest, the streaking field strength and other parameters, which would determine the partition
between inner and outer region and the range of short- vs. long-range parts of the potential. We do
not further investigate this option, since the full classical solutions are in excellent agreement withthe present quantum streaking results and therefore provide a good starting point for our further
analysis.

3.2.3

Physical interpretation of streaking time delay in single-photon ionization

We have just shown [81, 85, 86] that for the single-photon ionization process the streaking
time delay only arises from the propagation of the photoelectron in the combined potential of the72
streaking and the atomic fields, i.e.,∆ts = ∆tprop .

(3.40)

This implies that there is no contribution of the WS time delay ∆tWS in the streaking time delay∆ts , which however does not agree with previous interpretations that claim the streaking experimentis related to an intrinsic ∆tWS . In order to seek for an alternative physical interpretation for ∆ts
and also see whether or not ∆ts can be related to ∆tWS , we make use of the relation in Eq. (3.40)
and further simplify the accurate propagation time delay formula, Eq. (3.37). Noting that thesecond term of the square root in Eq. (3.37) is usually small, we can expand the square root to firstorder. Assuming α = 1, we obtain1∆ts = ∆tprop 'Es (ti )

Z Ts
ti



kz (t)dt,
Es (t) 1 −k0,z(3.41)

which still provides accurate results for liberation of the electron at the peak of Es (t) within its
central cycle as found in test calculations.
To avoid evaluating kz (t) by numerically solving the Newton’s equation, we approximate itusing
kz (t) ' kz(0) (t) − As (ti ) + As (t),(3.42)(0)where kz (t) is the momentum of the electron in the field-free atomic potential. This approximation
is based on the assumption that the coupling effect of the streaking field and the atomic potential
on the momentum is negligible and its accuracy can be seen from results in Fig. 3.14. Within thisapproximation, Eq. (3.41) then yields#(0)kz (t)As (ti )2
Es (t) 1 −
dt +
k0,z
2Es (ti )k0,z
ti
"
#
Z Ts
(0)
kz (t)
1'
Es (t) 1 −dt,
Es (ti ) ti
k0,z

1
∆ts '
Es (ti )Z Ts"(3.43)

where in the second step we have assumed that the electron is liberated at the peak of Es (t),i.e., at As (ti ) = 0. It is now instructive to further rewrite Eq. (3.43) as a sum by assuming thatthe streaking field and the electron momentum are approximately constant in the time interval

73
3k(0)(t)k(0)(t)−A(t )+A(t)i2.8momentum (a.u.)

k(t)
2.6
2.4
2.2
2
1.8
0

50

100
time (a.u.)

150

200

Figure 3.14: Comparison of momenta calculated using different methods: full numerical solution of
the Newton’s equation (red solid line) and Eq. (3.42) (green dashed line). We have considered anelectron propagating in the combined potential of the 1D Coulomb potential VC (z) and a 3-cycle
800 nm streaking pulse with Is = 1 × 1012 W/cm2 . As a reference the field-free momentum k (0) (t)of an electron propagating in VC (z) is also shown as blue solid line.

(0)

(0)

[tj , tj + δt], i.e., Es (t) ' Es (tj ) and kz (t) ' kz (tj ),
"
#
N
(0)
kz (tj )1 X
Es (tj ) 1 −δt∆ts 'Es (ti )
k0
j=1

=

N
X
Es (tj )j=1

Es (ti )(j)∆tfield-free ,(3.44)(j)where ∆tfield-free is a piecewise field-free time delay that the electron accumulates during its propagation in the time interval [tj , tj + δt] and over a related finite region [zj , zj + δz] of the potential
V (z) as compared to the propagation of a free particle over the same distance in space.
Eq. (3.44) provides us with an interesting interpretation of the streaking time delay: It isneither the WS time delay ∆tWS nor the simple sum of finite-range piecewise field-free time delays.
Instead, the piecewise field-free time delays are weighted by the streaking field strength present
when the electron wavepacket propagates over the corresponding part of the potential. Thus, we
can conclude that the streaking time delay strongly depends on the electron dynamics in the coupledatomic and time-varying streaking field potential.
Eq. (3.44) further indicates that in certain theoretical limits ∆ts can approach the WS timedelay ∆tWS . In particular, for short-range potentials V (r), ∆ts is approximately given by the first74

0
streaking delay
WS delay−50
time delay (as)time delay (as)

−20

−25

−30

zp=50

−100
−150
−200

−35
−250
−40

(a)
20

40
60
80
potential range zp (a.u.)time delay (as)

−18

100

−300
60

streaking delay
WS delay

(b)

80
100
120
XUV central frequency (eV)140

streaking delay w/ zp=20−20

streaking delay w/ zp=50−22

WS delay

−24
−26
−28
−30
−32
(c)
−34

500

1000
1500
2000
2500
streaking pulse wavelength (nm)3000

Figure 3.15: Results of numerical simulations for the streaking time delays ∆ts (solid line withcircles and diamonds) as a function of (a) the range xp of a short-range potential, (b) the frequencyof the ionizing XUV field, and (c) the wavelength of the streaking pulse are compared with thosefor the Wigner-Smith time delay [dashed lines with squares in (a) and (b) and solid circles in (c)].
Laser parameters are: IXUV = 1 × 1015 W/cm2 , TXUV = 600 as, ωXUV = 100 eV [(a) and (c)],φXUV = −π/2, Is = 1 × 1012 W/cm2 , Ns = 3 cycle [(a) and (b)], Ts = 32.02 fs [(c)], λs = 800 nm[(a) and (b)], and φs = −π/2.

75term (i.e., j = 1) of the sum in Eq. (3.44) if the electron wavepacket propagates over the full rangeof the potential during [ti , ti + δt] while Es (t) ' Es (ti ). This condition should be fulfilled in thefollowing theoretical limits: (a) the effective range of V (r) goes to zero, (b) the momentum of theelectron goes to infinity, or (c) the oscillation period of the streaking field goes to infinity. To test
our expectations we performed simulations for the following 1D potential:
V (z) = VC-WS (z) = − √Z1z 2 + a 1 + e(|z|−zp )/b,(3.45)

which is a product of Coulomb and Woods-Saxon potentials, where zp determines the effective
range of the potential. We have chosen b = 1 and zp > 10 such that VC-WS (z) ' VC (z) for |z| < zp ,while VC-WS (z) approaches zero quickly for |z| > zp . To obtain the WS time delays, we used the
back-propagation method introduced in the last section. In Fig. 3.15 we compare the results of
the numerical simulations for ∆ts (solid lines with circles and diamonds) as a function of (a) the
potential range zp , (b) the frequency of the ionizing XUV pulse, and (c) the wavelength of thestreaking pulse with the WS time delay ∆tWS (dashed lines with squares for (a) and (b) and solid
circles for (c)). All other parameters of the fields are kept the same in the simulations. As expected,∆ts approaches ∆tWS in each of the three theoretical limits listed above. One may notice that forlimit (a) a quantum proof based on the eikonal approximation was given in Ref. [38], in which theauthors obtained a similar conclusion that the propagation (or coupling) time delay is equal to the
WS time delay for a short-range potential.
As mentioned above, in earlier work (e.g, [38, 53, 58]) for the long-range Coulomb interaction
the streaking time delay was often separated into the sum of the field-free WS time delay anda contribution accounting for the coupling. Although we do not derive such a relation, we notethat in Eq. (3.44) Es (t1 ) ' Es (ti ) and therefore the first term of the sum is equal to the field-freetime delay, which the electron acquires during [ti , ti + δt]. The corresponding ’field-free’ delay
distance over which the electron travels during this time interval depends on the electron energy(or, XUV frequency) and the wavelength of the streaking field. In Fig. 3.16 we present classicalestimations for this distance by assuming that the electron is emitted at the peak of a 3-cycle 800

76
field−free delay distance (a.u.)

20
18

Es(t) changes 1%
Es(t) changes 5%

16
14
12
10
8
6
60

80
100
120
XUV central frequency (eV)140

Figure 3.16: Classical estimates for the ’field-free’ delay distance of an electron, released at thepeak of a three-cycle streaking pulse (wavelength of 800 nm) in the 1D Coulomb potential, as afunction of the XUV photon energy. This distance was determined by the position of the electron
at the time instant when the streaking field changed to 99% (solid line) and 95% (dashed line) of
the peak field strength.

nm laser pulse. We note that the calculated ’field-free’ delay distances are approximately equal toshort-range distances of atomic potentials.
The sum rule, ∆ts = ∆tWS + ∆tprop , used for the interpretation of the streaking time delays
appears to work for many systems and to some extent is related to our present interpretation,
Eq. (3.44). However, we want to emphasize that our interpretation is essentially different from theprevious one: We have demonstrated that the streaking time delay only arises from the propagationof the electron in the combined continuum of the streaking and ionic potentials, or alternatively
the coupling effect, while previously it was believed that an intrinsic WS time delay needs to be
added.
With our present interpretation, we are also able to answer the question whether or not fora short-range potential there is a contribution from the field coupling to the streaking time delay.
As mentioned above, consistent with our results, Ivanov et al. [38] showed that for a short-range
potential (e.g., the short-range part of the neon potential) there exists a propagation (coupling)time delay and it is equal to the WS time delay by applying the eikonal approximation. However,in Ref. [53], Nagele et al. pointed out that this conclusion cannot be supported by numerical

77
simulations when comparing ∆ts with ∆tWS + ∆tprop . The latter quantity was shown to be twiceas large as the streaking time delay ∆ts [53]. Based on the sum rule, it was claimed by Nagele
et al. [53] that ∆tprop = 0, i.e., no contribution due to the field coupling exists for a short-range
potential. From our analysis, it is obvious that this incorrect conclusion in Ref. [53] arises due tothe use of the sum rule, which is only an assumption and has never been rigorously proven. For
such a short-range potential considered in Refs. [38, 53], to reproduce the streaking time delay ∆ts ,one can use either the propagation time delay ∆tprop or the WS time delay ∆tWS , but never bothof them because the physical origin of ∆tWS is the propagation (or coupling) effect, which however
approaches the WS time delay for a short-range potential [see Fig. 3.15(a)].
Before proceding, we note that more subtle features in our results, such as the oscillations in
Fig. 3.11(b) and the change from a decrease to an increase in ∆ts as a function of zp in Fig. 3.15(a)can also be well understood with our classical interpretation, Eq. (3.44). The oscillations in Fig.
3.11(b) are caused by the Gaussian potential, which has a very short effective range. The streaking
delay contribution due to the Gaussian potential can be therefore written as∆tGs (z0 ) 'Es (tG ) G∆t ,Es (ti ) WS

(3.46)

where ti and tG are the instant of ionization and the time instant at which the electron reaches
the Gaussian potential respectively. Thus, the streaking time delay induced by the coupling of thestreaking field and the Gaussian potential should have an oscillation period proportional to thatof the streaking field, which can be easily confirmed from the results in Fig. 3.11(b) by notingthat tG ' z0 /k0 . On the other hand, the non-monotonic behavior of the streaking time delay inFig. 3.15(a) can be also explained by Eq. (3.43). Assuming that the electron is ionized at the peak of
the streaking field, the absolute value of the time delay is expected to increase over the first quarterof the streaking field cycle. However, then the field changes sign and the absolute value should start
decreasing. The turning point at zp = 50 agrees well with the classical estimate for the distancethe electron travels within the first quarter cycle after its release, namely zwp = k0 Ts /4 = 54 at thepresent parameters.

78
3.2.4

Finite-range time delay and its application to imaging

Up to now we have established a new physical interpretation for the time delay measured
in the streaking experiment, which explains the streaking time delay as a sum of field-weighted
field-free time delays. As pointed out above, one important consequence of this interpretation isthat, in general, the WS time delay does not present or contribute to the streaking time delay.
From Eq. (3.43), one clear sees that the coupling effect is present during the time interval of [ti , Ts ],i.e., between the ionization time and the end of the streaking field, which is in contrast to theoriginal definition of the WS time delay in which the propagation region extends to infinity. In
order to further support such a conclusion, we present results of numerical simulations of streaking
scenarios for the Coulomb-Gaussian potential, Eq. (3.38), in which the position of the Gaussian
potential can be changed to test the effective region that affects the streaking time delay.
We consider photoemission from the ground state of the 1D Coulomb potential VCG (z) with
Z = 3.0 and a = 2.0 by an XUV pulse with IXUV = 1 × 1015 W/cm2 , ωXUV = 100 eV, TXUV = 600
as, and φXUV = −π/2. The parameters of the streaking pulse are Is = 1 × 1012 W/cm2 , λs = 800
nm, and φs = −π/2. In Fig. 3.17(a), we show the extracted values for the streaking time delay∆ts (stars, squares, and circles) for different positions z0 of the Gaussian potential and a streaking
pulse having 8 cycles. The streaking time delay ∆ts strongly varies when the Gaussian potential islocated close to the center of the Coulomb potential and the amplitude of this variation increases
with an increase of the depth of the Gaussian potential. However, independent of the strengthof the Gaussian potential, ∆ts remains constant for z0 > zfinite ' 850, for which we obtain thesame numerical result with and without (diamond) the Gaussian potential. Our simulations showthat at the end of the streaking pulse the center of the outgoing electron wave packet is located atzfinite ' 850. Thus, ∆ts accounts for the presence of the Gaussian potential only when the electron
wave packet reaches the potential before the interaction with the streaking pulse ceases. This is, as
expected, in contrast with the WS time delay that accounts for the potential up to infinity basedon its original definition.

79

0

TDSE for V0=−0.5time delay (as)−5

TDSE for V0=−2.0−10TDSE for V0=−4.0−15

TDSE for V0=0−20numerical solution

−25
−30
−35
−40

(a)

−45
0

200
400
600
800
1000
Gaussian potential center z0 (a.u.)

1200

−24

time delay (as)

−24.5
−25
−25.5
−26
TDSE for V0=0

−26.5

TDSE for V0=−0.5
−27

4

6
8
10
12
14
streaking pulse cycle number(b)
16

Figure 3.17: (a) Streaking time delay ∆ts as a function of the position z0 of the Gaussian potential.
We compare results of TDSE simulations (stars, squares, and circles for V0 = −0.5, −2.0, and −4.0,
respectively) with those of classical calculations, Eq. (3.37) (red solid lines). Also shown is the
TDSE result without Gaussian potential (diamond on the right end). In the classical calculations
the parameter α varies between 0.985 and 1.084. (b) Comparison of ∆ts from TDSE calculations fordifferent cycle numbers of the streaking pulse with (green solid line with circles) and without (bluedashed line with squares) the Gaussian potential. The Gaussian potential is located at z0 = 650.
Other parameters are given in the text.

80
Table 3.3: Streaking time delays for different time intervals [t1 , t2 ) of the streaking trace. An 8-cycle
streaking pulse and two 1D potentials [VC (z) and VCG (z) (z0 = 650)], have been used for the TDSE
simulations. All values in this table are given in a.u..
[t1 , t2 )[−400, −200)[−200, 0)[0, 200)[200, 400]

Coulomb−0.9867−0.9985−1.0085−1.0250

Coulomb-Gaussian−1.3855−1.0765−1.0130−1.0250

To further confirm that it is the propagation distance of the electron during the presenceof the streaking field that determines the probing distance of the attosecond streak camera, we
determine ∆ts for streaking pulses having different durations while the Gaussian potential remains
located at z0 = 650. For short streaking pulses (Ns 6 5) we again obtain the same results for ∆ts
in simulations with [Fig. 3.17(b), green solid line with circles] and without (blue dashed line withsquares) the Gaussian potential. In each of these simulations the outgoing wave packet is locatedat z < z0 = 650 at the end of the streaking pulse. In contrast, as soon as the outgoing wave packet
reaches z0 for N > 5, ∆ts deviates from the result obtained in simulations without the Gaussian
potential.
According to our analysis the streaking time delay ∆ts depends on the time interval between
transition of the electron into the continuum, ti , and the end of the streaking pulse, Ts . We therefore
expect that for a long streaking pulse the time delay, extracted from different cycles of the streaking
pulse, varies, since for emission of the electron by the XUV pulse in the raising part of the streaking
pulse it will propagate over a larger distance until the streaking pulse ceases as compared to thecase when the electron is released in the trailing part of the streaking pulse.
To test our expectations we have performed quantum streaking simulations for the 1D
Coulomb potential VC (z), Eq. (3.20), and the combination of 1D Coulomb and Gaussian potentials
VCG (z), Eq. (3.38), using an 8-cycle streaking pulse. We then extract the time delay ∆ts by fitting
the streaked momentum over different intervals of the streaking pulse using Eq. (3.28). The results
are shown in Table 3.3, where negative times in the interval [t1 , t2 ) correspond to the raising part

81of the pulse, while positive times correspond to the trailing part of the pulse.
As expected, the extracted streaking time delay depends on the interval used for the fitting
procedure. In the case of the 1D Coulomb potential the variation of the time delays over the
streaking pulse is rather small. This indicates that the major part of the time delay is accumulated
shortly after the transition of the electron into the continuum near to the nucleus, where thepotential is strongest. In contrast for the combined Coulomb-Gaussian potential we see a strongchange of the time delay for application of the XUV pulse at the beginning of the streaking pulse
as compared to the other fitting intervals. This can be understood within our interpretation of thetime delay. If the electron is ionized early, it reaches the Gaussian potential during its propagation
at large streaking pulse field strengths. According to Eq. (3.44), this part of the potential therefore
contributes significantly to the accumulated time delay. On the other hand, if the electron is ionizedlater, it reaches the Gaussian potential either when the streaking pulse is weak or after the endof the streaking pulse (i.e., for the interval [200, 400]). In these cases the effect of the Gaussian
potential on ∆ts is small or not present at all, in agreement with Eq. (3.44). The results therefore
again confirm the finite-range property of the observed time delays. Furthermore, they actually
open a possibility to use the streaking technique to image the presence and even the position of anadditional potential (here, the Gaussian potential) within a long streaking pulse.
Previously, we have considered 1D examples which already reveal important aspects of thephysics behind the streaking measurement technique. We will now investigate the role of thepolarization direction in view of the finite-range property of the time delay and the detection of astatic potential at a distance from the location of the photoemission. To this end we consider thefollowing 2D potentials

V1 (x, y) = −√ Zr2 +a −√ Zr2 +aV2 (x, y) = − √where r =− V0,1 e−(r−r0 )2
σ2

for |x| 6 100,

(3.47)

for |x| > 100,(r−r0 )2 y 2
Z− V0,2 e− σ2,r2r2 + a(3.48)

px2 + y 2 . Both potentials are shown in Fig. 3.18 for Z = 2.0, a = 0.164, V0,1 = 0.5,

82

Figure 3.18: 2D model potentials, as defined in Eq. (3.47) [(a)] and Eq. (3.48) [(b)], plotted on alogarithmic scale as log[−V (x, y)].

V0,2 = 0.25, r0 = 140, and σ = 5.0. In our calculations the electron is initially bounded in the
ground state of the Coulomb potential with an eigenenergy of −2.0. We then used an XUV pulse
with IXUV = 1 × 1014 W/cm2 , ωXUV = 68 eV, NXUV = 10, and φXUV = −π/2 and a 3-cycle 400
nm streaking pulse with Is = 1 × 1012 W/cm2 and φs = −π/2 to streak the momentum of thephotoelectron. The TDSE is solved using a space-time grid with ∆x = ∆y = 0.3, Nx = Ny = 5000,
and ∆t = 0.05 in 2D Cartesian coordinates. The polarization of ionizing and streaking pulses
are kept parallel, but the polarization direction with respect to the orientation of the potentials isvaried.
To obtain the streaking trace as a single curve as in Fig. 1.7(a), we first perform a 2D
Fourier transform of the ionizing part of the wave function at the end of each simulation for agiven time delay between the XUV and streaking pulses. This 2D momentum distribution is thenintegrated over a small opening angle (±5◦ ) with respect to the polarization direction of the coaligned ionizing and streaking pulses, since the streaking effect is expected to be strongest along
the polarization direction of the streaking pulse [37]. Next, we obtain the expectation value of theresulting momentum distribution for a given time delay τ . Finally, by repeating the calculations
for application of the XUV ionizing pulse over the central cycle of the streaking pulse we obtain
the desired streaking trace as a function of the time delay between the two pulses and determine

83
the temporal offsets by comparison with the vector potential as before in the 1D cases. We havepropagated all wave packets for different XUV-streaking delays to the same distance in space,
therefore the conclusions based on the qualitative behavior of the results presented below do notdepend on the particular choice of wave function used for calculating the momentum distributions.
In Fig. 3.19(a) we present streaking time delays as a function of the number of cycles inthe streaking pulse for the 2D potential defined in Eq. (3.47). We compare results obtained forthe full potential in which the pulses are either polarized in x- (green line with stars) or in ydirection (red line with squares) with those in which we neglect the additional Gaussian potential
at a distance from the Coulomb potential (blue line with circles). Due to the spherical symmetryof the pure 2D Coulomb potential, in the latter case the results are independent of the choice ofthe polarization direction and we therefore show the results obtained for polarization along the xdirection only. In contrast, with the additional Gaussian potential the results strongly depend on thepolarization direction. While the streaking time delays with and without additional potential agreewith each other for polarization of the streaking field along the x-direction, they start to deviate
from N = 6 on for polarization in y-direction. The latter behavior is similar to the 1D example
studied before. The difference in the results is due to the fact that the time delay is determinedby the photoelectron dynamics in the coupled field of the atomic potential and the streaking fieldalong the polarization direction. For short streaking pulses the photoelectron does not reach thelocation of the additional Gaussian potential and the effect of the latter is therefore negligible,
while for longer pulses the coupling effect between the streaking field and the additional Gaussian
potential becomes present. The additional model potential is present along the y-polarization but
not in x-direction in agreement with our observations for the different polarization directions of the
field.
Next, in Fig. 3.19(b) we show the results obtained for the second 2D potential, Eq. (3.48).
Time delays for the full potential and different polarizations of the co-aligned ionizing and streaking
pulses (polarization in x-direction: green line with stars, at 45◦ : cyan line with crosses, and in ydirection: red line with squares) are compared with those obtained without additional potential

84

−44

−47
streak in x w/ Coulomb
streak in x w/ V1(x, y)streak in y w/ V1(x, y)time delay (as)time delay (as)

−46

−48

−50

−52

−54
3

−48

streak in x w/ Coulomb
streak in x w/ V2(x, y)−49

streak at 45o w/ V2(x, y)streak in y w/ V2(x, y)

−50
−51
−52

(a)
4

5
6
7
8
streaking pulse cycle number

9

10

−53
3

(b)
4

5
6
7
8
streaking pulse cycle number

9

10

time delay (as)

−47
−48

streak in x w/ Coulomb
streak in x w/ V2(x, y)−49

streak at 45o w/ V2(x, y)streak in y w/ V2(x, y)

−50
−51
−52
−53
3

(c)
4

5
6
7
8
streaking pulse cycle number

9

10

Figure 3.19: Time delay as a function of streaking pulse cycle number for the 2D potentials definedin Eq. (3.47) [(a)] and Eq. (3.48) [(b) and (c)] and polarizations of the streaking (and co-aligned
XUV ionizing) field in x-direction (green line with stars), at 45◦ (cyan line with crosses), and in
y-direction (red line with squares). The results are compared with those for the pure 2D Coulomb
potential without the additional potential (blue line with circles) and a streaking field polarized in
x-direction. In (c), the original streaking time delays for the 2D potential, Eq. (3.48) as shown in(b), have been shifted to match the result for the shortest streaking pulse for the pure 2D Coulomb
potential.

85(blue line with circles). We note that in this case the time delays do not coincide with eachother, even for short streaking pulses, for which the photoelectron wave packet does not reach theadditional potential when the streaking pulse ceases. Further analysis reveals that this discrepancy
can be explained as due to the scattering of the photoelectron at the additional potential, whichleads to different final momentum distributions as compared to that of the pure 2D Coulomb
potential case. In order to remove this scattering effect, we shift all time delay curves such thatthe results match for the shortest streaking pulse. The resulting modified time delays, shown in
Fig. 3.19(c), then reveal again the finite-range property since the curves deviate for N > 5. We
furthermore observe that the degree of deviation increases as the additional potential along thepolarization direction gets stronger.

3.2.5

Influences of pedestal and additional static field

So far, we have considered streaking pulses with a sin2 -envelope, which is only an ideal
assumption for theoretical convenience. In practice the shape of the streaking pulse is usually not
regular and sometimes has an oscillating tail. According to the classical analysis, Eq. (3.44), thestreaking time delay however depends on the instantaneous field strength during the propagationof the electron in the continuum. It is therefore interesting to study if and how the shape of thestreaking field envelope influences the observed time delays. In order to study this aspect we makeuse of the following pulse envelope
Eenv (βp , Tp ) =
E0 sin2 (πt/Ts )for 0 ≤ t ≤ Ts /2,






 (1 − βp )E0 sin2 (πt/Ts ) + βp E0 cos2 [π(t − Ts /2)/(Ts + 2Tp )] for Ts /2 ≤ t ≤ Ts ,
βp E0 cos2 [π(t − Ts /2)/(Ts + 2Tp )]





 0

(3.49)

for Ts ≤ t ≤ Ts + Tp ,else.

As shown in Fig. 3.20 the corresponding field has a pedestal in the trailing part of the pulse, which
depends on the parameters βp and Tp . For βp → 0 and/or Tp → 0, Eenv equals the commonly used86
sin2 -form of the envelope with duration Ts .
−3

streaking field strength (a.u.)

6

x 10

with pedestal
without pedestal

4
2
0
−2
−4
−6
0

200

400

600
time (a.u.)

800

1000

Figure 3.20: Streaking field with a pedestal (red solid line): E0 = 5.34 × 10−3 (i.e., Is = 1 × 1012
W/cm2 ), Ts = 331, βp = 0.2, and Tp = 750. As a comparison a basic 3-cycle 800 nm streaking
pulse (blue dashed line) is also shown.

We perform quantum streaking simulations for different values of the parameters βp and Tp
for the 1D Coulomb potential VC (z) with effective charges, Z = 1.0 and Z = 3.0 respectively, anda streaking pulse with Ts = 331 at a wavelength of 800 nm and a peak intensity of Is = 1 × 1012
W/cm2 . To ionize the electron from the ground state of each potential, we have used an XUV pulse
with IXUV = 1 × 1015 W/cm2 , ωXUV = 60 eV (Z = 1.0) or 100 eV (Z = 3.0), NXUV = 10 (Z = 1.0)
or 15 (Z = 3.0), and φXUV = −π/2. In each of these simulations we extract the time delay by
applying the XUV over the interval [0, Ts ] of the corresponding pulse. The results are shown in
Fig. 3.21 [blue solid line with circles in (a) and (c), blue solid line with squares in (b) and (d)].
We observe that the time delays obtained for the pedestals deviate from those for the sin2 pulse shape (Tp = 0 or βp = 0) for both potentials. This is in agreement with our interpretation ofa finite-range time delay and the dependence on the instantaneous streaking field strength during
the propagation of the electron [Eq. (3.44)], since the deviation increases as βp or Tp increases.
Quantitatively, our results however further show that the deviations due to the presence of thepedestal are rather small [e.g., about 1% for Z = 3.0 and βp = 1.0, see Fig. 3.21(d)]. This is due tothe fact that the effect of the pedestal comes into play when the electron is at rather large distances

87−11.36−24.28(b)−11.37−24.29−11.38−24.3
time delay (as)time delay (as)(a)−11.39−11.4−11.41−11.42−24.31−24.32−24.33−24.34

Z=1.0 w/ Estatic=0−11.43−11.44
0

200

400
600
800 1000
pedestal length Tp (a.u.)

Z=3.0 w/ Estatic=1e−5−24.36
0

1200

−10.8

200

400
600
800 1000
pedestal length Tp (a.u.)−24.05−11−24.1
time delay (as)

−10.9

−11.1
−11.2
−11.3
−11.4

(d)

−24.15−24.2−24.25
−24.3

−11.5

−24.35
Z=3.0

Z=1.0
0

1200

−24
(c)

time delay (as)Z=3.0 w/ Estatic=0−24.35Z=1.0 w/ Estatic=1e−5

0.2

0.4
0.6
pedestal strength βp

0.8

1

−24.4
0

0.2

0.4
0.6
pedestal strength βp

0.8

1Figure 3.21: Streaking time delay as a function of (a) and (b) the pedestal length Tp (βp = 0.2),
and (c) and (d) the pedestal strength βp (Tp = 750). The left column shows the results for Z = 1.0,
while the right column is for Z = 3.0. The influence of an additional static field on the streaking
time delay is also present (green dashed lines with asterisks and crosses) in (a) and (b).

88
−10

time delay (as)

−11.4
−15
−11.45−11.5−11.55 −8
10

−20
Z=1.0

−6

−4

10

10

Z=3.0
−25 −8
10

−7

10

−6

−5

10
10
static field strength (a.u.)−4

10

Figure 3.22: Streaking time delays as a function of the strength of the additional static field for abasic 3-cycle 800 nm streaking field for two potentials: blue solid line with circles for Z = 1.0 andgreen solid line with squares for Z = 3.0.

from the nucleus, at which the Coulomb potential is weak. Thus, we expect that for an atomic-likesystem pedestals and other deviations from a sin2 or Gaussian streaking pulse shape should nothave a large effect on the observed time delays.
Finally, we also investigate the influence of additional static electric fields on the attosecond
streaking time delay. In an experiment the presence of such additional fields may be used in order
to direct the photoelectrons towards a detector. We therefore include an additional static field term
Estatic in Eq. (2.32) and performed quantum streaking simulations for streaking pulses with andwithout pedestal. As before, we consider two 1D Coulomb potentials VC (z) with different effective
charges of Z = 1.0 and Z = 3.0 respectively. From the results presented in Fig. 3.22 we conclude
that the presence of additional static electric fields up to field strengths of Estatic = 10−5 do nothave a significant effect on the extracted time delays. The same conclusion can be also drawn
for streaking pulses with a pedestal, which can be seen by comparing the green dashed lines (forEstatic = 10−5 , stars for Z = 1.0 and crosses for Z = 3.0) with the blue solid lines (for Estatic = 0,
circles for Z = 1.0 and squares for Z = 3.0) in Fig. 3.21 (a) and (b).

Chapter 4Time Delays in Two-Photon Ionization

The advances in laser technology, in particular the ability to generate more and more intense
ultrashort XUV pulses, have provided us the possibilities to study high-order nonlinear processes,
such as the two-photon ionization (TPI) in XUV photon energy regime. The TPI, especially
its angular momentum distribution, has been extensively investigated both experimentally andtheoretically recently (e.g., [87–90]). These investigations provide a lot of information about the TPI
process, including the phase shifts and the dipole matrix elements for individual partial waves [88,89]. To continue and also extend our previous work, we however want to study this process froma different perspective, namely time resolving the TPI. Similar to what has been done for singlephoton ionization (SPI), we propose to use the attosecond streaking technique to study the TPI
and retrieve the instant of transition into the continuum, which would be possible based on ourprevious understanding of the propagation time delay. This chapter is organized as follows: First,
we present numerical results of streaking time delays obtained for TPI from numerical simulationsof streaking scenarios. By properly accounting for the measurement induced time delay, i.e., thepropagation time delay, we are able to retrieve the absorption time delay, i.e., the instant ofionization with respect to the center of the XUV pulse, for both resonant and nonresonant TPI
as well as its dependence on XUV laser parameters. Next, to further support our interpretationof the absorption time delay, we use two independent studies, namely the trajectory analysis andthe phase analysis, to reproduce our numerical results for the absorption time delay. Finally, we
review some of our previous work on two-photon coherent control that may relate to the potential

90
applications of the TPI delay.

Time delays extracted in two-photon ionization using streak camera1

4.1

The application of the attosecond streaking technique in studying SPI process has alreadygiven us interesting time information about this process and also initiated many fruitful discussions
about the origin of the time delay measured in the experiment. We therefore want to further apply
this technique to some more complex ultrafast processes, such as the resonant and nonresonant TPI,
in which “resonant” means that the central frequency of the XUV ionizing pulse is in resonance withthe energy difference of the initial ground state and one intermediate state. We will first report ournumerical results of time delays extracted in simulations of the streaking experiment, from which we
identify an apparent difference of the streaking time delay dependence on the XUV pulse duration
for these two types of TPI (Fig. 4.1). We then use our previously proposed classical approaches tocalculate the propagation time delay, which needs to be subtracted from the streaking time delay
to obtain the real photonionization instant with respect to the center of the XUV pulse. We arethen finally able to show that the nonresonant TPI occurs at the center of the XUV pulse while
the resonant TPI occurs with a positive delay with respect to the center and this delay changes
linearly with the duration of the XUV pulse.

4.1.1

Streaking time delay

As before, we use the standard Crank-Nicolson method to numerically solve the TDSE on agrid in space and time, 2p∂Ψ(r, t)=+ V (r) + (EXUV (t) + Es (t)) · r Ψ(r, t),
i
∂t
2

(4.1)

where p is the momentum operator and E(t) = E0 sin2 (πt/T ) cos(ωt + φ)ẑ for both laser fields
with peak amplitude E0 , pulse duration T , central frequency ω, and carrier-envelope phase φ of
the respective field. In all simulations we propagate the wave function on the grid for a sufficiently
1

The results of this section are presented in J. Su, H. Ni, A. Jaroń-Becker, and A. Becker, Time Delays in
Two-Photon Ionization, submitted.

91
long time until both laser pulses cease and the ionizing wave packet can be clearly separated fromthe bound wave function. We also confirm that the outgoing wave packets stay on the grid anddo not reach its boundaries. To obtain the momentum distribution, we either project the ionizing
wave packet onto the continuum eigenstates of the potential, or perform a Fourier transform. For
the latter choice, we have propagated the wave packet far enough from the core part such that theerror due to the planewave approximation is negligible.
We simulate the two-photon ionization process dressed by a streaking field for a helium (He)atom in the 3D Cylindrical coordinates as well as a 1D Coulomb model potential. The latter
one allows us to systematically study the dependence of the measured time delays on the laserparameters. For the 3D He atom, due to its azimuthal symmetry, we are able to solve this problem
in the 2D (ρ, z) coordinates. With spatial steps of ∆ρ = ∆z = 0.2 and a time step of ∆t = 0.02, we
obtain eigenenergies of the three lowest states as E1s = −0.90, E2s = −0.16, and E2p = −0.13, using
the single-active-electron model potential introduced in Ref. [71]. To ionize the electron initiallybound in 1s state, we choose the central frequency of the XUV field to be ωXUV = 15.63 and 21.07
eV, with the latter one corresponding to the energy difference between 1s and 2p states. The 1D√
Coulomb potential is chosen to have the form of V (x) = −Z/ x2 + a, with Z = 3.0 and a = 0.15,
which gives the energies of the three lowest states as E0 = −5.32, E1 = −2.31, and E2 = −1.30,
for ∆x = 0.05 and ∆t = 0.01. These three states are chosen to be well-separated such that usually
only the first excited state is involved in the resonant TPI. Accordingly, we choose ωXUV = 95.62
and 81.81 eV, corresponding to nonresonant and resonant TPI processes. Other parameters for
the XUV pulse are the same for the two potentials: IXUV = 1.0 × 1014 (nonresonant) or 1.0 × 1013
W/cm2 (resonant), and φXUV = −π/2. To streak the photoelectron in the continuum, we haveused a 3-cycle 2400 nm streaking field with Is = 1.0 × 1011 W/cm2 and φs = −π/2.
Fig. 4.2 presents our numerical results of time delays extracted from streaking traces measured
in the streaking experiment [91]. In each subfigure we show the change of the streaking time delay as
a function of the XUV pulse duration for resonant (green solid lines with squares) and nonresonant(blue solid lines with circles) ionization. We clearly see that the changing trends of the two processes

92

Figure 4.1: Schematic diagram for obtaining time delays in two-photon ionization (TPI) using
attosecond streaking technique. An electron initially bound in the ground state |gi is ionized byabsorbing two photons from the XUV pulse either through a resonant state |ei (resonant TPI)or not (nonresonant TPI), and then streaked by the long wavelength field in the continuum. Bychanging the relative delay between the XUV and streaking pulses, one obtains a streaking trace,from which a temporal shift (or streaking time delay ∆ts ) can be extracted. After accounting forthe propagation time delay ∆tprop that the electron accumulates in the continuum, we are able toaccess the absorption time delay ∆tabs of the TPI process.

nonresonant TPI
resonant TPI
SPI

15
streaking time delay (a.u.)

streaking time delay (a.u.)

15

10

5

0

nonresonant TPI
resonant TPI
SPI

10

5

0

−5
−5

(b)

(a)
50

100
XUV pulse duration (a.u.)

150

100
150
200
XUV pulse duration (a.u.)

Figure 4.2: Streaking time delay as a function of XUV pulse duration for (a) 1D Coulomb and(b) 3D He potentials. Three ionization processes are considered: nonresonant TPI (blue solid lines
with circles), resonant TPI (green solid lines with squares), and one-photon ionization (red solidlines with asterisks). Laser parameters are given in the text.

93
are different: for the resonant one the time delay increases with XUV duration, however, in the caseof the nonresonant one it remains almost unchanged. Even though the streaking time delay may
have some other contributions, the significant difference between the trends of the two delay curves
already strongly indicates that the intermediate state plays an important role in the temporal
dynamics of the TPI process. It appears that, for resonant TPI, after absorption of the first
photon, the electron remains in the intermediate state for some time, before it absorbs the second
photon and then transfers into the continuum. As we will show below, for long XUV pulses (e.g.,
TXUV > 150 for 3D He atom), the time the electron spends in the intermediate state is proportional
to the XUV duration. Moreover, this difference between resonant and nonresonant TPI appears to
be an universal phenomenon, since it occurs in both 3D He atom and 1D Coulomb potential.
To clearly distinguish our new observations from those in previous well-studied SPI processthat also uses attosecond streaking technique, we also show the dependence of the streaking timedelay (red solid lines with asterisks) on the XUV duration for SPI that is induced by a 2-ωXUV(nonresonant frequency) XUV pulse. The intensity of the XUV pulse is set to be IXUV = 1.0 × 1014
W/cm2 for 1D Coulomb, and IXUV = 1.0 × 1012 W/cm2 for 3D He, both of which guarantee thatthe photoionization process occurs in the perturbative intensity regime. Similar to the nonresonant
TPI process, the SPI time delay remains almost unchanged as the XUV pulse duration increases.
More interestingly, for the 1D Coulomb potential [Fig. 4.2(a)] the streaking time delay of SPI isexactly the same as that of the nonresonant TPI when the duration of the XUV pulse is longenough (e.g., TXUV > 50). This implies that the two-photon absorption occurs instantaneously as
the one-photon absorption, which occurs most likely at the center of the XUV pulse as it has beenshown in Chapter 3. However, in the 3D case [Fig. 4.2(b)], these two time delays are very similar
but not exactly equivalent. This may relate to the fact that the photoelectrons have different
angular momenta in the continuum (p wave for SPI, and a mixture of s and d waves for TPI). It
is well-known that different final states in the continuum result in different streaking time delays(e.g., [37, 58]).

94
4.1.2

Propagation time delay

It has been well-studied that the propagation of the photoelectron in the combined field of
the streaking pulse and the atomic potential can induce an extra shift in its final momentum, which
results in an temporal shift (or time delay) in the obtained streaking trace (e.g., [37, 38, 81]). Indeed,
one needs to account for this propagation time delay in order to obtain the instant of photoemission
as the instant at which the electron is set free in the continuum after the photon absorption.
Therefore, to extract a real time delay of a photoabsorption process from a streaking experiment, itis necessary to properly account for the propagation time delay, which is actually a measurementinduced delay. As discussed before [37, 81, 85, 92] for the SPI process, this propagation time delay
can be often well-reproduced by classical trajectory calculations. Here we can use the same idea toobtain the propagation time delays for TPI, because once the electron is in the continuum it does
not matter if it is ionized via a one- or two-photon process. For our 1D Coulomb case, we use theclassical method introduced in Ref. [81], in which we assume the electron is initially at x = 0 andsimulate its propagation in the combined Coulomb-streaking field by solving 1D Newton’s equation
directly. By fitting the classical streaking trace to the original streaking formula, we obtain the
propagation time delays, which are shown as green solid lines with squares in Fig. 4.3(a) and (c).
For the 3D He atom case, we use the classical-trajectory Monte Carlo method introduced in Ref. [37]in order to take into account both the distribution of the initial position of the electron in the boundstate and the final shape of the ionized wave packets in the continuum. To be clear, we have sampled
the initial position of the electron from the probability distribution of the initial eigenstate, andchosen the distribution of the initial momentum according to the angular shape of the wave packet(either s or d wave). The propagation time delays from different waves are added up coherently
and are also shown as green solid lines with squares in Fig. 4.3(b) and (d). One may note that,for both 1D and 3D models, the propagation time delays have very weak dependence on the XUVduration, which is expected, since the final electron energy that determines the propagation timedelay barely changes with the XUV duration.

95

1

1
(a)

0
−1
time delay (a.u.)time delay (a.u.)

0

(b)

−1

−2

−2
−3
−4
−5

−3
∆ ts
−4

40

20

∆ ts

∆ tprop60
80
XUV pulse duation (a.u.)∆ tprop−6∆ tabs100

20∆ tabs∆ ts∆ tprop100

150
XUV pulse duation (a.u.)∆ tabs200∆ tabs15
time delay (a.u.)time delay (a.u.)∆ tprop100
150
XUV pulse duation (a.u.)

15
10
5
0

10

5

0

−5
−10

∆ ts

−7

(c)
50

100
XUV pulse duation (a.u.)

150

−5

(d)
200

Figure 4.3: Absorption time delay (red solid lines with asterisks) as a function of XUV pulse
duration for nonresonant TPI in 1D Coulomb [(a)] and 3D He [(b)] potentials, and resonant TPI
in 1D Coulomb [(c)] and 3D He [(d)] potentials. Also shown are the streaking time delays (bluesolid lines with circles) and the propagation time delays (green solid lines with squares).

96
Table 4.1: Time delays for different parameters of the streaking field. We consider resonant TPI in
the 1D Coulomb potential streaked by a 3-cycle laser pulse. XUV laser parameters are: ωXUV =81.81 eV, NXUV = 45, IXUV = 1.0 × 1013 W/cm2 , and φXUV = −π/2.
λs (nm)
2400
3200
4800
4800

4.1.3

Is (W/cm2 )1.0 × 1011
1.0 × 1011
1.0 × 1011
1.0 × 1010

∆ts (a.u.)
5.65
5.15
4.41
4.45

∆tabs (a.u.)
10.34
10.34
10.33
10.36

Absorption time delay

Now we are able to retrieve the absorption time delay (red solid lines with asterisks in Fig. 4.3),
a real time delay in the photoabsorption process, after properly accounting for the propagation timedelay. For long XUV pulses used in the nonresonant TPI, one sees that the absorption time delay
is exactly zero for the 1D simulations and relatively close to zero for the 3D simulations. Thisis consistant with our previous finding that the nonresonant TPI occurs instantaneously as SPI,
which means both processes occur at the center of the XUV, i.e., without time delays. In contrast,
the electron appears to take a nonzero time to absorb two photons before it is set free to thecontinuum in resonant TPI. Moreover, the photon absorption time changes linearly with the XUVduration once the pulse is long enough. One may explain this difference by understanding the
TPI process in three very simple steps: First, the electron in the ground state is excited to theexcited state by absorbing one photon from the XUV field. Then it spends some time staying in
the excited state. Finally, it absorbs another photon and is set free to the continuum. Since there
is no intermediate state for the electron to stay for nonresonant TPI, it has to absorb two photons
simultaneously, which thus makes no difference to the case of SPI. Here we have to point out thatusing the attosecond streaking technique we only have the time information when the electron is set
free in the continuum, i.e., the instant when the second photon is absorbed. The time information
about the absorption of the first photon is not available based on the present technique even for
the nonresonant TPI.
In order to further support our interpretation, it is of our interest to confirm that the ab-9720
2400+1e11
3200+1e11
4800+1e11
4800+1e10

time delay (a.u.)15∆ tabs

10
5

∆ ts

0
∆ tprop
−5
−10

50

100
XUV pulse duation (a.u.)

150

Figure 4.4: Absorption time delay (red dash-dotted lines) as a function of XUV pulse duration forresonant TPI in 1D Coulomb potential for different streaking parameters (diamonds: λs = 2400
nm and Is = 1.0 × 1012 W/cm2 , asterisks: λs = 3200 nm and Is = 1.0 × 1012 W/cm2 , squares:
λs = 4800 nm and Is = 1.0 × 1012 W/cm2 , and crosses: λs = 4800 nm and Is = 1.0 × 1011 W/cm2 ).
Also shown are the streaking time delays (blue solid lines) and the propagation time delays (green
dashed lines).

−5

1.6

x 10

TDSE
perturbation theory

prob in 1st excited state

1.4
1.2
1
0.8
0.6
0.4
0.2
0

40

60
80
XUV pulse duration (a.u.)

100

Figure 4.5: Excitation probability as a function of the duration of the XUV pulse for the nonresonant
TPI in the 1D Coulomb potential. We have used two methods to calculate the excitation probability
in the first excited state: the numerical solution of TDSE (black open circles) and the perturbation
theory (red solid line).

98
sorption time delay we obtained for the TPI process are independent of the measurement, i.e.,
remains unchanged as the streaking pulse parameters change. In Fig. 4.4, as expected, we see thatthe absorption time delays (red dashed lines) are identically for different streaking parameters, although the streaking time delays measured directly from the experiment as well as the propagation
time delays are quite different in particular when the wavelength of the streaking field changes.
This again implies that correctly accounting for the propagation time delay is crucial for obtaining
a solid absorption time delay. To clearly show this measurement independence numerically, we
further present values of the absorption time delays in Table 4.1 for NXUV = 45 (TXUV = 94.05),
from which one sees that the deviation is within 1 as, i.e., a relative difference of 0.3%.
Before we further analyze our results, we want to explain some subtle details in our results,
such as the small oscillation in Fig. 4.3(a) and the nonlinear behavior in Fig. 4.3(d) when the
duration of the XUV pulse is small. These behaviors can be related to excitations to the excitedstates. When the XUV pulse is short in time domain, it has a broad bandwidth in frequency domain,which means certain excited states within the energy bandwidth can be involved in the nonresonant
TPI. Thus, the nonresonant TPI has a significant resonant contribution as in Fig. 4.3(a) or morethan one excited state is involved in the resonant TPI. Each of these may result in a deviationfrom the expected dependence of the absorption time delay on the XUV pulse duration. As an
example, a significant excitation probability in the first excited state of the 1D Coulomb potential
is observed when the XUV pulse is short as shown in Fig. 4.5(a), where we have used results fromboth the TDSE solution and the perturbation theory. We further note that the oscillations in thisfigure are due to the oscillating structure of the frequency spectrum of the sin2 XUV field.

4.2

Alternative studies of absorption time delays in two-photon ionization
In the last section, we have shown how one can extract the absorption time delay of the TPI

process using the attosecond streak camera and have also studied its dependence on the durationof the XUV pulse. In particular, the linear dependence of the absorption time delay on the XUVpulse duration for the resonant TPI appears to be very useful in potential control applications due

99
to this simple linear relation. In order to further confirm our new findings and also provide somemore understanding of these new results, we make use of two different methods in this section,
namely the trajectory analysis and the phase derivative analysis, both of which have already been
introduced in the study of the SPI process in section 3.1.3.

4.2.1

Trajectory analysis

It has been shown that the trajectory of the photoelectron in the continuum contains usefultime information about the photoabsorption process (e.g., [43, 44, 46, 52]). As introduced before,a typical way to obtain the WS time delay for the SPI process from the trajectory is to extrapolate
the linear part of the trajectory back to the time axis, where the trajectory can be obtained bysolving the TDSE numerically. The interception of this extrapolated line with the time axis is thedesired time delay. However, as emphasized before, this approach only works for ionization from ashort-range potential, while for a long-range potential such as the Coulomb potential the nonlinear
component in the trajectory makes it impossible to fit the trajectory to a linear line even for thepart that is far away from the origin.
For the systems considered in our TPI analysis, e.g., the 1D Coulomb potential, it therefore
appears on the first view impossible to make use of the trajectory analysis to obtain the time delay
because of the long-range property of the Coulomb potential. However, by noting the fact that each
TPI trajectory has a corresponding SPI trajectory that has the same long-range behavior if they
have the same final energy, we can obtain the time delay of the TPI trajectory with respect to the
SPI trajectory by treating the latter one as a reference instead of the free particle trajectory. Then,the previously mentioned divergence problem of the WS time delay for the Coulomb potential does
not occur since the divergent contributions of the two trajectories cancel out. A similar idea has
been used to understand the physical meaning of the energy derivative of the Coulomb phase shift
in Ref. [93], where the quantum and classical scattering trajectories are compared.
To illustrate this idea, we consider the TPI and SPI processes in the 1D Coulomb potential
with Z = 3.0 and a = 0.15, the same that we used in the last section. As usual the TDSE that

100

2000
(a)

1600
1540

1400

1520

1800

<x> (a.u.)<x> (a.u.)

1800

2000
TPI
SPI original
SPI adjusted

TPI
SPI original
SPI adjusted

(b)

1600

1400

∆t

1200

1200
1500
1260 1280 1300 1320

1000
800

1000

1200
1400
time (a.u.)

1600

1800

1000
500

16
15

700

800
900
time (a.u.)

1000

1100

7
original
adjusted

6(c)original
adjusted(d)5time delay (a.u.)time delay (a.u.)

600

14
13
12

4
3
2
1

11
10
1000

0
1200

1400
1600
<x> (a.u.)

1800

2000

−1
1000

1200

1400
1600
<x> (a.u.)

1800

2000

Figure 4.6: Trajectory analysis of time delays in TPI. We have shown the trajectories of thephotoelectron for three processes [(a) and (b), blue lines for TPI, green lines for original SPI, andred lines for adjusted SPI] and the time delay as a function of the position of the photoelectron [(c)and (d), green lines for original SPI, red lines for adjusted SPI]. The trajectory is calculated using
the expectation value hxi of the electron position as a function of time by solving TDSE. The timedelay is obtained as the time difference for the same hxi between the TPI and SPI trajectories. We
have considered both resonant (left column) and nonresonant (right column) TPI. Here the original
SPI means the SPI process initiated by an XUV pulse with ωSPI = 2ωTPI while the adjusted SPI
stands for using an XUV pulse with ωSPI = 2ωTPI + ∆ω, which makes the final energy of the SPI
wave packet exactly the same as that of the TPI wave packet. The inset in (a) shows how we obtainthe time delay between the TPI and SPI trajectories.

101
describes the interaction of the electron with the XUV pulse in length gauge is numerically solved
on a space-time grid. When the XUV pulse ceases and the ionizing wave packet can be clearly
separated from the bound wave function, we are able to obtain the trajectory of the photoelectron,i.e., the expectation value hxi of the position of the ionizing wave packet as a function of time. In
order to generate the corresponding SPI trajectory for a TPI process, we use an XUV pulse with adoubled central frequency to ionize the electron from the ground state, i.e., ωSPI = 2ωTPI . We choose
the pulse duration of the XUV pulses for TPI and SPI to be equal, which yields NSPI = 2NTPI ,and the intensities of the XUV pulses are the same as given in the last section for the individual
cases. Once we have obtained the trajectories, the time delay between the TPI and SPI processes
is given by the time difference between the two trajectories for the same electron position hxi.
The trajectories of the TPI (blue lines) and SPI (green lines) are presented in the upper panelof Fig. 4.6 for (a) for the nonresonant TPI with NTPI = 45 (TTPI = 94.05) and (b) the resonant
TPI with NTPI = 50 (TTPI = 89.40). Since the time delay is small (a few a.u.) as compared to thetime scale (a few thousands a.u.) on the horizontal axis, we use a small inset in (a) to show thedefinition of the time delay between TPI and SPI. To analyze this time delay, we show the timedelay (green lines) as a function of electron position hxi in the lower panel of Fig. 4.6 with (c) and(d) for the resonant and nonresonant TPI, respectively. The time delay does not remain constant
as hxi changes, which means that measuring the time delay at different positions of the electron
gives with different results. Of course, such a time delay is not well-defined and hence cannot be
useful to support our previous results obtained in streaking scenarios. However, the linear changeof the time delay with hxi is caused by the fact that the final energies of the TPI and SPI are notexactly equivalent for finite duration XUV pulses by choosing ωSPI = 2ωTPI . Indeed, the relation
for the expectation value of the final energy of the photoelectron hEf i = N ωXUV − Ip , is strictly
valid only in the limit of an infinite long pulse, for which its energy spectrum is a delta function.
For a typical attosecond XUV pulse this relation is only an approximation, which may sometimes
cause on the attosecond time scale large errors, e.g., for the present time delay study.
To clearly show the impact of this error on the time delay, let us recall the analytical formula

102
for the ionizing wave function of the photoelectron in a Coulomb potential (e.g., [92, 93]), which
for the present 1D Coulomb potential is
Z ln(2kx)
k2+ iϕc (k) + iϕabs (k)Ψ(x, t) ∝ exp ikx − i t + i
2
k

(4.2)

where k is the momentum of the electron, ϕc is the 1D Coulomb scattering phase shift, and ϕabs isthe phase acquired by the electron in the photoabsorption process that is in general not present for
a common scattering process. Using the stationary phase approximation, similar to our trajectory
analysis for SPI, we obtain the trajectory for the photoelectron in a photonionization process ast(x) =xZ ln(2k0 x)
dϕabs
dϕc
+
−
+
,
3
k0
dE E0
dE E0
k0

(4.3)

where k0 is the expectation value of the final asymptotic momentum of the photoelectron andE0 = k02 /2. Here the second term is a logarithmic term that changes slowly with x, which is
however responsible for the divergence property of the WS time delay for the long-range Coulomb
potential. As compared with the corresponding free-particle trajectory, the sum of the second
and third terms gives the time delay the photoelectron accumulates during its propagation in thecontinuum after the photoabsorption, while the last term is the time delay acquired by the electronin the photoabsorption process. We want to emphasize that Eq. (4.3) is general for both SPI and
TPI with the only difference is that the absorption time delay, i.e., the last term, is zero for SPI.
Considering the time delay between the TPI and SPI trajectories, we get
"
# "
#
x
Z ln(2k0TPI x)
dϕc
dϕabs
x
Z ln(2k0SPI x)
dϕc
∆t =
−
+
+
− SPI −+,SPI )3dE E TPI
dE E TPI
dE E SPI
k0TPI(k0TPI )3
k
(k
0
0
0
0
0
(4.4)
where we have dropped the absorption time delay term for SPI since it is zero. By assuming that
k0SPI = k0TPI + ∆k = k0 + ∆k and ∆k is sufficiently small, we can further simplify Eq. (4.4) as
∆t '

∆k
3∆k
dϕabs
x − 4 Z ln(2k0 x) +
,
2
dE E0
k0
k0

(4.5)

which contains a linear term in x that is responsible for the linear trend of the green lines inFig. 4.6(c) and (d). In this simplification step we have also assumed that the third terms in each

103

0.02
(a)

time delay (a.u.)

0
−0.02
−0.04
−0.06
streakingtrajectory original
trajectory adjusted

−0.08
−0.1
70

75

80
85
90
XUV pulse duration (a.u.)

95

20

100

(b)

streakingtrajectory original
trajectory adjusted
time delay (a.u.)

15

10

5

0

50

100
XUV pulse duration (a.u.)

150

Figure 4.7: Time delay obtained from trajectory analysis as a function of XUV pulse duration.
We have used two ways of calculating the time delays from the trajectories: one relies on fitting
the original time delay curve [e.g., green lines in Fig. 4.6(c) and (d)] to Eq. (4.5) (green lines withsquares), the other one relates the time delay to the constant delay line [e.g., red lines in Fig. 4.6(c)and (d)] (red lines with asterisks). Results for both nonresonant [(a)] and resonant [(b)] TPI are
shown. Also shown are the absorption time delays (blue lines with open circles) of TPI extracted
from streaking scenarios.

104
trajectory, i.e., the energy derivative of the scattering phase shift, cancel out. While such anassumption is valid for the 1D Coulomb potential, one needs to take special care for the moregeneral 3D case because different ionization processes will generally result in different final angularmomenta (p wave for SPI, s and d waves for TPI) and thus in different scattering phase shifts eventhough their final momenta are close. In the limit of ∆k → 0, Eq. (4.5) reveals that the time delay
between the TPI and SPI processes is given by
∆t =

dϕabs
.
dE E0

(4.6)

Now, it is apparent that there are two ways to obtain the time delay between TPI and SPIfrom the trajectory calculation: One is to extract the last constant term in Eq. (4.5) by fitting
the time delay curves [green lines in Fig. 4.6(c) and (d)] to Eq. (4.5). The other one is to satisfy
the condition of ∆k = 0 by changing the frequency ωSPI of the SPI XUV pulse. Then the timedelay is directly available from the constant time delay curve [red lines in Fig. 4.6(c) and (d)]. In
Fig. 4.7 we show results for the time delay obtained using these two methods for different XUVpulse durations. As expected, results from these two methods are in very good agreement for bothresonant and nonresonant TPI processes. More importantly, we find that both results also agree
well with the absorption time delays (blue lines with circles in Fig. 4.7) extracted from streaking
scenarios for TPI. This agreement is not surprising since our present time delay and the previous
extracted absorption time delay have the same reference point, namely the SPI process, which
occurs at the center of the XUV pulse that is the reference point of the absorption time delay. The
trajectory analysis presented here is a completely independent study of the time delay even without
the presence of the streaking field and therefore is another strong support of our new findings.

4.2.2

Phase derivative analysis

As widely discussed (e.g., [40, 41]), the time delay in a process can often be related tothe energy derivative of the phase that the particle acquires during this process (also see our
trajectory analysis in the last subsection). We therefore like to make use of this idea to seek

105
another independent way to confirm and understand our new findings. Based on second-ordertime-dependent perturbation theory, the complex amplitude cf of the ionizing wave packet in thecontinuum state |f i ionized from initial state |ii after the XUV pulse ceases can be written as [89],
cf =

X

Z T
µf m µmii∆f m t

e

Z t
E(t)

e

0

m

i∆mi t00

E(t )dt

0


dt,

(4.7)

0

where T is the duration of XUV pulse E(t), µjk and ∆jk = Ej − Ek are the dipole transition matrix
element and the energy difference between states |ji and |ki, respectively. For nonresonant TPI, if
all resonant states are out of the XUV bandwidth, Eq. (4.7) can be further simplified as [94, 95],
Z T
cf ∝ af =ei∆f i t E 2 (t)dt,

(4.8)

0

where we have dropped the dipole matrix element since its phase is related to the scattering (orpropagation) of the electron in the continuum, which is not of our interest. For resonant TPI, we
can further rewrite Eq. (4.7) as,Z T
cf ∝ af =i∆f r t

e
0

Z t
E(t)

i∆ri t0e0E(t )dt

0


dt,

(4.9)

0

by assuming only state |ri is within the field bandwidth. Now it is clear that the absorption timedelay can be related to the phase of the complex amplitude as (e.g., [38, 44]),∆tabs =

Td[arg(af )] − ,dE
2(4.10)

where we need to substract T /2 because our pulse starts from 0 instead of −T /2.
For a sin2 XUV profile considered in our simulations, Eq. (4.10) can be evaluated numerically.
As illustrated in Fig. 4.8, the agreement between absorption time delays extracted from streaking
measurement and those based on the phase analysis is acceptable. In particular, in the case of
nonresonant TPI, these two sets of time delays are exactly identical (all equal zero) for a long XUVpulse [e.g., T > 60 in Fig. 4.8(a)]. Furthermore, the perturbation theory also well predicts thelinear dependence of the time delay on the XUV duration for the resonant case [see Fig. 4.8(b)].
These successes confirm our expectation that the absorption time delay is directly related to theenergy derivative of the phase that the electron acquires during the photon absorption process.

106

0.5
absorption time delay (a.u.)

(a)
0

−0.5

−1
streakingphase derivative
−1.5

40

60
80
XUV pulse duation (a.u.)

100

20
absorption time delay (a.u.)

(b)
15

10

5
streakingphase derivative
0

50

100
XUV pulse duation (a.u.)

150

Figure 4.8: Comparison of absorption time delays from different methods: numerical streaking
experiment (blue dashed lines with circles) and phase analysis (red solid lines). As before, bothnonresonant [(a)] and resonant [(b)] TPI processes are studied.

107And, indeed, this phase information is accessible using the attosecond streaking technique oncethe propagation effect is correctly accounted for. More interestingly, for resonant TPI, the linear
dependence of the time delay on the XUV duration actually opens the possibility to use the durationof the ionizing pulse to control the emission time of the photoelectron, which may significantly
influence the dynamics of the subsequent processes in a chemical reaction.

4.3

Control of two-photon processes2
In the last two sections, we have studied the TPI process in atoms and reported and analyzed

new findings about the time delay in such a process. These findings, in particular the linear dependence of the absorption time delay on the XUV pulse duration, provide a new way of controllingelectron dynamics in physical and chemical interactions. Similar to this idea, at the beginning of
my Ph.D. studies, we have also done some work about coherent control of nonresonant two-photon
excitation (TPE) to bound and dissociative states in molecules [96]. Hence in this section we liketo briefly review the results of this project. The basic idea of this coherent control work is to utilize
a train of laser pulses that is phase modulated in its spectrum using the pulse shaping technique toeither enhance, or weaken, or even eliminate transitions between certain states in a molecule. We
want to emphasize that this work may not directly relate to the previous time delay work, however,it may give us some intuitive ideas and inspirations about the application of the new findings in
TPI. The section is organized as follows: First we like to give a brief introduction to coherent
control by introducing some recent work in this area. We then turn to our work by discussing thetheoretical models used for describing dissociative and bound states of a molecule in the second
subsection. Finally, as the major part of this work, the control of excitations to dissociative andbound states in the model potentials are investigated with the help of numerical simulations.
2

The results of this section are presented in J. Su, S. H. Chen, A. Jaroń-Becker, and A. Becker, Phys. Rev. A
84, 065402 (2011).

108
4.3.1

Short introduction to coherent control

Initial ideas of quantum coherent control were based on quantum path interferences using
phase-controlled laser fields [97] and pump-dump (pump-probe) schemes using sequences of laser
pulses with tunable delays [98, 99]. At first glance the former technique makes use of coherence
properties of light fields in the frequency domain while the latter takes advantage of the temporal
evolution of a process. However, converting the respective analysis from one domain to the otheroften adds a complementary view on a particular control scheme [100]. Recent development offemtosecond pulse shaping techniques extends the variety of schemes in quantum coherent control(for a review, see [101]). Since pulse shaping is usually achieved by manipulating amplitude, phase,and polarization of the frequency spectrum of a laser pulse, applications to coherent control oftransitions in atoms and molecules is often analyzed in the frequency domain based on the multipathway interference concept. For example, destructive and constructive interference among various
pathways can tune few-photon absorption probabilities to zero or maximum by modification of thespectral phase of the pulse (e.g., [94, 102–105]). In particular, the experiment in Ref. [94] hasdemonstrated a way to control the nonresonant TPE probability between two states in a caesium
atom with a femtosecond pulse train that is tailored using the pulse shaping technique.
Theoretically, Chen et al. [106] recently investigated some of these control schemes basedon spectral phase modulations from the complementary time-dependent perspective via numericalsolutions of the TDSE of an atom. The analysis shows that TPE probabilities in an (hydrogen)atom are controlled via destructive or constructive interferences between the amplitudes inducedby consecutive subpulses in a pulse train. The results also provide insights into the control of
atomic (2+1)-photon ionization processes. Here, we like to supplement these previous theoreticalstudies by the time-dependent analysis of the control of different types of two-photon transitions
in molecules. In particular, we consider two different types of nonresonant excitations, namely
bound-to-bound state and bound-to-dissociative state transitions.

109
−10

−5
−10

energy (eV)energy (eV)

−12

−14

−16

−15
−20
−25

−18

−30
(a)

−20
0

5
10
internuclear distance (a.u.)

(b)
15

−35
0

5
10
internuclear distance (a.u.)

15

Figure 4.9: Potential energy curves of the ground and lowest lying excited states of our 2D molecularmodel systems: (a) model with a dissociative excited state [Eq. (4.11)] and (b) model with boundexcited states [Eq. (4.12)]. We are interested in controlling excitations between states of the blue
and red curves.

4.3.2

Theoretical models

To study the control of excitations in molecules, we make use of two models describing theinteraction of a SAE diatomic molecule with an external field. It is well known that the first excitedstate (2pσu ) of the H+
2 molecule is a dissociative state. The following 2D model, which accounts
for the coupled electronic and nuclear dynamics, in a SAE diatomic molecule has been frequently
applied in simulations of H+
2 interacting with a laser field. The field-free Hamiltonian is given by
H1 (R, z) =

p2R
1
1
p2
1
−p
+√
+ z −p,2
22M2m
R2 + b(z + R/2) + a1
(z − R/2) + a2

(4.11)

where M and m are the reduced masses of the nuclei and the electron, respectively. R is theinternuclear distance, z is the electron position with respect to the center of mass of the nuclei,
and a1 , a2 as well as b are soft-core Coulomb parameters. In the case of H+
2 one chooses a1 = a2 .
However, the symmetry of the corresponding Hamiltonian prohibits a two-photon transition fromthe electronic ground state to the first dissociative state.
For our goal to study the coherent control of TPE processes, one could instead considertransitions to the second dissociative state. But, in view of the steep decline of the potential energy
curve of the first excited state as a function of the internuclear distance, (dominant) one-photon

110transitions to the first excited state cannot be neglected for an ultrashort laser pulse having abroad bandwidth. As an alternative, we choose to study a molecular model system with a1 6= a2 tobreak the symmetry of the Hamiltonian. In our studies, we arbitrarily choose a1 = 1.0, a2 = 2.0,
and b = 0.03. The adiabatic potential energy curves for the ground and first excited state of thissystem, calculated using ITP of the corresponding TDSE, are shown in Fig. 4.9(a). The equilibrium
distance of the two protons is R0 = 3.36, and the energy gap between the states at R0 is ∆E = 6.88
eV. The curve of the first excited state has a shallow well, which however does not influence thedissociative character of the state for wave packets pumped from the initial ground state close to
the equilibrium distance.
To study transitions between bound states, we use a molecular model system, in which the
Coulomb repulsion 1/R is replaced by a Morse potential in the Hamiltonian, i.e.,
H2 (R, z) =

h
i2
p2R
1
p2
1
+ z −p−p+ D 1 − e−λ(R−R0 ) − D, (4.12)
2M 2m(z + R/2)2 + a1
(z − R/2)2 + a2

where D is the depth of the well, R0 is the equilibrium distance, λ =pk/2D with k is the bondforce constant, and all the other parameters are the same defined as before. Using a1 = a2 = 4.0,
D = 0.4, λ = 0.5, and R0 = 2.0, the potential energy curves of the ground and first two excitedstates, shown in Fig. 4.9(b), are well separated. We consider nonresonant two-photon transitions tothe second excited state, which is possible in the present model since the energy difference between
the ground and the first excited states clearly exceeds the corresponding photon energy required
for the two-photon transition.
To investigate the time evolution of the TPE process, we consider the corresponding TDSE of
the molecular model systems interacting with a laser pulse linearly polarized along the internuclear
axis (j = 1, 2),i∂
Ψ(R, z; t) = [Hj (R, z) + E(t)z] Ψ(R, z; t),∂t(4.13)

which is solved with ∆R = 0.03, ∆z = 0.2, and ∆t = 0.01 using our common grid method. The
population in the ground and excited bound states are obtained by projection on the respective fieldfree states, which for the laser parameters used in the present study is a reasonable approximation

111
even during the interaction with the external field. In the case of the excited dissociative state thenPpopulation is calculated as Pexc = 1−
Pν , where ν denotes the vibrational mode in the electronic
ν=0ground state and n is the total number of the vibrational states considered. The ionization anddissociation probabilities are calculated as the outgoing probability flux at the respective boundariesof the grid. At the edges we use cos1/6 mask functions to suppress reflections. We have checkedthat the results are not influenced by remaining small inaccuracies due to the boundary conditions.
We use a grid of NR × Nz = 600 × 200 points (NR × Nz = 400 × 600 points) for the dissociative(bound) state model.

4.3.3

Control of excitation to dissociative state

First, we investigate control of TPE from a bound to a dissociative molecular state. To thisend, we consider the dissociative state model and a spectral phase modulated laser field of the form
E

ω
2



+ Ω = E0 sech



1.76Ω
∆ω


exp [iα cos(βΩ + φ)] .

(4.14)Such fields have been used recently by Silberberg and coworkers [94] to study the control of
two-photon transitions in atoms. For our numerical simulations we choose the central frequency
ω/2 = 3.44 eV (half of the energy gap of the desired transition), the bandwidth ∆ω = 0.2 eV,
the modulation depth α = 1.2024, and the modulation frequency β = 35 fs, while the modulation
phase φ is varied. A Fourier transform (FT) of these fields yields a pulse train in the temporal
domain [e.g., see Fig. 4.10(a)]. In our studies we fix the (overall) peak intensity in the pulse train
to I0 = 1 × 1012 W/cm2 , which for the present frequencies is within the perturbation regime.
In the atomic case it was shown [95] that using φ = 0 a dark pulse is formed and the final
population in the excited state at the end of the pulse train vanishes. On the other hand theexcitation probability was found to be maximized for pulse with φ = π/2, which is therefore calleda bright pulse. In contrast, we find that for the two-photon transition to a dissociative state in thepresent molecular model system that the final excitation probability is independent of the specificvalue of φ, as exemplified by the results for φ = 0 (left column) and π/2 (right column) in Fig. 4.10.

electric field (a.u.)

112
0.01

0.01

0

0

(a)

−0.01
−2

−2

probability

10

10

−6

−6

10

10
(c)

−10

10

(b)

−0.01

0

40

80

120

160

(d)

−10

10

0
time (fs)

40

80

120

160

electric field (a.u.)

Figure 4.10: Two-photon coherent control of excitation to a dissociative molecular state. The upper
panels show that the electric field distribution as a function of time [i.e., the Fourier transform ofEq. (4.14)] for (a) φ = 0 and (b) φ = π/2. The insets show the central field cycle and the carrierenvelope phase (CEP) of each subpulse. In the lower panels the time evolution of the probabilities
in the first excited state (blue solid line) is presented for (c) φ = 0 and (d) φ = π/2. We also show
the dissociation probabilities (red dashed lines).

0.01

0.01

0

0

(a)

−0.01
−2

−2

probability

10

10

−6

−6

10

10
(c)

−10

10

(b)

−0.01

0

40

80

120

160

(d)

−10

10

0
time (fs)

40

80

120

160

Figure 4.11: Same as Fig. 4.10, but for a molecular model system with fixed internuclear distanceR = R0 .

113
The excitation probability increases stepwise as a function of time with the interaction of every
subpulse in the train. Moreover, the final excitation probability is equal to the probability fordissociation at the end of the molecular model system (red dashed lines). This indicates it is thenuclear dynamics that causes the breakdown of the coherent control scheme.
In order to confirm this assumption, we perform another series of calculations in which we
fixed the internuclear distance R = R0 in Eq. (4.11) to suppress any nuclear dynamics but keep allother parameters unchanged except the central frequency ω/2 = 3.48 eV. The central frequency isslightly increased since the ground state energy of the fixed nuclei model is lower than that of the
2D model. The results of the fixed nuclei model for φ = 0 (Fig. 4.11, left column) and φ = π/2(Fig. 4.11, right column) agree with the findings for the atomic case. For φ = 0 the final excited
population is (close to) zero (dark pulse), while for φ = π/2 the population is maximized (brightpulse). Thus, from the temporal analysis of the process we conclude that coherent control of the
TPE to a dissociative molecular state using a train of pulses [or, a spectral phase modulated pulse of
the form given in Eq. (4.14)] fails since each of the wave packets pumped to the excited state by the
subpulses in the train quickly propagates to larger internuclear distances. Therefore, concerning thecontrol of the total population in the excited (dissociative) state there is no signature of destructive
nor constructive interferences between the wave packets generated from subsequent pulses.

4.3.4

Control of excitation to bound states

Next, we investigate the influence of nuclear dynamics on the control of two-photon excitations to bound molecular states. As shown above, this dynamics becomes effective for the present
control scheme over the time delay between two subsequent pulses in a pulse train, since thefundamental control mechanism is based on the interference between two electronic wave packets
induced by consecutive subpulses [106]. It is well known from studies in wave packet interferometry [107, 108] that, in the case of a two-pulse scenario, a control of the excitation probability canbe achieved via the pulse delay τ and the relative phase φrel of the two pulses. For example, the114effective coherent control of the excitation of a superposition of two vibrational states requires
τ =N

2π= N Trevival
E21 − E20

(4.15)and varying φrel . Here, E2i is the energy of the ith vibrational level of the second excited electronic
state of our model. The pulse delays have to coincide with the multiples of the revival period Trevivalof the vibrational wave packet created in the excited electronic state. This enables the control via
efficient interference of wave packets generated by subsequent pulses, since the contributions inboth vibrational levels of the second excited electronic state are in phase at these times.
The above analysis for the two-pulse sequence can be readily applied to the spectral phasemodulated fields used by Silberberg and coworkers by noting that the modulation frequency β
in Eq. (4.14) corresponds to the time delay τ between the consecutive subpulses in the train in
the temporal domain. Furthermore, the relative phase φrel between the consecutive subpulses isdetermined by the modulation phase φ in Eq. (4.14). Thus, setting β = N Trevival we expect toachieve an efficient control pattern by changing φ. To test our expectations, we do simulations bychoosing α = 1.2024, β = 8Trevival = 85.6682 fs, ω/2 = 5.83 eV, ∆ω = 0.1 eV, and I0 = 1 × 1011
W/cm2 . Note that we have chosen such a narrow bandwidth that only the lowest two vibrational
levels of the second electronic state can be excited. The results of the simulations (circles andasterisks) for excitation to these two vibrational states as a function of the modulation phase φ
are shown in Fig. 4.12. The same dependence of the two excitation probabilities on φ clearly
confirms our expectations. We also show results calculated from the second-order perturbation
theory (red dashed lines) in Fig. 4.12. Although the coefficient of the transition amplitude [94, 95]
is not taken into account in our perturbation calculations, by normalizing the values from theperturbation theory to the maximum probability from the TDSE, we can see a good agreement of
results from these two theories. The temporal analysis of the populations in the vibrational states(Fig. 4.13) clearly exhibits the destructive and constructive interference effects for the subsequent
pulses in the train for the spectral phase modulation of a dark (left column) and a bright pulse(right column).

115
−4

x 10

ν=0 TDSE
ν=0 pert.
ν=1 TDSE
ν=1 pert.

probability

2

1

0
0

0.2

0.4
0.6
φ (in units of π)

0.8

1Figure 4.12: Excitation probabilities to the ground vibrational (ν2nd = 0, circles and dashed line)
and the first excited vibrational states (ν2nd = 1, asterisks and dashed line) as functions of φ.
A spectral phase modulated field defined in Eq. (4.14) is used with parameters: α = 1.2024,
β = 8Trevival = 85.6682 fs, ω/2 = 5.83 eV, ∆ω = 0.1 eV, and I0 = 1 × 1011 W/cm2 . The symbols(circles and asterisks) are numerical results obtained by solving the TDSE, while the curves (dashedlines) are obtained using second-order perturbation theory.

−3

electric field (a.u.)

4

−3

x 10

4

0

0

(a)

−4
−3

(b)

−4
−3

10
probability

x 10

10

−6

−6

10

10
(c)

−9

(d)
−9

10

10
0

100

200

300

400

0
time (fs)

100

200

300

400

Figure 4.13: Two-photon coherent control of excitations to a superposition of two vibrational states.
The upper panels show that the electric field distribution as a function of time [i.e., the Fouriertransform of Eq. (4.14)] for (a) φ = 0.96π (dark pulse) and (b) φ = 0.46π (bright pulse). The insets
show the central field cycle and the CEP of each subpulse. In the lower panels the time evolutionof the probabilities in the ground (blue solid line) and first excited (red dashed line) vibrational
state is presented for the (c) dark and (d) bright pulses, respectively. Also shown is the ionization
probability (black dashed-dotted line).

Chapter 5

Conclusions

In this thesis we have systematically studied the attosecond time delays in single- and twophoton ionization using theoretical analyses and numerical simulations. For single-photon ionization(SPI) we have developed a new numerical method to calculate the Wigner-Smith-like time delay
from a time-dependent perspective and proposed an alternative interpretation for the time delay
measured in the streaking experiment. This new interpretation provides a deeper understanding of
the attosecond streaking technique and thus allows us to use it to study other ultrafast processes
such as the two-photon ionization (TPI) process.
In Chapter 2, we reviewed the numerical model and method used to describe laser-matterinteraction. We first reviewed the single-active-electron model, in which only one electron in a multielectron system is treated as active to reduce costs in numerical simulations. We then presented
our numerical approaches to solving the TDSE, including the methods of obtaining initial fieldfree states, propagating the wave function in time domain, and working with different absorbing
boundary conditions.
We carefully studied time delays in SPI in Chapter 3. This chapter consists of two major
aspects: On the one hand, we calculated the WS time delay for SPI from its fundamental definition using a new time-dependent method. Such a method provides a time-dependent analysisof the photoionization process and is more general for different potentials than the previous timeindependent definition. On the other hand, we investigated the time delay measured in SPI using
the attosecond streaking technique and proposed a new interpretation of such a measurement. Our

117
new interpretation suggests that the streaking time delay arises from the measurement, i.e., thepropagation of the electron in the coupled field of the streaking pulse and the ionic potential after
the absorption of the photon, and thus can be calculated using classical mechanics. The streaking
time delay can be further understood as the sum of piecewise field-free time delays weighted by the
relative instantaneous streaking field and only relates to a finite-range part of the ionic potential.
The latter finite-range property allows us to use the attosecond streaking technique as an imaging
tool.
We further extended our studies to present results on time resolving the TPI process in
Chapter 4. Based on our understanding of time delays in SPI, in particular, the new interpretation
for the streaking time delay, we are able to make use of the streaking technique to retrieve the timedelay of two-photon absorption with respect to the center of the XUV pulse. By accounting forthe propagation time delay obtained from classical simulations, we could retrieve the absorption
time delay from the streaking time delay extracted from the streaking trace. Our results revealthat the streaking time delay is zero for nonresonant TPI, but it is nonzero and changes linearlywith the XUV pulse duration for resonant TPI. We further confirmed our findings using two otherindependent methods: the trajectory analysis based on quantum simulations and the phase analysis
using second-order perturbation theory. Since the TPI time delay results may be further applied
to control the later dynamics followed by the photoionization process, we took the opportunity
here to review our previous work on coherent control of two-photon excitations in molecular model
systems.

Bibliography[1] F. Krausz and M. Ivanov. Review of Modern Physics, 81:163–234, 2009.
[2] A. L. Cavalieri, N. Muller, T. Uphues, V. S. Yakovlev, A. Baltuska, B. Horvath, B. Schmidt,
L. Blumel, R. Holzwarth, S. Hendel, M. Drescher, U. Kleineberg, P. M. Echenique, R. Kienberger, F. Krausz, and U. Heinzmann. Nature, 449:1029–1032, 2007.
[3] M. Schultze, M. Fiess, N. Karpowicz, J. Gagnon, M. Korbman, M. Hofstetter, S. Neppl, A. L.
Cavalieri, Y. Komninos, T. Mercouris, C. A. Nicolaides, R. Pazourek, S. Nagele, J. Feist,
J. Burgdörfer, A. M. Azzeer, R. Ernstorfer, R. Kienberger, U. Kleineberg, E. Goulielmakis,
F. Krausz, and V. S. Yakovlev. Science, 328:1658–1662, 2010.
[4] K. Klünder, J. M. Dahlström, M. Gisselbrecht, T. Fordell, M. Swoboda, D. Guenot, P. Johnsson, J. Caillat, J. Mauritsson, A. Maquet, R. Taı̈eb, and A. L’Huillier. Physical Review
Letters, 106:143002, 2011.
[5] M. Uiberacker, Th. Uphues, M. Schultze, A. J. Verhoef, V. Yakovlev, M. F. Kling,
J. Rauschenberger, N. M. Kabachnik, H. Schröder, M. Lezius, K. L. Kompa, H.-G. Muller,
M. J. J. Vrakking, S. Hendel, U. Kleineberg, U. Heinzmann, M. Drescher, and F. Krausz.
Nature, 446:627–632, 2007.
[6] P. Eckle, A. N. Pfeiffer, C. Cirelli, A. Staudte, R. Dörner, H. G. Muller, M. Büttiker, and
U. Keller. Science, 322:1525–1529, 2008.
[7] F. Remacle and R. D. Levine. Proceedings of the National Academy of Sciences of the UnitedStates of America, 103:6793–6798, 2006.
[8] P. B. Corkum and F. Krausz. Nature Physics, 3:381–387, 2007.
[9] L. Gallmann, C. Cirelli, and U. Keller. Annual Review of Physical Chemistry, 63:447–469,
2012.
[10] M. Ferray, A. L’Huillier, X. F. Li, L. A. Lompre, G. Mainfray, and C. Manus. Journal of
Physics B, 21:L31–L35, 1988.
[11] A. McPherson, G. Gibson, H. Jara, U. Johann, T. S. Luk, I. A. McIntyre, K. Boyer, andC. K. Rhodes. Journal of the Optical Society of America B, 4:595–601, 1987.
[12] J. L. Krause, K. J. Schafer, and K. C. Kulander. Physical Review A, 45:4998–5010, 1992.

119[13] J. L. Krause, K. J. Schafer, and K. C. Kulander. Physical Review Letters, 68:3535–3538,
1992.
[14] P. B. Corkum. Physical Review Letters, 71:1994–1997, 1993.
[15] M. Lewenstein, Ph. Balcou, M. Y. Ivanov, A. L’Huillier, and P. B. Corkum. Physical Review
A, 49:2117–2132, 1994.
[16] I. P. Christov, J. Zhou, J. Peatross, A. Rundquist, M. M. Murnane, and H. C. Kapteyn.
Physical Review Letters, 77:1743–1746, 1996.
[17] J. Zhou, J. Peatross, M. M. Murnane, H. C. Kapteyn, and I. P. Christov. Physical Review
Letters, 76:752–755, 1996.
[18] I. P. Christov, M. M. Murnane, and H. C. Kapteyn. Physical Review Letters, 78:1251–1254,
1997.
[19] T. Popmintchev, M.-C. Chen, P. Arpin, M. M. Murnane, and H. C. Kapteyn. Nature
Photonics, 4:822–832, 2010.
[20] T. Popmintchev, M.-C. Chen, D. Popmintchev, P. Arpin, S. Brown, S. Ališauskas, G. Andriukaitis, T. Balčiunas, O. D. Mcke, A. Pugzlys, A. Baltuška, B. Shim, S. E. Schrauth,A. Gaeta, C. Hernández-Garcı́a, L. Plaja, A. Becker, A. Jaroń-Becker, M. M. Murnane, and
H. C. Kapteyn. Science, 336:1287–1291, 2012.
[21] V. S. Yakovlev and A. Scrinzi. Physical Review Letters, 91:153901, 2003.
[22] M. Chini, K. Zhao, and Z. Chang. Nature Photonics, 8:178–186, 2014.
[23] M. Hentschel, R. Kienberger, Ch. Spielmann, G. A. Reider, N. Milosevic, T. Brabec,
P. Corkum, U. Heinzmann, M. Drescher, and F. Krausz. Nature, 414:509–513, 2001.
[24] E. Goulielmakis, M. Schultze, M. Hofstetter, V. S. Yakovlev, J. Gagnon, M. Uiberacker,
A. L. Aquila, E. M. Gullikson, D. T. Attwood, R. Kienberger, F. Krausz, and U. Kleineberg.
Science, 320:1614–1617, 2008.
[25] M.-C. Chen, C. Hernández-Garcı́a, C. Mancuso, F. Dollar, B. Galloway, D. Popmintchev, P.C. Huang, B. Walker, L. Plaja, A. Jaroń-Becker, A. Becker, T. Popmintchev, M. M. Murnane,and H. C. Kapteyn. Arxiv, http://arxiv.org/abs/1401.0240v1, 2014.
[26] J. A. Armstrong, N. Bloembergen, J. Ducuing, and P. S. Pershan. Physical Review, 127:1918–
1939, 1962.
[27] A. Rundquist, C. G. Durfee, Z. Chang, C. Herne, S. Backus, M. M. Murnane, and H. C.
Kapteyn. Science, 280:1412–1415, 1998.
[28] C. G. Durfee, A. R. Rundquist, S. Backus, C. Herne, M. M. Murnane, and H. C. Kapteyn.
Physical Review Letters, 83:2187–2190, 1999.
[29] E. Constant, D. Garzella, P. Breger, E. Mével, Ch. Dorrer, C. Le Blanc, F. Salin, and P. Agostini. Physical Review Letters, 82:1668–1671, 1999.

120[30] T. Popmintchev, M.-C. Chen, A. Bahabad, M. Gerrity, P. Sidorenko, O. Cohen, I. P. Christov,M. M. Murnane, and H. C. Kapteyn. Proceedings of the National Academy of Sciences,106:10516–10521, 2009.
[31] M. Uiberacker and F. Krausz. SPIE’s oemagazine, May:18–22, 2004.
[32] M. Ya. Schelev, M. C. Richardson, and A. J. Alcock. Applied Physics Letters, 18:354–357,
1971.
[33] D. J. Bradley, B. Liddy, and W. E. Sleat. Optics Communications, 2:391–395, 1971.
[34] J. Itatani, F. Quéré, G. L. Yudin, M. Y. Ivanov, F. Krausz, and P. B. Corkum. PhysicalReview Letter, 88:173903, 2002.
[35] M. Kitzler, N. Milosevic, A. Scrinzi, F. Krausz, and T. Brabec. Physical Review Letters,
88:173904, 2002.
[36] M. Drescher, M. Hentschel, R. Kienberger, G. Tempea, Ch. Spielmann, G. A. Reider, P. B.
Corkum, and F. Krausz. Science, 291:1923–1927, 2001.
[37] S. Nagele, R. Pazourek, J. Feist, K. Doblhoff-Dier, C. Lemell, K. Tökési, and J. Burgdörfer.
Journal of Physics B, 44:081001, 2011.
[38] M. Ivanov and O. Smirnova. Physical Review Letters, 107:213605, 2011.
[39] C. H. Zhang and U. Thumm. Physical Review A, 82:043405, 2010.
[40] E. P. Wigner. Physical Review, 98:145–147, 1955.
[41] F. T. Smith. Physical Review, 118:349–356, 1960.
[42] J. Su, H. Ni, A. Becker, and A. Jaroń-Becker. Physical Review A, 87:033420, 2013.
[43] J. Su, H. Ni, A. Becker, and A. Jaroń-Becker. Journal of Morden Optics, 60:1484–1491, 2013.
[44] A. S. Kheifets and I. A. Ivanov. Physical Review Letter, 105:233002, 2010.
[45] A. S. Kheifets, I. A. Ivanov, and I. Bray. Journal of Physics B, 44:101003, 2011.
[46] I. A. Ivanov. Physical Review A, 83:023421, 2011.
[47] I. A. Ivanov. Physical Review A, 86:023419, 2012.
[48] I. A. Ivanov, A. S. Kheifets, and V. V. Serov. Physical Review A, 86:063422, 2012.
[49] A. S. Kheifets. Physical Review A, 87:063404, 2013.
[50] V. V. Serov, V. L. Derbov, and T. A. Sergeeva. Physical Review A, 87:063414, 2013.
[51] G. Dixit, H. S. Chakraborty, and M. E. Madjet. Physical Review Letters, 111:203003, 2013.
[52] C. H. Zhang and U. Thumm. Physical Review A, 84:033401, 2011.
[53] S. Nagele, R. Pazourek, J. Feist, and J. Burgdörfer. Physical Review A, 85:033401, 2012.

121[54] J. C. Baggesen and L. B. Madsen. Physical Review Letters, 104:043602, 2010.
[55] L. R. Moore, M. A. Lysaght, J. S. Parker, H. W. van der Hart, and K. T. Taylor. PhysicalReview A, 84:061404, 2011.
[56] R. Pazourek, J. Feist, S. Nagele, and J. Burgdörfer. Physical Review Letters, 108:163001,
2012.
[57] D. Guénot, K. Klünder, C. L. Arnold, D. Kroon, J. M. Dahlström, M. Miranda, T. Fordell,M. Gisselbrecht, P. Johnsson, J. Mauritsson, E. Lindroth, A. Maquet, R. Taı̈eb, A. L’Huillier,and A. S. Kheifets. Physical Review A, 85:053424, 2012.
[58] J. M. Dahlström, A. L’Huillier, and A. Maquet. Journal of Physics B, 45:183001, 2012.
[59] J. M. Dahlström, D. Guenot, K. Klünder, M. Gisselbrecht, J. Mauritsson, A. L’Huillier,A. Maquet, and R. Taı̈eb. Chemical Physics, 414:53–64, 2013.
[60] M. Cardona and L. Ley. Photoemission in Solids I: General Principles. Springer-Verlag,
Berlin, 1978.
[61] E. Clementi and C. Roetti. Atomic Data and Nuclear Data Tables, 14:177–478, 1974.
[62] J. J. Yeh and I. Lindau. Atomic Data and Nuclear Data Tables, 32:1–155, 1985.
[63] I. N. Levine. Quantum Chemsitry. Prentice Hall, 2000.
[64] J. P. Perdew and A. Zunger. Physical Review B, 23:5048–5079, 1981.
[65] J. B. Krieger, Y. Li, and G. J. Iafrate. Physical Review A, 45:101–126, 1992.
[66] Y. Li, J. B. Krieger, and G. J. Iafrate. Physical Review A, 47:165–181, 1993.
[67] J. Chen, J. B. Krieger, Y. Li, and G. J. Iafrate. Physical Review A, 54:3939–3947, 1996.
[68] R. T. Sharp and G. K. Horton. Physical Review, 90:317–340, 1953.
[69] J. D. Talman and W. F. Shadwick. Physical Review A, 14:36–40, 1976.
[70] X. M. Tong and S. I. Chu. Physical Review A, 55:3406–3416, 1997.
[71] X. M. Tong and C. D. Lin. Journal of Physics B, 38:2593–2600, 2005.
[72] J. Crank and P. Nicolson. Advances in Computational Mathematics, 6:207–226, 1996.
[73] O. A. Sharafeddin, D. J. Kouri, and D. K. Hoffman. Canadian Journal of Chemistry, 70:686–
692, 1992.
[74] B. C. Hall. Lie Groups, Lie Algebras, and Representations. Springer New York, 2003.
[75] M. D. Feit, J. A. Fleck Jr., and A. Steiger. Journal of Computational Physics, 47:412–433,
1982.
[76] F. He, C. Ruiz, and A. Becker. Physical Review A, 75:053407, 2007.
[77] C. A. A. de Carvalhoa and H. M. Nussenzveig. Physics Reports, 364:83–174, 2002.

122[78] R. G. Newton. Scattering Theory of Waves and Particles. Springer-Verlag, New York, 1982.
[79] J. Itatani, J. Levesque, D. Zeidler, Hiromichi Niikura, H. Ppin, J. C. Kieffer, P. B. Corkum,and D. M. Villeneuve. Nature, 432:867–871, 2004.
[80] S. N. Yurchenko, S. Patchkovskii, I. V. Litvinyuk, P. B. Corkum, and G. L. Yudin. PhysicalReview Letters, 93:223003, 2004.
[81] J. Su, H. Ni, A. Becker, and A. Jaroń-Becker. Physical Review A, 88:023413, 2013.
[82] P. B. Corkum. Physical Review Letters, 71:1994–1998, 1993.
[83] O. Smirnova, A. S. Mouritzen, S. Patchkovskii, and M. Ivanov. Journal of Physics B, 40:F197–
F206, 2007.
[84] O. Smirnova, M. Spanner, and M. Ivanov. Journal of Physics B, 39:S307–S321, 2006.
[85] J. Su, H. Ni, A. Becker, and A. Jaroń-Becker. Physical Review A, 89:013404, 2014.
[86] J. Su, H. Ni, A. Becker, and A. Jaroń-Becker. Chinese Journal of Physics, 52:404–415, 2014.
[87] Y. Nabekawa, T. Shimizu, T. Okino, K. Furusawa, H. Hasegawa, K. Yamanouchi, and K. Midorikawa. Physical Review Letters, 97:153904, 2006.
[88] L. H. Haber, B. Doughty, and S. R. Leone. Physical Review A, 79:031401, 2009.
[89] K. L. Ishikawa and K. Ueda. Physical Review Letters, 108:033003, 2012.
[90] R. Moshammer, Th. Pfeifer, A. Rudenko, Y. H. Jiang, L. Foucar, M. Kurka, K. U. Kühnel,
C. D. Schröter, J. Ullrich, O. Herrwerth, M. F. Kling, X.-J. Liu, K. Motomura, H. Fukuzawa,
A. Yamada, K. Ueda, K. L. Ishikawa, K. Nagaya, H. Iwayama, A. Sugishima, Y. Mizoguchi,
S. Yase, M. Yao, N. Saito, A. Belkacem, M. Nagasono, A. Higashiya, M. Yabashi, T. Ishikawa,
H. Ohashi, H. Kimura, and T. Togashi. Optics Express, 19:21698–21706, 2011.
[91] J. Su, H. Ni, A. Jaroń-Becker, and A. Becker. Time Delays in Two-Photon Ionization.
submitted.
[92] R. Pazourek, S. Nagele, and J. Burgdörfer. Faraday Discussions, 163:353–376, 2013.
[93] C. W. Clark. American Journal of Physics, 47:683–684, 1979.
[94] D. Meshulach and Y. Silberberg. Nature, 396:239–242, 1998.
[95] D. Meshulach and Y. Silberberg. Physical Review A, 60:1287–1292, 1999.
[96] J. Su, S. H. Chen, A. Jaroń-Becker, and A. Becker. Physical Review A, 84:065402, 2011.
[97] P. Brumer and M. Shapiro. Chemical Physics Letters, 126:541–550, 1986.
[98] D. J. Tannor and S. A. Rice. Journal of Chemical Physics, 83:5013–5018, 1985.
[99] D. J. Tannor, R. Kosloff, and S. A. Rice. Journal of Chemical Physics, 85:5805–5820, 1986.
[100] S. A. Rice and M. Zhao. Optical Control of Molecular Dynamics. Wiley-Interscience, New
York, 2000.

123[101] A. M. Weiner. Optics Communication, 284:3669–3692, 2011.
[102] N. Dudovich, B. Dayan, S. M. Gallagher Faeder, and Y. Silberberg. Physical Review Letters,
86:47–50, 2001.
[103] P. Panek and A. Becker. Physical Review A, 74:023408, 2006.
[104] A. Gandman, L. Chuntonov, L. Rybak, and Z. Amitay. Physical Review A, 75:031401, 2007.
[105] T. Bayer, M. Wollenhaupt, C. Sarpe-Tudoran, and T. Baumert. Physical Review Letters,
102:023004, 2009.
[106] S. H. Chen, A. Jaroń-Becker, and A. Becker. Physical Review A, 82:013414, 2010.
[107] M. M. Salour and C. Cohen-Tannoudji. Physical Review Letters, 38:757–760, 1977.
[108] R. Teets, J. Eckstein, and T. W. Hansch. Physical Review Letters, 38:760–764, 1977.