# Quantum parameter estimation with optimal control

Jing Liu  Department of Mechanical and Automation Engineering, The Chinese
University of Hong Kong, Shatin, Hong Kong  Haidong Yuan
<EMAIL_ADDRESS>Department of Mechanical and Automation Engineering,
The Chinese University of Hong Kong, Shatin, Hong Kong

###### Abstract

A pivotal task in quantum metrology, and quantum parameter estimation in
general, is to design schemes that achieve the highest precision with given
resources. Standard models of quantum metrology usually assume the dynamics is
fixed, the highest precision is achieved by preparing the optimal probe states
and performing optimal measurements. However, in many practical experimental
settings, additional controls are usually available to alter the dynamics.
Here we propose to use optimal control methods for further improvement on the
precision limit of quantum parameter estimation. We show that by exploring the
additional degree of freedom offered by the controls higher precision limit
can be achieved. In particular we show that the precision limit under the
controlled schemes can go beyond the constraints put by the coherent time,
which is in contrast to the standard scheme where the precision limit is
always bounded by the coherent time.

###### pacs:

03.67.-a, 03.65.Yz, 03.65.-w.

##  I Introduction

Quantum metrology, which exploits quantum mechanical effects to achieve high
precision, has gained increased attention in recent years Giovannetti2011 ;
cavesprd ; rosetta ; VBRAU92-1 ; GIOV04 ; Fujiwara2008 ; Paris09 ; Escher2011
; Tsang2013 ; Rafal2012 ; durkin ; Demkowicz1 ; Alipour2014 ;
Interferometer1.1 ; Interferometer2 ; Lu ; YuanPRL ; Chaves ; Correa15 ; Toth
; Smerzi ; Berry ; Berni ; Datta2016 ; Yuan2016 . A typical metrological
procedure is to first encode the interested parameter $x$ on a probe state
$\rho_{0}$ via a parameter dependent dynamics $\mathcal{E}_{x}$, i.e.,
$\rho_{0}\xrightarrow{\mathcal{E}_{x}}\rho_{x}$, then perform a set of
Positive Operator Valued Measurements (POVM) on $\rho_{x}$. Based on the
measurement results an estimation $\hat{x}$ can then be obtained. It is known
that for unbiased estimation quantum Cramér-Rao bound sets a lower bound on
the precision Helstrom ; Holevo ; BRAU94 ; BRAU96  $\delta\hat{x}\geq
1/\sqrt{F},$where $\delta\hat{x}$ is the standard deviation and $F$ is the
quantum Fisher information (QFI). If the procedure is repeated $n$ times, then
$\delta\hat{x}\geq 1/\sqrt{nF}$ where the bound can be achieved in the
asymptotical limit.

In this standard procedure, the dynamics $\mathcal{E}_{x}$ is usually assumed
to be fixed, and the highest precision is achieved by preparing the optimal
probe state and performing the optimal POVM that saturates the quantum Cramér-
Rao bound. The obtained precision is often regarded as the ultimate precision.
However, in many experimental settings, additional controls are usually
available to alter the dynamics for further improvement of the precision
limit, this provides another degree of freedom for optimization.

The parallel scheme and the sequential scheme, as shown in Fig. 1, are two
standard schemes considered in quantum parameter estimation. It is known that
if the dynamics is unitary and the Hamiltonian takes the multiplication form
of the parameter, i.e., if $\mathcal{E}_{x}=e^{-ixH}$, then the two schemes
are equivalent GIOV06 ; while for general unitary dynamics
$\mathcal{E}_{x}=e^{-iH(x)}$, the parallel scheme is equivalent to the
controlled sequential scheme YuanPRL . For noisy quantum parameter estimation,
special controlled schemes, such as quantum error correction and dynamical
decoupling, have been used to improve the precision limit Kessler ; Dur ;
Arrad2014 ; Ozeri2013 ; Sekatski16 ; Tan ; lang ; Zhao ; Zheng ; Dur15 . The
controlled sequential scheme is more implementable on current experimental
settings than the parallel scheme, as high-fidelity controls on small systems
can now be routinely done while preparing large entangled states for the
parallel scheme is still very challenging. The controlled sequential scheme
thus starts to gain attention recently Kessler ; Dur ; Arrad2014 ; Ozeri2013 ;
Sekatski16 ; Adesso16 . Existing controlled schemes that use quantum error
correction or dynamical decoupling either need additional resources such as
auxiliary systems that are completely immune to noises or require the
underlying dynamics possessing certain symmetries, which restrict the scope of
the applications. Systematic methods that can design controls to improve the
precision limit for general dynamics are highly desired in practice.

Figure 1: (Color online) (a) Parallel scheme; (b) Sequential scheme; (c)
Controlled sequential scheme.

In this paper we propose to employ optimal quantum control methods, in
particular the GRadient Ascent Pulse Engineering (GRAPE) Khaneja05 , to design
controls for the improvement of the precision limit in quantum parameter
estimation. Such methods can be used to automatically obtain the optimal
controls for the improvement of the precision limit for general dynamics and
can easily incorporate practical constraints on the controls. It thus provides
a general method to design the controlled schemes in quantum metrology. With
this method we will show that the optimally controlled schemes can obtain
precision limits beyond the coherent time, which is in contrast to the
conventional schemes where the precision limit is always bounded by the
coherent time.

Figure 2: (Color online) The flow chart of the algorithm. The key for this
algorithm is to update the controls based on the gradient value.

##  II Methodology

In this article we consider the system whose dynamics can be described by the
master equation

$\partial_{t}{\rho}(t)=\mathcal{L}[\rho(t)],$ (1)

where $\mathcal{L}$ is a super-operator. For unitary evolution
$\mathcal{L}=-iH^{\times}$ where $H^{\times}(\rho)=\left[H,\rho\right]$; for
noisy evolution $\mathcal{L}=-iH^{\times}+\Gamma$ where $\Gamma$ denotes the
super-operator for the noisy process. The Hamiltonian of a controlled system
can be written as Khaneja05 ; Ernst1987

$H=H_{0}(x)+\sum_{k=1}^{p}V_{k}(t)H_{k},$ (2)

where $H_{0}(x)$ is the free evolution Hamiltonian, $x$ is the interested
parameter, $\sum_{k=1}^{p}V_{k}(t)H_{k}$ are control Hamiltonians with
$V_{k}(t)$ representing the amplitude of $k$th control field. Here we assume
the correlation in the environment decays much faster than the evolution of
the system under Eq. (1), and the Markovian approximation is still valid at
the presence of controls Breuer07 . For example, in Nuclear Magnetic
Resonance, the correlation time of the environment is around $10^{-6}$s and
the coherent time is around $0.1\sim 1$s Ernst1987 , if the time scale of the
control is around $10^{-3}$s, then the Markovian approximation is valid, and
the controls are fast enough to generate the desired operations. We also
assume the controls do not change the noisy operators, this holds under some
physical settings Sklarz2004 ; Tannor1999 ; note  but not in general. The
situations that noisy operators are affected by controls will be addressed in
another work.

To implement the GRAPE we will divide the evolution time $T$ into small time
steps, and within each time step $\Delta t$ the controls will be approximated
as constants. The final state at time $T$ can thus be written as
$\rho(T)=\Pi_{i=1}^{m}\exp(\Delta t\mathcal{L}_{i})\rho(0),$ here $m=T/\Delta
t$ is the number of time steps and $\mathcal{L}_{i}$ is the super-operator for
the $i$th time step. The multiplication in $\rho(T)$ is taken from right to
left.

GRAPE can obtain controls that optimize a given objective function. In this
article we focus on the local precision limit for the measurement of small
shifts around certain known values. Such local precision limit can be
quantified by the QFI, we will thus take the QFI as the objective function.
The QFI is defined as

$F(T)=\mathrm{Tr}\left[\rho(T)L_{\mathrm{s}}^{2}(T)\right],$ (3)

where $L_{\mathrm{s}}(T)$ denotes the symmetric logarithmic derivative (SLD)
which is the solution to the equation
$\partial_{x}\rho(T)=\left[\rho(T)L_{\mathrm{s}}(T)+L_{\mathrm{s}}(T)\rho(T)\right]/2$.
The flow of the algorithm is shown in Fig. 2 (detailed description is in
appendix A). Some steps of the algorithm may require the knowledge of $x$,
which is a-priori unknown, in that case an estimated value $\hat{x}$ will be
used and the controls will be updated adaptively. This, however, does not
affect the precision limit asymptotically.

Figure 3: (Color online) The QFI (normalized by $T^{2}$) as a function of
$\theta$. The dash-dotted green line, purple squares and circled red dots
represent the QFI with and without controls, respectively. The dashed black
line represent the value of QFI for unitary evolution. The target time $T=5$,
and decay rate $\gamma=0.1$. The states in the bracket in the legend represent
the corresponding initial state. The true values of $\omega_{0}$ is assumed to
be $1$.

In practical experiments, the measurements that can be taken are restricted.
It is thus also of practical importance to find the optimal controls that can
lead to the highest precision under a fixed measurement, which is quantified
by the classical Fisher information (CFI) $F_{\mathrm{cl}}$ under the
particular measurement, instead of the QFI. This can also be treated via
GRAPE. Given a set of POVM measurement $\\{E(y)\\}$ with
$\sum_{y}E(y)=\openone$, the probability of getting the measurement result $y$
is given by $p_{y|x}=\mathrm{Tr}(\rho(T)E(y))$, and the CFI is given by

$F_{\mathrm{cl}}(T)=\sum_{y}\frac{(\partial_{x}p_{y|x})^{2}}{p_{y|x}}.$ (4)

Figure 4:  (Color online) Transverse dephasing: (a) The evolution of QFI
(normalized by $T$) with and without controls. The red dots and black lines
represent the QFI with and without controls, respectively (the wiggles in the
black line is not numerical error but caused by some trigonometric functions
in the QFI in this case). The dashed blue line is the analytical solution for
the QFI (and CFI) under controls. The dash-dotted green line is the CFI
without controls. Decay rate $\gamma=0.1$, and the measurement for CFI is
$\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$. (b) Optimal controls obtained
from the GRAPE. (c) The enhanced QFI and CFI (normalized by $T^{2}$) as a
function of $\hat{\omega}_{0}$. The solid red and dash-dotted blue lines
represent the QFI and CFI under controls, respectively. The target time $T=20$
and $\gamma=0.2$. The dashed black line represents the value of QFI without
controls. The true values of $\omega_{0}$ in all panels are assumed to be $1$.

##  III Application

We first apply the algorithm to the phase estimation with a two-level system
under dephasing dynamics. The dynamics is given by Breuer07

$\partial_{t}\rho=-i\left[H,\rho\right]+\frac{\gamma}{2}\left(\sigma_{\vec{n}}\rho\sigma_{\vec{n}}-\rho\right),$
(5)

here the system Hamiltonian is
$H=\frac{1}{2}\omega_{0}\sigma_{3}+\vec{V}(t)\cdot\vec{\sigma}$ with
$\vec{V}(t)=(V_{1}(t),V_{2}(t),V_{3}(t))$,
$\vec{\sigma}=(\sigma_{1},\sigma_{2},\sigma_{3})$. $\sigma_{1}$, $\sigma_{2}$
and $\sigma_{3}$ are Pauli matrices. The dephasing is along
$\sigma_{\vec{n}}=\vec{n}\cdot\vec{\sigma}$ with
$\vec{n}=(\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta)$. Here
$\theta\in[0,\pi]$, $\phi\in[0,2\pi]$. $\omega_{0}$ is the parameter to be
estimated. Here we assume the controls can be performed along all three
directions, however the results hold as long as the controls span
$\mathfrak{su}(2)$.

In Fig. 3 we plotted the QFIs with different dephasing dynamics for $T=5$ (the
unit is taken in the order of $\omega^{-1}_{0}$). The different dephasing
dynamics are characterized by the angle $\theta$ ($\phi$ is taken as zero, as
we can always make a rotation along $\sigma_{3}$ direction to make $\phi$
equal to zero and such rotation does not affect the precision). From the
figure we can see that the highest enhancement, compared to the uncontrolled
schemes, occurs at $\theta=\pi/2$ where the noise is transverse to the
direction of the parameter, and the enhancement reduces when $\theta$ goes to
zero (parallel noises). We note that here no ancillary systems are used, which
is different from previous studies using quantum error correction where
ancillary systems are necessary Kessler ; Dur ; Arrad2014 ; Ozeri2013 ;
Sekatski16 . Besides, in this case, the highest precision does not strongly
dependent on the probe state, which can be seen in Fig. 3. $|0\rangle$ (purple
squares) and $|+\rangle=(|0\rangle+|1\rangle)/\sqrt{2}$ (dash-dotted green
line) provide almost the same precision under the optimal controlled scheme.

