Purification Dynamics in a Continuous-time Hybrid Quantum Circuit Model

We introduce a continuous time model of many-body quantum dynamics based on infinitesimal random unitary operations, combined with projective measurements. We consider purification dynamics in this model, where the system is initialized in a mixed state, which then purifies over time as a result of the measurements. By mapping our model to a family of effective 1D quantum Hamiltonians, we are able to derive analytic expressions that capture how the entropy of the system decays in time. Our results confirm the existence of two distinct dynamical phases, where purification occurs over a timescale that is exponential vs. constant in system size. We compare our analytic expressions for this microscopic model to results derived from field theories that are expected to capture such measurement-induced phase transitions, and find quantitative agreement between the two.


I. INTRODUCTION
The continuing development of programmable quantum devices with increasing numbers of degrees of freedom has led to a great deal of interest in addressing fundamental questions regarding the dynamics of information in many-body quantum systems [1][2][3][4][5][6][7][8][9][10][11][12].In recent years, there has been a particular focus on the competition between unitary operations, which generate entanglement, and local projective measurements, which are non-unitary processes that break entanglement.Models of dynamics that feature both of these ingredients are often referred to as hybrid quantum circuits, the study of which has led to the discovery of a sharp entanglement phase transition driven by the rate of measurements, separating regimes where many-body entanglement is either stable or fragile against these measurements [13][14][15][16][17][18][19][20][21][22][23][24].Typically, the studied geometry is that of a 1D chain of qudits, but similar transitions have also been found in more complex geometries such as random tensor networks [25][26][27][28].
The existence of this transition was first understood in terms of the entanglement structure of an ensemble of pure many-body states at equilibrium.Subsequent studies also revealed the existence of a simultaneous dynamical phase transition, which can be understood as the ability of the measurement protocol to learn an initially mixed state [29][30][31][32].The latter suggests a connection between the dynamics of hybrid quantum circuits and quantum error correcting codes [33], which by construction protect information against deleterious non-unitary processes.The transition was also shown to play an important role in the context of simulating the behaviour of open quantum systems [34][35][36][37].
These considerations have led to the notion of purification dynamics, where one studies how the entropy of an initially mixed state decreases over time as a result of the measurements.Away from the critical measurement rate we find two phases where the state purifies over a timescale that increases exponentially with system size ('mixed phase') or is independent of the system size ('purifying phase').To understand the phenomenology of these phases in a fully quantitative way, arguments based on capillary wave theory have been put forward [38].Using an effective field theory which is expected to capture the universal features of the transition, one can obtain concrete predictions of how the purity of the system will depend on time in each phase.However, direct verification of these predictions by means of a direct calculation from a microscopic model are as of yet lacking.
In this paper, we introduce and study a hybrid quantum circuit model of dynamics that is defined in continuous time, the properties of which we are able to calculate analytically.In particular, by means of a mapping onto an effective Hamiltonian, we are able to compute the time dependence of a particular family of operator-space entanglement measurements, which can be related to the purity of the system at a time t, starting from a maximally mixed initial state.We look in detail at both the mixed and purifying phases, as well as at the transition between them.Our results agree with those of capillary wave theory in both phases: the entropy decays exponentially with time in the purifying phase, and decreases as − log t in the mixed phase over an exponentially long time window [Eqs. (64,71)].In our calculations, we consider both periodic and open boundary conditions, and show that the two choices give rise to quantitatively different behaviour when in the mixed phase: in particular, a (1/2) log N contribution to the entropy appears when we impose periodic boundary conditions, but this is absent for open boundary conditions.We also look at the dynamics at criticality, where there exists a regime during which the entropy decays algebraically, Eq. (77).
The structure of our paper is as follows.In Section II we introduce a continuous time model of dynamics based on infinitesimal random unitary operations, and describe how one can calculate various measures of entanglement and information spreading in this model.We supplement the unitary dynamics with projective measurements in Section III, and in Section IV, we explain how the re-sultant unitary-projective dynamics can be mapped onto imaginary-time evolution under an effective 1D Hamiltonian.We then present our main quantitative results in Section V, giving analytic expressions that quantify how the purity of the system increases as a function of time in the purifying/mixed phase and at criticality.Finally, we discuss our results and conclude in Section VI.

