# Controlling quantum effects in enhanced strong-field ionisation with
machine-learning techniques

H Chomet, S Plesnik, C Nicolae, J Dunham, L Gover, T Weaving and C Figueira de
Morisson Faria  Department of Physics and Astronomy, University College
London, Gower Street, London WC1E 6BT, UK

###### Abstract

We study non-classical pathways and quantum interference in enhanced
ionisation of diatomic molecules in strong laser fields using machine learning
techniques. Quantum interference provides a ‘bridge’, which facilitates
intramolecular population transfer. Its frequency is higher than that of the
field, intrinsic to the system and depends on several factors, for instance
the state of the initial wavepacket or the internuclear separation. Using
dimensionality reduction techniques, namely t-distributed stochastic neighbour
embedding (t-SNE) and principal component analysis (PCA), we investigate the
effect of multiple parameters at once and find optimal conditions for enhanced
ionisation in static fields, and controlled ionisation release for two-colour
driving fields. This controlled ionisation manifests itself as a step-like
behaviour in the time-dependent autocorrelation function. We explain the
features encountered with phase-space arguments, and also establish a
hierarchy of parameters for controlling ionisation via phase-space Wigner
quasiprobability flows, such as specific coherent superpositions of states,
electron localisation and internuclear-distance ranges.

##  1 Introduction

Manipulating coherent superpositions of quantum states and nonclassical
pathways has been a central question to many areas of science, such as
coherent control [1, 2, 3] and quantum information [4]. Important applications
include controlling chemical reactions [5] or electron dynamics in ultrafast
molecular dissociation [6, 7, 8, 9], creating quantum switches [10, 11, 4] and
enhancing nonclassicality in extended systems such as light-harvesting
compounds [12, 13]. In this wide range of scenarios, decoherence must be kept
at bay, so that the timescales of interest are much shorter than the
decoherence times. This brings one’s attention to how coherent superpositions
of quantum states and nonclassical pathways may be controlled in attoscience.

Attoscience has emerged from the interaction of matter with very intense laser
fields, typically of the order of $10^{13}\mathrm{W}/\mathrm{cm}^{2}$ and
deals with some of the shortest time scales in nature, of the order of
$10^{-18}s$. These extremely short timescales bring about the possibility of
controlling real-time electron dynamics [14, 15], and have triggered many
applications. Examples are subfemtosecond imaging of matter (for reviews see,
e.g., [16, 17]), high-harmonic spectroscopy [18], ultrafast photoelectron
holography [19, 20], attosecond electron or hole migration [21, 22, 23, 24,
25], and, more recently, strong-field phenomena in solids [26, 27, 28, 29] and
nanostructures [30, 31]. It has even been speculated that, in the future,
laser-induced electron dynamics may lead to optoelectronic computers with
switching rates 100,000 times higher than existing digital electronic systems
[32, 33, 34].

Recently, it has been shown that quantum interference and nonclassical
pathways play a vital role in strong-field enhanced ionisation of diatomic
molecules [35]. Enhanced ionisation is known for over three decades [36, 37],
and consists in an increase of at least one order of magnitude in the
ionisation rate of a stretched molecule, in comparison with that of an atom
with a similar ionisation potential. It has been attributed to the narrowing
of the upfield potential barrier due to the presence of an adjacent, downfield
potential well, and to coupled charge-enhanced resonant states. Although the
phenomenon is widely known, phase-space tools, such as the Wigner quasi-
probability distributions [38] have shed new light on its behaviour [39, 40,
35]. Phase-space quasi-probability distributions have found enormous success
in many research areas such as quantum optics [41, 42, 43], quantum
information [44, 45], chemical physics [46, 47] and cold gases [48], but are
hitherto underused in attoscience (for a recent review on the overall
landscape within this research field see our article [49]). Nevertheless, they
provide valuable insight into the dynamics of these systems of interest and
have been used to study phenomena such as strong-field ionisation [50, 39, 40,
51], rescattering [52, 53] and entanglement [54]. For recent work
investigating different pathways in enhanced ionisation see [55, 56].

In 2008, momentum gates occurring in phase space under strong oscillating
laser fields were identified [39]. These gates facilitate the flow of quasi-
probability from one molecular centre to the other causing many ionisation
bursts within a field cycle. These bursts were later attributed to a non-
adiabatic response to the time-dependent field gradients [57, 40]. However,
our previous publication [35] dispelled this notion, by showing that these
momentum gates are present when there is no time dependence on the field [35],
or even no field at all [58], and are primarily caused by quantum
interference. The term quantum bridge was coined to describe this mechanism.

In the present article, we aim at controlling the pathways behind enhanced
ionisation, which have been identified in [35]. Some pathways may be
understood within a quasi-static, semiclassical picture, whereby the
quasiprobability flows along the field gradient, with tails following
equienergy curves [59, 50, 51]. Other pathways are enabled by quantum
interference, which creates a bridge and provides a passage for direct intra-
molecular quasi-probability flow. For a time-dependent field, one must bear in
mind that the phase-space configuration also changes as time progresses.
Therefore, the field can act both in favour of the cyclic motion of the
quantum bridges (and the enhanced ionisation) or against it. These different
temporal behaviours depend on a wealth of parameters and may be explored in
order to control strong-field ionisation. Recently, it has been shown that
quantum interference plays an important role even in the semiclassical
ionisation pathway [60].

Previous analyses focused on understanding the physical mechanisms behind the
quantum bridges, and were restricted to a homonuclear molecular potential in
the field-free setting, or under a strong static or monochromatic laser field.
The initial wavepacket was taken to be an upfield or downfield localised
Gaussian, or a symmetric coherent superposition thereof known as a cat state
[41]. Here, we focus on how the quantum bridges can be manipulated by using
different relative phases and wavepacket localisations in the initial
superposition state, molecular potentials of heteronuclear type, where each
well is weighted asymmetrically, and electric field configurations. This is a
complex task: With a wealth of tunable parameters, such as internuclear
separation, relative phase and wavepacket localisation in the initial wave
packet, differing nuclear charges, laser intensity and frequency as well as
the pulse shape, it necessitates optimisation techniques. Machine learning
methods have been explored in a wide range of scenarios in quantum physics
(see the reviews [61, 62, 63] and the recent tutorial [64]) and ultrafast
photonics [65]. Recently, they have been employed in the attosecond science
context to control attosecond pulses [66], retrieve the structure of a large
molecule from laser-induced electron diffraction patterns [67] and to track
quantum pathways within a spatio-temporal Feynman path integral framework
[68].

In this work, we employ machine-learning techniques for dimensionality
reduction, such as the T-Distributed Stochastic Neighbour Embedding (t-SNE)
[69] and Principal-Component Analysis (PCA) [70] to find overall trends and
patterns that allow us to establish a hierarchy of parameters and classify
different regimes for strong-field ionisation. Throughout, we are aided by a
qualitative phase-space analysis, which allows us to understand the physics
behind the patterns encountered. We first employ a model molecule in a static
field as a proof of concept, and, subsequently, use time-dependent two-colour
driving fields. The autocorrelation function exhibit several distinct
behaviours in time, including a stepwise profile in which ionisation is
switched on or off. We also provide a physical interpretation for these
features.

This article is organised as follows. In section 2 we discuss our model,
including its general features (section 2.1) and the physical quantities of
interest (section 2.2). Section 3 is devoted to the dimensionality reduction
methods and their implementation. Subsequently, in section 4 we provide a
qualitative analysis (4.1) and a proof of concept (4.2) for static fields, and
in section5 we show how strong-field ionisation in two-colour fields can be
optimised to produce a stepwise temporal behaviour. We commence by examining
the different clusters that result from the t-SNE (section 5.1) and analyse
them using physical arguments (section 5.2). Finally, in section 6 we
summarize our results and provide possible outlooks for the present studies.

##  2 Model

###  2.1 General features

Our system is a one-dimensional single-electron molecule whose wave function
dynamics $\Psi(x,t)$ are obtained by solving the time-dependent Schrödinger
equation (TDSE) in atomic units using the split operator method [71]

$i\frac{\partial}{\partial
t}\Psi(x,t)=\left(-\frac{1}{2}\frac{\partial^{2}}{\partial
x^{2}}+V(x)+xE(t)\right)\Psi(x,t).$ (1)

The quantity $E(t)$ is the external laser field and $V(x)$ is the molecular
binding potential given by

$V(x)=Z_{\mathrm{r}}V_{0}(x-R/2)+Z_{\mathrm{l}}V_{0}(x+R/2),$ (2)

where $R$ is the internuclear distance. The first and second term appearing on
the right-hand side of equation (2) shall be called the left and right
potential wells, respectively. Each potential well $V_{0}(x)$ is chosen to
have soft-core form

$V_{0}(x)=-\frac{1}{\sqrt{x^{2}+a}},$ (3)

where $a$ is known as the softening parameter. This removes the singularity at
$x=0$ of the true Coulomb potential, yet remains long-range [72, 73, 74]. This
parameter is typically chosen such that the ionisation energy coincides with
the ground state of the potential [75]; throughout the following analysis
$a=1$ is chosen. Using a Gaussian wavepacket (Eq. (6)) with $p_{0}=q_{0}=0.0$
and $\gamma=0.5$, this ground state energy is ${{E}_{{\rm sc}}}=-0.67$ a.u..

The homonuclear case is recovered by setting
$Z_{\mathrm{l}}=Z_{\mathrm{r}}=1$. The symbol $Z$ is chosen in reference to
nuclear charge; however, in this toy model, $Z_{i}$ is allowed to be a
continuous variable. Eq. (3) provides qualitative insight into the behaviour
of Coulomb potentials, which suffices for the objectives of the present work.
Recently, however, the soft-core potential has been modified to exhibit
quantitative agreement with realistic three dimensional models [76, 77].

The external laser field is either taken to be static, that is $E(t)=E_{0}$,
or a time-dependent linearly polarised polychromatic field such that

$E(t)=E_{0}\left[\mathrm{cos}(\omega t)+r_{t}\mathrm{cos}(b\omega
t+\phi)\right],$ (4)

where $b$ is the frequency ratio between the first and the second driving
wave, $r_{t}$ is the field-strength ratio, $\phi$ is the relative phase
between the two driving waves, and $\omega$ is the frequency of the
fundamental.

This molecular model is analysed in phase space. There, the bound and
continuum regions may be identified, and an example of is shown in Figure 1.
Fixed points (points that satisfy $\dot{x}=\dot{p}=0$) here are either centres
or saddles [78]. Centres correspond to minima of the effective potential, and
are surrounded by closed orbits. Saddles correspond to maxima of the effective
potential, and delimit different dynamical regions of phase space. For
instance, the central saddle separates the two molecular wells with energy
$E_{C}$, and the Stark saddle is located between the downfield well and the
continuum with energy $E_{S}$. The energy curves of those saddles in phase-
space are called separatrices. We define the energy difference between these
two saddles as

$\Delta E=E_{C}-E_{S}.$ (5)

$\Delta E$ is either negative (Fig. 1(a)) or positive (Fig. 1(c)) and
characterises two different phase space configurations: closed (Fig. 1(b)) and
open (Fig. 1(d)) separatrices.

Figure 1: Phase portraits for the one-dimensional homonuclear molecular models
described by the binding potentials (2), using internuclear separations of
$R=4$ and $R=6$ and a static field $E=0.07.$ a.u. [left panels], together with
the corresponding effective potentials [right panels]. The Stark and the
central saddles are indicated by the labels S and C in the figure, and their
energy by the labels $\mathrm{E_{S}}$ and $\mathrm{E_{C}}$. The field-free
separatrices and potentials are given by the dashed red lines.

