Exponential stabilization of quantum systems under continuous non-demolition measurements

TL;DR

This paper introduces a feedback control method using Brownian motion to exponentially stabilize a specific eigenstate in quantum systems under continuous non-demolition measurements, ensuring high-fidelity state preparation.

Contribution

It proposes a novel feedback strategy that employs Brownian noise to achieve exponential stabilization of a quantum eigenstate, advancing quantum control techniques.

Findings

Proves exponential convergence to the target eigenstate.

Demonstrates the effectiveness of Brownian-driven control in quantum stabilization.

Suggests simplified filters for practical implementation.

Abstract

We present a novel continuous-time control strategy to exponentially stabilize an eigenstate of a Quantum Non-Demolition (QND) measurement operator. In open-loop, the system converges to a random eigenstate of the measurement operator. The role of the feedback is to prepare a prescribed QND eigenstate with unit probability. To achieve this we introduce the use of Brownian motion to drive the unitary control actions; the feedback loop just adapts the amplitude of this Brownian noise input as a function of the system state. Essentially, it "shakes" the system away from undesired eigenstates by applying strong noise there, while relying on the open-loop dynamics to progressively reach the target. We prove exponential convergence towards the target eigenstate using standard stochastic Lyapunov methods. The feedback scheme and its stability analysis suggest the use of an approximate filter which only tracks the populations of the eigenstates of the measurement operator. Such reduced filters should play an increasing role towards advanced quantum technologies.

Equations

$$d{\rho }_t$$

$$d{Y}_t$$

$${\mathcal{D}}_L(\rho )$$

$${\mathcal{M}}_L(\rho )$$

$${\mathbf{p}}_k({\rho }_t):=\mathrm{Tr}\left({\rho }_t{\Pi }_k\right)\ge 0,$$

$$V_o(\rho )=\sum _{1\le k<k^{\prime }\le d}\sqrt{{\mathbf{p}}_k(\rho )}\sqrt{{\mathbf{p}}_{k^{\prime }}(\rho )},$$

$$\forall t\ge 0,\mathbb{E}[V_o({\rho }_t)]\le \exp (-rt)V_o({\rho }_0)$$

$$d{\mathbf{p}}_k({\rho }_t)=2\sqrt{\eta }\left({\lambda }_k-\sum _{k^{\prime }=1}^d{\lambda }_{k^{\prime }}{\mathbf{p}}_{k^{\prime }}({\rho }_t)\right){\mathbf{p}}_k({\rho }_t)d{W}_t.$$

$$d{\xi }_k=-\frac{1}{2}\eta ({\lambda }_k-\varpi (\xi ){)}^2{\xi }_{k}dt+\sqrt{\eta }({\lambda }_k-\varpi (\xi )){\xi }_{k}dW,$$

$$\begin{gathered} d({\xi }_k{\xi }_{k^{\prime }})=(d{\xi }_k){\xi }_{k^{\prime }}+{\xi }_k(d{\xi }_{k^{\prime }})+(d{\xi }_k)(d{\xi }_{k^{\prime }}) \\ =-\frac{1}{2}\eta ({\lambda }_k-\varpi (\xi ){)}^2{\xi }_k{\xi }_{k^{\prime }}dt+\sqrt{\eta }({\lambda }_k-\varpi (\xi )){\xi }_k{\xi }_{k^{\prime }}dW \\ -\frac{1}{2}\eta ({\lambda }_{k^{\prime }}-\varpi (\xi ){)}^2{\xi }_k{\xi }_{k^{\prime }}dt+\sqrt{\eta }({\lambda }_{k^{\prime }}-\varpi (\xi )){\xi }_k{\xi }_{k^{\prime }}dW \\ +\eta ({\lambda }_k-\varpi (\xi ))({\lambda }_{k^{\prime }}-\varpi (\xi )){\xi }_{k^{\prime }}{\xi }_{k}dt \\ =-\frac{1}{2}\eta ({\lambda }_k-{\lambda }_{k^{\prime }}{)}^2{\xi }_{k^{\prime }}{\xi }_{k}dt+\sqrt{\eta }({\lambda }_k+{\lambda }_{k^{\prime }}-2\varpi (\xi )){\xi }_k{\xi }_{k^{\prime }}dW. \end{gathered}$$

$$\mathcal{A}{V}_o=-\frac{\eta }{2}\sum _{k^{\prime }=1}^d\sum _{k^{\prime }<k}({\lambda }_k-{\lambda }_{k^{\prime }}{)}^2{\xi }_k{\xi }_{k^{\prime }}.$$

$$\mathcal{A}{V}_o\le -\frac{\eta }{2}\left(\underset{k^{\prime },k\ne k^{\prime }}{\min }({\lambda }_k-{\lambda }_{k^{\prime }}{)}^2\right)V_o.$$

$${\rho }_{t+dt}=e^{-iHudt}({\rho }_t+d{\rho }_t)e^{iHudt}$$

$$d\rho=-i\left(f[H~{},\rho]+\sqrt{\eta}\sigma[H~{},~{}L\rho+\rho L]\right)dt\\ +\big{(}\mathcal{D}_{L}(\rho)+\sigma^{2}\mathcal{D}_{H}(\rho)\big{)}dt+\ \left(\sqrt{\eta}\mathcal{M}_{L}(\rho)-i\sigma[H~{},~{}\rho]\right)dW.$$

$$d\rho =-if(\rho )[H,\rho ]dt+{\mathcal{D}}_L(\rho )dt+\sqrt{\eta }{\mathcal{M}}_L(\rho )dW.$$

$$d{v}_t=u_{t}dt=\sigma ({\rho }_t)d{B}_t.$$

$$d\rho =\left({\mathcal{D}}_L(\rho )+\sigma (\rho {)}^2{\mathcal{D}}_H(\rho )\right)dt+\sqrt{\eta }{\mathcal{M}}_L(\rho )dW+\sigma (\rho )i[\rho ,H]dB.$$

$$\sigma (\rho )=\bar{\sigma }\phi \left(\frac{{\max }_{k\ne \ell }{\mathbf{p}}_k(\rho )-p_{\text{min}}}{p_{\text{max}}-p_{\text{min}}}\right)$$

