Entropy Governed by the Absorbing State of Directed Percolation

TL;DR

This paper explores the entropy behavior in directed percolation, revealing universal relaxation at the phase transition and a new singularity in the active phase linked to the absorbing state, using tensor network methods.

Contribution

It uncovers a new singularity in the Renyi entropy away from the critical point and links the absorbing state to entropy dynamics in directed percolation.

Findings

Universal relaxation of Renyi entropy at the transition point

Discovery of a new entropy singularity in the active phase

Entanglement entropy behavior linked to the absorbing state

Abstract

We investigate the informational aspect of (1+1)-dimensional directed percolation, a canonical model of a nonequilibrium continuous transition to a phase dominated by a single special state called the "absorbing" state. Using a tensor network scheme, we numerically calculate the time evolution of state probability distribution of directed percolation. We find a universal relaxation of Renyi entropy at the absorbing phase transition point as well as a new singularity in the active phase, slightly but distinctly away from the absorbing transition point. At the new singular point, the second-order Renyi entropy has a clear cusp. There we also detect a singular behavior of "entanglement entropy," defined by regarding the probability distribution as a wave function. The entanglement entropy vanishes below the singular point and stays finite above. We confirm that the absorbing state, though its occurrence is exponentially rare in the active phase, is responsible for these phenomena. This interpretation provides us with a unified understanding of time evolution of the Renyi entropy at the critical point as well as in the active phase.

Equations

$$h_2(t)=H_2(t)/N\sim {\xi }_{\perp }(t{)}^{-1}\sim t^{-1/z},$$

$$P_N(S)={\delta }_{S=0}{P}_N(0)+{\delta }_{S\ne 0}(1-P_N(0))Q_N(S).$$

$$P_N(0)\sim \exp \left[-{\sigma }_0(N/{\xi }_{\perp })\right],$$

$$h_q=\min \left[\left\{q/(q-1)\right\}\left({\sigma }_0/{\xi }_{\perp }\right),h_q^{\prime }\right].$$

$$P_N(0,t)\sim \exp [-{\sigma }_{0}N/{\xi }_{\perp }(t)].$$

$$C_{ijk}$$

$$V_{ij}$$

$$=(1-j)+(2j-1)P[i]$$

$$W_{lmn}$$

$$=(1-n)+(2n-1)P[l+m].$$

$$\sum _{ij}{L}_{ijk}{L}_{ij{k}^{\prime }}^{\ast }={\delta }_{k{k}^{\prime }},\sum _{mn}{R}_{lmn}{R}_{l^{\prime }mn}^{\ast }={\delta }_{l{l}^{\prime }}.$$

Full text

Introduction.— We can define the phase of a macroscopic system out of equilibrium based on physical order parameters of a nonequilibrium steady-state distribution. The phase transition is confirmed theoretically and experimentally in many studiesDerrida et al. (1993); Sasamoto (1999); Domany and Kinzel (1984); Takeuchi et al. (2007, 2009); Sipos and Goldenfeld (2011); Sano and Tamai (2016); Derrida (2007); Henkel et al. (2008). However, we have little understanding of the informational aspect of a macroscopic nonequilibrium system. Recently, Wood et al. Wood et al. (2017) analytically calculated the statistical Rényi entropy of the asymmetric exclusion process Derrida et al. (1993); Sasamoto (1999); Derrida (2007). The phase boundary defined by the behavior of the Rényi entropy agrees with the conventional phase boundary of the asymmetric exclusion process for three nonequilibrium phases. It shows the potential power of informational quantities. However, computation of the entropy is technically difficult; no method for directly measuring entropy is known for the Monte Carlo simulation of a nonequilibrium system. Therefore, the knowledge about the informational aspect of a nonequilibrium system is limited than that of an equilibrium systemBarnett et al. (2013); Iaconis et al. (2013); Stéphan et al. (2014).

