Multi-parameter estimation with multi-mode Ramsey interferometry

One of the central objectives of quantum metrology concerns improving measurement precision with finite sized ensembles [1–4]. Most previous investigations have focused on single parameter estimation, of which the standard quantum limit (SQL) or the classical limit, $1/\sqrt{N}$, represents the minimal phase uncertainty achievable in an interferometric measurement using an ensemble of N uncorrelated particles [5]. Recently, the problem of estimating multiple parameters has attracted much interests [6–41], where the focus shifts to finding efficient strategies for estimating parameters corresponding to multiple commuting or non-commuting generators as precisely as possible. Potential applications of such studies include quantum imaging [8, 15, 28], sensor networks [37, 40], measurements of multidimensional fields [18], and joint measurements of multiple quadratures [10, 32–34], etc.

The main tasks in multi-parameter estimation are to generate an input quantum state capable of realizing the optimal precision limited by the laws of quantum mechanics, and to find a corresponding measurement scheme that achieves this precision. In the language of estimation theory, the former obtains a quantum state with the lowest quantum Cramér–Rao bound (QCRB) for a set of parameters to be estimated, while the latter provides measurement results of which the Cramér–Rao bound (CRB) equals the QCRB. For estimation of a single parameter, the latter can always be fulfilled using interferometry in which the second beam splitter acts as the inverse transformation (time reversed operation) to the first [42], such as in an atomic Ramsey interferometer and an optical Mach–Zehnder interferometer. However, the same does not apply in general for multi-parameter estimation.

This work considers the more specific case of multi-mode interferometry for estimating a set of parameters corresponding to mutually commuting generators using unentangled particles (as illustrated in figure 1). A probe state is generated by splitting a pure single-mode state $\left|i\right\rangle$ into multiple modes using a multi-mode beam splitter represented by a unitary transformation U₁. The prepared (D + 1)-mode probe state then undergoes phase accumulation, and is subsequently transformed by U₂ at the second beam splitter. The measured particle number distributions at the (D + 1)-outputs are used to estimate the D parameters in the end. Unlike the case of single parameter estimation, setting ${U}_{2}={U}_{1}^{\dagger }={U}_{1}^{-1}$ does not guarantee CRB will be equal to the QCRB in general when D > 1. Instead for a given U₁, U₂ has to be optimized numerically to reach QCRB. This becomes a cumbersome and tedious job particularly when the number of parameters to be estimated is large.

Zoom In Zoom Out Reset image size

As a main result to be reported in this paper, we prove that the Ramsey interferometric scheme can be straightforwardly generalized to estimation of multiple parameters (associated with commuting generators) using multi-mode pure states, if (but not iff) U is made orthogonal, i.e. when U is unitary and has only real matrix elements. We also illustrate how such orthogonal U can be constructed experimentally in a simple and scalable way. The influences of particle number fluctuation and detection noise will also be discussed.

So far, most measurement schemes which saturate the QCRB for multi-parameter estimation, if they exist, are found on a case by case basis. Important progresses have been made in this direction recently [8, 27, 32, 36, 39]. However many of the proposed measurement schemes are either not directly implementable or experimentally prohibitive, particularly when they involve measurements on entangled particles [6–41]. Therefore, generalization of Ramsey interferometry to multi-parameter estimation represents an interesting and timely advance.

This article is organized as follows: section 2 defines the problem we consider and gives the QCRB of a multi-mode probe state. Section 3 proves that for an orthogonal U₁, the CRB from setting ${U}_{2}={U}_{1}^{-1}$ is equal to the QCRB. Section 4 illustrates how to determine the optimal probe state that gives the lowest QCRB. In section 5, we show how an orthogonal U can always be constructed experimentally in a simple and scalable manner in an optical or atomic system. Finally, we consider the influence of particle number fluctuation and detection noise on the multi-mode Ramsey interferometer in sections 6 and 7, respectively. The article ends with appendices A and B containing further calculation details.

In this section,we define the problem we consider and give the QCRB for a given probe state. As shown in figure 1,the parameters we consider are encoded into quantum states with D + 1 modes,which can be implemented with photons split into multiple paths,or atoms with large spins. For unentangled particles,the interferometry can be discussed in terms of self-interference of individual particles [43]. Therefore,we consider an arbitrary single particle initial state $\left|i\right\rangle$,and a probe state $\left|{\psi }_{{\rm{p}}}\right\rangle ={U}_{1}\left|i\right\rangle ={\sum }_{k=0}^{D}{\alpha }_{k}\left|k\right\rangle$ after transformation U₁,with ${\alpha }_{k}$ being the probability amplitude in mode k. We assume that the probe state is pure for the time being and the interferometry is noiseless. The phase accumulation evolves the probe state into $\left|{\psi }_{{\boldsymbol{\phi }}}\right\rangle ={\sum }_{k=0}^{D}{\alpha }_{k}{{\rm{e}}}^{{\rm{i}}{\phi }_{k}}\left|k\right\rangle$. Interference from the first order coherence allows D (out of the D + 1) phases to be measured in the absence of an external reference. This is often carried out by choosing an arbitrary mode,say $\left|0\right\rangle$,as the reference,and measuring the relative phase shifts ${\theta }_{k}\equiv {\phi }_{k}-{\phi }_{0}$ (k = 1,2,..,D). However,if the parameters of interest ${\boldsymbol{\Theta }}\equiv \left\{{{\rm{\Theta }}}_{1},{{\rm{\Theta }}}_{2},\,\cdots ,{{\rm{\Theta }}}_{D}\right\}$ are not the same as ${\boldsymbol{\theta }}\equiv \{{\theta }_{1},{\theta }_{2},\,\ldots ,{\theta }_{D}\}$, optimized precision of ${\boldsymbol{\theta }}$ does not necessarily give the best precision for ${\boldsymbol{\Theta }}$ in general. We assume in the following that each parameter of interest, Θ_k, is a linear combination of $\left\{{\phi }_{0},{\phi }_{1},\cdots ,{\phi }_{D}\right\}$ in general, and the goal turns to finding a probe state that minimizes the total phase variance ${\left({\rm{\Delta }}{\boldsymbol{\Theta }}\right)}^{2}={\sum }_{k=1}^{D}{\left({\rm{\Delta }}{{\rm{\Theta }}}_{k}\right)}^{2}$.

With ${\boldsymbol{\Theta }}$ defined, the phases can in turn be written as ${\phi }_{k}={f}_{k}({\boldsymbol{\Theta }})$ (see appendix B for more details), and the probe state after phase accumulation becomes $\left|{\psi }_{{\boldsymbol{\phi }}}\right\rangle ={\sum }_{k=0}^{D}{\alpha }_{k}{{\rm{e}}}^{{{\rm{i}}{f}}_{k}\left({\boldsymbol{\Theta }}\right)}\left|k\right\rangle$. According to multi-parameter quantum estimation theory [44, 45], the lower bound of ${\left({\rm{\Delta }}{\boldsymbol{\Theta }}\right)}^{2}$ with an unbiased estimator is determined by the trace of the inverse of quantum Fisher information matrix (QFIM) ${{ \mathcal F }}^{Q}$:

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{\boldsymbol{\Theta }}\right)}^{2}\geqslant \mathrm{Tr}\left[{\left({N}_{\hat{M}}{{ \mathcal F }}^{Q}\right)}^{-1}\right],\end{eqnarray} \tag{ 1 }$$

where ${N}_{\hat{M}}$ is the number of experiments repeated (set to 1 hereafter for simplicity). Note that, the choice of a figure of merit for precision as the trace of the inverse of ${{ \mathcal F }}^{Q}$ in equation (1) is fully general, since the weight of each parameter ϕ_k can be adjusted by changing the coefficients of the linear combinations in ${\phi }_{k}={f}_{k}({\boldsymbol{\Theta }})$.

