Reaction fronts in stochastic exclusion models with three-site interactions

The microscopic structure and movement of reaction fronts in reaction–diffusion systems far from equilibrium are investigated. We show that some three-site interaction models exhibit exact diffusive shock measures, i.e. domains of different densities connected by a sharp wall without correlations. In all cases, fluctuating domains grow at the expense of ordered domains; the absence of growth is possible between ordered domains. It is shown that these models give rise to aspects not seen in nearest-neighbour models, namely double shocks and additional symmetries. A classification of the systems by their symmetries is given and the link of domain wall motion and a free fermion description is discussed.


Introduction
The emergence of patterns and fronts is a challenging problem in biology, chemistry and physics; for a review see [1]. In biology, for example, bacteria aggregate, building up regions with a high density in coexistence with regions with a low density of organisms [2]. A typical example in chemistry is the movement of reaction fronts. In physics, the movement of domain walls is directly related to the problem of coarsening, for example in magnetic systems. In [3] the phenomenon of hysteresis in a driven diffusive system is explained by the movement of shocks, i.e. a jump in the density profile. Various other phenomena in many-particle systems can be attributed to the emergence of shocks, for example, the first-order transition in the phase diagram of the asymmetric exclusion process (ASEP) [4] or phase coexistence in a driven diffusive system coupled to reaction kinetics [5,6]. Recently, shocks in quantum systems have also been discussed [7,8].
On a macroscopic level shock fronts are described by partial differential equations. The most prominent equations are the Fisher and the Burgers' equations. The Fisher equation [9] ∂ ∂t ρ(x, t) = D ∂ 2 ∂x 2 ρ(x, t) + aρ(x, t) − bρ 2 (x, t) (1) was originally proposed as a model for the propagation of a mutant gene. It shows travelling-wave solutions and may be used for modelling systems without conservation of the order parameter.
In this paper we demonstrate for some models how these macroscopic shocks originate from the microscopic dynamics. It is known that some driven diffusive systems can be described by the Fisher or Burgers' equation in the hydrodynamic limit. This limit is achieved by scaling the lattice constant to zero while keeping the overall length of the system constant; the time has to be rescaled appropriately. One of these models is the ASEP which is in the hydrodynamic limit described by the Burgers' equation. Another model is the branching and coalescing random walk, which is in the hydrodynamic limit described by the Fisher equation [11]. In these two models, the microscopic dynamics can be described by exact shock measures. An exact shock measure is a state where two product measures with different densities are connected. The time evolution of this state is given by a diffusion equation with respect to the position of the density step [12]- [14]. In particular, this implies that the microscopic structure of the system is known at all times.
The physical properties of large classes of one-species reaction-diffusion models with nearest-neighbour interactions have been widely studied [4], [15]- [17]. There are only four known models which show shocks without correlations [14,18,19]: the ASEP, the branching and coalescing random walk (BCRW), the asymmetric Kawasaki-Glauber process (AKGP) and the brick layer model. While the former three models are exclusion models where the number of particles on each lattice site is restricted to at most one, the particle number is not restricted in the latter one. Here, the investigation shall be extended to three-site interaction exclusion models [20,21] and it will be shown that three-site interactions give rise to models with exact shock measures which show aspects not seen in nearest-neighbour models, namely double shocks and additional symmetries even though no free fermion condition is met.
At this point we remark that the question of phase separation is directly linked to the movement of domain walls because coarsening is generic for this phenomenon. But although some of the models presented in this paper show growing domains we will argue that this mechanism cannot be used for constructing nonequilibrium models with two species (one type of particles and vacancies) showing phase separation in one dimension [22]- [24].
We present a survey of one-species models with three-site interactions and open boundaries whose time evolution is described by an exact diffusive shock measure, to be defined below. For the cases of shocks between two nonfluctuating phases (densities 0 and 1) and two fluctuating phases (both densities are different from 0 and 1) the list is complete. For the case of a shock from a nonfluctuating phase (density 0 or 1) to a fluctuating phase (density between 0 and 1) the variety of models is too large to give a complete survey-the number of free parameters rises from 12 to 56 when the interaction range is increased from 2 to 3. However, we present a classification of models with respect to their symmetries, where we considered models where at least two of the symmetries charge C, parity P and time T are valid independently.
We also address the question to what extent the description of the dynamics of the model by the movement of shock fronts is sufficient. To this end the interactions of the domain walls are determined and possible parallels to free fermion systems are discussed.