II. CONTINUOUS-TIME RANDOM CIRCUIT MODEL
In this section, we introduce a random unitary circuit (RUC) model of unitary dynamics, and describe how its entanglement properties can be analysed.We will later incorporate measurements into this model, which will allow us to study the dynamics of purification.
We consider a one-dimensional array of N qudits, each with a local Hilbert space of dimension d.The evolution is driven by a spatially local unitary circuit with a brickwork structure, illustrated in Fig. 1.In a given timestep τ = 1, 2, . .., two-site unitaries are applied to each pair of qudits on the odd bonds (2j − 1, 2j), followed by another layer of unitaries on the even bonds (2j, 2j + 1).These elementary two-site unitaries each have the same structure, also depicted in Fig. 1.First, single-site gates U ⊗V are applied, followed by evolution under some two-qudit Hamiltonian H for a time ∆t, and finally the change of basis is undone by applying the inverse single-site rotations U † ⊗V † .We denote the unitary operator describing the evolution from time 0 to τ as W (τ ).
Throughout this work, H will be treated as a free parameter of the model and it is kept fixed across both time and space.To simplify calculations, we will assume it is real, hermitian and symmetric under swapping the two qudits it acts on.The single-qudit unitaries will be sampled randomly and independently from the Haar ensemble for each unit cell.We will generally be interested in the limit where ∆t → 0, which we refer to as the continuous-time limit.Note that the state only evolves by an infinitesimal amount in each timestep, in contrast to discrete-time RUC models of quantum dynamics (e.g.Refs.[6,39]).A model of continuous-time dynamics was studied numerically in [40].Our method for constructing the unit cell is more general and more easily amenable to analytical treatment.
Our focus will be on the dynamics of entanglement and quantum information in these continuous-time models.For this purpose, it is useful to consider the Choi-Jamiolkowski (CJ) state |W (τ )⟩ corresponding to the unitary W (τ ).This state is defined on two copies of the system, which we can associate with the inputs and outputs of the time evolution operator.Formally, we have a=1 |a⟩ ⊗ |a⟩) consists of maximally entangled states between each input qudit and its corresponding output [41].Many important quantities that are used to diagnose the spreading of quantum information can be expressed as simple functions of this operator-state |W (τ )⟩ [5], and in particular we will find this representation useful when it comes to studying purification dynamics.
As is now common in studies of RUC dynamics, we use the Rényi entropies to quantify the entanglement properties of the state |W (τ )⟩ where n is some positive parameter.Here, ρ A is the reduced density matrix of |W (τ )⟩ corresponding to some subset A of inputs and outputs.Compared to the usual von Neumann entropy S vN , the Rényi entropies for n = 2, 3, . . .are more amenable to analytic studies, since they only involve integer moments of the density matrix and hence can be computed using a replica method.The von Neumann entanglement entropy S vN can be recovered by constructing an analytical continuation of the function and taking the limit n → 1 (see, e.g.Ref. [42]).
For the largest part of this work, we will only be concerned with the second Rényi entropy S (2) (ρ A ), which is the simplest to evaluate.This is a lower bound on the von Neumann entropy S vN ≥ S (2) , which in certain cases is known to be asymptotically tight [43].Since the purity Tr (ρ A ) 2 ≡ exp −S (2) it can be expressed using a fourfold copy of the evolution operator, which we denote Note here that the operator replicated in the expression differs from |W (τ )⟩ ⟨W (τ )| by a reshuffling of the indices.Henceforth, we will use this convention, but retain the essence of the CJ isomorphism by noting that we treat inputs and outputs on par when discussing Rényi entropies.
Define (unnormalized) states which live in the fourfold-replicated Hilbert space of each physical site j.In terms of these, we have where we denote the set of input (output) sites included in the region A as A in (A out ), and the states and similar for |Ψ Aout ⟩.To make progress, we look at the average of the Rényi entropy (5) over the random ensemble of unitary circuits.More precisely, we will evaluate the average purity as opposed to the average entropy, which is equivalent to performing averages inside the logarithm of Eq. ( 1).This simplification is common in analyses of RUCs [17], and still recovers the correctly-averaged von Neumann entropy if one takes the replica limit n → 1. Accordingly, averaging the purity amounts to replacing W (2) (τ ) with its ensemble average W (2) (τ ).As shown in App.A, W (2) (τ ) maps states spanned by tensor products of |I⟩ j , |S⟩ j to other such states, meaning we can focus on the restriction of this averaged operator to the subspace V (S 2 ) ⊗N , where V (S 2 ) = span(|I⟩ , |S⟩).
Because the single-site Haar-random unitaries appearing in each of the two-site elementary blocks of the circuit [Fig.1(b)] are sampled independently, we can consider the ensemble average of the evolution under a single one of these blocks, which we denote T : V (S 2 ) ⊗2 → V (S 2 ) ⊗2 .Using the Weingarten diagrammatic calculus as seen in App.A, we find that, for small ∆t, this map can be expressed as where i, j label the sites on which the unit cell acts, Ω(H) is a measure of the entangling power of the Hamiltonian and Wg is the Weingarten matrix corresponding to the symmetric group The induced evolution can be equivalently described using the effective imaginary-time Hamiltonian in terms of which the unit cell map is It is interesting to note that there are no contributions from odd powers of ∆t in the expansion of Eq. 7. Looking at the form of the Hamiltonian (10), we see that in the effective Hilbert space spanned by the states (3,4), the only mobile degrees of freedom are domain walls separating regions of |I⟩ from |S⟩, consistent with discrete-time RUCs discussed previously [6,39].The initial Hamiltonian H only enters the expression through its entangling rate Ω(H).This sets the overall timescale of quantum information transfer through the system.In App.B, we compute the transfer matrix for a higher number of replicas and show that this statement holds more generally.This result suggests that the qualitative behavior of entanglement dynamics derived from our model should be insensitive to most of the microscopic details, and hence applicable to a wide range of physical processes.
The propagator for the whole circuit T can be constructed by concatenating the two-site maps (11) according to the brickwork circuit structure illustrated in Fig. 1.We have 4 ).(12) We now define the effective time as t = τ ∆t 2 and take the limit τ → ∞, ∆t → 0 such that t is kept constant.Using the Suzuki-Trotter formula, we find the limit of the previous equation We reiterate here that this operator acts as the restriction of W (2) (τ ) to its invariant subspace V (S 2 ) ⊗N , and therefore may replace it in average entropy calculations (e.g.averaging Eq. 5).
In its current form, the effective Hamiltonian is not Hermitian, but can be made so through a local similarity transformation, a technique commonly encountered in the study of non-equilibrium dynamics [44].
If we define a new evolution operator by T = (Wg − 1 2 ) ⊗N T (Wg 2 ) ⊗2 , which gives us the Hermitian interaction where the overall strength is given by γ = 2Ω(H)/(d 2 − 1) 2 .This type of local interaction is found in the literature both as the quantum equivalent of the classical two-dimensional axial next-nearest neighbor Ising model (ANNNI) [45,46] or more recently as the Jordan-Wigner transform of the balanced interacting Kitaev chain [47,48].
In the limit of d → ∞ we are left with a simple ferromagnetic nearest-neighbor Hamiltonian, with each domain wall incurring an energy penalty of γ.The Hamiltonian is symmetric under the global spin-flip operator x , as can be seen through the commutation relation [C, H ij ] = 0.

III. INCLUDING MEASUREMENTS
In this section, we will introduce the formalism that can be used to incorporate measurements into the random circuit evolution.For the purpose of this work, we will consider projective measurements in the computational basis of each qudit that occur stochastically.The same framework can accommodate for weakmeasurement schemes as seen in Ref. [17].Due to the continuous nature of our circuits, the effective model will be identical in the two cases.
A projective measurement is a non-unitary stochastic process, where the wavefunction of the system |ψ⟩ collapses to a post-measurement state |m⟩ with probability Here, the set of wavefunctions {|m⟩} is the computational basis in which the measurement is performed and m = 1, 2, . . .d.For any fixed realisation of the random unitary circuit and positioning of the measurements, the final state of the system will depend on all the measurement outcomes m = (m 1 , m 2 , . ..).Thus, we can write the ensemble of final states as {(p m , |W m )⟩}, where p m is the joint probability of the measurement results, and |W m )⟩ is the (normalized) conditional state.As before, we will imagine the Choi-Jamiolkowski state, so |W m )⟩ is a state on a twofold copy of the system, and is constructed by preparing a maximally entangled state between the two copies in the computational basis, and evolving one of the copies under the evolution in question.
As is typical in the study of hybrid quantum circuits, our interest is on the statistics of the entanglement properties of individual conditional wavefunctions |W m )⟩; see, e.g.Refs.[14,29].The natural quantity to consider for this purpose is the von Neumann entropy, S vN (ρ A m ), where ρ A m is the reduced density matrix of |W m ⟩ over a subset of inputs and outputs A. Specifically, we would want to compute the average of this quantity over all realizations of the random circuit, measurement locations and measurement results, which we denote S vN (ρ A m ).However, this quantity is very difficult to compute directly in random circuit models.We follow Ref. [17] and introduce the series of related quantities