For a pure state $\left|{\psi }_{{\boldsymbol{\phi }}}\right\rangle$, the matrix elements of ${{ \mathcal F }}^{Q}$ are explicitly given by [44, 45]

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{Q}=4\mathrm{Re}\left[\left\langle {\partial }_{{{\rm{\Theta }}}_{l}}{\psi }_{{\boldsymbol{\phi }}}| {\partial }_{{{\rm{\Theta }}}_{n}}{\psi }_{{\boldsymbol{\phi }}}\right\rangle -\left\langle {\partial }_{{{\rm{\Theta }}}_{l}}{\psi }_{{\boldsymbol{\phi }}}| {\psi }_{{\boldsymbol{\phi }}}\right\rangle \left\langle {\psi }_{{\boldsymbol{\phi }}}| {\partial }_{{{\rm{\Theta }}}_{n}}{\psi }_{{\boldsymbol{\phi }}}\right\rangle \right],\end{eqnarray} \tag{ 2 }$$

where l, n = 1, 2, ⋯, D. The matrix elements of the D × D single-particle QFIM for $| {\psi }_{{\boldsymbol{\phi }}}\rangle$ are thus given by

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{Q}=4\left[\sum _{k=0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}-\sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial {f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}{\left|{\alpha }_{k^{\prime} }\right|}^{2}\right].\end{eqnarray} \tag{ 3 }$$

The QFIM is convex and additive [38]. For an uncorrelated but identically prepared N-particle product state, ${\left|{\psi }_{{\boldsymbol{\phi }}}\right\rangle }^{\otimes N}$, it is nothing but just the sum (N-times) of the single-particle QFIM. According to equation (1), the probe state that gives the best QCRB for estimating ${\boldsymbol{\Theta }}$ can be obtained by minimizing the trace of ${\left({{ \mathcal F }}^{Q}\right)}^{-1}$, via varying $| {\alpha }_{k}{| }^{2}$ under the normalization condition ${\sum }_{k=0}^{D}{\left|{\alpha }_{k}\right|}^{2}=1$. It is clear from equation (3) that the QCRB of a probe state depends only on the distribution of the particles $| {\alpha }_{k}{| }^{2}$ but not on the phase of α_k.

In the above section, we discuss how to calculate the ultimate sensitivity bound for a given beam splitting transformation U₁ and for a given set of parameters corresponding to mutually commuting generators. To saturate this bound, optimization of the second beam combining transformation (U₂) is required on a case by case basis, which is a tedious task for a large number of parameters. In this section, we prove that if (but not iff) the multi-mode unitary transformation U₁ = U is real (orthogonal), then ${U}_{2}={U}^{\dagger }={U}^{-1}$ followed by particle number detection afterwards gives the best precision allowed by the QCRB when ${\boldsymbol{\Theta }}\sim 0$.

For Ramsey interferometry, the state after the full interferometric protocol (before particle number detection) is represented by $\left|{\psi }_{\mathrm{out}}\right\rangle ={U}^{\dagger }{\prod }_{k=0}^{D}{{\rm{e}}}^{{\rm{i}}\left|k\right\rangle \left\langle k\right|{f}_{k}\left({\boldsymbol{\Theta }}\right)}U\left|i\right\rangle$. The CRB, which sets the minimal ${\left({\rm{\Delta }}{\boldsymbol{\Theta }}\right)}^{2}$ given a measurement scheme, can be calculated for any U using the classical Fisher information matrix (CFIM) [46]

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{C}\left({\boldsymbol{\Theta }}\right)=\sum _{m=0}^{D}\displaystyle \frac{1}{p\left(m| {\boldsymbol{\Theta }}\right)}\displaystyle \frac{\partial p\left(m| {\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial p\left(m| {\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}},\end{eqnarray} \tag{ 4 }$$

where $p(m| {\boldsymbol{\Theta }})={| \langle m| {U}^{\dagger }{\prod }_{k=0}^{D}{{\rm{e}}}^{{\rm{i}}| k\rangle \langle k| {f}_{k}({\boldsymbol{\Theta }})}U| i\rangle | }^{2}$ denotes the probability of finding a particle in mode $| m\rangle$ for a given ${\boldsymbol{\Theta }}$.

To show that a Ramsey interferometric scheme can be used to estimate multiple parameters to the SQL precision, we need to prove the CFIM given by equation (4) equals to QFIM given by equation (1) (since they correspond to the CRB and the QCRB, respectively). For small ${\boldsymbol{\Theta }}$, omitting the third order corrections, a Taylor series expansion around ${\boldsymbol{\Theta }}\sim 0$ gives

$$\begin{eqnarray}p\left(m| {\boldsymbol{\Theta }}\right)\simeq \left\{\begin{array}{ll}{\sum }_{k,k^{\prime} =0}^{D}{f}_{k}\left({\boldsymbol{\Theta }}\right){f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)\left\langle m\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|m\right\rangle & m\ne \ i,\\ 1+{\sum }_{k,k^{\prime} =0}^{D}{f}_{k}\left({\boldsymbol{\Theta }}\right){f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right){\left|\left\langle i\right|{U}^{\dagger }\left|k\right\rangle \right|}^{2}{\left|\left\langle k^{\prime} \right|U\left|i\right\rangle \right|}^{2}-{\sum }_{k=0}^{D}{f}_{k}{\left({\boldsymbol{\Theta }}\right)}^{2}{\left|\left\langle i\right|{U}^{\dagger }\left|k\right\rangle \right|}^{2} & m=i.\end{array}\right.\end{eqnarray} \tag{ 5 }$$

The derivatives of equation (5) with respect to any Θ_l can also be calculated directly (see appendix A). Substituting equation (5) and its derivatives into the definition of the CFIM (equation (4)) gives a complicated equation that looks vastly different from the QFIM of equation (2). Indeed, for an arbitrary U, the CFIM for an Ramsey interferometric scheme is not equal to the QFIM in most cases. We find, however, that when $\left\langle k\right|U\left|m\right\rangle =\left\langle m\right|{U}^{\dagger }\left|k\right\rangle$ for all k, m (meaning that all elements of U are real and U is orthogonal), the equation for the CFIM is simplified greatly and becomes (see appendix A)

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{C}\left({\boldsymbol{\Theta }}\sim 0\right)\simeq 4\sum _{m\ne i}\sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial {f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}}\left\langle i\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|m\right\rangle \left\langle m\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|i\right\rangle .\end{eqnarray} \tag{ 6 }$$

To further simplify the formula, we make use of the completeness of the basis $\left|i\right. \rangle \left. \langle i\right|+{\sum }_{m\ne i}\left|m\right. \rangle \left. \langle m\right|={\bf{1}}$. This gives

$$\begin{eqnarray}\begin{array}{c}\begin{array}{rcl}{{ \mathcal F }}_{l,n}^{C}\left({\boldsymbol{\Theta }}\sim 0\right) & \simeq & 4\displaystyle \sum _{m\ne i}\displaystyle \sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\mathop{\underbrace{\left|m\right. \rangle \left. \langle m\right|}}\limits_{m\ne i}{U}^{\dagger }\left|{k}^{{\rm{{\prime} }}}\right. \rangle \left. \langle {k}^{{\rm{{\prime} }}}\right|U\left|i\right. \rangle \\ & = & 4\left[\displaystyle \sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\left({\bf{1}}-\left|i\right. \rangle \left. \langle i\right|\right){U}^{\dagger }\left|{k}^{{\rm{{\prime} }}}\right. \rangle \left. \langle {k}^{{\rm{{\prime} }}}\right|U\left|i\right. \rangle \right]\\ & = & 4\left[\displaystyle \sum _{k=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\left|i\right. \rangle \right.\\ & & \left.-\displaystyle \sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\left|i\right. \rangle \left. \langle i\right|{U}^{\dagger }\left|{k}^{{\rm{{\prime} }}}\right. \rangle \left. \langle {k}^{{\rm{{\prime} }}}\right|U\left|i\right. \rangle \right]\\ & = & 4\left[\displaystyle \sum _{k=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}-\sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}{\left|{\alpha }_{{k}^{{\rm{{\prime} }}}}\right|}^{2}\right]={{ \mathcal F }}_{l,n}^{Q}.\end{array}\end{array}\end{eqnarray} \tag{ 7 }$$

In [27], Pezzè et al found the necessary and sufficient conditions (iff) for projective measurements which saturate the QCRB of a probe state. In their language, our measurement can be described by a set of projectors $\left\{\left|{{\rm{\Upsilon }}}_{k}\right\rangle \left\langle {{\rm{\Upsilon }}}_{k}\right|\right\}$, where $\left|{{\rm{\Upsilon }}}_{k}\right\rangle =U\left|k\right\rangle$. In the limit ${\boldsymbol{\Theta }}\sim 0$ and given that all elements of U are real, the projectors $\left\{\left|{{\rm{\Upsilon }}}_{k}\right\rangle \left\langle {{\rm{\Upsilon }}}_{k}\right|\right\}$ indeed satisfy the required condition given by their equation (7) in [27].

In this section, we demonstrate how to determine the optimal probe state. As an illustration, we consider the most common choice of ${{\rm{\Theta }}}_{k}={\theta }_{k}\equiv {\phi }_{k}-{\phi }_{0}$. In this case, the generator of parameter θ_k is proportional to $\left|k\right\rangle \left\langle k\right|$. Computing equation (3) and taking the trace of ${\left({{ \mathcal F }}^{Q}\right)}^{-1}$ gives (after dividing by particle number N) (see appendix B for more details)

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{{\boldsymbol{\theta }}}_{}\right)}^{2}\geqslant \displaystyle \frac{1}{N}\left[\displaystyle \frac{D}{4{\left|{\alpha }_{0}\right|}^{2}}+\sum _{k=1}^{D}\displaystyle \frac{1}{4{\left|{\alpha }_{k}\right|}^{2}}\right].\end{eqnarray} \tag{ 8 }$$

Minimizing equation (8) under the condition ${\sum }_{k=0}^{D}{\left|{\alpha }_{k}\right|}^{2}=1$ gives the optimal probe state described by

$$\begin{eqnarray}&&{\left|{\alpha }_{0}\right|}^{2}=\sqrt{D}/\left(D+\sqrt{D}\right),\end{eqnarray} \tag{ 9a }$$

$$\begin{eqnarray}&&{\left|{\alpha }_{k}\right|}^{2}=1/\left(D+\sqrt{D}\right),\quad \left(\mathrm{for}\,k\ne 0\right),\end{eqnarray} \tag{ 9b }$$

and the QCRB of

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{{\boldsymbol{\theta }}}_{\mathrm{opt}}\right)}^{2}={\left(D+\sqrt{D}\right)}^{2}/4N.\end{eqnarray} \tag{ 10 }$$

