Entanglement dynamics with string measurement operators

We explain how to apply a Gaussian-preserving operator to a fermionic Gaussian state. We use this method to study the evolution of the entanglement entropy of an Ising spin chain following a Lindblad dynamics with string measurement operators, focusing on the quantum-jump unraveling of such Lindbladian. We find that the asymptotic entanglement entropy obeys an area law for finite-range string operators and a volume law for ranges of the string which scale with the system size. The same behavior is observed for the measurement-only dynamics, suggesting that measurements can play a leading role in this context.

Understanding the role of entanglement in quantum many-body systems is a subject under intensive study [1,2].While the ground-state entanglement properties have been widely investigated and are relatively well understood, there are still several open questions regarding the stabilization of quantum correlations during the non-equilibrium dynamics.By exploiting a quasi-particle description it has been shown that in integrable systems, after a sudden quench, the entanglement entropy linearly increases in time, and eventually settles to a stationary value which linearly scales with the system size (volume-law scaling) [3,4].A similar volume-law scaling occurs in ergodic short-range systems [5], while in the presence of a strong disorder the situation is more involved, due to the appearance of many-body localization (see Ref. [6] for a review).
A natural extension of the above scenario is for a quantum system coupled to an external environment, a setup which is crucial to applications in quantum technologies and quantum computing, since any experimental platform is subject to noise-induced decoherence [7,8].In the hypothesis of weak and Markovian coupling with the bath, a good modeling for this setup is to assume the evolution of the reduced density matrix of the system to be ruled by a master equation in a Lindblad form [9].The system-bath framework can be also interpreted as a process in which the environment behaves as a classical stochastic process performing random measurements on the quantum system [10,11].The outcome of any realization of this process is a pure-state stochastic quantum trajectory and the density matrix obtained by averaging over such trajectories obeys a Lindblad-type dynamic evolution [12,13].
Recently, an increasing number of works have been focusing on MIPTs in fermionic systems described by quadratic Hamiltonians on a lattice, undergoing quantum trajectories under the action of random measurements of onsite quadratic operators [44][45][46][47][48][49][50][51][52][53][54][55][56][57][58][59][60][61].One of the reasons of this interest relies in the Gaussian structure of the states of such systems that, being entirely determined by two-point correlation functions, are suitable for a semi-analytical treatment up to large lattice sizes (of the order of hundreds of sites) [62][63][64][65][66][67][68].In particular, it has been found that the asymptotic volume-law scaling of the entanglement following a unitary evolution [69,70] is unstable under local measurements, independently of the measurement strength [44,46], and a MIPT may emerge.For example, the Ising chain exhibits a size-dependent crossover towards a subextensive or an area-law entanglement regime, depending on the system parameters, as well as on the measurement rate and the type of unraveling [47][48][49][50][51]71].In a version of this model with temporal noise also in the Hamiltonian, this crossover has been analytically proved to be a true transition by means of the replica method [72,73].Similarly, it has been found that entanglement transitions can arise, even in the absence of a coherent evolution, from the competition between different measurement protocols [74].In the same context of monitored fermionic models, a couple of recent works pointed out the importance of having long-range interactions to stabilize larger amounts of entanglement, if the system undergoes random onsite measurements [53,54] (see [75,76] for a similar phenomenon in quantum circuits).Similar results have been found also for shortranged Hamiltonian in the presence of long-range continuous measurement processes [77].
Inspired by these latter findings, here we follow a complementary strategy and address the role of non-local string measurement operators in the dynamics of entanglement for a quadratic Hamiltonian system.The purpose of this paper is twofold.On the one hand, we collect and review all the technical details to implement the dynamics of Gaussian states in the presence of Gaussian-preserving measuring operators.On the other hand, we discuss the effect of non-local measurement operators on the dynamics of the entanglement entropy, while keeping the Hamiltonian with short-range interactions.To this aim, we implement a quantum-jump unraveling of the Lindblad dynamics [9,10] for the Ising chain (mapped to the fermionic Kitaev chain) using strings of Pauli matrices as Lindblad operators (mapped to non-local two-point fermionic operators).The main result is that, when the string range r is comparable with the system size L, contrary to what happens in the presence of local jump operators, the measurement process can stabilize an asymptotic volume-law for the entanglement entropy.As a consequence, we cannot exclude that a MIPT can occur as a function of r.Quite interestingly, the same qualitative conclusion applies if we switch off the Hamiltonian and consider a measurement-only Lindblad dynamics induced by string operators, similarly to what has been studied in Ref. [42] for measurement-only models with nonlocal measurements of Pauli strings.
The paper is organized as follows.In Section 2 we briefly review the theory of quadratic fermionic Hamiltonians on a lattice.In Section 3 and Appendix A we present the technical details required to implement the dynamics of a Gaussian state under the effect of Gaussianpreserving exponentials of quadratic operators.We first set the general framework, and then shed light on the dynamical equations and on the normalization of the state.In Section 4 we apply such a machinery to string operators, a particular class of operators which generalize the one discussed in Ref. [49,50].Subsequently, in Section 5 we present and discuss numerical results for the Ising chain with non-local string measurement operators.Finally, our conclusions are drawn in Section 6.The other appendices contain some details on the quantum-jump measurement protocol (Appendix B) and an alternative approach to determine the evolved state with the string operator (Appendix C).
We consider a generic free-fermion system on a lattice with L sites, described by the quadratic Hamiltonian where ĉ( †) j denotes the anticommuting annihilation (creation) operator on the j-th site.To ensure the Hermiticity of Ĥ, the coefficient matrices obey Q † = Q and P T = −P.
In this context, a Gaussian state has the generic form and N a normalization constant.Here u and v are L × L matrices that obey the following symplectic unitarity conditions where I denotes the L × L identity matrix.These conditions are equivalent to say that the matrix is unitary and preserves the fermionic anticommutation relations: In fact it defines a new set of fermionic anticommuting operators that enjoy the property γκ |ψ〉 = 0 , ∀ κ ∈ {1, . . ., L} , meaning that |ψ〉 is the vacuum of the γ fermions, |ψ〉 ≡ |0〉 γ (the matrix is the so-called Bogoliubov transformation.For more details see, e.g., Refs.[66,78]).We remark that, because of the fermionic nature of the ĉ-operators, the matrix Z in Eq (3) is antisymmetric, as can be proved by exploiting Eqs.(4) [78].In this context, the Nambu spinor notation is useful.
The state in Eq. ( 2) is Gaussian, enjoys the Wick's theorem, and thus it is completely determined by two-point correlation functions.Using Eq. ( 7) one can show that these correlation functions have a simple expression in terms of the matrix: An important case where the state has the Gaussian form in Eq. ( 2) is the ground state of the Hamiltonian (1).The Hamiltonian reads The 2L × 2L Hermitian matrix is called the Bogoliubov-de Gennes Hamiltonian matrix, which can be diagonalized through the implementation of a Bogoliubov transformation of the form (5), such that The above transformation defines a new set of fermionic quasi-particles γκ [Eq.(6)], according to which the Hamiltonian of Eq. ( 1) can be written as Its ground state is a Gaussian state of the form in Eq. ( 2), that corresponds to the vacuum of the γκ operators [see Eq. ( 7)].
We note that the above formalism can be easily generalized to the case of time-varying coefficients Q i, j (t) and P i, j (t) in Eq. (1), by admitting a time dependence of the Bogoliubov transformation matrix = (t).Moreover, as we shall see below in Section 3, the form of a generic Gaussian state Eq. ( 2) is preserved by the application of the exponential of any operator quadratic in the ĉ fermions.

