††thanks: These two authors contribute equally††thanks: These two authors
contribute equally

# Variational principle for optimal quantum controls in quantum metrology

Jing Yang<EMAIL_ADDRESS>Department of Physics and Materials Science,
University of Luxembourg, L-1511 Luxembourg, Luxembourg  Shengshi Pang
<EMAIL_ADDRESS>School of Physics, Sun Yat-Sen University,
Guangzhou, Guangdong 510275, China  Zekai Chen Department of Physics and
Astronomy, University of Rochester, Rochester, New York 14627, USA  Andrew N.
Jordan  Institute for Quantum Studies, Chapman University, 1 University Drive,
Orange, CA 92866, USA  Department of Physics and Astronomy, University of
Rochester, Rochester, New York 14627, USA  Adolfo del Campo
<EMAIL_ADDRESS>Department of Physics and Materials Science, University
of Luxembourg, L-1511 Luxembourg, Luxembourg  Donostia International Physics
Center, E-20018 San Sebastián, Spain

###### Abstract

We develop a variational principle to determine the quantum controls and
initial state which optimizes the quantum Fisher information, the quantity
characterizing the precision in quantum metrology. When the set of available
controls is limited, the exact optimal initial state and the optimal controls
are in general dependent on the probe time, a feature missing in the
unrestricted case. Yet, for time-independent Hamiltonians with restricted
controls, the problem can be approximately reduced to the unconstrained case
via the Floquet engineering. In particular, we find for magnetometry with a
time-independent spin chain containing three-body interactions, even when the
controls are restricted to one and two-body interaction, that the Heisenberg
scaling can still be approximately achieved. Our results open the door to
investigate quantum metrology under a limited set of available controls, of
relevance to many-body quantum metrology in realistic scenarios.

Introduction. At the heart of quantum technology, quantum metrology aims at
improving the precision of parameter estimation Helstrom (1976); Holevo
(2011); Braunstein and Caves (1994); Paris (2009); Giovannetti  _et al._
(2011). Recent theoretical advances in quantum metrology Yuan and Fung (2015);
Yuan (2016); Pang and Jordan (2017); Yang  _et al._ (2017); Liu and Yuan
(2017); Yang  _et al._ (2021a) have led to the introduction of optimal control
protocols that maximizes the quantum Fisher information (QFI), the key
quantity in characterizing the precision in quantum metrology. The use of
control protocols for enhanced parameter estimation has been extensively
explored in many practical physical setups Pang and Jordan (2017); Yang  _et
al._ (2017, 2021a); Gefen  _et al._ (2019); Naghiloo  _et al._ (2017).

The identification of the optimal control protocol for parameter estimation
Pang and Jordan (2017); Yang  _et al._ (2017) can be done exploiting the
notion of adiabatic continuation of quantum states. It resembles the
engineering of shortcuts to adiabaticity (STA) by counterdiabatic driving
Demirplak and Rice (2003, 2005, 2008); Berry (2009); Guéry-Odelin  _et al._
(2019), whereby the system Hamiltonian is supplemented by some controls that
enforce adiabatic continuation along a prescribed trajectory. In the context
of STA, the auxiliary controls enforce parallel transport along the
eigenstates of the uncontrolled system Hamiltonian. By contrast, optimal
controls for quantum metrology enforce adiabatic continuation of the
eigenstates of an operator defined by the parametric derivative of the
estimation Hamiltonian Cabedo-Olaya  _et al._ (2020). This is intuitive as
such operator guarantees maximal distinguishability of states with a slight
change of the estimated parameter. The connection of optimal control for
quantum metrology and STA is not only fruitful in providing a geometric
justification of the required controls. As we show in this work, it suggests
solutions to common challenges.

In STA, the exact identification of the counterdiabatic control protocol is
not feasible in many-body systems in which spectral properties of the system
Hamiltonian cannot be found, e.g., in problems such as quantum optimization.
In addition, the set of available controls in a given experimental setup is
often restricted. This generally precludes the implementation of the exact STA
protocol. This is particularly important in the many-particle systems when
exact counterdiabatic controls involve non-local multiple-body interactions
beyond one-body and pair-wise potentials which are hard to realize in practice
del Campo  _et al._ (2012). In addition, the use of variational methods
provides then an alternative, by designing optimal protocols that are
realizable with a restricted set of controls Carlini  _et al._ (2007, 2006);
Rezakhani  _et al._ (2009); Opatrný and Mølmer (2014); Saberi  _et al._
(2014); Sels and Polkovnikov (2017); Claeys  _et al._ (2019); Chandarana  _et
al._ (2021); Wang  _et al._ (2015). Likewise, the exact construction of
optimal protocols in many-body quantum metrology Yang  _et al._ (2021a) (i)
requires access to the spectrum of the parametric derivative of the many-body
Hamiltonian which is hardly accessible in general, (ii) may require optimal
controls that are nonlocal, and hard to implement in the laboratory. Even in
single-qubit metrology using a NV center, superconducting qubit, or quantum
dot, certain control operations may be hard to implement. It is thus required
to develop a new formalism of optimal control for quantum metrology beyond the
state of the art Pang and Jordan (2017), by taking into account the ubiquitous
limitations on the controls and circumventing the requirement to access the
spectral properties of the Hamiltonian derivative. This is the problem solved
in this Letter by means of a novel variational approach that relies on a
metrological action. We note that our work is unrelated to other variational
approaches introduced in quantum metrology that do not involve quantum control
Kaubruegger  _et al._ (2021); Marciniak  _et al._ (2021).

Specifically, we assume that the space of control Hamiltonians is spanned by a
set of local basis operators and introduce a metrological action, which
includes contributions of the quantum Fisher information and the Schrödinger
equation, which is considered as a constraint to the optimization problem. We
derive the optimal control conditions and show that it reduces to the
unrestricted protocol in Ref. Pang and Jordan (2017) if the basis operators
generate the full space of the Hermitian operators. Furthermore, we show that
when the control Hamiltonian is restricted, the exact optimal solution for the
initial state and the control Hamiltonian depends on the probe time. However,
we find that for multiplicative time-independent Hamiltonians, even with
limited controls, one may identify approximate optimal controls, which reduce
to the unrestricted protocol and therefore become independent of the probe
time, using the high frequency expansion in Floquet engineering Rahav  _et
al._ (2003); Goldman and Dalibard (2014); Bukov  _et al._ (2015); Eckardt and
Anisimovas (2015); Goldman and Dalibard (2015); Eckardt (2017); Chen  _et al._
(2020). In particular, we show in an example in magnetometry that one can
avoid local three-body interaction terms making use only of local two-body
control Hamiltonians.

Optimal initial states and controls through variational principle. Consider
the quantum estimation of a parameter $\lambda$ in a general Hamiltonian
$H_{\lambda}(t)$. In this setting, one may find that the QFI may decrease as
the probe time increases Pang and Brun (2014). This motivates the introduction
of quantum controls $H_{\text{c}}(t)$ to enhance the QFI Yuan and Fung (2015);
Pang and Jordan (2017). When the quantum controls are introduced, the unitary
evolution $U(t)$ is generated by
$H_{\text{tot}}(t)=H_{\lambda}(t)+H_{\text{c}}(t)$. Let the initial time and
fixed final probe time be $0$ and $t_{f}$, respectively. The generator for
parameter estimation is given by
$G_{t_{f}}[U]=\int_{0}^{t_{f}}U^{\dagger}(\tau)\partial_{\lambda}H_{\lambda}(\tau)U(\tau)d\tau$
Pang and Jordan (2017) and the quantum Fisher information $I$ is given by the
variance of the generator Paris (2009), i.e., $I[U]=\text{Var}(G_{t_{f}}[U])$.

The optimization of the QFI $I[U]$ over the initial state yields the optimal
initial state $(\ket{\varphi_{+}(t_{f})}+\ket{\varphi_{-}(t_{f})})/\sqrt{2}$,
where $\ket{\varphi_{\pm}(t_{f})}$ are referred to as the maximum and minimum
eigenstates, as they are associated with the maximum and minimum eigenvalues
$\mu_{\pm}(t_{f})$ of $G_{t_{f}}[U]$ Yang  _et al._ (2021b). The maximum value
of the QFI over the initial states $\ket{\psi_{0}}$ is
$\max_{\ket{\psi_{0}}}I=\|G_{t_{f}}[U]\|^{2}/4$, where the norm of an operator
is defined by the difference between its maximum and minimum eigenvalues
Braunstein  _et al._ (1996); Boixo  _et al._ (2007); Giovannetti  _et al._
(2011). Our next goal is to maximize the quantum Fisher information $I[U]$
over all possible unitary dynamics under the condition that $U$ and
$H_{\text{c}}$ satisfy the Schrödinger equation
$\text{i}\partial_{t}U(t)=H_{\text{tot}}(t)U(t)$ as a constraint. We further
require the control Hamiltonian $H_{\text{c}}$ to be spanned by a limited set
of available linearly-independent terms is
$\\{\mathcal{X}_{i}\\}_{i=1}^{d_{\text{c}}}$. We denote
$\mathcal{V}_{\text{c}}\equiv\text{span}\\{\mathcal{X}_{i}\\}_{i=1}^{d_{\text{c}}}$
and expand $H_{\text{c}}\in\mathcal{V}_{\text{c}}$ as
$H_{\text{c}}(\tau)=\sum_{i=1}^{d_{\text{c}}}c_{i}(\tau)\mathcal{X}_{i}$. This
procedure amounts to removing the constraint on $H_{\text{c}}$, reducing the
optimization over $H_{\text{c}}$ to the optimization over
$c_{i}(\tau)\leavevmode\nobreak\ $ 111This is the approach pursued in e.g.,
the variational approach in shortcut-to-adiabaticity Sels and Polkovnikov
(2017); Claeys  _et al._ (2019). Another way of handling the constraints on
the control Hamiltonian is by introducing the following constraints
$f_{j}(H_{\text{c}}(\tau))=\text{Tr}\left\\{H_{\text{c}}(\tau)\mathcal{X}_{j}\right\\}=0,\,j=d_{\text{c}}+1,\,\cdots
N$ to disallow the terms $\mathcal{X}_{j}$, where $j=1,\,2\,\cdots
d_{\text{c}}$. This way of introducing the constraint is the one used in e.g.
quantum brachistochrone equation Carlini  _et al._ (2007, 2006). However, in
many-body quantum metrology, the number of disallowed nonlocal operators is
much more than the allowed local operators. Therefore the second approach may
introduce an intractable number of constraints and we shall pursue the first
approach of expanding $H_{\text{c}}$ in terms of basis operators in the main
text..

We shall denote $I_{0}=\max_{U,\,H_{\text{c}},\,\ket{\psi_{0}}}I[U]$. In
principle, $I_{0}$ can be computed through the variational principle by
constructing an appropriate action. We note that $I[U]$ is quartic in $U$
since it is quadratic in $G_{t_{f}}[U]$, which is itself quadratic in $U$.
This makes the variational calculus of $I[U]$ with respect to $U$ tedious. To
facilitate the calculation, we observe that

$\displaystyle I_{0}=$
$\displaystyle\max_{U,\,H_{\text{c}}}\max_{\ket{\varphi_{a,b}}}\left(\braket{\varphi_{a}\big{missing}}{G_{t_{f}}[U]\big{missing}}{\varphi_{a}}-\braket{\varphi_{b}\big{missing}}{G_{t_{f}}[U]\big{missing}}{\varphi_{b}}\right)^{2}$
(1)

under the constraint of Schrödinger equation and the condition that
$\ket{\varphi_{a,\,b}}$ are normalized. The introduction of two more
optimization variables $\ket{\varphi_{a,\,b}}$ allows us to remove the square
in Eq. (1) and transform the original optimization problem to the following
equivalent one,
$\max_{\ket{\varphi_{a,b}},\,U,\,H_{\text{c}}}S_{\text{I}}[\Delta\rho,\,U]$,
under aforementioned constraints, with the "information action" being defined
as

$\displaystyle S_{\text{I}}[\Delta\rho,\,U]$
$\displaystyle\equiv\braket{\varphi_{a}\big{missing}}{G_{t_{f}}[U]\big{missing}}{\varphi_{a}}-\braket{\varphi_{b}\big{missing}}{G_{t_{f}}[U]\big{missing}}{\varphi_{b}}$
$\displaystyle=\int_{0}^{t_{f}}\text{Tr}\left\\{\Delta\rho
U^{\dagger}(\tau)\partial_{\lambda}H_{\lambda}(\tau)U(\tau)\right\\}d\tau,$
(2)

and
$\Delta\rho\equiv\ket{\varphi_{a}}\bra{\varphi_{a}}-\ket{\varphi_{b}}\bra{\varphi_{b}}$.
The introduction of two additional auxiliary variables $\ket{\varphi_{a,\,b}}$
effectively renders $S_{\text{I}}[U]$ quadratic in $U$, unlike $I[U]$,
facilitating the variational calculus with respect to $U$. Upon introducing
the Lagrangian multipliers $\mu_{a,\,b}$ and $\Lambda(\tau)$, we obtain the
following "metrological action"