This precision can be reached in the asymptotic regime of large ${N}_{\hat{M}}$.

For comparison, we consider an individual estimation scheme which divides the N particles into D equal partitions, and uses each partition for measuring one ${\theta }_{k}$ through two-mode interferometry between $| 0\rangle$ and $| k\rangle$. Since the SQL of each θ_k in this case is $1/\sqrt{N/D}$, the lowest bound for the phase variance becomes

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{{\boldsymbol{\theta }}}_{\mathrm{ind}}\right)}^{2}={D}^{2}/N.\end{eqnarray} \tag{ 11 }$$

For D = 1 as in single parameter estimation, both equations (10) and (11) reduce to $1/N$ as expected (i.e. the SQL). For larger D, the simultaneous estimation scheme (equation (10)) always outperforms the individual estimation scheme (equation (11)).

We note that the results of equations (9) and (10) resemble an earlier study [8], where Humphreys et al considered a multi-mode entangled NOON state $\left|{\psi }_{\mathrm{in}}\right\rangle ={\alpha }_{0}\left|N,0,\,\ldots ,0\right\rangle \ +{\alpha }_{1}\left|0,N,\ldots ,\,0\right\rangle \ +\ldots +\ {\alpha }_{D}\left|0,0,\ldots ,N\right\rangle$. They found an optimal probe defined also by equation (9) and a QCRB N times smaller than equation (10), in agreement with the typical ratio between the SQL and the Heisenberg limit (HL). Their results reduce to ours when a multi-mode NOON state for N = 1 is considered.

The value of the reference mode ${\left|{\alpha }_{0}\right|}^{2}$ in equation (9) is $\sqrt{D}$ times larger than the other modes. Because ϕ₀ is referenced to by all ${\theta }_{k}\equiv {\phi }_{k}-{\phi }_{0}$, the measurement variance of ϕ₀ therefore contributes D times more to ${\left({\rm{\Delta }}{\boldsymbol{\theta }}\right)}^{2}$ than any other uncorrelated phases {ϕ_k}. Such a bias results from the choice of parameters. Consequently, equations (9) and (10) cannot be always optimal if the parameters of interest are different. For instance, should we consider a different set of parameters of interest, say, ${\varphi }_{k}\equiv {\phi }_{k}-{\phi }_{k-1}\,(k=1,2,\,..,D)$, i.e. the relative phase between the neighboring modes, repeating the above same procedures gives a minimal variance of ${\left({\rm{\Delta }}{{\boldsymbol{\varphi }}}_{\mathrm{opt}}\right)}^{2}=\displaystyle \frac{1}{4N}{\left[\sqrt{2}\left(D-1\right)+2\right]}^{2}$. If one measures {θ_k} instead of {φ_k} using the probe state given by equation (9) and then derives {φ_k} from {θ_k}, the resulting phase variance would be bounded by ${\left({\rm{\Delta }}{\boldsymbol{\varphi }}\right)}^{2}=\displaystyle \frac{1}{4N}\left[{(1+\sqrt{D})}^{2}+2(D-1)(\sqrt{D}+D)\right]$(for details see appendix B), a result always larger than ${\left({\rm{\Delta }}{{\boldsymbol{\varphi }}}_{\mathrm{opt}}\right)}^{2}$ for D > 1.

We now illustrate how a multi-mode Ramsey interferometer can be realized experimentally in a simple and scalable way. Here, the task reduces to designing an orthogonal U that generates the optimal probe state. We would illustrate our scheme first for an optical interferometer and then an atomic interferometer.

5.1. Multi-mode optical interferometer

We consider a design that employs a series of 2 × 2 non-polarizing beam splitters (BS_k, for k = 1, 2, 3, ⋯) for splitting particles into the optimal distributions $| {\alpha }_{k}{| }^{2}$ as shown in figure 2. In this case, each U^(k) which represents the transformation due to beam splitter BS_k acts only on two of the adjacent modes, leaving other modes untouched. The overall transformation $U={U}^{(D)}{U}^{(D-1)}\cdots {U}^{(1)}$ must be unitary since each lossless physical splitter U^(k) is unitary.

Zoom In Zoom Out Reset image size

To ensure that the resulting U and ${U}^{\dagger }$ constructed from these beam splitters are real (orthogonal), the most straightforward way is to make sure that each of the beam splitters behaves as a real 2 × 2 transformation. This criterion, which requires zero (or multiple of 2π) phase shifts for both the transmitted and reflected beams with respect to both input beams, is not automatically satisfied for any beam splitters. Fortunately, it is always possible to fulfil this criterion by adding respective phase compensating waveplate to each port of a beam splitter. After compensation, the matrix elements of a real U^(k) become ${U}_{k,k}^{(k)}={U}_{k+1,k+1}^{(k)}=\cos ({\eta }_{k}),{U}_{k,k+1}^{(k)}=-{U}_{k+1,k}^{(k)}\,=\sin ({\eta }_{k})$, and ${U}_{i,j\rlap{/}{\in }\{k,k+1\}}^{(k)}={\delta }_{i,j}$ (the Kronecker delta function), where ${\cos }^{2}{\eta }_{k}$ (${\sin }^{2}{\eta }_{k}$) represents the reflectance (transmittance) of ${\mathrm{BS}}_{k}$. Given an optimal distribution $| {\alpha }_{k}{| }^{2}$, the reflectance of BS_k should be chosen as

