A quantum open system model of molecular battery charged by excitons

TL;DR

This paper presents a quantum Markovian model for a molecular battery charged by excitons, providing insights into energy transfer efficiency and stability in biological and organic systems.

Contribution

It introduces an analytically solvable quantum open system model for molecular batteries charged by non-equilibrium excitons, grounded in quantum Markovian dynamics.

Findings

Relations between energy transfer efficiency and stability analyzed

Model applicable to biological and organic energy systems

Provides a mathematically rigorous framework for open quantum systems

Abstract

The analytically tractable model employing Quantum Markovian Master Equations, derived by weak coupling procedure and satisfying complete positivity, is proposed to describe a model of molecular battery charged by a non-equilibrium excitonic reservoir. The excitons are produced by non-equilibrium processes involving, e.g. light absorption and chemical reactions. Various relations concerning the efficiency of the involved processes of energy transfer and the stability of battery are discussed. The model can be treated as an initial step in applications of mathematically sound version of Quantum Theory of Open Systems to complex processes of energy transfer in biological system and man-made devices based on organic materials.

Equations

$$\Delta g=-\sum _{j=1}^K{\nu }_j{\mu }_j,$$

$$H_0=\sum _k{E}_a(k)a_k^{†}{a}_k+\sum _{\ell }{E}_b(\ell )b_{\ell }^{†}{b}_{\ell }.$$

$$H_{intra}=\sum _{k{k}^{\prime }}{a}_k^{†}{a}_{k^{\prime }}\otimes R_{k{k}^{\prime }}^{(a)}+\sum _{\ell {\ell }^{\prime }}{b}_{\ell }^{†}{b}_{{\ell }^{\prime }}\otimes R_{\ell {\ell }^{\prime }}^{(b)}$$

$$H_{inter}=\sum_{k\ell}\bigl{(}a^{\dagger}_{k}b_{\ell}+b^{\dagger}_{\ell}a_{k}\bigr{)}\otimes R^{(ab)}_{k\ell},$$

$$\frac{d}{dt}\rho =-\frac{i}{\hslash }[H_0,\rho ]+\mathcal{L}\rho ,$$

$$\mathcal{L}={\mathcal{L}}^{intra}+{\mathcal{L}}^{inter},$$

$${\mathcal{L}}^{intra}=\sum _{\{k{k}^{\prime }\}}{\mathcal{L}}_{k{k}^{\prime }}^{(a)}+\sum _{\{\ell {\ell }^{\prime }\}}{\mathcal{L}}_{\ell {\ell }^{\prime }}^{(b)},$$

$${\mathcal{L}}^{inter}=\sum _{\{k\ell \}}{\mathcal{L}}_{k\ell }^{(ab)},$$