$\displaystyle
S_{\text{M}}(\ket{\varphi_{a}},\,\ket{\varphi_{b}},\,U,\,H_{\text{c}})$
$\displaystyle\equiv
S_{\text{I}}[\Delta\rho,\,U]+S_{\text{S}}[U,\,H_{\text{c}}]$
$\displaystyle-\mu_{a}[\braket{\varphi_{a}\big{missing}}{\varphi_{a}}-1]+\mu_{b}[\braket{\varphi_{b}\big{missing}}{\varphi_{b}}-1]$
(3)

where the "Schrödinger action" is defined as

$S_{\text{S}}[U,\,H_{\text{c}}]\equiv\int_{0}^{t_{f}}\text{Tr}\left\\{\Lambda(\tau)[\text{i}\dot{U}(\tau)U^{\dagger}(\tau)-H_{\lambda}(\tau)-H_{\text{c}}(\tau)]\right\\}d\tau.$
(4)

We emphasize that $\ket{\varphi_{a,\,b}}$, $U$ and $H_{\text{c}}$ are
independent variables. The optimization over $\ket{\varphi_{a,\,b}}$ can be
easily implemented by differentiation with respect to them, which yields

$G_{t_{f}}[U]\ket{\varphi_{\alpha}}=\mu_{\alpha}\ket{\varphi_{\alpha}},\,\alpha=a,\,b.$
(5)

As shown in Sec. I in Yang  _et al._ (2021b)), in order for $I_{0}$ to take
the global maximum values over $\ket{\varphi_{a,\,b}}$, $\mu_{a,\,b}$ and
$\ket{\varphi_{a,\,b}}$ must be the maximum and minimum eigenvalues and
eigenvectors of $G_{t_{f}}[U]$, respectively.

Variation with respect to $H_{\text{c}}$ and $U$ gives the trace condition
$\text{Tr}\left\\{\Lambda(\tau)\mathcal{X}_{i}\right\\}=0$ for
$i\in[1,\,d_{\text{c}}]$ and the differential equation
$\dot{\Lambda}(\tau)-\text{i}[\Lambda(\tau),\,H_{\text{tot}}(\tau)]+\text{i}[\Delta\rho(\tau),\,\partial_{\lambda}H_{\lambda}(\tau)]=0$,
with the final condition $\Lambda(t_{f})=0$, where $\Delta\rho(\tau)\equiv
U(\tau)\Delta\rho U^{\dagger}(\tau)$. One can solve for $\Lambda(\tau)$ the
differential equation and substitute the result into the trace condition to
find

$\text{Tr}\left\\{\mathcal{X}_{i}\partial_{\lambda}[\Delta\rho(\tau)]\right\\}=0,\,i\in[1,\,d_{\text{c}}].$
(6)

Eq. (5) and Eq. (6) are our central results. The form of these equations in
the parameter-independent rotating frame can be found in Yang  _et al._
(2021b). They give the optimal initial state and optimal dynamics that
maximize the QFI when the quantum controls are restricted to the subspace
spanned by $\\{\mathcal{X}_{i}\\}$. In what follows, we discuss their
implications and applications.

The unrestricted control and the general feature of exact restricted controls.
As a first application of our results, let us assume there is no restriction
on the control Hamiltonians, that is, $\\{\mathcal{X}_{j}\\}$ spans the whole
space of traceless Hermitian operators. We shall see how the Pang-Jordan
protocol in Ref. Pang and Jordan (2017) is reproduced. In this case, Eq. (6)
is equivalent to $\partial_{\lambda}[\Delta\rho(\tau)]=0$. Taking the time
derivative on both sides yields

$\partial_{\tau}\partial_{\lambda}[\Delta\rho(\tau)]=-\text{i}[\partial_{\lambda}H_{\lambda}(\tau),\,\Delta\rho(\tau)]=0,\,\forall\tau\in[0,\,t_{f}],$
(7)

where we use the Liouville equation for $\Delta\rho(\tau)$. Since
$\ket{\varphi_{a,\,b}(\tau)}\equiv U(\tau)\ket{\varphi_{a,\,b}}$ are
associated with the non-degenerate eigenvalues $\pm 1$ of $\Delta\rho(\tau)$,
we conclude that $\ket{\varphi_{a,\,b}(\tau)}$ must also be eigenvectors of
$\partial_{\lambda}H_{\lambda}(\tau)$ at all times. That is,
$\partial_{\lambda}H_{\lambda}(\tau)\ket{\varphi_{a,\,b}(\tau)}=\nu_{a,\,b}(\tau)\ket{\varphi_{a,\,b}(\tau)},\,\forall\tau\in[0,\,t_{f}],$
where $\nu_{a,\,b}(\tau)$ is the eigenvalue of
$\partial_{\lambda}H_{\lambda}(\tau)$. It is then straightforward to check
that $\ket{\varphi_{a,\,b}}$ are eigenvectors of $G_{t_{f}}[U]$ with
eigenvalue $\int_{0}^{t_{f}}\nu_{a,\,b}(\tau)d\tau$. Thus any unitary dynamics
that preserves the adiabatic evolution of any pair of eigenstates of
$\partial_{\lambda}H_{\lambda}(\tau)$ satisfies Eq. (7) and is an extremal
solution satisfying $\delta S_{\text{M}}=0$. To further maximize the QFI among
the manifold of extremal solutions, one needs to further optimize the
difference between $\int_{0}^{t_{f}}\nu_{a}(\tau)$ and
$\int_{0}^{t_{f}}\nu_{b}(\tau)$. This requires $\nu_{a}(\tau)$ and
$\nu_{b}(\tau)$ to be the maximal and minimum eigenvalues of
$\partial_{\lambda}H_{\lambda}(\tau)$ at all times. When the unitary dynamics
preserves the adiabatic evolution of all eigenstates of
$\partial_{\lambda}H_{\lambda}(\tau)$, i.e.,
$U(\tau)=\sum_{\alpha}\ket{\varphi_{\alpha}(\tau)}\bra{\varphi_{\alpha}(0)}$,
where $\ket{\varphi_{\alpha}}$ denotes the eigenvectors of
$\partial_{\lambda}H_{\lambda}(\tau)$, one recovers the Pang-Jordan control
Hamiltonian Pang and Jordan (2017)
$H_{\text{c}}(\tau)=\text{i}\sum_{\alpha}\ket{\dot{\varphi}_{\alpha}(\tau)}\bra{\varphi_{\alpha}(\tau)}-H_{\lambda}(\tau).$
Note that if there is any level crossing in $\nu_{a}(\tau)$ at some instant
time in $[0,\,t_{f}]$, the way of labeling the eigenstates and eigenvalues is
not unique. Choosing $\nu_{+}(\tau)$ and $\nu_{-}(s)$ always as the maximum
and minimum eigenvalues of $\partial_{\lambda}H_{\lambda}(\tau)$ at all times,
the resulting QFI is the greatest among all the different ways of labeling the
eigenstates. With this labeling, the first-order time derivative of
$\ket{\varphi_{a}(\tau)}$ is discontinuous, which results in a $\delta$-pulse
in the control Hamiltonian $H_{\text{c}}$. This provides an alternative
understanding of the $\sigma_{x}$-like pulses in Ref. Yang  _et al._ (2017);
Pang and Jordan (2017); Naghiloo  _et al._ (2017).

In the general case, $\\{\mathcal{X}_{i}\\}_{i=1}^{d_{\text{c}}}$ do not span
the whole space of Hermitian operators, and the generic optimal solutions
$U(\tau)$ and $\ket{\varphi_{a,\,b}}$ implicitly depend on the probe time
$t_{f}$. This is due to the dependence of the generator $G_{t_{f}}[U]$ on
$t_{f}$, that may make the eigenvectors $\ket{\varphi_{a,\,b}}$ determining
the optimal initial state also dependent on $t_{f}$. The dependence on $t_{f}$
for the optimal unitary $U(\tau)$ can then be seen from Eq. (6). With these
observations, we take the derivative with respective to $t_{f}$ to obtain the
following differential-integral equation

$\partial_{t_{f}}G_{t_{f}}\ket{\varphi_{\alpha,\,t_{f}}}+G_{t_{f}}\ket{\partial_{t_{f}}\varphi_{\alpha,\,t_{f}}}=\partial_{t_{f}}\mu_{\alpha,\,t_{f}}\ket{\varphi_{\alpha,\,t_{f}}}+\mu_{\alpha,\,t_{f}}\ket{\partial_{t_{f}}\varphi_{\alpha,\,t_{f}}},$
(8)

where we have suppressed the dependence on $U$ in the generator $G_{t_{f}}$
for simplicity, and $\alpha=a,\,b$,
$\partial_{t_{f}}G_{t_{f}}=\allowbreak\text{i}\underline{\partial}_{t_{f}}U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}U_{t_{f}}(t_{f})$$+\text{i}U_{t_{f}}^{\dagger}(t_{f})\underline{\partial}_{t_{f}}\partial_{\lambda}U_{t_{f}}(t_{f})$$+U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}H_{\lambda}(t_{f})U_{t_{f}}(t_{f})$,
where $\underline{\partial}_{t_{f}}$ denotes the derivative with respect the
subscript $t_{f}$ instead of the one in the parenthesis Yang  _et al._
(2021b). Generally, Eq. (8) is difficult to solve analytically, while
numerical calculation is tractable. However, if $U_{t_{f}}(\tau)$ and
$\ket{\varphi_{a,\,t_{f}}}$ are independent of the subscript variable $t_{f}$,
Eq. (8) reduces to
$U^{\dagger}(t_{f})\partial_{\lambda}H_{\lambda}(t_{f})U(t_{f})\ket{\varphi_{\alpha}}=\partial_{t_{f}}\mu_{\alpha,\,t_{f}}\ket{\varphi_{\alpha}}$.
It then follows that $\ket{\varphi_{a}(t_{f})}=U(t_{f})\ket{\varphi_{a}}$ is
an eigenstates of $\partial_{\lambda}H_{\lambda}(t_{f})$ for all $t_{f}$, with
eigenvalue $\partial_{t_{f}}\mu_{\alpha,\,t_{f}}$. The solution reduces again
to the Pang-Jordan protocol Pang and Jordan (2017).

This suggests that when the control Hamiltonian is restricted to some non-
trivial subspace of the Hermitian operators, such that
$\partial_{\lambda}\Delta\rho(\tau)$ does not always vanish on $[0,\,t_{f}]$,
both the exact optimal controls and the exact optimal initial states depend on
the probing time $t_{f}$, making it challenging to find them analytically. In
particular, one can show that when only $U(\tau)$ depends on $t_{f}$ but
$\ket{\varphi_{\alpha}}$ is independent of $t_{f}$, this is the case as there
exists a time such that $\partial_{\lambda}\Delta\rho(\tau)\neq 0$ Yang  _et
al._ (2021b) .

Figure 1: QFI for various scenarios. (a) a single qubit with the sensing
Hamiltonian $H_{\lambda}=\lambda\sigma_{z}/2+\Delta\sigma_{x}/2$ (b-d) the
spin chain with the sensing Hamiltonian (13). For all figures (a-d)
$\lambda=\Delta=1$, the frequency of the drive $\omega=1826.67$. The red lines
satisfy the AFM condition (11) or (16) with
$c_{l}^{y}=\tilde{c}_{l}^{z}=c_{l}^{xy}=\tilde{c}_{l}^{zx}=10$ when $1\leq
l\leq 5$ and vanish for $l>5$. The total simulation time is $5000$ times the
fundamental periodic $2\pi/\omega$. For the case of lowest two harmonics in
green lines, $c_{l}^{y},\,\tilde{c}_{l}^{z},\,c_{l}^{xy},\,\tilde{c}_{l}^{zx}$
are non-vanishing only when $l=1,\,2$ and take values $10$.

