# Polarization in Strong-Field Ionization of Excited Helium

A. C. Bray  Department of Physics & Astronomy, University College London,
Gower St., London WC1E 6BT, UK  A. S. Maxwell  Department of Physics &
Astronomy, University College London, Gower St., London WC1E 6BT, UK  ICFO –
Institut de Ciencies Fotoniques, The Barcelona Institute of Science and
Technology, 08860 Castelldefels (Barcelona), Spain  Y. Kissin  Imperial
College London - Department of Physics, Prince Consort Rd, Kensington, London
SW7 2BW, UK  M. Ruberti  Imperial College London - Department of Physics,
Prince Consort Rd, Kensington, London SW7 2BW, UK  M. F. Ciappina  Physics
Program, Guangdong Technion – Israel Institute of Technology, Shantou,
Guangdong 515063, China  Technion – Israel Institute of Technology, Haifa,
32000, Israel  V. Averbukh  Imperial College London - Department of Physics,
Prince Consort Rd, Kensington, London SW7 2BW, UK  C. Figueira De Morisson
Faria  Department of Physics & Astronomy, University College London, Gower
St., London WC1E 6BT, UK

###### Abstract

We analyze how bound-state excitation, electron exchange and the residual
binding potential influence above-threshold ionization (ATI) in helium
prepared in an excited $p$ state, oriented parallel and perpendicular to a
linearly polarized mid-IR field. Using the ab initio B-spline algebraic
diagrammatic construction (ADC), and several one-electron methods with
effective potentials, including the Schrödinger solver Qprop, modified
versions of the strong-field approximation and the Coulomb quantum-orbit
strong-field approximation (CQSFA), we find that these specific physical
mechanisms leave significant imprints in ATI spectra and photoelectron
momentum distributions. Examples are changes of up to two orders of magnitude
in the high-energy photoelectron region, and ramp-like structures that can be
traced back to Coulomb-distorted trajectories. The present work also shows
that electron exchange renders rescattering less effective, causing
suppressions in the ATI plateau. Due to the long-range potential, the electron
continuum dynamics are no longer confined to the polarization axis, in
contrast to the predictions of traditional approaches. Thus, one may in
principle probe excited-state configurations perpendicular to the driving-
field polarization without the need for orthogonally polarized fields.

##  I Introduction

Above-threshold ionization (ATI) is a strong-field phenomenon in which an atom
absorbs more photons than are energetically required for it to ionize. Since
its first observation Agostini et al. (1979), ATI has been employed to provide
unprecedented insight into the interaction of intense laser radiation with
atoms or molecules. The high-energy photoelectrons released into the
continuum, together with subfemtosecond resolution, make ATI a powerful tool
for probing real-time electron dynamics. Applications of ATI directly
associated with imaging and target reconstruction are, for instance, laser-
induced electron diffraction (LIED) Zuo et al. (1996) and ultrafast
photoelectron holography Huismans et al. (2011). Such imaging methods combine
high photoelectron currents and the possibility of retrieving phase
differences without the need for elaborate interferometric schemes. See, e.g.,
de Morisson Faria and Maxwell (2020) for a review.

For that reason, orbit-based models have been used for over two decades due to
their huge predictive power. By relating quantum transition amplitudes to
electron paths, which may interfere for a given photoelectron energy, they
provide a physical picture of ATI as a laser-induced rescattering process
(Lewenstein et al. (1995); Lohr et al. (1997); Milošević and Ehlotzky (1998);
for a review see Becker et al. (2018)). Attosecond resolution is guaranteed in
the rescattering process and ensuing phenomena because it takes place within a
fraction of a field cycle, which, typically, is a few femtoseconds.
Established approaches such as the strong-field approximation (SFA) Keldysh
(1965); Faisal (1973); Reiss (1980) consider the target as a source term and
employ a structureless continuum, i.e., field-dressed plane waves. This has
allowed the transition amplitudes associated to strong-field phenomena to be
written as a laser-dressed Born-type series Lewenstein et al. (1995); Lohr et
al. (1997); Milošević and Ehlotzky (1998), which considerably reduces the
numerical effort involved. Besides the easy implementation, an obvious
advantage is that the interaction with the core is well-defined. This makes
rescattering processes straightforward to identify. Furthermore, the SFA can
be easily associated with interfering electron orbits, if formulated in
conjunction with saddle-point methods Becker et al. (2002); Kopold et al.
(2000); Salières et al. (2001).

For that reason, the insight provided by the SFA was vital for identifying the
overall structure and interpreting key features in ATI spectra. The ATI
spectrum consists of a low-energy region, up to $2U_{p}$, where $U_{p}$ is the
ponderomotive energy, and of a high energy region extending to approximately
$10U_{p}$. The $2U_{p}$ cutoff corresponds to the ‘classical’ kinetic energy
limit for ‘direct’ ATI electrons, which are freed into the continuum via
tunnel or multiphoton ionization, and reach the detector without further
interaction with the target. The ionization mechanism that predominates
depends on the Keldysh parameter for the system in question,
$\gamma=\sqrt{I_{p}/2U_{p}}$, where $I_{p}$ is the ionization potential. A
Keldysh parameter of greater than one implies multiphoton ionization, while
tunnelling is dominant if $\gamma\leq 1$ Keldysh (1965).

The high-energy region is mainly dominated by high-order ATI (HATI), which is
characterized by a plateau that does not significantly change across the
energy distribution, falling off at the cutoff point of approximately
$10U_{p}$. The widely accepted explanation of this feature is that the high-
energy photoelectrons result from the ionized electron being driven back by
the laser field to scatter elastically off its parent ion. The rescattered
electron is then accelerated further in the field to energies up to $10U_{p}$,
which can be derived from a classical model Paulus et al. (1994).
Alternatively, the freed electron may recombine with a bound state of the
parent ion, which leads to high-order harmonic generation (HHG), or recollide
inelastically with the core, releasing other electrons. The latter processes
give rise to laser-induced nonsequential double and multiple ionization (NSDI,
NSMI) de Morisson Faria and Liu (2011); Becker et al. (2012).

Nonetheless, the aforementioned picture is an oversimplification. For example,
ATI peaks have a substructure \- termed Freeman resonances - caused by the
ponderomotive shifts of states that produces resonant enhancements Freeman et
al. (1987). This effect can also cause broadening of the peaks and even large
energy shifts. Furthermore, recent discoveries have cast doubt on the
distinction between tunneling and multiphoton regime based on the Keldysh
parameter - for instance, photoelectron spectra obtained from solving the
time-dependent Schrödinger equation (TDSE) show features consistent with
multiphoton absorption even for laser intensities that correspond to the
tunnelling ionization regime Morishita et al. (2007). In addition to that,
specific features of the plateau can vary from system to system. Some atomic
targets, such as krypton, have a dropping plateau while others, including
xenon, have a flatter plateau. Specific laser intensities can also cause
resonance-like enhancements to the spectra. A variety of theoretical models
have been proposed in recent years to account for these surprising
observations, including Floquet theory Potvliege and Vučić (2006), channel-
closing theory Kopold et al. (2002); Milošević et al. (2007) and models based
on the analysis of TDSE solutions Muller (1999).

Moreover, dynamic effects such as charge migration Smirnova et al. (2009);
Mairesse et al. (2010); Lépine et al. (2014); Calegari et al. (2016); Kuleff
et al. (2016), bond breaking Wolter et al. (2016), nuclear degrees of freedom
Karamatskos et al. (2019) and multielectron resonances Shiner et al. (2011);
Faccialà et al. (2016) are expected to be important for larger systems.
Extended targets also imply a structured continuum, or a blurred bound-
continuum distinction. Modeling such dynamics is not an easy task even for the
simplest possible case, i.e., a single-electron continuum under the joint
influence of the laser field and a long-range potential. The residual
potential blurs the distinction between “direct” and “rescattered” as the
electron now may be lightly deflected, or undergo a soft scattering by the
core Maxwell and Figueira de Morisson Faria (2018). This goes beyond the
clear-cut distinction imposed by Born-type approaches such as the SFA Yan et
al. (2010); Lai et al. (2015). All the above suggests the driving-field
polarization as a tool to assess what the SFA leaves out. In fact,
experimental and theoretical studies have been carried out that demonstrate
the dependence of recollision outcomes - such as HATI or HHG - on the
polarization of the incident radiation Budil et al. (1993); Dietrich et al.
(1994). These studies show that the yield is increasingly suppressed as the
ellipticity of the laser pulse is increased. The generally accepted
explanation for this observation is that an elliptical laser will imbue the
ionized electron with a non-zero drift velocity that reduces the chances of
recollision with the parent ion by the electron Corkum (1993). Subsequently,
it was demonstrated that substituting laser polarization for atomic
polarization yielded the same order of magnitude suppression in the production
of HHG Averbukh (2004). The present work will extend this line of inquiry to
the dependence of the production of HATI on the atomic polarization. We will
focus for simplicity on the helium atom prepared in an excited state of
$ns^{1}(n+1)p^{1}$ interacting with a linearly polarized field. The
polarization of the _p_ orbital can then be oriented along the laser
polarization axis or in the plane perpendicular to it.

In case the atomic polarization is oriented parallel to the laser
polarization, we expect there to be no suppression to the process outlined
above, wherein a continuum electron is driven by the laser field to scatter
off the parent ion. However, if the atomic polarization is oriented
perpendicular to the laser polarization, ionization in the direction of the
laser polarization is precluded due to the atomic orbital symmetry. Therefore,
ionization is only possible if the electron acquires a nonzero drift velocity
in the plane perpendicular to the laser polarization. This should greatly
reduce the probability of a successful recollision with the core, as the
electron’s drift velocity, unimpeded by interaction with the laser, will
displace it away from the parent ion. Thus, in the high intensity, low
frequency regime for which the recollision model is applicable, we expect the
atomic polarization to significantly affect the HATI spectrum. On the other
hand, residual central potentials may favor ionization or rescattering off the
polarization axis, which deviates from the standard predictions. For instance,
in Mauger et al. (2010), rescattering in circularly polarized fields may occur
due to the interplay of the laser field and the long-range potential.

In this work, we consider excited helium parallel and perpendicular to the
driving-field polarization, for which we calculate photoelectron spectra and
momentum distributions. We use approaches that incorporate the residual
binding potential and/or the core dynamics: the B spline algebraic
diagrammatic construction (B-spline TD-ADC), the one-electron Schrödinger
solver Qprop Mosert and Bauer (2016); Tulsky and Bauer (2020), the SFA and the
Coulomb quantum-orbit strong-field approximation (CQSFA). They will allow for
an assessment of what the simple rescattering picture leaves out. The B-Spline
TD-ADC is an ab-initio method that solves the 3D atomic many electron time-
dependent Schrödinger equation for a neutral system Ruberti et al. (2014).
Specifically for the two-electron model studied in this work, it accounts for
electron exchange and recollision-induced excitations. In the strong-field
context, the B-Spline ADC has been used in high-order harmonic generation
Ruberti et al. (2018a) and pump-probe spectroscopy Ruberti et al. (2018b). The
CQSFA is an orbit-based method that incorporates the Coulomb potential and the
external laser field on equal footing. It has a high predictive power as it
allows for quantum mechanical pathways to be switched on and off at will. In
the context of ultrafast photoelectron holography, it led to possibly the most
thorough investigations of how holographic patterns form Maxwell et al.
(2017); Maxwell et al. (2018a); Maxwell and Figueira de Morisson Faria (2018);
Maxwell et al. (2018b). This work is organized as follows. In Sec. II, we
discuss the theoretical methods used in this work. Subsequently, in Sec. III,
we present our results for initial perpendicular and parallel polarized
states. Specifically, we investigate angle-integrated spectra (Sec. III), and
photoelectron momentum distributions (Sec. III.2). Finally, our main
conclusions are presented in Sec. IV.

##  II Theory

In order to calculate ATI spectra or photoelectron-momentum distributions, we
must compute the time-dependent wave function and project it onto an
asymptotic continuum state with a well-defined momentum. Below, we will
discuss several ways to do so. We will use atomic units throughout, unless
otherwise stated, and the dipole approximation.

###  II.1 B-Spline algebraic diagrammatic construction

We calculated the ATI spectra of excited, polarised, helium atoms using the
time-dependent (TD) B-spline algebraic diagrammatic construction (ADC)  _ab
initio_ method Ruberti et al. (2014); Simpson et al. (2016); Ruberti et al.
(2018a, b); You et al. (2019). Within the TD B-spline ADC approach, the 3D
many-electron time-dependent Schrödinger equations (TDSE) for the neutral
helium atom interacting with the intense mid-IR laser field, given by

$i\hbar\frac{\partial\left|\Psi_{N}\left(t\right)\right\rangle}{\partial t}\
=\ {\hat{H}}^{N}\left(t\right)\left|\Psi_{N}\left(t\right)\right\rangle,$ (1)

is solved by making the ansatz

$\left|\Psi_{N}\left(t\right)\right\rangle=\alpha_{\text{GS}}\left(t\right)\left|\Psi_{\text{GS}}^{N}\right\rangle+\sum_{j}{\alpha_{j}\left(t\right)\left|\Psi_{j}^{N}\right\rangle}$

for the time-dependent many-electron state of neutral helium. In where
$\left|\Psi_{\text{GS}}^{N}\right\rangle$ represents the ground state of
neutral helium, while the basis states $\left|\Psi_{j}^{N}\right\rangle$ refer
to the correlated many-electron configuration states of the ADC theory for
_N_ -electron neutral systems Ruberti (2019a).

The total time-dependent Hamiltonian in Eq. (1) for the time-evolution of the
system interacting with the pulse reads

$\hat{H^{N}(t)}=\hat{H_{0}^{N}}+\hat{D_{z}}E(t).$ (2)

Here $\hat{H_{0}^{N}}$ is the field-free many-electron ADC Hamiltonian
describing the neutral system. The laser-atom interaction is described in the
length gauge as $\hat{D_{z}}E(t)$, where the z component of the dipole
operator is $\hat{D_{z}}=\sum_{j=1}^{N}\hat{z_{j}}$ and the summation over the
$j$ index runs over all the $N$ electrons of the atom.

The single-particle basis set used in this approach consists of spherical
harmonics $Y_{l,m}\left(\theta,\varphi\right)$ for the angular part and
B-spline functions  _B i_(_r_) for the radial coordinate. The single particle
basis functions used in this calculation are therefore expressed as

$\psi_{i,l,m}\ =\frac{B_{i}\left(r\right)}{r}\
Y_{l,m}\left(\theta,\varphi\right).$ (3)