We next provide some analysis on the controlled scheme under different noises
to give some physical intuitions on how the controls actually helped improving
the precision limit.

Figure 5: (Color online) Parallel dephasing: (a) The evolution of QFI with and
without controls. The red dots and solid black lines represent the QFI with
and without controls, respectively. The solid yellow line represent the
maximum QFI with single $\pi/2$ pulse along $y$-axis at a proper time. The
green stars and dash-dotted green lines represent the CFI with and without
controls, respectively. The measurement for CFI is
$\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$ (b) QFI as a function of
$t_{0}$. The solid blue and triangular red lines represent the QFI with
rotation at $t_{0}$, for $T=15$ and $T=5$, respectively. The dotted black and
dashed green lines represent the QFIs without rotation.

###  III.1 Transverse dephasing

The improvement of the controlled scheme with transverse noises is shown in
Fig. 4(a). When the probe state is taken as $|+\rangle$, the obtained optimal
controls are $V_{1}(t)=V_{2}(t)=0$ and $V_{3}(t)=-0.5\omega_{0}$, shown in
Fig. 4(b). Such controls essentially keep the probe state at $|+\rangle$,
where it is not affected by the noises. The QFI under such controls is given
by (see appendix for detailed derivation)

$F(T)=\frac{2}{\gamma^{2}}\left(e^{-\gamma T}+\gamma T-1\right).$ (6)

It always increases with $T$ as $F^{\prime}(T)>0$. The precision limit thus is
not constrained by the coherent time. In contrast, without controls there is
an optimal time $T_{\mathrm{opt}}$ (which is determined by the decay rate, for
example when $\omega_{0}\gg\gamma$, $T_{\mathrm{opt}}\simeq 2/\gamma$), at
which the precision reaches the maximum, and beyond $T_{\mathrm{opt}}$ the QFI
starts to decrease with time, which can be seen in Fig. 4(a).

We note in this case the control $V_{3}(t)=-0.5\omega_{0}$ depends on the true
value which is a-priori unknown, in practice an estimated value
$\hat{\omega}_{0}$ need to be used and the controls need to be updated
adaptively according to the estimated value as
$V_{3}(t)=-0.5\hat{\omega}_{0}$. In Fig. 4(c) we plotted the improvement
provided by the controls with different estimation error, it can be seen that
the improvement is quite robust. For example assume the true value
$\omega_{0}=1$ and $T=20$, then as long as $\hat{\omega}_{0}\in[0.8,1.2]$ the
controlled scheme outperforms the uncontrolled scheme, and when
$\hat{\omega}_{0}\in[0.9,1.1]$, the QFI under the controlled scheme is more
than $10$ times larger than the value without controls, thus even with a
$10\%$ estimation error the controlled scheme still provides significant
improvement over the uncontrolled schemes.

If the measurement is fixed, for example the measurement is taken as
$\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$ (here
$|\pm\rangle=(|0\rangle\pm|1\rangle)/\sqrt{2}$), we can also use the optimal
control to improve the precision. From Fig. 4(a) we can see that with the
optimal controls the CFI can actually achieve the maximal QFI, indicating that
the precision limit under the optimal controlled scheme is insensitive to the
measurement performed on the final state, as long as it is projective. This is
because the optimal measurement(which is a projective measurement) can always
be rotated to the fixed measurement which corresponds to a counter rotation on
the probe state that can be achieved via controls. As a comparison the
precision without controls is also plotted, in this case the CFI oscillate
with time and can only reach the QFI for some specific time points, indicating
this measurement scheme is only optimal for some specific time points. From
Fig. 4(c) it can be seen that the precision obtained is also very robust
against the estimation error.

###  III.2 Parallel dephasing

We now provide some analysis for the case with parallel dephasing, which is
usually a more dominant noise for many physical systems NVcenter ; Register ,
and cannot be corrected by quantum error correction techniques even with
ancillary systems Kessler ; Dur ; Arrad2014 ; Ozeri2013 ; Sekatski16 .

The QFI under parallel dephasing are shown in Fig. 5(a). It can be seen that
the QFI under optimal control continues to increase beyond the coherent time,
while in contrast the QFI without control starts to decrease beyond the
coherent time.

To gain some intuition on how controls improved the precision we consider a
simple control strategy: we first prepare the probe state as $|+\rangle$ and
let it evolve under the natural evolution (without controls) for a period of
$t_{0}$, then apply a $\pi/2$-pulse along $y$-direction, after that let the
state evolve for another period of $T-t_{0}$ under the natural evolution. To
analyze the effect of this strategy we write the state with the Bloch
representation as $\rho=(\openone+\vec{r}\cdot\vec{\sigma})/2$, the initial
state $|+\rangle$ thus corresponds to $(r_{1}(0),r_{2}(0),r_{3}(0))=(1,0,0)$.
Under the free evolution the state evolves as

$\displaystyle r_{1}(t)$ $\displaystyle=$ $\displaystyle e^{-\gamma
t}\left[\sin(\omega_{0}t)r_{2}(0)+\cos(\omega_{0}t)r_{1}(0)\right],$ (7)
$\displaystyle r_{2}(t)$ $\displaystyle=$ $\displaystyle e^{-\gamma
t}\left[\cos(\omega_{0}t)r_{2}(0)-\sin(\omega_{0}t)r_{1}(0)\right],$ (8)
$\displaystyle r_{3}(t)$ $\displaystyle=$ $\displaystyle r_{3}(0),$ (9)

which gives $\vec{r}(t)=e^{-\gamma
t}(\cos(\omega_{0}t),-\sin(\omega_{0}t),0).$ If no controls are added, the QFI
for $\omega_{0}$ can be easily computed using the following formula zhong

$F(t)=|\partial_{\omega_{0}}\vec{r}(t)|^{2}+\frac{(\vec{r}(t)\cdot\partial_{\omega_{0}}\vec{r}(t))^{2}}{1-|\vec{r}(t)|^{2}},$
(10)

which gives $F(t)=t^{2}e^{-2\gamma t}$, the maximum is achieved at the
coherent time $T_{\mathrm{opt}}=1/\gamma$.

Figure 6: (Color online) Parallel dephasing: (a) The evolution of normalized
QFI (by $T$) with and without controls. The red dots and solid black line
represent the normalized QFI with and without controls, respectively. The
solid yellow line represent the maximum normalized QFI with single $\pi/2$
pulse along $y$-axis at a proper time. The green stars and dash-dotted green
lines represent the normalized CFI with and without controls, respectively.
The measurement for CFI is $\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$ (b)
The controls obtained from GRAPE for the dynamics with $T=20$. The initial
guessing is randomly. (c) The normalized QFI (by $T^{2}$) for different
$\omega_{0}$. the controls are obtained from the GRAPE for $\omega_{0}=1$. It
can be seen that the QFI of the controlled scheme is higher than the QFI of
the uncontrolled scheme as long as $|\omega_{0}-1|$ is not too big, i.e., as
long as the estimated value is reasonably good. The true values of
$\omega_{0}$ are assumed to be $1$ and decay rate $\gamma=0.1$ in all panels.

Now assume the target time is $T$ and we perform the rotation

$R_{y}=\left(\begin{array}[]{ccc}0&0&-1\\\ 0&-1&0\\\ 1&0&0\end{array}\right)$
(11)

at some time point $t_{0}<T$. The quantum state after the $R_{y}$-rotation is
$e^{-\gamma t_{0}}(0,\sin(\omega_{0}t_{0}),\cos(\omega_{0}t_{0}))$, which,
after another free evolution with a period of $\Delta t=T-t_{0}$, leads to the
final state $\vec{r}(T)=(r_{1}(T),r_{2}(T),r_{3}(T))$ where

$\displaystyle r_{1}(T)$ $\displaystyle=$ $\displaystyle e^{-\gamma
T}\sin(\omega_{0}\Delta t)\sin(\omega_{0}t_{0}),$ (12) $\displaystyle
r_{2}(T)$ $\displaystyle=$ $\displaystyle e^{-\gamma T}\cos(\omega_{0}\Delta
t)\sin(\omega_{0}t_{0}),$ (13) $\displaystyle r_{3}(T)$ $\displaystyle=$
$\displaystyle e^{-\gamma t_{0}}\cos(\omega_{0}t_{0}).$ (14)

The QFI can again be calculated from Eq. (10), with $t_{0}$ as a variable that
can be changed to maximize $F(T)$ (the explicit form of $F(T)$ is in the
appendix).

Figure 5(b) shows the QFI as a function of $t_{0}$, one can see that for
$T=5$, no matter when the rotation is performed, the QFI cannot be higher than
the QFI without rotations. However, for $T=15$, as long as the rotation is
performed at a proper time point, we can obtain an improved QFI. It can also
be seen that the QFI usually has multiple peaks with the variation of $t_{0}$,
the maximum peak may differ for different $T$. The maximum QFI thus may not be
smooth with respect to $T$.

In Fig. 5(a) we plotted the maximum QFI that can be achieved with this simple
control strategy. It can be seen that this strategy does not help improving
the QFI when $T$ is smaller than some time $T^{*}$ (which is approximately the
coherent time $\gamma^{-1}=10$ in this case, same behaviour are found for
other values of $\gamma$), however when $T$ gets big, a control pulse at a
proper time $t_{0}$ improves the QFI. The intuition of this simple strategy is
that although states in the $x$-$y$ plane have a fast rate of parametrization
under the Hamiltonian $\omega_{0}\sigma_{3}$, they are also affected most by
the parallel noise, when $T$ gets large the effect of noise overrides the
parametrization, applying pulses at proper time that rotate the states away
from $x$-$y$ plane help mitigate the noise effect thus improve the precision.
More rotations can further improve the precision and GRAPE essentially
provides a systematical way to find these rotations.

This is contrary to the conventional belief that coherent time sets the limit
on the achievable precision, which is particular useful for those systems
where the preparation of the probe states and the measurements are costly and
one would like to extract more information for each measurement. Note that for
the local precision limit which measures small shifts around certain known
value, the phase can still be distinguished even under a long evolution time.
For completely unknown phase, one needs to first evolve for a short time to
avoid the possible ambiguity as the phase may wrap around the $2\pi$ interval.
However after a rough estimation, the evolution time can get longer.

If the cost for the preparation and measurement is negligible, we should
compare the QFI per unit of time, which is called the normalized QFI. As shown
in Fig 6(a), with controls the maximum value of the normalized QFI is not
improved compared to the values without controls(dashed black line), which
indicates the normalized precision limit is still bounded by the coherent time
under the parallel dephasing. However, when a fixed measurement is considered,
for example the projective measurement
$\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$), the advantage of control shows
up. With optimal controls the CFI is very close to the maximum QFI, indicating
that the measurement $\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$ is almost
optimal under the controlled scheme, while without controls the CFI oscillates
with time and is usually far from the maximum QFI.

Figure 7: (Color online) Spontaneous emission: (a) normalized QFI (by $T$) as
a function of $T$. The red dots and solid black lines represent the QFI with
and without controls, respectively. And the dashed blue and dash-dotted green
lines represent the CFI with and without controls, respectively. Here
$\gamma_{+}=0$, $\gamma_{-}=0.1$. The true value of $\omega_{0}$ is 1. The
measurement for CFI is $\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$. (b) The
controls obtained from GRAPE for the dynamics with $T=20$. The initial
guessing is all zero. (c) The QFI and CFI (normalized by $T^{2}$) as a
function of $\omega_{0}$. The solid blue and dash-dotted green lines represent
the QFI and CFI under the controls given by GRPAE with $\hat{\omega}_{0}=1$,
reprectively. The dashed black line is the QFI without controls. (d) QFI as a
function of $t_{0}$. The solid blue and red lines represent the QFI with
rotation at $t_{0}$, for $T=16$ and $T=8$, respectively. The dashed and dash-
dotted black lines represent the QFIs without rotation.