Approximate solution via Floquet engineering.  We next show that for a time-
independent Hamiltonian $H_{\lambda}$ the restricted optimal controls can be
approximately engineered by the high-frequency control, known as the Floquet
engineering Rahav  _et al._ (2003); Goldman and Dalibard (2014); Bukov  _et
al._ (2015); Eckardt and Anisimovas (2015); Goldman and Dalibard (2015);
Eckardt (2017); Chen  _et al._ (2020). For time-independent Hamiltonians with
unrestricted controls, the minimum optimal control Hamiltonian is the one that
cancels the parts in $H_{\lambda}$ that do not commute with
$\partial_{\lambda}H_{\lambda}$ Pang and Jordan (2017). However, if any of the
required control Hamiltonians is not available, it is not clear how to reach
the Heisenberg scaling, $I_{0}^{\textrm{HS}}=4n^{2}t_{f}^{2}$, where $n$ is
the number of probes. Here, we search for a static control $H_{\text{c},\,0}$
and high-frequency driving controls $H_{\text{c}}^{(d)}(t)$, where
$H_{\text{c}}^{(d)}(t)\equiv\sum_{l\neq 0}H_{\text{c},\,l}e^{\text{i}l\omega
t}$, so that
$H_{\text{tot}}(t)=H_{\lambda}+H_{\text{c},\,0}+H_{\text{c}}^{(d)}(t)$, and
show that the unavailable controls in the minimum optimal control Hamiltonian
can be actually constructed approximately through a high-frequency expansion.
Let us first move to the rotating frame associated with
$\mathcal{U}(t)=e^{\text{i}K(t)}$, where $K(t)$ is the so-called kick operator
Goldman and Dalibard (2014). In it, the total Hamiltonian becomes is given by
the Floquet effective Hamiltonian $H_{\text{F}}$. When the driving frequency
is high enough, an expansion in orders of $1/\omega$ can be performed Goldman
and Dalibard (2014, 2015). To the first-order of $1/\omega$, one finds that
$K(t)=[1/(\text{i}\omega)]\sum_{l\neq
0}\frac{1}{l}\left(H_{\text{c},\,l}e^{\text{i}l\omega
t}-1\right)+O(1/\omega^{2})$ 222We choose the normalization $K(0)=0$ to all
orders of $1/\omega$, which is different from the normalization
$1/T\int_{0}^{T}K(t)dt=0$ used in Ref. Goldman and Dalibard (2014, 2015).
Therefore, the resulting expression of $K(t)$ is different from the one in
Ref. Goldman and Dalibard (2014, 2015), up to some irrelevant constant, which
does not affect the form of the Floquet effective Hamiltonian $H_{\text{F}}$
and

$H_{\text{F}}=H_{\lambda}+H_{\text{c},\,0}+1/\omega\sum_{l=1}^{\infty}[H_{\text{c},\,l},\,H_{\text{c},\,-l}]/l+O(1/\omega^{2}).$
(9)

One can explicitly show that the kick operator $K(t)$ consists of the
$H_{\text{c},\,l}$ with $l\neq 0$ and it is independent of the estimation
parameter. Thus, the QFI remains unchanged in the Floquet rotating frame. The
generator becomes
$G_{t_{f}}=\int_{0}^{t_{f}}e^{-\text{i}H_{\text{F}}\tau}\partial_{\lambda}H_{\lambda}e^{\text{i}H_{\text{F}}\tau}d\tau$
Yang  _et al._ (2021b). The key idea is that the unavailable controls in the
original static frame can be constructed through the commutator in Eq. (9).
This can be best illustrated using a simple qubit example with the Hamiltonian
$H_{\lambda}=\lambda\sigma_{z}/2+\Delta\sigma_{x}/2$ and
$\mathcal{V}_{\text{c}}=\\{\sigma_{y},\,\sigma_{z}\\}$. The term that does not
commute with $\partial_{\lambda}H_{\lambda}$ is $\Delta\sigma_{x}/2$, which is
not available in $\mathcal{V}_{\text{c}}$. Therefore, we consider
$H_{\text{c},\,0}=c_{0}^{y}\sigma_{y}+c_{0}^{z}\sigma_{z}$ and
$H_{\text{c}}^{(d)}(t)=\sum_{l\neq
0}(c_{l}^{y}\sigma_{y}+c_{l}^{z}\sigma_{z})e^{\text{i}l\omega t}$, where for
$l\neq 0$, $c_{l}^{y}=c_{-l}^{y*}$ and $c_{l}^{z}=c_{-l}^{z*}$ to guarantee
the Hermiticity of $H_{\text{c}}^{(d)}(t)$. According to Eq. (9), to the first
order of $1/\omega$, we find

$H_{F}=(\lambda/2+c_{0}^{z})\sigma_{z}+(\Delta/2-4/\omega\sum_{l=1}^{\infty}\text{Im}[c_{l}^{y}c_{l}^{z*}]/l)\sigma_{x}+c_{0}^{y}\sigma_{y}.$
(10)

The approximate optimal control fulfills $c_{0}^{y}=0$ and

$\omega=8/\Delta\sum_{l=1}^{\infty}\text{Im}[c_{l}^{y}c_{l}^{z*}]/l,$ (11)

which we call the amplitude-frequency matching (AFM) condition. The validity
of the high-frequency expansion requires that $\omega$ should be the largest
frequency in the original Hamiltonian, i.e., that
$\omega\gg\lambda,\,\Delta,c_{l}^{y},c_{l}^{z}$. Experimentally, for a laser
frequency that satisfies this condition, one can always tune the amplitude so
that Eq. (11) is satisfied. Conversely, one can also choose a proper laser
frequency for given amplitudes so that Eq. (11) is fulfilled. We emphasize
that when Eq. (11) is satisfied and the initial state in the Floquet rotating
frame is $(\ket{0}+\ket{1})/\sqrt{2}$, the optimal control conditions Eqs. (5,
6) are approximately satisfied up to the order of $1/\omega$ in the Floquet
rotating frame. The initial state in the lab and Floquet rotating frames is
the same as $K(0)=0$. Going back to the lab frame, we generate a solution to
Eqs. (5, 6) when the controls are restricted. This illustrates the power of
our approach beyond the Pang-Jordan protocol Pang and Jordan (2017).

A simple choice for the amplitudes involves taking for $l\geq 1$ both
$c_{l}^{y}$ and $\tilde{c}_{l}^{z}=\text{i}c_{l}^{z}$ to be real. This yields

$H_{\text{c}}^{(d)}(t)=2\sum_{l=1}^{\infty}\left[c_{l}^{y}\cos(l\omega
t)\sigma_{y}+\tilde{c}_{l}^{z}\sin(l\omega t)\sigma_{z}\right].$ (12)

As one can see from Fig. 1(a), for parameters satisfying the AFM condition Eq.
(11), the Heisenberg scaling is achieved. As a result, the first-order term in
the high-frequency expansion makes it possible to construct the $\sigma_{x}$
term by commuting the operators in $\mathcal{V}_{\text{c}}$. This idea can be
generalized to many-body systems, as shown in the following.

Restricted control in a quantum spin chain. Consider the sensing of magnetic
field using a spin chain

$H_{\lambda}=\frac{J}{2}\sum_{i=1}^{n}\sigma_{i}^{x}\sigma_{i+1}^{x}+\frac{\Delta}{2}\sum_{i=1}^{n}\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x}+\frac{\lambda}{2}\sum_{i=1}^{n}\sigma_{i}^{z},$
(13)

which contains both two-body and three-body interactions. We assume periodic
boundary conditions for simplicity and consider the set of allowed controls
$\mathcal{V}_{\text{c}}=\\{\sigma_{i}^{a},\,\sigma_{i}^{a}\sigma_{i+1}^{b}\\}$,
$(a,\,b=x,\,y,\,z)$, involving only one-body and nearest neighbor two-body
operators. When the controls are unrestricted, the minimum optimal control
Hamiltonian contains local two-body and local three-body terms. The first-part
in the control Hamiltonian would cancel the nearest neighbor two-body terms,
i.e., $H_{\text{c},\,0}=-J/2\sum_{i=1}^{n}\sigma_{i}^{x}\sigma_{i+1}^{x}$.
However, the term consisting of three-body operators cannot be canceled
directly through the allowed control set $\mathcal{V}_{\text{c}}$.

Figure 2: Numerical calculation of the QFI with respect to the number of qubit
$n$ for the sensing of magnetic field using the spin-chain Hamiltonian (13).
The initial state preparation and the values of the parameters $\lambda$,
$\Delta$, $\omega$, $c_{l}^{xy}$ and $\tilde{c}_{l}^{zx}$ are the same as Fig.
1 with $J=0$. (b-d) and $t_{f}=17.2$ is the total simulation time in Fig. 1.

Our goal is to construct the three-body operators using the first-order
correction in Eq. (9), that is, we expect to induce
$[H_{l},\,H_{-l}]\propto\sum_{i=1}^{n}\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x}$.
The commutator between one-body and two-body operators cannot produce a three-
body operator. To generate a three-body operator, we must commute the two-body
operators in $\mathcal{V}_{\text{c}}$ Yang  _et al._ (2021b). Therefore, one
can construct $H_{\text{c}}(t)=\sum_{l\neq
0}H_{\text{c},\,l}e^{\text{i}l\omega t}$, where
$H_{\text{c},\,l}=\sum_{i=1}^{n}\left(c_{li}^{xy}\sigma_{i}^{x}\sigma_{i+1}^{y}+c_{li,zx}\sigma_{i}^{z}\sigma_{i+1}^{x}\right)$,
$c_{li}^{xy}=c_{-li}^{xy*}$ and $c_{li}^{zx}=c_{-li}^{zx*}$ to ensure the
Hermiticity of $H_{\text{c}}(t)$ Yang  _et al._ (2021b). When $c_{li}^{xy}$
and $c_{li}^{zx}$ are respectively real and purely imaginary, the commutator
$[H_{\text{c},\,l},\,H_{\text{c},\,-l}]=-4\sum_{i=1}^{n}\text{Im}(c_{li}^{xy}c_{li+1}^{zx*})\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x}$.
The most general choice of the coefficients $c_{li}^{xy}$ and $c_{li}^{zx}$
which generates the $\sigma^{x}$-type three-body interaction is provided in
Sec. V in Yang  _et al._ (2021b). For simplicity, we further assume these
coefficients are homogeneous across the chain and take
$c_{li}^{xy}=c_{l}^{xy}$ and $\tilde{c}_{l}^{zx}=\text{i}c_{li}^{zx}$ to be
real for $l\geq 1$. This yields the high frequency driving control Hamiltonian

$\displaystyle H_{\text{c}}^{(d)}(t)$
$\displaystyle=2\sum_{l=1}^{\infty}\sum_{i=1}^{n}\left[c_{l}^{xy}\cos(l\omega
t)\sigma_{i}^{x}\sigma_{i+1}^{y}+\tilde{c}_{l}^{zx}\sin(l\omega
t)\sigma_{i}^{z}\sigma_{i+1}^{x}\right],$ (14)

while the effective Hamiltonian obtained as the leading term in the high-
frequency expansion becomes

$H_{F}=\frac{\lambda}{2}\sum_{i=1}^{n}\sigma_{i}^{z}+\left[\frac{\Delta}{2}-\frac{4}{\omega}\sum_{l=1}^{\infty}\frac{1}{l}(c_{l}^{xy}\tilde{c}_{l}^{zx})\right]\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x}.$
(15)

The three-body term is then canceled by tuning $c_{l}^{xy}$ and
$\tilde{c}_{l}^{zx}$ such that

$\omega=8/\Delta\sum_{l=1}^{\infty}(c_{l}^{xy}\tilde{c}_{l}^{zx}/l).$ (16)

As with the qubit case, the optimal initial state in both the lab and Floquet
rotating is the GHZ state. Equations (5, 6) are thus approximately satisfied
up to the order of $1/\omega$. As shown in Fig. 1(b-d), for parameters
satisfying the AFM condition (16), the Heisenberg scaling
$I_{0}^{\textrm{HS}}$ is achieved. Furthermore, even when one just takes the
lowest two harmonics in the driving, QFI is not very far below
$I_{0}^{\textrm{HS}}$. In Fig. 2, $I_{0}^{\textrm{HS}}$ is achieved with Eq.
(14) when the parameters satisfy Eq. (16). Again, the precision using only the
lowest two harmonics can approach $I_{0}^{\textrm{HS}}$ very closely.

Finally, we emphasize that a similar technique can be applied to design more
general driving protocols to cancel the effect of other types of three-body
interactions Yang  _et al._ (2021b). Experimental platforms where three-body
interactions either appear naturally or can be potentially engineered -such as
NMR system Tseng  _et al._ (1999); Peng  _et al._ (2009), Kitaev spin liquid
Takahashi  _et al._ (2021), superconducting circuit Petiziol  _et al._ (2021)
and quantum gas systems de Paz  _et al._ (2013); Signoles  _et al._ (2021);
Vishveshwara and Weld (2021); Büchler  _et al._ (2007)\- can be potentially
used to test the metrological protocol discussed here.

In summary, we have introduced a variational approach to quantum parameter
estimation and derived the optimal control equations under which the precision
is optimal when the available control Hamiltonians are limited. This approach
readily yields the optimal initial state and Hamiltonian controls, that are
generally dependent on the probe time, in contrast with the unconstrained
case. The implementation of the constrained optimal protocol in many-body
systems can be eased by Floquet engineering, as we have demonstrated in
applications to magnetometry. We hope that our results inspire new theoretical
and technological advances in quantum metrology with quantum many-body
systems. Many questions are open for further investigation, such as
determining the ultimate scaling bounds of the QFI under restricted local
controls and the application of our method to critically-enhanced quantum
metrology Rams  _et al._ (2018).