Formalism
On each lattice site k (k = 1, . . . , L) of a one-dimensional lattice there may be at most one particle (A) or a vacancy (∅). One can also consider these two-state systems as spin systems; a particle is represented by a down-spin (|↓ ) and a vacancy by an up-spin (|↑ ). The stochastic dynamics of the models are defined by a master equation [25], which is conveniently expressed in the quantum-Hamiltonian formalism for spin-1/2 chains as described in [4]. To each lattice configuration η we assign a basis vector |η which is given by the tensor product of the single-site states. The probability vector describing the system can then be written as where P(η, t) is the probability at time t to find the system in the state η. The time evolution of the system is described by a master equation which can be written as where H is the stochastic generator of the process. Owing to the analogy of equation (4) to the Schrödinger equation (in imaginary time), H is often called the Hamiltonian of the system. Conservation of probability requires that the sum of the entries of each column is zero, where s| = (1, 1, . . .) is the so-called summation vector. For finite interaction range and spatially homogeneous kinetics it is convenient to write the generator as where the local Hamiltonians h k contain the rates of the elementary local transitions and b 1 , b L−1 account for events at the left (resp. right) boundary. The operators h m include only operators acting on the sites m, m + 1 and m + 2; the operators b m include only operators acting on m and m + 1. They are formulated using the particle number operator (n k ), vacancy number operator (v k = 1 − n k ), the particle creation (s − k ) and annihilation (s + k ) operators. The lower index represents the lattice site on which the respective operator acts. The diagonal entries of the Hamiltonian have to be chosen such that conservation of probability is fulfilled which can be easily constructed by considering For example, diffusion to the right (A∅ → ∅A) with rate D r is written as h k = D r (s + k s − k+1 − n k v k+1 ).

Product and shock measures
In what follows, |· 1 denotes a probability vector for a single site and an operator without index a single-site operator. If no correlations are present the probability vectors are simply given by with In this case, the calculation of the action of the stochastic generator H can be simplified by using the following identities: A product measure |ρ is a stationary state of the system, if As a shock measure |ρ 1 , ρ 2 , k , we define the state which is a product measure with density ρ 1 up to and including site k and beginning from site k + 1 a product measure with density ρ 2 : We are interested in those systems, for which the time evolution of the shock measure is given by a diffusion equation (see figure 1) Starting from an initial density step the system will evolve into an exact diffusive shock measure defined by the number of sites is 500 and an average over 5000 systems was performed, the boundary conditions are periodic. As initial condition, we took a shock measure with steps at site 166 from 0.01 to 1/2 and at site 333 from 1/2 to 0. The left step is not stable because the product measure with ρ 1 = 0.01 is not a stationary state and the right step is not stable because there the diffusion equation is not fulfilled.
which is a time-dependent superposition of sharp shocks weighted with p l (t). This can be seen by solving equation (13) for an initial shock at position k: G(l, t|k, 0) is Green's function for this problem, for large arguments (k − l) and late times it approaches a Gaussian, as expected for a diffusion problem. From equation (15) we read off p l (t) = G(l, t|k, 0). If we choose a shock measure (12) as initial condition of a system which obeys equation (13) the form of the shock is preserved in time but due to the diffusion the state of the system will evolve into a superposition of shocks. Thus when performing an ensemble average the density profile is not a sharp step but smears out in time (see figure 2) as seen in Monte Carlo (MC) simulations. Nevertheless, a typical configuration of a single systems shows a sharp shock.
Equation (13) directly gives the diffusion coefficient D s and the shock velocity v s : The product measures to the left and to the right of the shock position are only possible for a special choice of boundary dynamics b 1 and b L−1 . A possible choice is always those boundary dynamics which is the effect of reservoirs with densities ρ 1 (resp. ρ 2 ) [26]: One imagines that the lattice is extended by two additional sites at each boundary and one calculates the action of h −1 and h −2 (h L−1 and h L ) assuming that the sites are occupied according to the product measures in the bulk and determines the effective rates for the boundary action. Thus, n −1 and s + −1 are substituted by ρ 1 , v −1 and s − −1 by (1 − ρ)-the other sites are handled accordingly. In this paper, we always consider the limit L → ∞ such that the shock is assumed not to hit the boundaries.
Given specific microscopic processes we test for the existence of diffusive shock measures by the following procedure: (a) Check whether |ρ 1 is a product measure of the system with periodic boundaries.
where c k−1 acts on sites k − 1 to k + 2. (c) Check whether |ρ 2 is a product measure of the system with periodic boundaries. (d) Set up the boundary processes b L−1 such that where d k+1 acts on sites k + 1 to k + 2. (e) Check whether In detail we perform the following: After having set up the local Hamiltonians h k the action on a product measure is brought into a diagonal form as described above, where h diag k contains only the operators n m and 1 (v can be eliminated using v = 1 − n) acting on sites m = {k, k + 1, k + 2}.
The requirement of a product measure (a) leads to the condition that the application of the Hamiltonian of the periodic system has to produce 'telescope'-sums of diagonal operators which results in five equations of the rates for each density.
Condition (e) can be checked using the identity and analogous identities for |ρ 1 , ρ 2 , k − 2 , |ρ 1 , ρ 2 , k + 1 and |ρ 1 , ρ 2 , k + 2 . A comparison of coefficients of the (only diagonal) operators on the left-and right-hand site of equation (19) then leads to another seven equations of the rates, ρ i (i = 1, 2) and δ j (j = 1, . . . , 5). Conservation of probability implies which is not an additional equation but simplifies the calculation. If one of the densities is either 0 or 1, the Hamiltonian cannot be brought to diagonal form using equations (10) and (21). However, in this case, the action of the Hamiltonian simplifies and an analogous comparison of coefficients of creation and annihilation operators is possible.

