an upper bound to the original problem, and the third inequality is a
consequence of $\mathrm{SO}(n)\subset\mathrm{O}(n)$. The task is then to
determine such an $\alpha$ which satisfies the first inequality of Eq. (974).
The core argument is a generalization of the Rietz method
[alon2004approximating], which proceeds by constructing a positive
semidefinite matrix $S\in\mathbb{R}^{mn\times mn}$ whose $(u,v)$th block is
defined as

$S_{uv}\coloneqq\mathopen{}\left(X_{u}Z-\alpha^{-1}\widetilde{\mathcal{P}}(X_{u}Z)\right)\mathclose{}\mathopen{}\left(X_{v}Z-\alpha^{-1}\widetilde{\mathcal{P}}(X_{v}Z)\right)\mathclose{}^{\mathsf{T}}.$
(975)

The expected value of this matrix is

$\begin{split}\operatorname*{\mathbb{E}}[S_{uv}]&=\operatorname*{\mathbb{E}}\mathopen{}\left[X_{u}Z(X_{v}Z)^{\mathsf{T}}-\alpha^{-1}\widetilde{\mathcal{P}}(X_{u}Z)(X_{v}Z)^{\mathsf{T}}\right.\mathclose{}\\\
&\quad\mathopen{}\left.-\alpha^{-1}(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}+\alpha^{-2}\widetilde{\mathcal{P}}(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}\\\
&=X_{u}\operatorname*{\mathbb{E}}\mathopen{}\left[ZZ^{\mathsf{T}}\right]\mathclose{}X_{v}-\alpha^{-1}\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(X_{u}Z)(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}\\\
&\quad-\alpha^{-1}\operatorname*{\mathbb{E}}\mathopen{}\left[(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}+\alpha^{-2}\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}.\end{split}$
(976)

Because $ZZ^{\mathsf{T}}$ is a Wishart matrix with covariance matrix
$\mathbb{I}_{n}/n$, we have
$\operatorname*{\mathbb{E}}[ZZ^{\mathsf{T}}]=\mathbb{I}_{n}$. Meanwhile,
$\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}$
is the quantity we wish to bound. To compute the expected values of the two
cross terms, we invoke Lemma 0.37 which holds because
$X_{u}X_{u}^{\mathsf{T}}=X_{v}X_{v}^{\mathsf{T}}=\mathbb{I}_{n}$:

$\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(X_{u}Z)(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}=\operatorname*{\mathbb{E}}\mathopen{}\left[(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}=\alpha_{\mathrm{SO}(n)}X_{u}X_{v}^{\mathsf{T}}.$
(977)

Thus, setting $\alpha=\alpha_{\mathrm{SO}(n)}$, we obtain

$\begin{split}\operatorname*{\mathbb{E}}[S_{uv}]&=X_{u}X_{v}^{\mathsf{T}}-X_{u}X_{v}^{\mathsf{T}}-X_{u}X_{v}^{\mathsf{T}}+\alpha_{\mathrm{SO}(n)}^{-2}\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}\\\
&=-X_{u}X_{v}^{\mathsf{T}}+\alpha_{\mathrm{SO}(n)}^{-2}\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\right]\mathclose{}.\end{split}$
(978)

Finally, using the fact that $C,S\succeq 0$, we have that $\langle
C,S\rangle\geq 0$ and so $\operatorname*{\mathbb{E}}\langle C,S\rangle\geq 0$,
which implies that

$\operatorname*{\mathbb{E}}\mathopen{}\left[\sum_{(u,v)\in E}\langle
C_{uv},\widetilde{\mathcal{P}}(X_{u}Z)\widetilde{\mathcal{P}}(X_{v}Z)^{\mathsf{T}}\rangle\right]\mathclose{}\geq\alpha_{\mathrm{SO}(n)}^{2}\sum_{(u,v)\in
E}\langle C_{uv},X_{u}X_{v}^{\mathsf{T}}\rangle.$ (979)

Then by Eq. (974) the claim follows. ∎

We now establish the value of

$\alpha_{\mathrm{SO}(n)}=\operatorname*{\mathbb{E}}\mathopen{}\left[\frac{1}{n}\sum_{i\in[n-1]}\sigma_{i}(Z_{1})\right]\mathclose{}$
(980)

from Lemmas 0.37 and 0.38. Because Lemma 0.37 is somewhat technical and the
argument is virtually unchanged by replacing $\mathcal{P}$ with
$\widetilde{\mathcal{P}}$, we refer the reader to Ref.
[bandeira2016approximating] for proof details. Instead, we simply note that
the only part of the proof for Lemma 0.37 which does depend on the change to
$\widetilde{\mathcal{P}}$ is the final result, wherein it is established that

$\operatorname*{\mathbb{E}}\mathopen{}\left[(MZ)\widetilde{\mathcal{P}}(NZ)^{\mathsf{T}}\right]\mathclose{}=\operatorname*{\mathbb{E}}\mathopen{}\left[Z_{1}\widetilde{\mathcal{P}}(Z_{1})^{\mathsf{T}}\right]\mathclose{}MN^{\mathsf{T}},$
(981)

where $Z_{1}\in\mathbb{R}^{n\times n}$ has entries i.i.d. from
$\mathcal{N}(0,1/n)$. Thus proving Lemma 0.38 is the key component in
establishing the value of the approximation ratio
$\alpha_{\mathrm{SO}(n)}^{2}$.

###### Proof (of Lemma 0.38).

Consider the singular value decomposition of $Z_{1}=U\Sigma
V^{\mathsf{T}}\in\mathbb{R}^{n\times n}$. Its special singular value
decomposition can be written as $Z_{1}=U(\Sigma
J_{UV^{\mathsf{T}}})(VJ_{UV^{\mathsf{T}}})^{\mathsf{T}}$, where
$J_{UV^{\mathsf{T}}}$ is the $n\times n$ diagonal matrix

$J_{UV^{\mathsf{T}}}\coloneqq\begin{bmatrix}I_{n-1}&0\\\ 0&\det
UV^{\mathsf{T}}\end{bmatrix}.$ (982)

Note that $J_{UV^{\mathsf{T}}}=J_{U}J_{V}$. Using the fact that the (special)
rounding operator returns

$\widetilde{\mathcal{P}}(Z_{1})=UJ_{UV^{\mathsf{T}}}V^{\mathsf{T}},$ (983)

we have

$\widetilde{\mathcal{P}}(Z_{1})Z_{1}^{\mathsf{T}}=UJ_{U}J_{V}\Sigma
U^{\mathsf{T}}.$ (984)

Because $Z_{1}$ is a random Gaussian matrix with i.i.d. entries, its singular
values and left- and right-singular vectors are distributed independently
[tulino2004random]. In particular, both $U$ and $V$ are distributed according
to the Haar measure on $\mathrm{O}(n)$. The expected value of Eq. (984) can
therefore be split into three independent averages:

$\operatorname*{\mathbb{E}}_{Z_{1}\sim\mathcal{N}(0,\mathbb{I}_{n}/n)}\mathopen{}\left[\widetilde{\mathcal{P}}(Z_{1})Z_{1}^{\mathsf{T}}\right]\mathclose{}=\operatorname*{\mathbb{E}}_{\Sigma\sim
D}\operatorname*{\mathbb{E}}_{U\sim\mathrm{O}(n)}\mathopen{}\left[UJ_{U}\operatorname*{\mathbb{E}}_{V\sim\mathrm{O}(n)}\mathopen{}\left[J_{V}\right]\mathclose{}\Sigma
U^{\mathsf{T}}\right]\mathclose{}.$ (985)

(We shall comment on the distribution $D$ of singular values later.) Because
$\mathrm{O}(n)$ is evenly divided into its unconnected $({+1})$\- and
$({-1})$-determinant components, the average determinant vanishes:
$\operatorname*{\mathbb{E}}_{V\sim\mathrm{O}(n)}[\det V]=0$. This leaves us
with

$\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(Z_{1})Z_{1}^{\mathsf{T}}\right]\mathclose{}=\operatorname*{\mathbb{E}}\mathopen{}\left[U\overline{\Sigma}U^{\mathsf{T}}\right]\mathclose{},$
(986)

where

$\overline{\Sigma}=\begin{bmatrix}\sigma_{1}(Z_{1})&&&\\\ &\ddots&&\\\
&&\sigma_{n-1}(Z_{1})&\\\ &&&0\end{bmatrix}.$ (987)

The Haar average over $U\sim\mathrm{O}(n)$ in Eq. (986) is well-known
[collins2006integration] to be proportional to the identity,
$\operatorname*{\mathbb{E}}[U\overline{\Sigma}U^{\mathsf{T}}]=\lambda\mathbb{I}_{n}$,
and the constant of proportionality can be determined by considering its
trace:

