Nonequilibrium phase transition of a one dimensional system reaches the absorbing state by two different ways

We study the nonequilibrium phase transitions from the absorbing phase to the active phase for the model of diseases spreading (Susceptible-Infected-Refractory-Susceptible (SIRS)) on a regular one-dimensional lattice. In this model, particles of three species (S, I, and R) on a lattice react as follows: \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$S+I\rightarrow 2I$$\end{document}S+I→2I with probability \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda $$\end{document}λ, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$I\rightarrow R$$\end{document}I→R after infection time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau _I$$\end{document}τI and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$R\rightarrow I$$\end{document}R→I after recovery time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau _R$$\end{document}τR. In the case of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau _R>\tau _I$$\end{document}τR>τI, this model has been found to have two critical thresholds separating the active phase from absorbing phases. The first critical threshold \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{c1}$$\end{document}λc1 corresponds to a low infection probability and the second critical threshold \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{c2}$$\end{document}λc2 corresponds to a high infection probability. At the first critical threshold \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{c1}$$\end{document}λc1, our Monte Carlo simulations of this model suggest the phase transition to be of directed percolation class (DP). However, at the second critical threshold \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{c2}$$\end{document}λc2 we observe that the system becomes so sensitive to initial values conditions which suggest the phase transition to be a discontinuous transition. We confirm this result using order parameter quasistationary probability distribution and finite-size analysis for this model at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{c2}$$\end{document}λc2.


Introduction
Nonequilibrium phase transition from the active states to the absorbing states has been attracted a lot of scientific community efforts recently [1,2,3].One of the most important efforts in this field is concerned by classify nonequilibrium systems into universal classes.In this sense, the directed percolation class (DP) is the most important class in the nonequilibrium phase transition to the absorbing states.Many models have been found to be belongings to this class, for example a contact process (CP), Domany-Kinzel cellular automaton (DK), Ziff-Gulari-Barshad (ZGB) model, paircontact process (PCP), threshold transfer process (TTP) and branching annihilating walks with an odd number of offspring [4,5,6,7,8,9,10].According to Janssen-Grassberger criterion [2,3,11] a model should belong to DP universality class if the model satisfies the following conditions, display a continuous transition into a unique absorbing state with a positive one-component order parameter, with short-range interactions and without quenched disorder or additional symmetries.In fact DP class seems to be even more general and has been found that, some systems belong to this universality class even if they are violate some of previous criterion for example in the longrange interactions [12] or certain models with many absorbing state [13,14,15,16] or fluctuating passive states [17].
Another important class of the nonequilibrium phase transitions to absorbing state is the voter universality class.This class has been observed in special case of models with a two symmetric (Z 2 symmetry) absorbing state [3,18,19,20].Models such as compact directed percolation (CDP), the 2A → ϕ and the 2A → A, the cellular automaton version of the nonequilibrium kinetic Ising model and Lévy-flight anomalous diffusion in annihilating random walks belong to this class.Parity-Conserving universality class (PC) [3,18,19,20] is another universality class to absorbing state.This class characterizes those models which conserve the number of particles modulo 2. Examples are the one-dimensional kinetic Ising models which combined finite temperature spin exchange dynamics and zero temperature spin-flip [21], branching and annihilating random walks with even number of offspring [22] and parity-conserving class of surface-catalytic models [23].
Nonequilibrium discontinuous phase transitions from an active state to an absorbing state have been also observed in the dimensions higher than one in many cases.For example in a two dimensional ZGB model and its modifications [4,6,7,24,25,26,27,28], a two dimensional reactiondiffusion contact-process-like model [29,30,31], a two lattice versions of the second Schlögl model (SSM) [32,33], a two dimensional naming game model [34,35,36], two and four dimensions deterministic conserved threshold transfer process [37,38] and the prisoner's dilemma with semisynchronous updates [39] on two dimension.However, discontinuous phase transition to absorbing states have been rarely seen in one dimension.This is due to the fact that, the fluctuations in low dimensions are strong which make the continuous phase transitions are likely to occur.Hinrichsen [1,2] argued that, first order phase transition can not occur in one dimensional nonequilibrium systems unless there are extra symmetries, conservation laws, long-range interactions or special boundary conditions.By any means in a one dimension, the first order phase transition has been observed in systems with conserved density [38,40], models with long-range interactions [41] and in the systems with multi-component [42,43].For a two-species reaction-diffusion process on a one dimension the renormalization group methods predict a first order phase transition [44,45], however the numerical simulations of that model have been yielded results in disagreement with the renormalization group prediction [46,47].The model candidate to violate Hinrichsen argumentations is the triplet creation model (TCM).This model is a single component and does not possess a conservation law or long-range interactions.Preceding study by Dickman and Tomé [31] had been suggested the first order phase transition of this model for a high value of diffusion rate (D ≥ 0.9).In sequence, Cardoso and Fontanari modified that value to be D ≥ 0.95 [48].Recently, the simulations results of TCM model by Park [49] have been shown that, the phase transition of this model is continuous for any value of D ≤ 0.98.More recent study of this model by Ódor and Dickman [50] suggests a continuous phase transition for any value of D < 1.
In this work we are going to study the phase transition from the absorbing state to the active state of the epidemic spreading model SIRS (Susceptible-Infected-Refractory-Susceptible) on a one dimensional regular network.This model has been proven to has a two critical threshold [51].We are interested to study the phase transition close to those critical thresholds.This work is organized as follows.In section 2, the model and simulations methods are described.Simulation results close to the first critical point of this model are presented and discussed in section 3. Simulation results close to the second critical point of this model are given and discussed in section 4. Conclusions are given in section 5.

