Adding dynamical generators in quantum master equations

TL;DR

This paper investigates the conditions under which dynamical generators in quantum master equations can be added without losing physical validity, highlighting limitations and providing criteria for justified addition.

Contribution

It establishes a general condition for the legitimate addition of generators in quantum master equations, especially under weak coupling and commuting Hamiltonian and interactions.

Findings

Generators can be added without losing physicality only under restrictive conditions.

Naively adding generators may not reflect the true system evolution in general.

Counterexamples show the exact dynamics often differ from the generator sum approximation.

Abstract

The quantum master equation is a widespread approach to describing open quantum system dynamics. In this approach, the effect of the environment on the system evolution is entirely captured by the dynamical generator, providing a compact and versatile description. However, care needs to be taken when several noise processes act simultaneously or the Hamiltonian evolution of the system is modified. Here, we show that generators can be added at the master equation level without compromising physicality only under restrictive conditions. Moreover, even when adding generators results in legitimate dynamics, this does not generally correspond to the true evolution of the system. We establish a general condition under which direct addition of dynamical generators is justified, showing that it is ensured under weak coupling and for settings where the free system Hamiltonian and all system-environment interactions commute. In all other cases, we demonstrate by counterexamples that the exact evolution derived microscopically cannot be guaranteed to coincide with the dynamics naively obtained by adding the generators.

Equations

$$\frac{d}{dt}{\rho }_S(t)={\mathcal{L}}_t[{\rho }_S(t)]={\mathcal{H}}_t[{\rho }_S(t)]+{\mathcal{D}}_t[{\rho }_S(t)],$$

$${\mathcal{L}}_t[\rho ]=-\mathrm{i}[H(t),\rho ]+\sum _{i,j=1}^{d^2-1}{\mathsf{D}}_{ij}(t)\left(F_i\rho F_j^{†}-\frac{1}{2}\{F_j^{†}{F}_i,\rho \}\right),$$

$${\int }_0^{t}d\tau \mathsf{D}(\tau )\ge 0.$$

$${\mathcal{L}}_t^{(1)}[\rho ]$$

$${\mathcal{L}}_t^{(2)}[\rho ]$$

$$\frac{d}{dt}{\rho }_S(t)=-\mathrm{i}{\mathrm{Tr}}_E\left[H_S+H_E+H_I,{\rho }_{SE}(t)\right],$$

$$\frac{d}{dt}{\bar{\rho }}_S(t)=-{\int }_0^{t}ds{\mathrm{Tr}}_E[{\bar{H}}_I(t),[{\bar{H}}_I(s),{\bar{\rho }}_{SE}(s)]],$$

$$\frac{d}{dt}{\rho }_S(t)=-\mathrm{i}[H_S,{\rho }_S(t)]+{\tilde{\mathcal{L}}}_t\left[{\rho }_S(t)\right].$$

$$\frac{d}{dt}{\rho }_S(t)=-\mathrm{i}\left[H_S,{\rho }_S(t)\right]+\sum _i{\tilde{\mathcal{L}}}_t^{(i)}\left[{\rho }_S(t)\right]-\sum _{i\ne j}{\int }_0^{t}ds{e}^{-\mathrm{i}{H}_S(t-s)}{\mathrm{Tr}}_{E_{ij}}\left[{\bar{H}}_{I_i}(t-s),\left[H_{I_j},{\rho }_{S{E}_{ij}}(s)\right]\right]e^{\mathrm{i}{H}_S(t-s)},$$

$${\rho }_{SE}(t)\approx {\rho }_S(t)\otimes \underset{i}{⨂}{\varrho }_{E_i}(t),$$

$$U_{SE}(t)=e^{-\mathrm{i}(H_S+{\sum }_i{H}_{E_i}+H_{I_i})t}=e^{-\mathrm{i}{H}_{S}t}\prod _i{e}^{-\mathrm{i}(H_{I_i}+H_{E_i})t}.$$

$$\hat{m}=\sum _{n=1}^N{\sigma }_z^{(n)}=\sum _{k=0}^N{m}_k{\Pi }_k,$$

$$m_k=-N+2k\text{for}k=0,...,N.$$

$${\rho }_E=\sum _{k=0}^N{q}_k{\Pi }_k.$$

$$p(m_k)=\mathrm{Tr}[{\rho }_E{\Pi }_k]=q_k\mathrm{Tr}{\Pi }_k.$$

$$\sum _{k=0}^N(m_k{)}^{s}p(m_k)\underset{N\to \infty }{⟶}{\int }_{-\infty }^{\infty }dm{m}^{s}p(m).$$

$$H_I=A\otimes \hat{m}$$

$${\mathrm{Tr}}_E\{H_I{\rho }_E\}=A\mathrm{Tr}\{\hat{m}{\rho }_E\}=A\sum _{k=0}^N{m}_{k}p(m_k)=0,$$

$$p(m)=\frac{1}{\sqrt{2\pi }\sigma }{e}^{-\frac{m^2}{2{\sigma }^2}},$$

$$p(m)=\frac{\lambda }{\pi ({\lambda }^2+m^2)},$$

$${\rho }_{SE}(t)=\sum _k{q}_k{\rho }_S^{(k)}(t)\otimes {\Pi }_k,$$

$${\rho }_S(t)={\mathrm{Tr}}_E{\rho }_{SE}(t)=\sum _{k}p(m_k){\rho }_S^{(k)}(t).$$

$${\dot{\rho }}_{SE}(t)$$

$$=-\mathrm{i}[H_S\otimes \text{openone}+A\otimes \hat{m},{\rho }_{SE}(t)]$$

$$=-\mathrm{i}[H_S\otimes \text{openone}+\sum _k{m}_{k}A\otimes {\Pi }_k,\sum _l{q}_l{\rho }_S^{(l)}(t)\otimes {\Pi }_l]$$

$$=-\mathrm{i}\sum _k{q}_k[H_S+m_{k}A,{\rho }_S^{(k)}(t)]\otimes {\Pi }_k.$$

$${\dot{\rho }}_S^{(k)}(t)=-\mathrm{i}[H_S+m_{k}A,{\rho }_S^{(k)}(t)],$$

$${\rho }_S(t)={\Lambda }_t[{\rho }_S(0)]=\sum _{k}p(m_k)U_S^{(k)}(t){\rho }_S(0)U_S^{(k)†}(t).$$

$${\dot{\mathbf{r}}}^{(k)}(t)\cdot \mathbf{\sigma }=-\mathrm{i}[H_S+m_{k}A,{\mathbf{r}}^{(k)}(t)\cdot \mathbf{\sigma }]$$

$$\mathbf{r}(t)={\mathbf{D}}_t\mathbf{r}(0)=\left[\sum _{k}p(m_k){\mathbf{R}}^{(k)}(t)\right]\mathbf{r}(0).$$

$${\mathbf{D}}_t:=\sum _{k}p(m_k){\mathbf{R}}^{(k)}(t)\underset{N\to \infty }{⟶}{\int }_{-\infty }^{\infty }dmp(m)\mathbf{R}(m,t),$$

$${\tilde{\mathbf{L}}}_t:={\mathbf{L}}_t-{\dot{\mathbf{R}}}_S(t){\mathbf{R}}_S(t{)}^{-1}={\mathbf{R}}_S(t){\bar{\mathbf{L}}}_t{\mathbf{R}}_S(t{)}^{-1},$$

Full text

I Introduction

It is generally impossible to completely isolate a small system of interest from the surrounding environment. Thus, dissipative effects caused by the environment are important in almost every quantum experiment, ranging from highly controlled settings, where much effort is invested in minimising them, to areas where the dissipation is the key object of interest. In many cases, exact modelling of the environment is not practical and its effect is instead accounted for by employing effective models describing the induced noise. Different approaches exist, e.g., quantum Langevin and stochastic Schrödinger equations Gardiner and Zoller (2004); Breuer and Petruccione (2002), quantum jump and state-diffusion models Percival (1998); Plenio and Knight (1998), or Hilbert-space averaging methods Gemmer and Michel (2006).

Arguably, the most widely applied approach is to use the quantum master equation (QME) description Gardiner and Zoller (2004); Breuer and Petruccione (2002). In this approach, the system evolution is given by a time-local differential equation, where the effect of the environment is captured by the dynamical generator. A master equation can be derived from a microscopic model of the system and environment, and their interaction, by tracing over the environment and applying appropriate approximations Gardiner and Zoller (2004); Breuer and Petruccione (2002). However, QMEs are also often applied directly, without explicit reference to an underlying model. In that case, care needs to be taken when several noise processes act in parallel, as simultaneous coupling to multiple baths in a microscopic model does not generally correspond to simple addition of noise generators. Moreover, when the Hamiltonian evolution of the system is modified, e.g., when controlling system dynamics by coherent driving Schmidt et al. (2011), the form of noise generators in a QME may significantly change. Additivity of noise at the QME level has been discussed recently for qubits when analysing dynamical effects of interference between different baths Chan et al. (2014); Mitchison and Plenio (2018), non-additivity of relaxation rates in multipartite systems Yu and Eberly (2006); Lankinen et al. (2016), as well as in the context of charge (excitation) transport Giusteri et al. (2017).

In this work, we address the questions of when:

• (i)

The naive addition of generators yields physically valid dynamics.

• (ii)

The corresponding evolution coincides with the true system dynamics derived from the underlying microscopic model.

First, we show that (i) is satisfied for generators which are commutative, semigroup-simulable (can be interpreted as a fictitious semigroup at each time instance), and preserve commutativity of the dynamics under addition. These reach beyond the case of Markovian generators for which (i) naturally holds. Outside of this class, we find examples of simple qubit QMEs which lead to unphysical dynamics. We observe that (ii) holds if and only if the cross-correlations between distinct environments can be ignored within a QME. We show this to be the case in the weak-coupling regime, extending previous results in this direction Cohen-Tannoudji et al. (1998); Chan et al. (2014); Schaller (2015). We also provide a sufficient condition for (ii) dictated by the commutativity of Hamiltonians at the microscopic level. We combine these generic considerations with a detailed study of a specific open system, namely a qubit interacting simultaneously with multiple spin baths, for which we provide examples where (ii) is not satisfied, while choosing the microscopic Hamiltonians to fulfil particular commutation relations.

Our results are of relevance to areas of quantum physics where careful description of dissipative dynamics plays a key role, e.g., in dissipative quantum state engineering Diehl et al. (2008); Verstraete et al. (2009); Metelmann and Clerk (2015); Reiter et al. (2016), dissipative coupling in optomechanics Aspelmeyer et al. (2014), or in dissipation-enhanced quantum transport scenarios Gurvitz and Prager (1996); Giusteri et al. (2017), including biological processes Lambert et al. (2013). In particular, they are of importance to situations in which QMEs are routinely employed to account for multiple sources of dissipation, e.g., in quantum thermodynamics Alicki (1979); Kosloff and Levy (2014); Vinjanampathy and Anders (2016); Goold et al. (2016) when dealing with multiple heat baths Skrzypczyk et al. (2011); Correa et al. (2013); Levy and Kosloff (2014); Mitchison and Plenio (2018) or in quantum metrology Maccone and Giovannetti (2011); Escher et al. (2011); Demkowicz-Dobrzański et al. (2012) where the relation between dissipation and Hamiltonian dynamics, encoding the estimated parameter, is crucial Chaves et al. (2013); Brask et al. (2015); Smirne et al. (2016); Haase et al. (2017).

The manuscript is structured as follows. In Sec. II, we discuss QMEs at an abstract level—as defined by families of dynamical generators whose important properties we summarise in Sec. II.1. We specify in Sec. II.2 conditions under which the addition of physically valid generators is guaranteed to yield legitimate dynamics. We demonstrate by explicit examples that even mild violation of these conditions may lead to unphysical evolutions.

In Sec. III, we view the validity of QMEs from the microscopic perspective. In particular, we briefly review in Sec. III.1 the canonical derivation of a QME based on an underlying microscopic model, in order to discuss the effect of changing the system Hamiltonian on the QME, as well as the generalisation to interactions with multiple environments. We then formulate a general criterion for the validity of generator addition in Sec. III.2, which we explicitly show to be ensured in the weak coupling regime, or when particular commutation relations of the microscopic Hamiltonians are fulfilled.

In Sec. IV, we develop an exactly solvable model of a qubit interacting with multiple spin baths, which allows us to explicitly construct counterexamples that disprove the microscopic validity of generator addition in all the regimes in which the aforementioned commutation relations do not hold. Finally, we conclude in Sec. V.

II Time-local quantum master equations

QMEs constitute a standard tool to describe reduced dynamics of open quantum systems. They provide a compact way of defining the effective system evolution at the level of its density matrix, $\rho_{S}(t)$, without need for explicit specification neither of environmental interactions nor the nature of the noise. Although a QME may be expressed in a generalised form as an integro-differential equation involving time-convolution Vacchini (2016), its equivalent (c.f. Chruściński and Kossakowski (2010)) and more transparent time-local formulation is typically favoured, providing a more direct connection to the underlying physical mechanisms responsible for the dissipation Gardiner and Zoller (2004); Breuer and Petruccione (2002). Given a time-local QME:

[TABLE]

all the information about the system evolution is contained within the dynamical generator, $\mathcal{L}_{t}$, that is uniquely defined at each moment of time $t$. Moreover, $\mathcal{L}_{t}$ can always be decomposed into its Hamiltonian and purely dissipative parts, i.e., $\mathcal{L}_{t}=\mathcal{H}_{t}+\mathcal{D}_{t}$ in Eq. (1) with $\mathcal{H}_{t}[\rho]=-\mathrm{i}[H(t),\rho]$ and some Hermitian $H(t)$ Gorini et al. (1976).

Although the QME (1) constitutes an ordinary differential equation, the system evolution may exhibit highly non-trivial memory features thanks to the arbitrary dependence of $\mathcal{L}_{t}$ on the local time-instance $t$, but also on the (fixed) initial time $t_{0}$ at which the evolution commences Chruściński and Kossakowski (2010)—which, without loss of generality, we choose to be zero ($t_{0}=0$) and drop throughout this work.

II.1 Physicality of dynamical generators

For the QME (1) to be physically valid, it must yield dynamics that is consistent with quantum theory. In particular, upon integration the QME must lead to a family of (dynamical) maps $\Lambda_{t}$ (parametrised by $t$) that satisfy ${\rho_{S}}(t)=\Lambda_{t}[{\rho_{S}}(0)]$ for any $t\geq 0$ and initial ${\rho_{S}}(0)$, with each $\Lambda_{t}$ being completely positive and trace preserving (CPTP) Jamiołkowski (1972); Choi (1975).

On the other hand, any QME (1) is unambiguously specified by the family of (dynamical) generators $\mathcal{L}_{t}$ appearing in Eq. (1). However, as discussed in App. A, although the CPTP condition can be straightforwardly checked for maps $\Lambda_{t}$, it does not directly translate onto the generators $\mathcal{L}_{t}$. As a result, for a generic QME its physicality cannot be easily inferred at the level of Eq. (1), unless its explicit integration is possible. Nevertheless, we formally call a family of dynamical generators $\mathcal{L}_{t}$ physical if the family of maps it generates consists only of CPTP transformations. In what follows (see also App. A.1), we describe properties of dynamical generators that ensure their physicality.

Any family of dynamical generators, whether physical or not, can be uniquely decomposed as Gorini et al. (1976)

[TABLE]