$\begin{split}n\lambda=\operatorname{tr}\mathopen{}\left(\lambda\mathbb{I}_{n}\right)\mathclose{}&=\operatorname{tr}\mathopen{}\left(\operatorname*{\mathbb{E}}[U\overline{\Sigma}U^{\mathsf{T}}]\right)\mathclose{}\\\
&=\operatorname*{\mathbb{E}}[\operatorname{tr}\overline{\Sigma}]\\\
&=\operatorname*{\mathbb{E}}\mathopen{}\left[\sum_{i\in[n-1]}\sigma_{i}(Z_{1})\right]\mathclose{}.\end{split}$
(988)

Hence $\lambda=\alpha_{\mathrm{SO}(n)}$ and so

$\operatorname*{\mathbb{E}}\mathopen{}\left[\widetilde{\mathcal{P}}(Z_{1})Z_{1}^{\mathsf{T}}\right]\mathclose{}=\alpha_{\mathrm{SO}(n)}\mathbb{I}_{n}.$
(989)

The corresponding statement for
$\operatorname*{\mathbb{E}}\mathopen{}\left[Z_{1}\widetilde{\mathcal{P}}(Z_{1})^{\mathsf{T}}\right]\mathclose{}$
follows completely analogously, essentially by interchanging the roles of $U$
and $V$. The entire argument is equivalent because $U$ and $V$ are i.i.d. ∎

To numerically evaluate $\alpha_{\mathrm{SO}(n)}$ we can use the linearity of
expectation,

$\alpha_{\mathrm{SO}(n)}=\alpha_{\mathrm{O}(n)}-\frac{1}{n}\operatorname*{\mathbb{E}}\mathopen{}\left[\sigma_{n}(Z_{1})\right]\mathclose{}.$
(990)

The distribution of singular values of random Gaussian matrices can be
analyzed from the theory of Wishart matrices. In particular,
$Z_{1}Z_{1}^{\mathsf{T}}=U\Sigma^{2}U^{\mathsf{T}}$ is a Wishart matrix with
covariance matrix $\mathbb{I}_{n}/n$, so the distribution of singular values
$\Sigma$ is the square root of the Wishart distribution of eigenvalues. The
quantity $\alpha_{\mathrm{O}(n)}$ in terms of the marginal distribution
$p_{n}^{\text{(avg)}}(x)$ of Wishart eigenvalues $x\in(0,\infty)$ was studied
in Ref. [bandeira2016approximating], yielding the expression

$\alpha_{\mathrm{O}(n)}=\frac{1}{\sqrt{n}}\int_{0}^{\infty}p_{n}^{\text{(avg)}}(x)\sqrt{x}\,dx.$
(991)

Note the factor of $n^{-1/2}$, which is introduced because the distribution
$p_{n}^{\text{(avg)}}(x)$ is normalized to have unit variance. An explicit
expression of $p_{n}^{\text{(avg)}}(x)$ can be found in Refs.
[bandeira2016approximating, Lemma 21] and [livan2011momets, Eq. (16)].

For our newly derived approximation ratio $\alpha_{\mathrm{SO}(n)}$, we need
to additionally evaluate the expected smallest singular value of this Wishart
distribution. This minimum-eigenvalue distribution was studied in Ref.
[edelman1988eigenvalues], wherein an analytical expression was derived (again
assuming unit variance):

$p_{n}^{\text{(min)}}(x)=\frac{n}{2^{n-1/2}}\frac{\Gamma(n)}{\Gamma(n/2)}\frac{e^{-xn/2}}{\sqrt{x}}U\mathopen{}\left(\frac{n-1}{2},-\frac{1}{2},\frac{x}{2}\right)\mathclose{}.$
(992)

Here, $U(a,b,z)$ with $a>0$ and $b<1$ is the Tricomi confluent hypergeometric
function, the unique solution to the differential equation

$z\frac{d^{2}U}{dz^{2}}+(b-z)\frac{dU}{dz}-aU=0$ (993)

with boundary conditions $U(a,b,0)=\Gamma(1-b)/\Gamma(1+a-b)$ and
$\lim_{z\to\infty}U(a,b,z)=0$. The expression for the average smallest
singular value is therefore

$\operatorname*{\mathbb{E}}\mathopen{}\left[\sigma_{n}(Z_{1})\right]\mathclose{}=\frac{1}{\sqrt{n}}\int_{0}^{\infty}p_{n}^{\text{(min)}}(x)\sqrt{x}\,dx.$
(994)

Altogether, we arrive at the integral expression for

$\alpha_{\mathrm{SO}(n)}=\frac{1}{\sqrt{n}}\int_{0}^{\infty}\mathopen{}\left[p_{n}^{\text{(avg)}}(x)-\frac{1}{n}p_{n}^{\text{(min)}}(x)\right]\mathclose{}\sqrt{x}\,dx.$
(995)

### LNCG Hamiltonian symmetries

Here we demonstrate the local $\mathrm{O}(n)$ symmetry discussed in Section
Vertex-marginal rounding. Consider an edge term

$H_{uv}=\sum_{i,j\in[n]}[C_{uv}]_{ij}\sum_{k\in[n]}P_{ik}^{(u)}\otimes
P_{jk}^{(v)}.$ (996)

Because $P_{ik}=\mathrm{i}\widetilde{\gamma}_{i}\gamma_{k}$ and the sum over
$k$ is independent of $C_{uv}$, we can factor out each
$\gamma_{k}^{(u)}\otimes\gamma_{k}^{(v)}$ and rewrite the Hamiltonian term as
as

$H_{uv}=-\sum_{i,j\in[n]}[C_{uv}]_{ij}\mathopen{}\left(\widetilde{\gamma}_{i}^{(u)}\otimes\widetilde{\gamma}_{j}^{(v)}\right)\mathclose{}\mathopen{}\left(\sum_{k\in[n]}\gamma_{k}^{(u)}\otimes\gamma_{k}^{(v)}\right)\mathclose{}.$
(997)

The operator $\sum_{k\in[n]}\gamma_{k}^{(u)}\otimes\gamma_{k}^{(v)}$ is
invariant to any orthogonal transformation $V\in\mathrm{O}(n)$ which acts
identically on both vertices:

$\begin{split}\mathcal{U}_{(\mathbb{I}_{n},V)}^{\otimes
2}\mathopen{}\left(\sum_{k\in[n]}\gamma_{k}^{(u)}\otimes\gamma_{k}^{(v)}\right)\mathclose{}(\mathcal{U}_{(\mathbb{I}_{n},V)}^{\otimes
2})^{\dagger}&=\sum_{k\in[n]}\mathopen{}\left(\sum_{\ell\in[n]}V_{k\ell}\gamma_{\ell}^{(u)}\right)\mathclose{}\otimes\mathopen{}\left(\sum_{\ell^{\prime}\in[n]}V_{k\ell^{\prime}}\gamma_{\ell^{\prime}}^{(v)}\right)\mathclose{}\\\
&=\sum_{\ell,\ell^{\prime}\in[n]}\mathopen{}\left(\sum_{k\in[n]}[V^{\mathsf{T}}]_{\ell
k}V_{k\ell^{\prime}}\right)\mathclose{}\gamma_{\ell}^{(u)}\otimes\gamma_{\ell^{\prime}}^{(v)}\\\
&=\sum_{\ell,\ell^{\prime}\in[n]}\delta_{\ell\ell^{\prime}}\gamma_{\ell}^{(u)}\otimes\gamma_{\ell^{\prime}}^{(v)}\\\
&=\sum_{\ell\in[n]}\gamma_{\ell}^{(u)}\otimes\gamma_{\ell}^{(v)}.\end{split}$
(998)

Because $\mathcal{U}_{(\mathbb{I}_{n},V)}$ acts trivially on all
$\widetilde{\gamma}_{i}$, it follows that

$\mathcal{U}_{(\mathbb{I}_{n},V)}^{\otimes
2}H_{uv}(\mathcal{U}_{(\mathbb{I}_{n},V)}^{\otimes 2})^{\dagger}=H_{uv}$ (999)

for each $(u,v)$. Finally, this symmetry can be straightforwardly extended to
all $m$ vertices:

$\mathcal{U}_{(\mathbb{I}_{n},V)}^{\otimes
m}H(\mathcal{U}_{(\mathbb{I}_{n},V)}^{\otimes m})^{\dagger}=H.$ (1000)

Now we investigate some consequences of this continuous symmetry. The
following lemma is particularly important, as it necessitates the use of the
one-body perturbation $\zeta H_{1}$ to break this symmetry when preparing of
eigenstates of $H$.

######  Lemma 0.39.

Let $|\psi\rangle$ be a nondegenerate eigenstate of $H$. Then for each single-
vertex marginal $\sigma_{v}\coloneqq\operatorname{tr}_{\neg
v}|\psi\rangle\\!\langle\psi|$, $v\in[m]$, we have

$Q(\sigma_{v})=0.$ (1001)

###### Proof.

Consider the expansion of its density matrix $|\psi\rangle\\!\langle\psi|$ in
the Majorana operator basis, up to the relevant one-body expectation values:

$|\psi\rangle\\!\langle\psi|=\frac{1}{d^{m}}\mathopen{}\left(\mathbb{I}^{\otimes
m}+\sum_{v\in[m]}\sum_{i,j\in[n]}[Q(\sigma_{v})]_{ij}\mathrm{i}\widetilde{\gamma}_{i}^{(v)}\gamma_{j}^{(v)}+\cdots\right)\mathclose{},$
(1002)

