text:
# Feedback-stabilization of an arbitrary pure state of a two-level atom

## I Introduction

Although classical models of feedback schemes have been used for a long time to control dynamical noise, an analogous quantum theory of feedback has been developed only in the last fifteen years . Recently there has been considerable interest in quantum feedback as a way to fight decoherence in isolated quantum systems, using the approach of Refs. . The central idea is to use a continuous measurement record, whose existence is due to the coupling of the system to a bath, to control the dynamics of the system so as to counteract the noise introduced by that bath and possibly other baths. For example, it has been suggested as a way to create optical squeezed states , to create micromaser number states , to correct errors in quantum bits , and to protect optical and microwave Schrödinger cat states against dissipation .

Decoherence in quantum systems can be loosely defined as loss of purity. Therefore the ultimate success in using feedback to fight decoherence would be to create an arbitrary stable pure state in the presence of dissipation. This goal was realized (better even that they realized) by Hoffman, Mahler and Hess (HMH)  for a very simple system: a resonantly driven two-level atom. They showed that by using the photocurrent derived from unit-efficiency homodyne detection of the atom’s fluorescence to control part of the driving field of the atom, it is possible to exactly cancel the noise introduced by the electromagnetic vacuum field when the atom is in a particular pure state. By choosing the driving strength and feedback strength appropriately, any pure state on the Bloch sphere may be picked out, although HMH claimed that only pure states on the lower half of the sphere would be stable under their scheme.

HMH chose to describe detection and feedback in their system in a way different from (but equivalent to) the standard approach in Refs. . In this paper we reformulate their theory using the latter approach. This has the advantage of enabling a number of generalizations of their results. First, we revisit the question of stability and find that, contrary to the claims of HMH, the states in the upper half of the Bloch sphere can be stabilized as well as those in the lower half (this is what was better than they realized). The only states which cannot be stabilized, in the sense that an arbitrary initial state would not always end up in the desired state, are those on the equator of the Bloch sphere; that is, those which are equal superpositions of excited and ground states.

Our second generalization is to consider how effective feedback is with $`\eta <1`$; that is, with non-unit-efficiency detection. In this case it is not possible to stabilize the atom at any pure state, except the ground state which is trivially stable by setting the driving and feedback to zero. Instead, we aim to produce a steady state which is as close as possible to a given pure state. For the two-level atom, this is equivalent to trying to create a state which is as pure as possible in a particular direction in Bloch space. Not surprisingly (given the above result), we find that states near the equator cannot be well-protected against decoherence. We also find an echo of the distinction HMH found between the upper and lower halves of the Bloch sphere, in that states in the upper half sphere are affected much more by loss of detection efficiency that those in the lower half.

The paper is organized as follows. In Sec. II we present the model of a driven two-level atom, including the stochastic Schrödinger equation for unit-efficiency homodyne detection. In Sec. III we use this equation to derive the driving and feedback required to stabilize the atom in an arbitrary pure state. These results agree with those of HMH. However, our stability analysis disagrees substantially with theirs. In Sec. IV we present entirely new analytical results relating to the effect of non-unit- efficiency detection. In Sec. V we give numerical simulations of the stochastic evolution equations, illustrating the issues discussed in the preceding two sections. In Sec. VI we summarize and interpret our results, explain their significance, and discuss the possibility of future work.

## II Homodyne Detection

### A Master Equation

Consider a atom, with two relevant levels $`\{|g,|e\}`$ and lowering operator $`\sigma =|ge|`$. Let the the decay rate be $`\gamma `$, and let it be driven by a resonant classical driving field with Rabi frequency $`2\alpha `$. This is as shown in Fig. 1, where for the moment we are omitting feedback by setting $`\lambda =0`$. This system is well-approximated by the master equation

$$\dot{\rho }=\gamma 𝒟[\sigma ]\rho i\alpha [\sigma _y,\rho ],$$
(1)
where the Lindblad  superoperator is defined as usual $`𝒟[A]BABA^{}\{A^{}A,B\}/2`$. In this master equation we have chosen to define the $`\sigma _x=\sigma +\sigma ^{}`$ and $`\sigma _y=i\sigma i\sigma ^{}`$ quadratures of the atomic dipole relative to the driving field. The effect of driving is to rotate the atom in Bloch space around the $`y`$-axis. The state of the atom in Bloch space is described by the three-vector $`(x,y,z)`$. It is related to the state matrix $`\rho `$ by

$$\rho =\frac{1}{2}\left(I+x\sigma _x+y\sigma _y+z\sigma _z\right).$$
(2)
It is easy to show that the stationary solution of the master equation (1) is

$`x_{\mathrm{ss}}`$ $`=`$ $`{\displaystyle \frac{4\alpha \gamma }{\gamma ^2+8\alpha ^2}},`$ (3)
$`y_{\mathrm{ss}}`$ $`=`$ $`0,`$ (4)
$`z_{\mathrm{ss}}`$ $`=`$ $`{\displaystyle \frac{\gamma ^2}{\gamma ^2+8\alpha ^2}}.`$ (5)
For $`\gamma `$ fixed, this is a family of solutions parameterized by the driving strength $`\alpha (\mathrm{},\mathrm{})`$. All members of the family are in the $`x`$$`z`$ plane on the Bloch sphere. Thus for this purpose we can reparametrize the relevant states using $`r`$ and $`\theta `$ by

$`x`$ $`=`$ $`r\mathrm{sin}\theta `$ (6)
$`z`$ $`=`$ $`r\mathrm{cos}\theta ,`$ (7)
where $`\theta [\pi ,\pi ]`$. Since

$$\mathrm{Tr}[\rho ^2]=\frac{1}{2}\left(1+x^2+y^2+z^2\right)$$
(8)
is a measure of the purity of the Bloch sphere, $`r=\sqrt{x^2+z^2}`$, the distance from the centre of the sphere, is also a measure of purity. Pure states correspond to $`r=1`$ and maximally mixed states to $`r=0`$.

The locus of solutions in this plane (an ellipse) is shown in Fig. 2. Since $`z_{\mathrm{ss}}<0`$, all solutions are in the lower half of the Bloch sphere. That is, we are restricted to $`|\theta |>\pi /2`$. Also, it is evident that the smaller $`|\theta |`$ is (that is, the more excited the atom is), the smaller $`r`$ is (that is, the less pure the atom is). At $`|\theta |=\pi `$ the stationary state is pure, but this is not surprising as it is simply the ground state of the atom with no driving. As $`|\theta |\pi /2`$ we have $`r0`$. This can only be approached asymptotically as $`|\alpha |\mathrm{}`$. In summary, the stationary states we can reach by driving the atom are limited, and generally far from pure.

### B Homodyne Detection

Now consider subjecting the atom to homodyne detection. As shown in Fig. 1, we assume that all of the fluorescence of the atom is collected and turned into a beam (represented in Fig. 1 by placing the atom at the focus of a mirror). Ignoring the vacuum fluctuations in the field, the annihilation operator for this beam is $`\sqrt{\gamma }\sigma `$, normalized so that the mean intensity $`\gamma \sigma ^{}\sigma `$ is equal to the number of photons per unit time in the beam. This beam then enters one port of a 50:50 beam splitter, while a strong local oscillator $`\beta `$ enters the other. To ensure that this local oscillator has a fixed phase relationship with the driving laser used in the measurement, it would be natural to utilize the same coherent light field source in the driving process and as the local oscillator in the homodyne detection. This is as shown in Fig. 1.

Again ignoring vacuum fluctuations, the two field operators exiting the beam splitter, $`b_1`$ and $`b_2`$, are

$$b_k=\left[\sqrt{\gamma }\sigma (1)^k\beta \right]/\sqrt{2}.$$
(9)
When these two fields are detected, the two photocurrents produced have means

$$\overline{I}_k=|\beta |^2(1)^k\left(\sqrt{\gamma }\beta \sigma ^{}+\sqrt{\gamma }\sigma \beta ^{}\right)+\gamma \sigma ^{}\sigma /2.$$
(10)
The middle two terms represent the interference between the system and the local oscillator.

The ideal limit of homodyne detection is when the local oscillator amplitude goes to infinity, which in practical terms means $`|\beta |^2\gamma `$. In this limit, the rate of the photodetections goes to infinity and thus it should be possible to change the point process of photocounts into a continuous photocurrent with white noise. Also, the only relevant quantity is the difference between the two photocurrents. Suitably normalized, this is

$$I(t)=\frac{I_1(t)I_2(t)}{|\beta |}=\sqrt{\gamma }e^{i\mathrm{\Phi }}\sigma ^{}+e^{i\mathrm{\Phi }}\sigma _\mathrm{c}(t)+\xi (t).$$
(11)
A number of aspects of Eq. (11) need to be explained. First, $`\mathrm{\Phi }=\mathrm{arg}\beta `$, the phase of the local oscillator (defined relative to the driving field). Second, the subscript c means conditioned and refers to the fact that if one is making a homodyne measurement then this yields information about the system. Hence, any system averages will be conditioned on the previous photocurrent record. Third, the final term $`\xi (t)`$ represents Gaussian white noise, so that

$$\xi (t)dt=dW(t),$$
(12)
an infinitesimal Wiener increment defined by

$`dW(t)^2=dt,`$ (13)
$`\mathrm{E}[dW(t)]=0.`$ (14)
Since the stationary solution of the master equation confines the state to the $`x`$$`z`$ plane, it makes sense to follow HMH by setting $`\mathrm{\Phi }=0`$. In that case,

$$I(t)=\sqrt{\gamma }\sigma _x_\mathrm{c}(t)+\xi (t).$$
(15)
That is, the deterministic part of the homodyne photocurrent is proportional to $`x_\mathrm{c}=\sigma _x_\mathrm{c}`$. This should be useful for controlling the dynamics of the state in the $`x`$$`z`$ plane by feedback, as we will consider in Sec. III. Of course, all that really matters here is the relationship between the driving phase and the local oscillator phase, not the absolute phase of either.

The conditioning process referred to above can be made explicit by calculating how the system state changes in response to the measured photocurrent. Assuming that the state at some point in time is pure (which will tend to happen because of the conditioning anyway), its future evolution can be described by the stochastic Schrödinger equation (SSE)

$$d|\psi _\mathrm{c}(t)=\widehat{A}_\mathrm{c}(t)|\psi _\mathrm{c}(t)dt+\widehat{B}_\mathrm{c}(t)|\psi _\mathrm{c}(t)dW(t).$$
(16)
This is an Itô stochastic equation  with a drift term and a diffusion term. The operator for the drift term is

$$\widehat{A}_\mathrm{c}(t)=\frac{\gamma }{2}\left[\sigma ^{}\sigma +\sigma _x_\mathrm{c}(t)\sigma \sigma _x_\mathrm{c}^2(t)/4\right]i\alpha \sigma _y,$$
(17)
while that for the diffusion is

$$\widehat{B}_\mathrm{c}(t)=\sqrt{\gamma }\left[\sigma \sigma _x_\mathrm{c}(t)/2\right].$$
(18)
Both of these operators are conditioned in that they depend on the system average

$$\sigma _x_\mathrm{c}(t)=\psi _\mathrm{c}(t)|\sigma _x|\psi _\mathrm{c}(t).$$
(19)
On average, the system still obeys the master equation (1). This is easiest to see from the stochastic master equation (SME), which allows for impure initial conditions. The SME can be derived from the SSE by constructing