$${\Delta }_{k,k^{\prime }}=\mathrm{Tr}\left({\Pi }_k{\mathcal{D}}_H({\Pi }_{k^{\prime }})\right),$$

$$\begin{gathered} {\Delta }_{k,k}=\mathrm{Tr}\left({\Pi }_{k}H{\Pi }_{k}H\right)-\mathrm{Tr}\left({\Pi }_k{H}^2\right) \\ =\mathrm{Tr}\left({\Pi }_{k}H{\Pi }_{k}H\right)-\mathrm{Tr}\left({\Pi }_{k}H({\Pi }_1+\dots +{\Pi }_d)H\right) \\ =-\sum _{k^{\prime }\ne k}\mathrm{Tr}\left({\Pi }_{k}H{\Pi }_{k^{\prime }}H\right). \end{gathered}$$

$$V_{\alpha }(\rho )=\sum _{s\in \{1,\dots ,d\}\setminus \{\ell \}}\sqrt{\sum _{k\in \{1,\dots ,d\}\setminus \{\ell \}}{\alpha }_{s,k}{\mathbf{p}}_k(\rho )}.$$

$$\sum _{k^{\prime }}{\Delta }_{k,k^{\prime }}{\alpha }_{s,k^{\prime }}=-{\beta }_{s,k}$$

$$c_{\ast }{V}_{\alpha }(\rho )\le \sqrt{1-{\mathbf{p}}_{\ell }(\rho )}\le c^{\ast }{V}_{\alpha }(\rho )$$

$${\alpha }_{\ast }=\underset{s,k\in \{1,\dots ,d\}\setminus \{\ell \}}{\min }{\alpha }_{s,k}\text{ and }{\alpha }^{\ast }=\underset{s,k\in \{1,\dots ,d\}\setminus \{\ell \}}{\max }{\alpha }_{s,k}.$$

$$\mathcal{A}{\mathcal{V}}_{\alpha }(\rho )=\frac{{\sigma }^2(\rho )}{2}{f}_{\alpha }(\rho )-\frac{\eta }{2}{g}_{\alpha }(\rho )-\frac{{\sigma }^2(\rho )}{8}{h}_{\alpha }(\rho )$$

$$f_{\alpha }(\rho )$$

$$g_{\alpha }(\rho )$$

$$h_{\alpha }(\rho )$$

$$d\sqrt{\sum _k{a}_k{\mathbf{p}}_k}=\frac{{\sum }_k{a}_{k}d{\mathbf{p}}_k}{2\sqrt{{\sum }_k{a}_k{\mathbf{p}}_k}}-\frac{{\left({\sum }_k{a}_{k}d{\mathbf{p}}_k\right)}^2}{8{\left({\sum }_k{a}_k{\mathbf{p}}_k\right)}^{3/2}}$$

$$\mathbb{E}[({\Sigma }_k{a}_{k}d{\mathbf{p}}_k{)}^2|\rho ]={\sigma }^2{\left({\Sigma }_k{a}_k\mathrm{Tr}\left(i[{\Pi }_k,H]\rho \right)\right)}^{2}dt+4\eta {\left({\Sigma }_k{a}_k({\lambda }_k-\mathrm{Tr}\left(L\rho \right)){\mathbf{p}}_k\right)}^{2}dt.$$

$${\mathcal{S}}_p\triangleq \left\{\rho \in \mathcal{S}|\exists j\ne \ell ,{\mathbf{p}}_j(\rho )\ge p\right\}.$$

$${\Sigma }_k{\alpha }_{s,k}\mathrm{Tr}\left({\Pi }_k{\mathcal{D}}_H(\rho )\right)=-{\beta }_{s,j}<0.$$

Full text

1 Introduction

The progress of methods for measuring and controlling quantum systems [32, 11] now allows the physics community to implement building blocks of quantum information processors [22, 29, 27] that pose challenging control problems. One of the elementary building blocks is the stabilization of a quantum system onto a target eigenstate of a measurement operator. In a quantum computer [26], such stabilization could be used e.g. for initializing the input states, or for providing auxiliary states that enable particular operations, like entangled states for quantum teleportation [7, 44, 30] and metrology [16] or magic states for performing T-gates [9]. Moreover, quantum error correction protects information by encoding it in a larger state space and rejecting deviations from the nominal code-space [20]. In this sense, strategies for stabilizing a state are a stepping stone towards stabilizing a code subspace and thus protecting information towards quantum information processing. An extension of the present setting towards quantum error correction can be found in [1, 13].

The defining element of our control setting is continuous quantum measurement, where a continuous-time signal provides weak information, associated to weak backaction, on the quantum state [4]. In general, this setting allows to weakly measure non-commuting observables in parallel and to have measurement channels associated to non-Hermitian operators like the one describing energy loss [10]. The most standard measurement for feedback control though is a Quantum non-demolition (QND) measurement [21]. As a close continuous-time counterpart to the projective measurement of a quantum observable, the measurement is characterized by a Hermitian operator; following stochastic dynamics, the system progressively converges towards one of the eigenstates of this measurement operator, in agreement with the stochastic measurement results. All eigenstates of the measurement operator are thus invariant under the dynamics, and as such the QND measurement itself can be considered as preparation tool [19]. Since the system converges to one of the eigenstates at random, an additional feedback control is necessary for preparing a particular target QND eigenstate. We here address how to add this feedback layer in continuous time.

Stabilization of QND eigenstates with continuous-time measurements has been investigated in the literature extensively [33, 39, 25, 38, 24]. Our proposal is meant to improve these results on three aspects. First and foremost, while previous proposals prove convergence in a proper probabilistic sense, their analysis technique does not provide an estimate of the rate of convergence. By using a Lyapunov technique and a novel control approach, we here prove exponential convergence. This is not unexpected since the the QND measurement process in open loop converges exponentially towards the set of QND eigenstates. Second, the feedback laws in existing work a priori depend on the full state $\rho$, which must be estimated in real-time. Our feedback law only depends on eigenstate populations which leads to a reduced filter. Third, existing feedback laws rely on conditional pulses, inspired by optimal discrete-time strategies, yet with the danger that such abrupt signals excite spurious dynamics in the fragile quantum system. Our strategy instead uses continuous control signals, consisting of Brownian noise whose gain is adapted as a function of the estimated state. We call this noise-assisted quantum feedback.