where we recall that
$[Q(\sigma_{v})]_{ij}=\langle\psi|\mathrm{i}\widetilde{\gamma}_{i}^{(v)}\gamma_{j}^{(v)}|\psi\rangle$.
Due to the symmetry [Eq. (1000)], for every $V\in\mathrm{O}(n)$ the state
$|\psi(V)\rangle={\mathcal{U}}_{(\mathbb{I}_{n},V)}^{\otimes m}|\psi\rangle$
is also an eigenvector of $H$ with the same eigenvalue. The one-body
expectation values of $|\psi\rangle$ are therefore transformed as

$\begin{split}{\mathcal{U}}_{(\mathbb{I}_{n},V)}^{\otimes
m}\mathopen{}\left(\sum_{v\in[m]}\sum_{i,j\in[n]}[Q(\sigma_{v})]_{ij}\mathrm{i}\widetilde{\gamma}_{i}^{(v)}\gamma_{j}^{(v)}\right)\mathclose{}({\mathcal{U}}_{(\mathbb{I}_{n},V)}^{\otimes
m})^{\dagger}&=\sum_{v\in[m]}\sum_{i,j\in[n]}[Q(\sigma_{v})]_{ij}\mathrm{i}\sum_{i^{\prime}\in[n]}V_{ii^{\prime}}\widetilde{\gamma}_{i^{\prime}}^{(v)}\gamma_{j}^{(v)}\\\
&=\sum_{v\in[m]}\sum_{i^{\prime},j\in[n]}[V^{\mathsf{T}}Q(\sigma_{v})]_{i^{\prime}j}\mathrm{i}\widetilde{\gamma}_{i^{\prime}}^{(v)}\gamma_{j}^{(v)}.\end{split}$
(1003)

Now suppose that $|\psi\rangle$ is nondegenerate. Then we have that
$|\psi\rangle\\!\langle\psi|=|\psi(V)\rangle\\!\langle\psi(V)|$ for all
$V\in\mathrm{O}(n)$, and in particular we can take the Haar integral over
$\mathrm{O}(n)$ of this identity:

$\begin{split}\int_{\mathrm{O}(n)}d\mu(V)|\psi(V)\rangle\\!\langle\psi(V)|=\int_{\mathrm{O}(n)}d\mu(V)|\psi\rangle\\!\langle\psi|=|\psi\rangle\\!\langle\psi|,\end{split}$
(1004)

where $\mu$ is the normalized Haar measure satisfying $\mu(\mathrm{O}(n))=1$.
Because the Haar integral over linear functions vanishes, i.e.,
$\int_{\mathrm{O}(n)}d\mu(V)\,V_{ij}=0$ [collins2006integration], it follows
that

$\int_{\mathrm{O}(n)}d\mu(V)\,V^{\mathsf{T}}Q(\sigma_{v})=0.$ (1005)

Furthermore, because $\mathrm{i}\widetilde{\gamma}_{i}^{(v)}\gamma_{j}^{(v)}$
are linearly independent (as elements of an operator basis), the equality
$|\psi\rangle\\!\langle\psi|=\int_{\mathrm{O}(n)}d\mu(V)|\psi(V)\rangle\\!\langle\psi(V)|$
implies that

$Q(\sigma_{v})=\int_{\mathrm{O}(n)}d\mu(V)\,V^{\mathsf{T}}Q(\sigma_{v})=0$
(1006)

for all $v\in[m]$. ∎

### The Pin group from quantum circuits

In the main text we showed that each $x\in\mathrm{Pin}(n)$ corresponds to the
eigenstates of a family of free-fermion Hamiltonians, which are (pure)
fermionic Gaussian states. Here we provide an alternative perspective of this
correspondence, using quantum circuits which prepare such states.

Recall that every $x\in\mathrm{Pin}(n)$ can be written as

$x=u_{1}\cdots u_{k}$ (1007)

for some $k\leq n$, where we may expand each $u_{j}\in S^{n-1}$ in the
standard basis as

$u_{j}=\sum_{i\in[n]}v_{i}^{(j)}e_{i}$ (1008)

for some unit vector $v^{(j)}\in\mathbb{R}^{n}$. The product of these unit
vectors can be expressed using the right-multiplication operator
$\rho_{u_{j}}$ acting on the identity element,

$\begin{split}x&=e_{\varnothing}x\\\ &=e_{\varnothing}u_{1}\cdots u_{k}\\\
&=(\rho_{u_{k}}\cdots\rho_{u_{1}})(e_{\varnothing})\\\ \end{split}$ (1009)

On the other hand, consider the so-called Clifford loader
[kerenidis2022quantum], a circuit primitive defined (in our notation) as

$\Gamma(v)=\sum_{i\in[n]}v_{i}\gamma_{i}$ (1010)

for any unit vector $v\in\mathbb{R}^{n}$. It is straightforward to check that
this operator is Hermitian and unitary, and Ref. [kerenidis2022quantum]
provides an explicit circuit constructions based on two-qubit Givens rotation
primitives.575757Givens rotations themselves are representations of fermionic
Gaussian transformations acting on two modes at a time. Using the relation
$\gamma_{i}=\rho_{i}\alpha$ and acting this circuit on the vacuum state
$|0^{n}\rangle\equiv|e_{\varnothing}\rangle$, we see that the state

$\begin{split}|x\rangle&=\Gamma(v^{(k)})\cdots\Gamma(v^{(1)})|0^{n}\rangle\\\
&=(\rho_{u_{k}}\alpha\cdots\rho_{u_{1}}\alpha)|e_{\varnothing}\rangle\\\
&=(-1)^{\binom{k}{2}}(\rho_{u_{k}}\cdots\rho_{u_{1}})|e_{\varnothing}\rangle\end{split}$
(1011)

indeed is equivalent to $x\in\mathrm{Pin}(n)$, up to a global sign (recall
that $\alpha^{2}=\mathbb{I}$ and
$\alpha|e_{\varnothing}\rangle=|e_{\varnothing}\rangle$). In other words, the
Clifford loader is precisely the quantum-circuit representation of generators
of the $\mathrm{Pin}$ group.

It is worth noting that in Ref. [kerenidis2022quantum] they construct
“subspace states” from this composition of Clifford loaders. In the language
of fermions, subspace states are Slater determinants: free-fermion states with
fixed particle number. Preparing Slater determinants in this fashion requires
that the unit vectors $v^{(1)},\ldots,v^{(k)}$ be linearly independent (and
thus, without loss of generality, they can be made orthonormal while
preserving the subspace that they span, hence the alternative name). However,
the definition of the Pin group demands all possible unit vectors in such
products, not just those which are linearly independent. Indeed, one can see
that if the state $|x\rangle$ is a Slater determinant, then its trace is an
integer, as

$\operatorname{tr}[Q(x)]=\langle
x|\sum_{i\in[n]}\mathrm{i}\widetilde{\gamma}_{i}\gamma_{i}|x\rangle=\langle
x|(n\mathbb{I}_{2^{n}}-2N)|x\rangle=n-2k\in\\{-n,\ldots,n\\},$ (1012)

where $N=\sum_{i\in[n]}a_{i}^{\dagger}a_{i}$ is the total number operator.
Clearly not all orthogonal matrices have integer trace, so Slater determinants
are insufficient to cover all of $\mathrm{Pin}(n)$. To reach the remaining
elements, we note that if the unit vectors are linearly dependent, then one
can show that the
$|x\rangle=\Gamma(v^{(k)})\cdots\Gamma(v^{(1)})|0^{n}\rangle$ does not have
fixed particle number, so $\langle x|N|x\rangle$ is not necessarily an
integer.

### Measurement schemes

In this section we comment on the efficient schemes available for measuring
the relevant expectation values. This is important even in the context of a
phase-estimation approach, as one needs to obtain the values of the decision
variables to perform the rounding procedure.

#### Tomography of edge marginals

To measure the energy (for variational approaches) or to perform edge
rounding, we require the expectation values of the two-body observables

$\displaystyle\Gamma_{ij}^{(u,v)}$
$\displaystyle=\sum_{k\in[n]}P_{ik}^{(u)}\otimes P_{jk}^{(v)},\quad
G=\mathrm{O}(n),$ (1013) $\displaystyle\widetilde{\Gamma}_{ij}^{(u,v)}$
$\displaystyle=\sum_{k\in[n]}\widetilde{P}_{ik}^{(u)}\otimes\widetilde{P}_{jk}^{(v)},\quad
G=\mathrm{SO}(n),$ (1014)

for each $(u,v)\in E$ and $i,j\in[n]$. When considering $G=\mathrm{O}(n)$,
because each $P_{ik}^{(u)}\otimes P_{jk}^{(v)}$ is a fermionic two-body
operator, we can straightforwardly apply the partial tomography schemes
developed for local fermionic systems, such as Majorana swap networks
[bonet2020nearly] or classical shadows [zhao2021fermionic]. In either case,
the measurement circuits required are fermionic Gaussian unitaries and the
sample complexity is $\mathcal{O}(N^{2}/\epsilon^{2})$, where $N=n|V|$ is the
total number of qubits and $\epsilon>0$ is the desired estimation precision of
each expectation value.