$$\begin{eqnarray}&&| \cos {\eta }_{1}{| }^{2}={\left|{\alpha }_{0}\right|}^{2},\end{eqnarray} \tag{ 12a }$$

$$\begin{eqnarray}&&| \cos {\eta }_{k}{| }^{2}=\displaystyle \frac{{\left|{\alpha }_{k-1}\right|}^{2}}{1-{\left|{\alpha }_{0}\right|}^{2}-{\left|{\alpha }_{1}\right|}^{2}\cdots -\ {\left|{\alpha }_{k-2}\right|}^{2}}\quad \left(k\geqslant 2\right).\end{eqnarray} \tag{ 12b }$$

In addition, extra phase compensators are needed in every arm of the interferometer to null out the difference in optical path lengths and to tune every phase shift ${\phi }_{k}$ to the region where ${\boldsymbol{\Theta }}$ can be measured most sensitively. If there is no detection noise, this region is ${\boldsymbol{\Theta }}$ close to zero, otherwise, it is shifted away from ${\boldsymbol{\Theta }}\sim 0$ (see section 7). We emphasize that compensating for ϕ_k to give the optimal sensitivity does not represent a flaw, in fact, as such tuning is needed in practically all real interferometric measurements near the SQL precision.

5.2. Multi-mode atomic interferometer

In interferometry of atoms with hyperfine spin $F,D=2F$ different parameters can be estimated. Analogous to the optical scheme, an arbitrary spin distribution can be constructed using a sequence of Rabi rotations between two adjacent Zeeman sublevels. Such rotations can be realized, for instance, using a two-photon Raman transition through an intermediate state as illustrated in figure 3. As long as the intermediate hyperfine levels have a different Landé g-factor from those involved in interferometry, one could perform Rabi rotations between any two adjacent sublevels by selectively detuned to a suitable intermediate states. To make sure that the individual transformation is orthogonal, every rotation should be performed along the ${\sigma }_{y}=\left(\begin{array}{cc}0 & -{\rm{i}}\\ {\rm{i}} & 0\end{array}\right)$ direction, such that ${U}^{(k)}=\exp (-{\rm{i}}{\sigma }_{y}{\beta }_{k})=\left[\begin{array}{cc}\cos ({\beta }_{k}) & -\sin ({\beta }_{k})\\ \sin ({\beta }_{k}) & \cos ({\beta }_{k})\end{array}\right]$, within the two-level subspace. However, as each of the Zeeman sublevel exhibits different shift inside a magnetic field and thus different phase accumulation rate, one would need to keep track of the phases of every levels and to account for them when performing individual Rabi rotations. While this is possible with current technologies in cold atom experiments, the process is perhaps too cumbersome to be practical, especially when atomic spin is large.

Zoom In Zoom Out Reset image size

For the aforementioned reasons, we restrict the transformation in the following to a single-pulse multi-mode Rabi rotation over an angle χ along the F_y direction (since the corresponding matrix $U=\exp (-{{\rm{i}}{F}}_{y}\chi )$ is always orthogonal for any atomic spin F), and study the performance of the Ramsey interferometric protocol for measuring ${\boldsymbol{\theta }}$. Experimentally, such a F_y rotation can be realized using a radio-frequency resonant with adjacent Zeeman sublevels, when the quadratic Zeeman shift is negligible. It transforms the initial state $\left|F,{m}_{i}\right\rangle$ into $\left|{\psi }_{{\rm{p}}}\right\rangle ={\sum }_{m=-F}^{F}{d}_{{m}_{i},m}^{F}\left(\chi \right)\left|F,m\right\rangle$ with the Wigner's (small) d-matrix. According to equation (8), the QCRB of this state is given by

$$\begin{eqnarray}&&\displaystyle \frac{1}{4N}\left[\displaystyle \frac{2F-1}{| {d}_{{m}_{0},{m}_{i}}^{F}\left(\chi \right){| }^{2}}+\sum _{m=-F}^{F}\displaystyle \frac{1}{| {d}_{{m}_{0},m}^{F}\left(\chi \right){| }^{2}}\right],\end{eqnarray} \tag{ 13 }$$

when $| F,{m}_{0}\rangle$ is chosen as the reference mode. Figures 4 (a)–(c) present the values of equation (13) for F = 1, 3, and 5, respectively (for m₀ = 0). This one-step-rotation scheme (OSRS), which employs the limited family of a single SU(2) transformation, is found to always outperform the individual measurement scheme (equation (11), grey dashed horizontal line) using a suitable initial state $\left|F,{m}_{i}\right\rangle$ and a rotation angle χ, at least up to F = 5 (figure 4(d)). The same conclusion is reached for parameters {φ_k}.

Zoom In Zoom Out Reset image size

Such multi-parameter estimation scheme can be useful when atoms are subjected to different sources of phase shifts simultaneously, as for example, with spin-1 ⁸⁷Rb atoms dressed by near-resonant microwaves while under a static magnetic field [47, 48], or spin-9/2 ⁸⁷Sr atoms placed in an optical lattice with polarization dependent light shifts, and collisions with background or non-condensed atoms.

In this section, we discuss the influence of particle number fluctuation of the probe state. Since quantum states with a definite large particle number are often difficult to prepare, we consider the situation when the particle number of the probe state fluctuates. Due to the superselection rule, such input state represents nothing but an incoherent superposition of different Fock state ${\rho }_{\mathrm{in}}={\oplus }_{N=0}^{+\infty }{Q}_{N}{\rho }^{\left(N\right)}$ in the absence of number coherences in the probe state and/or in the measurement strategy [49–51], where ${\rho }^{\left(N\right)}$ is the density matrix of the N-particle state and ${Q}_{N}$ the probability of having N particles. For coherent light of photons or an atomic Bose–Einstein condensate, the particle number obeys Poisson distribution with the probability ${Q}_{N}={{\rm{e}}}^{-\bar{N}}{\bar{N}}^{N}/N!$, where $\bar{N}$ denotes the mean particle number. Since QFIM is additive under a direct sum of density matrix ${\rho }^{\left(N\right)}$ in orthogonal subspaces [2], ${{ \mathcal F }}^{Q}({\rho }_{\mathrm{in}})={{ \mathcal F }}^{Q}\left[{\oplus }_{N}{Q}_{N}{\rho }^{\left(N\right)}\right]={\sum }_{N}{Q}_{N}{{ \mathcal F }}^{Q}\left[{\rho }^{\left(N\right)}\right]$. For unentangled N-particle states of the form ${\rho }^{(N)}={\rho }_{\mathrm{single}}^{\otimes N}={\left[| {\psi }_{{\rm{p}}}\rangle \langle {\psi }_{{\rm{p}}}| \right]}^{\otimes N},{{ \mathcal F }}^{Q}\left[{\rho }^{\left(N\right)}\right]=N{{ \mathcal F }}_{\mathrm{single}}^{Q}$, with ${{ \mathcal F }}_{\mathrm{single}}^{Q}$ being the QFIM of the single particle probe state ${\rho }_{\mathrm{single}}$. One has therefore ${{ \mathcal F }}^{Q}({\rho }_{\mathrm{in}})={\sum }_{N}{Q}_{N}N{{ \mathcal F }}_{\mathrm{single}}^{Q}={{ \mathcal F }}_{\mathrm{single}}^{Q}\bar{N}$. Similarly, it can be readily shown that ${{ \mathcal F }}^{C}={\sum }_{N}{Q}_{N}{{ \mathcal F }}^{C}[{\rho }_{\mathrm{single}}^{\otimes N}]={{ \mathcal F }}_{\mathrm{single}}^{C}\bar{N}$ for input state ${\rho }_{\mathrm{in}}$ [52] if particle numbers in all output ports are measured without detection noise. Since the CFIM of a single particle probe, ${{ \mathcal F }}_{\mathrm{single}}^{C}\approx {{ \mathcal F }}_{\mathrm{single}}^{Q}$ in our scheme, ${{ \mathcal F }}^{C}$ also equals approximately to ${{ \mathcal F }}^{Q}$ for probe state with fluctuating particle number and all our conclusions outlined above remain intact.