In Fig. 6(b), the optimal controls obtained from the GRAPE are plotted.
Generally the optimal control is not unique and the appearance of the controls
in Fig. 6(b) is due to the algorithm. Such kind of controls seem complicated,
but have been routinely implemented on physical systems, such as Nuclear
Magnetic Resonance Rancan ; Kocher ; Maximov ; Spindler ; Ryan . Various
techniques have also been developed to smooth the controls Bartels ; Li . And
as shown in Fig. 6(c) the controls obtained are again quite robust against the
estimation error (in this figure we first obtain the controls with
$\hat{\omega}_{0}=1$, then apply the controls to dynamics with different
$\omega_{0}$). It can be seen that the controlled scheme gains over the
uncontrolled scheme with a quite broad range ($\sim 10\%$) of estimation
error.

###  III.3 Spontaneous emission

We give some analysis for the controlled scheme at the presence of the
spontaneous emission, which is another major noise for many practical systems.
We consider the general master equation

$\displaystyle\partial_{t}\rho(t)$ $\displaystyle=$
$\displaystyle-i[H,\rho]+\gamma_{+}\left[\sigma_{+}\rho(t)\sigma_{-}-\frac{1}{2}\left\\{\sigma_{-}\sigma_{+},\rho(t)\right\\}\right]$
(15)
$\displaystyle+\gamma_{-}\left[\sigma_{-}\rho(t)\sigma_{+}-\frac{1}{2}\left\\{\sigma_{+}\sigma_{-},\rho(t)\right\\}\right],$

where $\sigma_{\pm}=(\sigma_{1}\pm i\sigma_{2})/2$ is a ladder operator,
$H=\frac{1}{2}\omega_{0}\sigma_{3}+\vec{V}(t)\cdot\vec{\sigma}$.

The effects of the controls are shown in Figure 7(a). In this case the
normalized QFI under the controlled scheme shows significant improvement over
the value without controls. And similar to the dephasing case, under the
measurement $\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$, the normazlied CFI
(dashed blue line) not only achieves the maximum value, but also ceases to
oscillate.

Again we use a simple control strategy with only one rotation to provide some
intuition on how controls helped improving the precision. For simplicity, we
assume $\gamma_{+}=0$ and $\gamma_{-}=\gamma$. In the Bloch representation,
without the controls the states evolves as

$\displaystyle r_{1}(t)$ $\displaystyle=$ $\displaystyle e^{-\frac{1}{2}\gamma
t}\left[\cos(\omega_{0}t)r_{1}(0)-\sin(\omega_{0}t)r_{2}(0)\right],$ (16)
$\displaystyle r_{2}(t)$ $\displaystyle=$ $\displaystyle e^{-\frac{1}{2}\gamma
t}\left[\cos(\omega_{0}t)r_{2}(0)+\sin(\omega_{0}t)r_{1}(0)\right],$ (17)
$\displaystyle r_{3}(t)$ $\displaystyle=$ $\displaystyle-1+e^{-\gamma
t}+e^{-\gamma t}r_{3}(0).$ (18)

If the initial state is taken as $|+\rangle$, then the QFI at time $t$ is
$F=e^{-\gamma t}t^{2}$ and the maximum is achieved at
$T_{\mathrm{opt}}=2/\gamma$.

Now consider a simple control strategy: we first let the initial state, which
is $|+\rangle$, evolve for some time $t_{0}$ under the free evolution(in this
case the free evolution will drive the state away from the $x-y$ plane of the
Bloch sphere), we then apply a control to rotate the state back to the $x-y$
plane, let it evolves for another period, $T-t_{0}$. The derivation of the QFI
under this simple control strategy is given in appendix F.

Figure 8: (Color online) Schematic for the effects of parallel dephasing and
spontaneous emission on the Bloch sphere. The blue spheres and green spheres
represent the initial and evolved state space in the noise.

In Fig. 7(d), we plotted the QFI under this simple control strategy as a
function of $t_{0}$ for $T=8,16$. Here the control strategy is different from
the parallel dephasing case. Under the parallel dephasing (the effects of the
parallel dephasing on the states are shown in Fig. 8(a)), the states in the
$x-y$ plane are affected most by the dephasing noise although they also
undergo the fastest parametrization. For a long evolution time, the dephasing
noise can override the parametrization, thus applying a control rotating the
state away from $x-y$ plane is beneficial under the parallel dephasing. While
under the spontaneous emission (with the effect on the states shown in Fig.
8(b)), the initial state is in the $x-y$ plane, which has the fastest
parameterizations, but the free evolution quickly drives the states away from
the $x-y$ plane before the noises override the parametrization. It is thus
beneficial to apply a control rotating the state back to the $x-y$ plane for a
fast parametrization. Also for the spontaneous emission the states in the
$x-y$ plane are not the states affected most by the noise.

###  III.4 Energy cost

Figure 9: (Color online) The energy cost $E(t)$ (in the scale of
$\omega_{0}=1$) for performing the controls in the case of transverse
dephasing (solid black line), parallel dephasing (dashed blue line) and
spontaneous emission (dash-dotted red line) with the target time $T=10$. The
dacay rates in all cases are $0.1$.

We provide some estimation on the energy cost for the optimal controls
implemented in some of the examples. In Fig. 9 we plot the energy cost
$E(t)=\sum_{k}\int^{t}_{0}V^{2}_{k}(\tau)d\tau$ of the optimal controls as a
function of time within the total time $T=10$. It can be seen that for the
examples, the implementations of optimal controls do not require much energy
and thus it will not be an obstacle in practice. Note that in the case of the
parallel dephasing, there is a sudden change near the coherent time. This is
due to the fact that under the parallel dephasing controls do not have much
effect far before the coherent time, therefore it needs little controls in
this regime. Strong controls only starts to appear near the coherent time.

##  IV Summary

For many current experimental settings the controlled sequential scheme is
more implementable than the parallel scheme since high-fidelity controls can
now be routinely done on many physical systems, such as nuclear magnetic
resonance Rancan ; Kocher ; Maximov ; Spindler ; Ryan , nitrogen-vacancy
centers Scheuer ; Dolde ; Rong  and cold atoms Bucker . GRAPE provides a
general method to obtain optimal controls for the improvement of the precision
limit, which is expected to find wide applications for many practical quantum
parameter estimation tasks.

As a demonstration we applied the method to the frequency estimation with
different noises and showed that GRAPE can improve the precision limit beyond
the limit set by the coherent time, which is contrary to the conventional
belief that coherent time sets the limit on the achievable precision. This is
particular useful for those systems where the measurements are costly and one
would like to extract more information for each measurement. For the dephasing
cases, we showed that the gain of the controlled scheme is most eminent when
the dephasing noise is orthogonal to the Hamiltonian, while when the noise is
parallel to the Hamiltonian, the controls do not increase the precision one
can obtain per unit of time. Future research includes characterizing the
dynamics and noises for which the controls are (un)useful.

The optimal control method can also be used for non-Markovian dynamics Reich ;
Rebentrost ; Schmidt  and easily incorporate various practical constraints on
pulse shape Dridi ; Machnes ; Koch ; Bartels ; Li , thus provides a versatile
tool for designing controlled schemes for various quantum parameter estimation
tasks.

###### Acknowledgements.

H.Yuan acknowledges partial financial support from RGC of Hong Kong with Grant
No. 538213.

##  Appendix A Algorithm description

GRAPE can obtain controls that optimize a given objective function. In this
article we focus on the local precision limit for the measurement of small
shifts around certain known values. Such local precision limit can be
quantified by the QFI, we will thus take the QFI as the objective function.
The QFI is given by
$F(T)=\mathrm{Tr}\left[\rho(T)L_{\mathrm{s}}^{2}(T)\right]$ where
$L_{\mathrm{s}}(T)$ denotes the symmetric logarithmic derivative (SLD) which
is the solution to the equation
$\partial_{x}\rho(T)=\left[\rho(T)L_{\mathrm{s}}(T)+L_{\mathrm{s}}(T)\rho(T)\right]/2$.
The flow of the algorithm is as following:

  1. 1.

guess initial values of $V_{k}(j)$ (here $V_{k}(j)$ denotes the $k$th control
at the $j$th time step);

  2. 2.

evolve the dynamics and obtain a trajectory of the system;

  3. 3.

calculate the QFI at the target time;

  4. 4.

calculate the gradient $\frac{\delta F(T)}{\delta V_{k}(j)}$;

  5. 5.

update $V_{k}(j)$ to $V_{k}(j)+\epsilon\frac{\delta F(T)}{\delta V_{k}(j)}$;

  6. 6.

restart from step 2 using the updated $V_{k}(j)$ until the QFI converges.

The detailed calculation of $\frac{\delta F(T)}{\delta V_{k}(j)}$ is in
appendix B. The gradient of the QFI can be written as

$\displaystyle\frac{\delta F(T)}{\delta V_{k}(j)}$ $\displaystyle=$
$\displaystyle\Delta
t\mathrm{Tr}\left[L_{\mathrm{s}}^{2}(T)\mathcal{M}_{j}^{(1)}\right]$ (19)
$\displaystyle-2\Delta^{2}t\mathrm{Tr}\left[L_{\mathrm{s}}(T)\\!\left(\mathcal{M}_{j}^{(2)}+\mathcal{M}_{j}^{(3)}\right)\right],$

where $\mathcal{M}^{(1)}_{j}$, $\mathcal{M}^{(2)}_{j}$ and
$\mathcal{M}^{(3)}_{j}$ are Hermitian operators and in the form below

$\displaystyle\mathcal{M}_{j}^{(1)}$ $\displaystyle=$ $\displaystyle
i\mathcal{D}_{j+1}^{m}H_{k}^{\times}(\rho_{j}),$
$\displaystyle\mathcal{M}_{j}^{(2)}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{j}\mathcal{D}_{j+1}^{m}H_{k}^{\times}\mathcal{D}_{i+1}^{j}\dot{H}_{0}^{\times}(\rho_{i}),$
(20) $\displaystyle\mathcal{M}_{j}^{(3)}$ $\displaystyle=$
$\displaystyle(1-\delta_{jm})\\!\\!\\!\sum_{i=j+1}^{m}\\!\\!\mathcal{D}_{i+1}^{m}\dot{H}_{0}^{\times}\mathcal{D}_{j+1}^{i}H_{k}^{\times}(\rho_{j}),$

here $\delta_{jm}$ is Kronecker delta function. $L_{\mathrm{s}}(T)$ is the SLD
of $\rho(T)$, $D_{j+1}^{m}\coloneqq\prod_{i=j+1}^{m}\exp(\Delta
t\mathcal{L}_{i})$ is the propagating superoperator from the $j$th time point
to the target time, $j<m$ (we will let $\mathcal{D}_{i}^{i^{\prime}}=\openone$
when $i>i^{\prime}$). $\rho_{j}=D_{1}^{j}\rho(0)$ is the state at the $j$th
time point. $H^{\times}_{k}=[H_{k},\cdot]$ and
$\dot{H}^{\times}_{0}=[\partial_{x}H_{0},\cdot]$.

In practical experiments, the measurements that can be taken are restricted.
It is thus also of practical importance to find the optimal controls that can
lead to the highest precision under a fixed measurement, which is quantified
by the CFI $F_{\mathrm{cl}}$ under the particular measurement, instead of the
QFI. This can also be treated via GRAPE. Given a set of POVM measurement
$\\{E(y)\\}$ with $\sum_{y}E(y)=\openone$, the probability of getting the
measurement result $y$ is given by $p_{y|x}=\mathrm{Tr}(\rho(T)E(y))$, and the
CFI is given by
$F_{\mathrm{cl}}(T)=\sum_{y}\frac{(\partial_{x}p_{y|x})^{2}}{p_{y|x}}.$ The
gradient of $F_{\mathrm{cl}}$ can be similarly obtained as Eq. (19), which is