$`d\left(|\psi _\mathrm{c}\psi _\mathrm{c}|\right)`$ $`=`$ $`\left(d|\psi _\mathrm{c}\right)\psi _\mathrm{c}|+|\psi _\mathrm{c}\left(d\psi _\mathrm{c}|\right)`$ (21)
$`+\left(d|\psi _\mathrm{c}\right)\left(d\psi _\mathrm{c}|\right),`$
using the Itô rule (13), and then identifying $`|\psi _\mathrm{c}\psi _\mathrm{c}|`$ with $`\rho _\mathrm{c}`$. The result is

$$d\rho _\mathrm{c}=dt\gamma 𝒟[\sigma ]\rho _\mathrm{c}idt\alpha [\sigma _y,\rho _\mathrm{c}]+dW(t)\sqrt{\gamma }[\sigma ]\rho _\mathrm{c},$$
(22)
where $`[A]BAB+BA^{}\mathrm{Tr}[AB+BA^{}]`$. Although this has been derived assuming pure initial conditions, it is valid for any initital conditions . This is also an Itô equation, which means the evolution for the ensemble average state matrix

$$\rho (t)=\mathrm{E}[\rho _\mathrm{c}(t)]$$
(23)
is found simply by averaging over the photocurrent noise term by using Eq. (14). This procedure yields the original master equation (1) again. The general term for the stochastic conditioned evolution of the system, be it described by a SSE or SME, is a quantum trajectory , and the quantum trajectory is said to unravel the master equation .

## III Feedback with unit-efficiency detection

### A SSE including feedback

We now include feedback onto the amplitude of the driving on the atom, proportional to the homodyne photocurrent, as done by HMH. This is as shown in Fig. 1, where the driving field passes through an electro-optic amplitude modulator controlled by the photocurrent, yielding a field proportional to $`\alpha +\lambda I(t)`$. This means that the feedback can be described by the Hamiltonian

$$H_{\mathrm{fb}}=\lambda \sigma _yI(t).$$
(24)
In this paper we are assuming instantaneous feedback, where the time delay in the feedback loop is negligible.

Since the homodyne photocurrent (11) is defined in terms of system averages and the noise $`dW(t)`$, the SSE including feedback can still be written as an equation of the form (16). The effect of the feedback Hamiltonian can be shown  to change the drift and diffusion operators to

$`\widehat{A}_\mathrm{c}(t)`$ $`=`$ $`{\displaystyle \frac{\gamma }{2}}\left[\sigma ^{}\sigma +\sigma _x_\mathrm{c}(t)\sigma \sigma _x_\mathrm{c}^2(t)/4\right]i\alpha \sigma _y`$ (26)
$`+\lambda \sqrt{\gamma }\left[i\sigma _x_\mathrm{c}(t)\sigma _y2\sigma ^{}\sigma \right]\lambda ^2/2,`$
$`\widehat{B}_\mathrm{c}(t)`$ $`=`$ $`\sqrt{\gamma }\left[\sigma \sigma _x_\mathrm{c}(t)/2\right]i\lambda \sigma _y.`$ (27)
Say we wish to stabilize the pure state with Bloch angle $`\theta `$, as defined in Eqs. (6) and (7), with $`r=1`$ of course. In terms of the ground and excited states, this state is

$$|\theta =\mathrm{cos}\frac{\theta }{2}|e+\mathrm{sin}\frac{\theta }{2}|g.$$
(28)
Now for this state to be stabilized we must have

$$\left[\widehat{A}_\mathrm{c}(t)dt+\widehat{B}_\mathrm{c}(t)dW(t)\right]|\theta |\theta .$$
(29)
We cannot say the left-hand-side should equal zero because a change in the overall phase still leaves the physical state unchanged. However, we can work with this equation, and simplify it by dropping all terms proportional to the identity operator in $`\widehat{A}_\mathrm{c}(t)`$ and $`\widehat{B}_\mathrm{c}(t)`$. We can also demand that it be satisfied for the deterministic and noise terms separately, because $`dW(t)`$ can take any value. This gives the two equations

$`\left(\sqrt{\gamma }\sigma i\lambda \sigma _y\right)|\theta `$ $``$ $`|\theta ,`$ (30)
$`[\gamma (\sigma ^{}\sigma +\mathrm{sin}\theta \sigma )i2\alpha \sigma _y`$ (31)
$`+\lambda \sqrt{\gamma }(i\mathrm{sin}\theta \sigma _y/2\sigma ^{}\sigma )]|\theta `$ $``$ $`|\theta ,`$ (32)
where we have put $`\sigma _x_\mathrm{c}(t)`$ equal to $`\mathrm{sin}\theta `$, its value for the state $`|\theta `$.

Solving the first equation easily yields the condition

$$\lambda =\frac{\sqrt{\gamma }}{2}(1+\mathrm{cos}\theta ).$$
(33)
This is equivalent to the feedback condition derived by HMH, stated as Eq. (35) of Ref. . Substituting this into the second equation gives, after some trigonometric manipulation, the second condition

$$\alpha =\frac{\gamma }{4}\mathrm{sin}\theta \mathrm{cos}\theta .$$
(34)
Again, this agrees with the driving strength of HMH, given as Eq. (44) of Ref. . It is worth emphasizing that the derivation given here is entirely different in detail from that of HMH, and so is an independent verification of their result. These functions are plotted in Fig. 3. Note that there are two points with the same values of both $`\lambda `$ and $`\alpha `$, at $`\theta =\pm \pi /2`$.

### B Stability

The preceding derivation seems to show that any pure state can be stabilized by a suitable choice of driving and feedback. Indeed our derivation proves that that if one prepares a state in exactly the pure state one desires, then the feedback scheme of HMH which we have analyzed will keep the system in that state. However, to discuss stability we need to know what will happen for states which are not initially in the desired state. To deal with this it is much more convenient to use the SME rather than the SSE, as will become apparent.

The SME can be constructed from the SSE in the same way as before. The result is

$`d\rho _\mathrm{c}`$ $`=`$ $`dt\gamma 𝒟[\sigma ]\rho _\mathrm{c}idt\alpha [\sigma _y,\rho _\mathrm{c}]`$ (37)
$`idt\lambda [\sigma _y,\sigma \rho _\mathrm{c}+\rho _\mathrm{c}\sigma ^{}]+dt(\lambda ^2/\gamma )𝒟[\sigma _y]\rho _\mathrm{c}`$
$`+dW(t)[\sqrt{\gamma }\sigma i\lambda \sigma _y]\rho _\mathrm{c}.`$
Also as before, this is an Itô stochastic equation, which means that the ensemble average can be found simply by dropping the stochastic terms. This time, the result is not the original master equation, but rather the feedback-modified master equation

$$\dot{\rho }=i[\alpha \sigma _y,\rho ]+𝒟[\sqrt{\gamma }\sigma i\lambda \sigma _y]\rho \rho .$$
(38)
Here we have put the Liouvillian superoperator $``$ in a manifestly Lindblad form.

Now we have shown already that the pure state $`\rho =|\theta \theta |`$ must be a solution of this master equation, for the appropriate choices of $`\lambda `$ (33) and $`\alpha `$ (34). But for it to be a stable solution we require all of the eigenvalues of the resulting $`_\theta `$ to have a negative real part (except for the one eigenvalue that is always zero, as required for $`_\theta `$ to be norm-preserving). It is not difficult to find these eigenvalues, and in terms of $`\theta `$ they are

$$\gamma /2,\gamma /2,\gamma \mathrm{cos}^2\theta .$$
(39)
Evidently the state $`|\theta `$ will be stable for all $`\theta `$ except $`\theta =\pm \pi /2`$. That is, all states are stable except those on the equator. This is contrary to the conclusion of HMH , based on a linearized stability analysis, that “long term stability of …inverted states \[i.e. states in the upper half plane\] cannot be achieved.” We emphasize that our stability analysis contains no approximations.

In hindsight, the lack of stability for pure states on the equator could have been predicted from the expressions (33) and (34). As discussed above and shown in Fig. 3, the values for driving and feedback for $`\theta =\pi /2`$ are the same as those for $`\theta =\pi /2`$. This means that both $`\rho =|\pi /2\pi /2|`$ and $`\rho =|\pi /2\pi /2|`$ are solutions of $`_\theta \rho =0`$ for $`\theta =\pi /2`$ or $`\pi /2`$. By linearity, any mixture of $`|\pi /2\pi /2|`$ and $`|\pi /2\pi /2|`$ will be a solution also. Hence any deviation away from one of these pure state will not necessarily be suppressed, and the system lacks stability. With random external perturbations, the system will eventually reach an equal mixture of $`|\pi /2\pi /2|`$ and $`|\pi /2\pi /2|`$, which is a state with $`r=0`$ (minimum purity). This is why we have plotted a value of $`r=0`$ in Fig. 3 for $`|\theta |=\pi /2`$. We also plot $`r`$ as a function of $`\theta `$ in Bloch space in Fig 2, giving the locus of states which can be stabilized by feedback. This can compared to the locus of the mixed states which are accessible by driving alone. We will return to the stability issue in the context of stochastic dynamics in Sec. V.

## IV Feedback with non-unit-efficiency detection

We have seen that the stochastic master equation is a very useful representation of a quantum trajectory, as it allows the unconditioned (deterministic) master equation to be derived easily, and this latter equation is all that is required for a completely rigorous stability analysis. The SME is also superior to the SSE in that allows inefficient detection to be described. In a real experiment this has to be taken into account. The effect of non-unit $`\eta `$ on feedback in the present system is of interest both in itself, and because of the extreme nonlinearity of the system dynamics as compared to other quantum optical feedback systems such as considered in Ref. .

As explained in Ref. , the homodyne photocurrent from a detection scheme with efficiency $`\eta `$ is

$$I(t)=\sqrt{\gamma }\sigma _x_\mathrm{c}(t)+\xi (t)/\sqrt{\eta }.$$
(40)
Here we have used a normalization such that the deterministic part does not depend on $`\eta `$. The effect of decreased efficiency is increased noise. This means that we can retain the same feedback Hamiltonian as above (24), without changing the significance of the feedback parameter $`\lambda `$. The SME with $`\eta <1`$, including feedback, is

$`d\rho _\mathrm{c}`$ $`=`$ $`dt\gamma 𝒟[\sigma ]\rho _\mathrm{c}idt\alpha [\sigma _y,\rho _\mathrm{c}]`$ (43)
$`idt\lambda [\sigma _y,\sigma \rho _\mathrm{c}+\rho _\mathrm{c}\sigma ^{}]+dt(\lambda ^2/\gamma \eta )𝒟[\sigma _y]\rho _\mathrm{c}`$
$`+dW(t)[\sqrt{\gamma \eta }\sigma i\lambda \sigma _y/\sqrt{\eta }]\rho _\mathrm{c}`$
The no-feedback SME, analogous to Eq. (22), can be obtained simply by setting $`\lambda =0`$, and was derived in Ref. .

Once again, it is easiest for the moment to just consider the ensemble average evolution by averaging $`dW`$ to zero. The Lindblad form of the resulting master equation is

$$\dot{\rho }=i[\alpha \sigma _y,\rho ]+𝒟[\sqrt{\gamma }\sigma i\lambda \sigma _y]\rho +(\lambda ^2/\eta )𝒟[\sigma _y]\rho .$$
(44)
We do not know a priori what values of $`\lambda `$ and $`\alpha `$ to choose to give the best results with inefficient detection, as the SSE analysis in Sec. III A obviously does not apply. Hence we simply solve for the stationary matrix in terms of $`\alpha `$ and $`\lambda `$. Using the Bloch representation we find