#### Tomography of vertex marginals

The vertex-rounding procedure requires the expectation values of only single-
qudit observables $P_{ij}^{(v)}$ or $\widetilde{P}_{ij}^{(v)}$ on each vertex
$v\in V$, $i,j\in[n]$. In this case the observables being measured commute
across vertices, so it suffices to talk about the tomography of a single
vertex, as the same process can be executed in parallel across all vertices.
Again, because these operators are fermionic one-body observables, the same
fermionic partial tomography technology [bonet2020nearly, zhao2021fermionic]
can be applied here, incurring a sampling cost of
$\mathcal{O}(n/\epsilon^{2})$. In fact, further constant-factor savings can be
achieved in the one-body setting by using the measurement scheme introduced in
Ref. [arute2020hartree]. This scheme requires only particle-conserving
fermionic Gaussian unitaries, which can be compiled with only half the depth
of the more general Gaussian unitaries required of the previous two methods.
Note that each operator $P_{ij}$ is of the form of either $XX$ or $YY$ when
$|i-j|=1$, and so they correspond precisely to the observables measured to
reconstruct the real part of the fermionic one-body reduced density matrix
[arute2020hartree].

#### Estimating observables via gradient method

Ref. [huggins2021nearly] introduces a quantum algorithm for estimating a large
collection of (generically noncommuting) $M$ observables $\\{O_{j}\mid
j\in[M]\\}$ to precision $\epsilon$ by encoding their expectation values into
the gradient of a function. This function is implemented as a quantum circuit
which prepares the state of interest and applies
$\tilde{\mathcal{O}}(\sqrt{M}/\epsilon)$ gates of the form
$c\text{-}e^{-\mathrm{i}\theta O_{j}}$, controlled on
$\mathcal{O}(M\log(1/\epsilon))$ ancilla qubits. Finally, using the algorithm
of Ref. [gilyen2019optimizing] for gradient estimation, one calls this circuit
$\tilde{\mathcal{O}}(\sqrt{M}/\epsilon)$ times to estimate the encoded
expectation values (the notation $\tilde{\mathcal{O}}(\cdot)$ suppresses
polylogarithmic factors). Although this approach demands additional qubits and
more complicated circuitry, it has the striking advantage of a quadratically
improved scaling in the number of state preparations with respect to
estimation error $\epsilon$, compared to the refinement of sampling error in
tomographic approaches. In our context, we have either $M=n^{3}|E|$ or
$M=n^{2}|V|$ observables of interest (satisfying a technical requirement of
having their spectral norms bounded by $1$), corresponding to the measurement
of edge or vertex terms respectively. The gates required are then simply
controlled Pauli rotations.

## Conclusion

The theory of quantum computation has sparked deep and profound questions
about the nature of computation and the computation of Nature. For computer
scientists, it has forced a reevaluation of what a “reasonable” computational
model means within the extended Church–Turing thesis [deutsch1985quantum,
bernstein1993quantum]. Physicists, meanwhile, have had to confront with the
possibility that their models of Nature might not properly reflect what she is
actually capable of [kempe2006complexity, poulin2011quantum].

Fortunately, quantum computers also offer to help us solve some difficult
problems that we face in concrete domains. This dissertation has addressed
this more down-to-earth affair, focusing on improvements to quantum algorithms
for simulating fermions so that we might realize that promise sooner rather
than later. As is always the case, we are left with more questions than when
we had started. In this chapter, I will opine on some of the most outstanding
questions that remain following the results presented within this
dissertation.

In Fermionic Partial Tomography via Classical Shadows, we developed an optimal
extension of classical shadows for fermionic systems, and in Group-theoretic
Error Mitigation Enabled by Classical Shadows and Symmetries we studied how
such protocols enjoy inherent robustness properties in the presence of quantum
symmetries. The motivation for these works is to squeeze out every drop of
performance from NISQ machines, from an algorithmic-design perspective. While
motivated for use in the NISQ era, I am interested in seeing how these results
can find applications more broadly. For example, in one collaboration that I
was involved with, we studied the efficiency of classical shadows within the
fault-tolerant regime [babbush2023quantum]. There, we considered simulating
the time evolution of interacting electrons, wherein partial-tomography
techniques are necessary to extract dynamical properties from the system (as
opposed to techniques like phase estimation, which primarily address static
properties). What can be said of enhancing classical shadows for fault-
tolerant applications and dynamics simulations? Impressive work on the former
aspect has already been accomplished by van Apeldoorn _et al._
[van2023quantum], showing how shadow tomography can be substantially improved
with access to state-preparation circuits. In particular, they achieve
Heisenberg-limited estimation of multiple observables, albeit with substantial
space complexity. To the latter aspect, Rall [rall2020quantum] showed how to
learn time-correlation functions using post-NISQ resources. Such properties
are both essential for studies of quantum many-body physics, yet highly
challenging to extract from a quantum computer. How much further can classical
shadows/shadow tomography be pushed, to offer practical solutions for these
problems?

More generally, the promise of NISQ appears to be wavering at the time of this
writing [obrien2023purification]. Nonetheless, this line of research continues
to provide us with interesting questions that we had not even thought to ask
before the advent of such machines. For example, a running theme of NISQ
algorithms has been to explore how much complexity can be offloaded onto a
classical processor before encountering exponential costs. The initial
formulation, VQE [peruzzo2014variational], supposed that the quantum computer
needs only to prepare a low-energy state; determining the circuit for doing so
can be offloaded to a classical processor. However, the hybrid
“quantum–classical quantum Monte Carlo” (QC-QMC) algorithm introduced by
Huggins _et al._ [huggins2022unbiasing] appears to challenge that notion.
Their experimental success is not necessarily attributed to how low energy of
a state that their noisy quantum computer produced, but instead how closely
its sign structure (relative phases between amplitudes of the wavefunction)
matched that of the true ground state. This is because a significant amount of
the computational hardness in QMC manifests as the so-called sign problem
[loh1990sign], whereby quantum negativity obstructs the convergence of Monte
Carlo simulations. With knowledge of the sign structure, the QMC algorithm
(run entirely classically) can be modified to avoid the sign problem.

This new perspective leads us to ask how we might further ease the job of the
NISQ machine in practice. To this end, we need to pinpoint precisely where the
hardness manifests within various classical simulation algorithms. While the
complexity-theoretic hardness is equivalent in all forms, it appears that some
tasks (perhaps like approximating the sign structure) might be robust than
others (like directly preparing an approximate ground state). As we move
toward larger and more accurate machines, can the hope for noisy quantum
information processing be rekindled? For example, suppose one is given a
device with below-threshold error rates. Such a machine could produce a
handful of perfect logical qubits, an impressive technological feat but still
easily classically simulated. On the other hand, using the physical qubits
_directly_ might grant access to ${\sim}10^{4}$ very-low-noise qubits. This
regime might not necessarily be “intermediate scale” anymore, but nonetheless
might be where the hidden promise of NISQ lies.

Finally, in Expanding the Reach of Quantum Optimization with Fermionic
Embeddings we investigated a hard classical optimization problem and showed
how to naturally embed it into a quantum (indeed, fermionic) problem, a
procedure referred to as “quantum relaxation.”585858To my knowledge, this
terminology/concept has only been considered by Fuller _et al._
[fuller2021approximate] and us [zhao2023expanding], although I believe it is a
very natural generalization of relaxations from classical approximation theory
that merits wider study. This work took the initial steps of this subject, and
we have already mentioned the many potential further directions to pursue at
the end of that chapter. Here, I will recap the vision I have for this line of
work with regards to approximation theory. For instance, Naor, Regev, and
Vidick [naor2014efficient] constructed an efficient algorithm for the
noncommutative Grothendieck problem which achieves an approximation ratio of
$1/2$ in the complex case, but only $1/2\sqrt{2}$ in the real case. On the
other hand, Briët, Regev, and Saket demonstrated an $\mathsf{NP}$-hard
threshold of approximation of $1/2$ for both cases. The real variant is
precisely what our fermionic representation addresses; can it then provide the
tools for designing a classical algorithm which improves the $1/2\sqrt{2}$
result to match the hardness threshold?

More broadly, what can be said about the hardness of the $\mathrm{SO}(N)$
problem, which has yet to be properly studied from a complexity-theoretic
perspective? Building off the work of Saunderson, Parrilo, and Willsky
[saunderson2014semidefinite], our results showed how to encode the determinant
condition into a linear constraint on an exponentially large vector space.
Ideally, this representation could elucidate rigorous statements about the
hardness of quadratic $\mathrm{SO}(N)$ optimization. In connection with
quantum algorithms, I would ideally like to see approximation thresholds in
terms of $\mathsf{QMA}$-hardness, which would establish the limits of
approximation even with quantum computers. Unfortunately, this avenue
currently seems intractable without an established quantum PCP theorem
[aharonov2013guest].


