Symmetric Logarithmic Derivative of Fermionic Gaussian States

Abstract

In this article, we derive a closed form expression for the symmetric logarithmic derivative of Fermionic Gaussian states. This provides a direct way of computing the quantum Fisher Information for Fermionic Gaussian states. Applications range from quantum Metrology with thermal states to non-equilibrium steady states with Fermionic many-body systems.

1. Introduction

Quantum metrology or quantum parameter estimation is the theory that studies the accuracy by which a physical parameter of a quantum system can be estimated through measurements and statistical inference. In many physical scenarios, quantities which are to be estimated may not be directly observable, either due to experimental limitations or on account of fundamental principles. When this is the case, one needs to infer the value of the variable after measurements on a given probe. This is essentially a parameter estimation problem whose solution may be found using methods from classical estimation theory [1] or, when quantum systems are involved, from its quantum counterpart [2]. Quantum metrology finds applications in a range of diverse fields, from fundamental physics, such as improving frequency and time standards [3,4,5], estimating parameters in quantum field theory [6,7], improving the accuracy of gravitational wave interferometry [8,9], to applied physics, such as thermometry [10,11], spectroscopy [12,13], imaging [14,15,16], magnetic field detection [17,18] navigation [19,20] and remote sensing [21,22]. Exploiting the quantum nature of physical systems provides a remarkable advantage in enhancing the accuracy of estimation problems, and exploring this possibility plays a pivotal role in the current swift development of quantum technology [23,24,25,26,27,28,29,30,31,32,33,34]. Simultaneous quantum estimation of multiple parameters provides better precision over individual estimation strategies with equivalent resources [35,36]. This fact has sparked interest in multiparameter quantum metrology in a variety of scenarios [35,36,37,38].

Recent advances in quantum metrology have shown that the accuracy in parameter estimation can be enhanced by employing peculiar quantum many-body state as a probe [34,39]. Conversely, quantum metrology tools may well be exploited in the characterisation of many-body systems. Noteworthy instances of many-body quantum systems are those experiencing quantum phase transitions. Indeed, quantum parameter estimation, with its intimate relation with geometric information, provides a novel and promising approach in the characterisation of equilibrium [40,41,42,43,44,45,46,47,48,49] and out-of-equilibrium quantum critical phenomena [50,51,52,53,54].

Fermionic Gaussian states play a major role in the derivation of exact and approximate solution of many-body problems of Fermionic and spin systems. Deriving closed form expressions of quantities involved in parameter estimation problems for many-body quantum systems is a major challenge. This work addresses this task in the special, yet relevant, case of arbitrary Fermionic Gaussian states.

The solution of a parameter estimation problem amounts to find an estimator, i.e., a mapping $\widehat{\mathit{\lambda }}=\widehat{\mathit{\lambda }}(x_1,x_2,\dots )$ from the set $\chi$ of measurement outcomes into the space of parameters $\mathit{\lambda }\in \mathcal{M}$. Optimal estimators in classical estimation theory are those saturating the Cramer-Rao (CR) inequality,

$${\mathrm{Cov}}_{\mathit{\lambda }}\left[\widehat{\mathit{\lambda }}\right]\ge J^c{\left(\mathit{\lambda }\right)}^{-1},$$

(1)

which poses a lower bound on the mean square error ${\mathrm{Cov}}_{\lambda }{\left[\widehat{\mathit{\lambda }}\right]}_{\mu \nu }=E_{\lambda }\left[{(\widehat{\lambda }-\lambda )}_{\mu }{(\widehat{\lambda }-\lambda )}_{\nu }\right]$ in terms of the Fisher information (FI)

$$J_{\mu \nu }^c\left(\mathit{\lambda }\right)={\int }_{\chi }d\widehat{\mathit{\lambda }}\left(x\right)p\left(\widehat{\mathit{\lambda }}\right|\lambda ){\partial }_{\mu }logp\left(\widehat{\mathit{\lambda }}\right|\lambda ){\partial }_{\nu }logp\left(\widehat{\mathit{\lambda }}\right|\lambda ).$$

(2)

