Satisfying fermionic statistics in the modeling of non-Markovian dynamics with one-electron reduced density matrices

Lindlbad derived a general equation to model environmental effects in open, time-dependent quantum systems while preserving the non-negativity of the density matrix.^1–3 However, when used in conjunction with the 1-particle reduced density matrix, or 1-RDM, fermion statistics are not inherently preserved.^4–7 Additional constraints have been derived on the Lindbladian operators to ensure that the 1-RDM obeys the Pauli exclusion principle.⁸ Here, we extend these constraints to treat non-Markovian dynamics.

Despite the prevalence of non-Markovian dynamics in a variety of physical systems such as information transfer in qubits and energy transfer in photosynthetic light harvesting complexes, general and accurate treatment of such dynamics remains a challenge.^9–20 Methodologies for non-Markovian dynamics include the quantum jump methods, quantum trajectory approaches, and master equation methods.^21–36 Recently we developed an extension of Lindblad’s formalism to treat non-Markovian dynamics by taking an ensemble of Lindbladian trajectories (ELT).³⁷ The advantage of the ELT method is that it provides a theoretically complete description of non-Markovian dynamics without violation of the complete positivity of the density matrix.³⁷ With a proper choice of the Lindbladian channels it can in principle match the system density matrix obtained from the exact system-bath density matrix.

Here we extend the ELT theory for non-Markovian dynamics to treat systems with multiple fermions. We show that the constraints that we previously derived for the Lindblad equation can be applied to each Lindblad trajectory of the ELT method to maintain the fermion statistics. In particular, we use the convexity of the constrained density matrices to prove that the ensemble inherits the constraints of the individual trajectories. These constraints are critical because without them the many-fermion system can collapse into a many-boson-like system in potentially extreme violation of the Pauli exclusion principle. The ELT with fermion constraints offers a theoretically complete description of non-Markovian dynamics without violation of the Pauli principle.

After developing the theoretical extension of the ELT method in Sec. II C, we apply the resulting ELT method in Sec. III. We apply the ELT method to three distinct systems of two fermions in three levels. While the electrons violate the fermion statistics without the constraints, correct fermion behavior is recovered with the constraints.

The constraints on the Lindbladian matrices to ensure fermionic statistics and the generalization of the Lindblad equation to the non-Markovian regime are reviewed in Secs. II A and II B, respectively. In Sec. II C we generalize the ELT method to treat systems with multiple fermions. While the derivations are presented in the notation of first quantization, all of the results can also be expressed without modification in the notation of second quantization.

The Lindblad equation for the 1-RDM is given by,

where ¹H is the one-body Hamiltonian, ¹D is the 1-RDM, $L(1D,1Cj)=1Cj1D1Cj†−12{1Cj†1Cj,1D}$ where ¹C_j are the one-body Lindbladians.^1–3 

The use of the 1-RDM requires additional constraints referred to as the N-representability conditions: (i) Hermiticity, (ii) normalization to N particles, and (iii) two linear matrix inequalities:

where ¹Q is the 1-hole RDM.^20,38–42 The nonnegativity of the 1-particle and 1-hole RDMs is equivalent to the Pauli exclusion principle that the eigenvalues of the 1-RDM must lie in the interval [0, 1].

In Ref. 8 we derived a constraint on the Lindbladian (bath) operators that is necessary and sufficient for the 1-RDM to obey the fermion statistics expressed in the Pauli exclusion principle. For a single dissipative channel, the Lindbladian matrix must be normal,⁸ 

or if there are multiple dissipative channels, the constraint can be written as a summation over multiple commutators. As previously shown, the simple constraint is a subset of the general constraint since if all Lindbladians are normal, then the summation is trivially zero. However, with multiple channels the Lindbladians have greater flexibility and are not required to be normal as long as cancellation among terms causes the sum to be zero.⁸ 

Taking the ensemble average of Lindbladian trajectories originating from different points in the system’s history produces a new density matrix,

where ω(τ_i) are statistical weights, τ_i are lag times, and $eL(τi)$ are the Lindbladian propagators.³⁷ Since this method takes positive trajectories from many different points in the system’s history, the system memory is built into the dynamics. Through density matrices from different times in the history the system has the opportunity to regain information that was previously lost to the surroundings. As shown in previous work, this method provides a full account of the system’s memory while maintaining the positivity of the density matrix.³⁷ 

The density matrix from the ELT method is an ensemble of individual Lindbladian trajectories. Since each trajectory produces a positive semidefinite density matrix, the new density matrix will also be positive semidefinite. Due to this relationship and the convexity of the sets of both the 1-electron and 1-hole RDMs, enforcing the Pauli principle on each trajectory is equivalent to enforcing the principle for the ensemble of trajectories. For completeness a formal proof of this result is given in the  Appendix. The constraint on the Lindbladians is thus similar to the constraint for the multiple dissipative channels case,