In a nutshell, the strategy that we propose is rather simple to summarize. We let the system evolve in open loop until it gets close to one of the QND eigenstates. If it is the good one, we are done. If it is another one, as we get closer to it, we increase noise input on the system, effectively shaking it away from the bad eigenstate such that it has another chance to converge to the good one. This strategy works with generic conditions on the control Hamiltonian and simple control logic. The detail of the convergence proof is somewhat challenging, but transparent to the end user. The use of a noisy input signal can also be understood as a necessity to induce global exponential convergence on a compact set.

The structure of this paper is as follows: Section 2 presents the dynamical model for QND measurements and introduces the control problem. Section 3 provides a brief overview of current feedback designs. In section 4 we introduce the use of Brownian motion to drive the controls and the standard control design that will be followed. Section 5 presents the main convergence theorem, proving exponential convergence towards a target eigenstate with a closed-loop Lyapunov function. In section 6 we present an approximated filter to estimate the eigenstate populations, inspired from the stability analysis of the previous section. Lastly, in section 7 we make numerical simulations on a spin $J$ system, illustrating the robustness of the control approach.

2 Open-loop QND dynamics and feedback goal

Consider a quantum system of finite dimension $n$. The state space is the set of density matrices $\mathcal{S}=\{\rho\in{\mathbb{C}}^{n\times n}:\rho=\rho^{\dagger},\rho\text{ positive semidefinite},\operatorname{Tr}\left(\rho\right)=1\}$. Here $\operatorname{Tr}\left(\ \cdot\ \right)$ denotes the trace, $A^{\dagger}$ is the complex conjugate transpose of $A$. The open-loop system for a continuous-time quantum measurement with a single measurement channel is governed by the Itō stochastic differential equations [5]:

[TABLE]

with the super-operators

[TABLE]

Here $L\in{\mathbb{C}}^{n\times n}$ is the measurement operator, $W$ is a standard Brownian motion, $Y\in{\mathbb{R}}$ corresponds to the measurement process and $\eta\in[0,1]$ to its efficiency. The second equation describes the stochastic measurement output signal, the first one describes the corresponding measurement backaction.

A QND measurement corresponds to a Hermitian measurement operator $L=L^{\dagger}$. Consider the spectral decomposition $L=\sum_{k=1}^{d}\lambda_{k}\Pi_{k}$ where $\lambda_{1},...\lambda_{d}$ are the distinct ($d\leq n$), real eigenvalues of $L$ with corresponding orthogonal projection operators $\Pi_{1},...,\Pi_{d}$ resolving the identity, i.e. $\sum_{k=1}^{d}\Pi_{k}=\boldsymbol{I}$. The population of the eigenspace $k$ is denoted by

[TABLE]

with the property $\sum_{k=1}^{d}\;\boldsymbol{p}_{k}(\rho_{t})=1$. The following Lemma summarizes the asymptotic behavior of (1) (without mentioning the corresponding measurement signal). It is based on an original exponential Lyapunov function $V_{o}$ providing an estimate of the convergence rate towards the set of stationary states.

Lemma 2.1

Consider the open-loop system (1)with initial condition $\rho_{0}\in\mathcal{S}$. Then any realization of (1) remains in $\mathcal{S}$. Moreover

1. (i)

For any $k$, the subspace population $\boldsymbol{p}_{k}(\rho_{t})$ is a martingale, i.e. $\mathbb{E}[\boldsymbol{p}_{k}(\rho_{t})]=\boldsymbol{p}_{k}(\rho_{0})$. 2. (ii)

If there exists $k\in\{1,2,...,d\}$ such that $\boldsymbol{p}_{k}(\rho_{0})=1$, then $\rho_{0}$ is a steady state of (1). 3. (iii)

The Lyapunov function

[TABLE]

decreases exponentially as

[TABLE]

*with rate *$r=\tfrac{\eta}{2}\;\min_{k,k^{\prime}}(\lambda_{k}-\lambda_{k^{\prime}})^{2}.$

In this sense, the open-loop system (1) converges, for all initial states, towards the set of invariant states described in point $(ii)$.

Proof: The fact that $\mathcal{S}$ is positively invariant for (1) is standard [5].

(i) For each $k$, the subspace populations $\boldsymbol{p}_{k}(\rho)$ follow the Itō SDE:

[TABLE]

Taking the expectation yields $\tfrac{d}{dt}\mathbb{E}[\boldsymbol{p}_{k}(\rho_{t})]=0$, so indeed $\mathbb{E}[\boldsymbol{p}_{k}(\rho_{t})]=\boldsymbol{p}_{k}(\rho_{0})$, $\forall t\geq 0$.

(ii) Take $\rho_{0}$ such that $\boldsymbol{p}_{k}(\rho_{0})=1$. Plugging into (1) we have $\mathcal{D}_{L}(\rho_{0})=L\rho_{0}L-\tfrac{1}{2}L^{2}\rho_{0}-\tfrac{1}{2}\rho_{0}L^{2}=\lambda_{k}^{2}\rho-\tfrac{1}{2}(2\lambda_{k}\rho_{0})=0$ and $\mathcal{M}_{L}(\rho_{0})=\sqrt{\eta}(L\rho_{0}+\rho_{0}L-\operatorname{Tr}\left(2L\,\rho_{0}\right)\rho_{0})=\sqrt{\eta}(2\lambda_{k}\rho_{0}-\operatorname{Tr}\left(2\lambda_{k}\rho_{0}\right)\rho_{0})=0$. Thus $\rho_{0}$ is a steady state of (1).

