Non-unitary dynamics of Sachdev-Ye-Kitaev chain

We construct a series of one-dimensional non-unitary dynamics consisting of both unitary and imaginary evolutions based on the Sachdev-Ye-Kitaev model. Starting from a short-range entangled state, we analyze the entanglement dynamics using the path integral formalism in the large $N$ limit. Among all the results that we obtain, two of them are particularly interesting: (1) By varying the strength of the imaginary evolution, the interacting model exhibits a first order phase transition from the highly entangled volume law phase to an area law phase; (2) The one-dimensional free fermion model displays an extensive critical regime with emergent two-dimensional conformal symmetry.

Recent years have witnessed tremendous breakthrough in many-body quantum dynamics. For a closed many-body quantum system decoupled from the environment, under the unitary dynamics, the interaction in the system can lead to chaos and thermalize all the small subsystems. The total wave function acts as its own heat bath and this phenomenon is referred as quantum thermalization [1,2].
The irreversible thermalization process can be avoided if we allow non-unitary evolution, which naturally arises in open quantum systems. Recently it is observed that a unitary dynamics subjected to repeated measurement can H p i < l a t e x i t s h a 1 _ b a s e 6 4 = " V v A t g 3 A O E k z Y x a i 1 K h x 4 a R S u G A k = " > A A A B 7 H i c b V B N S 8 N A E J 3 U r 1 q / q h 6 9 L B b B U 0 m k o N 6 K X n q s Y G q h j W W z n b R L N 5 u w u x F K 6 W / w 4 k E R r / 4 g b / 4 b t 2 0 O 2 v p g 4 P H e D D P z w l R w b V z 3

y m s r W 9 s b h W 3 S z u 7 e / s H 5 c O j l k 4 y x d B n i U h U O 6 Q a B Z f o G 4 E t l O F N A 4 F P o S j 5 n / 8 I R K 8 0 T e m 3 G K Q U w H k k e c U W M l v 9 H j j m v X H G r 7 h x k l X g 5 q U C O Z q / 8 1 e 0 n L I t R G i a o 1 h 3 P T U 0 w o c p w J n B a 6 m Y a U 8 p G d I A d S y W N U Q e T + b F T c m a V P o k S Z U s a M l d / T 0 x o r P U 4 D m 1 n T M 1 Q L 3 s z 8 T + v k 5 n o K p h w m W Y G J V s s i j J B T E J m n 5 M + V 8 i M G F t C m e L V s K G V F F m b D 4 l G 4 K 3 / P I q a V 1 U v V r 1 + q 5 W q d / k c R T h B E 7 h H D y 4 h D o 0 o A k + M O D w D K / w 5 k j n x X l 3 P h a t B S e f O Y Y / c D 5 / A J / n j p Y = < / l a t e x i t >
H q j,j+1 < l a t e x i t s h a 1 _ b a s e 6