Model and Methods
The model of epidemic spreading SIRS on the networks, can be described as follows [52]: consider a population of N particles residing on the sites of a lattice in which each particle is connected to k of its neighbors.The particles can exist in one stage of three stages, susceptible (S), infected (I) and refractory (R).The interaction between the particles on the lattice is as follows: the particles in state I on the network can infect any one of their neighbors which are in state S with probability λ at each time step (S + I → 2I).The particles in state I pass to the R state after an infection time τ I (I → R).The particles in state R return to the S state after a recovery time τ R (R → I).During the R phase, the particles are immune and do not infect.
For this model on the networks as it have been proven in Ref. [51], we have to distinguish between the following two cases: The first case happens when τ I ≥ τ R , where in this case SIRS model has only a one critical threshold λ c separates the active phase from absorbing phase.Second case happens when τ I ≤ τ R , in this case the SIRS model has a two critical threshold λ c1 and λ c2 in which the system is active in between them and die outside of them.In the second case, the first critical threshold λ c1 is corresponding to the situation where the infection probability λ is low.So in this case, the spreading of infection is limited and local therefore and when λ < λ c1 , the system will evolve to the absorbing state where all particles become susceptible (state S).In contrast the second critical threshold λ c2 corresponds to the situation where the infection probability λ is high.Therefore, in this case the infection will spread globally and quickly through the entire network.Now let us ask this question: What will happen for the system when λ is high enough such that, all the particles in the network become infected during a time which is less than or equal to τ I ?.In this case where τ R > τ I , so all the particles will approach the state R followed by S state during a time which is not longer than τ R .Then this case is also an absorbing state for this model.However this absorbing state is un-stationary where the particles will stay in this state only for a time which is not longer than τ R after that the system will approach the stable absorbing state S (see Fig. 1).Hence we can say that, when τ I ≤ τ R the SIRS model has a two absorbing states, even the second absorbing state is not stable but if the system reaches it, the system will evolve surely to the stable absorbing state (first absorbing state S).
As aforementioned, the absorbing state of this model is the state where the lattice becomes free of infected particles, i. e. the state S. SIRS model approaches this absorbing state by two different ways.At λ c1 the model reaches absorbing state due to of that, the strength of infection is very low hence, the average number of susceptible particles infected by an already infected one during the time τ I is less than one.Whereas at λ c2 the strength of infection is high such that, each infected particle infects all of its neighbors during the time τ I .Second critical threshold is equivalent to the state where all particles reach the state I during a time which is less than or equal τ I .In this case and where τ R > τ I then, instantaneously all the particles will approach the state R followed by the absorbing state S during a time which is not longer than τ R .we can consider the state where all particles on the lattice are infected (state I) as an absorbing state of this model however, this state is un-stationary will end up to the absorbing state S. In Fig. 1 we show the space-time evolution for a one dimensional lattice of 11 sites with periodic boundary condition.In this lattice each particle is connected to its first two neighbors.We set the infection probability λ = 1, infection time τ I = 2, recovery time τ R = 3 and the infection starts with one particle on the center of lattice.It is clear from the figure that, all particles on lattice will approach the absorbing state S after 11 time steps.
We simulate this model on a regular one dimensional lattice with periodic boundary conditions in which each particle on the lattice is connected to k = 3 of its nearest neighbors on each side.The system updates synchronically.In this work, we fix the values of infection time and recovery time to be τ I = 7 and τ R = 9 unless we state different.The order parameter ρ(t) is the density of infective particles I (active particles) where N is the total number of lattice sites and ⟨...⟩ stands for average over ensembles.Steady state of order parameter ρ s is the state when ρ s ≡ lim t → ∞ρ(t).In Fig. 2, we recreate the steady state of the density of active particle as function of λ given in Ref. [51].For each point in Fig. 2, we start the simulation from the initial density of active particle ρ(0) = 0.1 and averaged over 100 configurations after discarding 10 4 initial time steps.Figure shows clearly the two critical thresholds which we are interesting to study the phase transition at them.

