Integrable quantum many-body sensors for AC field sensing

Quantum sensing is inevitably an elegant example of the supremacy of quantum technologies over their classical counterparts. One of the desired endeavors of quantum metrology is AC field sensing. Here, by means of analytical and numerical analysis, we show that integrable many-body systems can be exploited efficiently for detecting the amplitude of an AC field. Unlike the conventional strategies in using the ground states in critical many-body probes for parameter estimation, we only consider partial access to a subsystem. Due to the periodicity of the dynamics, any local block of the system saturates to a steady state which allows achieving sensing precision well beyond the classical limit, almost reaching the Heisenberg bound. We associate the enhanced quantum precision to closing of the Floquet gap, resembling the features of quantum sensing in the ground state of critical systems. We show that the proposed protocol can also be realized in near-term quantum simulators, e.g. ion-traps, with a limited number of qubits. We show that in such systems a simple block magnetization measurement and a Bayesian inference estimator can achieve very high precision AC field sensing.


Supplementary Material for Integrable Quantum Many-Body Sensors for AC Field Sensing
Utkarsh Mishra 1 and Abolfazl Bayat 1 1 Institute of Fundamental and Frontier Sciences,

FLOQUET HAMILTONIAN
We present here the analytical description of the model where, J > 0 is the nearest-neighbor spin-spin interaction, h 0 is a DC external magnetic field which is tunable, σ x/y i are Pauli matrices at site i, and the periodic boundary conditions is assumed, i.e.,σ x/y N +1 =σ x/y 1 . The system represents interacting spin-1/2 particles with Ising type interaction on a 1D lattice of length N in a transverse field to serve as a many-body probe for sensing a time-periodic magnetic field, h(t). The Hamiltonian in Eq. (S1) can be solved analytically [S1, S2]. The first step is to map the spin operators,σ i , into fermionic operators,ĉ † i (ĉ i ), defined via the following transformations: where,σ ± j = (σ x j ±σ y j )/2. The Hamiltonian as a result transformed into a quadratic-fermionic form It is to be noted that the fermionic Hamiltonian in Eq. (S3) is transnational invariant and therefore by applying a Fourier transformationĉ j = 1 √ N k e ikjĉ k , it can be written in k−space. The full Hamiltonian can be decomposed into the sum of the Hamiltonian for each k−mode i.e., H = k H k with H k being the Hamiltonian of the k th subspace whereσ y p andσ z p are pseudospin operators in the pseudospin basis and k is termed as quasi-momentum which takes the values k = π N , 3π N ,. . . , (N −1)π N for even N . The Hamiltonian in Eq. (S1) can be decomposed into the sum of even and odd parity-conserving Hamiltonians and the ground state of the system belongs to the even parity subspace for every finite N and it assumes BCS like form, given by [S3] where u k = sin(θ/2), v k = cos(θ/2), and θ = tan −1 J sin(k) h0+J cos(k) . Thus, when u k = 1, v k = 0 for all k, it corresponds to a state with all spins in the eigenbasis ofσ z with eigenvalue +1. Under the dynamics given by the Hamiltonian in Eq. (S3) the states in the two parity sectors as well as those with different momentum evolves independently. Thus, for the unitary time dynamics, it is enough to consider states {|0 , c † k c † −k |0 } which leads to time evolved state in the form [S4] where u k (t) and v k (t) are the solutions of the Schrödinger equation The above analysis applies to any general time-dependent function h(t). For stroboscopic dynamics, i.e., time-evolution of the system monitored in the steps t = nτ , the state of the system at any time t = nτ , can be obtained by repeated application of the unitary operator as where , is the the time-evolution operator for single time-period τ and T denotes time ordered product. For the Dirac-delta function, it is possible to find the effective Hamiltonian, known as Floquet Hamiltonian H F , which generates equivalent dynamics. Therefore, the equation for the ensuing dynamics, namely Eq. (S8), can be simplified in term of the Floquet Hamiltonian as where Once the time-dependent state |Ψ(t = nτ ) is known, the time-dependent magnetization at the stroboscopic time, t = nτ , is expressed as m z (nτ ) = Ψ(nτ )| 1 N N i=1σ z i |Ψ(nτ ) . By employing the Jordan-Wigner transformation and expressingσ z in terms ofσ ± , we have where, u k (nτ ) and u k (nτ ) are solutions of the Schrödinger equation given in Eq. (S7).

