# Optimal control of molecular spin qudits

Alberto Castro  ARAID Foundation, Avenida de Ranillas 1, 50018 Zaragoza, Spain
Institute for Biocomputation and Physics of Complex Systems (BIFI) of the
University of Zaragoza, Edificio de Institutos Universitarios de
Investigación, Calle Mariano Esquillor, 50018 Zaragoza, Spain<EMAIL_ADDRESS>Adrián García Carrizo  Institute for Biocomputation and Physics of Complex
Systems (BIFI) of the University of Zaragoza, Edificio de Institutos
Universitarios de Investigación, Calle Mariano Esquillor, 50018 Zaragoza,
Spain  David Zueco  Instituto de Nanociencia y Materiales de Aragón (INMA),
CSIC-Universidad de Zaragoza, Pedro Cerbuna 12, 50009 Zaragoza, Spain
Departamento de Física de la Materia Condensada, Universidad de Zaragoza,
50009 Zaragoza, Spain  Fernando Luis  Instituto de Nanociencia y Materiales de
Aragón (INMA), CSIC-Universidad de Zaragoza, Pedro Cerbuna 12, 50009 Zaragoza,
Spain  Departamento de Física de la Materia Condensada, Universidad de
Zaragoza, 50009 Zaragoza, Spain

###### Abstract

We demonstrate, numerically, the possibility of manipulating the spin states
of molecular nanomagnets with shaped microwave pulses designed with quantum
optimal control theory techniques. The state-to-state or full gate
transformations can be performed in this way in shorter times than using
simple monochromatic resonant pulses. This enhancement in the operation rates
can therefore mitigate the effect of decoherence. The optimization protocols
and their potential for practical implementations are illustrated by
simulations performed for a simple molecular cluster hosting a single Gd3+
ion. Its eight accessible levels (corresponding to a total spin $S=7/2$) allow
encoding an $8$-level qudit or a system of three coupled qubits. All necessary
gates required for universal operation can be obtained with optimal pulses
using the intrinsic couplings present in this system. The application of
optimal control techniques can facilitate the implementation of quantum
technologies based on molecular spin qudits.

##  I  Introduction

A crucial challenge for the development of quantum simulation and quantum
computation is to scale up computational power while keeping the processor
robust against noise and limiting the complexity of control lines and
electronics [1, 2, 3]. A promising strategy is to replace qubits with
$d$-dimensional ($d>2$) quantum systems, or qudits [4, 5, 6], as the
elementary building blocks of the quantum architecture. The ability to
integrate nontrivial operations in a single physical system helps simplifying
some quantum algorithms [7, 8, 9]. In addition, it can also facilitate their
implementation by reducing the number of nonlocal operations, i.e. those
connecting different parts of the circuit.

Qudits have been realized with the multiple quantum states of diverse physical
systems, including photons [10], trapped ions [11], impurity nuclear spins in
semiconductors [12] and superconducting circuits [13]. Here, we focus on a
special class of systems, molecular nanomagnets [14, 15, 16, 17, 18] (see Fig.
1 for an illustrative example). These are coordination or supramolecular
complexes that consist of a magnetic core surrounded by a shell of organic
ligand molecules. Chemistry offers nearly unbound possibilities for the design
of spin qudits based on these molecules. The combination of one or several
$S>1/2$ transition metal or lanthanide ions with sufficiently weak magnetic
anisotropy and/or exchange interactions gives rise to a number of low-lying
magnetic levels. For suitably chosen molecular structures, these levels are
unequally spaced, making transitions between them addressable via microwave
resonant pulses. And nuclear spin states of the metal ions provide additional
resources [19, 20, 21, 22, 23].

Examples of molecular based electronic and electronuclear spin qudits, with
dimension $d$ ranging from $4$ up to $64$, have been reported recently [24,
25, 26, 27, 28, 19, 20, 21, 29, 30, 22, 23]. In addition, there have been
proposals for exploiting their multiple states to specific applications.
Relevant examples are the digital quantum simulation of spin-boson models
[31], where the qudit encodes boson states, and the implementation of quantum
error correction codes [32, 30, 33]. The latter is particularly promising, as
embedding in each basic unit, in this case a molecule, a suitably designed
protection against its specific decoherence sources might represent a huge
competitive advantage. Besides, the functionalities need not be defined a
priori. When the allowed transitions between different qudit states form a
universal set, these systems can be regarded as microscopic size universal
processors (or NISQs) [27, 28, 34, 29, 22, 18].

However, decoherence remains a serious limitation for fully unleashing the
potential of these otherwise very appealing systems. Even though some
specifically designed molecular spin qubits show remarkably long coherence
times $T_{2}$ [35, 36], decoherence tends to increase for higher spin or
higher nuclearity molecules. A sequence of necessarily imperfect gates might
become impractically long as compared to $T_{2}$, leading to large error
rates. In the paradigmatic example of a qudit-based quantum error protocol,
such effects can completely overcome the gain brought about by the code itself
[32, 30]. Clearly, this calls for more efficient methods to carry out such
operations.

Figure 1:  Top: Schematic image of a molecular spin subject to a static
magnetic field $\vec{B}$ and to an arbitrarily shaped microwave magnetic field
$\vec{b}$, here generated by a transmission line. Bottom: spin energy levels
of the GdW30 polyoxometalate cluster [37, 27], whose structure is shown in the
top. The $8$ spin states, associated with the $S=7/2$ spin of the central Gd3+
ion, enable encoding a $d=8$ qudit or three qubits. The coloured circular
arrows mark the seven transitions that can be implemented by means of resonant
monochromatic pulses.

In this work, we consider the possibility of applying quantum optimal control
theory (QOCT) techniques [38] to mitigate some of the limitations associated
with the use of monochromatic resonant pulses. This theory allows designing
more complex pulses in order to find the temporal shape of the external
perturbation that makes the evolution operator equal to a predefined gate, as
shown for the first time by Palao and Kosloff [39]. Numerous later
calculations employing similar schemes have been reported [40, 41, 42, 43, 44,
45, 46]. We illustrate its application to the efficient control of molecular
spin qudits by performing numerical simulations on a $d=8$ qudit encoded in
the $S=7/2$ spin states of a simple molecular cluster hosting a Gd3+ ion (Fig.
1) [37], for which a universal set of operations was realized experimentally
[27]. We have parameterized the pulses using a range of different frequencies,
and limiting the amplitude of each component (since otherwise the
optimizations tend to favour high intensity solutions, that may be
experimentally inaccessible, and moreover, would reduce the coherence times
via the excitation of unwanted levels). The purpose is to assess the expected
gains of using a multi-frequency setup, compared to the usual route to gate
construction through simple monochromatic rotations.

The manuscript is organized as follows. In Section II we introduce the spin
Hamiltonian describing the qudit, discuss standard control techniques based on
monochromatic resonant pulses and its limitations in terms of operations
speeds and or fidelities, and describe the optimization methodology used in
this work. Section III shows results obtained by these two methods and
discusses their differences. Finally, Section IV summarizes the conclusions
and future prospects derived from these results.

##  II  Model and methodology

###  II.1 Spin Hamiltonian: definition of qudit basis states

We consider a situation like that shown schematically in Fig. 1, where a
molecular spin qudit is tuned by an external magnetic field $\vec{B}$ and
controlled by a time-dependent perturbation able to induce transitions between
its different states. The underlying physics can be described by the following
spin Hamiltonian

$\hat{H}=\hat{H}_{\rm ZF}+\hat{H}_{\rm Zeeman}+\hat{H}_{\rm td}$ (1)

where $\hat{H}_{\rm ZF}$ is the zero-field spin Hamiltonian of the isolated
molecule, $\hat{H}_{\rm Zeeman}$ is the Zeeman interaction with $\vec{B}$ and
$\hat{H}_{\rm td}$ is a time-dependent control term. In general, $\hat{H}_{\rm
ZF}$ can be quite complex and include the single ion anisotropy of the
magnetic ions forming the molecular core, their mutual exchange and dipolar
interactions as well as the hyperfine couplings to the nuclear spins. Yet,
quite often this Hamiltonian can be simplified. This is the case when dealing
with mononuclear molecules, hosting one metal ion, or when exchange
interactions are sufficiently strong to ensure that only states with the
lowest energy total spin value $S$ are significantly populated at the relevant
temperatures. The spin response can then be well approximated with the help of
a “giant spin” approximation [14, 47]. In this work, we have used the
following expressions for $\hat{H}_{\rm ZF}$ and $\hat{H}_{\rm Zeeman}$

$\hat{H}_{\rm
ZF}=D\left[\hat{S}_{z}-\frac{1}{3}S(S+1)\right]+E\left(\hat{S}^{2}_{x}-\hat{S}^{2}_{y}\right)$
(2) $\hat{H}_{\rm Zeeman}=-g\mu_{\rm B}\vec{B}\hat{\vec{S}}$ (3)

where $S$ is the spin quantum number of the molecule,
$(\hat{S}_{x},\hat{S}_{y},\hat{S}_{z})$ are the spin operators, $D$ and $E$
are magnetic anisotropy constants, $g$ is the spin $g$-factor, and $\mu_{\rm
B}$ is the Bohr magneton. These expressions describe accurately the GdW30
molecular spin qudit, which we use in section III below to illustrate the
potential of QOCT techniques. In particular, the uniaxial magnetic anisotropy
$DS_{z}^{2}$ provides the level anharmonicity that is required to properly
address individual transitions between the qudit states. Yet, the methodology
is general and could be applied to more complex versions of the spin
Hamiltonian adapted to diverse implementations (e.g. those including weakly
coupled electronic spins or a combination of nuclear and electronic spin
states).

The last term in Eq. (1) provides the ability to control the quantum spin
states. As with the static terms, for the sake of simplicity we use the
following form

$\hat{H}_{\rm td}=-g\mu_{\rm B}f(t)\vec{b}\hat{\vec{S}}$ (4)

which corresponds to a time-dependent version of the Zeeman interaction term
(3). This Hamiltonian describes the most common spin control techniques, based
on electron paramagnetic resonance instrumentation, in which the absorption of
this signal is the crucial observable. In recent times, it has been shown that
spin qubits, including those in molecules, can also be manipulated by means of
electric field pulses [48, 49]. In this case, the time-dependent perturbation
introduces a modulation of the crystal field and the magnetic anisotropy terms
associated with it. Again, the methods described below are easily adaptable to
these situations.

###  II.2 Coherent control via monochromatic resonant pulses

Often, the temporal shape is a simple monochromatic term, e.g.

$f(t)=\lambda\cos(\omega t+\phi)\Pi_{t_{0}}^{t_{f}}(t)\,,$ (5)

where the amplitude is determined by $\lambda$, and $\Pi_{t_{0}}^{t_{f}}$ is
the rectangular function, that is equal to one if $t_{0}\leq t\leq t_{f}$,
zero otherwise (in reality, of course, the ramp up and down at $t_{0}$ and
$t_{f}$ are not abrupt). If the frequency is chosen to be close to one of the
resonances, and the amplitude is low enough, the rotating wave approximation
(RWA) may be applied, and the effect of these pulses can be worked out
analytically: if $j,k$ are the two levels linked by the resonance (let us
assume a perfect resonance, and that all other frequencies are well
separated), the evolution operator $\hat{U}(t)$ is:

$\hat{U}(t)=\hat{R}_{\vec{n}}^{(jk)}(\theta)\oplus\hat{I}^{(\overline{jk})}\,.$
(6)

This expression assumes the interaction representation, and, in order to
simplify the notation, we have set $t_{0}=0$. The superindex $(jk)$ on the
two-dimensional _rotation_ operator $\hat{R}_{\vec{n}}^{(jk)}(\theta)$ means
that it acts on the subspace spanned by the $j,k$ levels, whereas the rest of
the levels are unaffected ($\hat{I}^{(\overline{jk})}$ is the identity in all
but the $j,k$ levels). Within the basis spanned by the two states, $j,k$, and
the corresponding Pauli matrices $\sigma_{\alpha}$ (in that basis), the
rotation operator is given by:

$\hat{R}_{\vec{n}}(\theta)=\exp\left(-i\frac{\theta}{2}\vec{n}\cdot\hat{\vec{\sigma}}\right)\,.$
(7)

$\hat{\vec{\sigma}}$ is the vector of Pauli matrices. The rotation angle
$\theta$ is $\lambda g\mu_{B}|\mu_{jk}|t$, where $\mu_{jk}=\langle
j|\vec{b}\hat{\vec{S}}|k\rangle$ is the coupling matrix element, and $\vec{n}$
is the unit vector:

$\vec{n}=(\cos(\arg\mu_{jk}+\phi),-\sin(\arg\mu_{jk}+\phi),0)\,.$ (8)

The choice of $\phi$ determines the rotation axis: if $\phi=-\arg\mu_{jk}$, we
have a $R_{X}(\theta)$ rotation; if $\phi=-\arg\mu_{jk}-\pi/2$, we have a
$R_{Y}(\theta)$ rotation. We cannot have direct $R_{Z}(\theta)$ rotations, but
they can however be built by combinations of the former two. (we recall that
here, $X,Y,Z$ do not refer to any spatial direction, but to the Pauli matrices
defined in the basis spanned by the states $j,k$).

Let us consider $R_{X}(\theta)$ rotations in the following. By adjusting the
total pulse time $t_{f}$, one selects the angle $\theta$ and, for example,
performs a $\pi$-rotation, i.e. if

$t_{f}=t_{\pi}^{\lambda}=\frac{\pi}{\lambda g\mu_{B}|\mu_{jk}|}\,$ (9)

the rotation transforms the $j$ state into the $k$ state and viceversa:

$\hat{R}_{X}(\pi)=-i\left[\begin{array}[]{cc}0&1\\\ 1&0\end{array}\right]$
(10)

Note the presence of the $-i$ factor; it is an irrelevant global phase factor
if we consider the $(j,k)$ subspace as isolated, but it changes the phase with
respect to the rest of the levels.

By concatenation of several of these rotations, one may attempt to construct
arbitrary unitaries in any $2^{n}$-level system [27, 29, 22, 18]. Some
specific quantum gates cannot be constructed in this way, however, due
precisely to the presence of the phases mentioned above. This limitation can
be remedied: see [9] for a discussion on this issue, and for an easy solution
via the presence of an extra ancillary level.

The problem that cannot be remedied is the approximate character of all
previous expressions, that rely on the weakness of the perturbation amplitude
$\lambda$, thus avoiding the possibility of arbitrarily speeding up the
process by increasing that amplitude. If we require a minimum fidelity for the
$|j\rangle\to|k\rangle$ transformations, then a minimum amount of time is
necessary.

###  II.3 Quantum Optimal Control Theory (QOCT)

The previous arguments reveal an intrinsic limitation of monochromatic pulses
for the creation of fast quantum gates, which may be further complicated by
the presence of experimental constraints, _i.e._ the inability to increase the
magnetic field intensities. In order to create faster gates, one may attempt
to use more complex temporal shapes, _i.e._ combine various frequencies. We
therefore wonder how this possibility may help the control of molecular spin
qudits, i.e. whether operation times can be made substantially shorter than
the coherence times. For this purpose, we have applied QOCT. We summarize in
the following the method and the basic equations that we have employed.

We consider that the time-dependent pulse-shape function $f(t)$ alluded above
is parameterized with a set of values $u_{0},\dots,u_{M}\equiv u$: i.e.
$f=f(u;t)$. The evolution of the system is then determined by the Hamiltonian:

$\hat{H}(u;t)=\hat{H}_{0}+f(u;t)\hat{V}\,,$ (11)

where $\hat{H}_{0}=\hat{H}_{\rm ZF}+\hat{H}_{\rm Zeeman}$, and
$\hat{V}=-g\mu_{\rm B}\hat{\vec{S}}\vec{b}$.

The evolution operator $\hat{U}(u;t)$ is thus also determined by $u$. In the
interaction representation it evolves according to:

$\displaystyle i\frac{\partial}{\partial t}\hat{U}(u;t)$
$\displaystyle=f(u,t)\hat{\tilde{V}}(t)U(u;t)\,,$ (12)
$\displaystyle\hat{U}(u;0)$ $\displaystyle=\hat{I}\,,$ (13)

where $\hat{\tilde{V}}(t)=\exp(it\hat{H}_{0})\hat{V}\exp(-it\hat{H}_{0})$.

The goal is to find a set of parameters $u^{(0)}$ such that the evolution
operator is equal – or equivalent – to a given target gate $\hat{U}_{\rm G}$:
$\hat{U}(u^{(0)},T)=e^{ia}\hat{U}_{\rm G}$ for any irrelevant global phase
$a$. In the QOCT framework, this is achieved by finding a set of parameters
that leads to the maximization of a functional of the system evolution; in
this case this can be done by defining this functional as:

$F(\hat{U})=|\hat{U}\cdot\hat{U}_{\rm target}|^{2}$ (14)

where the dot product in the space of linear transformations that we have used
is the Fröbenius product:

$\hat{A}\cdot\hat{B}=\frac{1}{d}{\rm Tr}[\hat{A}^{\dagger}\hat{B}]\,.$ (15)

Here, $d$ is the space dimension. The functional thus defined acquires its
maximum value (one) when $\hat{U}$ is equal to the target gate, modulo a phase
factor. Since, as mentioned above, the evolution is determined by the
parameters $u$, the problem is reduced to the maximization of the function:

$G(u)=F\left[\hat{U}(u;t_{f})\right]\,.$ (16)

Many possible algorithms exist for finding the maxima of multivariate
functions such as $G$. Most of them require of a means to compute the gradient
of the function (in addition to the function itself). Optimal control theory
provides the mathematical tool to derive these gradients (essentially,
Pontryagin’s maximum principle [50]). For the case of our function $G$, the
gradient is given by:

$\frac{\partial G}{\partial u_{m}}(u)=2{\rm Im}\int_{0}^{t_{f}}\\!\\!\\!{\rm
d}t\;\frac{\partial f}{\partial
u_{m}}(u;t)\hat{B}(u;t)\cdot(\hat{\tilde{V}}(t)\hat{U}(u;t))\,,$ (17)

where the _costate_ $\hat{B}$ is defined by the following equations:

$\displaystyle i\frac{\partial}{\partial t}\hat{B}(u;t)$
$\displaystyle=f(u,t)\hat{\tilde{V}}(t)U(u;t)\,,$ (18)
$\displaystyle\hat{B}(u;t_{f})$ $\displaystyle=(\hat{U}_{\rm
target}\cdot\hat{U}(t_{f}))\hat{U}_{\rm target}\,.$ (19)

Note that it is an equation of motion similar to the one that determines the
evolution operator itself, except the boundary condition is given at the final
time of the propagation $t_{f}$ (it is a _final_ condition, instead of an
initial condition). In consequence, the computation of the gradient requires
two propagations: a forward propagation for $\hat{U}$, and a backward
propagation for $\hat{B}$.

It remains to specify the parameterization of $f$, an important task that
actually defines the set of allowed temporal shapes. This should be done with
the experimental capabilities in mind. In our case, we have opted for a simple
Fourier expansion:

$\displaystyle f(u,t)=\frac{1}{\sqrt{t_{f}}}u_{0}+\sum_{k=1}^{K}\left[\right.$
$\displaystyle u_{2k}\frac{2}{\sqrt{t_{f}}}\cos(\omega_{k}t)+$ $\displaystyle
u_{2k-1}\frac{2}{\sqrt{t_{f}}}\sin(\omega_{k}t)\left.\right]\,.$ (20)

The frequencies $\omega_{k}=2\pi k/t_{f}$, $k=1,\dots,K$ have a maximum
_cutoff_ value at $2\pi K/t_{f}$, that must be chosen big enough to include
the relevant natural frequencies of the spin qudit, but not so large that it
cannot be handled experimentally.

Some constraints have to be imposed on the allowed values for the parameters:
the pulse amplitude must start and end at zero: $f(u,0)=f(u,t_{f})=0$, which
translates into $\sum_{k=1}^{K}u_{2k}=0$. We have also imposed a zero value
for the average amplitude, $\int_{0}^{t_{f}}\\!\\!{\rm d}t\;f(u,t)=0$, which
translates into $u_{0}=0$. Finally, the generated magnetic field cannot have
arbitrary amplitudes. Therefore, in the calculations discussed below we have
set a maximum value. All these constraints have been added to the optimization
algorithm.

In order to implement these equations, we have used the qutip code as a base
[51, 52]. We have, however, not employed the QOCT algorithms provided by this
platform (at the time of writing of this article), but used the gradient,
computed as in Eq. (17), to feed our own QOCT code [53], that then utilizes a
general purpose function optimization algorithm: the Sequential Least-Squares
Quadratic Programming (SLSQP) algorithm [54] as implemented in the NLOPT
library [55].

##  III  Results

###  III.1 The GdW30 molecular spin qudit

In order to explore the potential of QOCT for the control of spin qudits, we
have chosen a system, GdW30 [37, 27], which is both well characterized and
relatively simple. The molecular structure of this polyoxometalate cluster is
shown in Fig. 1. It hosts a single Gd3+ ion and forms crystals with all
molecules oriented in the same manner. In addition, magnetically diluted
crystals can be grown by simply replacing Gd3+ with Y3+, which is chemically
equivalent but diamagnetic. This allows enhancing spin coherence times up to
$2-3$ $\mu$s [37] while keeping the possibility of orienting the magnetic
fields $\vec{B}$ and $\vec{b}$ along specific molecular axes. This molecule
shows a hard magnetic axis along the main molecular axis $z$. Its static spin
Hamiltonian can be well described by Eqs. (1), (2) and (3) with $S=7/2$,
$g=2,$ $D=1281$ MHz and $E=294$ MHz. The overall splitting of the $d=8$
multiplet is smaller than $1$ K, or $20.8$ GHz, thus ensuring that adjacent
level splittings lie within the reach of conventional EPR as well as of other
microwave technologies. In all calculations discussed below, we have chosen
the static magnetic field $\vec{B}$ to point in the $x$ direction (medium
axis) and set $B=0.15$ T. Under such conditions, the eigenstates of
$\hat{H}_{\rm ZF}+\hat{H}_{\rm Zeeman}$ become close approximations to pure
spin projections along $\vec{B}$. The time-dependent magnetic field $\vec{b}$
is perpendicular, and points in the $y$ direction (easy axis), thus inducing
transitions between adjacent energy levels.

###  III.2 Transition implementation via monochromatic pulses

In this section, we discuss the manipulation of the GdW30 spin using
monochromatic pulses resonant with the set of allowed transitions mentioned
above and shown in Fig. 1. As an illustration, we consider the application of
$\pi$ pulses linking every two of these states. Fig. 2 (top) displays the
transformation infidelities (i.e.
$1-|\langle\psi(t_{\pi}^{\lambda})|k\rangle|^{2}$) as a function of the time
$t_{\pi}^{\lambda}$ allocated to complete the operation, for the seven
$|j\rangle\to|k=j+1\rangle$ transitions. One can see that the error in the
outcome state is reduced as $t_{\pi}^{\lambda}\to\infty$. In fact, from the
logarithmic plot one may infer a quadratic behaviour:

$1-|\langle\psi(t_{\pi}^{\lambda})|k\rangle|^{2}=\mathcal{O}((1/t_{\pi}^{\lambda})^{2})\,.$
(21)

It becomes clear that, in order to ensure a given fidelity, a minimum time
(or, equivalently, a maximum amplitude) is required.

Figure 2:  Top: Infidelities of the seven main transitions in the GdW30
molecule, as function of the $\pi$-pulse time. Bottom: Infidelities of the
seven main transitions in the GdW30 molecule, as function of the total pulse
time, for pulses obtained with optimal control. Inset: Time-dependent shape of
the pulses used to generate two of those $|6\rangle$ to $|7\rangle$
transitions: one $\pi$-pulse (black), and one optimized pulse (red),
corresponding to the thick black and red dots, respectively, on the pink
lines.

###  III.3 Optimization of state-to-state transitions