Classification of models
To find the models with three-site interactions which exhibit shock diffusion one could, in principle, try to solve equation (13) for the general Hamiltonian. Since there are 56 microscopic processes (transitions from any of the 2 3 = 8 states to any different state) this task is tedious.
We facilitate the procedure by writing a computer program which does all the symbolical calculation (transformation of the Hamiltonian into diagonal form, gathering of coefficients) and sets up the constituting set of equations. These equations are then solved by standard mathematical software. Still, the general solution of the problem is far too complex to extract useful informations. It is hence useful to investigate physically motivated subclasses, as done in the following.

Shocks from
The completely empty lattice, ρ = 0, and the fully occupied one, ρ = 1, are nonfluctuating states. They are the two ground states of the zero-temperature Ising model. In this case only few processes play a role. All processes starting from ∅∅∅ or AAA are forbidden, because otherwise product measures with densities 0 and 1 would not be stationary solutions of the system. The processes which act on the configurations ∅A∅, A∅∅, A∅A and AA∅ are ineffectual because none of these configurations is possible, if the system is initialized with a shock measure. The remaining processes are The diffusion constants are then given by The solution for a shock between ρ 1 = 1 and ρ 2 = 0 is obtained by exchanging the roles of particles and holes (A ↔ ∅). This dynamics is a generalized zero-temperature Ising model where the domain wall between spin-up and spin-down regions diffuses freely. As the space symmetric processes have no influence on the upward shock they can be included giving rise to a model allowing for both, upward and downward shocks. All 10 rates are independent such that the drift of the domain walls can be chosen freely.

