aa r X i v : . [ qu a n t - ph ] F e b On the Markov evolution of the ρ -matrix of a subsystem. M.A.Braun
Dep. High-energy physics, Saint-Petersburg State University, Russia
February 4, 2021
Abstract
Evolution of the reduced density matrix for a subsystem is studied to determine deviationsfrom its Markov character for a system consisting of a closed chain of N oscillators with one ofthem serving as a subsystem. The dependence on N and on the coupling of the two subsystemsis investigated numerically. The found deviations strongly depend on N and the coupling. Inthe most beneficial case with N − If the whole system is split in two interacting subsystems, 0 and 1, the subsystem 0 under investi-gation (observable) and the other 1 as ”the bath”, one defines the reduced density matrix ρ as apartial trace of the total one ρ over the bath variables. The time evolution of the reduced densitymatrix ρ is of course uniquely determined by the evolution of the total density matrix and tosome degree should depend on the evolution of the ”bath”. It is customary to describe the timeevolution of the reduced density matrix by the Franke-Lindblad-Gorini-Kossakowsky-Sudarshan(FLGKS) equation [1, 2, 3, 4]. This equation is a consequence of the four properties assumed forthe evolution: conservation of unit trace, non-negativity of diagonal matrix elements, linearity andmarkovian character. It is the latter property which finally allows to derive the FLGKS equation(see [5, 6, 7, 8, 9]). And of course this is the least convincing assumption. In simple words itallows to predict the reduced matrix ρ ( t ) at time t provided it is known at some previous time t . Mathematically it implies that the derivative in time dρ ( t ) /dt is expressed by the value of ρ ( t ) itself and not by some integral of ρ ( t ′ ) over intermediate times t ′ , t < t ′ < t . In fact thisassumption looks rather doubtful from the start. Many authors have tried to derive this marko-vianity from certain specific properties of the bath (see e.g [6] and lecture notes [8]). All thesederivations actually include additional assumptions together with certain particular properties ofthe bath. In fact in these papers the attention was centered on deriving the FLKS equation orits certain generalizations for the condensed matter applications with the baths corresponding toa more or less realistic macroscopic and time dependent surroundings.In this note we abandon this purposeful approach and attempt a mini-investigation of themarkovianity itself for a bath of finite or even small dimensions. We study a more abstract problemof splitting a system in two without bothering about their degrees of freedom. As a simplest casewe study the exactly soluble model of N interacting one-dimensional oscillators one of which is ourobserved system 0 and the other N − ≡ N from 1 to N serve as the bath (system 1). In this1ase both ρ ( t ) and ρ ( t ) are explicitly known. For the Markov evolution, once we know the reduceddensity matrix ρ ( t ) at time t , its value at time t should be uniquely determined, that is ρ ( t )is fully determined by ρ ( t ). Any change in the initial total ρ matrix at time t which varies thestate of system 1 but preserves the initial value of the reduced matrix ρ ( t ) should not change ρ ( t ) at time t > t .Our study shows that of course generally this is not the case. The reduced density matrix attime t is not uniquely determined by its value at t < t but depends on the initial state of thebath, clearly manifesting non-markovianity of its evolution. This dependence depends on the timeelapsed in the course of evolution, on the properties of the bath, in particular on the number N − T is the average characteristic time of the bath (thatis the time at which the bath strongly varies) then with N = 100 at evolution times t of the order t/T ∼
50 the deviation from markovianity falls to 3%. However for other sort of couplings and ofcourse with a small bath of only 2 oscillators this is not the case and deviations of markovianitystay at the order 30% or even more.
We use the model borrowed from [10]. In this section we just recapitulate the main points necessaryfor the following.Take N oscillators on a one-dimensional circular lattice, H = 12 N = N − X a =0 (cid:2) p a + ω x a + λ ( x a − x a +1 ) (cid:3) , (1)with periodic boundary conditions x a + N = x a . Here ~ = 1 and ω has the units of energy. In thecalculations we put ω = 1. To solve the above system, one rewrites the Hamiltonian in terms ofthe normal modes, H = 12 N X k =0 (cid:2) | ˜ p k | + ˜ ω k | ˜ x k | (cid:3) , (2)where the transformation to the normal-mode basis is achieved by a (discrete) Fourier transform,˜ x k ≡ √ N N X a =0 exp( − πi kN a ) x a . (3)The normal mode momenta are defined as˜ p k ≡ √ N N X a =0 exp( 2 πi kN a ) p a . (4)They satisfy the standard commutation relations: [˜ x k , ˜ p k ′ ] = iδ kk ′ and [˜ x k , ˜ x k ′ ] = 0 = [˜ p k , ˜ p k ′ ].The normal-mode frequencies ˜ ω k are found in terms of the physical frequencies ω and λ in theHamiltonian (1) as follows [10]: ˜ ω k = ω + 4 λ sin πkN . (5)2quation (2) reduces the whole system to the set of N decoupled harmonic oscillators, with theground-state wave function ψ (˜ x , ˜ x , ˜ x , · · · ) = N Y k =0 (cid:18) ˜ ω k π (cid:19) / exp (cid:20) −
12 ˜ ω k | ˜ x k | (cid:21) . (6)To express the wave function in terms of the original variables x a in the position basis one writesthe Fourier transformation (3) between the position and normal-mode bases as ˜ x = S x , with S = 1 √ N . . . µ µ . . . µ N − µ µ . . . µ N − ... ... ... . . . ...1 µ N − µ N − . . . µ ( N − (7)where µ ≡ exp( − πi/N ). Since S is a unitary matrix, i.e. S † S = 1, the inverse transformation isgiven by x = S † ˜ x . In particular the ground state (6) is rewritten as ψ ( x a ) = N Y k =0 ( ω k π ) exp (cid:20) − x T A T x (cid:21) with A T = S † ˜ A T S and ˜ A T = diag(˜ ω , . . . , ˜ ω N − ) . (8)The relation ˜ ω k = ˜ ω N − k ensures that A T is real.Note that mode variables ˜ x are complex. Setting ˜ x = ξ + iη and ˜ x ∗ = ξ − iη one finds thatin terms of real variables our system consists of two identical copies of N oscillators with modecoordinates ξ and η . In terms of real variables we find H = 12 N X k =0 (cid:16) ∂ ∂ξ k + ∂ ∂η k + ˜ ω k ( ξ + η ) (cid:17) , (9)from which we conclude that the oscillators in ξ and η evolve independently. In particulary theevolution of the oscillators ξ is realized by the product of oscillator Green functions N Y k =0 (cid:16) g k πi (cid:17) / exp n i f k ( ξ k + ξ ′ k ) − g k ξ k ξ ′ k ] o (10)with f k = ˜ ω k cot(˜ ω k t ) and g j = ˜ ω k / sin(˜ ω k t ). The evolution of the oscillators η will be achieved bythe same product with ξ, ξ ′ → ηη ′ .We take our initial wave function as a Gaussian: ψ ( x ) = (cid:16) det(Re Ω) π (cid:17) / exp (cid:16) −
12 ( x Ω x ) (cid:17) , (11)3here ( x Ω x ) = P jk x j Ω jk x k and Ω is some matrix which characterizes the initial state of the wholesystem. Using ˜ x = Sx and so x = S − ˜ x where S is unitary we transform (11) to variables ˜ xψ (˜ x ) = (cid:16) det(Re ˜Ω) π (cid:17) / exp (cid:16) −
12 (˜ x † ˜Ω˜ x ) (cid:17) , (12)where ˜Ω = S Ω S − . In terms of real variables(˜ x † ˜Ω˜ x ) = ( ξ − iη | ˜Ω | ξ + iη ) = ( ξ ˜Ω ξ ) + ( η ˜Ω η ) − i ( η ˜Ω ξ ) + i ( ξ ˜Ω η ) . Due to symmetry of the matrix ˜Ω the two last terms cancel and we find ψ (˜ x ) = (cid:16) det(Re ˜Ω) π (cid:17) / exp (cid:16) −
12 ( ξ ˜Ω ξ ) −
12 ( η ˜Ω η ) (cid:17) , (13)which is the product of two independent wave functions for oscillators ξ and η . So they will involveindependently.Consider the evolution of oscillators ξ . We have ψ ( ξ, t ) = (cid:16) det(Re ˜Ω) π (cid:17) / N Y k =0 (cid:16) g k πi (cid:17) / exp (cid:16) i f k ξ k (cid:17) × Z dξ ′ k exp (cid:16) i f k ξ ′ − ig k ξ k ξ ′ k ] −
12 ( ξ ′ ˜Ω ξ ′ ) (cid:17) . The exponent P of the last exponential can be presented as P = X j,k (cid:16) − ξ ′ j ( ˜Ω jk − if k δ jk ) ξ ′ k (cid:17) ≡ (cid:16) ξ ′ ( ˜Ω − if ) ξ ′ (cid:17) , (14)where f jk ≡ f k δ jk . The Gaussian N dimensional integral over ξ ′ k gives Z Y k dξ ′ k e P = (cid:16) π N det( ˜Ω − if ) (cid:17) / exp h (cid:16) g ( ˜Ω − if ) − g (cid:17)i , (15)where g jk = g k δ jk . Taking into account the exponential outside the integral we find in the end ψ ( ξ ) = C exp h − (cid:16) ξ ˜Ω( t ) ξ (cid:17)i (16)where ˜Ω( t ) = g ( ˜Ω(0) − if ) − g − if (17)and we designated the initial ˜Ω in (13) as ˜Ω(0).The evolution of ψ ( η, t ) will give a similar expression with ξ → η . Combining these two resultswe shall find in variables ˜ xψ (˜ x, t ) = (cid:16) det(Re ˜Ω( t )) π (cid:17) / exp h − (cid:16) ˜ x † ˜Ω( t )˜ x (cid:17)i . (18)The normalization coefficient is obvious from the fact that the wave function remains correctlynormalized during the time evolution. 4ow we return to our initial real variables x . To do this it is necessary to rotate our matriceswith the transformation S . So we finally find ψ ( x, t ) = (cid:16) det(Re Ω( t )) π (cid:17) / exp h − (cid:16) x Ω( t ) x (cid:17) ) i , (19)where Ω( t ) = G (Ω(0) − iF ) G − iF (20)with F = S − f S, G = S − gS (21)and we designated the initial matrix Ω in (11) as Ω(0). In the position representation the ρ -matrix as an operator is constructed from the wave functionand its complex conjugate. ˆ ρ = Z N Y i =0 dx i dx ′ i ρ ( x, x ′ ) | x >< x ′ | (22)where x and x ′ are N -dimensional vectors x = { x , x , ...x N } and similarly for x ′ . At at the initialtime ρ ≡ ρ ( x, x ′ ) = ψ ( x )) ψ ∗ ( x ′ ) = (cid:18) det Re Ω) π (cid:19) / exp (cid:18) −
12 ( x Ω x ) −
12 ( x ′ Ω x ′ ) ∗ (cid:19) (23)and at time t the ρ -matrix will be ρ ′ = ρ ′ ( x, x ′ ) = ψ ( x, t ) ψ ∗ ( x ′ , t ) = (cid:18) det Re Ω ′ ) π (cid:19) / exp (cid:18) −
12 ( x Ω ′ x ) −
12 ( x ′ Ω ′ x ′ ) ∗ (cid:19) , (24)where Ω ′ ≡ Ω( t ) is determined from Ω ≡ Ω(0) by (20). Of course we have the standard propertyTr ˆ ρ = Tr ˆ ρ = 1 . We will be interested in analyzing the reduced density matrix. To that end, we partition thesystem into two subsystems, taking oscillator 0 to be our “system” and oscillators 1 . . . N − ≡ N to be the “bath”. We form the reduced density matrix of oscillator 0 as ρ = Tr ...N ρ , (25)where ρ is the density matrix of the entire system.In the position representationˆ ρ = Z dx dx ′ ρ ( x , x ′ ) | x ih x ′ | , (26)where ρ ( x , x ′ ) = Z (cid:16) N Y j =1 dx j (cid:17) ρ ( x , x , .., x N | x ′ , x ..., x N ) (27)5o do the integrations we separate variables x , ..., x N in the exponent of (23) with ¯ x k = x k for k = 1 , ...N . We introduce a vector v and matrix M of dimensions N and N × N respectively v k = Ω k , M jk = Ω jk , j, k = 1 , ..., N (28)( M is the minor 00 of Ω). We also denote the coordinates of the bath ¯ x k = x k , k = 1 , ..., N . Then −
12 ( x Ω x ) −
12 ( x ′ Ω x ′ ) ∗ = −
12 Ω x −
12 Ω ∗ x ′ − ¯ x ( vx + v ∗ x ′ ) − (cid:16) ¯ x ( M + M ∗ )¯ x (cid:17) . The integrations in (27) then give at the initial time ρ ( x , x ′ ) = r Re R − R π exp n − x R − x ′ R ∗ + x x ′ R o , (29)where R = Ω − a, R = | a | , a = (cid:16) v ( M + M ∗ ) − v (cid:17) . (30)The coefficient in(29) can be determined from the obvious property Tr ˆ ρ = 1.The purity µ = Tr ˆ ρ is generally less than unity. We find ρ ( x , x ′ ) = Z dxρ ( x , x ) ρ ( x, x ′ ) = Re R − R π r π Re R exp n − x (cid:16) R − R R (cid:17) − x ′ (cid:16) R ∗ − R R (cid:17) + xx ′ R R o . (31)Taking the trace we find µ = Tr ˆ ρ = r Re R − R Re R + R . (32)At time t we have the same formulas for the reduced matrix ρ ′ in which Ω → Ω ′ and R → R ′ : ρ ′ ( x , x ′ ) = r Re R ′ − R ′ π exp n − x R ′ − x ′ R ′ ∗ + x x ′ R ′ o , (33)where R ′ = Ω ′ − a ′ , R ′ = | a ′ | , a ′ = (cid:16) v ′ ( M ′ + M ′∗ ) − v ′ (cid:17) . (34)We recall that Ω ′ is determined via Ω by (20). Of course Tr ˆ ρ ′ = 1 and the purity µ ′ will be givenby µ ′ = Tr [ˆ ρ ′ ] = s Re R ′ − R ′ Re R ′ + R ′ . (35)Equations (29) – (34) describe evolution of the reduced density matrix in time. In the nextsection we shall study its properties and whether it is Markovian or not and if not to what degreeat different times and properties of system 1. 6 Time evolution of the reduced density matrix
The density matrix ρ is fully determined by the 2 × R = { R , R , R = R , R = R ∗ } .So our strategy will be to fix ˆ R at the initial time t and study if one can uniquely determine the ρ at some later t , that is to find in a unique manner its matrix ˆ R ′ . If this is possible then the processhas the Markovian character. However at time t matrix ˆ R ′ is actually expressed via matrix Ω ′ ,which is determined by Ω by the evolution equation for the whole system. We assume the initialsymmetric matrix Ω to be real. Then with fixed R and R we are left with N ( N − / M and (N-1) parameters in v constrained by the value of ˆ R , that is ( N +2)( N − / − R and R . So it may happen that withfixed ˆ R we shall not be able to find ˆ R ′ in a unique manner, since the latter will depend on extravariables in M and v left after fixing ˆ R and actually depending on the state of the ”bath”. In thiscase the evolution will be non-Markovian. The degree of the dependence on the mentioned extravariables will measure the scale of non-Markovianity. It will possibly depend on time, dynamicproperties of the bath and its coupling to system 0, that is on t and λ , having in mind that thevalue of ω just fixes the energy (and time) units.Our plan is to study this dependence numerically. We shall choose some initial state of thewhole system, that is matrix Ω, determine some of its matrix elements from fixed matrix ˆ R andvary the rest of them. Then we shall study the purity µ ′ at later time t and see if it changes whenˆ R is fixed but the rest of matrix elements in Ω change.We shall also study how this change depends on the value of λ choosing λ = 0 . λ = 1 . ω = 1As already mentioned we choose the initial matrix Ω to be real. We fix the matrix elementΩ = 1 at the initial moment t = 0. The rest N = N − kk k = 1 , N at the initial moment will be allowed to take different values, thus characterizing different states ofthe bath with bath parameters BP . We shall consider five possible sets of these diagonal matrixelements which preserve the average state of the bath, that is Tr Ω, and have the bath parameters BP ( i ), i = 1 , ... BP (1) : Ω kk = 1 , k = 1 , ...N . BP (2) : Ω = 1 . , Ω = 0 . , Ω kk = 1 , k = 3 , ...N . BP (3) : Ω N / ,N / = 1 . , Ω N / ,N / = 0 . , Ω kk = 1 for the rest k . BP (4) : Ω kk = 1 . , k = 1 , ...N / , Ω kk = 0 . k = N / , ...N . BP (5) : Ω kk = 1 . . k from k = 1 , ...N .Our aim will be to study the dependence of the evolution with these varying states of the bath.Note that these states of the bath are not completely equivalent. The subsystem 0 directly interactswith only oscillators 1 and N . So one expects that changing matrix elements Ω and Ω N N willhave greater influence of the bath state than changing other matrix elements or randomizing thechange. In this sense the idea of the ”bath” is better suited to BP (1), BP (3) and BP (5) ratherthan to BP (2) and BP (4).As a convenient signature we take the purity Tr ρ ( t ) during evolution. In the following westudy the bath dependence comparing the purity from the sets BP (1 , ,
5) and BP (1 , .. t = 0) R = 0, µ = 1 (case A) and R = 0 . µ = 0 .
592 (case B). In both cases at the initial moment R = Ω − R = 1 − R and non-diagonal elements Ω k , k = 1 , N are taken equal and adjusted to the fixed value of R .In particular for R = 0 they are all equal to zero.7 µ t 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 0 20 40 60 80 100 µ t Figure 1: The calculated purity µ at different times t for λ = 0 .
1, cases A(left panel) and B( rightpanel), bath parameter BP (1) and number of oscillators N = 101 λ = 0 . Note that as a function of t purity µ ( t ) in both cases oscillates with frequency ∼ .
5. So in theinterval from t = 0 up to t = 100 the plots of µ ( t ) show a lot of oscillations. This is illustrated inFig. 1 in which we show natural values of µ ( t ) for BP (1) and both cases A and B with N = 100Through the maze of oscillation however one can see the general trend of the purity, once aver-aged over the short range oscillations. To make the graphical presentation more convenient in thefollowing we show the purity averaged over these short range oscillation by the Bezier interpolation[11]. Case A
We first consider case A. Recall that then at t = 0 R = 0 and µ (0) = 1, so that the initialstate of system 0 is pure. In Fig. 2 we show in the left panel the averaged µ ( t ) for the set of BP ( i ), i = 1 , .. N = 11 and 101. In the right panel we show the averaged µ ( t ) for the restricted setof BP ( i ), i = 1 , , N .For illustration of the notion of the bath in Fig. 3 we also present the results including N = 3when the ”bath” is reduced to a pair of oscillators. In this case only bath parameters BP (1 ,
2) arepossible.To characterize the bath dependence one may use the dispersion of the values of µ obtained byusing different initial bath data. If we denote µ i as the purity obtained from the bath data BP ( i )then the width of the dispersion in the interval of times t ∈ ∆ t can be defined as w = max t ∈ ∆ t max i = j | µ i − µ j | , (36)where i, j belong to a set of initial conditions BP ( i ) used to find the width. In this way for each N we obtain three different widths w (2) ( N ), w (3) ( N ) and w (5) ( N ) from the sets BP (1 , BP (1 , , BP (1 , ... t in the interval0 < t < µ t 12345 0.94 0.95 0.96 0.97 0.98 0.99 1 1.01 1.02 0 20 40 60 80 100 µ t 135 0.94 0.95 0.96 0.97 0.98 0.99 1 1.01 1.02 0 20 40 60 80 100 µ t 12345 0.94 0.95 0.96 0.97 0.98 0.99 1 1.01 1.02 0 20 40 60 80 100 µ t 135 Figure 2: The calculated purity µ at different times t for λ = 0 .
1, case A, bath parameters BP (1 , ...
5) (left panels) and BP (1 , ,
5) (right panels) and number of oscillators N = 11 (upperpanels) and N = 101 (lower panels). The i -th curve corresponds to BP ( i ) µ t (3)1(11)1(101)1(3)2(11)2(101)2 Figure 3: The calculated purity µ at different times t for λ = 0 .
1, case A, bath parameters BP (1 , N = 3, 11 and 101. Curve ( N ) i corresponds to BP ( i ) with N oscillators.9 able 1 ( µ (0) = 1)∆ t w (2) (3) w (3) (11) w (3) (101) w (5) (11) w (5) (101)[0,20] 0.0615 0.0326 0.0381 0.0655 0.0655[20,40] 0.0613 0.0414 0.0387 0.0762 0.0793[40,60] 0.0591 0.0235 0.0206 0.0707 0.0718[60,80] 0.0618 0.0209 0.232 0.0306 0.0720[80,100] 0.0468 0.0171 0.00898 0.0563 0.0709 Case B
Now we consider case B when at t = 0 R = 0 .
325 and µ = 0 . µ ( t ) for the set of BP ( i ), i = 1 , .. n = 11 and 101. In the right panel we show the averaged µ ( t ) for the restricted set of BP ( i ), i = 1 , , N .For illustration of the notion of the bath in Fig. 5 we again present the results including N = 3when the ”bath” is reduced to a pair of oscillators. In this case only bath parameters BP (1 ,
2) arepossible.Similarly to case A we present the widths of the dispersion (36) in Table 2.
Table 2 ( µ (0) = 0 . t w (2) (3) w (3) (11) w (3) (101) w (5) (11) w (5) (101)[0,20] 0.247 0.132 0.131 0.182 0.155[20,40] 0.261 0.105 0.0577 0.125 0.121[40,60] 0.170 0.0930 0.0412 0.180 0.133[60,80] 0.275 0.0782 0.0364 0.157 0.124[80,100] 0.177 0.112 0.0278 0.128 0.128 λ = 1 . With λ = 1 . t = 0 up to t = 100. Comparing to the previous case with λ = 0 .
01 thistime interval roughly corresponds to a shorter one 0 < t <
67. This is illustrated in Fig. 6 inwhich we show natural values of µ ( t ) for BP (1) and both cases A and B with N = 100. As beforein the following we present the purity averaged over short range oscillation by the Bezier procedure. Case A
Again we start with case A. Recall that then at t = 0 R = 0 and µ (0) = 1, so that the initialstate of system 0 is pure. In Fig. 7 we show in the left panel the averaged µ ( t ) for the set of BP ( i ),10 µ t 12345 0.6 0.7 0.8 0.9 1 1.1 0 20 40 60 80 100 µ t 135 0.6 0.7 0.8 0.9 1 1.1 0 20 40 60 80 100 µ t 12345 0.6 0.7 0.8 0.9 1 1.1 0 20 40 60 80 100 µ t 135 Figure 4: The calculated purity µ at different times t for λ = 0 .
1, case B, bath parameters BP (1 , ... BP (1 , ,
5) (right panels) and number of oscillators N = 11 (upper panels) and N = 101 (lower panels). The i -th curve corresponds to BP ( i ) µ t (3)1(11)1(101)1(3)2(11)2(101)2 Figure 5: The calculated purity µ at different times t for λ = 0 .
1, case B, bath parameters BP (1 , N = 3, 11 and 101. Curve ( N ) i corresponds to BP ( i ) with N oscillators.11 µ t 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 0 20 40 60 80 100 µ t Figure 6: The calculated purity µ at different times t for λ = 1 .
0, cases A(left panel) and B( rightpanel), bath parameter BP (1) and number of oscillators N = 101 Table 3 ( µ (0) = 1)∆ t w (2) (3) w (3) (11) w (3) (101) w (5) (11) w (5) (101)[0,20] 0.0848 0.120 0.117 0.161 0.154[20,40] 0.0848 0.104 0.0239 0.141 0.120[40,60] 0.0846 0.215 0.0241 0.215 0.0946[60,80] 0.0848 0.0967 0.0123 0.115 0.0992[80,100] 0.0845 0.206 0.0226 0.224 0.0801 i = 1 , .. N = 11 and 101. In the right panel we show the averaged µ ( t ) for the restricted setof BP ( i ), i = 1 , , N .For illustration of the notion of the bath in Fig. 8 we once more present the results including N = 3 when the ”bath” is reduced to a pair of oscillators. In this case only bath parameters BP (1 ,
2) are possible.The width w defined by (36) is shown in Table 3. Recall that three different widths w (2) ( N ), w (3) ( N ) and w (5) ( N ) correspond to the sets BP (1 , BP (1 , ,
5) and BP (1 , ... t in the interval 0 < t < Case B
Now we consider case B when at t = 0 R = 0 .
325 and µ = 0 . µ ( t ) for the set of BP ( i ), i = 1 , .. n = 11and 101. In the right panel we show the averaged µ ( t ) for the restricted set of BP ( i ), i = 1 , , N .To illustrate the notion of the bath with its strong coupling λ = 1 . N = 3 when the ”bath” is reduced to a pair of oscillators. In this case onlybath parameters BP (1 ,
2) are possible.The widths of the dispersion (36) for case B and λ = 1 . µ t 12345 0.9 0.92 0.94 0.96 0.98 1 1.02 0 20 40 60 80 100 µ t 135 0.88 0.9 0.92 0.94 0.96 0.98 1 1.02 0 20 40 60 80 100 µ t 12345 0.9 0.92 0.94 0.96 0.98 1 1.02 0 20 40 60 80 100 µ t 135 Figure 7: The calculated purity µ at different times t for λ = 1 .
0, case A, bath parameters BP (1 , ...
5) (left panels) and BP (1 , ,
5) (right panels) and number of oscillators N = 11 (upperpanels) and N = 101 (lower panels). The i -th curve corresponds to BP ( i ) µ t (3)1(11)1(101)1(3)2(11)2(101)2 Figure 8: The calculated purity µ at different times t for λ = 1 .
0, case A(left panel), bath parameters BP (1 ,
2) and number of oscillators N = 3, 11 and 101. Curve ( N ) i corresponds to BP ( i ) with N oscillators. 13 µ t 12345 0.6 0.7 0.8 0.9 1 1.1 0 20 40 60 80 100 µ t 135 0.6 0.7 0.8 0.9 1 1.1 0 20 40 60 80 100 µ t 12345 0.6 0.7 0.8 0.9 1 1.1 0 20 40 60 80 100 µ t 135 Figure 9: The calculated purity µ at different times t for λ = 1 .
0, case B, bath parameters BP (1 , ... BP (1 , ,
5) (right panels) and number of oscillators N = 11 (upper panels) and N = 101 (lower panels). The i -th curve corresponds to BP ( i ) µ t (3)1(11)1(101)1(3)2(11)2(101)2 Figure 10: The calculated purity µ at different times t for λ = 1 .
0, case B, bath parameters BP (1 , N = 3, 11 and 101. Curve ( N ) i corresponds to BP ( i ) with N oscillators.14 able 4 ( µ (0) = 0 . t w (2) (3) w (3) (11) w (3) (101) w (5) (11) w (5) (101)[0,20] 0.201 0.198 0.136 0.234 0.183[20,40] 0.201 0.146 0.0602 0.198 0.170[40,60] 0.201 0.161 0.0543 0.200 0.156[60,80] 0.201 0.170 0.0450 0.187 0.160[80,100] 0.200 0.206 0.0440 0.235 0.151 We start with the weak coupling to the bath usually assumed for the evolution (”the Born-Markovapproximation” [6, 8]). Inspection of our results demonstrated in Figs. 1 -5 and Tables 1. and 2.leads to the following conclusions.In all cases evolution of the purity for the ρ -matrix of the subsystem (system 0) is not Markovian.However deviations from the Markov evolution substantially depend on the property of the bathand elapsed time.If the bath is composed of only two other oscillators 1 and 2 then deviations from the Markovevolutions are quite strong and practically persist at all times from t = 0 to t = 100. The scale ofdeviations measured by the width w is of the order 50%.With the growth of the number of oscillators composing the bath the deviation from the Markovbehavior diminish and they also diminish with time. This is especially visible if the bath changedoes not directly involves the interaction with the subsystem, that is for bath parameters BP (1), BP (3) and BP (5). In this case for the bath with 10 oscillators the deviations from the Markovevolution fall to 10% . With the bath composed of 100 oscillators the deviations fall from 10% at0 < t <
20 to 3% at 80 < t < t = 30 the state becomes practically pure and remains such at latertimes.With a raised coupling constant, that is with a strong bath-system interaction, from Figs.6-10and Tables 3 and 4 we may conclude that on the general the evolution is of the same pattern15s with a weak interaction, with raised widths. So the deviations from the markovian evolutionbecome stronger. Also these deviations do not substantially diminish with time. The purity of thesystem 0 in this case does not generally become restored at large times. With the initial purity µ (0) = 0 .
592 it freezes at the order 0.8. A remarkable exception is the behavior with the bath of101 oscillators without direct contact with the system 0. This behavior is not very different fromthe weak coupling case and leads to the nearly markovian evolution at long times with somewhatraised widths. Also in this case the system tends to the pure one, as with a weak coupling.So in conclusion one finds that the assumption of the Markovian behavior is not a bad ap-proximation provided the bath has a ”macroscopic” character, is randomly connected with thesubsystem and the time of evolution is long enough. This is true both for weak and strong couplingof the system with the bath.Note that in [8, 12] certain criteria for the Markov evolution of the ρ matrix of the open systemwere presented in the form of estimations of the error of the corresponding Lindblad equation.Unfortunately they refer to the case when the bath is described by a bounded operator and socannot be applied to our bath of a set of oscillators. As an alternative in [8] a new coarse grainedLindblad-like equation was proposed for the ρ matrix averaged over some prescribed intervals oftime together with the corresponding error bounds. This time averaging probably is equivalentto our Bezier procedure used for plotting our results but it has no relation to our widths readfrom the results without any averaging. So again it is very difficult to compare the degree of non-Markovianity obtained from this pure theoretical derivation and our numerical results. Note thatthe errors found in [8] in all cases grow exponentially with the time of evolution. In contrary ourwidths in the most beneficial case (with the bath better corresponding to its standardly assumedproperties) steadily diminish with time indicating restoration of markovianity at large enoughtimes as stressed above. So probably the error bounds presented in [8] were too stringent and theapplicability of the Lindblad equation has a wider time limitation.We have to stress that our conclusions have been derived only from the behavior of the purity,which of course is only one of the properties of the ρ -matrix for the subsystem. Still we think thatthe purity gives a very conclusive manifestation of the global behavior of this ρ matrix. We cannotexclude that there exist some observables which are more sensitive to the bath parameters and sonot described by the Markovian evolution. It is not easy to pinpoint such observables apriori. Theycan be found only in the study of some concrete problems, which we postpone for future studies. References [1] Franke V A, Theor. Math. Phys.
406 (1976).[2] Lindblad G, Commun. Math. Phys.
119 (1976).[3] Lindblad G, Rep. Math. Phys.
393 (1976).[4] Gorini V, Kossakowski A and Sudarshan E C G, J. Math. Phys.
821 (1976).[5] P.Pearle, Eur. J.Phys.
805 (2012).[6] C.A.Brasil, F.F.Fanchili, R.de J.Napolitano, arxiv: 1110.2122 [quant-ph], Revista Brasilleirade Ensino de Fisica, (2013) 1303.[7] D.A.Lidar, Z.Bihary, K.B.Whaley, Chem. Phys.,
35 (2001).168] D.A.Lidar, lecture notes , arXiv: 1902.00967 (2019).[9] E.Mozgunov, D.Lidar, Quantum , 227 (2020).[10] R.Jefferson, R.C.Myers, arXiv:1707.08570 [hep-ph].[11] Bartels, R. H., Beatty, J. C., Barsky, B. A. ”Bezier Curves.” in An Introduction to Splines forUse in Computer Graphics and Geometric Modelling . Ch. 10, pp. 211-245, San Francisco, CA:Morgan Kaufmann, 1998.[12] F.Nathan, M.S.Rudner, Phys. Rev.