(iii) $\;V_{o}$ is a positive definite function on $\mathcal{S}$ and it equals [math] only when $\boldsymbol{p}_{\ell}(\rho)=1$ for some $\ell$. It remains to check that it is a supermartingale with exponential decay. By Itō’s formula (14), the variable $\xi_{k}:=\sqrt{\boldsymbol{p}_{k}}$ satisfies

[TABLE]

with $\varpi(\xi)=\sum_{k=1}^{d}\lambda_{k}\xi_{k}^{2}$.

Since $V_{o}(\rho)=\sum_{1\leq k<k^{\prime}\leq d}^{d}\xi_{k}\xi_{k^{\prime}}$, consider the computation with Itō’s formula:

[TABLE]

The Markov generator $\mathcal{A}$ (see (15)) following from the expectation of this equation thus yields:

[TABLE]

Since each $\xi_{k}(t)$ remains non-negative for all $t$, this readily yields

[TABLE]

By Theorem A.1 in Appendix, $V_{o}$ decays exponentially towards zero, concluding the proof. $\square$

The open-loop system corresponding to a QND measurement (1),(2) thus converges towards an eigenstate of $L$, satisfying $\boldsymbol{p}_{k}(\rho)=1$ for some $k$, for each realization. However, the particular eigenstate $k$ will be random, with correlated measurement results indicating which state has been chosen; from (i), the probability to converge towards the particular eigenstate $k$ is equal to $\boldsymbol{p}_{k}(\rho_{0})$.

The control objective is to ensure convergence to a target QND eigenstate, indexed by $\ell\in\{1,\ldots,d\}$, for all realizations. More precisely, we will design a continuous stochastic real feedback process $v$, depending on the state $\rho$, such that $\lim_{t\rightarrow\infty}\mathbb{E}[\boldsymbol{p}_{\ell}(\rho_{t})]=1$ with exponential convergence rate, for any initial condition $\rho_{0}\in\mathcal{S}$. The feedback action is modeled by means of a unitary control operation $U_{t}=e^{-iHudt}$ during the infinitesimal interval $[t,t+dt]$, where $H=H^{\dagger}$ is the actuator Hamiltonian and $dv=udt$ is the feedback input signal. Applied on (1), the closed-loop system reads then:

[TABLE]

with measurement process $Y$ still given by (2). Regarding exponential convergence, we aim to provide a global Lyapunov function $V(\rho)$ such that in closed-loop, $V(\rho_{t})$ is a supermartingale with exponential decay for all $t\geq 0$ and all $\rho_{0}\in\mathcal{S}$. Providing feedback controls that ensure decay of $V$ in this sense has remained so far an open issue.

3 Existing feedback designs

Measurement-based quantum feedback has been investigated thoroughly in the literature, considering a wide array of applications such as state preparation [33, 25], state purification [15, 43] or continous-time quantum error correction [1, 31]. Nevertheless, to our knowledge, control design has been done mainly on two feedback architectures: static output feedback [42] (so-called Markovian feedback in the physics literature) and Bayesian feedback [17, 18, 25] involving a full state estimate. Intermediate approaches have attracted much less attention and promise an interesting research area.

3.1 Static output feedback

Static output feedback corresponds to quantum feedback of the form $dv=udt=fdt+\sigma dY$, where $f$ and $\sigma$ are constant. From (5) the dynamics in closed-loop read (see e.g. [41, 40]):

[TABLE]

The simplicity of the feedback scheme makes it attractive for experimental implementations, since it avoids any overhead associated to dynamical computations in the feedback loop; in particular it needs no quantum state observer. With proper tuning of the constants $f,\sigma$, Markovian feedback allows to exponentially stabilize a range of target states, with direct algebraic proofs [42, 40, 35, 36, 34].

Unfortunately, an invariance argument shows that precisely the QND eigenstates, which are of particular interest in quantum engineering, are not stabilizable using this static output feedback with $f,\sigma$ just constants ([12]).

It seems that a more involved controller is needed in order to bias the stochastic evolution towards a prescribed QND measurement eigenstate. One solution is to turn to state feedback, assuming an underlying quantum state observer.

3.2 State feedback

The standard state feedback takes the form $\;dv=udt=f(\rho)dt\;$. Thus the control signal is a deterministic scalar function of the state $\rho$. The closed-loop model for this quantum state feedback just takes the form:

[TABLE]

Stabilization of QND eigenstates under this controlled dynamics has been treated extensively in the literature [1, 39, 44, 25, 37, 38, 24]. These results have succeeded in proving asymptotic convergence, in proper probabilistic settings. There is however room for improvement on a few aspects.

A first aspect is the convergence speed. Most papers do not provide any convergence rate; as the strongest result so far to our knowledge, [24] provide an estimate of the Lyapunov exponent for a qubit, valid for the final approach of the target state after an unspecified final initial transient. By Lemma 2.1, the open-loop system under QND measurements converges towards the set of its steady states at global exponential speed. The absence of a proven similar property for the selection of one target QND steady state thus appears as an avoidable gap.

A second aspect is that full state feedback, as is used in the above proposals, may appear unnecessary for this application. Running a quantum state estimate in real-time does pose experimental challenges, given the very short timescales involved (nanoseconds) and the plan to ultimately scale quantum computers to high-dimensional systems. In the present task, the asymptotic behavior is directly visible on the measurement signal. Indeed, the measurement signal corresponding to $\boldsymbol{p}_{k}(\rho_{t})=1$ for all $t\geq 0$ is $Y_{t}=Y_{0}+2\lambda_{k}t+W_{t}$, with thus an expectation that directly informs on the eigenstate via the drift $\lambda_{k}\,t$, and a standard deviation in $\sqrt{t}$. This suggests that simple filtering should allow to essentially solve the task too. The necessity of keeping a full state observer, including quantum coherences, comes from the actuation strategy. This could be improved.

A third aspect, more related to model uncertainties, is that many existing proposals work with short pulses, inspired from the discrete-time counterpart. In actual implementations, it may be more cautious to use smoother control signals in order to avoid exciting spurious dynamics.

With respect to these two main approaches, we thus aim for an intermediate solution using a reduced estimator and smoother controls, while providing an exponential convergence guarantee. For this we resort to a control signal with a novel structure, called noise-assisted feedback.