For unbiased estimators, the mean square error is equal to the covariance matrix ${\mathrm{Cov}}_{\mathit{\lambda }}{\left[\widehat{\mathit{\lambda }}\right]}_{\mu \nu }=E_{\mathit{\lambda }}\left[{\widehat{\mathit{\lambda }}}_{\mu }{\widehat{\mathit{\lambda }}}_{\nu }\right]-E_{\mathit{\lambda }}\left[{\widehat{\mathit{\lambda }}}_{\mu }\right]E_{\mathit{\lambda }}\left[{\widehat{\mathit{\lambda }}}_{\nu }\right].$ The expression (1) should be understood as a matrix inequality. In general, one writes

$$\mathrm{tr}\left(W{\mathrm{Cov}}_{\mathit{\lambda }}\left[\widehat{\mathit{\lambda }}\right]\right)\ge \mathrm{tr}\left(W{J}^c{\left(\mathit{\lambda }\right)}^{-1}\right),$$

where W is a given positive definite cost matrix, which allows the uncertainty cost of different parameters to be weighed unevenly.

In the classical estimation problem, both in the single parameter case, and in the multiparameter one, the bound is saturable in the limit of an infinite number of repetitions of an experiment using the maximum likelihood estimator [55]. However, an interesting difference between multiparameter and single parameter metrology arises due to the correlation between parameters. Indeed, it may well happen that the resulting Fisher information matrix is non-diagonal. This means that the estimators for the parameters will not be independent. In a separate scheme in which all parameters except the ${\lambda }_{\mu }$ are perfectly known, the single parameter CR bound implies that the uncertainty of estimating ${\lambda }_{\mu }$ is lower bounded by $\mathrm{Var}\left(\widehat{\mathit{\lambda }}\right)\ge 1/J_{\mu \mu }^c$. On the other hand, in the simultaneous scenario in which all parameters are estimated at the same time, one finds $\mathrm{Var}\left(\widehat{\mathit{\lambda }}\right)\ge {\left(J^c{\left(\mathit{\lambda }\right)}^{-1}\right)}_{\mu \mu }$. From basic algebra of positive-definite matrices, we have that ${\left(J^c{\left(\mathit{\lambda }\right)}^{-1}\right)}_{\mu \mu }\ge 1/J^c{\left(\mathit{\lambda }\right)}_{\mu \mu }$, with equality holding only in the case when all off-diagonal elements vanish. Since asymptotically the CR bound is saturable, it implies that the equivalence between the simultaneous and separate scheme in the limit of a large number of experiment repetitions can only hold if F is a diagonal matrix, and hence there are no statistical correlations between the estimators [56].

Clearly, for any real positive definite matrix, one can perform an orthogonal rotation to a new basis in which the matrix is diagonal. This simply means that there are always linear combinations of the parameters for which the diagonality conditions hold. This choice should be, however, contrasted with the physical opportunity of performing such a rotation, as the choice of the parameters we are interested in may arise as a result of physical considerations and in this sense determine a preference in a specific basis.

While the fundamental objects in classical Fisher information are parameter-dependent probability-distribution of the data, the fundamental objects involved in the quantum estimation problem are the density matrices $\rho \left(\mathit{\lambda }\right)$ labelled by $\mathit{\lambda }\in \mathcal{M}$. In the quantum scenario, we therefore face the additional challenge of determining the optimal measurement for extracting most of the information on the parameters of interest from the quantum states. In the single parameter case, the situation is relatively simple. Maximization of the classical Fisher information over all quantum measurements yields the quantity referred to as the quantum Fisher information (QFI). The key object involved in the calculation of the QFI is the so-called symmetric logarithmic derivative (SLD), L, which is implicitly defined as the Hermitian operator satisfying the equation

$$\frac{d\rho }{d\lambda }=\frac{1}{2}\left(\rho L+L\rho \right).$$

(3)

The QFI can be calculated using the formula:

$$J=\mathrm{Tr}\left(\rho L^2\right).$$

(4)

One can always choose the projective measurement in the eigenbasis of the SLD that yields FI equal to the QFI. Hence, the QFI determines the ultimate achievable precision of estimating the parameter on density matrices $\rho \left(\mathit{\lambda }\right)$ in the asymptotic limit of an infinite number of experiment repetitions. In a multiparameter scenario, a direct generalization of single parameter CR bound leads to the multiparameter quantum Cramer-Rao (QCR) bound [2,26,57], which reads