= " q z c e T K V O z S 7 b 3 0 h C d R B 5 F h Y Z U l 8 = " > A A A B 8 n i c b V B N S 8 N A E N 3 U r 1 q / q h 6 9 L B Z B U E o i B f V W 9 N J j B f s B a S y b 7 a b d d p O N u x O h h P M L x U 8 e q v 8 e a / c d v m o K 0 P B h 7 v z T A z z 8 F 1 2 D b 3 1 Z u Z X V t f S O / W d j a 3 t n d K + f N L V M F G U N K o V U b Z 9 o J n j E G s B B s H a s G
A l 9 w V r + 6 H b q t 5 6 Y 0 l x G 9 z C O m R e S f s Q D T g k Y y a 1 1 0 + H 5 8 M y Z P D x 2 i y W 7 b M + A l 4 m T k R L K U O 8 W v z o 9 S Z O Q R U A F 0 d p 1 7 B i 8 l C j g V L B J o Z N o F h M 6 I n 3 m G h q R k G k v n Z 0 8 w S d G 6 e F A K l M R 4 J n 6 e y I l o d b j 0 D e d I Y G B X v S m 4 n + e m 0 B w 5 a U 8 i h N g E Z 0 v C h K B Q e L p / 7 j H F a M g x o Y Q q r i 5 F d M B U Y S C S a l g Q n A W X 1 4 m z Y u y U y l f 3 1 V K 1 Z s s j j w 6 Q s f o F D n o E l V R D d V R A 1 E k 0 T N 6 R W 8 W W C / W u / U x b 8 1 Z 2 c w h + g P r 8 w d 2 o p C + < / l a t e x i t > represents the inter-cluster SYK coupling term. We compute Rényi entropy of a single interval in the box.
exhibit non-thermal phases if we follow the quantum trajectory of the many-body wave function. More strikingly, by varying the measurement rate, there is a continuous entanglement phase transition [3][4][5][6][7][8][9][10][11][12][13][14][15][16][17]. In the phase with slow measurement rate, the state remains highly entangled and the entanglement entropy obeys volume law scaling, while in the phase with fast measurement rate, the entanglement entropy obeys area law scaling. Notice that the explicit form of the measurement is not important and can be either a projective measurement or a more generalized weak measurement [5,15].
Motivated by these findings, we consider the following question: For a non-unitary dynamics U ∼ exp(−iHt) governed by a non-Hermitian Hamiltonian can we realize an entanglement phase transition by varying g? In the above equation, both H 1 and H 2 are Hermitian Hamiltonians. More specifically, in this paper, we consider the following interaction: H 1 is a Hamiltonian describing interaction between different sites and H 2 is a Hamiltonian defined at each site. H 2 for example can describe the coupling of the onsite degrees of freedom to an external field. For such non-unitary dynamics, in the limit g = 0, we expect that the steady state will typically saturate to a highly entangled state with a volume law scaling, while in the limit g → ∞, this becomes a purely imaginary evolution and the steady state is a trivial product state with zero entanglement entropy. In a strongly interacting system, it is not obvious if there is a phase transition occurring at finite g. To address the above questions, we consider a one-dimensional (1D) non-unitary dynamics constructed from Sachdev-Ye-Kitaev (SYK) model [18,19] and explore the possible phase transition in it.
The SYK model is a fermionic system with random all-to-all interaction [18][19][20]. This model can be analytically solved in the large N limit. Due to this special property, this model has vast applications in different fields including high energy, condensed matter physics and quantum information theory. Many variants of the SYK model have been constructed to study quantum chaos, quantum gravity and non-Fermi liquid analytically [21][22][23][24][25][26][27][28][29][30]. In particular, there are studies on the entanglement entropy of the SYK model [31][32][33][34][35][36][37], where transitions to the replica wormhole [38,39] solution are found in the long time limit.
In this paper, we will use the SYK model to construct a set of 1D chain models and explore the non-unitary dynamics harbored in them. We study the entanglement dynamics by using the path integral formalism, which can be obtained self-consistently by virtue of the large N nature of the SYK model. By varying g, we observe different entanglement scaling behaviors which correspond to different saddle point solutions. We hope these results in the large N limit could shed light on the more generalized phase transition in interacting systems at finite N where analytical tools are lacking.
The rest of the paper is organized as follows. In Sec. 2 we define the SYK chain models, derive the path integral formalism for the Rényi entropy of these models, and write down the corresponding saddle point equations which are amenable to numerical study. Then in Sec. 3 we apply the formalism to two models: the interacting model with inter-cluster SYK 4 coupling and intra-cluster SYK 2 coupling, and the non-interacting model with inter-cluster SYK 2 coupling and intra-cluster SYK 2 coupling. We find that the steady state of the interacting model exhibits either a volume-law or an area-law phase as the coupling λ is varied, and the two phases are separated by a first order transition, while that of the non-interacting model exhibits critical behavior for any finite coupling λ. We finally summarize and discuss these results in Sec. 4. Some derivation details, as well as results for another model with both inter-and intra-cluster SYK 4 interactions, are given in the Appendices.