The initial excited states used in the present simulations are
$\ket{\Psi_{N}\left(t=0\right)}_{\bot}=\ket{1s2p_{x}}$ and
$\ket{\Psi_{N}\left(t=0\right)}_{\parallel}=\ket{1s2p_{z}}$ in the
perpendicular and parallel set-ups, respectively. In this work, we have used
the lowest level of the ADC-hierarchy compatible with a correct description of
the strong-field ionization of excited helium by the mid-IR laser pulses, i.e.
ADC(1). Within ADC(1), the configuration manifold used to describe the
ionisation of helium by the laser pulses, via solving the TDSE, consists of
the singly excited one-hole–one-particle $(1h-1p)$ configurations. The $1h-1p$
manifold allows one to describe ionization as well as excitation of helium
atom in any of the singly-excited bound states. The time propagation of the
unknown coefficients $\alpha_{j}\left(t\right)$ of the B-spline ADC(1) many-
electron wave-function is performed by means of the Arnoldi–Lanczos, algorithm
Ruberti et al. (2018a, b); Ruberti (2019b). The B Spline ADC within ADC(1) has
been successfully used in the strong-field regime to model the intensity-
dependent interference minimum that is present in the high-order harmonic
spectra of $CO_{2}$ Ruberti et al. (2018a).

###  II.2 One-electron models

The methods below approximate the system by a one-electron system, whose
evolution is given, in atomic units, by the single-electron time-dependent
Schrödinger equation (TDSE)

$i\partial_{t}|\psi(t)\rangle=H(t)|\psi(t)\rangle,\,$ (4)

where the Hamiltonian

$H(t)=H_{a}+H_{I}(t)$ (5)

describes the joint influence of the binding potential and the external field.
Thereby,

$H_{a}=\frac{\hat{\mathbf{p}}^{2}}{2}+V(\hat{\mathbf{r}})$ (6)

gives the field-free one-electron atomic Hamiltonian and $\hat{\mathbf{r}}$
and $\hat{\mathbf{p}}$ denote the position and momentum operators,
respectively. The coupling with the field is given by the interaction
Hamiltonian $H_{I}(t)=\hat{\mathbf{r}}\cdot\mathbf{E}(t)$ and
$H_{I}(t)=\hat{\mathbf{p}}\cdot\mathbf{A}(t)+\mathbf{A}^{2}/2$ in the length
and velocity gauges, respectively, where $\mathbf{E}(t)=-d\mathbf{A}(t)/dt$ is
the electric field of the external laser field and $\mathbf{A}(t)$ the
corresponding vector potential.

Eq. (5) equation is either solved numerically using the freely available
software Qprop Mosert and Bauer (2016); Tulsky and Bauer (2020) or semi-
analytically using the SFA or the CQSFA. Throughout, we employ the effective
potential

$V(\mathbf{r}(\tau))=-\frac{1+f(r(\tau))}{r(\tau)},$ (7)

where

$f(r)=a_{1}e^{-a_{2}r}+a_{3}re^{-a_{4}r}+a_{5}e^{-a_{6}r}$ (8)

and $r(\tau)=\sqrt{\mathbf{r}({\tau})\cdot\mathbf{r}({\tau})}$. The
coefficients are chosen as $a_{1}=1.231$ a.u., $a_{2}=0.662$ a.u.,
$a_{3}=-1.325$ a.u., $a_{4}=-1.236$ a.u., $a_{5}=-0.231$ a.u. and
$a_{6}=0.480$ a.u. Maxwell et al. (2020). These parameters were obtained by
fitting to a numerical potential computed by self-interaction free density
functional theory Tong and Lin (2005); Tong and Chu (1997).

For Qprop, the TDSE is solved in the velocity gauge, while for the semi-
analytic approaches the length-gauge Hamiltonian is used. This is mainly for
practical reasons: the numerical solution of the TDSE converges faster for the
velocity gauge, while the length gauge gives better results for ATI in the SFA
Bauer et al. (2005) and CQSFA de Morisson Faria and Maxwell (2020). For Qprop
we refer to Mosert and Bauer (2016); Tulsky and Bauer (2020), while the semi-
analytic models will be summarized below. In all one-electron models, for the
perpendicular initial $p$ state, the valued real $p_{x}$ orbital is used via
the corresponding coherent superposition of states with $m=\pm 1$.

In Qprop it is possible to compute strong field ionization for the initial
states $\ket{\psi_{2p_{\pm 1}}}$, where the subscript $\pm 1$ refers to the
quantum number $m$. The ‘x’ oriented state may be computed by the following
superposition,

$\ket{\psi_{2p_{x}}}=\frac{1}{\sqrt{2}}\left(\ket{\psi_{2p_{-1}}}-\ket{\psi_{2p_{+1}}}\right)$
(9)

One can simply take the amplitude components (both real and imaginary) for the
initial states $\ket{\psi_{2p_{\pm 1}}}$ and compute the coherent
superposition above (9) to receive $\ket{\psi_{2p_{x}}}$ amplitude, which is
needed to provide the yield for the photoelectron momentum distributions.

A convenient starting point for the semi-analytic expressions is the
Schrödinger equation (5) in integral form, namely

$U(t,t_{0})=U_{a}(t,t_{0})-i\int^{t}_{t_{0}}U(t,t^{\prime})H_{I}(t^{\prime})U_{a}(t^{\prime},t_{0})dt^{\prime}\,,$
(10)

where $U_{a}(t,t_{0})=\exp[iH_{a}(t-t_{0})]$ is the time-evolution operator
associated with the field-free Hamiltonian, and the time evolution operator

$U(t,t_{0})=\mathcal{T}\exp\bigg{[}i\int^{t}_{t_{0}}H(t^{\prime})dt^{\prime}\bigg{]},$
(11)

where $\mathcal{T}$ denotes time-ordering, relates to the full Hamiltonian
$H(t)$ evolving from an initial time $t_{0}$ to a final time $t$.
Alternatively, one may construct the integral equation as

$U(t,t_{0})=U^{(V)}(t,t_{0})-i\int^{t}_{t_{0}}U(t,t^{\prime})VU^{(V)}(t^{\prime},t_{0})dt^{\prime},$
(12)

where $U^{(V)}(t,t_{0})$ is the Volkov time-evolution operator, associated
with the Hamiltonian

$H^{(V)}=\frac{\hat{\mathbf{p}}^{2}}{2}+H_{I}(t).$ (13)

We wish to calculate the transition amplitude
$\left\langle\psi_{\textbf{p}}(t)|U(t,t_{0})|\psi_{0}\right\rangle$ from a
bound state $\left|\psi_{0}\right\rangle$ to a final continuum state
$|\psi_{\textbf{p}}(t)\rangle$ with momentum $\mathbf{p}$. Using Eq. (10) and
the orthogonality between continuum and bound states, it can be written as

$M(\mathbf{p})=-i\lim_{t\rightarrow\infty}\int_{-\infty}^{t}dt^{\prime}\left\langle\psi_{\textbf{p}}(t)|U(t,t^{\prime})H_{I}(t^{\prime})|\psi_{0}(t^{\prime})\right\rangle\,,$
(14)

with
$\left|\psi_{0}(t^{\prime})\right\rangle=\exp[iI_{p}t^{\prime}]\left|\psi_{0}\right\rangle$,
where $I_{p}$ is the ionization potential and we have set $t_{0}=0$. One
should note that, apart from considering a single active electron, no
approximation has been made in the time propagation described by Eq. (14).
Some approximations and/or further assumptions will be employed next.

####  II.2.1 Strong-field approximation

The Strong-Field Approximation (SFA), also known as the Keldysh-Faisal-Reiss
(KFR) theory Keldysh (1965); Faisal (1973); Reiss (1980); Amini et al. (2019),
can be obtained by constructing a field-dressed Born series using Eq. (12),
and inserting it into Eq. (10).

Physically, the SFA relies upon the assumption that all bound states apart
from the initial state of the system do not need to be considered. That is,
the state of the system at a time $t$ can be written as

$\left|\psi(t)\right\rangle=\left|\psi_{0}(t)\right\rangle+\int
d\mathbf{p}\,b(\mathbf{p},t)\left|\mathbf{\psi_{p}}\right\rangle,$ (15)

where $\left|\mathbf{\psi_{p}}\right\rangle$ are the continuum states, and
$\left|\psi_{0}(t)\right\rangle$ is the population in the ground state, often
written $\left|\psi_{0}(t)\right\rangle=a(t)\ket{\psi_{0}}$. If the laser
field is not intense enough for the ionization to reach the saturation level,
a widespread approximation is $a(t)=e^{iI_{p}t}$, which means that the ground-
state depletion can be neglected.

#### Transition amplitudes

For ATI, the lowest nonvanishing term of the Born-type series is obtained if
the full time evolution operator is replaced by the Volkov time evolution
operator $U^{(V)}(t,t^{\prime})$ in Eq. (14). The subsequent term will include
a single interaction in Eq. (12). This will lead to two terms that give the
amplitude for producing a photoelectron with a particular well-defined
momentum,

$M(\mathbf{p})=M_{\mathrm{Dir}}(\mathbf{p})+M_{\mathrm{Resc}}(\mathbf{p}).$
(16)

The term $M_{\mathrm{Dir}}(\mathbf{p})$ can be physically interpreted as a
direct contribution from electrons that reach the detector with no further
interaction with the core, while $M_{\mathrm{Resc}}(\mathbf{p})$ is associated
with electrons that are driven back to the parent ion by the laser field and
are subsequently elastically rescattered.

The expression for the direct term is

$M_{\mathrm{Dir}}(\mathbf{p})=\frac{-i}{(2\pi)^{3/2}}\int_{\infty}^{\infty}dtE(t)e^{iS(\mathbf{p},t)}e^{iI_{p}t}\int
d\mathbf{r}e^{-i\left[\mathbf{p}+\mathbf{A}(t)\right]\cdot\mathbf{r}}r\cos\theta\psi_{0}(\mathbf{r})$
(17)

where $\theta$ is the polar angle, $H_{I}(t)$ is the length-gauge atom field
interaction and the action reads

$S(\mathbf{p},t)=\frac{1}{2}\int_{-\infty}^{t}d\tau[\mathbf{p}+\mathbf{A}(\tau)]^{2}.$
(18)

The rescattered term is given by

$\displaystyle M_{\mathrm{Resc}}(\mathbf{p})$
$\displaystyle=-\int_{-\infty}^{\infty}dt\int_{-\infty}^{t}dt^{\prime}E(t^{\prime})e^{iI_{p}t^{\prime}}e^{-iS(\mathbf{k},t)}e^{iS(\mathbf{k},t^{\prime})}e^{iS(\mathbf{p},t)}$
$\displaystyle\times\int d\mathbf{k}\frac{1}{(2\pi)^{3}}\int
d\mathbf{r^{\prime}}e^{i(\mathbf{k}-\mathbf{p})\cdot\mathbf{r^{\prime}}}V(\mathbf{r^{\prime}})\times\frac{1}{(2\pi)^{3/2}}\int
d\mathbf{r}e^{-i(\mathbf{k}+\mathbf{A(t^{\prime})})\cdot\mathbf{r}}r\cos\theta\psi_{0}(\mathbf{r})$
(19)

The three steps in the rescattering process delineated in Eq. (19) are: (i)
ionization at time $t$ to a Volkov state with canonical momentum $\mathbf{k}$;
(ii) propagation in the laser field from $t^{\prime}$ to $t$, where
$t>t^{\prime}$; (iii) scattering off the ion to a final Volkov state with
canonical momentum $\mathbf{p}$ at time $t$, where $V$ can be a model
potential for the ionic system in question. It is very common to calculate the
SFA transition amplitudes employing saddle-point methods, as they simplify the
numerical effort and provide a very intuitive, clear picture in terms of
rescattering electron trajectories (for details see, e.g., Figueira de
Morisson Faria et al. (2002)). However, the standard saddle-point method
selects the rescattering processes that are exactly along the driving-field
polarization. This will be problematic for initial bound states aligned
perpendicularly to the field polarization. A discussion of how to compute
these matrix elements explicitly and how to implement saddle-point methods in
the present context is provided in the appendix.

For early discussions of how to incorporate rescattering in ATI see Refs.
Lewenstein et al. (1995) and Lohr et al. (1997). Further detail about the SFA
and strong-field ionization is provided in the tutorial Ivanov et al. (2005)
and in the review articles Popruzhenko (2014); Becker et al. (2018).
Introducing the term $M_{\mathrm{Resc}}(\mathbf{p})$ is also known as the
“improved strong-field approximation” Hasović et al. (2007).

#### Bound-state transitions

The standard formulation of the SFA method assumes that the active electron
occupies a specific ground state and that it can only transition into
continuum states through interaction with the laser. The influence of other
bound states on the dynamics of the system is typically not considered.

However, if the initial state of the system is excited to a point where it is
within energetic reach of other excited states, this assumption may become
problematic. In the case of the helium atom excited into the 21 _P_ state,
particular attention should be paid to the transition dipole between the 21
_P_ state and the 21 _S_ state. Here, the energy gap is relatively small,
approximately 0.0221 a.u (0.6 eV), and the transition dipole moment is large,
approximately 2.9 a.u. In principle, other possible bound-state transitions,
such as 21 _P_ – 11 _S_ or 21 _P_ – 31 _S_ can be more safely ignored as the
energy difference that characterizes these transitions is significantly
larger. Excitation has also been incorporated in the SFA in the context of
NSDI Shaaran et al. (2010); Shaaran and de Morisson Faria (2010); Wang et al.
(2012); Hao et al. (2014); Maxwell and Faria (2015, 2016). For coherent
superpositions of states in the SFA which, however, are only coupled via the
continuum see, e.g., Figueira de Morisson Faria and Augstein (2010).

To incorporate this into the SFA model above, we modify the system’s state
space vector (see Eq. (15)) as follows:

$\left|\Psi(t)\right\rangle=a_{2^{1}S}(t)\left|\psi_{2^{1}S}\right\rangle+a_{2^{1}P}(t)\left|\psi_{2^{1}P}\right\rangle+\int
d\mathbf{p}\hskip
2.84544ptb(\mathbf{p},t)\left|\psi_{\mathbf{p}}\right\rangle.$ (20)

This would modify the first term in Eq. (16) to

$M_{\mathrm{Dir}}(\mathbf{p})=\frac{-i}{(2\pi)^{3/2}}\sum_{j}\int_{\infty}^{\infty}dtE(t)e^{iS(\mathbf{p},t)}a_{j}(t)\int
d\mathbf{r}e^{-i\left[\mathbf{p}+\mathbf{A}(t)\right]\cdot\mathbf{r}}r\cos\theta\psi_{j}(\mathbf{r}),$
(21)