The initial state of the electron $\Psi(x,0)$ is approximated as a Gaussian
wavepacket,

$\Psi(x,0)=\langle
x\mid\Psi(0)\rangle=\left(\frac{\gamma}{\pi}\right)^{\frac{1}{4}}\exp\left\\{-\frac{\gamma}{2}\left(x-q_{0}\right)^{2}+\mathrm{i}p_{0}\left(x-q_{0}\right)\right\\},$
(6)

of width $\gamma$ and initial momentum $p_{0}=0$. The state shall be called
localised in the left well when $q_{0}=-R/2$, denoted
$\Psi_{\mathrm{l}}(x,0)$, or in the right well when $q_{0}=R/2$, denoted
$\Psi_{\mathrm{r}}(x,0)$. Alternatively, we may form delocalised coherent
superpositions. We introduce some asymmetry between the initial left and right
wavepackets in a delocalised setting. As such, the initial state in [35] is
generalised by including a relative phase $\theta\in[-\pi,\pi]$ and a
localisation parameter $\alpha\in[0,1]$, giving

$\Psi_{\alpha,\theta}(x,0)\coloneqq\frac{\sqrt{\alpha}\Psi_{\mathrm{l}}(x,0)+\sqrt{1-\alpha}e^{i\theta}\Psi_{\mathrm{r}}(x,0)}{\sqrt{1+2\sqrt{\alpha(1-\alpha)}\cos(\theta)\mathcal{I}_{o}}}.$
(7)

Under some circumstances, we also consider different widths $\gamma_{l}$ and
$\gamma_{r}$ for the wavefunctions $\Psi_{\mathrm{l}}$ and
$\Psi_{\mathrm{r}}$, respectively. The constant of normalisation in the
denominator depends on the initial state overlap $\mathcal{I}_{0}$:

$\mathcal{I}_{0}=\int\Psi_{\text{l}}^{*}(x,0)\Psi_{\text{r}}(x,0)\mathrm{d}x.$
(8)

Note that $\Psi_{1,\theta}(x,0)=\Psi_{\mathrm{l}}(x,0)$ and
$\Psi_{0,\theta}(x,0)=\Psi_{\mathrm{r}}(x,0)$ (up to a global phase), whilst
intermediate values of $\alpha$ produce coherent superpositions of varying
weight. If one considers $\Psi_{\frac{1}{2},0}(x,0)$ and
$\Psi_{\frac{1}{2},\pi}(x,0)$, this will lead to even and odd cat states, for
which the wavepacket localisation in each centre are equally weighted.

###  2.2 Relevant quantities

Of particular interest is the phenomenon of enhanced ionisation, characterised
by anomalous peaks in the ionisation rate for relatively large internuclear
separations. The effect is understood to arise from the presence of strongly
coupled charge resonant states in combination with a narrowing of the
effective potential for the upfield charge centre [36, 37]. The latter is a
result of the neighbouring downfield centre enabling efficient population
transfer from upfield to the continuum. The ionisation rate $\Gamma$ from time
$t=0$ to $t=T$ is quantified using

$\Gamma=-\ln\left(\frac{|\mathcal{P}(T)|^{2}}{|\mathcal{P}(0)|^{2}}\right)\frac{1}{T},$
(9)

where $\mathcal{P}(t)$ is the probability function at $t$, given by

$\mathcal{P}(t)=\int^{+\infty}_{-\infty}\Psi^{*}(x,t)\Psi(x,t)dx.$ (10)

This definition of ionisation rate was used in the seminal paper [36] in the
context of enhanced ionisation of molecules. Numerically, the integral limits
of Eq. (10) are finite, taken throughout the box of size $x=-100$ to $x=100$
a.u. Because of irreversible ionisation Eq. (10) will decrease with time and
be less than unity. No absorber was used, and reflections are minimised by
using a total grid size twice as large as the box size. The final time is
chosen to be $T=150~{}\mathrm{a.u.}$ as it has been used in previous
publications [49, 35, 58]. Throughout, $\Gamma$ will be written in arbitrary
units.

Eq. (9) will be employed throughout in order to test our optimisation methods
for static fields. The static field in this work is chosen in such a way that
the left well is located downfield and the right well upfield.

Another relevant quantity, which provides insight into the resulting wave
function dynamics, is the autocorrelation function

$a(t)=\int\Psi^{*}(x,t)\Psi(x,0)\mathrm{d}x,$ (11)

which is given by the overlap integral between the initial and time propagated
wave function. The autocorrelation function will be employed to assess the
behaviour of the system in time-dependent electric fields.

The Wigner quasi-probability distribution [79] is a useful tool in the
analysis of phase-space dynamics (for a review see [80]). It is defined as

$W(x,p,t)=\frac{1}{\pi}\int_{\infty}^{-\infty}\mathrm{d}\xi\Psi^{*}(x+\xi,t)\Psi(x-\xi,t)\mathrm{e}^{2\mathrm{i}p\xi}.$
(12)

The Wigner function is real and normalised, with the property that its
marginals correspond to physical probability distributions for each conjugate
variable, respectively. W(x,p,t) may assume negative values, and this can
provide an indication of nonclassicality [81].

##  3 Dimensionality Reduction

The present problem depends on many parameters: The initial wavepacket
localisation $\alpha$, its widths $\gamma_{l}$ and $\gamma_{r}$ and relative
phase $\theta$; molecular parameters like the weights $Z_{l}$ and $Z_{r}$ and
the inter-nuclear distance $R$; and many external electric-field parameters
such as $\phi$, $r_{t}$, $b$ and $E_{0}$. The main challenge to overcome when
studying the effect of a large number of parameters simultaneously is how to
visualise the results. Indeed, each data point will exist in a high-
dimensional space (the number of dimensions equal to the number of parameter
used), and ideally we would like to project our results down to a two-
dimensional space.

To this end we will use T-Distributed Stochastic Neighbour Embedding (t-SNE),
an unsupervised dimensionality reduction technique, meaning it will find
overall trends and patterns without any prior knowledge on the origin of the
data. It can be briefly summarised in two steps:

  * •

Define $G_{xy}$ and $C_{xy}$ as the probability of picking a pair of points
(x,y) in high dimensional and two dimensional space respectively. Both
probabilities are characterised such that they are larger if the neighbouring
points are close. $G_{xy}$ is a Gaussian distribution and $C_{xy}$ a Cauchy
distribution.

  * •

Minimise the Kullback–Leibler divergence [82] between the two distributions
with respect to their projection in two-dimensional space using the gradient
descent method [83].

A more detailed explanation of the workings of the t-SNE method can be found
in [69]. As a reality check, we will compare all t-SNE projections to
Principal Component Analysis (PCA) [70] dimensional reduction. Indeed, t-SNE
is stochastic, so each run can yield different results. It focuses on
preserving the local structure of the data and is a non-linear technique. On
the other hand, PCA is linear, deterministic, and preserves global properties
while potentially losing low-variance deviations between neighbours. Therefore
by comparing both projections we can guarantee the accuracy of the t-SNE
results. For clarity, some key results obtained with the PCA are given in the
Appendix A. To make use of the t-SNE capabilities we need to include
information about the wavefunction in the static field case or its evolution
in the time-dependent case.

###  3.1 Time evolution characterisation

In order to understand the system dynamics through the use the dimensionality
reduction methods, we wish to characterise the electronic dynamics by a single
value. In the static-field case we will employ the ionisation rate $\Gamma$
(9). For the time-dependent case we choose to quantify the shape of the
autocorrelation function (11), as it carries more temporal information. A
preliminary investigation showed that the autocorrelation function can take
three distinct shapes (Figure 2): ‘steps’, ‘constant’, and ‘other’, which
represent different behaviours.

The classification of autocorrelation functions into the defined categories is
performed by a simple deterministic decision tree algorithm. First, we check
whether a function is ‘constant’, that is, if it does not fall below a given
threshold value. Next, we look for the steps. Steps are characterised by the
intervals where the autocorrelation function stays approximately constant,
followed by a steep decline. Hence, if we find a constant interval longer than
some ‘critical length’, we categorise it as a ‘step’. We define a constant
interval as an interval where the difference of any two values is smaller than
a critical threshold, ‘$\epsilon$’. If the function is not caught by these
checks, it is classified as ‘other’. In the results presented in section 5 the
threshold value is set to 0.8, the critical length is set to 100, while the
critical $\epsilon$ is 0.03. All of the autocorrelation function arrays used
have 3000 elements.

Figure 2: Examples of the three potential outputs of the autocorrelation
function sorting algorithm. The parameters used vary greatly but are all
within the range shown in Table. 1.

###  3.2 Method

To obtain a data point, $N-1$ parameters are randomly uniformly distributed
within the range shown in Table. 1. The electronic wavefunction is then
evolved using the time-dependent Schrödinger equation and we characterise the
time evolution information by a single value as described in the above section
In case of a static external field studied in section 4 this is done by
computing the ionisation rate $\Gamma$ (Eq. (9)). For the time-dependent field
in section 5 we quantify the shape of the electron autocorrelation function
(see section 3.1). The $N-1$ parameters chosen randomly as well as the single
value characterising the time evolution form a data point of size $N$.

In certain situations specific data points are removed from the final data set
in order to focus the visualisation on the matter at hand. For example in
section 5.1, only 1332 data points are used from the 1 million originally
computed using random parameters. Indeed, in order to focus on the ‘step’
output of the autocorrelation function, approximately 99.9% of data points
with outputs ‘constant’ and 99.97% of data points with output ‘other’ are
removed.

We then use the t-SNE to reduce our N-dimensional data points to 2 dimensions,
projecting the data set onto the two axis $\mathrm{y}_{1}$ and
$\mathrm{y}_{2}$. The data points will have clustered into different groups,
the nature of which will be understood by plotting the projection as a
function of one of the input parameters. This will allow us to determine which
parameter or which combination of parameters influence the ionisation rate or
the autocorrelation function step, and to what degree. The overall method
described here and used throughout this paper is summarised schematically in
Fig. 3. Finally, the PCA method is used to ensure that the observed behaviour
is neither a remnant of the stochastic nature of the algorithm, nor due to the
removing of data points from the data set. Examples of PCA results are shown
in Appendix A.

Static Field | $\alpha$ | $\theta$ | $\gamma_{r}$ | $\gamma_{l}$ | R | $Z_{r}$ | $Z_{l}$ | $E_{0}$
---|---|---|---|---|---|---|---|---
min | 0.0 | 0.0 | 0.2 | 0.2 | 3.0 | 1.0 | 1.0 | 0.04
max | 1.0 | $\pi$ | 1.0 | 1.0 | 12.0 | 2.0 | 2.0 | 0.08

Time-dependent Field | $\alpha$ | $\theta$ | $\gamma_{r}$ | $\gamma_{l}$ | R | $Z_{r}$ | $Z_{l}$ | $E_{0}$ | $b$ | $r_{t}$ | $\phi$
---|---|---|---|---|---|---|---|---|---|---|---
min | 0.0 | 0.0 | 0.2 | 0.2 | 4.0 | 1.0 | 1.0 | 0.04 | 0.5 | 0.0 | -$\pi$
max | 1.0 | $\pi$ | 1.0 | 1.0 | 10.0 | 1.5 | 1.5 | 0.08 | 2.0 | 1.0 | $\pi$

Table 1: Table of all parameters used throughout this paper with the minimum
and maximum value of the range they are randomly chosen from.  Figure 3:
Schematic outline of the methodology used throughout the paper.

##  4 Static fields

