# On-the-fly ab initio semiclassical evaluation of electronic coherences in
polyatomic molecules reveals a simple mechanism of decoherence

Nikolay V. Golubev<EMAIL_ADDRESS>Tomislav Begušić  Jiří Vaníček
<EMAIL_ADDRESS>Laboratory of Theoretical Physical Chemistry, Institut
des Sciences et Ingénierie Chimiques, Ecole Polytechnique Fédérale de Lausanne
(EPFL), CH-1015 Lausanne, Switzerland

###### Abstract

Irradiation of a molecular system by an intense laser field can trigger
dynamics of both electronic and nuclear subsystems. The lighter electrons
usually move on much faster, attosecond time scale but the slow nuclear
rearrangement damps ultrafast electronic oscillations, leading to the
decoherence of the electronic dynamics within a few femtoseconds. We show that
a simple, single-trajectory semiclassical scheme can evaluate the electronic
coherence time in polyatomic molecules accurately by demonstrating an
excellent agreement with full-dimensional quantum calculations. In contrast to
numerical quantum methods, the semiclassical one reveals the physical
mechanism of decoherence beyond the general blame on nuclear motion. In the
propiolic acid, the rate of decoherence and the large deviation from the
static frequency of electronic oscillations are quantitatively described with
just two semiclassical parameters—the phase space distance and signed area
between the trajectories moving on two electronic surfaces. Because it
evaluates the electronic structure on the fly, the semiclassical technique
avoids the “curse of dimensionality” and should be useful for preselecting
molecules for experimental studies.

###### pacs:

Valid PACS appear here

††preprint:

Recent progress in laser technologies Krausz and Ivanov (2009); Young  _et
al._ (2018); Mauritsson  _et al._ (2018) has revolutionized the field of
atomic and molecular physics. In particular, tremendous developments of
coherent light sources enabled the creation of sub-femtosecond laser pulses
with remarkably well controlled parameters Krausz (2016). Using state-of-the-
art lasers, one is able to initiate and probe processes that are driven solely
by the electron correlation, i.e., to study and manipulate electron dynamics
on its natural time scale Kraus  _et al._ (2015).

Experimental measurements of the electron motion in isolated atoms were
reported Goulielmakis  _et al._ (2010); Xie  _et al._ (2012), whereas a direct
evidence of ultrafast electron dynamics in molecules remains a point of debate
Lépine  _et al._ (2014). In particular, there are contradictions between
recent experimental studies Calegari  _et al._ (2014); Lara-Astiaso  _et al._
(2018) claiming to have observed the ultrafast electronic processes in
molecules and theoretical investigations Arnold  _et al._ (2017); Vacher  _et
al._ (2017) performed on systems of similar complexity. The disagreement is
centered around the question on how strong is the influence of the slow
nuclear motion on the dynamics of electronic density. Extensive ab initio
calculations for several molecules Arnold  _et al._ (2017); Vacher  _et al._
(2017) demonstrated that the nuclear dynamics leads to the decoherence of the
electronic wave packet on the time scale of a few femtoseconds which can make
experimental observations of the electronic motion impossible. At the same
time, long-lasting electronic coherences were reported for the ionized
propiolic acid Despré  _et al._ (2018) and iodoacetylene Jia  _et al._ (2019),
suggesting that the influence of nuclear motion on the electronic dynamics is
very case-specific and requires careful investigation.

Understanding the interplay between the nuclear rearrangement and ultrafast
electronic motion requires a concerted description of the electron-nuclear
dynamics. Being one of the most powerful approaches for this purpose, the
multi-configurational time-dependent Hartree (MCTDH) method Meyer  _et al._
(1990); Raab  _et al._ (1999); Beck  _et al._ (2000) was recently applied for
describing electronic coherence Arnold  _et al._ (2017); Despré  _et al._
(2018). Although this rigorous technique takes into account all quantum
effects, such as tunneling and non-adiabatic transitions, it suffers from
exponential scaling and also requires constructing global potential energy
surfaces (PESs).

An alternative strategy for simulating coupled electron-nuclear dynamics
employs a trajectory-guided Gaussian basis to represent the evolving
wavepacket and an “on-the-fly” evaluation of the electronic structure. These
“direct dynamics” approaches calculate the PESs along trajectories only, thus
avoiding the precomputation of globally fitted surfaces, and sample only the
relevant regions of the configuration space. Among these methods, the closest
in spirit to MCTDH are the variational multi-configurational Gaussians (vMCG)
Burghardt  _et al._ (1999); Worth and Burghardt (2003); Worth  _et al._
(2008); Šulc  _et al._ (2013), but many others exist, ranging from ab initio
multiple spawning Curchod and Martínez (2018) and other Gaussian basis methods
Shalashilin and Child (2004); Saita and Shalashilin (2012) to more approximate
mixed quantum-classical Tully (1990); Li  _et al._ (2005) and semiclassical
Herman and Kluk (1984); Tatchen and Pollak (2009); Ceotto  _et al._ (2009)
approaches.

Treating electronic coherence with trajectory-based techniques was pioneered
by Bearpark, Robb and co-workers Vacher  _et al._ (2015a, b); Jenkins  _et
al._ (2016a, b) who propagated multiple trajectories representing an initially
delocalized wave packet using Ehrenfest approximation. Although this technique
captures decoherence due to a superposition of coherent oscillations with
different frequencies appearing at the respective nuclear geometries, it
completely ignores the decoherence due to the quantum motion of the wave
packet resulting in the accumulation of phase along the propagated trajectory.
The latter mechanism is referred to as the phase jitter Fiete and Heller
(2003) in general literature on quantum decoherence or, more specifically, as
dephasing in the case of electron-nuclear processes. Allowing the wave packet
to evolve quantum mechanically, dephasing mechanism was taken into account and
the electronic coherence upon ionization of a system was simulated in several
molecules using the direct dynamics version of vMCG scheme Vacher  _et al._
(2017) and its Ehrenfest-based variant Jenkins  _et al._ (2018).

Here, we use the thawed Gaussian approximation (TGA) Heller (1975); Begušić
_et al._ (2018); Lasser and Lubich (2020), one of the simplest semiclassical
approaches for molecular dynamics, to evaluate the influence of nuclear motion
on the ultrafast electronic dynamics and to find a simple, yet detailed
mechanism of decoherence, which is not available in basis set methods such as
MCTDH or vMCG. Within the TGA, the nuclear wavefunction is described by a
single Gaussian wave packet whose center follows Hamilton’s equations of
motion and whose time-dependent width and phase are propagated using the local
harmonic approximation of the PES.

We simulate coupled electron-nuclear dynamics taking place after outer-valence
ionization of two polyatomic molecules: propiolic acid (HC2COOH) and its amide
derivative propiolamide (HC2CONH2). Propiolic acid provides us a perfect
system for validating the semiclassical TGA because the electronic coherences
in this molecule were recently calculated using a full quantum MCTDH approach
Despré  _et al._ (2018). The propiolamide molecule, in turn, is studied here
for the first time.

The starting point of our investigations is a neutral molecule in its ground
electronic and vibrational states. The ionization of the system performed by
the ultrashort laser pulse can bring the molecule to a non-stationary
superposition of ionic states, thus launching a coupled dynamics of electronic
and nuclear wave packets. A coherent superposition of multiple electronic
states triggers oscillations of the charge along a molecular chain. This
purely electronic mechanism was termed charge migration Cederbaum and Zobeley
(1999); Kuleff and Cederbaum (2014) to distinguish it from a more common
charge transfer driven by nuclei Sun  _et al._ (2017). Although the charge
migration is governed by the electronic motion, it is strongly coupled with
the nuclear dynamics and therefore can crucially affect the behavior of the
whole molecule.

Previous calculations of the ionization spectrum of propiolic acid showed
Golubev  _et al._ (2017); Despré  _et al._ (2018) that, due to the electron
correlation, the ground and the second excited ionic states of the molecule
are a strong mixture of two one-hole configurations: an electron missing in
the highest occupied molecular orbital (HOMO) and an electron missing in the
HOMO–2. Therefore, a sudden removal of an electron either from HOMO or from
HOMO–2 will create an electronic wave packet, which will initiate charge
migration oscillations between the carbon triple bond and the carbonyl oxygen
with a period of about 6.2 fs, determined by the energy gap between the first
and the third cationic states Golubev  _et al._ (2017). Due to the planar
geometry of the propiolic acid and a large energy gap between remaining ionic
states, the indicated superposition can be obtained in an experiment by an
appropriate orientation of the molecule with respect to the laser
polarization.

Figure 1: Electronic coherence measured by the time-dependent overlap
$\chi_{13}(t)$ of the nuclear wave packets propagated on the first and third
cationic states of propiolic acid after removal of a HOMO electron. Dynamics
was performed with the full quantum MCTDH method (“quantum”) or with the
semiclassical vertical-Hessian thawed Gaussian approximation (“semiclass.”),
and either with the vibronic coupling (“VC”) or on-the-fly (“OTF”)
Hamiltonian, both obtained at the ab initio many-body Green’s function ADC(3)
level of theory. All simulations except the “quantum VC” calculation employed
the diabatic approximation (“diab.”), which neglects nonadiabatic couplings
between the diabatic states, or the adiabatic approximation (“adiab.”), which
neglects nonadiabatic couplings between adiabatic states.

To describe the ionization process, we used the non-Dyson algebraic-
diagrammatic-construction (ADC) scheme Schirmer  _et al._ (1998) to represent
the one-particle Green’s function. We chose a rather computationally expensive
ADC method because it allows treating the ionization process explicitly
starting from the neutral state of a system. This important advantage over
other, more conventional electronic structure approaches allows us to estimate
populations of the resulting ionic states created after removal of an electron
from a molecule. [See Sec. I Supplementary Material (SM) for details.]

Within the sudden and Franck–Condon approximations, the ionization is modelled
by projecting the ground (electronic and nuclear) neutral state of a molecule
onto the ionic subspace of the system. After ionization, a single nuclear
Gaussian wave packet on each involved ionic surface is propagated
independently from the others (non-adiabatic effects are neglected). The
center of each Gaussian is guided by a single classical trajectory, while the
width and phase are propagated using the single-Hessian Begušić  _et al._
(2019) variant of the TGA.

Within the Born–Huang representation Born and Huang (1954) of the molecular
wavefunction, the expectation value of an electronic operator
$\hat{O}(\mathbf{r},\mathbf{R})$ can be expressed as

$\langle\hat{O}\rangle(t)=\sum_{i,j}\int\chi_{i}^{\ast}(\mathbf{R},t)O_{ij}(\mathbf{R})\chi_{j}(\mathbf{R},t)d\mathbf{R},$
(1)

where
$O_{ij}(\mathbf{R})=\int\Phi_{i}^{\ast}(\mathbf{r},\mathbf{R})\hat{O}(\mathbf{r},\mathbf{R})\Phi_{j}(\mathbf{r},\mathbf{R})d\mathbf{r}$
denote the $\mathbf{R}$-dependent matrix elements of the electronic operator
between electronic states $i$ and $j$, and the quantities
$\chi_{i}(\mathbf{R},t)$ are the time-dependent nuclear wave packets
propagated on the corresponding PESs.

If both the operator $\hat{O}(\mathbf{r},\mathbf{R})$ and the electronic
states $\\{\Phi_{i}(\mathbf{r},\mathbf{R})\\}$ depend weakly on nuclear
coordinates $\mathbf{R}$, Eq. (1) can be further simplified as

$\langle\hat{O}\rangle(t)\approx\sum_{i,j}O_{ij}\chi_{ij}(t),$ (2)