S(n)
where M labels a particular configuration of measurement locations in spacetime, which occurs with probability p M , and m runs over all measurement results for the given configuration M .These quantities are related to measurement-averaged Rényi entropies, with the main difference that each outcome is weighted by p n i M .The additional factor of d |M |(n−1) ensures that the correct order of magnitude, in powers of d, of the correct weight is preserved, and only deviations from it are amplified by the number of replicas.The renormalization is also performed on average, i.e. we compute the average of the numerator and the denominator independently.Knowledge of this quantity for all integers n ≥ 2 can be used to recover the average entanglement entropy SA of subsystem A using the replica limit Each term in the sums over M appearing in the numerator and the denominator in Eq. ( 15) is a scalar that depends linearly on the tensor which is analogous to the duplicated state in Eq. (2) defined earlier.The angled brackets are a short-hand notation for the weighted sum on the RHS.For simplicity, we once again revert to the normal operator indices, but keep in mind that the probabilities are obtained from expectation values of the appropriate projectors in the CJ state |W m ⟩.
To see how this tensor evolves as the circuit progresses, let us consider how ⟨W (n) ⟩ is updated when a new measurement is performed on site i, the outcome of which we denote m.Each of the n replicas transforms via the action of a projector P (i) m , which corresponds to the mth computational basis state for qudit i.Since all measurement outcomes m are summed over in Eq. ( 17), we find The effect of post-selection is included by assuming a perfect correlation of the measurement results in all n replicas.
Since adding an infinitesimal time evolution to the averaged tensor only leads to linear transformations by left multiplication due to both the chaotic dynamics and the measurements, we can proceed again by mapping the evolution of ⟨W (n) ⟩ to a reduced system with an effective Hamiltonian.If we focus again on the twofold replica n = 2, we see that the action of the measurement opera-tor on the reduced Hilbert space at each site is Therefore, we find that the new vector space V M (S 2 ) = span(|I⟩ , |S⟩ , |O⟩) is closed under measurements.If we promote this to the reduced Hilbert space of the entire chain V ⊗N M (S 2 ), we find that this is also closed under the action of the Haar averaged unit cell between any pair of sites.To show this, we can consider the properties of the following linear combination It is straightforward to show that this becomes null under any contraction between a normal and a complex conjugate leg.Due to the rules of the Weingarten calculus, this means that such local states are preserved by averaged unit cells.This is summarised in the following equation In Appendix C we give an explicit representation of the new operator T M ij , acting on V M (S 2 ) ⊗2 .We find that evolution in subspaces that contain |X⟩ states happen at a different rate Γ, independent of the rate of information propagation γ.This is defined by and can be qualitatively understood as an energy cost associated with |X⟩ states.In App.D we derive a more explicit relation between the rates Γ, γ and the microscopic Hamiltonian H.The new state |X⟩, which appears after a measurement, ensures that we obtain the correct correlations between measurements performed consecutively at short time intervals on the same qudit.The timescale 1/Γ represents the time it takes a qudit to relax before we can obtain new information by measuring again in the same basis.For the rest of this work, we set Γ → ∞, such that no measurement inertia can be observed.In App.C, we show that doing so is effectively equivalent to projecting out the |X⟩ state and working in the previous 2-dimensional reduced Hilbert space V (S 2 ).The action of the measurements is also projected onto this subspace and can be expressed as It can be shown that this same operator is obtained in the reduced subspace if we consider instead measurements in random bases.
In the following, the measurements are distributed through the circuit according to an independent Poisson process for each site, at some uniform rate f (in the natural time units of the continuous model).The transfer matrix at time t under both random dynamics and measurements is then given by an effective imaginary-time evolution T eff (t) = exp(−tH eff ), with H eff given by From Eq. 5 and Eq. 15 we see that we can express the second moment of the entanglement entropy of some subregion A at time t using matrix elements of the transfer matrix The denominator acts as a normalization factor, so using the expression above allows us to safely drop constant terms in the effective Hamiltonian.
We can perform a similar analysis for the case of multiple replicas.Using the results in App.B and the limits d, Γ → ∞ we show that the effective Hamiltonian of the n'th replica theory is given by where M (n) is the generalization of the operator in Eq. ( 23) that acts as and D is a diagonal two-site operator with entries given by with D(σ, τ ) the bi-invariant metric on S n given by the Hamming distance between σ and τ , i.e. the number of elements that are not mapped onto themselves under τ −1 σ.This form is manifestly consistent with the expected symmetry group S n × S n .It is interesting to note that the d → ∞ limit does not result in the fine tuned S n! -symmetric Potts model observed in circuits with fully Haar random unit cells [17].

IV. FERMIONIC MAPPING
If we take the limit of large local dimension d and keep only the leading contributions, we obtain dynamics driven by where g = 2f /γ.This is easily recognized as the transverse field Ising model (TFIM) in 1D, subject to open boundary conditions.It is well-known that this can be mapped to a system of non-interacting fermions using the Jordan-Wigner transformation [49].In this section, we will introduce the general formalism used to compute quantities of the form shown in Eq. 25.
We start by constructing a set of non-local Majorana operators as defined on all sites i = 1, 2 . . .N [50].These operators are Hermitian (γ µ ) † = γ µ and obey the standard anticommutation relations where the indices µ, ν are understood to run over all 2N previously defined operators.From the definition, we get the additional relation such that the product of all Clifford operators is This operator anti-commutes with all the Majorana fermions {C, γ µ } = 0 and it is a conserved quantity, since it commutes with the full Hamiltonian [C, H] = 0. We can couple Majorana fermions living on adjacent sites into domain wall creation and annihilation operators for i = 0, 1 . . .N − 1, where we assume periodic boundary condition N = 0.These obey the typical anticommutation relations A simple calculation shows that such that the number operator of the fermionic mode at some site i ̸ = 0 is a projector onto configurations that have a domain wall between sites i and i + 1.With this convention, the Hamiltonian becomes a quadratic form This can be more succinctly expressed using the Bogoliubov-de Gennes notation where a = (a 0 , a 1 , . . ., a N −1 , a † 0 , a † 1 , . . ., a † N −1 ) T .The matrix D is called the grand-dynamical matrix and it obeys the particle-hole symmetry equation .
The exponentials of such Hamiltonians are most easily treated using the algebra of fermionic Gaussian states, as worked out in Ref. [51].We will briefly outline some of the results relevant to our calculation.It is convenient to define a Gaussian state through its generating quadratic form as with a normalisation constant Z(W ) chosen such that Tr ρ[W ] = 1.By Wick's theorem, such many-body states are fully characterized by their two-body correlation matrix, defined as The correlation matrix is related to the generator of the quadratic form through the useful relations where it is assumed that 1 − Γ is invertible.Using a special case of the Baker-Campbell-Hausdorff formula, it is shown that fermionic Gaussian states are closed under multiplication and we have with the new generating matrix Ω given by If we denote the correlation matrix of Ω by Γ × Γ ′ , with Γ and Γ ′ the correlation matrices of W and W ′ respectively, the following formula is proven in Ref. [51] Γ Inner products can be easily computed using the following trace formula where the ambiguity of the sign is in general a complex issue, but this will not be a problem for our purposes.