where $d$ is the Hilbert space dimension and $\{F_{i}\}_{i=1}^{d^{2}}$ is any orthonormal operator basis with $\operatorname{Tr}\{F_{i}^{\dagger}F_{j}\}=\delta_{ij}$ and all $F_{i}$ traceless except $F_{d^{2}}=\openone/\sqrt{d}$. The Hamiltonian part, $\mathcal{H}_{t}$, of the generator in Eq. (1) is then determined by $H(t)$ of Eq. (2), while the dissipative part $\mathcal{D}_{t}$ is defined by the Hermitian matrix $\mathsf{D}(t)$. Although general criteria for physicality of dynamical-generator families are not known, two natural classes of physical dynamics can be identified based on the above decomposition.

In particular, when $\mathsf{D}(t)$ is positive semidefinite, $\mathcal{L}_{t}$ is said to be of Gorini-Kossakowski-Sudarshan-Lindblad (GKSL) form Gorini et al. (1976); Lindblad (1976). If this is the case for all $t\geq 0$, then the corresponding evolution is not only physical but also CP-divisible, i.e., the corresponding family of maps can be decomposed as $\Lambda_{t}=\tilde{\Lambda}_{t,s}\Lambda_{s}$, where $\tilde{\Lambda}_{t,s}$ is CPTP for all $0\leq s\leq t$. This property is typically associated with Markovianity of the evolution Rivas et al. (2014); Breuer et al. (2016); de Vega and Alonso (2017).

Furthermore, when, in addition, $H$ and $\mathsf{D}$ in Eq. (2) are time-independent, the dynamics forms a semigroup, such that the generator and map families are directly related via $\Lambda_{t}=\exp[t\mathcal{L}]$ with all $\mathcal{L}_{t}=\mathcal{L}$ Alicki (2002). See App. A.1 for a more detailed discussion of different types of evolutions.

For the purpose of this work, we also identify another important class of physical dynamics:

Definition 1**.**

A given dynamical family $\Lambda_{t}$ is semigroup-simulable (SS) if for any $t\geq 0$ the map $\mathcal{Z}_{t}=\log\Lambda_{t}$ is of the GKSL form (2) with some $\mathsf{D}(t)\geq 0$.

Formally, $\mathcal{Z}_{t}$ constitutes the instantaneous exponent of the dynamics, satisfying $\Lambda_{t}=e^{\mathcal{Z}_{t}}$ (see Ref. Chruściński (2014) and App. A.1). Physicality of the evolution is then guaranteed by the GKLS form of $\mathcal{Z}_{t}$, because at any $t$ the dynamical map $\Lambda_{t}$ can be interpreted as a fictitious semigroup $\Lambda_{t}=\left.e^{\mathcal{Z}_{t}\tau}\right|_{\tau=1}$ generated by $\mathcal{Z}_{t}$ (at this particular time instance). $\Lambda_{t}$ must therefore be CPTP at $t$.

In general, it is not straightforward to verify whether a given QME (1) yields SS dynamics Chruściński (2014), even after decomposing its dynamical generators according to Eq. (2). However, in the special case of commutative dynamics, for which $[\mathcal{L}_{s},\mathcal{L}_{t}]=0$ (or equivalently $[\Lambda_{s},\Lambda_{t}]=0$) for all $s,t\geq 0$, one may directly identify the SS subclass, because (see App. A.1):

Lemma 1**.**

Any commutative dynamics is SS iff (if and only if) for any $t\geq 0$ the decomposition of its dynamical generators (2) fulfills

[TABLE]

In short, we term any such semigroup-simulable and commutative dynamics SSC.

Note that the condition (3) is clearly weaker than positive semi-definiteness, $\mathsf{D}(t)\geq 0$, at all times. Hence, there exist commutative dynamics which are SS but not CP-divisible. However, let us emphasise that there also exist commutative dynamics which are physical but not even SS. An explicit example is provided by the eternally non-Markovian model of Ref. Hall et al. (2014), as well as by other instances of random unitary Andersson et al. (2007); Chruściński and Wudarski (2013) and phase covariant Smirne et al. (2016) qubit dynamics, which we discuss in detail in App. B.

II.2 Additivity of dynamical generators

We define the notion of additivity for families of dynamical generators as follows:

Definition 2**.**

Two physical families of generators $\mathcal{L}^{(1)}_{t}$ and $\mathcal{L}^{(2)}_{t}$ are additive if all their non-negative linear combinations, $\alpha\mathcal{L}^{(1)}_{t}+\beta\mathcal{L}^{(2)}_{t}$ with $\alpha,\beta\geq 0$, are also physical.

Note that according to Def. 2 a pair of generator families can be additive only if each of them is individually rescalable—remains physical when multiplied by an non-negative scalar, i.e., $\mathcal{L}_{t}\to\alpha\mathcal{L}_{t}$ remains physical for any $\alpha\geq 0$. However, as such a multiplication does not invalidate the GKSL form of the decomposition (2) or the condition (3), it follows that any generator family which is CP-divisible or SSC must be rescalable.

From the linear algebra perspective Rockafellar (1970), one may formally define the vector space containing families of dynamical generators. Physical generators then form its particular subset. Rescalability of a given $\mathcal{L}_{t}$ states then that the whole ray $\{\alpha\mathcal{L}_{t}\}_{\alpha\geq 0}$ lies within the physical set. Additivity of $\mathcal{L}^{(1)}_{t}$ and $\mathcal{L}^{(2)}_{t}$, on the other hand, means that all the elements of the convex cone, $\{\alpha\mathcal{L}^{(1)}_{t}+\beta\mathcal{L}^{(2)}_{t}\}_{\alpha,\beta\geq 0}$, are physical.

For CP-divisible dynamics, we observe that when both $\mathcal{L}^{(1)}_{t}$ and $\mathcal{L}^{(2)}_{t}$ are of the GKLS form or even form a semigroup, so must any non-negative linear combination of them. Hence, it naturally follows that generator families describing CP-divisible evolutions constitute a convex cone contained in the physical set, with semigroups forming a subcone. Furthermore, as we demonstrate in App. A.3:

Lemma 2**.**

Any SSC generator families $\mathcal{L}^{(1)}_{t}$ and $\mathcal{L}^{(2)}_{t}$ are additive if upon addition, $\alpha\mathcal{L}^{(1)}_{t}+\beta\mathcal{L}^{(2)}_{t}$ with any $\alpha,\beta\geq 0$, they yield commutative dynamics.

Hence, the non-negative linear span of any such SSC pair forms a convex cone contained in the physical set. Moreover, one may then naturally expand such a cone by considering more than two, in particular, a complete set of SSC generator families whose non-negative linear combinations are all commutative. We term the convex cone so-constructed a particular SSC class.

In Fig. 1, we schematically depict the cross-section of the set of physical generator families, which then also cuts through the convex cones containing generator families of the aforementioned dynamical subclasses. Importantly, as all physical dynamics do not form a convex cone in the vector space, the ones lying within such a hyperplane are described by a non-convex set that, in turn, contains the convex sets of: CP-divisible dynamics, its semigroup subset, as well as ones representing particular SSC classes.

Now, as indicated in Fig. 1 by the dashed line, by adding a generator family that is SSC but not CP-divisible (i.e., non-Markovian Rivas et al. (2014); Breuer et al. (2016); de Vega and Alonso (2017)) and another physical family, even a semigroup, which does not commute with the first—i.e., is not contained within the corresponding SSC class—one may obtain unphysical dynamics. Consider an example of two purely dissipative ($\mathcal{L}_{t}=\mathcal{D}_{t}$ in Eq. (1)) qubit generators:

[TABLE]

where $\sigma_{\pm}=(\sigma_{x}\pm\mathrm{i}\sigma_{y})/2$ and $\sigma_{x}$, $\sigma_{y}$, $\sigma_{z}$ are the Pauli operators, and $\gamma_{1}(t)$, $\gamma_{2}(t)$ are chosen such that the generators are physical. Importantly, the families $\mathcal{L}^{(1)}_{t}$ and $\mathcal{L}^{(2)}_{t}$ despite being commutative, do not commute between one another. They belong to different SSC classes of qubit dynamics (see App. B), namely, random-unitary Andersson et al. (2007); Chruściński and Wudarski (2013) and phase-covariant Smirne et al. (2016) evolutions, respectively.

In order to prove the situation indicated in Fig. 1, we construct examples in which both $\mathcal{L}^{(1)}_{t}$ and $\mathcal{L}^{(2)}_{t}$ are physical, but their sum is not. We take instances of $\gamma_{1}(t)$ and $\gamma_{2}(t)$ with one rate being constant (semigroup), and the other taking negative values for some times (non-Markovian) while fulfilling $\int_{0}^{t}d\tau\gamma(\tau)\geq 0$ of Eq. (3) (SSC). Two simple examples are provided by choosing $\gamma_{1}(t)=\sin(\omega t)$ and $\gamma_{2}(t)=\gamma$ and vice versa, with $\gamma$ and $\omega$ being positive constants. We consider then the generator family $\mathcal{L}^{(1)}_{t}+\mathcal{L}^{(2)}_{t}$ and solve analytically in App. B.2.2 for the families of maps $\Lambda_{t}$ that arise in both cases. For each $\Lambda_{t}$, we compute the eigenvalues of its Choi-Jamiołkowski (CJ) matrix—all of which must be non-negative at all times for the map to be CPTP (see App. A.1). We depict them as a function of time in Fig. 2 for a choice of parameters which clearly demonstrates that the physicality is, indeed, invalidated at finite times. Their negativity, as demonstrated in App. B.2.2, can also be verified analytically.

In App. B.2.2, we also consider additional choices of $\gamma_{1}(t)$ and $\gamma_{2}(t)$ for the generators (4), in order to show that the same conclusion holds when both semigroup and non-Markovian contributions come from explicit microscopic derivations. In particular, as the generators describe dephasing (4a) and spontaneous-emission (4b) processes, we consider their non-Markovian forms derived (see App. B.1) from spin-boson and Jaynes-Cummings models, respectively Breuer and Petruccione (2002).

Note that it follows from the above observations that physicality of a non-Markovian QME can be easily broken by addition of even a time-invariant (semigroup) dissipative term. Moreover, one should be extremely careful when dealing with dynamics described by generator families that are not even rescalable, e.g., see App. A.2: ones that exhibit singularities at finite times Andersson et al. (2007), are derived assuming weak-coupling interactions Gaspard and Nagaoka (1999), or lead to physical (even commutative) but non-SS dynamics Hall et al. (2014).

III Microscopic approach to QMEs

Let us recall that the QME (1) constitutes an effective description of the reduced dynamics, whose form must always originate from an underlying physical mechanism responsible for both free (noiseless) and dissipative parts of the system evolution. In particular, given a microscopic model one should arrive at Eq. (1) starting from a closed dynamics describing the evolution of: the system, its environment, as well as their interaction; after tracing out the environmental degrees of freedom Gardiner and Zoller (2004); Breuer and Petruccione (2002).

III.1 Microscopic derivation of a QME

In a microscopic model of an evolving open quantum system, as illustrated in Fig. 3(a), one considers a system of interest $S$, coupled to an environment $E$ that is taken sufficiently large for the total system to be closed. The global evolution is then unitary, $U_{SE}(t)=\exp[-\mathrm{i}(H_{S}+H_{E}+H_{I})t]$, being determined by the free Hamiltonians $H_{S}$ and $H_{E}$, and the system-environment interaction $H_{I}$. In the Schrödinger picture, the reduced state of the system, ${\rho_{S}}(t)=\operatorname{Tr}_{E}{\rho_{SE}}(t)$, evolves as

[TABLE]

where ${\rho_{SE}}$ is the total system-environment state.

If the environment and the system are initially uncorrelated, so that ${\rho_{SE}}(0)={\rho_{S}}(0)\otimes{\rho_{E}}$, and ${\rho_{E}}$ is stationary, i.e., $\left[H_{E},{\rho_{E}}\right]=0$, Eq. (5) can be conveniently rewritten as (see also App. C.2) Breuer and Petruccione (2002); Rivas and Huelga (2012):

[TABLE]

where by the bar, $\bar{\bullet}:=e^{\mathrm{i}\left(H_{S}+H_{E}\right)t}\bullet\,e^{-\mathrm{i}\left(H_{S}+H_{E}\right)t}$, we denote the interaction picture with respect to the free system-environment Hamiltonian $H_{S}+H_{E}$. Eq. (6) constitutes the integro-differential QME discussed at the beginning of Sec. II that, in practice, is typically recast into the time-local form (1), which after returning to the Schrödinger picture (see App. C.3 and C.4) reads:

[TABLE]

Importantly, $\tilde{\mathcal{L}}_{t}$ above can be unambiguously identified as the dynamical generator—containing both Hamiltonian and dissipative parts as in Eq. (1)—that arises purely due to the interaction with the environment; with the system free evolution (dictated by the system Hamiltonian $H_{S}$) being explicitly separated.

III.1.1 Dependence on the system Hamiltonian

However, as detailed in App. C.5, despite the separation of terms in Eq. (7) the form of the dynamical generator, $\tilde{\mathcal{L}}_{t}$, may in general strongly depend on the system Hamiltonian $H_{S}$. Crucially, this means that the evolution of systems with different $H_{S}$, which interact with the same type of environment, cannot generally be modelled with the same QME after simply changing the $H_{S}$ in Eq. (7). However, under certain circumstances this can be justified.

In App. C.6, we discuss in detail the natural cases when variations of $H_{S}$ do not affect the form of $\tilde{\mathcal{L}}_{t}$ in Eq. (7), yet we summarise them here by the following lemma:

Lemma 3**.**

Consider a change $H_{S}\rightarrow H^{\prime}_{S}(t)=H_{S}+V(t)$ in Eq. (5). The corresponding time-local QME can be obtained by just replacing $H_{S}$ with $H^{\prime}_{S}(t)$ in Eq. (7) while keeping $\tilde{\mathcal{L}}_{t}$ unchanged, if at all times $[V(t),H_{I}]=0$ and either $[H_{S},H_{I}]=0$ or $[V(t),H_{S}]=0$ (or both).

Unfortunately, if the above sufficient condition cannot be met, one must, in principle, rederive the QME (7) and the corresponding generator $\tilde{\mathcal{L}}_{t}$ for $H^{\prime}_{S}(t)$. Moreover, such treatment is required independently of the interaction strength, i.e., also in the weak-coupling regime discussed below. A prominent physical example is provided by coherently driven systems, for which $V(t)$ represents the externally applied force. In their case, it is common that the time-dependence of $V(t)$ is naturally carried over onto, and significantly amends, the dynamical generator irrespectively of the coupling strength Rivas et al. (2010); Rivas and Huelga (2012).

III.1.2 Generalisation to multiple environments

Another important question one should pose is under what conditions the full derivation of the QME (7) can also be bypassed when dealing with a system that simultaneously interacts with multiple environments—as depicted in Fig. 3(b). Motivated by the analysis of Sec. II.2, one may then naively expect that, given multiple additive generator families describing each separate interaction, $\tilde{\mathcal{L}}^{(i)}_{t}$, they should be simply added to construct the overall QME of the form (7) with $\tilde{\mathcal{L}}_{t}=\sum_{i}\tilde{\mathcal{L}}^{(i)}_{t}$ 111Note that we are not concerned with the internal structure of the system that is crucial when discussing, e.g., additivity of decay rates for a bipartite system with each of its parts coupled to a different reservoir Yu and Eberly (2006)..

Such a procedure may, however, lead to incorrect dynamics, as may be demonstrated by considering explicitly the microscopic model that incorporates interactions with multiple environments—with now $H_{E}=\sum_{i}H_{E_{i}}$ and $H_{I}=\sum_{i}H_{I_{i}}$ in Eq. (5). Following the derivation steps of the time-local QME (7), while assuming its existence both in the presence of each single environment and all of them, one arrives at a generalised QME (see also App. D):