$$\mathrm{tr}\left(W\mathrm{Cov}\left(\widehat{\mathit{\lambda }}\right)\right)\ge \mathrm{tr}\left(W{J}^{-1}\right),$$

(5)

where

$$J_{\mu \nu }=\frac{1}{2}\mathrm{Tr}\rho \{L_{\mu },L_{\nu }\}$$

(6)

is the quantum Fisher information matrix (QFIM), W is the cost matrix, and $L_{\mu }$ is the SLD implicitly defined by Equation (3), with $\rho$ derived with respect to the parameter ${\lambda }_{\mu }$.

Unlike the single parameter case, in the multiparameter scenario, the QCR bound cannot always be saturated. Intuitively, this is due to the incompatibility of the optimal measurements for different parameters. A sufficient condition for the saturation is indeed $[L_{\mu },L_{\nu }]=0$, which is however not a necessary condition. Within the comprehensive framework of quantum local asymptotic normality (QLAN) [58,59,60,61], a necessary and sufficient condition for the saturation of the multiparameter QCR bound is given by [62]

$${\mathcal{U}}_{\mu \nu }=-\frac{i}{4}\mathrm{Tr}\rho [L_{\mu },L_{\nu }]=0\forall \mu ,\nu ,$$

(7)

and it is known as compatibility condition [62]. In the context of quantum information geometry, and quantum holonomies of mixed states, ${\mathcal{U}}_{\mu \nu }$ is known as mean Uhlmann curvature [53].

In this paper, we derive a closed form expression of the SLD of Fermionic Gaussian states, which are of fundamental importance in the analysis of steady-states of both equilibrium and non-equilibrium quantum many-body systems, and their applications to quantum metrology.

The paper is organised as follows. In the next section, we shortly review the main properties of Fermionic Gaussian states. In Section 3, an explicit form for the calculation of the SLD is derived. In the last section, the conclusions are drawn.

2. Fermionic Gaussian State

We review here the main properties of Fermionic Gaussian states (FGSs). Let us consider a system of n fermionic particles described by creation and annihilation operators $c_j^{†}$ and $c_j$. These operators obey the canonical anticommutation relations,

$$\{c_j,c_k\}=0\{c_j,c_k^{†}\}={\delta }_{jk}.$$

(8)

Let us define the Hermitian Majorana operators as

$${\omega }_{2j-1}:=c_j+c_j^{†},{\omega }_{2j}:=i(c_j-c_j^{†}),$$

(9)

which are generators of a Clifford algebra, and satisfy the following anti-commutation relations

$$\{{\omega }_j,{\omega }_k\}=2{\delta }_{jk}.$$

(10)

Fermionic Gaussian states are defined as states that can be expressed as

$$\rho =\frac{e^{-\frac{i}{4}{\mathit{\omega }}^T\Omega \mathit{\omega }}}{Z},Z:=\mathrm{Tr}\left[e^{-\frac{i}{4}{\mathit{\omega }}^T\Omega \mathit{\omega }}\right],$$

(11)

where $\Omega$ is a $2n\times 2n$ real antisymmetric matrix and $\mathit{\omega }:={({\omega }_1\cdots {\omega }_{2n})}^T$ is a 2n-dimensional array of Majorana Fermions. As any antisymmetric real matrix, $\Omega$ can be cast in the following canonical form by an orthogonal matrix Q, i.e.,

$$\begin{array}{ccc}\hfill \Omega & =Q^T\underset{k=1}{\overset{n}{⨁}}\left(\begin{array}{cc}0& {\Omega }_k\\ -{\Omega }_k& 0\end{array}\right)Q\hfill & \hfill Q^T=Q^{-1},\end{array}$$

(12)

where $\pm i{\Omega }_k$ are $\Omega$’s eigenvalues. Let

$$\mathbf{z}={(z_1,\cdots ,z_{2n})}^T:=Q\mathit{\omega }$$

(13)

be the vector of Majorana fermions in the eigenmode representation. Hence,

$$\begin{array}{cc}\hfill \rho & =\frac{1}{Z}\prod _k\left[cosh\left(\frac{{\Omega }_k}{2}\right)-isinh\left(\frac{{\Omega }_k}{2}\right)z_{2k-1}{z}_{2k}\right],\hfill \end{array}$$

(14)