In this study, overcoming the technical difficulty by the tensor network method, we will focus on the informational aspect of a $(1+1)$-dimensional directed percolation (DP) through the Rényi entropy and the entanglement entropy. While the Rényi entropy is discussed frequently in conventional statistics, the introduction of the entanglement entropy to the analysis of the classical model is somehow unconventional. Its utility in the context of the classical models was recently discussed Bañuls and Garrahan (2019); Helms et al. (2019). If we regard a percolating direction of objects as a time direction, the $(d+1)$-dimensional DP corresponds to the $d$-dimensional reaction-diffusion process. There are generally two phases in a reaction-diffusion process, the active phase with finite density of objects and the inactive phase with zero density in the long-time limit. As shown in many numerical studies for DP(see a review Henkel et al. (2008)), the transition between the active phase and the inactive phase is critical, and the concept of the universality is also extended as the DP universality class. Because of the extreme simplicity, the universality class is expected to be ubiquitousTakeuchi et al. (2007, 2009); Sipos and Goldenfeld (2011); Vazquez et al. (2011); Sano and Tamai (2016); Sampaio Filho et al. (2018); Deng and Ziff (2018). In the following, after we first briefly introduce the model of DP and the numerical method for the Rényi entropy and the entanglement entropy, we will report the results of both entropies as well as their interpretations.

Directed percolation.— A (1+1)-dimensional DP is defined on a square lattice rotated by 45∘. A site can be either active or inactive. An active state can percolate to the nearest-neighbor sites with a probability $p$. However, the direction of percolation is limited. If it is downward, there are two nearest-neighbor sites from which a site may be percolated. We can regard the row of sites as a one-dimensional system at a time. Then, the DP describes a reaction-diffusion process in a one-dimensional system. This DP is called a bond DP, and it is a special case of Domany-Kinzel (DK) automatonDomany and Kinzel (1984). In the DK automaton, a probability of an active state depends on the number of active states $n$ in the nearest-neighbor sites at the previous time as $P[n]$. If there is no active nearest neighbor just before the current time, the probability of having an active state is zero, i.e., $P[0]=0$. The bond DP is defined as $P[1]=p,P[2]=p(2-p)$. There is a special state in a DK automaton in which there is no active site. It is called the absorbing state because a system cannot escape from it. If $P[1]$ and $P[2]$ are small, the state of the system eventually arrives at the absorbing state after a long time. Then, the evolution of the automaton is stopped. If $P[1]$ and $P[2]$ are large enough, the steady-state distribution with other states can exist in the thermodynamic limit. Thus, there are an active phase and an inactive phase in a DK automaton, and the nonequilibrium phase transition between them is critical. It is called an absorbing phase transition. The behavior at the critical point of DK automatons is universal. It is called the DP universality class.

Numerical method for Rényi and entanglement entropies.— A state probability distribution, $P(S)$, can be mapped to a normalized “wave function” Doi (1976) of which amplitude is proportional to a state probability distribution as $|\psi(t)\rangle={\exp[H_{2}(t)/2]}{\sum_{S}P(S)\ |S\rangle}$. Here, $H_{2}$ is the second-order Rényi entropyRényi (1961). It is a generalization of the Shannon entropy and is defined as $H_{q}=-(q-1)^{-1}\log\sum_{S}[P(S)]^{q}$, where $q$ is called the order. The Shannon entropy is recovered when $q\to 1$.

Recently, various new techniques have been developed in the calculation of a wave function. The most promising approach is the use of a tensor network representation which is a composite tensor defined by tensor contractions. Here, we use a one-dimensional tensor network, i.e., matrix product states (MPS), to represent a wave function of a one-dimensional system. Since the transfer matrix in a master equation is an operator to the wave function which is written as a tensor network(see also App. 1), we can calculate the time step evolution of a wave function as a tensor contraction between themJohnson et al. (2010, 2015); Hotta (2016). Holding a canonical form of MPSShi et al. (2006), we can efficiently calculate an approximated MPS at the next time step by the time-evolving block decimation methodVidal (2003, 2004)(see also App. 2). We can control the precision of MPS representation by a bond dimension of the virtual index between two tensors in MPS(see also App. 3). Since the second-order Rényi entropy appears in a normalization factor of a wave function, we can efficiently calculate it from a tensor contraction of two MPSs.