[TABLE]

where $\bar{H}_{I_{i}}(\tau)=e^{\mathrm{i}\left(H_{S}+H_{E_{i}}\right)\tau}H_{I_{i}}e^{-\mathrm{i}\left(H_{S}+H_{E_{i}}\right)\tau}$, $\tilde{\mathcal{L}}^{(i)}_{t}$ is the generator arising when only the $i$th environment is present, $\rho_{SE_{ij}}$ denotes the joint-reduced state of the system and environments $i$ and $j$, while $\operatorname{Tr}_{E_{ij}}$ stands for the trace over these environments.

Crucially, the naive addition of generators would lead to a QME that contains only the first two terms in Eq. (8). In particular, it would completely ignore the last term, which we here name the cross-term, as it accounts for the cross-correlations that may emerge between each two environments due to their indirect interaction being mediated by the system.

III.2 Microscopic validity of generator addition

The generalisation of the QME to multiple environments (8) allows one to unambiguously identify when the true dynamics derived microscopically coincides with the evolution obtained by naively adding the generators.

Observation 1**.**

A dynamical generator corresponding to a system simultaneously interacting with multiple environments can be constructed by simple addition of the generators associated with each individual environment iff the cross-term in Eq. (8) identically vanishes.

In what follows, we show that Obs. 1 allows one to prove the validity of generator addition in the weak-coupling regime. However, as the cross-term in Eq. (8) involves a time-convolution integral, in order to prove that it identically vanishes given the microscopic Hamiltonians satisfy particular commutation relations, we must consider the dynamics in its integrated form—at the level of the corresponding dynamical (CPTP) map.

III.2.1 Weak-coupling regime

Lemma 4**.**

The dynamical generator of the evolution of a system simultaneously interacting with multiple environments in the weak coupling regime can be constructed by simple addition of the generators associated with each individual environment.

Here, we summarise the proof of Lemma 4 that can be found in App. D.1, and generalises the argumentation of Ref. Cohen-Tannoudji et al. (1998); Schaller (2015) applicable to the more restricted regime in which the Born-Markov approximation is valid.

In particular, it applies to any QME valid in the weak-coupling regime (c.f. de Vega and Alonso (2017); Breuer and Petruccione (2002); Rivas and Huelga (2012)) which is derived using the ansatz:

[TABLE]

where ${\rho_{S}}(t)$ is the reduced system state at $t$, while $\varrho_{E_{i}}(t)$ can be arbitrarily chosen for $t>0$—it does not need to represent the reduced state of the $i$th environment, $\rho_{E_{i}}(t)=\operatorname{Tr}_{\neg E_{i}}{\rho_{SE}}(t)$, as long as it initially coincides with its stationary state, i.e., $\varrho_{E_{i}}(0)=\rho_{E_{i}}$. Let us emphasise that the assumption (9) is employed only at the derivation stage, so that the QME so-obtained, despite correctly reproducing the reduced dynamics of the system under weak coupling, may yield upon integration closed dynamics with overall system-environment states that significantly deviate from the ansatz (9) and, in particular, its tensor-product structure Rivas et al. (2010).

In App. D.1, we first employ operator Schmidt decomposition Bengtsson and Życzkowski (2006) to reexpress each of the interaction Hamiltonians in Eq. (8) as $H_{I_{i}}=\sum_{k}A_{i;k}\otimes B_{k}^{E_{i}}$, i.e, as a sum of operators that act separately on the system and corresponding environment. This decomposition, together with the ansatz (9) allows us to rewrite the overall QME (8) in terms of correlation functions involving only pairs of baths. Furthermore, the tensor-product structure of Eq. (9) ensures that each of these reduces to a product of single-bath correlation functions. Hence, as any single-bath (one-time) correlation function can always be assumed to be zero, every summand in the cross-term of the QME (8) must independently vanish.

Note that, in particular, this holds for all QMEs derived using the time-convolutionless approach Breuer et al. (2001) up to second order in all the interaction parameters representing coupling strengths for each environment.

III.2.2 Commutativity of microscopic Hamiltonians

Next, we investigate the implications that commutativity of the system, environment, and interaction Hamiltonians has on the validity of generator addition. We consider the cases when all $H_{I_{i}}$ commute with each other (II), with $H_{S}$ (IS), or with all the $H_{E_{i}}$ (IE), and summarise the results in Fig. 4. We find that:

Lemma 5**.**

Only when the interaction Hamiltonians commute among themselves and with the system free Hamiltonian, i.e., $[H_{I_{i}},H_{I_{j}}]=0$ and $[H_{I_{i}},H_{S}]=0$ for all $i,j$; can the overall QME be constructed by adding dynamical generators associated individually with each environment—ignoring the cross-term in Eq. (8).

Again, we summarise here the proof of Lemma 5 that can be found in App. D.2. However, in contrast to the discussion of the weak-coupling regime, we are required to return to the microscopic derivation of the QME (8).

Crucially, the commutativity of interaction Hamiltonians with one another as well as with $H_{S}$—the region $\text{II}\cap\text{IS}$ marked ‘Yes’ in Fig. 4—assures that the unitary of the global von Neumann equation (5) factorises, i.e.:

[TABLE]

As a result, the system dynamics is described by a product of commuting CPTP maps, $\bar{\rho}_{S}(t)=\prod_{i}\tilde{\Lambda}_{t}^{(i)}[\bar{\rho}_{S}(0)]$, associated with each individual environment and given by $\tilde{\Lambda}_{t}^{(i)}[\bullet]=\operatorname{Tr}_{E_{i}}\!\{e^{-\mathrm{i}(H_{I_{i}}+H_{E_{i}})t}\,(\bullet\otimes\rho_{E_{i}})\,e^{\mathrm{i}(H_{I_{i}}+H_{E_{i}})t}\}$. By differentiating the dynamics with respect to $t$, it is then evident that the QME takes the form (8) with each $\tilde{\mathcal{L}}^{(i)}_{t}=\dot{\tilde{\Lambda}}_{t}^{(i)}\circ(\tilde{\Lambda}_{t}^{(i)})^{-1}$ and the cross-term being, indeed, absent. Note that, as all $\tilde{\mathcal{L}}^{(i)}_{t}$ must then represent generator families belonging to a common commutative class, if each of them yields dynamics that is also SS, they all must belong to the same SSC class in Fig. 1.

In all other cases marked ‘No’ in Fig. 4, the commutativity does not ensure the generators to simply add. We demonstrate this by providing explicit counterexamples based on a concrete microscopic model, for which the evolution of the system interacting with each environment separately, as well as all simultaneously, can be explicitly solved. It is sufficient to do so for the settings in which either all $H_{I_{i}}$ commute with all $H_{E_{i}}$ and $H_{S}$ (intersection $\text{IS}\cap\text{IE}$ in Fig. 4), or all $H_{I_{i}}$ commute with each other and all $H_{E_{i}}$ (intersection $\text{II}\cap\text{IE}$), since it then follows that neither II, IS, nor IE alone can ensure the validity of generator addition.

Note that it is known that $[H_{I},H_{S}]=0$ implies the evolution to be CP-divisible Khalil and Richert (2014). Thus, as families of CP-divisible generators are additive (c.f. Fig. 1), our counterexample for $\text{IS}\cap\text{IE}$ below corresponds to a case where generator addition results in dynamics which is physical but does not agree with the microscopic derivation. For a setting in which $H_{I_{i}}$ do not commute neither among each other nor with $H_{S}$, validity of adding generators has been discussed in Ref. Chan et al. (2014).

IV Spin-magnet model

In this section, we construct counterexamples to the validity of adding generators at the QME level for the relevant cases summarised in Fig. 4. Inspired by Ref. Allahverdyan et al. (2013); Perarnau-Llobet and Nieuwenhuizen (2017), we consider a single qubit (spin-$1/2$ particle) in contact with multiple magnets—environments consisting of many spin-$1/2$ systems. Within this model, the closed dynamics of the global qubit-magnets system can be solved and, after tracing out the magnets’ degrees of freedom, the exact open dynamics of the qubit can be obtained. As a result, we can determine the QME generators describing dynamics of the qubit when coupled to one or more magnets, so that comparison with evolutions obtained by adding the corresponding generators can be explicitly made.

IV.1 Magnet as an environment

Within our model, we allow the system free Hamiltonian $H_{S}$ to be chosen arbitrarily, yet, for simplicity, we take the free Hamiltonian of the environment to vanish, $H_{E}=0$. As a result, any initial environment state is stationary, with $[\rho_{E},H_{E}]=0$ trivially for any $\rho_{E}$.

The environment is represented by a magnet that consists of $N$ spin-1/2 particles, for which we introduce the magnetisation operator:

[TABLE]

where $\Pi_{k}$ is the projector onto the subspaces with magnetisation $m_{k}$ (i.e., with $k$ spins pointing up). The magnetisation $m_{k}$ takes $N+1$ equally spaced values between $-N$ and $N$:

[TABLE]

Consistently with Sec. III, the initial state of the spin-magnet system reads ${\rho_{SE}}(0)={\rho_{S}}(0)\otimes{\rho_{E}}$, where we take the initial magnet state to be a classical mixture of different magnetisations, i.e.,

[TABLE]

The initial probability for an observation of the magnetisation to yield $m_{k}$ is then

[TABLE]

In the limit of large $N$, $p(m_{k})$ approaches a continuous distribution $p(m)$, whose moments can then be computed as follows:

[TABLE]

We consider only interaction Hamiltonians which couple the system qubit to the magnet’s magnetisation, i.e,

[TABLE]

with $A$ being an arbitrary qubit observable. However, in order to be consistent with Sec. III, we must impose that

[TABLE]

which implies that we must restrict to distributions $p(m_{k})$ (and $p(m)$ in the $N\to\infty$ limit) with zero mean.

In the examples discussed below, we consider two initial magnetisation distributions for the magnet in the asymptotic $N$ limit. In particular, we consider a Gaussian distribution:

[TABLE]

which formally corresponds to the asymptotic limit of a magnet being described by a microcanonical ensemble Allahverdyan et al. (2013)—its every spin configuration being equally probable, with $q_{k}=1/2^{N}$ in Eq. (14), yielding a binomial distribution of magnetisation with variance equal to the number of spins ($\sigma^{2}=N$). We also consider the case when the magnetisation follows a Lorentzian distribution in the $N\to\infty$ limit, i.e.:

[TABLE]

parametrised by the scale parameter $\lambda$ (specifying the half width at half maximum).

Given the above initial magnet state (13) and the interaction Hamiltonian (16), the global system-magnet state constitutes at all times a mixture of states with different magnet magnetisation. In particular, it can be decomposed at any $t\geq 0$ as

[TABLE]

where every $\rho_{S}^{(k)}$ can be understood as the (normalised) state of the system conditioned on the magnet possessing the magnetisation $m_{k}$. Consequently, the full reduced system state at time $t$ reads

[TABLE]

Crucially, within the model each of the conditional states $\rho_{S}^{(k)}$ in Eq. (21) evolves independently. In order to show this, we substitute the system-environment state (20) and the microscopic Hamiltonians into the global von Neumann equation (5) to obtain

[TABLE]

As no coupling between different magnetisation subspaces (labelled by $k$) is present, after rewriting the l.h.s. above using Eq. (20), one obtains a set of uncoupled differential equations for each conditional state:

[TABLE]

with $\rho_{S}^{(k)}(t)=\rho_{S}(0)$ for each $k$.

Hence, every $\rho_{S}^{(k)}$ evolves unitarily within our model with $U_{S}^{(k)}(t):=\exp[-\mathrm{i}(H_{S}+m_{k}A)t]$, while the overall evolution of the qubit (21) is given by the dynamical map, $\Lambda_{t}$, corresponding to a mixture of such (conditional) unitary transformations distributed according to the initial magnetisation distribution of the magnet, $p(m_{k})$, i.e.:

[TABLE]

Furthermore, as $\Lambda_{t}$ constitutes a mixture of unitaries in the model, it must be unital, i.e., for all $t\geq 0$: $\Lambda_{t}[\openone]=\openone$.

IV.1.1 Bloch ball representation

We rewrite the above qubit dynamics employing the Bloch ball representation Bengtsson and Życzkowski (2006), i.e., $\rho\equiv\frac{1}{2}(\openone+\mathbf{r}\cdot{\boldsymbol{\sigma}})$ with the Bloch vector $\mathbf{r}$ unambiguously specifying a qubit state $\rho$. Then, Eqs. (23) and (24) read, respectively:

[TABLE]

and

[TABLE]

The rotation matrices above, $\mathbf{{R}}^{(k)}$, constitute the SO(3) representations of the unitaries $U_{S}^{(k)}\in\text{SU(2)}$ in Eq. (24), and are thus similarly mixed according to $p(m_{k})$.

The qubit dynamical map, $\Lambda_{t}$ of Eq. (24), is represented by an affine transformation of the Bloch vector:

[TABLE]

which is linear due to $\Lambda_{t}$ being unital within the magnet model, i.e., does not contain a translation.

Now, as the spaces of physical $\Lambda_{t}$ (dynamical maps) and $\mathbf{{D}}_{t}$ (affine transformations) are isomorphic Bengtsson and Życzkowski (2006), the dynamical generators of the former $\mathcal{L}_{t}:=\dot{\Lambda}_{t}\circ\Lambda_{t}^{-1}$ (see App. A.1) directly translate onto $\mathbf{{L}}_{t}:=\dot{\mathbf{{D}}}_{t}\,\mathbf{{D}}_{t}^{-1}$ of the latter, with the map composition and inversion replaced by matrix multiplication and inversion, respectively. Moreover, as the vector spaces containing families of generators defined in this manner must also be isomorphic, all the notions described in Sec. II.2—in particular, rescalability and additivity—naturally carry over.

However, in order to define the Bloch ball representation of the environment-induced generator $\tilde{\mathcal{L}}_{t}$ in the QME (7), we must correctly relate it to the interaction and Schrödinger pictures of the dynamics, summarised in App. C.1. In general (see App. C.4 for the derivation from the dynamical maps perspective), the Bloch ball representation of $\tilde{\mathcal{L}}_{t}$ reads

[TABLE]

where $\mathbf{{R}}_{S}(t)\in\text{SO(3)}$ is the rotation matrix of the Bloch vector that represents the qubit unitary map, $U_{S}(t):=\exp[-\mathrm{i}H_{S}t]\in\text{SU(2)}$, induced by the system free Hamiltonian $H_{S}$.

As stated in Eq. (28), $\tilde{\mathbf{{L}}}_{t}$ may be equivalently specified with help of $\bar{\mathbf{{L}}}_{t}:=\dot{\bar{\mathbf{{D}}}}_{t}\,\bar{\mathbf{{D}}}^{-1}_{t}$, i.e, the Bloch ball representation of the dynamical generator defined in the interaction picture, $\bar{\mathcal{L}}_{t}:=\dot{\bar{\Lambda}}_{t}\circ\bar{\Lambda}_{t}^{-1}$—see also App. C.3 for its formal microscopic definition. Importantly, $\bar{\mathbf{{L}}}_{t}$ may be directly computed for a given $\mathbf{{D}}_{t}$ of Eq. (27) by first transforming it to the interaction picture, i.e., determining $\bar{\mathbf{{D}}}_{t}:=\mathbf{{R}}_{S}^{-1}(t)\,\mathbf{{D}}_{t}$ that is the Bloch-equivalent of $\bar{\Lambda}_{t}[\bullet]:=U_{S}^{\dagger}(t)\,\Lambda_{t}\!\left[\bullet\right]U_{S}(t)$ discussed in App. C.1.