We will start our study with an analysis of enhanced ionisation in static
fields. While the ultimate goal is to understand and control enhanced
ionisation in time-dependent fields, the purpose of this static field section
is two-fold:

  * •

Obtain an in depth understanding of how certain parameters influence the
system, and expand on the conclusions drawn in [49, 35].

  * •

Establish the effectiveness of the t-SNE dimensionality reduction technique as
a method of simultaneous multi-parameter analysis. This will allow us to use
it as a primary tool in section 5.

###  4.1 Qualitative analysis

As a starting point, we will perform a qualitative investigation of how
specific parameters influence the quantum bridges and the ionisation rate in a
static field using the generalised initial wavepacket given by Eq. (6). We
employ phase-space arguments and Wigner functions to facilitate the
interpretation.

####  4.1.1 Initial electron wavepacket.

Our earlier work on the subject [35] led to the conclusion that a major
predictor of enhanced ionisation was the localisation of the initial electron
wavepacket. However, only three configurations were studied: localised around
the right well, localised around the left well, and delocalised equally around
both centres. Here we start with a more thorough analysis of this parameter by
considering it a continuous variable, $\alpha$, see Eq. (7). As shown in Fig.
4, the ionisation rate is linearly proportional to the electron localisation,
with $\alpha=0.0$ (initially localised upfield) leading to the highest
ionisation rate. This can be understood as the downfield centre narrowing the
effective-potential barrier for the upfield centre, hence enhancing upfield
ionisation.

Figure 4: Ionisation rate $\Gamma$ for an initial wave packet
$\Psi_{\alpha,\theta}(x,0)$ with respect to the localisation parameter
$\alpha$, comparing maximum $\theta=2.4$ and minimum $\theta=-0.8$. The
external static field has strength $E_{0}=0.06$, and the molecular parameter
used are internuclear distance $R=7$, $Z_{r}=1.0$ and $Z_{l}=1.0$.

In addition to a continuous localisation parameter, the state of the initial
electron wavepacket is varied by a relative phase factor $\theta$, see Eq.
(7), between the superposition of $\Psi_{r}$ and $\Psi_{l}$. When using an
equally delocalised wavepacket, $\alpha=0.5$, the ionisation rate as a
function of $\theta$ (in rad) peaks around $\theta=2.4$ and is minimal for
$\theta=-0.8$. This has little overall influence on the ionisation rate as
seen in Fig. 4. Finally, the initial wavepacket widths $\gamma_{r}$ and
$\gamma_{l}$, see Eq. (6), are set to differ from $\gamma=0.5$, which is the
minimal ground-state energy of a field free single-centre soft-core potential
with softening parameter $a=1.0$. This has no major influence on the
ionisation rate and is not shown here.

####  4.1.2 Heteronuclear molecule $Z_{l,r}$.

Finally, we expand our study to include heteronuclear molecules. We analyse
the ionisation rate while varying the value of molecular weights $Z_{r}$ and
$Z_{l}$, of Eq. 2. The study of heteronuclear molecules is motivated by
enabling a change in the critical internuclear distance where the phase-space
separatrices change configuration and in the dynamics of the quantum bridges.
From Fig. 5 studying the ionisation rate as a function of both charges $Z_{r}$
and $Z_{l}$ simultaneously shows a large region of suppressed ionisation as
$Z_{r}>1.3$. The effect of $Z_{l}$ is a lot less pronounced, leading to high
ionisation rates even at $Z_{l}=2.0$. While open separatrices do not guarantee
a high ionisation rate, closed separatrices always suppress ionisation. To
have a better understanding of the sharp decrease in ionisation rate with
respect to $Z_{r}$, we look at the Wigner quasi-probability distribution in
Fig. 6. Comparing the bound region around the right molecular centre in (a)
and (b), as $Z_{r}$ increases, the right molecular well deepens greatly. This
leads the upfield population to stay trapped and greatly suppresses
ionisation. On the other hand, the increase of the left bound region in (c)
due to an increase in $Z_{l}$ does not greatly influence the ionisation rate.
This is because in enhanced ionisation the majority of the ionised population
comes from the upfield centre, as shown in [35]. There is very little change
to the upfield (right-hand-side here) centre when increasing $Z_{l}$.

Figure 5: Contour plot of the ionisation rate $\Gamma$ as a function of
$Z_{l}$ and $Z_{r}$ using an external static field $E_{0}=0.06$ with
internuclear distance $R=7$. The white lines indicate the separatrix energy
difference $\Delta E$, with the 0.0 line indicating the transition from open
separatrices ($\Delta E$ > 0) to closed separatrices ($\Delta E$ < 0). Figure
6: Wigner quasi-probability distribution of the ionisation from a static field
$E_{0}=0.06$ with internuclear distance $R=7$ for (a) $Z_{r}=1.0$ and
$Z_{l}=1.0$, (b) $Z_{r}=1.5$ and $Z_{l}=1.0$ and (c) $Z_{r}=1.0$ and
$Z_{l}=1.5$. The snapshots are taken at $t=25~{}\mathrm{a.u.}$.

Now, a more complete overview would go through the effect of $Z_{l}$ and
$Z_{r}$ with respect to various parameters, which would be slow and
exhaustive. Instead we propose in the following section to use dimensionality
reduction techniques to analyse the effect of multiple parameters at once.

###  4.2 Proof of concept - machine learning techniques

We now aim to study the effect of multiple parameters on the ionisation rate,
starting with a limited number in order to fully understand the results and to
focus on the effect of the chosen parameters. Those are the internuclear
distance R, the width of the initial superposed gaussian upfield and downfield
wavepackets $\gamma_{l}$ and $\gamma_{r}$, the phase difference $\theta$
between them, and the localisation parameter of the wavepacket $\alpha$. The
value of these parameters is chosen from a predetermined range described in
section 3.2.

Data points comprised of these five parameters as well as the ionisation rate
are projected onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the
t-SNE algorithm, shown in Figure 7. This method is described in section 3.2.

Figure 7: Plot of the 3000 6-dimensional data points that have been projected
onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the t-SNE
dimensional reduction technique described in section 2.2. Each subplot has the
projected data plotted as a function an initial input parameter. The initial
input parameters include (a) the ionisation rate $\Gamma$ (given in arbitrary
units), (b) the internuclear distance $R$, (c) the electron localisation
$\alpha$. Also included but not shown here are the left initial wavepacket
width $\gamma_{l}$, the right wavepacket width $\gamma_{r}$ and the phase
difference between the left and right initial wavepackets $\theta$. PCA
distributions using the same data points have been analysed for accuracy (not
shown). The external field is taken to be static with $E_{0}=0.0534$ and the
molecular weights $Z_{r}=1.0$ and $Z_{l}=1.0$, matching the results from [35].

From Figure 7, we see that the projected data bring into focus the connection
between the ionisation rate, the internuclear distance and the electron
localisation. These match the results from [35] and section 4.1. That is, peak
ionisation yield is present at and around the peak internuclear distance of
$R=6.8$ a.u., shown in Figure 7(b), and when the initial wavepacket is
localised upfield, shown in Figure 7(c). The low ionisation rate regions, seen
in Figure 7(a), split into two distinct parts seen in Fig. 7(b) which
correspond respectively to regions of too low and too high internuclear
distance, studied in depth in [35]. The ionisation rate falls linearly as the
electron localisation changes from localised upfield to localised downfield,
as seen in Fig. 4. The effect of the gaussian widths $\gamma$ and the phase
difference $\theta$ is a lot weaker than those discussed above, and are not
shown here.

We will now look at the full range of parameters. In addition to the previous
six, these following results will include the effects of the nuclear weight
number $Z_{r}$ and $Z_{l}$, and the electric field strength $E_{0}$. Because
of these additional parameters, we posit that there is no longer a clear
internuclear distance peak, but instead a larger range of optimal internuclear
distances, dependent on $E_{0}$, $Z_{r}$ and $Z_{l}$. For that reason the
difference in energy between the central and Stark separatrix $\Delta E$ will
also be analysed, but will not be an input parameter.

As seen in Figure 8(a), the algorithm successfully separates two clusters of
low ($\Gamma<0.002$) and high ($\Gamma>0.05$) ionisation yield. Within the
high ionisation cluster, the data points are further distinguished by a very
high ($\Gamma>0.3$) ionisation yield region. From the distribution of other
parameters in these two clusters, we can determine which lead to a high
ionisation yield.

Figure 8: Plot of the 16 995 9-dimensional data points that have been
projected onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the t-SNE
dimensional reduction technique described in section 2.2. Each subplot has the
projected data plotted as a function an initial input parameter. The initial
input parameters include (a) the ionisation rate $\Gamma$, (b) the right
molecular well depth $Z_{r}$, (c) the electron localisation $\alpha$, (d) the
field strength $E_{0}$, (e) the internuclear distance $R$. Also included but
not shown here are the left molecular well depth $Z_{l}$, the left initial
wavepacket width $\gamma_{l}$, the right wavepacket width $\gamma_{r}$ and the
phase difference between the left and right initial wavepackets $\theta$. In
subplot (f) the separatrix energy difference $\Delta E$ is shown for each data
point, however it is not one of the input parameters. The input data set has
been reduced from 100 000 data points to only have either $\Gamma<0.002$ or
$\Gamma>0.05$. PCA distributions using the original 100 000 data points are
shown in Appendix A.

The strongest predictor of high ionisation rate is the nuclear depth factor
$Z_{r}$, [see Fig. 8(b)], which in this configuration is associated to the
upfield molecular well. As shown in Figure 5, the ionisation rate falls of
sharply as soon as $Z_{r}>1.0$. While the value of $Z_{r}$ is originally
chosen randomly in the range [1:2], the average value of $Z_{r}$ in the high
ionisation yield cluster is 1.16, and 1.12 in the $\Gamma>0.3$ region. Indeed,
90% of the $\Gamma>0.3$ ionisation yield data points have $Z_{r}<1.25$.

In Fig. 8(c), echoing the conclusions from Fig. 7 and Fig. 4, only 4% of the
high ionisation yield data points have the initial wavepacket of the electron
localised downfield, with $\alpha>0.8$. Similarly, 30% of the high ionisation
yield cluster data points are for an initial wavepacket localised upfield,
with $\alpha<0.2$. This proportion goes up to 72% for $\Gamma>0.3$.

From Fig. 8(d), the median initial field strength $E_{0}$ of the high
ionisation yield cluster is 0.07, despite the parameter being originally
randomly chosen from a [0.04:0.08] range. This means the majority of high
ionisation yield data points have a high field strength. Indeed, 82% of data
points have $E_{0}>0.06~{}\mathrm{a.u.}$ in that cluster and 88% in the
$\Gamma>0.3$ region.

As hypothesised earlier, the clear internuclear distance peak from Fig. 7 is
now less pronounced in Fig. 8 (e) as we are now also varying the molecular
weights $Z$ and the field amplitude $E_{0}$. Still, 72% of data points within
the very high, $\Gamma>0.3$, ionisation yield region have internuclear
distances $R$ in the range $6<R<9$. This is surprising as the initial range
for $R$ is randomly taken from [3:12]. Therefore, one would expect that one
third of the data points would be in this range.

Finally, in Figure 8(f), we can deduce the effect of the separatrices on the
ionisation yield. Indeed if the separatrices are nested (see Fig. 1), meaning
that $E_{S}>E_{C}$, or $\Delta E<0$, then the ionisation rate is suppressed.
In contrast, 99.5% data points in the high ionisation yield cluster have
$\Delta E>-0.5$, and 86% of them have $\Delta E>0$. This separation is present
despite $\Delta E$ not being used as an input parameter.