$`x_{\mathrm{ss}}`$ $`=`$ $`4\alpha \eta ^2(\gamma +2\sqrt{\gamma }\lambda )/D,`$ (45)
$`y_{\mathrm{ss}}`$ $`=`$ $`0,`$ (46)
$`z_{\mathrm{ss}}`$ $`=`$ $`\sqrt{\gamma }\eta (\sqrt{\gamma }+2\lambda )(\gamma \eta +4\sqrt{\gamma }\eta \lambda +4\lambda ^2)/D,`$ (47)
where

$`D`$ $`=`$ $`\gamma ^2\eta ^2+6\gamma ^{3/2}\eta ^2\lambda +2\gamma \eta (3+4\eta )\lambda ^2+16\sqrt{\gamma }\eta \lambda ^3`$ (49)
$`+\mathrm{\hspace{0.17em}8}(\alpha ^2\eta ^2+\lambda ^4).`$
The question now arises, what do we mean by “best results” for the feedback system. We cannot hope anymore to produce stable pure states anywhere on the Bloch sphere. However, we can pick a direction $`\theta `$ on the Bloch sphere and ask how close can we get to a pure state? That is, we use the radius $`r`$ in Eq. (6) and Eq. (7) as the quantity to be maximized, for each $`\theta `$. From these two equations we have

$$\mathrm{tan}\theta =x_{\mathrm{ss}}/z_{\mathrm{ss}}.$$
(50)
¿From Eqs. (45) and (47) we can immediately find the desired driving in terms of $`\lambda `$ and $`\theta `$ as

$$\alpha =\frac{\sqrt{\gamma }\eta (\sqrt{\gamma }+2\lambda )(\gamma \eta +4\sqrt{\gamma }\eta \lambda +4\lambda ^2)\mathrm{tan}\theta }{4\eta ^2(\gamma +2\sqrt{\gamma }\lambda )}.$$
(51)
The aim is then, for each $`\theta `$, to find the feedback $`\lambda `$ which maximizes

$$r^2=x_{\mathrm{ss}}^2+z_{\mathrm{ss}}^2.$$
(52)
This was done numerically using matlab.

The results of our calculations are shown in Fig. 4, where we plot the locus in Bloch space of the best (most pure) stationary states which can be achieved by feedback from non-unit-efficiency detection. We use a variety of values of $`\eta `$. A number of points are worth noting. First, and most obviously, the degree of purity (measured by the $`r`$, the distance from the origin) decreases with $`\eta `$. Second, the gap at the equator for $`\eta =1`$ quickly widens, so that the purity of the best states with $`\theta `$ close to $`\pi /2`$ is small. Third, the purity of the best states in the upper half of the Bloch sphere is affected much more by loss of detection efficiency than those in the lower half. Fourth, in the limit $`\eta =0`$, the best solutions correspond to the no-feedback solutions shown in Fig. 2. This is not surprising, since with $`\eta =0`$ the photocurrent contains no information about the system (as the noise is infinitely large) and hence there is no point doing feedback. Since the stationary states with no feedback are confined to the lower half of the Bloch sphere, this explains why the best states with feedback in the lower half are less affected as $`\eta `$ decreases than those in the upper half.

For the particular value $`\eta =0.8`$ we plot in Fig. 5 the values of $`\alpha `$ and $`\lambda `$ (as well as purity, quantified as $`r^2`$) versus $`\theta `$. By comparing this plot with Fig. 3 one obtains some idea of the effect of inefficient detection. A number of features remain the same. First, $`\alpha `$ is antisymmetric in $`\theta `$, while $`\lambda `$ is symmetric. Recalling that the deterministic part of the feedback is proportional to $`\lambda \sigma _x_\mathrm{c}`$, the feedback itself is actually antisymmetric as well as the driving. Second, the magnitude of the feedback is zero for $`|\theta |=\pi `$ (the ground state) and increases monotonically to a maximum of $`\sqrt{\gamma }`$ at $`\theta =0`$ (in the direction of the excited state). Third, the driving is zero at the ground state and at $`\theta =0`$, and also changes sign as one passes through the equatorial place. The most obvious difference between the parameters for $`\eta =1`$ and those for $`\eta =0.8`$ parameter is that the latter have a discontinuity at $`|\theta |=\pi /2`$. The feedback parameter $`\lambda `$ jumps as one crosses the equatorial plane, while the driving $`\alpha `$ asymptotes to $`+\mathrm{}`$ on one side and $`\mathrm{}`$ on the other. These extreme variations in the driving do not prevent the best purity from approaching zero in the equatorial plane.

## V Stochastic Dynamics

### A Stochastic Bloch Equations

So far we have considered the stochastic conditioned dynamics for the system state in order to derive the parameters $`\lambda `$ and $`\alpha `$ such that for $`\eta =1`$ those dynamics are banished in the steady state. In this section we will consider them in more detail, highlighting the difference between the $`\eta =1`$ case and the $`\eta <1`$ case, and also looking in more detail at the special case of $`|\theta |=\pi /2`$. The most convenient way to treat the stochastic dynamics in general is through the stochastic Bloch equations (SBE). These are simply the stochastic equations for the conditioned Bloch vector, defined by

$$\rho _\mathrm{c}=\frac{1}{2}\left(I+x_\mathrm{c}\sigma _x+y_\mathrm{c}\sigma _y+z_\mathrm{c}\sigma _z\right).$$
(53)
¿From the SME including feedback (43), we find

$`\left(\begin{array}{c}dx_\mathrm{c}\\ dy_\mathrm{c}\\ dz_\mathrm{c}\end{array}\right)`$ $`=`$ $`dt\left(\begin{array}{ccc}\gamma /22\kappa & 0& 2\alpha \\ 0& \gamma /2& 0\\ 2\alpha & 0& \gamma 2\kappa \end{array}\right)\left(\begin{array}{c}x_\mathrm{c}\\ y_\mathrm{c}\\ z_\mathrm{c}\end{array}\right)`$ (71)
$`dt\left(\begin{array}{c}0\\ 0\\ 2\lambda \sqrt{\gamma }+\gamma \end{array}\right)+dW(t)\times `$
$`\left(\begin{array}{c}\sqrt{\gamma \eta }x_\mathrm{c}^2+(\sqrt{\gamma \eta }+2\lambda /\sqrt{\eta })z_\mathrm{c}+\sqrt{\gamma \eta }\\ \sqrt{\gamma \eta }x_\mathrm{c}y_\mathrm{c}\\ (\sqrt{\gamma \eta }+2\lambda /\sqrt{\eta })x_\mathrm{c}\sqrt{\gamma \eta }x_\mathrm{c}z_\mathrm{c}\end{array}\right),`$
where $`\kappa =\lambda ^2/\eta +\lambda \sqrt{\gamma }`$. If we ignore the final (noise) term, we get the Bloch equations from the master equation (44).

### B unit-Efficiency

In the case $`\eta =1`$, considered in Sec. III, both the deterministic and stochastic dynamics disappeared in the steady state for the appropriate choice of $`\alpha `$ and $`\lambda `$. Because the stationary solution of the SSE was a unique pure state, that was necessarily also the stationary solution of the master equation found by averaging over the noise in the equivalent SME. Thus there was no distinction between the unconditioned and conditioned states. There are two exceptions to this lack of distinction. The first is in the transients, before the system reaches its steady state. The second is for the special case $`|\theta |=\pi /2`$. In this subsection we investigate these exceptions.

First, the transient behaviour was simulated using the SBE with $`\eta =1`$. We chose the initial state to be the ground state, and evolved the system stochastically from $`t=0`$ to $`t=10\gamma ^1`$. With this choice of initial condition, $`y_\mathrm{c}=0`$ for all time. We verified that in each trajectory $`x_\mathrm{c}^2+z_\mathrm{c}^2=1`$ to a good approximation (indicating a pure state), but that the ensemble averages over many trajectories

$$x=\mathrm{E}[x_\mathrm{c}],z=\mathrm{E}[z_\mathrm{c}]$$
(73)
obey the deterministic Bloch equations. A typical trajectory for $`\theta =\pi /6`$ is shown in Fig. 6(a). We see that the initial evolution is very erratic, but that on a time scale of a few $`\gamma ^1`$ the system relaxes towards a steady state which is pure and stationary. By $`t=10\gamma ^1`$ the system is locked in a stable pure state for all intents. We have also illustrated another typical trajectory in Bloch space, in Fig. 7.

It is easy to verify that by putting $`\eta =1`$ and

$$(x_\mathrm{c},y_\mathrm{c},z_\mathrm{c})=(\mathrm{sin}\theta ,0,\mathrm{cos}\theta )$$
(74)
in the right-hand-side of the SBE (LABEL:SBE1), we obtain complete cancellation. If we wish, we can follow HMH and separate the noise term into the contribution from feedback (proportional to $`\lambda `$) and that contribution present even without feedback (the rest). We interpret the latter stochasticity as being due to the quantum measurement we are making, with its underlying probabilistic nature. Obviously the fluctuation due to measurement is canceled by the feedback, as HMH point out. It is equally important that the deterministic dynamics are also canceled at this point.

The story for the special case $`\theta =\pi /2`$ is quite different. For this case the SBEs are

$`\left(\begin{array}{c}dx_\mathrm{c}\\ dy_\mathrm{c}\\ dz_\mathrm{c}\end{array}\right)`$ $`=`$ $`dt\left(\begin{array}{ccc}0& 0& 0\\ 0& \gamma /2& 0\\ 0& 0& \gamma /2\end{array}\right)\left(\begin{array}{c}x_\mathrm{c}\\ y_\mathrm{c}\\ z_\mathrm{c}\end{array}\right)`$ (88)
$`+\sqrt{\gamma }dW(t)\left(\begin{array}{c}1x_\mathrm{c}^2\\ x_\mathrm{c}y_\mathrm{c}\\ x_\mathrm{c}z_\mathrm{c}\end{array}\right).`$
Here the three eigenvalues in Eq. (39) are clearly evident. Both $`z_c`$ and $`y_\mathrm{c}`$ will decay to zero (as required for $`\theta =\pi /2`$), and their noise terms vanish at that point. By contrast, the equation for $`x_\mathrm{c}`$ is independent of the others, and is purely stochastic:

$$dx_\mathrm{c}=\sqrt{\gamma }dW(t)(1x_\mathrm{c}^2).$$
(89)
Clearly the equatorial pure states with $`x_\mathrm{c}=\pm 1`$ are stationary solutions to this problem. Also, the system will tend to one of these states. We can see this by calculating

$$d\mathrm{E}[x_\mathrm{c}^2]=\gamma dt\mathrm{E}\left[(1x_\mathrm{c}^2)^2\right],$$
(90)
which is always positive. That is, on average $`x_\mathrm{c}^2`$ always increases. But it is also clear that $`x_\mathrm{c}`$ has no preference to go to either of these states. Hence they are not stable. The ensemble average $`x`$ is unchanging under this evolution. Thus a perturbation which moves the state from $`x_\mathrm{c}=1`$ to $`x_\mathrm{c}=1ϵ`$ say, will result in a proportion $`ϵ/2`$ of the states ending up at $`x_\mathrm{c}=1`$, and a proportion $`1ϵ/2`$ ending up at $`x_\mathrm{c}=1`$.

We have illustrated these features by showing two typical trajectories in Fig. 6(b). Once again, the initial evolution is highly erratic, but the system reaches a fixed point on a time scale of a few $`\gamma ^1`$. However, with the same initial condition (the ground state), one trajectory ends up at $`x_\mathrm{c}=1`$ and the other at $`x_\mathrm{c}=1`$.

### C Non-unit-efficiency