# Real-time quantum state estimation in circuit QED via Bayesian approach

Yang Yang  Beili Gong  Wei Cui<EMAIL_ADDRESS>School of Automation
Science and Engineering, South China University of Technology, Guangzhou
510641, China

###### Abstract

Using a circuit QED device, we present a theoretical study of real-time
quantum state estimation via quantum Bayesian approach. Suitable conditions
under which the Bayesian approach can accurately update the density matrix of
the qubit are analyzed. We also consider the correlation between some basic
and physically meaningful parameters of the circuit QED and the performance of
the Bayesian approach. Our results advance the understanding of quantum
Bayesian approach and pave the way to study quantum feedback control and
adaptive control.

###### pacs:

42.50.Pq, 42.50.Dv, 85.25.-j

##  I Introduction

Quantum state estimation that aims at determining the states of the quantum
systems using the measurement records is instrumental to quantum physics
experiments and has attracted much attentions Bisio2009 ; Doherty2000 ;
Ladd2010 . It can be divided into two categories: quantum state tomography and
real-time quantum state estimation. Quantum state tomography completes the
task of reconstructing the initial states of the systems Bisio2009 . There
have already been lots of methods applied to achieve efficient quantum state
tomography, including maximum likelihood method Six2016 , linear regression
method Qi2013  and compressive sensing method Gross2009 ; Flammia2012 ;
Kalev2015 . Different with the classical worlds, even though we adopt
continuous weak measurements, the quantum states after measurements are still
different with the original states due to the back-action of the measurements,
the spontaneous decay and dephasing of the systems Braginsky1992 ; Wiseman2009
. Thus, the real-time monitoring of the states is widely demanded as well, in
which the filter equation Wiseman2009 ; Tanaka2012  and the Bayesian approach
Korotkov1999 ; Korotkov2016  are the two most widely used methods.

Comparing with other methods, the Bayesian approach mainly advantages in its
great rapidity and low computational complexity. Many quantum experiments
based on the Bayesian approach have been realized Weber2014 ; Vijay2012 ;
Gillett2010 . The core content of the quantum Bayesian estimation is about
updating the diagonal elements of the qubitKorotkov1999 , and the methods to
update the non-diagonal elements were also raised recently Wang2014 ; Feng2016
; Qin2017 ; Korotkov2016 . Until now, the Bayesian approach has already been
enriched to be an effective method that can be used as a useful tool in
quantum information technology and quantum control Zhang2017 . However, there
is no strict proof of its establishment. Our manuscript is going to consider
the conditions under which the Bayesian approach can real-time precisely
monitor the quantum states, that is, to analyze the correlation between the
values of the system parameters and the performance of the Bayesian approach
numerically.

This paper is organized as follows. In Sec. II, we briefly review the
dispersive measurement in circuit QED. In Sec. III, we present the quantum
Bayesian approach as well as examining its estimation performance. In Sec. IV,
we explore the conditions under which the method can real-time precisely track
the evolution of the qubit numerically. We plot the error curves of the system
parameters and find out all the existing critical points. We summarize our
conclusions in Sec. V.

##  II The model

Superconducting circuit QED has been regarded as a promising candidate for the
realization of scalable quantum computing because of its controllability and
design flexibility Blais2004 . Moreover, it is also a great platform to test
quantum feedback control and quantum estimation Xiang2013 ; Slichter2012 ;
Hatridge2013 . Shortly after being proposed in Ref. Blais2004 , it was
implemented experimentally Wallraff2004 ; Chiorescu2004 . It provides us with
several simple high-fidelity readout mechanisms. For example one can realize
the continuous weak measurement of the system by performing homodyne
measurement on the cavity field.

The circuit-QED system consists of a superconducting qubit and a microwave
cavity Blais2004 . One can describe the qubit, the cavity and their coupling
with the Jaynes-Cummings Hamiltonian. In order to realize the measurement of
the qubit, the cavity is driven with a tone of amplitude $\varepsilon_{d}$ and
frequency $\omega_{d}$. The Hamiltonian of the whole system is given by
($\hbar=1$).

$\displaystyle H=$
$\displaystyle\omega_{r}a^{\dagger}a+\omega_{q}\sigma_{z}/2+g(a\sigma_{+}+a^{\dagger}\sigma_{-})$
$\displaystyle+\varepsilon_{d}(t)e^{-i\omega_{d}t}a^{\dagger}+\varepsilon^{*}_{d}(t)e^{i\omega_{d}t}a,$
(1)

where $\omega_{r}$ is the cavity frequency, $\omega_{q}$ is the qubit
transition frequency, and $g$ is the coupling strength. We consider the system
in a dispersive regime, $|\Delta|=|\omega_{r}-\omega_{q}|\gg|g|$. Under this
condition and in the rotating frame, we can use the effective Hamiltonian to
describe the system Gambetta2006 .

$\displaystyle H_{\rm eff}=\Delta_{r}a^{{\dagger}}a+{{\tilde{\omega}_{q}}\over
2}\sigma_{z}+\chi
a^{\dagger}a\sigma_{z}+\varepsilon_{d}(t)a^{\dagger}+\varepsilon^{*}_{d}(t)a,$
(2)

where $\Delta_{r}=\omega_{r}-\omega_{d}$ and
$\tilde{\omega}_{q}=\omega_{q}+\chi$ with $\chi=g^{2}/\Delta$ representing the
dispersive coupling strength between the cavity photon number and the qubit.

We use the following quantum stochastic master equation to describe the
evolution of the qubit-cavity state $\rho_{T}(t)$ Wiseman2009 ,

$\displaystyle\dot{\rho}_{T}(t)=$ $\displaystyle-{i\over\hbar}[H_{\rm
eff},\rho_{T}(t)]+\kappa D[a]\rho_{T}(t)+\gamma_{1}D[\sigma]\rho_{T}(t)$
$\displaystyle+\gamma_{\varphi}D[\sigma_{z}]\rho_{T}(t)/2+\sqrt{\kappa}H[ae^{-i\phi}]\rho_{T}(t)\xi(t),$
(3)

where $\kappa$, $\gamma_{1}$ and $\gamma_{\varphi}$ represent the damping rate
of the resonator, qubit decay and dephasing of the qubit, respectively.
$\xi(t)$ is a Gaussian white noise coming from the stochastic quantum-jump,
and satisfies $E\left[\xi(t)\right]=0$ and
$E\left[\xi(t)\xi(t^{\prime})\right]=\delta(t-t^{\prime})$. The superoperators
$D[c]\rho=c\rho c^{\dagger}-c^{\dagger}c\rho/2-\rho c^{\dagger}c/2$ and
$H[c]\rho=c\rho+\rho c^{\dagger}-{\rm
Tr}\left\\{\left[c+c^{\dagger}\right]\rho\right\\}\rho$. We simply consider
the single quadrature measurement, through which what we actually measure is a
photocurrent Gambetta2008 . The photocurrent measured can be expressed as:

$\displaystyle I_{\phi}(t)=\sqrt{\kappa\eta}{\left\langle
ae^{-i\phi}+a^{\dagger}e^{i\phi}\right\rangle}_{t}+\xi(t).$ (4)

With the initial state prepared by $\left|e(g)\right\rangle$, we have the
coherent state of the cavity field $\left|\alpha_{e(g)}\right\rangle$ whose
amplitude evolution can be described with the following equation Gambetta2006
.

$\displaystyle\dot{\alpha}_{e(g)}=-i\varepsilon_{d}(t)-i(\Delta_{r}\pm\chi)\alpha_{e(g)}-\kappa\alpha_{e(g)}/2.$
(5)

Ref. Gambetta2008  defined
$\beta(t)=\alpha_{e}(t)-\alpha_{g}(t)=|\beta|e^{-i\theta_{\beta}}$ to
represent the information about the qubit involved in the measurement
outcomes.

Figure 1: (Color online) The measured signal (a) and the filtered signal (b)
of a single measurement in circuit QED with $t_{m}=0.01{\rm\mu s}$,
$\omega_{q}=0$ Hz, $\Delta_{r}=0$ Hz, $\chi=1$ MHz, $\varepsilon_{d}=15$ MHz,
$\kappa=20$ MHz, $\alpha_{e}=\alpha_{g}=0$, $\gamma_{1}=\gamma_{\psi}=0$ and
$\phi=\pi$.

What one finally want to control is the qubit, so it’s necessary to find an
equation that gives the evolution of the qubit itself instead of the qubit-
cavity system. Through applying a displacement transformation on the Eq. (3)
and tracing over the cavity state, one can obtain an effective stochastic
master equation (SME) for the qubit Gambetta2008 :

$\displaystyle\dot{\rho}(t)=$ $\displaystyle-i{\tilde{\omega}_{q}+B(t)\over
2}[\sigma_{z},\rho(t)]+{\Gamma_{d}(t)+\gamma_{\varphi}\over
2}D[\sigma_{z}]\rho(t)$
$\displaystyle+\gamma_{1}D[\sigma_{-}]\rho(t)-\sqrt{\Gamma_{ci}(t)}M[\sigma_{z}]\rho\xi(t)$
$\displaystyle+i{\sqrt{\Gamma_{ab}(t)}\over 2}[\sigma_{z},\rho]\xi(t),$ (6)