Acknowledgement. We are grateful to Aurélia Chenu for critical reading of the
manuscript and helpful feedback. We thank Nicolas P. Bigelow and Fernando J.
Gómez-Ruiz for useful discussions. Support from National Natural Science
Foundation of China (NSFC) grant no. 12075323, NSF grant PHY 1708008, NASA/JPL
RSA 1656126, and US Army Research Office grant no. W911NF-18-10178 is greatly
acknowledged.

## References

  * Helstrom (1976) C. W. Helstrom,  _Quantum Detection and Estimation Theory_ (Academic Press, 1976).
  * Holevo (2011) A. S. Holevo,  _Probabilistic and Statistical Aspects of Quantum Theory_ (Springer Science & Business Media, 2011).
  * Braunstein and Caves (1994) S. L. Braunstein and C. M. Caves, Physical Review Letters 72, 3439 (1994).
  * Paris (2009) M. G. A. Paris, International Journal of Quantum Information 07, 125 (2009).
  * Giovannetti  _et al._ (2011) V. Giovannetti, S. Lloyd, and L. Maccone, Nature Photonics 5, 222 (2011).
  * Yuan and Fung (2015) H. Yuan and C.-H. F. Fung, Physical Review Letters 115, 110401 (2015).
  * Yuan (2016) H. Yuan, Physical Review Letters 117, 160801 (2016).
  * Pang and Jordan (2017) S. Pang and A. N. Jordan, Nature Communications 8, 14695 (2017).
  * Yang  _et al._ (2017) J. Yang, S. Pang, and A. N. Jordan, Physical Review A 96, 020301 (2017).
  * Liu and Yuan (2017) J. Liu and H. Yuan, Physical Review A 96, 012117 (2017).
  * Yang  _et al._ (2021a) J. Yang, S. Pang, A. del Campo, and A. N. Jordan, arXiv:2104.07120 [cond-mat, physics:quant-ph] (2021a), arXiv:2104.07120 [cond-mat, physics:quant-ph] .
  * Gefen  _et al._ (2019) T. Gefen, A. Rotem, and A. Retzker, Nature Communications 10, 4992 (2019).
  * Naghiloo  _et al._ (2017) M. Naghiloo, A. N. Jordan, and K. W. Murch, Phys. Rev. Lett. 119, 180801 (2017).
  * Demirplak and Rice (2003) M. Demirplak and S. A. Rice, The Journal of Physical Chemistry A 107, 9937 (2003).
  * Demirplak and Rice (2005) M. Demirplak and S. A. Rice, The Journal of Physical Chemistry B 109, 6838 (2005).
  * Demirplak and Rice (2008) M. Demirplak and S. A. Rice, The Journal of Chemical Physics 129, 154111 (2008).
  * Berry (2009) M. V. Berry, Journal of Physics A: Mathematical and Theoretical 42, 365303 (2009).
  * Guéry-Odelin  _et al._ (2019) D. Guéry-Odelin, A. Ruschhaupt, A. Kiely, E. Torrontegui, S. Martínez-Garaot, and J. G. Muga, Rev. Mod. Phys. 91, 045001 (2019).
  * Cabedo-Olaya  _et al._ (2020) M. Cabedo-Olaya, J. G. Muga, and S. Martínez-Garaot, Entropy 22, 1251 (2020).
  * del Campo  _et al._ (2012) A. del Campo, M. M. Rams, and W. H. Zurek, Phys. Rev. Lett. 109, 115703 (2012).
  * Carlini  _et al._ (2007) A. Carlini, A. Hosoya, T. Koike, and Y. Okudaira, Physical Review A 75, 042308 (2007).
  * Carlini  _et al._ (2006) A. Carlini, A. Hosoya, T. Koike, and Y. Okudaira, Physical Review Letters 96, 060503 (2006).
  * Rezakhani  _et al._ (2009) A. T. Rezakhani, W.-J. Kuo, A. Hamma, D. A. Lidar, and P. Zanardi, Physical Review Letters 103, 080502 (2009).
  * Opatrný and Mølmer (2014) T. Opatrný and K. Mølmer, New J. Phys. 16, 015025 (2014).
  * Saberi  _et al._ (2014) H. Saberi, K. Opatrný, T.and Mølmer, and A. del Campo, Phys. Rev. A 90, 060301 (2014).
  * Sels and Polkovnikov (2017) D. Sels and A. Polkovnikov, Proceedings of the National Academy of Sciences 114, E3909 (2017).
  * Claeys  _et al._ (2019) P. W. Claeys, M. Pandey, D. Sels, and A. Polkovnikov, Physical Review Letters 123, 090602 (2019).
  * Chandarana  _et al._ (2021) P. Chandarana, N. N. Hegade, K. Paul, F. AlbarrÃ¡n-Arriagada, E. Solano, A. del Campo, and X. Chen, “Digitized-counterdiabatic quantum approximate optimization algorithm,” (2021), arXiv:2107.02789 [quant-ph] .
  * Wang  _et al._ (2015) X. Wang, M. Allegra, K. Jacobs, S. Lloyd, C. Lupo, and M. Mohseni, Physical Review Letters 114, 170501 (2015).
  * Kaubruegger  _et al._ (2021) R. Kaubruegger, D. V. Vasilyev, M. Schulte, K. Hammerer, and P. Zoller, Phys. Rev. X 11, 041045 (2021).
  * Marciniak  _et al._ (2021) C. D. Marciniak, T. Feldker, I. Pogorelov, R. Kaubruegger, D. V. Vasilyev, R. van Bijnen, P. Schindler, P. Zoller, R. Blatt, and T. Monz, arXiv preprint arXiv:2107.01860 (2021).
  * Rahav  _et al._ (2003) S. Rahav, I. Gilary, and S. Fishman, Physical Review Letters 91, 110404 (2003).
  * Goldman and Dalibard (2014) N. Goldman and J. Dalibard, Physical Review X 4, 031027 (2014).
  * Bukov  _et al._ (2015) M. Bukov, L. D’Alessio, and A. Polkovnikov, Advances in Physics 64, 139 (2015).
  * Eckardt and Anisimovas (2015) A. Eckardt and E. Anisimovas, New Journal of Physics 17, 093039 (2015).
  * Goldman and Dalibard (2015) N. Goldman and J. Dalibard, Physical Review X 5, 029902 (2015).
  * Eckardt (2017) A. Eckardt, Reviews of Modern Physics 89, 011004 (2017).
  * Chen  _et al._ (2020) Z. Chen, J. D. Murphree, and N. P. Bigelow, Physical Review A 101, 013606 (2020).
  * Pang and Brun (2014) S. Pang and T. A. Brun, Physical Review A 90 (2014), 10.1103/PhysRevA.90.022117.
  * Yang  _et al._ (2021b) J. Yang, S. Pang, Z. Chen, A. Jordan, and A. del Campo,  _Supplemental Material_ (2021b).
  * Braunstein  _et al._ (1996) S. L. Braunstein, C. M. Caves, and G. J. Milburn, annals of physics 247, 135 (1996).
  * Boixo  _et al._ (2007) S. Boixo, S. T. Flammia, C. M. Caves, and J. Geremia, Physical Review Letters 98, 090401 (2007).
  * Note (1) This is the approach pursued in e.g., the variational approach in shortcut-to-adiabaticity Sels and Polkovnikov (2017); Claeys  _et al._ (2019). Another way of handling the constraints on the control Hamiltonian is by introducing the following constraints $f_{j}(H_{\text{c}}(\tau))=\text{Tr}\left\\{H_{\text{c}}(\tau)\mathcal{X}_{j}\right\\}=0,\,j=d_{\text{c}}+1,\,\cdots N$ to disallow the terms $\mathcal{X}_{j}$, where $j=1,\,2\,\cdots d_{\text{c}}$. This way of introducing the constraint is the one used in e.g. quantum brachistochrone equation Carlini  _et al._ (2007, 2006). However, in many-body quantum metrology, the number of disallowed nonlocal operators is much more than the allowed local operators. Therefore the second approach may introduce an intractable number of constraints and we shall pursue the first approach of expanding $H_{\text{c}}$ in terms of basis operators in the main text.
  * Note (2) We choose the normalization $K(0)=0$ to all orders of $1/\omega$, which is different from the normalization $1/T\intop\ilimits@_{0}^{T}K(t)dt=0$ used in Ref. Goldman and Dalibard (2014, 2015). Therefore, the resulting expression of $K(t)$ is different from the one in Ref. Goldman and Dalibard (2014, 2015), up to some irrelevant constant, which does not affect the form of the Floquet effective Hamiltonian $H_{\text{F}}$.
  * Tseng  _et al._ (1999) C. H. Tseng, S. Somaroo, Y. Sharf, E. Knill, R. Laflamme, T. F. Havel, and D. G. Cory, Phys. Rev. A 61, 012302 (1999).
  * Peng  _et al._ (2009) X. Peng, J. Zhang, J. Du, and D. Suter, Phys. Rev. Lett. 103, 140501 (2009).
  * Takahashi  _et al._ (2021) M. O. Takahashi, M. G. Yamada, D. Takikawa, T. Mizushima, and S. Fujimoto, Phys. Rev. Research 3, 023189 (2021).
  * Petiziol  _et al._ (2021) F. Petiziol, M. Sameti, S. Carretta, S. Wimberger, and F. Mintert, Phys. Rev. Lett. 126, 250504 (2021).
  * de Paz  _et al._ (2013) A. de Paz, A. Sharma, A. Chotia, E. Maréchal, J. H. Huckans, P. Pedri, L. Santos, O. Gorceix, L. Vernac, and B. Laburthe-Tolra, Phys. Rev. Lett. 111, 185305 (2013).
  * Signoles  _et al._ (2021) A. Signoles, T. Franz, R. Ferracini Alves, M. Gärttner, S. Whitlock, G. Zürn, and M. Weidemüller, Phys. Rev. X 11, 011011 (2021).
  * Vishveshwara and Weld (2021) S. Vishveshwara and D. M. Weld, Phys. Rev. A 103, L051301 (2021).
  * Büchler  _et al._ (2007) H. P. Büchler, A. Micheli, and P. Zoller, Nature Physics 3, 726 (2007).
  * Rams  _et al._ (2018) M. M. Rams, P. Sierant, O. Dutta, P. Horodecki, and J. Zakrzewski, Phys. Rev. X 8, 021022 (2018).

Supplemental Materials

###### Contents

  1. I Optimization over the initial state
  2. II Variational Calculus
  3. III  Dependence of $U_{t_{f}}(\tau)$ on the probe time $t_{f}$
  4. IV Generator and optimal control in parameter-independent rotating frame
  5. V Commutators in many-spin systems
  6. VI More general driving protocols

##  I Optimization over the initial state

Our goal is to maximize $\text{Var}[\mathcal{O}]\big{|}_{\ket{\varphi_{0}}}$
over all the quantum states $\ket{\varphi_{0}}$, where $\mathcal{O}$ is some
Hermitian operator. We note the following gauge invariance

$\displaystyle\mathcal{O}$
$\displaystyle\to\tilde{\mathcal{O}}\equiv\mathcal{O}-c,$ (S1)
$\displaystyle\text{Var}[\mathcal{O}]\big{|}_{\ket{\varphi_{0}}}$
$\displaystyle\to\text{Var}[\tilde{\mathcal{O}}]\big{|}_{\ket{\varphi_{0}}}=\text{Var}[\mathcal{O}]\big{|}_{\ket{\varphi_{0}}},$
(S2)

where $c$ is a constant. Introducing a Lagrange multiplier to account for the
normalization of $\ket{\varphi_{0}}$, we aim at optimizing the following
function

$\displaystyle f$
$\displaystyle\equiv\text{Var}[\tilde{\mathcal{O}}]\big{|}_{\ket{\varphi_{0}}}-\mu\braket{\varphi_{0}\big{missing}}{\varphi_{0}}$
$\displaystyle=\braket{\varphi_{0}\big{missing}}{\tilde{\mathcal{O}}^{2}\big{missing}}{\varphi_{0}}-\braket{\varphi_{0}\big{missing}}{\tilde{\mathcal{O}}\big{missing}}{\varphi_{0}}^{2}-\mu\braket{\varphi_{0}\big{missing}}{\varphi_{0}}.$
(S3)

One can choose
$c=\braket{\varphi_{0}\big{missing}}{\mathcal{O}\big{missing}}{\varphi_{0}}$,
so that $f$ becomes

$f=\braket{\varphi_{0}\big{missing}}{\tilde{\mathcal{O}}^{2}\big{missing}}{\varphi_{0}}-\mu\braket{\varphi_{0}\big{missing}}{\varphi_{0}}.$
(S4)

Taking derivatives with respect to $\ket{\varphi_{0}}$, one can easily obtain