In the unit-efficiency case the stationary solution of the master equation is (except for $`|\theta |=\pi /2`$), a pure state. This is very special in that in means that every unraveling of the master equation as a SSE or SME must end up in this same pure state also. For non-unit-efficiency we have found the most pure stable state for each $`\theta `$. In this case we must use a SME to unravel the master equation, since the conditioned state will not be pure in general, because of the inefficient detection. Since the deterministic steady state is not pure (except for $`|\theta |=\pi `$), the quantum trajectories need not end up in this state. Instead, the quantum state in an individual trajectory may continue to evolve stochastically even when the system is in steady state, and the equivalence to the deterministic evolution may hold only on average. On the other hand, it is also possible that the quantum trajectories do all end up in the deterministic steady state, since we expect the conditioned state to be mixed anyway.

It turns out that with the optimal values of $`\alpha `$ and $`\lambda `$ defined in Sec. IV, the actual behaviour is the first option described above. That is, the system state continues to vary stochastically in the long time limit, but is constrained so that the time-averaged state is equal to the solution of the deterministic master equation. We show this in Bloch space Fig. 7 for $`\eta =0.8`$ and $`\theta =\pi /6`$. We see that the amount of randomness in the system state in the long-time limit is quite large even for fairly high $`\eta `$.

This result suggests another question: would a different choice for $`\lambda `$ be able to reduce, or even eliminate, the randomness in the steady-state quantum trajectory, even though it would necessarily be at the expense of the purity of the deterministic stationary solution. \[Recall that for a given $`\lambda `$, $`\alpha `$ is still necessarily fixed by Eq. (51).\] To test this idea we tried choosing $`\lambda `$ based not on maximizing $`r^2`$ as in Eq. (52), but on minimizing

$$N_\theta (\lambda )=\left|\left(\begin{array}{c}\sqrt{\gamma \eta }x_{\mathrm{ss}}^2+(\sqrt{\gamma \eta }+2\lambda /\sqrt{\eta })z_{\mathrm{ss}}+\sqrt{\gamma \eta }\\ \sqrt{\gamma \eta }x_{\mathrm{ss}}y_{\mathrm{ss}}\\ (\sqrt{\gamma \eta }+2\lambda /\sqrt{\eta })x_{\mathrm{ss}}\sqrt{\gamma \eta }x_{\mathrm{ss}}z_{\mathrm{ss}}\end{array}\right)\right|^2.$$
(91)
That is, we minimize the noise terms in the SBE Eq. (LABEL:SBE1). Note that we have replaced the conditioned Bloch variables $`x_\mathrm{c}`$ etc. with the deterministic stationary solutions $`x_{\mathrm{ss}}`$ etc., and that the dependence of these stationary solutions on $`\alpha `$ and $`\lambda `$ add a further, implicit, dependence on $`\lambda `$ to $`N_\theta (\lambda )`$. This is a sensible procedure if the aim is realizable, and the noise in the solutions is reduced or eliminated so that the conditioned states are approximately or exactly equal to the deterministic stationary solution.

It turns out that this procedure cannot significantly reduce the amount of steady-state randomness in the quantum trajectories below that resulting from minimizing the deterministic stationary purity. In fact, for all values of $`\eta `$ we considered, the variation of $`\lambda `$ (as a function of $`\theta `$) based on minimizing the noise was indistinguishable by eye from that based on maximizing the purity. This is not too surprising, but could not have been predicted a priori.

## VI Conclusion

We have given a rigorous analysis of the anti-decoherence feedback scheme proposed by Hofmann, Mahler and Hess . They proposed modulating the driving of a two-level atom using the instantaneous homodyne photocurrent, in order to stabilize the atom in an arbitrary known pure state. We have shown that, for detection efficiency $`\eta =1`$, the pure states thus produced are stable. This is contrary to the conclusion of HMH, that only pure states in the lower half of the Bloch sphere would be stable. The one exception we found is for pure states on the equator. Although they are fixed points of the dynamics, they are not stable. A small perturbation away from one fixed pure state leads to a proportionally small fraction of the ensemble ending up in the diametrically opposite pure state.

It is nevertheless possible to obtain an asymmetry between the upper and lower halves of the Bloch sphere, reminiscent of the conclusion of HMH, if one considers detection efficiencies less than one. In this case, it is no longer possible to stabilize the system in a given pure state, so we choose the feedback and driving so that the solution of the master equation (including feedback) is as close as possible to a given pure state. We find that the purity (which measures this closeness) of states thus produced decays to zero as $`\eta `$ decreases to zero, for states in the upper half of the Bloch sphere. By contrast, those in the lower half do not decay to zero. This is readily understandable since in the absence of feedback (which is the situation which must prevail when the detection efficiency goes to zero), the master equation with driving alone has stationary solutions in the lower-half plane with non-zero purity. The purity decays most rapidly with $`\eta `$ for states near the equator, which is unsurprising given the instability of states on the equator even for $`\eta =1`$.

In the non-unit-efficiency case, the state of the system conditioned on the homodyne measurement results continues to evolve stochastically even in the long-time limit, where the ensemble average evolution has reached the desired most-pure state. Moreover, it seems that any other choice of driving and feedback will result in more, not less, randomness in the steady-state quantum trajectory.

Our results are significant in a number of ways. First, they show the power of the quantum trajectory and master equation techniques developed in Refs. . Those techniques were particularly useful for illuminating subtle questions regarding the stability of pure states, and for treating inefficient detection. Second, the physical system (the two-level atom) may one day find application as a quantum bit in quantum information technology . In that eventuality, the ability to stabilize the atom against dissipation in an arbitrary (known) pure state may be useful. Third, the system is a simple but non-trivial example of quantum feedback in a nonlinear system (the two-level atom). Thus the effectiveness of feedback, and in particular the influence of non-unit- efficiency detection on this effectiveness, is of interest for what it may tell us about other more complicated nonlinear systems.

In this last context, it would be of interest to also consider the effect of non-Markovian feedback; that is, feedback with a time delay or non-flat loop response function. This is much more difficult to treat than Markovian feedback because the Lindblad master equations derived in Refs.  do not apply. Analytical solutions for non-Markovian feedback are possible for linear systems . For a nonlinear system like the two-level atom, numerical simulations, or novel analytical approaches, are necessary. This is an issue we plan to explore in future work.

## Acknowledgments

This work has been supported by the Australian Research Council, the University of Queensland, and the Department of Employment, Education and Training, Australia. We would like to acknowledge useful discussions with Gerard Milburn.


format: text

text:
# The mass distribution in the innermost regions of Spiral Galaxies.

## 1 Introduction

There is increasing direct evidence that, at the centers of bulge-dominated galaxies there reside Massive Dark Objects (MDO), which are probably the remnants of the engines that once powered the QSO phenomenon (Ho, 1998; Kormendy & Richstone, 1995). In fact, virtually every hot galaxy hosts a MDO/BH with a mass ranging from $`10^8M_{}`$ to $`2\times 10^{10}M_{}`$, similar to those related to the QSO phenomenon. For disk galaxies, the situation is different and much more uncertain. A direct determination of the central mass has been obtained only in very few cases which include our own Galaxy where a black hole with a mass $`2\times 10^6M_{}`$ has been discovered (Ghez et al.1999, Genzel 1998; see also Salucci et al 1999 for other few cases). Remarkably, these masses do not exceed $`10^7M_{}`$; however, the lack of detections of very massive objects ($`M_{MDO}>10^8M_{}`$) cannot be ascribed to observational biases. In fact, the rotation curves (RC), in great number available down to $`r_{in}100pc`$, could easily expose central bodies with masses of the order of $`10^8M_{}`$, given that usually $`V(r_{in})10km/s`$, or equivalently, the stellar mass inside $`100pc`$ barely reaches $`10^7M_{}`$. On the other hand, the same RC analysis determining the MDO mass, obtain also the mass distribution of the innermost regions. This is particularly important in the case of low luminosity DM-dominated spirals (see Persic and Salucci, 1990) for which it is generally difficult to disentangle the disk component from the whole bulk of gravitating matter. In any cases, let us stress that to proper model the region where the luminous matter dominates is indispensable to infer the structural properties of dark matter.

The aim of this work is to study the innermost kinematics of spirals in order to a) investigate the distribution of matter in relation to that of light b) set stringent upper limits on the MDO masses. The plan of the paper is the following: in section 2 we analyse a sample drawn from the 967 rotation curves (RC’s) of Persic and Salucci (1995), in section 3 we derive and discuss the disk properties and in sections 4 and 5 we derive upper limits on the MDO/BH masses and discuss the results. In this paper we assume $`H_0=75km/sMpc^1`$, and take, as the reference magnitude in the B-band, $`M_{}=20.5`$, which translates to $`M_{}21.9`$ in the I-band (e.g. Rhee, 1997). All luminosities correspond to the I-band unless otherwise specified.

## 2 Inner Rotation Curves of Late Types Spirals and mass modeling

Recently, about a thousand rotation curves of spiral galaxies (PS95), tracing the kinematics inside the central kpc, have been available. Among these, about one hundred have at least one measurement $`<200pc`$ and several $`<5001000pc`$. For these objects the detection of a MDO of $`10^8M_{}`$ or more, if present, is guaranteed and so is the determination of the stellar disk mass, provided that in this region the stellar disk is the major mass component.

Therefore, we select from PS95 the largest sample of rotation curves equally distributed in magnitude interval and of sufficient high-quality for our purposes, by setting the following criteria: i) each RC has at least 4 measurements inside a radius of $`340(1\frac{(M_I+15)}{9})pc(350500)pc`$, ii) for each RC, the innermost data is situated at a radius $`r_{in}<100(1\frac{(M_I+15)}{18})(100150)pc`$. Let us notice that to restrict the criteria reduces the number of RC’s without improving their already high quality that however rapidly decreases as they are relaxed. The rotation curves resulting from the selection are generally smooth, axisymmetric and with negligible non circular perturbations. It is evident that there is a very good agreement between the rotation fields on the receding and the approaching arm (see PS95). The sample (Sample B in Salucci et al., 2000) has 83 objects, well mixed in luminosity and almost equally distributed between Sb-Sc and Sd-Im Hubble types.

We aim to reproduce the RCs by a mass model featuring: (i) a Freeman disk of length-scale $`R_D`$, (derived in PS95) and mass $`M_D`$, which contributes to the circular velocity as:

$$V_d^2(R)=\frac{1}{2}GM_D/R_Dx^2(I_0K_0I_1K_1)_{x/2}$$
$`(1)`$
where $`I_n,K_n`$ are modified Bessel functions and $`x=R/R_D`$, and (ii) a MDO which contributes to $`V(R)`$ as:

$$V_{MDO}^2(R)=GM_{MDO}/R$$
$`(2a)`$
where the MDO mass is:

$$M_{MDO}=fG^1V^2(r_{in})r_{in}$$
$`(2b)`$
with $`0<f1`$.

In Fig.(1) we compare the velocity data with the best-fitting mass models in the cases of 1) $`f=0`$ (only disk, OD) mass model and of 2) $`f>0`$ (disk+ black hole) mass model (in practice: $`f1`$). As result, we find no central-body-dominated (CBD) rotation curve, i.e. no RC shows the Keplerian fall-off expected for a CBD RC. On the contrary, rotation curves are strikingly close to those predicted by the OD model, i.e. by a self-gravitating exponential thin disk with radially constant mass-to-light ratio. The results of the mass modeling are: in 57 objects the OD velocity curve accounts for the rotation data in a excellent way and it lies within the data error bars (for $`R<R_{IBD}`$ see below). In 14 cases111349-G6, 38-G23, 35-G18, 157-G20, 436-G3, 468-G23, 547-G24, 550-G7, N4348, 358G-63, 582-G12, 547-G31, I96099, N4705the luminous matter accounts for the rotation curve in a satisfactory way: data and model differ utmost by $`2\sigma `$ and no evidence for an additional component emerges. In 9 cases222347-G84, 328-G43, 249-G16, 84-G10, N755, 79-G14, M222025, N7218, 249-G16 the OD fits are reasonable expecially if we consider that the corresponding RC’s have non-negligible internal dispersion and/or some asymmetry. The low-luminosity galaxy 545-G3 is dominated by dark matter at any radius. Finally, in 2 cases; 346-G26, 410-G19, the OD model reproduces the data with some difficulty.