In this and the following section, we apply the QOCT methods described in
section II to quantum operations on GdW30 having different targets. First, we
optimize resonant transitions and sequences of these. Then, we tackle the
optimization of quantum gates. The goal here is to see how OCT permits to
increase the fidelities shown in Fig. 2 (top), by allowing for the presence in
the pulse of other frequencies, besides the resonant one. We have performed
QOCT calculations 111 Note that the OCT for the simpler problem of state
population – in contrast to the the harder problem of the creation of a given
evolution operator or gate – requires equations that are slightly different to
the ones described in Section II. considering, for each transition, the same
total propagation time used to create Fig. 2 (top). Each of these propagation
times $t_{\pi}^{\lambda}$ corresponds to a $\pi$-pulse amplitude $\lambda$
[Eq. (9)], that we have used now to set a bound for the Fourier expansion
coefficients in the OCT maximization: the temporal shape of the microwave
field is given by Eq. (20), where
$|\frac{2u_{j}}{\sqrt{t_{\pi}^{\lambda}}}|\leq\lambda$. Having in mind the
level splitting present in GdW30 and typical experimental capabilities, we
have set the frequency cutoff at 8 GHz. The resulting (in)-fidelities are
displayed in Fig. 2 (bottom). The shaped pulses permit to decrease those
infidelities with respect to the simple $\pi-$pulse values, down to negligible
values for all but the shortest total transition times (we have set a 10-7
threshold to stop the search algorithm, and hence the flat curves for the
longer times). In Fig. 2 (inset), we also compare the shape of control pulses
corresponding to a resonant transition (monochromatic $\pi$-pulse) and to the
optimized one from states $|6\rangle$ to $|7\rangle$. These examples
correspond to duration values marked by (respectively) black (top) and red
dots (bottom). It is clear that the optimized pulse achieves a much better
fidelity in a much shorter time.

Figure 3:  Infidelities for the $|0\rangle\to|7\rangle$ and
$|0\rangle\to\frac{1}{\sqrt{2}}(|0\rangle-i|7\rangle)$ transitions, as a
function of operation time, for the pulses obtained as a sequence of
monochromatic resonant pulses ($\pi$ or $\pi/2$), and for the optimal pulses
obtained with QOCT.

Let us consider now the more general case of optimizing transitions between
states, say $|n\rangle$ and $|m\rangle$, that are not directly coupled by the
external field, i.e. $\langle n|S_{y}|m\rangle=0$. A possible solution is to
concatenate a series of $\pi$-pulses, between intermediate states. In fact,
this is a criterion for universality: if any two states can be connected
through others, the qudit can perform any unitary and, in this sense, can be
regarded as a universal quantum processor. However, the time needed is
proportional to the number of pulses, which in practice is a limitation. Thus,
in this case, QOCT offers a clear advantage by replacing a single shaped pulse
to achieve, and accelerate, a process that would otherwise require a sequence
of $7$ monochromatic pulses.

In Fig. 3 we have tested this by comparing the infidelity in the transition
$|0\rangle\to|7\rangle$, using a sequence of $\pi$-pulses between adjacent
levels, $|k\rangle\to|k+1\rangle$, and using QOCT. The improvement is quite
significant. In particular, we have run over a range of amplitudes $\lambda$,
that determine the $\pi$-pulse length for each transition,
$t_{\pi}^{\lambda}(k\to k+1)$. The full $|0\rangle\to|7\rangle$ process then
requires $t_{f}=\sum_{k=0}^{6}t_{\pi}^{\lambda}(k\to k+1)$. The plot displays
the fidelity achieved with these pulse sequences. Then, for each of those
times, we have performed QOCT calculations, once again setting a bound for the
amplitudes of the individual Fourier terms equal to $\lambda$. The plot shows
how, even at very short operation times, the fidelities achieved by the
optimized pulses are almost equal to one. In terms of time scales, a $99$ %
fidelity can be achieved in less than $10$ ns, thus much shorter than
$T_{2}\simeq 2$ $\mu$s, while reaching the same result with a sequence of
monochromatic pulses would take more than $1$ $\mu$s.

In Fig 3 we also show results of a similar calculation, but using the
$|0\rangle\to\frac{1}{\sqrt{2}}(|0\rangle-i|7\rangle)$ state as target, a
superposition state that can be reached with a $\pi/2$-pulse corresponding to
the $|0\rangle\to|1\rangle$ transition, followed by the same previous sequence
of $\pi$-pulses that raises the state through the next adjacent levels. The
results are qualitatively similar to the ones obtained for the full
$|0\rangle$ to $|7\rangle$ transition, thus showing that the speed enhancement
achieved by the application of QOCT is not restricted to any particular class
of transitions. This allows targeting the optimization of complex gates, which
is discussed next.

Figure 4:  Optimal pulse in the time domain (top), and its power spectrum in
the frequency domain (bottom), of the optimal pulse obtained for the
$\theta=\pi/4$ Deutsch gate.

###  III.4 Quantum gate optimization

Finally, we proceed to our true objective: the search for non-trivial pulse
shapes that realize quantum gates, with high fidelities, in short times. As
target gates, we have chosen the family of Deutsch gates [57, 58]:

$D(\theta)=\left[\begin{array}[]{ccc}I_{6}&0&0\\\
0&i\cos(\theta)&\sin(\theta)\\\
0&\sin(\theta)&i\cos(\theta)\end{array}\right]\,,$ (22)

where $I_{6}$ is the $6\times 6$ identity matrix. Note that this family
includes the Toffoli gate, for $\theta=\pi/2$. The reason for focusing on this
set of gates is that it is universal: any circuit can be constructed by
combination of these components.

The total propagation time $t_{\rm f}$ has been set to $20$ times the maximum
natural period of the field-free Hamiltonian, i.e. the period corresponding to
the smallest transition frequency. For the choice of static magnetic field
used here ($\vec{B}=0.15T\vec{e}_{x}$), $t_{\rm f}\sim 10$ ns.

As discussed above, the optimizations are performed constraining the allowed
parameter set, such that each sinusoidal (or cosinusoidal) term in the
expansion (20) has a maximum amplitude, i.e.

$|\frac{2u_{k}}{\sqrt{t_{f}}}|b\leq b_{\rm max}\,.$ (23)

In the calculations shown here, we have set $b_{\rm max}$ to a relatively high
value of $20$ mT in order to make the operation, and therefore also the
computational, times manageably short. The influence of $b_{\rm max}$ on
$t_{\rm f}$ is discussed below. The optimization algorithm is an iterative
process that we stop when the quality of the gate, measured as
$F(\hat{U}(u;t_{f}))=|\hat{U}(u,t_{f})\cdot\hat{D}(\theta)|^{2}$, reaches a
certain threshold, that for these calculations we have set to $0.99$.

Figure 4 displays the results obtained for $\theta=\pi/4$ as an example (the
results obtained for other angles are qualitatively similar). The top panel
shows $f(u^{(0)};t)$ in real time, whereas the bottom panels displays its
power spectrum. Both plots demonstrate the complexity of the pulses, that do
not have dominant frequencies.