V. DYNAMICS OF PURIFICATION
In the preceding sections, we developed a formalism that allows us to study entanglement dynamics in our continuous-time random quantum circuit models.Here, we focus specifically on purification dynamics in these models.Namely, starting from an initial mixed state, we are interested in how fast the state of the system is purified by measurements.We will be particularly interested in the purification transition that occurs as a function of measurement frequency f [29], which is thought to be concomitant with the measurement-induced entanglement transition separating area-and volume-law phases [13][14][15][16]18].Thanks to the exact solvability of our model in the d → ∞ limit, we are able to compute analytical expressions for the order parameters of this dynamical phase transition, and infer the key critical exponents.

A. Setup and phase diagram
The setup we study is as in Ref. [29]: The system is initialized in a maximally mixed state, which is represented in the above formalism by the input state |I⟩ ⊗N .After some evolution time t, the purity of the state of the system will have increased from its initial value due to the measurements.As explained previously, we will use the quantity (15) as a measure of the the typical entropy of the ensemble of states.Since we are looking at the purity of the entire state after a time t, the set A that appears in Eq. ( 25) will contain all of the output qubits.Accordingly, we can express the quantity in question in terms of the transfer matrix T eff (t) The purification transition that occurs in our model is associated with a quantum phase transition in the effective Hamiltonian H eff , which generates the time evolution operator T eff (t).Based on the phase diagram of the TFIM, we can deduce that such a transition must occur at the critical measurement rate g c = 1, i.e. f c = γ/2.In the spin basis (29), the two phases correspond to the for g > 1, and a spontaneous symmetry-broken phase for g < 1.
For the problem in hand, the relevant order parameter that we use to distinguish these phases is not a correlation function, as is usually the case, but rather the many-body overlap appearing inside the logarithm in Eq. (49).To provide intuition into how this quantity behaves either side of the transition, we can reformulate our expression for S(2) (t) as follows.Since |S⟩ ⊗N = C |I⟩ ⊗N , the above fraction becomes equal to the expectation value of C in the state where | Ψ(t)⟩ := |Ψ(t)⟩ / ⟨Ψ(t)|Ψ(t)⟩ is the wavefunction after imaginary time evolution under H eff /2, appropriately normalized.
If the measurement rate is sufficiently high such that the Hamiltonian ( 29) is in a symmetry-unbroken phase, then the ground state is non-degenerate and thus | Ψ(t)⟩ inherits the symmetry of the Hamiltonian.Since the Hamiltonian is also gapped, we see that the (accordingly normalized) state |Ψ(t)⟩ = exp(−tH eff /2) |I⟩ ⊗N converges to the ground state exponentially quickly.The ground state must be an eigenstate of C, whose eigenvalues are ±1, so we can then conclude that ⟨C⟩ Ψ(t) → 1 exponentially quickly as t → ∞, and hence S(2) → 0 at a rate independent of the system size, as expected in this regime.When g < 1 the symmetry is spontaneously broken.In this case, the ground eigenspace is doubly degenerate in the thermodynamic limit N → ∞, and the effect of the transfer matrix at long times is to project onto this subspace.The projected state may no longer be an eigenstate of C, so we can have a non-zero residual entropy.As we will see, this residual entropy is extensive, with a log(N ) correction [Eq.( 69)].
While this picture allows us to understand the transition at a qualitative level, to obtain an analytic expression for the residual entropy we will instead use the fermionic mapping detailed in the previous section.We will find it convenient to work with states of definite fermion parity, and hence we define the density matrices ρ ± = |±⟩ ⟨±|, where |±⟩ := (|I⟩ ⊗N ± |S⟩ ⊗N )/ √ 2, which are eigenstates of C. We can then write where the parameter Θ is defined by We note that Eqs.(51,52) are quite general, and could be applied even if we didn't take the d → ∞ limit.
The above considerations help us anticipate the existence of two distinct dynamical phases, consistent with previous work on purification dynamics in discrete time random circuit models, which we refer to as 'mixed' (g < 1) and 'purifying' phases (g > 1), following Ref.[29].In the following, we derive analytical expressions for the time dependence of Θ, which in turn determines the Rényi entropy S(2) (t) via Eq.( 51).We will use these expressions later to understand the nature of the two phases and the transition between them at a quantitative level.

B. Expressions for Θ(t)
While we have so far left the boundary conditions unspecified, in computing Θ(t) we will consider open and periodic boundary conditions separately in our calculations.The conventional timescale γ = 1 is employed throughout this chapter.

Periodic boundary conditions
We start by considering periodic boundary conditions, which can be realised by introducing additional random unitary gates that act between sites 1 and N in the original circuit model.In this case, a standard calculation shows that the Jordan-Wigner-transformed Hamiltonian (39) acquires an additional term which imposes either periodic or antiperiodic boundary conditions depending on the fermion parity sector one works in (see, e.g.Ref. [49]).Taking N to be even from hereon for simplicity, the even (odd) parity sector features antiperiodic (periodic) boundary conditions.These are sometimes referred to as Ramond and Neveu-Schwarz sectors, respectively.
Thanks to the translation invariance of the system, the single particle Hamiltonian D can be block diagonalized using momentum eigenstates, whose wavevector is quantized to k l = lπ/N , with l ∈ {1, . . ., N − 1}.In the even parity sector l must be odd to be compatible with the boundary conditions, and likewise vice-versa.Recognizing that the states |±⟩ are the ground states of the Hamiltonian in the g → 0 limit in each parity sector, we can express Θ as a ratio of products over k l modes, with even l in the numerator and odd l in the denominator.
In Appendix E, we show that where and the energy eigenvalues are given by The factor e t in Eq. ( 54) accounts for the modes k = 0, k = π, which are only present in the odd parity sector.

Open boundary conditions
We can also consider open boundary conditions (OBCs), where the Hamiltonian H eff no longer features the term connecting sites 1 and N .While the change of boundary condition makes little difference in the purifying phase, we will find that in the mixed phase, quantitative differences between OBCs and PBCs can be seen in the behaviour of the Rényi entropy.As such, will focus mainly on mixed phase g < 1, although of the following holds true throughout the phase diagram.
With OBCs, the JW-transformed does not contain a term that manifestly depends on the fermion parity sector.This leaves us with the problem of diagonalizing the single-particle matrix D, which is now the same in both parity sectors.Since the system is no longer translation invariant, we cannot treat momentum eigenmodes separately as we did before.Instead, we must explicitly compute the correlation matrices Γ[−tD] and Γ ± , and use the more general expression (53).To do so, we calculate the single-particle eigenstates, which form the columns of a real orthogonal eigenvector matrix O.In terms of these, we obtain a spectral decomposition of the grand dynamical matrix D = OΛO −1 .The eigenvalues Λ come in pairs due to the particle-hole symmetry (41), which we arrange as Λ = diag(λ 0 , λ 1 , λ 2 , . . .λ N −1 , −λ 0 , −λ 1 , −λ 2 , . . .− λ N −1 ) with λ i > 0 and in non-decreasing order.The corresponding pairs of eigenvectors are also related through |−λ i ⟩ = η |λ i ⟩, with η defined in Eq. (41).In terms of these eigenvalues, the correlation matrix Γ[−tD] becomes which can be substituted directly into Eq.( 53).In Appendix F, we show that the single particle eigenstates take the form of sinusoids, whose wavevectors k l can be found as the solutions to the equation lying in the interval [0, π) and labeled in increasing order.In terms of these wavevectors, the eigenvalues λ l themselves again follow the well-known dispersion for the TFIM, Eq. (56).In the mixed phase, there is also a single imaginary solution to the above, which we label k 0 = iK, corresponding to a Majorana edge mode localized at the two boundaries of the chain.The energy of this edge mode λ 0 can be shown to be exponentially small in the system size N , and this energy is associated with a long timescale λ −1 0 which we will show is responsible for the slow decay of purity in this phase.