$\tilde{\mathcal{O}}^{2}\ket{\varphi_{0}}=\mu\ket{\varphi_{0}}.$ (S5)

All eigenstates of $\mathcal{\tilde{\mathcal{O}}}$ must be an eigenstate of
$\tilde{\mathcal{O}}^{2}$, but not vice versa. This is because if there is a
pair of opposite eigenvalues of $\mathcal{\tilde{\mathcal{O}}}$, the
eigenstate of $\tilde{\mathcal{O}}^{2}$ could be the superposition of the
eigenstates of $\tilde{\mathcal{O}}$ corresponding to these two opposite
eigenvalues. Therefore, either $\ket{\varphi_{0}}$ is an eigenstate of
$\tilde{\mathcal{O}}$ or $\ket{\varphi_{0}}$ is the superposition of two
eigenstates of $\tilde{\mathcal{O}}$ with opposite eigenvalues. In the former
case $f=0$ corresponds to the extremal minimum value. In the latter case,
i.e.,

$\ket{\varphi_{0}}=\alpha\ket{\varphi_{a}}+\beta\ket{\varphi_{b}},$ (S6)

where

$\displaystyle\mathcal{O}\ket{\varphi_{a}}$
$\displaystyle=\mu_{a}\ket{\varphi_{a}},$ (S7)
$\displaystyle\mathcal{O}\ket{\varphi_{b}}$
$\displaystyle=\mu_{b}\ket{\varphi_{b}},$ (S8)

$\mu_{a},\,\mu_{b}$ are not equal to each other and satisfy

$\mu_{a}+\mu_{b}-2\braket{\varphi_{0}\big{missing}}{\mathcal{O}\big{missing}}{\varphi_{0}}=0.$
(S9)

Then, the function $f$ becomes

$f(a,\,b)=a\mu_{a}^{2}+b\mu_{b}^{2}-(a\mu_{a}+b\mu_{b})^{2},$ (S10)

where $a\equiv|\alpha|^{2}$ and $b\equiv|\beta|^{2}$. We proceed to introduce
a Lagrangian multiplier to perform further maximization over $a$ and $b$

$g(a,\,b)=a\mu_{a}^{2}+b\mu_{b}^{2}-(a\mu_{a}+b\mu_{b})^{2}+\xi(a+b-1).$ (S11)

Setting $\partial_{a}g(a,\,b)=\partial_{b}g(a,\,b)=0$, one find $a=b=1/2$.
Therefore, we find $f(a,\,b)=(\mu_{a}-\mu_{b})^{2}/4$. Finally, optimizing
over all the possible distinct eigenvalues of $\mathcal{O}$, we find the
global maximum values of $f$. That is, $f_{\max}=(\mu_{+}-\mu_{-})^{2}/4$,
with the corresponding optimal initial state
$\ket{\varphi_{0}}=(\ket{\varphi_{+}}+e^{\text{i}\theta}\ket{\varphi_{-}})/\sqrt{2}$,
where $\ket{\varphi_{\pm}}$ correspond to the maximum and minimum eigenvalues
$\mu_{\pm}$, respectively, and $\theta$ is an arbitrary phase.

##  II Variational Calculus

For the variation with respect to $U$ we have

$\delta S_{\textrm{M}}=\delta S_{\textrm{I}}+\delta S_{\textrm{S}}=0,$ (S12)

where

$\displaystyle\delta S_{\textrm{I}}$
$\displaystyle=-\int_{0}^{t_{f}}\text{Tr}\left\\{\Delta\rho
U^{\dagger}(\tau)\delta
U(\tau)U^{\dagger}(\tau)\partial_{\lambda}H_{\lambda}(\tau)U(\tau)\right\\}d\tau+\int_{0}^{t_{f}}\text{Tr}\left\\{\Delta\rho
U^{\dagger}(\tau)\partial_{\lambda}H_{\lambda}(\tau)\delta
U(\tau)\right\\}d\tau$
$\displaystyle=\int_{0}^{t_{f}}\text{Tr}\left\\{U^{\dagger}[\Delta\rho(\tau),\,\partial_{\lambda}H_{\lambda}(\tau)]\delta
U(\tau)\right\\}d\tau,$ (S13) $\displaystyle\delta S_{\textrm{S}}$
$\displaystyle=\text{i}\text{Tr}\left\\{U^{\dagger}(\tau)\Lambda(\tau)\delta
U(\tau)\right\\}\big{|}_{\tau=0}^{\tau=t_{f}}-\text{i}\int_{0}^{t_{f}}d\tau\text{Tr}\left\\{\left[U^{\dagger}(\tau)\Lambda(\tau)\dot{U}(\tau)U^{\dagger}(\tau)+\frac{d}{d\tau}\left(U^{\dagger}(\tau)\Lambda(\tau)\right)\right]\delta
U(\tau)\right\\}$
$\displaystyle=\text{i}\text{Tr}\left\\{U^{\dagger}(\tau)\Lambda(\tau)\delta
U(\tau)\right\\}\big{|}_{\tau=0}^{\tau=t_{f}}-\text{i}\int
dt\text{Tr}\left\\{U^{\dagger}\left[\Lambda\dot{U}U^{\dagger}+\dot{\Lambda}+U\dot{U}^{\dagger}\Lambda\right]\delta
U(\tau)\right\\},$ (S14)

where we have used the fact that $\delta U^{\dagger}=-U^{\dagger}\delta
UU^{\dagger}$ in Eq. (S13). Eq. (S12) yields

$\dot{\Lambda}(\tau)-\text{i}[\Lambda(\tau),\,H_{\text{tot}}(\tau)]+\text{i}[\Delta\rho(\tau),\,\partial_{\lambda}H_{\lambda}(\tau)]=0,$
(S15)

where

$\displaystyle\Delta\rho(\tau)$ $\displaystyle\equiv U(\tau)\Delta\rho
U^{\dagger}(\tau),$ (S16) $\displaystyle H_{\text{tot}}(\tau)$
$\displaystyle\equiv\text{i}\dot{U}(\tau)U^{\dagger}(\tau)=H_{\lambda}(\tau)+H_{\text{c}}(\tau),$
(S17)

with the essential boundary condition

$\Lambda(t_{f})=0.$ (S18)

We observe that Eq. (S15) may be solved in terms of the evolution operator
because Eq. (S15) is a driven Liouville equation and the Green’s function is
the product of unitary evolution operators for the forward and backward
directions multiplied by the Heaviside function. To see this in a clear
manner, let us note that

$\Lambda_{\text{hom}}(\tau)=U(\tau,\,t_{f})\Lambda(t_{f})U^{\dagger}(\tau,\,t_{f})=0,$
(S19)

where $U(\tau,\,s)$ is the evolution operator from time $s$ to time $\tau$,
which satisfies

$\text{i}\dot{U}(\tau,\,s)=H_{\text{tot}}(\tau)U(\tau,\,s),$ (S20)

where the time derivative is applied to the variable $\tau$. One can
explicitly check that for the impulse applied at time $s$, the Green’s
function is

$\mathscr{G}(\tau,\,s)=-U(\tau,\,s)\Theta(s-\tau)U^{\dagger}(\tau,\,s),\,s\in[0,\,t_{f}],$
(S21)

which satisfies

$\partial_{\tau}\mathscr{G}(\tau,\,s)=\text{i}[\mathscr{G}(\tau,\,s),\,H_{\text{tot}}(\tau)]+\text{$\delta$}(\tau-s),$
(S22)

and the homogeneous initial condition $\mathscr{G}(t_{f},\,s)=0$. Note that we
do not contract $U(\tau,\,s)$ and $U^{\dagger}(\tau,\,s)$ in Eq. (S21).
Therefore, we find

$\displaystyle\Lambda(\tau)$
$\displaystyle=\Lambda_{\text{inhom}}(\tau)=\text{i}\int_{0}^{t_{f}}U(\tau,\,s)[\Delta\rho(s),\,\partial_{\lambda}H_{\lambda}(s)]\Theta(s-\tau)U^{\dagger}(\tau,\,s)ds$
$\displaystyle=\text{i}U(\tau)[\Delta\rho,\,\int_{\tau}^{t_{f}}U^{\dagger}(s)\partial_{\lambda}H_{\lambda}(s)U(s)ds]U^{\dagger}(\tau)$
$\displaystyle=\text{i}U(\tau)[\Delta\rho,\,G_{t_{f}}[U]-G_{\tau}[U]]U^{\dagger}(\tau)$
$\displaystyle=-\text{i}U(\tau)[\Delta\rho,\,G_{\tau}[U]]U^{\dagger}(\tau),$
(S23)

where we have used the fact that $[\Delta\rho,\,G_{t_{f}}[U]]=0$. We note the
alternative concise expression for the generator Pang and Jordan (2017)

$G_{\tau}[U]=\text{i}U^{\dagger}(\tau)\partial_{\lambda}U(\tau).$ (S24)

Substituting Eq. (S24) into Eq. (S23) yields,

$\displaystyle\Lambda(\tau)$ $\displaystyle=U(\tau)\Delta\rho
U^{\dagger}(\tau)\partial_{\lambda}U(\tau)U^{\dagger}(\tau)-\partial_{\lambda}U(\tau)\Delta\rho
U^{\dagger}(\tau)$
$\displaystyle=-[U(\tau)\Delta\rho\partial_{\lambda}U^{\dagger}(\tau)+\partial_{\lambda}U(\tau)\Delta\rho
U^{\dagger}(\tau)]=-\partial_{\lambda}[\Delta\rho(\tau)],$ (S25)

where we have used
$\partial_{\lambda}U(\tau)U^{\dagger}(\tau)=-U(\tau)\partial_{\lambda}U^{\dagger}(\tau)$.
Therefore, the trace condition
$\text{Tr}\left\\{\Lambda(\tau)\mathcal{X}_{i}\right\\}=0$ becomes Eq. (6) in
the main text.

Mathematically, we note that $\Delta\rho$ is independent of the estimation of
$\lambda$ from the very beginning when we construct the metrological action
Eq. (3) in the main text. Practically, after solving the optimization problem,
one may find that the optimal value of $\Delta\rho$ depends on $\lambda$,
since $\Delta\rho$ is related the maximum and minimum eigenvectors of
$G_{t_{f}}[U]$, which depends on $\lambda$ in general. If this is the case,
the optimal initial state should be prepared with the prior knowledge of the
estimation parameter, which is very close to the true value of $\lambda$, as
we focus on the ultra-sensitive estimation regime.

##  III  Dependence of $U_{t_{f}}(\tau)$ on the probe time $t_{f}$

Let us first calculate explicitly $\partial_{t_{f}}G_{t_{f}}$. We note the
generator $G_{t_{f}}$ has the following two alternative forms

$G_{t_{f}}=\text{i}U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}U_{t_{f}}(t_{f})=\int_{0}^{t_{f}}U_{t_{f}}^{\dagger}(\tau)\partial_{\lambda}H_{\lambda}(\tau)U_{t_{f}}(\tau)d\tau.$
(S26)

Then it is straightforward to calculate

$\displaystyle\partial_{t_{f}}G_{t_{f}}$
$\displaystyle=\int_{0}^{t_{f}}\underline{\partial}_{t_{f}}\left[U_{t_{f}}^{\dagger}(\tau)\partial_{\lambda}H_{\lambda}(\tau)U_{t_{f}}(\tau)\right]d\tau+U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}H_{\lambda}(t_{f})U_{t_{f}}(t_{f})$
$\displaystyle=\underline{\partial}_{t_{f}}\left[\text{i}U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}U_{t_{f}}(t_{f})\right]+U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}H_{\lambda}(t_{f})U_{t_{f}}(t_{f})$
$\displaystyle=\text{i}\underline{\partial}_{t_{f}}U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}U_{t_{f}}(t_{f})+\text{i}U_{t_{f}}^{\dagger}(t_{f})\underline{\partial}_{t_{f}}\partial_{\lambda}U_{t_{f}}(t_{f})+U_{t_{f}}^{\dagger}(t_{f})\partial_{\lambda}H_{\lambda}(t_{f})U_{t_{f}}(t_{f})$
(S27)

where $\underline{\partial}_{t_{f}}$ denotes the derivative with respect the
subscript $t_{f}$ instead of the one in the parenthesis.

When $\ket{\varphi_{a,\,t_{f}}}$ is independent of $t_{f}$, after multiplying
both sides by $U_{t_{f}}(t_{f})$ from the left, Eq. (8) in the main text
becomes

$\partial_{t_{f}}\tilde{G}_{t_{f}}\ket{\varphi_{a,\,t_{f}}(t_{f})}=\partial_{t_{f}}\mu_{a,\,t_{f}}\ket{\varphi_{a,\,t_{f}}(t_{f})},\,\forall
t_{f},\,a=\pm,$ (S28)