MODEL AND METHOD
We consider non-unitary time evolution of SYK chain generated by the following non-Hermitian Hamiltonian: where H p inter-cluster is the inter-site interaction with coupling strength J and H q intra-cluster is the onsite interaction with coupling strength V . From now on, we drop their subscript for conciseness. At each site, there are N Majorana fermions. H p x,x+1 denotes p-body random interaction between neighboring sites x and x + 1 and H q x is q-body random interaction at site x (See the cartoon in Fig. 1). Their explicit forms are where the subscript i m , j n = 1, 2, . . . , N is the Majorana flavor index at each site and the superscript x = 1, 2, ..., L is the site of the SYK chain; we assume periodic boundary condition for the chain so that L + 1 ≡ 1. J x,x+1 i1···i p/2 j1···j p/2 and V x i1i2···iq are time-independent Gaussian random variables with vanishing mean and the following variance Note that both H p x,x+1 and H q x are Hermitian Hamiltonian. Therefore in the time evolution governed by exp(−iHt), H p x,x+1 provides the unitary evolution while H q x is responsible for the imaginary evolution. We describe the path integral formalism for this non-unitary evolution in the next subsection.

Rényi entropy under non-unitary evolution
We are interested in the entanglement dynamics of certain pure state. We prepare an initial state |ψ(t = 0) = ⊗ L x=1 |{0} x which is the state annihilated by all the complex fermions defined by c x,j = χ x,2j−1 +iχ x,2j as c x,j |{0} x = 0, for all x ∈ 1, 2, ..., L and j = 1, 2, ..., N . Since our initial state is a product state between different sites, there is no spatial entanglement.
We then evolve |ψ(t = 0) under the Hamiltonian (2) to entangle different sites. At time t, the state reads where is introduced to ensure normalization of the state |ψ(t) . From now on, we keep the dependence of random parameters implicit. It is useful to write Z(t)|ψ(t) and Z(t) as path-integrals: Here we introduce f (τ ) to specify whether the real-time evolution is forward or backward: f (τ ) = i for forward evolutions e −iH p t and f (τ ) = −i for backward evolutions e iH p t 1 . We use b.c. to denote proper boundary conditions specified by the pictorial representations: where we have separated out the Majorana modes with even/odd indices as χ x e/o and x ∈ [1, L]. The solid lines represent the evolution and the black dots denote the initial state. The dotted lines represent interaction between fermions, which contains coupling between different sites. For Z(t)|ψ(t) , additional quantum state is attached to the free ends.
We are interested in the Rényi entropy of |ψ(t) . We divide the SYK chain into subsystems A and B, with A containing sites x = 1, 2, ..., L A . The reduced density matrix ρ A = Tr B ρ(t) is obtained by tracing over the degrees of freedom in B. A pictorial representation of ρ A reads: Here the red line represents the A subsystem and the blue line represents the B subsystem. The nth Rényi entropy is then defined as where Z n (t, L A ) is given by sewing n copies of reduced density matrix ρ A , as shown in Fig. 2, where we have labeled contour C anti-clockwise by a single time variable s ∈ [0, 4nt] (See Appendix A for details). The contour itself again defines the path-integral for computing the Rényi entropy, with boundary conditions indicated by dashed lines and black dots. From the contour, it is obvious that there is a symmetry by interchanging A and B: S (n) (t, L A ) = S (n) (t, L − L A ). In addition, we have Z n (t, 0) = Z n (t, L) = Z n (t), which is consistent with the fact that the full system is in a pure state.
To proceed, we need to average over random variables J x,x+1 i1...i p/2 j1...j p/2 and V x i1...iq after computing log Tr A ρ n A for each realization of J x,x+1 i1...i p/2 j1...j p/2 and V x i1...iq . Additional disorder replicas are required to accomplish this. In the SYK models that we consider, however, the computation can be simplified by working within a disorder replica diagonal ansatz: both numerical analysis and analytical arguments [40][41][42][43] suggest that at leading order of 1/N , the Rényi entropy can be approximated by and the result becomes exact in the large N limit. After disorder average, to compute Z n (t) one could further introduce the bilocal fields G x and Σ x [20,36]. Leaving details for Appendix A, we find As aforementioned, the time parameters s, s ∈ [0, 4nt] are parameterized on the contour C. P (s, s ) is defined as P (s, s ) = 1, s, s both parameterizing even or odd fields, 0, s parameterizing odd (even) fields and s parameterizing even (odd) fields. In the large N limit, since the action (14) is proportional to N , we obtain the saddle point solution Therefore the Rényi entropy is where S n,saddle is the on-shell action (18) obtained from the solution of the Eqs. (16). Explicitly, we have