4 Noise-assisted feedback scheme

In our feedback scheme, the control signal $dv=udt$ still depends on (part of) the state like in the Bayesian feedback approach, but instead of just involving a deterministic function we drive $dv$ by an exogenous Brownian noise $dB$ independent of $W_{t}$:

[TABLE]

The feedback control occurs by making the noise gain $\sigma(\rho)$ a continuously differentiable function of $\rho$. Using in (5) the Baker-Campbell-Hausdorff (BCH) formula $e^{xH}\rho e^{-xH}=\sum_{j}T_{j}\frac{x^{j}}{j!},$ where $T_{0}=\rho$ and $T_{j+1}=[H,T_{j}]$ for $j\geq 0$, and applying the Itō rules $(dv)^{2}=\sigma^{2}(\rho)dt$, $dtdv=dtdB=dWdB=dt^{2}=0$, yields the following closed loop dynamics with two independent Wiener processes $W$ and $B$ :

[TABLE]

The main idea comes down to noise discouraging the system to converge towards an eigenstate different from $\boldsymbol{p}_{\ell}(\rho)=1$. Indeed, as the open-loop QND dynamics stochastically converges to one of the QND eigenstates, but on the average does not move closer to any particular one, it is sufficient to activate noise only when the state is close to a bad equilibrium in order to “shake it away” and induce global convergence to the target. Accordingly, we consider the feedback (6) with the following gain law:

[TABLE]

where $\varphi\geq 0$ is a smooth saturating function on $[0,1]$, i.e. $\varphi(]-\infty,0])=\{0\}$ and $\varphi([1,+\infty[)=\{1\}$, with parameters $\bar{\sigma}>0$ and $1>p_{\text{\tiny max}}>p_{\text{\tiny min}}>\frac{1}{2}$. Since $\sum_{k}\boldsymbol{p}_{k}=1$, each $\boldsymbol{p}_{k}\geq 0$ and $p_{\text{\tiny min}}>1/2$, the argument of the max can only change when $\max_{k\neq\ell}\boldsymbol{p}_{k}(\rho)-p_{\text{\tiny min}}<0$. Therefore, the function $\rho\mapsto\sigma(\rho)$ is smooth despite the use of a max in its definition.

5 Exponential stabilization via noise-assisted feedback

We now construct a Lyapunov function and prove that our feedback design ensures its exponential convergence. A particular point for $n$-level systems, compared to $2$-level systems like in [12, 24], is to take the limited actuation into account.

Inspired by [3], we consider the $d\times d$ real symmetric matrix $\Delta$ with components

[TABLE]

combining the spectral decomposition $L=\sum_{k}\lambda_{k}\Pi_{k}$ with the actuator Hamiltonian $H$. Its off-diagonal elements are non negative since $\Delta_{k,k^{\prime}}=\operatorname{Tr}\left(\Pi_{k}H\Pi_{k^{\prime}}H\right)\geq 0$ for $k\neq k^{\prime}$. Its diagonal elements are non positive and

[TABLE]

Thus $\Delta$ is a Laplacian matrix.

Theorem 5.1

Assume that $L$ is nondegenerate, i.e., $d=n$ and each projector $\Pi_{k}$ is a rank one projector. Consider the closed-loop system (7) with feedback gain $\sigma(\rho)$ given by (8) for a given projector $\Pi_{\ell}$ with $\ell\in\{1,\ldots,d\}$. Assume that the graph associated to the Laplacian matrix $\Delta$ defined in (9) is connected. Then there exists $\overline{p}\in]\tfrac{1}{2},1[$ such that for any choice of parameter $\bar{\sigma}>0$ and parameters $1>p_{\text{\tiny max}}>p_{\text{\tiny min}}\geq\overline{p}$, the closed-loop trajectories converge exponentially to $\rho=\Pi_{\ell}$, in the sense that: there exist constants $\nu>0$ and $C>0$ (depending on $\bar{\sigma}$, $p_{\text{\tiny max}}$, $p_{\text{\tiny min}}$) for which $\;\mathbb{E}\left[\sqrt{1-\boldsymbol{p}_{\ell}(\rho_{t})}\right]\leq Ce^{-\nu t}\sqrt{1-\boldsymbol{p}_{\ell}(\rho_{0})}$ for any initial state $\rho_{0}\in\mathcal{S}$.

If the graph associated to $\Delta$ is not fully connected then there exists a partition of $\{1,\ldots,d\}=I\cup J$ ($I,J\neq\emptyset$, $I\cap J=\emptyset$) such that $H=\Pi_{I}H\Pi_{I}+\Pi_{J}H\Pi_{J}$ with $\Pi_{I}=\sum_{k\in I}\Pi_{k}$ and $\Pi_{J}=\sum_{k\in J}\Pi_{k}$. Then any trajectory $\rho_{t}$ of (5) with any feedback scheme starting from $\operatorname{Tr}\left(\rho_{0}\Pi_{I}\right)=0$ (resp. $\operatorname{Tr}\left(\rho_{0}\Pi_{J}\right)=0$), satisfies $\operatorname{Tr}\left(\rho_{t}\Pi_{I}\right)=0$ (resp. $\operatorname{Tr}\left(\rho_{t}\Pi_{J}\right)=0$) for all $t>0$. Thus, closed-loop convergence to $\boldsymbol{p}_{\ell}=1$ with $\ell\in I$ is impossible when $\operatorname{Tr}\left(\rho_{0}\Pi_{I}\right)=\sum_{k\in I}\boldsymbol{p}_{k}(\rho_{0})<1$. In this sense the above connectivity condition on the graph of $\Delta$ cannot be weakened.

Proof:

(– Lyapunov function construction –) We do not use directly $\sqrt{1-\boldsymbol{p}_{\ell}(\rho)}$ as a closed-loop Lyapunov function. Instead we construct a closed-loop Lyapunov function $V(\rho)$ equivalent to $\sqrt{1-\boldsymbol{p}_{\ell}}$ (i.e. $c_{*}V(\rho)\leq\sqrt{1-\boldsymbol{p}_{\ell}(\rho)}\leq c^{*}V(\rho)$ with $0<c_{*}<c^{*}$) such that $\mathcal{A}V\leq-rV$. More precisely, we use

[TABLE]

The positive parameters $\alpha_{s,k}$ will be given by solving $d-1$ linear systems, indexed by $s$:

[TABLE]

with $\beta_{s,k}>0$ for $k\neq\ell$ and $\beta_{s,\ell}=-\sum_{k\neq\ell}\beta_{s,k}$. Standard arguments used in [3] guarantee under the connectivity assumption that there exists, for each $s$, a unique solution $(\alpha_{s,k})$ such that $\alpha_{s,k}>0$ for $k\neq\ell$ and $\alpha_{s,\ell}=0.$ (see, e.g. [6, Chapter 4]). When the $(d-1)\times d$ matrix $\beta$ is chosen of maximal rank $d-1$, the obtained matrix $\alpha$ is also of maximal rank $d-1$. Since $\sum_{k\in\{1,\ldots,d\}\setminus\{\ell\}}\boldsymbol{p}_{k}=1-\boldsymbol{p}_{\ell}$, one has

[TABLE]

where $c_{*}=\frac{1}{(d-1)\sqrt{\alpha^{*}}}$ and $c^{*}=\frac{1}{(d-1)\sqrt{\alpha_{*}}}$ with

[TABLE]

(– expression of the stochastic Markov generator –) The rest of the proof consists in showing that for any such choice of maximal rank matrix $\beta$, the resulting $V_{\alpha}$ becomes an exponential Lyapunov function as soon as $\overline{p}$ is close enough to $1$ and $p_{\text{\tiny min}}>\overline{p}$. This is based on the following simple but slightly tedious computations of $\mathcal{A}V_{\alpha}$:

[TABLE]

with

[TABLE]

These expressions are obtained with the following general formula based on Itō rules ($a_{k}>0$ constant)

[TABLE]

with $\mathbb{E}[d\boldsymbol{p}_{k}|\rho]=\sigma^{2}(\rho)\operatorname{Tr}\left(\Pi_{k}\mathcal{D}_{H}(\rho)\right)~{}dt$ and

[TABLE]

(– closed-loop essential contribution –) For $p\in[0,1]$ let

[TABLE]

Take $\rho\in\mathcal{S}_{p}$ with $p=1$. Then there exists $j\in\{1,\ldots,d\}\setminus\{\ell\}$ such that $\rho=\Pi_{j}$, thus $\sum_{k}\alpha_{s,k}\boldsymbol{p}_{k}(\rho)=\alpha_{s,j}>0$ and

[TABLE]

Consequently $f_{\alpha}(\rho)=-\sum_{s\in\{1,\ldots,d\}\setminus\{\ell\}}\frac{\beta_{s,j}}{\sqrt{\alpha_{s,j}}}<0$ with $V_{\alpha}(\rho)=\sum_{s}\sqrt{\alpha_{s,j}}>0$.

By continuity of $f_{\alpha}/V_{\alpha}$ on $\mathcal{S}_{1/2}$, there exist $\epsilon_{f}>0$ and $\overline{p}\in]1/2,1[$ such that

[TABLE]

Taking $p_{\text{\tiny min}}>\overline{p}$, we ensure that the feedback will only be turned on when it contributes a negative term to $\mathcal{AV}_{\alpha}$.

(– open-loop essential contribution –) For all $\rho\in\mathcal{S}\setminus\{\Pi_{\ell}\}$, $\chi(\rho)=g_{\alpha}(\rho)/V_{\alpha}(\rho)$ is well defined. Since $\operatorname{Tr}\left(\rho L\right)=\sum_{k}\lambda_{k}\boldsymbol{p}_{k}(\rho)$, the function $\chi(\rho)$ depends only on the populations $\boldsymbol{p}_{k}$. Consider the following parametrization exploiting the degree 0 homogeneity of $\chi$ in the populations:

[TABLE]

For $\rho\in\mathcal{S}\setminus\{\Pi_{\ell}\}$, the function $\chi$ admits the following smooth expression with the variables $r\in]0,1]$ and $x_{k}\in[0,1]$ satisfying $\sum_{k\neq\ell}x_{k}=1$:

[TABLE]

with $\varpi(r,x)=(1-r)\lambda_{\ell}+r\left(\sum_{k\neq\ell}\lambda_{k}x_{k}\right)$. We study its extension to the compact set with $r\in[0,1]$. Clearly $\chi(r,x)\geq 0$. Consider the solutions $(r,x)$ of $\chi(r,x)=0$. Necessarily, they satisfy

[TABLE]

Since the $(d-1)\times d$ matrix $\alpha_{s,k}$ is of maximal rank $d-1$, $\forall k\neq\ell$ we have $(\lambda_{k}-\varpi(r,x))x_{k}=0.$ Taking the sum over $k\neq\ell$, one gets $(1-r)\left(\lambda_{\ell}-\sum_{k\neq\ell}\lambda_{k}x_{k}\right)=0$ implying two possibilities:

• •

if $\lambda_{\ell}=\sum_{k\neq\ell}\lambda_{k}x_{k}$, then $\varpi(r,x)=\lambda_{\ell}$ and the condition before summing requires $(\lambda_{k}-\lambda_{\ell})x_{k}=0$ for $k\neq\ell$. Since $\lambda_{k}\neq\lambda_{\ell}$ this implies that $x_{k}=0$ $\forall k\neq\ell$. But this is not possible since $\sum_{k\neq\ell}x_{k}=1$.

• •

if $r=1$, then $\varpi(r,x)=\sum_{k^{\prime}\neq\ell}\lambda_{k^{\prime}}x_{k^{\prime}}$. For $k\neq\ell$ the conditions become $\left(\lambda_{k}-\sum_{k^{\prime}\neq\ell}\lambda_{k^{\prime}}x_{k^{\prime}}\right)x_{k}=0$ before summing. Since $\lambda_{k}\neq\lambda_{k^{\prime}}$ for $k\neq k^{\prime}$ and $x_{k}\in[0,1]$ with $\sum_{k\neq\ell}x_{k}=1$, there necessarily exists $j\neq\ell$ such that $x_{j}=1$ and $x_{k}=0$ for $k\notin\{\ell,j\}$.