On the other hand, as $\tilde{\mathbf{{L}}}_{t}$ and $\bar{\mathbf{{L}}}_{t}$ are linearly related via Eq. (28), their vector spaces must be isomorphic. Hence, in what follows, we may equivalently stick to the interaction picture and consider $\bar{\mathbf{{L}}}_{t}$ instead, in particular, $\bar{\mathbf{{L}}}^{(1)}_{t}+\bar{\mathbf{{L}}}^{(2)}_{t}\;\Leftrightarrow\;\tilde{\mathbf{{L}}}^{(1)}_{t}+\tilde{\mathbf{{L}}}^{(2)}_{t}$, when verifying the validity of generator addition.

IV.1.2 Example: magnet-induced dephasing

To illustrate the model, we first consider the case when it is employed to provide a simple microscopic derivation of the qubit dephasing dynamics. We take:

[TABLE]

with the system Hamiltonian being absent, so that all the generators, $\mathbf{{L}}_{t}=\bar{\mathbf{{L}}}_{t}=\tilde{\mathbf{{L}}}_{t}$ in Eq. (28), become equivalent. From Eq. (25) we get

[TABLE]

which just yields rotations of the Bloch ball around the $z$ axis with angular speed depending on the magnetisation $m_{k}$, i.e., Eq. (27) with (in Cartesian coordinates)

[TABLE]

for the $N\to\infty$ limit.

Integrating Eq. (31) over the initial magnetisation distribution $p(m)$, we obtain the affine transformation (27) to asymptotically read

[TABLE]

with $f(t)=\frac{1}{2}\sigma^{2}g^{2}t^{2}$ and $f(t)=\gamma gt$ in case of the Gaussian (18) and Lorentzian (19) distributions $p(m)$, respectively. Hence, the corresponding generators (28) in Bloch ball representation take a simple form:

[TABLE]

which corresponds to the standard dephasing generator with a time-dependent rate (as defined in Eq. (75) of App. B), i.e.:

[TABLE]

with $\gamma(t)=\sigma^{2}g^{2}t$ in the Gaussian case, and constant (semigroup) $\gamma(t)=\lambda g$ in the Lorentzian case.

IV.2 Counterexamples to sufficiency of the commutativity assumptions

We now prove the regions marked ‘No’ in Fig. 4. In particular, we provide explicit counterexamples which assure that the commutativity assumption—associated with the particular region of the Venn diagram—is generally not sufficient for the system dynamics to be recoverable by simple addition of the generators attributed to each of the environments. In order to do so, it is enough to consider the scenario in which the qubit is independently coupled to just two magnets via the mechanism described above.

IV.2.1 $\text{IS}\cap\text{IE}$ commutativity assumption

We start with an example of dynamics, in which the interaction Hamiltonians (trivially) commute with the free system and all the environmental Hamiltonians, but not among each other. In particular, we simply set:

[TABLE]

with subscripts $\{1,2\}$ labelling to the first and the second magnet. As in the case of the dephasing-noise derivation above, the system Hamiltonian is absent, so the generators in Eq. (28) coincide with $\mathbf{{L}}_{t}=\bar{\mathbf{{L}}}_{t}=\tilde{\mathbf{{L}}}_{t}$.

In the case of simultaneous coupling to two magnets, Eq. (20) naturally generalises to

[TABLE]

where $q_{k,k^{\prime}}$ now represents the joint probability of finding the first and the second magnet in magnetisations $m_{k}$ and $m_{k^{\prime}}$, respectively, while $\rho_{S}^{(k,k^{\prime})}(t)$ stands for the corresponding conditional reduced state of the system.

Consequently, the (conditional) von Neumann equation (25), which now must be derived for $H_{I}=H_{I_{1}}+H_{I_{2}}$, describes the dynamics of Bloch-vectors that represent each conditional state, $\rho_{S}^{(k,k^{\prime})}(t)$, being also parametrised by the two indices $k$ and $k^{\prime}$, i.e.:

[TABLE]

Eq. (37) leads to coupled equations in the Cartesian basis, i.e. (dropping the indices $k,k^{\prime}$ and the explicit time-dependence for simplicity):

[TABLE]

which can be analytically solved to obtain the $\mathbf{{R}}$-matrix in Eq. (27)—labelled $\mathbf{{R}}_{12}$ to indicate that both magnets are involved. $\mathbf{{R}}_{12}(m_{1},m_{2},t)$ possesses now two magnetisation parameters associated with each of the magnets, and we state its explicit form in App. E.1.

Furthermore, we can straightforwardly obtain the solution of the equations of motion when only one of the magnets is present by simply setting either $g_{1}=0$ or $g_{2}=0$ in $\mathbf{{R}}_{12}$. In presence of only the first magnet ($g_{2}=0$), we recover the magnet-induced dephasing noise described above—with $\mathbf{{R}}_{12}(m_{1},m_{2},t)$ simplifying to $\mathbf{{R}}_{1}(m_{1},t)$ that takes exactly the form (31). On the other hand, when only the second magnet ($g_{1}=0$) is present, which couples to the system via $\sigma_{x}$ rather than $\sigma_{z}$, see Eq. (35), we obtain $\mathbf{{R}}_{2}(m_{2},t)$ as in Eq. (31) but with coordinates cyclically exchanged—see App. E.1 for explicit expressions.

We then average each $\mathbf{{R}}_{\mathsf{x}}$, where $\mathsf{x}=\{12,1,2\}$ denotes the magnet(s) being present, over the initial magnetisations in the $N\!\to\!\infty$ limit. This way, we obtain the affine maps (27) representing the corresponding qubit dynamics for all the three cases in the asymptotic $N$ limit as:

[TABLE]

where in presence of both magnets $m_{12}\equiv(m_{1},m_{2})$ and $p(m_{12})\equiv p(m_{1})\,p(m_{2})$; and we take each $p(m_{i})$ to follow a Gaussian distribution (18) with variance $\sigma_{i}$. Similarly, we obtain the integral expressions for the time-derivatives of the affine maps, $\dot{\mathbf{{D}}}_{t}^{(\mathsf{x})}=\int\!dm_{\mathsf{x}}\,p(m_{\mathsf{x}})\,\dot{\mathbf{{R}}}_{\mathsf{x}}(m_{\mathsf{x}},t)$ after also computing analytically all the corresponding $\dot{\mathbf{{R}}}_{\mathsf{x}}$.

Finally, we choose particular values of $g_{1}$, $g_{2}$, $\sigma_{1}$, $\sigma_{2}$ and time $t$, in order to numerically compute the integrals over the magnetisation parameters and obtain all $\mathbf{{D}}_{t}^{(\mathsf{x})}$ and $\dot{\mathbf{{D}}}_{t}^{(\mathsf{x})}$. Our choice, allows us then to explicitly construct dynamical generators $\mathbf{{L}}^{(\mathsf{x})}_{t}=\dot{\mathbf{{D}}}_{t}^{(\mathsf{x})}(\mathbf{{D}}_{t}^{(\mathsf{x})})^{-1}$, which importantly exhibit $\mathbf{{L}}^{(12)}_{t}\neq\mathbf{{L}}^{(1)}_{t}+\mathbf{{L}}^{(2)}_{t}$, see App. E.1.

Hence, we conclude that the commutativity of the interaction Hamiltonians with both system and environment Hamiltonians, but not with each other, cannot assure the generators to simply add at the level of the QME—as denoted by the ‘No’ label in the region of Fig. 4 representing the $\text{IS}\cap\text{IE}$ commutativity assumption.

IV.2.2 $\text{II}\cap\text{IE}$ commutativity assumption

In order to construct an example of reduced dynamics in which, at the microscopic level, the interaction Hamiltonians commute with each other and all free environmental Hamiltonians but not the system Hamiltonian, we consider again a two-magnet model but this time set:

[TABLE]

In contrast to the previous example, since $H_{S}\neq 0$, in order to investigate the validity of generator addition we must either consider the environment-induced generators $\tilde{\mathbf{{L}}}_{t}$ defined in Eq. (28) or the generators $\bar{\mathbf{{L}}}_{t}$ directly computed in the interaction picture.

We do the latter and compute both $H_{I_{i}}$ in the interaction picture, i.e., $\bar{H}_{I_{i}}(t)$, which are obtained by replacing $\sigma_{z}$ Pauli operators in Eq. (40) with

[TABLE]

Importantly, ${\bar{\sigma}}_{z}(t)$ should be interpreted as a (time-dependent) operator $A$ in the general expression (16) for $H_{I}$ that, in contrast to the previous case, is now identical for both magnets. Thus, inspecting the general expression for the dynamics (25), we obtain the equation of motion for the Bloch vector in the interaction picture, $\bar{\mathbf{r}}^{(k,k^{\prime})}(t):=\mathbf{{R}}_{S}^{-1}(t)\,\mathbf{r}^{(k,k^{\prime})}(t)$, that represents the qubit state conditioned on the first and second magnet possessing magnetisations $m_{k}$ and $m_{k^{\prime}}$, respectively, as

[TABLE]

which leads to coupled equations (again, dropping the indices $k,k^{\prime}$ and the explicit time-dependence):

[TABLE]

As before, see App. E.2, we solve the above equations of motion in order to obtain the $\bar{\mathbf{{R}}}$-matrix of Eq. (27) in the interaction picture, i.e., $\bar{\mathbf{{R}}}_{12}(m_{1},m_{2},t)$. Again, by setting either $g_{2}=0$ or $g_{1}=0$, we obtain expressions for $\bar{\mathbf{{R}}}_{1}$ and $\bar{\mathbf{{R}}}_{2}$, respectively, corresponding to the cases when only first or second magnet is present. We then also compute all $\dot{\bar{\mathbf{{R}}}}_{\mathsf{x}}$ with $\mathsf{x}=\{12,1,2\}$, in order to arrive at integral expressions for both the affine maps and their time-derivatives, i.e., $\bar{\mathbf{{D}}}_{t}^{(\mathsf{x})}$ and $\dot{\bar{\mathbf{{D}}}}_{t}^{(\mathsf{x})}$, respectively, computed now in the interaction picture.

As in the previous example, we take initial magnetisation distributions of both magnets to be Gaussian and fix all the model parameters (i.e., $\omega$ and $g_{1}$, $g_{2}$ for the system and interaction, $\sigma_{1}$, $\sigma_{2}$ for the magnets, as well as the time $t$) in order to numerically perform the integration over magnetisations $m_{\mathsf{x}}$. We then find, see App. E.2, a choice of parameters for which it is clear that the dynamical generators, $\bar{\mathbf{{L}}}^{(\mathsf{x})}_{t}=\dot{\bar{\mathbf{{D}}}}_{t}^{(\mathsf{x})}(\bar{\mathbf{{D}}}_{t}^{(\mathsf{x})})^{-1}$, fulfil $\bar{\mathbf{{L}}}^{(12)}_{t}\neq\bar{\mathbf{{L}}}^{(1)}_{t}+\bar{\mathbf{{L}}}^{(2)}_{t}$.

Hence, we similarly conclude that the commutativity of all the interaction Hamiltonians with each other, and all the free Hamiltonians of environments also cannot assure the generators to simply add at the level of the QME—proving the ‘No’ label in the region of Fig. 4 representing the $\text{II}\cap\text{IE}$ commutativity assumption.

V Conclusions

We have investigated under what circumstances modifications to open system dynamics can be effectively dealt with at the master equation level by adding dynamical generators. We have identified a condition—semigroup simulability and commutativity preservation—applicable beyond Markovian (CP-divisible) dynamics which guarantees generator addition to yield physical evolutions. We have also demonstrated by considering simple qubit generators that even mild violation of this condition may yield unphysical dynamics under generator addition.

Moreover, even when physically valid, generator addition does not generally correspond to the real evolution derived from a microscopic model describing interactions with multiple environments. We have formulated a general criterion under which the addition of generators associated with each individual environment yields the correct dynamics.

We have then shown that this condition is generally satisfied in the weak-coupling regime, whenever it is correct to use a master equation derived assuming a tensor-product ansatz for the global state describing the system and environments. Finally, we have demonstrated that, at the microscopic level, the commutativity of interaction Hamiltonians among each other and with the system Hamiltonian also ensures addition of dynamical generators to give the correct dynamics.

We believe that our results may prove useful in areas where the master equation description of open quantum systems is a common workhorse, including quantum metrology, thermodynamics, transport, and engineered dissipation.

Acknowledgements.

We would like to thank L. Aolita and N. Bernades for interesting discussions on non-Markovianity that sparked this work, as well as A. Smirne, M. Lostaglio, S. Huelga, L. Correa and A. S. Sørensen for helpful exchanges. J.K. and B.B. acknowledge support from the Spanish MINECO (Grant QIBEQI FIS2016-80773-P and Severo Ochoa SEV-2015-0522), Fundació Privada Cellex, Generalitat de Catalunya (SGR875 and CERCA Program). J.K. is also supported by the EU Horizon 2020 programme under the MSCA Fellowship Q-METAPP (no. 655161), while B.B. by the ICFO-MPQ Fellowship. M.P.-L. acknowledges also support from the Alexander von Humboldt Foundation.

Appendix A QMEs as families of dynamical generators

A.1 Describing open system dynamics

A.1.1 Physically valid quantum dynamics

A particular evolution of an open quantum system is formally represented by a continuous family of density matrices $\left\{\rho_{t}\in\mathcal{B}(\mathcal{H}_{d})\right\}_{t\geq 0}$ that describe the system state at each time $t\geq 0$. The system evolution is then defined by a family of dynamical maps (quantum channels Bengtsson and Życzkowski (2006)), $\left\{\Lambda_{t}\right\}_{t\geq 0}$ with $\Lambda_{0}=\mathcal{I}$ being the identity map, such that for any initial system state, $\rho_{0}$, the state at time $t\geq 0$ is given by

[TABLE]

Importantly, for a given dynamics to be physical the family $\left\{\Lambda_{t}\right\}_{t\geq 0}$ must consist of completely-positive trace preserving (CPTP) maps. Only then, for any given enlarged initial state $\varrho(0)\in\mathcal{B}(\mathcal{H}_{d}\otimes\mathcal{H}_{d}^{\prime})$ with arbitrary $d^{\prime}=\dim\mathcal{H}_{d}^{\prime}$, the state at every time $t\geq 0$, i.e., $\varrho(t)=\Lambda_{t}\otimes\mathcal{I}[\varrho(0)]$, is guaranteed to be correctly described by a positive semidefinite matrix.

In practice, any linear map $\Lambda:\,\mathcal{B}(\mathcal{H}_{d})\to\mathcal{B}(\mathcal{H}_{d^{\prime}})$ may be verified to be CPTP by constructing its corresponding Choi-Jamiołkowski (CJ) matrix $\Omega_{\Lambda}\in\mathcal{B}(\mathcal{H}_{d^{\prime}}\otimes\mathcal{H}_{d})$ defined as (Jamiołkowski, 1972; Choi, 1975):

[TABLE]

with $\left|\psi\right\rangle=\sum_{i=1}^{d}\left|i\right\rangle\!\left|i\right\rangle$ and $\left\{\left|i\right\rangle\right\}_{i=1}^{d}$ being some orthonormal basis spanning $\mathcal{H}_{d}$. In particular, a map $\Lambda$ is CP and TP iff its CJ matrix is positive semi-definite, i.e., $\Omega_{\Lambda}\geq 0$, and satisfies $\mathrm{Tr}_{\mathcal{H}_{d^{\prime}}}\!\left\{\Omega_{\Lambda}\right\}=\openone_{d}$, respectively.

A.1.2 Dynamical generators