Shocks from
For each Hamiltonian with three-site interactions there is either none, exactly one product measure with a density in the open interval ρ ∈ (0; 1), or infinitely many. The latter is only possible for particle-conserving Hamiltonians. The proof of this assertion is possible for the case of three-site interactions, but is rather technical and is therefore omitted here. It is based on analysing the conditions on the rates and densities for the existence of a product measure, equation (11), for the general three-site interaction Hamiltonian.
Instead, we present a general consideration why the existence of two fluctuating stationary product measures is impossible for non-particle-conserving Hamiltonians. This is not intended to be a mathematically rigorous proof, it is rather a physical consistency check. This consideration is related to the positive rates conjecture; see [22,23].
Let us assume that there are two fluctuating stationary densities ρ 1 = ρ 2 for a non-particleconserving Hamiltonian. Then due to fluctuations there is a finite probability that in the stationary state with density ρ 1 a region of density ρ 2 emerges. As ρ 1 is a stationary density this region has to be suppressed and vanishes again. Consequently, the domain walls ρ 1 |ρ 2 and ρ 2 |ρ 1 are moving towards each other shortening the domain of ρ 2 . However, if this is the case, ρ 2 cannot be a stationary state, as a region of ρ 1 emerging in a phase of ρ 2 is growing. Hence the assumption of the existence of two fluctuating stationary densities was wrong.
Note that Gacs' error-correcting model [22,23] does not represent a counter-example for this statement. It is crucial that product measures are considered: For this case no dynamics is able to determine the length of a one-dimensional domain-the domain growth is independent of the domain size, while phase separation requires a faster growth of larger domains.
In figure 2 the time evolution of a non-particle-conserving system is shown. The particle density in the low-density region is not stationary and increases until the first shock vanishes. The second shock disappears due to another mechanism; here the diffusion equation of the shock is not fulfilled and consequently the shock dissolves.
The situation is basically different if particle-conserving Hamiltonians are considered. In this case, the system is not ergodic and hence a region of different density cannot evolve in the bulk of the system. Consequently several stationary densities are possible.
As nonconserving Hamiltonians do not allow for the existence of two fluctuating phases we only need to investigate particle-conserving Hamiltonians and recover the solution of the investigation of nearest-neighbour interactions, the asymmetric simple exclusion process (ASEP) [14]: The inclusion of next-nearest-neighbour interaction does not lead to further models in this case.

Shocks from
While the existence of two fluctuating stationary states is impossible, the existence of one fluctuating and one nonfluctuating phase is easy to construct. The density ρ 1 = 0 is stationary if all processes from the empty lattice, ∅∅∅, are prohibited. This is a violation of the assumption underlying the positive rates conjecture as certain types of fluctuations are absent. Therefore, two stationary states are possible; however, the nonfluctuating phase is unstable and the fluctuating phase will always enlarge at the expense of the nonfluctuating one.
For the existence of a product measure ρ 2 ∈ (0; 1) it is then also necessary that no process to the empty lattice is present. Apart from these constraints no further processes can be excluded a priori. We restrict ourselves to the case of a shock between ρ 1 = 0 to ρ 2 ∈ (0; 1). The cases of a shock from ρ 1 = 1 to ρ 2 ∈ (0; 1) can be obtained by exchanging particles and holes, and the case of the fluctuating phase to the left can be obtained by a parity transformation.
As argued above, the degrees of freedom are too many for a complete investigation and we classify the systems by the symmetries charge C, parity P and time T . In this context, charge symmetry is the invariance of the microscopic processes under the exchange of particles and holes (A ↔ ∅), i.e. for each process in the model there exists the C-symmetric one with the same rate. Although this picture is quite artificial when applied to particles it is natural in the language of spins where in the absence of an external field the symmetry between the up-and down-spin is obvious. Parity symmetry is the invariance of the microscopic processes under the exchange of left and right, i.e. for each process in the model there exists the P-symmetric one with the same rate. Time symmetry is invariant under the transformation of t → −t. A stochastic model is T -invariant if detailed balance with respect to its stationary distribution is fulfilled [25], i.e. these models are able to reach an equilibrium steady state. This is the case if for all states η 1 ,η 2 the transition rates w η 1 →η 2 , w η 2 →η 1 and probabilities of finding the system in the configurations P(η 1 ), P(η 2 ), obey the equation w η 1 →η 2 P(η 1 ) = w η 2 →η 1 P(η 2 ), i.e. there is no net current between states. If the stationary state of a system is a product measure, the validity of detailed balance is easily checked or refuted, because the calculation of the probabilities of the system to be in a specific state is trivial and additionally the absence of correlations permits one to investigate only the local rates instead of the configurations of the whole lattice.
One has to distinguish between a combined symmetry, for example PT , from an independent symmetry, here expressed by the symbol '∧', for example P ∧ T . While in the former case the system is invariant after applying the symmetries one after each other, in the latter case the system is invariant under the symmetries applied each by themselves.

