Eternal non-Markovianity: from random unitary to Markov chain realisations

The theoretical description of quantum dynamics in an intriguing way does not necessarily imply the underlying dynamics is indeed intriguing. Here we show how a known very interesting master equation with an always negative decay rate [eternal non-Markovianity (ENM)] arises from simple stochastic Schrödinger dynamics (random unitary dynamics). Equivalently, it may be seen as arising from a mixture of Markov (semi-group) open system dynamics. Both these approaches lead to a more general family of CPT maps, characterized by a point within a parameter triangle. Our results show how ENM quantum dynamics can be realised easily in the laboratory. Moreover, we find a quantum time-continuously measured (quantum trajectory) realisation of the dynamics of the ENM master equation based on unitary transformations and projective measurements in an extended Hilbert space, guided by a classical Markov process. Furthermore, a Gorini-Kossakowski-Sudarshan-Lindblad (GKSL) representation of the dynamics in an extended Hilbert space can be found, with a remarkable property: there is no dynamics in the ancilla state. Finally, analogous constructions for two qubits extend these results from non-CP-divisible to non-P-divisible dynamics.

A realistic modelling of many quantum phenomena inevitably needs to take into account the interaction of our system of interest with environmental degrees of freedom. Thus, in order to describe the quantum system dynamics appropriately, one is often forced to deal with open quantum systems. A very relevant and well understood class of such open quantum system dynamics follows from a Markov master equation of GKSL form 1,2 . Non-Markovian behaviour may arise from a structured environment or strong system-environment interaction 3 . Non-Markovian systems are very challenging: in the often-employed projection operator formalism their dynamics involves memory kernels 4,5 . Other approaches range from path integrals 6,7 , over hierarchical equations of motion (HEOM) for the reduced density matrix 8,9 , to hierarchies of stochastic pure states (HOPS) 10,11 . Sometimes time-convolutionless master equations can be used 12 . During the last few years, due to tremendous experimental progress in quantum technologies in many different areas and more and more refined measurement schemes, specific investigations of non-Markovian quantum dynamics, where GKSL is no longer applicable, have become possible [13][14][15][16][17][18] . Recent experiments also demonstrate how to use non-Markovianity for entanglement preservation 19 and for a quantum information protocol 20 .
The theory of non-Markovian quantum dynamics is much less developed than the GKSL class and subject of tremendous research over the last decade and more. A very valid point of view would be to call any dynamics other than GKSL semigroup evolution "non-Markovian". A more detailed analysis, however, reveals an astonishing variety of possible definitions of what constitutes non-Markovian dynamics [21][22][23] , and therefore a large number of definitions and measures of non-Markovianity have been proposed [24][25][26][27][28][29][30][31][32] . So far, most studies are based on the effective dynamics of the reduced density operator, other consider the full dynamics of system and environment 33,34 .
As mentioned earlier, in some cases of interest, the open system dynamics may be written in terms of a time-local master equation involving time-dependent functions as prefactors with otherwise GKSL form. Then, for some periods of time negative decay rates may show up, which according to some measures indicates non-Markovian dynamics 35,36 . Recently, a remarkable master equation for a qubit was presented involving an 1 Institut für Theoretische Physik, Technische Universität Dresden, D-01062, Dresden, Germany. 2 Institute of Physics, Faculty of Physics, Astronomy and Informatics, Nicolaus Copernicus University, Grudziadzka 5/7, 87-100, Toruń, Poland. 3  always negative decay rate of an otherwise GKSL-type-looking master equation. It was termed the master equation of eternal non-Markovianity (ENM master equation) 35,37 .
We expect non-Markovian dynamics to be related to some form of memory-dependence arising from the dynamics of the environmental degrees of freedom. This is why non-Markovianity is associated to a "backflow of information" 25,[38][39][40] or to the occurrence of quantum memory 23 , or simply to a joint complex system-environment dynamics 41 . In such cases, the measures detect non-Markovianity. In this contribution we want to emphasize, however, that the reverse need not be true: there are non-Markovian master equations (according to one of the definitions), whose physical realisation does not support any notion of such "memory effects". Instead, either there is no dynamical environment at all, the dynamics can be realised by a classical Markov process or, when embedded in a larger Hilbert space, there is no dynamics of the environmental state.
In this paper we derive the ENM master equation from an appropriate mixture of Markov dynamics in two (related) ways: one is based on random unitary evolution, the second approach uses a mixture of Markov GKSL maps. By highlighting the equivalence of all these dynamics on the reduced level, we show explicitly how ENM evolution of a qubit could be realised in a laboratory either with a white noise or with a classical jump process with time independent jump probabilities. Moreover, also the bipartite GKSL representation, for which the ancilla state is frozen, is possible. Nonetheless, we may choose to describe the dynamics in terms of a negative-rate time-local master equation, or, involving a non-trivial memory integral. These findings support the point of view that the interpretation of non-Markovianity is elusive and great care has to be taken when talking about memory effects based solely on a reduced (master equation) description.