In conclusion, the t-SNE enables us to visualise the effect of multiple
parameters simultaneously. All the conclusions drawn in section 4.1 are
visualised in Fig. 8. Parameters with little to no effect are swiftly
identified and separated from relevant parameters. Not only that, but a
hierarchy of the effect of different parameters is also established much more
rapidly. For example the effect of $Z_{l}$ on the ionisation rate being
greater than that of $\alpha$. Because of this we will be able to use the
t-SNE in the following section to rapidly assert what parameters are relevant
and to what degree.

##  5 Time dependent fields - control of enhanced ionisation - optimising
autocorrelation step function

When using a polychromatic field (see Eq. 4) and specific parameters, we can
obtain a ‘controlled’ ionisation release, which translates to a step function
in the autocorrelation function plot. In the following section we will be
looking at what parameters lead to step functions and why specific
configurations lead to a controlled ionisation release.

###  5.1 t-SNE application to step functions

####  5.1.1 The 3 clusters and electron localisation.

Figure 9: Plot of the 1332 12-dimensional data points that have been projected
onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the t-SNE
dimensionality reduction technique described in section 2.2. Each subplot has
the projected data plotted as a function of an initial input parameter (shown
by the colour bar). All units are in [a.u.]. unless stated otherwise. The
input parameters shown in this figure are (a) The step function shape obtained
by the sorting algorithm in section 3.1 that outputs ‘constant’, ‘other’, or
the length of the step in the autocorrelation function and (b) The initial
electron wavepacket localisation represented by the parameter $\alpha$ in Eq.
(7). The full list of input parameters as well as details on the initial data
set used can be found in section 2.2. PCA distributions using 500 000 data
points are shown in Appendix A.

The t-SNE algorithm separates the data into three different clusters, whose
nature can be seen in the two plots of Fig. 9. Therein, the data points are
coloured with respect to (a) their autocorrelation function type and (b) their
initial electron localisation $\alpha$. Indeed, one cluster groups all
autocorrelation functions classified as ‘constant’ or ‘other’, while the other
two contain the autocorrelation functions classified as ‘steps’. Those two
clusters are understood by looking at the distribution with respect to the
initial electron localisation. All autocorrelation functions that are
classified as ‘steps’ have their initial electron wavepacket localised
completely at one centre. The clusters are separated between having the
initial electron wavepacket localised at the left centre ($\alpha>0.95$, or
dark red) and localised at the right centre ($\alpha<0.05$, or dark blue).
From this it is clear that the initial wavepacket must be localised either to
the left or to the right for there to be steps. This is explained when looking
at Fig. 14, in section 5.2.1.

####  5.1.2 Difference between the two ‘step’ clusters and the ‘other’ and
‘constant’ cluster.

Figure 10: Plot of the 1332 12-dimensional data points that have been
projected onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the t-SNE
dimensional reduction technique described in section 2.2. Each subplot has the
projected data plotted as a function an initial input parameter (shown by the
colour bar). All units are in [a.u.]. unless stated otherwise. The input
parameter shown in this figure are (a) The ratio $r_{t}$ between the two
colour field amplitudes, (b) the overall field amplitude $E_{0}$ in Eq. (4)
and (c) the internuclear distance $R$. The full list of input parameters as
well as details on the initial data points used can be found in section 2.2.
PCA distributions using 500 000 data points have been checked for accuracy
(not shown).

From looking at the different clusters, we can find other correlations, shown
in Fig. 10, between certain parameters and the presence of step functions.
most notably, the field ratio $r_{t}$ between the two driving waves, in
subplot (a). Indeed, while originally the value of $r_{t}$ is chosen randomly
from a range [0:1], the two step clusters have a median $r_{t}$ value of 0.81
and an average of 0.75. Moreover, longer step lengths (more than 150 a.u., or
orange/red) have a $r_{t}$ median value of 0.85 and an average of 0.8.

Next, the overall electric field strength $E_{0}$, chosen at random from a
range of [0.04:0.08], has an median and average value of
$E_{0}=0.07~{}\mathrm{a.u.}$. These can be understood as by seeing that both a
high field ratio $r_{t}$ and a high field strength $E_{0}$ lead to a higher
field peak during which more of the population will escape, causing a deeper
drop in the autocorrelation function.

Finally, the internuclear distance is clearly closer to the peak internuclear
distance of $R=7~{}\mathrm{a.u.}$, which is also confirmed by the standard
deviation being only 1.05 instead of 1.3. $R$ is also for $92.6\%$ of all data
in the step clusters above a minimum internuclear distance of
$R=6~{}\mathrm{a.u.}$. The cause of this effect is more complex, and is
discussed in detail in section 5.2.2.

####  5.1.3 Symmetry between the ‘localised to the left’ and ‘localised to the
right’ step clusters.

Figure 11: Plot of the 1332 12-dimensional data points that have been
projected onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the t-SNE
dimensionality reduction technique described in section 2.2. Each subplot has
the projected data plotted as a function an initial input parameter (shown by
the colour bar). All units are in [a.u.]. unless stated otherwise. The input
parameters shown in this figure are (a) the phase difference between the left
and right initial wavepackets $\theta$, (b) the left initial wavepacket width
$\gamma_{l}$, (c) the right wavepacket width $\gamma_{r}$, (d) the left
molecular well depth $Z_{l}$ and (e) the right molecular well depth $Z_{r}$.
The full list of input parameters as well as details on the initial data
points used can be found in section 2.2. PCA distributions using 500 000 data
points have been checked for accuracy (not shown).

The value of $\alpha$ is not the only difference between the two step
clusters. The next series of parameters, shown in Fig. 11, are connected to
the symmetry between the left and right molecular wells. First, there is quite
naturally the distribution of $\gamma_{l}$ and $\gamma_{r}$. From their
definition in Eq. (6), $\gamma_{l}$ only influences the initial wavepacket
localised to the left. Therefore, if the initial wavepacket was localised to
the left, only the initial left wavepacket width $\gamma_{l}$ is a relevant
parameter, and vice versa. Since, as shown in Fig. 9, both step clusters use
initially localised wavepackets exclusively, only one $\gamma$ parameter
influences each cluster and the $\theta$ parameter (subplot (a)) has no effect
on them. As seen in Fig. 11 (b) and (c), 22.8% of $\gamma_{l}$ in the
localised right step cluster are above 0.8, but only 6.8% are above 0.8 in the
localised left step cluster. Similarly, 27.5% of $\gamma_{r}$ in the localised
left step cluster are above 0.8, while 5.3% are above 0.8 in the localised
right step cluster. The width $\gamma=0.5$ corresponds to the minimal ground-
state energy of a field free single-centre soft-core potential with our
current parameters, and the vast majority of step functions have a $\gamma$
value that is within the range [0.2 : 0.8]. As $\gamma$ deviates from 0.5,
part of the initial wavepacket bleeds out of the bound region despite being
initially localised. This leads to oscillations during the flat ‘step’ portion
of the autocorrelation function and explains the absence of initial wavepacket
widths $\gamma$ above 0.8 in the step clusters.

The other set of parameters symmetric with respect to the left/right direction
are the molecular weights $Z_{l}$ and $Z_{r}$. Their effect on one cluster
should be mirrored on the other. From Figs. 11(d) and (e), while there is no
clear connection between the localised left cluster and $Z_{r}$, $Z_{l}$ is
more likely to be close to 1.0 (and vice-versa). Indeed there are 42% of data
points that have $Z_{l}<1.2$ in the localised right cluster, but 70% in the
localised left cluster. Similarly, 39% of data points have $Z_{r}<1.2$ in the
localised left cluster, but 67% in the localised right cluster. This is
further investigated in section 5.2.2.

####  5.1.4 Link between different parameters and the step length.

Figure 12: Plot of the 1332 12-dimensional data points that have been
projected onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the t-SNE
dimensionality reduction technique described in section 2.2. Each subplot has
the projected data plotted as a function of an initial input parameter (shown
by the colour bar). All units are in [a.u.] unless stated otherwise. The input
parameters shown in this figure are (a) The field-frequency ratio $b$ and (b)
the field offset $\phi$ in Eq. (4). The full list of input parameters as well
as details on the initial data points used can be found in section 2.2. PCA
distributions using 500 000 data points have been checked for accuracy (not
shown).

Going back to the distribution with respect to the autocorrelation function
type in Fig. 9, those two ‘step’ clusters are then also distinguished into
regions of ‘short’ (less than 150 a.u., or lime yellow) and ‘long’ (more than
150 a.u., or orange/red) step length. From Fig. 12, we can already glean that
the field-frequency ratio $b$ affects the step length. However, there is no
clear link with the field offset $\phi$. This is further expanded upon in
section 5.2.3.

###  5.2 Analysis of the step functions

In order to obtain a controlled burst of ionisation, two pre-requisites are
needed: First, a time interval in which the population stays bounded, and,
second, a short burst of enhanced ionisation. Following the conclusions drawn
from section 5.1, we already know what parameter range is required or
preferred. In this section, we will provide an intuitive, physics-based
analysis for the controlled-ionisation burst conditions.

####  5.2.1 Electron localisation.

From Fig. 9, it is clear that the initial electron wavepacket must be
localised, either around the left or right molecular well. Fig. 13 illustrates
the difference between using an initial localised wavepacket and a delocalised
wavepacket. This, along with the Wigner function shown in Fig. 14, allows us
to understand the mechanism behind the controlled ionisation release.

To obtain a burst of enhanced ionisation, the system must be in an optimal
configuration. This was initially found in a study using static fields in
[35], and further expanded upon in section 4. A key point is that the initial
wavepacket must be localised upfield. This allows the population to escape
directly through the quantum bridge to the semiclassical escape pathway and
ultimately the continuum. It also stops the quantum bridge from cycling
through the momentum space and bringing the population back to the upfield
centre. Similarly, the population must stay bounded for an interval of time,
meaning that the system must be in the configuration with lowest ionisation
rate. As seen in Fig. 4, for static fields, this happens when the initial
wavepacket is localised downfield. For the time dependent polychromatic field,
the terms ‘upfield’ and ‘downfield’ become relative as the configuration
changes with the direction of the field. However in either case the initial
wavepacket must be localised on one of the molecular centres. As seen in the
behaviour of the Wigner function in Fig. 14 [right column] these two
conditions are not possible if the initial wavepacket is delocalised. Because
of the quantum interference in between the two molecular centres, there is
always the cyclic motion of the quantum bridge.

Figure 13: Comparison of auto-correlation functions at R = 6.8 a.u. with
initial wavepacket localised upfield or delocalised. Direction of the field is
denoted by the colour of the background, blue regions are positive values,
pink are negative. The field is shown as the orange dashed lines with values
shown on the rightmost y-axis and has form $E\propto\cos{\omega
t}+\cos{2\omega t}$ where $\omega=0.057$ a.u. The stars denote autocorrelation
function values at times $t=10~{}\mathrm{a.u.}$, $t=25~{}\mathrm{a.u.}$,
$t=50~{}\mathrm{a.u.}$, and $t=85~{}\mathrm{a.u.}$.

Figure 14:  Wigner quasiprobability distributions computed for internuclear
distances R = 6.8 a.u. with initial wavepacket localised at the right
potential well ($\alpha=0.0$) [left column] or delocalised between the two
molecular wells ($\alpha=0.5$)[right column]. The external field has frequency
$\omega=0.057$ a.u. and strength $E_{0}=0.0534~{}\mathrm{a.u.}$. The snap
shots are taken at times $t=10~{}\mathrm{a.u.}$ for (a) and (a’),
$t=25~{}\mathrm{a.u.}$ for (b) and (b’), $t=50~{}\mathrm{a.u.}$ for (c) and
(c’), and $t=85~{}\mathrm{a.u.}$ for (d) and (d’). These times are indicated
by stars in Fig. 13.