$$\begin{array}{cc}\hfill Z& =\prod _{k}2cosh\left(\frac{{\Omega }_k}{2}\right).\hfill \end{array}$$

(15)

Gaussian states are completely specified by the two-point correlation matrix

$${\mathrm{\Gamma }}_{jk}:=1/2\mathrm{Tr}\left(\rho [{\omega }_j,{\omega }_k]\right),\mathrm{\Gamma }={\mathrm{\Gamma }}^{†}=-{\mathrm{\Gamma }}^T,$$

(16)

which is an imaginary antisymmetric matrix. Let’s recall some basic properties of the correlation function. As for any Fermionic State, all odd-order correlation functions of FGSs are zero, due to the parity super-selection rule. In FGSs, all even-order correlations, higher than two, can be obtained from $\mathrm{\Gamma }$ by Wick’s theorem [63], i.e.,

$$\mathrm{Tr}(\rho {\omega }_{k_1}{\omega }_{k_2}\dots {\omega }_{k_{2p}})=\mathrm{Pf}\left({\mathrm{\Gamma }}_{k_1{k}_2\cdots k_{2p}}\right),1\le k_1<\dots <k_{2p}\le 2n$$

(17)

and ${\mathrm{\Gamma }}_{k_1{k}_2\cdots k_{2p}}$ is the corresponding $2p\times 2p$ submatrix of $\mathrm{\Gamma }$. $\mathrm{Pf}{\left({\mathrm{\Gamma }}_{k_1{k}_2\cdots k_{2p}}\right)}^2=det\left({\mathrm{\Gamma }}_{k_1{k}_2\cdots k_{2p}}\right)$ is the Pfaffian. An especially useful case is the four-point correlation function

$$\mathrm{Tr}\left(\rho {\omega }_j{\omega }_k{\omega }_l{\omega }_m\right)=a_{jk}{a}_{lm}-a_{jl}{a}_{km}+a_{jm}{a}_{kl},$$

(18)

where $a_{jk}:={\mathrm{\Gamma }}_{jk}+{\delta }_{jk}$. As

$${\mathrm{\Gamma }}_{jk}=\frac{2i}{Z}\frac{\partial Z}{\partial {\Omega }_{jk}},$$

(19)

one can show that

$$\begin{array}{c}\hfill \mathrm{\Gamma }=tanh\left(i\frac{\Omega }{2}\right).\end{array}$$

(20)

The correlation matrix is diagonal in the same basis of $\Omega$ and its eigenvalues read ${\gamma }_k=tanh({\Omega }_k/2)$. Hence

$$\begin{array}{cc}\hfill \rho & =\prod _{k=1}^n\frac{1-i|{\gamma }_k|z_{2k-1}{z}_{2k}}{2},\hfill \end{array}$$

(21)

where $|{\gamma }_k|\le 1$. Hence, the Gaussian fermionic state can be factorised into a tensor product $\rho ={⨂}_k{\rho }_k$ of density matrices of the eigenmodes ${\rho }_k:=\frac{1-i|{\gamma }_k|z_{2k-1}{z}_{2k}}{2}$. Note that, for ${\gamma }_k=\pm 1$, one has ${\Omega }_k=\pm \infty$, making the definition (11) of Gaussian states not well defined, unlike Equation (21), showing that the latter offer an appropriate parameterisation even in those extremal points. Notice that $|{\gamma }_k|=1$ corresponds to a fermionic mode ${\tilde{c}}_k=1/2(z_{2k-1}+z_{2k})$ being in a pure state, as it is clear from the following explicit expression for the purity of the states ${\rho }_k$:

$$\begin{array}{c}\hfill \mathrm{Tr}\left[{\rho }_k^2\right]=\frac{det{\left[2cosh\left({\Omega }_k\right)\right]}^{\frac{1}{2}}}{det\left[2cosh\left(\frac{{\Omega }_k}{2}\right)\right]},\end{array}$$

(22)

$$\begin{array}{c}\hfill \mathrm{Tr}\left[{\rho }^2\right]=\frac{det{\left[2cosh\left(i\Omega \right)\right]}^{\frac{1}{2}}}{det\left[2cosh\left(i\frac{\Omega }{2}\right)\right]}=\sqrt{det\left(\frac{1+{\mathrm{\Gamma }}^2}{2}\right)}.\end{array}$$