One may associate with any given dynamics (44) the family of dynamical generators, $\left\{\mathcal{L}_{t}\right\}_{t\geq 0}$, that specify for each $t\geq 0$ the time-local QME stated in Eq. (1) of the main text (Andersson et al., 2007; Chruściński, 2014; Hall et al., 2014), i.e.:

[TABLE]

with $\dot{\bullet}\equiv\frac{{d}}{{d}t}\bullet$ and the dynamical generator being then formally defined at each $t\geq 0$ as

[TABLE]

where $\Lambda_{t}^{-1}$ is the inverse ($\Lambda_{t}^{-1}\circ\Lambda_{t}=\mathcal{I}$) of the dynamical map at time $t$, and is not necessarily CPTP. However, $\Lambda_{t}^{-1}$ in general may cease to exist at certain time-instances, at which the corresponding generators $\mathcal{L}_{t}$ then become singular, even though the family of maps is perfectly smooth. Nevertheless, under particular conditions (Andersson et al., 2007), the resulting QME (46) can still be integrated and yield correctly the original dynamics (44).

In the other direction, given a family of dynamical generators $\left\{\mathcal{L}_{t}\right\}_{t\geq 0}$, one may write the corresponding dynamical map at any time $t$ with help of a time-ordered exponential, expressible in the Dyson-series form (Andersson et al., 2007; Chruściński, 2014):

[TABLE]

where $\mathcal{S}_{t}^{(0)}\left[\bullet\right]=\mathcal{I}$ and for all $i\geq 1$:

[TABLE]

A.1.3 Instantaneous generators

Given a family of dynamical maps $\left\{\Lambda_{t}\right\}_{t\geq 0}$ defining the evolution (44), one can also construct the corresponding family of instantaneous generators $\left\{\mathcal{X}_{t}\right\}_{t\geq 0}$, defined for each $t\geq 0$ as (Chruściński, 2014):

[TABLE]

being the so-called instantaneous exponent, such that:

[TABLE]

where the above expression, in contrast to Eq. (48), does not involve the time-ordering operator $\mathcal{T}_{\leftarrow}$.

Although the instantaneous generators $\mathcal{X}_{t}$ cannot be directly used to construct the QME (46), the family of dynamical generators $\{\mathcal{L}_{t}\}_{t\geq 0}$ can be formally related to the family of instantaneous ones. In particular, by substituting into Eq. (47)

[TABLE]

and $\Lambda_{t}^{-1}=e^{-\mathcal{Z}_{t}}$, one observes that

[TABLE]

A.1.4 Commutative dynamics

A dynamics is defined to be commutative if the maps describing the evolution in Eq. (44), or equivalently—as follows from Eq. (48)—all the dynamical generators defining the QME (46) commute between one another, i.e., for all $s,t\geq 0$:

[TABLE]

Moreover, Eq. (51) implies then that also all instantaneous exponents must commute with each another, $\left[\mathcal{Z}_{s},\mathcal{Z}_{t}\right]=0$, but also with instantaneous generators, $0=\partial_{s}\left[\mathcal{Z}_{s},\mathcal{Z}_{t}\right]=\left[\mathcal{X}_{s},\mathcal{Z}_{t}\right]$. Hence, in case of commutative dynamics the dynamical and instantaneous generators must coincide at all times, as

[TABLE]

so that Eqs. (48) and (51) both, respectively, read

[TABLE]

A.1.5 CP-divisible dynamics

The dynamics is said to be divisible into CPTP maps—CP-divisible (or Markovian (Rivas et al., 2014; Breuer et al., 2016; de Vega and Alonso, 2017))—if its corresponding family of maps, $\left\{\Lambda_{t}\right\}_{t\geq 0}$, satisfies for all $0\leq s\leq t$:

[TABLE]

where $\tilde{\Lambda}_{t,s}$ is a CPTP map itself.

At the level of dynamical generators, this is equivalent to the statement that all $\left\{\mathcal{L}_{t}\right\}_{t\geq 0}$ are of the Gorini-Kosakowski-Sudarshan-Linblad (GKSL) form (Gorini et al., 1976; Lindblad, 1976):

[TABLE]

where $H_{t}$ is a time-dependent Hermitian operator, $\Phi_{t}$ is a completely-positive (CP) map, and they, respectively, represent the Hamiltonian, $\mathcal{H}_{t}$, and dissipative, $\mathcal{D}_{t}$, parts in the QME (1) of the main text. $\Phi_{t}^{\star}$ is the dual map of $\Phi_{t}$ that—given a Kraus representation of the CP map $\Phi_{t}$, i.e., a set of operators $\left\{V_{i}(t)\right\}_{i}$ satisfying $\sum_{i}V_{i}(t)^{\dagger}V_{i}(t)=0$ for all $t\geq 0$ such that $\Phi_{t}\!\left[\bullet\right]=\sum_{i}V_{i}(t)\bullet V_{i}(t)^{\dagger}$ (Sudarshan et al., 1961)—is defined as $\Phi_{t}^{\star}\!\left[\bullet\right]:=\sum_{i}V_{i}^{\dagger}(t)\bullet V_{i}(t)$.

Thus, one may rewrite Eq. (58) also as (Breuer and Petruccione, 2002):

[TABLE]

which—after fixing a particular orthonormal basis of matrices $\left\{F_{j}\right\}_{j}$ satisfying $\operatorname{Tr}\!\left\{F_{i}^{\dagger}F_{j}\right\}=\delta_{ij}$ in which each $V_{i}(t)=\sum_{j}\mathsf{V}_{ij}(t)F_{j}$—can be further rewritten as in Eq. (2) of the main text:

[TABLE]

with the time dependence of the dissipative part being now fully contained within the matrix $\mathsf{D}(t)$.

Although any dynamical generator $\mathcal{L}_{t}$, constituting a traceless and Hermiticity-preserving operator, can be decomposed as above, the GKSL form (58) ensures that for all $t\geq 0$ there exists a matrix $\mathsf{V}(t)$ such that $\mathsf{D}(t)=\mathsf{V}(t)^{\dagger}\mathsf{V}(t)$. Hence, it follows that any dynamics is CP-divisible iff one may at all times decompose its corresponding dynamical generators according to Eq. (60) with some positive semi-definite $\mathsf{D}(t)\geq 0$.

A.1.6 Semigroup dynamics

An important subclass of commutative and CP-divisible dynamics are semigroups, for which the whole evolution is determined by a single fixed generator $\mathcal{L}$,

[TABLE]

which in order to describe physical dynamics (so that all $\Lambda_{t}$ are CPTP) must be of the GKSL form (58) with both the Hamiltonian $H$ and the positive semi-definite matrix $\mathsf{D}\geq 0$ in Eq. (60) being now time-independent.

A.1.7 Semigroup-simulable (SS) dynamics

Definition 3**.**

We define a map $\Lambda_{t}$ to be (instantaneously) semigroup-simulable (SS) at time $t$ if its corresponding instantaneous exponent $\mathcal{Z}_{t}=\log\Lambda_{t}$ in Eq. (50) is of the GKSL form, i.e.:

[TABLE]

where, similarly to Eq. (58), $\tilde{H}_{t}$ and $\tilde{\Phi}_{t}$ are some Hermitian operator and CP map, respectively. If all instantaneous exponents, $\{\mathcal{Z}_{t}\}_{t\geq 0}$, can be decomposed according to Eq. (62), we term the whole dynamics to be SS.

Importantly, given that $\Lambda_{t}$ is SS at time $t$, a semigroup parametrised by $\tau\geq 0$ with physical (of GKSL form) generator $\mathcal{L}=\mathcal{Z}_{t}$ may be defined:

[TABLE]

so that it coincides with the original map at $\tau=1$, $\tilde{\Lambda}(t)_{\tau=1}=\Lambda_{t}$, or, in other words, “simulates” its action at this particular instance of “fictitious time” $\tau$.

Observation 2**.**

The SS property provides a sufficient but not necessary condition for physicality of dynamics.

If for a dynamical family $\left\{\Lambda_{t}=e^{\mathcal{Z}_{t}}\right\}_{t\geq 0}$ all its instantaneous exponents are of the GKSL form (62), it must consist of maps which coincide with semigroups at all $t\geq 0$ and, hence, all must be CPTP. In the other direction, however, there exist dynamics that are not SS but nonetheless physical. Examples may be found by considering instances of, e.g.,random unitary and phase covariant, qubit evolutions, as shown below in App. B.1.

In case of commutative dynamics, it follows from Eq. (56) that $\mathcal{Z}_{t}=\int_{0}^{t}{d}\tau\mathcal{L}_{\tau}$, so that one may explicitly connect the decomposition (60) of the dynamical generator at time $t$ with the one of the instantaneous exponent in Eq. (62), as follows

[TABLE]

where $\int_{0}^{t}{d}\tau H_{\tau}$ constitutes then $\tilde{H}_{t}$ in Eq. (62).

Crucially, the decomposition (64) proves Lemma 1 of the main text, as it becomes clear that the GKSL form of $\mathcal{Z}_{t}$ is then fully ensured by the condition

[TABLE]

stated in Eq. (3) of the main text. Hence, given a commutative dynamics for which condition (65) holds, it must also be SS—constitute a semigroup-simulable and commutative (SSC) evolution.

In some previous works (Wolf et al., 2008; Cubitt et al., 2012), the SS property has been identified as Markovianity of the dynamics. Let us emphasise that such a notion is non-trivially related to the concept of CP-divisibility introduced in App. A.1.5, which is more commonly associated with Markovianity (Rivas et al., 2014; Breuer et al., 2016; de Vega and Alonso, 2017). The CP-divisibility ensures $\mathsf{D}(t)\geq 0$ in Eq. (60) at all times, so that (in case of commutative dynamics) the SSC condition (65) is trivially fulfilled. However, as $\mathsf{D}(t)\geq 0$ is a stronger requirement, there must exist (also commutative) evolutions that are SS but not CP-divisible, e.g., instances of qubit dynamics discussed below in App. B. This fact can also be understood by inspecting Eq. (53), from which it is clear that the GKSL form (62) of the instantaneous exponent $\mathcal{Z}_{t}$ (and, hence, of $\mathcal{X}_{t}=\dot{\mathcal{Z}}_{t}$) does not generally ensure the corresponding dynamical generator $\mathcal{L}_{t}$ to also be of GKSL form (58).

A.2 Rescalability of dynamical generators

Definition 4**.**

We define a physical family of dynamical generators $\left\{\mathcal{L}_{t}\right\}_{t\geq 0}$ to be rescalable if by multiplying all its elements by any non-negative constant, $\alpha\geq 0$, one obtains a generator family,

[TABLE]

that also yields physical dynamics.

Given a family of dynamical maps $\left\{\Lambda_{t}\right\}_{t\geq 0}$, by rescaling its corresponding dynamic generators $\left\{\mathcal{L}_{t}\right\}_{t\geq 0}$, as in Eq. (66), we obtain a family of maps, $\left\{\Lambda_{t}^{\prime}\right\}_{t\geq 0}$, that according to Eq. (48) reads

[TABLE]

Crucially, as the above Dyson series includes now the factor $\alpha\geq 0$, it is non-trivial to determine whether the resulting map $\Lambda_{t}^{\prime}$ is CPTP; even in case of commutative dynamics for which time-ordering, $\mathcal{T}_{\leftarrow}$, can be dropped.

The rescalability, however, is naturally ensured in case of CP-divisible evolutions (and, hence, semigroups), as the GKSL form (58) of any dynamical generator is then trivially carried over to $\mathcal{L}^{\prime}_{t}$ in Eq. (66) for any $\alpha\geq 0$.

On the other hand, any family of dynamical generators yielding SSC dynamics must also be rescalable. As the instantaneous exponents are then related to the dynamical generators via $\mathcal{Z}_{t}=\int_{0}^{t}{d}\tau\mathcal{L}_{\tau}$, they transform similarly to Eq. (66) with $\mathcal{Z}^{\prime}_{t}:=\alpha\mathcal{Z}_{t}$. Thus, the condition (65) ensuring their GKSL form (62) is fulfilled for any $\alpha\geq 0$.

Nevertheless, non-rescalability of generators also naturally emerges in some particular situations, e.g., when dealing with:

Dynamical generators with singularities,

which emerge in case of evolutions whose family of dynamical maps, $\{\Lambda_{t}\}_{t\geq 0}$, contains non-invertible CPTP maps. In this case, the dynamics can be unambiguously recovered from the dynamical generators only for times smaller than $T$, denoting the occurrence of the (first) singularity Andersson et al. (2007). As a result, even though the dynamics is physical despite $\{\mathcal{L}_{t}\}_{t\geq 0}$ containing singular generators, as soon as $\alpha\neq 1$ in Eq. (66) the integrability of the corresponding QME (46)—and, hence, the physicality—is lost for times $t\geq T$. We provide an explicit example of such a phenomenon below in App. B.2.1, where we discuss the Jaynes-Cummings model describing a qubit that undergoes spontaneous emission Breuer and Petruccione (2002).

Weak-coupling-based generators,

which are approximate and only valid for a particular timescale $T$ ($0\leq t\leq T$). Consider a family $\left\{\lambda^{2}\mathcal{L}_{t}\right\}_{0\leq t\leq T}$ of generators derived by employing a microscopic model and assuming the system-environment coupling constant, $\lambda$, small enough, so that the weak-coupling approximation to $O(\lambda^{2})$ holds and $\mathcal{L}^{real}_{t}\approx\lambda^{2}\mathcal{L}_{t}$ (Gaspard and Nagaoka, 1999) (e.g., by assuming the Redfield form of the QME (Redfield, 1957)). One may then simply interpret the rescaling factor as the square of the coupling constant, $\alpha=\lambda^{2}$. Importantly, such a generator family is guaranteed to yield physical dynamics—a family of CPTP maps—only on timescales with $T\ll\lambda^{-2}$ (Dümcke and Spohn, 1979). Hence, by rescaling the generators with large enough $\alpha=\lambda^{2}$ or, in other words, by choosing strong enough coupling, one must at some point invalidate the weak-coupling approximation and, eventually, the physicality.

Commutative but not SS dynamics.

Although all families of dynamical generators that lead to SSC dynamics must be rescalable, the commutativity property alone is not enough. A direct example is provided by the eternally non-Markovian model introduced in Ref. Hall et al. (2014) and discussed below in App. B.1. In particular, when rescaling its generators according to Eq. (66), one obtains dynamics that is not physical for any $0\leq\alpha<1$ (Benatti et al., 2017).

A.3 Additivity of dynamical generators

The Def. 2 of the main text may be restated in a more detailed form as:

Definition 5**.**

Two families of physical and rescalable dynamical generators $\left\{\mathcal{L}_{t}^{(1)}\right\}_{t\geq 0}$ and $\left\{\mathcal{L}_{t}^{(2)}\right\}_{t\geq 0}$ are additive, if all their non-negative linear combinations,

[TABLE]

with $\alpha,\beta\geq 0$, yield families of dynamical generators, $\left\{\mathcal{L}_{t}^{\prime}\right\}_{t\geq 0}$, that are physical.

Firstly, we realise that (as for rescalability) all pairs of generator families describing CP-divisible evolutions must be additive, as by adding families of CP-divisible dynamics according to Eq. (68) one obtains generators that are also of the GKSL form (58).

On the other hand, by considering generator families describing SSC dynamics, we observe that:

Lemma 6**.**