C ∧ P ∧ T symmetric systems with shock measures.
The 56 possible three-site interaction transitions can be arranged in 11 minimal models, which obey C ∧ P ∧ T symmetry each (but are not necessarily described by exact diffusive shock measures). If we exclude those which have transitions involving the configurations ∅∅∅ or AAA only six remain. There are 63 combinations which can be built out of six elements. We checked all combinations and found that there is only one model which is C ∧ P ∧ T -symmetric and has an exact diffusive shock measure as solution:   Monte Carlo simulation of the C ∧ P ∧ T symmetric model, the number of sites is 500 and an average over 10 000 systems was performed. As initial condition we took a shock measure with steps at site 166 from 0 to 1/2 and at site 333 from 1/2 to 1.
In this model diffusion is combined with branching processes which are only possible in a neighbouring lattice site if the subsequent site is empty and its reversal, the coalescence process.
Owing to the C ∧ P symmetry this model has the property that both the shock from 0 or 1 to 1/2 and the shock from 1/2 to 0 or 1 is stable. Hence, double shocks 0| 1 2 |0 and 0| 1 2 |1 are possible. We calculated the time evolution of an initial double shock (0| 1 2 |1) in a periodic system by a Monte Carlo simulation as shown in figure 3. As calculated, the shock moves with a drift of 1 and such that always the fluctuating phase penetrates the nonfluctuating phase. The shock 1|0 between site 500 and 0 evolves into two steps, 1| 1 2 and 1 2 |0, since these shocks drift into the nonfluctuating phase and split up. As a consequence of the periodic boundary conditions the two shock fronts moving into the nonfluctuating phase coalesce after finite time and the nonfluctuating phase vanishes. This is a consequence of the stability of the fluctuating phase as argued in section 4.2. A detailed discussion of the reaction of shock fronts will be given in section 5.

C ∧ P symmetric systems with shock measures.
There are 17 minimal models obeying C ∧ P symmetry; if we exclude those which have transitions involving the configurations ∅∅∅ and AAA only eight remain. Thus, in this case, 255 combinations have to be checked, but no additional model besides the C ∧ P ∧ T model can be found.  Thus, in this case, 8191 combinations have to be checked. We found 14 models, as presented in the appendix. Owing to the P symmetry in each of these models downward shocks are also stable.
In figure 4 we show a Monte Carlo simulation of the model D (see the appendix) with ω = 1/2 on a ring. As predicted two aspects can be observed: Firstly, both the upward and the downward shock are stable. Secondly, the fluctuating phase spreads into the nonfluctuating one until the latter finally vanishes. We remind that due to the superposition of shocks the ensemble average does not exhibit a sharp step although each single realization does.

C ∧ T symmetric systems.
There are 16 minimal models obeying C ∧ T symmetry; if we exclude those which include the configurations ∅∅∅ or AAA only nine remain. Thus, in this case, 511 combinations have to be checked. Besides the C ∧ P ∧ T model we find the following two models: In this model only branching to the right is possible in the presence of another zero at the nearest-neighbouring site and the corresponding coalescence processes. The domain wall 0| 1 2 is not fluctuating even though the domain is.

The second model is a combination of branching and coalescence processes in both directions
The domain wall performs an unbiased diffusion with diffusion constant D = 2.
For these two C ∧ T models, a downward shock has not necessarily to be stable as it is the case for the models which are P-symmetric. Indeed, a downward shock is not stable in the two models because otherwise the space-reflected versions of the processes would show a stable upward shock and would constitute additional C ∧ T models.

Some further models.
Here we present some models which do not belong to the classes presented above.
The following model is the totally asymmetric exclusion process combined with activated Langmuir kinetics used in [3] to show hysteresis in driven diffusive systems: The simplest model with a fluctuating shock front not included in nearest-neighbour interaction models is where particles are only allowed to hop over a vacancy to the right. This model is PT -symmetric, because reversing the direction of the process is equal to exchanging left and right. Another simple model is This model is again PT -symmetric by the same argument as in the model presented above. In this case one has to pay attention to the rates when reversing the time direction; the time reversed process is For ω = 1 the density is ρ 2 = 1 2 .