(23)

3. Symmetric Logarithmic, Derivative of Fermionic Gaussian States

We will derive here an explicit formula for the calculation of the SLD for Fermionic Gaussian states. To this end, we review a useful expression adapted from [64] needed for the derivation of the symmetric logarithmic derivative of density matrices in the exponential form

$$\rho =e^{D\left(\mathit{\lambda }\right)}.$$

(24)

Clearly, a Gaussian Fermionic state can be expressed in the exponential form (24) by identifying

$$D=-\frac{i}{4}{\mathit{\omega }}^T·\Omega ·\mathit{\omega }-11lnZ.$$

(25)

Notice that the above parameterisation is well defined only in the case of full-rank density matrices. As usual, for the case of extremal conditions $|{\gamma }_k|=1$, an eigenvalue of the correlation function should be carried out as a limiting procedure.

The starting point is the expression derived in equation (2.1) of [65] for derivative of density operators

$$\begin{array}{c}\hfill \dot{\rho }={\int }_0^1{e}^{sD}\dot{D}{e}^{(1-s)D}ds,\end{array}$$

(26)

where dots represent derivatives with respect to a parameter $\lambda$. One can use the nested- commutator relation

$$\begin{array}{cc}\hfill e^{D}A{e}^{-D}& =A+[D,A]+\frac{1}{2!}\left[D,[D,A]\right]+\cdots \hfill \\ & =\sum _{n=0}^{\infty }\frac{1}{n!}{\mathcal{C}}^n\left(A\right)=e^{\mathcal{C}}\left(A\right),\hfill \end{array}$$

(27)

where ${\mathcal{C}}^n\left(A\right)$, a linear operation on A, denotes the nth-order nested commutator $\left[D,\dots ,[D,A]\right]$, with ${\mathcal{C}}^0\left(A\right)=A$. Applying this relation to the expression (26) leads to

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \dot{\rho }{\rho }^{-1}& =\dot{D}+\frac{1}{2!}[D,\dot{D}]+\frac{1}{3!}\left[D,[D,\dot{D}]\right]+\cdots \hfill \\ & =\sum _{n=0}^{\infty }\frac{1}{(n+1)!}{\mathcal{C}}^n\left(\dot{D}\right)=h\left(\mathcal{C}\right)\left(\dot{D}\right),\hfill \end{array}\end{array}$$

(28)

where h is the generating function of the expansion coefficients in Equation (28),

$$\begin{array}{c}\hfill h\left(t\right)=1+\frac{t}{2!}+\frac{t^2}{3!}+\cdots =\frac{e^t-1}{t}.\end{array}$$

(29)

Using the definition of symmetric logarithmic derivative, i.e.,

$$\dot{\rho }=\frac{1}{2}\left(L\rho +\rho L\right),$$

(30)

and that of density matrix in exponential form (24), one gets

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \dot{\rho }{\rho }^{-1}& =\frac{1}{2}\left(L+e^{D}L{e}^{-D}\right)\hfill \\ & =\frac{1}{2}\left(L+\sum _{n=0}^{\infty }\frac{1}{n!}{\mathcal{C}}^n\left(L\right)\right)=r\left(\mathcal{C}\right)\left(L\right),\hfill \end{array}\end{array}$$

(31)

where the generating function is $r\left(t\right)=(e^t+1)/2$. Suppose that the SLD adopts the form,

$$\begin{array}{cc}\hfill L& =\sum _{n=0}^{\infty }{f}_n{\mathcal{C}}^n\left(\dot{D}\right)=f\left(\mathcal{C}\right)\left(\dot{D}\right),\hfill \end{array}$$

(32)

with the generating function

$$\begin{array}{c}\hfill f\left(t\right)=f_0+f_{1}t+f_2{t}^2+\cdots \end{array}$$

(33)

to be determined. Plugging Equation (32) into Equation (31) yields

$$\begin{array}{c}\hfill \dot{\rho }{\rho }^{-1}=r\left(\mathcal{C}\right)\left[f\left(\mathcal{C}\right)\left(\dot{D}\right)\right]=r\circ f\left(\mathcal{C}\right)\left(\dot{D}\right)=r·f\left(\mathcal{C}\right)\left(\dot{D}\right),\end{array}$$

(34)