$\displaystyle\frac{\delta F_{\mathrm{cl}}(T)}{\delta V_{k}(j)}$
$\displaystyle=$ $\displaystyle\Delta
t\mathrm{Tr}\left(\tilde{L}_{2}\mathcal{M}_{j}^{(1)}\right)$ (21)
$\displaystyle-2\Delta^{2}t\mathrm{Tr}\left[\tilde{L}_{1}\left(\\!\mathcal{M}_{j}^{(2)}\\!+\\!\mathcal{M}_{j}^{(3)}\\!\right)\right],$

where

$\displaystyle\tilde{L}_{1}$ $\displaystyle=$
$\displaystyle\sum_{y}(\partial_{x}\ln p_{y})E(y),$ (22)
$\displaystyle\tilde{L}_{2}$ $\displaystyle=$
$\displaystyle\sum_{y}(\partial_{x}\ln p_{y})^{2}E(y).$ (23)

Here $\mathcal{M}^{(1,2,3)}_{j}$ takes the same form as in Eq. (20) (see
appendix C for detailed derivation).

##  Appendix B Gradient for QFI

The dynamics for the density matrix of a system can be described by the
following general master equation

$\partial_{t}\rho(t)=\mathcal{L}[\rho(t)],$ (24)

where $\mathcal{L}$ is a super-operator and $\mathcal{L}=-iH^{\times}+\Gamma$.
$H^{\times}=[H,\cdot]$ and $\Gamma$ is the super-operator for the noise part.
The Hamiltonian here is

$H=H_{0}(x)+\sum_{k=1}^{p}V_{k}(t)H_{k},$ (25)

where $H_{0}(x)$ is the free evolution Hamiltonian, $x$ is the interested
parameter, $\sum_{k=1}^{p}V_{k}(t)H_{k}$ are control Hamiltonians with
$V_{k}(t)$ representing the amplitude of $k$th control field.

Our objective function is the QFI and the goal is to find the optimal control
to obtain the maximum QFI. The QFI is defined as

$F(T)=\mathrm{Tr}(L^{2}_{\mathrm{s}}(T)\rho(T)),$ (26)

where $L_{\mathrm{s}}(T)$ is the SLD operator at target time $T$ and is
determined by the equation
$2\partial_{x}\rho(T)=\rho(T)L_{\mathrm{s}}(T)+L_{\mathrm{s}}(T)\rho(T)$.

Before utilizing GRAPE to obtain the optimal control, it is necessary to know
the corresponding gradient for QFI on the control coefficients, i.e.,
$\frac{\delta F(T)}{\delta V_{k}(j)}$ based on the general master equation. In
this section we will show the detailed calculation for the general dynamics
given by Eq. (24).

Since
$F(T)=\mathrm{Tr}\left[L_{\mathrm{s}}^{2}(T)D_{j+1}^{m}\rho_{j}\right]=\mathrm{Tr}\left[\lambda_{j}\rho_{j}\right]$,
here $L_{\mathrm{s}}(T)$ is the symmetric logarithm derivative of $\rho(T)$,
$D_{j+1}^{m}\coloneqq\prod_{i=j+1}^{m}\exp(\Delta t\mathcal{L}_{i})$ is the
propagating superoperator from the $j$th time point to the target time,
$j<m$(we will let $\mathcal{D}_{i}^{i^{\prime}}=\openone$ when
$i>i^{\prime}$), $\rho_{j}=D_{1}^{j}\rho(0)$ is the state at the $j$th time
point and $\lambda_{j}=L_{\mathrm{s}}^{2}(T)D_{j+1}^{m}$. The gradient of the
QFI with respect to the controls at the $j$th time step $\frac{\delta
F(T)}{\delta V_{k}(j)}$ can then be computed

$\frac{\delta F(T)}{\delta
V_{k}(j)}=\mathrm{Tr}\left(\frac{\delta\lambda_{j}}{\delta
V_{k}(j)}\rho_{j}\right)+\mathrm{Tr}\left(\lambda_{j}\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right).$ (27)

we calculate both terms in the following.

1) First we calculate $\delta\rho_{j}/\delta V_{k}(j)$. Here the only term
contains $V_{k}(j)$ is the propagator at the $j$th time point, $e^{\Delta
tL_{j}}$. Since $\rho_{j}=e^{\Delta tL_{j}}\rho_{j-1}$, we have

$\frac{\delta\rho_{j}}{\delta V_{k}(j)}=\frac{\delta e^{\Delta
t\mathcal{L}_{j}}}{\delta V_{k}(j)}\rho_{j-1}.$ (28)

It is known that the derivative of an exponential operator is
$\partial_{x}e^{A(x)}=\int_{0}^{1}e^{sA}\left(\partial_{x}A\right)e^{(1-s)A}ds.$
Thus,

$\frac{\delta e^{\Delta t\mathcal{L}_{j}}}{\delta
V_{k}(j)}=\int_{0}^{1}e^{\tau\Delta t\mathcal{L}_{j}}\left(\Delta
t\frac{\delta\mathcal{L}_{j}}{\delta V_{k}(j)}\right)e^{-\tau\Delta
t\mathcal{L}_{j}}d\tau e^{\Delta t\mathcal{L}_{j}}.$ (29)

Since
$\mathcal{L}_{j}(\cdot)=-i[H_{0}+\sum_{k}V_{k}(j)H_{k},\cdot]+\Gamma(\cdot)$,
we have $\frac{\delta\mathcal{L}_{j}}{\delta V_{k}(j)}=-iH_{k}^{\times},$
where $H_{k}^{\times}$ represents the commutation superoperator, i.e.,
$H_{k}^{\times}A=[H_{k},A]$. We thus have

$\frac{\delta e^{\Delta t\mathcal{L}_{j}}}{\delta V_{k}(j)}=-i\Delta
t\int_{0}^{1}e^{\tau\Delta t\mathcal{L}_{j}}H_{k}^{\times}e^{-\tau\Delta
t\mathcal{L}_{j}}d\tau e^{\Delta t\mathcal{L}_{j}},$ (30)

which can be rewritten as

$\frac{\delta e^{\Delta t\mathcal{L}_{j}}}{\delta V_{k}(j)}=-i\Delta
t\int_{0}^{1}e^{\tau\Delta t\mathcal{L}_{j}^{\times}}H_{k}^{\times}d\tau
e^{\Delta t\mathcal{L}_{j}}.$ (31)

Expand it with the Taylor series,

$\displaystyle\frac{\delta e^{\Delta t\mathcal{L}_{j}}}{\delta V_{k}(j)}$
$\displaystyle=$ $\displaystyle-i\Delta
t\int_{0}^{1}\sum_{n=0}^{\infty}\frac{\left(\tau\Delta
t\right)^{n}d\tau}{n!}\left(\mathcal{L}_{j}^{\times}\right)^{n}H_{k}^{\times}e^{\Delta
t\mathcal{L}_{j}}$ (32) $\displaystyle=$
$\displaystyle-i\sum_{n=0}^{\infty}\frac{\left(\Delta
t\right)^{n+1}}{\left(n+1\right)!}\left(\mathcal{L}_{j}^{\times}\right)^{n}H_{k}^{\times}e^{\Delta
t\mathcal{L}_{j}}$ $\displaystyle=$ $\displaystyle-i\Delta
tH_{k}^{\times}e^{\Delta t\mathcal{L}_{j}},$

where the last equation we used the first order approximation. Thus

$\frac{\delta\rho_{j}}{\delta V_{k}(j)}=-i\Delta tH_{k}^{\times}\rho_{j}.$
(33)

2) Next we calculate $\mathrm{Tr}[\rho_{j}\delta\lambda_{j}/\delta V_{k}(j)]$.
We first consider the cases when $j<m$. Since
$\lambda_{j}=L_{\mathrm{s}}^{2}(T)D_{j+1}^{m}$, and $D_{j+1}^{m}$ does not
contain $V_{k}(j)$, thus $\frac{\delta\lambda_{j}}{\delta
V_{k}(j)}=\frac{\delta L_{\mathrm{s}}^{2}(T)}{\delta V_{k}(j)}D_{j+1}^{m},$ we
then have

$\displaystyle\mathrm{Tr}\left(\frac{\delta\lambda_{j}}{\delta
V_{k}(j)}\rho_{j}\right)$ $\displaystyle=$
$\displaystyle\mathrm{Tr}\left(\frac{\delta L_{\mathrm{s}}(T)}{\delta
V_{k}(j)}L_{\mathrm{s}}(T)\rho(T)\right)$ (34)
$\displaystyle+\mathrm{Tr}\left(L_{\mathrm{s}}(T)\frac{\delta
L_{\mathrm{s}}(T)}{\delta V_{k}(j)}\rho(T)\right).$

where we used the fact that $D_{j+1}^{m}\rho_{j}=\rho(T)$.

Now take the functional derivative at both sides of the equation
$\partial_{x}\rho(T)=[\rho(T)L_{\mathrm{s}}(T)+L_{\mathrm{s}}(T)\rho(T)]/2$,
then multiply $L_{\mathrm{s}}(T)$ and take the trace, we get

$\displaystyle\mathrm{Tr}\left[\frac{\delta\left(\partial_{x}\rho(T)\right)}{\delta
V_{k}(j)}L_{\mathrm{s}}(T)\right]$ (35) $\displaystyle=$
$\displaystyle\mathrm{Tr}\left[\frac{\delta\rho(T)}{\delta
V_{k}(j)}L_{\mathrm{s}}^{2}(T)\right]+\frac{1}{2}\mathrm{Tr}\left[\rho(T)\frac{\delta
L_{\mathrm{s}}(T)}{\delta V_{k}(j)}L_{\mathrm{s}}(T)\right]$
$\displaystyle+\frac{1}{2}\mathrm{Tr}\left[\frac{\delta
L_{\mathrm{s}}(T)}{\delta V_{k}(j)}\rho(T)L_{\mathrm{s}}(T)\right].$

Compare with Eq. (34), we then have

$\displaystyle\mathrm{Tr}\left[\frac{\delta\left(\partial_{x}\rho(T)\right)}{\delta
V_{k}(j)}L_{\mathrm{s}}(T)\right]$ (36) $\displaystyle=$
$\displaystyle\mathrm{Tr}\left[\frac{\delta\rho(T)}{\delta
V_{k}(j)}L_{\mathrm{s}}^{2}(T)\right]+\frac{1}{2}\mathrm{Tr}\left(\frac{\delta\lambda_{j}}{\delta
V_{k}(j)}\rho_{j}\right),$

which means

$\displaystyle\mathrm{Tr}\left(\frac{\delta\lambda_{j}}{\delta
V_{k}(j)}\rho_{j}\right)$ (37) $\displaystyle=$ $\displaystyle
2\mathrm{Tr}\\!\\!\left[\frac{\delta\left(\partial_{x}\rho(T)\right)}{\delta
V_{k}(j)}L_{\mathrm{s}}(T)\right]\\!\\!-\\!\\!2\mathrm{Tr}\\!\\!\left[\frac{\delta\rho(T)}{\delta
V_{k}(j)}\\!L_{\mathrm{s}}^{2}(T)\right]\\!\\!.$

Since $\mathrm{Tr}\left[\frac{\delta\rho(T)}{\delta
V_{k}(j)}L_{\mathrm{s}}^{2}(T)\right]=\mathrm{Tr}\left[L_{\mathrm{s}}^{2}(T)D_{j+1}^{m}\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right]=\mathrm{Tr}\left(\lambda_{j}\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right),$ we thus have

$\displaystyle\mathrm{Tr}\left(\frac{\delta\lambda_{j}}{\delta
V_{k}(j)}\rho_{j}\right)$ (38) $\displaystyle=$ $\displaystyle
2\mathrm{Tr}\left[\left(\frac{\delta(\partial_{x}\rho(T))}{\delta
V_{k}(j)}\right)\\!\\!L_{\mathrm{s}}(T)\right]\\!\\!-\\!\\!2\mathrm{Tr}\left(\lambda_{j}\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right)$ $\displaystyle=$ $\displaystyle
2\mathrm{Tr}\left[\partial_{x}\\!\\!\left(\frac{\delta\rho(T)}{\delta
V_{k}(j)}\right)\\!\\!L_{\mathrm{s}}(T)\right]\\!\\!-\\!\\!2\mathrm{Tr}\left(\lambda_{j}\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right)\\!\\!,$

where the last equality we assume the functional derivative and partial
differentiation can be exchanged. Substitute Eq. (38) into Eq. (27), we can
obtain the expression for the gradient, which is