Evolving Gaussian states
We now discuss the evolution of a generic Gaussian state |ψ〉 under the effect of an operator of the form being s ∈ , ξ = {1, i}, and a generic Hermitian quadratic operator.For simplicity of notations, we assume all coefficients in Â to be real, but our analysis can be easily generalized to the complex case.Due to the fermionic anticommutation rules, we have D T = D, O T = −O.Notice that, if ξ = i, the operator in Eqs.(12) describes the real-time evolution according to a (pseudo)-Hamiltonian Â.On the other hand, if ξ = 1, it describes an imaginary-time evolution (related to the construction of the thermal ensemble at inverse temperature β = 1/s), that can be thought of as the analytic continuation of the purely imaginary evolution.As discussed in Section 4, in some particular cases, the action played by such kind of operator on a given state can be intended as part of a measurement process.
In what follows we discuss how to evaluate, for any (real) symmetric Â operator, the sevolved state |ψ(s)〉 = M (s) |ψ〉 , taking as initial state |ψ〉 a generic Gaussian state of the form (2). In Appendix A we show that, because of the quadratic structure of Â, the state |ψ(s)〉 keeps the Gaussian form of Eq. (2): and u(s), v(s) submatrices of the s-evolved Bogoliubov matrix (s).In the next Subsection we are going to show how to evaluate the matrix (s) (an alternative derivation is provided in Appendix A).