where the index $j$ runs over the excited states to be included in the
calculation. The same modification is made to the expression for the term
$M_{\mathrm{Resc}}(\mathbf{p})$ in Eq. (16). The model now describes
ionization from a series of strongly coupled bound states. We determine the
time-dependent behaviour of the bound states in a laser field with the system
of coupled equations:

$\frac{da_{j}}{dt}=I_{p_{j}}a_{j}(t)-\sum_{i\neq
j}\mathbf{E}(t)\cdot\mathbf{d}_{j-i}a_{i}(t),$ (22)

where $I_{p_{j}}$ is the ionization potential for the excited state $j$ and
$\mathbf{d}_{j-i}$ is the matrix element of the dipole transition from one
excited state to another with the orientation of this element determined by
the spatial polarization of the excited orbitals. If the number of excited
states under consideration is two and the quantity
$2\mathbf{E}_{0}\cdot\mathbf{d}_{j-i}$ (where $\mathbf{E}_{0}$ is the maximum
field strength) is much greater than the bound-state energy difference between
the states $I_{P_{j}}-I_{P_{i}}$, then the system of equations in (22) can be
solved perturbatively Ivanov and Corkum (1993).

The solution for our helium system is

$\displaystyle a_{2^{1}S}^{0}(t)$
$\displaystyle=ie^{-i\tilde{I}_{p}t}\sin\left[\int_{0}^{t}\mathbf{E}(t^{\prime})\cdot\mathbf{d}_{2^{1}S-2^{1}P}dt^{\prime}\right]$
(23) $\displaystyle a_{2^{1}P}^{0}(t)$
$\displaystyle=e^{-i\tilde{I}_{p}t}\cos\left[\int_{0}^{t}\mathbf{E}(t^{\prime})\cdot\mathbf{d}_{2^{1}S-2^{1}P}dt^{\prime}\right],$
(24)

where $\tilde{I}_{p}$ is the average energy of the two excited states and the
0 in the superscript indicates that this solution is to zeroth order. This
means that the energy difference is neglected and the states are considered as
degenerate. The condition for this solution
$2\mathbf{E}_{0}\cdot\mathbf{d}_{2^{1}S-2^{1}P}\gg
I_{P_{2^{1}P}}-I_{P_{2^{1}S}}$ is not valid in the case of a finite pulse at
the edges of the envelope, but remains applicable provided an alternative
condition is met: that $\omega\gg I_{P_{2^{1}P}}-I_{P_{2^{1}S}}$ where
$\omega$ is the field frequency. It can be readily seen from Eqs. (23) and
(24) that the solution for the perpendicular configuration is

$\displaystyle a_{2^{1}S}^{0}(t)$ $\displaystyle=0$ $\displaystyle
a_{2^{1}P}^{0}(t)$ $\displaystyle=e^{-i\tilde{I}_{p}t}.$

This is because the 21 _P_ – 21 _S_ transition is not possible due to the
symmetry constraints of the dipole transition in the case of the perpendicular
configuration.

####  II.2.2 The Coulomb quantum-orbit strong-field approximation

Instead of constructing iterative schemes around Eq. (10), one of which is the
SFA, one may also employ path-integral methods and construct approaches that
incorporate the field and the binding potential on equal footing. One of such
approaches is the Coulomb quantum-orbit strong-field approximation (CQSFA) Lai
et al. (2015); Maxwell et al. (2017). It is constructed from the transition
amplitude

$\displaystyle M(\mathbf{p}_{f})=$
$\displaystyle-i\lim_{t\rightarrow\infty}\int_{-\infty}^{t}dt^{\prime}\int
d\mathbf{\tilde{p}}_{0}\left\langle\mathbf{\tilde{p}}_{f}(t)|U(t,t^{\prime})|\mathbf{\tilde{p}}_{0}\right\rangle\left\langle\mathbf{\tilde{p}}_{0}|H_{I}(t^{\prime})|\psi_{0}(t^{\prime})\right\rangle,\,$
(25)

which is obtained by inserting a closure relation in the initial field-dressed
momentum $\mathbf{\tilde{p}}_{0}=\mathbf{p}_{0}+\mathbf{A}(t^{\prime})$ in Eq.
(10) and is exact for a one-electron system. Eq. (25) takes the system from
the initial photoelectron bound state
$\left|\psi_{0}(t^{\prime})\right\rangle\,$ to its final momentum state
$\ket{\psi_{\mathbf{p}}(t)}=\left|\mathbf{\tilde{p}}_{f}(t)\right\rangle\,$,
where the tilde indicates field-dressing. The
$\mathbf{\tilde{p}}_{0}=\mathbf{p}_{0}+\mathbf{A}(t^{\prime})$ and
$\mathbf{\tilde{p}}_{f}(t)=\mathbf{p}_{f}+\mathbf{A}(t)$, where
$\mathbf{A}(\tau)$, with $\tau=t,t^{\prime}$ is the vector potential, denote
the electron’s initial and final field-dressed momenta, respectively. The
time-evolution operator $U(t,t^{\prime})$ is associated to the full
Hamiltonian, and thus includes the laser field and  the binding potential. We
take the interaction Hamiltonian $H_{I}(t^{\prime})$ to be in the length,
gauge. Applying time-slicing techniques Kleinert (2009); Milošević (2013a)
gives

$\displaystyle M(\mathbf{p}_{f})$
$\displaystyle=-i\lim_{t\rightarrow\infty}\int_{-\infty}^{t}dt^{\prime}\int
d\mathbf{\tilde{p}}_{0}\int_{\mathbf{\tilde{p}}_{0}}^{\mathbf{\tilde{p}}_{f}(t)}\mathcal{D}^{\prime}\mathbf{\tilde{p}}\int\frac{\mathcal{D}\mathbf{r}}{(2\pi)^{3}}$
$\displaystyle\hskip 56.9055pt\times
e^{iS(\mathbf{\tilde{p}},\mathbf{r},t,t^{\prime})}\langle\mathbf{\tilde{p}}_{0}|H_{I}(t^{\prime})|\psi_{0}\rangle\,,$
(26)

where $\mathcal{D}^{\prime}\mathbf{p}$ and $\mathcal{D}\mathbf{r}$ are the
integration measures for the path integrals Kleinert (2009); Lai et al.
(2015), and the prime denotes a restriction. The action and the Hamiltonian in
Eq. (26) reads as