The conclusions in section 3 and appendix A are reached assuming noiseless particle number detections. When detection noise is taken into consideration, the optimal sensitivity typically shifts away from ${\boldsymbol{\Theta }}\sim 0$. For example, for single parameter estimation using Ramsey interferometry in an atomic clock, the measurement is usually performed near θ ∼ π/2, a region least sensitive to detection noise.

Here, we study numerically the effects of detection noise to the multi-parameter Ramsey interferometry using the example of two parameter estimation. We consider estimation of θ₁ and θ₂ using the optimal probe state given by equation (9). Starting from the initial state $| 0\rangle ={(0,1,0)}^{\dagger }$, we choose an orthogonal U given by the SU(2) rotation of a spin-1 system along F_y-direction

$$\begin{eqnarray}U=\exp \left(-{{\rm{i}}{F}}_{y}\chi \right)=\left(\begin{array}{ccc}\tfrac{1}{2}+\tfrac{\cos \chi }{2} & -\tfrac{\sin \chi }{\sqrt{2}} & \tfrac{1}{2}-\tfrac{\cos \chi }{2}\\ \tfrac{\sin \chi }{\sqrt{2}} & \cos \chi & -\tfrac{\sin \chi }{\sqrt{2}}\\ \tfrac{1}{2}-\tfrac{\cos \chi }{2} & \tfrac{\sin \chi }{\sqrt{2}} & \tfrac{1}{2}+\tfrac{\cos \chi }{2}\end{array}\right).\end{eqnarray} \tag{ 14 }$$

Here, χ is set to 0.2774π to give the optimal probe state. The simulated procedure consists of applying U, phase accumulation $\exp \left[{\rm{i}}\left({\theta }_{1}| 1\rangle \langle 1| +{\theta }_{2}| -1\rangle \langle -1| \right)\right]$ ($| 1\rangle ={\left(1,0,0\right)}^{\dagger }$, $| -1\rangle ={\left(0,0,1\right)}^{\dagger }$), and ${U}^{\dagger }$, followed by population detection with or without including noise.

When there exists no detection noise, the CFIM (equation (4)) of the aforementioned protocol is directly computed and the trace of its inverse is used to obtain the CRB of ${\left({\rm{\Delta }}{\boldsymbol{\theta }}\right)}^{2}$. Figure 5(a) compares the value of the corresponding result to the QCRB of the individual measurement scheme (equation (11)) for {θ₁, θ₂}∈(0, π), illustrated by the parameter $\zeta =-10{\mathrm{log}}_{10}\left[{({\rm{\Delta }}{\boldsymbol{\theta }})}^{2}/{({\rm{\Delta }}{{\boldsymbol{\theta }}}_{\mathrm{ind}})}^{2}\right]$. The region surrounded by the white dashed curve represents the {θ₁, θ₂}-space where the proposed scheme outperforms the individual estimation scheme. It shows that the proposed scheme works well even for ${\boldsymbol{\theta }}$ far away from zero. We emphasize that the probe state defined by equation (9) gives always the best QCRB for any ${\boldsymbol{\theta }}$. However, application of the reversed transformation ${U}^{\dagger }$ followed by a population measurement is not necessary the optimal measurement scheme when ${\boldsymbol{\theta }}$ is away from zero, which explains the deficiency of the scheme over some parameter space.

Zoom In Zoom Out Reset image size

When detection noise is present, we numerically simulate the estimation process of the two parameters $\{{\theta }_{1},{\theta }_{2}\}$ using 10⁴ three-mode (spin-1) atoms with a detection resolution (noise) of 14 atoms (typical numbers achievable in cold-atom experiments [47, 48]). For each pair of {θ₁, θ₂}, we first compute the probability of detecting an atom in the output mode $m,p\left(m| {\theta }_{1},{\theta }_{2}\right)={\left|\left\langle m\right|{U}^{\dagger }\exp \left[{\rm{i}}\left({\theta }_{1}\left|1\right\rangle \left\langle 1\right|+{\theta }_{2}\left|-1\right\rangle \left\langle -1\right|\right)\right]U\left|0\right\rangle \right|}^{2}$. For each run, we perform Monte-Carlo simulation on the outcome for each of the 10⁴ atoms according to the distribution of $p\left(m| {\theta }_{1},{\theta }_{2}\right)$ and obtain ${N}_{m}^{{\prime} }$ (the total number of particles in mode m without detection noise). We then add to ${N}_{m}^{{\prime} }$ a random detection noise featuring a normal distribution with an average of zero and a standard deviation of 14 to obtain N_m. The maximal likelihood method (which can saturate the CRB [46] in the asymptotic limit and is unbiased) is then used to estimate {θ₁, θ₂}. The likelihood function given by $L\left({\theta }_{1},{\theta }_{2}| {N}_{1},{N}_{0},{N}_{-1}\right)={\prod }_{m=\mathrm{1,0},-1}p{\left(m| {\theta }_{1},{\theta }_{2}\right)}^{{N}_{m}}$ is maximized by varying {θ₁, θ₂} to obtain the estimated $\{{\theta }_{1}^{{\prime} },{\theta }_{2}^{{\prime} }\}$. The simulation is repeated independently over ${N}_{\hat{M}}=1000$ times. The estimated results from the 1000 simulations are then used to obtain ${\left({\rm{\Delta }}{\boldsymbol{\theta }}\right)}^{2}$, whose ratio to the QCRB of the individual measurement scheme using the same number of particles, is shown in figure 5(b). Although detection noise degrades the sensitivity of the proposed scheme, the discussed scheme is seen to maintain its advantage over the individual measurement scheme over a large parameter space. Similar to typical single parameter estimation scenario, the position of the minimum ${\left({\rm{\Delta }}{\boldsymbol{\theta }}\right)}^{2}$ is seen to shift away from zero. The star in figure 5(b) (near ${\theta }_{1}={\theta }_{2}=0.3\pi$) denotes the position of the maximum precision for the scenario we consider, where it is $\zeta \simeq 0.6\,\mathrm{dB}$ more sensitive than the individual measurement scheme, but is 0.77 dB less than the optimal $\zeta \simeq 1.37\,\mathrm{dB}$ for noiseless detection.

In short, the influence of detection noise to a multi-mode Ramsey interferometer is similar to that to a single-mode Ramsey interferometer.

In summary, we show that the Ramsey interferometric scheme can be extended to estimation of multiple parameters (associated with commuting generators) using multi-mode pure states, if (but not iff) the multi-mode beam splitter U is orthogonal, i.e. all matrix elements of U are real and ${{UU}}^{\dagger }={\bf{1}}$. We then discuss how to obtain the optimal probe state, and how to construct U experimentally in a simple and scalable manner. We find that the proposed scheme remains intact even under particle number fluctuation and detection noise. The results of this study can be useful to applications in multi-mode optical sensing and quantum phase imaging.

This work is supported by the National Key R&D Program of China (Grant No. 2018YFA0306503 and No. 2018YFA0306504) and the NSFC (Grant No. 91 636 213, No. 91 736 311, No. 11 574 177, No. 91 836 302, and No. 11 654 001).