$$\displaystyle\mathcal{L}^{(a)}_{kk^{\prime}}\rho=\frac{1}{2}\Gamma^{(a)}_{kk^{\prime}}\Bigl{(}[a_{k}a_{k^{\prime}}^{\dagger}\,\rho,\,a_{k^{\prime}}a_{k}^{\dagger}]+[a_{k}a_{k^{\prime}}^{\dagger},\,\rho\,a_{k^{\prime}}a_{k}^{\dagger}]$$

$$\displaystyle+e^{-(E_{a}(k)-E_{a}(k^{\prime})/k_{B}T}\bigl{(}[a_{k}^{\dagger}a_{k^{\prime}}\,\rho,\,a_{k^{\prime}}^{\dagger}a_{k}]+[a_{k}^{\dagger}a_{k^{\prime}},\,\rho\,a_{k^{\prime}}^{\dagger}a_{k}]\bigr{)}\Bigr{)},$$

$$\displaystyle\mathcal{L}^{(b)}_{\ell\ell^{\prime}}\rho=\frac{1}{2}\Gamma^{(b)}_{\ell\ell^{\prime}}\Bigl{(}[b_{\ell}b_{\ell^{\prime}}^{\dagger}\,\rho,\,b_{\ell^{\prime}}b_{\ell}^{\dagger}]+[b_{\ell}b_{\ell^{\prime}}^{\dagger},\,\rho\,b_{\ell^{\prime}}b_{\ell}^{\dagger}]$$

$$\displaystyle+e^{-(E_{b}(\ell)-E_{b}(\ell^{\prime})/k_{B}T}\bigl{(}[b_{\ell}^{\dagger}b_{\ell^{\prime}}\,\rho,\,b_{\ell^{\prime}}^{\dagger}b_{\ell}]+[b_{\ell}^{\dagger}b_{\ell^{\prime}},\,\rho\,b_{\ell^{\prime}}^{\dagger}b_{\ell}]\bigr{)}\Bigr{)},$$

$$\displaystyle\mathcal{L}^{(ab)}_{k\ell}\rho=\frac{1}{2}\gamma_{k\ell}\Bigl{(}[a_{k}b_{\ell}^{\dagger}\,\rho,\,b_{\ell}a_{k}^{\dagger}]+[a_{k}b_{\ell}^{\dagger},\,\rho\,b_{\ell}a_{k}^{\dagger}]$$

$$\displaystyle+e^{-(\epsilon_{kl}-\Delta g)/k_{B}T[\epsilon_{kl}]}\bigl{(}[a_{k}^{\dagger}b_{\ell}\,\rho,\,b_{\ell}^{\dagger}a_{k}]+[a_{k}^{\dagger}b_{\ell},\,\rho\,b_{\ell}^{\dagger}a_{k}]\bigr{)}\Bigr{)}.$$

$$\bar{\rho}[\mu_{a},\mu_{b}]={Z}^{-1}\exp\Bigl{\{}-\frac{1}{k_{B}T}\Bigl{[}\sum_{k}\bigl{(}E_{a}(k)-\mu_{a}\bigr{)}a^{\dagger}_{k}a_{k}+\sum_{l}\bigl{(}E_{b}(l)-\mu_{b}\bigr{)}b^{\dagger}_{l}b_{l}\Bigr{]}\Bigr{\}}.$$

$$\Delta \mu ={\mu }_a-{\mu }_b,$$

$$\mathcal{L}^{inter}\bar{\rho}=\sum_{k,\ell}\gamma_{k\ell}\Bigl{\{}\bigl{[}e^{-X_{k\ell}}-e^{-Y_{k\ell}}\bigr{]}\bar{\rho}\,a_{k}a_{k}^{\dagger}b_{\ell}^{\dagger}b_{\ell}+\bigl{[}e^{X_{k\ell}}e^{-Y_{k\ell}}-1\bigr{]}\bar{\rho}\,a_{k}^{\dagger}a_{k}b_{\ell}b_{\ell}^{\dagger}\Bigr{\}},$$

$$X_{k\ell }=\beta ({ϵ}_{k\ell }-\Delta \mu ),Y_{k\ell }=\beta [{ϵ}_{k\ell }]({ϵ}_{k\ell }-\Delta g).$$

$$f_a(k)=\frac{1}{e^{\beta (E_a(k)-{\mu }_a)}+1},f_b(\ell )=\frac{1}{e^{\beta (E_b(\ell )-{\mu }_b)}+1}$$

$$\|\mathcal{L}^{inter}\bar{\rho}\|_{1}\leq\sum_{k,\ell}\gamma_{k\ell}|1-e^{X_{k\ell}-Y_{k\ell}}|\Bigl{\{}e^{-X_{k\ell}}[1-f_{a}(k)]f_{b}(\ell)+f_{a}(k)[1-f_{b}(\ell)]\Bigr{\}}$$

$$X_{k\ell }\simeq Y_{k\ell }$$

$$\Delta\mu\simeq\Bigl{(}1-\frac{T}{T[\epsilon_{k\ell}]}\Bigr{)}\epsilon_{k\ell}+\frac{T}{T[\epsilon_{k\ell}]}\Delta g.$$

$$\Delta\mu=\Bigl{(}1-\frac{T}{T[\bar{E}_{g}]}\Bigr{)}\bar{E}_{g}+\frac{T}{T[\bar{E}_{g}]}\Delta g.$$

$$H_m=\hslash {\omega }_0(A^{†}-{\xi }_0|1⟩⟨1|)(A-{\xi }_0|1⟩⟨1|)+E_{el}|1⟩⟨1|,$$

$$U=e^{{\xi }_0(A^{†}-A)|1⟩⟨1|},$$

$$A\mapsto B=UA{U}^{†}=A-{\xi }_0|1⟩⟨1|,H_m=\hslash {\omega }_0{B}^{†}B+E_{el}|1⟩⟨1|.$$

$$|0⟩⟨0|\mapsto |0⟩⟨0|,|1⟩⟨1|\mapsto |0⟩⟨0|,|1⟩⟨0|\mapsto E_{10}=e^{{\xi }_0(B^{†}-B)}|1⟩⟨0|$$

$$H_{int}^{(1)}$$

$$H_{int}^{(2)}$$

$$H_{int}^{(3)}$$

$$e^{\frac{i}{\hslash }{H}_{m}t}(A+A^{†})e^{-\frac{i}{\hslash }{H}_{m}t}$$

$$e^{\frac{i}{\hslash }{H}_{m}t}|1⟩⟨0|e^{-\frac{i}{\hslash }{H}_{m}t}$$

$$W(\alpha )=e^{(\alpha A^{†}-\bar{\alpha }A)}=e^{-|\alpha {|}^2/2}{e}^{\alpha A^{†}}{e}^{-\bar{\alpha }A}.$$

Full text

I Introduction

Energy transduction processes on the molecular level are responsible for functioning of living organisms and operation of various man-made devices like photovoltaic, thermoelectric or chemical cells. The ultimate microscopic theory of these phenomena should be quantum and derivable from the first principles, i.e. from the fully quantum Hamiltonian models. The natural mathematical formalism to achieve this task is provided by the Quantum Theory of Open Systems Alicki:07 ,Breuer ,Huelga , in particular by the theory of Quantum Markovian Master Equations (QMME). Those equations of motion for the reduced density matrices of the open systems are obtained by elimination of the environmental degrees of freedom. From the plethora of various approaches and mathematical tools we select the method introduced by Davies Davies combining the weak coupling (Born), Markovian and secular approximations. The obtained QMMEs satisfy the requirement of complete positivity,GKS ,Lindblad possess a clear mathematical structure and are consistent with thermodynamics SpohnLeb ,Alicki:79 ,QT-review . In particular the last property is important because the universal laws of thermodynamics provide the independent validity control of the often involved approximation schemes.

The aim of this paper is to present an analytically tractable model based on QMMEs which are derivable from the full Hamiltonian dynamics and describe a certain sequence of energy transduction processes. Namely, in the first step the electronic system (solid state or organic material) is weakly coupled to a stationary but non-equilibrium environment. The thermal and/or chemical energy supplied by the reservoir is used to produce a gas of excitons which can be treated as a “chemical fuel” used, for example, to power various types of engines (e.g. photovoltaic cells, photosynthetic reaction centers). The models of such engines based on QMMEs and employing self-oscillation mechanisms were studied for example in Solar , Universal . Here, the different process of energy transduction is discussed – charging of “molecular battery” by non-equilibrium exciton gas. The presented model of molecular battery, involving a single “reaction coordinate” and two electronic states, can be seen, for example, as a toy model of ATP/ADP complex while the whole discussed sequence of energy transfers can mimic the complicated processes in photosynthetic light harvesting systems. Finally, the process of spontaneous discharging of molecular battery at ambient temperature is studied in order to estimate stability of battery.

The main advantage of the applied version of QMMEs is their particular mathematical structure with clear thermodynamical interpretation. This allows to derive various relations characterizing the thermodynamical efficiency of energy transduction processes and the corresponding trade-offs between physical parameters of the model.

For the reader’s convenience the Appendix with a brief derivation of QMME is added. An intermediate step, called recently refined weak coupling Rivas and derived much earlier in Alicki:89 as completely positive non-Markovian Born approximation, is used. This formalism provides a mathematically consistent and apparently numerically accurate interpolation between a slow highly non-Markovian short-time dynamics and long-time Markovian one governed by QMME.

II Exciton factory

Exciton is a neutral elementary excitation of an electronic system composed of electron and hole forming a bound state. In semiconductors the binding energy is usually small with respect to ambient thermal energy $k_{B}T$ and hence one can treat exciton as an essentially free electron excited from the valence to conduction band leaving behind a hole in the valence band. The bands are assumed to be separated by the gap $E_{g}>>k_{B}T$ and the valence band is completely filled at zero temperature. Therefore each excited electron can be identified with the created exciton. For organic materials the HOMO (highest occupied molecular orbital) level replaces the valence band maximum and the LUMO (lowest unoccupied molecular orbital) level the conduction band minimum. Here, binding energy of the exciton is higher but again the picture of independent electrons with a certain self-consistent free Hamiltonian provides a reasonable approximation.

The exciton production is the first step of the discussed energy transfer process involving also charging of molecular battery and its controlled or spontaneous discharging, as presented schematically on Fig.1.

To describe creation of excitons we use a fully quantum mechanical model of a system with relevant degrees of freedom corresponding to electrons distributed in two bands $A$ (higher) and $B$ (lower), separated by the energy gap $E_{g}$, and interacting with a non-equilibrium environment. The environment consists of two components: the cold, typically phonon bath at the ambient temperature $T$ and the driving hot bath which can be in a non-equilibrium but stationary state characterized by the local temperature $T[\epsilon]$ depending on the energy scale $\epsilon$ localT and fixed chemical potentials $\mu_{j}$ corresponding to different chemical species being subjected to a chemical reaction. The example of such a bath is a stationary but generally non-thermal light source parametrized by photon populations $n(\omega)=[e^{\hbar\omega/k_{B}T[\hbar\omega]}-1]^{-1}$. A chemical bath can produce a single exciton at the expense of the free energy $\Delta g$ per single reaction step given by

[TABLE]

where $\nu_{j}$ are stoichiometric coefficients, negative for reactants and positive for reaction products (see Appendix).

II.1 Model Hamiltonians

The electrons distributed in two bands $A$ and $B$ are described by two sets of fermionic annihilation and creation operators $a_{k}$, $a^{\dagger}_{k}$ and $b_{l}$, $b^{\dagger}_{l}$, respectively. The electrons are treated as non-interacting fermions moving in a self-consistent potential with the Hamiltonian

[TABLE]

Here, $k$ and $\ell$ is a short-hand notation for a set of quantum numbers characterizing a single-electron energy eigenstate. Notice, that the spatial structure of the sample is hidden in the corresponding wave functions $\phi_{k}(\mathbf{r})$, $\phi_{\ell}(\mathbf{r})$ written in the position representation. The often used semi-classical picture of “band bending” at the interfaces between different materials is not consistent with the fully quantum description.

We can assume that the energy gap $E_{g}$ is higher than the typical Debye energy of phonons. Then, the cold phonon bath can cause only intraband transitions and the corresponding interaction Hamiltonian can be written as

[TABLE]

Here, $R^{(a)}_{kk^{\prime}},R^{(b)}_{\ell\ell^{\prime}}$ are hermitian matrices consisting of operators acting on the Fock space of phonons. The detailed form of those operators is not relevant for the further derivations.

Because the intraband transitions are usually dominated by the interaction with cold phonon bath we can assume that the driving hot bath contributes to interband transitions only. In this case the interaction Hamiltonian reads

[TABLE]

where, again, the exact form of bath operators is not relevant.

We are interested in the situation when at zero temperature the lower band $B$ is completely filled, hence the operator $a^{\dagger}_{k}b_{\ell}$ ($b^{\dagger}_{\ell}a_{k}$) creates (annihilates) an electron-hole pair (exciton) from the ground state of electronic systems.

II.2 Quantum Markovian Master Equation

The derivations of QMME for the similar models have been presented in previous papers (e.g. Solar , Universal ), therefore we reproduce here only the final results and discuss their physical interpretation. However, for reader’s convenience the basis steps of the derivation for a general QMME are presented in the Appendix together with selected properties of them.

The QMME for the reduced density matrix $\rho(t)$ acting on the fermionic Fock space of electrons reads

[TABLE]

where $H_{0}$ is the Hamiltonian given by (2) with renormalized single-electron energies including lowest order corrections due to the interaction with the baths in the self-consistent approximation. The irreversible processes are governed by the generator $\mathcal{L}$ with the following structure

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

The summations in the formulas of above are performed over the sets of indices corresponding to positive energy differences between pairs of states: $\{kk^{\prime}\}=\{kk^{\prime};E_{a}(k)-E_{a}(k^{\prime})\geq 0\},\{\ell\ell^{\prime}\}=\{\ell\ell^{\prime};E_{b}(\ell)-E_{b}(\ell^{\prime})\geq 0\},\{k\ell\}=\{k\ell;E_{a}(k)-E_{b}(\ell)>\geq 0\}$, $\epsilon_{kl}=(E_{a}(k)-E_{b}(\ell))$.

Excitons are created by the operators $a_{k}^{\dagger}b_{\ell}$ and annihilated by $b_{\ell}^{\dagger}a_{k}$. As explained in the Appendix the Boltzmann factor in (11) contains the contribution $\Delta g$ to energy difference which accounts for chemical energy provided by the hot bath.

The term (9) ((10)) describes the following intraband processes:

a) electron relaxation from the state $k$ ($l$) to the state $k^{\prime}$ ($l^{\prime}$ ) accompanied by a positive energy release to the heat bath at the ambient temperature ,

b) the inverse process of electron transfer from the lower energy state to the higher energy one with the probability suppressed by the Boltzmann factor.