Consequently, the nonnegative smooth function $\chi(r,x)$ vanishes only at $d-1$ isolated points, where $r=1$ and $x_{k}=\delta_{kj}$ for some $j\neq\ell$. By continuity, for any $p\in]0,1[$, there exists $\theta_{p}>0$ such that $\forall\rho\in\mathcal{S}\setminus(\mathcal{S}_{p}\cup\{\Pi_{\ell}\})$, we have $\chi(\rho)\geq\theta_{p}$. This proves that

[TABLE]

(– bringing all pieces together –) To conclude, consider $\mathcal{A}V_{\alpha}$ given in (10). Since $h_{\alpha}\geq 0$, we have

[TABLE]

Consider the feedback gain $\sigma(\rho)$ with $p_{\text{\tiny min}}>\overline{p}$. Then from the closed-loop essential contribution,

[TABLE]

and from the open-loop essential contribution,

[TABLE]

Thus for all $\rho$ in $\mathcal{S}$, one has $\mathcal{AV_{\alpha}}(\rho)\leq-\nu V_{\alpha}(\rho)$ where $\nu$ is equal to $\min\left(\tfrac{\bar{\sigma}^{2}\epsilon_{f}}{2},\tfrac{\eta\theta_{p_{\text{\tiny max}}}}{2}\right)$. A direct application of Theorem A.1 recalled in Appendix ensures $\mathbb{E}(V_{\alpha}(\rho_{t})]\leq V_{\alpha}(\rho_{0})e^{-\nu t}$. $\square$

6 Observer and approximate quantum filtering

The feedback is based on the value $\rho_{t}$ of the quantum state, which is not directly measured. One has to reconstruct in real-time this quantum state via a quantum filter, a quantum state-observer. For (7), it reads (see e.g. [8]):

[TABLE]

A standard result [33, 39, 25, 2] ensures that the resulting observer/controller closed-loop system

[TABLE]

converges almost surely towards the stationary state $\rho=\widehat{\rho}=\Pi_{\ell}$ as soon as $\operatorname{Tr}\left(\rho_{0}\widehat{\rho}_{0}\right)>0$ and the convergence conditions of Theorem 5.1 are satisfied.

Practical implementation of such quantum filter could be a problematic issue: when the dimension $n$ of the system is large, it requires to store and update in real-time the $n(n-1)/2$ components of the operator $\widehat{\rho}$. Since the feedback depends only on the populations, i.e., the diagonal of $\rho$ in the eigenbasis of the measurement operator $L$, the development of a reduced-order (possibly approximate) quantum filter depending only on these populations is suggested.

A first reduction consists in replacing $\rho_{t}$ in the feedback law by $\widehat{\varrho}_{t}$ corresponding to the Bayesian estimate of $\rho_{t}$ knowing $\rho_{0}$ and $Y_{\tau}$ for $\tau\in[0,t]$. One thus discards here the knowledge of $B_{t}$. Then, one can prove that $\widehat{\varrho}_{t}$ obeys to the following stochastic differential equation:

[TABLE]

This filter involves less computations but still the full matrix $\widehat{\varrho}$. Since the feedback law only depends on the $\boldsymbol{p}_{k}(\widehat{\varrho})=\operatorname{Tr}\left(\widehat{\varrho}\Pi_{k}\right)$, it would be ideal to have a reduced filter involving only those variables. In cases like [13], the filter (11) can be directly and exactly reduced to an autonomous system on the $\boldsymbol{p}_{k}(\widehat{\varrho})$, as we discard the random coherences among eigenstates induced by $dB_{t}$. However, when a single noise process $dB_{t}$ drives many levels at once, the filter (11) has to model correlations among the various random coherences and this precludes an exact reduction. Nevertheless, by neglecting those feedback-induced correlations (recalling that feedback actuation is often turned off), we can set coherences $\left\langle k|\widehat{\varrho}|k^{\prime}\right\rangle$, $k\neq k^{\prime}$ to zero, and propose an approximate filter for populations $\widehat{\boldsymbol{p}}_{k}$ to estimate $\boldsymbol{p}_{k}(\widehat{\varrho})$. This amounts to replacing $\mathcal{D}_{H}(\widehat{\varrho})$ by a population transfer via the Laplacian matrix $\Delta$ defined in (9). It yields:

[TABLE]

where $\varpi(\widehat{\boldsymbol{p}})=\sum_{k^{\prime}=1}^{d}\lambda_{k^{\prime}}\widehat{\boldsymbol{p}}_{k^{\prime}}$. This approximate filter requires to store and update in real-time only $d$ real numbers. For any measurement trajectory $Y_{t}$, the components of $\widehat{\boldsymbol{p}}$ remain nonnegative and their sum equal to one. In open-loop ($\sigma\equiv 0$), this population filter is exact.

7 Simulation and robustness issues

Theorem 5.1 ensures for $p_{\text{\tiny min}}$ close enough to $1$ exponential closed-loop convergence of (5) with noise-assisted feedback (6),(8). This section is devoted to numerical estimation of closed-loop convergence rates and investigation of the related robustness on a specific quantum system already considered in [25, 33]: a spin $J$ system of dimension $n=2J+1$ where the measurement operator is

[TABLE]

and the actuator Hamiltonian is a tridiagonal matrix

[TABLE]

The Hilbert space is spanned by the $2J+1$ orthonormal vectors $|J\text{-}m\rangle$ for $m=0,\ldots,2J$.

All the simulations below correspond to detection efficiency $\eta=0.8$ and $J=2$ ($n=5$), for which:

[TABLE]

According to Lemma 2.1, the open-loop convergence rate is $\eta/2=0.4$. As a generic control goal we choose to stabilize the state $\Pi_{\ell}=|0\rangle\langle 0|$ associated to the zero eigenvalue of $L$. For the noise-assisted feedback gain (8), we take $\bar{\sigma}=\sqrt{5\eta}$, $p_{\text{\tiny max}}=p_{\text{\tiny min}}+0.05$, $p_{\text{\tiny min}}=0.9$ or $p_{\text{\tiny min}}=0.6$ with the saturation function $\varphi(s)=\min(1,\max(0,s))$. For each case we simulate a set of 1000 realizations starting from the fully depolarized state $\rho_{0}=\boldsymbol{I}/5$. We estimate the evolution of $\mathbb{E}[\sqrt{1-\operatorname{Tr}\left(\Pi_{\ell}\rho_{t}\right)}]\equiv\mathbb{E}[\sqrt{1-\left\langle 0|\rho_{t}|0\right\rangle}]$ by taking the ensemble average over these 1000 realizations. The Laplacian matrix

[TABLE]

inherits the tridiagonal structure of $H$ and thus admits a connected graph, so Theorem 5.1 predicts exponential convergence with our controller using a sufficiently high value of $p_{\text{\tiny min}}$ and the filter (11).

We first illustrate the exponential convergence rate $\nu$. The proof of Theorem 5.1 only provides a very loose bound on both $\nu$ and the necessary $p_{\text{\tiny min}}$, so we here stick to numerical simulations. From the system analysis, the qualitative trend should be that lower values of $p_{\text{\tiny min}}$ imply more frequent feedback corrections and thus faster convergence. Figure 1 illustrates a simulation set with $p_{\text{\tiny min}}=0.9$, which fits an exponential convergence at rate $\nu\approx 0.04$. In a second simulation set, see Figure 2, we have pushed this to $p_{\text{\tiny min}}=0.6$ and observed a much faster convergence rate around $\nu\approx 0.2$, i.e. one half of the open-loop convergence rate. This suggests that such noise-assisted feedback can be tuned to achieve convergence rates similar to the open-loop one.

We next investigate the behavior with the approximate reduced filter on the population vector $\widehat{\boldsymbol{p}}$. The simulations of Fig. 3 differ from the ones of Fig. 2 just by replacing in the feedback law the ideal populations $\boldsymbol{p}$ by the approximated ones $\widehat{\boldsymbol{p}}$ solutions of (12). One observes still an exponential convergence, with a reasonable decrease of the convergence rate from $0.19$ to $0.12$, illustrating the practical interest of such low-dimensional filters.

Finally, on Fig. 4 we further investigate robustness of the control strategy by using the reduced filter (12) with a feedback delay of $0.5$ time units in the closed-loop simulations. Despite a decrease by a factor two of the convergence rate, the fact a feedback latency of $1/4$ of the open-loop convergence time does not destabilize this feedback scheme suggests promising robustness properties.

8 Concluding remarks

We have approached the problem of stabilizing a QND measurement eigenstate in continuous-time by introducing Brownian noise to drive the control. The use of noise provides a simple controller that shakes away spurious steady states in closed loop and thus achieves exponential stabilization of a target eigenstate.

The present work still leaves room for improvement in, at least, the following directions:

• •

While our proof of exponential convergence can provide an estimate of the convergence rate, we did not try to optimize the Lyapunov function nor the other control parameters in order to maximize the speed of convergence. We have shown numerically that the closed-loop convergence rate apparently can be made similar to the open-loop convergence rate with our approach, suggesting promising results in a precise analysis of convergence rate.

• •

In numerical simulations, the reduced approximate filter (12) appears good enough to achieve global exponential stabilization. We conjecture that this can be proven, for this filter and possibly for even simpler filters based on direct output signal filtering or sparse $\widehat{\boldsymbol{p}}$.

• •

The capability of performing several quantum measurements and of applying different unitary feedback controls on a single quantum system motivates the study of multi-input multi-output (MIMO) quantum feedback schemes. A MIMO version of static output feedback was introduced in [14] and an implementation was made in [11]. In addition, we have presented, in the context of quantum error correction [13], a MIMO scheme of this noise-assisted feedback to stabilize a manifold of quantum states. A general MIMO extension of Theorem 5.1 should be feasible along the same lines as the present work.

Appendix A Tools from stochastic stability

We refer the reader to [23] for further reference on these fundamental results. We consider concrete instances of Itō stochastic differential equations (SDEs) on $\mathbb{R}^{n}$ of the form

[TABLE]

where $W_{t}$ is a standard Brownian motion on $\mathbb{R}^{k}$, and $\mu,\theta$ are regular functions of $x$ with image in $\mathbb{R}^{n}$ and $\mathbb{R}^{n\times k}$ respectively, satisfying the usual conditions for existence and uniqueness of solutions [23] on $\mathcal{S}$, a compact and positively invariant subset of $\mathbb{R}^{n}$.

Let $\mathcal{I}:=\{x\in\mathcal{S}:\mu(x)=\theta(x)=0\}$ be an invariant set of (13). Let $V(x)$, a nonnegative real-valued twice continuously differentiable function on $\mathcal{S}\setminus\mathcal{I}$, with $V(x)=0$ implying $x\in\mathcal{I}$. Itō’s formula on $V$ yields [28]

[TABLE]

The Markov generator $\mathcal{A}$ of the SDE (13) is defined for any function $V$ in its domain by

[TABLE]

It is related to the SDE (14) by [28, Chapter 7]

[TABLE]

The stochastic counterpart of Lyapunov’s second method provides sufficient conditions for stochastic stability by analyzing the generator $\mathcal{A}V$, e.g.:

Theorem A.1** (Khasminskii [23])**

If for some $r>0$ the Markov generator satisfies $\mathcal{A}V(x)\leq-rV(x)$ $\forall x\in\mathcal{S}$, then $\mathbb{E}[V(x(t))]\leq V(x(0))\;\exp(-r\,t)$ $\forall t\geq 0$, i.e., $V(x(t))$ is a supermartingale with exponential decay. Thus $\lim_{t\rightarrow\infty}\mathbb{E}[V(x(t))]=0$; since $V(x)=0$ implies $x\in\mathcal{I}$, this implies convergence towards $\mathcal{I}$ of all the solutions of Eq. (13) starting in $\mathcal{S}$.