Hayden-Preskill Recovery in Hamiltonian Systems

Information scrambling refers to the unitary dynamics that quickly spreads and encodes localized quantum information over an entire many-body system and makes the information accessible from any small subsystem. While information scrambling is the key to understanding complex quantum many-body dynamics and is well-understood in random unitary models, it has been hardly explored in Hamiltonian systems. In this Letter, we investigate the information recovery in various time-independent Hamiltonian systems, including chaotic spin chains and Sachdev-Ye-Kitaev (SYK) models. We show that information recovery is possible in certain, but not all, chaotic models, which highlights the difference between information recovery and quantum chaos based on the energy spectrum or the out-of-time-ordered correlators. We also show that information recovery probes transitions caused by the change of information-theoretic features of the dynamics.

Introduction.A central challenge in modern physics is to characterize the dynamics in far-from-equilibrium quantum systems.The Hayden-Preskill protocol [1] offers an operational approach toward this goal and has been attracting much attention .The protocol addresses if the information initially localized in a small subsystem can be recovered from other subsystems after unitary time evolution.If the unitary dynamics is sufficiently random, the information is rapidly encoded into the whole system, and information can be recovered from any small subsystem [1].This phenomenon is can be thought of an information-theoretic manifestation of complex quantum dynamics and is called the Hayden-Preskill recovery.The unitary dynamics that leads to the Hayden-Preskill recovery is referred to as information scrambling [1].
Despite these progresses, information recovery in Hamiltonian systems has been rarely explored [40].It is widely believed that the Hayden-Preskill recovery is possible in quantum chaotic systems, but the original analysis strongly relies on the random unitary assumption, which is unlikely to be satisfied even approximately by time-independent Hamiltonian dynamics [8,41].Furthermore, quantum chaos is commonly characterized by eigenenergy statistics, which is a static property, but information recovery is about the dynamical properties.Thus, the relation between quantum chaos and information recovery is not a priori trivial.
In this Letter, we investigate in detail the information recovery in various Hamiltonian systems.We first provide a class of Hamiltonians that do not lead to information scrambling.This includes chaotic spin-1/2 chains, such as the Heisenberg model with random magnetic field and the mixed field Ising model.Notably, they saturate OTOCs for local observables but do not achieve the Hayden-Preskill recovery, demonstrating the difference between the saturation of OTOCs for local observables and information scrambling.
We then confirm information scrambling in the Sachdev-Ye-Kitaev (SYK) Hamiltonian [11,12,42,43], which is a canonical holographic dual to quantum gravity [13,14,[44][45][46].The SYK model is known to have 'scrambling' features in many senses, such as saturation of OTOCs [11,12] for local observables, the maximum quantum Lyapunov exponent [47], and an RMT-like energy statistics [44,48,49].Our result adds another scrambling feature to the model, that is, it achieves the Hayden-Preskill recovery.We also show that sparse variants [50] achieve the information recovery as well, possibly helping experimental realizations of the protocol.
We finally address the question whether the information recovery can reveal novel quantum many-body phenomena.Using a variant of SYK models, we affirmatively answer to this question: the information recovery can capture a transition that was previously unknown.The transition is caused by a drastic change of informationtheoretic structures of the Hamiltonian dynamics.This is of interest as it characterizes complex quantum manybody dynamics from the quantum information viewpoint.
The Hayden-Preskill protocol.Given quantum manybody system S of N qubits, we encode quantum informa-qubits qubits qubits qubits FIG. 1.A diagram of the Hayden-Preskill protocol.Time flows from bottom to top.Horizontal lines imply that the qubits connected by the line may be entangled.The initial states on AR and BB ′ are given by a maximally entangled state of k ebits, i.e., k EPR pairs, that keeps track of quantum information in A, and a purified state |ξ(β)⟩ of a thermal state on B at the inverse temperature β, respectively.The system S := AB undergoes Hamiltonian dynamics by ĤS, and then, is split into two subsystems, C of ℓ qubits and D of N − ℓ qubits.By applying a quantum channel D onto B ′ C, one aims to decode the quantum information in A, that is, to recover the k EPR pairs between Â and R.This protocol has a natural interpretation in the context of information paradox [1].See also a tutorial [51].
tion into a local subsystem A ⊆ S of k qubits (k ≪ N ).We then let the system S undergo the Hamiltonian timeevolution U Ĥ (t) := e −i Ĥt for some time t, where Ĥ is the Hamiltonian in S.After the time-evolution, the information is tried to be recovered from an ℓ-qubit subsystem C ⊆ S. Throughout our analysis, we assume C ⊆ B := S \ A as far as ℓ ≤ N − k.The question is how large ℓ should be for a successful recovery.
The answer depends on the initial state in B as well as available resources in the recovery process.Here, we assume that the initial state in B is a thermal state ξ B (β) at inverse temperature β.Based on its eigendecomposition ξ B (β) = j p j (β) |ψ j ⟩⟨ψ j | B , we introduce a purified state |ξ(β)⟩ on the system BB ′ .When the subsystem B ′ is traced out, the marginal state on B is the original thermal state ξ B (β).We consider the scenario in which the subsystem B ′ can be used in the recovery process.This has a natural interpretation in the black hole information paradox [1] and is of considerable interest.See Fig. 1.
The recovery of quantum information is formally defined by introducing a virtual reference system R that keeps track of the quantum information.Denoting k Einstein-Podolsky-Rosen (EPR) pairs bertween A and R by |Φ⟩ AR , we set the initial state on the system . The subsystem S undergoes the Hamiltonian time-evolution by ĤS , resulting in the state ) be the marginal state on CB ′ R, which is given by taking the partial trace over D of |Ψ(t, β)⟩ SB ′ R .In the following, we indicate the subsystem over which the partial trace is taken by omitting the subsystem from the superscript.
Following the standard convention [1], the recovery error is defined by (1) Here, the minimum is taken over all possible quantum operations D, namely, all completely-positive and tracepreserving (CPTP) maps, from CB ′ to A, and ∥ρ∥ 1 = Tr ρ † ρ is the trace norm.Due to the Holevo-Helstrom theorem [52], the trace norm between two quantum states characterizes how well they can be distinguished and is suitable to quantify the recovery error.We normalize ∆ Ĥ (t, β) so that 0 ≤ ∆ Ĥ (t, β) ≤ 1. See S1 [53].
Computing ∆ Ĥ (t, β) is in general intractable due to the minimization over the CPTP maps.The decoupling condition provides a necessary and sufficient condition for the recovery in terms of the state Ψ DR (t, β) [54][55][56].Using the condition, a calculable, and typically good, upper bound on ∆ Ĥ (t, β) is obtained: where Σ Ĥ (t, , and π R = I R /2 k is the completely mixed state on R. See also S1 [53].Note that the quantity Σ Ĥ (t, β) is closely related to the mutual information between R and D, which is equivalent to the OMI [19], but the OMI leads to a worse bound than Eq. ( 2).The recovery error ∆ Ĥ (t, β) is also related to OTOCs.If one can compute OTOCs for all observables on the kqubit subsystem A and the ℓ-qubit subsystem C, or all the 4 k+ℓ operators that form an operator basis on AC, the recovery error could be evaluated [57,58].However, this is computationally intractable as OTOCs for at least 4 k+ℓ operators are needed.Note that the existing studies of OTOCs in Hamiltonian systems are mostly about the cases with k = ℓ = 1, and hence, do not provide much insight into the recovery error.
Random unitary model and information scrambling.The Hayden-Preskill protocol was understood well in a random unitary model, where the time evolution e −i ĤS t is replaced with a Haar random unitary.The model does not have a parameter corresponding to time t, and its recovery error ∆ Haar (β) satisfies, with high probability, ∆ Haar (β) ≤ ∆Haar (β) := min{1, 2 2 ] is the Renyi-2 entropy [1,24,56].
From Eq. (3), ∆ Haar (β) ≪ 1 if ℓ ≫ ℓ Haar,th (β).In particular, ℓ Haar,th (0) = k.Hence, if the system B is initially at infinite temperature, the k-qubit quantum information 2. Semilogarithmic plot of the late-time values of ∆ for ĤXXZ, ĤIsing, and ĤSYK 4 .To compare, the values for the Haar random unitary are also plotted.k = 1 for all models, and N = 12 (Nq = 12) is chosen for the XXZ and Ising (SYK4) models.The dimension of the Haar random unitary is 2 12 .Note that the conservation of the z component of spins (parity) is considered for the Heisenberg spin chain (SYK model).The values of (g, h) are those discussed in [61].The averages for t = (1, 2, . . . ., 10) × 10 6 are plotted as the late-time value.For random-field average, 16 samples are taken.
in A is recoverable with exponential precision from any subsystem of size larger than k, which is independent of N .This phenomena is referred to as the Hayden-Preskill recovery.Following the original proposal [1], we refer to the dynamics achieving the Hayden-Preskill recovery as information scrambling.
Hamiltonians without information scrambling.Due to the facts that the information scrambling occurs in the random unitary model and that quantum chaos can be characterized by RMT, information scrambling has been commonly studied in relation with quantum chaos.However, information scrambling is not necessarily related to the quantum chaos in terms of the RMT-like energy statistics.For instance, we can analytically show that the dynamics of any commuting Hamiltonians is not information scrambling (see S2 [53]), while they can have RMT-like features [59,60] in the energy spectrum.
More illuminating instances are the spin-1/2 chains such as the Heisenberg with random magnetic field, ĤXXZ = , where h j are independently sampled from a uniform distribution in [−W, W ], and the mixedfield Ising with constant magnetic field, ĤIsing = − Here, X j , Y j , and Z j denote the Pauli matrices on site j.In contrast to the fact that both have integrable-chaotic transitions by varying parameters , our numerical analysis reveals that the recovery errors are not as small as the random unitary model for any values of parameters  at any time t.This is demonstrated in Fig. 2, where the late-time values of ∆ Ĥ (t, β = 0) are plotted for these Hamiltonians.Hereafter, we set k to 1 in all numerics throughout the paper for the sake of computational tractability.It is observed that, by increasing ℓ, ∆ Ĥ decays inversepolynomially or more slowly.We also provide in S3 and S4 [53] evidences that these values do not depend on the system size N .This implies that, although Fig. 2 is for N = 12 and ℓ ≤ N − 1 = 11, we can infer ∆ Ĥ for larger N and ℓ by extrapolation.By doing so, we may observe that ∆ Ĥ for W ≳ 1 may possibly stay nearly constant in the large-N limit unless ℓ ≈ N .
We can also investigate lower bounds on the recovery errors based on the mutual information, which we denote by ∆ Ĥ (see S1 B [53] for the derivation).The bound is not tight, but we show in Fig. 3 that the lower bounds for ĤXXZ and ĤIsing scale similarly to those of the upper bounds.That is, they decay inverse-polynomially or more slowly as ℓ increases and possibly stay almost constant if W ≳ 1 unless ℓ ≈ N .
As both upper and lower bounds scale similarly, we reasonably conclude ∆ Ĥ = Ω(1/poly(ℓ)) in the large-N limit.This is in sharp contrast to the exponential decay of the recovery error in the random unitary model and implies that the dynamics of these Hamiltonians is not information scrambling in any parameter region.
More closely looking at Figs. 2 and 3, ∆ Ĥ is likely to be dependent on the parameters of the Hamiltonians.It is known that ĤXXZ shows integrable-chaotic-MBL transitions as W increases and that the system is chaotic for W ≈ 0.5 [66,70].However, this chaotic transition does not seem to have strong consequence to information scrambling as both upper and lower bounds on ∆ Ĥ for W = 0.5 are only slightly smaller than that in the integrable case with W = 0.For ĤIsing , while the parameter (g, h) = (1.08,0.3) leads to the most chaotic feature in the entanglement structure [61], both upper and lower bounds on ∆ Ĥ for that value can be worse than other parameters.This also indicates that information scrambling differs from quantum chaos and may not be able to be inferred from static features of Hamiltonians.
The fact that information scrambling is not observed in these systems does not contradict to the saturation of OTOCs for local, typically single-qubit, observables at late time when the parameters are appropriately set [83][84][85][86][87][88].Our numerital results rather indicate that OTOCs for multi-qubit observables are not saturated in such cases [57,58], which may be of independent interest.
Original and sparse SYK Hamiltonians.From these results, it is likely that more drastic Hamiltonians are needed to achieve the information scrambling.We next consider the SYK model, SYK 4 , of 2N q Majorana fermions: ĤSYK4 = 1≤a1<a2<a3<a4≤2Nq J a1a2a3a4 ψa1 ψa2 ψa3 ψa4 , (4) with ψj being Majorana fermion operators.The couplings J a1a2a3a4 are independently chosen at random from the Gaussian with average zero and σ Since the parity symmetry of SYK 4 leads to deviations in the information recovery [9,10,24,89], we focus on the even-parity sector and set N = N q − 1.The recovery error of the corresponding random unitary model is given by ∆′ Haar (β) = min{1, 2 (ℓ Haar,th (β)−ℓ)− 1 4 }.See S5 [53] for details.We have also checked that the effect by the periodicity, characterized by N q mod 4, is negligible.
In Fig. 4, we numerically plot the upper bound on the recovery error, ∆SYK4 (t, β = 0), against time t for various ℓ.It clearly shows that ∆SYK4 quickly approaches ∆′ Haar .This is also the case for β > 0. We estimate that ∆SYK4 converges to ∆′ Haar before time t = O( N q ), which qualitatively supports the fast scrambling conjecture [2][3][4].Hence, the SYK 4 dynamics, while differs from Haar random, has an excellent agreement with the prediction by RMT and achieves the Hayden-Preskill recovery.See also Figs. 3 and 2.
The situation remains the same even for a sparse simplification of SYK 4 , spSYK 4 .In spSYK 4 , the number of non-zero random coupling constant is fixed to K cpl .It recovers SYK 4 when K cpl = 2Nq 4 , but K cpl = O(N q ) is known to suffice to have chaotic features and to reproduce holographic properties [90,91].
In Fig. 5, we plot the upper bound on the recovery error for a further simplified sparse SYK model (±spSYK 4 ), in which a half of the non-zero couplings is set to 1/ K cpl and the other half to −1/ K cpl [92].We observe that, when K cpl ≳ 30 = O(N q ), this simplification does not change the upper bound on the recovery error from that of the Haar value.Hence, ±spSYK 4 with K cpl = O(N q ) suffices to reproduce information-theoretic properties of SYK 4 as well as its chaotic features.As this number of non-zero couplings is substantially smaller than the original SYK 4 , which has K cpl = O(N q 4 ), this would help experimental realizations of the Hayden-Preskill protocol in many-body systems.See S6 [53] for details.
Probing transitions by the Hayden-Preskill protocol.Yet another SYK model attracting much attention is the SYK 4+2 model [93].The Hamiltonian is where ĤSYK2 = i 1≤b1<b2≤2Nq K b1b2 ψb1 ψb2 , and θ ∈ [0, π/2] is a mixing parameter.The coupling constants {K b1b2 } satisfy K b2b1 = −K b1b2 and are normalized for the variance of eigenenergies of ĤSYK4+2 (θ) to be unity.The SYK 4+2 model has a peculiar energy-shell structure in the sense of the local density of states in Fock space, which shows drastic changes by varying θ.Accordingly, the range of θ ∈ [0, π/2] is divided into four regimes I, II, III, and IV [94,95].In I, only one energyshell is dominant in the whole Hilbert space, and it is quantum chaotic.As θ increases the size of the energyshell becomes diminished, and O(poly(N q )) energy-shells appear in II and III.The energy statistics remains RMTlike in these two regimes.Characterizing physics in II and III has been under intense investigations [96].In IV, the number of energy-shell approaches O(exp(N q )), and Fock-space localization is observed.
For sufficiently small and large θ, the behaviour of ∆SYK4+2 can be naturally understood.For small θ, the model is approximately SYK 4 .As we have observed above, the dynamics of SYK 4 quickly achieves the Hayden-Preskill recovery.Hence, this should also be the case in the regime I.In contrast, for sufficiently large θ, the model is almost SYK 2 and the Fock-space localization occurs.Thus, information recovery should not be possible, resulting in the absence of information scrambling in the regime IV.In contrast, ∆SYK4+2 is smoothly changing for the intermediate values of θ, which is seemingly in tension with the division of the regimes II and III in terms of the energy-shell structure.
To understand the intermediate plateau, we shall recall that, in II and III, transitions from one energy-shell to the other are strongly suppressed, which effectively results in the division of the whole Hilbert space into O(poly(N q )) energy-shells [94].Additionally, it is known that the dynamics in each energy-shell seems to be approximately Haar random within the subspace [95].The intermediate plateau of ∆SYK4+2 can be explained from these common features in II and III.Since the whole Hilbert space is effectively divided into smaller ones, within which the dynamics remains still Haar random, the unitary dynamics in II and III induces partial decoupling [97] rather than decoupling.In this case, the recovery error is given in the form of 2 ℓ ′ th −ℓ + ∆ rem [24].Here, ℓ ′ th ≈ ℓ Haar,th + O( √ k) and ∆ rem quantifies the amount of information that cannot be recovered unless ℓ ≈ N q .As we set k = 1 in our numerics, ℓ ′ th is hardly observed in our analysis.In contrast, ∆ rem is clearly observed as an intermediate plateau.
It is known that ∆ rem is inverse proportional to the standard deviation of energy in D. As the standard deviation of energy in B shall be O( N q ), that in D is at least O( N q ).Hence, we can qualitatively estimate that ∆ rem = O(1/ N q ), which remains nearly constant unless ℓ ≈ N q .
From this perspective, we can understand the two transitions as reflections of the changes of decoupling properties.In I, the combined regime over II and III, and IV, the SYK 4+2 dynamics leads to full, partial, and no decoupling, respectively.Accordingly, each regime has qualitatively different behaviours in the information recovery.The emerging difference between II and III in the energy-shell picture should be an artifact due to the fact that the energy-shell is viewed in the Fock basis, which is not necessarily physically intrinsic to the system.
Summary and discussions.In this Letter, we have studied the information recovery in various Hamiltonian systems and have shown that information scrambling in the sense of information recovery does not always coincide with quantum chaos.Spin chains are unlikely to be information scrambling, while they are quantum chaotic in energy spectrum and saturate OTOCs for local observables.In contrast, the (sparse) SYK models are information scrambling and have the latter two properties.
We have also demonstrated a potential use of the information recovery protocol to find new transitions caused by a information-theoretic mechanics.
It is open if any local spin models can be information scrambling since the family of SYK models, the only models that are information scrambling in our analysis, does not have spatially local interactions.It will be also of interest to further explore the direction of characterizing various quantum phases in the information-theoretic manner.

The work was supported by Grants-in-Aid for
We apply this relation to Eq. (S4).A purification of Ψ DR and π R are |Ψ⟩ SB ′ R and |Φ⟩ ÂR , respectively.We denote by |σ⟩ DE a purification of σ D in Eq. (S4) by some system E.Then, there is a partial isometry V DB ′ → ÂE such that By tracing out DE, we obtain where D is a CPTP map obtained from the isometry V and the partial trace.As the left-hand side is exactly ∆ Ĥ by the identification of the system Â with A, this provides an upper bound on the recovery error.The lower bound in Eq. ( S3) is obtained using Eqs.(S5) and (S6): Here, the second and the second last lines follow from Eq. (S5), the third line from the Uhlmann's theorem, and the fourth line from the monotonicity of the fidelity under the partial trace.While the bounds in Eq. ( S3) replace the analysis of the recovery error with that of decoupling, the degree of decoupling Σ opt Ĥ is, in general, computationally intractable due to the minimization over all states in the system D. To circumvent this issue, the state σ to be minimized in Eq. ( S4) is commonly set to Ψ D .By doing so, we obtain a computationally tractable upper bound as where This is the bound we have used throughout the paper.Note that a computationally tractable lower bound cannot be obtained in this approach.

B. Upper and lower bounds from the mutual information
Instead of the decoupling approach, we may use the mutual information of the state |Ψ(t, β)⟩ SB ′ R to bound the recovery error, leading to Here, I(D : R) To derive the upper bound in Eq. (S16), we simply use the relation between the quantity Σ Ĥ , defined just below Eq. (S15), and the mutual information I(D : R) Ψ , namely, Σ Ĥ ≤ ln 2 2 I(D : R) Ψ (see, e.g., [98]).From this and the upper bound on the recovery error in Eq. (S15), we obtain As for the lower bound, let σD be argmin 1 2 ∥Ψ DR − σ D ⊗ π R ∥ 1 , and ΨDR be σD ⊗ π R .Using the fact that I(D : R) Ψ = 0, the mutual information can be rewritten as where we used that H(R) Ψ = H(R) Ψ as Ψ R = ΨR = π R .Applying the Alicki-Fannes-Winter inequality [99,100], the last expression is bounded from above by a function of the trace norm between Ψ DR and ΨDR , which is nothing but Σ opt Ĥ .We then have I(D : R) Ψ ≤ f (Σ opt Ĥ ), where f (x) = 2kx + (1 + x)h x 1+x , and h(x) = −x log x − (1 − x) log(1 − x) for 0 ≤ x ≤ 1 is the binary entropy.Since f (x) is monotonically increasing function for any x ≥ 0, it has the inverse f −1 , and we have f −1 I(D : R) Ψ ≤ Σ opt Ĥ .Together with the lower bound on ∆ Ĥ (t, β) in terms of Σ opt Ĥ (S3), we have It is obvious that Eq. (S16) is in general worse than Eq.(S3) because the former is obtained by substituting into Eq.(S3) the bounds on the degree of decoupling, Σ opt Ĥ , in terms of the mutual information I(D : R) Ψ .However, the lower bound in Eq. (S16) has an advantage from a computational viewpoint since it is calculable if the mutual information is.Based on this consideration, we used in the main text an upper and a lower bound on the recovery error such as Before we leave this section, we point out that I(D : R) Ψ has been studied as the OMI, especially for the infinite temperature (β = 0) [19,[35][36][37].Using Eq. (S16), the results in the literature can be rephrased in terms of the recovery error.

S2. COMMUTING HAMILTONIAN MODELS
We provide an in-depth analysis of the Hayden-Preskill protocol for commuting Hamiltonians.We mean by a commuting Hamiltonian the one in the form of where ĤA and ĤB are Hamiltonians non-trivially acting only on A and B, respectively, and ĤA:B is an interacting Hamiltonian between A and B that commutes with the other two.Below, we assume that A ⊆ D, which implies C ⊆ B as S = AB = CD.Due to the commuting condition, the time-evolving operator generated by ĤS is decomposed to e −i ĤS t = e −i ĤA:B t e −i ĤA t e −i ĤB t .Since a thermal state ξ(β) B ∝ e −β ĤB is invariant under the time-evolution by ĤB , we have Ψ DR (t, β) = Tr C e −i ĤA:B t e −i ĤA t Φ AR ⊗ ξ(β) B e i ĤA t e i ĤA:B t .We then use the property of maximally entangled state that, for any unitary U A on the system A, U A |Φ⟩ AR = Ū R |Φ⟩ AR , where ¯indicates the complex conjugate.
Denoting e −i ĤA t by U A (t), it follows that Ψ DR (t, where we used the unitary invariance of the trace norm and π R = I R /2 k in the last line.We now introduce an extended region B of B in terms of ĤA:B .Namely, B is the union of the system B and the set of qubits on which ĤA:B acts non-trivially.As we assumed that C ⊆ B, it holds that C ⊆ B. We denote by ∂ B the boundary of B in terms of ĤA:B , i.e., ∂ B = B \ B. By taking the trace over B \ C in the right-hand side Eq.(S25) and using the monotonicity of the trace norm under the partial trace, we have where the last line holds since the interaction Hamiltonian ĤA:B nontrivially acts only within B. Let A − be A\∂ B. As ∂ B ⊆ A, A − is not empty in general.Using this notation, we have Tr is a maximally entangled state between A − and R − , and π R+ is the completely mixed state in R + .Using κ, defined by the number of qubits in A that interact with B by ĤA:B Using Eq. (S5) and the fact that , we obtain Σ opt Ĥ (t, β) ≥ 1 − 2 κ−k .Substituting this into the lower bound in Eq. (S3), we have Thus, when κ < k, the recovery error ∆ Ĥ is bounded from below by a constant.
Although the lower bound in Eq. (S28) becomes trivial when κ = k, this is just due to the analytical derivation.To illustrate this, let us particularly consider the Sherrington-Kirkpatrick (SK) Hamiltonian given by ĤSK := − n<m J nm Z n ⊗ Z m , where J nm is chosen from a given distribution.Similarly to Eq. (S23), we first divide the Hamiltonian into three, such as ĤSK = ĤA + ĤB + ĤA:B , where ĤA and ĤB non-trivially act only on A and B, respectively.Following the above argument, it holds that We further decompose ĤA:B into ĤA:B = Ĥ(A:B)∩D + ĤA:C , where the former term non-trivially acts only on D and the latter is the rest.The unitary invariance of the trace norm leads to By applying a CPTP map onto R that maps ρ R to j ⟨j| ρ R |j⟩ |j⟩⟨j| R , where {|j⟩} j (j = 0, . . ., 2 k − 1) is the Pauli-Z basis in R, we have where Ω AR := 2 −k j |j⟩⟨j| A ⊗ |j⟩⟨j| R .Note that π R is invariant under the above CPTP map, and that the second inequality follows as C ⊆ B. We now use the relation that, for n ∈ A and m ∈ C, e iJnmZn⊗Zmt (|j⟩⟨j| A ⊗ ξ(β) B )e −iJnmZn⊗Zmt = |j⟩⟨j| A ⊗ e (−1) jn iJnmZm ξ(β) B e −(−1) jn iJnmZm , (S32) where we express j in binary as j 1 . . .j k (j n = 0, 1).The terms such as e (−1) jn iJnmZm then disappears when Tr B is taken.Hence, it follows that Hence, we have where we have used the lower bound given in Eq. (S5) and the last line follows from the direct calculation.This leads to

S3. THE HEISENBERG MODEL WITH RANDOM MAGNETIC FIELD
We consider a one-dimensional quantum spin chain with site-dependent random magnetic field to the z direction.The Hamiltonian is in which N is the number of S = 1/2 spins, h j are independently sampled from a uniform distribution in [−W, W ], J z is the ratio of the coupling in the z direction to that in the xy plane.The model has different features depending on the choice of J z and W .For J z = 0 and h j = 0 (j = 1, 2, . . ., N ), the Hamiltonian becomes ĤXY which is integrable.For J z = 1 and W > 0, the model has been extensively studied as a prototypical model of MBL in one spatial dimension [62][63][64][65][66][67][68][69].For system sizes accessible by numerical diagonalization, various measures of localization point to the MBL transition at finite critical W [70]. Note, however, that the location of the transition to the genuine MBL phase in the thermodynamic limit have been heavily debated in more recent studies.[71][72][73][74][75][76][77][78][79][80][81] In the following, we compute the upper and lower bounds on the recovery error for this model based on Eq. (S22), which are denoted by ∆XXZ and ∆ XXZ , respectively.When J z = 0, we denote the bounds by ∆XY and ∆ XY .
In Fig. S1, we plot the time dependence of ∆XY and ∆ XY for the XY model (J = 0) without random magnetic field.We set k = 1.Although there exists a gap between ∆XY and ∆ XY , both of them varies in a similar manner as time t increases.The plots at β = 0 for different N = 8, 10, 12 are qualitatively similar to each other, and they do not converge to a constant because we are here working with a single, integrable Hamiltonian for each N .Introducing finite temperature, β = 10, increases the value of ∆, reflecting the decrease of the effective dimension of the Hilbert space.
The time dependence of ∆XXZ and ∆ XXZ for the XXZ model with J = 1 are given in Fig. S2.The two quantities decay similarly.For W = 0, the Hamiltonian reduces to the Heisenberg model.As it is integrable, the values are not converging.For W = 1, the average over 16 samples is plotted.The late-time behavior of these averaged values of ∆XXZ and ∆ XXZ are smoother compared to the W = 0 case.The value decreases as ℓ increases.For ℓ = 9 = N − 1, we always obtain ∆XXZ = 0, which is because the system is in the S z = 0 subspace.
In Fig. S3, the sample averages of ∆XXZ and ∆ XXZ for various values of ℓ and W are provided.The upper and lower bounds again show similar time dependence.We observe that ∆XXZ monotonically increases and approaches unity as W increases.This is naturally expected since the system shows MBL for large W , which should prevent the quantum information initially localized in the subsystem A from spreading over the whole system S.It is worth noting that, even though a small W introduces a chaotic behavior to the system [66,70], ∆XXZ remains at high values, indicating the failure of information recovery.
In Fig. S4, we plot the late-time values of ∆XXZ and ∆ XXZ against ℓ for β = 0, W = 0.5, 2 and for various N .The plots almost overlap for all N except when ℓ = N − 1.Although our numerical analysis is for N up to 12, which forces ℓ to be at most 12 as ℓ ≤ N , this result indicates that we can infer the values of ∆XXZ and ∆ XXZ for larger ℓ simply by extrapolation.As both in Fig. S4 seem to decay inverse-polynomially or more slowly, we may reasonably conclude that the recovery error decays similarly as ℓ increases, that is, ∆ XXZ = Ω(1/poly(ℓ)) unless ℓ ≈ N .

S4. ISING MODEL WITH UNIFORM MAGNETIC FIELD
A translationally invariant spin chain with nearest-neighbor Ising-type interaction and magnetic field, is exactly solvable if g = 0 or h = 0.For other choices of (g, h), the model is non-integrable, and the case with (g, h) = (1.05,0.5) [82] has been often studied as a prototypical model of chaotic spin chain.In Fig. S5, we plot the values of ∆Ising and ∆ Ising for β = 0 and various ℓ.We observe that both decrease slowly as time increases and are likely to converge at the late time.We have also checked finite β.The values of ∆Ising and ∆ Ising are generally larger, indicating less efficient error correction.
For the late-time values of ∆Ising and ∆ Ising against ℓ for several values of N , see Fig. S4.Similarly to the XXZ case, the plots almost overlap except for ℓ = N − 1, and the values decay inverse-polynomially or more slowly as ℓ increases, except ℓ = N − 1.Thus, we may again reasonably conclude that ∆ Ising = Ω(1/poly(ℓ)) unless ℓ ≈ N .

C. Convergence time
To check the time scale needed for ∆SYK4 converging to ∆Haar , we plot in Fig. S8 ∆(t, β = 0) for ℓ = 4, 6, 8, 10 against re-scaled time t/ N q .The plot indicates that the convergence time is likely to be N q .

FIG. 6 .
FIG. 6.The late-time value of ∆SYK 4+2 (t, β = 0) plotted for various ℓ against the value of δ.Nq is set to 13.The number of samplings is 64.The lines connecting the data points are guide to the eye.The horizontal dashed lines indicate ∆′ Haar the von Neumann entropy of the state Ψ R , and H(R|D) Ψ = H(DR) Ψ − H(D) Ψ is the conditional entropy.For simplicity, we omitted (t, β).Note that the mutual information I(D : R) Ψ satisfies I(D : R) Ψ + I(CB ′ : R) Ψ = 2k since |Ψ⟩ DCB ′ R is pure.Hence, the mutual information I(D : R) Ψ in the bounds (Eq.(S16)) can be replaced with 2k − I(R : CB ′ ) Ψ .