where $B(t)=2\chi\rm{Re}[\alpha_{g}(t)\alpha_{e}(t)^{*}]$ describes the
generalized ac-Stark shift of the qubit energy due to the dynamic fluctuation
of cavity field. The superoperator $M[c]\rho=\left[c-{\rm
Tr}(c\rho)\right]\rho/2+\rho\left[c-{\rm Tr}(c\rho)\right]/2$. The effect of
the coupling with a cavity on the qubit is transformed to an additional
dephasing rate $\Gamma_{d}=2\chi\rm{Im}[\alpha_{g}(t)\alpha_{e}(t)^{*}]$.
$\Gamma_{ci}(t)=\kappa|\beta(t)|^{2}\rm{cos}^{2}(\phi-\theta_{\beta})$ and
$\Gamma_{ab}(t)=\kappa|\beta(t)|^{2}\rm{sin}^{2}(\phi-\theta_{\beta})$
respectively represent the measurement efficiency and the measurement back-
action rate for a single measurement Gambetta2008 . In this reduced
representation, the measured photocurrent can also be expressed as

$\displaystyle I(t)$
$\displaystyle=\sqrt{\Gamma_{ci}}{\left\langle\sigma_{z}\right\rangle}_{t}+\xi(t)+\sqrt{\kappa}|\mu(t)|\rm{cos}(\theta_{\mu}-\phi)$
$\displaystyle=\bar{I}(t)+\sqrt{\kappa}|\mu(t)|\rm{cos}(\theta_{\mu}-\phi),$
(7)

where $\mu(t)=\alpha_{e}(t)+\alpha_{g}(t)=|\mu|e^{-i\theta_{\mu}}$. Once we
obtain the measurement records during a time interval and want to estimate the
system state in the next time with the initial qubit state known, the common
method is the filter equation

$\displaystyle\dot{\varrho}(t)=$ $\displaystyle-i{\tilde{\omega}_{q}+B(t)\over
2}[\sigma_{z},\varrho(t)]+{\Gamma_{d}(t)+\gamma_{\varphi}\over
2}D[\sigma_{z}]\varrho(t)$
$\displaystyle-\sqrt{\Gamma_{ci}(t)}M[\sigma_{z}]\varrho(\bar{I}(t)-\sqrt{\Gamma_{ci}(t)}{\left\langle\sigma_{z}\right\rangle}_{t})$
$\displaystyle+i{\sqrt{\Gamma_{ab}(t)}\over
2}[\sigma_{z},\varrho](\bar{I}(t)-\sqrt{\Gamma_{ci}(t)}{\left\langle\sigma_{z}\right\rangle}_{t})~{}$
$\displaystyle+\gamma_{1}D[\sigma_{-}]\rho(t).$ (8)

However, it is time-consuming and computationally complex to calculate the
evolution of the filter equation, while the time consumed to real-time
estimation of the quantum states is hoped to be as short as possible in
quantum feedback control. Thus, Korotkov put forward the quantum Bayesian
approach Korotkov1999 ; Korotkov2016  to deal with this problem.

Figure 2: (Color online) Real-time quantum state estimation via quantum
Bayesian approach. Figs. 2(a)-2(c) detailly depict the estimation of the three
Bloch components $\langle\hat{\sigma}_{i}\rangle$ with quantum Bayesian
approach (solid blue curves) and the original evolution of
$\langle\sigma_{i}\rangle$ (dash red curves) under the standard SME (6). Fig.
2(d) depicts the trajectories in the Bloch sphere, where the blue curve is the
estimated trajectory and the red curve is the evolution one under the standard
SME (6). To make it easier to observe, we only intercept the first $0.75\rm\mu
s$ of the trajectory. Here $\eta=0.1$ and the rest parameters are the same as
those in Fig. 1.

##  III Bayesian approach for single quadrature measurement

As we mentioned in the previous section, what we actually measure in the
circuit-QED is a photocurrent that satisfies Eq. (4). Due to the noise effect,
we apply a low-pass filter process to the photocurrent and average the
measured output signal over a time interval $t_{m}$,

$\displaystyle{I_{m}}=\frac{1}{{{t_{m}}}}\int_{0}^{{t_{m}}}{{I_{\phi}}(t)dt}.$
(9)

When the initial state is determined, based on Eq. (5) we can calculate the
evolution of the amplitude of the cavity field. Thus we have the ensemble-
average filtered output

$\displaystyle{{\bar{I}}_{g(e)}}=\frac{1}{{{t_{m}}}}\int_{0}^{{t_{m}}}{(2\sqrt{\kappa})\operatorname{Re}[{\alpha_{g(e)}}(t){e^{-i\phi}}]}dt.$
(10)

The measurement output can be proved to satisfy a determined Gaussian
distributionKorotkov1999 ; Wang2014 :

$\displaystyle{P_{g(e)}}({I_{m}})=\frac{1}{{\sqrt{2\pi
D}}}\exp[-{({I_{m}}-{\bar{I}_{g(e)}})^{2}}/2D_{I}],$ (11)

where $D_{I}=1/t_{m}$ characterizes the variance of the distribution. Here we
plot the measured signal and filtered signal of an individual measurement in
Fig.1.

As the diagonal elements of the qubit density matrix describe the probability
distribution of the system state, they can be regarded as the priori
knowledge. With the determined distribution Eq. (11) and the measurement
output, the classical Bayesian formula is proposed to update the diagonal
elements of the quantum state.

$\displaystyle{\hat{\rho}_{gg}}({t_{m}})={\rho_{gg}}(0){P_{g}}({I_{m}})/N,$
$\displaystyle{\hat{\rho}_{ee}}({t_{m}})={\rho_{ee}}(0){P_{e}}({I_{m}})/N,$
$\displaystyle
N={\rho_{gg}}(0){P_{g}}({I_{m}})+{\rho_{ee}}(0){P_{e}}({I_{m}}).$ (12)

If the qubit is a pure state, one can approximate $\rho_{ge}(t_{m})$ as
Korotkov1999

$\displaystyle{\tilde{\rho}_{ge}}({t_{m}})={\rho_{ge}}(0){e^{-i{{\tilde{\omega}}_{q}}{t_{m}}}}\sqrt{{P_{g}}({I_{m}}){P_{e}}({I_{m}})}/N.$
(13)

However this result was rather imprecise and quite different with the filter
equation (8) due to the purity degradation and stochastic fluctuations of the
systems. Thus, in order to further correct the estimation of the non-diagonal
elements, Ref. Korotkov2016  detailly analyzed the measurement-induced
backaction and found that besides the diagonal part the phase backaction, the
dephasing due to nonideality of the measurement and ac Stark shift of the
qubit frequency may play more significant role in producing evolution of the
non-diagonal elements. Based on the stochastic master equation (6), Refs.
Wang2014 ; Feng2016  further considered the cavity-field-fluctuation effects
on the qubit and proposed a different method to accurately estimate the non-
diagonal matrix elements. The updated quantum Bayesian approach can be
summarized as follows

$\displaystyle{\hat{\rho}_{ge}}({t_{m}})={{\tilde{\rho}}_{ge}}({t_{m}})De(t_{m})\exp\\{-i[{\Phi_{1}}({t_{m}})+{\Phi_{2}}({t_{m}})]\\},$
$\displaystyle
De(t_{m})=\left|\langle{{\alpha_{e}}({t_{m}})}|{{\alpha_{g}}({t_{m}})}\rangle\right|,$
$\displaystyle{\Phi_{1}}({t_{m}})=\int_{0}^{{t_{m}}}{B(t)dt},$
$\displaystyle{\Phi_{2}}({t_{m}})=-\int_{0}^{{t_{m}}}{\sqrt{{\Gamma_{ba}}(t)}\bar{I}_{\phi}(t)dt},$
$\displaystyle\bar{I}_{\phi}(t)=I_{\phi}(t)-\sqrt{\kappa}|\mu|\cos(\theta_{\mu}-\phi),$
(14)

where $De$ is the purity degradation factor, $\Phi_{1}$ and $\Phi_{2}$ are the
two additional phase factors.

With Eqs. (12) and Eqs. (14) we can real-time estimate the quantum state. In
Fig. 2, we plot the the estimation of the three Bloch components
$\langle\hat{\sigma}_{i}\rangle$ with Bayesian approach and the evolution of
$\langle\sigma_{i}\rangle$. Moreover the estimated trajectory and the
evolution trajectory are depicted in Fig. 2(d), from which it is clear that
the Bayesian approach can real-time accurately estimate the quantum state of
the circuit QED system. To make it easier to observe, we only intercept the
first 0.75$\mu s$ of the trajectories. Here $\eta=0.1$ and the rest parameters
are the same as those in Fig. 1.

##  IV Suitable conditions of the Bayesian approach