Any pair of SSC dynamics with generator families $\left\{\mathcal{L}_{t}^{(1)}\right\}_{t\geq 0}$ and $\left\{\mathcal{L}_{t}^{(2)}\right\}_{t\geq 0}$ whose addition (68) yields commutative dynamics $\left\{\mathcal{L}_{t}^{\prime}\right\}_{t\geq 0}$ for any $\alpha,\beta\geq 0$ must be additive.

Proof.

As all the families $\left\{\mathcal{L}_{t}^{(1)}\right\}_{t\geq 0}$, $\left\{\mathcal{L}_{t}^{(2)}\right\}_{t\geq 0}$ and $\left\{\mathcal{L}_{t}^{\prime}\right\}_{t\geq 0}$ are commutative, their instantaneous exponents also add according to Eq. (68), i.e., $\mathcal{Z}^{\prime}_{t}=\alpha\mathcal{Z}_{t}^{(1)}+\beta\mathcal{Z}_{t}^{(2)}$. Moreover, as $\left\{\mathcal{L}_{t}^{(1)}\right\}_{t\geq 0}$, $\left\{\mathcal{L}_{t}^{(2)}\right\}_{t\geq 0}$ are SS, both $\mathcal{Z}_{t}^{(1)}$ and $\mathcal{Z}_{t}^{(2)}$ must satisfy Eq. (65) with $\mathsf{\Gamma}^{(1)}(t)\geq 0$ and $\mathsf{\Gamma}{}^{(2)}(t)\geq 0$. Hence, any family $\left\{\mathcal{L}^{\prime}_{t}\right\}_{t\geq 0}$ must also be SS (and, hence, physical), as recomputing condition (65) for $\mathcal{Z}^{\prime}_{t}$ with help of Eq. (64) it reads

[TABLE]

and is trivially fulfilled for any $\alpha,\beta\geq 0$. ∎

Appendix B Rescalability and additivity of qubit dynamical generators

B.1 RU and PC classes of qubit dynamics

We consider two important classes of commutative qubit dynamics, namely, random unitary (RU) Andersson et al. (2007); Chruściński and Wudarski (2013) and phase-covariant (PC) Smirne et al. (2016) evolutions. In order to provide their physically motivated instances, we explicitly discuss exemplary microscopic derivations for the (generalised) dephasing and amplitude damping models that fall into the RU and PC classes, respectively.

B.1.1 Random unitary (RU) dynamics

Random unitary (RU) qubit dynamics are formed by considering smooth families of Pauli channels, which up to unitary transformations represent the most general qubit unital ($\Lambda[\openone]=\openone$) maps Vacchini (2012). In particular, any RU evolution is described by a qubit QME (46) with (Andersson et al., 2007; Chruściński and Wudarski, 2013):

[TABLE]

so that the RU generator family $\left\{\mathcal{L}_{t}\right\}_{t\geq 0}$ is fully specified by the three rates $\gamma_{k}(t)$ defining a diagonal form of the general $\mathsf{D}$-matrix in Eq. (60), with Pauli operators $\left\{\openone,\sigma_{x},\sigma_{y},\sigma_{z}\right\}$ constituting a basis for two-dimensional Hermitian matrices. Hence, it directly follows that the GKSL form (58) of RU generators, and hence the CP-divisibility (57) of the dynamics, is ensured iff at all times all $\gamma_{k}(t)\geq 0$ are non-negative in Eq. (70).

One may straightforwardly verify that any RU dynamics (70) is commutative, so that by Eq. (55) dynamical and instantenuous generators coincide, and the instantaneous exponents according to Eq. (64) read:

[TABLE]

with $\Gamma_{k}(t):=\int_{0}^{t}{d}\tau\gamma_{k}(\tau)$. Hence, it directly follows from the condition (65) that any RU dynamics is SS iff

[TABLE]

Note that, as one may easily construct families of dynamical generators (70) that satisfy Eq. (72) without requiring $\gamma_{k}(t)\geq 0$ for all $t$, there exist RU dynamics that are SS but not CP-divisible.

On the other hand, an explicit condition for the physicality of RU dynamics is known (Andersson et al., 2007; Chruściński and Wudarski, 2015). In particular, RU dynamics is physical iff for all the cyclic permutations of $i,j,k\in\{x,y,z\}$ (i.e., such that $\epsilon_{ijk}=1$):

[TABLE]

where each $\mu_{i}(t):=\exp\!\left[-2\left(\Gamma_{j}(t)+\Gamma_{k}(t)\right)\right]$. It is easy to verify that the physicality condition (73) is less restrictive than the SS condition (72). Hence, there exist RU dynamics that are physical but not SS, despite being commutative.

Eternally non-Markovian model.

An example of RU dynamics that is physical but not SS is also provided by the eternally non-Markovian model introduced in Ref. Hall et al. (2014), which corresponds to the following choice of rates in Eq. (70):

[TABLE]

for which the physicality condition (73) holds, even though $\gamma_{z}(t)<0$ (and hence $\Gamma_{z}(t)<0$) for all $t\geq 0$.

Dephasing dynamics.

The simplest example of RU dynamics (70) is provided by the dephasing model:

[TABLE]

where $\sigma_{\boldsymbol{n}}=\boldsymbol{n}\cdot{\boldsymbol{\sigma}}=\sum_{i}n_{i}\sigma_{i}$, and $\sigma_{\boldsymbol{n}}^{2}=\openone$ implies $\left\|\boldsymbol{n}\right\|=1$. The unit vector $\boldsymbol{n}$ should be interpreted as a choice (a passive rotation in the Bloch-ball picture) of the Pauli-operator basis, in which then Eq. (75) corresponds to (rank-one Pauli) RU dynamics (70) with only a single term present in the sum. One may easily verify that for the dephasing model to be physical $\Gamma(t)=\int_{0}^{t}{d}\tau\gamma(\tau)\geq 0$, with the notions of physicality and SS then trivially coinciding.

The dephasing dynamics (75) can be explicitly obtained by considering various microscopic derivations, in which a qubit is coupled to a large environment via some $H_{\mathrm{int}}\propto\sigma_{\boldsymbol{n}}\otimes O_{\mathrm{env}}$. In Sec. IV.1.2 of the main text, we provide a compact example by using a toy-model of a qubit coupled to a large magnet. The most common microscopic derivation, however, is constructed by considering a qubit coupled to a large, thermal bosonic bath Breuer and Petruccione (2002). The interaction is then modelled by $H_{\mathrm{int}}\propto\sigma_{\boldsymbol{n}}\otimes\left(\sum_{k}g_{k}\hat{a}_{k}+g_{k}^{\star}\hat{a}_{k}^{\dagger}\right)$ which couples the qubit to a bosonic reservoir of an Ohmic-like spectral density (Leggett et al., 1987):

[TABLE]

where $\omega_{c}$ represents the reservoir cutoff frequency, while $s\geq 0$ is the so-called Ohmicity parameter.

Assuming further the reservoir to be at zero temperature, the dynamical generators describing the qubit evolution take then exactly the form (75) with the dephasing rate reading (Haikka et al., 2013):

[TABLE]

where $\Gamma[s]$ above represents the Euler gamma function. Moreover, one may show that the dephasing rate temporarily takes negative values iff $s>2$, so that the dynamics ceases then to be CP-divisible (Haikka et al., 2013).

B.1.2 Phase-covariant (PC) dynamics

A phase-covariant (PC) qubit evolution corresponds to a family of dynamical maps that possess azimuthal symmetry with respect to rotations about the $z$ axis in the Bloch-ball representation. The most general PC dynamics is described by a qubit QME (46) with (Smirne et al., 2016):

[TABLE]

which represents a combination of relaxation, excitation and dephasing processes occurring with rates: $\gamma_{-}(t)$, $\gamma_{+}(t)$ and $\gamma_{z}(t)$, respectively; while $\sigma_{\pm}:=\frac{1}{\sqrt{2}}\left(\sigma_{x}\pm\mathrm{i}\sigma_{y}\right)$ are the transition operators.

Although dynamical generators commute within each of the RU (70) and PC (78) classes of dynamics, they do not generally commute in between the two. In fact, their common commutative subset corresponds to all unital PC evolutions for which $\mathcal{L}_{t}$ in Eqs. (70) and (78) coincide with $\gamma_{x}(t)=\gamma_{y}(t)=\frac{1}{2}\gamma_{+}(t)=\frac{1}{2}\gamma_{-}(t)$. Note that, the eternally non-Markovian model with decay rates specified in Eq. (74) is, in fact, both RU and PC, while the dephasing (RU) dynamics belongs to the PC class only when aligned along the $z$ direction, i.e., when $\boldsymbol{n}=\{0,0,1\}$ in Eq. (75).

As $\left\{\openone,\sigma_{+},\sigma_{-},\sigma_{z}\right\}$ equivalently constitute a basis for two-dimensional Hermitian matrices, the properties of PC dynamics can be determined analogously to the RU case. In particular, considering now $k=\{+,-,z\}$, a given PC evolution is CP-divisible iff all $\gamma_{k}(t)\geq 0$ at all times, while it is SS iff all $\Gamma_{k}(t)=\int_{0}^{t}\gamma_{k}(\tau){d}\tau\geq 0$ for any $t\geq 0$. Moreover, it is not hard to verify that the family of PC generators (78) is physical iff for all $t\geq 0$ (Smirne et al., 2016):

[TABLE]

where $\eta_{||}(t)\!:=\!e^{-\delta(t)}$, $\eta_{\perp}(t)\!:=\!e^{-\frac{1}{2}\left(\delta(t)-4\Gamma_{z}(t)\right)}$, $\kappa(t)\!:=\!e^{-\delta(t)}\int_{0}^{t}\!d\tau\,e^{\delta(\tau)}\!\left[\gamma_{+}(\tau)-\gamma_{-}(\tau)\right]$ and $\delta(t)\!:=\!\Gamma_{+}(t)+\Gamma_{-}(t)$.

Hence, similarly to the case of RU dynamics, as the physicality condition (79) is less restrictive than the SS requirement ($\Gamma_{k}(t)\geq 0$), there exist (commutative) PC dynamics that are not SS but still physically legitimate. The eternally non-Markovian model Hall et al. (2014), being both RU and PC, provides again an appropriate example.

Amplitude damping dynamics.

The most common example of PC dynamics, which is not RU, is the amplitude damping evolution that represents the pure relaxation process, i.e., spontaneous emission of a two-level (qubit) system (Breuer and Petruccione, 2002); and corresponds to the choice $\gamma_{+}(t)=\gamma_{z}(t)=0$ in Eq. (78), i.e.,

[TABLE]

Its canonical microscopic derivation stems from the Jaynes-Cummings interaction model (Breuer and Petruccione, 2002), $H_{\mathrm{int}}\propto\sigma_{+}\otimes\left(\sum_{k}g_{k}\hat{a}_{k}+g_{k}^{\star}\hat{a}_{k}^{\dagger}\right)$, in which the qubit is coupled to a cavity possessing Lorentzian spectral density (Li et al., 2010):

[TABLE]

where $\Delta$ describes the difference between qubit transition, $\omega_{0}$, and cavity central frequencies, while $\lambda$ represents the cavity spectral width. Crucially, such a model—after tracing out degrees of freedom of the cavity—leads to a qubit QME (46) with the dynamical generator (80), whose relaxation rate reads (Li et al., 2010):

[TABLE]

and does not exhibit a singular behaviour as long as the real part of the complex parameter $d:=\sqrt{(\lambda-\mathrm{i}\Delta)^{2}-2\gamma_{0}\lambda}$ is positive.

However, this is not the case in the on-resonance ($\Delta=0$), strong-coupling ($\gamma_{0}>\frac{\lambda}{2}$) regime, in which $d$ becomes purely imaginary, $d=\mathrm{i}\left|d\right|$ with $\left|d\right|=\sqrt{2\gamma_{0}\lambda-\lambda^{2}}$, so that the relaxation rate (82) simplifies to

[TABLE]

and diverges at every $t=\frac{2}{\left|d\right|}\!\left(\mathrm{arccot}\!\left(\frac{-\lambda}{\left|d\right|}\right)+n\pi\right)$ with $n\in\mathbb{N}^{+}$ (Breuer and Petruccione, 2002).

B.2 Counterexamples to rescalability and additivity

B.2.1 Non-rescalable qubit dynamics

An explicit example of non-rescalable qubit dynamics is provided by the Jaynes-Cummings model of spontaneous emission described just above, considered in the on-resonance, strong-coupling regime. It leads to an example of dynamics described in App. A.2 with dynamical generators being singular due to the damping rate (83) being divergent periodically in $t$.

Considering then the rescaled version of the amplitude-damping generator (80), i.e., $\mathcal{L}^{\prime}_{t}$ of Eq. (66), one may simply integrate the resulting QME (46) for any $t$ in order to explicitly determine the form of the rescaled dynamical map $\Lambda^{\prime}_{t}$ in Eq. (67). For instance, when setting $\gamma_{0}=3/2$ and $\lambda=1$ for simplicity in Eq. (83), the corresponding CJ-matrix (45), $\Omega_{\Lambda^{\prime}_{t}}$, may be explicitly computed and its non-zero eigenvalues read

[TABLE]

Crucially, although in the case of original dynamics (when $\alpha=1$) both $\lambda_{\pm}^{\mathrm{vals}}\geq 0$ for any $t\geq 0$, only for times before the occurrence of the first singularity, i.e., when $t<T=\frac{2}{\left|d\right|}\left(\mathrm{arccot}\!\left(\frac{-\lambda}{\left|d\right|}\right)+\pi\right)=\sqrt{2}\left[\pi-\arctan(\sqrt{2})\right]$, the eigenvalues are guaranteed to be real and non-negative independently of $\alpha$. In particular, for any $t\geq T$ one may easily find $\alpha\geq 0$ ($\alpha\neq 1$) such that the eigenvalues (84) take complex values with the QME being, in fact, not even integrable.

B.2.2 Non-additive qubit dynamics

In order to construct counterexamples to additivity, we consider dynamical generators given in Eq. (4) of the main text and make specific choices for their dissipation rates $\gamma_{1}(t)$ and $\gamma_{2}(t)$. We then solve the QME obtained after adding the generators as in Eq. (68) with some $\alpha,\beta\geq 0$, i.e.,

[TABLE]

in order to explicitly compute the corresponding family of maps $\{\Lambda^{\prime}_{t}\}_{t\geq 0}$. Crucially, we find in this way families containing maps that cease to be CPTP—with their CJ matrices, $\{\Omega_{\Lambda^{\prime}_{t}}\}_{t\geq 0}$ as defined in Eq. (45), exhibiting negative eigenvalues at some time instances.

In order to solve the QME (85), we choose a qubit operator basis:

[TABLE]

which allows us to use the matrix and vector representations for generators and states, respectively. As $\operatorname{Tr}[\hat{\mu}_{i}\hat{\mu}_{j}]=\delta_{ij}$, any generator $\mathcal{L}$ may then be represented by a matrix $\mathbf{{M}}$ with entries $M_{ij}=\operatorname{Tr}[\hat{\mu}_{i}\mathcal{L}[\hat{\mu}_{j}]]$, while any state $\rho$ by a vector $\mathbf{x}$ with components $x_{i}=\operatorname{Tr}[\rho\hat{\mu}_{i}]$ ($x_{0}=\operatorname{Tr}[\rho]/\sqrt{2}=\sqrt{2}/2$ by definition). The QME (85) is then equivalent to the set of linear, coupled differential equations:

[TABLE]

where $\mathbf{{M}}_{t}^{\prime}$ is the matrix representation of $\mathcal{L}_{t}^{\prime}$.

For our first example, we take dephasing and amplitude damping generators of Eq. (4) in the main text to be non-Markovian and semigroup, respectively, with

[TABLE]