QUANTUM FISHER INFORMATION OF A BLOCK OF SIZE L
To calculate the quantum Fisher information of a subsystem of size L, we need to calculate the reduced density matrix of L sites (we consider L contiguous sites). The calculation of reduced density matrix requires partial tracing of complimentary degrees of freedoms. For quadratic Hamiltonians of the form given in Eq. (S3), this is accomplished by noting the relation between the density matrix and the matrix of single-particle correlations of L sites [S5]. The form of the reduced density matrix which reproduces correct correlation matrix on L sites is given by ρ L = e −Ω /Z, where Ω is quadratic in fermionic operators with energies i i.e., Ω = L i=1 iĉ † iĉ i and Z is the normalization constant. By breaking the complex fermions into Majorana basis defined asĉ i = 1 2 (a 2i−1 + ia 2i ) andĉ † i = 1 2 (a 2i−1 − ia 2i ), the density matrix of blocks of size L can be represented as free fermionic Gaussian state, given by where Here Ω is a 2L × 2L real antisymmetric matrix and a ≡ (a 1 , . . . , a 2L ) T is a 2L-dimensional array of Majorana fermions. The Gaussian state in Eq. (S11) is completely characterized by a two-point correlation matrix Γ whose elements are given by The Γ matrix is an antisymmetric matrix. The elements of the Γ matrix can be obtained in terms of C and I matrix where, In terms of these matrices, we have [S6] where, [·] ( [·]) represents the real (imaginary) part, i, j = 1, . . . , L, δ i,j is discrete Kronecker delta function, and the elements of the time-dependent correlation matrix, C(t), and anomalous correlation matrix, F(t), are given as and Once the Γ matrix is known, the quantities of interest can be expressed in terms of the Γ matrix. For example, the symmetric logarithmic derivative has been obtained for the Gaussian states, Eq. (S11), and in the Majorana basis it has the following form where K is a Hermitian anti-symmetric matrix of dimension 2L × 2L, ζ is a real vector and Λ is a real number. The matrix elements of the K in the eigenbasis of Γ = 2L r=1 γ r |γ r γ r | are given by [S7] with (∂ h1 Γ) rs = r|∂ h1 Γ|s and the partial derivative is taken with respect to the parameter to be estimated. By substitutingL into Tr[ρ LL 2 ], the QFI in the eigenbasis of Γ is expressed as [S7, S8] We have used this expression to obtain the Fisher information of a block of size L in our time-dependent model. The expression of QFI is valid for all parameters except when γ r = γ s = ±1. At these values of γ s, the density matrix ρ L becomes singular and is not well defined.