Bogoliubov matrix evolution
To characterize the evolved state |ψ(s)〉, one only needs to know how the Bogoliubov matrix (s) evolves with M (s).Exploiting the identity M (s) −1 = M (−s), we can write where |ψ〉 is a generic Gaussian state given by Eq. ( 2).The operator is defined through the evolved fermions obeying anticommutation relations ĉq (s), ĉ † p (s) = δ p,q .Note that, for ξ = 1, the transformation M (s) is not unitary and the conjugation operation does not commute with the evolution one, hence ĉ † q (s) to determine the evolved state, one has to calculate Ẑ(s).This can be done by evaluating Given the quadratic structure of the operator Â, the above commutators read Therefore At this point, we make the Ansatz that the fermionic operators depend on s only through a Bogoliubov-like matrix The expression for the submatrices of the ev (s) is given in Eq. ( 29).Note that, while for ξ = i one has ūev (s) = [u ev (s)] * and vev (s) = [v ev (s)] * (so that ev (s) is a true Bogoliubov transformation), in general this is not true.The above Ansatz is equivalent to assume that the s-evolved state |ψ(s)〉 is a Bogoliubov vacuum γµ (s) |ψ(s)〉 = 0, ∀µ.By substituting this Ansatz in Eq. ( 20), we obtain a set of differential equations for the elements of the evolved matrix: with the initial condition: u ev q,µ (0) = [ū ev q,µ (0)] * = u q,µ , and v ev q,µ (0) = [v ev q,µ (0)] * = v q,µ .In a compact form, these equations read being The matrix ev (s) is not a priori characterizing a fermionic state in that it should be unitary, to guarantee the normalization of the state, and should obey Eqs. ( 4), to impose the anticommutation relations for the Ψ(s).While the former requirement is automatically fulfilled for ξ = i (because of the unitarity of the evolution operator), in the case ξ = 1 it must be forced by performing a QR-decomposition, as discussed in details in the next subsection.Once the unitarity of the matrix is imposed, Eqs. ( 4) are satisfied if the Bogoliubov matrix preserves the symplectic structure of Eq. ( 6), and this happens for Â being Hermitian.By substituting the normalized s-evolved Bogoliubov matrix in Z(s) [cf.Eq. ( 14)], we obtain the expression for |ψ(s)〉.If the above conditions are satisfied, this is a Gaussian fermionic state, as the one in Eq. (13).A more detailed discussion about this point is provided in Appendix A.

QR decomposition
We now discuss how to restore the unitarity of the s-evolved Bogoliubov matrix ev (s), when considering a norm-non preserving evolution with ξ = 1, that is in general obtained by QRdecomposing the s-evolved Bogoliubov matrix ev (s).This procedure keeps the Z(s) matrix in Eq. ( 14) -and then the Gaussian state -invariant, as we are going to see below.Nevertheless, it is necessary to apply the known formulae [78] for evaluating the entanglement entropy and the local observables.
The QR-decomposition is a procedure that allows one to decompose a generic matrix as K = U Q R, with U Q a unitary matrix and R an upper triangular one.In what follows, we provide an argument to show that, if Â = Â † , the unitary matrix obtained by QR-decomposing the sevolved Bogoliubov matrix ev (s) describes a fermionic state, i.e., it has the symplectic form of Eq. ( 6).The straightforward way would be to check it by construction, a procedure that requires some involved calculations.Here we choose a different strategy.Since R is positive definite (we have checked it numerically), the QR decomposition must be unique.Therefore, we assume U Q to have the right structure and then we show that it exists a decomposition that is compatible with this assumption.
Let us consider a symmetric operator Â defined in Eq. (12b) and the associated coefficient matrix defined in Eq. (23).For simplicity let us assume Â to have nearest-neighbors interactions, i.e.D i, j = O i, j = 0 if j ̸ = i + 1, and let as assume D i,i+1 = O i,i+1 (e.g. the Kitaev Hamiltonian).The most general case does not have an analytical solution, but the main argument of this section should hold independently.Because of the symmetries of D and O, we have that {D, O} = 0. 3 We notice that By exploiting the block diagonal form of Eq. ( 25) and the identities we can write where, for shortness of notation, we have posed Notice that [ch ± (s)] † = ch ± (s), while sh † (s) = −sh(s).Moreover, since {D, O} = 0, we have We can QR-decompose ev (s) making the Ansatz that the (s) preserves the symplectic form of Eq. ( 6), so that (30) can be rewritten as In what follows, we show that there exists a decomposition compatible with this Ansatz, by imposing Eq. ( 4) and that, because of the unicity, it must be the only possible.Let us assume to have derived r 1 (s) by construction.We can explicitly write (s) = ev (s)R(s) −1 , exploiting the fact that the inverse of an upper block triangular matrix remains an upper block triangular matrix: We can substitute these expressions in Eq. ( 4) to obtain Since [ū ev ] T (s) u ev (s) + [v ev ] T (s) v ev (s) = I, the above equation provides a relation between two of the blocks of the matrix R(s): By imposing the second constraint in Eq. ( 4), we find the relation for r 3 (s) As a last comment we notice that the matrix Z(s) in Eq. ( 14) -and then the Gaussian stateare left invariant by the application of the QR decomposition:

Number-preserving operators
We now consider a generic number-preserving operator where I, J are two sets of indices labeling certain sites on the system, and apply MI,J = e β ÂI,J , ( to a Gaussian state.We focus on this particular example of the theory presented in Section 3.1, since it generalizes the construction discussed in Ref. [49,50].By exploiting the anticommutation relations, it is possible to show that Ân where the last equality holds, since j and m run in the same set.Therefore, if I ∩ J = ∅, we have Â2 I,J = 0 and MI,J = 1 + β ÂI,J , where 1 denotes the identity operator.Otherwise In what follows, without loss of generality, we set β ≡ log(α + 1)/α, so that MI,J = 1 + ÂI,J reduces to that studied in Ref. [49,50] for I = J = {i}.As discussed in Appendix B, the application of such a MI,J can be thought as the measurement step of a quantum-jump dynamics.We want to evaluate the effect of this operator on a Gaussian state Following the procedure discussed in Section 3.1, we write the evolved state as where we introduced the transformed fermionic operators k( †) q = MI,J ĉ( †) q M −1 I,J .By exploiting the fact that M −1 I,J = 1 − (α + 1) −1 ÂI,J , we derive the exact expression of k( †) This leads to the following expressions for the evolved Bogoliubov matrix: As discussed in Appendix C, an alternative derivation can be obtained from the evolved antisymmetric matrix Z.

Two-site operator
We now consider a non-local (nl) operator with support on two sites (I = J = {i, j}): Ânl = 1+ Ânl I,J .Without losing in generality, in the following we will focus on one dimensional systems and set j = i + r (with r > 0).As we shall see in the next Section, a measurement operator as the one in Eq. ( 46) may unveil a richer physics than local ones.In fact, in the spin-1/2 language, it corresponds to a string operator.The evolved Bogoliubov matrix can be derived from Eqs. (44) which, in this case, read ūev, nl q,µ = u * q,µ + (δ q,i + δ q,i+r )(u * i,µ + u * i+r,µ ) , (47c) v ev, nl q,µ = v q,µ + (δ q,i + δ q,i+r )(v i,µ + v i+r,µ ) .(47d)

Quantum jumps with string operators
In this section we use the machinery presented above to describe the evolution of a quantum Ising chain under the effect of a quantum-jump dynamics generated by M nl I,J [cf.Eq. ( 46) and below], with emphasis on the evolution of the entanglement entropy.
We consider the Hamiltonian where σα j (α = x, y, z) are the usual spin-1/2 Pauli matrices acting on the j-th site, while J > 0 denotes the ferromagnetic spin-spin coupling and h the transverse magnetic field.In the thermodynamic limit, the ground state exhibits a zero-temperature continuous transition at |h/J| = 1, separating a paramagnetic (symmetric) phase, for |h/J| > 1, from a ferromagnetic (symmetry-broken) phase, for |h/J| < 1 [79,80].Through a Jordan-Wigner transformation [81,82], the Ising model ( 48) can be mapped into the so-called Kitaev chain This Hamiltonian is of the form (1), with Hereafter we work in units of ħ h = 1, set J = 1 as the energy scale, and assume periodic boundary conditions for fermions.We implement the quantum-jump protocol described in Appendix B on the fermionic model above, by initializing the system in the ground state of the Hamiltonian (50) with a given value of h and then Trotterizing the time evolution in steps δt ∝ (4Lγ) −1 , to ensure convergence, being γ the strength of the measurements.The emerging dynamics is one of the possible unravelings [10] of the Lindblad master equation [9] being Âj (r) = (ĉ † j + ĉ † j+r )(ĉ j + ĉj+r ) the same operator as in Eq. ( 46), with i ≡ j + r. 4 In the spin language, this operator is a sum of strings of Pauli matrices, of length r: with σ± j = ( σx j ± i σ y j )/2.Hereafter, we refer to it as the string operator.We point out that our framework is similar to the one presented in Ref. [42], although with some important qualitative differences.The former is described by a random sequence of stabilizer-formalism measurements of Pauli string operators.These operators have a fixed range r and are taken from a random distribution of strings, located around different sites, and made of different sets of onsite Pauli operators.
In this work, in contrast, we are following an unraveling of the quantum master equation Eq. ( 52).As detailed in Appendix B, we apply random quantum-jump measurements of string operators 1 + Âj (r), located at different sites, but with the structure fixed as the one in Eq. ( 53).This is sufficient to achieve a measurement-induced nontrivial dynamics, with important consequences on the entanglement generation.
It is important to emphasize that the operators 1 + Âj (r) have the exponential form given by Eqs.(38) and (40) (see also Sec. 4).In this way, the state keeps its Gaussian form under the application of the quantum jumps, and a numerical analysis of quite large system sizes is possible.

Entanglement entropy
Let us focus our attention on the entanglement entropy.Given a pure state |ψ〉, there is an unambiguous way to quantify its quantum correlations.We divide the system in two subchains, one with length ℓ and the other with length L − ℓ.After evaluating the density matrix ρ ℓ = Tr L−ℓ |ψ〉 〈ψ| reduced to the subsystem ℓ, the entanglement entropy is quantified by its von Neumann entropy [7,8]: For the Ising model we are discussing herewith, it is known that after, a sudden quench in a purely Hamiltonian dynamics, the entanglement entropy approaches a stationary value that increases linearly with ℓ, thus obeying a volume-law behavior [3,4].As discussed in Refs.[44][45][46][47][48][49][50][51], this dynamical volume law is generally destroyed by any arbitrary small local weak measurement process.In this scenario, the asymptotic entanglement entropy is conjectured to exhibit either a logarithmic scaling or an area-law behavior, depending on the strength of measurements and of the Hamiltonian parameters.We are going to show that the situation is quite different for the case of the nonlocal operators defined in Eq. ( 53).
In practice, to avoid spurious contributions from classical correlations, we evaluate Eq. ( 54) on each stochastic quantum trajectory |ψ t 〉 (defined in Appendix B), and then average over different stochastic realizations. 5For simplicity of notation, we denote the ensemble-averaged entanglement entropy as and where the symbol • • • denotes the average over different stochastic realizations.We also define the asymptotic entanglement entropy as We remark that, in all the simulations shown hereafter, we keep the size of the partition fixed to ℓ = L/4 and study the scaling with L.
Below we present two different scenarios for fixed γ = 0.5 and for r = 1, 4, L/2: the measurement-only dynamics (Section 5.2) and the full quantum-jump dynamics (Section 5.3).Finally, in Section 5.4, we discuss how these results may change when varying the measurement rate γ, as well as the range r and the partition size ℓ.

Measurement-only dynamics
We first neglect the Hamiltonian contribution, i.e., we implement the quantum-jump protocol described in Appendix B with Ĥ = 0.The Hamiltonian is only used to set the initial state as its corresponding ground state with h = 0.5.However, we have carefully checked that the choice of the initial state does not affect the asymptotic value reached by S ℓ (t), but only its transient dynamics. 7The advantage of studying a measurement-only dynamics is to isolate the effects of measuring the string operators, without considering the interplay with the unitary dynamics.
The three upper panels of Figure 1 report the ensemble-averaged entanglement entropy S ℓ (t) versus the rescaled time γt, for three ranges of the string measurement operator: r = 1 [panel (a)], r = 4 [panel (b)], and r = L/2 [panel (c)], and for different system sizes L. We note that, for any finite value of r (i.e., which does not depend on L), the entanglement entropy displays metastable plateaus, whose lifetime generally increases with some power of L (see the discussion below, following the rescaled data of Fig. 3).This may cause some difficulties in estimating the long-time asymptotic behavior, especially for large system size.However we checked that the qualitative scaling of the entanglement entropy with L does not change if one considers the value reached in one of such steady values.Moreover, since their lifetime grows with the size, in the thermodynamic limit we can treat the first emerging plateau as the effective stationary state.
Going into the details, for r = 1 we observe that the measurement-only dynamics generally destroys entanglement in time and the system stabilizes on a low-entangled asymptotic state [Fig.1(a)].This is confirmed by the data reported in Fig. 1(d), which shows the behavior of the asymptotic entanglement entropy S ℓ with the system size, as extrapolated by the second plateau stabilizing around the rescaled time γt ≳ 10 −2 L 2 [see also Fig. 3(a)].A decreasing behavior in L is clearly visible: S ℓ should saturate to a constant value for small enough ratio r/L, i.e., in the limit of local measurements, thus obeying an area-law behavior, as expected by a comparison with the case where onsite operators are measured [49,50].
When the correlation range of the dissipative string operator is increased, as for r = 4 [Fig.1(b)], the measurements still produce a low entangled state.In fact, the asymptotic value of the entanglement entropy approaches a constant value with L, with superimposed staggered oscillations, probably due to finite-size commensurability effects between r, ℓ, L [Fig.1(e)].It is worth noticing that, also in this case, we extrapolated the asymptotic value S ℓ over the first observable metastable plateau itself.When going to even later times, we find that the subsequent plateau, corresponding to a lower value of the entropy, is still obeying an area-law scaling (at least for L ≤ 48, so that we were able to reach such plateau).
The above scenario dramatically changes when r is chosen in such a way to be a thermodynamic fraction of L [see Fig. 1(c), for r = L/2].In that case, the measurement-only dynamics is projecting over Bell pairs at arbitrary distances, thus quickly leading the bipartite entanglement to an asymptotic stationary value which undergoes a volume-law behavior.This is clarified by Fig. 1(f), which shows S ℓ /L versus L −1 : the saturation of S ℓ /L to a finite value for L −1 → 0 (see the dashed line, which extrapolates to the asymptotic value in the thermodynamic limit) clearly evidences a linear growth S ℓ ∝ L for large system sizes.

Interplay with the Hamiltonian dynamics
We now switch on the Kitaev Hamiltonian Ĥ and follow the full quantum-jump dynamics for h = 0.1, starting from the ground state with h i = 0.5.As for the case with Ĥ = 0, we have checked that the asymptotic values reached by the entanglement entropy are independent of h i , while they may depend on h (although the qualitative behaviors with L are unaffected by this choice).
Figure 2: Same plots of Figure 1, but in the simultaneous presence of the Kitaev Hamiltonian Ĥ and of the string measurement process.We observe qualitatively similar behaviors, although with different time scales.Here we fix γ = 0.5, J = 1, h = 0.1, and start from the ground state of Ĥ with h i = 0.5.
].The data shown have been obtained for γ = 0.5.In Fig. 2(a) we focus on the local case r = 1.The long-time averaged entanglement entropy attains a stationary value that, as for the measurement-only case, slightly decreases with the system size, thus reflecting an area-law behavior [see Fig. 2(d), which shows S ℓ corresponding to the plateau reached after a transient time γt ≳ 10 −2 L 2 , as visible in Fig. 3(c)].The value of the asymptotic entanglement entropy is slightly larger than the one obtained for the measurement-only dynamics, without the entangling action of Ĥ [compare with Fig. 1(a)].
Figure 2(b) reports the case r = 4, where we observe convergence in time to a single plateau, up to the times we were able to achieve numerically.Comparing the time behavior of the various curves with the corresponding ones for Ĥ = 0 [see Fig. 1(b)], it is reasonable to associate this only plateau with the late-time second plateau emerging with Ĥ = 0.In the presence of the Hamiltonian dynamics, the transient time appears to be reduced, with respect to that for the measurement-only dynamics.The asymptotic data shown in Fig. 2(e) (x-axis in logarithmic scale) suggest that S ℓ may grow sub-extensively with L, for small system sizes, while we expect it to eventually attain an area-law regime for large enough L. Physically, this can be understood by the fact that, for L much larger than r, the dissipation becomes effectively local and the entanglement growth is thus eventually suppressed.
In contrast with that and analogously as for the case without the Hamiltonian, if the range r of the jump operator is comparable with the total size L, the measurements can stabilize a volume-law scaling of the entanglement, as we can see in Fig. 2(c), referring to r = L/2.Here S ℓ quickly attains, in time, an asymptotic value that follows a volume-law scaling with L [see the data in Fig. 2(f) for S ℓ /L versus L −1 , which display convergence to a finite value in the limit L −1 → 0, and compare them with Fig. 1(f)].
It is finally worth commenting on the metastable plateaus displayed by the entanglement entropy, when considering ranges not extensively scaling with the system size.The curves reported in Fig. 3 are for r = 1 and correspond to those of Fig. 1(a) and Fig. 2(a), without and with the Hamiltonian, respectively.Here we changed the time axis according to t → γt/L 2 and noticed a fairly good convergence to an asymptotic scaling behavior with L. In fact, as highlighted in the two insets, the time t ⋆ at which the entanglement entropy S ℓ (t ⋆ ) reaches a given threshold value scales as t * ∼ L 2 .
Figure 3: Same plots as in Fig. 1(a) and Fig. 2(a) for r = 1 (i.e., without and with the Hamiltonian Ĥ, respectively), such that the time (x axis) has been rescaled as t → γt/L 2 .The two insets show the time t ⋆ at which the entanglement entropy S ℓ (t ⋆ ) reaches the threshold value 0.75 (dashed line in the two main frames), as a function of L 2 .
A less clear situation emerges at larger values of r.For example, when r = 4, with the numerical data at our disposal we could not find a clear signature of rescaling with L. Nonetheless, it emerges that the metastable plateau for Ĥ = 0 has a lifetime that grows superlinearly with L.