Brownian dynamics
We can also generalize the above formalism to a non-unitary Brownian SYK model in which J x,x+1 and/or V x are/is independent Gaussian random variables in time, with vanishing mean and the following variance i.e. the correlation is non-vanishing only at equal-time. A variant of this model was studied in Ref. [44]. Previously, there has been extensive study on the Brownian unitary circuit models [44][45][46][47][48][49][50]. These models can provide analytical solutions for the quantum dynamics even with a small onsite Hilbert space and can sometimes show different dynamics from systems without randomness [51]. Here we plan to explore the influence of the temporal correlation on the nonunitary quantum dynamics.
One can apply the entire procedure of the last subsection to computing the Rényi entropy of the Brownian SYK model with minimal change. Note that in the Brownian model, a time-dependent random variable admits only equaltime correlation t = t , which happens at the redefined times s = 2jt + s and s = 2(j + 1)t − s on the contour for Z n (t). Therefore we need to introduce an additional projection to the Brownian term As an example, when the inter-site interaction becomes Brownian, the saddle point equation (16b) becomes Similar modification applies when V x becomes Brownian. The on-shell action can still be computed as in (18).

Numerical details
The complex form of the saddle point equations (16) hampers analytical treatment for generic parameters of J, V , t, and L A , and numerical methods are exploited instead. We will hereby focus on the second Rényi entropy with n = 2. We discretize s ∈ [0, 8t) to 4L points, and represent G x and Σ x as 4L × 4L matrices. Numerically, for moderate J and V , we take L = O(1) × t and finite size scaling for L is performed as the final procedure to obtain the L → ∞ result.
One solves the numerical version of Eqs. (16) iteratively: at the nth step, one plugs the Green's function obtained from the previous step, G x } converge within numerical precision. The numerical version of (16) is where ∆t = 2t/L. P and f are the straight-forward discretization of their definitions (see Eq. (15) and the paragraph below Eq. (8)): where we have pictorially presented the 4L × 4L matrix of P , with the gray (white) representing unit (vanishing) entries . The partial operator ∂ τ in the continuum version (16) is now discretized in Eqs. (22) to the non-interacting Green's function (G 0 where sgn(x) is the sign function with sgn(0) = 0, the row and column indices are understood in the modulo sense i + 4L ≡ i, and the unmentioned entries are identically zero. Note this substitution ∂ τ → (G 0 x ) −1 is crucial to ensure numerically accurate result with small number of points. We introduce the following terminology: a solution is replica diagonal if for all x, the Green's function G 0 x is non-zero only for elements that are non-zero in G 0 x 2 , and a solution is replica quasi-diagonal if this property is satisfied for all x in the subsystem bulks (i.e. interior)Å,B and is violated on the subsystem boundaries ∂A and ∂B.
We consider two different types of initial conditions:  By definition, the initial condition Type 2(D) is replica diagonal, and the initial condition Type 2(QA/QB) is replica quasi-diagonal. When q ≥ 4, the replica diagonal/quasi-diagonal property is preserved under the iteration and the initial condition Type 2(D) (Type 2(QA/QB)) always leads to replica diagonal (quasi-diagonal) self-consistent solution. We have also checked that other initial conditions would not lead to new saddle points.
After solving Eqs. (22), we compute S 2,saddle (L A ) as a function of the subsystem size L A using the discretized version of Eq. (18). The final result is obtained by using Eq. (17)

PHYSICAL PROPERTIES OF THE MODEL
In this section we present numerical results along with analytical understandings for our model with different (p, q). We focus on models with (p, q) = (4, 2) and (2,2) in the next two subsections. For the (4, 2) model, we find that the steady state is a volume-law entanglement state for V J and is an area-law entanglement state for J V . These two phases are separated by a first order transition. On the other hand, for the (2, 2) model, we find the steady state is always in a critical phase for finite V /J. We will briefly mention properties of the (p, q) = (4, 4) model in Appendix C.

Inter-cluster SYK4 with intra-cluster SYK2
The (p, q) = (4, 2) model consists of Hermitian, inter-cluster SYK 4 interaction with strength J and non-Hermitian, intra-cluster SYK 2 interaction with strength V . We prepare a product state and let it evolve under this Hamiltonian described by Eq. (6). We study the scaling of the second Rényi entropy S (2) of the state as a function of J, V , subsystem size L A , and evolution time t.

V J
Before discussing the non-unitary dynamics, we first consider the unitary dynamics at V = 0. Since the system has local interaction, we expect the entanglement entropy for a subsystem grows linearly in time and saturates to a constant proportional to L A at late time. We compute S (2) by numerically solving the saddle point solution and present the results in in Fig. 3(a). In summary, the entanglement dynamics satisfies where the entanglement velocity v E depends on the magnitude of J and the transition occurs at t = s 0 L A /v E . Since the initial state we consider has energy expectation value E = 0, the final steady state is indeed a maximally entangled thermal state at infinite temperature with s 0 = log 2/2. Notice that the entanglement dynamics in the large N model is distinct from that in a chaotic system with finite or small N , in which the energy conservation law leads to a diffusive dynamics for S (2) (t), i.e., S (2) (t) ∼ √ t [52,53]. This slow dynamics disappears in the SYK chain model in the large N limit.
The above behavior still holds when V takes a small finite value, although with s 0 < log 2/2. For instance, we present the V = 0.3 result in Fig. 3(b), which is analogous to the V = 0 result. This indicates that when V J, there exists a phase in which the steady state entanglement entropy has volume law scaling.
Here we provide an explain for this entanglement scaling from the saddle point solution perspective. The S (2) in Fig. 3(a) and Fig. 3 Fig. 4(a), starting from an initial condition of Type 1, S (2) quickly saturates to a volume law scaling, including a branch with S (2) ∝ L A (Type 1(A)) and a branch with S (2) ∝ L B (Type 1(B)). On the other hand, initial condition Type 2 always leads to a constant S (2) that is insensitive to subsystem size. Furthermore, for small V 0.3, we observe that the replica quasi-diagonal saddle obtained using Type 2(QA/QB) and the replica diagonal saddle obtained using Type 2(D) give slightly different results. Both solutions give S (2) that grows linearly in time, however the replica quasi-diagonal saddle has a lower S (2) value (see Fig. 5(a)). Interestingly, as V increases, the replica quasi-diagonal saddle gradually merges with the replica diagonal saddle, and as V > 0.3, both initial conditions Type 2(Q) and Type 2(D) give the same replica diagonal solution. We will henceforth refer to the saddle point solutions by the initial conditions that lead to them.
The above analysis suggests that early time S (2) is determined by the Type 2 saddle which grows linearly in time, and that late time S (2) is determined by the Type 1 saddle which eventually saturates to a constant that is proportional to L A . The time dynamics is presented in Fig. 3(c) and Fig. 3(d) for V = 0 and V = 0.3 respectively and is consistent with the result in Eq. (26). In addition, in Fig. 3(d), we verify that no finite size effect exists in both saddle solutions.
We can further understand the Type 1 and Type 2 saddles by analyzing the Green's functions. Here we take V = 0.2 as an example. The S (2) computed from different saddles are plotted in Fig. 4(a), and the subsystem bulk Green's functions for different saddles are shown in Figs. 4(b)-4(d). The Green's function corresponding to Type 2 saddles has a replica diagonal form (Fig. 4(b)) inside the subsystem bulk. This means there exists self-consistent replica diagonal/quasi-diagonal solution of the Green's functions when computing Z n (t). Compared to the Green's functions used for computing the normalization Z(t) (with proper time labeling), the Green's functions for Z n (t) are modified only in the neighborhood of subsystem boundaries, since only near the boundaries is there a change of the Schwinger-Dyson equation. Consequently, terms in Z n (t) and nZ(t) that involve only Green's functions in the subsystem bulkÅ andB cancel, and the Rényi entropy receives contribution only from the Green's functions on the subsystem boundaries ∂A and ∂B. This leads to a S (2) which is independent of the subsystem size. On the contrary, if we take a replica non-diagonal saddle as in Fig. 4(c) and Fig. 4(d), the Green's functions become highly replica non-diagonal in one of the subsystems (A or B, depending on the initial conditions being Type 1(A) or Type 1(B)). Consequently, there is non-vanishing contribution to the entropy even deeply inside the subsystem bulk. This may give rise to a S (2) proportional to the subsystem size. We mention in passing that the Green's functions also distinguish replica diagonal and replica quasi-diagonal saddles due to their difference at the subsystem boundaries.
As an example, we show in Fig. 6  in Fig. 4(b). The boundary Green's function profiles confirm that as V increases, the replica quasi-diagonal solution converges to the replica diagonal solution. Interestingly, the simultaneous existence of such replica non-diagonal and replica diagonal/quasi-diagonal solutions is directly related to the recent resolution of the black hole information paradox. In the language of gravity theories, such highly replica non-diagonal solution is known as a "replica wormhole" [38,39]. In the story of the information paradox, without such replica wormhole, the entropy of the system would grow unboundedly, which ultimately violates unitarity. We find similar behavior exists in our model: without finding the Type 1 saddles, we would conclude that the entropy of our system grows linearly, finally exceeding the maximal entropy L A log 2/2. When the "replica wormhole" saddle point is taken into account, the entropy becomes well-behaved.

V J
Now we turn to the other limit J/V → 0. First, the J = 0 case defines purely imaginary evolution governed by the intra-cluster SYK 2 interaction. The system evolves into its ground state after long time evolution. Since the SYK 2 model is a free fermion model and has vanishing zero temperature entropy [20], the ground state is a simple product state, i.e., where |ψ x GS denotes the SYK 2 ground state at site x and there is no entanglement between different sites. When J = 0, we expect that the steady state entropy S (2) still satisfies area law and confirm this result in Fig. 7(a). This scaling behavior can be understood as follows: As long as J V , we can perturbatively study S (2) by taking the Type 2(D) saddle solution. The leading behavior of S (2) receives corrections only from the subsystem boundaries and can be written in closed form (see Appendix B for derivation): where G c is the conformally invariant Green's function for the ground state of SYK 2 valid for s − s ≈ V −13 . Plugging G c into Eq. (28), we get where the prefactor of the logarithmic term b = J 2 3π 2 V 2 is cutoff-independent, while a, c are cutoff-dependent finite constant. The first term of S (2) is responsible for the area law scaling behavior of the steady state and the J 2 /V 2 scaling form is confirmed numerically in Fig. 5(b). The second term contains the leading time dependence at small J/V which shows a decaying behavior that is linear in log t/t 2 , see the right inset of Fig. 5(b). We also verify the explicit form of the coefficient b, see the left inset of Fig. 5(b).

The transition at the intermediate regime
The above analysis suggests that there exists a volume law phase when V J and an area law phase when V J. We further explore the intermediate regime and find that there is a phase transition separating these two phases. In Fig. 7, we analyze the Type 2 saddle solution for S (2) as a function of time and we find that it changes dramatically across V /J = 0.39. When V /J = 0.38, it grows linearly in time while when V /J = 0.39 or 0.4, it decays with time and saturates to a constant. This result indicates that the transition is first order and occurs close to 0.39. In the volume law phase, the steady state entanglement entropy is determined by the Type 1 saddle solution, which scales linearly with the subsystem size up to L/2. The slope decreases as we increase V /J. As we enter the area law phase, S (2) is bounded by the Type 2 saddle solution which decays with time and saturates to a finite constant.
We further analyze the volume law phase close to the transition point, V /J 0.39, and find that the behavior of S (2) is slightly different from that of the V ≤ 0.3 region at early time. We take V /J = 0.36 as an example, see Fig. 8(a). When the subsystem size is small, S (2) ∼ αL A , with the coefficient α = α(t) a decreasing function in time. This decaying behavior originates from the Type 1 saddle solution and is responsible for the spike observed in Fig. 8(b). At late time, this coefficient saturates to a finite constant and the steady state still exhibits volume law scaling.
We summarize the main result for the (4, 2) model in the phase diagram shown in Fig. 9. Three distinct early time behaviors of S (2) (L A ) (dashed and dotted lines) as well as two distinct late time behaviors (red lines) are sketched.  correlation relation (19).
• The first order transition occurs at V /J ∼ 0.07 for the (4 B , 2) model, while it occurs at V /J 0.39 for the (4,2) model.
Unlike the (4 B , 2) model, the (4 B , 2 B ) model does not exhibit phase transitions. The entanglement dynamics in the entire parameter range of V /J is described by Eq. (26), and the steady state is a volume-law entanglement state 4 .

Inter-cluster SYK2 with intra-cluster SYK2
Another variant of the model is the (p, q) = (2, 2) model. This free fermion model exhibits critical phases distinct from the previous interacting model. We again first consider the limit V = 0. The entropy S (2) as a function of subsystem size and time is shown in Figs. 10(a) and 10(b). We observe that S (2) at early time increases linearly with √ t, and at late time exponentially saturates to a steady state value (see the inset of Fig. 10(b)). The subsystem scaling of this steady state value agrees well with the following expression (see Fig. 10(a)) where S SYK2 (λ) is the eigenstate Rényi entanglement entropy for a subsystem with N λ Majorana modes of a singlecluster SYK 2 model with N Majorana modes [31,34,54,55]: This result is much smaller than the (log 2/2)L A scaling in the (4, 2) model with V = 0. The entanglement entropy takes a different form when V = 0. When V = J = 1, as shown in Fig. 10(c), the steady state entropy has the form which is the same as that for the ground state of a 1D critical system with periodic boundary condition. Furthermore, starting from a product state, S (2) increases as at early time with the same coefficient α. Such scaling behavior implies an emergent two-dimensional conformal symmetry with dynamical exponent z = 1. Similar behavior has also been found in the random non-unitary free fermion dynamics in the small N limit [56]. We then tune away from the ratio V /J = 1 and find that the critical scaling behavior of S (2) is retained in a large parameter range. More interestingly, we observe that the coefficient α is linearly proportional to the ratio J/V (see inset of Fig. 10(c)), indicating that this model remains critical as long as V /J is finite.
Here we take periodic boundary condition and partition the system into four connected intervals with end points x 1,2,3,4 . We present I (2) in Fig. 10(d) and find that it is only a function of cross-ratio, which is defined as This numerical result provides further evidence that this model has conformal symmetry. In particular, we find that when η → 0, I (2) ∼ η, the same as that for non-unitary random free fermion dynamics in the small N limit [56]. The very short time behavior of S (2) in the limit J V can be treated perturbatively, see Appendix C.

DISCUSSION AND CONCLUSION
In this work we study from the entanglement perspective the non-unitary dynamics of a 1D SYK chain H = J x H p x,x+1 − iV x H q x . We derive the large N self-consistent saddle point solution to the nth Rényi entropy using the path integral method, and numerically study the second Rényi entropy as a function of time, subsystem size, and the dimensionless coupling strength. We find that for the (p, q) = (4, 2) model, as we vary the ratio V /J, the steady state can be either in a volume-law phase or an area-law phase. The steady states of these two phases correspond to distinct saddle point solutions, and the two phases are separated by a first order transition. For the non-interacting (p, q) = (2, 2) model, the steady state exhibits critical behavior for any finite ratio V /J, indicating the emergent two-dimensional conformal symmetry in this model. There are a lot of possible interesting extensions of our work, which we briefly mention. Firstly, in the (4, 2) model, the transition is first order when inter-cluster coupling is described by the regular SYK or Brownian SYK interaction. It would be interesting to study the 1/N correction in this model and check if it can round this first order transition to a second order transition. Secondly, the emergent conformal symmetry observed in (2, 2) model may need further investigation. We plan to analytically study the critical saddle point in this large N model and find the connection with the critical phenomena observed in the small N model. Lastly, we can use the method developed in this paper to construct other non-unitary dynamics and explore the exotic phases in them.