$\frac{\delta F(T)}{\delta
V_{k}(j)}=2\mathrm{Tr}\left[\partial_{x}\left(\frac{\delta\rho(T)}{\delta
V_{k}(j)}\right)L_{\mathrm{s}}(T)\right]+i\Delta
t\mathrm{Tr}\left(\lambda_{j}H_{k}^{\times}\rho_{j}\right).$ (39)

We now derive $\partial_{x}\left(\frac{\delta\rho(T)}{\delta
V_{k}(j)}\right)$. Since $\frac{\delta\rho(T)}{\delta
V_{k}(j)}=\mathcal{D}_{j+1}^{m}\frac{\delta\rho_{j}}{\delta V_{k}(j)},$ we
have

$\displaystyle\partial_{x}\left(\frac{\delta\rho(T)}{\delta V_{k}(j)}\right)$
$\displaystyle=$ $\displaystyle\Delta
t\\!\\!\sum_{i=j+1}^{m}\mathcal{D}_{i+1}^{m}\left(\partial_{x}\mathcal{L}_{i}\right)\mathcal{D}_{j+1}^{i}\frac{\delta\rho_{j}}{\delta
V_{k}(j)}$ (40)
$\displaystyle+\mathcal{D}_{j+1}^{m}\partial_{x}\left(\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right),$

where in the first term we used the fact that in the first order

$\partial_{x}e^{\Delta t\mathcal{L}_{i}}=\Delta
t\left(\partial_{x}\mathcal{L}_{i}\right)e^{\Delta t\mathcal{L}_{i}}$ (41)

From Eq. (33), we then have $\partial_{x}\left(\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right)=\partial_{x}(-i\Delta tH_{k}^{\times}\rho_{j})=-i\Delta
tH_{k}^{\times}\partial_{x}\rho_{j}$. Now as

$\displaystyle\partial_{x}\rho_{j}$ $\displaystyle=$ $\displaystyle\Delta
t\sum_{i=1}^{j}\mathcal{D}_{i+1}^{j}\left(\partial_{x}\mathcal{L}_{i}\right)\mathcal{D}_{1}^{i}\rho_{0}$
(42) $\displaystyle=$ $\displaystyle\Delta
t\sum_{i=1}^{j}\mathcal{D}_{i+1}^{j}\left(\partial_{x}\mathcal{L}_{i}\right)\rho_{i},$

one have $\partial_{x}\left(\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right)=-i\Delta^{2}t\sum_{i=1}^{j}H_{k}^{\times}\mathcal{D}_{i+1}^{j}\left(\partial_{x}\mathcal{L}_{i}\right)\rho_{i}.$
With this expression, we have

$\mathcal{D}_{j+1}^{m}\partial_{x}\left(\frac{\delta\rho_{j}}{\delta
V_{k}(j)}\right)=-i\Delta^{2}t\sum_{i=1}^{j}\mathcal{D}_{j+1}^{m}H_{k}^{\times}\mathcal{D}_{i+1}^{j}\left(\partial_{x}\mathcal{L}_{i}\right)\rho_{i}.$

Thus

$\displaystyle\partial_{x}\left(\frac{\delta\rho(T)}{\delta
V_{k}(j)}\right)\\!\\!$ $\displaystyle=$
$\displaystyle\\!\\!-i\Delta^{2}t\\!\sum_{i=j+1}^{m}\mathcal{D}_{i+1}^{m}\left(\partial_{x}\mathcal{L}_{i}\right)\mathcal{D}_{j+1}^{i}H_{k}^{\times}\rho_{j}$
(43)
$\displaystyle\\!\\!-i\Delta^{2}t\sum_{i=1}^{j}\mathcal{D}_{j+1}^{m}H_{k}^{\times}\mathcal{D}_{i+1}^{j}\\!\\!\left(\partial_{x}\mathcal{L}_{i}\right)\\!\rho_{i}.$

Note that here we cannot discard $\Delta^{2}t$ as the second order, as we have
a summation which can effectively add up to cancel one order of $\Delta t$
(for example $\Delta^{2}t\sum_{i=1}^{m}1=m\Delta^{2}t=T\Delta t$).

Furthermore, as
$\mathcal{L}_{i}(\cdot)=-i\left[H_{0}+\sum_{k}V_{k}(i)H_{k},\cdot\right]+\Gamma(\cdot),$
and only the free Hamiltonian $H_{0}$ contains $x$, we have
$\partial_{x}\mathcal{L}_{i}=-i\left[\partial_{x}H_{0},\cdot\right]=-i(\partial_{x}H_{0})^{\times}$.
Thus, Eq. (43) can be expressed by

$\displaystyle\partial_{x}\left(\frac{\delta\rho(T)}{\delta
V_{k}(j)}\right)\\!\\!$ $\displaystyle=$
$\displaystyle\\!\\!-\Delta^{2}t\\!\\!\sum_{i=j+1}^{m}\mathcal{D}_{i+1}^{m}\left(\partial_{x}H_{0}\right)^{\times}\mathcal{D}_{j+1}^{i}H_{k}^{\times}\rho_{j}$
(44)
$\displaystyle\\!\\!-\Delta^{2}t\sum_{i=1}^{j}\mathcal{D}_{j+1}^{m}H_{k}^{\times}\mathcal{D}_{i+1}^{j}\\!\\!\left(\partial_{x}H_{0}\right)^{\times}\\!\\!\rho_{i}\\!.$

Multiplying $L_{\mathrm{s}}(T)$ on both sides of the equation above and taking
the trace gives

$\displaystyle\mathrm{Tr}\left[\partial_{x}\left(\frac{\delta\rho(T)}{\delta
V_{k}(j)}\right)L_{\mathrm{s}}(T)\right]$ (45) $\displaystyle=$
$\displaystyle-\Delta^{2}t\\!\\!\sum_{i=j+1}^{m}\mathrm{Tr}\left[L_{\mathrm{s}}(T)\mathcal{D}_{i+1}^{m}\left(\partial_{x}H_{0}\right)^{\times}\mathcal{D}_{j+1}^{i}H_{k}^{\times}\rho_{j}\right]$
$\displaystyle-\Delta^{2}t\\!\\!\sum_{i=1}^{j}\mathrm{Tr}\left[L_{\mathrm{s}}(T)\mathcal{D}_{j+1}^{m}H_{k}^{\times}\mathcal{D}_{i+1}^{j}\left(\partial_{x}H_{0}\right)^{\times}\\!\\!\rho_{i}\right]\\!\\!.$

Utilizing above equation, one can obtain the final expression for the
gradient, which is

$\displaystyle\frac{\delta F(T)}{\delta V_{k}(j)}$ $\displaystyle=$
$\displaystyle-2\Delta^{2}t\sum_{i=j+1}^{m}\mathrm{Tr}\left[L_{\mathrm{s}}(T)\mathcal{D}_{i+1}^{m}\left(\partial_{x}H_{0}\right)^{\times}\mathcal{D}_{j+1}^{i}H_{k}^{\times}\rho_{j}\right]$
(46)
$\displaystyle-2\Delta^{2}t\sum_{i=1}^{j}\mathrm{Tr}\left[L_{\mathrm{s}}(T)\mathcal{D}_{j+1}^{m}H_{k}^{\times}\mathcal{D}_{i+1}^{j}\left(\partial_{x}H_{0}\right)^{\times}\rho_{i}\right]$
$\displaystyle+i\Delta
t\mathrm{Tr}\left[L_{\mathrm{s}}^{2}(T)\mathcal{D}_{j+1}^{m}H_{k}^{\times}\rho_{j}\right].$

For the case when $j=m$, the gradient is $\frac{\delta F(T)}{\delta
V_{k}(m)}=2\mathrm{Tr}\left[\left(\partial_{x}\frac{\delta\rho_{m}}{\delta
V_{k}(m)}\right)L_{\mathrm{s}}(T)\right]+i\mathrm{Tr}\left[L_{\mathrm{s}}^{2}(T)\Delta
tH_{k}^{\times}\rho_{m}\right].$ In this case, $\delta\rho_{m}/\delta
V_{k}(m)=-i\Delta tH_{k}^{\times}\rho_{m}$, thus

$\displaystyle\mathrm{Tr}\left[\left(\partial_{x}\frac{\delta\rho_{m}}{\delta
V_{k}(m)}\right)L_{\mathrm{s}}(T)\right]$ (47) $\displaystyle=$
$\displaystyle-\Delta^{2}t\sum_{i=1}^{m}\mathrm{Tr}\left[L_{\mathrm{s}}(T)H_{k}^{\times}D_{i+1}^{m}\left(\partial_{x}H_{0}\right)^{\times}\rho_{i}\right].$

the gradient is then

$\displaystyle\frac{\delta F(T)}{\delta V_{k}(m)}$ $\displaystyle=$
$\displaystyle-2\Delta^{2}t\sum_{i=1}^{m}\mathrm{Tr}\left[L_{\mathrm{s}}(T)H_{k}^{\times}D_{i+1}^{m}\left(\partial_{x}H_{0}\right)^{\times}\rho_{i}\right]$
(48) $\displaystyle+i\Delta
t\mathrm{Tr}\left[L_{\mathrm{s}}^{2}(T)H_{k}^{\times}\rho_{m}\right].$

Combine this equation with Eq. (46), the gradient of the QFI can be written
compactly as the from in the main text.

##  Appendix C Gradient for CFI

It is known that CFI is

$F_{\mathrm{cl}}(T)=\sum_{y}\frac{\left(\partial_{\theta}p_{y}\right)^{2}}{p_{y}},$
(49)

where $p_{y}=\mathrm{Tr}(\rho(x,T)E(y))$. Here $E(y)$ is a POVM measurement
which satisfying $\sum_{y}E(y)=\openone$. To calculate the gradient, we need
to know

$\displaystyle\frac{\delta p_{y}}{\delta V_{k}(j)}$ $\displaystyle=$
$\displaystyle\mathrm{Tr}\left[\frac{\delta\rho(T)}{\delta
V_{k}(j)}E(y)\right]$ (50) $\displaystyle=$
$\displaystyle\mathrm{Tr}\left[D_{j+1}^{m}\frac{\delta\rho_{j}}{\delta
V_{k}(j)}E(y)\right]$ $\displaystyle=$ $\displaystyle-i\Delta
t\mathrm{Tr}\left[E(y)D_{j+1}^{m}H_{k}^{\times}\rho_{j}\right]$
$\displaystyle=$ $\displaystyle-\Delta
t\mathrm{Tr}\left[E(y)\mathcal{M}_{j}^{(1)}\right].$

Then we have

$\displaystyle\frac{\delta\left(\partial_{x}p_{y}\right)}{\delta V_{k}(j)}$
$\displaystyle=$ $\displaystyle-i\Delta
t\mathrm{Tr}\left[E(y)\partial_{x}\left(D_{j+1}^{m}H_{k}^{\times}\rho_{j}\right)\right]$
(51) $\displaystyle=$ $\displaystyle-i\Delta
t\mathrm{Tr}\Big{\\{}E(y)\Big{[}\left(\partial_{x}D_{j+1}^{m}\right)H_{k}^{\times}\rho_{j}$
$\displaystyle+D_{j+1}^{m}H_{k}^{\times}\partial_{x}\rho_{j}\Big{]}\Big{\\}}.$

From previous calculations, we know

$\displaystyle\partial_{x}D_{j+1}^{m}$ $\displaystyle=$ $\displaystyle\Delta
t\sum_{i=j+1}^{m}D_{i+1}^{m}\left(\partial_{x}\mathcal{L}_{i}\right)D_{j+1}^{i},$
(52) $\displaystyle\partial_{x}\rho_{j}$ $\displaystyle=$ $\displaystyle\Delta
t\sum_{i=1}^{j}D_{i+1}^{j}\left(\partial_{x}\mathcal{L}_{i}\right)\rho_{i},$
(53)

then for $j\neq m$,