OPTIMAL VERSUS SUB-OPTIMAL MEASUREMENTS
In order to obtain the exact form of the optimal POVM, we first define the symmetric logarithmic derivative operatorL to satisfy By inserting the spectral decomposition form of ρ L in Eq. (S18), one can easily obtain the following expression forL where, |∂ h1 λ p = ∂|λp ∂h1 . The QFI then can be written as F Q (h 1 ) = Tr[ρ LL 2 ] and the eigenvectors ofL provide the optimal POVM projectors {Π r } [S9]. In general, the symmetric logarithmic derivative operatorL and thus its eigenvectors, which are the optimal measurement basis, depend on the unknown parameter h 1 . Therefore, finding the optimal measurement basis is one of the big challenges in quantum estimation theory which often hinders saturating the quantum Cramér-Rao bound. In fact, in most of the cases, the bound is only achievable when sophisticated adaptive methods for updating the measurement basis are employed [S10-S14].
Here we outline the calculation of symmetric logarithmic derivative for a single and two-qubit density matrix. A general single-qubit state of the system is written as ρ 1 = 1 2 (I + m. σ), where m = Tr(ρ 1 σ) and σ = (σ x ,σ y ,σ z ). It is to be noted that the Hamiltonian H(t) is invariant under fermionic parity transformation, z . This implies, as shown in [S15, S16], that m x , m y = 0. Thus, we get a single-site density matrix which is diagonal in the eigenbasis ofσ z i.e., in {| ↑ , | ↓ } basis. The symmetric logarithmic derivative in this basis is given byL = On the other hand, the eigenvectors ofL are {| ↑ , | ↓ }. If the set of POVM is constructed using the projections onto the eigenvectors ofL, then the classical Fisher information F C is given by where, p ↑ = ↑ |ρ 1 | ↑ and p ↓ = ↓ |ρ 1 | ↓ . A further simplification of F C gives F C = F Q . For the density matrix of two nearest-neighbor sites, one also needs to calculate the two-point correlatorsσ s i ⊗ σ s i+1 (s, s = x, y, z). It can be seen that the correlators such asσ x i ⊗σ z i+1 andσ y i ⊗σ z i+1 vanishes due to invariance of the Hamiltonian under parity transformation [S17, S18]. Since periodic boundary conditions are assumed, the nearest-neighbor state is independent of which two neighboring sites are chosen for constructing the density matrix. The two-site density matrix of the system, therefore, is given by where, The other non-zero elements are given as u 32 = u * 23 , and u 41 = u * 14 . The non-zero correlatorsσ s i ⊗σ s i+1 can be obtained using the formalism presented in [S19]. Once the two-site density matrix is obtained, the symmetric logarithmic derivative for the two-qubit state can be calculated. We find that the symmetric logarithmic derivative with respect to the state ρ 2 (t) is given byL The eigenvectors of symmetric logarithmic derivative are given by where . As two of the eigenvectors namley | 1 and | 2 depends on the parameter to be estimated i.e., on h 1 , the classical Fisher information obtained from the measurement of Eq. (S24) may not be equal to the quantum Fisher information [S20].
The QFI is only a bound in the Cramér-Rao inequality and it is not guaranteed to be saturated unless an optimal measurement basis, as well as an optimal estimator, are chosen. As mentioned before, the measurement basis is determined by eigenvectors of the symmetric logarithmic derivative operatorL. In a general case, the optimal measurement basis depends on the sensing parameter h 1 which by definition is unknown making the optimal sensing impractical. Normally, to overcome this problem, one has to update the measurement basis adaptively by extracting partial information about the unknown parameter using a sequence of non-optimal measurement basis [S10-S14]. Due to practical constraints, even if the optimal measurement basis is known (e.g. from an adaptive strategy) its implementation may not be feasible. Therefore, one of the desired issues in quantum metrology problems is to find a suitable measurement basis that is close to the optimal one and is independent of h 1 .
We consider a simple, though sub-optimal, measurement which is independent of h 1 . The measurement is the block magnetization, which for a block of size L takes L + 1 outcomes from O 1 = +L (when all the qubits are | ↑ ), O 2 = L − 2 (when except one qubit the rest are in the state | ↑ ) until O L+1 = −L (when all the qubits are | ↓ ). Each of the outcomes O r appears with the probability p r . Then one can get the corresponding classical Fisher information F C . Note that the block magnetization is easily measurable in ion traps [S21-S23] and superconducting quantum devices [S24-S26]. In Figs. S1(a)-(h) we plot both the CFI, computed for the block magnetization, and the QFI as a function of time t = nτ for various block sizes. In all these plots the control field is fixed to be h 0 /J = 1 while the left and the right panels represent the results for h 1 /J = 0.1 and h 1 /J = 0.2, respectively. For the sake of clarity, the right panels are only shown for the later times. Interestingly, although the block magnetization is not the optimal measurement the resulted F C takes values greater than unity. This suggests that such a simple measurement can be used for efficient sensing.

SENSING SQUARE PULSES
To see the generality of our approach, we also consider square pulsed form of the periodic field, given by the following equation where w characterizes the width of the pulse over an interval of τ . The Floquet evolution operator over a time-period is given by The local density matrix of the system reaches to a steady state under the AC field of the form given in Eq. (S25). Similar to the case of Dirac-delta kick pulse, the Floquet gap takes its minimum along a straight line on the h 0 −h 1 plane and F ss Q exhibits peaks along the same line. In Fig. S2(a) we plot the steady state quantum Fisher information F ss Q with respect to h 0 and h 1 for a block of size L = 4. The quantum Fisher information shows its peak exactly a straight line. In Fig. S2(b), we plot both max(F ss Q ) and min(∆ F ). The lines of max(F ss Q ) and min(∆ F ) coincides. This shows the generality of the fact that the vanishing Floquet gap results in a higher sensitivity in steady state quantum metrology.

COMPARISON WITH OTHER PROTOCOLS
Now in this section, we outline some of the key points about our protocol addressing its efficacy as compared to other existing protocols for AC-field sensing. We made this comparison with the two main existing schemes, namely spin echo and GHZ-based schemes.
Spin echo utilizes a coherent superposition of spin states and a sequence of external pulses. The typical platform for realizing the spin echo mechanism for AC-field sensing is the nitrogen vacancy centers [S27-S29]. However, the scheme is limited by the maximum time interval needed to accumulate phase and the quality of the coherent superposition of the spin states. To further enhance the precision, one can increase the number of spins. Although, it is crucial that the spins remain non-interacting as any interaction between them acts like decoherence and decreases the precision. In fact, in the case of spin ensembles, the interaction is inevitable and one has to utilize a sophisticated pulse sequence and dynamical decoupling scheme [S30]. Our proposal, however, takes a fundamentally different route as it exploits the interaction between the particles to drive the subsystems into a steady state. This naturally spares us from any