Efficiency of quantum controlled non-Markovian thermalization

We study optimal control strategies to optimize the relaxation rate towards the fixed point of a quantum system in the presence of a non-Markovian dissipative bath. Contrary to naive expectations that suggest that memory effects might be exploited to improve optimal control effectiveness, non-Markovian effects influence the optimal strategy in a non trivial way: we present a necessary condition to be satisfied so that the effectiveness of optimal control is enhanced by non-Markovianity subject to suitable unitary controls. For illustration, we specialize our findings for the case of the dynamics of single qubit amplitude damping channels. The optimal control strategy presented here can be used to implement optimal cooling processes in quantum technologies and may have implications in quantum thermodynamics when assessing the efficiency of thermal micro-machines.


I. INTRODUCTION
Controlling quantum systems by using time-dependent fields [1] is of primary importance in different branches of science, ranging from chemical reactions [2,3], NMR [4], molecular physics [5] to the emergent quantum technologies [6][7][8]. Investigations on optimal control of open quantum systems mostly focus on memoryless environments [9][10][11] and specifically on those situations where the reduced dynamics can be described by a Markovian master equation of the Lindblad form [12]. In this context, optimal control applications to open quantum systems have been explored in different settings [5,[13][14][15][16][17][18] and recently the ultimate limits to optimal control dictated by quantum mechanics in closed and open systems [19][20][21][22] and the complexity of dealing with manybody systems [23][24][25] have been determined. Timeoptimal quantum control has been extensively discussed for one qubit systems in a dissipative environment [9,26] and the optimal relaxation times determined in [27]. These studies might have both fundamental and practical applications, for example in assessing the ultimate efficiency of quantum thermal machines [28], or to implement fast cooling schemes which have already proven to be advantageous [29,30]. However, introducing a Markovian approximation requires some constraints on system and environment, which may not be valid in general [31,32]. Consequently, incorporating non-Markovian (NM) effects of the environment, in a sense that will be defined more precisely below, might be a necessity in a many experimental situations. Recently, the possible influence of memory effects in the quantum speed of evolution [33] and on quantum control [34] have been analyzed. Here, we present a study of the optimal control strategies to manipulate quantum systems in the presence of NM dissipative baths and compare the performance of optimal control with the case of operating subject to a Markovian (M) environment.
Intuitively, the absence of memory effects in the dynamics of open quantum systems is linked to the possi-bility of identifying well separated time scales in the evolution of system and environment. Recently, a number of proposals have been put forward to quantitatively characterize this effect in terms of explicit non-Markovianity measures [35,36]. In this light, one can define an evolution to be Markovian if described by a quantum dynamical semigroup (time-homogeneous Lindbladian evolution) [37], which would be the traditional Markovianity considered in most previous work on open system control. However, other definitions encompass this as a special case while allowing for more general, non-homogeneous generators, albeit still ensuring the unidirectionality of the system-environment information backflow and therefore the absence of memory effects in the dynamics of the system. Relevant for our analysis is the definition of Markovian evolution in terms of the divisibility of the associated dynamical map [38]. When the dynamics is parametrized using a time-local master equation, the requirement of trace and hermiticity preservation, yields a generator of the forṁ whereL(t) is a time dependent Lindblad superoperator, γ k (t) are generalized (i.e. not necessarily positive) decay rates, the A k (t)'s form an orthonormal basis for the operators for the system, see e.g. Ref. [40] (hereafter has been set equal to one for convenience) and H s (t) is the effective Hamiltonian acting on the system. Equation (1) generalizes the familiar Lindbladian structure to include NM effects while maintaining a time-local structure. However, apart from same special cases, it not known which are the conditions which H s (t), A k (t), and γ k (t) have to satisfy in order to guarantee Complete Positivity (CP) [39][40][41][42][43], i.e. the fundamental prerequisite which under fairly general assumptions is needed to describe a proper quantum evolution [31,32]. In what follows we will focus on a simplified scenario where the γ k (t)'s either are null or coincide with an assigned function γ(t), and where the A k (t)'s are explicitly timeindependent. Accordingly, in the absence of any control Hamiltonian applied during the course of the evolution, we assume a dynamical evolution described by the equatioṅ where L is a (time-independent) Lindblad generator characterized by having a unique fixed point ρ f p (i.e. L(ρ) = 0 iff ρ = ρ f p ). For this model, in the absence of any Hamiltonian term (i.e. H s (t) = 0) CP over a time interval [0, T ] is guaranteed when [39] t 0 while divisibility (i.e, Markovianity) is tantamount to the positivity of the single decay rate at all times [40]: if there exists a time interval where γ(t) becomes negative, the ensuing dynamics is no longer divisible and the evolution is NM. In this context we will assume a control Hamiltonian H s (t) to represent time-localized infinitely strong pulses, which induce instantaneous unitary transformations at specific control times. This corresponds to are time independent operators which act impulsively on the system at t = t j (δ(t) being the Dirac delta-function), at which instants one can neglect the contribution from the non-unitary part, and represent the master equation bẏ . Therefore the resulting dynamics is described by a sequence of free evolutions induced by the noise over the intervals t ∈ [t j , t j+1 ] interweaved with where D j = exp tj+1 tjL (t) , U(· · ·) = U (· · ·)U † and "•" is the composition of super-operators. When only two control pulses are applied (the first U in at the very beginning and the second U f in at the very end of the temporal evolution), the non-unitary evolution is described by Eq.
(2) and CP of the trajectory (4) is automatically guaranteed by Eq. (3), the scenario corresponding to the realistic case where one acts on the system with very strong control pulses at the state preparation stage and immediately before the measuring stage. When more U j 's are present, the situation however becomes more complex. There is no clear physical prescription which one can follow to impose the associated dynamics on the system at least when the dissipative evolution is assumed to be NM. Consider for instance the case U f in •D 1 •U 1 •D 0 •U in (ρ(0)) where U 1 is a non trivial unitary. Even admitting that the latter is enforced by applying at time t 1 a strong instantaneous control pulse, there is absolutely no clear evidence that the open dynamics for t ≥ t 1 should be still described by the same generalized Lindbladian γ(t)L, the system environment being highly sensitive to whatever the system itself has experienced in its previous history.
Keeping in mind the above limitations, in this work we attempt for the first time a systematic study of optimal control protocols which would allow one to speed up the driving of a generic (but known) initial state ρ(0) toward the fixed point ρ f p of the bare dissipative evolution for the model of Eq. (2) which explicitly includes NM effects. We arrive at the quantum speed limit times when application of only two control pulses U in and U f in , at initial and final times respectively, is enough to follow the optimal trajectory. On the other hand, we present lower bounds to the same when optimal control strategy demands unitary pulses at intermediate times as well. We show that the efficiency of optimal control protocols is not determined by the M/NM divide alone but it depends drastically on the behaviour of the NM channel: if the system displays NM behavior before reaching the fixed point for the first time, NM effects might be exploited to obtain an increased optimal control efficiency as compared to the M scenario. On the contrary, NM effects are detrimental to the optimal control effectiveness if information back-flow occurs only after the system reaches the fixed point (see Fig. 1). These results are valid irrespective of the detailed description of the system, i.e. its dimension, Hamiltonian, control field, or the explicit form of the dissipative bath.

II. THE MODEL
The divisibility measure for the model Eq. (2) is equivalent to the characterization of memory effects by means of the time evolution of the trace distance [44]. This provides an intuitive characterization of the presence of memory effects in terms of a temporary increase in the distinguishability of quantum states as a result of an information back-flow from the system and into the environment that is absent when the evolution is divisible [45]. As a result, a divisible evolution for which the single decay rate γ(t) ≥ 0 at all times will exhibit a monotonic decrease of the trace distance of any input state towards a (assumed to be unique) fixed point ρ f p of the Lindblad generator L [46]. On the contrary, as illustrated in Fig. 1, the behaviour of the trace distance can be non-monotonic when the dynamics is NM. In this case, there exist time interval(s) where γ(t) becomes negative. Denoting by d(t) = ||ρ(t) − ρ f p || the trace distance between ρ(t) and the fixed point ρ f p , it straightforwardly follows thatḋ(t) ≤ 0 ∀t in the M limit. Looking at this quantity one can classify NM dynamics into two distinct classes (see Fig. 1): the first one (Class A) is defined by those dynamics where the system reaches the fixed point at time T F before γ(t) changes sign i.e., γ(t) ≥ 0 anḋ d(t) ≤ 0 for 0 ≤ t < T F . In this case, the NM dynamics reaches the fixed point ρ f p and then start to oscillate. On the other hand, Class B dynamics is characterized by γ(t) that changes sign (and correspondinglyḋ(t) > 0) at some time t < T F , that is the solutions of the equation γ(t s ) = 0 are such that t s < T F for at least one s. In contrast, in the M dynamics d(t) always decreases monotonically and asymptotically to d(t → ∞) = 0.
NM channels of class A/B arise from different physical implementations. As an illustration, the damped Jaynes Cummings model exemplifies a Class A dynamics. Here a qubit is coupled to a single cavity mode which in turn is coupled to a reservoir consisting of harmonic oscillators in the vacuum state (see Eq. (6)) [31,[47][48][49]. On the other hand, dynamics similar to Class B can arise for example in a two level system in contact with an environment made of another two level system, as realized recently in an experimental demonstration of NM dynamics [50].
As we will see hereafter, the difference between Class A and B appears to drastically affect the performance of any possible optimal control strategy to improve the speed of relaxation of the system towards the fixed point.
We assume full knowledge of the initial state and we allow for an error tolerance of 0 < ǫ ≪ 1, considering that the target is reached whenever the condition |d(t)| ≤ ǫ is satisfied. To obtain a lower bound on the minimum time T QSL needed to fulfill such constraint we restrict our analysis to the ideal limit of infinite control which allows us to carry out any unitary transformation instantaneously along the lines of (and with all the limitations associated with) the formalism detailed in Eq. (4). In the limit of infinite control an important role is played by the Casimir invariants Γ j (j = 2, 3, ..., N for a N level system). The Casimir invariants of a state ρ are related to the trace invariants Tr(ρ j ) (j = 2, 3, ..., N ) and they cannot be altered by unitary transformation alone [10,51]. For example, a two level system has a single Casimir invariant -its purity P = Tr ρ 2 , which remains unchanged under any unitary transformation. Consequently, any optimal strategy with the controls restricted to unitary transformations only, would be to reach a state ρ characterired by all Casimir invariants same as those of ρ f p in the minimum possible time. Following this we can apply a unitary pulse to reach the fixed point instantaneously. Clearly, any constrained control will at most be as efficient as the results we present hereafter, based on the analysis we have presented previously for the case of M dynamics [5,27].
In what follows, we will analyse Class A and B channels independently.
Class A: As shown in Fig. 1, in the NM regime d(t) goes to zero at t = T F when ρ = ρ f p and L(ρ f p ) = 0. At the same time we expect |γ(t)| → ∞ at t ≈ T F in order to have finiteρ(t) = γ(t)L(ρ(t)) even for L(ρ(t)) ≈ L(ρ f p ) = 0, as is required for a non-monotonic d of the form shown in Fig. 1. Notice that γ(t) and hence the time t = T F at which γ(t) → ∞ are in general independent of ρ. Consequently any optimal control protocol which involves unitary transformation of ρ(t) generated by H s (t) at earlier times t < T F followed by non-unitary relaxation to ρ f p is expected to be ineffective in this case and we have T QSL = T F . That is, the gain (or efficiency) of optimal trajectory in the NM class is R A N M = T F /T QSL = 1. One can easily see T F /T QSL = 1 implies absence of any speed up, whereas any advantage one gains by optimal control can be quantified by T F /T QSL > 1. On the other hand in the M limit γ(t) = γ 0 is finite and constant, and the system relaxes asymptotically to the fixed point in the absence of any control. In this case we introduce an error tolerance ǫ ≪ 1, such that we say the target state is reached if |d(T F )| ≤ ǫ. Clearly, T F increases with decreasing ǫ diverging to T F → ∞ in the limit ǫ → 0, as can be expected for finite γ 0 . Therefore the above argument of |γ(t)| → ∞ at t ≈ T F does not apply in this case and in general one can expect the time of evolution to depend on the initial state. Consequently the quantum speed up ratio R A M can exceed R A N M ≈ 1, as is explicitly derived below in the case of a two level system in presence of an amplitude damping channel. Similar arguments apply also in the case when an additional unitary transformation is needed at the end of the evolution to reach ρ f p , where R A M → ∞ for ǫ → 0 [27]. Class B: Here we focus on systems of Class B where as already mentioned γ(t) changes sign for t s < T F with s = 1, . . . , N s . Clearly, in this case |γ(t)| does not necessarily diverge for any t. Consequently the arguments presented above for class A fails to hold any longer and the time of relaxation to the fixed point can in general be expected to depend on ρ(t) (and hence on H s ). Furthermore, it might be possible to exploit the NM effects such that even thoughḋ(t) > 0 for t 1 ≤ t ≤ t 2 one can, by ap-plication of optimal control, make sure thatΓ j (t) > 0 and maximum ∀ t, j (where we have assumed Γ j (t = 0) < Γ jf ∀ j and Γ jf denotes the jth Casimir invariants for the fixed point ρ f p ). This presents the possibility of exploiting NM effects to achieve better control as opposed to the M dynamics, as is presented below for the case of a two level system in the presence of an amplitude damping channel. However we stress that this is not a general result and explicit examples can be constructed where this is actually not true.

A. Generalized amplitude damping channel
Let us now analyze in detail the generic formalism outlined above for the specific case of a two level system in contact with NM amplitude damping channels of the two classes introduced before.
We consider the non-unitary dissipative dynamics described by the time local master equation acting on a 2×2 reduced density matrix ρ(t) of a qubit and we consider the time independent Linbladian L given by with σ ± being the raising/lowering qubit operators. The system evolution is given by Eq. 2, and we will focus on two different functional dependence of the parameter γ(t) corresponding to the Class A and B dynamics. We will analyze the system evolution following the Bloch vector r representing the state ρ = (I + r. σ) /2 inside the Bloch sphere, where the unitary part of the dynamics generated by H s induce rotations, thus preserving the purity P = (1+| r| 2 )/2. In contrast, in general the action of the noise is expected to modify the purity as well. Class A: An example of this class of dynamics is obtained under the assumption In the above expression λ and γ 0 are two positive constants whose ratio determines the bath behavior: λ > 2γ 0 corresponds to a M bath, whereas g becomes imaginary in the limit γ 0 > λ/2 resulting in NM dynamics. In the NM limit of γ 0 ≫ λ, 1 the bath time scale is determined by the product λγ 0 and is independent of the specific form of the super-operator L. It can be easily seen that γ(t) increases monotonically from 0 to γ(t) → ∞ at lim γ0/λ→∞ where T F is independent of the initial state and the system reaches the fixed point when γ(t) diverges. With this choice of γ(t) in Eq. (6) the time scale is given by √ 2λγ 0 in the NM limit γ 0 ≫ λ, while γ(t) ≈ γ 0 sets the time scale in the M limit λ ≫ γ 0 . Therefore the time taken to reach the fixed point can be expected to decrease as 1/ √ λγ 0 in the NM limit while it scales as ∼ 1/γ 0 in the M limit.
As mentioned above, this form of γ(t) arises in the damped Jaynes-Cummings model. Here we restrict ourselves to a single excitation in the qubit-cavity system. In this context the parameter λ in γ(t) denotes the spectral width of the coupling to the reservoir, while γ 0 characterizes the strength of the coupling.
In the absence of any control the qubit relaxes to a fixed point ρ f p characterized by the Bloch vector 0, 0, 1−e β 1+e β and the optimal control we analyze here aims to accelerate the relaxation towards this state with unconstrained unitary control. Following the strategy proposed above, we look for the extremal speed v = ∂P/∂t of purity change for every r. For this model, the speed of change of purity is given by v(r, θ, t) = −γ(t)(e β − 1)r cos θ + r 2r f p 1 + cos 2 θ : note that a positive (negative) v denotes increasing (decreasing) purity. The two strategies differ slightly in case of cooling or heating (i.e. the final purity is lower or higher than the initial one); but both cases correspond to applying unitary rotations at the beginning and at the end (for heating) of the dynamical evolution, thus yielding a trajectory of the form Eq. (4) which is fully compatible with the CP requirement and which doesn't pose any problem in terms of physical implementation (see discussion in Sec. II). Specifically we need to apply unitary control so that the system evolves along θ = π till the final purity is reached in the case of cooling, while θ = 0 is the optimal path in the case of heating [27], in agreement with a recent work on quantum speed limit in open quantum systems [33] (see Appendix IV for details). Our analysis clearly shows that T cool QSL decreases with increasing λ as ∼ 1/ √ λ for small λ, finally saturating to λ independent constant values in the M limit (λ ≫ 1). However, it would be misleading to conclude about the role of Markovianity on T QSL from this alone, since both T cool QSL and T heat QSL decrease with increasing γ 0 as well. In particular, they scale as ∼ 1/γ 0 in the M limit of small γ 0 , while the scaling changes to 1/ √ γ 0 for λ/γ 0 → 0. Indeed, the behavior depends more on the specific path in the (γ 0 , λ) plane rather than whether the system is M or NM (see Fig. 2a). However, if one analyzes the speedup obtained by means of optimal control strategies, the scenario changes: in this case (as shown in Fig. 2b) the ratio R A = T F /T cool QSL clearly distinguishes between M and NM dynamics with typical limiting values given by lim ǫ→0,λ/γ0→∞ R A M → 2 and lim ǫ→0,λ/γ0→0 R A N M → 1 (see Appendix for details). Class B: Finally, we investigate a particular case belonging to the class B dynamics and compare it to the previous case. As in the previous case we consider a time evolution described by a master equation Eq. (2) with L given by Eq. (5); however for our present purpose we formulate a γ(t) given by with ζ, Ω being two positive constants satisfying the CP condition Eq. (3). With this choice, in the absence of the control Hamiltonian, NM effects manifest themselves for (2n + 1)π/2 < Ωt < (2n + 3)π/2 for integer n ≥ 0 as γ(t) changes sign at Ωt = (2n + 1)π/2, simultaneously altering the sign ofḋ(t) toḋ(t) > 0. With a proper choice of parameters one can make γ(t) (and hence d(t)) exhibit oscillatory dynamics for 1 > d > 0. As for the previous example, also in this case the extremals of v are independent of γ(t) and determined by L(ρ(t)) only. Therefore they occur at exactly the same points as for class A (6), i.e., at θ = 0, π and arccos(r/r f ). The unconstrained optimal strategies are now modified as follows. In the case of cooling, the optimal strategy is to follow the path θ = π for 0 ≤ t < T cool QSL,N M , during which time the purity increases monotonically, where we have assumed the system reaches the target at t = T cool QSL,N M < π/(2Ω) for simplicity. Consequently the optimal strategy demands a single pulse at time t = 0 only to make θ = π, which corresponds to an evolution operator of the form U f in • D • U in . Clearly, this evolution is CPT as already discussed in section I with D depending on γ t and L (see Eqs. (5) and (8)) and is thus possible to implement physically. In Fig. 3(a) we report the speedup obtained by such optimal strategy for different values of Ω in Eq. (8), where we have taken ǫ = 0.01 large enough so that ΩT cool QSL,N M < π/2. In this case we arrive at the M limit by setting Ω = 0; as can be clearly seen the speed up in the NM limit is such that R B N M > R B M , ∀ Ω > 0 showing that there exist scenarios where NM effects can be exploited to improve the control effectiveness. On the contrary, if the system does not reach the target for ΩT cool QSL,N M < π/2, the optimal strategy changes: at Ωt = π/2, γ(t) and hence v change sign, leading to decrease (increase) of purity for θ = π (θ = 0). Interestingly, we can take advantage of this effect by making θ = 0 at t = π/(2Ω), where v exhibits a maximum for γ t < 0. As mentioned before, application of a unitary pulse during the course of an evolution may change the form of γ t and L (see Eq. (4)). However, for simplicity let us assume v changes sign at t = (2n + 1)π/2 and assumes extremum values at θ = 0, π and arccos (−r/r f p ) even in presence of unitary control. One can easily extend our analysis to a more generic case where the simplifying assumptions do not hold by following the path of maximum (minimum) v for cooling (heating). Let us first consider the case where the system reaches the target r z = r f − ǫ at t = T cool QSL,N M < 3π/(2Ω). In such a scenario, as depicted in Fig. 3 (b), we let the system evolve freely for π/(2Ω) < t ≤ T cool QSL,N M , following which we take the system to θ = π and r z = −r f + ǫ, thus obtaining the desired goal. Clearly, an optimal path exists in case of the class B non-Markovian channel, which if possible to be followed by application of suitable unitary controls, helps in cooling and in particular might make it possible to reach the fixed point in finite time (if ǫ = 0). Generalization to the case ΩT cool QSL,N M > 3π/2 where multiple π rotations are needed is straightforward. However, we emphasize that the strategy presented above for ΩT cool QSL,N M > π/2 follows an evolution operator of the form Eq. (4) with unitary pulses applied at intermediate times. Consequently our analysis gives a lower bound to T cool QSL,N M only, achieved by following the optimal path shown in Fig. 3(b), for which at present we do not have any implementation strategy.
Finally, we address the problem of heating the system in the shortest possible time, which amounts to minimizing v ∀ t. Therefore, the optimal path dictates to set θ = 0 at t = 0 and then let it evolve freely till Ωt = π/2, where γ(t) and hence v change sign. Unfortunately in contrast to the cooling problem, now v > 0 ∀θ making it impossible to take advantage of the NM effects to accelerate the evolution. However, even in this case one can always minimize the unwanted effect of information backflow for γ(t) < 0 by increasing θ to θ = arccos(−r f p /r) where v has a minimum.

III. CONCLUSION
We have studied the effectiveness of unconstrained optimal control of a generic quantum system in the presence of a non-Markovian dissipative bath. Contrary to common expectations, the speedup does not crucially depend on the Markovian versus non-Markovian divide, but rather on the specific details of the non-Markovian evolution. We showed that the speed up drastically depends on whether the system dynamics is monotonic or not before reaching the fixed point for the first time, as determined by the trace distance to the fixed point (class A and class B dynamics respectively). Indeed, in the former case, the speed up obtained via optimal control is always higher in presence of a Markovian bath as compared to a non-Markovian one, while the reverse can be true in the latter case. Finally, we have presented some specific examples of these findings for the case of a two level system subject to an amplitude damping channel. Note that, in the more realistic scenario where constraints are present, the presented results serve as theoretical bounds to the optimal control effectiveness.