Time-local master equations and negative decay rates
For any total Hamiltonian of system and environment and for any product initial state, the dynamics of an open quantum system can be expressed in terms of the dynamical map ρ(t) = Λ t [ρ(0)], with Λ t completely positive and trace preserving (CPT). If Λ t is an invertible map then one finds the corresponding time-local generator follows. Assuming the semi-group property Λ t+s = Λ t Λ s , the generator takes the GKSL form 1,2 (ħ = 1): with decay rates γ i (t) ≥ 0, for all i. Equation (2) defines a reasonable dynamics if applied to any state at any time and therefore defines a CP-divisible dynamical map Λ t 43 , i.e. the dynamical map Λ t satisfies the following property Λ t = Λ t,s Λ s and the family of maps (propagators) Λ t,s is CPT for any t > s. It seems natural to regard dynamical maps Λ t with master equations of type (2) for which γ i (t) < 0 for some i and some t as candidates for non-Markovian quantum dynamics. In these cases, the dynamical map is no longer CP-divisible. Indeed, some authors 35 propose to use the negativity of decoherence rates as a definition of non-Markovianity of the dynamics. This approach is based on the fact that the canonical form of the master equation, defined in analogy to the Markov case (so the time dependent Lindblad operators are traceless, normalized and mutually orthogonal), is unique. Consequently, to all CPT maps generated by a master equation of form (2) one can uniquely assign a set of γ i (t).
Actually, one also considers Λ t,s which is not necessarily CP. If Λ t,s is positive for all t > s then one calls the evolution P-divisible. Recently, this notion was refined in ref. 44 as follows: the evolution is k-divisible if Λ t,s is k-positive. CP-divisibility is fully characterised by the corresponding time-local generator t  -all local decoherence rates γ i (t) are always non-negative. P-divisibility is more difficult to characterise on the level of the generator. One has the following property: if Λ t is P-divisible, then for all Hermitian operators X, where ||⋅|| 1 is a trace norm. Actually, when Λ t is invertible then (3) implies P-divisibility. This property is very close to the so-called BLP condition 25 which says that Λ t defines Markovian evolution if where σ k are the Pauli spin operators.
, and {i, j, k} run over the cyclic permutations of {1, 2, 3}, one has the following conditions which guarantee that the evolution Λ t is CPT: Clearly, the corresponding dynamical map is CP-divisible iff γ k (t) ≥ 0. Interestingly, the dynamical map is P-divisible iff the weaker conditions are satisfied 46,47 i j given the validity of (6). Actually, in this case P-divisibility is equivalent to the BLP condition (4).
Using the geometric measure of non-Markovianity based on the volume of admissible states 31 , our one qubit dynamics is classified as Markov, too, as for all times γ 1 (t) + γ 2 (t) + γ 3 (t) > 0 46 is satisfied.
An interesting example of the generator was proposed in ref. 35 -the ENM master equation, with where one rate is always negative: γ 3 (t) < 0 for all t > 0. One easily checks that (6) are satisfied and hence the dynamical map is CPT. Clearly, the corresponding dynamical map is not CP-divisible because of the negativity of γ 3 (t). Moreover, conditions (7) are also satisfied which implies that the map is P-divisible 46,47 . Is this evolution non-Markovian? Based on the concept of CP-divisibility it is clearly non-Markovian. However, it satisfies condition (4), hence it is Markovian according to BLP. In the following we want to argue that the meaning of non-Markovianity for non-CP-divisible maps like those generated by (2) with an always negative rate (8) needs to be discussed carefully. In particular, it can be highly misleading here to relate the formal property of "non-Markovianity" according to one of its definitions to some notion of "complex system-environment dynamics" or "backflow of information" from environment to system as will be exemplified in this paper.
We show that there is a whole family of master equations of type (2) with γ i (t) < 0, for some i and times t that (i) turn out to arise from random unitary Schrödinger dynamics, (ii) are mere mixtures of Markovian semi-group dynamics, (iii) allow for a physical realisation based on a classical Markov process. With these observations in mind, it is obvious, that the ENM master equation (or its two-qubits extension, see the "From one to two qubits dynamics and breaking also P-divisibility" and the "Bipartite GKSL representation" sections) needs not be related to any information backflow from dynamical environment. The particular choice (8) turns out to be a special case of this more general family of evolutions.