$\displaystyle\frac{\delta\left(\partial_{x}p_{y}\right)}{\delta V_{k}(j)}$
$\displaystyle=$
$\displaystyle-\Delta^{2}t\mathrm{Tr}\Big{[}\Big{(}E(y)\sum_{i=j+1}^{m}D_{i+1}^{m}\dot{H}_{0}^{\times}D_{j+1}^{i}H_{k}^{\times}\rho_{j}$
(54)
$\displaystyle+\sum_{i=1}^{j}D_{j+1}^{m}H_{k}^{\times}D_{i+1}^{j}\dot{H}_{0}^{\times}\rho_{i}\Big{)}\Big{]}.$

for $j=m$, there is

$\displaystyle\frac{\delta\left(\partial_{x}p_{y}\right)}{\delta V_{k}(m)}$
$\displaystyle=$ $\displaystyle-i\Delta
t\mathrm{Tr}\left[E(y)H_{k}^{\times}\partial_{x}\rho_{m}\right]$
$\displaystyle=$
$\displaystyle-\Delta^{2}t\mathrm{Tr}\left[E(y)H_{k}^{\times}\sum_{i=1}^{m}D_{i+1}^{j}\dot{H}_{0}^{\times}\rho_{i}\right].$

Thus, combined above equations, we have

$\frac{\delta\left(\partial_{x}p_{y}\right)}{\delta
V_{k}(j)}=-\Delta^{2}t\mathrm{Tr}\left\\{E(y)\left[\mathcal{M}_{j}^{(2)}+\mathcal{M}_{j}^{(3)}\right]\right\\}.$

Finally, the gradient is

$\displaystyle\frac{\delta F_{\mathrm{cl}}(T)}{\delta V_{k}(j)}$
$\displaystyle=$ $\displaystyle\sum_{y}\frac{\delta}{\delta
V_{k}(j)}\left(\frac{(\partial_{x}p_{y})^{2}}{p_{\mathrm{T}}(y|x)}\right)$
(55) $\displaystyle=$
$\displaystyle\sum_{y}2\frac{\partial_{x}p_{y}}{p_{y}}\left[\frac{\delta\left(\partial_{x}p_{y}\right)}{\delta
V_{k}(j)}\right]-\left(\frac{\partial_{x}p_{y}}{p_{y}}\right)^{2}\frac{\delta
p_{y}}{\delta V_{k}(j)}$ $\displaystyle=$
$\displaystyle\sum_{y}-2\Delta^{2}t\frac{\partial_{x}p_{y}}{p_{y}}\mathrm{Tr}\left[E(y)\left(\mathcal{M}_{j}^{(2)}+\mathcal{M}_{j}^{(3)}\right)\right]$
$\displaystyle+\Delta
t\left(\frac{\partial_{x}p_{y}}{p_{y}}\right)^{2}\mathrm{Tr}\left(E(y)\mathcal{M}_{j}^{(1)}\right).$

The gradient for CFI is then obtained.

##  Appendix D Analytical solution for transverse dephasing noise

For the dynamics with transverse dephasing noises, the controls obtained from
the GRAPE are shown in Fig. 4(b), which are $V_{x}=0$, $V_{y}=0$ and
$V_{z}=-\omega_{0}/2$. The QFI under such control is the same as the QFI under
free evolution for $\omega_{0}=0$, which can actually be computed
analytically. In the following we give a detailed calculation.

Under the Bloch representation
$\rho=\frac{1}{2}\left(\openone+\vec{r}\cdot\vec{\sigma}\right),$ the initial
state $|+\rangle$ can be expressed as $\vec{r}(0)=(1,0,0)$. From the master
equation we can obtain the differential equations for the Bloch vector as

$\displaystyle\partial_{t}r_{1}(t)$ $\displaystyle=$
$\displaystyle\omega_{0}r_{2}(t),$ (56) $\displaystyle\partial_{t}r_{2}(t)$
$\displaystyle=$ $\displaystyle-\gamma r_{2}(t)-\omega_{0}r_{1}(t),$ (57)
$\displaystyle\partial_{t}r_{3}(t)$ $\displaystyle=$ $\displaystyle-\gamma
r_{3}(t).$ (58)

The solution of these equations are

$\displaystyle r_{1}(t)$ $\displaystyle=$ $\displaystyle e^{-\frac{1}{2}\gamma
t}\left[\frac{\gamma}{a}\sinh\left(\frac{1}{2}at\right)+\cosh\left(\frac{1}{2}at\right)\right],$
(59) $\displaystyle r_{2}(t)$ $\displaystyle=$
$\displaystyle-\frac{2\omega_{0}}{a}e^{-\frac{1}{2}\gamma
t}\sinh\left(\frac{1}{2}at\right),$ (60) $\displaystyle r_{3}(t)$
$\displaystyle=$ $\displaystyle 0,$ (61)

where $a=\sqrt{\gamma^{2}-4\omega^{2}_{0}}$. We will now compute
$f(\rho_{\omega_{0}},\rho_{\omega_{0}+\delta\omega_{0}})$ where
$f(\rho_{1},\rho_{2})=\mathrm{Tr}\sqrt{\sqrt{\rho_{1}}\rho_{2}\sqrt{\rho_{1}}}$
denotes the fidelity, the QFI can then be obtained from the second order
expansion of $f(\rho_{\omega_{0}},\rho_{\omega_{0}+\delta\omega_{0}})$ for
$\omega_{0}=0$. It is easy to see that $\vec{r}(t)|_{\omega_{0}=0}=(1,0,0)$,
which is $|+\rangle=\frac{1}{\sqrt{2}}(|0\rangle+|1\rangle)$ and the fidelity
between $|+\rangle$ and an evolved state with a general $\omega_{0}$ is

$f=\sqrt{\langle+|\frac{1}{2}\left(\openone+\vec{r}(t)\cdot\vec{\sigma}\right)|+\rangle}=\sqrt{\frac{1}{2}+\frac{1}{2}r_{1}(t)}.$
(62)

For a small $\delta\omega_{0}$, up to the second order we have
$a=\gamma-\frac{2}{\gamma}\delta^{2}\omega_{0}$ and
$\frac{\gamma}{a}=1+\frac{2}{\gamma^{2}}\delta^{2}\omega_{0}$, then

$\displaystyle r_{1}(t)$ $\displaystyle=$ $\displaystyle e^{-\frac{1}{2}\gamma
t}\left[\frac{\gamma}{a}\sinh\left(\frac{1}{2}at\right)+\cosh\left(\frac{1}{2}at\right)\right]$
(63) $\displaystyle=$ $\displaystyle
e^{-\frac{t}{\gamma}\delta^{2}\omega_{0}}+\frac{2}{\gamma^{2}}\delta^{2}\omega_{0}\frac{e^{-\frac{t}{\gamma}\delta^{2}\omega_{0}}-e^{-\gamma
t}e^{\frac{t}{\gamma}\delta^{2}\omega_{0}}}{2}$ $\displaystyle=$
$\displaystyle 1-\frac{1}{\gamma^{2}}\left(e^{-\gamma t}+\gamma
t-1\right)\delta^{2}\omega_{0}.$

Thus the fidelity is

$\displaystyle f(\rho_{0},\rho_{\delta\omega_{0}})$ $\displaystyle=$
$\displaystyle\sqrt{1-\frac{1}{2}\frac{1}{\gamma^{2}}\left(e^{-\gamma
t}+\gamma t-1\right)\delta^{2}\omega_{0}}$ (64) $\displaystyle=$
$\displaystyle 1-\frac{1}{4}\frac{1}{\gamma^{2}}\left(e^{-\gamma t}+\gamma
t-1\right)\delta^{2}\omega_{0}.$

The QFI can then be obtained from the second order term as

$F(t)=\frac{2}{\gamma^{2}}\left(e^{-\gamma t}+\gamma t-1\right).$ (65)

Now we consider the classical Fisher information under optimal controls.
Taking the measurement as $\\{|+\rangle\langle+|,|-\rangle\langle-|\\}$, then
the probabilities are $p_{+}(t)=(1+r_{1}(t))/2$ and $p_{-}=(1-r_{1}(t))/2$.
The corresponding classical Fisher information is

$F_{\mathrm{class}}(t)=\frac{(\partial_{\omega_{0}}p_{+})^{2}}{p_{+}p_{-}}=\frac{(\partial_{\omega_{0}}r_{1})^{2}}{1-r^{2}_{1}(t)}.$
(66)

For the the controls $\hat{\omega}_{0}\sigma_{3}/2$, where $\hat{\omega}_{0}$
is very close to $\omega_{0}$, based on Eq. (63), one can see
$1-r^{2}_{1}(t)=F(t)\delta^{2}\omega_{0}$ and
$\partial_{\omega_{0}}r_{1}(t)=-F(t)\delta\omega_{0}$, which indicates
$F_{\mathrm{class}}(t)=F(t)$, i.e., the measurement
$\\{|+\rangle\langle+\rangle,|-\rangle\langle-|\\}$ is the optimal measurement
to access the quantum Fisher information.

##  Appendix E Parallel dephasing

Here we consider the simple control strategy for the dynamics with parallel
dephasing noises. Recall that the strategy is to first prepare the probe state
at $|+\rangle$ and let it evolve under the natural evolution (without
controls) for a period of $t_{0}$, then apply a $\pi/2$-pulse along
$y$-direction and let it evolve for another period of $T-t_{0}$. As shown in
the main text the final state at $T$ in the Bloch representation is given by
$\vec{r}(T)=(r_{1}(T),r_{2}(T),r_{3}(T))$ with

$\displaystyle r_{1}(T)$ $\displaystyle=$ $\displaystyle e^{-\gamma
T}\sin(\omega_{0}\Delta t)\sin(\omega_{0}t_{0}),$ (67) $\displaystyle
r_{2}(T)$ $\displaystyle=$ $\displaystyle e^{-\gamma T}\cos(\omega_{0}\Delta
t)\sin(\omega_{0}t_{0}),$ (68) $\displaystyle r_{3}(T)$ $\displaystyle=$
$\displaystyle e^{-\gamma t_{0}}\cos(\omega_{0}t_{0}),$ (69)

from which we can obtain the QFI using the following formula zhong

$F(T)=|\partial_{\omega_{0}}\vec{r}(T)|^{2}+\frac{(\vec{r}(T)\cdot\partial_{\omega_{0}}\vec{r}(T))^{2}}{1-|\vec{r}(T)|^{2}},$
(70)

specifically,

$\displaystyle F(T)$ $\displaystyle=$ $\displaystyle e^{-2\gamma
t_{0}}t_{0}^{2}\sin^{2}(\omega_{0}t_{0})$ $\displaystyle+e^{-2\gamma
T}[t_{0}^{2}+T(T-2t_{0})\sin^{2}(\omega_{0}t_{0})]$
$\displaystyle+\frac{t_{0}^{2}(e^{-2\gamma T}-e^{-2\gamma
t_{0}})^{2}\sin^{2}(\omega_{0}t_{0})\cos^{2}(\omega_{t}t_{0})}{1-e^{-2\gamma
T}\sin^{2}(\omega_{0}t_{0})-e^{-2\gamma t_{0}}\cos^{2}(\omega_{0}t_{0})}.$

##  Appendix F Spontaneous emission

For the non-controlled scheme, recall the Hamiltonian of this example is
$H=\frac{1}{2}\omega_{0}\sigma_{3}$ and the master equation for the
spontaneous emission is

$\displaystyle\partial_{t}\rho(t)$ $\displaystyle=$
$\displaystyle-i[H,\rho]+\gamma_{+}\left[\sigma_{+}\rho(t)\sigma_{-}-\frac{1}{2}\left\\{\sigma_{-}\sigma_{+},\rho(t)\right\\}\right]$
(71)
$\displaystyle+\gamma_{-}\left[\sigma_{-}\rho(t)\sigma_{+}-\frac{1}{2}\left\\{\sigma_{+}\sigma_{-},\rho(t)\right\\}\right],$

The solution for the master equation is

$\displaystyle\rho_{00}(t)$ $\displaystyle=$ $\displaystyle
e^{-(\gamma_{+}+\gamma_{-})t}\rho_{00}(0)\\!+\\!\frac{\gamma_{+}}{\gamma_{+}+\gamma-}\\!\left[1-e^{-(\gamma_{+}+\gamma_{-})t}\right],$
$\displaystyle\rho_{01}(t)$ $\displaystyle=$ $\displaystyle
e^{-i\omega_{0}t-\frac{1}{2}(\gamma_{+}+\gamma_{-})t}\rho_{01}(0).$