The term (11) describes the following interband processes:

c) electron transfer from the higher energy state $k$ in the band $A$ to the lower energy state $l$ in the band $B$,

d) the inverse process of electron transfer with the probability controlled by suitable Boltzmann factors including energy dependent effective temperature $T[\cdot]$ and the chemical energy $\Delta g$.

In principle all relaxation parameters can be derived from the underlying Hamiltonian models with properly chosen stationary states of the baths. In practice phenomenological values obtained e.g. from the measurements of light absorption and electron thermalisation rates can be used.

II.3 Stationary state of electronic system

The analysis of the stationary state for the QMME (5) - (11) is based on the following facts:

i) the intraband thermalization processes governed by (9), (10) preserve the number of electrons in each band and for a fixed initial numbers of electrons drive the system into a product of two canonical Gibbs states at the temperature $T$,

ii) the interband transitions described by (11) preserve the total number of electrons but modify the relative occupation of both band,

iii) because the number of electrons in a macroscopic system is very large we can use the equivalence of canonical and grand canonical ensembles to derive the convenient approximate form of the stationary state for $\mathcal{L}$.

Taking i)-iv) into account one concludes that the approximate stationary state of the generator $\mathcal{L}$ has a form of the product of grand canonical ensembles at the ambient temperature $T$, but with two different electrochemical potentials $\mu_{a},\mu_{b}$