Through the above mapping between a distribution and a wave function, we can formally introduce the concept of entanglement of a quantum state into a state probability distribution. The entanglement entropy of a quantum state is the Shannon entropy of a quantum subsystem as $H_{E}(A)=-\mathrm{Tr}\ \rho_{A}\log\rho_{A}$, where $\rho_{A}$ is a reduced density operator of a subsystem $A$. The entanglement entropy quantifies a quantum correlation between a subsystem $A$ and its complement $\bar{A}$. In general, we can define the Schmidt decomposition of a quantum state as $|\psi\rangle=\sum_{k}e^{-(1/2)\xi_{k}}|\psi_{A}^{k}\rangle\otimes|\psi_{\bar{A}}^{k}\rangle$, where $|\psi_{A}^{k}\rangle$ and $|\psi_{\bar{A}}^{k}\rangle$ are orthogonal states in subsystem $A$ and its complement $\bar{A}$, respectively. Then, the entanglement entropy is rewritten as $H_{E}(A)=\sum_{k}\xi_{k}\exp(-\xi_{k})$ and $\{\xi_{k}\}$ is called an entanglement spectrum. Since the entanglement spectrum is a key component in a canonical form of MPS, we can efficiently calculate it from the canonical form of MPS.

Critical relaxation of Rényi entropy.— The density of active sites is an order parameter of DP. If we start from a homogeneously full active initial state, the density of active sites shows a power-law decay at the critical point. For example, critical points of the bond DP and the site DP are at $p_{c}=0.644700185(5)$ and $0.70548515(20)$, respectively Jensen (2004). The decay exponent is universal in the DP universality class. Figure 1 shows the critical relaxation of the density and the Rényi entropy per site at the critical points of the bond DP and the site DP of $8192$ sites with open boundary condition. The tensor network with the bond dimension 250 was used for the calculation. Here, to reduce the boundary effect, the density is measured at the central site, whereas the Rényi entropy is defined on the central half of the system. The site DP is a special case in DK automaton as $P[0]=0,P[1]=P[2]=p$, and it also belongs to the DP universality class. Thus, the decay exponents for both DPs should be equal. Our estimates of the decay exponent based on Fig. 1, $\delta_{\textrm{bond}}=0.1593(1)$ for the bond DP and $\delta_{\textrm{site}}=0.158(1)$ for the site DP. They are consistent with a previously estimated value $\delta=0.159464$ Jensen (1999) within the error. Thus, we confirm the reliability of MPS calculation for criticality. The Rényi entropy per site is also an order parameter of the absorbing phase transition because it takes a finite value only in the active phase. As shown in Fig. 1, Rényi entropy at the critical point also shows a power-law decay. Since the decay exponents in Fig. 1 are $0.632(5)$ and $0.625(7)$ for bond DP and site DP, respectively, we confirm the universal relaxation in the Rényi entropy. As we see below, this universal relaxation can be expected as

[TABLE]

where $\xi_{\bot}$ is the spatial correlation length. Indeed, the value of decay exponent in Fig. 1 is consistent with the reciprocal of the dynamical critical exponent, $1/z=0.632613(6)$Jensen (1999).

Rényi entropy and entanglement entropy in active phase.— If the initial probability distribution consists of a fully active state, the Rényi entropy per site starts from zero. In the early time region, it rapidly grows. After that, it quickly converges to a finite value in the active phase. Figure 2 shows the Rényi entropy per site extrapolated to the thermodynamic limit by a $1/N$ quadratic fitting for a system size $N\leq 4096$(see also App. 4). As shown in the inset, there is a cusp at $p_{2}^{*}=0.6785(5)$ in the active phase.