where the nuclear overlaps
$\chi_{ij}(t)=\int\chi_{i}^{\ast}(\mathbf{R},t)\chi_{j}(\mathbf{R},t)d\mathbf{R}$
represent the populations of electronic states when $i=j$ and the electronic
coherences Arnold  _et al._ (2017); Vacher  _et al._ (2017); Despré  _et al._
(2018) when $i\neq j$. Equation (2) shows that the factors $\chi_{ij}(t)$
provide the only source of time dependence in the expectation value of the
electronic operator and thus can serve as convenient properties to quantify
the decoherence time. Electronic coherence evaluated as the nuclear overlap
$\chi_{ij}(t)$ is a special case of _fidelity amplitude_ (also called
Loschmidt echo), a quantity measuring sensitivity of quantum dynamics to
perturbations and defined as the overlap at time $t$ of two wavepackets,
initially the same but propagated with different Hamiltonians Peres (1984);
Jalabert and Pastawski (2001); Vaníček (2006); Gorin  _et al._ (2006); Jacquod
and Petitjean (2009); Mollica and Vaníček (2011); Zambrano and Ozorio de
Almeida (2011); Goussev  _et al._ (2012).

We computed the electronic coherences $\chi_{ij}(t)$ generated after ionizing
an electron from the HOMO of the propiolic acid; the initial molecular wave
packet was an equally weighted and phase-synchronized superposition of the
first and third cationic states. Figure 1 shows the electronic coherence
evaluated by five different schemes. We adopted the vibronic-coupling (VC)
Hamiltonian from Ref. Despré  _et al._ (2018) to perform MCTDH simulations
taking into account all nuclear degrees of freedom. The full quantum-
mechanical calculations show (red solid line in Fig. 1) that the electronic
oscillations are strongly influenced by nuclear motion—the coherences are
completely suppressed within first 15 fs Despré  _et al._ (2018).

To validate the applicability of the TGA, we performed semiclassical
calculations using adiabatic version of the VC Hamiltonian, where the PESs
were obtained by diagonalizing the four-state VC model used in the MCTDH
calculations. Our simulation shows (blue dashed line in Fig. 1) that on the
short time scale the TGA gives results almost identical to the full quantum
MCTDH calculations. The small deviations start to appear at longer times due
to the nonadiabatic effects, neglected in the TGA. Note that TGA is exact
within the VC model when the diabatic PESs are not coupled to each other, and
therefore provides results identical to the MCTDH simulations performed on
such VC Hamiltonian (yellow solid line in Fig. 1).

Figure 2: Geometric interpretation of Eq. (5): The reduced action
$S_{\text{red}}$ is the gray phase-space area enclosed by the curve
$\mathbf{C}$ consisting of the classical trajectory propagated in state $i$
with potential energy $V_{i}$ and connecting $\mathbf{X}_{i}(t)$ with
$\mathbf{X}_{0}$ (blue curve), classical trajectory propagated in state $j$
with potential energy $V_{j}$ and connecting $\mathbf{X}_{0}$ with
$\mathbf{X}_{j}(t)$ (red curve), and a straight (black) line connecting
$\mathbf{X}_{j}(t)$ with $\mathbf{X}_{i}(t)$. The correct sign of
$S_{\text{red}}$ is obtained by taking the curve integral along $\mathbf{C}$
in the direction indicated by the arrows. The phase space distance betwen the
final points $\mathbf{X}_{j}(t)$ and $\mathbf{X}_{i}(t)$, given by the length
$d$ of the black line segment, determines the coherence decay [i.e., the decay
of $|\chi_{ij}(t)|$], while the reduced action $S_{\text{red}}$, equal to the
gray area, affects the frequency of oscillations of the coherence
$\chi_{ij}(t)$.

We also performed semiclassical calculations with on-the-fly evaluation of the
electronic structure at the same ab initio level of theory as that used in the
construction of the VC Hamiltonian (green dash-dotted line in Fig. 1). In this
case, the wave packet can potentially evolve beyond a simple model used for
fitting PESs. In particular, VC Hamiltonian typically uses a rather primitive
approximation of PESs for nuclear configurations formed by superposition of
normal modes. Allowing the wave packet to evolve according to the exact
Hamiltonian computed on the fly makes it possible to visit nuclear regions
inaccessible within the VC Hamiltonian and thus to take anharmonicity effects
into account (see Sec. II of SM for details). This is reflected in our on-the-
fly calculations, which predict the electronic motion with a similar
oscillation period, but a slightly faster decay of the electronic coherence
than within the VC model. Remarkably, because the effect of using the on-the-
fly potential is much larger than the effect of including the nonadiabatic
couplings, the semiclassical on-the-fly result of the TGA is most likely more
accurate than the quantum result of the MCTDH calculation with the VC model!

In Sec. III of SM, we derive an analytical expression

$\chi_{ij}(t)=e^{-d(t)^{2}/4\hbar}e^{iS(t)/\hbar}$ (3)

for the semiclassically evaluated coherence in case the two Gaussians have
fixed widths. Here $d(t)$ is phase-space distance (in mass- and frequency-
scaled phase-space coordinates $\mathbf{R}$ and $\mathbf{P}$) between the
centers of the two Gaussian wavepackets at time $t$,

$S(t)=S_{\text{red}}(t)-\Delta E\,t,$ (4)

is the classical action, $\Delta
E=V_{j}(\mathbf{R}_{0})-V_{i}(\mathbf{R}_{0})$ is the energy gap at the
initial point $\mathbf{R}_{0}$, and

$S_{\text{red}}(t)=\oint_{\mathbf{C}}\mathbf{P}^{T}\cdot d\mathbf{R}$ (5)

is the reduced action equal to the signed area within the closed curve
$\mathbf{C}$ shown in Fig. 2 Zambrano and Ozorio de Almeida (2011).

Figure 3: Semiclassical analysis of the electronic coherence $\chi_{13}(t)$
from Fig. 1. Comparison of the coherence computed with the on-the-fly version
of TGA (green dash-dotted line), analytical semiclassical expression (3) (red
solid line), factors describing the decay (blue solid line) and oscillations
(gray solid line) of coherence in the presence of nuclear motion, and the
undamped coherence in the absence of nuclear motion (gray dashed line).

The analytical expression (3) provides a simple, semiclassical interpretation
of the effect of nuclear dynamics on electronic coherence (see Fig. 3): The
diverging nuclear trajectories affect not only the absolute value of
$\chi_{ij}(t)$, which, as expected, decays as a Gaussian function of the
phase-space distance $d(t)$, but also frequency of electronic oscillations. In
the absence of nuclear motion, the electronic coherence would oscillate with
frequency $\Delta E/\hbar$, but now, due to nuclear dynamics, the phase of
electronic oscillations at time $t$ is modified by the area
$S_{\text{red}}(t)$ divided by $\hbar$ [see Eqs. (3), (4), and (5)]. It is
easy to see that if the potential energy surfaces are simply vertically
shifted, i.e., if $V_{j}=V_{i}+\Delta E$, then $\mathbf{X}_{i}=\mathbf{X}_{j}$
and $d=S_{\text{red}}=0$, implying that the electronic coherence
$\chi_{ij}(t)$ is not affected by nuclear dynamics.

Let us turn to the electron-nuclear dynamics driven by the ionization of the
propiolamide. Similar to the spectrum of the propiolic acid Golubev  _et al._
(2017), in the energy range 10–14 eV only the four states shown in Fig. 2 in
Sec. IV of the SM are present. The strong electron correlation between valence
orbitals in the neutral propiolamide leads to appearance of the almost equal
in weights one-hole configurations in the ionic states. Therefore, an
ultrashort (sudden) ionization of the molecule will inevitably create an
electronic wave packet and trigger dynamics of the electron density between
the carbon triple bond and amide moiety. The molecule is planar and thus
belongs to the $C_{s}$ symmetry group which allows assignment of the ionic
states to two irreducible representations: the second and fourth states belong
to the A′, while the first and third states to A′′. As for propiolic acid,
orientation of the propiolamide with respect to polarization of the ionizing
laser field can be used to populate only the states of interest.

Figure 4: Top panel: Time evolution of the electronic coherences
$\chi_{ij}(t)$ created by the equally weighted coherent superposition of the
second and fourth cationic states (ionization from the HOMO–1, A′ symmetry),
and the first and third cationic states (ionization from the HOMO, A′′
symmetry). Bottom panel: Comparison of the electronic purity function
$\text{Tr}(\rho^{2})$ for the propiolamide and propiolic acid molecules.

Figure 4 shows the evolution of the electronic coherence in the propiolamide.
Taking advantage of the molecular symmetry, we simulate dynamics occurring
after removal of the HOMO (populating the first and third cationic states) and
HOMO–1 (populating the second and fourth cationic states) electrons. Starting
from the equally weighted and phase-matched superposition of the electronic
states, in both cases the nuclear motion perturbs electronic oscillations and
leads to decoherence. Coherence times of different molecules can be compared
using the purity $\operatorname{Tr}(\rho(t)^{2})$ Arnold  _et al._ (2017),
where the electronic density matrix $\rho(t)$ is related to the matrix of
nuclear overlaps from Eq. (2) by transposition: $\rho_{ij}(t)=\chi_{ji}(t)$.
Due to decoherence, the purity decays from the value
$\text{Tr}(\rho(0)^{2})=1$ for the initially pure state to the value $1/n$ for
the equally weighted mixture of $n$ states.

Our simulations demonstrate that, contrary to the propiolic acid molecule, for
which long-lasting coherences were observed (see bottom panel of Fig. 4 and
also Ref. Despré  _et al._ (2018)), the initially pure superpositions in the
propiolamide evolve to mixed states in just a few femtoseconds. Importantly,
the energy gaps between the involved electronic states of propiolamide are
larger than those for the propiolic acid (see Fig. 2 of SM and Ref. Golubev
_et al._ (2017)), which leads to faster oscillations of electronic density
along the molecular chain. Despite the rather short coherence time, due to the
faster charge migration the electronic density in the propiolamide has enough
time to perform one clear oscillation (see Fig. 4 above and Fig. 3 of SM).
Moreover, the existence of the strong hole-mixing in both symmetries of the
propiolamide can be used to induce oscillations of the charge along different
directions in the molecule. Dependence of the charge migration on the
molecular orientation provides an important advantage for experimental
measurements utilizing the time-resolved high-harmonic generation (HHG)
spectroscopy employed recently by Wörner and co-workers Kraus  _et al._
(2015). Alignment of the molecule with respect to the pump pulse should be
reflected in the resulting HHG spectra and thus can be used as a direct
evidence of the ultrafast electron dynamics.

Although trajectory-based direct dynamics methods were previously used to
estimate electronic coherences in various polyatomic molecules Vacher  _et
al._ (2017); Jenkins  _et al._ (2018), converged results typically required a
large number of trajectories. The vMCG used in these studies take into account
non-adiabatic transitions and tunneling effects, but require solving rather
complicated equations of motion and make interpreting obtained results less
intuitive than the simple picture provided here in Eq. (3) and Figs. 2 and 3.
Our approach based on the TGA can be viewed as a very special case of far more
general vMCG—namely, it can be classified as a single-Gaussian, non-
variational, multi-set, non-frozen, adiabatic, single-Hessian version of vMCG.
In contrast to the application of vMCG in Refs. Vacher  _et al._ (2017);
Jenkins  _et al._ (2018), our implementation of the TGA uses a multi-set
approach, where a single Gaussian function with relaxed parameters is used for
every involved electronic state. By approximately taking into account quantum
properties of the wave packet, such as its width and phase, the TGA captures
the dephasing mechanism while maintaining sufficient accuracy, especially at
short time scales. A detailed comparison of the TGA and vMCG is provided in
Sec. V of SM—remarkably, the TGA, which uses only a single classical
trajectory per electronic state, yields in the propiolic acid better results
than the single-set version of vMCG with 31 variational trajectories (see Fig.
4 in SM). While the multi-set version of vMCG gives similar results to TGA
also with 1 Gaussian, no improvement is seen by using 8 Gaussians (Fig. 4 of
SM). In the TGA, using only two classical trajectories was not only sufficient
but also crucial for revealing the simple physical mechanism of decoherence.