Although the quantum Bayesian approach provides us with great rapidity and low
computational complexity, there is no theoretical proof about its
establishment, and it still remains an assumption. So exploring the suitable
conditions of the quantum Bayesian approach might be a critical issue. For
simplicity, we consider some basic and physical meaningful parameters involved
in the experiment stated above and concentrate on finding out the ranges of
these parameters in which the Bayesian approach can work properly. The reduced
master equation Eq. (6) is regarded as the standard evolution of the qubit. It
should be noted that the system works in the dispersive regime,
$|\Delta|=|\omega_{r}-\omega_{q}|\gg|g|$. This regime limits $\chi$ to be much
smaller than $|\Delta|$ (usually approximately equal to $\omega_{q}$).

Figure 3: (Color online) Evolution of the error functions $E_{\gamma_{1}}$,
$Ez_{\gamma_{1}}$, $Ey_{\gamma_{1}}$ and $Ex_{\gamma_{1}}$ for a total time of
10$\mu s$ and with 5000 trajectories, as functions of the parameter
$\gamma_{1}$. In Fig. 3(a) the blue solid curve, the red dashed curve, the
yellow pecked curve and the purple dotted curve represent $E_{\gamma_{1}}$
$Ez_{\gamma_{1}}$, $Ey_{\gamma_{1}}$ and $Ex_{\gamma_{1}}$, respectively. The
red filled dot is the critical point of $\gamma_{1}$ at which $E_{\gamma_{1}}$
reaches 0.1. For clarity, the time evolution of the estimated
$\langle\hat{\sigma}_{z}\rangle$ (solid blue curve) and the evolution of
$\langle\sigma_{z}\rangle$ (dash red curve) at this point are plotted in Fig.
3(b).  Figure 4: (Color online) Evolutions of the error functions $E_{a}$,
$Ex_{a}$, $Ey_{a}$ and $Ez_{a}$ for a total time of 10$\mu s$ and with 5000
trajectories, as functions of the parameters $\chi$, $\Delta_{r}$,
$\varepsilon_{d}$, $\eta$, $\gamma_{\varphi}$ and $\kappa$ in Figs. 4(a)-4(f),
respectively. The blue solid curves, the red dashed lines, the yellow pecked
curves and the purple dotted curves represent $E_{a}$ $Ez_{a}$, $Ey_{a}$ and
$Ex_{a}$ respectively. The red filled dot in each subfigure represents at
which the overall error $E_{a}$ reaches 0.1.

The basic and physical meaningful parameters in the circuit-QED system we
considered are $\gamma_{1}$, $\gamma_{\varphi}$, $\chi$, $\varepsilon_{d}$,
$\Delta_{r}$, $\kappa$ and $\eta$. For a specific parameter $a$, we
characterize the overall performance of the Bayesian approach by defining an
error $E_{a}$, and the performances of the three Bloch components
($\langle\sigma_{x}\rangle$, $\langle\sigma_{y}\rangle$,
$\langle\sigma_{z}\rangle$ ) estimation are defined as $Ex_{a}$, $Ey_{a}$, and
$Ez_{a}$, respectively. Due to the dissipation property of the system, the
ensemble evolution of the master equation (6) would tend to be a steady state
$\rho^{\ast}$, or ($\langle\sigma_{x}^{\ast}\rangle$,
$\langle\sigma_{y}^{\ast}\rangle$, $\langle\sigma_{z}^{\ast}\rangle$). To make
the performance analysis more accurate, we consider the error functions
$Ex_{a}$, $Ey_{a}$, and $Ez_{a}$ are defined in an unstable time interval

$\displaystyle S_{i}$ $\displaystyle=$
$\displaystyle\bigg{\\{}t:\bigg{|}\langle\sigma_{i}(t)\rangle-\langle\sigma_{i}^{\ast}(t)\rangle\bigg{|}>\epsilon_{1}\bigg{\\}}$
$\displaystyle\cup\bigg{\\{}t:\bigg{|}\langle\hat{\sigma}_{i}(t)\rangle-\langle\sigma_{i}^{\ast}(t)\rangle\bigg{|}>\epsilon_{1}\bigg{\\}},~{}~{}~{}~{}~{}~{}i=x,~{}y,~{}z,$

where $\epsilon_{1}$ is a steady-state precision. Throughout this manuscript
we set $\epsilon_{1}=1\%$. For a given parameter $a$, the error functions
$Ex_{a}$, $Ey_{a}$, and $Ez_{a}$ are defined as

$\displaystyle
Ex_{a}=\int_{S_{x}}\frac{[\langle\sigma_{x}(t)\rangle-\langle\hat{\sigma}_{x}(t)\rangle]^{2}}{\parallel
S_{x}\parallel}dt,$ $\displaystyle
Ey_{a}=\int_{S_{y}}\frac{[\langle\sigma_{y}(t)\rangle-\langle\hat{\sigma}_{y}(t)\rangle]^{2}}{\parallel
S_{y}\parallel}dt,$ (16) $\displaystyle
Ez_{a}=\int_{S_{z}}\frac{[\langle\sigma_{z}(t)\rangle-\langle\hat{\sigma}_{z}(t)\rangle]^{2}}{\parallel
S_{z}\parallel}dt.$

Thus, we can calculate the overall error for a specific parameter $a$

$\displaystyle
E_{a}=\alpha_{x}{Ex}_{a}+\alpha_{y}{Ey}_{a}+\alpha_{z}{Ez}_{a},$ (17)

where $\alpha_{x}$, $\alpha_{y}$ and $\alpha_{z}$ are the weighting factors
used to achieve a balance between the three Bloch components. It is clear that
the diagonal elements play a more important role in the quantum Bayesian
approach. Thus, throughout this manuscript we set
$\alpha_{x}=\alpha_{y}=0.25$, and $\alpha_{z}=0.5$.

Table 1: The suitable conditions of the system parameters ($E_{a}\leq 0.1$) Parameter |  $\chi$ |  $\Delta_{r}$ |  $\varepsilon_{d}$ |  $\eta$ |  $\gamma_{1}$ |  $\gamma_{\varphi}$ |  $\kappa$
---|---|---|---|---|---|---|---
Suitable conditions |  $\leq 19.57$MHz |  $\leq 0.5949$MHz |  $\leq 70.11$MHz |  $\leq 0.4817$ |  $\leq 5.841$kHz |  $\leq 0.2036$MHz |  -
Table 2: The suitable conditions of the system parameters ($E_{a}\leq 0.05$) Parameter |  $\chi$ |  $\Delta_{r}$ |  $\varepsilon_{d}$ |  $\eta$ |  $\gamma_{1}$ |  $\gamma_{\varphi}$ |  $\kappa$
---|---|---|---|---|---|---|---
Suitable conditions |  $\leq 12.1$MHz |  $\leq 0.1959$MHz |  $\leq 53.4$MHz |  $\leq 0.2372$ |  $\leq 1.821$kHz |  $\leq 0.0467$MHz |  -

In order to maintain the estimation accuracy one can set an error upper-bound
$\epsilon_{2}$, and make the overall error $E_{a}\leq\epsilon_{2}$. Once
$E_{a}$ exceeds $\epsilon_{2}$, the Bayesian approach is considered unable to
work properly. Without loss of generality, we set $\epsilon_{2}=0.1$ in our
numerical experiments. Moreover, we find that $E_{a}$ exceeds $\epsilon_{2}$
only when certain parameter becomes too large, thus we define the point at
which $E_{a}$ reaches $\epsilon_{2}$ as the critical point. When certain
parameter exceeds its critical point, the Bayesian approach is regarded as
inapplicable. In the following we numerically find the critical value of a
specific parameter $a$, below or above which the Bayesian approach works well.
We keep the other parameter constant and make the parameter $a$ change within
a certain range. Based on the ensemble evolution of the system state, we
calculate the evolution of the error function $E_{a}$ with respect to the
parameter $a$. In this way, we can obtain the error curves of these
parameters, which describe the relationships between the values of the
parameters and the performance of the Bayesian approach.

In Fig. 3, we take $\gamma_{1}$ for example. The error curves of $\gamma_{1}$
for a total time of 10$\mu$s and with $5000$ single evolutions are illustrated
in Fig. 3(a). The blue solid curve, the red dashed curve, the yellow pecked
curve and the purple dotted curve represent $E_{\gamma_{1}}$,
$Ez_{\gamma_{1}}$, $Ey_{\gamma_{1}}$ and $Ex_{\gamma_{1}}$, respectively. The
red filled dot is the critical point of $\gamma_{1}$ at which $E_{\gamma_{1}}$
reaches $\epsilon_{2}$. It is clear that the increasing of $\gamma_{1}$ mainly
results in increasing of $Ez_{\gamma_{1}}$ and the critical point of
$\gamma_{1}$ is around $5.841$ kHz. For clarity, the time evolution of the
estimated $\langle\hat{\sigma}_{z}\rangle$ (solid blue curve) and the
evolution of the original $\langle\sigma_{z}\rangle$ (dash red curve) at the
critical point are plotted in Fig. 3(b) from which one can see that the
Bayesian approach doesn’t work well under such a parameter setting. This
result demonstrates the effectiveness of the proposed error definition.