Reactions of domain walls
Up to now we have only investigated whether the shock fronts are stable and, if so, how they move. The movement of the shocks can be used to describe the dynamics of the systems as their positions are sufficient to characterize the state of the system. This has been used, for example, in [3] where the dynamics with many degrees of freedom could be reduced to an effective one-particle system. To describe a system completely by the position of the domain walls we additionally have to investigate how the domain walls affect each other.
A special case would be if the domain walls do not interact at all. This means that their rates do not change in the presence of another wall-certainly the possibility of mutual annihilation has to be included as a phase may vanish if two boundaries meet. In this case, one could interpret the dynamics as the motion of annihilating random walkers. For those systems a direct link to a free fermion system has been discussed in [4] and, thus, for the unaffected movement of domain walls a description by free fermions could be possible. To this end one may apply the Jordan-Wigner [27] transformation which converts spin-1/2 operators into fermionic creation and annihilation operators.
The significance of the Jordan-Wigner transformation in this context is that some Hamiltonians of spin-1/2 systems transform into fermion Hamiltonians which include only bilinear expressions of the fermionic operators-this can be regarded as a system of free fermions. For these Hamiltonians, additional techniques for calculating the dynamical properties are available. Interestingly the dynamics of all particle models showing stable shock fronts with two-site interactions can be represented by the free motion of domain walls without interactions [14]. For the BCRW and the AKGP this is directly related to the free fermion character of these systems.
It is the purpose of this section to investigate potential relations of the three-site interactions models found above to fermion systems, since a link between the domain wall motion and free fermion behaviour would be interesting. We first discuss the dynamics of domain walls in detail and then turn to the transformation into fermion systems.

Dynamics of domain walls
When investigating the interaction of domain walls, the situation simplifies again by the fact that within the domains the probabilities are given by a product measure. By this, two domain walls may only influence each other if the distance is smaller than three lattice sites and thus only operators acting on a small range have to be included.
The first model to be considered is the C ∧ P ∧ T model. In figure 3, on the one hand, it is shown how a domain wall 1|0 splits into two domain walls 1| 1 2 and 1 2 |0 and, on the other, how two domain walls 1 2 |0 and 0| 1 2 coalesce. These two cases are now studied in detail analytically. If the system is characterized simply by the position of the domain walls without further correlations it will be sufficient to describe the dynamics by states |k, l; ρ 1 , ρ 2 , ρ 3 = · · · ⊗|ρ 1 ⊗ · · · | ρ 1 To describe the evolution of the step 1|0 we define |k, l ≡ |k, l; 1, 1 2 , 0 . In the following, the densities are omitted for the sake of simplicity. By applying the Hamiltonian of the C ∧ P ∧ T model on this state one gets Thus, the time evolution of the system can be completely described by the movement of the domain walls, no additional correlations evolve. Starting from 1|0 both domain walls simultaneously move by one lattice site creating a fluctuating domain of two sites with density 1 2 . The time evolution of a state in which the domain walls are separated by a lattice site is given by the separate movement of both domain walls by one lattice site. A state in which the domain walls are separated by two or more lattice sites evolves simply by the rates δ 2 = 2 and δ 3 = 1 as calculated before. The time evolution of the distance of the domain walls l − k is illustrated in figure 5. Once the domain walls are separated it is not possible that they coalesce again as the l − k = 1 is an isolated point. Thus, the movements of the domain walls are not independent; a repulsive interaction is present.
Next, the interaction of the shock fronts 1 2 |0 and 0| 1 2 shall be investigated. To this end we now define |k, l ≡ |k, l; 1 2 and test again whether the dynamics can be described in terms of the |k, l . Applying the Hamiltonian on the state |k − 1, k + 1 shows that additional correlations appear and that consequently the time evolution cannot be described simply by the location of the domain walls. Nevertheless, it is instructive to analyse the appearance of the correlations in detail to reveal the link to free fermion systems. The correlations can be compensated by including the process A∅A AAA. By this the C symmetry is broken (the ρ = 1 phase is not stable anymore), but it is still a P ∧ T model.
Choosing the rate for the forward and backward reactions to be 2 one gets |k − 1, l + |k, l + 1 + 2|k + 1, l + 2|k, l − 1 − 6|k, l , l − k 3, where | 1 2 = |k, l | k−l=1 is the product measure with density ρ = 1 2 . This state is stationary which is recovered by the vanishing time derivative for k − l = 1. The time evolution of the system can be described completely by the movement of the shocks and their dynamics is independent of each other until they meet, then both are annihilated. Consequently, this system is a candidate for the description by free fermions. However, by including the process A∅A AAA, we get the model C of the appendix which is the BCRW-only a two-site interaction model. The transformation of the BCRW into a free fermion system is known [4,14].
As a second example the model J of the appendix is chosen as it is one of the simplest models. Again we define |k, l ≡ |k, l; ρ, 0, ρ , and apply the Hamiltonian. The description by the states |k, l is only closed if we set the parameter of the model w = 1 for which ρ = 1/2 and some of the rates vanish. One gets Next the time evolution of a state |k, l ≡ |k, l; 0, 1 2 , 0 is investigated. We find For k − l = 2 it turns out that the time evolution of this state cannot be described by a superposition of shock measures, i.e. additional correlations emerge. Thus, no independent movement of the shock fronts is possible in this case. This suggests that this is not a free fermion model.