Markov dephasing dynamics
To start with, consider simple dephasing dynamics of a qubit given by a master equation of GKSL type 48 Since this CPT map is unital, the dynamics is of random unitary or random external field type [49][50][51][52][53] . In fact, a physical realisation of Eq. (9) for pure initial states is easily obtained from a fluctuating field ξ(t) driving the unitary Schrödinger dynamics: we find for an arbitrary initial condition: As shown in Supplemetary Information, the noise average can easily be performed analytically to give the solution of (9) in Kraus form

Mixture of Markov dephasing dynamics
Now we allow the direction α n of the dephasing to be random with probability distribution α p n ( ). From (12) we see that with α α α = α n n n n ( ( ), ( ), ( )) 1 2 3 , the averaged dynamics depends on the second order correlations only. Due to an overall orthogonal freedom of the whole problem, we may assume a diagonal (x kl ) and will from now on use the notation k k 2 assuming that x kl = 0 for k ≠ l. As final result, the dynamical map arising from averaging over noise ξ(t) and direction α n is again a map given in Kraus form by: The three positive parameters x 1 , x 2 , x 3 , with x 1 + x 2 + x 3 = 1 (the Cartesian variances of the distribution) are the only quantities of α p n ( ) that determine the dynamics. In Bloch representation this corresponds to a monotonic and (in general) anisotropic shrinking of the Bloch sphere (see Supplementary Information).
Therefore, it also follows that (15) can be obtained from a mixture of just three orthogonal dephasing directions along the Cartesian axes. Accordingly, the underlying dynamical map may be written as a mixture of three Markov (semigroup) dynamical maps according to (9). The variances may thus be seen as probabilities x k of choosing either of three semigroup evolutions  e t k for the dynamics. We conclude that dephasing dynamics in random directions can be written in two ways as a mixture of CP-divisible maps. Representation (15) is a continuous mixture of unitary (Schrödinger) time evolutions, while in (16) we have a discrete, finite sum of irreversible Markov GKSL dynamics. As we will show next, the corresponding master equation is just (5), with possibly negative rates.
Master equation and negativity of decay rates. As shown in Supplementary Information, we find that the map Λ t from (16) satisfies the time-local master equation with the generator of the dynamics acting on density operators according to t k k k k 1 3  as in (5). The time dependent decoherence rates can be expressed as As we will work out in detail, these rates need not be positive. Thus, the random mixture of Markovian dephasing leads to a time-local master equation with possibly negative decay rates.
Discussion of the negativity of the rates. The parameter set of variances (or probabilities) x 1 , x 2 , x 3 with x 1 + x 2 + x 3 = 1 and x k positive represents a triangular area in 3-dimensional space spanned by the vectors  r = (x 1 , x 2 , x 3 ), see Fig. 1. We refer to that set as the parameter triangle. We display in Fig. 1 (hatched) that subset of parameters, for which all γ k (t) are positive at that particular time t: (a) t = 0, (b) some intermediate time t > 0, and (c) t → ∞. Clearly, initially for t = 0, all γ k (0) = 2x k ≥ 0 are non-negative. Later, only a symmetric triangular-star shaped region near the centre reaching out to the tips of the parameter triangle corresponds to choices of parameters for which all γ k (t) are non-negative. Regions near the edges of the parameter triangle but away from the lines connecting the vertices with the centre of the triangle correspond to choices of the x k that lead to a negative γ k (t) for some t > t * . As t → ∞, an asymptotic finite area of that shape remains (we call it asymptotic area) for which all γ k (t) ≥ 0 for all times. We will investigate the shape and size of that area in more detail later.
The rates have the following seven properties: (i) All rates start off non-negatively, γ k (0) = 2x k ≥ 0. (ii) At most one γ k (t) can turn negative. (iii) Once a γ k (t) turns negative at t = t * , it remains negative ever after: γ k (t) < 0 for all and x 3 = 0. In those cases the dephasing is complete in that direction, and the corresponding probability distribution has a form p( α n ) = p(n i , n j )δ(n k ). (vi) For all parameters outside the asymptotic parameter area there exists some time t * > 0, so that for all t < t * all γ k (t) are positive, and for all t > t * one of the γ k (t) is negative. (vii) We have γ 1 (t) + γ 2 (t) ≥ 0 for all times (and cyclic) and thus, the dynamics is P-divisible for all times and all choices of parameters 44,46,54 .
We thus see that (quasi-)GKSL dynamics is only realised for our dephasing in random directions-process for choices of (x 1 , x 2 , x 3 ) within the asymptotic parameter area. Outside that area one of the rates turns negative eventually (or immediately, for values at the border) and thus, the corresponding CPT map is not CP-divisible. Remarkably, for all possible choices of parameters, the map is P-divisible 46 .
If the possible parameters (x 1 , x 2 , x 3 ) are uniformly distributed over the parameter triangle, the probability for the corresponding dephasing process in random directions to be of quasi-GKSL type is just the area of the asymptotic area relative to the full parameter triangle.
As expanded in detail in Supplementary Information, in an appropriate parametrization, the shape of the asymptotic area is determined by one of Newton's cubic curves 55 , x y x y 0 (20) 2 For the relative area of parameters outside the asymptotic (hatched) area, we find