In this section, we show in detail that the proposed Ramsey-like multi-mode interferometric scheme with particle number measurement can always saturate the QCRB for small phase shift ${\boldsymbol{\Theta }}$, given that the matrix elements of the beam-splitting unitary transformation U are real (or U is orthogonal). The proposed scheme starts with splitting an initial state $\left|i\right\rangle$ by a unitary transformation U, followed by a phase accumulation process and a reversed transformation ${U}^{\dagger }$, and finally ends with measuring the projection probability in mode $| m\rangle$. The projection probability in mode $| m\rangle$ after the Ramsey interferometer can be explicitly written as

$$\begin{eqnarray}&&p\left(m| {\boldsymbol{\Theta }}\right)={\left|\left\langle m\right|{U}^{\dagger }\displaystyle \prod _{k=0}^{D}{{\rm{e}}}^{{\rm{i}}\left|k\right\rangle \left\langle k\right|{f}_{k}\left({\boldsymbol{\Theta }}\right)}U\left|i\right\rangle \right|}^{2}.\end{eqnarray} \tag{ A1 }$$

Omitting the third order corrections, a Taylor series expansion around ${\boldsymbol{\Theta }}\sim 0$ gives for $m\ne i$

$$\begin{eqnarray}&&p\left(m| {\boldsymbol{\Theta }}\right)\simeq \sum _{k,k^{\prime} =0}^{D}{f}_{k}\left({\boldsymbol{\Theta }}\right){f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)\left\langle m\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|m\right\rangle ,\end{eqnarray} \tag{ A2 }$$

and for m = i

$$\begin{eqnarray}&&p\left(i| {\boldsymbol{\Theta }}\right)\simeq 1+\sum _{k,k^{\prime} =0}^{D}{f}_{k}\left({\boldsymbol{\Theta }}\right){f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right){\left|\left\langle i\right|{U}^{\dagger }\left|k\right\rangle \right|}^{2}{\left|\left\langle k^{\prime} \right|U\left|i\right\rangle \right|}^{2}-\sum _{k=0}^{D}{f}_{k}{\left({\boldsymbol{\Theta }}\right)}^{2}{\left|\left\langle i\right|{U}^{\dagger }\left|k\right\rangle \right|}^{2}.\end{eqnarray} \tag{ A3 }$$

The derivatives $p\left(i| {\boldsymbol{\Theta }}\right)$ with respect to any ${{\rm{\Theta }}}_{l}$ is given by

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{\partial }}p\left(i| {\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}=2\left[\displaystyle \sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)| \left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle {| }^{2}\left. \langle i\right|{U}^{\dagger }\left|{k}^{{\rm{{\prime} }}}\right. \rangle {| }^{2}-\sum _{k=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}{f}_{k}\left({\boldsymbol{\Theta }}\right){\left|\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \right|}^{2}\right],\end{eqnarray} \tag{ A4 }$$

and the derivatives with respect to any ${{\rm{\Theta }}}_{l}$ for $m\ne i$ can be calculated as

$$\begin{eqnarray}\begin{array}{rcl}\displaystyle \frac{\partial p\left(m| {\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}} & = & \sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}{f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)\left\langle m\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|m\right\rangle \\ & & +\sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}{f}_{k}\left({\boldsymbol{\Theta }}\right)\left\langle m\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|m\right\rangle .\end{array}\end{eqnarray} \tag{ A5 }$$

The two terms in equation (A5) are equivalent if all the matrix elements of U ($\left\langle k\right|U\left|m\right\rangle$) are real numbers. In this case, equation (A5) can be simplified by summing up two terms as