where $\ket{\varphi_{a,\,t_{f}}(t_{f})}=U_{t_{f}}(t_{f})\ket{\varphi_{a}}$,
and

$\displaystyle\partial_{t_{f}}\tilde{G}_{t_{f}}$ $\displaystyle\equiv
U_{t_{f}}(t_{f})\partial_{t_{f}}G_{t_{f}}U_{t_{f}}^{\dagger}(t_{f})$
$\displaystyle=\text{i}\underline{\partial}_{t_{f}}U_{t_{f}}(t_{f})\partial_{\lambda}U_{t_{f}}^{\dagger}(t_{f})+\text{i}\underline{\partial}_{t_{f}}\partial_{\lambda}U_{t_{f}}(t_{f})U_{t_{f}}^{\dagger}(t_{f})+\partial_{\lambda}H_{\lambda}(t_{f})$
$\displaystyle=\text{i}\partial_{\lambda}\left[\underline{\partial}_{t_{f}}U_{t_{f}}(t_{f})U_{t_{f}}^{\dagger}(t_{f})\right]+\partial_{\lambda}H_{\lambda}(t_{f})$
$\displaystyle=\partial_{\lambda}\left[H_{\lambda}(t_{f})+A_{\lambda,\,t_{f}}(t_{f})\right].$
(S29)

Here, $A_{\lambda,\,t_{f}}(t_{f})$ is defined as

$A_{\lambda,\,t_{f}}(\tau)\equiv-\text{i}U_{t_{f}}(\tau)\underline{\partial}_{t_{f}}U_{t_{f}}^{\dagger}(\tau),$
(S30)

which bears the same form as $G_{t_{f}}$ and therefore has the integral
representation

$A_{\lambda,\,t_{f}}(\tau)=-\int_{0}^{\tau}U_{t_{f}}^{\dagger}(s)\partial_{t_{f}}H_{\text{c},\,t_{f}}(s)U_{t_{f}}(s)ds.$
(S31)

The physical meaning of $A_{\lambda,\,t_{f}}(\tau)$ is the following: When
$A_{\lambda,\,t_{f}}(\tau)$ vanishes as in the Pang-Jordan protocol Pang and
Jordan (2017), Eq. (S28) implies that $\ket{\varphi_{\alpha,\,t_{f}}(\tau)}$
is always an eigenstates of $\partial_{\lambda}H_{\lambda}(\tau)$, as can be
seen from Eq. (7) in the main text, given that $\Delta\rho_{t_{f}}(\tau)$
always commutes with $\partial_{\lambda}H_{\lambda}(\tau)$ for any
$\tau\in[0,\,t_{f}]$. Eq. (S28) indicates that, as long as long as
$[\partial_{\lambda}A_{\lambda,\,t_{f}}(t_{f}),\,\partial_{\lambda}H(t_{f})]\neq
0$, $\Delta\rho_{t_{f}}(t_{f})$ does not commute with
$\partial_{\lambda}H_{\lambda}(t_{f})$ in general, which in turn implies that
$\partial_{\lambda}[\Delta\rho_{t_{f}}(\tau)]$ does not vanishes at all times.
Therefore Eq. (6) in the main text becomes non-trivial in the sense that
$\mathcal{X}_{j}$ must be orthogonal to the non-vanishing values of
$\partial_{\lambda}\Delta\rho_{t_{f}}(\tau)$ at least for some time $\tau$.

##  IV Generator and optimal control in parameter-independent rotating frame

For a parameter-independent unitary operator $\mathcal{U}(t)$, consider the
interaction frame associated with $\mathcal{U}(t)$ and make the transformation

$\tilde{U}(t)=\mathcal{U}(t)U(t).$ (S32)

Then, the Schrödinger equation becomes

$\text{i}\partial_{t}\tilde{U}(t)=\tilde{H}_{\text{tot}}(t)\tilde{U}(t),$
(S33)

where

$\tilde{H}_{\text{tot}}(t)=\mathcal{U}(t)H_{\text{tot}}(t)\mathcal{U}^{\dagger}(t)-\text{i}\mathcal{U}(t)\partial_{t}\mathcal{U}^{\dagger}(t).$
(S34)

Using Eqs. (S32, S34), the generator can be rewritten as

$\displaystyle G_{t_{f}}$
$\displaystyle=\int_{0}^{t_{f}}\tilde{U}^{\dagger}(\tau)\mathcal{U}(\tau)\partial_{\lambda}H_{\lambda}\mathcal{U}^{\dagger}(\tau)\tilde{U}(\tau)d\tau$
$\displaystyle=\int_{0}^{t_{f}}\tilde{U}^{\dagger}(\tau)\partial_{\lambda}\left[\mathcal{U}(\tau)H_{\text{tot}}(\tau)\mathcal{U}^{\dagger}(\tau)-\text{i}\mathcal{U}(\tau)\partial_{t}\mathcal{U}^{\dagger}(\tau)\right]\tilde{U}(\tau)d\tau$
$\displaystyle=\int_{0}^{t_{f}}\tilde{U}^{\dagger}(\tau)\partial_{\lambda}\tilde{H}_{\text{tot}}(\tau)\tilde{U}(\tau)d\tau,$
(S35)

where we have used the important fact that $\mathcal{U}(t)$ is independent of
the parameter $\lambda$. Eq. (6) in the interaction frame can be rewritten as

$\text{Tr}\left\\{\mathcal{X}_{i}\partial_{\lambda}[\Delta\rho(\tau)]\right\\}=\text{Tr}\left\\{\tilde{\mathcal{X}}_{i}(\tau)\partial_{\lambda}[\Delta\tilde{\rho}(\tau)]\right\\},$
(S36)

where

$\displaystyle\tilde{\mathcal{X}}_{i}(\tau)$
$\displaystyle\equiv\mathcal{U}(\tau)\mathcal{X}_{i}\mathcal{U}^{\dagger}(\tau),$
(S37) $\displaystyle\Delta\tilde{\rho}(\tau)$
$\displaystyle\equiv\tilde{U}(\tau)\Delta\rho\tilde{U}^{\dagger}(\tau).$ (S38)

##  V Commutators in many-spin systems

In this section, we discuss several properties for many-spin systems. We
define the $k$-body Hermitian basis operator

$\mathcal{X}_{k}^{\alpha_{k}}=\sigma_{i_{1}}^{\alpha_{i_{1}}}\sigma_{i_{2}}^{\alpha_{i_{2}}}\cdots\sigma_{i_{k}}^{\alpha_{i_{k}}},$
(S39)

where the distinct indices $i_{1}<i_{2},\,\cdots<i_{k}$ take values in
$[1,\,n]$ , $\alpha_{i_{k}}\in\\{x,\,y,\,z\\}$ and
$\alpha_{k}\equiv(\alpha_{i_{1}},\,\cdots,\alpha_{i_{k}})$. The $0-$body
operator $\mathcal{X}_{0}$ is defined to be a constant. We first discuss a
lemma concerning the general property of the commutator between a $k-$body
operator and a $l-$body operator. We denote by $p$ the number of identical
subscripts for two many-body operators $\mathcal{X}_{k}^{\alpha_{k}}$ and
$\mathcal{X}_{l}^{\alpha_{l}}$. For $p=0$, we know
$[\mathcal{X}_{k}^{\alpha_{k}},\,\mathcal{X}_{l}^{\alpha_{l}}]=0$. For
positive $p$, we have the following lemma:

######  Lemma 1.

Consider two many-body spin operators, where
$\mathcal{X}_{k}=\sigma_{i_{1}}^{\alpha_{i_{1}}}\sigma_{i_{2}}^{\alpha_{i_{2}}}\cdots\sigma_{i_{k}}^{\alpha_{i_{k}}}$
and
$\mathcal{X}_{l}=\sigma_{j_{1}}^{\alpha_{j_{1}}}\sigma_{j_{2}}^{\alpha_{j_{2}}}\cdots\sigma_{j_{l}}^{\alpha_{j_{l}}}$.
Without loss of generality, assume $i_{1}=j_{1},\cdots,\,$ $i_{p}=j_{p}$ with
the remaining subscripts being distinct and
$\alpha_{i_{1}}=\alpha_{j_{1}},\,\cdots,\,$ $\alpha_{i_{q}}=\alpha_{j_{q}}$
with the remaining $p-q$ pairs of superscripts being distinct, i.e.,
$\alpha_{i_{q+1}}\neq\alpha_{j_{q+1}},\,\cdots,\,$
$\alpha_{i_{p}}\neq\alpha_{j_{p}}$. Then,

$[\mathcal{X}_{k},\,\mathcal{X}_{l}]=\begin{cases}0&p-q=\text{even}\\\
\mathcal{X}_{k+l-(p-q)}&p-q=\text{odd}\end{cases},$ (S40)

where $p\in[0,\,\min\\{k,\,l\\}]$ and $q\in[0,\,p]$.

###### Proof.

It is straightforward calculate

$[\mathcal{X}_{k},\,\mathcal{X}_{l}]=\left[\prod_{r=1}^{p}\left(\sigma_{i_{r}}^{\alpha_{i_{r}}}\sigma_{j_{r}}^{\alpha_{j_{r}}}\right)-\prod_{r=1}^{p}\left(\sigma_{j_{r}}^{\alpha_{j_{r}}}\sigma_{i_{r}}^{\alpha_{i_{r}}}\right)\right]\sigma_{i_{p+1}}^{\alpha_{i_{p+1}}}\dots\sigma_{i_{k}}^{\alpha_{i_{k}}}\sigma_{i_{p+1}}^{\alpha_{i_{p+1}}}\dots\sigma_{i_{l}}^{\alpha_{i_{l}}}.$
(S41)

With the identity
$\sigma_{i_{r}}^{\alpha_{i_{r}}}\sigma_{i_{r}}^{\alpha_{j_{r}}}=\delta^{\alpha_{i_{r}}\alpha_{j_{r}}}+\text{i}\epsilon_{\alpha_{i_{r}}\alpha_{j_{r}}\beta}\sigma_{i_{r}}^{\beta}$,
where $\epsilon_{\alpha_{i_{r}}\alpha_{j_{r}}\beta}$ is the Levi-Civita
symbol, we find

$\sigma_{i_{r}}^{\alpha_{i_{r}}}\sigma_{i_{r}}^{\alpha_{j_{r}}}=\begin{cases}1&\alpha_{i_{r}}=\alpha_{j_{r}}\\\
-\sigma_{i_{r}}^{\alpha_{j_{r}}}\sigma_{i_{r}}^{\alpha_{i_{r}}}\sim\sigma_{i_{r}}^{\beta_{i_{r}}}&a_{i_{r}}\neq
a_{j_{r}}\end{cases}.$ (S42)

Therefore,

$\displaystyle\prod_{r=1}^{p}\left(\sigma_{i_{r}}^{\alpha_{i_{r}}}\sigma_{j_{r}}^{\alpha_{j_{r}}}\right)$
$\displaystyle=\prod_{r=q+1}^{p}\left(\sigma_{i_{r}}^{a_{i_{r}}}\sigma_{j_{r}}^{\alpha_{j_{r}}}\right)=\prod_{r=q+1}^{p}\left(\sigma_{i_{r}}^{\alpha_{j_{r}}}\sigma_{i_{r}}^{\alpha_{i_{r}}}\right)(-1)^{p-q},$
(S43)

$\prod_{r=1}^{p}\left(\sigma_{j_{r}}^{\alpha_{j_{r}}}\sigma_{i_{r}}^{\alpha_{i_{r}}}\right)=\prod_{r=q+1}^{p}\left(\sigma_{j_{r}}^{\alpha_{j_{r}}}\sigma_{i_{r}}^{\alpha_{i_{r}}}\right),$
(S44)

and we conclude that

$\prod_{r=1}^{p}\left(\sigma_{i_{r}}^{\alpha_{i_{r}}}\sigma_{j_{r}}^{a_{j_{r}}}\right)-\prod_{r=1}^{p}\left(\sigma_{j_{r}}^{\alpha_{j_{r}}}\sigma_{i_{r}}^{\alpha_{i_{r}}}\right)=\begin{cases}0&p-q=\text{even}\\\
\prod_{r=q+1}^{p}\sigma_{i_{r}}^{\beta_{i_{r}}}&p-q=\text{odd}\end{cases},$
(S45)

which completes the proof. ∎

We next consider the following commutator