In the Bloch representation, we have

$\displaystyle r_{1}(t)$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}(\gamma_{+}+\gamma_{-})t}\left[\cos(\omega_{0}t)r_{1}(0)-\sin(\omega_{0}t)r_{2}(0)\right],$
$\displaystyle r_{2}(t)$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}(\gamma_{+}+\gamma_{-})t}\left[\cos(\omega_{0}t)r_{2}(0)+\sin(\omega_{0}t)r_{1}(0)\right],$
$\displaystyle r_{3}(t)$ $\displaystyle=$
$\displaystyle\frac{\gamma_{+}-\gamma_{-}}{\gamma_{+}+\gamma_{-}}\left[1-e^{-(\gamma_{+}+\gamma_{-})t}\right]+e^{-(\gamma_{+}+\gamma_{-})t}r_{3}(0).$

For the initial state $|+\rangle$, the evolved Bloch vector reads

$\displaystyle r_{1}(t)$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}(\gamma_{+}+\gamma_{-})t}\cos(\omega_{0}t),$ (72)
$\displaystyle r_{2}(t)$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}(\gamma_{+}+\gamma_{-})t}\sin(\omega_{0}t),$ (73)
$\displaystyle r_{3}(t)$ $\displaystyle=$
$\displaystyle\frac{\gamma_{+}-\gamma_{-}}{\gamma_{+}+\gamma_{-}}\left[1-e^{-(\gamma_{+}+\gamma_{-})t}\right].$
(74)

From these expressions, we have
$|\partial_{\omega_{0}}\vec{r}(T)|^{2}=e^{-(\gamma_{+}+\gamma_{-})T}T^{2}$ and
$\vec{r}(T)\cdot\partial_{\omega_{0}}\vec{r}(T)=0$, thus, the QFI at target
time for non-controlled scheme is

$F=e^{-(\gamma_{+}+\gamma_{-})T}T^{2}.$ (75)

Next we perform a single rotation strategy as an intuitive mechanism for the
effec of control: the Bloch vector is rotated by the control to $x-y$ plane
along the $y$ axis at time $t_{0}$. Before the rotation, the Bloch vector is

$\displaystyle r_{1}(t_{0})$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}(\gamma_{+}+\gamma_{-})t_{0}}\cos(\omega_{0}t_{0}),$ (76)
$\displaystyle r_{2}(t_{0})$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}(\gamma_{+}+\gamma_{-})t_{0}}\sin(\omega_{0}t_{0}),$ (77)
$\displaystyle r_{3}(t_{0})$ $\displaystyle=$
$\displaystyle\frac{\gamma_{+}-\gamma_{-}}{\gamma_{+}+\gamma_{-}}\left[1-e^{-(\gamma_{+}+\gamma_{-})t_{0}}\right].$
(78)

In the following we assume $\gamma_{+}=0$ and $\gamma_{-}=\gamma$, above
expressions then reduce to

$\displaystyle r_{1}(t_{0})$ $\displaystyle=$ $\displaystyle e^{-\gamma
t_{0}}\cos(\omega_{0}t_{0}),$ (79) $\displaystyle r_{2}(t_{0})$
$\displaystyle=$ $\displaystyle e^{-\gamma t_{0}}\sin(\omega_{0}t_{0}),$ (80)
$\displaystyle r_{3}(t_{0})$ $\displaystyle=$ $\displaystyle-1+e^{-\gamma
t_{0}}.$ (81)

Now we perform the rotation $R_{y,1}$ in the form below

$\left(\begin{array}[]{ccc}\frac{r_{1}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}&0&\frac{r_{3}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}\\\
0&1&0\\\
-\frac{r_{3}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}&0&\frac{r_{1}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}\end{array}\right),$

where $\bar{\omega}_{0}$ is the true value of $\omega_{0}$. After the
rotation, we have

$R_{y,1}\vec{r}(t_{0})=\left(\begin{array}[]{c}\frac{r_{1}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{1}(t_{0})+r_{3}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{3}(t_{0})}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}\\\
r_{2}(t_{0})\\\
\frac{-r_{3}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{1}(t_{0})+r_{1}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{3}(t_{0})}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}\end{array}\right).$
(82)

Then the Bloch vector at target time $T$ reads

$\displaystyle r_{1}(T)$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}\gamma(T-t_{0})}\left\\{\cos\left[\omega_{0}(T-t_{0})\right]\frac{r_{1}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{1}(t_{0})+r_{3}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{3}(t_{0})}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}-\sin\left[\omega_{0}(T-t_{0})\right]r_{2}(t_{0})\right\\},$
(83) $\displaystyle r_{2}(T)$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}\gamma(T-t_{0})}\left\\{\cos\left[\omega_{0}(T-t_{0})\right]r_{2}(t_{0})+\sin\left[\omega_{0}(T-t_{0})\right]\frac{r_{1}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{1}(t_{0})+r_{3}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{3}(t_{0})}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}\right\\},$
(84) $\displaystyle r_{3}(T)$ $\displaystyle=$
$\displaystyle-1+e^{-\gamma(T-t_{0})}+e^{-\gamma(T-t_{0})}\frac{-r_{3}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{1}(t_{0})+r_{1}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}r_{3}(t_{0})}{\sqrt{r_{1}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}+r_{3}^{2}(t_{0})|_{\omega_{0}=\bar{\omega}_{0}}}}.$
(85)

At the point $\omega_{0}=\bar{\omega}_{0}$, the Bloch vector is

$\displaystyle r_{1}(T)|_{\omega_{0}=\bar{\omega}_{0}}$ $\displaystyle=$
$\displaystyle e^{-\frac{1}{2}\gamma
T}\left\\{\cos\left[\bar{\omega}_{0}(T-t_{0})\right]\sqrt{\cos^{2}(\bar{\omega}_{0}t_{0})+e^{\gamma
t_{0}}-2+e^{-\gamma
t_{0}}}-\sin\left[\bar{\omega}_{0}(T-t_{0})\right]\sin(\bar{\omega}_{0}t_{0})\right\\},$
(86) $\displaystyle r_{2}(T)|_{\omega_{0}=\bar{\omega}_{0}}$ $\displaystyle=$
$\displaystyle e^{-\frac{1}{2}\gamma
T}\left\\{\cos\left[\bar{\omega}_{0}(T-t_{0})\right]\sin(\omega_{0}t_{0})+\sin\left[\bar{\omega}_{0}(T-t_{0})\right]\sqrt{\cos^{2}(\bar{\omega}_{0}t_{0})+e^{\gamma
t_{0}}-2+e^{-\gamma t_{0}}}\right\\},$ (87) $\displaystyle
r_{3}(T)|_{\omega_{0}=\bar{\omega}_{0}}$ $\displaystyle=$
$\displaystyle-1+e^{-\gamma(T-t_{0})}.$ (88)

The derivative of Bloch vector at $\omega_{0}=\bar{\omega}_{0}$ is

$\displaystyle[\partial_{\omega_{0}}r_{1}(T)]|_{\omega_{0}=\bar{\omega}_{0}}$
$\displaystyle=$ $\displaystyle e^{-\frac{1}{2}\gamma
T}\Big{\\{}-(T-t_{0})\sin\left[\bar{\omega}_{0}(T-t_{0})\right]\sqrt{\cos^{2}(\bar{\omega}_{0}t_{0})+e^{\gamma
t_{0}}-2+e^{-\gamma t_{0}}}+t_{0}\sin\left[\bar{\omega}_{0}(2t_{0}-T)\right]$
(89)
$\displaystyle-\cos\left[\bar{\omega}_{0}(T-t_{0})\right]\sin(\bar{\omega}_{0}t_{0})\left(\frac{t_{0}\cos(\bar{\omega}_{0}t_{0})}{\sqrt{\cos^{2}(\bar{\omega}_{0}t_{0})+e^{\gamma
t_{0}}-2+e^{-\gamma t_{0}}}}-T\right)\Big{\\}},$

and

$\displaystyle[\partial_{\omega_{0}}r_{2}(T)]|_{\omega_{0}=\bar{\omega}_{0}}$
$\displaystyle=$ $\displaystyle e^{-\frac{1}{2}\gamma
T}\Big{\\{}(T-t_{0})\cos\left[\bar{\omega}_{0}(T-t_{0})\right]\sqrt{\cos^{2}(\bar{\omega}_{0}t_{0})+e^{\gamma
t_{0}}-2+e^{-\gamma t_{0}}}+t_{0}\cos\left[\bar{\omega}_{0}(2t_{0}-T)\right]$
(90)
$\displaystyle-\sin\left[\bar{\omega}_{0}(T-t_{0})\right]\sin(\bar{\omega}_{0}t_{0})\left(\frac{t_{0}\cos(\bar{\omega}_{0}t_{0})}{\sqrt{\cos^{2}(\bar{\omega}_{0}t_{0})+e^{\gamma
t_{0}}-2+e^{-\gamma t_{0}}}}+T\right)\Big{\\}},$

and

$[\partial_{\omega_{0}}r_{3}(T)]|_{\omega_{0}=\bar{\omega}_{0}}=\frac{\left[e^{-\gamma
T}-e^{-\gamma(T-t_{0})}\right]t_{0}\sin(\bar{\omega}_{0}t_{0})}{\sqrt{\cos^{2}(\bar{\omega}_{0}t_{0})+e^{\gamma
t_{0}}-2+e^{-\gamma t_{0}}}}.$ (91)

The QFI can then be obtained via Eq. (70).