$$\begin{eqnarray}&&\displaystyle \frac{\partial p\left(m| {\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}=2\sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}{f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)\left\langle m\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|m\right\rangle ,\quad \quad m\ne i.\end{eqnarray} \tag{ A6 }$$

By substituting equations (A2), (A3) and (A6) into the classical Fisher information matrix (CFIM) [46],

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{C}\left({\boldsymbol{\Theta }}\right)=\sum _{m}\displaystyle \frac{1}{p\left(m| {\boldsymbol{\Theta }}\right)}\displaystyle \frac{\partial p\left(m| {\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial p\left(m| {\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}},\end{eqnarray} \tag{ A7 }$$

we obtain (when ${\boldsymbol{\Theta }}\sim 0$ and $\mathrm{Im}\left[\langle k| U\left|m\right\rangle \right]=0$)

$$\begin{eqnarray}&&\displaystyle \sum _{m\ne i}\left\{\displaystyle \frac{4\left[\displaystyle \sum _{k,k^{\prime} ,k^{\prime\prime} ,k\prime\prime\prime =0}^{D}\tfrac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\tfrac{\partial {f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}{f}_{k^{\prime\prime} }\left({\boldsymbol{\Theta }}\right){f}_{k\prime\prime\prime }\left({\boldsymbol{\Theta }}\right)\left\langle m\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k^{\prime\prime} \right\rangle \left\langle k^{\prime\prime} \right|U\left|m\right\rangle \left\langle m\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k\prime\prime\prime \right\rangle \left\langle k\prime\prime\prime \right|U\left|m\right\rangle \right]}{\displaystyle \sum _{k^{\prime\prime} ,k\prime\prime\prime =0}^{D}{f}_{k^{\prime\prime} }\left({\boldsymbol{\Theta }}\right){f}_{k\prime\prime\prime }\left({\boldsymbol{\Theta }}\right)\left\langle m\right|{U}^{\dagger }\left|k^{\prime\prime} \right\rangle \left\langle k^{\prime\prime} \right|U\left|i\right\rangle \left\langle i\right|{U}^{\dagger }\left|k\prime\prime\prime \right\rangle \left\langle k\prime\prime\prime \right|U\left|m\right\rangle }\right\}.\end{eqnarray} \tag{ A8 }$$

The term m = i is missing from equation (A8) because $\displaystyle \frac{\partial p\left(i| {\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\propto {{\rm{\Theta }}}^{2}\sim 0$ while $p\left(i| {\boldsymbol{\Theta }}\right)\approx 1$. The numerator of equation (A8) can be factorized as

$$\begin{eqnarray}&&\begin{array}{l}4\left[\displaystyle \sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial {f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}}\langle m| {U}^{\dagger }\left|k\right\rangle \langle k| U\left|i\right\rangle \langle m| {U}^{\dagger }\left|k^{\prime} \right\rangle \langle k^{\prime} | U\left|i\right\rangle \right.\\ \quad \times \,\left.\displaystyle \sum _{k^{\prime\prime} ,k\prime\prime\prime =0}^{D}{f}_{k^{\prime\prime} }\left({\boldsymbol{\Theta }}\right){f}_{k\prime\prime\prime }\left({\boldsymbol{\Theta }}\right)\langle i| {U}^{\dagger }\left|k^{\prime\prime} \right\rangle \langle k^{\prime\prime} | U\left|m\right\rangle \langle i| {U}^{\dagger }\left|k\prime\prime\prime \right\rangle \langle k\prime\prime\prime | U\left|m\right\rangle \right],\end{array}\end{eqnarray} \tag{ A9 }$$

where the second summation cancels the denominator of equation (A8), giving a CFIM of the form

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{C}\left({\boldsymbol{\Theta }}\sim 0\right)\simeq 4\sum _{m\ne i}\sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial {f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}}\left\langle m\right|{U}^{\dagger }\left|k\right\rangle \left\langle k\right|U\left|i\right\rangle \left\langle m\right|{U}^{\dagger }\left|k^{\prime} \right\rangle \left\langle k^{\prime} \right|U\left|i\right\rangle .\end{eqnarray} \tag{ A10 }$$

To further simplify the formula, we use again the condition that the matrix elements of U ($\left\langle k\right|U\left|m\right\rangle$) are real numbers, thus $\left\langle m\right|{U}^{\dagger }\left|k\right\rangle =\left\langle k\right|U\left|m\right\rangle$ and $\left\langle k\right|U\left|i\right\rangle =\left\langle i\right|{U}^{\dagger }\left|k\right\rangle$. This gives

$$\begin{eqnarray}\begin{array}{c}\begin{array}{rcl}{{ \mathcal F }}_{l,n}^{C}\left({\boldsymbol{\Theta }}\sim 0\right) & \simeq & 4\displaystyle \sum _{m\ne i}\displaystyle \sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\mathop{\underbrace{\left|m\right. \rangle \left. \langle m\right|}}\limits_{m\ne i}{U}^{\dagger }\left|{k}^{{\rm{{\prime} }}}\right. \rangle \left. \langle {k}^{{\rm{{\prime} }}}\right|U\left|i\right. \rangle \\ & = & 4\left[\displaystyle \sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\left({\bf{1}}-\left|i\right. \rangle \left. \langle i\right|\right){U}^{\dagger }\left|{k}^{{\rm{{\prime} }}}\right. \rangle \left. \langle {k}^{{\rm{{\prime} }}}\right|U\left|i\right. \rangle \right]\\ & = & 4\left[\displaystyle \sum _{k=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\left|i\right. \rangle \right.\\ & & \left.-\displaystyle \sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}\left. \langle i\right|{U}^{\dagger }\left|k\right. \rangle \left. \langle k\right|U\left|i\right. \rangle \left. \langle i\right|{U}^{\dagger }\left|{k}^{{\rm{{\prime} }}}\right. \rangle \left. \langle {k}^{{\rm{{\prime} }}}\right|U\left|i\right. \rangle \right]\\ & = & 4\left[\displaystyle \sum _{k=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}-\sum _{k,{k}^{{\rm{{\prime} }}}=0}^{D}\displaystyle \frac{{\rm{\partial }}{f}_{k}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{l}}\displaystyle \frac{{\rm{\partial }}{f}_{{k}^{{\rm{{\prime} }}}}\left({\boldsymbol{\Theta }}\right)}{{\rm{\partial }}{{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}{\left|{\alpha }_{{k}^{{\rm{{\prime} }}}}\right|}^{2}\right]={{ \mathcal F }}_{l,n}^{Q}.\end{array}\end{array}\end{eqnarray} \tag{ A11 }$$

In the first line of equation (A11), the completeness of the basis $\left|i\right. \rangle \left. \langle i\right|+{\sum }_{m\ne i}\left|m\right. \rangle \left. \langle m\right|={\bf{1}}$ is invoked. The final result is identical to ${{ \mathcal F }}_{l,n}^{Q}$ given by equation (3). Thus this proves that the multi-mode Ramsey interferometer we consider here can always saturate the QCRB.

The above proof also explains why we limit the beam splitter to orthogonal matrix. This is crucial for the steps from equation (A5) to equation (A6) and from equation (A10) to (A11). Our proof by no means excludes the existence of U with non-real matrix elements which saturates the QCRB. But the general structures of such U are beyond our current knowledge.

In addition, the above proof also requires all elements $\left\langle m\right|U\left|k\right\rangle$ to be real even when $m,k\ne i$, which is the reason why it is insufficient to only require that the transformations from the input light in figure 2 to the outputs of all BS_k are real. Instead, all beam splitters involved must act as real 2 × 2 transformations.

As discussed in the main text, the probe state after phase accumulation takes the form $\left|{\psi }_{{\boldsymbol{\phi }}}\right\rangle ={\sum }_{k=0}^{D}{\alpha }_{k}{{\rm{e}}}^{{\rm{i}}{\phi }_{k}}\left|k\right\rangle$. If the parameters of interest ${\boldsymbol{\Theta }}\equiv \left\{{{\rm{\Theta }}}_{1},{{\rm{\Theta }}}_{2},\ \cdots ,{{\rm{\Theta }}}_{D}\right\}$ are linear combinations of ${\phi }_{k}$, the probe state can be expressed as $\left|{\psi }_{{\boldsymbol{\phi }}}\right\rangle ={\sum }_{k=0}^{D}{\alpha }_{k}{{\rm{e}}}^{{\rm{i}}\cdot {f}_{k}\left({\boldsymbol{\Theta }}\right)}\left|k\right\rangle$, where ${f}_{k}({\boldsymbol{\Theta }})$ are linear functions of ${\boldsymbol{\Theta }}$. The derivative of the state above w.r.t. ${{\rm{\Theta }}}_{l}$ is

$$\begin{eqnarray}&&\left|{\partial }_{{{\rm{\Theta }}}_{l}}{\psi }_{{\boldsymbol{\phi }}}\right\rangle ={\rm{i}}\displaystyle \sum _{k=0}^{D}\frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}{\alpha }_{k}{{\rm{e}}}^{{\rm{i}}\cdot {f}_{k}\left({\boldsymbol{\Theta }}\right)}\left|k\right\rangle .\end{eqnarray} \tag{ B1 }$$

For a pure state $\left|{\psi }_{{\boldsymbol{\phi }}}\right\rangle$, the matrix elements of ${{ \mathcal F }}^{Q}$ are explicitly given by [44, 45]

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{Q}=4\mathrm{Re}\left[\left\langle {\partial }_{{{\rm{\Theta }}}_{l}}{\psi }_{{\boldsymbol{\phi }}}\left|{\partial }_{{{\rm{\Theta }}}_{n}}{\psi }_{{\boldsymbol{\phi }}}\right.\right\rangle -\left\langle {\partial }_{{{\rm{\Theta }}}_{l}}{\psi }_{{\boldsymbol{\phi }}}\,| \,{\psi }_{{\boldsymbol{\phi }}}\right\rangle \left\langle {\psi }_{{\boldsymbol{\phi }}}\,| \,{\partial }_{{{\rm{\Theta }}}_{n}}{\psi }_{{\boldsymbol{\phi }}}\right\rangle \right],\end{eqnarray} \tag{ B2 }$$

where $l,n=1,2,\,\cdots ,$ and D. Substituting equation (B1) into equation (B2) gives the matrix elements of quantum Fisher information matrix (QFIM) of $| {\psi }_{{\boldsymbol{\phi }}}\rangle$

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{Q}=4\left[\sum _{k=0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}-\sum _{k,k^{\prime} =0}^{D}\displaystyle \frac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}\displaystyle \frac{\partial {f}_{k^{\prime} }\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{n}}{\left|{\alpha }_{k}\right|}^{2}{\left|{\alpha }_{k^{\prime} }\right|}^{2}\right].\end{eqnarray} \tag{ B3 }$$

In the case of ${{\rm{\Theta }}}_{k}={\theta }_{k}={\phi }_{k}-{\phi }_{0}$ $(k=1,2,\cdots ,D)$, one can choose ${f}_{0}({\boldsymbol{\Theta }})={\phi }_{0},{f}_{k}({\boldsymbol{\Theta }})={{\rm{\Theta }}}_{k}+{\phi }_{0}$ for $k=1,2,\cdot \cdot \cdot ,D$. The matrix elements of N-particle QFIM can be calculated with $\tfrac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}={\delta }_{k,l}$, leading to the result

$$\begin{eqnarray}&&{{ \mathcal F }}_{n,l}^{Q}=4N\left[{\left|{\alpha }_{l}\right|}^{2}{\delta }_{l,n}-{\left|{\alpha }_{l}\right|}^{2}{\left|{\alpha }_{n}\right|}^{2}\right].\end{eqnarray} \tag{ B4 }$$

Note that since there are D + 1 ϕ_k but only D Θ_k to be estimated, one of the ${f}_{k}({\boldsymbol{\Theta }})$ can chosen at will without affecting the final results. For example, for ${{\rm{\Theta }}}_{k}={\phi }_{k}-{\phi }_{0}$, one can also choose ${f}_{1}({\boldsymbol{\Theta }})={\phi }_{1},{f}_{0}({\boldsymbol{\Theta }})={\phi }_{1}-{{\rm{\Theta }}}_{1},{f}_{k}({\boldsymbol{\Theta }})={{\rm{\Theta }}}_{k}+{{\rm{\Theta }}}_{1}-{\phi }_{1}$ for $k=2,\cdot \cdot \cdot ,D$. Substituting the so-chosen ${f}_{k}({\boldsymbol{\Theta }})$ into equation (B3) gives the same results as equation (B4).

The inverse of equation (B4) can be obtained analytically as

$$\begin{eqnarray}&&\displaystyle \frac{1}{N}\left[\mathrm{diag}\left(\displaystyle \frac{1}{4{\left|{\alpha }_{1}\right|}^{2}},\displaystyle \frac{1}{4{\left|{\alpha }_{2}\right|}^{2}}\cdots \displaystyle \frac{1}{4{\left|{\alpha }_{D}\right|}^{2}}\right)\ +\ \displaystyle \frac{{\bf{G}}}{4{\left|{\alpha }_{0}\right|}^{2}}\right],\end{eqnarray} \tag{ B5 }$$

where G is a d × d all-ones matrix. Taking the trace of equation (B5) gives

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{{\boldsymbol{\theta }}}_{}\right)}^{2}\geqslant \displaystyle \frac{1}{N}\left[\displaystyle \frac{D}{4{\left|{\alpha }_{0}\right|}^{2}}+\sum _{k=1}^{D}\displaystyle \frac{1}{4{\left|{\alpha }_{k}\right|}^{2}}\right].\end{eqnarray} \tag{ B6 }$$

To find the optimal probe state and the corresponding total phase variance (optimal QCRB), we minimize equation (B6) under the normalization condition ${\sum }_{k=0}^{D}{\left|{\alpha }_{k}\right|}^{2}=1$. Setting the derivatives $\partial \left[{\left({\rm{\Delta }}{\boldsymbol{\theta }}\right)}^{2}\right]/\partial {\left|{\alpha }_{k}\right|}^{2}$ to zero for any k = 1, 2, ⋯, and D gives a set of equations

$$\begin{eqnarray}&&-\displaystyle \frac{1}{{\left|{\alpha }_{k}\right|}^{4}}+\displaystyle \frac{D}{{\left(1-{\left|{\alpha }_{1}\right|}^{2}-{\left|{\alpha }_{2}\right|}^{2}-\cdots -{\left|{\alpha }_{D}\right|}^{2}\right)}^{2}}=0.\end{eqnarray} \tag{ B7 }$$

Solving the equations above gives the optimal probe state described by

$$\begin{eqnarray}&&{\left|{\alpha }_{0}\right|}^{2}=\sqrt{D}/\left(D+\sqrt{D}\right),\end{eqnarray} \tag{ B8a }$$

$$\begin{eqnarray}&&{\left|{\alpha }_{k}\right|}^{2}=1/\left(D+\sqrt{D}\right),\quad \left(\mathrm{for}\,k\ne 0\right),\end{eqnarray} \tag{ B8b }$$

and a QCRB of

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{{\boldsymbol{\theta }}}_{\mathrm{opt}}\right)}^{2}={\left(D+\sqrt{D}\right)}^{2}/4N.\end{eqnarray} \tag{ B9 }$$

In the scenario where the parameters are defined as ${{\rm{\Theta }}}_{1}={\varphi }_{1}={\phi }_{1}-{\phi }_{0}$, ..., ${{\rm{\Theta }}}_{D}={\varphi }_{D}={\phi }_{D}-{\phi }_{D-1}$, or the phase difference between two neighboring modes, $\tfrac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}=1$ for k ≥ l and $\tfrac{\partial {f}_{k}\left({\boldsymbol{\Theta }}\right)}{\partial {{\rm{\Theta }}}_{l}}=0$ for k < l. The matrix elements of ${{ \mathcal F }}^{Q}$ therefore become