It should be noted that due to the multiple trajectories in the ensemble, this constraint is analogous to the more general constraint previously derived. It is not necessary that each Lindbladian matrix be normal, only that the sum of all the commutators vanishes. For example, two Lindbladian channels could be paired so that each channel has a bath matrix that is the adjoint of the other channel’s bath matrix, causing the sum in Eq. (6) to vanish.

The ELT method was previously applied to the Jaynes-Cummings model, which is a two-level system interacting with the cavity field.^37,43,44 To verify the extended version of the ELT method to capture fermion statistics, we employ here a three-level system with two fermions.²¹ There are three primary ways of arranging three different energy levels, as shown in Fig. 1. The ELT method with and without constraints on the Lindbladian is employed to calculate occupation numbers of the three states for each of the three systems.

For the three systems shown in Fig. 1 the following Lindbladian matrices are used,

and

respectively for the unconstrained calculation. For the constrained calculation, the number of trajectories is doubled to include an additional trajectory for the adjoint of each dissipative operator such that $∑i[Ci,Ci†]=0$⁠.

The Lindblad equation in Eq. (1) is used to propagate both the 1-particle and 1-hole density matrices from each time in the history of the ELT method. The density matrices are obtained through a fourth-fifth-order Runge-Kutta method for solving differential equations in Maple.^45,46 For all calculations the initial conditions consist of one fermion in state |2⟩, one fermion in state |1⟩, and one hole in state |0⟩.

To obtain a set of weights that contains non-Markovian dynamics in the 3-level system, we use weights that are optimized to reproduce the complete non-Markovian time dynamics of the two-level detuned Jaynes-Cummings model with parameters reported in Ref. 37. The optimization of the weights is performed in Maple^45,46 by a sequential quadratic programming algorithm.⁴⁷ To recap the parameters, γ₀ = 1.091 is inversely proportional to the timescale of system decay, λ = 0.3γ₀ is inversely proportional to the timescale of bath decay, and Δ = 2.4γ₀ is the amount of detuning. For all calculations, $1λ$ is used as the unit of time. Two sets of weights are calculated, one set for the 1-particle RDM and one set for the 1-hole RDM, as shown in Fig. 2.

Although these weights are obtained from numerical optimization, it is interesting to note that when accessing memory further back, for τ > 10 the weights of the 1-particle RDM and 1-hole RDM trajectories alternate in their ebbs and flows. This alternation corresponds to the flow of energy into the bath (1-particle RDM) and from the bath (1-hole RDM). This intuitively makes sense as the 1-particle RDM and 1-hole RDM are complimentary entities related by the Pauli exclusion principle.

The first system is the 3-level ladder system, where the occupied state |2⟩ decays to the occupied state |1⟩, which in turn decays into the unoccupied state |0⟩ as shown in Fig. 3. The second system is the 3-level ∧ system where the occupied state |2⟩ decays to the occupied state |1⟩ and the unoccupied state |0⟩, as shown in Fig. 4. The third system is the 3-level ∨ system where the two occupied states |2⟩ and |1⟩ decay into the unoccupied state |0⟩, as shown in Fig. 5.

In all three configurations, it can be seen that when the Lindbladians are unconstrained, the fermions behave as bosons and pile up in a given state, violating the Pauli exclusion principle. However, when the constraint on the Lindbladian matrices is applied, the fermions obey the Pauli exclusion principle and each of the state occupations remains physically reasonable.

Many physical systems where non-Markovian dynamics are prevalent are also many-electron systems. For example, a method which models electron transport through a photosynthetic light harvesting complex needs to preserve fermion statistics. Moreover, if this system is embedded in a complex protein environment, it also must capture non-Markovian dynamics. Here, we generalized a non-Markovian method based on an ensemble of Lindbladian trajectories to treat systems of multiple fermions. Lindbladian operators from each time in the system’s history are constrained to preserve fermion statistics. These constraints were verified to maintain fermion statistics through a variety of three-level systems with two electrons.

Many non-Markovian methods rely on the Lindbladian formalism, such as Monte-Carlo quantum jump methods or several master equation methods.^21,48 The derived constraints on the Lindbladian bath operators are potentially useful in not only the ELT method but other non-Markovian methods. The present work offers an initial but important step towards improving the treatment of non-Markovian dynamics in many-fermion quantum systems.

D.A.M. gratefully acknowledges the National Science Foundation (NSF) Grant No. CHE-1152425, and the United States Army Research Office (ARO) Grant No. W911NF-16-1-0152.

The set of density matrices with eigenvalues in the interval [0, 1] is convex, and hence, the ensemble average of matrices in the set produces another matrix in the set. Because this result is central to the present work, we present a proof of the set’s convexity in the following theorem.