$\displaystyle\big{[}\sum_{i=1}^{n}z_{i}^{\alpha\beta}\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta},\,\sum_{j=1}^{n}w_{j}^{\gamma\delta*}\sigma_{j}^{\gamma}\sigma_{j+1}^{\delta}\big{]}=\sum_{i=1}^{n}\sum_{j=1}^{n}z_{j\alpha\beta}w_{j}^{\gamma\delta*}[\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta},\,\sigma_{j}^{\gamma}\sigma_{j+1}^{\delta}]$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{n}z_{i}^{\alpha\beta}w_{i}^{\gamma\delta*}[\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta},\,\sigma_{i}^{\gamma}\sigma_{i+1}^{\delta}]+\sum_{i=1}^{n}z_{i}^{\alpha\beta}w_{i+1}^{\gamma\delta*}\sigma_{i}^{\alpha}[\sigma_{i+1}^{\beta},\,\sigma_{i+1}^{\gamma}]\sigma_{i+2}^{\delta}$
$\displaystyle+$
$\displaystyle\sum_{i=1}^{n}z_{i+1}^{\alpha\beta}w_{i}^{\gamma\delta*}\sigma_{i}^{\gamma}[\sigma_{i+1}^{\alpha},\,\,\sigma_{i+1}^{\delta}]\sigma_{i+2}^{\beta},$
(S46)

where $z_{i}^{\alpha\beta}$ and $w_{j\gamma\delta}$ are complex numbers and we
have used

$\sum_{i=1}^{n}z_{i}^{\alpha\beta}w_{i-1}^{\gamma\delta}[\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta},\,\sigma_{i-1}^{\gamma}\sigma_{i}^{\delta}]=\sum_{i=1}^{n}z_{i+1}^{\alpha\beta}w_{i}^{\gamma\delta}\sigma_{i}^{\gamma}[\sigma_{i+1}^{\alpha},\,\,\sigma_{i+1}^{\delta}]\sigma_{i+2}^{\beta},$
(S47)

thanks to the periodic boundary condition. Using Eq. (S49), it is
straightforward to obtain

$\displaystyle\big{[}\sum_{i=1}^{n}z_{i}^{\alpha\beta}\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta},\,\sum_{j=1}^{n}z_{j}^{\alpha\beta*}\sigma_{j}^{\alpha}\sigma_{j+1}^{\beta}\big{]}$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{n}z_{i}^{\alpha\beta}z_{i+1}^{\alpha\beta*}\sigma_{i}^{\alpha}[\sigma_{i+1}^{\beta},\,\sigma_{i+1}^{\alpha}]\sigma_{i+2}^{\beta}+\sum_{i=1}^{n}z_{i+1}^{\alpha\beta}z_{i}^{\alpha\beta*}\sigma_{i}^{\alpha}[\sigma_{i+1}^{\alpha},\,\sigma_{i+1}^{\beta}]\sigma_{i+2}^{\beta}$
$\displaystyle=$
$\displaystyle-2\text{i}\epsilon_{\alpha\beta\gamma}\sum_{i=1}^{n}(z_{i}^{\alpha\beta}z_{i+1}^{\alpha\beta*}-z_{i}^{\alpha\beta*}z_{i+1}^{\alpha\beta})\sigma_{i}^{\alpha}\sigma_{i+1}^{\gamma}\sigma_{i+2}^{\beta}$
$\displaystyle=$ $\displaystyle
4\epsilon_{\alpha\beta\gamma}\sum_{i=1}^{n}\text{Im}(z_{i}^{\alpha\beta}z_{i+1}^{\alpha\beta*})\sigma_{i}^{\alpha}\sigma_{i+1}^{\gamma}\sigma_{i+2}^{\beta}.$
(S48)

With Eqs. (S46, S48), one finds

$\displaystyle[\sum_{i=1}^{n}z_{i}^{\alpha\beta}\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta}+\sum_{j=1}^{n}w_{j}^{\gamma\delta}\sigma_{j}^{\gamma}\sigma_{j+1}^{\delta},\,\sum_{i=1}^{n}z_{i}^{\alpha\beta*}\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta}+\sum_{j=1}^{n}w_{j}^{\gamma\delta*}\sigma_{j}^{\gamma}\sigma_{j+1}^{\delta}]$
$\displaystyle=$ $\displaystyle
4\epsilon_{\alpha\beta\mu}\sum_{i=1}^{n}\text{Im}(z_{i}^{\alpha\beta}z_{i+1}^{\alpha\beta*})\sigma_{i}^{\alpha}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\beta}+4\epsilon_{\gamma\delta\lambda}\sum_{i=1}^{n}\text{Im}(w_{i}^{\alpha\beta}w_{i+1}^{\alpha\beta*})\sigma_{i}^{\gamma}\sigma_{i+1}^{\lambda}\sigma_{i+2}^{\delta}$
$\displaystyle-$ $\displaystyle
4\epsilon_{\beta\gamma\mu}\sum_{i=1}^{n}\text{Im}(z_{i}^{\alpha\beta}w_{i+1}^{\gamma\delta*})\sigma_{i}^{\alpha}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\delta}-4\epsilon_{\alpha\delta\lambda}\sum_{i=1}^{n}\text{Im}(z_{i+1}^{\alpha\beta}w_{i}^{\gamma\delta*})\sigma_{i}^{\gamma}\sigma_{i+1}^{\lambda}\sigma_{i+2}^{\beta}+\mathcal{X}_{\leq
2}^{\alpha\beta\gamma\delta},$ (S49)

where we have used the fact that
$[A+B,\,A^{\dagger}+B^{\dagger}]=[A,\,A^{\dagger}]+[B,\,B^{\dagger}]+([A,\,B^{\dagger}]+\text{h.c.})$
and that

$\mathcal{X}_{\leq 2}^{\alpha\beta\gamma\delta}\equiv
2\text{i}\sum_{i=1}^{n}\text{Im}(z_{i}^{\alpha\beta}w_{i}^{\gamma\delta*})[\sigma_{i}^{\alpha}\sigma_{i+1}^{\beta},\,\sigma_{i}^{\gamma}\sigma_{i+1}^{\delta}]$
(S50)

is at most two-body.

To cancel the three-body interaction
$\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x}$ in Eq. (13) in the main text,
we would like the second term of Eq. (S46) to be the three body operator. One
may take $\alpha=x$, $\beta=y$, $\gamma=z$ and $\delta=x$ in Eq. (S46). This
results in

$[\sum_{i=1}^{n}z_{i}^{xy}\sigma_{i}^{x}\sigma_{i+1}^{y},\,\sum_{j=1}^{n}w_{j}^{zx*}\sigma_{j}^{z}\sigma_{j+1}^{x}]=2\text{i}\sum_{i=1}^{n}z_{i}^{xy}w_{i+1}^{zx*}\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x},$
(S51)

where have used the fact that

$\displaystyle[\sigma_{i}^{x}\sigma_{i+1}^{y},\,\sigma_{i}^{z}\sigma_{i+1}^{x}]$
$\displaystyle=0,$ (S52)
$\displaystyle\sigma_{i}^{z}[\sigma_{i+1}^{x},\,\sigma_{i+2}^{x}]\sigma_{i+2}^{y}$
$\displaystyle=0,$ (S53)

according to Eq. (S45). Thus, if we take

$H_{\text{c},\,l}=\sum_{i=1}^{n}\left(c_{li}^{xy}\sigma_{i}^{x}\sigma_{i+1}^{y}+c_{li}^{zx}\sigma_{i}^{z}\sigma_{i+1}^{x}\right)$
(S54)

discussed in the main text, according to Eq. (S49), one can easily find

$\displaystyle[H_{\text{c},\,l},\,H_{\text{c},\,-l}]$
$\displaystyle=4\sum_{i=1}^{n}\text{Im}(c_{li}^{xy}c_{li+1}^{xy*})\sigma_{i}^{x}\sigma_{i+1}^{z}\sigma_{i+2}^{y}$
$\displaystyle+4\sum_{i=1}^{n}\text{Im}(c_{li}^{zx}c_{li+1}^{zx*})\sigma_{i}^{z}\sigma_{i+1}^{y}\sigma_{i+2}^{x}$
$\displaystyle-4\sum_{i=1}^{n}\text{Im}(c_{li}^{xy}c_{li+1}^{zx*})\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x}.$
(S55)

To cancel the three-body interactions in Eq. (13), we would like the
coefficients in the first and second terms on the r.h.s. of Eq. (S55) to
vanish, while the coefficient in last term on the r.h.s. does not vanish. To
determine the conditions for this to be the case, we first introduce the
following lemma

######  Lemma 2.

Given a set of complex number $\\{z_{i}\\}_{i=1}^{n}$ and periodic boundary
condition $z_{n+1}=z_{1}$, the condition

$\text{Im}(z_{i}z_{i+1}^{*})=0,\,\forall i$ (S56)

is satisfied if and only if

$\frac{\text{Re}(z_{i})}{\text{Im}(z_{i})}=\alpha,$ (S57)

with $\alpha$ independent of the index $i$.

###### Proof.

Assuming Eq. (S57), Eq. (S56) is obvious. The other direction can be proved by
separating $z_{i}$ into real and imaginary part, i.e.,
$z_{i}=u_{i}+\text{i}v_{i}$. Eq. (S56) leads to $u_{i}/v_{i}=u_{i+1}/v_{i+1}$,
$\forall i$. We may thus set $u_{i}/v_{i}=\alpha$, where $\alpha$ is real
constant and independent of the index $i$. This concludes the proof. ∎

Therefore, according Lemma 2, in order to obtain vanishing coefficients in the
first and second terms on the r.h.s. of Eq. (S55), the coefficients should
take the following form

$\displaystyle c_{li}^{xy}$
$\displaystyle=(\alpha_{l}^{xy}+\text{i})v_{li}^{xy},$ (S58) $\displaystyle
c_{li}^{zx}$ $\displaystyle=(\alpha_{l}^{zx}+\text{i})v_{li}^{zx}.$ (S59)

As a result, we obtain

$\text{Im}(c_{li}^{xy}c_{li+1}^{zx*})=(\alpha_{l}^{xy}-\alpha_{l}^{zx})v_{li}^{xy}v_{li+1}^{zx},$
(S60)

and

$[H_{\text{c},\,l},\,H_{\text{c},\,-l}]=4(\alpha_{l}^{zx}-\alpha_{l}^{xy})\sum_{i=1}^{n}v_{li}^{xy}v_{li+1}^{zx}\sigma_{i}^{x}\sigma_{i+1}^{x}\sigma_{i+2}^{x}.$
(S61)

Thus, it is necessary that $\alpha_{l}^{zx}\neq\alpha_{l}^{xy}$ to have non-
vanishing $[H_{\text{c},\,l},\,H_{\text{c},\,-l}]$. The choice in the main
text corresponds to the case where $\alpha_{l}^{xy}=\infty$ and
$\alpha_{l}^{zx}=0$, so that $c_{li}^{xy}$ is real and $c_{li}^{zx}$ is purely
imaginary.

##  VI More general driving protocols

In this section, we generalize the mechanism of the cancellation of the the
$\sigma_{x}$-type interactions discussed in the main text to arbitrary type of
three-body interactions,
$J_{\kappa\mu\lambda}/2\sum_{i=1}^{n}\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda}$.
We employ Eqs. (S46, S48, S49) and Lemma 2 to construct the high frequency
drive $H^{(d)}(t)$ in the main text. Note that when assuming two-body controls
are accessible, any two-body interaction out of the commuting
$H_{\text{c},\,l}$ and $H_{\text{c},\,-l}$ $(l\geq 1)$ can be cancel by adding
control to the static control $H_{\text{c},\,0}$. Therefore, one should not
concern about the last term $\mathcal{X}_{\leq 2}^{\alpha\beta\gamma\delta}$
on the r.h.s. of Eq. (S49) , which is at most two-body interaction according
to Lemma 1. Now let us discuss case by case according to the index degeneracy.

  1. 1.

For three-body interaction
$J_{\kappa\kappa\kappa}/2\sum_{i=1}^{n}\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\kappa}$,
we take $\alpha=\delta=\kappa$, $\beta=\lambda$ and $\gamma=\mu$ in Eq. (S46),
where $\kappa,\,\mu$ and $\lambda$ are distinct and
$[\sigma_{i}^{\lambda},\,\sigma_{i}^{\mu}]=2\text{i}\sigma_{i}^{\kappa}$, so
that the second term on the r.h.s generates
$\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\kappa}$ . Then Eq.
(S46) becomes

$[\sum_{i=1}^{n}z_{i}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda},\,\sum_{j=1}^{n}w_{j}^{\mu\kappa*}\sigma_{j}^{\mu}\sigma_{j+1}^{\kappa}]=2\text{i}\sum_{i=1}^{n}z_{i}^{\kappa\lambda}w_{i+1}^{\mu\kappa*}\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\kappa}$
(S62)

We construct

$H_{\text{c},\,l}=\sum_{i=1}^{n}\left(c_{li}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda}+c_{li}^{\mu\kappa}\sigma_{i}^{\mu}\sigma_{i+1}^{\kappa}\right)$
(S63)

and assume the coefficients $c_{li}^{\kappa\lambda}$ and $c_{li}^{\mu\kappa}$
satisfy Eq. (S57).Then according to Eq. (S49) and Lemma 2, one can easily find