[TABLE]

The states $\bar{\rho}[\mu_{a},\mu_{b}]\equiv\bar{\rho}$ are stationary for the generator $\mathcal{L}^{intra}$ for any choice of chemical potentials. Even their restrictions to subspaces of the Fock space with fixed electron numbers are stationary. To find the best approximation to the stationary state with respect to the full dynamics given by QMME (5) we have to minimize the trace norm of $\mathcal{L}^{inter}\bar{\rho}$. This allows to obtain the optimal parameter

[TABLE]

while the chemical potential $\mu_{b}$ is fixed by the averaged electron number in the system. A series of simple computations yields the following expression for $\mathcal{L}^{inter}\bar{\rho}$

[TABLE]

where

[TABLE]

Introducing the Fermi-Dirac distributions

[TABLE]

one obtains the following estimation of the trace norm

[TABLE]

The trace norm is small if for those energy differences $\epsilon_{k\ell}$ which enter the RHS of (17) with the substantial weights, the approximate equality holds

[TABLE]

what implies

[TABLE]

Remembering that $\epsilon_{k\ell}\geq E_{g}>>k_{B}T$ and that the intraband thermalization is much faster than the interband one we expect that only the energy differences close to the gap, i.e. $\epsilon_{k\ell}\simeq E_{g}+\mathcal{O}(k_{B}T)$ matter. Hence, in the formula (19) $\epsilon_{k\ell}$ can be replaced by the effective, weakly dependent on the temperature, energy gap $\bar{E}_{g}\simeq E_{g}+\mathcal{O}(k_{B}T)$. As a consequence the final formula for the chemical potential of the non-equilibrium excitonic bath fed by the external temperature gradient and chemical reactions reads

[TABLE]

This relation appears in various contexts, for example describes the open circuit voltage ($eV_{oc}\equiv\mu_{a}-\mu_{b}$) in photovoltaic, thermoelectric and chemical cells hot ,Universal . Its thermodynamical meaning is quite clear. The chemical energy of an exciton $\Delta\mu$ can be in principle completely converted into work by exciton recombination in a certain external circuit attached to the excitonic reservoir. The formula shows that if the exciton is created by thermal energy the heat - chemical potential conversion efficiency is limited by the Carnot bound while in the case of isothermal conditions ($T=T[\bar{E}_{g}]$) the whole energy released in the chemical reaction can be, in principle, transferred into work ($\Delta\mu=\Delta g$).

III Molecular battery

We consider a model of a molecule with two relevant electronic states $|0\rangle,|1\rangle$ , a ground state and the first excited state, respectively, separated by the energy gap $E_{el}$. The other relevant degree of freedom can be called reaction coordinate and is modeled by the harmonic oscillator with a single, independent of the electronic state, frequency $\omega_{0}$ and the annihilation and creation operators $A,A^{\dagger}$ ( a similar model of a “quantum switch” has been studied in Alicki:2013 ). It is not difficult to generalize this model to many electronic states with different frequencies for the reaction coordinate at the price of much more complicated formulas.

The Hamiltonian of the model which is illustrated on Fig.1 reads

[TABLE]