where the identity $r\circ f=r·f$, between the combination, $r\circ f$, of the two functions and their simple product, $r·f$, arises from ${\mathcal{C}}^n\left({\mathcal{C}}^m\left(A\right)\right)={\mathcal{C}}^{n+m}\left(A\right)$. Comparing Equation (34) with Equation (28) leads to the following relation between generating functions,

$$\begin{array}{c}\hfill f\left(t\right)=\frac{h\left(t\right)}{r\left(t\right)}=\frac{tanh(t/2)}{t/2}=\sum _{n=0}^{\infty }\frac{4(4^{n+1}-1)B_{2n+2}}{(2n+2)!}{t}^{2n},\end{array}$$

(35)

where $B_{2n+2}$ is the $(2n+2)$th Bernoulli number. Comparing Equations (33) with (35), we have

$$\begin{array}{c}\hfill f_n=\left\{\begin{array}{cc}{\displaystyle \frac{4(4^{n/2+1}-1)B_{n+2}}{(n+2)!},}\hfill & \mathrm{for} \mathrm{even} n,\hfill \\ 0,\hfill & \mathrm{for} \mathrm{odd} n.\hfill \end{array}\right.\end{array}$$

(36)

The vanishing of the odd-order of $f_n$s is a consequence of the Hermiticity of L, which makes $f\left(t\right)$ an even function. From the definition (25) of D, one straightforwardly finds

$$\dot{D}=-\frac{1}{2}{\mathit{\omega }}^T·\dot{\Omega }·\mathit{\omega }-\frac{\dot{Z}}{Z}11,$$

(37)

which shows that $\dot{D}$ is itself a quadratic function of the Majorana Fermion operators, where $\dot{\Omega }$ is an antisymmetric real matrix given by

$$\begin{array}{cc}\hfill \dot{\Omega }& =Q^T\underset{k=1}{\overset{n}{⨁}}\left(\begin{array}{cc}0& {\dot{\Omega }}_k\\ -{\dot{\Omega }}_k& 0\end{array}\right)Q+i[R,\Omega ],\hfill \end{array}$$

(38)

with $R:=i{Q}^T\dot{Q}$. Therefore, for a Gaussian Fermionic state, the operator $\dot{D}$ can be written in a canonical form in terms of the Majorana Fermions of its own eigenmodes, as

$$\dot{D}=-\frac{i}{4}\sum _k{\tilde{\Omega }}_k[{\tilde{z}}_{2k-1},{\tilde{z}}_{2k}]-\frac{\dot{Z}}{Z}11=\sum _k{\tilde{\Omega }}_k\left({\tilde{c}}_k^{†}{\tilde{c}}_k-\frac{1}{2}\right)-\frac{\dot{Z}}{Z}11,$$

(39)

where ${\tilde{c}}_k:=\frac{1}{2}({\tilde{z}}_{2k-1}+i{\tilde{z}}_{2k})$, ${\tilde{c}}_k^{†}:=\frac{1}{2}({\tilde{z}}_{2k-1}-i{\tilde{z}}_{2k})$ are the ordinary annihilation and creation operators of the eigenmodes of $\dot{D}$, and ${\tilde{\Omega }}_k$ the corresponding eigenvalues. It is straightforward to derive the commutation relations between $\dot{D}$ and Fermionic operators,

$$\begin{array}{c}\hfill \left[\dot{D},{\tilde{c}}_k\right]=-{\tilde{\Omega }}_k{\tilde{c}}_k,\left[\dot{D},{\tilde{c}}_k^{†}\right]={\tilde{\Omega }}_k{\tilde{c}}_k^{†},\end{array}$$

(40)

and for quadratic operators, also

$$\begin{array}{c}\hfill \left[\dot{D},{\tilde{c}}_j^{†}{\tilde{c}}_k\right]=({\tilde{\Omega }}_j-{\tilde{\Omega }}_k){\tilde{c}}_j^{†}{\tilde{c}}_k,\left[\dot{D},{\tilde{c}}_j^{†}{\tilde{c}}_k^{†}\right]=({\tilde{\Omega }}_j+{\tilde{\Omega }}_k){\tilde{c}}_j^{†}{\tilde{c}}_k^{†}.\end{array}$$

(41)

Consequently, one finds

$$\begin{array}{c}\hfill f\left(\mathcal{C}\right)\left({\tilde{c}}_j^{†}{\tilde{c}}_k\right)=f({\tilde{\Omega }}_j-{\tilde{\Omega }}_k){\tilde{c}}_j^{†}{\tilde{c}}_k,\end{array}$$

(42)

$$\begin{array}{c}\hfill f\left(\mathcal{C}\right)\left({\tilde{c}}_j^{†}{\tilde{c}}_k^{†}\right)=f({\tilde{\Omega }}_j+{\tilde{\Omega }}_k){\tilde{c}}_j^{†}{\tilde{c}}_k^{†}.\end{array}$$

(43)

The above expression, plugged into formula (32), shows that L is at most quadratic in the Fermionic operators.

Due to the quadratic dependence of L on the Fermionic operator, clearly L can be expressed as a quadratic polynomial in the Majorana Fermions in the following form

$$L_{}=:\frac{1}{2}{\mathit{\omega }}^T·K·\mathit{\omega }+{\mathit{\zeta }}_{}^T·\mathit{\omega }+\eta ,$$

(44)

where $K:={\left\{K_{jk}\right\}}_{jk=1}^{2n}$ is $2n\times 2n$ Hermitian anti-symmetric matrix, $\mathit{\zeta }:={\left\{{\zeta }^k\right\}}_{k=1}^{2n}$ is a $2n$ real vector, and $\eta$ is a real number. Note that any odd-order correlation function for a Gaussian Fermionic state vanishes identically, then

$$\langle {\omega }_k\rangle =\mathrm{Tr}\left(\rho {\omega }_k\right)=0\forall k=1\cdots 2n.$$

(45)

By differentiating the above equation, one readily shows that the linear term in (44) is identically zero

$$0=\frac{d}{d\lambda }\mathrm{Tr}\left(\rho {\omega }_k\right)=\mathrm{Tr}\left({\omega }_k\dot{\rho }\right)=\frac{1}{2}\mathrm{Tr}\left({\omega }_k\{L,\rho \}\right)=\mathrm{Tr}\left(\rho \{{\mathit{\zeta }}^T\mathit{\omega },{\omega }_k\}\right)={\zeta }^k,$$

where ${\zeta }^k$ is the k-th component of $\mathit{\zeta }$, and, in the fourth equality, one takes into account that the third order correlations vanish. The quantity $\eta$ can be determined from the trace preserving condition, i.e.,

$$0=\frac{d}{d\lambda }\mathrm{Tr}\rho =\mathrm{Tr}\left(\dot{\rho }\right)=\mathrm{Tr}\left(\rho L\right),$$

(46)

which, after plugging in Equation (44), leads to

$$\eta =-\frac{1}{2}\mathrm{Tr}\left(\rho {\mathit{\omega }}^{T}K\mathit{\omega }\right)=\frac{1}{2}\mathrm{Tr}\left(K_{}\mathrm{\Gamma }\right).$$

(47)

In order to determine K, let’s take the differential of ${\mathrm{\Gamma }}_{jk}=1/2\mathrm{Tr}\left(\rho [{\omega }_j,{\omega }_k]\right)$, then

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\dot{\mathrm{\Gamma }}}_{jk}=\frac{1}{2}\mathrm{Tr}\left(\dot{\rho }[{\omega }_j,{\omega }_k]\right)& =\frac{1}{4}\mathrm{Tr}\left(\{\rho ,L_{}\}[{\omega }_j,{\omega }_k]\right)\hfill \\ & =\frac{1}{8}\mathrm{Tr}\left(\{\rho ,{\mathit{\omega }}^T{K}_{}\mathit{\omega }\}[{\omega }_j,{\omega }_k]\right)+\frac{{\eta }_{}}{4}\mathrm{Tr}\left(\rho [{\omega }_j,{\omega }_k]\right)\hfill \\ & =\frac{1}{8}\sum _{lm}{K}_{}^{lm}\mathrm{Tr}\left(\{\rho ,[{\omega }_l,{\omega }_m]\}[{\omega }_j,{\omega }_k]\right)+\frac{{\eta }_{}}{2}{\mathrm{\Gamma }}_{jk}\hfill \\ & ={(\mathrm{\Gamma }{K}_{}\mathrm{\Gamma }-K_{})}_{jk}+\frac{1}{2}\left[\eta -\frac{1}{2}\mathrm{Tr}\left(K_{}\mathrm{\Gamma }\right)\right]{\mathrm{\Gamma }}_{jk},\hfill \end{array}\end{array}$$