where $\omega,\gamma>0$ are some fixed constants. Crucially, since for all $t\geq 0$:

[TABLE]

and similarly in the case of semigroup $\gamma_{2}(t)$, both generator families are SSC and hence rescalable, so that their additivity may be unambiguously considered.

Considering their non-negative linear combinations, $\mathcal{L}_{t}^{\prime}=\alpha\mathcal{L}^{(1)}_{t}+\beta\mathcal{L}^{(2)}_{t}$, one obtains generator families with

[TABLE]

Solving Eq. (87), one finds

[TABLE]

where $x_{0}(0)=\sqrt{2}/2$ and

[TABLE]

The four eigenvalues of the CJ matrices, $\Omega_{\Lambda^{\prime}_{t}}$, for the corresponding family of maps $\{\Lambda^{\prime}_{t}\}_{t\geq 0}$ read

[TABLE]

and are plotted in Fig. 5(a) for $\alpha=\beta=1$, $\omega=2$, and $\gamma=1$. For $t=\pi$, the integral in Eq. (92) evaluates to $I\approx 23.36$, and it is easy to check that two of the eigenvalues are negative. Hence, the evolution is clearly unphysical.

For the second example, we consider the symmetric case with the dissipation rates exchanged, i.e.:

[TABLE]

By the same argumentation as before both generators are SSC (and thus rescalable) and upon addition yield

[TABLE]

which after solving Eq. (87) leads to

[TABLE]

with again $x_{0}(0)=\sqrt{2}/2$ and now

[TABLE]

The four CJ eigenvalues this time read:

[TABLE]

and are plotted in Fig. 5(b) for $\alpha=\beta=1$, $\omega=1$, and $\gamma=1/2$. For $t=2\pi$ the integral (97) yields $I\approx-204.81$, and again two of the eigenvalues are negative, proving the evolution to be unphysical.

We repeat the above analysis, but considering this time the dissipation rates of either dephasing or amplitude-damping in Eq. (4) of the main text to have a functional form derived explicitly from an underlying microscopic model yielding non-Markovian dynamics.

Firstly, in an analogy to Eq. (88), we consider the dephasing rate to be specified by Eq. (77)—as if the qubit were coupled to a bosonic reservoir with an Ohmic-like spectrum—while the damping rate to be constant. In this case, we can solve Eq. (87) numerically at each $t$ for given parameter settings, in order to compute the corresponding CJ eigenvalues. These are plotted in Fig. 5(c) for $\alpha=\beta=1$, $\omega_{c}=1$, $s=4.5$, which corresponds to a super-Ohmic spectrum Haikka et al. (2013), and $\gamma=1$. We observe that the evolution becomes unphysical around $t=\pi/2$.

Secondly, we consider the symmetric case in an analogy to Eq. (94), this time setting the dephasing to be constant, while the damping rate to the one of Eq. (82)—derived basing on the Jaynes-Cummings microscopic model in which the qubit is coupled to a cavity with a Lorentzian frequency spectrum—in the off-resonant ($\Delta\neq 0$) regime. Again, we find the CJ eigenvalues by solving Eq. (87) numerically for fixed parameter values. These are plotted in Fig. 5(d) for $\alpha=\beta=1$, $\Delta=3$, $\lambda=0.05$, $\gamma_{0}=150$, and $\gamma=1$. We observe again that the evolution becomes unphysical, this time a bit before $t=\pi/2$.

Appendix C Microscopic derivations of QMEs

We consider the situation depicted in Fig. 3(a) of the main text, in which a system of interest and its environment evolve under closed dynamics determined by a time-invariant total (T) Hamiltonian—consisting of Hamiltonians associated with the system (S), the environment (E) and their interaction (I):

[TABLE]

C.1 Interaction and Schrödinger pictures

The interaction picture (IP), which we denote here with an over-bar, is then defined in the same manner for all operators and states acting on the system-environment Hilbert space, i.e., as

[TABLE]

for any given $O\in\mathcal{B}(\mathcal{H}_{S}\otimes\mathcal{H}_{E})$ that is specified in the Schrödinger picture (SP).

In contrast, a general dynamical map, $\Lambda_{t,t_{0}}$, that describes the evolution of solely the system between the initial time $t_{0}$ and some later $t$ transforms from SP to IP (and vice versa) as:

[TABLE]

where by

[TABLE]

we denote the unitary transformation induced by the system free Hamiltonian, $H_{S}$. However, as we consider throughout this work dynamical maps that commence at zero time ($t_{0}=0$), see Eq. (44), Eq. (101) simplifies to

[TABLE]

which allows us to explicitly compute how the corresponding dynamical generators of $\Lambda_{t}$ and $\bar{\Lambda}_{t}$ transform between the SP and IP.

In particular, defining the IP-based dynamical generator in accordance with Eq. (47) as

[TABLE]

and substituting for $\bar{\Lambda}_{t}$ according to (103), we obtain

[TABLE]

where in the last line we have used the definition of $\mathcal{U}_{t}^{S}$ (102), and accordingly defined the SP-based dynamical generator, i.e., as in Eq. (47):

[TABLE]

C.2 QME in the integro-differential form

The von Neumann equation describing the unitary evolution of the closed system-enviroment (SE) system, i.e., Eq. (5) of the main text, in the IP reads:

[TABLE]

Assuming the SE to initially be in a product state,

[TABLE]

with $\rho_{E}$ being a stationary state of the environment that satisfies $\left[\bar{H}_{E}(t),\rho_{E}\right]=\left[H_{E},\rho_{E}\right]=0$, one may write the integral of Eq. (109) as:

[TABLE]

Tracing out the environment in Eq. (109), so that its l.h.s. reduces to $d\bar{\rho}_{S}(t)/dt$ and substituting into its r.h.s. for $\bar{\rho}_{SE}(t)$ according to Eq. (111), one arrives at the integro-differential equation describing the system density matrix in the IP at time $t$:

[TABLE]

The first term in Eq. (C.2) may be dropped, as without loss of generality one may impose

[TABLE]

by shifting the zero point energy of Hamiltonians, i.e., by changing $H_{I}$ and $H_{S}$ as follows

[TABLE]

so that condition (113) is ensured, given $\left[H_{E},\rho_{E}\right]=0$, without affecting the total Hamiltonian $H_{T}$ in Eq. (99).

As a result, we obtain the QME in its integro-differential form that does not involve any approximations, but only assumes Eq. (110) with $\left[H_{E},\rho_{E}\right]=0$,

[TABLE]

and constitutes Eq. (6) of the main text.

C.3 QME in the time-local form

The QME (115) despite being compact and exact is typically not of much use, as it involves the full system-environment state and a time-convoluted integral. Nevertheless, one may always formally rewrite it as a function of the system state at a given time.

After integrating the closed von Neumann dynamics (109), one should arrive at

[TABLE]

with the unitary rotation being formally defined as a time-ordered exponential:

[TABLE]

Now, as the reduced state of the system is obtained at any time by tracing out the environment, the dynamical map, $\bar{\Lambda}_{t}$, associated solely with the system evolution in the IP may be identified as

[TABLE]

Hence, if at a given $t$ one can compute the inverse of the dynamical map, i.e., $\bar{\Lambda}_{t}^{-1}$ such that $\rho_{S}(0)=\bar{\Lambda}_{t}^{-1}\!\left[\bar{\rho}_{S}(t)\right]$, as well as its the time-differential $\dot{\bar{\Lambda}}_{t}$, which is now formally determined by Eq. (115) as

[TABLE]

one may equivalently rewrite the exact QME (115) into its time-local form (in the IP):

[TABLE]

where $\bar{\mathcal{L}}_{t}=\dot{\bar{\Lambda}}_{t}\circ\bar{\Lambda}_{t}^{-1}$ is the IP-based dynamical generator defined in Eq. (104).

C.4 QME in the Schrödinger picture and the environment-induced generator

We rewrite the QME (120) in the SP as

[TABLE]

where in accordance with Eq. (107) the IP-based dynamical generator must be transformed to

[TABLE]

As a result, we arrive at the time-local QME as stated in Eq. (7) of the main text:

[TABLE]

and identify with help of Eq. (122):

[TABLE]

as the environment-induced dynamical generator, which can be then associated solely with the impact of the environment on the system.

Still, as $\tilde{\mathcal{L}}_{t}$ generally contains both Hamiltonian and dissipative parts, i.e., $\tilde{\mathcal{L}}_{t}\left[\bullet\right]=-\mathrm{i}\left[H(t),\bullet\right]+\mathcal{D}_{t}\left[\bullet\right]$, one may conveniently rewrite the QME (121) as Breuer and Petruccione (2002); Rivas and Huelga (2012):

[TABLE]

which allows to explicitly identify $H(t)$ as the environment-induced Hamiltonian correction to the system free evolution, e.g., representing the Lamb shifts when describing atom-light interactions Breuer and Petruccione (2002). Whereas, $\mathcal{D}_{t}$ in Eq. (125) may then be entirely associated with the dissipative impact of the environment.

Lastly, let us emphasise that when investigating whether by adding generators associated with different environments—i.e., the generators $\tilde{\mathcal{L}}_{t}$ in Eq. (121) obtained by considering the impact of each environment separately—one reproduces the correct dynamics, it is equivalent to consider all the generators in the IP.

As the IP-based generators, $\bar{\mathcal{L}}_{t}$, are linearly related to the environment-induced ones, $\tilde{\mathcal{L}}_{t}$, by Eq. (124), the vector spaces formed by their families must be isomorphic. Hence, the notions of rescalability and additivity of generator families, discussed in Sec. II.2 of the main text and App. A above, are naturally carried over between the two, e.g., for any $\alpha,\beta,t\geq 0$:

[TABLE]

with Eq. (124) relating all $\tilde{\mathcal{L}}^{\mathsf{x}}_{t}=\mathcal{U}_{t}^{S}\circ\bar{\mathcal{L}}^{\mathsf{x}}_{t}\circ\mathcal{U}_{t}^{S\dagger}$ for each $\mathsf{x}=\{\prime,(1),(2)\}$.

C.5 $H_{S}$-covariant dynamics

Although the dynamical generators $\tilde{\mathcal{L}}_{t}$ defined in Eq. (124) arise due to the presence of the environment, their form may still strongly depend on the system free Hamiltonian $H_{S}$. Thus, $\tilde{\mathcal{L}}_{t}$ and, in particular, both its Hamiltonian $H(t)$ and dissipative parts $\mathcal{D}_{t}$ in Eq. (125) cannot be generally associated with the properties of just the environment and the interactions. In fact, only in very special cases the form of $\mathcal{L}_{t}$ can be derived from $H_{I}$, $H_{E}$, and $\rho_{E}$.

An important example is provided when the system and interaction Hamiltonians alone commute:

[TABLE]

As the global unitary $\bar{U}_{SE}$ in Eq. (117) then also commutes with $H_{S}$, $\left[\bar{U}_{SE},H_{S}\right]=0$, the IP-based map $\bar{\Lambda}_{t}$ in Eq. (118) is assured to be $H_{S}$-covariant, i.e., to commute with any $H_{S}$-induced unitary (102) (and so must trivially the SP-based $\Lambda_{t}$), so that for any $s\geq 0$ (Holevo, 1993):

[TABLE]

As noted above, the $H_{S}$-covariance must be naturally inherited by the IP-based dynamical generators $\bar{\mathcal{L}}_{t}$ (Holevo, 1993; Vacchini, 2010), which, in turn, must then coincide with the environment-induced ones, with $\tilde{\mathcal{L}}_{t}=\bar{\mathcal{L}}_{t}$ in Eq. (124). As a result, the form of $\tilde{\mathcal{L}}_{t}$ must then, indeed, be independent of $H_{S}$.

C.6 Externally modifying the system Hamiltonian

In general, a modification of the system free Hamiltonian $H_{S}$ may affect both the Hamiltonian and the dissipative parts of the generator $\tilde{\mathcal{L}}_{t}$ in Eq. (125). Nevertheless, let us consider a transformation:

[TABLE]

with $V(t)$ being an arbitrary (potentially time-dependent) Hermitian operator. Crucially, by considering particular commutation relations satisfied by the microscopic Hamiltonians of Eq. (99) and $V(t)$, one may identify two important cases for which the microscopic rederivation of the QME (125) can be bypassed—with the impact of $V(t)$ being directly accountable for at the level of the QME:

$\left[H_{S},H_{I}\right]=0$, $\forall_{t\geq 0}\!\!:\left[V(t),H_{I}\right]=0$.

If the system Hamiltonian commutes with the interaction Hamiltonian—so that the dynamics is $H_{S}$-covariant—and so does the perturbation $V(t)$ for all $t$, then the modified dynamics must also be $H^{\prime}_{S}$-covariant, as $\left[H^{\prime}_{S}(t),H_{I}\right]=0$ at all times. Hence, the form of $\tilde{\mathcal{L}}_{t}$ in Eq. (123) is unaffected by the modification of $H_{S}$, remaining fully determined by $H_{E}$, $H_{I}$ and ${\rho_{E}}$. Moreover, the new dynamics is then correctly described by simply replacing $H_{S}$ with $H^{\prime}_{S}(t)$ in Eq. (123) (or Eq. (125)).

$\left[H_{S},H_{I}\right]\neq 0$, $\forall_{t\geq 0}\!\!:\left[V(t),H_{I}\right]=\left[V(t),H_{S}\right]=0$.

The above conclusion also holds when dealing with non-$H_{S}$-covariant dynamics, given that $V(t)$ commutes with both the interaction and the system Hamiltonian. As then $\left[H_{S},H_{I}\right]\neq 0$, the generator $\tilde{\mathcal{L}}_{t}$ in principle depends on $H_{S}$. However, without affecting the total Hamiltonian $H_{T}$ in Eq. (99) and hence the dynamics, we may redefine the interaction Hamiltonian as $H_{I}^{\prime}:=H_{I}+H_{S}$, pretending the system Hamiltonian to be absent. In such a fictitious picture, the QME (123) possesses just the second term with $\tilde{\mathcal{L}}^{\prime}_{t}$ now being derived based on $H_{I}^{\prime}$. As importantly $\left[V(t),{H}_{I}^{\prime}\right]=0$ is fulfilled at all times, it becomes clear that the dynamics must be $V(t)$-covariant. Hence, the perturbation must lead to a QME that may be equivalently obtained by simply adding $V(t)$ to the Hamiltonians in Eq. (125)—even though $H(t)$ and $\mathcal{D}_{t}$ non-trivially depend on the original $H_{S}$ (but not on $V(t)$).

Lastly, let us note that in case of $H_{S}$-covariant dynamics and $\left[H_{S},H_{I}\right]=0$, one may play a similar trick in order to deal with the case when $\left[V(t),H_{S}\right]=0$ but $\left[V(t),H_{I}\right]\neq 0$, so that the modified dynamics is no longer guaranteed to be $H^{\prime}_{S}$-covariant. By redefining the interaction Hamiltonian this time as $H_{I}^{\prime}:=V(t)+H_{I}$, which importantly commutes with $H_{S}$, it becomes clear that the $H_{S}$-covariance must be preserved. Nonetheless, although the dynamical generator $\tilde{\mathcal{L}}_{t}$ in Eq. (123) remains then independent of $H_{S}$, the form of $\tilde{\mathcal{L}}_{t}$ may depend on $V(t)$ and must thus be rederived, i.e., based now on $H_{I}^{\prime}$.

Appendix D Microscopic validity of generator addition

D.1 Weak-coupling regime

Below, we prove Lemma 4 stated in the main text, in particular, we show that under weak coupling the cross-term in Eq. (8) can always be assumed to vanish. Hence, in accordance with Obs. 1, it is then valid to add dynamical generators corresponding to each individual environment, without necessity to rederive the overall QME.