C. Behaviour of the Rényi entropy
With the above expressions in place, we are now ready to study the behaviour of the Rényi entropy in the two phases, as well as the critical point which separates them.

Mixed phase-periodic boundary conditions
In the mixed phase, the Rényi entropy shows quantitative differences depending on the boundary conditions, and thus we will consider both PBCs and OBCs in the regime g < 1, starting with the former.Beginning with Eq. (54), we can use complex integration methods to transform the alternating product over even and odd k modes into an infinite product with x q 's found as the solutions of the equation in the interval (K, ∞), where K = − log g is the point in the complex plane where the dispersion function λ(iK) has a zero.This form makes it manifestly clear that Θ < 1.The details of this calculation are given in App.E. Except at criticality, for sufficiently large N we can have N K ≫ 1, which in turn implies N x q ≫ 1 for all q.In this case, by virtue of the approximation log tanh(y) ≈ −2e −2y for y ≫ 1, we can approximate log Θ by an integral where dq/dx is the density of solutions and can be found by differentiating Eq. (60).The result to highest order in powers of To relate this expression to the Rényi entropy (15), we focus on the regime where t scales no faster than polynomial in N , such that the small factor e −N K dominates, making − log Θ(t) itself small.Then, if we make the approximations 1 − Θ ≈ − log Θ and 1 + Θ ≈ 2, we can deduce the following form for the Rényi entropy, valid for a broad window of times e N K ≫ t ≫ 2/(1 − g 2 ) (restoring the original units of time) where the term o(1) represents terms that tend to zero in the limit of large t or N .We see that the entropy decreases very slowly in time in this window, which is a defining feature of the mixed phase.Finally, for times t that scale exponentially with system size γt ≳ e N K , such that | log Θ(t)| ≫ 1, we find that the entropy decays as

Mixed phase-open boundary conditions
We now wish to compute the same quantity with open boundary conditions in the mixed phase g < 1.In particular, we are interested in the regime during which the entropy decays very slowly.As such, we can separate out the bulk single-particle eigenstates, whose energies lie above the bulk gap ∆ = (1 − g), from the Majorana edge mode, which is exponentially small in N .In particular, as long as one is not too close to criticality (1−g) ≫ 1/N , the Majorana eigenvalue can be approximated as This indicates that there is a regime of times ∆ −1 ≪ t ≪ λ −1 0 during which the transient bulk modes have decayed away exp(−tλ i ) ≈ 0, while the Majorana mode has not decayed.The approximate correlation matrix in this regime will then take the form In App.F we show that the form of the eigenvectors O can be found explicitly and the parameter Θ can be expressed using Vandermonde determinants.The factorization of the latter is known, leading to the exact final expression where k l for 1 ≤ l ≤ N −1 are the wavevectors of the bulk modes, defined in Eq. (58), and K is the spatial decay rate of the edge mode.As with the analogous expression for periodic boundary conditions (54), this product can be evaluated with the help of complex integration techniques, which we describe in Appendix G.We find that the large N asymptotic expression of Θ takes the form FIG. 2: Qualitative plot showing the behaviour of the system entropy on different timescales in the mixed(red) and purifying(black) phase.In the mixing phase, the entropy undergoes a period of very slow logarithmic decay up to exponentially long times in the system size.
Again, we focus on the regime where t scales polynomially with N , in which case the right hand side of the above is small.Moreover, noting that Θ should be no greater than unity, we find that the above expression should only be trusted in the regime We view the above constraint as a condition for validity of the approximation (67) made earlier.Then, in this regime we can make the same series of approximations as before to relate Θ(t) to the Rényi entropy.Thus, for t c ≲ t ≪ e N K , we find At very long times, when t scales exponentially with N such that |log Θ(t)| ≫ 1, we find that the entropy decays as Together, Eqs.(71, 72) characterize the salient features of purification dynamics in our model in the mixed phase sufficiently far from criticality (e −N K ≪ 1).

Purifying phase
The purifying case corresponds to the regime g > 1, where measurements occur so often that they overcome the scrambling and an initially mixed state quickly becomes pure.Looking at the spectrum of the singleparticle matrix D, one finds that all eigenvalues are at least as large as the bulk gap ∆ = (g − 1), which sets a timescale t ≳ (g −1) −1 after which the correlation matrix Γ[−tD] will have converged close to its t → ∞ limit.
Using Eq. 59 we see that there is exactly one root x 0 in the interval 0 < x 0 < log g.For large enough N we have that N log g ≫ 1, so for the other roots tanh N x q /2 ≈ 1.This means that for t ∼ poly(N ) we can make the approximation S(2) (t) ≈ N x 0 , (73) so x 0 is the entropy density in the chain at time t.The equation determining x 0 can be written as (74) Taking t ≫ 1 we see that the solution is approximately and follows the expected decay rate set by the spectral gap ∆.Since the system is disordered in this regime, the result is expected to hold in the thermodynamic limit, irrespective of the boundary conditions.

Critical point
We now address dynamics at the critical point g = 1.While we are no longer able to reliably make an approximation of the kind (67) in the open boundary condition case, we find that the purification dynamics for periodic boundary conditions is amenable to analytical treatment in this regime.In particular, Eq. (59) continues to hold at g = 1, with the the equation for the roots x q now given as solutions to the equation 2t sinh This equation is still transcendental, making it difficult to find a universal expression for its solutions.However, we can study the three regimes t ≪ 1, 1 ≪ t ≪ N and t ≫ N separately.
Early times t ≪ 1.-In the initial time frame t ≪ 1 it can be shown that the entropy is approximately given by S The logarithmic divergence at the origin has a simple intuitive explanation: since the averaging is performed over the purities rather than the entropies and we work in a d → ∞ system, the only scenarios that contribute to the average purity at early times are those where the entire chain is measured.If this happens, the system is immediately purified, since we can neglect the unitary evolution at t ≪ 1.The entropy is then simply the logarithm of the probability that all qudits are measured within the time t, which is p ≈ (f t) N = (gt/2) N .This intuitive picture matches the exact answer we found above at criticality, and is expected to hold for all values of g in both the periodic and open boundary conditions.Intermediate times 1 ≪ t ≪ N .-Inthis regime we can assume N x q ≫ 1 but x q ≪ 1 for all solutions x q that contribute meaningfully to the value of Θ.Using these approximations we find that the x q are equally spaced and the formula for log Θ is calculated as a geometric sum.The final expression for the entropy is The algebraic relationship S(2) ∝ t −1 is an important feature and only occurs exactly at criticality.
Late times t ≫ N .-Finally, in the long time limit, we see that we can approximate Θ by where (a, q) ∞ is the q-Pochammer symbol.This leads to an exponential decay of the entropy