The semiclassical vertical-Hessians TGA used in this paper can be further
improved by calculating Hessians along the propagated trajectory and thus
taking into account more complicated situations, e.g., dissociation of a
molecule. Extensions of the TGA, such as the extended thawed Gaussian
approximation Lee and Heller (1982); Patoz  _et al._ (2018) or Hagedorn
wavepackets Hagedorn (1998); Lasser and Lubich (2020), which propagate a
Gaussian multiplied by a linear or general polynomial, can make on-the-fly
semiclassical simulations even more accurate.

In conclusion, we implemented a simple and efficient on-the-fly semiclassical
approach to understand the effects of nuclear motion on electronic coherence
in molecules. Although the propiolic acid and propiolamide have very similar
ionization spectra, our calculations predict that their electronic coherence
times differ substantially. The simple method was validated by comparison with
the full-dimensional quantum calculations performed using the MCTDH and vMCG
methods. As suggested by Fig. 1, neglecting the nonadiabatic couplings by the
TGA may not be a severe approximation because, even in systems with strong
couplings, the nonadiabatic effects typically start playing a significant role
only at times longer than the ultrashort decoherence time scale. If one
suspects an exceptional importance of nonadiabatic effects even at times
before the electronic coherence has decayed to zero, it is possible to
validate the applicability of the TGA without performing expensive quantum
simulations by verifying adiabaticity with on-the-fly semiclassical
calculations based, e.g., on surface hopping Tully (1990) or multiple-surface
dephasing representation Zimmermann and Vaníček (2010, 2012, 2012); Prlj  _et
al._ (2020).

Despite its limitations, the presented technique for evaluating coherence can
help breaking the “curse of dimensionality” appearing in the quantum treatment
of large molecules, which can be crucial for full-dimensional simulations of
ultrafast electronic processes in biologically relevant systems Calegari  _et
al._ (2016). Being able to treat molecules with a few hundred atoms, this
technique can help shed light on the continuing debates on the role of quantum
coherence in biology Lambert  _et al._ (2013), quickly preselect molecules
suitable for further experimental investigations, and support theoretically
recent experimental observations of attosecond electron dynamics in realistic
molecular systems.

###### Acknowledgements.

Authors thank Alexander Kuleff for many valuable discussions and the Swiss
National Science Foundation for financial support through the NCCR MUST
(Molecular Ultrafast Science and Technology) Network. N. V. G. acknowledges
the support by the Branco Weiss Fellowship—Society in Science, administered by
the ETH Zürich.

## References

  * Krausz and Ivanov (2009) F. Krausz and M. Ivanov, Rev. Mod. Phys. 81, 163 (2009).
  * Young  _et al._ (2018) L. Young, K. Ueda, M. Gühr, P. H. Bucksbaum, M. Simon, S. Mukamel, N. Rohringer, K. C. Prince, C. Masciovecchio, M. Meyer, A. Rudenko, D. Rolles, C. Bostedt, M. Fuchs, D. A. Reis, R. Santra, H. Kapteyn, M. Murnane, H. Ibrahim, F. Légaré, M. Vrakking, M. Isinger, D. Kroon, M. Gisselbrecht, A. L’Huillier, H. J. Wörner, and S. R. Leone, J. Phys. B: At., Mol. Opt. Phys. 51, 032003 (2018).
  * Mauritsson  _et al._ (2018) J. Mauritsson, G. Vampa, and C. Vozzi, J. Opt. 20, 110201 (2018).
  * Krausz (2016) F. Krausz, Phys. Scr. 91, 063011 (2016).
  * Kraus  _et al._ (2015) P. M. Kraus, B. Mignolet, D. Baykusheva, A. Rupenyan, L. Horny, E. F. Penka, G. Grassi, O. I. Tolstikhin, J. Schneider, F. Jensen, L. B. Madsen, A. D. Bandrauk, F. Remacle, and H. J. Wörner, Science 350, 790 (2015).
  * Goulielmakis  _et al._ (2010) E. Goulielmakis, Z.-H. Loh, A. Wirth, R. Santra, N. Rohringer, V. S. Yakovlev, S. Zherebtsov, T. Pfeifer, A. M. Azzeer, M. F. Kling, S. R. Leone, and F. Krausz, Nature 466, 739 (2010).
  * Xie  _et al._ (2012) X. Xie, S. Roither, D. Kartashov, E. Persson, D. G. Arbó, L. Zhang, S. Gräfe, M. S. Schöffler, J. Burgdörfer, A. Baltuška, and M. Kitzler, Phys. Rev. Lett. 108, 193004 (2012).
  * Lépine  _et al._ (2014) F. Lépine, M. Y. Ivanov, and M. J. J. Vrakking, Nat. Photonics 8, 195 (2014).
  * Calegari  _et al._ (2014) F. Calegari, D. Ayuso, A. Trabattoni, L. Belshaw, S. De Camillis, S. Anumula, F. Frassetto, L. Poletto, A. Palacios, P. Decleva, J. B. Greenwood, F. Martin, and M. Nisoli, Science 346, 336 (2014).
  * Lara-Astiaso  _et al._ (2018) M. Lara-Astiaso, M. Galli, A. Trabattoni, A. Palacios, D. Ayuso, F. Frassetto, L. Poletto, S. De Camillis, J. Greenwood, P. Decleva, I. Tavernelli, F. Calegari, M. Nisoli, and F. Martín, J. Phys. Chem. Lett. 9, 4570 (2018).
  * Arnold  _et al._ (2017) C. Arnold, O. Vendrell, and R. Santra, Phys. Rev. A 95, 033425 (2017).
  * Vacher  _et al._ (2017) M. Vacher, M. J. Bearpark, M. A. Robb, and J. a. P. Malhado, Phys. Rev. Lett. 118, 083001 (2017).
  * Despré  _et al._ (2018) V. Despré, N. V. Golubev, and A. I. Kuleff, Phys. Rev. Lett. 121, 203002 (2018).
  * Jia  _et al._ (2019) D. Jia, J. Manz, and Y. Yang, J. Phys. Chem. Lett. 10, 4273 (2019).
  * Meyer  _et al._ (1990) H.-D. Meyer, U. Manthe, and L. S. Cederbaum, Chem. Phys. Lett. 165, 73 (1990).
  * Raab  _et al._ (1999) A. Raab, G. A. Worth, H.-D. Meyer, and L. S. Cederbaum, J. Chem. Phys. 110, 936 (1999).
  * Beck  _et al._ (2000) M. Beck, A. Jäckle, G. Worth, and H.-D. Meyer, Phys. Rep. 324, 1 (2000).
  * Burghardt  _et al._ (1999) I. Burghardt, H.-D. Meyer, and L. S. Cederbaum, J. Chem. Phys. 111, 2927 (1999).
  * Worth and Burghardt (2003) G. A. Worth and I. Burghardt, Chem. Phys. Lett. 368, 502 (2003).
  * Worth  _et al._ (2008) G. Worth, M. Robb, and B. Lasorne, Mol. Phys. 106, 2077 (2008).
  * Šulc  _et al._ (2013) M. Šulc, H. Hernández, T. J. Martínez, and J. Vaníček, J. Chem. Phys. 139, 034112 (2013).
  * Curchod and Martínez (2018) B. F. E. Curchod and T. J. Martínez, Chem. Rev. 118, 3305 (2018).
  * Shalashilin and Child (2004) D. V. Shalashilin and M. S. Child, J. Chem. Phys. 121, 3563 (2004).
  * Saita and Shalashilin (2012) K. Saita and D. V. Shalashilin, J. Chem. Phys. 137, 22A506 (2012).
  * Tully (1990) J. C. Tully, J. Chem. Phys. 93, 1061 (1990).
  * Li  _et al._ (2005) X. Li, J. C. Tully, H. B. Schlegel, and M. J. Frisch, J. Chem. Phys. 123, 084106 (2005).
  * Herman and Kluk (1984) M. F. Herman and E. Kluk, Chem. Phys. 91, 27 (1984).
  * Tatchen and Pollak (2009) J. Tatchen and E. Pollak, J. Chem. Phys. 130, 041103 (2009).
  * Ceotto  _et al._ (2009) M. Ceotto, S. Atahan, G. F. Tantardini, and A. Aspuru-Guzik, J. Chem. Phys. 130, 234113 (2009).
  * Vacher  _et al._ (2015a) M. Vacher, D. Mendive-Tapia, M. J. Bearpark, and M. A. Robb, J. Chem. Phys. 142, 094105 (2015a).
  * Vacher  _et al._ (2015b) M. Vacher, L. Steinberg, A. J. Jenkins, M. J. Bearpark, and M. A. Robb, Phys. Rev. A 92, 040502 (2015b).
  * Jenkins  _et al._ (2016a) A. J. Jenkins, M. Vacher, M. J. Bearpark, and M. A. Robb, J. Chem. Phys. 144, 104110 (2016a).
  * Jenkins  _et al._ (2016b) A. J. Jenkins, M. Vacher, R. M. Twidale, M. J. Bearpark, and M. A. Robb, J. Chem. Phys. 145, 164103 (2016b).
  * Fiete and Heller (2003) G. A. Fiete and E. J. Heller, Phys. Rev. A 68, 022112 (2003).
  * Jenkins  _et al._ (2018) A. J. Jenkins, K. E. Spinlove, M. Vacher, G. A. Worth, and M. A. Robb, J. Chem. Phys. 149, 094108 (2018).
  * Heller (1975) E. J. Heller, J. Chem. Phys. 62, 1544 (1975).
  * Begušić  _et al._ (2018) T. Begušić, J. Roulet, and J. Vaníček, J. Chem. Phys. 149, 244115 (2018).
  * Lasser and Lubich (2020) C. Lasser and C. Lubich,  (2020), arXiv:2002.00624 [math.NA] .
  * Cederbaum and Zobeley (1999) L. Cederbaum and J. Zobeley, Chem. Phys. Lett. 307, 205 (1999).
  * Kuleff and Cederbaum (2014) A. I. Kuleff and L. S. Cederbaum, J. Phys. B: At., Mol. Opt. Phys. 47, 124002 (2014).
  * Sun  _et al._ (2017) S. Sun, B. Mignolet, L. Fan, W. Li, R. D. Levine, and F. Remacle, J. Phys. Chem. A 121, 1442 (2017).
  * Golubev  _et al._ (2017) N. V. Golubev, V. Despré, and A. I. Kuleff, J. Mod. Opt. 64, 1031 (2017).
  * Schirmer  _et al._ (1998) J. Schirmer, A. B. Trofimov, and G. Stelter, J. Chem. Phys. 109, 4734 (1998).
  * Begušić  _et al._ (2019) T. Begušić, M. Cordova, and J. Vaníček, J. Chem. Phys. 150, 154117 (2019).
  * Born and Huang (1954) M. Born and K. Huang,  _Dynamical Theory of Crystal Lattices_ , Oxford Classic Texts in the Physical Sciences (Oxford University Press, Oxford, 1954).
  * Peres (1984) A. Peres, Phys. Rev. A 30, 1610 (1984).
  * Jalabert and Pastawski (2001) R. A. Jalabert and H. M. Pastawski, Phys. Rev. Lett. 86, 2490 (2001).
  * Vaníček (2006) J. Vaníček, Phys. Rev. E 73, 046204 (2006).
  * Gorin  _et al._ (2006) T. Gorin, T. Prosen, T. H. Seligman, and M. Žnidarič, Phys. Rep. 435, 33 (2006).
  * Jacquod and Petitjean (2009) P. Jacquod and C. Petitjean, Adv. Phys. 58, 67 (2009).
  * Mollica and Vaníček (2011) C. Mollica and J. Vaníček, Phys. Rev. Lett. 107, 214101 (2011).
  * Zambrano and Ozorio de Almeida (2011) E. Zambrano and A. M. Ozorio de Almeida, Phys. Rev. E 84, 045201(R) (2011).
  * Goussev  _et al._ (2012) A. Goussev, R. A. Jalabert, H. M. Pastawski, and D. Wisniacki, Scholarpedia 7, 11687 (2012).
  * Lee and Heller (1982) S.-Y. Lee and E. J. Heller, J. Chem. Phys. 76, 3035 (1982).
  * Patoz  _et al._ (2018) A. Patoz, T. Begušić, and J. Vaníček, J. Phys. Chem. Lett. 9, 2367 (2018).
  * Hagedorn (1998) G. A. Hagedorn, Ann. Phys. (N. Y.) 269, 77 (1998).
  * Zimmermann and Vaníček (2010) T. Zimmermann and J. Vaníček, J. Chem. Phys. 132, 241101 (2010).
  * Zimmermann and Vaníček (2012) T. Zimmermann and J. Vaníček, J. Chem. Phys. 136, 094106 (2012).
  * Zimmermann and Vaníček (2012) T. Zimmermann and J. Vaníček, J. Chem. Phys. 137, 22A516 (2012).
  * Prlj  _et al._ (2020) A. Prlj, T. Begušić, Z. T. Zhang, G. C. Fish, M. Wehrle, T. Zimmermann, S. Choi, J. Roulet, J.-E. Moser, and J. Vaníček, J. Chem. Theory Comput. 16, 2617 (2020).
  * Calegari  _et al._ (2016) F. Calegari, A. Trabattoni, A. Palacios, D. Ayuso, M. C. Castrovilli, J. B. Greenwood, P. Decleva, F. Martín, and M. Nisoli, J. Phys. B: At., Mol. Opt. Phys. 49, 142001 (2016).
  * Lambert  _et al._ (2013) N. Lambert, Y.-N. Chen, Y.-C. Cheng, C.-M. Li, G.-Y. Chen, and F. Nori, Nat. Phys. 9, 10 (2013).