Phase Translation at the first critical point λ c1
For a general view about the kind of the phase transition at the first critical point, we start the simulation of this model with the typical space-time evolution beside the λ c1 .In Fig. 3 we show the typical space-time evolution of this model when the simulation starts initially from a single active seed located at the center of lattice for the values of the parameters λ = 0.090 and λ = 0.094.Fig. 3 seems to be similar to the typical space-time evolution of the systems which undergo DP phase transition [2].
To confirm if the phase transition in this case is in the DP universality class, we are going to calculate some values of the critical exponents of this model, and before that we will determine the value of critical point λ c1 .In Fig. 4 we plot the average steady state of density of active particles ρ s at various values of the infection probability λ.In the simulation, we use a lattice of size N = 10 4 averaged over 200 realizations after discarding the first 10 4 time steps.Fig. 4 shows clearly that, the system crosses from the absorbing phase to the active phase at specific value of the parameter λ.For the best estimations, the value of the critical point seems to converge to λ c1 = 0.906 ± 0.004.Using this result we can determine one of the critical exponents related to this model where, it is known that, for the continuous phase transitions and as the control parameter λ approaches the critical point λ c , the stationary value of the order parameter ρ s vanishes asymptotically according to a power law as follows [1,2,3]: Inset of Fig. 4 shows the logarithmic plot of ρ s as function of the distance from the critical point (λ − λ c1 ), which shows clearly the power law behavior.Estimated value of the critical exponent β from the inset of Fig. 4 gives us β = 0.281 ± 0.005 which consists very well with the value of β = 0.276 for the (1 + 1) DP universality class [1,2,3].
To extract a further critical exponent, we perform time dependent Monte Carlo simulations of this model starting from a fully occupied lattice.As we know for continuous phase transitions, at the critical point λ c , the time evolution of the order parameter ρ(t) decays asymptotically according to the following power law [1,2,3] Where δ is the critical exponent which equal to 0.159464 (6) for DP universality class in the 1 + 1 dimension [2].Fig. 5 shows the density of active particle ρ(t) as function of time on a logarithmic scale.At the critical point the system clearly shows a power law decay of the active particles.For the best Figure 3: Typical space-time evolutions for λ = 0.090 (right) and λ = 0.094 (left).Simulation starts from a single active particle (black) and time increases downward.
fitting, the value of the critical exponent we find to be δ = 0.159 ± 0.005 which is again consistent with the value of the critical exponent for DP universality class in the 1 + 1 dimension.
Hence we can confirm that, for the case when τ I < τ R the phase transition from the absorbing state to active state for SIRS model at the first critical point λ c1 is of kind DP universality class in (1 + 1) dimension.Here we should mention to that, close to λ c1 , the values of λ are low enough in which the stable absorbing state S is dominated state on the system i.e., it is impossible for the system to approach the un-stationary absorbing state I in this case.Therefore the accessible absorbing state for the system close to the first critical point only the state S, consequently in this case the system satisfies Janssen Grassberger criterion expect this model is multi-component system.