D. Comparison to result from field theory
Having derived expressions for the time-dependence of the Rényi entropy of our model of hybrid quantum dynamics, it is instructive to compare our findings to the approach introduced in Ref. [38].There, the authors invoke an effective field theory known as capillarywave theory, which was first developed to model the dynamics of domain walls in the low-temperature phase of the Ising model.The correspondence between the two is rooted in mapping between discrete-time hybrid quantum circuits and two-dimensional ferromagnets, see e.g.Refs.[14,27].The parameters of the theory are a phenomenological surface tension σ and inverse temperature β, and once these are fixed, it is possible to find approximations for the time-dependence of the Rényi entropy starting from a mixed initial state in the associated discrete-time monitored quantum circuit model.
Upon comparing their expression to our results, we find that the same universal features hold.In particular, for both cases, there is a marked regime of times t ∼ poly(N ) in the mixed phase during which the entropy decays as an extensive constant with a − log t contribution.The sensitivity to boundary conditions we see [Eq.(64) vs. (71)] can also be understood in the capillary-wave picture as a consequence of the difference in configurational entropies of the endpoints of a domain wall for periodic versus open boundary conditions.Moreover, by looking at the prefactor of the term proportional to N , we can relate the microscopic parameters of our model to the phenomenological parameters of the field theory; in particular, we can fix the code rate βσ = K = − log g, which vanishes non-analytically at the transition g = 1.

VI. DISCUSSION
Our work introduces a class of random unitary circuits following a brickwork geometry, where each unit cell performs an infinitesimally small unitary transformation.We show that the limiting case of the construction above leads to a continuous stochastic process through the many-body Hilbert space.We show that the non-equilibrium behaviour of statistical averages of a large class of operator-space entanglement measures (the Rényi entropies) of this dynamical process can be obtained as equilibrium partition functions in an effective quantum spin system, governed by a universal, timeindependent Hamiltonian.The construction relies on an initial microscopic Hamiltonian describing local interactions, but we prove that this only enters the effective quantum information dynamics by setting the overall timescale.
We only perform a thorough investigation of the second Rényi entropy, where the effective theory is the spin-1/2 ferromagnetic TFIM, with an integrability breaking term that becomes quadratically small in the local dimension d.The ground state becomes degenerate in the thermodynamic limit and it is ferromagnetically ordered.Taking a phenomenological perspective, the two types of stable ordering roughly correspond to the measuring agent having full knowledge or no knowledge about the state of the system.The lowest energy excitations are topological domain walls, and roughly represent the geometric boundaries of our knowledge.We show that local measurements can also be studied within the same framework by adding an extra state to the spins of the effective system.When the local tumbling rate of the microscopic Hamiltonian is sufficiently strong, this extra state is adiabatically eliminated, and the effect of measurements is to introduce a transverse magnetic field whose strength is proportional to the measurement frequency.When this exceeds a critical threshold, the system undergoes an Ising-type phase transition into a disordered phase.This is recognized as the purification transition observed in numerical studies of similar models [29,31].The signature of the transition is a logarithmically decaying residual uncertainty in the state of the system after a purification procedure using uncorrelated local measurements, which is present only if the system is in the ordered phase.
We identify the order parameter corresponding to the residual second Rényi entropy in the effective model and prove exact product expansion formulae that can be used to calculate it in both open and closed boundary conditions.Complex integration techniques are used to find thermodynamic limit approximations on various timescales, and we see that their scaling agrees with field theoretic arguments.The method is not restricted to the residual entropy of the whole chain, and could be adapted to calculations of other second Rényi entropies.Universal characteristics of the transition such as the critical exponents must be the same as for the effective 1D quantum Ising theory.The transition in the von Neumann entropy requires higher replica analysis and may be of a different universality class, but we expect a qualitatively similar behaviour away from criticality.We begin the investigation of higher replica calculations by proving a formula for the matrix elements of the effective Hamiltonian in App.B. In contrast with similar models studied in literature, our circuits are not expected to lead to a percolation transition of the von Neumann entropy, even in the d → ∞ limit.This is because the small gate action limit ∆t → 0 is taken first, making the Hartley entropy S 0 undefined for any d.
To simplify our calculations, we have set the local tumbling rate Γ to infinity, but it may be interesting to investigate how it affects the transition.This introduces measurement inertia, wherein less information is gained by consecutively measuring the same qudit at intervals less than ∼ 1/Γ.If the measurement frequency grows beyond this, the qudits become effectively Zeno-locked.To the best of our knowledge, the growth of entanglement in this regime has not been previously investigated.FIG.3: Tensor contraction diagrams corresponding to the 2 non-trivial matrix elements found for the case of k = 2 replicas.Note that the complex conjugate operators were transposed to better fit in the plane, so that the top operators become the inverses of the bottom ones.The convention of having the bottom legs as inputs and top legs as outputs is also reversed for these blocks.
where the matrix T µν,στ defined in the last line is the matrix we are looking for, which tells us how vectors in V (S 2 ) ⊗ V (S 2 ) evolve under the action of a unit cell.If we let T and T ⋆ denote the 4 by 4 matrices with entries defined above, and let Wg be the 2 by 2 matrix with elements Wg στ = Wg(στ −1 ), then we can compute the transfer matrix knowing the results of all contraction diagrams from the equation Since the unit cell is just the identity for t = 0, the matrix T is simply the identity at zeroth order when expanding in powers of t.The correction to T ⋆ at order t 2 was computed above and can be written as If we plug this correction into Eq.A8 we get the result presented in Eq. 7 of the main text.