Varying γ and r
In the numerical analysis presented so far, the measurement rate has been always kept fixed.It is however quite natural that the parameter γ should have an impact on the stabilization of the entanglement, when the measurement-induced dynamics competes with the (unitary) Hamiltonian dynamics.In fact, for fixed and finite values of r, we observe a change of behavior in the entanglement scaling with L, from a logarithmic growth to an area-law scaling, when varying γ and in the presence of Ĥ.An example of this is reported in Fig. 4(a) where we show, for r = 1, the asymptotic entanglement entropy S ℓ versus the system size L for different measurement rates.In this case, a clear distinction between the two above mentioned regimes is observable for γ c ≈ 0.15.This behavior agrees with what has been already observed for a dynamics with onsite measurement operators [49,50].To support this result in Fig. 4(b) we plot the rescaled entanglement entropy S/ log(L) vs γ for different L. As expected, within the error bars, the curves collapse for γ ≲ 0.15, indicating a logarithmic scaling.For larger γ, instead, the various curves drop, suggesting that the asymptotic entanglement entropy grows slower than log(L).On the other hand, for r growing extensively with L, we expect that the volume-law scaling of the asymptotic entanglement entropy is robust.
Coming to the dependence of the entanglement dynamics on the range r of the string measurement operators, differently from the onsite case, this should crucially depend on the system size L, since, as discussed above, the entanglement scaling is affected by the ratio L/r.Unfortunately, the numerical effort required to derive the results does not allow us to reach quantitative conclusions, because of the strong finite-size effects we have to deal with.However, our numerics suggests the existence of a critical value of n = L/r (n > 1), separating a regime in which the asymptotic entanglement entropy does not change with the system size (i.e., n ≫ 1), to a regime where it exhibits a linear growth with L (i.e., n ∼ 1).
Figure 4: (a): Asymptotic entanglement entropy S ℓ versus the system size L, for r = 1 and different measurement rates (color scale).Here we fix J = 1 and quench the field from h i = 0.5 to h = 0.1, as in Fig. 2. For γ ≲ 0.15 we observe a logarithmic scaling, while for γ ≳ 0.15 an area-law behavior sets in.Note the logarithmic scale on the x axis.(b): Rescaled asymptotic entanglement S/ log(L) vs γ.For γ ≲ 0.15 the rescaled curves (are overlapped within the error bars), while they separate for larger measurement rates.
From our simulations, we cannot rule out the possibility that, in the thermodynamic limit, a volume-law scaling emerges for any non-diverging value of n.A qualitative argument supporting this possibility is the following.Assuming that two subsequent measurements with overlapping operators of range r (which denotes the length of the corresponding string in spin-1/2 language) occur, these may strongly influence each other and thus generate large entanglement.Roughly speaking, the ratio between the probability p ov that the operators in subsequent measurements overlap and the probability p nov that they do not is p ov /p nov ∼ 2r/(L − 2r).In the thermodynamic limit, this ratio remains nonvanishing only when r = O(L) (i.e., n ∼ 1).Thus it is reasonable to assume a volume-law scaling of the entanglement entropy whenever L/r is not diverging.
We finally point out that, although the numerical analysis described above has been performed for a partition size fixed to ℓ = L/4, we have done further checks to test the impact of this choice on our results.We found that the general scenario is qualitatively robust, at least for small and for large ranges r of the string measurement operator, compared with ℓ.Different behaviors, however, may appear when r becomes of the order of the partition size ℓ, suggesting that the important ratio to distinguish between the different entanglement scalings is r/ℓ, rather than r/L.