In the previous example, we have set a total operation time $t_{\rm f}$ and an
amplitude bound $b_{\rm max}$. Those magnitudes are of course related: if, for
a given $t_{f}$, we set a too low amplitude bound, the quality of the gate (as
measured by Eq. (14) will also be too low. In fact, if we fix a threshold for
acceptable gate quality (say, $0.99$), for each given propagation time there
will be a minimum amplitude bound necessary for the QOCT algorithm to return a
successful pulse. Even with optimal pulse shapes, we need a minimum of field
amplitudes in order to get a high quality gate. We have therefore studied this
issue, computing the minimum amplitude bound that can be used to constrain the
QOCT calculation in order to get a given gate, as a function of operation
time. The results are shown in Fig. 5.

This type of plot helps to ascertain whether or not the gate operations are
experimentally feasible, as in practice there is a hardware bound on the field
amplitudes that can be used. Given this limit, one may learn from the plot
what operation times are feasible, even with shaped optimal pulses. Obviously,
the lower the available amplitudes, the longer the operations times must be.

Figure 5:  Minimal bound on the amplitudes ($y$-axis) that allows to obtain a
$\theta=\pi/2$ Deutsch gate for a given total operation time ($x$-axis).

##  IV  Conclusions

The previous results show that the use of complex pulses, engineered by
optimal control techniques, provides a method for improving the speed of
operations performed on spin qudits. The advantages are already noticeable for
the realization of elementary transitions, when the application of
monochromatic pulses is limited by the need of keeping the excitation
amplitudes sufficiently low. Yet, they become even more important when dealing
with more complex operations. Then, QOCT allows reaching any given fidelity of
the outcome wave function with a single control pulse replacing the often long
sequence of resonant pulses required by standard techniques. This possibility
is especially relevant for algorithms that involve transitions between
relatively disconnected states. For a necessarily limited coherence spin time,
this difference can represent a big gain in the performance of such protocols.

The molecular qudit design and the way information is encoded on its spin
states can be adapted to suit best the requirements of specific quantum
protocols [18, 59], thus offering a vast choice of possible molecular
platforms and implementations. Optimal control techniques described in this
work are flexible enough to be made compatible with almost any of them.
Although we have here considered a specific molecule for illustrative
purposes, QOCT can deal with any spin model and with diverse control
interactions (e.g. magnetic or electric field pulses). Therefore, it can be
adapted to boost the implementation of diverse algorithms. We feel that it
will be of special relevance to quantum error correction, because reaching a
fidelity improvement with such protocols critically depends on the ratio
between the implementation time and $T_{2}$. With an additional computational
cost, one can even consider optimizing the control pulses to best compensate
for the actual sources of decoherence in each molecule, mainly dephasing by
nuclear spins located in the ligand shell surrounding the magnetic core.

The use of more sophisticated control pulses represents a challenge to
experimental implementations. Most commercial EPR systems work with relatively
narrow excitation bands and a reduced choice of pulse shapes. In recent years
such systems have been complemented with waveguide generators able to
arbitrarily design the excitation pulses [60, 61, 62]. Still, these set-ups
are limited to frequencies lying sufficiently close to the cavity resonance
frequency. In order to expand the frequency window, one can resort to on-chip
circuits, with the excitation being driven by an open transmission line. This
scheme, illustrated by Fig. 1, has been applied to investigate coherent
control of NV centers in diamond [63] and to perform broadband spectroscopy of
GdW30 and other molecular spin qudits [27, 22]. It can also be used to read
out the outcome, either by looking at the frequency dependent absorption, in a
projective measurement, or by coupling it to a superconducting resonator that
can perform nondemolition, dispersive measurements of the qudit states [64].
Some experimental systems combining superconducting resonators and broadband
control lines have been recently reported [65]. Another promising
implementation is based on the combination of single-molecule electronics with
gates or coils able to locally generate arbitrarily shaped electric or
magnetic microwave pulses. Experiments performed on molecules trapped between
point contacts or between a metal substrate and a STM tip have provided the
first measurements of spin coherence in individual molecules [66] and achieved
the realization of Grover’s search algorithm using three nuclear spin states
in a Tb-based molecule [28].

In summary, the application of quantum optimal control theory to operate the
states of molecular spin qudits offers remarkable prospects to improve their
performance, compensating for their not too long coherence times. Equipped
with these techniques, many more molecular systems and applications can become
feasible, thus contributing to an alternative and promising path towards large
scale quantum computation and simulation.

###### Acknowledgements.

We acknowledge the financial support from grants FIS2017-82426-P,
PGC2018-094792-B-I00, RTI2018-096075-B-C21 and PCI2018-093116 funded by
MCIN/AEI/ 10.13039/501100011033 and “ERDF A way of making Europe”, grant
PID2020-115221GB-C41/AEI/10.13039/501100011033, the European Union’s Horizon
2020 research and innovation programme (QUANTERA project SUMO, FET-OPEN grant
862893 FATMOLS), the Gobierno de Aragón grant E09-17R-Q-MAD and the CSIC
Quantum Technology Platform PT-001.

## References

  * Arute and et al [2019] F. Arute and et al, Quantum supremacy using a programmable superconducting processor, Nature 574, 505– (2019).
  * Wu  _et al._ [2021] Y. Wu, W.-S. Bao, S. Cao, F. Chen, M.-C. Chen, X. Chen, T.-H. Chung, H. Deng, Y. Du, D. Fan, M. Gong, C. Guo, C. Guo, S. Guo, L. Han, L. Hong, H.-L. Huang, Y.-H. Huo, L. Li, N. Li, S. Li, Y. Li, F. Liang, C. Lin, J. Lin, H. Qian, D. Qiao, H. Rong, H. Su, L. Sun, L. Wang, S. Wang, D. Wu, Y. Xu, K. Yan, W. Yang, Y. Yang, Y. Ye, J. Yin, C. Ying, J. Yu, C. Zha, C. Zhang, H. Zhang, K. Zhang, Y. Zhang, H. Zhao, Y. Zhao, L. Zhou, Q. Zhu, C.-Y. Lu, C.-Z. Peng, X. Zhu, and J.-W. Pan, Strong quantum computational advantage using a superconducting quantum processor, Phys. Rev. Lett. 127, 180501 (2021).
  * Xao  _et al._ [2021] X. Xao _et al._ , CMOS-based cryogenic control of silicon quantum circuits, Nature 593, 205 (2021).
  * Gottesman  _et al._ [2001] D. Gottesman, A. Kitaev, and J. Preskill, Encoding a qubit in an oscillator, Phys. Rev. A 64, 012310 (2001).
  * Leuenberger and Loss [2001] M. N. Leuenberger and D. Loss, Quantum computing in molecular magnets, Nature 410, 789 (2001).
  * Brennen  _et al._ [2005] G. K. Brennen, D. P. O’Leary, and S. S. Bullock, Criteria for exact qudit universality, Phys. Rev. A 71, 052318 (2005).
  * Pirandola  _et al._ [2008] S. Pirandola, S. Mancini, S. L. Braunstein, and D. Vitali, Minimal qudit code for a qubit in the phase-damping channel, Phys. Rev. A 77, 032309 (2008).
  * Lanyon  _et al._ [2009] B. P. Lanyon, M. Barbieri, M. P. Almeida, T. Jennewein, T. C. Ralph, K. J. Resch, G. J. Pryde, J. L. O’Brien, G. Alexei, and A. G. White, Simplifying quantum logic using higher-dimensional hilbert spaces, Nature Physics 5, 134 (2009).
  * Kiktenko  _et al._ [2015] E. Kiktenko, A. Fedorov, A. Strakhov, and V. Man’ko, Single qudit realization of the deutsch algorithm using superconducting many-level quantum circuits, Physics Letters A 379, 1409 (2015).
  * Lapkiewicz  _et al._ [2011] R. Lapkiewicz, P. Li, C. Schaeff, N. K. Langford, S. Ramelow, M. Wieśniak, and A. Zeilinger, Experimental non-classicality of an indivisible quantum system, Nature 474, 490 (2011).
  * Ringbauer  _et al._ [2021] M. Ringbauer, M. Meth, L. Postler, R. Stricker, R. Blatt, P. Schindler, and T. Monz, A universal qudit quantum processor with trapped ions (2021), arXiv:2109.06903 [quant-ph] .
  * Asaad  _et al._ [2020] S. Asaad, V. Mourik, B. Joecker, M. A. I. Johnson, A. D. Baczewski, H. R. Firgau, M. T. Madzik, V. Schmitt, J. J. Pla, F. E. Hudson, K. M. Itoh, J. C. McCallum, A. S. Dzurak, A. Laucht, and A. Morello, Coherent electrical control of a single high-spin nucleus in silicon, Nature 579, 205 (2020).
  * Neeley  _et al._ [2009] M. Neeley, M. Ansmann, R. C. Bialczak, M. Hofheinz, E. Lucero, A. D. O’Connell, D. Sank, H. Wang, J. Wenner, A. N. Cleland, M. R. Geller, and J. M. Martinis, Emulation of a quantum spin with a superconducting phase qudit, Science 325, 722 (2009).
  * Gatteschi and Sessoli [2003] D. Gatteschi and R. Sessoli, Quantum tunneling of magnetization and related phenomena in molecular materials, Angewandte Chemie International Edition 42, 268 (2003).
  * Aromí  _et al._ [2012] G. Aromí, D. Aguilà, F. Luis, S. Hill, and E. Coronado, Design of magnetic coordination complexes for quantum computing., Chem. Soc. Rev. 41, 537 (2012).
  * Atzori and Sessoli [2019] M. Atzori and R. Sessoli, The second quantum revolution: Role and challenges of molecular chemistry, J. Am. Chem. Soc. 141, 11339 (2019).
  * Gaita-Ari$\tilde{\text{n}}$o  _et al._ [2019] A. Gaita-Ari$\tilde{\text{n}}$o, F. Luis, S. Hill, and E. Coronado, Molecular spins for quantum computation., Nature Chem. 11, 301 (2019).
  * Carretta  _et al._ [2021] S. Carretta, D. Zueco, A. Chiesa, Á. Gómez-León, and F. Luis, A perspective on scaling up quantum computation with molecular spins, Applied Physics Letters 118, 240501 (2021).
  * Moreno-Pineda  _et al._ [2017] E. Moreno-Pineda, M. Damjanović, O. Fuhr, W. Wernsdorfer, and M. Ruben, Nuclear spin isomers: Engineering a et4n[dypc2] spin qudit, Angewandte Chemie International Edition 56, 9915 (2017).
  * Moreno-Pineda  _et al._ [2018] E. Moreno-Pineda, C. Godfrin, F. Balestro, W. Wernsdorfer, and M. Ruben, Molecular spin qudits for quantum algorithms, Chem. Soc. Rev. 47, 501 (2018).
  * Hussain  _et al._ [2018] R. Hussain, G. Allodi, A. Chiesa, E. Garlatti, D. Mitcov, A. Konstantatos, K. S. Pedersen, R. De Renzi, S. Piligkos, and S. Carretta, Coherent manipulation of a molecular ln-based nuclear qudit coupled to an electron qubit, Journal of the American Chemical Society 140, 9814 (2018).
  * Gimeno  _et al._ [2021] I. Gimeno, A. Urtizberea, J. Román-Roche, D. Zueco, A. Camón, P. J. Alonso, O. Roubeau, and F. Luis, Broad-band spectroscopy of a vanadyl porphyrin: a model electronuclear spin qudit, Chem. Sci. , (2021).
  * Chicco  _et al._ [2021] S. Chicco, A. Chiesa, G. Allodi, E. Garlatti, M. Atzori, L. Sorace, R. De Renzi, R. Sessoli, and S. Carretta, Controlled coherent dynamics of [vo(tpp)], a prototype molecular nuclear qudit with an electronic ancilla, Chem. Sci. 12, 12046 (2021).
  * Luis  _et al._ [2011] F. Luis, A. Repollés, M. J. Martínez-Pérez, D. Aguilà, O. Roubeau, D. Zueco, P. J. Alonso, M. Evangelisti, A. Camón, J. Sesé, L. A. Barrios, and G. Aromí, Molecular prototypes for spin-based cnot and swap quantum gates., Phys. Rev. Lett. 107, 117203 (2011).
  * Aguilà  _et al._ [2014] D. Aguilà, D. Barrios, V. Velasco, O. Roubeau, A. Repollés, P. Alonso, J. Sesé, S. Teat, F. Luis, and G. Aromí, Heterodimetallic [lnln’] lanthanide complexes: Toward a chemical design of two-qubit molecular spin quantum gates., J. Am. Chem. Soc. 136, 14215 (2014).
  * Ferrando-Soria  _et al._ [2016] J. Ferrando-Soria, E. Moreno-Pineda, A. Chiesa, A. Fernández, S. A. Magee, S. Carretta, P. Santini, I. Vitorica-Yrezabal, F. Tuna, E. J. L. McInness, and R. E. P. Winpenny, A modular design of molecular qubits to implement universal quantum gates., Nat. Commun. 7, 11377 (2016).
  * Jenkins  _et al._ [2017] M. D. Jenkins, Y. Duan, B. Diosdado, J. J. García-Ripoll, A. Gaita-Ariño, C. Giménez-Saiz, P. J. Alonso, E. Coronado, and F. Luis, Coherent manipulation of three-qubit states in a molecular single-ion magnet, Phys. Rev. B 95, 064423 (2017).
  * Godfrin  _et al._ [2017] C. Godfrin, A. Ferhat, R. Ballou, S. Klyatskaya, M. Ruben, W. Wernsdorfer, and F. Balestro, Operating quantum states in single magnetic molecules: Implementation of grover’s quantum algorithm, Phys. Rev. Lett. 119, 187702 (2017).
  * Luis  _et al._ [2020] F. Luis, P. J. Alonso, O. Roubeau, V. Velasco, D. Zueco, D. Aguilà, J. I. Martínez, L. A. Barrios, and G. Aromí, A dissymmetric [gd2] coordination molecular dimer hosting six addressable spin qubits, Communications Chemistry 3, 176 (2020).
  * Macaluso  _et al._ [2020] E. Macaluso, M. Rubín, D. Aguilà, A. Chiesa, J. I. M. L. A. Barrios, P. J. Alonso, O. Roubeau, F. Luis, G. Aromí, and S. Carretta, A heterometallic [lnln’ln] lanthanide complex as a qubit with embedded quantum error correction, Chem. Sci. 11, 10337 (2020).
  * Tacchino  _et al._ [2021] F. Tacchino, A. Chiesa, R. Sessoli, I. Tavernelli, and S. Carretta, A proposal for using molecular spin qudits as quantum simulators of light–matter interactions, J. Mater. Chem. C 9, 10266 (2021).
  * Chiesa  _et al._ [2020] A. Chiesa, E. Macaluso, F. Petiziol, S. Wimberger, P. Santini, and S. Carretta, Molecular nanomagnets as qubits with embedded quantum-error correction, The Journal of Physical Chemistry Letters 11, 8610 (2020).
  * Chiesa  _et al._ [2021] A. Chiesa, F. Petiziol, E. Macaluso, S. Wimberger, P. Santini, and S. Carretta, Embedded quantum-error correction and controlled-phase gate for molecular spin qubits, AIP Advances 11, 025134 (2021).
  * Godfrin  _et al._ [2018] C. Godfrin, R. Ballou, E. Bonet, M. Ruben, S. Klyatskaya, W. W., and F. Balestro, Generalized Ramsey interferometry explored with a single nuclear spin qudit, npj Quantum Inf. 4, 53 (2018).
  * Bader  _et al._ [2014] K. Bader, D. Dengler, S. Lenz, B. Endeward, S.-D. Jiang, P. Neugebauer, and J. van Slageren, Room temperature quantum coherence in a potential molecular qubit, Nat. Commun 5, 5304 (2014).
  * Zadrozny  _et al._ [2015] J. M. Zadrozny, J. Niklas, O. G. Poluektov, and D. E. Freedman, Millisecond coherence time in a tunable molecular electronic spin qubit., ACS Cent. Sci. 1, 488 (2015).
  * Martínez-Pérez  _et al._ [2012] M. J. Martínez-Pérez, S. Cardona-Serra, C. Schlegel, F. Moro, P. J. Alonso, H. Prima-García, J. M. Clemente-Juan, M. Evangelisti, A. Gaita-Ariño, J. Sesé, J. van Slageren, E. Coronado, and F. Luis, Gd-based single-ion magnets with tunable magnetic anisotropy: Molecular design of spin qubits, Phys. Rev. Lett. 108, 247213 (2012).
  * Brif  _et al._ [2010] C. Brif, R. Chakrabarti, and H. Rabitz, Control of quantum phenomena: past present and future, New Journal of Physics 12, 075008 (2010).
  * Palao and Kosloff [2002] J. P. Palao and R. Kosloff, Quantum computing by an optimal control algorithm for unitary transformations, Phys. Rev. Lett. 89, 188301 (2002).
  * Brion  _et al._ [2006] E. Brion, D. Comparat, and G. Harel, Implementation of a cnot gate in two cold rydberg atoms by the nonholonomic control technique, The European Physical Journal D - Atomic, Molecular, Optical and Plasma Physics 38, 381 (2006).
  * Schirmer [2009] S. Schirmer, Implementation of quantum gates via optimal control, Journal of Modern Optics 56, 831 (2009).
  * Reich  _et al._ [2012] D. M. Reich, M. Ndong, and C. P. Koch, Monotonically convergent optimization in quantum control using krotov’s method, The Journal of Chemical Physics 136, 104103 (2012).
  * Hou  _et al._ [2014] S. Hou, L. Wang, and X. Yi, Realization of quantum gates by lyapunov control, Physics Letters A 378, 699 (2014).
  * Chou  _et al._ [2015] Y. Chou, S.-Y. Huang, and H.-S. Goan, Optimal control of fast and high-fidelity quantum gates with electron and nuclear spins of a nitrogen-vacancy center in diamond, Phys. Rev. A 91, 052315 (2015).
  * Arai and Ohtsuki [2015] K. Arai and Y. Ohtsuki, Reduced-dynamics approach for optimally designing unitary transformations, Phys. Rev. A 92, 062302 (2015).
  * Dong  _et al._ [2016] D. Dong, C. Wu, C. Chen, B. Qi, I. R. Petersen, and F. Nori, Learning robust pulses for generating universal quantum gates, Scientific Reports 6, 36090 (2016).
  * Wilson  _et al._ [2007] A. Wilson, E.-C. Yang, D. N. Hendrickson, and S. Hill, On the validity of the giant spin approximation and its application to single-molecule magnets, Polyhedron 26, 2065 (2007), proceedings of the 10th International Conference on Molecule-based Magnets (ICMM 2006).
  * Thiele  _et al._ [2014] S. Thiele, F. Balestro, R. Ballou, S. Klyatskaya, M. Ruben, and W. Wernsdorfer, Electrically driven nuclear spin resonance in single-molecule magnets, Science 344, 1135 (2014).
  * Liu  _et al._ [2019] J. Liu, J. Mrozek, W. K. Myers, G. A. Timco, R. E. P. Winpenny, B. Kintzel, W. Plass, and A. Ardavan, Electric field control of spins in molecular magnets, Phys. Rev. Lett. 122, 037202 (2019).
  * Pontryagin  _et al._ [1962] L. S. Pontryagin, V. G. Boltyanskii, R. V. Gamkrelidze, and E. F. Mishechenko,  _The Mathematical Theory of Optimal Processes_ (John Wiley & Sons, New York/London, 1962).
  * Johansson  _et al._ [2012] J. Johansson, P. Nation, and F. Nori, Qutip: An open-source python framework for the dynamics of open quantum systems, Computer Physics Communications 183, 1760 (2012).
  * Johansson  _et al._ [2013] J. Johansson, P. Nation, and F. Nori, Qutip 2: A python framework for the dynamics of open quantum systems, Computer Physics Communications 184, 1234 (2013).
  * [53] A. Castro et. al., qocttools, https://acbarrigon.gitlab.io/qocttools/.
  * Kraft [1994] D. Kraft, Algorithm 733: Tomp–fortran modules for optimal control calculations, ACM Trans. Math. Softw. 20, 262–281 (1994).
  * [55] S. G. Johnson, The nlopt nonlinear-optimization package, http://github.com/stevengj/nlopt.
  * Note [1] Note that the OCT for the simpler problem of state population – in contrast to the the harder problem of the creation of a given evolution operator or gate – requires equations that are slightly different to the ones described in Section II.
  * Deutsch and Penrose [1989] D. E. Deutsch and R. Penrose, Quantum computational networks, Proceedings of the Royal Society of London. A. Mathematical and Physical Sciences 425, 73 (1989).
  * Shi [2018] X.-F. Shi, Deutsch, toffoli, and cnot gates via rydberg blockade of neutral atoms, Phys. Rev. Applied 9, 051001 (2018).
  * Petiziol  _et al._ [2021] F. Petiziol, A. Chiesa, S. Wimberger, P. Santini, and S. Carretta, Counteracting dephasing in molecular nanomagnets by optimized qudit encodings, npj Quantum Information 7, 10.1038/s41534-021-00466-3 (2021).
  * Tseitlin  _et al._ [2011] M. Tseitlin, R. W. Quine, G. A. Rinard, S. S. Eaton, and G. R. Eaton, Digital epr with an arbitrary waveform generator and direct detection at the carrier frequency, J. Mag. Res. 213, 119 (2011).
  * Spindler  _et al._ [2016] P. E. Spindler, P. Schöps, A. M. Bowen, B. Endeward, and T. F. Prisner, Shaped pulses in epr, in  _eMagRes_ (John Wiley & Sons, Ltd, 2016) pp. 1477–1492, https://onlinelibrary.wiley.com/doi/pdf/10.1002/9780470034590.emrstm1520 .
  * Prisner [2019] T. F. Prisner, Shaping epr: Phase and amplitude modulated microwave pulses, Journal of Magnetic Resonance 306, 98 (2019).
  * de Lange  _et al._ [2010] G. de Lange, Z. H. Wang, D. Ristè, V. V. Dobrovitski, and R. Hanson, Universal dynamical decoupling of a single solid-state spin from a spin bath, Science 330, 60–63 (2010).
  * Álvaro Gómez-León  _et al._ [2021] Álvaro Gómez-León, F. Luis, and D. Zueco, Dispersive readout of molecular spin qudits (2021), arXiv:2109.14639 [quant-ph] .
  * Keyser  _et al._ [2020] A. K. Keyser, J. J. Burnett, S. E. Kubatkin, A. V. Danilov, M. Oxborrow, S. E. de Graaf, and T. Lindström, Pulsed electron spin resonance of an organic microcrystal by dispersive readout, Journal of Magnetic Resonance 321, 106853 (2020).
  * Zhang  _et al._ [2021] X. Zhang, C. Wolf, Y. Wang, H. Aubin, T. Bilgeri, P. Willke, A. J. Heinrich, and T. Choi, Electron spin resonance of single iron phthalocyanine molecules and role of their non-localized spins in magnetic interactions, Nature Chemistry 10.1038/s41557-021-00827-7 (2021).


#  Matrix discrepancy and the log-rank conjecture

Benny Sudakov, István Tomon ETH Zurich, _e-mail_ :
<EMAIL_ADDRESS>Research supported in part by SNSF grant
200021_196965.Umeå University, _e-mail_<EMAIL_ADDRESS>
###### Abstract

Given an $m\times n$ binary matrix $M$ with $|M|=p\cdot mn$ (where $|M|$
denotes the number of 1 entries), define the _discrepancy_ of $M$ as
$\operatorname{disc}(M)=\displaystyle\max_{X\subset[m],Y\subset[n]}\big{|}|M[X\times
Y]|-p|X|\cdot|Y|\big{|}$. Using semidefinite programming and spectral
techniques, we prove that if $\operatorname{rank}(M)\leq r$ and $p\leq 1/2$,
then

$\operatorname{disc}(M)\geq\Omega(mn)\cdot\min\left\\{p,\frac{p^{1/2}}{\sqrt{r}}\right\\}.$

We use this result to obtain a modest improvement of Lovett’s best known upper
bound on the log-rank conjecture. We prove that any $m\times n$ binary matrix
$M$ of rank at most $r$ contains an $(m\cdot 2^{-O(\sqrt{r})})\times(n\cdot
2^{-O(\sqrt{r})})$ sized all-1 or all-0 submatrix, which implies that the
deterministic communication complexity of any Boolean function of rank $r$ is
at most $O(\sqrt{r})$.

##  1 Introduction

The log-rank conjecture, proposed by Lovász and Saks [11] in 1988, is one of
the fundamental open problems in communication complexity. It states that for
any Boolean function $f:X\times Y\rightarrow\\{-1,1\\}$ of rank $r$, its
deterministic communication complexity $\mbox{CC}^{det}(f)$ is bounded by
$\operatorname{polylog}(r)$. We refer the reader to Kushilevitz and Nisan [6]
for exact definitions, or the recent survey of Lee and Shraibman [9] for a
detailed overview of the problem. The log-rank conjecture has a number of
combinatorial interpretations, which will be the main focus of our paper.

Lovász and Saks [11] observed (see also [12]) that this problem is closely
related to bounding the chromatic number of graphs, whose adjacency matrix has
bounded rank. Indeed, if $c(r)$ denotes the maximum of $\mbox{CC}^{det}(f)$
among every $f$ of rank $r$, then $\chi(G)\leq 2^{c(r+1)}$ for every graph
$G$, whose adjacency matrix has rank at most $r$. Furthermore, the log-rank
conjecture is related to finding monochromatic submatrices in low-rank binary
matrices. To this end, let $\alpha(r)>0$ be such that any binary matrix
$M\in\\{0,1\\}^{m\times n}$ of rank at most $r$ contains an all-0 or all-1
submatrix of area at least $mn/2^{\alpha(r)}$ (here, the _area_ of a matrix
refers to the product of the number of rows and columns). It was proved by
Nisan and Wigderson [13] that $c(f)=O((\log
r)^{2}+\sum_{i=0}^{\log_{2}r}\alpha(r/2^{i}))$.

The best known lower bound $c(r)=\tilde{\Omega}((\log r)^{2})$ is due to Göös,
Pitassi, and Watson [3]. On the other hand, the best known upper bound was due
to Lovett [12], who showed that $\alpha(r)\leq O(\sqrt{r}\cdot\log r)$, and
thus $c(r)=O(\sqrt{r}\cdot\log r)$. Here, we provide the modest improvement
that $c(r)=O(\sqrt{r})$ by establishing the following improvement on
$\alpha(r)$.

######  Theorem 1.1.

There exists a constant $c>0$ such that the following holds. Let
$M\in\\{0,1\\}^{m\times n}$ such that $\operatorname{rank}(M)\leq r$. Then
there exists $X\subset[m]$ and $Y\subset[n]$ such that

$|X|\geq\frac{m}{2^{c\sqrt{r}}}\mbox{ and }|Y|\geq\frac{n}{2^{c\sqrt{r}}},$

and $M[X\times Y]$ contains only 0 or only 1 entries.

One of the key components of the proof of Lovett [12] is a result of Linial,
Mendelson, Schechtman, and Shraibman [8] (see also [9]) on the _discrepancy_
of low-rank matrices. We examine an alternative definition of discrepancy,
inspired by the one commonly used in graph theory [1, 14]. This notion of
discrepancy provides stronger bounds for sparse matrices, which allows us to
remove the log-factor from Lovett’s bound. Given a binary matrix
$M\in\\{0,1\\}^{m\times n}$, let $|M|$ denote the number of 1 entries of $M$.
Writing $p=\frac{|M|}{mn}$, define the _discrepancy_ of $M$ as

$\operatorname{disc}(M)=\max_{X\subset[m],Y\subset[n]}\big{|}|M[X\times
Y]|-p|X|\cdot|Y|\big{|}.$

Viewing $M$ as the bipartite graph, whose vertex classes are the rows and
columns, and whose edges correspond to the 1 entries of $M$, this definition
of discrepancy is a bipartite analogue of the graph discrepancy introduced in
the 80’s by Erdős, Goldberg, Pach, and Spencer [1] and extensively studied
since then. Our main contribution is the following lower bound on the
discrepancy of low-rank matrices.

######  Theorem 1.2.

There exists a constant $c>0$ such that the following holds. Let
$M\in\\{0,1\\}^{m\times n}$ such that $\operatorname{rank}(M)\leq r$. If
$p=\frac{|M|}{mn}\leq 1/2$, then $\operatorname{disc}(M)\geq
cmn\cdot\min\left\\{p,\frac{p^{1/2}}{\sqrt{r}}\right\\}$.

Note that this bound is the best possible. In case $p\gg 1/r$, i.e. when the
minimum is attained by $p^{1/2}/\sqrt{r}$, consider a random $r\times r$
binary matrix $R$ in which each entry is $1$ with probability $p$. It is easy
to show that $\operatorname{disc}(R)=\Theta(\sqrt{p}r^{3/2})$ with high
probability. For every $m,n$, then one can construct the $m\times n$ matrix
$M$ by repeating each row of $R$ $m/r$-times, and repeating each column
$n/r$-times. Then $\operatorname{rank}(M)\leq r$ and
$\operatorname{disc}(M)=\frac{mn}{r^{2}}\cdot\operatorname{disc}(R)=\Theta(mnp^{1/2}/\sqrt{r})$.
On the other hand, if $p\ll 1/r$, then the minimum is attained by $p$.
Trivially, $|M|$ is always an upper bound on the discrepancy, so our theorem
is tight in this case as well.

In order to prove Theorem 1.2, we follow the recent framework of Räty,
Sudakov, and Tomon [14]. While their results are concerned about the usual
graph theoretic notion of discrepancy, working in a bipartite setting
simplifies many of their ideas. We consider a semidefinitie relaxation of the
discrepancy, and then apply spectral techniques to prove our lower bounds. We
remark that this approach is fairly different from (and perhaps more
elementary than) that of Linial, Mendelson, Schechtman, and Schraibman [9],
which relies on John’s ellipsoid theorem [5].

##  2 Discrepancy

In this section, we present basic notions and results about matrices and
discrepancy, and introduce a semidefinite relaxation. We omit floors and
ceilings whenever they are not crucial.

An $m\times n$ binary matrix $M$ naturally corresponds to the bipartite graph,
whose two vertex classes are identified with $[m]=\\{1,\dots,m\\}$ and $[n]$,
and there is an edge between $i$ and $j$ if $M_{i,j}=1$. While we will use the
language of binary matrices, this correspondence is good to keep in mind for
many arguments.

Given $M\in\\{0,1\\}^{m\times n}$, let $p(M)=\frac{|M|}{mn}$ be the _density_
of $M$, and let $d(M)=\frac{2|M|}{m+n}$ (which is the average degree of the
corresponding bipartite graph). In what follows, fix some
$M\in\\{0,1\\}^{m\times n}$ and write $p=p(M)$. Given $X\subset[m]$ and
$Y\subset[n]$, $M[X\times Y]$ is the submatrix of $M$ induced by the rows $X$
and columns $Y$. Define

$\operatorname{disc}_{M}(X,Y)=\operatorname{disc}(X,Y)=|M[X,Y]|-p|X||Y|.$

Furthermore, define the _positive discrepancy_ of $M$ as

$\operatorname{disc}^{+}(M)=\max_{X\subset[m],Y\subset[n]}\operatorname{disc}(X,Y)$

and the _negative discrepancy_ as

$\operatorname{disc}^{-}(M)=\max_{X\subset[m],Y\subset[n]}-\operatorname{disc}(X,Y).$

Then
$\operatorname{disc}(M)=\max\\{\operatorname{disc}^{+}(M),\operatorname{disc}^{-}(M)\\}$.
Note that as $\operatorname{disc}(\emptyset,\emptyset)=0$, we have
$\operatorname{disc}^{+}(G)\geq 0$ and $\operatorname{disc}^{-}(G)\geq 0$.
Also, if $X$ and $X^{\prime}$ are disjoint, then $\operatorname{disc}(X\cup
X^{\prime},Y)=\operatorname{disc}(X,Y)+\operatorname{disc}(X^{\prime},Y)$.
Next, we show that $\operatorname{disc}^{-}(M)$ and
$\operatorname{disc}^{+}(M)$ cannot be too different.

######  Claim 2.1.

$\operatorname{disc}^{+}(M)=\Theta(\operatorname{disc}^{-}(M))$.

###### Proof.

Let $X\subset[m]$ and $Y\subset[n]$ such that
$\operatorname{disc}^{+}(M)=\operatorname{disc}(X,Y).$ Let
$X^{\prime}=[m]\setminus X$ and $Y^{\prime}=[n]\setminus Y$. Then,
$0=\operatorname{disc}([m],[n])=\operatorname{disc}(X,Y)+\operatorname{disc}(X^{\prime},Y)+\operatorname{disc}(X,Y^{\prime})+\operatorname{disc}(X^{\prime},Y^{\prime})$.
Hence, at least one of the terms is at most $-\operatorname{disc}(X,Y)/3$,
which shows that $\operatorname{disc}^{-}(M)\geq\operatorname{disc}^{+}(M)/3$.
Similarly, $\operatorname{disc}^{+}(M)\geq\operatorname{disc}^{-}(M)/3$. ∎

Next, we show that it is enough to consider submatrices of $M$ that intersect
exactly half of the rows and half of the columns to get the discrepancy up to
a constant.

######  Claim 2.2.

There exists $X_{0}\subset[m]$ and $Y_{0}\subset[n]$ such that $|X_{0}|=m/2$
and $|Y_{0}|=n/2$, and

$\operatorname{disc}(X_{0},Y_{0})\leq-\operatorname{disc}^{-}(M)/12.$

###### Proof.

Let $X\subset[m]$ and $Y\subset[n]$ such that
$|\operatorname{disc}(X,Y)|=\operatorname{disc}(M)$. Assume that
$\operatorname{disc}(X,Y)\geq 0$, the other case can be handled similarly.

First, let us assume that $|X|\leq m/2$. Let $X_{1}$ be a random $m/2$ element
subset of $[m]$ containing $X$. Then

$\displaystyle\mathbb{E}(\operatorname{disc}(X_{1},Y))$
$\displaystyle=\operatorname{disc}(X,Y)+\mathbb{E}(\operatorname{disc}(X_{1}\setminus
X,Y))=\operatorname{disc}(X,Y)+\frac{m/2-|X|}{m-|X|}\operatorname{disc}(A\setminus
X,Y)$
$\displaystyle\geq\operatorname{disc}(X,Y)-\frac{m/2-|X|}{m-|X|}\operatorname{disc}(X,Y)\geq\frac{1}{2}\operatorname{disc}(X,Y).$

Here, the first inequality holds by the maximality of $(X,Y)$. Therefore,
there exists a choice for $X_{1}$ such that
$\operatorname{disc}(X_{1},Y)\geq\frac{1}{2}\operatorname{disc}(X,Y)$. In case
$|X|\geq m/2$, let $X_{1}$ be a random $m/2$ element subset of $X$. Similarly
as before, by considering the expectation, there exists a choice for $X_{1}$
such that
$\operatorname{disc}(X_{1},Y)\geq\frac{1}{2}\operatorname{disc}(X,Y)$.

Now let $X_{2}\subset[m]$ and $Y_{1}\subset[n]$ such that $|X_{2}|=m/2$, and
$|\operatorname{disc}(X_{2},Y_{1})|$ is maximal under these conditions. Assume
that $\operatorname{disc}(X_{2},Y_{1})\geq 0$, the other case can be handled
similarly. Then
$\operatorname{disc}(X_{2},Y_{1})\geq\operatorname{disc}(X_{2},Y)$. If
$|Y_{1}|\leq n/2$, let $Y_{2}$ be a random $n/2$ element subset of $[n]$ that
contains $Y_{1}$. As before, by considering the expectation of
$\operatorname{disc}(X_{2},Y_{2})$, we get that there is a choice for $Y_{2}$
such that
$\operatorname{disc}(X_{2},Y_{2})\geq\frac{1}{2}\operatorname{disc}(X_{2},Y_{1})\geq\frac{1}{4}\operatorname{disc}(X,Y)$.
If $|Y_{1}|\geq n/2$, then we can proceed as above by taking $Y_{2}$ to be a
random $m/2$ element subset of $Y_{1}$.

In conclusion, we get $X_{2}\subset[m]$ and $Y_{2}\subset[n]$ such that
$|X_{2}|=m/2$, $|Y_{2}|=n/2$, and
$|\operatorname{disc}(X_{2},Y_{2})|\geq\frac{1}{4}|\operatorname{disc}(X,Y)|$.
If $\operatorname{disc}(X_{2},Y_{2})\leq 0$, we are done, so we may assume
that $\operatorname{disc}(X_{2},Y_{2})>0$. Let $X_{3}=[m]\setminus X_{2}$, and
$Y_{3}=[n]\setminus Y_{2}$. Writing the equality
$0=\operatorname{disc}(X_{2},Y_{2})+\operatorname{disc}(X_{3},Y_{2})+\operatorname{disc}(X_{2},Y_{3})+\operatorname{disc}(X_{3},Y_{3})$,
at least one of the terms is at most
$-\operatorname{disc}(X_{2},Y_{2})/3\leq-|\operatorname{disc}(X,Y)|/12$.
Hence, we can choose $X_{0}$ and $Y_{0}$ such that $|X_{0}|=n/2$,
$|Y_{0}|=m/2$ and
$\operatorname{disc}(X_{0},Y_{0})\leq-\operatorname{disc}(M)/12\leq-\operatorname{disc}^{-}(M)/12$.
∎

Consider the following relaxation of the positive discrepancy. Given
$\mathbf{x}\in\mathbb{R}^{m}$ and $\mathbf{y}\in\mathbb{R}^{n}$, let

$\operatorname{disc}(\mathbf{x},\mathbf{y})=\mathbf{x}^{T}(M-pJ)\mathbf{y},$

where $J$ is the $m\times n$ sized all-1 matrix. Also, define

$\operatorname{disc}_{0}^{+}(M)=\max_{\mathbf{x}\in[-1,1]^{m},\mathbf{y}\in[-1,1]^{n}}\operatorname{disc}(\mathbf{x},\mathbf{y}).$

######  Claim 2.3.

$\operatorname{disc}^{+}(M)=\Theta(\operatorname{disc}_{0}^{+}(M)).$

###### Proof.

Let $\mathbf{x}\in[-1,1]^{m},\mathbf{y}\in[-1,1]^{n}$ be such that
$\operatorname{disc}_{0}^{+}(M)=\operatorname{disc}(\mathbf{x},\mathbf{y})$.
As $\operatorname{disc}(\mathbf{x},\mathbf{y})$ is a linear function in every
variable, we can find $\mathbf{x}\in\\{-1,1\\}^{m}$ and
$\mathbf{y}\in\\{-1,1\\}^{n}$ achieving the maximum. Let
$X_{1}=\\{i\in[m]:\mathbf{x}(i)=1\\}$, $X_{2}=[m]\setminus X_{1}$,
$Y_{1}=\\{i\in[n]:\mathbf{y}(i)=1\\}$ and $Y_{2}=[n]\setminus Y_{1}$. Then

$\operatorname{disc}(\mathbf{x},\mathbf{y})=\operatorname{disc}(X_{1},Y_{1})-\operatorname{disc}(X_{1},Y_{2})-\operatorname{disc}(X_{2},Y_{1})+\operatorname{disc}(X_{2},Y_{2}).$

Hence, at least one of the four terms has absolute value at least
$\frac{1}{4}\operatorname{disc}_{0}^{+}(M)$. But this shows that either
$\operatorname{disc}^{+}(M)$ or $\operatorname{disc}^{-}(M)$ is at least
$\frac{1}{4}\operatorname{disc}_{0}^{+}(M)$, and we are done by Claim 2.1. ∎

Let us further relax the definition of discrepancy by assigning vector values
to the rows and columns. That is, let

$\operatorname{pdisc}_{0}(M)=\max\sum_{i,j}M_{i,j}\langle
v_{i},w_{j}\rangle-p\sum_{i,j}\langle v_{i},w_{j}\rangle,$

where the maximum is taken over all
$v_{1},\dots,v_{m},w_{1},\dots,w_{n}\in\mathbb{R}^{m+n}$ such that
$||v_{i}||_{2},||w_{j}||_{2}\leq 1$. Clearly,
$\operatorname{pdisc}_{0}(M)\geq\operatorname{disc}^{+}(M)$. However, it
follows from Grothendieck’s inequality [4] that
$\operatorname{pdisc}_{0}(M)=O(\operatorname{disc}^{+}(M))$ also holds.

######  Lemma 2.4 (Grothendieck’s inequality).

There exists a universal constant $K>0$ such that the following holds. For a
matrix $Q\in\mathbb{R}^{n\times m}$ let
$\beta=\sup\left\\{\sum_{i,j}Q_{i,j}x_{i}y_{j}:x_{i},y_{j}\in[-1,1]\right\\}$
and $\beta^{*}=\sup\left\\{\sum_{i,j}Q_{i,j}\langle
v_{i},w_{j}\rangle:v_{i},w_{j},\in\mathbb{R}^{n+m},||v_{i}||_{2}\leq
1,||w_{i}||_{2}\leq 1\right\\}$. Then $\beta\leq\beta^{*}\leq K\beta$.

######  Claim 2.5.

$\operatorname{pdisc}_{0}(M)=\Theta(\operatorname{disc}^{+}(M))=\Theta(\operatorname{disc}^{-}(M))$

###### Proof.

Apply Groethendieck’s inequality with the matrix $Q\in\mathbb{R}^{m\times n}$
defined as $Q_{i,j}=M_{i,j}-p$. ∎

Let $N=m+n$. Let $A\in\mathbb{R}^{N\times N}$ be the adjacency matrix of the
bipartite graph corresponding to $M$, and let $L$ be the adjacency matrix of
the complete bipartite graph on the same vertex set. Formally, $A$ is the
symmetric matrix defined as $A_{i,j+m}=A_{j+m,i}=M_{i,j}$ for
$(i,j)\in[m]\times[n]$, and $A_{i,j}=0$ for $(i,j)\in[m]\times[m]$ and
$(i,j)\in[m+1,N]\times[m+1,N]$. Note that
$\operatorname{rank}(A)=2\operatorname{rank}(M)$. On the other hand, we can
write $L=\frac{1}{2}(e\cdot e^{T}-f\cdot f^{T})$, where $e$ is the all 1
vector, and $f\in\mathbb{R}^{N}$ is defined $f(i)=1$ if $i\in[m]$ and
$f(i)=-1$ if $i\in[m+1,N]$. Call $A$ the _symmetrization_ of $M$.

Given a matrix $X\in\mathbb{R}^{N\times N}$, define its discrepancy as

$\operatorname{disc}_{M}(X)=\operatorname{disc}(X)=\langle X,A\rangle-p\langle
X,L\rangle,$

where $\langle.,.\rangle$ denotes the Frobenius (i.e. entry-wise) inner
product. Furthermore, define

$\operatorname{pdisc}(M)=\max\left\\{\operatorname{disc}(X):X\mbox{ is
symmetric, positive semidefinite, }\forall i\in[N],X_{i,i}\leq 1\right\\}.$

Then $\operatorname{pdisc}(M)=2\operatorname{pdisc}_{0}(M)$, as every $X$
satisfying the conditions in the definition of $\operatorname{pdisc}(M)$
corresponds to the Gramm-matrix of vectors
$v_{1},\dots,v_{m},w_{1},\dots,w_{n}$ of length at most $1$. In conclusion, we
arrive to the following.

######  Claim 2.6.

$\operatorname{pdisc}(M)=\Theta(\operatorname{disc}^{+}(M))=\Theta(\operatorname{disc}^{-}(M)).$

Finally, let us argue that considering only those matrices with equal number
of rows and columns does not restrict us too much. Let $M^{\otimes}$ be the
$(mn)\times(mn)$ sized matrix we get by repeating each row of $M$ $n$-times,
and then repeating every column $m$ times. It is easy to see that

$\frac{1}{(mn)^{2}}\operatorname{disc}^{+}_{0}(M^{\otimes})=\frac{1}{mn}\operatorname{disc}^{+}_{0}(M),$

and $\operatorname{rank}(M)=\operatorname{rank}(M^{\otimes})$.

##  3 Low rank matrices

Now let us turn to bounding $\operatorname{pdisc}(M)$. Let $A$ be the
adjacency matrix of the bipartite graph corresponding to $M$, defined above.
Let $\lambda_{1}\geq\dots\geq\lambda_{N}$ be the eigenvalues of $A$, and let
$v_{1},\dots,v_{N}$ be a corresponding orthonormal set of eigenvectors. Note
that as $A$ is the adjacency matrix of a bipartite graph, we have
$\lambda_{i}=-\lambda_{N+1-i}$ for $i\in[N]$, and the vectors $v_{i}(j)$ and
$v_{N+1-i}(j)$ agree on $[m]$, and are opposite on $[m+1,N]$, see Chapter 11
of [10]. Formally, $v_{i}(j)=f(j)\cdot v_{N+1-i}(j)$ for every $i,j\in[N]$,
where $f(j)=1$ if $j\in[m]$ and $-1$ if $j\in[m+1,N]$. Note that
$\lambda_{1}\geq\dots\geq\lambda_{n}\geq
0\geq\lambda_{n+1}\geq\dots\geq\lambda_{N}$.

######  Lemma 3.1.

Let $a_{1},\dots,a_{N}\geq 0$ and let $X=\sum_{i=1}^{N}a_{i}\cdot v_{i}\cdot
v_{i}^{T}$. Then

$\operatorname{disc}(X)\geq\sum_{i=1}^{N}a_{i}\lambda_{i}-\frac{pN}{2}\cdot\max_{i}(a_{i}+a_{N+1-i}).$

###### Proof.

Write $A=\sum_{i=1}^{N}\lambda_{i}v_{i}\cdot v_{i}^{T}$, and observe the
identity $\langle v\cdot v^{T},w\cdot w^{T}\rangle=\langle v,w\rangle^{2}$.
Then

$\langle X,A\rangle=\sum_{i=1}^{N}\sum_{j=1}^{N}a_{i}\lambda_{j}\langle
v_{i},v_{j}\rangle^{2}=\sum_{i=1}^{N}\lambda_{i}a_{i}.$

One the other hand,

$\langle X,L\rangle=\frac{1}{2}\sum_{i=1}^{N}a_{i}(\langle
v_{i},e\rangle^{2}+\langle
v_{i},f\rangle^{2})=\frac{1}{2}\sum_{i=1}^{N}(a_{i}+a_{N+1-i})\langle
v_{i},e\rangle^{2}\leq\frac{N}{2}\max_{i}(a_{i}+a_{N+1-i}).$

Here, the second equality holds by noting that $\langle v_{i},e\rangle=\langle
v_{N+1-i},f\rangle$, and the last inequality as $\sum_{i=1}^{N}\langle
v_{i},e\rangle^{2}=||e||_{2}^{2}=N.$ ∎

In what follows, let us assume that $m=n$, so $N=2n$. Let
$d=d(M)=\frac{|M|}{n}=\frac{pN}{2}$, which is the average degree of the
corresponding bipartite graph $G$, and let $\Delta=\Delta(M)$ denote the
maximum number of 1 entries in a row or a column, that is, the maximum degree
of $G$.

######  Lemma 3.2.

$\operatorname{pdisc}(M)\geq\frac{1}{\Delta}\sum_{i=2}^{n}\lambda_{i}^{3}.$

###### Proof.

Let $X=\frac{1}{\Delta}\sum_{i=1}^{n}\lambda_{i}^{2}\cdot v_{i}\cdot
v_{i}^{T}$, which is positive semidefinite. We also claim that $X_{i,i}\leq 1$
for every $i\in[N]$. Indeed, note that
$\frac{1}{\Delta}A^{2}-X=\frac{1}{\Delta}\sum_{i=n+1}^{N}\lambda_{i}^{2}\cdot
v_{i}\cdot v_{i}^{T}$ is positive semidefinite. Therefore all diagonal entries
of $\frac{1}{\Delta}A^{2}-X$ are non-negative. Moreover, the $i$-th diagonal
entry of $A^{2}$ counts the number of 1 in the $i$-th row of $A$. Therefore
every diagonal entry of $\frac{1}{\Delta}A^{2}$ is at most $1$. Thus,
$X_{i,i}\leq 1$ for every $i\in[N]$, which implies that
$\operatorname{pdisc}(M)\geq\operatorname{disc}(X)$. Then, by Lemma 3.1 (with
$a_{i}=\lambda_{i}^{2}$ for $1\leq i\leq n$ and $a_{i}=0$ for $n+1\leq i\leq
N$), we have

$\operatorname{pdisc}(M)\geq\operatorname{disc}(X)\geq\frac{1}{\Delta}\left(\sum_{i=1}^{n}\lambda_{i}^{3}-\frac{pN}{2}\cdot\lambda_{1}^{2}\right)=\frac{1}{\Delta}\left(\sum_{i=1}^{n}\lambda_{i}^{3}-d\lambda_{1}^{2}\right)\geq\frac{1}{\Delta}\sum_{i=2}^{n}\lambda_{i}^{3}.$

In the last inequality, we used the well known result that the maximal
eigenvalue $\lambda_{1}$ of a graph $G$ is always at least its average degree
$d$, see e.g. [10]. ∎

Next, we show that if $M$ has low rank and not too large maximum degree, then
it has large discrepancy.

######  Lemma 3.3.

Let $\operatorname{rank}(M)\leq r$, $d\leq n/2$ and $\Delta\leq 1.1d$. Then
$\operatorname{pdisc}(M)\geq\frac{d^{1/2}n^{3/2}}{7\sqrt{r}}$.

###### Proof.

Note that $\sum_{i=1}^{2n}\lambda_{i}^{2}=\mbox{tr}(A)=2dn$, so
$\sum_{i=2}^{n}\lambda_{i}^{2}\geq dn-\Delta^{2}\geq dn-1.21d^{2}\geq dn/3$.
Here, we used that the maximum degree $\Delta$ of a graph is always at least
the largest eigenvalue $\lambda_{1}$, see e.g. [10]. Moreover, as
$\operatorname{rank}(M)\leq r$, we have $\lambda_{r+1}=\dots=\lambda_{n}=0$,
so $\sum_{i=2}^{r}\lambda_{i}^{2}\geq dn/3$. Applying the inequality between
the quadratic and cubic mean, we get

$\left(\frac{\sum_{i=2}^{r}\lambda_{i}^{3}}{r-1}\right)^{1/3}\geq\left(\frac{\sum_{i=2}^{r}\lambda_{i}^{2}}{r-1}\right)^{1/2}\geq\frac{d^{1/2}n^{1/2}}{\sqrt{3(r-1)}}.$

Hence, applying Lemma 3.2,
$\operatorname{pdisc}(M)\geq\frac{1}{1.1d}\sum_{i=2}^{r}\lambda_{i}^{3}\geq\frac{d^{1/2}n^{3/2}}{7\sqrt{r}}$.
∎

We improve the previous result by showing that the maximum degree condition
can be omitted.

######  Lemma 3.4.

Let $\operatorname{rank}(M)\leq r$, $d\leq n/2$. Then
$\operatorname{pdisc}(M)\geq
c\cdot\min\left\\{dn,\frac{d^{1/2}n^{3/2}}{\sqrt{r}}\right\\}$ for some $c>0$.

###### Proof.

Let $D=\min\left\\{dn,\frac{d^{1/2}n^{3/2}}{7\sqrt{r}}\right\\}$ and
$\delta=0.01$. For $i\in[n]$, write $\deg_{r}(i)$ for the number of 1 entries
in row $i$, and let $\deg_{c}(i)$ be the number of 1 entries in column $i$.
Let $U_{r}=\\{i\in[n]:\deg_{r}(i)\geq(1+\delta)d\\}$, and similarly, let
$U_{c}=\\{i\in[n]:\deg_{c}(i)\geq(1+\delta)d\\}$. Finally, write $t_{r}$ for
the number of 1 entries in $M[U_{r}\times[n]]$, and define $t_{c}$
analogously. We have

$\operatorname{disc}(U_{r},[n])=t_{r}-p|U_{r}|n=\sum_{i\in
U_{r}}(\deg_{r}(i)-d)\geq\sum_{i\in
U_{r}}\frac{\delta}{2}\cdot\deg_{r}(i)=\frac{\delta}{2}t_{r}.$

Similarly, $\operatorname{disc}([n],U_{c})\geq\frac{\delta}{2}t_{c}$.
Therefore, if $t_{r}+t_{c}\geq 0.01D$, then
$\operatorname{disc}^{+}(M)\geq\frac{\delta}{2}\cdot\max\\{t_{r},t_{c}\\}\geq\frac{1}{4}10^{-4}D$.
Hence, as $\operatorname{pdisc}(M)=\Theta(\operatorname{disc}^{+}(M))$, we may
assume that $t_{r}+t_{c}<0.01D$.

Change all 1 entries of $M$ that are contained in either $U_{r}\times[n]$ or
$[n]\times U_{c}$ into a 0 entry, and let $M^{\prime}$ be the resulting
matrix. Then

$|M^{\prime}|\geq|M|-t_{r}-t_{c}\geq|M|-0.01D\geq|M|-0.01dn,$

so $d^{\prime}=d(M^{\prime})\geq 0.99d$ and $d^{\prime}\leq d$. Furthermore,
$\Delta(M^{\prime})\leq(1+\delta)d\leq 1.1d^{\prime}$, and
$\operatorname{rank}(M^{\prime})\leq r$. Therefore, we can apply Lemma 3.3 to
conclude that $\operatorname{pdisc}(M^{\prime})\geq\frac{d^{\prime
1/2}n^{3/2}}{7\sqrt{r}}\geq D/2$. Let $X$ be matrix such that
$\operatorname{disc}_{M^{\prime}}(X)=\operatorname{pdisc}(M^{\prime})$, and
let $A^{\prime}$ be the symmetrization of $M^{\prime}$. Then

$\displaystyle|\operatorname{disc}_{M}(X)-\operatorname{disc}_{M^{\prime}}(X)|$
$\displaystyle=|\langle X,A-A^{\prime}\rangle-\frac{d-d^{\prime}}{n}\langle
X,L\rangle|$ $\displaystyle\leq\langle
J,A-A^{\prime}\rangle+2(d-d^{\prime})n\leq 4(|M|-|M^{\prime}|)\leq 0.04D.$

Hence,
$\operatorname{disc}_{M}(X)\geq\operatorname{disc}_{M^{\prime}}(X)-0.04D\geq
D/2-0.04D\geq D/4$, finishing the proof. ∎

Note that the previous lemma immediately implies Theorem 1.2 by considering
$M^{\otimes}$ instead of $M$, and using that
$\operatorname{pdisc}(M)=\Theta(\operatorname{disc}(M))$. It will be
convenient to rewrite Lemma 3.4 into a slightly more convenient form.

######  Corollary 3.5.

There exists $c>0$ such that the following holds. Assume that
$\operatorname{rank}(M)\leq r$ and $1/(8r)\leq p\leq 1/2$. Then there exists
$X,Y\subset[n]$ such that $|X|=|Y|=n/2$, and

$p(M[X\times Y])\leq p-c\frac{p^{1/2}}{\sqrt{r}}.$

###### Proof.

Note that $p=d/n$, so for $1/(8r)\leq p$, we have
$\min\\{dn,\frac{d^{1/2}n^{3/2}}{\sqrt{r}}\\}=\Omega(p^{1/2}n^{2}/\sqrt{r})$.
Therefore,
$\operatorname{disc}^{-}(M)=\Omega(\operatorname{pdisc}(M))=\Omega(p^{1/2}n^{2}/\sqrt{r})$
by Lemma 3.4. By Claim 2.2, there exist $X\subset[n]$ and $Y\subset[n]$ such
that $|X|=|Y|=n/2$ and
$-\operatorname{disc}(X,Y)=\Omega(\operatorname{disc}^{-}(M))$. But then if
$M^{\prime}=M[X\times Y]$, we have

$|M^{\prime}|=p|X||Y|+\operatorname{disc}(X,Y)\leq
p|X||Y|-\Omega(p^{1/2}n^{2}/\sqrt{r}),$

which implies

$p(M^{\prime})=\frac{|M^{\prime}|}{|X||Y|}\leq
p-\Omega\left(\frac{p^{1/2}}{\sqrt{r}}\right).$

∎

##  4 Log-rank conjecture

Let $z(M)$ denote the largest $z$ for which there exist $X\subset[n]$,
$Y\subset[n]$, such that $|X|=|Y|=z$ and $M[X\times Y]$ has only 0 entries.
Based on a lemma of Gavinsky and Lovett [2], we show that sufficiently sparse
low-rank matrices $M$ satisfy $z(M)=\Omega(n)$.

######  Lemma 4.1.

Let $M\in\\{0,1\\}^{n\times n}$ such that $\operatorname{rank}(M)\leq r$ and
$p(M)\leq\frac{1}{8r}$. Then $z(M)\geq n/4$.

###### Proof.

Let $X,Y\subset[n]$ be the set of rows and columns containing more than
$n/(4r)$ entries equal to 1, respectively. As $p(M)\leq\frac{1}{8r}$, we have
$|X|,|Y|\leq n/2$. Let $X^{\prime}=[n]\setminus X$ and
$Y^{\prime}=[n]\setminus Y$. Assume that $z(M[X^{\prime}\times
Y^{\prime}])<n/4$, then we find a permutation matrix of size $r+1$ in
$M[X^{\prime}\times Y^{\prime}]$, contradicting that
$\operatorname{rank}(M)\leq r$.

For $k=1,\dots,r+1$, we find a $k\times k$ sized permutation matrix greedily.
Suppose that we have already found $A\subset X^{\prime}$ and $B\subset
Y^{\prime}$ such that $|A|=|B|=k$, and $M[A\times B]$ is a permutation matrix.
Let $X_{0}\subset X^{\prime}$ be the set of all rows that intersect a column
in $B$ in a 1 entry, and let $Y_{0}\subset Y^{\prime}$ be the set of all
columns that intersect a row in $A$ in a 1 entry. Then $|X_{0}|,|Y_{0}|\leq
kn/(4r)\leq n/4$, so $|X^{\prime}\setminus X_{0}|\geq n/4$ and
$|Y^{\prime}\setminus Y_{0}|\geq n/4$. Since we assumed that
$z(M[X^{\prime}\times Y^{\prime}])<n/4$, there exists $i\in
X^{\prime}\setminus X_{0}$ and $j\in Y^{\prime}\setminus Y_{0}$ such that
$M_{i,j}=1$, which means that $(A\cup\\{i\\})\times(B\cup\\{j\\})$ induces a
$(k+1)\times(k+1)$ sized permutation matrix in $M$. ∎

######  Theorem 4.2.

There exists a constant $c>0$ such that the following holds. Let
$M\in\\{0,1\\}^{n\times n}$ such that $\operatorname{rank}(M)\leq r$ and
$p(M)\leq 1/2$. Then $z(M)\geq n/2^{c\sqrt{r}}$.

###### Proof.

We will proceed by a density decrement argument. Let $n_{0}=n$ and $M_{0}=M$,
and let $c_{0}$ be the constant given by Corollary 3.5. We define a sequence
$M_{0},M_{1},\dots$ of submatrices of $M$ with decreasing density. If
$M_{i}\in\\{0,1\\}^{n_{i}\times n_{i}}$ is already defined for $i\geq 0$ with
$p_{i}:=p(M_{i})\leq p(M)$, and we have $p_{i}\geq 1/(8r)$, we define
$M_{i+1}$ by setting $M_{i+1}=M_{i}[X\cup Y]$, where $X,Y\subset[n_{i}]$,
$|X|=|Y|=n_{i}/2$, and

$p(M_{i}[X\times Y])\leq p_{i}-c_{0}\frac{p_{i}^{1/2}}{\sqrt{r}}.$

Such $X$ and $Y$ exist by Corollary 3.5. On the other hand, if
$p(M_{i})<1/(8r)$, then set $I=i$ and stop.

Observe that $n_{i}=n/2^{i}$ and $p_{i+1}\leq
p_{i}-c_{0}\frac{p_{i}^{1/2}}{\sqrt{r}}$ for $i=0,\dots,I-1$. This implies
that if $C$ is a sufficiently large constant, then for every
$x\in[1/(8r),1/2]$, the number of indices $i$ such that $x\leq p_{i}\leq 2x$
is at most $C\sqrt{r}x^{1/2}$. Applying this for every $x=2^{-k}$, we get that

$I\leq C\sqrt{r}\sum_{i=1}^{\infty}2^{-i/2}\leq 10C\sqrt{r}.$

As $p(M_{I})\leq 1/(8r)$ and $n_{I}\geq n/2^{I}\geq n/2^{10C\sqrt{r}}$, we
have $z(M_{I})\geq n_{I}/2=n/2^{c\sqrt{r}}$ by Lemma 4.1 with a suitable
constant $c>0$, finishing the proof. ∎

From this, the proof of Theorem 1.1 is immediate.

###### Proof of Theorem 1.1.

Let $M\in\\{0,1\\}^{m\times n}$ such that $\operatorname{rank}(M)\leq r$.
Without loss of generality, we may assume that $p(M)\leq 1/2$, otherwise we
can consider $J-M$. Note that $\operatorname{rank}(J-M)\leq r+1$. Writing
$n^{\prime}=mn$, $M^{\otimes}\in\\{0,1\\}^{n^{\prime}\times n^{\prime}}$
satisfies that $p(M^{\otimes})\leq 1/2$. Hence, by the previous theorem,
$z(M^{\otimes})\geq n^{\prime}/2^{c\sqrt{r}}$ for some constant $c>0$. Note
that each row of $M$ is repeated $n$ times in $M^{\otimes}$, and each column
is at most $m$ times, so we can find $X\subset[m]$ and $Y\subset[n]$ such that
$|X|\geq z(M^{\otimes})/n$, $|Y|\geq z(M^{\otimes})/m$, and $M[X\times Y]$ has
only 0 entries. But then $|X|\geq m/2^{c\sqrt{r}}$ and $|Y|\geq
n/2^{c\sqrt{r}}$, finishing the proof. ∎

## References

  * [1] P. Erdős, M. Goldberg, J. Pach, and J. Spencer. _Cutting a graph into two dissimilar halves._ J. Graph Theory 12 (1988): 121–131.
  * [2] D. Gavinsky and S. Lovett. _En route to the log-rank conjecture: New reductions and equivalent formulations._ Electronic Colloquium on Computational Complexity (ECCC’13) (2013) 20, 80.
  * [3] M. Göös, T. Pitassi, T. Watson. _Deterministic Communication vs. Partition Number._ SIAM Journal on Computing, 47 (6) (2018): 2435–2450.
  * [4] A. Grothendieck. _Résumé de la théorie métrique des produits tensoriels topologiques._ Bol. Soc. Mat. Sao Paulo, 8 (1953): 1–79.
  * [5] F. John. _Extremum problems with inequalities as subsidiary conditions._ In Studies and Essays Presented to R. Courant on his 60th Birthday, January 8, 1948, pages 187–204. Interscience Publishers, Inc., New York, N. Y., 1948.
  * [6] E. Kushilevitz and N. Nisan. Communication complexity. Cambridge University Press, New York, NY (1997).
  * [7] T. Lee and A. Shraibman. _Around the log-rank conjecture._ Israel Journal of Mathematics 256 (2023): 441–477.
  * [8] N. Linial, S. Mendelson, G. Schechtman, and A. Shraibman. _Complexity measures of sign matrices._ Combinatorica 27 (4) (2007): 439–463.
  * [9] N. Linial and A. Shraibman. _Learning complexity vs. communication complexity._ Combinatorics, Probability & Computing 18, 1–2 (2009): 227–245.
  * [10] L. Lovász. Combinatorial Problems and Exercises, 2nd Ed., AMS Chelsea Publishing, 2007.
  * [11] L. Lovász and M. Saks. _Lattices, Möbius functions and communication complexity._ Annual Symposium on Foundations of Computer Science (1988): 81–90.
  * [12] S. Lovett. _Communication is Bounded by Root of Rank._ Journal of the ACM, 63 (1) (2016): 1:1–1:9.
  * [13] N. Nisan and A. Wigderson. _On rank vs. communication complexity._ Proceedings of the 35rd Annual Symposium on Foundations of Computer Science (1994): 831–836.
  * [14] E. Räty, B. Sudakov, and I. Tomon. _Positive discrepancy, MaxCut, and eigenvalues of graphs._ preprint, arxiv:2311.02070 (2023)