$[H_{\text{c},\,l},\,H_{\text{c},\,-l}]=-4\sum_{i=1}^{n}\text{Im}(c_{li}^{\kappa\lambda}c_{li+1}^{\mu\kappa*})\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\kappa}.$
(S64)

Assuming the coefficients are independent of $i$, the three-body interaction
$J_{\kappa\kappa\kappa}\sum_{i=1}^{n}\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\kappa}$
can be cancelled if

$\frac{J_{\kappa\kappa\kappa}}{2}=\sum_{l=1}^{\infty}\frac{1}{l}\text{Im}(c_{l}^{\kappa\lambda}c_{l}^{\mu\kappa*}).$
(S65)

  2. 2.

For three-body interaction
$J_{\kappa\mu\kappa}/2\sum_{i=1}^{n}\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\kappa}$,
where $\mu\neq\kappa$, we take $\alpha=\delta=\kappa$, $\beta=\lambda$ and
$\gamma=\kappa$ in Eq. (S46), where $\kappa,\,\mu$ and $\lambda$ are distinct
and
$[\sigma_{i}^{\lambda},\,\sigma_{i}^{\kappa}]=2\text{i}\epsilon_{\lambda\kappa\mu}\sigma_{i}^{\mu}$.
Eq. (S46) becomes

$[\sum_{i=1}^{n}z_{i}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda},\,\sum_{j=1}^{n}w_{j}^{\kappa\kappa*}\sigma_{j}^{\kappa}\sigma_{j+1}^{\kappa}]=\mathcal{X}_{2}^{\kappa\lambda\kappa\kappa}+2\text{i}\epsilon_{\lambda\kappa\mu}\sum_{i=1}^{n}z_{i}^{\kappa\lambda}w_{i+1}^{\kappa\kappa*}\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\kappa},$
(S66)

where the two-body interaction $\mathcal{X}_{2}^{\kappa\lambda\kappa\kappa}$
can be always cancelled through $H_{\text{c},\,0}$. One can then construct

$H_{\text{c},\,l}=\sum_{i=1}^{n}\left(c_{li}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda}+c_{li}^{\kappa\kappa}\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\right).$
(S67)

According to Eq. (S49), one can readily find

$[H_{\text{c},\,l},\,H_{\text{c},\,-l}]=4\epsilon_{\lambda\kappa\mu}\sum_{i=1}^{n}\text{Im}(z_{i}^{\kappa\lambda}w_{i+1}^{\kappa\kappa*})\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\kappa}+\mathcal{X}_{\leq
2}^{\alpha\beta\gamma\delta}$ (S68)

with a similar AFM condition (S65).

  3. 3.

For three-body interaction
$J_{\kappa\mu\lambda}/2\sum_{i=1}^{n}\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda}$,
where $\kappa,\,\mu$ and $\lambda$ are distinct, we take
$\alpha=\gamma=\kappa$, $\beta=\delta=\lambda$ in Eq. (S46) and obtain

$\displaystyle[\sum_{i=1}^{n}z_{i}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda},\,\sum_{j=1}^{n}w_{j}^{\mu\kappa*}\sigma_{j}^{\kappa}\sigma_{j+1}^{\lambda}]$
$\displaystyle=2\text{i}\epsilon_{\lambda\kappa\mu}\sum_{i=1}^{n}z_{i}^{\kappa\lambda}w_{i+1}^{\kappa\lambda*}\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda}+2\text{i}\epsilon_{\kappa\lambda\mu}\sum_{i=1}^{n}z_{i+1}^{\kappa\lambda}w_{i}^{\kappa\lambda*}\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda}$
$\displaystyle=-2\text{i}\epsilon_{\kappa\lambda\mu}\sum_{i=1}^{n}[z_{i}^{\kappa\lambda}w_{i+1\kappa\lambda}^{*}-z_{i+1}^{\kappa\lambda}w_{i}^{\kappa\lambda*}]\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda},$
(S69)

which suggests us to use the two-body interaction
$\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda}$ as the only drive. Indeed, if we
take

$H_{\text{c},\,l}=\sum_{i=1}^{n}c_{li}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda},$
(S70)

for $l\geq 1$ then according to Eq. (S48), we find

$[H_{\text{c},\,l},\,H_{\text{c},\,-l}]=4\epsilon_{\kappa\mu\lambda}\sum_{i=1}^{n}\text{Im}(c_{li}^{\kappa\lambda}c_{li+1}^{\kappa\lambda*})\sigma_{i}^{\kappa}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda}.$
(S71)

According to Lemma 2, as long as the ratio
$\text{Re}(c_{li}^{\kappa\lambda})/\text{Im}(c_{li}^{\kappa\lambda})\equiv\alpha_{li}^{\kappa\lambda}$
varies across the chain, the three-body interaction is non-zero. Denoting

$c_{li}^{\kappa\lambda}=(\alpha_{li}^{\kappa\lambda}+\text{i})v_{li}^{\kappa\lambda},$
(S72)

we find

$\text{Im}(c_{li}^{\kappa\lambda}c_{li+1}^{\kappa\lambda*})=(\alpha_{li+1}^{\kappa\lambda}-\alpha_{li}^{\kappa\lambda})v_{li+1}^{\kappa\lambda}v_{li}^{\kappa\lambda}.$
(S73)

A simple example could be
$\alpha_{li+1}^{\kappa\lambda}-\alpha_{li}^{\kappa\lambda}=\Delta\alpha_{l}^{\kappa\lambda}$
and $v_{li}^{\kappa\lambda}=v_{l}^{\kappa\lambda}$ for $i=1,\,\cdots n-1$.
Therefore
$\alpha_{li+1}^{\kappa\lambda}=\alpha_{l1}^{\kappa\lambda}+i\Delta\alpha_{l}^{\kappa\lambda}$
and

$\text{Im}(c_{li}^{\kappa\lambda}c_{li+1}^{\kappa\lambda*})=\Delta\alpha_{l}^{\kappa\lambda}[v_{l}^{\kappa\lambda}]^{2},$
(S74)

for $i=1,\,\cdots n-2$. Now we shall choose $\alpha_{ln}^{\kappa\lambda}$ and
$v_{ln,\,\kappa\lambda}$ such that

$\displaystyle\text{Im}(c_{ln-1,\,\kappa\lambda}c_{ln,\,\kappa\lambda}^{*})$
$\displaystyle=\left[\alpha_{ln,\,\kappa\lambda}-\alpha_{l1}^{\kappa\lambda}-(n-2)\Delta\alpha_{l}^{\kappa\lambda}\right]v_{ln,\,\kappa\lambda}v_{l}^{\kappa\lambda}=\Delta\alpha_{l}^{\kappa\lambda}[v_{l}^{\kappa\lambda}]^{2},$
(S75)
$\displaystyle\text{Im}(c_{ln,\,\kappa\lambda}c_{l1,\,\kappa\lambda}^{*})$
$\displaystyle=\left[\alpha_{l1}^{\kappa\lambda}-\alpha_{ln,\,\kappa\lambda}\right]v_{l}^{\kappa\lambda}v_{ln,\,\kappa\lambda}=\Delta\alpha_{l}^{\kappa\lambda}[v_{l}^{\kappa\lambda}]^{2},$
(S76)

from which one finds

$\displaystyle\alpha_{ln}^{\kappa\lambda}$
$\displaystyle=\alpha_{1}-\left(1-\frac{n}{2}\right)v_{l}^{\kappa\lambda},$
(S77) $\displaystyle v_{ln,\,\kappa\lambda}$
$\displaystyle=\frac{2\Delta\alpha_{l}^{\kappa\lambda}}{2-n}.$ (S78)

So the AFM condition become

$4\epsilon_{\kappa\mu\lambda}\sum_{l=1}^{\infty}\frac{\Delta\alpha_{l}^{\kappa\lambda}[v_{l}^{\kappa\lambda}]^{2}}{l}+\frac{J_{\kappa\mu\lambda}}{2}=0.$
(S79)

  4. 4.

Finally, it turns out, up to the first-order of the high frequency expansion,
there is no simple scheme to cancel the interaction of the type
$\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\lambda}$. However the
combination
$1/2\sum_{i=1}^{n}(J_{\kappa\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}+J_{\mu\mu\lambda}\sigma_{i}^{\mu}\sigma_{i+1}^{\mu})\sigma_{i+2}^{\lambda}$,
where $\kappa,\,\mu$ and $\lambda$ are distinct and
$J_{\kappa\kappa\lambda}=J_{\mu\mu\lambda}$, , can be cancel in the following
simple way. Taking $\alpha=\kappa$, $\beta=\lambda$, $\gamma=\mu$ and
$\delta=\lambda$, one finds

$[\sum_{i=1}^{n}z_{i}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda},\,\sum_{j=1}^{n}w_{j}^{\mu\lambda*}\sigma_{j}^{\mu}\sigma_{j+1}^{\lambda}]=2\text{i}\epsilon_{\lambda\mu\kappa}\sum_{i=1}^{n}(z_{i}^{\kappa\lambda}w_{i+1}^{\mu\lambda*}\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\lambda}+z_{i+1}^{\kappa\lambda}w_{i}^{\mu\lambda*}\sigma_{i}^{\mu}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda})+\mathcal{X}_{2}^{\kappa\lambda\mu\lambda}.$
(S80)

So if we take

$H_{\text{c},\,l}=\sum_{i=1}^{n}\left(c_{li}^{\kappa\lambda}\sigma_{i}^{\kappa}\sigma_{i+1}^{\lambda}+c_{li}^{\mu\lambda}\sigma_{i}^{\mu}\sigma_{i+1}^{\lambda}\right),$
(S81)

and assume the coefficients $c_{li}^{\kappa\lambda}$ and $c_{li}^{\mu\lambda}$
satisfy Eq. (S57), the according to Eq. (S49), one find

$[H_{\text{c},\,l},\,H_{\text{c},-l}]=-4\epsilon_{\lambda\mu\kappa}\sum_{i=1}^{n}\sum_{i=1}^{n}[\text{Im}(c_{li}^{\kappa\lambda}c_{li+1}^{\mu\lambda*})\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}\sigma_{i+2}^{\lambda}+\text{Im}(c_{li+1}^{\kappa\lambda}c_{li}^{\mu\lambda*})\sigma_{i}^{\mu}\sigma_{i+1}^{\mu}\sigma_{i+2}^{\lambda}]+\mathcal{X}_{2}^{\kappa\lambda\mu\lambda}$
(S82)

Assuming $c_{li}^{\kappa\lambda}$ and $c_{li}^{\mu\lambda}$ are homogeneous,
since satisfy Eq. (S57), they can take the form
$c_{li}^{\kappa\lambda}=(\alpha_{l}^{\kappa\lambda}+\text{i})u_{li}^{\kappa\lambda}$
and
$c_{li}^{\mu\lambda}=(\alpha_{l}^{\mu\lambda}+\text{i})v_{li}^{\mu\lambda}$
respectively. Therefore

$\displaystyle\text{Im}(c_{li}^{\kappa\lambda}c_{li+1}^{\mu\lambda*})$
$\displaystyle=(\alpha_{l}^{\mu\lambda}-\alpha_{l}^{\kappa\lambda})u_{li}^{\kappa\lambda}v_{li+1}^{\mu\lambda},$
(S83) $\displaystyle\text{Im}(c_{li+1}^{\kappa\lambda}c_{li}^{\mu\lambda*})$
$\displaystyle=(\alpha_{l}^{\mu\lambda}-\alpha_{l}^{\kappa\lambda})u_{li+1}^{\kappa\lambda}v_{li}^{\mu\lambda}.$
(S84)

Then one can take $u_{i}=\delta_{1,\,i}$ and $v_{i}\neq 0$ for $i\neq n$ and
$v_{n}=0$. Then one can easily show

$\text{Im}(c_{li}^{\kappa\lambda}c_{li+1}^{\mu\lambda*})=(\alpha_{l}^{\mu\lambda}-\alpha_{l}^{\kappa\lambda})u_{li}^{\kappa\lambda}v_{li+1}^{\mu\lambda},$
(S85)
$4\epsilon_{\lambda\mu\kappa}\sum_{l=1}^{\infty}\frac{(\alpha_{l}^{\mu\lambda}-\alpha_{l}^{\kappa\lambda})u_{l}^{\kappa\lambda}v_{l}^{\mu\lambda}}{l}=\frac{J_{\kappa\kappa\lambda}}{2}.$
(S86)

Similar construction also holds for the three-body interaction
$1/2\sum_{i=1}^{n}\sigma_{i+2}^{\lambda}(J_{\lambda\kappa\kappa}\sigma_{i}^{\kappa}\sigma_{i+1}^{\kappa}+J_{\lambda\mu\mu}\sigma_{i}^{\mu}\sigma_{i+1}^{\mu})$
with $J_{\lambda\kappa\kappa}=J_{\lambda\mu\mu}$, which will not be discussed
here.