# Riemannian Geometry Approach for Minimizing Distortion and its Applications

Dror Ozeri
Rafael–Advanced Defense Systems Israel
email<EMAIL_ADDRESS>
###### Abstract

Given an affine transformation $T$, we define its Fisher distortion
$Dist_{F}(T)$. We show that the Fisher distortion has Riemannian metric
structure and provide an algorithm for finding mean distorting transformation
- namely - for a given set $\\{T_{i}\\}_{i=1}^{N}$ of affine transformations,
find an affine transformation $T$ that minimize the overall distortion
$\sum_{i=1}^{N}Dist_{F}^{2}(T^{-1}T_{i}).$ The mean distorting transformation
can be useful in some fields -in particular we apply it for rendering affine
panoramas.

##  1 Introduction

Registration algorithms between different objects are common in computer
vision and engineering. Most important usages are planar panoramas [6], point
cloud and point set registration [17, 10].

A common practice of the registration process is choosing one of the objects
to be the reference object (i.e., the identity transformation). However,
unlike rigid transformation, non-rigid transformation can visually distort the
original object; choosing an arbitrary or wrong reference object may lead to
significant distortions.

Figure 1:  The space of lower triangular matrices with positive diagonal
$L_{2}^{+}$. The Identity matrix $e$ is the origin and does not distort the
image. The matrix $h$ is located somewhere on the manifold and distorts the
image. The matrix $g$ is located further on the manifold and further distorts
the image. The distortion is proportional to the distance symbolized by the
red curve (see Theorem 2).

(a)

(b)

Figure 2:  Tissot’s indicatrix of the Mercator projection (a) and The Wagner
IV projection (b) . The Mercator projection preserves angles but distorts
areas, the Wagner IV projection preserves areas but distorts angles.
Projection images were taken from [12]

###  1.1 Related Work

####  1.1.1 Choosing the Best Reference on Planar Panorama

The problem of estimating the best reference image in a panorama is
surprisingly understudied. In [1], Fabio Bellavia et al. investigated this
issue with a different distortion measure. In [7], Capel suggests taking the
central image in the panorama. Indeed in many implementations of panorama this
method is chosen (for example: OpenCV implementation [4]).

In this study, we show that a natural geometric notion of distortion can be
defined and minimized, yielding a new reference plane that is not necessarily
one of the initial planes. Moreover, the method is indifferent to the
dimension of the space involved and can be applied not only for images. We
focus only on the space of affine transformations.

As far as we know, reducing distortion using this method is a novel method and
was never studied before.

####  1.1.2 Map Projection Distortion

For natural reasons, the concept of distortion has roots in the field of
cartography (see [14]). First Introduced by Nicolas Auguste Tissot as Tissot’s
indicatrix (ellipse of distortion) in 1878 [18] (see Fig 2). Airy, Kavrayskiy
[11, 13] and others developed several distortion measures; one of those
criteria is the popular Airy-Kavrayskiy criterion:

$E_{AK}=\sqrt{\int\epsilon_{AK}^{2}d\sigma}=\sqrt{\int[log^{2}{\frac{a}{b}}+log^{2}{ab}]d\sigma}=\sqrt{\int[log^{2}{a}+log^{2}{b}]d\sigma},$
(1)

where a $a$ and $b$ are the infinitesimal linear scales of the ellipse. The
first term $\epsilon_{a}^{2}=\log^{2}{\frac{a}{b}}$ measures the angular
distortion and the second term $\epsilon_{p}^{2}=\log^{2}{ab}$ measures the
area distortion (see [14]). Laskowski [15] has linked between $a$ and $b$ and
the singular values of the local Jacobi matrix of the projection, yielding
precisely the Fisher distortion that we define later in this paper. In the
following sections, we will present an alternative geometric derivation of the
Fisher distortion.

##  2 Preliminaries and Definitions

###  2.1 Riemannian manifolds and the space of positive definite matrices

The space of positive definite matrices $P_{n}^{+}$, has been the focus of
much research in recent years. It has a rich Riemannian structure. Interesting
metrics on $P_{n}^{+}$ include Euclidean, Inverse Euclidean, Cholesky-
Euclidean, log-Cholesky [16] and the Fisher metric:

$d_{F}(A,B)=\sqrt{||log(A^{-\frac{1}{2}}BA^{-\frac{1}{2}})||^{2}}=\sqrt{\sum_{i=1}^{N}\log^{2}(\lambda_{i}(A^{-1}B))},$
(2)

where $\lambda_{i}$ are the matrix eigenvalues [2]. The Fisher metric is
fundamental, holds a variety of important invariant properties and will be the
main focus of this paper.

###  2.2 Fréchet mean

The Fréchet mean on a metric space is a generalization of the arithmetic mean
(or center of mass) in Euclidean space. It is well known that for
$\\{x_{i}\in\mathbb{R}^{n}\\}_{i=1}^{N}$ the arithmetic mean of the set
$\overline{x}=\frac{\sum_{i=1}^{N}x_{i}}{N}$ minimize the expression
$\sum_{i=1}^{N}||x-x_{i}||^{2}$. A Fréchet mean of a set of points
$\\{\nu_{i}\in M\\}_{i=1}^{N}$ in a metric space $M$ is an element of the
metric space that globally minimize the expression
$\sum_{i=1}^{N}d^{2}(x,\nu_{i})$, where $d$ is the metric. Fréchet means are
used in a variety of scientific and mathematical applications.

The space $P_{n}^{+}$ endowed with the Fisher metric holds an important
property relevant for our case - Every finite subset of $P_{n}^{+}$ has a
unique Fréchet mean (see [3]).

###  2.3 Cholesky Decomposition

It is known that every positive definite matrix $P$ has a decomposition of
$P=LL^{t},$ where $L$ is a lower triangular matrix; This is called the
Cholesky decomposition of the matrix. This decomposition is unique if we
restrict the diagonal elements of $L$ to be positive, namely the Lie group
$L_{n}^{+}$. Furthermore, the map $\Phi:x\to xx^{t}$ is a diffeomorphism
between the group $L_{n}^{+}$ and the $P_{n}^{+}$ as Riemannian manifolds. We
call the map $\Phi$ \- the Cholesky map, and $\Phi^{-1}(P)$ is called the
Cholesky factor of $P$.

###  2.4 Pullback of a metric

A common notion in Riemannian geometry is a pullback of a metric. The pullback
of a metric is a new metric defined on the preimage of a map. In particular we
are interested in the pullback of the Fisher metric into the space $L_{n}^{+}$
using the Cholesky map:

$\Phi^{*}d_{F}(l_{i},l_{j})=d_{F}(\Phi(l_{i}),\Phi(l_{j})).$ (3)

###  2.5 QR Decomposition

Throughout this paper we will use the following version of the QR
decomposition:

######  Theorem 1.

Any matrix $A\in GL_{n}(\mathbb{R})$ can be decomposed uniquely to $A=LQ$
where $L\in L_{n}^{+}$ and $Q\in O_{n}$ and $O_{n}$ is the set of orthogonal
transformation of $\mathbb{R}^{n}$.

###  2.6 Definition of the Fisher Distortion

Let $f(\textbf{x}):\mathbb{R}^{n}\to\mathbb{R}^{n}$ be an affine
transformation, then $f(\textbf{x})=A\textbf{x}+\textbf{b}$, for some $A\in
GL_{n}(\mathbb{R})$ and $\textbf{b}\in\mathbb{R}^{n}$. The Fisher distortion
of $f$ is

$Dist_{F}(f)=Dist_{F}(A)=\sqrt{\sum_{i=1}^{N}\log^{2}(\sigma_{i})},$ (4)

where $\sigma_{i}$ are $A$ singular values.

##  3 Mean Distorting Transformation Problem

The mean distorting transformation [MDT] problem is the following: Let
$\\{A_{i}\in GL_{n}(\mathbb{R})\\}_{i=1}^{N}$. Find
$\Lambda(A_{1},A_{2},...,A_{n})\in GL_{n}(\mathbb{R})$ that satisfies:

$\Lambda(A_{1},A_{2},...,A_{N})=\operatorname*{argmin}_{A\in
GL_{n}(\mathbb{R})}\sum_{i=1}^{N}Dist_{F}^{2}(A^{-1}A_{i}).$ (5)

.

In section 3.3 we will show that $\Lambda(A_{1},A_{2},...,A_{N})$ is unique up
to multiplication by the right with orthogonal transformation.

###  3.1 Motivation for the Fisher Distortion

Besides being studied in the map projection domain, two main properties make
the Fisher distortion a natural candidate for the Mean Distorting
Transformation Problem.
1\. Rigid Invariant \- The Fisher Distortion is invariant under left and right
multiplication by orthogonal transformations. i.e

$Dist_{F}(q_{1}\circ T\circ q_{2})=Dist_{F}(T),\forall q_{1},q_{2}\in O_{n}.$
(6)

This property is simply a consequence of the invariance of the singular values
under orthogonal matrix multiplication.
2\. Geometric Mean Property \- Let $\\{D_{i}\in D_{n}^{+}\\}_{i=1}^{N}$ be a
set of diagonal matrices with positive entries. Then it is easy to show that
the MDT of this set is simply the geometric mean of the matrices:

$T=\left(\prod_{i=1}^{N}D_{i}\right)^{\frac{1}{N}}.$ (7)

.

###  3.2 Relation between the Fisher metric and Fisher Distortion

Our main result is the following Theorem:

######  Theorem 2.

Let $\Phi^{*}d_{F}$ be the pullback metric on $L_{n}^{+}$ of the Fisher metric
$d_{F}$ defined on $P_{n}^{+}$, then

$\ \Phi^{*}d_{F}(l_{i},l_{j})=2Dist_{F}(l_{i}^{-1}l_{j})\;\forall
l_{i},l_{j}\in L_{n}^{+}.$ (8)

###### Proof.

First we show that $\Phi^{*}d_{F}$ is left invariant metric:
Let $l_{k}\in L_{n}^{+}$ , then by definition of $\Phi^{*}d_{F}$:

$\displaystyle\Phi^{*}d_{F}(l_{k}l_{i},l_{k}l_{j})=d_{F}(\Phi(l_{k}l_{i}),\Phi(l_{k}l_{j})).$
(9)

But,

$\displaystyle
d_{F}(\Phi(l_{k}l_{i}),\Phi(l_{k}l_{j}))=\sqrt{\sum_{m=1}^{N}\log^{2}(\lambda_{m}(\Phi(l_{k}l_{i})^{-1}\Phi(l_{k}l_{j}))}.$
(10)

Observe that

$\displaystyle\Phi(l_{k}l_{i})^{-1}\Phi(l_{k}l_{j})=(l_{k}^{t})^{-1}\Phi(l_{i})^{-1}\Phi(l_{j})l_{k}^{t},$
(11)

the matrices similarity means that both admit the same eigenvalues and by (2)
we have

$\Phi^{*}d_{F}(l_{k}l_{i},l_{k}l_{j})=\Phi^{*}d_{F}(l_{i},l_{j}),$ (12)

which means the metric $\Phi^{*}d_{F}$ is left invariant.
Now,

$\displaystyle\Phi^{*}d_{F}(l_{i},l_{j})=\Phi^{*}d_{F}(I,l_{i}^{-1}l_{j})=d_{F}(I,\Phi(l_{i}^{-1}l_{j}))=$
(13)
$\displaystyle\sqrt{\sum_{m=1}^{N}\log^{2}(\lambda_{m}(\Phi(l_{i}^{-1}l_{j})))}.$

By simple proprieties of the singular values, we have:
$\sigma_{m}^{2}(l_{i}^{-1}l_{j})=\lambda_{m}(\Phi(l_{i}^{-1}l_{j})$), Finally,

$\displaystyle\Phi^{*}d_{F}(l_{i},l_{j})=\sqrt{\sum_{m=1}^{N}\log^{2}(\lambda_{m}(\Phi(l_{i}^{-1}l_{j}))}=$
(14) $\displaystyle
2\sqrt{\sum_{m=1}^{N}\log^{2}(\sigma_{m}(l_{i}^{-1}l_{j}))}=2Dist_{F}(l_{i}^{-1}l_{j}).$

∎

A conclusion of Theorem 2 is that Fisher distortion of a matrix in $L_{n}^{+}$
is proportional to its distance from the identity matrix (see Fig 1).

###  3.3 Algorithm for finding MDT

Let $\\{A_{i}\in GL_{n}(\mathbb{R}\\}_{i=1}^{N}$ and let $T\in
GL_{n}(\mathbb{R})$. The Fisher distortion is invariant under orthogonal
transformation, therefore we can use QR decomposition of Theorem 1 and
decompose $A_{i}=L_{i}Q_{i}$ , $T=LQ$ where $L,L_{i}\in L_{n}^{+}$ and
$Q,Q_{i}\in O_{n}$.
Note that
$\sum_{k=1}^{N}Dist_{F}^{2}(A^{-1}A_{i})=\sum_{k=1}^{N}Dist_{F}^{2}(L^{-1}L_{i})$.
From Theorem 2, we derive that finding an MDT is equivalent to finding a
Fréchet mean on the space $L_{n}^{+}$ w.r.t the metric $\Phi^{*}d_{F}$. For
simplicity notice that

$\displaystyle Mean_{\Phi^{*}d_{F}}(L_{1},...,L_{n})=$ (15)
$\displaystyle\Phi^{-1}[Mean_{d_{F}}(\Phi(L_{1}),..,\Phi(L_{n})]=$
$\displaystyle\Phi^{-1}[Mean_{d_{F}}(A_{1}A_{1}^{t},..,A_{n}A_{n}^{t})].$

Therefore calculation of the Fréchet mean can be performed on the space
$P_{n}^{+}$ \- resulting in the following simple algorithm:

Input : A set of linear transformations $\\{A_{i}\in
GL_{n}(\mathbb{R})\\}_{i=1}^{N}$

Output : $T\in GL_{n}(\mathbb{R})$, where $T$ minimize
$\sum_{i=1}^{N}{Dist_{F}^{2}(T^{-1}A_{i})}$.

1 Calculate $P_{i}=A_{i}A_{i}^{t}$.

2 Find the Fréchet mean $\overline{P}$ of $\\{P_{i}\\}_{i=1}^{N}$ w.r.t Fisher
metric.

3 $T=\Phi^{-1}(\overline{P})$.

Algorithm 1 Mean Distorting Transformation Algorithm

In step 2, note that the Fréchet mean on the space of $P_{n}^{+}$ has no
explicit formula; however, there are several algorithms for obtaining a
numerical solution. For example - gradient based algorithms [3] and
multiplicative power means (MPM) algorithm [8].

While our algorithm supplies a unique solution, notice that any matrix $A$
satisfying $AA^{t}=\overline{P}$ is also a solution. The latter can only
happen if $A=TQ$, where $Q$ is any orthogonal matrix.

##  4 Results and Applications

###  4.1 Affine Panorama

One of the immediate applications of Algorithm 1 is in the affine panorama
flow. A complete affine panorama flow is described in [6]. In this flow,
images $\\{I_{i}\\}_{i=1}^{N}$ of a particular scene are registered to each
other by affine transformations $\\{T_{i}\\}_{i=1}^{N}$, where $T_{i}$ is
affine transformation mapping $I_{i}$ to a common plane. Classically, the
transformations are obtained by solving the system of equations
$T_{i}\vec{p}_{ik}=T_{j}\vec{q}_{jk}$ (in the least squares sense), where
$\vec{p}_{ik}$ and $\vec{q}_{jk}$ are corresponding points in the $I_{i}$ and
$I_{j}$ images.

When solving this system of equations, we have the degree of freedom of
multiplying by global affine transformation $T$. Particular choices of $T$ are
taking $T=T_{1}^{-1}$ or $T=T_{c}^{-1}$ where $c$ is the index of the the
central frame in the panorama. A more natural and symmetric choice is using
Algorithm 1 on the linear part of the affine transformations
$\\{A_{i}\\}_{i=1}^{N}$ where $T_{i}\vec{v}=A_{i}\vec{v}+\vec{b}$. After
finding the MDT solution, a rigid global correction can be applied to the
images in the following way:
First we apply a global rotation $q=exp[-\frac{\sum_{i=1}^{N}log(q_{i})}{N}]$,
where $r_{i}q_{i}$ is the QR-factorization of $A_{i}$. Afterward, we apply the
unique global shift transformation $S$ that aligns the panorama with the
overall image axis, see Fig 3 for panorama comparison between panorama
generated with a reference fixed image and MDT method.

###  4.2 3D Point cloud registration

Less researched than its rigid counterpart, non-rigid point cloud registration
appears in the literature, and various affine point cloud registration
algorithms were proposed (see [10, 9]). Similarly to the 2D case, an MDT on
the affine transformations can reduce the overall distortion of both point
sets. Note that our method has no dimensionality assumption and that Algorithm
1 should work well on any dimension.

##  5 Further Work

The natural extension of this work is to find a mean distorting transformation
for non-affine transformation. If $f:\mathbb{R}^{n}\to\mathbb{R}^{n}$ is
$C^{1}$, we can use the affine definition of the Fisher distortion of the
linearization of $f$ (Jacobi matrix $J_{f}$):

$Dist_{F}(f)=\sqrt{\int Dist_{F}^{2}(J_{f})d\mu}$ (16)

In this case, the above definition is similar to the Airy-Kavrayskiy criterion
distortion measure of map projection.

The generalized MDT statement of the problem:
Let $\\{f_{i}:D_{i}\to D\\}_{i=1}^{N}$ be a set of $C^{1}$ functions. Find a
$C^{1}$ invertible function $f:D\to D$ that minimizes the overall distortion:

$\sum_{i=1}^{N}Dist_{F}^{2}(f^{-1}\circ f_{i}).$ (17)

A compelling case to study is the case where $f_{i}\in PGL(3,\mathbb{R})$ are
2D-homographies. Those studies will enable us to apply an MDT method on planar
panoramas, which are a generalization of affine panoramas.

##  6 Acknowledgments

The author would like to thank Adi Shasha for fruitful discussions on this
topic.

---
(a)
(b)
(c)
(d)
Figure 3: Visual comparison of panoramas between the fixed reference image
method and Mean Distorting Transformation (MDT) correction applied on the
fixed reference method: (a) Park panorama (2 images), left image is fixed; (b)
MDT correction. (c) Yard Panorama (9 images), rightmost image is fixed. (d)
MDT correction. Raw images of (c) and (d) where taken from Adobe Panorama
dataset (see [5]).

## References

  * [1] Fabio Bellavia and Carlo Colombo.  Estimating the best reference homography for planar mosaics from videos.  In VISAPP (3), pages 512–519, 2015.
  * [2] Rajendra Bhatia.  On the exponential metric increasing property.  Linear Algebra and its applications, 375:211–220, 2003.
  * [3] Dario A. Bini and Bruno Iannazzo.  Computing the karcher mean of symmetric positive definite matrices.  Linear Algebra and its Applications, 438(4):1700–1710, 2013.  16th ILAS Conference Proceedings, Pisa 2010.
  * [4] Gary Bradski and Adrian Kaehler.  Learning OpenCV: Computer vision with the OpenCV library.  ” O’Reilly Media, Inc.”, 2008.
  * [5] Jonathan Brandt.  Transform coding for fast approximate nearest neighbor search in high dimensions.  In 2010 IEEE Computer Society Conference on Computer Vision and Pattern Recognition, pages 1815–1822. IEEE, 2010.
  * [6] Matthew Brown and David G Lowe.  Automatic panoramic image stitching using invariant features.  International journal of computer vision, 74(1):59–73, 2007.
  * [7] David Capel.  Image mosaicing.  In Image Mosaicing and super-resolution, pages 47–79. Springer, 2004.
  * [8] Marco Congedo, Alexandre Barachant, and Ehsan Kharati Koopaei.  Fixed point algorithms for estimating power means of positive definite matrices.  IEEE Transactions on Signal Processing, 65(9):2211–2220, 2017.
  * [9] Shaoyi Du, Nanning Zheng, Gaofeng Meng, and Zejian Yuan.  Affine registration of point sets using icp and ica.  IEEE Signal Processing Letters, 15:689–692, 2008.
  * [10] Shaoyi Du, Nanning Zheng, Shihui Ying, and Jianyi Liu.  Affine iterative closest point algorithm for point set registration.  Pattern Recognition Letters, 31(9):791–799, 2010.
  * [11] G. B. Airy Esq.  Liii. explanation of a projection by balance of errors for maps applying to a very large extent of the earth’s surface; and comparison of this projection with other projections.  The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science, 22(149):409–421, 1861.
  * [12] Tobias Jung.  Cc by-sa 4.0 license.  https://map-projections.net.
  * [13] VV Kavrajskij.  Matematiceskaja Kartografija.  Redbaza Goskartotresta, 1934.
  * [14] Kerkovits Krisztián.  Comparing finite and infinitesimal map distortion measures.  International Journal of Cartography, 5(1):3–22, 2019.
  * [15] Piotr H Laskowski.  The traditional and modern look at tissot’s indicatrix.  The American Cartographer, 16(2):123–133, 1989.
  * [16] Zhenhua Lin.  Riemannian geometry of symmetric positive definite matrices via cholesky decomposition.  SIAM Journal on Matrix Analysis and Applications, 40(4):1353–1370, 2019.
  * [17] A Makovetskii, S Voronin, V Kober, and D Tihonkih.  Affine registration of point clouds based on point-to-plane approach.  Procedia engineering, 201:322–330, 2017.
  * [18] Auguste Tissot.  Mémoire sur la représentation des surfaces et les projections des cartes géographiques.  Nouvelles annales de mathématiques: journal des candidats aux écoles polytechnique et normale, 17:145–163, 1878.


# Transverse mode-encoded quantum gate on a silicon photonic chip

Lan-Tian Feng111These authors contributed equally to this work. Key Laboratory
of Quantum Information, CAS, University of Science and Technology of China,
Hefei 230026, China.  Synergetic Innovation Center of Quantum Information $\&$
Quantum Physics, University of Science and Technology of China, Hefei 230026,
China.  Ming Zhang111missing State Key Laboratory for Modern Optical
Instrumentation, Centre for Optical and Electromagnetic Research, Zhejiang
Provincial Key Laboratory for Sensing Technologies, Zhejiang University,
Zijingang Campus, Hangzhou 310058, China.  Ningbo Research Institute, Zhejiang
University, Ningbo 315100, China.  Xiao Xiong  Di Liu  Yu-Jie Cheng  Fang-Ming
Jing  Xiao-Zhuo Qi  Yang Chen  De-Yong He  Guo-Ping Guo  Guang-Can Guo  Key
Laboratory of Quantum Information, CAS, University of Science and Technology
of China, Hefei 230026, China.  Synergetic Innovation Center of Quantum
Information $\&$ Quantum Physics, University of Science and Technology of
China, Hefei 230026, China.  Dao-Xin<EMAIL_ADDRESS>State Key
Laboratory for Modern Optical Instrumentation, Centre for Optical and
Electromagnetic Research, Zhejiang Provincial Key Laboratory for Sensing
Technologies, Zhejiang University, Zijingang Campus, Hangzhou 310058, China.
Ningbo Research Institute, Zhejiang University, Ningbo 315100, China.  Xi-Feng
<EMAIL_ADDRESS>Key Laboratory of Quantum Information, CAS, University
of Science and Technology of China, Hefei 230026, China.  Synergetic
Innovation Center of Quantum Information $\&$ Quantum Physics, University of
Science and Technology of China, Hefei 230026, China.

###### Abstract

As an important degree of freedom (DoF) in integrated photonic circuits, the
orthogonal transverse mode provides a promising and flexible way to increasing
communication capability, for both classical and quantum information
processing. To construct large-scale on-chip multimode multi-DoF quantum
systems, a transverse mode-encoded controlled-NOT (CNOT) gate is necessary.
Here, through design and integrate transverse mode-dependent directional
coupler and attenuators on a silicon photonic chip, we demonstrate the first
multimode implementation of a two-qubit quantum gate. With the aid of state
preparation and analysis parts, we show the ability of the gate to entangle
two separated transverse mode qubits with an average fidelity of $0.89\pm
0.02$ and the achievement of 10 standard deviations of violations in the
quantum nonlocality verification. In addition, a fidelity of $0.82\pm 0.01$
was obtained from quantum process tomography used to completely characterize
the CNOT gate. Our work paves the way for universal transverse mode-encoded
quantum operations and large-scale multimode multi-DoF quantum systems.

###### pacs:

In integrated optics, a significant enhancement of information processing and
communication capability is an ultimate goal due to the exponential growth in
optical interconnection and communication requirements for both classical
Dai2013  and quantum Wang2018  information applications. To achieve large-
scale photonic quantum systems, multiphoton, multiple degrees of freedom
(DoFs) and high-dimensional encoding is becoming attractive and essential
Wang2019 . In recent years, a lot of technical and experimental efforts to
achieve universal quantum operations have been done with different DoFs, such
as single qubit operation and two-qubit controlled-NOT (CNOT) gate with path
or polarization encoding O'Brien2009 ; O'Brien2003 ; Kiesel2005 ; Politi2008 ;
Crespi2011 ; Shadbolt2012 , as well as two-qubit quantum gate for time-bin
qubits Humphreys2013 ; Lo2020  and frequency-bin qubits Lu2019 . In addition,
multi-DoF hybrid encoding has also shown impressive applications such as
super-dense coding Barreiro2008  and multi-dimensional teleportation Wang2015
.

More recently, the orthogonal transverse modes supported by one multimode
waveguide are valued and developed rapidly for their high-dimensional
scalability, compactness and coherent conversion with other DoFs Li2019 ;
Luo2014 ; Dai2018 ; Feng2016 . As a new degree of freedom, transverse mode
provides the flexibility and scalability for a wide scope of novel
applications, such as dense optical interconnection Dai2013 ; Luo2014 ;
Dai2018 , quantum information science Feng2016 ; Mohanty2017 ; Feng2019 ;
Paesani2020 , nonlinear photonics Wang2017 ; Kittlaus2017 , and so on.
Meanwhile, a variety of multimode devices such as mode (de)multiplexer
Ding2013 ; Wang20142 , grating coupler Lai2018 , switch Stern2015 ; Jia2017 ;
Xiong2017 , filter Huang2016 , sharp bending and mode-independent crossing
Liu2019 ; Li2018 , and many other building blocks Dai2021  have been developed
intensively. However, as a potential on-chip DoF, the full manipulation of
transverse mode qubits, especially a two-qubit CNOT gate, is still not studied
extensively and challenging.

Figure 1: Transverse mode-encoded quantum CNOT gate. (a) Common circuit
implementation of the path-encoded quantum CNOT gate. MIPS: measurement-
induced controlled phase; EC: energy compensator. (b) The transverse mode in a
multimode waveguide can be directly used for quantum information encoding in a
compact manner. (c) Simulated light propagation in the designed TMDDC. As the
two identical multimode waveguides gradually merge, the TE0 mode remains in
the top waveguide without any coupling, while the TE1 mode is cross-coupled
with a ratio of 2/3 from the upper waveguide to the lower one. (d) Simulated
light propagation in the designed MMA. (e) Schematic of the entire silicon
photonic integrated circuit. Photon pairs generated in free space are coupled
to the chip through a fiber array. After state evolution in the chip, the
photon pairs are coupled to another fiber array and received by single-photon
detectors. The dashed black frame marks out the multimode section, where both
modes are supported and function together as a transverse mode encoded CNOT
gate. The other parts of the circuits are single-mode and used for state
preparation and analysis. MMI BS: multimode interference beam splitter; MMA:
multimode attenuator; TMDDC: transverse mode dependent directional coupler.

In this paper, we report the on-chip generation, manipulation, and measurement
of transverse mode-encoded quantum entangled states based on a compact CNOT
gate. The CNOT gate relies on two newly developed multimode devices that can
selectively manipulate different transverse modes, called unbalanced
transverse mode-dependent directional coupler (TMDDC) and multimode
attenuators (MMAs). This multimode waveguide CNOT gate is very compact and has
a size of approximately 10$\times$160 $\upmu$m2. We characterize it by
performing quantum state tomography and quantum process tomography
measurements. The gate realizes the entanglement of two separable qubits with
an average state fidelity of $\bar{F}_{s}=0.89\pm 0.02$ and a process fidelity
of $F_{p}=0.82\pm 0.01$. Because transverse mode-encoding is compatible with
other encoding protocols and facilitates the construction of quantum photonic
integrated circuits (QPICs) with higher dimensions and capacity, our work
paves the way for universal transverse mode-encoded quantum operations and
practical applications in compact and scalable photonic quantum processing
with multiple DoFs.

A quantum processor requires a cascade of various unitary operations, which
consist of a finite series of basic logic gates. It has been proven that a
two-qubit CNOT gate, which flips the target qubit $\left|t\right\rangle$
depending on the state of the control qubit $\left|c\right\rangle$, together
with single-qubit logic gates, is sufficient for universal quantum operation
Nielsen2000 . In QPICs, the universal quantum operation is always implemented
using beam splitters and phase controllers O'Brien2009 ; Politi2008 ;
Crespi2011 ; Shadbolt2012 . The common circuit implementation of the quantum
CNOT gate is shown in Fig. 1a. It is mainly composed of three unbalanced
directional couplers, one for quantum interference to induce the controlled-
phase (CPhase) operation, and the other two for energy compensation. Two 50/50
beam splitters are used to transform the CPhase operation to the CNOT
operation by applying the relation $\mathrm{CNOT}=(I\otimes
H)\mathrm{CPhase}(I\otimes H)$, where $I$ and $H$ denote the identity and
Hadamard gates. This transformation can also be easily achieved by selecting
the computation basis of the target qubit, as demonstrated in Ref. Crespi2011
.

Among the various material systems for QPICs, silicon-on-insulator was chosen
here because of its excellent optical properties and CMOS compatibility
Bogaerts2005 . The designed multimode waveguide has a cross-section of
$760\times 220~{}\textrm{nm}^{2}$ and supports the two lowest-order TE
polarization transverse modes, that is, the TE0 and TE1 modes, as shown in
Fig. 1b. The modes are orthogonal to each other and no cross-coupling between
them, thus they can be directly used for quantum information encoding. We have
newly designed two mode dependent devices, TMDDC and MMA, which can achieve
mode selective coupling and loss, as shown in Fig. 1c and 1d. By cascading a
TMDDC and two MMAs, we constructed one multimode two-qubit CNOT gate, Fig. 1e.

The TMMDC was based on a symmetric directional coupler, in which the incident
modes (i.e., TE0, TE1) can be partially transferred to the adjacent waveguide
by evanescent coupling. Here, to achieve a mode selection beam splitter, we
chose a wide waveguide width $w_{0}=760~{}\textrm{nm}$ and a large gap
$g=400~{}\textrm{nm}$. In this situation, the evanescent field of TE0 mode
will be restricted in the waveguide with barely no power transfer, but the TE1
mode still has enough coupling strength and a coupling length $L=32.5~{}$µm
was used to realize a cross-coupling ratio of 2/3. Fig. 1(c) shows the
simulated light propagation in the designed structure when the TE0 and TE1
modes were launched. The MMA was designed to attenuate 2/3 power of the TE0
mode and guarantee lossless propagation of TE1 mode for energy compensation.
In order to realize such an energy attenuator, a multimode interferometer
(MMI) was introduced due to its compact size ( <3$\times$10 µm2) and low
insertion losses (IL <0.5 dB). For the lossless propagation of TE1, we set the
MMI length determined by the interference-induced self-imaging of TE1 mode. In
this case, the TE0 mode shows imperfect self-imaging and the output energy
depends on the device width. Here, we chose a device width $w=2.5~{}$ µm and
length $l=7.5~{}$ µm, which will maintain lossless propagation of TE1 mode
while attenuating 2/3 of the power of TE0 mode. In order to radiate away the
unwanted energy left in the attenuator which would cause unnecessary
interference and negatively affect the output, we introduced a butterfly-shape
MMI to reduce reflection of TE0 mode by adding 4 acute angles to the MMI. Fig.
1(d) shows the simulated light propagation in the designed attenuator when the
TE0 and TE1 modes were launched. For more information on these two structures,
see the Supplementary Materials (Section I). The clever design of controlling
transverse modes separately here enables us to further cascading them for
realizing complex functions with mode encoding such as two-qubit CNOT gate.

In the experiment, a PPKTP crystal was used to generate frequency-degenerate
photon pairs centered at wavelength of 1540 nm. More details are provided in
the Supplementary Materials (Section II). The photon pairs were collected with
a single-mode fiber array and then coupled into the TE0 mode in the on-chip
single-mode waveguides via grating couplers. To encode photons with different
transverse modes, mode (de)multiplexers Dai2013  were introduced before
(after) the photons have passed through the CNOT gate. After state evolution
in the chip, the photons were coupled using another fiber array and detected
with superconducting nanowire single-photon detectors for the coincidence
measurement. To prepare and detect the different quantum states, path encoding
Mach-Zehnder (MZ) interferometers were performed. Micro-heaters controlled by
a homemade multi-channel direct-current regulated power supply were used to
thermally tune the phases of the photons. All the heaters were individually
characterized and displayed high quality. See the Supplementary Materials
(Section IV) for further details.

To characterize the chip with different input qubits, we defined logical and
physical qubits with the following states:
$\left|0\right\rangle_{c}=\left|\mathrm{TE}_{0}\right\rangle$,
$\left|1\right\rangle_{c}=\left|\mathrm{TE}_{1}\right\rangle$,
$\left|0\right\rangle_{t}=(\left|\mathrm{TE}_{0}\right\rangle+\left|\mathrm{TE}_{1}\right\rangle)/\sqrt{2}$,
$\left|1\right\rangle_{t}=(\left|\mathrm{TE}_{0}\right\rangle-\left|\mathrm{TE}_{1}\right\rangle)/\sqrt{2}$,
where the subscripts $c$ and $t$ denote the control and target qubits. This
configuration has a success probability of $P=1/9$, and it can be rendered
deterministic by making use of ancillary resources, measurements, and feed-
forward Knill2001 .

Figure 2: Quantum interference in the TMDDC. Both the control and target
photons are in the TE1 mode and have an adjustable time delay between them.
The dashed blue line denotes the lower limit for an ideal quantum
interferometer with 1:2 split ratio. The error bars represent $\pm\sqrt{N}$
for the raw coincidence data (black rectangles) $N$. The dots denote the
experimental data, and the curves the corresponding Gaussian fit.

Two-photon quantum interference is essential in the CNOT gate and the gate
operation succeeds when two photons arrive at the directional coupler at the
same time, thus the performance of the CNOT gate is directly related to the
interference visibility. We first characterized the quantum interference in
the multimode directional coupler. The input/output control and target photons
were $\left|c\right\rangle=\left|\mathrm{TE}_{1}\right\rangle$ and
$\left|t\right\rangle=\left|\mathrm{TE}_{1}\right\rangle$, respectively. By
adjusting the time delay between the two photons, a quantum interference
fringe was observed, as shown in Fig. 2. An interference visibility of
$0.82\pm 0.03$ was derived from Gaussian fitting before background subtraction
(the raw data used, the same below) where the error was the fitting standard
error. Here, the visibility is defined as
$V=(C_{\rm{max}}-C_{\rm{min}})/C_{\rm{max}}$, where $C_{\rm{max}}$ and
$C_{\rm{min}}$ are the maximum and minimum of fitted data, respectively. For
the configuration in Fig. 1(a), the quantum interference should ideally give
$C_{\rm{min}}=0.2C_{\rm{max}}$, and thus the ideal interference visibility is
0.80, as indicated by the dashed blue line in Fig. 2. The excellent agreement
between the measured and ideal values for the quantum interference visibility
proves that the TMDDC works well. In principle, the CNOT gate operation
between the control and target qubits can be realized at zero time delay (dip
position).

Figure 3: Quantum state tomography of the generated Bell states. The real
parts of the measured density matrices of the output states with the different
input states of (a) $\left|+\right\rangle_{c}\left|0\right\rangle_{t}$, (b)
$\left|-\right\rangle_{c}\left|0\right\rangle_{t}$, (c)
$\left|+\right\rangle_{c}\left|1\right\rangle_{t}$, and (d)
$\left|-\right\rangle_{c}\left|1\right\rangle_{t}$. The imaginary parts of the
measured density matrices are negligible. The clear bars represent the ideal
density matrices for the maximally entangled Bell states $|\Phi^{\pm}\rangle$
and $|\Psi^{\pm}\rangle$.

The basic function of a CNOT gate is entangling two separate qubits. This is
one fundamental operation in quantum information processing and represents the
most non-classical implication of quantum mechanics. Here, we use this
operation to characterize the CNOT gate with an input of
$\\{\left|\pm\right\rangle_{c}\left|0\right\rangle_{t},\left|\pm\right\rangle_{c}\left|1\right\rangle_{t}\\}$
where
$\left|+\right\rangle_{c}=(\left|0\right\rangle_{c}+\left|1\right\rangle_{c})/\sqrt{2}$
and
$\left|-\right\rangle_{c}=(\left|0\right\rangle_{c}-\left|1\right\rangle_{c})/\sqrt{2}$.
The output states were collected on-chip and reconstructed by quantum state
tomography James2001 .

Figure 4: Fig. 4. Quantum process tomography measurement results. Real parts
of the density matrices for (a) an ideal CNOT gate $\chi_{\rm{ideal}}$ and (b)
the measured $\chi_{\rm{exp}}$. The imaginary part of $\chi_{\rm{exp}}$ is
negligible. Here, $\bm{X}$, $\bm{Y}$, and $\rm{Z}$ represent the Pauli
matrices $\sigma_{\rm{X}}$, $\sigma_{\rm{Y}}$, and $\sigma_{\rm{Z}}$,
respectively.

With a perfect CNOT gate, the input states are expected to be converted into
the maximally entangled Bell states, that is,

$\displaystyle|\Phi^{\pm}\rangle=\frac{1}{\sqrt{2}}(\left|0\right\rangle_{c}\left|0\right\rangle_{t}\pm\left|1\right\rangle_{c}\left|1\right\rangle_{t}),$
(1)
$\displaystyle|\Psi^{\pm}\rangle=\frac{1}{\sqrt{2}}(\left|0\right\rangle_{c}\left|1\right\rangle_{t}\pm\left|1\right\rangle_{c}\left|0\right\rangle_{t}).$

The corresponding density matrix is denoted as
$\hat{\rho}_{\rm{ideal}}=|\phi\rangle\langle\phi|$, where
$|\phi\rangle=|\Phi^{\pm}\rangle,|\Psi^{\pm}\rangle$. Fig. 3 shows the real
parts of the reconstructed density matrices $\hat{\rho}_{\rm{mea}}$ and
$\hat{\rho}_{\rm{ideal}}$ (crystal clear bars). The raw fidelity of the
generated state, defined as
$F_{s}=[\rm{Tr}(\sqrt{\sqrt{\hat{\rho}_{\rm{mea}}}\hat{\rho}_{\rm{ideal}}\sqrt{\hat{\rho}_{\rm{mea}}}})]^{2}$,
were calculated for the different output states as
$F_{|\Phi^{+}\rangle}=0.89\pm 0.03$, $F_{|\Phi^{-}\rangle}=0.88\pm 0.03$,
$F_{|\Psi^{+}\rangle}=0.90\pm 0.02$, $F_{|\Psi^{-}\rangle}=0.87\pm 0.02$. The
averaged fidelity was obtained as $\bar{F}_{s}=0.89\pm 0.02$, which confirms
the entangling function of the CNOT gate. The fidelity errors were obtained by
a 100-times Monto Carlo calculation with the experimental data subjected to
Gaussian statistics. In addition, the linear entropy and tangle were
calculated to quantify the mixture and amount of entanglement for these states
James2001 ; Wei2003 ; Coffman2000 . The details are provided in the
Supplementary Materials (Section V).

The nonlocality nature of the output states was further verified using the
Clauser-Horne-Shimony-Holt (CHSH) inequality Clauser1969 ,

$\hat{S}(\hat{\rho})=E(a,b)-E(a,b^{\prime})+E(a^{\prime},b)+E(a^{\prime},b^{\prime}).$
(2)

Here, $a$, $a^{\prime}$ and $b$, $b^{\prime}$ are two sets of measurement
direction settings for different qubits. For a state that can be described by
a local classical theory, the non-locality parameter $\hat{S}$ should be less
than or equal to two. In other words, if the output state $\hat{\rho}$ cannot
be described by a local classical theory, the inequality is violated. For
entangled Bell states, quantum mechanics predicts a maximal value of
$\hat{S}_{\rm{max}}=2\sqrt{2}$ with a selected measurement basis. In this
study, the phases in the MZIs at the measuring part were individually set as
$\\{0^{\circ},180^{\circ},45^{\circ},225^{\circ}\\}$ for the control qubit and
$\\{22.5^{\circ},202.5^{\circ},67.5^{\circ},247.5^{\circ}\\}$ for the target
qubit. In this way, $\hat{S}$ values for the different output states were
obtained as $\hat{S}_{|\Phi^{+}\rangle}=2.54\pm 0.05$,
$\hat{S}_{|\Phi^{-}\rangle}=2.50\pm 0.05$, $\hat{S}_{|\Psi^{+}\rangle}=2.51\pm
0.05$, $\hat{S}_{|\Psi^{+}\rangle}=2.43\pm 0.05$. The average value of the
non-locality parameter is $\hat{S}_{\rm{Bell}}=2.50\pm 0.05$, which violates
the inequality by 10 standard deviations.

Finally, we performed quantum process tomography to fully characterize the
quantum device. We prepared 16 different input states from the set
$\\{\left|0\right\rangle,\left|1\right\rangle,\left|0\right\rangle+\left|1\right\rangle,\left|0\right\rangle+i\left|1\right\rangle\\}^{\otimes
2}$and analyzed the output states. For each output state, 16 analysis bases
were projected for the quantum state tomography measurement. Therefore, 256
projections were performed to fully characterize the two-qubit CNOT gate. The
experimentally reconstructed process matrix $\chi_{\rm{exp}}$ is presented
together with the ideal $\chi_{\rm{ideal}}$ in Fig. 4. Using the definition of
the process fidelity
$F=[\rm{Tr}(\sqrt{\sqrt{\chi_{\rm{exp}}}\chi_{\rm{ideal}}\sqrt{\chi_{\rm{exp}}}})]^{2}$
Crespi2011 , we obtained $F=0.82\pm 0.01$, which is similar to or better than
those reported in recent works Crespi2011 ; Politi2008 ; Shadbolt2012 .

The technological progress we have made here is multifaceted. First, as new
multimode devices, the beam splitter and attenuator we designed and
manufactured can also be used as basic components for various optical
information processing applications. Second, the CNOT gate we achieved is the
first and indispensable demonstration of a two-qubit logic gate with
transverse mode-encoding. Together with the single qubit rotation operation
realized previously Mohanty2017 , the CNOT gate makes universal quantum
computation with transverse mode-encoded qubits possible. Compared to the chip
size of path- and polarization-encoding CNOT gates based on dielectric
waveguides (approximately 105-107 $\upmu$m2) Crespi2011 ; Politi2008 ,
transverse-mode encoding enables us to further shrink the chip size to
approximately 10$\times$160 $\upmu$m2 and improve the chip stability when it
becomes necessary to cascade multiple logic gates in large-scale circuits. The
gate is promising for high-dimensional quantum information applications by
introducing more transverse-modes. The gate is also compatible with the use of
other degrees of freedom for encoding Feng2016  and offers a platform for
hybrid encoding with diverse functions. Fully on-chip manipulation of
transverse mode-encoded qubits is readily achievable with the integration of
an entangled photon pair source using a multimode waveguide Feng2019 . It
should be noted that the imperfect gate fidelity reported here is mainly due
to the fabrication error, which resulted in insertion loss for the higher-
order transverse mode. In principle, the higher-order transverse modes
supported by a multimode waveguide can be transmitted with negligible loss
through the use of expected advanced fabrication technologies.

In conclusion, we have demonstrated the manipulation of transverse mode-
encoded qubits using a CNOT gate, which consists of newly developed multimode
directional coupler and attenuator. The measured quantum state tomography
shows that the entangled Bell states were generated by the gate operation with
an average fidelity of $\bar{F}_{s}=0.89\pm 0.02$, and that quantum locality
was violated by 10 standard deviations. Quantum process tomography was also
performed to fully characterize the device, which exhibited a process fidelity
of $F_{p}=0.82\pm 0.01$. The device is compact and extendable to higher
dimensions, and paves the way for on-chip quantum information processing with
multiple degrees of freedom and higher capacity.

This work was supported by the National Natural Science Foundation of China
(NSFC) (Nos. 61590932, 11774333, 62061160487, 12004373, 91950205,
61961146003), National Science Fund for Distinguished Young Scholars
(61725503), the Strategic Priority Research Program of the Chinese Academy of
Sciences (No. XDB24030601), the National Key R & D Program (No.
2016YFA0301700), Anhui Initiative in Quantum Information Technologies (No.
AHY130300), Zhejiang Provincial Natural Science Foundation (LZ18F050001,
LD19F050001), the Postdoctoral Science Foundation of China (No. 2020M671860)
and the Fundamental Research Funds for the Central Universities. This work was
partially carried out at the USTC Centre for Micro and Nanoscale Research and
Fabrication.

## References

  * (1) D. Dai, J. Wang, and Y. Shi. Y. Silicon mode (de)multiplexer enabling high capacity photonic networks-on-chip with a single-wavelength-carrier light. _Opt. Lett._ 38, 1422-1424 (2013).
  * (2) J. Wang, S. Paesani, Y. Ding, R. Santagati, P. Skrzypczyk, A. Salavrakos, J. Tura, R. Augusiak, L. Mancinska, D. Bacco, D. Bonneau, J. W. Silverstone, Q. Gong, A. Acin, K. Rottwitt, L. K. Oxenlowe, J. L. O’Brien, A. Laing, and M. G. Thompson. Multidimensional quantum entanglement with large-scale integrated optics. _Science_ 360, 285–291 (2018).
  * (3) J. Wang, F. Sciarrino, A. Laing, and M.G. Thompson. Integrated photonic quantum technologies. _Nature Photon._ https://doi.org/10.1038/s41566-019-0532-1 (2019).
  * (4) J. L. O’Brien, A. Furusawa, and J. Vučković. Photonic quantum technologies. _Nature Photon._ 3, 687-695 (2009).
  * (5) J. L. O’Brien, G. J. Pryde, A. G. White, T. C. Ralph, and D. Branning. Demonstration of an all-optical quantum controlled-NOT gate. _Nature_ 426, 264–267 (2003).
  * (6) N. Kiesel, C. Schmid, U. Weber, R. Ursin, and H. Weinfurter. Linear Optics Controlled-Phase Gate Made Simple. _Phys. Rev. Lett._ 95, 210505 (2005).
  * (7) A. Politi, M. J. Cryan, J. G. Rarity, S. Yu, and J. L. O’brien. Silica-on-silicon waveguide quantum circuits. _Science_ 409, 646-649 (2008).
  * (8) A. Crespi, R. Ramponi, R. Osellame, L. Sansoni, I. Bongioanni, F. Sciarrino, G. Vallone, and P. Mataloni. Integrated photonic quantum gates for polarization qubits. _Nat. Commun._ 2, 566 (2011).
  * (9) P. J. Shadbolt, M. R. Verde, A. Peruzzo, A. Politi, A. Laing, M. Lobino, J. C. F. Matthews, M. G. Thompson, and J. L. O’Brien. Generating, manipulating and measuring entanglement and mixture with a reconfigurable photonic circuit. _Nature Photon._ 6, 45-49 (2012).
  * (10) P. C. Humphreys, B. J. Metcalf, J. B. Spring, M. Moore, X. Jin, M. Barbieri, W. S. Kolthammer, and I. A. Walmsley. Linear Optical Quantum Computing in a Single Spatial Mode. _Phys. Rev. Lett._ 111, 150501 (2013).
  * (11) H. P. Lo, T. Ikuta, N. Matsuda, T. Honjo, W. J. Munro, and H. Takesue. Quantum Process Tomography of a Controlled-Phase Gate for Time-Bin Qubits. _Phys. Rev. Applied_ 13, 034013 (2020).
  * (12) H. H. Lu, J. M. Lukens, B. P. Williams, P. Imany, N. A. Peters, A. M. Weiner, and P. Lougovski. A controlled-NOT gate for frequency-bin qubits. _npj Quantum Information_ 5, 24 (2019).
  * (13) J. T. Barreiro, T. C. Wei, and P. G. Kwiat. Beating the channel capacity limit for linear photonic superdense coding. _Nature Phys._ 4, 282–286 (2008).
  * (14) X. L. Wang, X. D. Cai, Z. E. Su, M. C. Chen, D. Wu, L. Li, N. L. Liu, C. Y. Lu, and J. W. Pan. Quantum teleportation of multiple degrees of freedom of a single photon. _Nature_ 518, 516–519 (2015).
  * (15) C. Li, D. Liu, D. Dai, Multimode silicon photonics. _Nanophotonics_ 8, 227-247 (2019).
  * (16) L. W. Luo, N. Ophir, C. P. Chen, L. H. Gabrielli, C. B. Poitras, K. Bergmen, and M. Lipson. WDM-compatible mode-division multiplexing on a silicon chip. _Nat. Commun._ 5, 3069 (2014).
  * (17) D. Dai, C. Li, S. Wang, H. Wu, Y. Shi, Z. Wu, S. Gao, T. Dai, H. Yu, and H. K. Tsang. $10$-channel mode (de) multiplexer with dual polarizations. _Laser Photonics Rev._ 12, 1700109 (2018).
  * (18) L. T. Feng, M. Zhang, Z. Y. Zhou, M. Li, X. Xiong, L. Yu, B. S. Shi, G. P. Guo, D. X. Dai, X. F. Ren, and G. C. Guo. On-chip coherent conversion of photonic quantum entanglement between different degrees of freedom. _Nat. Commun._ 7, 11985 (2016).
  * (19) A. Mohanty, M. Zhang, A. Dutt, S. Ramelow, P. Nussenzveig, and M. Lipson. Quantum interference between transverse spatial waveguide modes. _Nat. Commun._ 8, 14010 (2017).
  * (20) L. T. Feng, M. Zhang, X. Xiong, Y. Chen, H. Wu, M. Li, G. P. Guo, G. C. Guo, D. X. Dai, and X. F. Ren. On-chip transverse-mode entangled photon pair source. _npj Quantum Information_ 5, 2 (2019).
  * (21) S. Paesani, M. Borghi, S. Signorini, A. Maïnos, L. Pavesi, and A. Laing. Near-ideal spontaneous photon sources in silicon quantum photonics. _Nat. Commun._ 11, 2505 (2020).
  * (22) C. Wang, X. Xiong, N. Andrade, V. Venkataraman, X. F. Ren, G. C. Guo, and M. Lončar. Second harmonic generation in nano-structured thin-film lithium niobate waveguides. _Opt. Express_ 25, 6963-6973 (2017).
  * (23) E. A. Kittlaus, N. T. Otterstrom, and P. T. Rakich, On-chip inter-modal Brillouin scattering. _Nat. Commun._ 8, 15819 (2017).
  * (24) Y. Ding, J. Xu, F. D. Ros, B. Huang, H. Ou, and C. Peucheret. On-chip two-mode division multiplexing using tapered directional coupler-based mode multiplexer and demultiplexer. _Opt. Express_ 21, 10376–10382 (2013).
  * (25) J. Wang, P. Chen, S. Chen, Y. Shi, and D. Dai. Improved 8-channel silicon mode demultiplexer. _Opt. Express_ 22, 12799–12807 (2014).
  * (26) Y. Lai, Y. Yu, S. Fu, J. Xu, P. P. Shum, and X. Zhang. Compact double-part grating coupler for higher-order mode coupling. _Opt. Lett._ 43, 3172–3175 (2018).
  * (27) B. Stern, X. Zhu, C. P. Chen, L. D. Tzuang, J. Cardenas, K. Bergman, and M. Lipson. On-chip mode-division multiplexing switch. _Optica_ 2, 530-535 (2015).
  * (28) H. Jia, T. Zhou, L. Zhang, J. Ding, X. Fu, and L. Yang. Optical switch compatible with wavelength division multiplexing and mode division multiplexing for photonic networks-on-chip. _Opt. Express_ 25, 20698–20707 (2017).
  * (29) Y. Xiong, R. Priti, and O. Ladouceur. High-speed two-mode switch for mode-division multiplexing optical networks. _Optica_ 4, 1098–1102 (2017).
  * (30) Q. Huang, K. Jie, Q. Liu, Y. Huang, Y. Wang, and J. Xia. Ultra-compact, broadband tunable optical bandstop filters based on a multimode one-dimensional photonic crystal waveguide. _Opt. Express_ 24, 20542-20553 (2016).
  * (31) Y. Liu, K. Xu, S. Wang, W. Shen, H. Xie, Y. Wang, S. Xiao, Y. Yao, J. Du, Z. He and Q. Song. Arbitrarily routed mode-division multiplexed photonic circuits for dense integration. _Nat. Commun._ 10, 3263 (2019).
  * (32) S. Li, Y. Zhou, J. Dong, X. Zhang, E. Cassan, J. Hou, C. Yang, S. Chen, D. Gao, and H. Chen. Universal multimode waveguide crossing based on transformation optics. _Optica_ 5, 1549-1556 (2018).
  * (33) D. Dai, W. Zhao, J. Guo, and D. Liu. Multimode silicon photonic devices. _Integrated Optics: Devices, Materials, and Technologies XXV_ 11689, 116890N (2021).
  * (34) M. Nielsen, and I. Chuang. Quantum Computation and Quantum Information, Cambridge University Press, Cambridge (2000).
  * (35) W. Bogaerts, R. Baets, P. Dumon, V. Wiaux, S. Beckx, D. Taillaert, B. Luyssaert, J. V. Campenhout, P. Bienstman, and D. V. Thourhout. Nanophotonic waveguides in silicon-on-insulator fabricated with CMOS technology. _J. Lightwave Technol._ 23, 401-412 (2005).
  * (36) E. Knill, R. Laflamme, and G. J. Milburn. A scheme for efficient quantum computation with linear optics. _Nature_ 409, 46-52 (2001).
  * (37) D. F. V. James, P. G. Kwiat, W. J. Munro, and A. G. White. Measurement of qubits. _Phys. Rev. A_ 64, 052312 (2001).
  * (38) T. C. Wei, K. Nemoto, P. M. Goldbart, P. G. Kwiat, W. J. Munro, and F. Verstraete. Maximal entanglement versus entropy for mixed quantum states. _Phys. Rev. A_ 67, 022110 (2003).
  * (39) V. Coffman, J. Kundu, and W. K. Wootters. Distributed entanglement. _Phys. Rev. A_ 61, 052306 (2000).
  * (40) J. F. Clauser, M. A. Horne, A. Shimony, and R. A. Holt. Proposed experiment to test local hidden-variable theories. _Phys. Rev. Lett._ 23, 880 (1969).