We calculate the time evolution of entanglement entropy between the left and the right halves of the system from the fully active initial state. Figure 3 shows entanglement entropy from the system size 128 to 8192. The bond dimension in MPS is 50111Results in Fig. 3 do not change when the bond dimension is larger than 50.. Solid and dashed lines represent entanglement entropy of $p=0.678$ and $0.679$, respectively. The entanglement entropy has a peak for both cases. The peak time of $p=0.678$ almost converges for more than 4096 sites. However, the peak time of $p=0.679$ diverges with the system size. Solid and dashed lines in the inset represent the levels in the entanglement spectrum from $\xi_{0}$ to $\xi_{2}$ for $p=0.678$ and $0.679$ with 512 sites, respectively. The crossover time at which two levels $\xi_{0}$ and $\xi_{1}$ cross agrees with the peak time of entanglement entropy. Figure 4 shows the crossover time, $\tau_{0}$, from $p=0.676$ to $0.681$ for less than 8192 sites. When we increase the system size, $\tau_{0}$ diverges for $p\geq 0.6786$. However, $\tau_{0}$ converges or its increase becomes weaker for $p\leq 0.6785$. In the thermodynamic limit, the entanglement entropy has a finite peak time for $p<p_{E}^{*}=0.67855(5)$, but no peak for $p>p_{E}^{*}$. The asymptotic behavior of entanglement entropy also changes at $p_{E}^{*}$. In a long time region, the entanglement entropy eventually decreases and converges to zero for $p<p_{E}^{*}$. On the other hand, it converges to a finite value for $p>p_{E}^{*}$. As shown in the inset, the asymptotic behavior of a gap $\Delta_{0}$ between $\xi_{0}(\tau_{0})$ and $\xi_{1}(\tau_{0})$ also changes at $p_{E}^{*}$. The gap is open for $p<p_{E}^{*}$, but it is closed for $p>p_{E}^{*}$. The position of a singularity of entanglement entropy $p_{E}^{*}$ is equal to the position of a singularity of the second-order Rényi entropy $p_{2}^{*}$ within the error. As we see below, these singularities can be explained by the effect of the extremely rare absorbing state.

Effect of absorbing state.— The probability of the absorbing state monotonically increases, because the system cannot escape from it. Therefore, the absorbing state is the only steady state for any finite system. To define a nontrivial steady state, we must consider an infinite system. We then consider a marginal distribution for a subsystem that consists of sequential $N$ sites. In the long time limit, it safely converges to a steady-marginal distribution $P_{N}(S)$. From our calculation, we observe that this steady-state marginal distribution consists of two components, the delta-function at the zero-density state and the continuous component around some nonzero density. Namely,

[TABLE]

The delta-function component can be interpreted as a “quasiabsorbing” state, which is a consequence of the fact that the zero-density state plays a special role even in the marginal distribution because of the absence of the nucleation, i.e., active sites occur only by the influence from the outside of the boundary. Here, since the occurrence of the absorbing state is exponentially rare in the active phase, we assume an asymptotic form of the probability of quasiabsorbing state as

[TABLE]

where $\exp(-\sigma_{0})$ is the probability of having the zero-density state in a correlated region where $\xi_{\bot}$ is the region’s size. In other words, $\sigma_{0}$ is the renormalized chemical potential (of holes).

In the active phase, the Rényi entropy per site of a continuous-part $Q_{N}$ converges to $h_{q}^{\prime}\equiv-(q-1)^{-1}\log\sum_{S}[Q_{N}(S)^{q}]/N$. Using (3), Rényi entropy per site for a steady-marginal distribution in the large $N$ limit is written as

[TABLE]

The singularity at $p_{q}^{*}$ can be understood as the crossover from the first-term dominant region to the second term dominant one, so that $\{q/(q-1)\}\{\sigma_{0}(p_{q}^{*})/\xi_{\bot}(p_{q}^{*})\}=h_{q}^{\prime}(p_{q}^{*})$222Based on (4), we can expect the singular behavior of the Rényi entropy near $p_{c}$ with the exponent $\nu_{\bot}$. The value of the exponent estimated based on all points of Fig. 2 below $p_{2}^{*}$ is $1.093(5)$. Our estimation of the exponent is consistent with the value of $\nu_{\bot}=1.096854(4)$Jensen (1999)..

The entanglement entropy is the Shannon entropy for the singular value part of the amplitude matrix $\mathbb{P}$ as $\mathbb{P}_{S_{L},S_{R}}=\langle S_{L},S_{R}|\psi\rangle$ where $S_{L}$ and $S_{R}$ denote states of left and right half systems. Based on (2), $\mathbb{P}=\mathbb{P}^{0}+\mathbb{Q}$, where ${\mathbb{P}^{0}}_{S_{L},S_{R}}=\exp(Nh_{2}/2)P_{N}(0)\delta_{S_{L}=S_{R}=0}$. From (4), the Frobenius norm333The Frobenius norm of matrix $\mathbb{A}$ is $\sqrt{\textrm{Tr}\mathbb{A}^{\dagger}\mathbb{A}}$. It is the square root of the summation of squared singular values. $|\mathbb{P}^{0}|\gg|\mathbb{Q}|$ for $p<p_{2}^{*}$, and vice versa in the large $N$ limit. The quasiabsorbing state is a direct product state. Therefore, the entanglement entropy is zero for $p<p_{2}^{*}$ and may be finite for $p>p_{2}^{*}$. This result is consistent with the common singular point for the entanglement entropy and the second-order Rényi entropy in the bond DP.