Conclusions
In this paper we have investigated exact diffusive shock measures in one-dimensional reactiondiffusion systems with next-nearest-neighbour interactions and open boundaries. We distinguish the following three cases: 1. The connection of two nonfluctuating phases, the two densities are 0 and 1. The conditions that both states are stable exclude many models and we find many next-nearest-neighbour models as solution generalizing the Glauber Ising model at zero temperature.
We restricted ourselves to the case of completely ordered connected phases as initial conditions. It would be interesting to investigate how a system evolves out of random initial conditions. In this scenario, coarsening of the ordered domains or the emergence of a third stationary state which is fluctuating is possible.
2. The connection of two fluctuating phases, both densities are between 0 and 1. It is argued that, in general, nonconservative models cannot have two fluctuating product measures as solution, in agreement with the positive rates conjecture. Consequently, only conservative models have to be investigated and we recover the ASEP as the most general solution. Hence, the inclusion of three-site interactions does not lead to models not known from the investigation of nearest-neighbour interactions.
3. The connection of a nonfluctuating phase to a fluctuating phase, one density 0 or 1 and one between 0 and 1. In this case, numerous models exist and we classify the systems with respect to their symmetry. There is only one model which is C ∧ P ∧ T -invariant, and this is as well the only model which is C ∧ P-invariant. Two additional models are found that are C ∧ T -invariant and 14 models are found that are P ∧ T -invariant.
We stress that the mechanisms of exact shock measures are not suitable to construct a (one-species) model which shows phase separation on a ring. On the one hand, although models with one fluctuating phase and one nonfluctuating phase allow for stable up and downward shocks, the nonfluctuating phase will always vanish on a ring because the two shock fronts always enlarge the active region. On the other hand, conservative models, which in principle, allow for shocks between two fluctuating phases, are unable to show both up and downward shocks. This is a consequence of the collective velocity v c (ρ), which describes the movement of the centre of mass of a disturbance in a region of a certain density ρ. In order that a shock is stable a disturbance has to tend towards the shock, v c (ρ 1 ) > v s > v c (ρ 2 ), where v s is the shock velocity. Obviously, this equation can only hold either for the upward or the downward shock. Note that this consideration only holds for short-ranged, homogeneous one-species models; it is known that phase separation is possible in models with defects [28] or several species of particles [29,30].
Although it is shown that double shocks 0|ρ|0 are possible we argue that this cannot be used to construct models which show phase separation in one-dimensional, short-ranged periodic systems with a single species.
We have also investigated the influence of shock fronts on each other. A case of special interest is when the fronts move independently, except for the possibility of mutual annihilation. In the C ∧ P ∧ T model this cannot be observed; it turns out that it has to be combined with an additional process-violating C symmetry. However, by this, the BCRW which is a two-site interaction model is recovered for which the independence of shock fronts is known.
We conclude that there is no direct connection of models whose time evolution is given by exact diffusive shock measures and free fermion systems.