(48)

where the last equality is obtained with the help of Equation (18) and using the antisymmetry of $\mathrm{\Gamma }$ and K under the exchange of j and k. Finally, according to Equation (47), the last term vanishes and we obtain the following (discrete time) Lyapunov equation

$$\dot{\mathrm{\Gamma }}=\mathrm{\Gamma }{K}_{}\mathrm{\Gamma }-K_{}.$$

(49)

The above equation can be formally solved by

$$K_{}={({\mathrm{Ad}}_{\mathrm{\Gamma }}-11)}^{-1}\left(\dot{\mathrm{\Gamma }}\right),$$

(50)

where ${\mathrm{Ad}}_{\mathrm{\Gamma }}\left(X\right):=\mathrm{\Gamma }X{\mathrm{\Gamma }}^{†}$ is the adjoint action. In the eigenbasis of $\mathrm{\Gamma }$, (i.e., $\mathrm{\Gamma }|j\rangle ={\gamma }_j|j\rangle$), it reads

$$\langle j|K_{}|k\rangle ={\left(K_{}\right)}_{jk}=\frac{{\left(\dot{\mathrm{\Gamma }}\right)}_{jk}}{{\gamma }_j{\gamma }_k-1}=-\frac{{\dot{\Omega }}_k}{2}{\delta }_{jk}+tanh\frac{{\Omega }_j-{\Omega }_k}{2}\langle j|\dot{k}\rangle ,$$