Appendix B: Effective dynamics for multiple replicas
Here we give an overview of the method used to compute the transfer matrix using a number of replicas k > 2. This generalizes the calculation discussed in App.A and it can be reduced to a simple exercise in combinatorics.Suppose we have T the unit cell operator and denote ⟨κϵ| T |στ ⟩ the result of contracting its legs according to permutation operators κ, ϵ, σ, τ ∈ S k .This will lead to diagrams similar to those seen in Fig. 3, but with k exponential operators of each type.To compute this, we must once again look at the Taylor expansions of the operators.Contributions at order t 2 come from the second order term in the expansion of each box and from all pairings of operators expanded to first order.All of these will reduce to one of three possible diagrams, with coefficients of d raised to some power depending on the number of loops that appear.These three are tr H 2 , tr tr 1 (H) 2 and tr(H) 2 , the same terms that appear in the definition of Ω(H).
To find the result of the calculation, consider the following auxiliary construction.Let i = 1, . . ., k and i = 1, . . ., k label a series of boxes and let σ, κ ∈ S k .Denote the first set of boxes by B and the second by B. Now imagine that we connect the bottom of the i boxes to the top of the i boxes according to the map i → σ(i), and similarly the top of the i boxes to the bottom of the i boxes according to i → κ(i).The result is that the boxes are now tied together in chains of varying lengths, each chain containing an equal number of i and i boxes.If we call l the number of chains, then we let µ = 1, . . ., l index the chains and form sets S µ as These sets exactly partition the full set of boxes with no overlaps We now return to our initial problem and consider the sets S L µ and S R ν formed by the left (σ, κ) and right (τ, ϵ) pairs of permutations.The indices are allowed to go from 1 to the number of chains formed by the respective pair.Denote by M µν the number of elements from B that are found in both S L µ and S R ν , and by M µν the number of elements from B found in the same 2 partitions.Remarkably, counting the number of diagrams of the three types mentioned above that show up in the transfer matrix reveals that the correction to second order in t is given by the simple expression (B3) This result shows that even for a higher number of replicas, the characteristics of the evolution are still independent of the microscopic Hamiltonian H, which only sets the overall timescale through the same function Ω(H) as in the k = 2 case.The interactions produced by this transfer matrix may have an interesting interpretation in terms of domain walls, but a discussion of this is beyond the scope of the current work.
Appendix C: Evolution of the X state Following a similar calculation of the matrix elements in the transfer matrix as shown in App.A we can also derive the unitary evolution of the X states that appear after a measurement.If we introduce a projection operator onto this state P X , such that P X |X⟩ = |X⟩ and P X V (S 2 ) = 0 and denote by Q X its complement so that P X + Q X = I, the total effective Hamiltonian is given by where H ij is the restriction to V (S 2 ) of the effective Hamiltonian in Eq. 10 of the main text, the operator X is given by and the XX interaction energy is .
(C3) We see that the only effect of Γ is to raise the energy of configurations that include X states.It should then be intuitively clear that taking a very large Γ will result in the dynamics being projected into the X-free subspace of the Hilbert space.The rest of the appendix constitutes a formal proof of this fact.
Denote by P the projector onto the X-free subspace of the Hilbert space and Q = I −P its complement.To simplify notation let H denote the full effective Hamiltonian over the course of this proof.This includes interactions of the type H M ij above between all nearest-neighbors in the chain and selective measurements M at some frequency f on all sites (defined in Eq. 19).Then we can separate H as where the diagonal parts become H 0 and the coupling between sectors is ∆H.If we let U = exp(−tH) be the imaginary time propagator then its evolution equation is We move into the interaction picture of H 0 by letting U(t) = exp(−tH 0 )U I (t).The evolution equation of U I (t) is The boundary conditions on all quantities of interest lie in the X-free subspace of the Hilbert space and configurations that include X states have 0 overlap with this sector.Therefore, we are only interested in the reduced propagator PU I (t)P.This evolves according to We see that the evolution equation requires knowledge about QU I (t)P.This evolves according to We integrate the two equations from 0 to t and substitute the form of QU I (t)P into the PU I (t)P equation to arrive at the integral form Let us now look at the operator norm of the correction term on the RHS.Since this is both sub-additive and sub-multiplicative we have Our aim is now to find a bound for each norm in the expression.The connection between the different sectors ∆H is bound and independent of Γ, as it appears only due to the measurements and is not a result of the random unitaries.The PH 0 P drives evolution in the X-free subspace, which was also shown to be independent of Γ. Therefore we define m(t) = sup t ′ ∈(0,t) ∥Q∆HP∥ Pe t ′ H0 P∆HQ , (C12) a real and continuous function.Note that m(t) is strictly positive and non-decreasing by construction.Furthermore, it is independent of Γ by the previous arguments.For the QH 0 Q part of the Hamiltonian, the energy of all configurations is raised by some constant proportional to Γ, as it is clear from Eq. C2 and Eq.C3.Therefore there must exist positive constants a, b such that the following inequality holds Qe −(t ′ −t ′′ )H0 Q ≤ e −(t ′ −t ′′ )(aΓ−b) , (C13) for all t ′ > t ′′ .The energy E(Γ) = aΓ − b can be interpreted as a lower bound on the ground state energy of the restricted Hamiltonian QH 0 Q.We are only interested in the large Γ limit, so we will assume E(Γ) > 0. For the final part we have since the norm of a projector is always 1.If we introduce the notation ψ(t) = ∥δU(t)∥/m(t) and combine all of the previous inequalities, we have (C15) To proceed, we invert the order of the integrals, so that

(C16)
We can now perform the integration over t ′ .A quick calculation shows that It is convenient to separate the Hamiltonian into an entangling and a non-entangling sector.We define the non-entangling sector of the basis to be spanned by those operators which act trivially as identity on at least one of the qudits, and the entangling sector is everything else.We can write this split as and equivalently for the matrices In a microscopic physical picture the H int stores information about the nearest neighbor interactions in the chain, while H loc generates local dynamics such as interactions with magnetic fields.
For simplicity, we can assume tr(H) = 0 such that h 00;00 = 0. Additionally, since we assumed H is symmetric under swapping the two qudits, we get (D13) so from Eq. 8 we have showing that only the interaction part of the Hamiltonian leads to transfer of information, as expected.
Appendix E: Evaluation of Θ(t) for periodic boundary conditions With periodic boundary conditions, the noninteracting fermionic Hamiltonian introduced in Eq. ( 39) can be block diagonalized in terms of states with definite quasimomentum.Working in units where γ = 1, a standard computation gives (E1) where a k = j e −ikj a j are annihilation operators for fermions with wavevector k.The set of k-space points K p in the above sum depends on the fermion parity sector p [49].Assuming an even number of sites N , the odd parity sector p = 1 has periodic boundaries (PBCs) e ikN = +1 while the even parity sector p = 0 has antiperiodic boundaries (ABCs) e ikN = −1.Since the anomalous terms pair +k and −k modes, it is helpful to single out the positive k-modes Note that the p = 1 case has one fewer k value, which is resolved by accounting for modes k = 0 and k = π, for which the anomalous terms vanish Overall, this allows us to write the Hamiltonian in matrix form We can write H k = ⃗ h(k)•⃗ σ, where ⃗ σ is a 3-vector of Pauli matrices, and we have When this k-space Hamiltonian H k is diagonalized, we obtain energies consistent with the dispersion relation (56) quoted in the main text.As a reminder, the object we wish to calculate is Now, we note that the state |+⟩ (|−⟩) is the ground state of the Hamiltonian H eff in the even (odd) parity sector, for the specific case g = 0.Both objects in each of the traces in Eq. (E8) are k-diagonal, and so both the numerator and denominator of the above become a product over k.Each factor of the product is given by an overlap between two states of a two-level system spanned by |VAC⟩ and a † k a † −k |VAC⟩, namely the overlap between the ground state with g = 0 and the thermal state with temperature t for nonzero g.We should be careful to include the modes k = 0, π, that are present in the denominator, separately.Denoting ρ +,k as the (pure) 2 × 2 density matrix for mode k in the state |+⟩, we have ρ Meanwhile, using the diagonal representation of H k we get where is the 3-vector that specifies the ground state of the 2level Hamiltonian H k .Now, since we can compute the dot product ⃗ n +,k • ⃗ n H,k = (sin 2 k + cos k(cos k − g))/ϵ k = (1 − g cos k)/ϵ k , and the relevant factor for each k becomes Finally, the relevant factor for the k = 0, π modes comes from recognizing that in the state |−⟩, we have a † 0 a 0 = 1 and a † π a π = 0.Because H p=0 conserves these occupations, we get a factor of e t in the numerator.Overall, we obtain the expression which was quoted in the main text, Eq. (54).Now, to evaluate the products in the above, we define the function where z is now a complex variable which corresponds to the wavevctor k on the real axis.If one considers the integral of f (z) around a thin contour Γ encircling the real interval [−π, π], cutting the line at the edges at a straight angle, then using the residue theorem we can show that By deforming the integration contour, the above can be re-expressed as which integrates on both sides of the branch cut chosen along the imaginary line.We do not include it explicitly but we keep in mind that the contour has a small gap such as to not cross the real line.Consider the rotated θ function θ r (x) = −iθ(ix) (we leave the t-dependence implicit), and equivalently for the energy function ϵ r .We have where we introduced tan ϕ(x) = (g cosh(x) − 1)/ϵ r (x).
By symmetry we can focus on the integral above the real line and have where K = − log g is the point where ϵ r (x) = 0. Consider the Weierstrass factorization of the cosine and apply it to θ r .The integral will exactly cancel for all terms that are analytic in the top half plane, so we only need to keep the factors that have a 0 above K on the imaginary line.The 0 at K is exactly cancelled by the denominator in the expression of θ r .Then the integral is .