Phase Transition at the second critical point λ c2
As we increase the value of λ toward the second critical point λ c2 we observe that, at a specific value of λ (which is λ > 0.15 for the system of size N = 10 4 other parameters as same as the parameters we have used in the previous section) the steady state of average density of infected particles ρ(t) becomes strongly dependent on its initial values ρ(0) as the Fig. 6 shows.In Fig. 6 we plot the average values of ρ(t) as function of time for various values of initial conditions ρ(0).Fig. 6b, shows three plots of ρ(t) as function of time for the case when λ = 0.16, in which we start simulations from different initial values ρ(0).The figure shows clearly that, for the initial values ρ(0) = 0.1 and ρ(0) = 0.3 the system attains the same steady state, however when ρ(0) = 0.6 the system saturates at a different steady state.Fig. 6c shows the situation for λ = 0.20 with initial conditions same as in Fig. 6b.In this case the system dependence on its initial values becomes more apparent, where the three plots approach different steady states for the three values of initial conditions.For λ = 0.12, Fig. 6a is given for comparison, where in this case the curves reach same steady state and the system does not depend on its initial conditions.Here we should mention to that, in Monte Carlo simulations of Fig. 6, we take the average of ρ(t) over all configurations those survive or not.Another point we mention here is that, the sensitivity of the steady state of the order parameter to initial conditions in the model of disease spreading has been observed in Figure 4: The steady state of density of particles at various values of the λ for the same parameters in Fig. 2. The inset is Log-log plot of the ρ s as function of the distance to the critical point.
By careful inspection of the time evolution of ρ(t) we observe that, there are some configurations quickly become trapped in absorbing state during a time which is not longer than τ I + τ R (in our case it is 17 time steps).This time is on the average is exactly the time need it the particles to go through one infection cycle.Whereas some other configurations survive for a longer time.We also observe that, for fixed values of ρ(0) the density of trapped configurations increases as the value of λ becomes higher.In contrast increasing in the values of ρ(0) causes increasing in the density of trapped configurations for fixed values of λ.Fig. 7 shows the density of trapped configurations (DTCO) as function of density of initial active particles ρ(0) for the value of λ = 0.16.In this figure, the trapped configurations are those configurations reach the absorbing state during the time which is less than or equal 17 time steps.Figure shows clearly that, whenever the initial density ρ(0) of active particles is ρ(0) > 0.4, all the configurations reach the absorbing state.However, when the values of ρ(0) < 0.4, we can see increasing in the number of surviving configurations as the value of ρ(0) decreases.
We can deduce from the previous results that, the chance for the system to approach the state where all the particles in the network are infected (absorbing state I) becomes possible as the value of λ goes toward the higher values.This possibility enhances as the density of initial values are increased.We can understand the relation between the values of infection probability λ and the initial values density ρ(0) of this model from the mean field approximation (see Ref. [56]).Whereas the increasing on the number of I particles at the next time will proportional to the number of particles of kind I and S at this time, i. e. I(t+1) ∝ I(t)+λS(t)I(t), then, if the value of λS(I)I(t) becomes high enough such that I(t + 1) = N the system will reach the un-stationary absorbing state I. Therefore, we can say that, for high values of λ un-stationary absorbing state becomes the dominating state in system.
Because of that dependence for the system on its initial values conditions near λ c2 , we have faced difficulty in determining the kind of the phase transition or even accurately determine the critical point close to λ c2 .We should mention also to that, we could not get any kind of power law behavior near the expected value of the critical point using the time dependent dynamics starting from a fully occupied lattice or form a single active seed.However, and according to Refs.[46,57], the initial values dependence is an indicator of a discontinuous phase transition.Therefore, for better understanding about the mechanism of the phase transition in this case, we perform the order parameter quasistationary probability distribution for this model.As claimed by Refs.[27,29,47,58,59,60] the order parameter quasistationary probability distribution is bimodal in the neighborhood of a discontinuous phase transition in contrast to a continuous phase transition where there will be only a single pick.Fig. 8 shows the results of our Monte Carlo simulations for the cell-occupancy histogram distribution (P) of our model for cells of 100 sites at the center of lattice of N = 10 3 particles, at a various values of λ.The variable n in Fig. 8 is the number of active particles.The quasistationary distribution is clearly bimodal distribution.This result enhances assumption the discontinuous phase transition of this model at λ c2 .For comparison the inset of Fig. 8 shows the cell-occupancy histogram distribution in the vicinity of λ c1 (λ = 0.1), where the system undergoes a continuous phase transition.It is clear in this case there is only a single pick.
The results we find previously suggest strongly that, the phase transition at the second critical point λ c2 is discontinuous.However, for a more clarification we study the quasistationary behavior of this system beside the critical point λ c2 .For that we achieve a finite-size analysis, which is a more reliable procedure as recently proposed in Refs.[33,36].According to that procedure, the difference between the pseudotransition point λ N (where N denotes the system volume) and the transition point λ c2 scales with N −1 according to the relation λ N = λ c2 + aN −1 .Therefore, in order to determine accurately the value of λ N we use the system order parameter variance χ = N ( ρ 2 − ⟨ρ⟩ 2 ).This quantity has been proven to has a peak at the value of pseudotransition point λ N [33,36].We restrict our simulation here, only on surviving configurations.Fig. 9a shows the order parameter variance χ as function of λ.Whereas Fig. 9b shows how the values of pseudotransition point λ L scale with the values of system size N −1 .Extrapolation of N → ∞ yields the critical point for this model to be λ c2 = 0.24 ± 0.01.
Final point we have studied for this model is the space-time evolution of the active particles close to λ c2 .Fig. 10 shows the infected particles with red color during the time of 10 3 time steps for a system of N = 10 3 particles at the value of λ = 0.231 and λ = 0.241.In both figures simulation starts at t = 0 with all particles are in the state S expect for a one active particle I at the center of lattice.Figures clearly show that, the spreading of active particles are compact in a behavior remind us the spreading behavior of CDP models.However, the difference here is that, the particles have a finite time to stay in the active phase.The coexistence of small compact isolated islands of active particles with high regions of inactive ones is again indeed consistent with a discontinuous phase transition.
Finally we mention to that, the models of disease spreading such as a minimal vaccinationepidemic model and Susceptible-Infected-Susceptible (SIS) model have been also found to show either a continuous or a discontinuous active to absorbing phase transition [53,54,55].Additionally, we can deduce some similarities between the phase transition in this model and the phase transition of ZGB model on two dimensional lattice where both models have two critical thresholds.In both models the first critical threshold is corresponding to the continuous DP class and the second critical threshold is corresponding to the discontinuous phase transition.ZGB model has two absorbing states, the first one is at small values of adsorption rate (beside the first critical point) and the second one (beside second critical point) is at high values of adsorption rate [4].SIRS also has a one absorbing state at a low infection rate (beside the first critical point) and unstable absorbing state at a high infection rate (beside second critical point).