where $\xi_{0}$ describes the dimensionless displacement of the minimum of the potential curve for excited electronic state with respect to the ground state and is equal to the square root of the Huang-Rhys factor.

It is convenient to introduce the polaron transformation in the form of the following unitary operator

[TABLE]

which yields the new set of operators for the harmonic oscillator part and a new form of the Hamiltonian

[TABLE]

For electronic states

[TABLE]

We treat the molecule as an open system subject to the weak interaction with the stationary environment including also the other molecular degrees of freedom. The lowest order interaction Hamiltonians can be chosen as

[TABLE]

where $F_{j},j=1,2,2$ are statistically independent environment observables.

The interaction Hamiltonian ${H}^{(1)}_{int}$ is used to account for the thermalization process of the reaction coordinate to the ambient temperature. The second Hamiltonian ${H}^{(2)}_{int}$ commutes with the molecular Hamiltonian, hence describes various pure decoherence processes of electronic states. The last Hamiltonian ${H}^{(3)}_{int}$ governs transitions between electronic states which in our model are dominated by the interaction with the non-equilibrium excitonic reservoir discussed in the previous Section.

In order to construct the proper QMME using again the weak coupling approach we have to calculate the Heisenberg dynamics of relevant molecular operators entering the interaction Hamiltonians (25) – (27)

[TABLE]

where $W(\alpha)$ denotes the Weyl operator

[TABLE]

The transition operators entering the QMME for the reduced density matrix of the molecule are obtained from the decomposition of the RHS of eqs.(28), (29) into terms oscillating with relevant Bohr frequencies.

III.1 Thermalization of reaction coordinate

To simplify the analysis of the molecular dynamics we consider first only the thermalization to ambient temperature cause by the coupling of the type (25), (26). The QMME is constructed using the decomposition into oscillating terms (28) and coupling spectra for the thermal bath at the temperature $T$

[TABLE]

where $\langle\cdot\rangle_{T}$ denotes the average with respect to Gibbs state of the bath. The following Kubo-Martin-Schwinger (KMS) relation holds

[TABLE]

Standard derivation procedure leads to the following QMME equation for the reduced density matrix of the molecule

[TABLE]

where the dissipation rate $\gamma=G_{1}(\omega_{0})$, and the pure decoherence rate

[TABLE]

The master equation (33) is a well known exactly solvable one which describes thermalization of the quantum harmonic oscillator (with equilibrium position conditioned on the electronic state) accompanied by the independent decoherence process for a 2-level electron system. Any initial state evolves into the mixture of the conditioned Gibbs states

[TABLE]

with weights given by the initial occupation probabilities of electronic states $|0\rangle,|1\rangle$.

III.2 Coupling to exciton bath

The interaction between the exciton bath and the molecule which is responsible for the process of “charging”’ of molecular battery can be derived from the fundamental (screened) Coulomb interaction between electrons in the whole system. Within the (nonrelativistic) formalism of the quantum field theory we describe the total electronic systems in terms of position-dependent (we omit spin variables for simplicity) fermionic quantum fields $\Psi(\mathbf{r})$, $\Psi^{\dagger}(\mathbf{r})$. Then the Coulomb interaction is described by the following Hamiltonian

[TABLE]

with the effective electron-electron interaction potential $V(\mathbf{r},\mathbf{r}^{\prime})$.

The electronic quantum field can be decomposed using the single-electron wave functions : $\phi^{a}_{k}(\mathbf{r})$, $\phi^{b}_{\ell}(\mathbf{r})$ corresponding to bands A and B, and $\chi_{0}(\mathbf{r})$ , $\chi_{1}(\mathbf{r})$ to electronic states of the molecule. Then

[TABLE]

where $c_{0}$ and $c_{1}$ are new fermionic annihilation and creation operators for an electron localized in the molecule.

Inserting (37) into (36) one obtains a sum of different terms each containing two annihilation and two creation operators. Most of the terms describe the electron-electron interaction in the excitonic bath which were already taken into account in the self-consistent single-electron Hamiltonian (2), the other describe some residual potential between the bath and the molecule. The only terms which can account for the energy transfer between the bath and the molecule are represented by the following effective interaction Hamiltonian

[TABLE]

where the coupling coefficients $g_{kl}$ might be in principle computed from the fundamental interaction Hamiltonian. The energy transfer is effective under the condition $E_{g}\simeq E_{el}$ which implies that only the terms containing the operators $a_{k}^{\dagger}b_{\ell}c_{1}^{\dagger}c_{0}$ and $b_{\ell}^{\dagger}a_{k}c_{0}^{\dagger}c_{1}$ describe real resonant exchange process depicted on Fig.1. Notice that because we consider only a single electron occupying the molecule the operator $c_{1}^{\dagger}c_{0}$ is equivalent to $|1\rangle\langle 0|$ and hence the Hamiltonians (38) is a special case of (27).

III.3 Charging the battery

The interaction between the exciton bath and the molecular battery yields irreversible processes in the battery which can be described by the additional term added to the QMME (33), derived again using weak coupling formalism. The first step is the computation of the coupling spectrum given by the formula (31) with $F_{3}=\sum_{k\ell}g_{kl}(a_{k}^{\dagger}b_{\ell}+b_{\ell}^{\dagger}a_{k})$ and for the thermal equilibrium state replaced by the stationary non-equilibrium one (12). The standard computation yields the following results, for $\omega\geq 0$,

[TABLE]

and for negative frequencies

[TABLE]

with the Fermi-Dirac distributions (16).

Inserting (16) into (39) and (40) yields the KMS relation

[TABLE]