So far, we have been discussing singularity of the entropies due to the crossover from the “absorbing-state dominating” regime to the “continuous-part dominating” regime. While the transitions discussed above are fictitious, the “real” critical phenomena can also be understood by this crossover. To see this, we extend (3) to the probability $P_{N}(0,t)$ allowing the time dependence of the correlation length as

[TABLE]

As shown in Fig. 5, the absorbing-state entropy $-2\log P_{N}(0,t)/N$, equal to the first term of (4), is a decreasing function of time, because $\xi_{\bot}(t)$ monotonically increases. For $p=0.68>p_{2}^{*}=0.6785(5)$, it never completely catches up with $h_{2}(t)$. When $p_{c}<p=0.67<p_{2}^{*}$, it catches up with it at a certain point of time and once it does it dominates the Rényi entropy; $h_{2}$ and the absorbing-state entropy agree with each other afterward. These are the behaviors that are naturally expected from (4) and (5). However, this behavior remains essentially the same even at the critical point $p=p_{c}$, with the only difference being that the asymptotic curve is power-law decay in (1). As shown in Fig. 5, the agreement between the slope of the asymptotic behavior and $z$ confirms the validity of the interpretation with the simple assumption (5).

Conclusion.— We focus on the informational aspect of the (1+1)-dimensional DP through a statistical Rényi entropy and an entanglement entropy. We numerically find a universal relaxation with a slope $1/z$ at the critical point $p_{c}$, and a cusp at $p_{2}^{*}=0.6785(5)$ in the active phase of the bond DP for the second-order Rényi entropy. As we increase $p$ passing $p^{*}_{2}$, the dynamical behavior of the entanglement entropy also changes. The entanglement entropy asymptotically goes to zero below $p^{*}_{2}$ and converges to a finite value above $p^{*}_{2}$. Based on the assumption (5) for the asymptotic form of the probability of quasiabsorbing state, we can understand all of their behavior as the domination of the absorbing state of which the occurrence is exponentially rare in the active phase. Since the existence of the absorbing state is the unique character of DP, the universal relaxation and the cusp in the Rényi entropy are hallmarks of DP(see also App. 5). Our scenario can be extended to the higher-dimensional DP and the DP-like models.

In this study, we show an application of a tensor network scheme to DP as a canonical model of a nonequilibrium continuous transition. The tensor network scheme gives us detailed information of the state probability distribution as the Rényi entropy. We hope for the further use of tensor networks to understand the characterization of nonequilibrium systems from the informational point of view.

Acknowledgment.— K. H. appreciates fruitful comments from T. Okubo and K. Tamai, and K.H. acknowledges many discussions with R. Sato and K. Nataochi about the numerical schemes in the early stage of this study. This work was supported by JSPS KAKENHI Grant No. 17K05576, and by MEXT as Exploratory Challenge on Post-K computer (Frontiers of Basic Science: Challenging the Limits). The computation in this work has been done using the facilities of the Supercomputer Center, the Institute for Solid State Physics, the University of Tokyo and the facilities of the Supercomputer Center at Kyoto University.

Appendix 1 Exact tensor network representation of a transfer matrix of

the Domany-Kinzel automaton

The Domany-Kinzel (DK) automatonDomany and Kinzel (1984) is defined on a square lattice rotated 45 degrees in Fig. 6(a). A state on a site can be either active or not. An active state can propagate to the nearest neighbor sites at the next time. A probability of an active state depends on the number of active states, $n$, in the nearest neighboring sites at the previous time as $P[n]$. Always, $P[0]=0$. If the rotated square lattice of the DK automaton is transformed to a regular lattice from Fig. 6(a) to Fig. 6(b), the DK automaton is a stochastic process on a chain.

Since an active state in the DK automaton is locally propagated, the transfer matrix of DK automaton can be exactly rewritten as a tensor network as shown in Fig. 7. This is a transfer matrix of a two time-step evolution for four sites with two open boundaries which corresponds to Fig. 6(a). The tensor $C$, $V$, and $W$ are defined as

[TABLE]