Conclusions
In the summary, we have studied the phase transition from the absorbing phase to active phase for the model of infection spearing SIRS on the one dimensional network.This model has been found to has a two critical points where the infection survives in between those critical points and dies out outside of them.The two critical points correspond to low infection rate and high infection rate.Using Monte Carlo simulations we have found that, whereas the the phase transition at the first critical point is of kind the DP universal class, the phase transition at the second critical point is of kind first order phase transition.In this manner, the presence of continuous and discontinuous phase transitions has been also confirmed in the models of disease spreading suchas a minimal vaccination-epidemic model and SIS model [53,54,55].We can also compare the phase transition in this model with the phase transition in ZGB model.Both models have two critical points, where the phase transition at the first critical point is of kind DP class and the phase transition at the second critical point is discontinuous.However, we should mention here to that, the system we have studied here is a one dimensional system whereas ZGB is a two dimensional system.

Figure 2 :
Figure 2: The steady state of the density of particles ρ s at various values of the λ when N = 10 4 , k = 3, τ I = 7 and τ R = 9.

Figure 6 :
Figure 6: The time evolution of average value of density of particles as function of time at different values of initial conditions for N = 10 4 , k = 3, τ I = 7 and τ R = 9: a) for λ = 0.12 b) for λ = 0.16 c) for λ = 0.20.For each curve we averaged over 200 realizations.

Figure 7 :
Figure 7: Density of trapped configurations (DTCO) as function of time at different values of initial values conditions when N = 10 4 , k = 3, τ I = 7, τ R = 9, and λ = 0.16, each point in the figure averaged over 2000 configurations.

Figure 9 :
Figure 9: (a) Order parameter variance χ versus the parameter λ.(b) Value of λ N for which χ is maximum vs N −1