Conclusions
We have unveiled some aspects of the entanglement-entropy dynamics of fermionic many-body Gaussian states, in the presence of a Gaussian-preserving evolution, with a special focus on the possible emergence of measurement-induced phase transitions in such systems.The first part of the paper presents a detailed technical discussion on how to maintain and treat Gaussianity, when exponentials of quadratic fermion operators are applied to a Gaussian state.Our purpose is to provide the reader with a comprehensive guide to implement fermionic Gaussian evolutions and study the dynamics (not necessarily unitary) of such systems.Among all the possible treatable evolutions, we focused on the quantum-jump unraveling of a Lindblad dynamics induced by a particular class of quadratic fermion operators, having a simple definition in the fermionic language, but a non-trivial string structure in terms of the Pauli matrices.We called them string operators.
In the second part of the paper, after having derived the equation of motions for a quantumjump dynamics induced by the string operators, we have solved it for the quantum Ising chain, focusing on the time behavior of the entanglement entropy.Our main result is that the scaling of the entanglement entropy with the system size L strongly depends on the range r of the string operators.Because of numerical limitations we cannot investigate system sizes large enough to claim that this is an entanglement transition rather than a simple crossover.For short strings (r finite), we see that the asymptotic entanglement entropy obeys a non-extensive scaling behavior.In this case, the measurement operator is effectively local in space and thus it asymptotically destroys quantum correlations that may be generated by the unitary dynamics.For string ranges comparable with the system size (r ∼ L), the measurement process may become highly entangling and is able to stabilize a volume-law scaling in the asymptotic entanglement entropy.
Remarkably, the picture above is valid both for a measurement-only dynamics and for a dynamics where also the effect of a short-range Hamiltonian is considered.So, the measurement dynamics of the string operators is strong enough to induce different entanglement scalings by itself.
It would be tempting to address dynamical behaviors induced by more than a single type of string measurement operators (e.g., with different r), in such a way to explore the possibility of having frustration of the measurements [42,74].A further natural extension of our work could be to understand how the MIPT scenario would change in the presence of long-range measurement operators (e.g., by considering some power-law decay of dissipators with the distance).Finally, novel dynamical mechanisms may emerge in systems beyond the free-fermion paradigm, due to the role of interactions in the entanglement dynamics under variable-range monitoring.
In this form, the Lindblad dynamics ruled by Eq. (B.1) can be interpreted as a deterministic non-Hermitian time evolution of a pure quantum state |ψ t 〉 generated by Ĥeff , plus a stochastic part given by the possibility of applying mj to the state |ψ t 〉 [cf. the second term in the righthand-side of Eq. (B.3)].For mj Hermitian, this dynamics can be thought of as an occasional, yet abrupt, measurement process in which, at any time, the system has a chance to be POVMmeasured [7,8], i.e., to undergo the action of one of the operators mj , taken from a distribution (see more details below).
In this framework, Eq. (B.1) can be solved by implementing the following stochastic process.The time evolution is discretized in time slices of step d t.Then, at each step: 1.With probability π j = γ 〈 m † j mj 〉 d t, the jump operator is applied: Note that the Lindblad master equation in Eq. (B.1) is invariant under a shift of the jump operators mj → 1 + mj .In practice, one can implement the protocol described above with the Lindblad operators 1 + mj (as we actually do in Sec.5), instead of mj , and obtain the same averaged density matrix, although the single trajectories |ψ t 〉 may behave differently.We remark that, since the von Neumann entropy of Eq. ( 54) quantifies the amount of genuine quantum correlations of a pure state, we first evaluate it over the single pure-state quantum trajectory and then perform the ensemble averaging.It is important to fix the order of the two operations of measuring and ensemble averaging, since they are commuting only for linear functions of the state ρ (as for expectation values of observables), but not for quantities as the one in Eq. (54).For a more detailed discussion see, e.g., Refs.[44,49,50].