The other parameters in the circuit-QED system we considered are $\chi$,
$\Delta_{r}$, $\varepsilon_{d}$, $\eta$, $\gamma_{\varphi}$, and $\kappa$. We
find that the values of all parameters except $\kappa$ have the ability to
result in rejecting of the Bayesian approach. The phenomenon of the errors
evolution with respect to $\kappa$ is plotted in Fig. 4(f). This is reasonable
because $\kappa$ represents the damping rate of the resonator and only has an
effect on the photocurrent Eq. (4). The rest error curves of the parameters
considered are also plotted in Fig. 4, from which the critical points can be
easily found to be around 19.57MHz, 0.5949MHz, 70.11MHz, 0.4817 and 0.2036MHz,
corresponding to $\chi$, $\Delta_{r}$, $\varepsilon_{d}$, $\eta$ and
$\gamma_{\varphi}$, respectively. All of these five parameters have
significant effects on the evolution of the error, where $\chi$ and
$\Delta_{r}$ mainly affect the estimation of the diagonal elements of the
density matrix. On the contrary, $\varepsilon_{d}$, $\eta$ and
$\gamma_{\varphi}$ mainly affect the estimation of the non-diagonal elements.
For convenience of reference, we conclude the suitable conditions that allow
the overall error $E_{a}\leq 0.1$ in Tab. 1. The definition of the error
$E_{a}$ as well as the values of the weighting factors can be adjusted
according to the specific precision requirement. In Tab. 2 we display the
suitable conditions of various systems parameters that allows the overall
error $E_{a}\leq 0.05$. There is no doubt that the selection becomes much
narrower. As a preliminary work, this paper only considers the situation under
the constraint of a single parameter, while the research about the
collaborative effects of the multi-parameters might be of interest for further
studies.

##  V Conclusion

In conclusion, this paper theoretically studies the Bayesian approach which
completes the task of real-time quantum state estimation in circuit QED. The
main content of our work is to analyse the suitable conditions under which the
Bayesian approach can real-time accurately estimate the quantum state. In
detail we explore the correlation between some basic and physically meaningful
parameters of the circuit QED and the performance of the Bayesian approach.
Our work helps to gain a deeper understanding of the quantum Bayesian approach
and promotes the application of the quantum Bayesian approach in quantum
control. Moreover, through the analysis of the correlation between the system
parameters and the performance of the method, the roles that these parameters
play in the system can be understood more vividly.

###### Acknowledgements.

This work is supported by National Natural Science Foundation of China under
Grant 11404113, and the Guangzhou Key Laboratory of Brain Computer Interaction
and Applications under Grant 201509010006.

## References

  * (1) A. Bisio, G. Chiribella, G. M. D’Ariano, S. Facchini, and P. Perinotti, Optimal Quantum Tomography, IEEE Journal of Selected Topics in Quantum Electronics 15, 1646-1660 (2009).
  * (2) A. C. Doherty, S. Habib, K. Jacobs, H. Mabuchi, and S. M. Tan, Quantum feedback control and classical control theory, Phys. Rev. A 62, 012105 (2000).
  * (3) T. D. Ladd, F. Jelezko, R. Laflamme, Y. Nakamura, C. Monroe, and J. L. O’Brien, Quantum Computers, Nature 464, 45-53 (2010).
  * (4) P. Six, P. Campagneibarcq, I. Dotsenko, A. Sarlette, B. Huard, and P. Rouchon1, Quantum state tomography with non-instantaneous measurements, imperfections and decoherence, Phys. Rev. A, 93, 012109 (2016).
  * (5) B. Qi, Z. Hou, L. Li, D. Dong, G. Xiang, and G. Guo, Quantum State Tomography via Linear Regression Estimation, Sci. Rep. 3, 3496 (2013).
  * (6) D. Gross, Y. K. Liu, S. T. Flammia, S. Becker, and J. Eisert, Quantum state tomography via compressed sensing, Phys. Rev. Lett. 105, 150401 (2009).
  * (7) S. T. Flammia, D. Gross, Y. K. Liu, and J. Eisert, Quantum Tomography via Compressed Sensing: Error Bounds, Sample Complexity, and Efficient Estimators, New J. Phys. 14, 095022 (2012).
  * (8) A. Kalev, R. L. Kosut, and I. H. Deutsch, Quantum tomography protocols with positivity are compressed sensing protocols, NPJ Quan. Info. 1, 15018 (2015).
  * (9) V. Braginsky, F. Y. Khalili, and E. Thorne,_Quantum measurement_(Cambridge University Press, Cambridge, 1992).
  * (10) H. M Wiseman and G. J. Milburn, _Quantum Measurement and Control_ (Cambridge University Press, Cambridge, 2009).
  * (11) S. Tanaka and N. Yamamoto, Robust adaptive measurement scheme for qubit-state preparation, Phys. Rev. A 86, 062331 (2012).
  * (12) A. N. Korotkov, Continous quantum measurement of a double dot, Phys. Rev. B 60, 5737 (1999).
  * (13) A. N. Korotkov, Quantum Bayesian approach to circuit QED measurement with moderate bandwidth, Phys. Rev. A, 94, 042326 (2016).
  * (14) S. J. Weber, A. Chantasri, J. Dressel, A. N. Jordan, K. W. Murch, and I. Siddiqi, Mapping the optimal route between two quantum states, nature 511 ,570-573 (2014).
  * (15) R. Vijay, C. Macklin, D. H. Slichter, S. J. Weber, K. W. Murch, R. Naik1, A. N. Korotkov, and I. Siddiqi, Stabilizing Rabi oscillations in a superconducting qubit using quantum feedback, Nature 490, 77-80 (2012).
  * (16) G. G. Gillett, R. B. Dalton, B. P. Lanyon, M. P. Almeida, M. Barbieri, G. J. Pryde, J. L. O’Brien, K. J. Resch, S. D. Bartlett, and A. G. White, Experimental feedback control of quantum systems using weak measurements, Phys. Rev. Lett. 104, 080503 (2010).
  * (17) P. Wang, L. Qin, and X. Q. Li, Quantum Bayesian rule for weak measurements of qubits in superconducting circuit QED, New J. Phys. 16, 123047 (2014).
  * (18) W. Feng, P. F. Liang, L. P. Qin and X. Q. Li, Exact quantum Bayesian rule for qubit measurements in circuit QED, Sci. Rep. 6, 20492 (2016).
  * (19) L. Qin, L. Xu, W. Feng, and X. Q. Li, Qubit state tomography in a superconducting circuit via weak measurements, New J. Phys. 19, 033036 (2017).
  * (20) J. Zhang, Y. Liu, R. B. Wu, K. Jacobs, and F. Nori, Quantum feedback: theory, experiments, and applications, Phys. Rep. 679,1-60 (2017).
  * (21) A. Blais, R. S. Huang, A. Wallraff, S. M. Girvin, and R. J. Schoelkopf, Cavity quantum electrodynamics for superconducting electrical circuits: An architecture for quantum computation, Phys. Rev. A 69, 062320 (2004).
  * (22) Z. L. Xiang, S. Ashhab, J. Q. You, and F. Nori, Hybrid quantum circuits: Superconducting circuits interacting with other quantum systems, Rev. Mod. Phys. 85, 623-654 (2013).
  * (23) D. H. Slichter, R. Vijay, S. J. Weber, S. Boutin, M. Boissonneault, J. M. Gambetta, A. Blais, and I. Siddiqi, Measurement-induced qubit state mixing in circuit QED from up-converted dephasing noise, Phys. Rev. Lett. 109, 153601 (2012).
  * (24) M. Hatridge, S. Shankar, M. Mirrahimi, F. Schackert, K. Geerlings, T. Brecht, K. M. Sliwa, B. Abdo, L. Frunzio, S. M. Girvin, R. J. Schoelkopf, and M. H. Devoret, Quantum back-action of an individual variable-strength measurement, Science 339, 178-181 (2013).
  * (25) A. Wallraff, D. I. Schuster, A. Blais, and L. Frunzio, Strong coupling of a single photon to a superconducting qubit usingcircuit quantum electrodynamics, Nature 431, 162-167 (2004).
  * (26) I. Chiorescu, P. Bertet, K. Semba, Y. Nakamura, C. J. P. M. Harmans, and J. E. Mooij, Coherent dynamics of a flux qubit coupled to a harmonic oscillator, Nature 431, 159-162 (2004).
  * (27) J. Gambetta, A. Blais, D. I. Schuster, A. Wallraff, L. Frunzio, J. Majer, M. H. Devoret, S. M. Girvin, and R. J. Schoelkopf, Qubit-photon interactions in a cavity: Measurement-induced dephasing and number splitting, Phys. Rev. A 74, 042318 (2006).
  * (28) J. Gambetta, A. Blais, M. Boissonneault, A. A. Houck, D. I. Schuster, and S. M. Girvin, Quantum trajectory approach to circuit QED: Quantum jumps and the Zeno effect, Phys. Rev. A 77, 012112 (2008).