Here each index is a binary variable, 0 or 1. The composite tensor of $C$ and $W$ enclosed by a dashed line is a two-site operator. At an edge of a chain, it also includes the tensor $V$. Therefore, there are four different composite tensors. The transfer matrix of a two time-step evolution is a cascade of these composite tensors (two-site operators) from left to right and the opposite.

Appendix 2 Canonical form of a matrix product state

A canonical form of a matrix product state(MPS)Shi et al. (2006) consists of three types of tensors as $L$, $R$, and $\Lambda$ shown in Fig. 8(a). $L$ and $R$ are isometries as follows:

[TABLE]

$\Lambda$ is a positive diagonal matrix. Based on a singular value decomposition(SVD)Shi et al. (2006), we can transform an MPS to a canonical form. If we apply a two-site operator to a $\Lambda$-neighboring site in a canonical form(see Fig. 8(b)), we can apply an SVD to a composite tensor from $R^{\prime}$, $R^{\prime\prime}$, $\Lambda$, and a two-site operator to obtain a new approximated MPS. This method is called the time-evolving block decimation methodVidal (2003, 2004). Then, the new MPS again takes a new canonical form in which the position of $\Lambda$ is shifted to right. This SVD approximation means a global minimization of the norm as $||\ |\psi^{\prime}\rangle-T_{i,i+1}|\psi\rangle\ ||^{2}$, where $T_{i,i+1}$ is a two-site operator. The exact tensor network representation of a transfer matrix of the DK automaton is a sequence of two-site operators as shown in Fig. 7. Using a canonical form of MPS, we can keep a total accuracy in the time evolution of the DK automaton. Unlike in Johnson et al. (2010, 2015); Hotta (2016), if we do not use a canonical form, it is hard to calculate the long time behavior of DK automaton with high accuracy.

Appendix 3 Time evolution of a density in the bond DP

In our numerical simulation, we use a canonical form of MPS which approximates a state probability distribution. The bond dimension of a canonical form should be large enough for a high precision. Fig. 9 shows the time evolution of a density in the bond DP by using a canonical form with a bond dimension $D=120$. Here, in order to reduce the effect of open boundaries, the density is defined as the average density of two sites at the center of a system. The system size is $N=512$ except for $N=256$ at $p=0.7$. Results of the canonical form are consistent with those of Monte Carlo simulations for the same system size. There is no visible dependence on the finite bond dimension and the finite system size for these results.

Appendix 4 Rényi entropy of the steady-state distribution of the

bond DP in the thermodynamic limit

In order to estimate a Rényi entropy, we need to check the dependence on the system size and the bond dimension $D$ of MPS. Fig. 10 shows the time evolution of the Rényi entropy of the bond DP by $D=50$, and $100$ for $N=4096$. There is no visible difference between $D=50$ and $D=100$.

The Rényi entropy quickly converges to the value of the steady-state distribution before $t=2000$. Fig. 11 shows the values of Rényi entropy near $p_{2}^{*}=0.6785(5)$ at $t=2000$. There is the strong finite-size effect of a system size $N$. Using $1/N$ quadratic fitting curves, we estimated the value of the Rényi entropy in the thermodynamic limit. Fig. 12 shows the Rényi entropy in the thermodynamic limit. The error shown in Fig. 12 (the inset of Fig.2 of the main text) is the extrapolation error with respect to the system size. The error due to the finiteness of $D$ is much smaller. As shown in Fig. 12, there is a cusp slightly above the bond-DP threshold, $p_{c}=0.644700185(5)$ Jensen (2004). The position of the cusp is between $p=0.678$ and $0.679$. Therefore, we read the error in our estimate of $p_{2}^{*}$ is smaller than $0.0005$ as $p_{2}^{*}=0.6785(5)$.

Appendix 5 Rényi entropy of the steady-state

distribution of the site DP in the thermodynamic limit

Using the same procedure in the bond-DP case, we estimated the Rényi entropy of the site DP in the thermodynamic limit from the data of $N\leq 1024$. As shown in Fig. 13, there is a similar cusp in the Rényi entropy of the site DP in the thermodynamic limit. The position of cusp is at $p_{2}^{*}=0.756(1)$ slightly above the site-DP threshold, $p_{c}=0.70548515(20)$ Jensen (2004).