The localisation of the wavepacket must work in conjunction with the sign of
the external field for the stepwise behaviour to be present in the
autocorrelation function. If the wavepacket is localised to the left and the
external field $E(t)$ is positive, or if the wavepacket is localised to the
right and $E(t)$ is negative, this will result in downfield localisation.
Alternatively, upfield localisation occurs for a wavepacket localised to the
left and negative $E(t)$, or a wavepacket localised to the right and positive
$E(t)$.

In order to obtain a step function, the localisation must be upfield during
the drop and downfield during the flat part of the step. In Fig. 13, the
different field signs are shown by the different coloured shaded areas.
Indeed, during the initial drop, the external field is positive and the
population localised to the right, and therefore upfield. It escapes rapidly
through the quantum bridge as seen by the Wigner function in Fig. 14(a) and
the separatrices are open. As soon as the field changes sign, the
approximately constant behaviour in the autocorrelation function begins. As
seen in Figs. 14(b), (c) and (d), during that time the wavepacket localisation
is now downfield. Because of that, despite the separatrices being sometimes
open, like in subplots (b) and (d), ionisation is suppressed. When the field
is no longer strongly negative around $t=60$, the separatrices are closed, as
shown in subplot (c). This also stops the population from escaping. The
autocorrelation function therefore has a step until the field becomes positive
again, with the separatrices strongly open. From this analysis and the
clustering obtained by the t-SNE distribution, it is clear the electron
localisation is the most important parameter in obtaining a controlled
ionisation release.

####  5.2.2 Separatrices.

Figure 15: Different plots of autocorrelation functions [Top row] and $\Delta
E$ [Bottom row] for $R=4~{}\mathrm{a.u.}$, $R=7~{}\mathrm{a.u.}$ and
$R=10~{}\mathrm{a.u.}$ with [Left column] $E_{0}=0.04~{}\mathrm{a.u.}$,
[Middle column] $E_{0}=0.06~{}\mathrm{a.u.}$ and [Right column]
$E_{0}=0.08~{}\mathrm{a.u.}$. The external laser field is represented by the
red dotted line (not to scale). We consider an initial wavepacket localised to
the right, with $\gamma_{r}$ = 0.5, and equal molecular weights
$Z_{r}=Z_{l}=1.0$. The external field amplitude and frequency ratios are
$r_{t}=1.0$ and $b=0.5$, respectively.  Figure 16: Evolution of the Wigner
function with internuclear distance $R=4~{}\mathrm{a.u.}$ and external field
parameters $r_{t}=1.0$, $b=0.5$ and $E_{0}=0.04~{}\mathrm{a.u.}$. The snap
shots are taken at times (a) $t=65~{}\mathrm{a.u.}$, (b)
$t=85~{}\mathrm{a.u.}$ and (c) $t=100~{}\mathrm{a.u.}$. Initial wavepacket
localised to the right, with $\gamma_{r}$ = 0.5. Molecular weights
$Z_{r}=Z_{l}=1.0$. Figure 17: Different plots of autocorrelation functions and
$\Delta E$ for [Left column] $Z_{r}=1.0~{}\mathrm{a.u.}$,
$Z_{r}=1.25~{}\mathrm{a.u.}$ and $Z_{r}=1.5~{}\mathrm{a.u.}$ with
$Z_{l}=1.0~{}\mathrm{a.u.}$. and [Right column] $Z_{l}=1.0~{}\mathrm{a.u.}$,
$Z_{l}=1.25~{}\mathrm{a.u.}$ and $Z_{l}=1.5~{}\mathrm{a.u.}$ with
$Z_{r}=1.0~{}\mathrm{a.u.}$. The initial wavepacket is localised to the right
and other parameters are set to $R=7~{}\mathrm{a.u.}$,
$E=0.06~{}\mathrm{a.u.}$, $r_{t}=1.0$ and $b=0.5$. The external laser field is
represented by the red dotted line (not to scale).

On top of the electron localisation, the presence of step functions is
strongly dependent on the value of $E_{0}$, $r_{t}$, $Z_{l}$, $Z_{r}$ and $R$,
as seen in Fig. 10 and Fig. 11. Both $E_{0}$ and $r_{t}$ have a similar effect
on the external laser field: when $r_{t}$ is close to 1.0 and when $E_{0}$ is
large, the field peak increases in absolute value. As shown in Fig. 8, when
using a static field, increasing $E_{0}$ leads to higher ionisation yield. In
our time-dependent case, when $E_{0}$ or $r_{t}$ is high, the value of $\Delta
E$ at the field peak is much higher, meaning the separatrices are wide open.
This phenomenon is best seen in Fig. 15 for a fixed value of $R$. In the
second row it is clear that the increase in $E_{0}$ leads to higher maximum
values of $\Delta E$. For example, at internuclear distance $R=10$,
$E_{0}=0.04$ has $\Delta E=0.7$ at the field maximum, while $E_{0}=0.08$ has
$\Delta E=1.6$. Compared to the corresponding autocorrelation functions in the
first row, we see that when matched by the correct electron localisation, this
will lead too a much steeper drop during the decreasing part of the
autocorrelation function step. Using our previous example, at $E_{0}=0.04$ the
first step plateaus at 0.9, while $E_{0}=0.08$ at 0.5.

The effect of the internuclear distance $R$ is at a glance not clear. When
focusing on the initial drop, the connection between $\Delta E$ and the
initial drop of the correlation seem contradictory. Intuitively one would
expect a that a higher maximum $\Delta E$ lead to a steeper drop of the
autocorrelation function, as was the case for varying $r_{t}$ and $E_{0}$. But
the opposite appears true with varying $R$. This is understood by results in
[35], and illustrated in Fig. 15. As $R$ is becomes too large, the quantum
bridge between the two molecular wells weakens and the ionisation rate falls.
As consequence the initial drop of the autocorrelation function is low. With
$E_{0}$ being low as well, this leads to a ‘constant’ autocorrelation function
as seen in Fig. 15 for $R=10$ and $E=0.04$. When $R$ is too small, the quantum
bridges between the molecular wells are at their strongest, and initially
facilitate population transfer, leading to a steeper drop. However, as seen in
section 4, a small $R$ does not lead to a high ionisation rate. That is
understood in results from [35], and illustrated by the Wigner functions in
Fig.16: the quantum bridge brings population back to the right-side core. This
cyclic motion does not follow the direction of the external field. As a
consequence, the quasiprobability distribution does not stay within one centre
during the ‘flat’ portion of the step autocorrelation function, as illustrated
in Fig.15.

While the effect of $Z_{l}$ and $Z_{r}$ differs greatly depending on the step
cluster analysed, it is mirrored, so conclusions on one can be expanded to
both. Focusing on the initially localised to the right cluster from Fig.11, we
see that while $Z_{l}$ has little influence, the vast majority of step
autocorrelation functions has $Z_{r}<1.2$. This is illustrated in Fig.17: The
initial drop of the step autocorrelation function greatly reduces as $Z_{r}$
increases. From the corresponding $\Delta E$ plot, the peak separatrix energy
difference goes from 0.7 for $Z_{r}$ = 1.0 to 0.6 for $Z_{r}$ = 1.5. Compared
to the effect of the separatrix energy difference of $E_{0}$ and its effect on
the initial autocorrelation function drop, it appears something else must be
at play. A better understanding is obtained when looking at the ionisation
rate in static fields as a function of $Z_{r}$ and $Z_{l}$. As explained in
section 4.1, the biggest difference between Fig. 6 is not the separatrix
energy difference, but the range of the bound region, especially in the
momentum space. Back to the time dependent situation, as $Z_{r}$ increases the
population (here initially localised to the right) stays trapped at the right
molecular well, and the initial step drop will be much shorter.

####  5.2.3 Step length and two colour field frequencies.

From Fig. 12 we can already observe that the step length is connected to the
value of the two colour field frequency ratio $b$. Indeed, while step
functions of step length less than $110~{}\mathrm{[a.u.]}$ have on average
$b=1.5$, step function of step length greater than $150~{}\mathrm{[a.u.]}$
have on average $b=0.9$.

Figure 18: Plot of the 27 194 8-dimensional data points that have been
projected onto two axis $\mathrm{y}_{1}$ and $\mathrm{y}_{2}$ using the t-SNE
dimensional reduction technique described in section 2.2. All units are in
[a.u.]. unless stated otherwise. Each subplot has the projected data plotted
as a function a initial input parameter. The initial input parameters include
(a) The output of the step function algorithm from section 3.1 that is the
length of the step in the autocorrelation function and (b) the frequency-ratio
of the second colour field $b$ as well as not shown here the offset of the
second field $\phi$, the field ratio $r_{t}$, the internuclear distance $R$,
the left molecular well depth $Z_{l}$, the right molecular well depth $Z_{r}$,
the right initial wavepacket width $\gamma_{r}$. All data points use an
initially localised to the right wavepacket with $E_{0}=0.07$. The input data
set has been reduced from 1 000 000 data points to only look at step
autocorrelation functions. PCA distributions using the original 1 000 000 data
points have been analysed for accuracy (not shown).

Figure 19: Plot of the corresponding autocorrelation function and external
field in each of the four different clusters in Fig. 18, in order of
increasing field frequency ration $b$. Constant parameters are $E_{0}=0.07$,
$R=7$, $\alpha=0.0$, $\gamma_{r}=0.5$ . The parameters used in each subplot
are (a) from C1; $b=0.5$, $r_{t}=1.0$, $\phi=-0.2$, $Z_{r}=1.2$ and
$Z_{l}=1.2$, (b) from C2; $b=0.7$, $r_{t}=1.0$, $\phi=-0.2$, $Z_{r}=1.2$ and
$Z_{l}=1.2$, (c) from C3; $b=1.3$, $r_{t}=1.0$, $\phi=1.0$, $Z_{r}=1.1$ and
$Z_{l}=1.1$ and (d) from C4; $b=2.0$, $r_{t}=1.0$, $\phi=-0.13$, $Z_{r}=1.05$
and $Z_{l}=1.5$.

In order to view the relationship between the step length and other parameters
in greater detail, we will now fix certain parameters in order to get a
clearer picture of the effect the others, as well as only conserve data points
with ‘step’ autocorrelation functions. First, we will start with the initial
wavepacket localised at the right centre. This is equivalent to focusing only
on the ‘localised to the right’ cluster of Fig. 9. Because of this, we also
remove the initial wavepacket width $\gamma_{l}$ and $\theta$ parameters.
Finally, we will fix the field strength to $E_{0}=0.07~{}\mathrm{a.u.}$. These
results are shown in Fig. 18.

The t-SNE algorithm separates the data into four different clusters. The C1
and C3 clusters represent step lengths between $110$ and $260$. They are
distinguished by the value of the frequency ratio $b$ of the two-color fields:
the largest C1 having an average of $b=0.65$ and the smallest C3 $b=1.32$.
Next, the C2 cluster is the smallest and corresponds to all autocorrelation
functions with a very long (greater than $260$) step length. Finally, the C4
cluster is the largest cluster and corresponds to all autocorrelation
functions with a short step length. Here the average step length is $1.09$. An
example of autocorrelation function and external field within each cluster is
shown in Fig. 19.

As seen in the previous sections, the temporal location of the field peak,
when matched with the population localisation, determines the step drop. From
this, we can assume that the specific times for which the field peaks occur
match the step length. This can be seen in Fig. 19 (a) and (b), when
increasing the second field frequency-ratio from $b=0.5$ to $b=0.7$, the
distance between two field peaks increases and the step length goes from 187
to 285. This is also confirmed when looking at an example step function from
clusters (c) and (d): The distance between field peaks matches the step
length.