$$\begin{eqnarray}&&{{ \mathcal F }}_{l,n}^{Q}=4N\left[\sum _{k\geqslant \max \left(n,l\right)}^{D}{\left|{\alpha }_{k}\right|}^{2}-\left(\sum _{k^{\prime} \geqslant l}^{D}{\left|{\alpha }_{k^{\prime} }\right|}^{2}\right)\left(\sum _{k\geqslant n}^{D}{\left|{\alpha }_{k}\right|}^{2}\right)\right].\end{eqnarray} \tag{ B10 }$$

Taking the trace of the inverse of equation (B10) gives the lower bound of ${\left({\rm{\Delta }}{\boldsymbol{\varphi }}\right)}^{2}$

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{\boldsymbol{\varphi }}\right)}^{2}\geqslant \displaystyle \frac{1}{N}\left(\displaystyle \frac{1}{4{\left|{\alpha }_{D}\right|}^{2}}+\displaystyle \frac{1}{4{\left|{\alpha }_{0}\right|}^{2}}+\sum _{k=1}^{D-1}\displaystyle \frac{1}{2{\left|{\alpha }_{k}\right|}^{2}}\right).\end{eqnarray} \tag{ B11 }$$

Similarly, by minimizing equation (B11) under the normalization condition, the optimal probe reads

$$\begin{eqnarray}&&{\left|{\alpha }_{0}\right|}^{2}={\left|{\alpha }_{D}\right|}^{2}=1/\left[\sqrt{2}\left(D-1\right)+2\right],\end{eqnarray} \tag{ B12a }$$

$$\begin{eqnarray}&&{\left|{\alpha }_{k}\right|}^{2}=\sqrt{2}/\left[\sqrt{2}\left(D-1\right)+2\right],\quad \left(\mathrm{for}\,k\ne 0,D\right),\end{eqnarray} \tag{ B12b }$$

and the corresponding QCRB is found to be

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{{\boldsymbol{\varphi }}}_{\mathrm{opt}}\right)}^{2}=\displaystyle \frac{1}{4N}{\left[\sqrt{2}\left(D-1\right)+2\right]}^{2}.\end{eqnarray} \tag{ B13 }$$

If one measures $\left\{{\theta }_{1},{\theta }_{2},\ \ldots ,{\theta }_{D}\right\}$ with the input state given by equation (B8) and estimates $\left\{{\varphi }_{1},{\varphi }_{2},\ \cdots ,{\varphi }_{D}\right\}$ from the measured ${\theta }_{k}$, the ${\left({\rm{\Delta }}{\boldsymbol{\varphi }}\right)}^{2}$ is bounded by [46]

$$\begin{eqnarray}&&{\left({\rm{\Delta }}{\boldsymbol{\varphi }}\right)}^{2}\geqslant \mathrm{Tr}\left[{\boldsymbol{J}}{\left({{ \mathcal F }}_{{\boldsymbol{\theta }}}^{Q}\right)}^{-1}{{\boldsymbol{J}}}^{T}\right]=\displaystyle \frac{1}{4N}\left[{\left(1+\sqrt{D}\right)}^{2}+2\left(D-1\right)\left(\sqrt{D}+D\right)\right],\end{eqnarray} \tag{ B14 }$$

where J is the Jacobian matrix defined as ${{\boldsymbol{J}}}_{k,l}=\displaystyle \frac{\partial {\varphi }_{k}\left({\boldsymbol{\theta }}\right)}{\partial {\theta }_{l}}$. The result of equation (B14) is larger than that of equation (B13) for D > 1. Thus it is always better to estimate $\left\{{\varphi }_{1},{\varphi }_{2},\ \ldots ,{\varphi }_{D}\right\}$ directly using the probe state given by equation (B12).