A simple inspection of the results shows that in every galaxy only for $`RR_{IBD}`$333Inner Baryon Dominancethe dark component begins to significantly contribute to the total gravitational potential. We then confirm the picture of Salucci & Persic (1999) and Salucci etal, (2000) according to which, the luminous matter dominates an innermost region of spirals of size $`R_{IBD}`$. This truly baryonic scale turns out to be a function of galaxy luminosity by ranging from $`0.5R_D0.5kpc`$ to $`2R_D30kpc`$ along the luminosity sequence Such a transition can be seen, e.g. in the RC of 55-G4: OD curve reproduces the data only out to $`R=R_{IBD}0.5kpc`$. More in general, the best-fit models (see Fig 1) clearly show that, with the exception of the least luminous galaxy, at $`R0.5kpc`$, the DM contribution to the RC lies below the detection threshold. On the other hand, in order to become the major component at a $`23`$ disk scale-lenghts, the dark mass, negligible at $`0.5kpc`$ must strongly increase with radius.

We conclude by stating that the OD mass models reproduce the innermost kpc of the rotation curves extremely well: on this scale, there is it no hint of a dark component.

## 3 Disk properties

Disk masses are derived from eq (1); Although the uncertainties on the disk mass-to-light ratio values must include model-fit uncertainties, photometry errors, and uncertainties on the assumed inclination and distance (the latter present also for Tully-Fisher distances). The uncertainty budget can be roughly estimated as: $`\delta M_D/L_I0.4M_D/L_I`$.

The present sample of DM-free mass models is a most suitable one for investigating the disk masses in relation with the photometric properties of galaxies. In detail, this is done with a sub-sample of Sample B comprising 28 objects for which $`I`$-band luminosities and $`BI`$ colors are available. In Table 1 we report best-fit disk masses, color, mass-to-light ratios $`M_D/L_I`$ and luminosity for these objects. B-I color are good indicators of stellar mass-to-light ratio (Vazdekis et al 1996). In Fig.2 we the show $`𝒴_I`$ vs. $`BI`$ colors alongside with with the same relation predicted by population synthesis models. We find

$$𝒴_I0.63\times \left(\frac{1.6L_B}{L_I}\right)^2$$
$`(3a)`$
in good agreement with Vadzekis et al