(51)

where, in the second equality, we made use of the relation ${\gamma }_k=tanh({\Omega }_k/2)$, which yields the following diagonal ${\left(\dot{\mathrm{\Gamma }}\right)}_{jj}=(1-{\gamma }_j^2){\dot{\Omega }}_j$ and off-diagonal ${\left(\dot{\mathrm{\Gamma }}\right)}_{jk}=({\gamma }_k-{\gamma }_j)\langle j|\dot{k}\rangle$ terms. This expression is well defined everywhere except for ${\gamma }_j={\gamma }_k=\pm 1$, where the Gaussian state $\rho$ becomes singular (i.e., it is not full rank). In this condition, the expression (51) for the SLD L may become singular. Nevertheless, the boundness of the function $|tanh\frac{{\Omega }_j-{\Omega }_k}{2}|\le 1$ in (51) shows that such a singularity is relatively benign. One can show that the condition ${\gamma }_j={\gamma }_k=\pm 1$ produces, at most, removable singularities in the quantum Fisher Information Matrix (cf. [66]). This allows the quantum Fisher information matrix to be extended by continuity from the set of full-rank density matrices to the subset with ${\gamma }_j={\gamma }_k=\pm 1$.

4. Conclusions

In this work, we derived a general expression for the symmetric logarithmic derivative of an arbitrary Fermionic Gaussian state. We obtained a compact expression in terms of a correlation matrix of a FGS, which allows for the calculation of the quantum Fisher information. This provides a way of assessing the ultimate precision of parameter estimation problems in many-body systems involving Fermionic Gaussian states, through the quantum Cramer-Rao bound. Moreover, the expression of the SLD allows for the explicit derivation of the eigenbasis associated with the optimal quantum measurement to be performed for the estimation of a parameter of interest. The generality of the method also offers a way of evaluating the so-called compatibility condition in multiparameter quantum estimation problems [62]. Indeed, due to the quantum nature of the underlying probe, the multiparameter estimation problem may not saturate the multiparameter quantum Cramer-Rao bound. Unlike classical estimation problems, the non-commutativity of the observables involved in the optimal quantum measurements associated with independent parameters may prevent the quantum CR bound from being saturated [29]. An explicit quantitative condition [62] for such a compatibility can be easily derived, once the formula for the SLD is given. The general framework presented provides a way of easily assessing the above-mentioned quantities. Moreover, the explicit expression of the SLD, analogous to Bosonic Gaussian estimation problems [31], can be exploited, in combination with efficient numerical algorithms, to find an optimal Fermionic Gaussian probe that minimises the overall multiparameter estimation problem [67].