##  6 Conclusions

In this work we use dimensionality reduction techniques, namely the the
t-distributed stochastic neighbor embedding (t-SNE) method and Principal
component analysis (PCA) to investigate the effect of multiple parameters at
once in enhanced ionisation. We show how quantum effects in strong-field
ionisation of stretched diatomic molecules may be understood, classified and
synchronised with the external laser field in order to create tailored
ionisation bursts. Thereby, the t-SNE was crucial to establish a hierarchy of
parameters, and manipulate the relevant time scales. These time scales are
either dictated by the field, or by the molecular system. The latter are
associated to quantum-interference structures that provide a direct, intra-
molecular population transfer via a quantum pathway in phase space, the
quantum bridges [35, 58]. The results presented here have been mainly obtained
with the t-SNE, while the PCA was used to rule out possible algorithm-
dependent artefacts.

The effectiveness and accuracy of the t-SNE is first showcased on the
ionisation rate using a static field. The expected results are found using
phase-space arguments and Wigner quasiprobability distributions. They
completely match the conclusions drawn from the t-SNE distributions.
Furthermore, a hierarchy in the importance of different parameters is quickly
and effectively demonstrated. Because of this, we confidently used the t-SNE
to project data of autocorrelation functions and their parameters for time-
dependent fields. Using a two-colour field, we obtained a controlled
ionisation release, which translates to a stepwise behaviour in the
autocorrelation function. By using the t-SNE to project multi-dimensional data
sets into 2 dimensions, ionisation profiles obtained while varying multiple
parameters at once are organised and understood.

A conclusion from previous work [35] was that a key to understand enhanced
ionisation is the interplay between the intra-molecular quantum bridges and
the external field. This also holds for controlled ionisation release
(represented by a step function in the autocorrelation function), where the
system requires a short burst of enhanced ionisation followed by a time
interval for which the population stays bounded.

Quantum bridges are tied to many predictors for a step function in the
autocorrelation function, including the strongest one: electron localisation.
A delocalised initial wavepacket creates an interference pattern in both
positive and negative momentum space. This leads to ionisation bursts that do
not follow the external field. Therefore, we cannot obtain a time interval for
which the population stays bounded. Moreover, the quantum bridges are highly
tied to the internuclear distance. For that reason, configurations that
suppress momentum gates do not lead to the burst of ionisation needed. This is
the case when the internuclear distance is too large. Alternatively, a too
small internuclear distance leads to cyclic motion in the momentum space,
regardless of the initial state [35]. This also stops the population from
staying bounded at one centre.

On the other hand, a higher field peak maximum leads to stronger ionisation
bursts. A stronger field peak is achieved by increasing the field amplitude
$E_{0}$ or adding a second field. This is neatly quantifiable by the
separatrix energy difference $\Delta E$. From this we determine bounded
regions with closed separatrices and ionisation burst regions when the
localisation is upfield with open separatrices. The timing of these field
peaks matching with the upfield population lead to the start of the
autocorrelation function drops. The field frequency determines the distance
between the field peaks and is therefore the key to controlling the length of
the steps. Finally, changing the molecular weights $Z_{r,l}$ affects the
presence of step functions. As the upfield molecular weight increases, the
upfield bound region increases and this greatly suppresses ionisation.

In summary, we can understand the physical cause and requirements for
controlled ionisation release by separating the phenomenon into two steps, a
short burst of enhanced ionisation and a time interval where the population
stays bounded. For static fields, the optimal configuration for enhanced
ionisation is pinpointed to an initial wavepacket localised upfield at optimal
internuclear distance with a high intensity external field and a low upfield
molecular weight. The separatrices are open with a high separatrix energy
difference. This also holds during the time-dependent ionisation burst time
interval, and will be periodic as the upfield and downfield configurations
will change along with the sign of the external field, matching the field
frequency. For the population to stay bounded, either the field strength must
be low, separatrices closed or barely open, or the population must be
localised in the downfield centre and therefore not follow the changes in the
field strength and direction. In addition, the population must not cycle
around through momentum gates, meaning the internuclear distance cannot be too
short. The t-SNE was greatly effective at separating results into clusters and
therefore visualising results due to specific combinations of parameters. For
example, an overall analysis of the effect of $Z_{r}$ on the step function
autocorrelation function reveals only negligible effects. However, the step
function results were separated into two clusters $\alpha<0.2$ and
$\alpha>0.8$. The effect of the increased bound region and ionisation
suppression from $Z_{r}$ only affecting the $\alpha<0.2$ cluster is brought
into focus. Moreover, the t-SNE presented a hierarchy of the correlation
between different parameters and the quantity of study. From that, the role of
electron localisation was pinpointed as the key parameter.

The importance of electron localisation is in agreement with the findings by
several groups, reported over many years. In fact, resolving and ultimately
controlling electron localisation in extended or dissociating molecules holds
the promise of steering chemical reactions. This has led to many studies, both
theoretical and experimental, employing, for instance, pump-probe schemes [84,
85, 7, 8, 86, 9], CEP stabilised few-cycle pulses [87, 6], long-wavelength
fields [88, 89] or synthesised wave forms [90]. Moreover, orthogonally
polarised fields have also been applied to trace or control the ionisation
site [91, 92, 93]. Regarding enhanced ionisation, electron localisation in the
upfield well, together with coupled charge-resonant states, were widely
mentioned in the literature as major contributing factors. This holds both for
seminal enhanced ionisation papers [36, 37] and more recent work in which
models for time-dependent ionisation bursts have been developed [94, 57, 40,
95].

Since potential wells occur in a wide range of physical systems, the present
techniques can be applied not only to more realistic molecular models, but
also to solids and nanostructures. Thereby, a crucial issue would be to
incorporate other degrees of freedom, multielectron dynamics and assess the
role of decoherence. Loss of coherence can be caused by many physical
mechanisms, such as coupling with additional degrees of freedom, intensity
fluctuations, and incoherent emission from across the focal volume (for a
brief discussion of some decoherence mechanisms in the context of strong-field
quantum sensing see [96]). Recent studies using pump-probe schemes in
$\mathrm{H}_{2}$ have shown that coupling of electronic and vibrational
degrees of freedom affects coherence and may hinder electron localisation
[97]. Still, studies in large molecules indicate that phase relations may be
preserved even in systems with many degrees of freedom [98], and could even be
controlled leading to tailored ionisation enhancement or suppression [99]. Our
simplified model is meant as a proof of concept and a first step towards more
realistic scenarios.

## Acknowledgements

We would like to thank Eduard Oravkin, Pavol Drotar and Samo Hromadka for
useful discussions. This work was funded by grant No. EP/J019143/1, from the
UK Engineering and Physical Sciences Research Council (EPSRC).

##  Appendix A Principal Component Analysis

In this Appendix we will present a few of the PCA results obtained throughout
this study. The aim of PCA is to compute the eigenvectors of the covariance
matrix (principal components) of the multidimensional data set and recast the
data along the two largest eigenvectors. In other words, the first principal
component can be defined as the direction that maximises the variance of the
projected data. Throughout this paper all t-SNE distributions were computed
along with PCA distributions for comparison. In this appendix we present two
examples, one from the static field analysis in section 4 and one from the
time-dependent step autocorrelation function analysis in section 5, Fig. 20
and Fig. 21 respectively.

Figure 20: Plot of the 100 000 9-dimensional data points that have been
projected onto their two principal component axis using PCA. The initial
parameters used are equal to that of Figure. 8. Each subplot has the projected
data plotted as a function an initial input parameter. The initial input
parameters include (a) the ionisation rate $\Gamma$, (b) the electron
localisation $\alpha$ and (c) the right molecular well depth $Z_{r}$. Figure
21:  Plot of the 500 000 12-dimensional data points that have been projected
onto their two principal component axis using PCA. The initial parameters used
are equal to that of Figure. 9. Each subplot has the projected data plotted as
a function an initial input parameter. The initial input parameters include
(a) The step function shape obtained by the sorting algorithm in section 3.1
that outputs ‘constant’, ‘other’, or the length of the step in the
autocorrelation function and (b) the electron localisation $\alpha$.

The PCA distribution in Fig. 20, obtained for the ionisation rate data set
calculated using the static field, does not exhibit individual clusters
corresponding to different ionisation rates. However, the distribution is
organised following the different ionisation rates. From this we can see that
both upfield electron localisation ($\alpha=0$) and a right molecular weight
$Z_{r}=1.0$ lead to the highest ionisation rates.

When it comes to the time-dependent autocorrelation function, the PCA
distribution in Fig. 21 does break into three distinct clusters corresponding
to the different autocorrelation outputs: constant, step, and other. This
allows us to draw similar conclusions to those obtained with the t-SNE, for
instance that the electron must be localised to the right or left in order for
the autocorrelation function to have a stepwise shape. However, contrary to
the t-SNE distribution shown in Fig. 9, the step autocorrelation initially
localised to the right data points and those initially localised to the left
are not separated. That means that parameters that affect those two groups
differently will not have their effect visually represented, for example
$Z_{r}$ and $Z_{l}$.