Finally, we compare the kinematical disk masses $`M_D`$ with those derived from the galaxy specto-photometry: $`M_{phot}=𝒴_I(BI)L_I`$ (let us remind that $`𝒴_I`$ is the average stellar mass-to-light ratio for a stellar population of color $`BI`$. As a quite general outcome of the past history of galaxies, Vazdekis et al. 1996 found

$$𝒴_I=b(BI)+c$$
$`(3b)`$
where $`b0.8`$ and the value of the constant $`c`$ is irrelevant. We find:

$$M_D=(0.96\pm 0.1)M_{phot}$$
$`(4)`$
(28 d.o.f., see Fig 4), i.e. dynamical and photometric mass estimates statistically coincide. Notice that the disk mass and I-band luminosity are not directly proportional:

$$logM_D=(0.79\pm 0.1)logL_I+const$$
$`(5)`$
with the log slope being significantly different from unity. Then, the existence of a relationship between the maximum rotation velocity and the galaxy luminosity is far from being a trivial one, in that the $`I`$-band luminosity is not a straightforward measure of the stellar mass. It is worth to note that in the 3-D space defined by the ($`logL_I`$, $`log(M_D/L_I)`$, $`BI`$) coordinate vector, spirals are not randomly distributed, but occupy a very thin plane (see Fig. 2).

## 4 MDO/BH masses. Upper limits

The very good spatial resolution of the RC’s of our sample plays a crucial role in the determination of strict and reliable MDO upper limits. In fact, since at small radii $`VR`$ and the mass inside $`r_{in}`$ scales as $`r_{in}^3`$, a spatial resolution of $`r_{in}100pc`$ yields to a mass resolution of $`10^6M_{}`$, given that $`V(r_{in})10km/s`$. In detail, the ratio of the rotational velocities $`(V_{MDO}/V_d)^2(M_{MDO}/r_{in})^3`$, implies that a higher spatial resolutions lead to stronger limits on $`M_{MDO}`$ (or easier detection).

We now force the presence of a central MDO/BH by adopting the maximum possible MDO mass compatible with the rotation curve. We determine the upper mass limits by assuming $`f=1`$ or $`f=f_{max}`$ with $`f=f_{max}`$ the value that brings the model RC $`1\sigma `$ higher than the innermost data. Note that for $`f=1`$ all the mass inside $`r_{in}`$ is in the black hole. We typically find: $`f_{max}=0.70.8`$ and that, in spite of having an additional parameter, the disk+black hole model performs significantly worse than the OD model.

Notice that only in a few cases $`f_{max}0.2`$, furthermore, no result changes if we assume$`f_{max}=1`$ for all the objects. $`M_{BH}`$ is then obtained by substituting in equation (2b) the corresponding values of $`f`$. These are shown in Table 1: it is evident that typically $`M_{MDO}<<10^8M_{}`$ especially at low luminosity. It is clear that these limits are short of the mass residing in a accreting BH which is required to power high redshifts quasars, thus, spirals do not host QSO remnants (see also Salucci etal 2000). More specifically, if central MDO’s of masses $`2M_{MDO}`$ (and still $`<10^8M_{}`$) were actually present at the center of the late type spirals of our sample, they would have affected the inner kpc of the available RC’s more strongly than MDO’s of $`10^9M_{}`$ affect the inner kinematics of ellipticals (Magorrian et al., 1998).

In Fig. 5 we show the relationship $`M_{BH}`$ vs  the stellar disk mass $`M_D`$, and, as a reference, we plot the extrapolation to spirals of the elliptical’s stellar mass-MDO/BH mass relationship $`M_{MDO}=0.005\times M_{stars}`$ (Magorrian et al.,1998, Kormendy & Richstone 1995 ). We immediately realize, that in spirals, the MDO masses, if different from zero, are however much smaller than those detected at the centers of ellipticals of same stellar luminosity; whether this is related to their less prominent bulge content or does reflect also a morphological mass segregation, will be considered elsewhere (Salucci et al 2000). Notice that in the same paper we perform a detailed investigation on possible biases occurring in the present estimate of the MDO upper limits.

## 5 DISCUSSION

Accurate mass models of the innermost kpc of 83 spirals derived from high spatial resolution RC’s, reveal that, inside this region, the luminous matter fully accounts for the observed kinematics: a stellar disk of constant mass-to-light ratio (in some cases in conjunction with a spheroid) is virtually the only mass component. In this region, with the exception of very high luminosity galaxy with complex dynamics, we find

$$V^2(R)\mathrm{𝑙𝑖𝑔ℎ𝑡}(<R)/R$$
$`(6)`$
which is indicative of the absence of a sizable dark component. Therefore, the only possibility for the presence of a dark halo in the central regions of spirals is that it conceals itself below the detection threshold.

More in detail, inside the innermost kpc or so, the contribution to the gravitational potential from a dark component cannot exceed 10%; this is in disagreement with claims and scenarios in which a dark halo has a major role at any radii. Moreover, in the same objects studied here, the DM, undetected at $`R<1/3R_D`$, begins to dominate the mass distribution at $`R3R_D`$ (PSS). Then, the dark mass inside a giv en radius $`R`$ must increase with $`R`$ very steeply, e.g. as in the case of a constant density distribution.

In a DM-free environment, we have been able to estimate the disk masses also for objects in which the global influence of the dark matter is relevant and to relate disk colors and masses across a large range of galaxy luminosities. The next step will be to pursue a coordinated kinematical and photometric study to tackle open cosmological issues.

## 6 Appendix

In the appendix we produce the two Tables of galaxy and MDO/BH properties.


text:
# Error correction for mutually interacting qubits

## I Introduction and model

In previous papers  I have discussed how unwanted mutual interactions between qubits can be a source of error in a quantum computer (or quantum memory register) and studied some of the characteristic features of this type of error. It was found in  that, if the characteristic interaction energy is $`\mathrm{}\delta `$, the fidelity $``$ of the register decreases with time $`t`$ and with the total number of qubits $`N`$ as

$$e^{N(\delta t)^2}$$
(1)
Equation (1) suggests that $`(\delta t)^2`$ can be interpreted as a sort of failure probability per qubit. In  it was shown that the general scaling (1) with $`t`$ holds even when the register is being acted on by “gates” such as one might use in the course of a quantum computation: that is, the error accumulates quasi-coherently in time, and there is no “quantum Zeno effect” associated with normal gate operation. (This suggests, in particular, that it is not necessary to simulate a full computation to estimate the magnitude of this effect; instead, it should be enough to calculate the survival probability of an appropriate initial state of the register under the action of the appropriate Hamiltonian .)

Mutual interactions are Hamiltonian and hence, in principle, reversible, and there exists, in fact, a variety of special techniques for eliminating, or greatly reducing, these errors: in particular, techniques similar to, and inspired by, the use of refocusing pulses in magnetic resonance . Nonetheless, it is worthwhile to look at the question of how well these errors could be eliminated by the use of ordinary error-correction methods based on quantum error-correcting codes , and in particular on concatenated coding, if only because one would like to be able to make an informed choice between various available methods. Moreover, although pulse methods could in principle be used to deal with certain kinds of environmental noise, it seems very likely that, in large-scale implementations of quantum computers, error-correcting codes will have to be used, in any event, to deal with at least some of these random environmental errors. Assuming, then, that these codes would already be in place, it is natural to ask how much one could get from them.

While there is no question that the present theory of quantum error-correction codes covers mutually interacting qubits as a special case, this is not a case that is often explicitly discussed, and, in particular, it is not immediately clear whether, in order to obtain a threshold estimate for these kinds of errors one should add them coherently (adding error amplitudes) or incoherently (adding error probabilities). This is one of the issues explored in this paper, for the two best-known quantum error-correcting codes: the 5-qubit and 7-qubit codes, and their concatenated forms. The conclusion, at least for the kind of interaction Hamiltonian assumed, is that these codes could fully handle the errors arising from mutual interactions, almost certainly at no extra cost in terms of resources needed (e.g., depth of concatenation). This is an encouraging result which may not have been immediately obvious before.

The interaction to be considered here is one that results in pure phase shifts, conditioned on the relative state of neighboring qubits, according to a Hamiltonian

$$H=\underset{n=1}{\overset{N1}{}}\mathrm{}\delta \sigma _{zn}\sigma _{z,n+1}$$
(2)
Examples of physical systems leading to interactions of this form have been discussed in previous publications ; in  it has also been shown that, for instance, changing $`\sigma _z`$ to $`\sigma _x`$ (bit errors instead or phase errors) does not modify substantially the scaling (1).

In order to simplify the calculation, the (unrealistic) assumption is made that the qubits are divided into physical blocks, each of length $`N_c`$, which do not interact with each other; that is, the interaction Hamiltonian is applied separately to each such block. The blocks are intended to represent the “logical” qubits in a computer encoded for error correction. When the code is concatenated, the decoupling assumption is applied at the lowest (deepest) level: For instance, for a code concatenated once, a logical qubit is made of $`N_c\times N_c`$ physical qubits, and it is assumed that the $`N_c`$ lowest-level blocks do not interact with each other. While unrealistic, it is expected that this approximation should not modify the results very much, based on the numerical calculations made using a similar “decoupling of registers” approximation in ref. .

In the following Section, the action of this Hamiltonian on a register encoded for error correction is analyzed, and then in Section III the analysis is extended to a concatenated code. The goal is to derive an expression for the “failure probability” per qubit, that is to say, the probability that the evolution under the Hamiltonian (2) may lead to a state which cannot be recovered by the error-correction protocol. Section II shows a few general features of this failure probability, such as its dependence on the initial state considered, and the way some error amplitudes add coherently and some do not, depending on the code considered. These results pave the way for the more complicated situation discussed in Section III. It is to be noted that the failure probability is calculated here for arbitrarily large interaction strength, for up to two-level-deep concatenated codes.

## II Results for simple distance 3 codes

Distance 3 codes can correct any of the three basic one-qubit errors represented by the three Pauli matrices $`\sigma _x`$, $`\sigma _y`$, and $`\sigma _z`$. They will typically fail whenever two errors happen simultaneously in the same logical block.

The Hamiltonian (2) causes two phase errors in each physical block every time it is applied, but the actual damage it does to an encoded state turns out to depend on the code used, and on the state it acts upon. Two examples are considered in this section.

The evolution operator resulting from (2) is, for a single logical block,

$$U(t)=\underset{n=1}{\overset{N_c1}{}}\left(\mathrm{cos}\delta t+i\sigma _{zn}\sigma _{z,n+1}\mathrm{sin}\delta t\right)$$
(3)
since the Pauli matrices for different qubits commute. Thus, a simple perturbative estimate of the failure probability would be as follows: if it is assumed that each term of the form $`\sigma _{zn}\sigma _{z,n+1}`$ leads to a distinct (orthogonal) uncorrectible error, there are $`N_c1`$ such terms in (3), and the sum of the squares of their amplitudes is

$`P_{pert}(t)`$ $`=(N_c1)\left(\mathrm{cos}^2\delta t\right)^{N_c1}\mathrm{sin}^2\delta t`$ (5)
$`(N_c1)(\delta t)^2`$
for sufficiently small interaction strength and time $`\delta t`$.

As it turns out, this expression is not, in general, correct, even for small $`\delta t`$, because of two assumptions: that all terms of the form $`\sigma _{zn}\sigma _{z,n+1}`$ lead invariably to fatal (uncorrectible) errors, and that the different errors are orthogonal. This is explored in detail in what follows.

### A The 7-qubit code

Consider the 7-qubit code introduced by Steane :

$`|0_L={\displaystyle \frac{1}{\sqrt{8}}}(`$ $`|0000000+|0001111+|0110011`$ (9)
$`+|0111100+|1010101+|1011010`$
$`+|1100110+|1101001)`$
$`|1_L={\displaystyle \frac{1}{\sqrt{8}}}(`$ $`|1111111+|1110000+|1001100`$ (12)
$`+|1000011+|0101010+|0100101`$
$`+|0011001+|0010110)`$
Because the code is non-degenerate, all 7 states of the form $`\sigma _{zn}|0_L`$ (with $`n=1,\mathrm{},7`$) are orthogonal to each other and to $`|0_L`$. Altogether, then, these eight states form an orthogonal basis of the space spanned by the eight kets appearing in (9). But then, since any combination of products of $`\sigma _{zn}`$ acting on $`|0_L`$ does nothing but to change the relative signs in the superposition of these eight kets, it follows that it can always be written as a superposition of “one-error” states:

$$U(t)|0_L=C_0|0_L+\underset{n=1}{\overset{7}{}}\sigma _{zn}|0_L$$
(13)
and this is always correctible. A similar reasoning applies to the state $`|1_L`$, which means that, if the initial state of the system is either $`|0_L`$ or $`|1_L`$, the evolution operator (3) does not introduce any fatal errors at all.

On the other hand, the initial state of the system is much more likely to be a coherent superposition of $`|0_L`$ and $`|1_L`$, and for such states the situation is different, because, although error correction will restore $`|0_L`$ and $`|1_L`$ separately, it will also introduce a phase difference between them. The following useful relations can be easily established from the explicit forms (II A):

$`\sigma _{z1}\sigma _{z2}|0_L=\sigma _{z5}\sigma _{z6}|0_L`$ $`=\sigma _{z3}|0_L`$ (15)
$`\sigma _{z2}\sigma _{z3}|0_L=\sigma _{z4}\sigma _{z5}|0_L=\sigma _{z6}\sigma _{z7}|0_L`$ $`=\sigma _{z1}|0_L`$ (16)
$`\sigma _{z3}\sigma _{z4}|0_L`$ $`=\sigma _{z7}|0_L`$ (17)
$`\sigma _{z1}\sigma _{z2}|1_L=\sigma _{z5}\sigma _{z6}|1_L`$ $`=\sigma _{z3}|1_L`$ (18)
$`\sigma _{z2}\sigma _{z3}|1_L=\sigma _{z4}\sigma _{z5}|1_L=\sigma _{z6}\sigma _{z7}|1_L`$ $`=\sigma _{z1}|1_L`$ (19)
$`\sigma _{z3}\sigma _{z4}|1_L`$ $`=\sigma _{z7}|1_L`$ (20)
and it is easily seen from this that when error correction is applied to the state $`\sigma _{zn}\sigma _{z,n+1}(\alpha |0_L+\beta |1_L)`$, the result will be the state $`(\alpha |0_L+\beta |1_L)`$, that is, a state with the wrong relative phase. Equations (II A) also show that the assumption that different $`\sigma _{zn}\sigma _{z,n+1}`$ lead to orthogonal error states is incorrect: for instance, $`\sigma _{z1}\sigma _{z2}`$ and $`\sigma _{z5}\sigma _{z6}`$ lead to the same error state.

Using Eqs. (II A) it is straightforward to derive the correct form of the approximate failure probability for small $`\delta t`$. When only the terms of first order in $`\delta t`$ are kept in the expression (3), acting on a general state $`\alpha |0_L+\beta |1_L`$, one obtains

$$\delta t(3\sigma _{z1}+2\sigma _{z3}+\sigma _{z7})(\alpha |0_L+\beta |1_L)$$
(21)
which means that the probability that after error correction one might end up with the state $`\alpha |0_L+\beta |1_L`$ is $`(3^2+2^2+1^2)(\delta t)^2=14(\delta t)^2`$. One might directly call this the failure probability (since the state $`\alpha |0_L+\beta |1_L`$ is generally not equal to the initial state $`\alpha |0_L+\beta |1_L`$) or define the failure probability as including explicitly the square of the overlap between the initial and the final state (see Eq. (23), below, for the precise definition), in which case one finds

$$P14(\delta t)^2\left(1\left(|\alpha |^2|\beta |^2\right)^2\right)$$
(22)
where the $``$ sign is a reminder that this is only valid for small $`\delta t`$. The failure probability (22) is largest for the symmetric linear superpositions $`(|0_L\pm |1_L)/\sqrt{2}`$. Note how it is different from the result (5) (for $`N_c=7`$), obtained under the assumption that the errors generated by different $`\sigma _{zn}\sigma _{z,n+1}`$ were all orthogonal, and it is also not equal to the same expression with $`N_c1`$ replaced by $`(N_c1)^2`$, which would correspond to a “worst case scenario” in which all the different errors added coherently. The coefficient $`14`$ falls instead in between the two extremes of 6 $`(N_c1)`$ and 36 ($`(N_c1)^2`$).

The definition of the failure probability $`P`$ used above is

$$P=1\left|\mathrm{\Psi }_0|U(t)|\mathrm{\Psi }_0\right|^2\underset{n=1}{\overset{N_c}{}}\left|\mathrm{\Psi }_0|\sigma _{z,n}U(t)|\mathrm{\Psi }_0\right|^2$$
(23)
that is, one minus the total probability that, after a time $`t`$, the system may be found in either the original state $`|\mathrm{\Psi }_0`$, or one of the $`N_c`$ correctible states $`\sigma _{z,n}|\mathrm{\Psi }_0`$ (there are, of course, other correctible states, of the form $`\sigma _{x,n}|\mathrm{\Psi }_0`$ and $`\sigma _{y,n}|\mathrm{\Psi }_0`$, but they cannot be accessed by this interaction operator).

Using this definition it is relatively straightforward to calculate $`P`$ exactly for arbitrarily large $`\delta t`$. The following relation is useful:

$$0_L|\sigma _{z,n_1},\mathrm{},\sigma _{z,n_k}U|0_L=(1)^k1_L|\sigma _{z,n_1},\mathrm{},\sigma _{z,n_k}U|1_L$$
(25)
$$0_L|\sigma _{z,n_1},\mathrm{},\sigma _{z,n_k}U|1_L=0$$
(26)
where $`k=0,\mathrm{},7`$. This follows directly from the fact that $`U`$ involves always the product of an even number of $`\sigma _z`$ operators, and that the code words used in $`|0_L`$ are the complementaries of the ones used in $`|1_L`$.

For future use it is convenient to define

$$f_0=\left|0_L|U|0_L\right|^2=\left|1_L|U|1_L\right|^2$$
(28)
$$f_1=\underset{n=1}{\overset{7}{}}\left|0_L|\sigma _{z,n}U|0_L\right|^2=\underset{n=1}{\overset{7}{}}\left|1_L|U|1_L\right|^2$$
(29)
Direct calculation then yields

$$f_0=1f_1=\mathrm{cos}^2(\delta t)\mathrm{cos}^2(2\delta t)\mathrm{cos}^2(3\delta t)$$
(30)
and

$$P(t)=f_1\left(1\left(|\alpha |^2|\beta |^2\right)^2\right)$$
(31)
and it is easy to verify that (22) is the lowest-order term in the expansion of (31) in powers of $`\delta t`$.

The failure probability (31) for the case $`|\alpha |=|\beta |`$ is plotted in Figure 1 (solid line) as a function of $`\delta t`$, from $`\delta t=0`$ to $`\pi /2`$. From $`0`$ to $`\pi `$, the curve is symmetric about $`\pi /2`$, as is to be expected from the form of the evolution operator (3).

### B The 5-qubit code

In 1996 Laflamme et al.  introduced a 5-qubit code which can correct general one-qubit errors. A possible encoding is

$`|0_L={\displaystyle \frac{1}{\sqrt{8}}}(`$ $`|00000+|01111|10011+|11100`$ (34)
$`+|10101+|11010+|00110+|01001)`$
$`|1_L={\displaystyle \frac{1}{\sqrt{8}}}(`$ $`|11111+|10000+|01100|00011`$ (36)
$`|01010|00101+|11001+|10110)`$
The codewords can be obtained from the 7-qubit code by deleting any two qubits (in this case, the first two), but the pattern of signs is nontrivial.

Note that the logical 0 and 1 states are still superpositions of eight basis vectors each (where the basis vectors are simultaneous eigenstates of all the $`\sigma _{z,n}`$); hence, the ten one-error states $`\sigma _{z,n}|0_L`$ and $`\sigma _{z,n}|1_L`$, plus the two no-error states $`|0_L`$ and $`|1_L`$ are not enough to span the space resulting from all the possible combinations of the 16 basis kets. As a result, there are now no states which exhibit fully correctible evolution (unlike the $`|0_L`$ and $`|1_L`$ did for the 7-qubit code). However, we still have the phenomenon of non-orthogonal errors. In particular, we find that

$$\sigma _{z2}\sigma _{z3}|0_L=\sigma _{z4}\sigma _{z5}|0_L$$
(38)
$$\sigma _{z2}\sigma _{z3}|1_L=\sigma _{z4}\sigma _{z5}|1_L$$
(39)
and neither of these expressions is equivalent to a linear combination of single-qubit errors acting on $`|0_L`$ or $`|1_L`$; in fact, they are orthogonal to all such states, and to the no-error states $`|0_L`$ and $`|1_L`$ as well.

The calculation of $`P(t)`$, as given by the definition (23) is a little more involved, although the results (II A) still hold. Defining as before

$$f_0=\left|0_L|U|0_L\right|^2=\left|1_L|U|1_L\right|^2$$
(41)
$$f_1=\underset{n=1}{\overset{5}{}}\left|0_L|\sigma _{z,n}U|0_L\right|^2=\underset{n=1}{\overset{5}{}}\left|1_L|U|1_L\right|^2$$
(42)
one finds, thanks to Eqs. (II A), that $`P(t)`$ can still be written as

$$P(t)=1f_0f_1\left(|\alpha |^2|\beta |^2\right)^2$$
(43)
but this time, unlike for the 7-qubit code, one does not have $`f_0=1f_1`$. Instead,

$$f_0=\mathrm{cos}^4(\delta t)\mathrm{cos}^2(2\delta t)$$
(45)
$$f_1=\mathrm{sin}^2(\delta t)\left(1+\mathrm{cos}^2(\delta t)\mathrm{cos}^2(2\delta t)\right)$$
(46)
The failure probability (43) for the case $`|\alpha |=|\beta |`$ is plotted in Figure 2 (solid line) as a function of $`\delta t`$, from $`\delta t=0`$ to $`\pi /2`$. For small values of $`\delta t`$ Eq. (43) reduces to

$$P(\delta t)^2\left[62\left(|\alpha |^2|\beta |^2\right)^2\right]$$
(47)
The reason $`f_01f_1`$ for the 5-qubit code is that, as mentioned above, the states $`|0_L`$ and $`\sigma _{z,n}|0_L`$ do not form a complete basis of the space spanned by the eight kets which make up $`|0_L`$. For the 7-qubit code, instead, they do, and, since the operator $`U`$ acts entirely within that space (when the initial state is $`|0_L`$), any part of $`U|0_L`$ which is not accounted for by the original state $`|0_L`$ has to show up in one of the $`\sigma _{z,n}|0_L`$ (and the same for the initial state $`|1_L`$).

## III Concatenated codes

A relatively straightforward way to obtain a code that can correct more than one-qubit errors is to concatenate a one-error correcting code, in principle, as many times (or levels deep) as necessary . In a code concatenated once, each logical qubit consits of $`N_c`$ blocks of $`N_c`$ (lower-level) qubits each. This code can protect against general errors in any three lower-level qubits, and also against many other combinations of multiple lower-level errors. The minimum number of lower-level errors required to cause the code to fail is four, two in one upper level block (causing that whole block to fail) and two in a different one.

In this section I look first at the 7-qubit and 5-qubit codes concatenated once. A physical (lowest-level) qubit will be identified by two indices, $`n`$ and $`m`$, running from $`1`$ to $`N_c`$: the first index denotes the upper-level block it belongs to, and the second index the qubit’s position within that block. With the assumption that different blocks do not interact, the time-evolution operator (3) becomes

$$U(t)=\underset{n=1}{\overset{N_c}{}}\underset{m=1}{\overset{N_c1}{}}\left(\mathrm{cos}\delta t+i\sigma _{nm}\sigma _{n,m+1}\mathrm{sin}\delta t\right)$$
(48)
From now on, for simplicity, the subscript $`z`$ on $`\sigma _z`$ will be dropped, as only $`\sigma _z`$ is used.

### A The 7-qubit code

If all the errors generated by the individual terms in (48) were orthogonal, it would be a relatively simple matter to calculate the probability of an uncorrectible error to lowest order. It would be given by the terms of the form $`\sigma _{n_1m_1}\sigma _{n_1,m_1+1}\sigma _{n_2m_2}\sigma _{n_2,m_2+1}`$, with $`n_1n_2`$ and $`m_1,m_2`$ arbitrary. There are $`C(N_c,2)`$ ways to choose the indices $`n_1`$ and $`n_2`$, and for each choice there are $`(N_c1)^2`$ different ways to choose $`m_1,m_2`$, so there is a total of

$$\left(\genfrac{}{}{0pt}{}{N_c}{2}\right)\left(N_c1\right)^2$$
(49)
terms in (48) which contribute to the failure probability in lowest order, and all these terms are proportional to $`(\delta t)^4`$ (after squaring the probability amplitudes).

However, as has been shown already in the previous section, not all the terms enumerated in (49) are mutually orthogonal. For the 7-qubit code, according to (II A), there are two sets of values for $`m_1`$ (one set of three values, and one set of two values) which result in identical states. When the equivalent choices for $`m_2`$ are factored in and added, and the probability amplitudes of the corresponding states are squared, one finds, instead of the factor $`(N_c1)^2=36`$ in (49), a factor

$$(3^2+2^2+1^2)(3^2+2^2+1^2)=14^2=196$$
(50)
As a result, to lowest order, the failure probability for this two-level code is given by

$`P^{(2)}`$ $`196\left({\displaystyle \genfrac{}{}{0pt}{}{7}{2}}\right)(\delta t)^4\left(1\left(|\alpha |^2|\beta |^2\right)^2\right)`$ (52)
$`=4116(\delta t)^4\left(1\left(|\alpha |^2|\beta |^2\right)^2\right)`$
for a state $`\alpha |0_L+\beta |1_L`$. The dependence on $`\alpha `$ and $`\beta `$ exhibited by (52) may not be not immediately obvious at this point (it will be derived later), but it shows that, for the concatenated code also, the individual states $`|0_L`$ and $`|1_L`$ are fully correctible.

Note that the maximum error probability $`P_{max}^{(2)}`$ for the once-concatenated code, according to (52), is related to the corresponding one for the plain code, (22), by the formula

$$P_{max}^{(2)}=\left(\genfrac{}{}{0pt}{}{7}{2}\right)\left(P_{max}^{(1)}\right)^2$$
(53)
and the right-hand side of (53) is simply the probability that two blocks might fail independently, times the number of ways to choose two blocks. This is because, to this order (recall that these are only approximate results for small $`\delta t`$), errors across different blocks are orthogonal, ie., they add incoherently. This can also be seen from Equation (50): within each block there is some coherent addition of error amplitudes (which results in the terms $`3^2`$ and $`2^2`$), but then the total error probabilities for each block are simply multiplied.

This is an important result which generalizes to higher orders, as discussed later in this Section (subsection C). The key point is that all the errors which would cause a particular block to fail are orthogonal to all the errors that would cause a different block to fail, to lowest order. This is because to lowest order all it takes is two errors in one block to cause that block to fail, and because the original code is nondegenerate, which means that all expectation values of the form $`\xi _L|\sigma _i\sigma _j|\eta _L`$ (where $`\eta ,\xi =0,1`$) within the same block are zero.

Equation (53) does not, however, hold for arbitrarily large $`\delta t`$, when multiple (more than two) errors within the same block need to be taken into account. For a code concatenated only once it is still possible, albeit very cumbersome, to derive this failure probability for arbitrary $`\delta t`$. The definition (23) has to be changed to

$$P=1\left|\mathrm{\Psi }_0|U(t)|\mathrm{\Psi }_0\right|^2\underset{n}{}\left|\mathrm{\Psi }_0|A_nU(t)|\mathrm{\Psi }_0\right|^2$$
(54)
where the $`\{A_n\}`$ represent a maximal set of correctible error operators leading to orthogonal states. The difficulty is that now, unlike in the previous Section, not all possible correctible errors lead to orthogonal states. For instance, any number of errors within one single block are now correctible, but certainly not all of these lead to orthogonal states. Thus, even though we started out with a nondegenerate code, the concatenated code is degenerate .

As far as I can tell, the set $`\{A_n\}`$ in (54) must be constructed by inspection. For the 7-qubit code, concatenated once, I find that the maximum number of correctible errors acting on a single block and leading to orthogonal states is 15, which can be chosen to be

$`\sigma _{nm}(m=1,\mathrm{},7)`$ (55)
$`\sigma _{n1}\sigma _{nm}(m=2,\mathrm{},7)`$ (56)
$`\sigma _{n2}\sigma _{n3}`$ (57)
$`\sigma _{n1}\sigma _{n4}\sigma _{n5}`$ (58)
Of these, the first 7 are, of course, fully correctible at the lower level, whereas the others, upon correction at the lower level, lead to what look like “block” phase errors at the higher level.

The single-block errors (58) can be combined in a large number of ways among the different blocks to yield other, orthogonal, correctible errors. Fortunately, however, the actual number of such combinations which make a nonvanishing contribution to (54) is relatively small. To see this, divide the error operators acting on a single block into those involving an even number of $`\sigma `$’s (denoted as $`E_n`$ if they act on the $`n`$-th block) and those involving an odd number of $`\sigma `$’s (denoted as $`O_n`$). Then, for a general multiblock error of the form $`A=E_{n_1}\mathrm{}E_{n_k}O_{m_1}\mathrm{}O_{m_k}`$, Eqs. (II A) for a single block imply that the expectation value

$$0_L|O_{n_1}\mathrm{}O_{n_k}E_{m_1}\mathrm{}E_{m_l}U|0_L$$
(59)
will only be nonzero if the corresponding combination $`0_L|\sigma _{n_1}\mathrm{}\sigma _{n_k}|0_L`$ is nonzero for the single, non-concatenated, code. The same is true for the state $`|1_L`$, and there are no cross-terms between $`|0_L`$ and $`|1_L`$ (also as a direct consequence of the single-block Eqs. (II A)). Thus, to evaluate (54), it is sufficient to identify the combinations of numbers $`\{n_1,\mathrm{},n_k\}`$ which lead to a nonzero $`0_L|\sigma _{n_1}\mathrm{}\sigma _{n_k}|0_L`$ for the 7-qubit code, and put odd errors, chosen from the list (58), in the blocks $`n_1,\mathrm{},n_k`$, augmented perhaps with even errors, always chosen from the list (58), in any other blocks, provided, of course, that the total error be correctible. This quickly reduces the total of different possibilities to a manageable number: for instance, since an even error in one block would cause that block to fail, one cannot have even errors in more than one block at a time.

For reference, for the 7-qubit code, the sets of values of $`n_1,\mathrm{},n_k`$ to be used are: $`(1,2,3)`$, $`(1,4,5)`$, $`(1,6,7)`$, $`(2,4,6)`$, $`(2,5,7)`$, $`(3,4,7)`$, $`(3,5,6)`$, the complementary set ($`(4,5,6,7)`$, $`(2,3,6,7)`$, and so on), and the set $`(1,2,3,4,5,6,7)`$. Then, possible error operators $`A_n`$ to be used in (54) might be products such as such as $`\sigma _{1,m_1}\sigma _{2,m_2}\sigma _{3,m_3}`$, or $`\sigma _{1,m_1}\sigma _{2,m_2}\sigma _{3,m_3}\sigma _{4,1}\sigma _{4,m_4}`$, or $`\sigma _{n,1}\sigma _{n,m}`$. In terms of the functions $`f_0`$, $`f_1`$ defined in (II A), the term(s) $`_{m_1,m_2,m_3}|0_L|\sigma _{1,m_1}\sigma _{2,m_2}\sigma _{3,m_3}|0_L|^2`$ contribute a term $`f_0^4f_1^3`$ to the expression (54); the term $`_{m_1,m_2,m_3,m_4}|0_L|\sigma _{1,m_1}\sigma _{2,m_2}\sigma _{3,m_3}\sigma _{4,1}\sigma _{4,m_4}|0_L|^2`$ contributes a term $`f_0^3f_1^4`$ (since the product of two $`\sigma `$’s in the same block is equivalent, except for a sign, to a single $`\sigma `$ acting on the appropriate state, cf. Eqs. (II A)), and so on.

Collecting all the terms, and keeping track carefully of which ones have opposite signs in the states $`|0_L`$ and the states $`|1_L`$, one finds

$$P^{(2)}(t)=\left(1f_0^77f_0^6f_128f_0^4f_1^37f_0^3f_1^421f_0^2f_1^5\right)\left(1\left(|\alpha |^2|\beta |^2\right)^2\right)$$
(60)
with $`f_0`$ and $`f_1`$ given by (30), which does reduce to (52) for small $`\delta t`$.

The exact result (60) is plotted as a function of $`\delta t`$ in Figure 1 (dashed line), for $`|\alpha |=|\beta |`$. It is, perhaps, remarkable that it does not look all that different from the single-encoding result (31), except in the very small $`\delta t`$ region. It is also clear that Eq. (53) does not hold for large $`\delta t`$, which means that for large $`\delta t`$ we do see the effects of errors across different blocks adding coherently.

### B The 5-qubit code

For the 5-qubit code, one again has to identify the appropriate maximal set of error operators leading to orthogonal errors. For a single block, one finds the following set of 15:

$`\sigma _{nm}(m=1,\mathrm{},5)`$ (61)
$`\sigma _{n1}\sigma _{nm}(m=2,\mathrm{},5)`$ (62)
$`\sigma _{n2}\sigma _{nm}(m=3,\mathrm{},5)`$ (63)
$`\sigma _{n1}\sigma _{n2}\sigma _{nm}(m=3,\mathrm{},5)`$ (64)
As in the previous subsection, it is found that only certain combinations of odd errors across different blocks lead to nonzero contributions to (54). These turn out to be much fewer than before, being limited to $`n_1,\mathrm{},n_k=(1,2,5),(1,3,4)`$, and $`(2,3,4,5)`$.

The result, valid for arbitrarily large $`\delta t`$, turns out to be now

$`P^{(2)}(t)=1`$ $`(f_0^5+5f_0^4(1f_0)+f_1^4+4f_0f_1^3(1f_1))`$ (66)
$`(4f_0f_1^3+6f_0^2f_1^28f_0^2f_1^3)\left(|\alpha |^2|\beta |^2\right)^2`$
with $`f_0`$ and $`f_1`$ given by Eqs. (43). The factors of $`1f_0`$ appearing in (66) come from the sum of the single-block expectation values squared of $`\sigma _1\sigma _m`$, $`m=2,\mathrm{},5`$ and $`\sigma _2\sigma _m`$, $`m=3,\mathrm{},5`$; these seven error operators, acting on the singly encoded $`|0_L`$ or $`|1_L`$, generate the seven orthogonal states necessary to complete a basis of the space spanned by the 8 kets appearing in (II B). The factors of $`1f_1`$ come from the sum of the single-block expectation values squared of $`\sigma _1\sigma _2\sigma _m`$, $`m=3,\mathrm{},5`$: these 3 error operators turn out to be the complement of the 5 ones appearing in the definition (39) of $`f_1`$.

The limit of (66) for small $`\delta t`$ is

$$P^{(2)}(t)(\delta t)^4\left[36024\left(|\alpha |^2|\beta |^2\right)^2\right]$$
(67)
As usual, the failure probability is largest for the symmetric superposition state $`(|0_L\pm |1_L)/\sqrt{2}`$. The corresponding coefficient, $`360`$, can be obtained also from the simple arguments used to derive the result (52) in the previous subsection: it equals $`C(5,2)(2^2+1^2+1^2)^2`$, where the $`2^2`$ stands for the amplitude square of the sum of the two error terms which lead to identical errors, as shown in Eq. (II B) (and the two $`1^2`$ stand for the other two orthogonal errors).

The exact result (66) is plotted in Figure 2 (dashed line) as a function of $`\delta t`$, for $`|\alpha |=|\beta |`$. In this case, unlike for the 7-qubit code, the difference with the single-encoding formula (solid line, equation (43)) is substantial; the twice-encoded qubit has a substantially smaller failure probability for most values of $`\delta t`$. The failure probability actually goes to zero at $`\delta =\pi /2`$, where Eq. (48) yields $`U=i_{n=1}^5\sigma _{n1}\sigma _{n5}`$. Upon inspection, it is found to be a peculiarity of the concatenated (two-levels deep) 5-qubit code that

$$\underset{n=1}{\overset{5}{}}\sigma _{n1}\sigma _{n5}\left(\alpha |0_L+\beta |1_L\right)=\sigma _{11}\sigma _{15}\sigma _{22}\sigma _{32}\sigma _{42}\sigma _{52}\left(\alpha |0_L+\beta |1_L\right)$$
(68)
and the state on the right-hand side of (68) is fully correctible, since it amounts to only one block failing.

### C Results for arbitrarily deep concatenation

Clearly, the exact calculation of $`P^{(n)}`$ for a code concatenated $`n`$ levels deep, with $`n>2`$, becomes much too cumbersome to be feasible. However, the results in the previous subsections show that the leading term in $`\delta t`$ can be calculated using relatively simple arguments.

Consider, for instance, a code concatenated three levels deep. It will fail if at least two of the high-level qubits fail, but these high-level qubits are themselves two-level deep codes. One can choose the two highest-level qubits that fail in any of $`C(N_c,2)`$ ways, and for each of them the number of terms in the expansion of the evolution operator which lead to a qubit failure are those computed for the previous level. So, for the 7-qubit code, we can immediately iterate (52) to get

$`P_{max}^{(3)}`$ $`\left({\displaystyle \genfrac{}{}{0pt}{}{7}{2}}\right)\left(4116(\delta t)^4\right)^2`$ (70)
$`\left({\displaystyle \genfrac{}{}{0pt}{}{7}{2}}\right)\left(P_{max}^{(2)}\right)^2`$
where the subscript “max” indicates that we are looking at the failure probability for the state that maximizes it (the linear combination $`(|0_L\pm |1_L)/\sqrt{2}`$), and Eq. (22) has been used for $`P^{(1)}`$.

Equation (70) contains the essence of a recursion relation which yields the failure probability at the $`n`$-th level of concatenation ($`n2`$) as

$$P_{max}^{(n)}\left(\genfrac{}{}{0pt}{}{7}{2}\right)^1\left(\left(\genfrac{}{}{0pt}{}{7}{2}\right)P_{max}^{(1)}\right)^{2^{n1}}$$
(71)
Figure 3 shows the approximate $`P^{(n)}`$, calculated from (71), for $`n=1,\mathrm{},4`$ (solid lines), and also the exact results (31) and (60), for comparison. Based on these examples, it would seem that Equation (71) overestimates the failure probability, which would make it a somewhat conservative estimate

As encoding depth $`n`$ increases, the failure probability (71) will, for a given $`\delta t`$ increase or decrease depending on whether $`C(7,2)P_{max}^{(1)}`$ is greater than or less than one. This yields the “threshold” value for $`(\delta t)^2`$, below which virtually error-free operation can be achieved for sufficiently deep encoding. Using (22), one finds the condition:

$$(\delta t)^2<\frac{1}{2114}=3.4\times 10^3$$
(72)
(Figure 3 shows this as the point where all the curves cross.) This threshold (72) is comparable to some estimates for the correction of independent, random errors by concatenate codes, if one thinks of $`(\delta t)^2`$ as a sort of failure probability for a single physical qubit. As pointed out in the Introduction, this is a natural interpretation (in spite of the fact that two physical qubits at a time are modified by this interaction), because of the result (1).

A similar calculation yields the threshold result for the 5-qubit code (using (47) and C(5,2)=10):

$$(\delta t)^2<\frac{1}{106}=1.7\times 10^2$$
(73)
## IV Discussion

The results in the previous section suggest that, in principle, the correction of even rather large undesired interactions between the qubits does not require special methods or resources beyond those needed to control independent, random errors arising from the interaction of the qubits with the environment. This observation is based, first, on the fact that estimates of the threshold for fault-tolerant computation with independent environmental errors tend to be rather more restrictive than either (72) or (73) ; and, second, on the fact that for most proposed quantum computing systems the unwanted interactions betwen qubits would easily satisfy the threshold conditions (72) or (73), provided that the time in between error corrections is not excessive.

As an example, consider a system of nuclear spins, in a solid, for instance (such as in the proposal by Kane ), separated by $`150`$ Å, and interacting via the magnetic dipole-dipole interaction, which would naturally lead to a Hamiltonian like (2). For such a system $`\delta `$ may be estimated as

$$\delta \frac{\mu _0\mu _N^2}{4\pi \mathrm{}d^3}10^2\text{s}^1$$
(74)
Here $`\mu _N`$ is the nuclear magneton.

Without error correction $`\delta t`$ becomes already of the order of 1 after 100 s, and long calculations are impossible. But assume that one uses error correction and that complete error correction of an $`n`$-order encoded qubit requires a time $`t=c_n\tau _{gate}`$, where the time needed to perform an elementary logic gate for this system is about $`\tau _{gate}=10^5`$ s. Then the threshold is easily reached unless $`c_n`$ is exceptionally large, which would make error correction generally impossible anyway (a large $`c_n`$ would mean essentially that the computer architecture does not allow for a large degree of parallelism when correcting errors across different qubits).

Of course, the results presented here have been obtained only for a specific kind of interaction Hamiltonian and a specific geometry, and would need to be modified for other cases (e.g., two-dimensional geometries, with more nearest-neighbors per qubit, or to include non-nearest neighbor interactions). There is also the simplifying assumption made from the outset, of neglecting interactions beween qubits in different blocks. These complications might require recalculations for specific physical systems and hardware configurations, but I do not expect them to change substantially the threshold estimates of the previous Section.

For instance, for the kind of “linear” geometry considered here, including interactions between the last qubit of one block and the first one of the next block would not modify the failure probability to lowest order, since a single two-qubit error which straddles two blocks does not cause either block to fail; hence, to make any one block fail, more of these errors are required than of the kind of “internal” errors considered in this paper.

The key result of Section III concerning the incoherent addition of failure probabilities for different blocks, to lowest order, should be quite general: it depends only on the base code being a nondegenerate distance 3 code and the interaction Hamiltonian being a sum of pairwise products of elementary error operators for individual qubits.

The most important “real-life” complication which has been ignored here is, rather, the fact that, in practice, error correction will not be a clean and purely mathematical operation, but a physical one, involving interacting systems and therefore itself subject to error. The theory of fault-tolerant error correction has shown that even when this is accounted for it should still be possible to achieve “almost” perfect error correction, provided that the error rate per gate is below a certain threshold, but it is here where one finds the (possibly) very small error thresholds which constitute the main challenge for large-scale quantum computing.

I am grateful to R. Laflamme and D. A. Meyer for comments. This research has been sponsored by the National Science Foundation and by the Army Research Office.