## References

  * (1) V. Giovannetti, S. Lloyd and L. Maccone, Advances in quantum metrology, Nature Photonics 5, 222-229 (2011).
  * (2) V. Giovannetti, S. Lloyd, and L. Maccone, Quantum-enhanced measurements: Beating the standard quantum limit. Science 306, 1330-1336 (2004).
  * (3) C. M. Caves, Quantum-mechanical noise in an interferometer, Phys. Rev. D 23, 1693-1708 (1981).
  * (4) P. M. Anisimov, G. M. Raterman, A. Chiruvelli, W. N. Plick, S. D. Huver, H. Lee and J. P. Dowling, Quantum Metrology with Two-Mode Squeezed Vacuum: Parity Detection Beats the Heisenberg Limit, Phys. Rev. Lett. 104, 103602 (2010).
  * (5) S. L. Braunstein, Quantum limits on precision measurements of phase, Phys. Rev. Lett. 69, 3598-3601 (1992).
  * (6) A. Fujiwara and H. Imai, A fibre bundle over manifolds of quantum channels and its application to quantum statistics, J. Phys. A 41, 255304 (2008).
  * (7) M. G. A. Paris, Quantum estimation for quantum technology, Int. J. Quant. Inf. 7 125 (2009).
  * (8) B. M. Escher, R. L. de Matos Filho and L. Davidovich, General framework for estimating the ultimate precision limit in noisy quantum-enhanced metrology, Nat. Phys. 7, 406-411 (2011).
  * (9) R. Demkowicz-Dobrzanski and L. Maccone, Using Entanglement Against Noise in Quantum Metrology, Phys. Rev. Lett. 113, 250801 (2014).
  * (10) M. Tsang, Quantum metrology with open dynamical systems, New J. Phys. 15, 073005 (2013).
  * (11) R. Demkowicz-Dobrzanski, J. Kolodynski and M. Guta, The elusive Heisenberg limit in quantum-enhanced metrology, Nat. Commun. 3, 1063 (2012).
  * (12) S. Knysh, V. N. Smelyanskiy and G. A. Durkin, Scaling laws for precision in quantum interferometry and the bifurcation landscape of the optimal state, Phys. Rev. A 83, 021804 (2011).
  * (13) S. Alipour, M. Mehboudi and A. T. Rezakhani, Quantum metrology in open systems: dissipative Cramér-Rao bound, Phys. Rev. Lett. 112, 120405 (2014).
  * (14) J. Joo, W. J. Munro and T. P. Spiller, Quantum metrology with entangled coherent states, Phys. Rev. Lett. 107, 083601 (2011).
  * (15) J. Liu, X.-X. Jing and X. Wang, Phase-matching condition for enhancement of phase sensitivity in quantum metrology, Phys. Rev. A 88, 042316 (2013).
  * (16) X.-M. Lu, S. Yu and C. H. Oh, Robust quantum metrological schemes based on protection of quantum Fisher information, Nat. Commun. 6, 7282 (2015).
  * (17) H. Yuan and C.-H. F. Fung, Optimal Feedback Scheme and Universal Time Scaling for Hamiltonian Parameter Estimation, Phys. Rev. Lett. 115, 110401 (2015).
  * (18) R. Chaves, J. B. Brask, M. Markiewicz, J. Kolodynski and A. Acin, Noisy Metrology beyond the Standard Quantum Limit, Phys. Rev. Lett. 111, 120401 (2013).
  * (19) L. A. Correa, M. Mehboudi, G. Adesso and A. Sanpera, Individual Quantum Probes for Optimal Thermometry, Phys. Rev. Lett. 114, 220405 (2015).
  * (20) G. Toth and I. Apellaniz, Quantum metrology from a quantum information science perspective, J. Phys. A 47, 424006 (2014).
  * (21) L. Pezze and Augusto Smerzi, Mach-Zehnder Interferometry at the Heisenberg Limit with Coherent and Squeezed-Vacuum Light, Phys. Rev. Lett. 100, 073601 (2008).
  * (22) D. W. Berry, M. Tsang, M. J. W. Hall, and H. M. Wiseman, Quantum Bell-Ziv-Zakai Bounds and Heisenberg Limits for Waveform Estimation, Phys. Rev. X 5, 031018 (2015).
  * (23) A. A. Berni, T. Gehring, B. M. Nielsen, V. Handchen, M. G. A. Paris and U. L. Andersen, Ab initio quantum-enhanced optical phase estimation using real-time feedback control, Nature Photonics 9, 577-581 (2015).
  * (24) T. Baumgratz and A. Datta, Quantum Enhanced Estimation of a Multidimensional Field, Phys. Rev. Lett. 116, 030801 (2016).
  * (25) H. Yuan, Sequential feedback scheme outperforms the parallel scheme for Hamiltonian parameter estimation, Phys. Rev. Lett. 117, 160801 (2016).
  * (26) C. W. Helstrom, _Quantum Detection and Estimation Theory_ (Academic, New York, 1976).
  * (27) A. S. Holevo, _Probabilistic and Statistical Aspects of Quantum Theory_ (North-Holland, Amsterdam, 1982).
  * (28) S. L. Braunstein and C. M. Caves, Statistical distance and the geometry of quantum states, Phys. Rev. Lett. 72, 3439-3443 (1994).
  * (29) S. L. Braunstein, C. M. Caves and G. J. Milburn, Generalized uncertainty relations: Theory, examples, and Lorentz invariance, Ann. of Phys. 247, 135-173 (1996).
  * (30) V. Giovannetti, S. Lloyd and L. Maccone, Quantum Metrology, Phys. Rev. Lett. 96, 010401 (2006).
  * (31) Q.-S. Tan, Y. Huang, X. Yin, L.-M. Kuang and X. Wang, Enhancement of parameter-estimation precision in noisy systems by dynamical decoupling pulses, Phys. Rev. A 87, 032102 (2013).
  * (32) J. E. Lang, R. B. Liu and T. S. Monteiro, Dynamical-Decoupling-Based Quantum Sensing: Floquet Spectroscopy, Phys. Rev. X 5, 041016 (2015).
  * (33) N. Zhao and Z.-q. Yin, Room-temperature ultrasensitive mass spectrometer via dynamical decoupling, Phys. Rev. A 90, 042118 (2014).
  * (34) Q. Zheng, L. Ge, Y. Yao and Q.-j. Zhi, Enhancing parameter precision of optimal quantum estimation by direct quantum feedback, Phys. Rev. A 91, 033805 (2015).
  * (35) P. Sekatski, M. Skotiniotis and W. Dur, Dynamical decoupling leads to improved scaling in noisy quantum metrology, New J. Phys. 18, 073034 (2016).
  * (36) P. Sekatski, M. Skotiniotis, J. Kolodynski and W. Dur, Quantum metrology with full and fast quantum control, arXiv: 1603.08944.
  * (37) E. M. Kessler, I. Lovchinsky, A. O. Sushkov and M. D. Lukin, Quantum Error Correction for Metrology, Phys. Rev. Lett. 112, 150802 (2014).
  * (38) W. Dur, M. Skotiniotis, F. Frowis and B. Kraus, Improved Quantum Metrology Using Quantum Error Correction, Phys. Rev. Lett. 112, 080801 (2014).
  * (39) G. Arrad, Y. Vinkler, D. Aharonov and A. Retzker, Increasing Sensing Resolution with Error Correction, Phys. Rev. Lett. 112, 150801 (2014).
  * (40) R. Ozeri, Heisenberg limited metrology using Quantum Error-Correction Codes, arXiv: 1310.3432.
  * (41) R. Nichols, T. R. Bromley, L. A. Correa and G. Adesso, Practical quantum metrology in noisy environments, Phys. Rev. A 94, 042101 (2016).
  * (42) N. Khaneja, T. Reiss, C. Hehlet, T. Schulte-Herbruggen and S. J. Glaser, Optimal control of coupled spin dynamics: design of NMR pulse sequences by gradient ascent algorithms, J. Magn. Res. 172, 296-305 (2005).
  * (43) H.-P. Breuer and F. Petruccione, _The Theory of Open Quantum Systems_ (Oxford University Press, Oxford, 2007).
  * (44) R. R. Ernst, G. Bodenhausen and A. Wokaun, _Principles of Nuclear Magnetic Resonance in One and Two Dimensions_ (Clarendon Press, Oxford, 1987).
  * (45) S. E. Sklarz, D. J. Tannor, and N. Khaneja, Optimal control of quantum dissipative dynamics: Analytic solution for cooling the three-level $\Lambda$ system, Phys. Rev. A 69, 053408 (2004).
  * (46) D. J. Tannor and A. Bartana, On the Interplay of Control Fields and Spontaneous Emission in Laser Cooling, J. Phys. Chem. A 103, 10359-10363 (1999).
  * (47) Another way to understand this assumption is to take it as the limit of discrete pulses, where the evolution with the controls within $\delta t$ is $e^{(\mathcal{L}+\mathcal{H}_{c})\delta t}$, here $\mathcal{L}$ is the original operator and $\mathcal{H}_{c}$ comes from the added control, using the Trotter formula this can be thought of as the original evolution interspersed by fast controls, i.e., $e^{(\mathcal{L}+\mathcal{H}_{c})\delta t}\approx(e^{\mathcal{L}\delta t/n}e^{\mathcal{H}_{c}\delta t/n})^{n}$. The continuous controls can thus be taken as the limit of the discrete pulses where the noises are not affected as long as the controls can be made fast enough (fast controls mean that the discrete pulses $e^{\mathcal{H}_{c}\delta t/n}$ can be realized by a stronger control field with much less time so the evolution of $\mathcal{L}$ can be neglected).
  * (48) S. Machnes, U. Sander, S. J. Glaser, P. de Fouquieres, A. Gruslys, S. Schirmer and T. Schulte-Herbrüggen, Comparing, optimizing, and benchmarking quantum-control algorithms in a unifying programming framework, Phys. Rev. A 84, 022305 (2011).
  * (49) L. Childress, M. V. Gurudev Dutt, J. M. Taylor, A. S. Zibrov, F. Jelezko, J. Wrachtrup, P. R. Hemmer and M. D. Lukin, Coherent Dynamics of Coupled Electron and Nuclear Spin Qubits in Diamond, Science 314, 281-285 (2006).
  * (50) M. V. Gurudev Dutt, L. Childress, L. Jiang, E. Togan, J. Maze, F. Jelezko, A. S. Zibrov, P. R. Hemmer and M. D. Lukin, Quantum Register Based on Individual Electronic and Nuclear Spin Qubits in Diamond, Science 316, 1312-1316 (2007).
  * (51) W. Zhong, Z. Sun, J. Ma, X. Wang and F. Nori, Fisher information under decoherence in Bloch representation, Phys. Rev. A 87, 022337 (2013).
  * (52) G. Rancan, T. T. Nguyen and S. J. Glaser, Gradient ascent pulse engineering for rapid exchange saturation transfer, J. Magn. Res. 252, 1-9 (2015).
  * (53) S. S. Kocher, T. Heydenreich, S.J. Glaser, Visualization and analysis of modulated pulses in magnetic resonance by joint time-frequency representations, J. Magn. Res. 249, 63-71 (2014).
  * (54) van I. Maximov, Z. Tosner and Niels Chr. Nielsen, Optimal control design of NMR and dynamic nuclear polarization experiments using monotonically convergent algorithms, J. Chem. Phys. 128, 184505 (2008).
  * (55) P. E. Spindler, Y. Zhang, B. Endeward, N. Gershernzon, T. E. Skinner, S. J. Glaser and T. F. Prisner, Shaped optimal control pulses for increased excitation bandwidth in EPR, J. Magn. Res. 218, 49-58 (2012).
  * (56) C. A. Ryan, C. Negrevergne, M. Laforest, E. Knill and R. Laflamme, Liquid-state nuclear magnetic resonance as a testbed for developing quantum control methods, Phys. Rev. A 78, 012328 (2008).
  * (57) B. Bartels, Smooth Optimal Control of Coherent Quantum Dynamics, Ph.D. Dissertation, School of Soft Matter, Freiburg Institute for Advanced Studies, Albert-Ludwigs-Universitat Freiburg (2015).
  * (58) J.-S. Li, J. Ruthsa, T.-Y. Yub, H. Arthanarib and G. Wagnerb, Optimal pulse design in quantum control: A unified computational method, Proceedings of the National Academy of Sciences 108 (5), 1879-1884 (2011).
  * (59) J. Scheuer, X. Kong, R. S. Said, J. Chen, A. Kurz, L. Marseglia, J. Du, P. R. Hemmer, S. Montangero, T. Calarco, B. Naydenov, and F. Jelezko, Precise qubit control beyond the rotating wave approximation, New J. Phys. 16, 093022 (2014).
  * (60) F. Dolde, V. Bergholm, Y. Wang, I. Jakobi, B. Naydenov, S. Pezzagna, J. Meijer, F. Jelezko, P. Neumann, T. Schulte-Herbruggen, J. Biamonte, and J. Wrachtrup, High-fidelity spin entanglement using optimal control, Nat. Commun. 5, 3371 (2014).
  * (61) X. Rong, J. Geng, Z. Wang, Q. Zhang, C. Ju, F. Shi, C.-K. Duan and J. Du, Implementation of Dynamically Corrected Gates on a Single Electron Spin in Diamond, Phys. Rev. Lett. 112, 050503 (2014)
  * (62) R. Bucker, T. Berrada, S. van Frank, J.-F. Schaff, T. Schumm, J. Schmiedmayer, G. Jager, J. Grond, and U. Hohenester, Vibrational state inversion of a Bose-Einstein condensate: optimal control and state tomography, J. Phys. B 46, 104012 (2013).
  * (63) D. M. Reich, N. Katz and C. P. Koch, Exploiting Non-Markovianity for Quantum Control. Sci. Rep. 5, 12430 (2015).
  * (64) P. Rebentrost, I. Serban, T. Schulte-Herbrüggen and F. K. Wilhelm, Optimal Control of a Qubit Coupled to a Non-Markovian Environment, Phys. Rev. Lett. 102, 090401 (2009).
  * (65) R. Schmidt, A. Negretti, J. Ankerhold, T. Calarco and J. T. Stockburger, Optimal Control of Open Quantum Systems: Cooperative Effects of Driving and Dissipation, Phys. Rev. Lett. 107, 130404 (2011).
  * (66) G. Dridi, M. Lapert, J. Salomon, S. J. Glaser and D. Sugny, Discrete-valued-pulse optimal control algorithms: Application to spin systems, Phys. Rev. A 92, 043417 (2015).
  * (67) S. Machnes, U. Sander, S. J. Glaser, P. de Fouquieres, A. Gruslys, S. Schirmer and T. Schulte-Herbruggen, Comparing, optimizing, and benchmarking quantum-control algorithms in a unifying programming framework, Phys. Rev. A 84, 022305 (2011).
  * (68) C. P. Koch, Controlling open quantum systems: tools, achievements, and limitations, J. Phys.: Condens. Matter 28 213001 (2016).