(E21)
Since the function tϵ r (x) + ϕ(x) is strictly increasing from −π/2 to ∞ as x goes from 0 to ∞ each of the q terms will have exactly one node on the integration line.The value of the logarithm differs by exactly 2πi starting from the node up to ∞.If we denote by x q the unique solution of the equation tϵ(x) + ϕ(x) = qπ/2, the value of the integral is (E22) Then we see that as claimed in the main text.
The normalization constants a j are chosen such that for all j (and the additional j = 0 mode) although their explicit values are not important in our calculation.
The wavevectors k j are in [0, π) (except for k 0 , which becomes imaginary at the transition) and are the solutions of the equation The corresponding eigenvalues are given by the wellknown dispersion relation for the TFIM λ 2 = 1 + g 2 − 2g cos k. (F11) Since the matrices u and v are split into even and odd sectors one can show there must be some linear dependence, such that det u = det v = 0.This problem disappears when looking at the i, j ≥ 1 sector, such that the inverses u −1 ij and v −1 ij exist.We can then proceed to calculate the determinants in the expression for Θ.We saw that at long enough times, we can make the approximation which becomes the same as Eq.68 in the main text with the rewriting of the complex wavenumber K = ik 0 .
Appendix G: Large-N asymptotic expression of entropy In this appendix, we will derive the large-N asymptotic behavior of the exact product expansion of Θ in Eq. 68 of the main text.We start by introducing the auxiliary complex functions f (z) = g sin z 1 − g cos z , ϕ(z) = arctan f (z) = 1 2i log 1 − ge −iz 1 − ge iz .

(G1)
In terms of these, we can write the quantization equation as sin(kN − ϕ(k)) = 0, where This function is not analytic at the simple poles z = k l and z = ±i log g and at the branch points z = ±iK.Note that K becomes exponentially close to − log g in the limit of KN → ∞.If we integrate the function around a contour that encircles (but remains close to) the real axis interval [−π, π] and apply the residue theorem we obtain where we note that only half the residue at the edge poles is included.Comparing this to the product expansion of Θ we find that log Θ = log tanh We are free to deform the integration contour without changing the result of the integral, so long as we do not cross any singularity.The integrals over the sides of the strip are equal to We will focus on the contour integration around the singularity in the top half of the complex plane z ∼ iK.The integration on the bottom half follows by analogy, but we find that it has no contribution to the order of approximation we are interested in, so we simply discard it.If we make the change of coordinate z → iz ′ = z − i log 1 g and consider points sufficiently far from the origin that along with the additional approximation that −N log g ≫ 1, then g ′ (z ′ ) = g(z) is simplified to (G10) In this coordinate, the integration is performed along a dumbbell contour, going along the branch cut on the positive half of the real axis and encircling the singularity at the origin.The radius ϵ is always understood to obey the condition in Eq.G9.Since the function decays very fast along the positive real axis for large N , we can make a further approximation on our contour (G11) where the error we incur is of order O( 1 KN ) in the integral.We can then express the integral contribution to log Θ as ) log z 2 .

(G12)
The circular region has a contribution of (G14)

FIG. 1 :
FIG. 1: Schematic representation of the random circuit geometry for open boundary conditions.The construction of the unit cells is illustrated on the right.The Hamiltonian H and evolution time ∆t are kept fixed, but the random unitaries U, V are sampled independently at each spacetime location in the circuit.

d 2 d− 1 a 2 = 2d 3 (d 2 − 1 ) 2 d− 1 ahh 2 = d 2 d− 1 a 2 d− 1 a
D7) and the hermiticity condition H = H † of the Hamiltonian implies h −a−b;−s−t = h * ab;st ω ab+st , (D8) where the inverse −a is with respect to modulo d addition a ⊕ (−a) = 0. We can start by calculating Γ for this Hamiltonian.If we notice that all elements in our basis other than I are traceless, we have tr 1 (H) = d−1 a,b,s,t=0 h ab;st X s d Z t d tr X a d Z square this result and apply the braiding equation D3 we get tr 1 (H) 2 = d 2 d−1 a,b,s,t=0h 00;ab h 00;st X a d Z b d X s d Z t d = ,b,s,t=0h 00;ab h 00;st ω at X a⊕s d Z b⊕t d .(D10)Finally,we take the trace of this and use Eq.D7 to get Γ = 2d (d 2 − 1) 2 tr tr 1 (H) ,b,s,t=0 h 00;ab h 00;st ω at tr X a⊕s d 00;ab h 00;st ω at δ a⊕s,0 δ b⊕t,0 00;ab h 00;−a−b ω −ab can write this more compactly in terms of the Frobenius norm of the local part of the HamiltonianΓ = d 4 (d 2 − 1) 2 ∥h loc ∥ 2 F ,(D12)with a coefficient that goes to 1 in the d → ∞ limit.It is clear that this part of the Hamiltonian is only relevant in randomizing the state of individual qudits after measurement and does not play any role in the transfer of information across the chain.For the information transfer rate γ, we can compute tr H ,b,s,t=0 h ab;st h −a−b;−s−t ω −ab−st = d ,b,s,t=0 |h ab;st | 2 ,