## References

  * [1] Moshe Shapiro and Paul Brumer  “Coherent Control of Atomic, Molecular, and Electronic Processes” 42, Advances In Atomic, Molecular, and Optical Physics  Academic Press, 2000, pp. 287–345  DOI: https://doi.org/10.1016/S1049-250X(08)60189-5
  * [2] Christopher Bäuerle et al.  “Coherent control of single electrons: a review of current progress”  In _Reports on Progress in Physics_ 81.5 IOP Publishing, 2018, pp. 056503  DOI: 10.1088/1361-6633/aaa98a
  * [3] Y. Hikosaka et al.  “Coherent control in the extreme ultraviolet and attosecond regime by synchrotron radiation”  In _Nature Communications_ 10.1, 2019, pp. 4988  DOI: 10.1038/s41467-019-12978-w
  * [4] Andreas Reiserer and Gerhard Rempe  “Cavity-based quantum networks with single atoms and optical photons”  In _Rev. Mod. Phys._ 87 American Physical Society, 2015, pp. 1379–1418  DOI: 10.1103/RevModPhys.87.1379
  * [5] Dominik Peller et al.  “Sub-cycle atomic-scale forces coherently control a single-molecule switch”  In _Nature_ 585.7823, 2020, pp. 58–62  DOI: 10.1038/s41586-020-2620-2
  * [6] M.. Kling et al.  “Control of Electron Localization in Molecular Dissociation”  In _Science_ 312.5771, 2006, pp. 246–248  DOI: 10.1126/science.1126259
  * [7] G. Sansone et al.  “Electron localization following attosecond molecular photoionization”  In _Nature_ 465.7299, 2010, pp. 763–766  DOI: 10.1038/nature09084
  * [8] K.. Singh et al.  “Control of Electron Localization in Deuterium Molecular Ions using an Attosecond Pulse Train and a Many-Cycle Infrared Pulse”  In _Phys. Rev. Lett._ 104 American Physical Society, 2010, pp. 023001  DOI: 10.1103/PhysRevLett.104.023001
  * [9] H. Xu et al.  “Observing electron localization in a dissociating H2+ molecule in real time”  In _Nature Communications_ 8.1, 2017, pp. 15849  DOI: 10.1038/ncomms15849
  * [10] J.. Cirac, P. Zoller, H.. Kimble and H. Mabuchi  “Quantum State Transfer and Entanglement Distribution among Distant Nodes in a Quantum Network”  In _Phys. Rev. Lett._ 78 American Physical Society, 1997, pp. 3221–3224  DOI: 10.1103/PhysRevLett.78.3221
  * [11] H.. Kimble  “The quantum internet”  In _Nature_ 453.7198, 2008, pp. 1023–1030  DOI: 10.1038/nature07127
  * [12] Mark M. Wilde, James M. McCracken and Ari Mizel  “Could light harvesting complexes exhibit non-classical effects at room temperature?”  In _Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences_ 466.2117, 2010, pp. 1347–1363  DOI: 10.1098/rspa.2009.0575
  * [13] Edward J. O’Reilly and Alexandra Olaya-Castro  “Non-classicality of the molecular vibrations assisting exciton energy transfer at room temperature”  In _Nature Communications_ 5.1, 2014, pp. 3012  DOI: 10.1038/ncomms4012
  * [14] Ferenc Krausz and Misha Ivanov  “Attosecond physics”  In _Rev. Mod. Phys._ 81 American Physical Society, 2009, pp. 163–234  DOI: 10.1103/RevModPhys.81.163
  * [15] Franck Lépine, Misha Y. Ivanov and Marc J.. Vrakking  “Attosecond molecular dynamics: fact or fiction?”  In _Nature Photonics_ 8 Nature Publishing Group, a division of Macmillan Publishers Limited. All Rights Reserved., 2014, pp. 195  DOI: 10.1038/nphoton.2014.25
  * [16] Manfred Lein  “Molecular imaging using recolliding electrons”  In _Journal of Physics B: Atomic, Molecular and Optical Physics_ 40.16, 2007, pp. R135  URL: http://stacks.iop.org/0953-4075/40/i=16/a=R01
  * [17] P Salières et al.  “Imaging orbitals with attosecond and Ångström resolutions: toward attochemistry?”  In _Reports on Progress in Physics_ 75.6, 2012, pp. 062401  URL: http://stacks.iop.org/0034-4885/75/i=6/a=062401
  * [18] J P Marangos  “Development of high harmonic generation spectroscopy of organic molecules and biomolecules”  In _Journal of Physics B: Atomic, Molecular and Optical Physics_ 49.13, 2016, pp. 132001  URL: http://stacks.iop.org/0953-4075/49/i=13/a=132001
  * [19] Y Huismans et al.  “Time-Resolved Holography”  In _Science_ 331.6013, 2011, pp. 61–65  DOI: 10.1126/science.1198450
  * [20] C Figueira Morisson Faria and A S Maxwell  “It is all about phases: ultrafast holographic photoelectron imaging”  In _Reports on Progress in Physics_ 83.3 IOP Publishing, 2020, pp. 034401  DOI: 10.1088/1361-6633/ab5c91
  * [21] Olga Smirnova et al.  “High harmonic interferometry of multi-electron dynamics in molecules”  In _Nature_ 460.7258, 2009, pp. 972–977  DOI: 10.1038/nature08253
  * [22] Y. Mairesse et al.  “High Harmonic Spectroscopy of Multichannel Dynamics in Strong-Field Ionization”  In _Phys. Rev. Lett._ 104 American Physical Society, 2010, pp. 213601  DOI: 10.1103/PhysRevLett.104.213601
  * [23] F. Calegari et al.  “Ultrafast electron dynamics in phenylalanine initiated by attosecond pulses”  In _Science_ 346.6207 American Association for the Advancement of Science, 2014, pp. 336–339  DOI: 10.1126/science.1254061
  * [24] Alexander I. Kuleff, Nikolai V. Kryzhevoi, Markus Pernpointner and Lorenz S. Cederbaum  “Core Ionization Initiates Subfemtosecond Charge Migration in the Valence Shell of Molecules”  In _Phys. Rev. Lett._ 117 American Physical Society, 2016, pp. 093002  DOI: 10.1103/PhysRevLett.117.093002
  * [25] Francesca Calegari et al.  “Charge migration induced by attosecond pulses in bio-relevant molecules”  In _Journal of Physics B: Atomic, Molecular and Optical Physics_ 49.14, 2016, pp. 142001  URL: http://stacks.iop.org/0953-4075/49/i=14/a=142001
  * [26] Shambhu Ghimire et al.  “Strong-field and attosecond physics in solids”  In _Journal of Physics B: Atomic, Molecular and Optical Physics_ 47.20 IOP Publishing, 2014, pp. 204030  DOI: 10.1088/0953-4075/47/20/204030
  * [27] G. Vampa et al.  “Strong-field optoelectronics in solids”  In _Nature Photonics_ 12.8, 2018, pp. 465–468  DOI: 10.1038/s41566-018-0193-5
  * [28] Chao Yu, Shicheng Jiang and Ruifeng Lu  “High order harmonic generation in solids: a review on recent numerical methods”  In _Advances in Physics: X_ 4.1 Taylor & Francis, 2019, pp. 1562982  DOI: 10.1080/23746149.2018.1562982
  * [29] Lisa Ortmann and Alexandra S. Landsman  “Chapter Two - High-harmonic generation in solids” 70, Advances In Atomic, Molecular, and Optical Physics  Academic Press, 2021, pp. 103–156  DOI: https://doi.org/10.1016/bs.aamop.2021.04.002
  * [30] M F Ciappina et al.  “Attosecond physics at the nanoscale”  In _Reports on Progress in Physics_ 80.5 IOP Publishing, 2017, pp. 054401  DOI: 10.1088/1361-6633/aa574e
  * [31] M F Ciappina and M Lewenstein  “21st Century Nanoscience – A Handbook”  CRC Press, 2019  DOI: 10.1201/9780367333003
  * [32] Agustin Schiffrin et al.  “Optical-field-induced current in dielectrics”  In _Nature_ 493.7430, 2013, pp. 70–74  DOI: 10.1038/nature11567
  * [33] Martin Schultze et al.  “Controlling dielectrics with the electric field of light”  In _Nature_ 493.7430, 2013, pp. 75–78  DOI: 10.1038/nature11720
  * [34] Annkatrin Sommer et al.  “Attosecond nonlinear polarization and light–matter energy transfer in solids”  In _Nature_ 534.7605 Nature Publishing Group, 2016, pp. 86–90
  * [35] Heloise Chomet, Dhruva Sarkar and C Figueira Morisson Faria  “Quantum bridges in phase space: interference and nonclassicality in strong-field enhanced ionisation”  In _New. J. Phys._ 21.12 IOP Publishing, 2019, pp. 123004  DOI: 10.1088/1367-2630/ab55fa
  * [36] T. Zuo and A.. Bandrauk  “Charge-resonance-enhanced ionization of diatomic molecular ions by intense lasers”  In _Phys. Rev. A_ 52 American Physical Society, 1995, pp. R2511–R2514  DOI: 10.1103/PhysRevA.52.R2511
  * [37] Tamar Seideman, M.. Ivanov and P.. Corkum  “Role of Electron Localization in Intense-Field Molecular Ionization”  In _Phys. Rev. Lett._ 75 American Physical Society, 1995, pp. 2819–2822  DOI: 10.1103/PhysRevLett.75.2819
  * [38] E. Wigner  “On the Quantum Correction For Thermodynamic Equilibrium”  In _Phys. Rev._ 40 American Physical Society, 1932, pp. 749–759  DOI: 10.1103/PhysRev.40.749
  * [39] Feng He, Andreas Becker and Uwe Thumm  “Strong-field modulated diffraction effects in the correlated electron-nuclear motion in dissociating H 2+”  In _Phys. Rev. Lett._ 101.21 APS, 2008, pp. 213002  DOI: 10.1103/PhysRevLett.101.213002
  * [40] Norio Takemoto and Andreas Becker  “Time-resolved view on charge-resonance-enhanced ionization”  In _Phys. Rev. A_ 84.2 APS, 2011, pp. 023401  DOI: 10.1103/PhysRevA.84.023401
  * [41] Wolfgang P Schleich  “Quantum optics in phase space”  John Wiley & Sons, 2011  DOI: 10.1002/3527602976
  * [42] Stephen M. Barnett and Paul M. Radmore  “Methods in theoretical quantum optics”  Clarendon Press, 2005
  * [43] Ulf Leonhardt  “Essential quantum optics: from quantum measurements to black holes”  Cambridge University Press, 2010  DOI: 10.1017/CBO9780511806117
  * [44] Samuel L. Braunstein and Peter Loock  “Quantum information with continuous variables”  In _Rev. Mod. Phys._ 77 American Physical Society, 2005, pp. 513–577  DOI: 10.1103/RevModPhys.77.513
  * [45] Alessio Serafini  “Quantum continuous variables: a primer of theoretical methods”  CRC Press, TaylorFrancis Group, 2017
  * [46] William H. Miller  “The Semiclassical Initial Value Representation: A Potentially Practical Way for Adding Quantum Effects to Classical Molecular Dynamics Simulations”  In _The Journal of Physical Chemistry A_ 105.13 American Chemical Society, 2001, pp. 2942–2955  DOI: 10.1021/jp003712k
  * [47] William H. Miller  “Quantum dynamics of complex molecular systems”  In _Proceedings of the National Academy of Sciences_ 102.19 National Academy of Sciences, 2005, pp. 6660–6664  DOI: 10.1073/pnas.0408043102
  * [48] P.B. Blakie et al.  “Dynamics and statistical mechanics of ultra-cold Bose gases using c-field techniques”  In _Advances in Physics_ 57.5 Taylor & Francis, 2008, pp. 363–455  DOI: 10.1080/00018730802564254
  * [49] H. Chomet and C. Figueira de Morisson Faria  “Attoscience in phase space”  In _The European Physical Journal D_ 75.7, 2021, pp. 201  DOI: 10.1140/epjd/s10053-021-00199-0
  * [50] Attila Czirják et al.  “The Wigner function for tunneling in a uniform static electric field”  In _Opt. Commun._ 179.1-6 Elsevier, 2000, pp. 29–38  DOI: 10.1016/S0030-4018(99)00591-X
  * [51] C Zagoya et al.  “Quantum and semiclassical phase-space dynamics of a wave packet in strong fields using initial-value representations”  In _New. J. Phys._ 16.10 IOP Publishing, 2014, pp. 103040  DOI: 10.1088/1367-2630/16/10/103040
  * [52] Stefanie Gräfe, Jens Doose and Joachim Burgdörfer  “Quantum phase-space analysis of electronic rescattering dynamics in intense few-cycle laser fields”  In _Journal of Physics B: Atomic, Molecular and Optical Physics_ 45.5 IOP Publishing, 2012, pp. 055002  DOI: 10.1088/0953-4075/45/5/055002
  * [53] C. Baumann, H.-J. Kull and G.. Fraiman  “Wigner representation of ionization and scattering in strong laser fields”  In _Phys. Rev. A_ 92 American Physical Society, 2015, pp. 063420  DOI: 10.1103/PhysRevA.92.063420
  * [54] Attila Czirják, Szilárd Majorosi, Judit Kovács and Mihály G Benedict  “Emergence of oscillations in quantum entanglement during rescattering”  In _Physica Scripta_ T153 IOP Publishing, 2013, pp. 014013  DOI: 10.1088/0031-8949/2013/t153/014013
  * [55] Han Xu et al.  “Experimental observation of the elusive double-peak structure in R-dependent strong-field ionization rate of H2+”  In _Scientific Reports_ 5.1, 2015, pp. 13527  DOI: 10.1038/srep13527
  * [56] Kunlong Liu and Ingo Barth  “Distinguishing two mechanisms for enhanced ionization of ${\mathrm{H}}_{2}{}^{+}$ using orthogonal two-color laser fields”  In _Phys. Rev. A_ 103 American Physical Society, 2021, pp. 013103  DOI: 10.1103/PhysRevA.103.013103
  * [57] Norio Takemoto and Andreas Becker  “Multiple Ionization Bursts in Laser-Driven Hydrogen Molecular Ion”  In _Phys. Rev. Lett._ 105 American Physical Society, 2010, pp. 203004  DOI: 10.1103/PhysRevLett.105.203004
  * [58] D Kufel, Heloise Chomet and C Figueira Morisson Faria  “Alternative quantisation condition for wavepacket dynamics in a hyperbolic double well”  In _J. Phys. A Math. Theor._ 54.3 IOP Publishing, 2020, pp. 035304  DOI: 10.1088/1751-8121/abd267
  * [59] NL Balazs and A Voros  “Wigner’s function and tunneling”  In _Ann. Phys._ 199.1 Elsevier, 1990, pp. 123–140  DOI: 10.1016/0003-4916(90)90370-4
  * [60] Szabolcs Hack et al.  “Quantum interference in strong-field ionization by a linearly polarized laser pulse and its relevance to tunnel exit time and momentum”  In _Physical Review A_ 104.3 APS, 2021, pp. L031102
  * [61] Vedran Dunjko and Hans J Briegel  “Machine learning and artificial intelligence in the quantum domain: a review of recent progress”  In _Reports on Progress in Physics_ 81.7 IOP Publishing, 2018, pp. 074001  DOI: 10.1088/1361-6633/aab406
  * [62] Giuseppe Carleo et al.  “Machine learning and the physical sciences”  In _Rev. Mod. Phys._ 91 American Physical Society, 2019, pp. 045002  DOI: 10.1103/RevModPhys.91.045002
  * [63] Juan Carrasquilla  “Machine learning for quantum matter”  In _Advances in Physics: X_ 5.1 Taylor & Francis, 2020, pp. 1797528  DOI: 10.1080/23746149.2020.1797528
  * [64] Anna Dawid et al.  “Modern applications of machine learning in quantum sciences”  arXiv, 2022  DOI: 10.48550/ARXIV.2204.04198
  * [65] Goëry Genty et al.  “Machine learning and applications in ultrafast photonics”  In _Nature Photonics_ 15.2, 2021, pp. 91–101  DOI: 10.1038/s41566-020-00716-4
  * [66] Chon-Teng Belmiro Chu, Yae-Lin Sheu and Shih-I Chu  “Bayesian optimal control of the ultrashort circularly polarized attosecond pulse generation by two-color polarization gating”  In _Opt. Express_ 29.21 OSA, 2021, pp. 32900–32909  DOI: 10.1364/OE.438212
  * [67] Xinyao Liu et al.  “Machine learning for laser-induced electron diffraction imaging of molecular structures”  In _Communications Chemistry_ 4.1, 2021, pp. 154  DOI: 10.1038/s42004-021-00594-z
  * [68] Xiwang Liu et al.  “Deep Learning for Feynman’s Path Integral in Strong-Field Time-Dependent Dynamics”  In _Phys. Rev. Lett._ 124 American Physical Society, 2020, pp. 113202  DOI: 10.1103/PhysRevLett.124.113202
  * [69] Laurens Maaten and Geoffrey Hinton  “Visualizing Data using t-SNE”  In _Journal of Machine Learning Research_ 9.86, 2008, pp. 2579–2605  URL: http://jmlr.org/papers/v9/vandermaaten08a.html
  * [70] Karl Pearson F.R.S.  “LIII. On lines and planes of closest fit to systems of points in space”  In _The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science_ 2.11 Taylor & Francis, 1901, pp. 559–572  DOI: 10.1080/14786440109462720
  * [71] MD Feit, JA Fleck Jr and A Steiger  “Solution of the Schrödinger equation by a spectral method”  In _J. Comput. Phys._ 47.3 Elsevier, 1982, pp. 412–433  DOI: 10.1063/1.444501
  * [72] Juha Javanainen, Joseph H Eberly and Qichang Su  “Numerical simulations of multiphoton ionization and above-threshold electron spectra”  In _Phys. Rev. A_ 38.7 APS, 1988, pp. 3430  DOI: 10.1103/PhysRevA.38.3430
  * [73] JH Eberly  “Scale variation in a one-dimensional model of an atom interacting with a strong laser field”  In _Phys. Rev. A_ 42.9 APS, 1990, pp. 5750  DOI: 10.1103/PhysRevA.42.5750
  * [74] Q Su and JH Eberly  “Model atom for multiphoton physics”  In _Phys. Rev. A_ 44.9 APS, 1991, pp. 5997  DOI: 10.1103/PhysRevA.44.5997
  * [75] SC Rae, X Chen and K Burnett  “Saturation of harmonic generation in one-and three-dimensional atoms”  In _Phys. Rev. A_ 50.2 APS, 1994, pp. 1946  DOI: 10.1103/PhysRevA.50.1946
  * [76] Szilárd Majorosi, Mihály G. Benedict and Attila Czirják  “Improved one-dimensional model potentials for strong-field simulations”  In _Phys. Rev. A_ 98 American Physical Society, 2018, pp. 023401  DOI: 10.1103/PhysRevA.98.023401
  * [77] Szilárd Majorosi et al.  “Density-based one-dimensional model potentials for strong-field simulations in $\mathrm{He}$, $\mathrm{H}_{2}^{+}$, and $\mathrm{H}_{2}$”  In _Phys. Rev. A_ 101 American Physical Society, 2020, pp. 023405  DOI: 10.1103/PhysRevA.101.023405
  * [78] D.K. Arrowsmith and C.M. Place  “Dynamical Systems: Differential Equations, Maps, and Chaotic Behaviour”  New York: Routledge, 1992
  * [79] E. Wigner  “On the Quantum Correction For Thermodynamic Equilibrium”  In _Phys. Rev._ 40 American Physical Society, 1932, pp. 749–759  DOI: 10.1103/PhysRev.40.749
  * [80] J. Weinbub and D.. Ferry  “Recent advances in Wigner function approaches”  In _Applied Physics Reviews_ 5.4, 2018, pp. 041104  DOI: 10.1063/1.5046663
  * [81] MG Benedict and A Czirják  “Wigner functions, squeezing properties, and slow decoherence of a mesoscopic superposition of two-level atoms”  In _Phys. Rev. A_ 60.5 APS, 1999, pp. 4034  DOI: 10.1103/PhysRevA.60.4034
  * [82] S. Kullback and R.. Leibler  “On Information and Sufficiency”  In _The Annals of Mathematical Statistics_ 22.1 Institute of Mathematical Statistics, 1951, pp. 79–86  DOI: 10.1214/aoms/1177729694
  * [83] Claude Lemaréchal  “Cauchy and the Gradient Method”, 2012
  * [84] G.. Yudin et al.  “Attosecond photoionization of coherently coupled electronic states”  In _Phys. Rev. A_ 72 American Physical Society, 2005, pp. 051401  DOI: 10.1103/PhysRevA.72.051401
  * [85] Feng He, Camilo Ruiz and Andreas Becker  “Control of Electron Excitation and Localization in the Dissociation of $\mathrm{H}_{2}{}^{+}$ and Its Isotopes Using Two Sequential Ultrashort Laser Pulses”  In _Phys. Rev. Lett._ 99 American Physical Society, 2007, pp. 083002  DOI: 10.1103/PhysRevLett.99.083002
  * [86] Zhuo Wang, Kunlong Liu, Pengfei Lan and Peixiang Lu  “Control of electron localization in highly excited states with two ultraviolet laser pulses”  In _Journal of Physics B: Atomic, Molecular and Optical Physics_ 48.1 IOP Publishing, 2014, pp. 015601  DOI: 10.1088/0953-4075/48/1/015601
  * [87] Vladimir Roudnev, B.. Esry and I. Ben-Itzhak  “Controlling ${\mathrm{H}\mathrm{D}}^{+}$ and ${\mathrm{H}}_{2}^{+}$ Dissociation with the Carrier-Envelope Phase Difference of an Intense Ultrashort Laser Pulse”  In _Phys. Rev. Lett._ 93 American Physical Society, 2004, pp. 163601  DOI: 10.1103/PhysRevLett.93.163601
  * [88] Kunlong Liu, Weiyi Hong, Qingbin Zhang and Peixiang Lu  “Wavelength dependence of electron localization in the laser-driven dissociation of H2$+$”  In _Opt. Express_ 19.27 OSA, 2011, pp. 26359–26369  DOI: 10.1364/OE.19.026359
  * [89] Kunlong Liu, Qingbin Zhang, Pengfei Lan and Peixiang Lu  “Anomalous isotopic effect on electron-directed reactivity by a 3-$\mu$m midinfrared pulse”  In _Opt. Express_ 21.4 OSA, 2013, pp. 5107–5116  DOI: 10.1364/OE.21.005107
  * [90] Pengfei Lan, Eiji J. Takahashi and Katsumi Midorikawa  “Efficient control of electron localization by subcycle waveform synthesis”  In _Phys. Rev. A_ 86 American Physical Society, 2012, pp. 013418  DOI: 10.1103/PhysRevA.86.013418
  * [91] J. Wu et al.  “Probing the tunnelling site of electrons in strong field enhanced ionization of molecules”  In _Nature Communications_ 3.1, 2012, pp. 1113  DOI: 10.1038/ncomms2130
  * [92] J. Wu et al.  “Coulomb Asymmetry in Strong Field Multielectron Ionization of Diatomic Molecules”  In _Phys. Rev. Lett._ 108 American Physical Society, 2012, pp. 043002  DOI: 10.1103/PhysRevLett.108.043002
  * [93] Kunlong Liu and Ingo Barth  “Identifying the Tunneling Site in Strong-Field Ionization of ${{\mathrm{H}}_{2}}^{+}$”  In _Phys. Rev. Lett._ 119 American Physical Society, 2017, pp. 243204  DOI: 10.1103/PhysRevLett.119.243204
  * [94] Freek Kelkensberg, Giuseppe Sansone, Misha Y. Ivanov and Marc Vrakking  “A semi-classical model of attosecond electron localization in dissociative ionization of hydrogen”  In _Phys. Chem. Chem. Phys._ 13 The Royal Society of Chemistry, 2011, pp. 8647–8652  DOI: 10.1039/C1CP20058E
  * [95] Norio Takemoto and Andreas Becker  “Visualization and interpretation of attosecond electron dynamics in laser-driven hydrogen molecular ion using Bohmian trajectories”  In _The Journal of Chemical Physics_ 134.7, 2011, pp. 074309  DOI: 10.1063/1.3553178
  * [96] A.. Maxwell, A. Serafini, S. Bose and C. de Morisson Faria  “Quantum estimation in strong fields: In situ ponderomotive sensing”  In _Phys. Rev. A_ 103 American Physical Society, 2021, pp. 043519  DOI: 10.1103/PhysRevA.103.043519
  * [97] Marc J.. Vrakking  “Control of Attosecond Entanglement and Coherence”  In _Phys. Rev. Lett._ 126 American Physical Society, 2021, pp. 113203  DOI: 10.1103/PhysRevLett.126.113203
  * [98] Jorge Delgado et al.  “Molecular fragmentation as a way to reveal early electron dynamics induced by attosecond pulses”  In _Faraday Discuss._ The Royal Society of Chemistry, 2021  DOI: 10.1039/D0FD00121J
  * [99] Shuai Li et al.  “Control of electron recollision and molecular nonsequential double ionization”  In _Communications Physics_ 3.1 Nature Publishing Group, 2020, pp. 1–8