which describes the ratio of transition probabilities per unit time of excitation and relaxation processes. This relation concerns a pair of energy levels of any system weakly coupled to the excitonic bath.

In order to discuss the shape of the coupling spectrum $G_{3}(\omega)$ we need some estimations for the magnitudes of relevant parameters of the model. The typical parameters satisfy the conditions

[TABLE]

The existence of the gap and energy conservation hidden in Dirac deltas in (39), (40) imply that $G_{3}(\omega)$ can be strictly positive only for Bohr frequencies $\omega>E_{g}$. Notice that because the states of the molecule which are coupled through the interaction with the exciton bath must involve different electronic states the set of relevant (and positive) Bohr frequencies is the following

[TABLE]

where $\omega_{el}=E_{el}/\hbar$.

Any initial state of the battery weakly coupled to the excitonic baths and thermalized by the mechanism described in the previous section tends to a final stationary state which according to (35), (41) and (43) has form

[TABLE]

and generally differs from the equilibrium state of battery given by (44) with $\Delta\mu=0$.

Remark The computation of the coupling spectrum $G_{3}(\omega)$ was done using as dynamics of the exciton bath the unitary evolution governed by the Hamiltonian (2) and not the full irreversible dynamics including the interaction with hot and cold baths. However, the dissipative effects caused by the later produce only a certain indeterminacy in the energy conservation which is expressed here by the Dirac deltas in formulas (39), (40). The coupling spectrum $G_{3}(\omega)$ is essentially different from zero for the value $\omega>E_{g}/\hbar$ which are high in comparison to the relaxation rates induced by cold and hot bath (weak coupling assumption). Therefore, the dissipative effects lead only to replacement of Dirac deltas by their regularized models with the widths of the order of ($\hbar\times$ relaxation rates) what does not change essentially the coupling spectrum.

III.4 Work extraction from the battery

The battery in the stationary state (44) is a source of useful work. The maximal extractable work is given by ergotropy, the name introduced in ergo and defined for an arbitrary state $\rho$ of a quantum system equipped with a Hamiltonian $H$ as

[TABLE]

where $\sigma_{\rho}$ is a unique passive state Pusz associated with $\rho$. This state is diagonal in the Hamiltonian basis, and possesses the same eigenvalues as $\rho$, ordered in the non-increasing order with energy.

For a system composed of $n$ identical, and non-interacting subsystems one can define the asymptotic maximal work, per subsystem, extractable by application of arbitrary many-particle external Hamiltonian perturbations

[TABLE]

where $H_{n}$ is a sum of identical single-subsystem Hamiltonians.

It was shown in AlFan that

[TABLE]

where $\rho_{\bar{\beta}}=Z^{-1}e^{-\bar{\beta}H}$ is a unique Gibbs state with the same von Neumann entropy as the state $\rho$. Generally, $\bar{W}_{max}\geq W_{max}$, is much easier to compute and the formula (47) can be interpreted in terms of non-equilibrium free energy.

In the case of molecular battery the basic ingredients of the relation (47) can be computed analytically. Namely, the von Neumann entropies of the stationary state (44) and the Gibbs state of the battery at the inverse temperature $\bar{\beta}$ read

[TABLE]

Now from the equation $S(\rho_{st})=S(\rho^{(m)}_{\bar{\beta}})$ one can numerically compute $\bar{\beta}$ and then insert into the final formula for the extractable work

[TABLE]

For low temperature regime one can obtain from (47)

[TABLE]

where $\Theta(x)$ is the Heaviside function.

III.5 Stability of molecular battery

The molecular battery in the excited electronic state $\rho^{(1)}$ detached from the excitonic bath still remains in contact with the heat bath at the ambient temperature $T$ corresponding to various environmental degrees of freedom. Such interaction leads to spontaneous discharging of the battery as well as to much less probable (suppressed by the Boltzmann factor $e^{-E_{el}/k_{B}T}$) spontaneous recharging. Of course, the lower is the discharging rate the more useful is the battery.

The dynamics of the molecule in the environment at ambient temperature is given by the Master equation of the form

[TABLE]

Here, $\mathcal{L}_{0}$ describes thermalization process of the reaction coordinate and pure decoherence of the electronic degree of freedom given by the non-Hamiltonian part of the RHS of eq.(33). The generator $\mathcal{L}_{1}$ accounts for transitions between different electronic states caused by thermal fluctuations and derived using the interaction Hamiltonian of the type (27). This additive structure can be justified assuming that the leading interaction of the molecule with an environment is mediated by the dipole moment of the molecule which is a sum of contributions from nuclei and electrons. The first one enters the interaction Hamiltonian of the form (25) leading to reaction coordinate thermalization and the second one leads to (27) and accounts for spontaneous discharging/charging processes.

As we shall see the structure of $\mathcal{L}_{1}$ is quite involved and therefore we cannot solve analytically the full QMME (53). However, the stability of battery can be characterized by the initial rate of the probability flow out of the “excited Gibbs state” $\rho^{(1)}$ (see (35))

[TABLE]

Indeed, the state $\rho^{(1)}$ commutes with $H_{m}$ and is stationary with respect to $\mathcal{L}_{0}$ hence only $\mathcal{L}_{1}$ contributes to the initial dynamics of $\rho^{(1)}$.

The generator $\mathcal{L}_{1}$ can be obtained again using the weak coupling approach and possesses the following form (we use a different ordering of terms than in the formula (33))

[TABLE]

The transition operators ${V}_{m}$ can be computed by decomposition of of (29) into terms oscillating with the frequencies $-(E_{el}/\hbar+m\omega_{0})$ and using (30) what gives