$S(\mathbf{\tilde{p}},\mathbf{r},t,t^{\prime})=I_{p}t^{\prime}-\int^{t}_{t^{\prime}}[\dot{\mathbf{p}}(\tau)\cdot\mathbf{r}(\tau)+H(\mathbf{r}(\tau),\mathbf{p}(\tau),\tau]d\tau,$
(27)

and

$H(\mathbf{r}(\tau),\mathbf{p}(\tau),\tau)=\frac{1}{2}\left[\mathbf{p}(\tau)+\mathbf{A}(\tau)\right]^{2}+V(\mathbf{r}(\tau)),$
(28)

respectively. In Eqs. (27) and (28), the intermediate momentum $\mathbf{p}$
and coordinate $\mathbf{r}$ have been parameterized with regard to the time
$\tau$ and $\mathbf{\tilde{p}}=\mathbf{p}(\tau)+\mathbf{A}(\tau)$.

One should note that no Born-type expansion was used in the derivation of Eq.
(26). This implies that the CQSFA cannot be viewed as a field-dressed
perturbative series with regard to the binding potential, and that the
distinction between direct and rescattered electrons is blurred. This is even
clearer when Eq. (26) is solved using saddle-point methods, i.e., we search
for values of $t^{\prime}$, $\mathbf{r}$ and $\mathbf{p}$ that render the
action (27) stationary. This gives

$\frac{\left[\mathbf{p}(t^{\prime})+\mathbf{A}(t^{\prime})\right]^{2}}{2}+V(\mathbf{r}(t^{\prime}))+\dot{\mathbf{p}}(t^{\prime})\cdot\mathbf{r}(t^{\prime})=-I_{p},$
(29) $\mathbf{\dot{p}}=-\nabla_{r}V(\mathbf{r}(\tau))$ (30)

and

$\mathbf{\dot{r}}=\mathbf{p}+\mathbf{A}(\tau).$ (31)

Eq. (29) yields the energy conservation at tunnelling, and Eqs. (30) and (31)
give the subsequent electron propagation.

We compute Eq. (25) along a two-pronged contour that starts in at the complex
time $t^{\prime}=t^{\prime}_{r}+it^{\prime}_{i}$, goes parallel to the
imaginary time axis, i.e., $t^{\prime}$ to $t^{\prime}_{r}$, and,
subsequently, along the real time axis from the real ionization time
$t^{\prime}_{r}$ to the final time $t\rightarrow\infty$ Popruzhenko and Bauer
(2008); Yan and Bauer (2012); Torlina and Smirnova (2012); Torlina et al.
(2013). The first and second arms of the contour are associated with under-
the-barrier dynamics and continuum propagation, respectively.

In the under-the-barrier part of the contour, we kept the momentum constant
and equal to $\mathbf{p}_{0}$ Maxwell et al. (2018b); de Morisson Faria and
Maxwell (2020). This reduces Eq. (29) to

$\frac{\left[\mathbf{p}_{0}+\mathbf{A}(t^{\prime})\right]^{2}}{2}=-I_{p}.$
(32)

In the continuum propagation, the action can be simplified according to the
procedure discussed in Lai et al. (2015); Maxwell et al. (2017); Shvetsov-
Shilovski et al. (2016); Maxwell et al. (2020).

The CQSFA transition amplitude obtained using the saddle-point approximation
reads

$M(\mathbf{p}_{f})\propto-i\lim_{t\rightarrow\infty}\sum_{s}\bigg{\\{}\det\bigg{[}\frac{\partial\mathbf{p}_{s}(t)}{\partial\mathbf{r}_{s}(t_{s})}\bigg{]}\bigg{\\}}^{-1/2}\mathcal{C}(t_{s})e^{iS(\mathbf{\tilde{p}}_{s},\textbf{r}_{s},t,t_{s})},$
(33)

where $\mathbf{p}_{s}$ and $\mathbf{r}_{s}$ and $t_{s}$ are determined by Eqs.
(30),(31) and (32), respectively. The term in brackets gives the stability of
a specific trajectory, and

$\mathcal{C}(t_{s})=\sqrt{\frac{2\pi
i}{\partial^{2}S(\mathbf{\tilde{p}}_{s},\textbf{r}_{s},t,t_{s})/\partial
t^{2}_{s}}}\langle\mathbf{p}+\mathbf{A}(t_{s})|H_{I}(t_{s})|\psi_{0}\rangle.$
(34)

In practice, we use $\partial\mathbf{p}_{s}(t)/\partial\mathbf{p}_{s}(t_{s})$
instead of the stability factor in Eq. (33). This may be obtained employing a
Legendre transformation and will not affect the action if the electron starts
from the origin. The continuum propagation is performed by solving the inverse
problem, i.e., given a final momentum $\mathbf{p}_{f}(t)$, we find the initial
momentum $p_{0}$ at the tunnel exit. This is defined as
$\mathbf{r}_{0}(t^{\prime}_{r})$, where

$\mathbf{r}_{0}(\tau)=\int_{t^{\prime}}^{\tau}(\mathbf{p}_{0}+\mathbf{A}(\tau^{\prime}))d\tau^{\prime}.$
(35)

A further approximation used here is to take a real tunnel exit

$z_{0}=\mathrm{Re}[r_{0z}(t^{\prime}_{r})],$ (36)

where $r_{0z}$ indicates the component of the tunnel trajectory along the
field polarization direction. This will simplify the problem and lead to real
trajectories in the continuum. For discussions of how to solve the full
complex problem see Pisanty and Ivanov (2016); Maxwell et al. (2018b).

The term $\langle\mathbf{p}+\mathbf{A}(t_{s})|H_{I}(t_{s})|\psi_{0}\rangle$
contains the influence of the initial bound state, i.e., its geometry, and
will be important for the present work. Explicitly,

$\langle\mathbf{p}|H_{I}(t_{s})|\psi_{0}\rangle=\int
d^{3}\mathbf{r}\frac{\exp(-i\mathbf{p}\cdot\mathbf{r})}{(2\pi)^{3/2}}\mathbf{r}\cdot\mathbf{E}(t_{s})\psi_{0}(\mathbf{r}),$
(37)

where

$\psi_{0}(\mathbf{r})=\sum_{\alpha}\frac{c_{\alpha}}{\sqrt{2}(2\pi)^{3/2}}\phi_{\alpha}(\mathbf{r}),$
(38)

and $\phi_{\alpha}(\mathbf{r})$ is given by a Gaussian basis set

$\phi_{\alpha}(\mathbf{r})=\sum_{v}b_{v}x^{\beta_{x}}y^{\beta_{y}}z^{\beta_{z}}\exp\left(-\xi_{v}(x^{2}+y^{2}+z^{2})\right).$
(39)

For details see Figueira de Morisson Faria and Augstein (2010).

In the CQSFA, one may identify four main types of orbits:

  1. 1.

Type 1 orbits behave like short, direct SFA trajectories that are slowed down
by the Coulomb potential. They leave in the direction of the detector and
their transverse momentum does not change direction.

  2. 2.

Type 2 orbits are field-dressed Kepler hyperbolae whose start times are
displaced in half a cycle from those in type 1 orbits. They start on the
“wrong” side and are turned towards the detector. They are lightly deflected
by the potential and can be loosely associated with the long SFA orbits. Their
transverse momentum does not change its sign during continuum propagation.

  3. 3.

Type 3 orbits are also field-dressed Kepler hyperbolae starting in the same
half cycle as type 2 orbits, but interact more strongly with the core. Their
transverse momentum changes direction during continuum propagation. They have
no counterpart in the SFA as they depend on the interplay between the driving
field and the binding potential and are not embedded in a Born-type method. In
previous work, we have approximated orbit 3 within an analytic piecewise model
in which rescattering is incorporated. However, mimicking the effect of orbit
3 and accurately reproducing the spider was only possible by including the
acceleration caused by the Coulomb potential during the continuum propagation
and the Coulomb phase Maxwell et al. (2018a). This is a substantially
different scenario than that encountered with the SFA, be it direct or
rescattered, or in approaches focused on short-range potentials Zhou et al.
(2016). In the limit of vanishing residual potential for the continuum
propagation, orbits 2 and 3 become degenerate and tend to the long direct SFA
orbit Lai et al. (2015).

  4. 4.

Type 4 orbits start in the same half cycle as type 1 orbits, but go around the
core before reaching the continuum. Their behaviour is close to that of SFA
rescattered trajectories, and they may reach energies of around $10U_{p}$.
However, one should notice that the process they undergo is different and
their propagation in the continuum is not restricted to the field-polarisation
axis.

One should note that the energies of orbits 2 and 3 may go beyond the direct
ATI cutoff of $2U_{p}$, due to the influence of the residual potential. This
contributes to blur the distinction between direct and rescattered, as
mentioned above. A limitation of the CQSFA is that it does not properly
account for processes involving excited bound states. The basis chosen in Eq.
(25) by using the closure relation in $\tilde{\mathbf{p}}_{0}$ works well for
continuum states but does not reproduce excitation processes accurately.
Still, it does take into account highly excited states which are strongly
coupled with the field and behave similarly to the continuum. The four types
of orbits stated above were first identified in Yan et al. (2010).

##  III Results

In the following, we will analyse how different polarisations in the initial
bound states influence angle-integrated photoelectron spectra and
photoelectron momentum distributions, with emphasis on several
counterintuitive features. We used a bandwidth-limited mid-IR pulse of
intensity $I=3.2\times 10^{13}$ W/cm2, frequency $\omega=0.775$ eV, that
corresponds to a wavelength $\lambda=1600$ nm, linear polarisation along the z
direction, with a ponderomotive energy ($U_{p}$) of $7.65$ eV. Unless
otherwise stated, the binding energy used for this initial state is $3.286$ eV
($I_{p}=0.1238$ a.u.). For this parameter range, we have verified that the
electron in the ground state orbital 1s is not appreciably affected by the IR
pulse in the TDSE solvers.

For the B-Spline ADC method and the SFA, the pulse used is described by a
cosine squared envelope for the vector potential

$\mathbf{A}(t)=-\frac{A_{0}}{\omega}\cos^{2}\left(\frac{\pi
t}{T}\right)\sin(\omega t)\hat{e}_{z},\;\\{-T/2\leq t\leq T/2\\}$ (40)

of amplitude $A_{0}$ and polarization vector $\hat{e}_{z}$, and with a pulse
duration $T$ set to 10 cycles (ca. 50 fs).

The calculations have been performed using a linear B-spline knot sequence
Ruberti et al. (2014) with a radial box radius of Rmax = 500 a.u. and  _N b_ =
625 radial B-splines. The maximum angular momentum employed in the monocentric
expansion of Eq. (3) was  _l_ max = 30. Convergence of the results with
respect to the basis set parameters has been checked. A complex absorbing
potential (CAP) with starting radius RCAP = 400 a.u. has been used to absorb
the wavefunction and avoid its reflections from the grid boundary. Each and
every B-spline ADC(1) energy-dependent photoelectron spectra presented here
were calculated by coupling the TD B-spline ADC method to the t-surff
technique.

For the CQSFA, for simplicity we consider a linearly polarized monochromatic
wave, whose vector potential is described by

$\mathbf{A}(t)=A_{0}\cos{\omega t}\hat{e}_{z}.$ (41)

This is a good approximation for long enough pulses. In practice, we take the
field to be over 20 cycles long, and include ionization events within up to
four laser cycles, that is, the range of ionization times $t^{\prime}$ lies
within $0\leq t^{\prime}\leq 8\pi/\omega$. A restricted temporal unit cell is
necessary for practical purposes, but may introduce some boundary effects. For
a complete discussion see Werby et al. (2021). For Qprop, we have employed a
four-cycle flat-top pulse with additional linear half-cycle turn on and off
for the comparison to the CQSFA and SFA models. In the comparison of Qprop
with the B-spline ADC model, we have matched the appropriate pulse of
$\sin^{2}$, whose vector potential is given by Eq. (40).

###  III.1 ATI Spectra

####  III.1.1 One and two-electron Schrödinger solvers

We will start by discussing the angle-integrated ATI spectra obtained from the
numerical solutions of the TDSE employed in this work, calculated with the
diagrammatic B-Spline ADC method and with the freely available one-electron
Schrödinger solver Qprop Mosert and Bauer (2016); Tulsky and Bauer (2020).
These results are presented in Fig. 1, for the parallel and perpendicular
initial-state configurations (upper and lower panel, respectively). For the
parallel configuration, the spectra obtained with the B-Spline ADC method
exhibits a long ramp-like structure extending from near the ionization
threshold up to the high photoelectron energy of $10U_{p}$, while the Qprop
computations exhibit a low-energy region followed by a ramp-like decrease from
$2U_{p}$ to $4U_{p}$ (see pink shaded region), with a flat plateau up to
$10U_{p}$. The high-energy suppression, in the case of TD B-Spline ADC, is
caused by the “potential” that the excited/ionised electron experiences, which
has contributions from both the Coulomb and exchange terms calculated with the
Hartree-Fock orbitals of the helium atom. These exchange terms effectively
repel the rescattered electron if it gets too close to the core, and thus
render rescattering less effective. This influence increases with the
photoelectron final energy, as higher-energy electrons must approach the core
more closely. These ramp-like features are not altered significantly if the
bound-state transition $2^{1}P-2^{1}S$ is removed from the B-Spline ADC
computation (see blue and green curves).

Figure 1: Angle-integrated photoelectron spectra computed for excited helium
($I_{p}=0.1238$ a.u.) with numerical solutions of the time-dependent
Schrödinger equation, for pulses with peak intensity $I_{0}=3.2\times 10^{13}$
W/cm2, ($U_{p}=7.65$ eV ), and wavelength $\lambda=1600$ nm. The orange
(dashed) curve gives the one-electron computation, performed with the free
software Qprop Mosert and Bauer (2016); Tulsky and Bauer (2020), while the
green (solid) and the blue (dotted) curves have been computed with the
diagrammatic B-Spline ADC method including and excluding the bound-state
transition $2^{1}P-2^{1}S$, respectively. For the Qprop calculation a 4-cycle
$\sin^{2}$ pulse was used, while for the B-Spline ADC a 10-cycle $\sin^{2}$
pulse was taken and pulse length of 53 fs. The upper and lower panels have
been calculated for $z$ (parallel) and $x$ (perpendicular) orientation of the
initial $2p$ bound states of helium with regard to the driving-field
polarization. The blue curve in the lower panel has been computed for parallel
polarization, but removing the bound-state transitions, the most important
among them being $2^{1}P-2^{1}S$. Each curve has been normalized to its
largest value to facilitate a qualitative comparison. The pink and violet
shaded areas highlight the photoelectron energy ranges $2U_{p}\leq E_{k}\leq
4U_{p}$ and $4U_{p}\leq E_{k}\leq 6U_{p}$, respectively.

On the other hand, if the initial $2p$ state is aligned perpendicular to the
driving-field polarization (lower panel in Fig. 1), the computations lead to
strikingly different results. For clarity, we have included the parallel-
polarized spectrum without the $2^{1}P-2^{1}S$ transition. One should note
that, due to selection rules, this transition is forbidden in the
perpendicular configuration. Up to the photoelectron energy of roughly
$2.5U_{p}$, there are at most subtle differences, with the Qprop spectrum
being slightly larger than the B-Spline ADC spectra. Subsequently, one sees a
steep ramp-like suppression for the (one-electron) Qprop computation up to
$4U_{p}$ succeeded by a second ramp of much gentler slope. The B-Spline ADC
spectra are significantly less suppressed for energies higher than $2U_{p}$,
regardless of configuration. Thus, in contrast to the parallel-aligned case,
electron exchange appears to enhance rescattering. The results suggest that
exchange effects, which are included in the B-spline ADC simulation, lead to a
decrease of the dependence of the yield on the direction of the driving field
with respect to the initial atomic polarization, especially in the
(approximately) $2U_{p}-6U_{p}$ energy region. Interestingly, for the
perpendicularly polarized initial state the B-Spline ADC computation exhibits
multiple ramps instead of the monotonic behavior identified for parallel
alignment. Thereby, one may identify two ramps, for photoelectron energy
ranges $2U_{p}\leq E_{k}\leq 4U_{p}$ and $4U_{p}\leq E_{k}\leq 6U_{p}$,
respectively. Furthermore, in comparison to the parallel-aligned case (blue
curve), the B-Spline ADC spectrum computed for the perpendicular configuration
(green curve) exhibits a decrease of one order of magnitude for energies
higher than $4U_{p}$ even in the absence of the $2^{1}P-2^{1}S$ transition.
Still, this behaviour occurs for much larger energies than those in the Qprop
spectrum. The remarkably different behaviour from Qprop suggests a
multielectron character for the second ramp as well.

####  III.1.2 Comparison with orbit-based methods

Figure 2: Angle-integrated photoelectron spectra computed for excited helium
($I_{p}=0.1238$ a.u.) using Qprop (orange and brown dotted curves for parallel
and perpendicular cases, respectively) and the strong-field approximation
(SFA) (remaining curves), for the same peak-field intensities and wavelengths
as in the previous figures. For the Qprop calculation a flat-top four cycle
pulse with an additional half a cycle turn on and off was used, while for the
SFA we have taken only a single cycle. Each curve has been normalized to its
largest value to facilitate a qualitative comparison. For the perpendicular
configuration, the SFA spectrum is one order of magnitude smaller than its
parallel counterpart in the direct region (energies smaller than $2U_{p}$ and
the rescattered plateau is two orders of magnitude smaller. The top panel
(parallel alignment) also shows all scattering (Direct and Rescattered)
contributions in the SFA for the bound coupled state, $2s2p$ (purple wide
dashed curve). The pink and violet shaded areas highlight the photoelectron
energy ranges $2U_{p}\leq E_{k}\leq 4U_{p}$ and $4U_{p}\leq E_{k}\leq 6U_{p}$,
respectively.

For that reason, we will now compare the outcome of Qprop with that of semi-
analytical, orbit-based methods. We will start by discussing the spectra
obtained with the strong-field approximation (SFA), which is shown in Fig. 2.
In all cases, the energy regions in the SFA spectra can be clearly associated
with direct and rescattered ATI, with a distinct cutoff at photoelectron
energy $2U_{p}$ followed by long plateaus up to the rescattered ATI cutoff of
$10U_{p}$. The direct contributions dominate up to approximately $4U_{p}$,
with the plateau prevailing for higher energies. For perpendicular alignment,
the plateau is much more suppressed than in the parallel-aligned case. This is
consistent with the picture put across by the SFA, that rescattering occurs
for a small angular range around the polarization axis. Because there is not
much probability density with which the returning electronic wave packet can
overlap, the resulting signal will be suppressed. This adds up to an overall
suppression that occurs for similar reasons at the instant of ionization.
Since tunneling is highly directional and occurs predominantly along and in
the vicinity of the polarization axis, the overall signal will be around two
orders of magnitude weaker in the perpendicular case. We have, however,
normalized the yield in the lower panels of the figure in order to perform a
qualitative comparison.

In comparison with the Qprop calculation, the SFA plateau is strongly
suppressed for parallel orientation, although they follow the same trend (see
upper panel in Fig. 2). In fact, it is noteworthy that there are only two
orders of magnitude difference between the direct and rescattered signal in
the Qprop case, while for the SFA this difference amounts to four orders of
magnitude. Moreover, the behavior in the direct region is much flatter than
for the SFA, with a ramp for photoelectron energies between $2U_{p}$ and
$4U_{p}$ (see pink shaded region). Interestingly, including the
$2^{1}P-2^{1}S$ bound-state transition in the SFA as discussed in Sec. II.2
leads to subtle changes in the plateau, but does not alter its overall shape
or intensity.

For perpendicular-aligned initial states (lower panel in Fig. 2), the Qprop
and SFA computations are much more similar, with a steep decay in several
orders of magnitude for $2U_{p}\leq E_{k}\leq 4U_{p}$ (dubbed ‘mid-energy
ramp’) followed by a plateau. Still, in that energy region the TDSE result
always follows that of the SFA from above, and exhibit a ramp-like structure.
One should note that this mid-energy ramp appears in a single-active electron
setting. This is a key difference with regard to the behavior seen for the
second ramp in the B-Spline ADC computations; see discussion of Fig. 1. The
mid-energy ramp is associated with hybrid orbits, which are missed by the SFA,
and will be analyzed next.

Figure 3: Angle-integrated photoelectron spectra computed for excited helium
($I_{p}=0.1238$ a.u.) in initial $2p$ states oriented parallel (z) and
perpendicular (x) to the laser-field polarization, using Qprop (dotted orange
and brown lines) and the coulomb quantum-orbit strong-field approximation
(CQSFA) (purple and red lines, dashed and solid respectively), for the same
peak-field intensities and wavelengths as in Fig. 1. For the Qprop
calculation, a flat-top four cycle pulse with an additional half a cycle turn
on and off was used, while for the CQSFA, we have considered ionization events
within one cycle of a monochromatic field. Each curve has been normalized to
its largest value to facilitate a qualitative comparison. The pink and violet
shaded areas highlight the photoelectron energy ranges $2U_{p}\leq E_{k}\leq
4U_{p}$ and $4U_{p}\leq E_{k}\leq 6U_{p}$, respectively.

In Fig. 3, we plot the spectra obtained with Qprop against those from the
CQSFA. The spectra obtained with the CQSFA exhibit an almost constant, much
higher plateau, whose onset occurs for energies immediately higher than
$2U_{p}$. This is in stark contrast with the ramp-like behavior obtained in
the previous cases. Furthermore, a direct comparison of the parallel and
perpendicular configurations shows a suppression of at least one order of
magnitude for the plateau in the latter case. This is much less than that
observed in all other methods, and it is also curious that the parallel and
perpendicular spectra exhibit very similar shapes. This suggests that, for the
parameter range employed, and in particular the ionization potential $I_{p}$ =
0.1238 a.u. of the excited states, the CQSFA orbits are behaving in an unusual
way, which is much more dictated by the binding potential than by the laser
field.

A further investigation is presented in Fig. 4, in which, in addition to the
previous curves, calculated for excited helium, we have included a CQSFA
computation for a slightly larger ionization potential ($I_{p2}$ = 0.1750
a.u.) and $2p$ initial states parallel and perpendicular to the laser-field
polarization. This is quite different from the behavior we have observed if
the correct ionization potential ($I_{p1}$ = 0.1238 a.u.) is taken. Thereby,
one can see that, for the perpendicular and parallel-oriented cases, the
spectra practically overlap up to energies up to $2U_{p}$. However, for the
perpendicular case there is a long ramp in the photoelectron energy region
$2U_{p}\leq E_{k}\leq 4U_{p}$ followed by a strongly suppressed ATI plateau.
This behavior is in better qualitative agreement with the Qprop outcome than
that obtained for the correct ionization potential. Next, we will analyze
these behaviors in more detail using photoelectron momentum distributions.

Figure 4: Angle-integrated photoelectron spectra computed for excited helium
using Qprop (orange and brown dotted curves) and the Coulomb quantum-orbit
strong-field approximation (CQSFA) (purple and red curves, wide dashed and
solid respectively) computed for excited helium ($I_{p1}=0.1238$ a.u.) using
the same parameters as in Fig. 3, but with two additional (blue and green,
dot-dahsed and small dashes respectively) curves computed using the CQSFA for
a slightly larger ionization potential ($I_{p2}=0.1750$ a.u.). Each curve has
been normalized to its largest value to facilitate a qualitative comparison,
except for those computed for the larger ionization potential $I_{p2}$. The
pink and violet shaded areas highlight the photoelectron energy ranges
$2U_{p}\leq E_{k}\leq 4U_{p}$ and $4U_{p}\leq E_{k}\leq 6U_{p}$, respectively.

###  III.2 Photoelectron momentum distributions

In the following, we will have a closer look at the electron dynamics by
inspecting photoelectron momentum distributions (PMDs). For simplicity, we
restrict ourselves to a one-electron scenario employing the CQSFA and using
Qprop Mosert and Bauer (2016); Tulsky and Bauer (2020) as a benchmark. These
results are presented in Fig. 5, for parallel and perpendicularly polarized
helium (upper and lower rows, respectively) in a four-cycle trapezoidal pulse.
There is an overall suppression of the signal in the perpendicular
configuration, but in the present discussion we will focus on qualitative
features and the physics behind them. All distributions exhibit well-defined
ATI rings and several holographic structures, such as a spider-like pattern
around the polarization axis, fan-shaped fringes close to the ionization
threshold and spiral-like patterns perpendicular to the driving-field
polarization. These structures have been discussed in previous publications
Maxwell et al. (2017); Maxwell et al. (2018a); Maxwell and Figueira de
Morisson Faria (2018); Maxwell et al. (2018b) and have been observed in
numerical computations and experiments (for details see the recent review
article de Morisson Faria and Maxwell (2020)). There are, however, some
qualitative differences between the CQSFA and the TDSE computations. First,
the holographic patterns in the TDSE outcome are asymmetric with regard to the
$p_{\perp}$ axis and the ATI rings are less sharp. This asymmetry is very
likely caused by a strong bound-state depletion, which is present in the TDSE
case. Depletion will affect the probabilities associated with different
ionization events such that the later events will be suppressed, and this will
cause distortions in the patterns, making the distributions more similar to
those from a few-cycle pulse. As such, we have verified that the initial
excited-state population is quickly reduced for the parameter range employed
(not shown). This effect has not been incorporated in the CQSFA. The CQSFA
also has some asymmetry for a different reason. Despite using a monochromatic
field, we must restrict the ionization events to a finite temporal range,
which introduce some artificial asymmetry depending on how this range is
chosen. This can be eliminated by the incoherent averaging procedure described
in Werby et al. (2021).

Figure 5: Photoelectron momentum distributions computed using the one-electron
Schrödinger solver Qprop [panels (a) and (d)] and the Coulomb quantum-orbit
strong-field approximation (CQSFA) [panels (b) and (e)] with all relevant
orbit types for excited helium at the given $I_{p}$ = 0.1238 a.u. In panels
(c) and (f), we have considered the CQSFA and a slightly larger ionization
potential ($I_{p}$ = 0.175 a.u.) as a test case. The upper and lower rows have
been computed for initial bound states aligned parallel and perpendicular to
the driving-field polarization, respectively. For Qprop, we have used a four-
cycle flat-top pulse with an additional half cycle turn on and off. For the
CQSFA, we have approximated the driving field by a monochromatic wave. In all
cases, we considered a driving-field peak intensity $I_{0}=3.2\times 10^{13}$
W/cm2 ($U_{p}=7.65$ eV), and wavelength $\lambda=1600$ nm. The thick circular
dashed line indicates the direct ATI cutoff $2U_{p}$ and the yield in each
panel has been normalized to its maximum value, but there are roughly between
one and two orders of magnitude difference between the parallel and
perpendicular configuration. Furthermore, the yield for the perpendicular case
is suppressed by a factor of $10^{2}$ for CQSFA and $10$ for Qprop,
respectively. The green dashed rectangle on panels (d) and (f) indicate the
suppressed rescattering ridge for a particular perpendicular momentum, whereas
the dashed triangle shows off axis suppression in panels (b), (c), (e) and
(f). Finally, the red squares seen on all panels show the spilling beyond the
$2U_{p}$ cutoff as mentioned above.

Another noteworthy aspect in the PMDs are superimposed circular regions of
radius $10U_{p}$ centered at $(p_{\parallel},p_{\perp})=(\pm
2\sqrt{U_{p}},0)$, which overlap at
$(p_{\parallel},p_{\perp})=(0,2\sqrt{U_{p}})$ along the perpendicular momentum
axis. The boundaries of such regions determine a ridge, which is a well known
indicator of rescattering. The rescattering ridge is present in the CQSFA
computation, regardless of the orientation of the initial orbital [see middle
panels in the figure], but is strongly suppressed for the TDSE in the case of
perpendicular orientation [see lower left panel]. This suppression is not
isotropic, and remnants of the ridge can be seen around the axis $p_{\perp}$
perpendicular to the driving-field polarization. If, on the other hand, we
consider a slightly larger ionization potential for the CQSFA, we also observe
a very faint rescattering ridge for the perpendicular case (see right column
in the figure). In addition to that, there are suppressions of the
photoelectron signal around the $p_{\perp}$ axis, which are stronger for
perpendicular alignment. Subtler features are a spilling of the spider-like
structure beyond the direct $2U_{p}$ cutoff and a wedge-like suppression
(shown on Fig. 5 by the dashed square and triangle, respectively.) along the
$p_{\parallel}$ axis for the perpendicular bound-state configuration in all
approaches.

We will next analyze the above features in terms of CQSFA orbits. With that
aim in mind, in Fig. 6, we plot PMDs obtained using individual CQSFA orbits.
The plots illustrate the momentum region they occupy, and how they are
influenced by the parallel and perpendicular configurations. We will commence
by looking at the parameters in Figs. 5(b) and (e), which, apart from an
overall suppression of roughly two orders of magnitude in the perpendicular
case, look qualitatively similar.

Figure 6: Single-orbit electron momentum distributions for excited helium
($I_{p}=0.1238$ a.u.) using the same field parameters as in Fig. 5 and a
single cycle, for parallel and perpendicular alignment with regard to the
laser-field polarization. From left to right, we present the contributions of
orbit 1 [panels (a) and (b)], orbit 2 [panels (c) and (d)], orbit 3 [panels
(e) and (f)] and orbit 4 [panels (g) and (h)], where the first column is the
parallel configuration and the second is the perpendicular configuration. The
thick circular dashed lines indicate the direct-ATI cutoff energy $2U_{p}$,
and the yield on each panel has been normalized to its maximum value, but
there is an overall suppression of over one order of magnitude for the
perpendicular case. Orbits 3 and 4 are not defined in the grey regions
displayed in panels (e), (f), (g) and (h).

Orbits 1 and 2 lead to approximately elliptical PMDs, whose major axis is
along $p_{\parallel}$ and which are mostly constrained by the $2U_{p}$ direct
ATI cutoff. This is expected as both orbits behave primarily as direct SFA
orbits. Spilling beyond this cutoff occurs mostly along the driving-field
polarization and is caused by the presence of the Coulomb potential. For
perpendicular orientation, there is a wedge-like suppression along the
$p_{\parallel}$ axis for the contributions of orbit 1 [see Fig. 6(b) in
comparison with Fig. 6(a)], whose angular range broadens as $|p_{\parallel}|$
increases. This effect can be understood as follows: if an electron propagates
along orbit 1, its final and initial momentum will be similar, and, for high
enough final momenta, close to the standard SFA Maxwell et al. (2017). Thus,
if the electron is freed along the field-polarization axis, it will continue
to propagate along this direction. For small angles, its propagation direction
will not change significantly. The $2p$ orbital is strongly directional, so
that if it is oriented perpendicular to the laser-field polarization,
tunneling along the $p_{\parallel}$ axis will be hindered. The wedge-like
shape comes from the fact that this is not a perfect mapping, and the final
and initial momenta will only coincide for large values of $p_{\parallel}$.
Still, the above discussion will approximately hold. In contrast, as orbit 2
is essentially a field-dressed hyperbola, the binding potential plays a
stronger role and the transverse momentum mapping is non-trivial. In this
case, the wedge-like feature is absent.

The PMDs associated with orbits 3 and 4, shown in the bottom rows in Fig. 6,
occupy momentum regions within a broader angular range than the previous
single-orbit distributions. This is due to them undergoing a stronger
deflection and interaction with the core. Previous studies have shown that
orbit 3 has a hybrid character Maxwell et al. (2018a); Maxwell and Figueira de
Morisson Faria (2018), with no counterpart in the SFA, while orbit 4 behaves
essentially like a rescattered orbit. The PMD associated with orbit 3 mainly
occupies the momentum region defined by the overlapping rescattering ridges,
with a slight spilling near the polarization axes. In both parallel and
perpendicular configurations, there is a suppression along $p_{\parallel}=0$,
which, however, is more pronounced in the latter case. Orbit 4 leads to a PMD
with clear rescattering ridges for both parallel and perpendicular
configurations and a caustic structure at the boundary. This indicates that,
for the present parameter range, the dynamics associated with orbit 4 are not
predominantly located along the polarization axis, as would be expected in an
SFA treatment of rescattered trajectories. On the other hand, if we consider a
slightly larger ionization potential, there is a much stronger suppression of
the ridge for orbit 4 and of the signal around the $p_{\perp}$ axis for orbit
3. This is illustrated in Fig. 7.

Figure 7: Single-orbit distributions calculated for orbits 3 [first row] and 4
[second row] using the CQSFA for the same field as in Fig. 5 and a single
cycle, but considering a slightly higher ionization potential ($I_{p}=0.175$
a.u.). The first columns panels (a) and (c) have been computed for parallel
case. Whereas perpendicular alignment with regard to the driving-field
polarization, is the second column, panels (b) and (d) and the circular lines
indicate the direct ATI cutoff at energy $2U_{p}$ across all panels. The two
upper panels have been normalized to their highest value, but there is a
difference of one to two orders of magnitude between the parallel and
perpendicular cases. In the lower panels the same scale was used.

These are great examples of how rescattering may no longer be restricted to
the polarization axis due to the Coulomb potential. For the correct value of
the ionization potential associated with the excited state, the CQSFA tunnel
exit is located in a region for which the Coulomb potential is dominant. This
suggests that orbit 4 will behave in a much less directional way than what one
expects within the SFA framework. It will be able to probe orbitals oriented
perpendicular to the driving-field polarization, and the rescattering ridge
will be present. On the other hand, for higher ionization potentials, the
tunnel exit will be located in the region for which the field is dominant.
This will reduce the angular range for which rescattering may occur, and will
cause a strong suppression for photoelectron energies higher than $2U_{p}$.

Figure 8: The behavior of orbit 4 for different ionization potentials
($I_{p}=0.150$ a.u. and $I_{p}=0.175$ a.u.), plotted over one cycle of the
driving field and calculated for the same driving-field parameters as in Fig.
6. The tunnel exit for the orange curve occurs at $5.46$ a.u in parallel
direction, while for the blue curve it is located at $5.89$ a.u. The final
momentum components associated with the orange and blue curves are $1.90$ a.u
in the parallel and $0.53$ a.u. in the perpendicular direction.

This issue is illustrated in Fig. 8, which shows that this transition depends
quite critically on the interplay between the ionization potential and the
ponderomotive energy. For the ionization potential $I_{p}=0.15$ a.u., which is
slightly larger than that of excited helium (orange curve), the binding
potential dominates. This causes orbit 4 to leave in a perpendicular direction
to the field, until approximately 2 a.u. Subsequently, the orbit leaves and
returns to the core region propagating in the expected direction, that is,
parallel to the driving-field polarization. Still, the large initial
transverse momentum acquired by orbit 4 makes it much less directional than
what is predicted by the SFA. This means that, for that parameter range, this
orbit will also probe orbitals with perpendicular alignment. This explains the
similar CQSFA spectra in the previous section, and the prominent rescattering
ridges in Fig. 6(h) and Fig. 5(e). In contrast, if the tunnel exit is located
in a region for which the field is dominant (blue curve), orbit 4 will behave
in a much more directional way, practically along the polarization axis. This
will cause an overall suppression for initial bound states oriented
perpendicular to the driving-field polarization. Nonetheless, there are still
momentum components in the perpendicular direction and the orbit does not
leave parallel to the field. This, together with the presence of intermediate
orbits, leads to remnants of the rescattering ridge in the PMDs.

Therefore, even in a single-electron setting, the residual binding potential
influences rescattering. It gives rise to orbits of hybrid character, and, in
some instances, makes ionization and/or recollision deviate from the driving-
field polarization axis. Hybrid orbits are absent in Born-type approaches such
as the SFA, and will contribute for the anisotropic suppression in the
rescattering ridge reported in this work. Ionization, and in some cases
propagation, perpendicular to the laser-polarization axis may occur and will
partially hinder the suppression of the rescattered ATI signal. One should
note, however, that the CQSFA results for excited helium strongly deviate from
the TDSE computations. Interestingly, the TDSE agrees more with the CQSFA if a
higher (artificial) ionization potential is taken in the CQSFA.

##  IV Conclusions

In this work, we have probed excited helium in initial $2p$-state
configurations parallel and perpendicular to a mid-IR, linearly polarized
field. We use different numerical methods, namely the B-spline algebraic
diagrammatic construction (ADC) Ruberti et al. (2014) and the one-electron
Schrödinger solver Qprop Mosert and Bauer (2016); Tulsky and Bauer (2020),
whose outcomes were compared to those of two semi-analytic methods: the
Strong-Field Approximation (SFA) and the Coulomb quantum-orbit strong-field
approximation (CQSFA). Our main goal was to understand the limitations of the
rescattering model in its simplest form, i.e., that dictated by the SFA, and
what subtleties must be taken into consideration. Overall, it is noteworthy
that different computations may lead to spectra whose shapes and intensities
are quite distinct.

The SFA predicts a sharp decrease in the ATI signal for energies above the
direct ATI cutoff $2U_{p}$ for all cases, and a strong suppression in the
high-energy ATI plateau, which extends up to the energy of 10$U_{p}$, for the
perpendicular initial state. This is due to tunnel ionization and recollision
being located around the field polarization axis. This implies that both
processes are suppressed due to the initial-state geometry. First, there will
not be a substantial probability density to undergo tunneling and, upon the
act of rescattering, there will not be a significant overlap of the returning
wavepacket with the orthogonally oriented bound state. Furthermore, in the SFA
no changes that may occur in the core and alter its geometry are incorporated.
This sharp decrease is similar to that resulting from changing the driving-
field polarization, instead of that of the initial bound state. Any deviations
from those patterns mean that the dynamics are more intricate, either because
ionization or recollision may occur off axis, or because processes involving
the core must be taken into consideration.

The results presented in this paper reveal a subtler picture, with the angle-
integrated ATI spectra exhibiting ramp-like structures in some computations.
This is most striking for the B-spline ADC computation, which goes beyond a
single-active electron and has the least degree of physical approximations.
The B-spline ADC spectra exhibit multiple ramps and a smooth transition from
direct to rescattered ATI, both for parallel and perpendicular orientation. A
ramp-like structure is expected if there are many superimposed processes
and/or substantial depletion. This will tend to weaken effects caused by the
initial bound-state geometry or suppress high-energy features by removing
parts of the bound-state population via several channels. It could also be
caused by mechanisms rendering rescattering less effective.

It is noteworthy that the B-Spline ADC spectra behave in strikingly different
ways from those obtained with the one-electron computations. The spectra from
Qprop and the analytical models exhibit a much more obvious transition from
the direct to the rescattered regimes, while the B Spline ADC spectra present
a continuous ramp. Removal of a specific bound-state excitation channel
($2^{1}P-2^{1}S$) in the B-Spline ADC influences the spectra. However, the
changes are too subtle to justify all the discrepancies. The current results
suggest that in the parallel-aligned case the presence of electron exchange
renders rescattering less effective by introducing repulsive effects, while
for the perpendicular-aligned configuration multielectron effects enhance
tunneling and low-energy rescattering. Because exchange is non-local, it
renders both effects less directional, which decreases the suppression in the
perpendicular case. Interestingly, a suppression of one order of magnitude
between parallel and perpendicular polarised initial $p$ states has also been
observed for ATI spectra computed for argon using the time-dependent density
functional theory (TDDFT) Tong et al. (2015). At the ADC(1) level at which our
simulations were done, we are close in accuracy to the TDDFT, and both
computations incorporates electron exchange. This suggests that rescattering
in ATI is a delicate effect that probes the fine details of the core structure
and core-photoelectron interaction including its non-local (exchange)
component. The influence of electron exchange will also be energy dependent.
The closer to the core the rescattering electron gets, the higher the
photoelectron energy will be. This also means that exchange effects will
become increasingly important, leading to a ramp.

Still, even for single-electron computations, the spectra exhibit quite
different features, depending on what is included or left out. For instance,
in the parallel-aligned case, the SFA spectrum shows a much stronger drop in
signal for the plateau than that of Qprop. In the perpendicular-aligned
scenario, the Qprop outcome follows the SFA spectrum from above in the direct
region, and decays sharply after the direct ATI cutoff. In contrast, in the
CQSFA the plateau is overestimated unless an artificially large ionization
potential is taken. Moreover, it varies much less dramatically with the
initial states’ orientation than all other methods considered in this work.

The features mentioned above can be understood in greater depth by looking at
angle-resolved photoelectron momentum distributions (PMDs) and by analyzing
the electron orbits involved. This analysis has been performed for the CQSFA
and Qprop. The agreement is good for the main features, such as the
rescattering ridge and a suppression near the polarization axis that occurs
for the perpendicular configuration. However, the holographic patterns are
much more irregular for Qprop than for the CQSFA. Furthermore, in the
perpendicular configuration the rescattered ATI signal is much more suppressed
for Qprop. This includes most of the signal beyond the energy of $2U_{p}$ and
the rescattering ridge. A better agreement is obtained if the ionization
potential is artificially increased in the CQSFA computation.

The irregularities in the Qprop holographic patterns are explained by the huge
amount of bound-state depletion, which favours certain events over others and
compromises the contrast in the quantum interference patterns. Furthermore,
there may be other mechanisms through which the electron is freed in the
continuum, such as over-the-barrier ionization and coupling with highly
excited states. These features are not incorporated in the CQSFA. The over-
enhancement of the rescattered contributions in the perpendicular case is due
to the Coulomb potential being dominant upon ionization, for the CQSFA orbit
type that leads to the rescattering ridge. This introduces a large
perpendicular momentum component in the orbits, so that they initially tunnel
in a perpendicular direction and may return farther away from the polarization
axis.

Further studies indicate that a tunnel exit in a region for which the field
dominates, which can be achieved with a slightly larger ionization potential,
will alter this behavior. In this case, the orbit will be more localized
around the polarization axis, so that tunneling and rescattering will be
suppressed for perpendicular configurations. Interestingly, this suppression
is closer to the behavior observed for the TDSE computation. The TDSE outcome
will be influenced by ac Stark shifts, additional ionization channels, bound-
state broadening and over-the-barrier ionization, which are not included in
the CQSFA. Support for this assumption is provided by further TDSE
computations using weaker binding energies, which exhibit a rescattering ridge
for the perpendicular-aligned initial states (not shown). Some combination of
these effects likely leads to the rescattering behaviour observed in the CQSFA
when a larger ionization potential is used. Discussions of how to include
Stark shifts in the SFA are provided in Dimitrovski et al. (2010); Etches and
Madsen (2010). For a proposal of how to measure Stark shifts in excited
helium, see He et al. (2011); Chini et al. (2012).

If, on the one hand, this suggests that the CQSFA may need modifications for
highly excited, loosely bound states, on the other hand this is a nice
illustration of how the interplay of the binding potential and the external
laser field may lead to counterintuitive behaviors and move rescattering or
ionization away from the driving-field direction. This, together with the core
dynamics effects observed for the B-spline ADC spectra, shows that there may
be other possibilities for imaging than those based exclusively on the field
exploring bound-state geometry or driving-field shapes. Finally, one should
note that preparing targets in excited and oriented states is within
experimental reach. For instance, one may use the Magneto-optical trap recoil
ion momentum spectroscopy (MOTRIMS) technique (see, e.g. Blieck et al. (2008);
Fischer et al. (2012) in the context of electron collisions and Yuan et al.
(2020) for multiply ionized Rubidium in a strong elliptically polarized
field). Moreover, one may prepare the 2 1P state of helium either by resonant
excitation from the metastable 2 1S state Williams and Fry (1968); Fry and
Williams (1969), or by a direct extraction from the ground state Larsson et
al. (1995); Gisselbrecht et al. (1999). The polarization of the excited state
can be controlled as described in Johansson et al. (2003), and the lifetime of
the 21P state is several orders of magnitude larger than the timescales
involved in the present work. Usual excitation results in a superposition of
the ground and excited state helium. However, for the parameters employed in
this work, the ground state is too strongly bound to result in any appreciable
ATI. Therefore, we are probing the excited-state population and whatever
population stays in the ground state is simply lost for us.

##  V Acknowledgements

We would like to thank Dr Bridgette Cooper for her insights and tutorials in
using GAMESS and Dr Emilio Pisanty for his knowledge on colour schemes to
enhance plots without the loss of key features. ASM acknowledges grant
EP/P510270/1 and CFMF grant No. EP/J019143/1, both funded by the UK
Engineering and Physical Sciences Research Council (EPSRC). ASM and MC
acknowledge support from ERC AdG NOQIA, Agencia Estatal de Investigación
(“Severo Ochoa” Center of Excellence CEX2019-000910-S, Plan National FIDEUA
PID2019-106901GB-I00/10.13039 / 501100011033, FPI), Fundació Privada Cellex,
Fundació Mir-Puig, and from Generalitat de Catalunya (AGAUR Grant No. 2017 SGR
1341, CERCA program, QuantumCAT _U16-011424 , co-funded by ERDF Operational
Program of Catalonia 2014-2020), MINECO-EU QUANTERA MAQS (funded by State
Research Agency (AEI) PCI2019-111828-2 / 10.13039/501100011033), EU Horizon
2020 FET-OPEN OPTOLogic (Grant No 899794), and the National Science Centre,
Poland-Symfonia Grant No. 2016/20/W/ST4/00314, Marie Sklodowska-Curie grant
STRETCH No. 101029393, “La Caixa” Junior Leaders fellowships (ID100010434),
and EU Horizon 2020 under Marie Skłodowska-Curie grant agreement No 847648
(LCF/BQ/PI19/11690013, LCF/BQ/PI20/11760031, LCF/BQ/PR20/11770012). MR and VA
acknowledge support from EPSRC/DSTL via MURI award No. EP/N018680/1.

## Appendix. SFA Matrix elements

In this Appendix, we provide the main steps for calculating the integrals

$\mathcal{I}_{1}=\int d\mathbf{r}e^{i\mathbf{k}\cdot\mathbf{r}}V(\mathbf{r})$
(42)

and

$\mathcal{I}_{2}=\int
d\mathbf{r}e^{-i\mathbf{k}\cdot\mathbf{r}}r\cos\theta\psi_{0}(\mathbf{r})$
(43)

in Eqs. (17) and (19), for the effective binding potential (7). To avoid
singularities in the Fourier transform of this potential, it is multiplied
with a damping factor: $V(r)\rightarrow V(r)\times e^{-\epsilon r}$. Numerical
simulations show that this factor does not distort the shape of the ATI
spectrum, although it can shift its overall magnitude; for this study we have
chosen $\epsilon=1.0$. For other methods to treat this singularity see, e.g.,
Ref. de Morisson Faria and Lewenstein (2005).

Using the spatial representation

$\Psi_{0}(\mathbf{r})=R_{n_{0}l_{0}}(r)Y_{l_{0}m_{0}}(\Omega)$

of the initial state, where $\Omega$ is the solid angle, Eqs. (42) and (43)
can be simplified as

$\displaystyle\mathcal{I}_{1}$ $\displaystyle=$
$\displaystyle-4\pi\left[\frac{1}{\epsilon^{2}+k^{2}}+\frac{a_{1}}{\epsilon^{2}+a_{2}^{2}+k^{2}}+\frac{2a_{3}a_{4}}{(\epsilon^{2}+a_{4}^{2}+k^{2})^{2}}+\frac{a_{5}}{\epsilon^{2}+a_{6}^{2}+k^{2}}\right]$
(44) $\displaystyle\mathcal{I}_{2}$ $\displaystyle=$
$\displaystyle\frac{(4\pi)^{3/2}}{\sqrt{3}}\sum_{lm}i^{-l}Y_{lm}(\Omega_{k})\int
drr^{3}R_{n_{0}l_{0}}(r)j_{l}(kr)\int d\Omega
Y_{lm}^{*}(\Omega)Y_{10}(\Omega)Y_{l_{0}m_{0}}(\Omega,)$ (45)

where we have inserted the damping factor into Eq. (7), substituted
$\cos\theta$ by a spherical harmonic and used the multipole expansion for the
exponential function in three dimensions. In Eq. (45), the solid angle
$\Omega$ is defined in coordinate space and the solid angle $\Omega_{k}$ with
regard to the intermediate momentum $\mathbf{k}$.

The spatial radial integral can be solved analytically for an initial state
described by a hydrogenic wave function. For the ground state (1 _s_) we use

$R_{10}(r)=2Z^{3/2}e^{-Zr},$

for the excited state (2 _s_) we have

$R_{20}(r)=\frac{Z^{3/2}}{\sqrt{2}}(1-Zr/2)e^{-Zr/2},$

and for the polarizable excited state (2 _p_) we use

$R_{21}(r)=\frac{Z^{5/2}}{2\sqrt{6}}re^{-Zr/2},$

where $Z$ is the nuclear charge.

For the ground state 1 _s_ , where $l_{0}=0$, the sum in Eq. (45) is
restricted to $l=1$ and the spatial integral reads as

$2Z^{3/2}\int
drr^{3}e^{-Zr}\left(\frac{\sin(kr)}{(kr)^{2}}-\frac{\cos(kr)}{kr}\right)=\frac{16Z^{5/2}k}{\left(Z^{2}+k^{2}\right)^{3}}.$

Similarly, for 2 _s_ we obtain

$\frac{Z^{3/2}}{\sqrt{2}}\int
drr^{3}e^{-Zr/2}(1-Zr/2)\left(\frac{\sin(kr)}{(kr)^{2}}-\frac{\cos(kr)}{kr}\right)=\frac{2048kZ^{5/2}(2k^{2}-Z^{2})}{\sqrt{2}\left(Z^{2}+4k^{2}\right)^{4}}.$

For the excited state 2 _p_ , where $l_{0}=1$, the sum over $l$ runs from 0 to
2, giving the following three expressions in ascending order of $l$

$\frac{Z^{5/2}}{2\sqrt{6}}\int
drr^{3}e^{-Zr/2}\frac{\sin(kr)}{k}=\frac{Z^{5/2}}{2\sqrt{6}}\frac{12Zk(Z/2-k)(Z/2+k)}{k\left(Z^{2}/4+k^{2}\right)^{4}}$
$\frac{Z^{5/2}}{2\sqrt{6}}\int
drr^{3}e^{-Zr/2}\left(\frac{\sin(kr)}{k^{2}r}-\frac{\cos(kr)}{k}\right)=\frac{Z^{5/2}}{2\sqrt{6}}\left(\frac{2(3(Z/2)^{2}-k^{2})}{k\left((Z/2)^{2}+k^{2}\right)^{3}}-\frac{6((Z/2)^{4}-6(Z/2)^{2}k^{2}+k^{4})}{k\left((Z/2)^{2}+k^{2}\right)^{4}}\right)$

$\frac{Z^{5/2}}{2\sqrt{6}}\int
drr^{3}e^{-Zr/2}\left(\frac{3\sin(kr)}{k^{3}r^{2}}-\frac{\sin(kr)}{k}-\frac{3\cos(kr)}{k^{2}r}\right)=\frac{Z^{5/2}}{2\sqrt{6}}\left[\mathcal{A}_{1}(Z,k)-\mathcal{A}_{2}(Z,k)-\mathcal{A}_{3}(Z,k)\right],$

with

$\mathcal{A}_{1}(Z,k)=\frac{3Z}{k^{2}\left((Z/2)^{2}+k^{2}\right)^{2}}$
$\mathcal{A}_{2}(Z,k)=\frac{12Zk(Z/2-k)(Z/2+k)}{k\left((Z/2)^{2}+k^{2}\right)^{4}}$

and

$\mathcal{A}_{3}(Z,k)=\frac{3Z((Z/2)^{2}-3k^{2})}{k^{2}\left((Z/2)^{2}+k^{2}\right)^{3}}.$

The angular integral of the product of three spherical harmonics is given by
the product of Clebsch-Gordan coefficients ($CG$), whose general form reads as

$\int d\Omega
Y_{l_{1}m_{1}}^{*}(\Omega)Y_{l_{2}m_{2}}(\Omega)Y_{l_{3}m_{3}}(\Omega)=\sqrt{\frac{(2l_{2}+1)(2l_{3}+1)}{4\pi(2l_{1}+1)}}CG(l_{2},l_{3},l_{1};m_{2},m_{3},m_{1})CG(l_{2},l_{3},l_{1};0,0,0).$

If the initial state is perpendicularly polarized, then the angular part
$Y_{l_{0}m_{0}}(\Omega)$ in Eq. (45) will have some combination of
$m_{0}=1,-1$. Therefore, the only nonzero contribution to the rescattering
term will come from the element in the summation for which $m=1,-1$. One
should note that, so far, the rescattered ATI transition amplitude (19)
incorporates all possible paths for this three-step process are integrated
over, all possible ionization times for the ionization (when the laser is
active) and rescattering processes, and all the possible intermediate Volkov
states with canonical momentum $\mathbf{k}$. The integrals in Eq. (19) can be
calculated using saddle-point methods, which simplify the numerical effort
involved and are easily relatable to an intuitive, orbit-based picture
Salières et al. (2001); Nayak et al. (2019); Figueira de Morisson Faria et al.
(2002). However, depending on the problem at hand care must be taken. If the
saddle approximation is used in all three dimensions, then $k_{x}=0$ and
$k_{y}=0$, this constrains the vector $\mathbf{k}$ to be fixed in the
z-direction, leaving no contribution at all and not enabling meaningful
comparisons to the numerical simulation results to be made. Physically, this
implies that only those rescattering processes that occur exactly along the
laser polarization axis are included in this approximation, and the remaining
processes are selected out. For discussions in the context of molecules see
Chirilă and Lein (2009).

We therefore use the saddle point approximation for just two dimensions in the
integral and numerically integrate over the dimension of the atomic
polarization (chosen in the $x$ direction). The action can be expanded in
three dimensions, so that

$S(\mathbf{k};t,t^{\prime})=\frac{1}{2}\intop_{t^{\prime}}^{t}dt^{\prime\prime}\left[(k_{z}+A_{z}(t^{\prime\prime}))^{2}+k_{x}^{2}+k_{y}^{2}\right]$

Each variable in the integrand that we approximate is substituted by its
stationary value obtained with the saddle-point approximation, obtained
setting the derivative of the action is equal to zero. This gives

$\displaystyle k_{z}^{st}$
$\displaystyle=-\frac{1}{t-t^{\prime}}\intop_{t^{\prime}}^{t}dt^{\prime\prime}A_{z}(t^{\prime\prime})$
(46) $\displaystyle k_{y}^{st}$ $\displaystyle=0.$ (47)

for $\partial S(\mathbf{k};t,t^{\prime})/\partial k_{z}=0$ and $\partial
S(\mathbf{k};t,t^{\prime})/\partial k_{y}=0$, respectively. For each dimension
in which we use the approximation, we make the following substitution

$\intop d\textbf{k}\left[\frac{2\pi}{\alpha+i(t-t^{\prime})}\right]^{1/2},$
(48)

where $\alpha$ is a small parameter to avoid a singularity when $t\rightarrow
t^{\prime}$. If the saddle point approximation is used in three dimensions,
then this parameter is not necessary unless the initial state has _p_
symmetry. This is because the exponent in the last integral of Eq. (19) goes
to zero in the limit when $t\rightarrow t^{\prime}$ and the spherical Bessel
function of such an argument tends to zero for all orders apart from for
$j_{0}$ \- this function only appears in the summation for an initial _p_
state. Thus, modifying (45) with (46), (47) and (48), we obtain

$\displaystyle M_{\mathrm{Resc}}(\mathbf{p})$
$\displaystyle=-\intop_{-\infty}^{\infty}dt\intop_{-\infty}^{t}dt^{\prime}\left[\frac{2\pi}{\alpha+i(t-t^{\prime})}\right]E(t^{\prime})e^{iI_{P}t^{\prime}}e^{-iS(k_{z}^{sp};t,t^{\prime})}e^{iS(\mathbf{p},t)}$
$\displaystyle\int
dk_{x}e^{-iS(k_{x};t,t^{\prime})}\times\frac{1}{(2\pi)^{3}}\int
d\mathbf{r^{\prime}}e^{i(\mathbf{k}-\mathbf{p}).\mathbf{r^{\prime}}}V(\mathbf{r^{\prime}})\times\frac{1}{(2\pi)^{3/2}}\int
d\mathbf{r}e^{-i(\mathbf{k}+\mathbf{A(t^{\prime})}).\mathbf{r}}\cos\theta\psi_{0}(\mathbf{r})$

where $\mathbf{k}=k_{z}^{sp}\hat{z}+k_{y}^{sp}\hat{y}+k_{x}\hat{x}$. This
specific approach is being used both for parallel and perpendicular bound-
state orientation. This integral was also performed in Milošević (2013b).

## References

  * Agostini et al. (1979) P. Agostini, F. Fabre, G. Mainfray, G. Petite, and N. K. Rahman, Phys. Rev. Lett. 42, 1127 (1979), URL https://link.aps.org/doi/10.1103/PhysRevLett.42.1127.
  * Zuo et al. (1996) T. Zuo, A. Bandrauk, and P. Corkum, Chemical Physics Letters 259, 313 (1996), ISSN 0009-2614, URL http://www.sciencedirect.com/science/article/pii/0009261496007865.
  * Huismans et al. (2011) Y. Huismans, A. Rouzée, A. Gijsbertsen, J. H. Jungmann, A. S. Smolkowska, P. S. W. M. Logman, F. Lépine, C. Cauchy, S. Zamith, T. Marchenko, et al., Science 331, 61 (2011), URL http://science.sciencemag.org/content/331/6013/61.
  * de Morisson Faria and Maxwell (2020) C. F. de Morisson Faria and A. S. Maxwell, Reports on Progress in Physics 83, 034401 (2020), URL https://doi.org/10.1088%2F1361-6633%2Fab5c91.
  * Lewenstein et al. (1995) M. Lewenstein, K. C. Kulander, K. J. Schafer, and P. H. Bucksbaum, Phys. Rev. A 51, 1495 (1995), URL https://link.aps.org/doi/10.1103/PhysRevA.51.1495.
  * Lohr et al. (1997) A. Lohr, M. Kleber, R. Kopold, and W. Becker, Phys. Rev. A 55, R4003 (1997), URL https://link.aps.org/doi/10.1103/PhysRevA.55.R4003.
  * Milošević and Ehlotzky (1998) D. B. Milošević and F. Ehlotzky, Phys. Rev. A 58, 3124 (1998), URL https://link.aps.org/doi/10.1103/PhysRevA.58.3124.
  * Becker et al. (2018) W. Becker, S. P. Goreslavski, D. B. Milošević, and G. G. Paulus, Journal of Physics B: Atomic, Molecular and Optical Physics 51, 162002 (2018), URL https://doi.org/10.1088%2F1361-6455%2Faad150.
  * Keldysh (1965) L. Keldysh, Sov. Phys. JETP 20, 1307 (1965).
  * Faisal (1973) F. H. M. Faisal, Journal of Physics B: Atomic and Molecular Physics 6, L89 (1973), URL https://doi.or/10.1088/0022-3700/6/4/011.
  * Reiss (1980) H. R. Reiss, Phys. Rev. A 22, 1786 (1980), URL https://link.aps.org/doi/10.1103/PhysRevA.22.1786.
  * Becker et al. (2002) W. Becker, F. Grasbon, R. Kopold, D. Milošević, G. Paulus, and H. Walther (Academic Press, 2002), vol. 48 of _Advances In Atomic, Molecular, and Optical Physics_ , pp. 35–98, URL https://www.sciencedirect.com/science/article/pii/S1049250X02800064.
  * Kopold et al. (2000) R. Kopold, D. B. Milošević, and W. Becker, Phys. Rev. Lett. 84, 3831 (2000), URL https://link.aps.org/doi/10.1103/PhysRevLett.84.3831.
  * Salières et al. (2001) P. Salières, B. Carré, L. Le Déroff, F. Grasbon, G. G. Paulus, H. Walther, R. Kopold, W. Becker, D. B. Milošević, A. Sanpera, et al., Science 292, 902 (2001), ISSN 0036-8075, eprint https://science.sciencemag.org/content/292/5518/902.full.pdf, URL https://science.sciencemag.org/content/292/5518/902.
  * Paulus et al. (1994) G. G. Paulus, W. Becker, W. Nicklich, and H. Walther, Journal of Physics B: Atomic, Molecular and Optical Physics 27, L703 (1994), URL https://doi.org/10.1088%2F0953-4075%2F27%2F21%2F003.
  * de Morisson Faria and Liu (2011) C. F. de Morisson Faria and X. Liu, Journal of Modern Optics 58, 1076 (2011), eprint https://doi.org/10.1080/09500340.2010.543958, URL https://doi.org/10.1080/09500340.2010.543958.
  * Becker et al. (2012) W. Becker, X. Liu, P. J. Ho, and J. H. Eberly, Rev. Mod. Phys. 84, 1011 (2012), URL https://link.aps.org/doi/10.1103/RevModPhys.84.1011.
  * Freeman et al. (1987) R. R. Freeman, P. H. Bucksbaum, H. Milchberg, S. Darack, D. Schumacher, and M. E. Geusic, Phys. Rev. Lett. 59, 1092 (1987), URL https://link.aps.org/doi/10.1103/PhysRevLett.59.1092.
  * Morishita et al. (2007) T. Morishita, Z. Chen, S. Watanabe, and C. D. Lin, Phys. Rev. A 75, 023407 (2007), URL https://link.aps.org/doi/10.1103/PhysRevA.75.023407.
  * Potvliege and Vučić (2006) R. M. Potvliege and S. Vučić, Phys. Rev. A 74, 023412 (2006), URL https://link.aps.org/doi/10.1103/PhysRevA.74.023412.
  * Kopold et al. (2002) R. Kopold, W. Becker, M. Kleber, and G. G. Paulus, Journal of Physics B: Atomic, Molecular and Optical Physics 35, 217 (2002), URL https://doi.org/10.1088%2F0953-4075%2F35%2F2%2F302.
  * Milošević et al. (2007) D. B. Milošević, E. Hasović, M. Busuladžić, A. Gazibegović-Busuladžić, and W. Becker, Phys. Rev. A 76, 053410 (2007), URL https://link.aps.org/doi/10.1103/PhysRevA.76.053410.
  * Muller (1999) H. G. Muller, Phys. Rev. Lett. 83, 3158 (1999), URL https://link.aps.org/doi/10.1103/PhysRevLett.83.3158.
  * Smirnova et al. (2009) O. Smirnova, Y. Mairesse, S. Patchkovskii, N. Dudovich, D. Villeneuve, P. Corkum, and M. Y. Ivanov, Nature 460, 972 (2009), ISSN 1476-4687, URL https://doi.org/10.1038/nature08253.
  * Mairesse et al. (2010) Y. Mairesse, J. Higuet, N. Dudovich, D. Shafir, B. Fabre, E. Mével, E. Constant, S. Patchkovskii, Z. Walters, M. Y. Ivanov, et al., Phys. Rev. Lett. 104, 213601 (2010), URL https://link.aps.org/doi/10.1103/PhysRevLett.104.213601.
  * Lépine et al. (2014) F. Lépine, M. Y. Ivanov, and M. J. J. Vrakking, Nature Photonics 8, 195 (2014), URL https://doi.org/10.1038/nphoton.2014.25.
  * Calegari et al. (2016) F. Calegari, A. Trabattoni, A. Palacios, D. Ayuso, M. C. Castrovilli, J. B. Greenwood, P. Decleva, F. Martín, and M. Nisoli, Journal of Physics B: Atomic, Molecular and Optical Physics 49, 142001 (2016), URL http://stacks.iop.org/0953-4075/49/i=14/a=142001.
  * Kuleff et al. (2016) A. I. Kuleff, N. V. Kryzhevoi, M. Pernpointner, and L. S. Cederbaum, Phys. Rev. Lett. 117, 093002 (2016), URL https://link.aps.org/doi/10.1103/PhysRevLett.117.093002.
  * Wolter et al. (2016) B. Wolter, M. G. Pullen, A.-T. Le, M. Baudisch, K. Doblhoff-Dier, A. Senftleben, M. Hemmer, C. D. Schröter, J. Ullrich, R. Moshammer, et al., Science 6310, 308 (2016).
  * Karamatskos et al. (2019) E. T. Karamatskos, G. Goldsztejn, S. Raabe, P. Stammer, T. Mullins, A. Trabattoni, R. R. Johansen, H. Stapelfeldt, S. Trippel, M. J. J. Vrakking, et al., The Journal of Chemical Physics 150, 244301 (2019), eprint https://doi.org/10.1063/1.5093959, URL https://doi.org/10.1063/1.5093959.
  * Shiner et al. (2011) A. D. Shiner, B. E. Schmidt, C. Trallero-Herrero, H. J. Wörner, S. Patchkovskii, P. B. Corkum, J.-C. Kieffer, F. Légaré, and D. M. Villeneuve, Nature Physics 7, 464 (2011), ISSN 1745-2481, URL https://doi.org/10.1038/nphys1940.
  * Faccialà et al. (2016) D. Faccialà, S. Pabst, B. D. Bruner, A. G. Ciriolo, S. De Silvestri, M. Devetta, M. Negro, H. Soifer, S. Stagira, N. Dudovich, et al., Phys. Rev. Lett. 117, 093902 (2016), URL https://link.aps.org/doi/10.1103/PhysRevLett.117.093902.
  * Maxwell and Figueira de Morisson Faria (2018) A. S. Maxwell and C. Figueira de Morisson Faria, J. Phys. B At. Mol. Phys. 51, 124001 (2018), eprint 1802.00789, URL http://arxiv.org/abs/1802.00789.
  * Yan et al. (2010) T.-M. Yan, S. V. Popruzhenko, M. J. J. Vrakking, and D. Bauer, Phys. Rev. Lett. 105, 253002 (2010), URL https://link.aps.org/doi/10.1103/PhysRevLett.105.253002.
  * Lai et al. (2015) X. Y. Lai, C. Poli, H. Schomerus, and C. Figueira De Morisson Faria, Phys. Rev. A 92, 1 (2015), ISSN 10941622.
  * Budil et al. (1993) K. S. Budil, P. Salières, A. L’Huillier, T. Ditmire, and M. D. Perry, Phys. Rev. A 48, R3437 (1993), URL https://link.aps.org/doi/10.1103/PhysRevA.48.R3437.
  * Dietrich et al. (1994) P. Dietrich, N. H. Burnett, M. Ivanov, and P. B. Corkum, Phys. Rev. A 50, R3585 (1994), URL https://link.aps.org/doi/10.1103/PhysRevA.50.R3585.
  * Corkum (1993) P. B. Corkum, Phys. Rev. Lett. 71, 1994 (1993), URL https://link.aps.org/doi/10.1103/PhysRevLett.71.1994.
  * Averbukh (2004) V. Averbukh, Phys. Rev. A 69, 043406 (2004), URL https://link.aps.org/doi/10.1103/PhysRevA.69.043406.
  * Mauger et al. (2010) F. Mauger, C. Chandre, and T. Uzer, Phys. Rev. Lett. 105, 083002 (2010), URL https://link.aps.org/doi/10.1103/PhysRevLett.105.083002.
  * Mosert and Bauer (2016) V. Mosert and D. Bauer, Comp. Phys. Commun. 207, 452 (2016).
  * Tulsky and Bauer (2020) V. Tulsky and D. Bauer, Comp. Phys. Commun. 251, 107098 (2020).
  * Ruberti et al. (2014) M. Ruberti, V. Averbukh, and P. Decleva, The Journal of Chemical Physics 141, 164126 (2014), eprint https://doi.org/10.1063/1.4900444, URL https://doi.org/10.1063/1.4900444.
  * Ruberti et al. (2018a) M. Ruberti, P. Decleva, and V. Averbukh, Phys. Chem. Chem. Phys. 20, 8311 (2018a), URL http://dx.doi.org/10.1039/C7CP07849H.
  * Ruberti et al. (2018b) M. Ruberti, P. Decleva, and V. Averbukh, Journal of Chemical Theory and Computation 14, 4991 (2018b), ISSN 1549-9618, URL https://doi.org/10.1021/acs.jctc.8b00479.
  * Maxwell et al. (2017) A. S. Maxwell, A. Al-Jawahiry, T. Das, and C. F. d. M. Faria, Phys. Rev. A 96, 023420 (2017), eprint 1705.01518.
  * Maxwell et al. (2018a) A. S. Maxwell, A. Al-Jawahiry, X. Y. Lai, and C. Figueira de Morisson Faria, J. Phys. B At. Mol. Opt. Phys. 51, 044004 (2018a).
  * Maxwell et al. (2018b) A. S. Maxwell, S. V. Popruzhenko, and C. F. d. M. Faria, Phys. Rev. A 98, 063423 (2018b), URL https://link.aps.org/doi/10.1103/PhysRevA.98.063423.
  * Simpson et al. (2016) E. R. Simpson, A. Sanchez-Gonzalez, D. R. Austin, Z. Diveki, S. E. E. Hutchinson, T. Siegel, M. Ruberti, V. Averbukh, L. Miseikis, C. S. Strüber, et al., New Journal of Physics 18, 083032 (2016), URL https://doi.org/10.1088/1367-2630/18/8/083032.
  * You et al. (2019) D. You, K. Ueda, M. Ruberti, K. L. Ishikawa, P. A. Carpeggiani, T. Csizmadia, L. G. Oldal, H. N. G, G. Sansone, P. K. Maroju, et al., New Journal of Physics 21, 113036 (2019), URL https://doi.org/10.1088/1367-2630/ab520d.
  * Ruberti (2019a) M. Ruberti, Journal of Chemical Theory and Computation 15, 3635 (2019a), ISSN 1549-9618, URL https://doi.org/10.1021/acs.jctc.9b00288.
  * Ruberti (2019b) M. Ruberti, Phys. Chem. Chem. Phys. 21, 17584 (2019b), URL http://dx.doi.org/10.1039/C9CP03074C.
  * Maxwell et al. (2020) A. S. Maxwell, C. F. d. M. Faria, X. Lai, R. Sun, and X. Liu, Phys. Rev. A 102, 033111 (2020), URL https://link.aps.org/doi/10.1103/PhysRevA.102.033111.
  * Tong and Lin (2005) X. M. Tong and C. D. Lin, Journal of Physics B: Atomic, Molecular and Optical Physics 38, 2593 (2005), URL https://doi.org/10.1088%2F0953-4075%2F38%2F15%2F001.
  * Tong and Chu (1997) X.-M. Tong and S.-I. Chu, Phys. Rev. A 55, 3406 (1997), URL https://link.aps.org/doi/10.1103/PhysRevA.55.3406.
  * Bauer et al. (2005) D. Bauer, D. B. Milo šević, and W. Becker, Phys. Rev. A 72, 023415 (2005), URL https://link.aps.org/doi/10.1103/PhysRevA.72.023415.
  * Amini et al. (2019) K. Amini, J. Biegert, F. Calegari, A. Chacón, M. F. Ciappina, A. Dauphin, D. K. Efimov, C. F. de Morisson Faria, K. Giergiel, P. Gniewek, et al., Reports on Progress in Physics 82, 116001 (2019), URL https://doi.org/10.1088/1361-6633/ab2bb1.
  * Figueira de Morisson Faria et al. (2002) C. Figueira de Morisson Faria, H. Schomerus, and W. Becker, Phys. Rev. A 66, 043413 (2002), URL https://link.aps.org/doi/10.1103/PhysRevA.66.043413.
  * Ivanov et al. (2005) M. Y. Ivanov, M. Spanner, and O. Smirnova, Journal of Modern Optics 52, 165 (2005), eprint https://doi.org/10.1080/0950034042000275360, URL https://doi.org/10.1080/0950034042000275360.
  * Popruzhenko (2014) S. V. Popruzhenko, Journal of Physics B: Atomic, Molecular and Optical Physics 47, 204001 (2014), URL https://doi.org/10.1088/0953-4075/47/20/204001.
  * Hasović et al. (2007) E. Hasović, M. Busuladžić, A. Gazibegović-Busuladžić, D. B. Milošević, and W. Becker, Laser Physics 17, 376 (2007), ISSN 1555-6611, URL https://doi.org/10.1134/S1054660X07040135.
  * Shaaran et al. (2010) T. Shaaran, M. T. Nygren, and C. Figueira de Morisson Faria, Phys. Rev. A 81, 063413 (2010), URL https://link.aps.org/doi/10.1103/PhysRevA.81.063413.
  * Shaaran and de Morisson Faria (2010) T. Shaaran and C. F. de Morisson Faria, Journal of Modern Optics 57, 984 (2010), eprint https://doi.org/10.1080/09500340903414619, URL https://doi.org/10.1080/09500340903414619.
  * Wang et al. (2012) B. Wang, Y. Guo, J. Chen, Z.-C. Yan, and P. Fu, Phys. Rev. A 85, 023402 (2012), URL https://link.aps.org/doi/10.1103/PhysRevA.85.023402.
  * Hao et al. (2014) X. Hao, J. Chen, W. Li, B. Wang, X. Wang, and W. Becker, Phys. Rev. Lett. 112, 073002 (2014), URL https://link.aps.org/doi/10.1103/PhysRevLett.112.073002.
  * Maxwell and Faria (2015) A. S. Maxwell and C. F. d. M. Faria, Phys. Rev. A 92, 023421 (2015), URL https://link.aps.org/doi/10.1103/PhysRevA.92.023421.
  * Maxwell and Faria (2016) A. S. Maxwell and C. F. d. M. Faria, Phys. Rev. Lett. 116, 143001 (2016), URL https://link.aps.org/doi/10.1103/PhysRevLett.116.143001.
  * Figueira de Morisson Faria and Augstein (2010) C. Figueira de Morisson Faria and B. B. Augstein, Phys. Rev. A 81, 043409 (2010), URL https://link.aps.org/doi/10.1103/PhysRevA.81.043409.
  * Ivanov and Corkum (1993) M. Y. Ivanov and P. B. Corkum, Phys. Rev. A 48, 580 (1993), URL https://link.aps.org/doi/10.1103/PhysRevA.48.580.
  * Kleinert (2009) H. Kleinert, _Path Integrals in Quantum Mechanics, Statistics, Polymer Physics, and Financial Markets_ (World Scientific, 2009).
  * Milošević (2013a) D. B. Milošević, Journal of Mathematical Physics 54, 042101 (2013a), URL https://doi.org/10.1063/1.4797476.
  * Popruzhenko and Bauer (2008) S. V. Popruzhenko and D. Bauer, J. Mod. Opt. 55, 2573 (2008), ISSN 09500340, eprint 0803.1972.
  * Yan and Bauer (2012) T.-M. Yan and D. Bauer, Phys. Rev. A 86, 053403 (2012), ISSN 1050-2947, URL http://link.aps.org/doi/10.1103/PhysRevA.86.053403.
  * Torlina and Smirnova (2012) L. Torlina and O. Smirnova, Phys. Rev. A 86, 043408 (2012), URL https://link.aps.org/doi/10.1103/PhysRevA.86.043408.
  * Torlina et al. (2013) L. Torlina, J. Kaushal, and O. Smirnova, Phys. Rev. A 88, 053403 (2013), URL https://link.aps.org/doi/10.1103/PhysRevA.88.053403.
  * Shvetsov-Shilovski et al. (2016) N. I. Shvetsov-Shilovski, M. Lein, L. B. Madsen, E. Räsänen, C. Lemell, J. Burgdörfer, D. G. Arbó, and K. Tőkési, Phys. Rev. A 94, 013415 (2016), URL https://link.aps.org/doi/10.1103/PhysRevA.94.013415.
  * Pisanty and Ivanov (2016) E. Pisanty and M. Ivanov, Phys. Rev. A 93, 1 (2016), ISSN 10941622, eprint 1507.00011, URL https://journals.aps.org/pra/pdf/10.1103/PhysRevA.93.043408.
  * Zhou et al. (2016) Y. Zhou, O. I. Tolstikhin, and T. Morishita, Phys. Rev. Lett. 116, 173001 (2016), URL https://link.aps.org/doi/10.1103/PhysRevLett.116.173001.
  * Werby et al. (2021) N. Werby, A. S. Maxwell, R. Forbes, P. H. Bucksbaum, and C. F. d. M. Faria, Phys. Rev. A 104, 013109 (2021), URL https://link.aps.org/doi/10.1103/PhysRevA.104.013109.
  * Tong et al. (2015) X.-M. Tong, G. Wachter, S. A. Sato, C. Lemell, K. Yabana, and J. Burgdörfer, Phys. Rev. A 92, 043422 (2015), URL https://link.aps.org/doi/10.1103/PhysRevA.92.043422.
  * Dimitrovski et al. (2010) D. Dimitrovski, C. P. J. Martiny, and L. B. Madsen, Phys. Rev. A 82, 053404 (2010), URL https://link.aps.org/doi/10.1103/PhysRevA.82.053404.
  * Etches and Madsen (2010) A. Etches and L. B. Madsen, Journal of Physics B: Atomic, Molecular and Optical Physics 43, 155602 (2010), URL https://doi.org/10.1088/0953-4075/43/15/155602.
  * He et al. (2011) F. He, C. Ruiz, A. Becker, and U. Thumm, Journal of Physics B: Atomic, Molecular and Optical Physics 44, 211001 (2011), URL https://doi.org/10.1088/0953-4075/44/21/211001.
  * Chini et al. (2012) M. Chini, B. Zhao, H. Wang, Y. Cheng, S. X. Hu, and Z. Chang, Phys. Rev. Lett. 109, 073601 (2012), URL https://link.aps.org/doi/10.1103/PhysRevLett.109.073601.
  * Blieck et al. (2008) J. Blieck, X. Fléchard, A. Cassimi, H. Gilles, S. Girard, and D. Hennecart, Review of Scientific Instruments 79, 103102 (2008), eprint https://doi.org/10.1063/1.2994151, URL https://doi.org/10.1063/1.2994151.
  * Fischer et al. (2012) D. Fischer, D. Globig, J. Goullon, M. Grieser, R. Hubele, V. L. B. de Jesus, A. Kelkar, A. LaForge, H. Lindenblatt, D. Misra, et al., Phys. Rev. Lett. 109, 113202 (2012), URL https://link.aps.org/doi/10.1103/PhysRevLett.109.113202.
  * Yuan et al. (2020) J. Yuan, Y. Ma, R. Li, H. Ma, Y. Zhang, D. Ye, Z. Shen, T. Yan, X. Wang, M. Weidemüller, et al., Chinese Physics Letters 37, 053201 (2020), ISSN 1741-3540, URL http://dx.doi.org/10.1088/0256-307X/37/5/053201.
  * Williams and Fry (1968) W. L. Williams and E. S. Fry, Phys. Rev. Lett. 20, 1335 (1968), URL https://link.aps.org/doi/10.1103/PhysRevLett.20.1335.
  * Fry and Williams (1969) E. S. Fry and W. L. Williams, Phys. Rev. 183, 81 (1969), URL https://link.aps.org/doi/10.1103/PhysRev.183.81.
  * Larsson et al. (1995) J. Larsson, E. Mevel, R. Zerne, A. L’uillier, C. G. Wahlstrom, and S. Svanberg, Journal of Physics B: Atomic, Molecular and Optical Physics 28, L53 (1995), URL https://doi.org/10.1088/0953-4075/28/3/003.
  * Gisselbrecht et al. (1999) M. Gisselbrecht, D. Descamps, C. Lyngå, A. L’Huillier, C.-G. Wahlström, and M. Meyer, Phys. Rev. Lett. 82, 4607 (1999), URL https://link.aps.org/doi/10.1103/PhysRevLett.82.4607.
  * Johansson et al. (2003) A. Johansson, M. K. Raarup, Z. S. Li, V. Lokhnygin, D. Descamps, C. Lyngå, E. Mevel, J. Larsson, C.-G. Wahlström, S. Aloise, et al., The European Physical Journal D - Atomic, Molecular, Optical and Plasma Physics 22, 3 (2003), ISSN 1434-6079, URL https://doi.org/10.1140/epjd/e2002-00213-1.
  * de Morisson Faria and Lewenstein (2005) C. F. de Morisson Faria and M. Lewenstein, Journal of Physics B: Atomic, Molecular and Optical Physics 38, 3251 (2005), URL https://doi.org/10.1088/0953-4075/38/17/014.
  * Nayak et al. (2019) A. Nayak, M. Dumergue, S. Kühn, S. Mondal, T. Csizmadia, N. Harshitha, M. Füle, M. Upadhyay Kahaly, B. Farkas, B. Major, et al., Physics Reports 833, 1 (2019), ISSN 0370-1573, saddle point approaches in strong field physics and generation of attosecond pulses, URL https://www.sciencedirect.com/science/article/pii/S0370157319303096.
  * Chirilă and Lein (2009) C. C. Chirilă and M. Lein, Phys. Rev. A 80, 013405 (2009), URL https://link.aps.org/doi/10.1103/PhysRevA.80.013405.
  * Milošević (2013b) D. B. Milošević, Phys. Rev. A 88, 023417 (2013b), URL https://link.aps.org/doi/10.1103/PhysRevA.88.023417.