C Evolution of the state: A different derivation
In this appendix we provide a different derivation of the analytic expression for the evolved state |ψ〉 MI,J , by evaluating the operator Ẑ(1) = 1 2 p,q Z p,q ĉ † p (1) ĉ † q (1), in Eq. ( 16).To this purpose we notice that ÂI,J , Ẑ = 1 2 p,q i∈I Z p,J δ q,i + Z q,J δ p,i ĉ † p ĉ † q , (C.1) where we have defined Ẑ = 1 2 p,q Z p,q ĉ † p ĉ † q and Z p,J = j∈J Z p, j .Consequently Ẑ, [ MI,J , Ẑ] = 0.This allows to write MI,J e Ẑ = e Ẑ MI,J + MI,J , Ẑ , (C. being Ẑ′ = 1 2 i∈I p,q Z p,q + Z p,J δ q,i + Z q,J δ p,i ĉ † p ĉ † q . (C.5) For the local case M l I,J in Subsec 4.1, the jumped state |ψ〉 M l I,J is described by the antisymmetric matrix Z l ′ , with Z l ′ p,q = (1 + δ q,i + δ p,i )Z p,q .(C.6) On the other hand, for the non-local case M nl I,J in Subsec.4.2, the antisymmetric matrix Z nl ′ reads Z nl ′ p,q = Z p,q + (δ p,i + δ p,i+r )(Z i,q + Z i+r,q ) .(C.7)

Figure 2 displays
Figure 2 displays S ℓ (t) versus γt (upper panels) and the asymptotic value S ℓ versus L, for the three ranges of r = 1 [panels (a),(d)], r = 4 [panels (b),(e)], and r = L/2 [panels (c),(f)].The data shown have been obtained for γ = 0.5.In Fig.2(a) we focus on the local case r = 1.The long-time averaged entanglement entropy attains a stationary value that, as for the measurement-only case, slightly decreases with the system size, thus reflecting an area-law behavior [see Fig.2(d), which shows S ℓ corresponding to the plateau reached after a transient time γt ≳ 10 −2 L 2 , as visible in Fig.3(c)].The value of the asymptotic entanglement entropy is slightly larger than the one obtained for the measurement-only dynamics, without the entangling action of Ĥ [compare with Fig.1(a)].Figure2(b) reports the case r = 4, where we observe convergence in time to a single plateau, up to the times we were able to achieve numerically.Comparing the time behavior of the various curves with the corresponding ones for Ĥ = 0 [see Fig.1(b)], it is reasonable to associate this only plateau with the late-time second plateau emerging with Ĥ = 0.In the presence of the Hamiltonian dynamics, the transient time appears to be reduced, with respect to that for the measurement-only dynamics.The asymptotic data shown in Fig.2(e) (x-axis in logarithmic scale) suggest that S ℓ may grow sub-extensively with L, for small system sizes, while we expect it to eventually attain an area-law regime for large enough L. Physically, this can be understood by the fact that, for L much larger than r, the dissipation becomes effectively local and the entanglement growth is thus eventually suppressed.In contrast with that and analogously as for the case without the Hamiltonian, if the range r of the jump operator is comparable with the total size L, the measurements can stabilize a volume-law scaling of the entanglement, as we can see in Fig.2(c), referring to r = L/2.Here S ℓ quickly attains, in time, an asymptotic value that follows a volume-law scaling with L [see the data in Fig.2(f) for S ℓ /L versus L −1 , which display convergence to a finite value in the limit L −1 → 0, and compare them with Fig.1(f)].It is finally worth commenting on the metastable plateaus displayed by the entanglement entropy, when considering ranges not extensively scaling with the system size.The curves reported in Fig.3are for r = 1 and correspond to those of Fig.1(a) and Fig.2(a), without and with the Hamiltonian, respectively.Here we changed the time axis according to t → γt/L 2 and noticed a fairly good convergence to an asymptotic scaling behavior with L. In fact, as highlighted in the two insets, the time t ⋆ at which the entanglement entropy S ℓ (t ⋆ ) reaches a given threshold value scales as t * ∼ L 2 .