[TABLE]

Combining the eqs. (35), (54), (III.5) and (56) one obtains

[TABLE]

where

[TABLE]

As typically $\hbar\omega_{el}>>\hbar\omega_{0}>>k_{B}T$ it is instructive to compute $\Gamma_{\downarrow}$ in the limit case $T=0$. Then

[TABLE]

For $S>>1$ the Poisson distribution on the RHS of (59) is concentrated around the mean value $\bar{m}=S$ what, for smooth enough coupling spectrum $G(\omega)$ leads to the final formula

[TABLE]

Due to the KMS condition (32), for $T\to 0$, the coupling spectrum $G(\omega)$ vanishes for $\omega<0$. Therefore, the charged state of the molecular battery becomes stable for $S>\frac{\omega_{el}}{\omega_{0}}$. One can expect that in the case of finite temperatures the rapid suppression of spontaneous discharging is replaced by the Boltzmann factor $\sim e^{-(S\hbar\omega_{0}-E_{el})}$ (compare with the computation in Alicki:2013 performed for the model with $E_{el}=0$).

One should notice that the strategy of increasing Huang-Rhys factor to improve the performance of battery possesses its limitations. Namely, the Huang-Rhys factor appears also in the decoherence rate (34) in the QMME (33). To extract work from the battery one has to apply external time-dependent perturbation of the Hamiltonian which does not commute with the Hamiltonian of the battery. Under such perturbation, the previous pure decoherence process causes now energy dissipation reducing the efficiency of work extraction.

IV Concluding remarks

Although quantum theory of open systems is an active research topic in the field of chemical, biological, and condensed-phase physics, in most cases little attention is payed to the mathematical consistency of the applied formalism. The typical examples are the various versions of Redfield equation or approximations involving path integrals which neither satisfy the requirements of transforming density matrices into density matrices, nor the general restrictions imposed by thermodynamics. The presented method based on systematic and careful weak coupling approximation allows to construct the evolution equations of the GKLS type satisfying the mentioned above natural requirements. Moreover, those equations possess a highly symmetric mathematical structure which allows to obtain analytical results for quite complex systems as illustrated by the model of molecular battery charged by excitons. Those analytical results can be extrapolated to more complicated situations which require numerical methods.

Quite often it is argued that the QMME yield much less accurate solutions than, for instance, Redfield equation, in particular for the short-time evolution Redfield . In the Appendix the basis approximations, used in the derivation of QMME, are discussed. It is shown, indeed, that a clear separation of various time scales is necessary to justify the Markovian approximation. Generically, solutions of QMME describe accurately the evolution of open systems for long enough times what is often sufficient to discuss the approach to asymptotic stationary state. In the case of various not well-separated time scales a careful choice of the reference physical Hamiltonian used in the first step of the derivation – transition to the interaction picture – can improve the accuracy of QMME (see e.g. Vsystem ). Finally, the “refined weak coupling” formalism Alicki:89 , Rivas allows to construct completely positive dynamical maps interpolating with a good accuracy between the strongly non-Markovian short-time dynamics and asymptotic long-time evolution governed by QMME.

Acknowledgement

The work is a part of the ICTQTIRAP project of FNP, financed by structural funds of EU.

Appendix. Quantum Markovian Master Equations

The brief derivation of QMME based on cumulant expansion and refined weak coupling approach is presented and basic properties of QMME are discussed, in particular for heat and chemical baths Alicki:89 Rivas .

IV.1 From refined weak coupling to MME

A quantum system with a ”bare” system Hamiltonians $H^{0}$ interacts with a reservoir equipped with the Hamiltonian $H_{R}$, by means of the interaction Hamiltonian $H_{int}=\lambda S\otimes R$ with $S=S^{\dagger}$, $R=R^{\dagger}$, with $\lambda$ being a “small” coupling strength ( generalization to more complicated $H_{int}$ is straightforward). The initial state of the baths satisfies

[TABLE]

The reduced, system-only dynamics in the interaction picture is defined as a partial trace

[TABLE]

where the unitary propagator in the interaction picture is given by the ordered exponential

[TABLE]

with

[TABLE]

Notice, that $S(t)$ is defined with respect to the renormalized, physical Hamiltonian $H$ which can be expressed as

[TABLE]

The renormalizing counterterms containing powers of $\lambda$ are often called * Lamb-shift* corrections and are due to the interaction with the bath.

The cumulant expansion for the reduced dynamics can be defined as

[TABLE]

Comparing terms with the same power of $\lambda$ in (66) and in the Dyson expansion of (62) one finds that $\mathcal{K}^{(1)}=0$ and the Born approximation to the exact dynamics in the interaction picture has form:

[TABLE]

with

[TABLE]

Here $F(s)=\mathrm{Tr}(\rho_{R}R(s)R)$ and the Hamiltonian term

[TABLE]

is often called Lamb-shift Hamiltonian. In the Schroedinger picture and for the Markovian limit $H_{L}(t)$ can be omitted because its contribution is compensated by the counterterm in the renormalized Hamiltonian (65).

Notice, that the mathematical structure of the superoperator $\mathcal{K}^{(2)}(t)$ allows to recast it into the standard Gorini-Kossakowski-Lindblad-Sudarshan (GKLS) form and hence its exponential (67) defines a family of completely positive trace preserving dynamical maps. This mathematically consistent approximation was proposed first in Alicki:89 and rediscovered recently in Rivas where numerical tests shown its usefulness.

We show now that for long enough time the following Markovian approximation is valid

[TABLE]

where $\mathcal{L}$ is a GKLS generator. To find its form we first decompose $S(t)$ into its Fourier components