Memory kernel master equation
It is worth noting that the dynamics (15) can also be described with a master equation involving a memory ker- t 0 For our dynamics, we find a kernel of the following form: Interestingly, the memory kernel K(t) has the following structure    is just the weighted sum of the three Cartesian GKSL dephasing generators. The non-local part depends on three smooth functions η k (t).
As observed in ref. 56 and confirmed here, a local in time master equation description of the dynamics has complementary properties to a memory kernel master equation, in the sense that a "nice" functional form in one formulation may lead to a more singular description in the other.
We see that the mixture of Markovian dephasing dynamics studied in this paper "    3 " may well be written in a form involving a "memory integral", that is, apart from the more or less clear local term K loc it contains a truly non-local part K nloc (t). In open quantum system dynamics, non-local master equations of type (22) appear naturally from a dynamical environment, as, for instance, in the Nakajima-Zwanzig approach 4, 5 . Obviously, no dynamical environment exists in our constructions.

Classical Markov process representation of dynamics
So far we have acknowledged that the simple mixture of Markovian dynamics may well lead to a master equation involving negative rates. Remarkably, as we will explain in this section, that latter master equation may easily be simulated using a classical Markov process.
We start with the Kraus representation of the dephasing dynamics in random directions, Eq. (15). We introduce the unitarily transformed states ρ k := σ k ρ(0)σ k , k = 0, …, 3 (with σ 0 = 1) and corresponding probabilities p k (t) such that the state at time t reads ρ ρ = ∑ = t p t ( ) ( ) can be read off from the Kraus representation (15). As elaborated upon in Supplementary Information, these probabilities are solutions of the rate equations k j j k j k j k with positive an time-independent rates Γ 0→k = x k , Γ k→0 = 1, and all other rates being zero. The corresponding transitions are displayed in Fig. 2. Most remarkably, despite the negativity of the rates of the underlying quantum master equation, its solution ρ(t) can be obtained from the classical Markov master equation (29) according to the following construction. Take a classical process between four classical states {r 0 , r 1 , r 2 , r 3 } as determined from the classical master equation (29). For a transition from state r 0 to some r k (with k = 1, 2, 3), apply the unitary transformation σ k to the state, so that ρ 0 → ρ k occurs with rate Γ 0→k = x k . Equally, if a jump from r k (k = 1, 2, 3) back to r 0 occurs, again apply the unitary σ k to the current state so that ρ k → ρ 0 with rate Γ k→0 = 1. No other jumps can take place, see Fig. 2.
By construction, is the solution (15). Consequently, one can also simply generate the probability distribution p k (t) simulating the classical Markov process and afterwards accordingly mix the final density matrix using the four ρ k .
We have managed to describe the process (5) t j j j j 0 3 we obtain the following equation for the probability 4-vector P t ( ) (for clarity we suppress the time dependence of γ j (t)): where γ 0 (t) := γ 1 (t) + γ 2 (t) + γ 3 (t). It can be rewritten in the form of a Pauli master equation k j j k j k j k with γ 0→j (t) = γ j→0 (t) = γ j (t), γ k→j (t) = γ l (t), for k ≠ j ≠ l (k, j, l = 1, 2, 3). As for the quantum master equation the transition rates can turn negative, equation (32) does not define a proper Markov process. The solution can be obtained from the propagator T(t) as given in Supplementary Information. For the initial condition = P (0) (1, 0, 0, 0) T we find positive P k (t) for all t. Thus, despite the negative rates, the master equation (32) defines a proper evolution for a probability distribution for that particular choice of P (0). Similarly, for that initial condition only, we have Λ 0 = 1.
Due to the negative rates one is tempted to think of (32) as representing a non-Markovian jump process. Yet, it is clear that = P t p t ( ) ( ) for initial condition (1, 0, 0, 0) T , and P t ( ) is therefore also a solution of a Markovian jump process (29). Hence, one should also be careful with the interpretation of classical master equations involving negative rates.
ScIentIfIc REPORTS | 7: 6379 | DOI:10.1038/s41598-017-06059-5 = . P t ( ) 0 (35) 3 Hence P 3 (t) is irrelevant and the solution is generated from (28) via the simplified Markov chain: Realisation with orthogonal states. Note that the ρ j = σ j ρ(0)σ j are not mutually orthogonal, so they cannot be distinguished faithfully by a measurement. However, a truly classical implementation involving four classical (i.e. orthogonal) states can be found by expanding the dimension of our system to four qubits ( ). For this construction, we define the following extended dynamics involving the three ancilla qubits: where U k are unitary operators, specified below.
Tracing out the ancilla (B) degrees of freedom, this dynamics reduces to (5).
To construct the four orthonormal states we write the initial density operator in diagonal form: with orthonormal vectors |φ 1 〉, |φ 2 〉 and non-negative probabilities p 1 , p 2 = 1 − p 1 . Our four-qubit states are defined in the following way: where the U k are chosen, such that U k |ψ i 〉 = |ψ 2k+i 〉 and  ψ ∈ l 8 , l = 1, …, 8 are mutually orthogonal and normalized.
These four vectors of course don't build a basis of the  16 . Nonetheless, if we set |Ψ 0 〉〈Ψ 0 | as the initial state of our four-qubit system and let it evolve according to (37), the output state is always a mixture of these four states.
Consequently, we get a realisation of the dynamics (5) with distinguishable states. In a lab, therefore, one might choose to measure in a time-continuous fashion the actual four-qubit state such as to have a time-continuous (Markov) realisation of the classical process described in Fig. 2. By construction, the ensemble mean of the corresponding reduced states, at all times of continuous monitoring, is a solution of the original negative-rate master equation.
From one to two qubits dynamics and breaking also P-divisibility. From the non-CP-divisibility of the one qubit dynamical map (16) one can conclude that the corresponding map for two qubits, where the first qubit undergoes the dynamics (16) and the second one is frozen, is not P-divisible. Nonetheless, also in this case we can find a classic Markov process representation, which can be realised with orthogonal states.
To show this we expand the initial state of the two-qubit system in a following form: AB ikmn ikmn i k m n 1 2 where |φ 1 〉, |φ 2 〉 are the eigenstates of the first qubit A (with the corresponding eigenvalues p 1 , p 2 , p 1 + p 2 = 1) and |ψ 1 〉, |ψ 2 〉 are two orthogonal states of the second qubit B. Equation (38) represents a general initial state, also entangled states are included.
The coefficients a ikmn are some complex numbers, which have to satisfy , following from ρ AB (0) = Tr C (ρ 0 ) = Tr C (|ξ 0 〉〈ξ 0 |). Per construction, Eq. (40) is fulfilled, also the positivity of the p 1 , p 2 is guaranteed for all c ikl . The other conditions for a ikmn put some constraints on the possible choice of c ikl .
The other jump states take the form (j = 1, 2, 3): where the unitary V j are chosen, such that V j |χ l 〉 = |χ 4j+l 〉 and  χ ∈ i 16 , i = 1, …, 16 are mutually orthogonal and normalized. Consequently, |ξ 0 〉, …, |ξ 3 〉 are mutually orthogonal. To achieve this we have to extend our Hilbert space by four qubits, so overall our system consists of six qubits.
Accordingly, also the dynamics of two qubits, where the first undergoes (16) and the second is frozen, can be mapped on the (time-continuous limit of the) Markov jump process graphically represented in Fig. 2, where the states ρ k are redefined. From this we conclude, that there are non-P-divisible maps, for which a classical Markov process description is possible. Therefore, both non-CP-divisibility 46 , but also the weaker non-P-divisibility 38 , are questionable indicators for the occurrence of memory effects associated with dynamics of environmental degrees of freedom.
Bipartite GKSL representation. Interestingly, the dynamics defined by (15) may be represented via where  denotes a time independent bipartite GKSL generator. This construction is based on the correlated projection method 58,59 : one defines the initial state of the bipartite system to be the following quantum-classical state where |i〉 are orthonormal vectors in   = E 3 . Suppose now that the generator  gives rise to , that is, the bipartite evolution preserves the structure (42). Then the partial trace Note that in general this prescription does not define a dynamical map 58 and hence one arrives at the legitimate map (41).
Let us define  by Such a bipartite Markovian dynamics, which potentially gives rise to the non-Markovian evolution on the reduced level, was already widely described in the recent literature, e.g. in refs 59-63. Notice however the qualitative difference of our description to the cited one: as is apparent from (43) in our case the dynamics of the ancilla state is frozen (the reduced density matrix of the ancilla does not change) and there is never any entanglement between the system and an ancilla. That means that the ancilla is only a "casual bystander" during the whole dynamics t > 0. Consequently, it is hard to see any information backflow in this construction.
The corresponding GKSL master equation also exists in the extended two qubits case: . The dynamics of the first qubit is defined by (15), the second one and the ancilla state are frozen. Notice, that the initial state of the two qubits can be chosen arbitrarily.
Also here the ancilla is only a "casual bystander" during the whole evolution t > 0. Actually, as can be easily seen from the above construction, such an embedding in a bipartite GKSL equation with a "casual bystander" ancilla is possible for all dynamics, which can be written as a time-independent mixture of GKSL evolutions.  may contain exactly one decoherence rate γ k (t) which is negative for t > t * . Based on the concept of CP-divisibility such evolution is immediately classified as non-Markovian. Interestingly, within this class the evolution is P-divisible and hence Markovian according to the concept of information flow 25 . This is, therefore, another example showing that these two concepts do not coincide. Equivalently, Λ t satisfies memory kernel master equation with the memory kernel K(t) possessing apart from the local part δ + + x x x t ( ) ( ) 1 1 2 2 3 3    a non-trivial non-local term suggesting the presence of memory effects.

Conclusions
More interestingly, however, we showed that Λ t may be easily realised as stochastic averaging of the purely unitary evolution governed by dephasing dynamics in random directions. Alternatively, there is a realisation based on a classical Markov process, where the probabilities p k (t) are governed by a classical Pauli master equation. Such a classical Markov representation exists also for a non-P-divisible dynamics of an extended two qubit system. In both cases a description with a bipartite GKSL equation, where the ancilla state is frozen, is possible, too. These realisations show that actually there is no room for physical memory effects. This proves that the interpretation of both time-local and memory kernel master equations with respect to memory effects is a delicate issue. A reduced description may not suffice to study the physics of memory in terms of information flow.