The following proof can be regarded as an extension of the argumentation found in Ref. Cohen-Tannoudji et al. (1998); Schaller (2015), which applies to the more stringent regime in which the Born-Markov approximation holds.

Firstly, we perform the operator Schmidt decomposition of each interaction Hamiltonian (indexed by $i$) Bengtsson and Życzkowski (2006):

[TABLE]

where $\{A_{i;k}\}_{k}$ and $\{B_{k}^{E_{i}}\}_{k}$ form then sets of Hermitian operators that act separately on the system and corresponding environment subspaces, i.e., $\mathcal{H}_{S}$ and $\mathcal{H}_{E_{i}}$, respectively.

As noted above, this decomposition preserves its tensor-product structure in the IP, which is now defined according to Eq. (100) with the free system-environment Hamiltonian incorporating multiple environments, $H_{S}+\sum_{i}H_{E_{i}}$. Hence, carrying out here the analysis in the IP for compactness of the expressions, we rewrite the general and exact QME (8) of the main text, which describes a system interacting with multiple environments, as

[TABLE]

where the second term above is the crucial, inter-environment cross-term whose absence assures the validity of generator addition.

There exist various approaches to obtain simplified forms of QMEs for the weak-coupling regime de Vega and Alonso (2017); Breuer and Petruccione (2002); Schaller (2015); Rivas and Huelga (2012). Here, in order keep the derivation general and emphasise necessary requirements for our arguments to apply, we assume that the appropriate QME under weak coupling is derived after approximating the global system-environments state at every time $t\geq 0$ as

[TABLE]

where $\rho_{S}(t)=\operatorname{Tr}_{E}\rho_{SE}(t)$ is the reduced state of the system at time $t$. Although in the weak-coupling approximations Breuer and Petruccione (2002); Rivas and Huelga (2012) the separable state of each environment in Eq. (133) is frequently taken to be its reduced state at time $t$, i.e., $\varrho_{E_{i}}(t)\equiv\rho_{E_{i}}(t):=\operatorname{Tr}_{\neg E_{i}}\rho_{SE}(t)$, in what follows it can be chosen arbitrarily—as long as for all environments (labelled by $i$) $\varrho_{E_{i}}(0)=\rho_{E_{i}}$ to maintain consistency with the derivation in App. C.2.

Let us stress that the tensor-product ansatz of Eq. (133) for the system-environments state is employed only to obtain the form of the QME valid under the weak coupling, and does not force the solutions of the QME to actually be separable states. In particular, the resulting QME, before tracing out environmental degrees of freedom, can yield upon integration states $\bar{\rho}_{SE}(t)$ that strongly deviate from the form (133) already at moderate times $t$, even though the validity of the dynamics—and, hence, the QME employed—is still assured by weak coupling Rivas et al. (2010).

Thanks to the condition (133), the crucial cross-term within the exact QME (132) can be reexpressed as

[TABLE]

where now

[TABLE]

is the two-bath correlation function that is independent of the reduced system state, being evaluated only on $\bar{\varrho}_{E_{ij}}(t):=\bar{\varrho}_{E_{i}}(t)\otimes\bar{\varrho}_{E_{j}}(t)$. Note that, as $i$ and $j$ are just labels of distinct environments, the correlation function (135) is symmetric with $\mathscr{C}_{[i;k][j;l]}(t,s;s^{\prime})=\mathscr{C}_{[j;l][i;k]}(s,t;s^{\prime})$.

Moreover, Eq. (133) assures all the correlation functions (135) to factorise, so that for any $s,s^{\prime},t\geq 0$:

[TABLE]

reducing to a product of single-bath correlation functions:

[TABLE]

Furthermore, as the two-bath correlation function in Eq. (134) factorises to Eq. (137) with $s=s^{\prime}$, the whole cross-term (134) is guaranteed to vanish whenever at all times $t\geq 0$ for each environment (labelled by $i$) and its each operator $B_{k}^{E_{i}}$ (labelled by $k$):

[TABLE]

Crucially, the condition Eq. (139) can always be ensured by shifting adequately the interaction and the system Hamiltonians without affecting the total Hamiltonian (99) (similarly to Eq. (114) of App. C.2). In particular, one can redefine the system and each interaction Hamiltonian to be generally time-dependent and read:

[TABLE]

so that the decomposition (130) of the interaction Hamiltonian for each environment becomes

[TABLE]

with the new correlation function (139) identically vanishing by construction, as for any $t\geq 0$:

[TABLE]

For consistency, le us also note that the necessary requirement $\operatorname{Tr}_{E}\!\left\{\bar{H}_{I_{i}}^{\prime}(t)\rho_{E_{i}}\right\}=0$, introduced in App. C.2, is then trivially fulfilled for every environment. Decomposing $H^{\prime}_{I_{i}}(t)$ according to Eq. (142) and remembering that $[H_{E_{i}},\rho_{E_{i}}]=0$ for each $i$, one gets (in the IP):

[TABLE]

as each $\mathscr{C}_{i;k}^{\prime}(0)=0$ is zero by Eq. (143).

Note that, in particular, the above argumentation holds for all QMEs derived using the time-convolutionless approach Breuer et al. (2001) up to the second order in all the interaction parameters—in which case the QME (131) is from the start assumed to exhibit a time-local form, rather than involve a time-convolution integral.

On the other hand, the most conservative Born-Markov approximation discussed in Ref. Cohen-Tannoudji et al. (1998); Schaller (2015) enforces every $\varrho_{E_{i}}(t)$ in Eq. (133) to be at all times the initial, stationary state $\rho_{E_{i}}$ of each environment. As a result, all the single-bath correlation functions, $\mathscr{C}_{i;k}(t,s)$ in Eq. (138), become then $t$- and $s$-independent due to $[H_{E_{i}},\rho_{E_{i}}]=0$, and identically vanish by Eq. (143). Hence, then trivially $\mathscr{C}_{[i;k][j;l]}=\mathscr{C}_{i;k}\mathscr{C}_{j;l}=0$ at all times.

D.2 Commutativity of microscopic Hamiltonians

Here, we provide the proof of Lemma 5 stated in the main text, which assures that dynamical generators associated with each individual environment can be simply added at the QME level, if the interaction Hamiltonians commute between each other and with the system Hamiltonian. This condition corresponds to the $\text{II}\cap\text{IS}$ region in the Venn diagram of Fig. 4—marked ‘Yes’ to indicate the validity of generator addition.

Let us note that whenever for all $i$ and $j$:

[TABLE]

with $H_{I}:=\sum_{i}H_{I_{i}}$ being the full interaction Hamiltonian, the global unitary dynamical operator (117) can be decomposed in the SP, as follows

[TABLE]

so that at the level of the corresponding unitary maps:

[TABLE]

where we have defined $\mathcal{U}_{t}^{IE_{i}}[\bullet]:=e^{-\mathrm{i}(H_{E_{i}}+H_{I_{i}})t}\bullet e^{\mathrm{i}(H_{E_{i}}+H_{I_{i}})t}$, and by $\prod$ we denote also the conjugation of multiple maps, i.e., for a given set of maps $\{\Lambda_{i}\}_{i}$:

[TABLE]

As a result, after straightforwardly generalising Eq. (118) to multiple environments and transforming it to the SP, we can generally write the system reduced state at a given time $t$ as

[TABLE]

where in (151) by $\cdots$ we mean repeating the procedure for the $2$nd environment, and similarly in (152) for all the other environments. The overall dynamical (SP-based) map $\tilde{\Lambda}_{t}=\prod_{i}\tilde{\Lambda}_{t}^{(i)}$ is $H_{S}$-independent and constitutes a composition of maps $\tilde{\Lambda}_{t}^{(i)}[\bullet]:=\operatorname{Tr}_{E_{i}}\!\left\{\mathcal{U}_{t}^{SE_{i}}[\bullet\otimes\rho_{E_{i}}]\right\}$, each of which describing the impact of the $i$th environment.

As discussed in App. C.5 above, the condition $[H_{I},H_{S}]=0$ in Eq. (145) ensures the dynamics to be $H_{S}$-covariant—as explicitly manifested in Eq. (153) in which $\mathcal{U}_{t}^{S}$, thanks to commuting with all $\mathcal{U}_{t}^{SE_{i}}$, commutes also with all the $\tilde{\Lambda}_{t}^{(i)}$ maps. On the other hand, as all $\mathcal{U}_{t}^{SE_{i}}$ commute between one another due to $[H_{I_{i}},H_{I_{j}}]=0$ in Eq. (145), the overall map $\tilde{\Lambda}_{t}$ in Eq. (153) could have been constructed by composing the maps $\tilde{\Lambda}_{t}^{(i)}$ in any order. Hence, all maps originating from interactions with different reservoirs commute, i.e., for all $i$ and $j$:

[TABLE]

As a consequence, while due to the $H_{S}$-covariance all the dynamical generators induced by separate environments, i.e., $\tilde{\mathcal{L}}^{(i)}_{t}$ in Eq. (124) indexed now by $i$, coincide with their corresponding IP-based generators, they also add at the level of the QME (123) thanks to Eq. (154). Computing explicitly the dynamical generator induced by all the environments together, i.e., $\tilde{\mathcal{L}}_{t}$ associated with the overall map $\tilde{\Lambda}_{t}$ in (153), we have

[TABLE]

where $\tilde{\mathcal{L}}^{(i)}_{t}$ is the generator corresponding to the interaction with the $i$th environment alone, and we have used the commutativity of the maps (154).

Appendix E Spin-magnet model

We provide here the details and explicit form of relevant quantities for the calculations presented in Sec. IV.2 of the main text, discussing the counterexamples to the commutativity assumptions based on the spin-magnet model.

E.1 $\text{IS}\cap\text{IE}$ commutativity assumption

We solve the equations of motion (38) that describe the Bloch vector dynamics in order to obtain an explicit form of the $\mathbf{{R}}$-matrix in Eq. (27), which is then parametrised by magnetisations of the two magnets, $m_{1}$ and $m_{2}$, and reads:

[TABLE]

When only the first magnet is present ($g_{2}=0$), the above expression reduces to

[TABLE]

while, when in contact with only the second magnet ($g_{1}=0$), it becomes

[TABLE]

We also analytically compute the time-derivative $\dot{\mathbf{{R}}}_{12}$ (as well as $\dot{\mathbf{{R}}}_{1}$ and $\dot{\mathbf{{R}}}_{2}$, after setting $g_{2}=0$ and $g_{1}=0$, respectively) which we, however, do not include here due its cumbersome form.

With the exact expressions (156-158) at hand, we can explicitly write each affine map $\mathbf{{D}}_{t}^{(\mathsf{x})}$ with $\mathsf{x}=\{12,1,2\}$ according to Eq. (39) of the main text, i.e., as an average of the corresponding $\mathbf{{R}}_{\mathsf{x}}$ over the Gaussian distributions $p(m_{i})$ of fixed variance $\sigma_{i}$ in Eq. (18); and similarly in case of the time-derivatives $\dot{\mathbf{{D}}}_{t}^{(\mathsf{x})}$ by averaging $\dot{\mathbf{{R}}}_{\mathsf{x}}$.

We perform the averaging integrals numerically after fixing the parameters $g_{1}$, $g_{2}$, $\sigma_{1}$, $\sigma_{2}$, and the time $t$. As a result, we obtain the expressions of dynamical generators $\mathbf{{L}}^{(1)}_{t}$, $\mathbf{{L}}^{(2)}_{t}$, and $\mathbf{{L}}^{(12)}_{t}$ by substituting into $\mathbf{{L}}^{(\mathsf{x})}_{t}=\dot{\mathbf{{D}}}_{t}^{(\mathsf{x})}(\mathbf{{D}}_{t}^{(\mathsf{x})})^{-1}$ the relevant affine maps and their time-derivatives.

For instance, when taking $g_{1}\!=\!g_{2}=2$ and $\sigma_{1}\!=\!\sigma_{2}=1$, we obtain at $t=0.5$:

[TABLE]

and

[TABLE]

which provides the desired example of $\mathbf{{L}}^{(12)}_{t}\neq\mathbf{{L}}^{(1)}_{t}+\mathbf{{L}}^{(2)}_{t}$.

E.2 $\text{II}\cap\text{IE}$ commutativity assumption

We observe that, in order to solve the equations of motion (43) stated in the main text, which describe the dynamics of the Bloch vector in the IP, i.e., $\bar{\mathbf{r}}(t)=\mathbf{{R}}_{S}^{-1}(t)\,\mathbf{r}(t)$, it is convenient to move to a rotating frame defined as $\check{\mathbf{r}}(t):=\mathbf{{V}}(t)\,\bar{\mathbf{r}}(t)$, where

[TABLE]

is an orthogonal matrix such that $\mathbf{{V}}(t)=\mathbf{{P}}\,\mathbf{{R}}_{S}(t)$ with $\mathbf{{P}}=\textrm{diag}\{1,-1,1\}$ and $\mathbf{{R}}_{S}(t)$ being the SO(3) representation of qubit unitary $U_{S}(t)=e^{-\mathrm{i}H_{S}t}$ induced by the system free Hamiltonian $H_{S}=\frac{1}{2}\omega\sigma_{x}$ of Eq. (40). Hence, $\check{\mathbf{r}}(t)=\mathbf{{P}}\,\mathbf{r}(t)$ can be interpreted as the Bloch vector in the SP with the $y\to-y$ coordinate inverted.

Defining also $\gamma=g_{1}m_{1}+g_{2}m_{2}$ for compactness, we obtain a simpler set of equations of motion:

[TABLE]

which can be explicitly solved, yielding

[TABLE]

We then construct the $\mathbf{{R}}$-matrix determining the affine map $\mathbf{{D}}_{t}$ in Eq. (27) by transforming back the above $\check{\mathbf{{R}}}$-matrix to the IP, so that

[TABLE]

We do not enclose here the explicit forms, but, as in the previous example, we also compute all the relevant $\bar{\mathbf{{R}}}_{\mathsf{x}}$ and $\dot{\bar{\mathbf{{R}}}}_{\mathsf{x}}$ in the IP, with $\mathsf{x}=\{12,1,2\}$, which allow us to obtain the integral expressions for the corresponding affine maps, $\bar{\mathbf{{D}}}_{t}^{(\mathsf{x})}$, and their time-derivatives, $\dot{\bar{\mathbf{{D}}}}_{t}^{(\mathsf{x})}$. Again, we choose the initial magnetisations of both magnets to be Gaussian distributed with both $p(m_{i})$ as in Eq. (18) of fixed variance $\sigma_{i}$.

As before, we perform the averaging integrals numerically, after fixing the model parameters—now: $\omega$, $g_{1}$, $g_{2}$, $\sigma_{1}$, $\sigma_{2}$, and the time $t$—in order to obtain numerical expressions for all $\bar{\mathbf{{L}}}^{(\mathsf{x})}_{t}=\dot{\mathbf{{D}}}_{t}^{(\mathsf{x})}(\mathbf{{D}}_{t}^{(\mathsf{x})})^{-1}$.

For example, when choosing $g_{1}=g_{2}=2$, $\sigma_{1}=\sigma_{2}=1$ and $\omega=2$, we obtain at $t=0.5$:

[TABLE]

giving

[TABLE]

On the other hand, we find in presence of both magnets:

[TABLE]

what, thus, provides an instance of $\bar{\mathbf{{L}}}^{(12)}_{t}\neq\bar{\mathbf{{L}}}^{(1)}_{t}+\bar{\mathbf{{L}}}^{(2)}_{t}$.