[TABLE]

where the set $\{\omega\}$ contains Bohr frequencies of the physical Hamiltonian with the spectral decomposition

[TABLE]

Then we can rewrite the expression (68) (with neglected Hamiltonian term) as

[TABLE]

and use two crucial approximations:

[TABLE]

which are valid for $t\gg\max\{1/(\omega-\omega^{\prime})\}$. Applying these two approximation we obtain $K(t)\rho_{S}=(t/\hbar^{2})\sum_{\omega}S_{\omega}\rho_{S}S_{\omega}^{{\dagger}}{G}(\omega)+(\mathrm{similar}$ $\mathrm{terms})$, and hence it follows from Eq. (70) that $\mathcal{L}$ is a special case of the GKLS generator derived rigorously for the first time by Davies Davies . Returning to the Schroedinger picture one obtains the following QMME:

[TABLE]

Remarks:

(i) The absence of off-resonant terms in Eq. (75) (so-called secular approximation) is the crucial property of the Davies generator which can be interpreted as averaging in time of fast oscillating terms. It implies also the commutation of $\mathcal{L}$ with the Hamiltonian part $[H,\cdot]$. In the case of almost degenerated Bohr frequencies one can use in the interaction picture a modified Hamiltonian with certain strictly degenerated Bohr frequencies to improve the approximation Vsystem .

(ii) The positivity $G(\omega)\geq 0$ follows from the Bochner’s theorem and is a necessary condition for the complete positivity of the QMME.

(iii) The presented derivation shows implicitly that the notion of bath’s correlation time, often used in the literature, is generally not well-defined – Markovian behavior involves a rather complicated cooperation between system and bath dynamics.

It follows from (iii) that, contrary to what is often done in phenomenological treatments, one cannot combine arbitrary $H$’s with a given GKLS generator. This is particularly important in the context of thermodynamics of controlled quantum open system, where it is common to fix Markovian, irreversible part of dynamics and apply arbitrary control Hamiltonians. Uncontrolled derivations of QMMEs can easily lead to violation of the laws of thermodynamics Benatti Levy ,Hofer .

Remark All the results including the expressions for (68) and (75) can be easily extended to more complicated interaction Hamiltonians

[TABLE]

In particular the most general form of QMME obtained by the weak coupling procedure reads

[TABLE]

where for all $\{\omega\}$, the matrix $[G_{\beta\alpha}(\omega)]$ is positively defined.

)

IV.2 Thermal and chemical reservoirs

If the reservoir is a quantum system at thermal equilibrium state the additional Kubo-Martin-Schwinger (KMS) condition holds

[TABLE]

where $T$ is the bath’s temperature. As a consequence of (78) the Gibbs state

[TABLE]

is a stationary solution of (75). Under mild conditions (e.g : ”the only system operators commuting with $H$ and $S$ are scalars”) the Gibbs state is a unique stationary state and any initial state relaxes towards equilibrium (”0-th law of thermodynamics”). A convenient parametrization of the corresponding thermal generator reads

[TABLE]

where

[TABLE]

A chemical bath Universal is a collection of $K$ separated baths each consisting of many molecules of a given type. The whole bath is described by the density matrix corresponding to the grand canonical ensemble at the temperature $T_{1}$ and the chemical potentials $\mu_{j}$

[TABLE]

Here, $N_{j}$ is the operator describing the number of $j$-type molecules in the bath and $H_{R}$ is the total Hamiltonian of the chemical bath satisfying

[TABLE]

The condition (83) means that the molecules in the chemical bath do not undergo chemical reactions unless they became coupled by the quantum system. This coupling can be described by the rotating-wave approximation version of the interaction Hamiltonian (76)

[TABLE]

Here, $S^{+}$ increases the energy of the system and the operator $R^{-}$ describes the possible transitions from higher to lower energy levels of the bath including chemical reactions. One can think about $R^{-}$ as a sum of products of operators describing three possible processes: a) transitions preserving all numbers of molecules, b) “annihilation operators” corresponding to reactants, c) “creation operators” corresponding to reaction products. Hence, $R^{-}$ satisfies the following relation

[TABLE]

where $\nu_{j}$ are stoichiometric coefficients, negative for reactants and positive for reaction products. The expression

[TABLE]

can be interpreted as the Gibbs free energy transferred during a single step of chemical reaction to the system. Obviously, the operator $S^{-}\otimes R^{+}$ describes the reverse processes consuming energy from the system. If the operators $R^{\pm}$ commute with all number operators $N_{j}$, then $\Delta g=0$, no chemical reactions occur, and the bath acts as a purely thermal bath at the temperature $T_{1}$.

Applying again the standard derivation of the QMME one obtains the equation (80) with the modified relaxation coefficient

[TABLE]

and $e^{-\beta\omega}\gamma(\omega)$ replaced by the excitation coefficient

[TABLE]

The fictitious dynamics governed by the extended Hamiltonian $\tilde{H}_{R}=H_{R}-\sum_{j=1}^{K}\mu_{j}N_{j}$ and applied to the operators $R^{\pm}$ yields

[TABLE]

Using (89), and again the KMS condition for the grand canonical ensemble $\tilde{\rho}_{R}$ treated as an “extended” Gibbs state and the “extended” operators $\tilde{R}^{\pm}(t)$ one obtains the following relation

[TABLE]

which generalizes the detailed balance condition (78) to chemical baths. One can show that such structure of QMME implies the validity of the second law of thermodynamics SpohnLeb - QT-review and Onsager relations in the close to equilibrium regimeAlicki:76 .