Magnetically-Driven Quantum Heat Engines: The Quasi-Static Limit of Their Efficiency

Abstract

The concept of a quantum heat engine (QHEN) has been discussed in the literature, not only due to its intrinsic scientific interest, but also as an alternative to efficiently recover, on a nanoscale device, thermal energy in the form of useful work. The quantum character of a QHEN relies, for instance, on the fact that any of its intermediate states is determined by a density matrix operator. In particular, this matrix can represent a mixed state. For a classical heat engine, a theoretical upper bound for its efficiency is obtained by analyzing its quasi-static operation along a cycle drawn by a sequence of quasi-equilibrium states. A similar analysis can be carried out for a quantum engine, where quasi-static processes are driven by the evolution of ensemble-averaged observables, via variation of the corresponding operators or of the density matrix itself on a tunable physical parameter. We recently proposed two new conceptual designs for a magnetically-driven quantum engine, where the tunable parameter is the intensity of an external magnetic field. Along this article, we shall present the general quantum thermodynamics formalism developed in order to analyze this type of QHEN, and moreover, we shall apply it to describe the theoretical efficiency of two different practical implementations of this concept: an array of semiconductor quantum dots and an ensemble of graphene flakes submitted to mechanical tension.

1. Introduction

The interesting subject of quantum thermodynamics [1], which is largely based on the theory of quantum open systems [1,2], provides a theoretical framework to study quantum heat engines (QHEN). The analysis of a QHEN operating far from thermal equilibrium is an interesting and, to a great extent, an open problem. The dynamics of such systems, whose states are determined by a reduced density matrix operator out of equilibrium, is described by different approximations to the master equation [1,2,3].

As an alternative to recover thermal energy in the form of useful work on a nanoscale device, QHENs have been proposed in the literature [4,5,6,7,8,9,10,11,12,13,14,15,16,17,18]. Within the general definition of a QHEN, whose working fluid is of a quantum mechanical nature, it is important to distinguish those that have a reciprocating operation [3,19]. It has been proven [19] that under not too restrictive conditions, a reciprocating QHEN converges to a stationary limit cycle. Such a limit cycle, in direct analogy to a classical engine, can be pictured as a sequence of equilibrium states, with a trajectory determined by a minimal set of variables [1,19,20]. Interesting examples of this concept are constituted by photosynthesis in plants [21], as well as human-designed photocells [11,12], where the working substance are thermalized photons. Moreover, it has been recently proposed that if the reservoirs are also of a quantum mechanical nature, these could be prepared into quantum coherent states [11,12] or into squeezed thermal states [11,12,22]. These examples by no means constitute the only possible configurations, since a number of different designs based on alternative principles have been proposed in the literature, such as entangled states in a qubit [23] and quantum mechanical versions of the Diesel [15] and the Otto cycle [16,22,24,25]. Conceptually, a statistical ensemble of confined single-particle systems can undergo a cycle of reversible transformations driven by a generalized external field. The driving field can be a mechanical force [4,5,6,7,9] that, by tuning a confinement length-scale, modifies the inter-level spacing of the single-particle spectrum, thus inducing a sequence of transitions on the statistical population of the single-particle states, in close analogy to a classical gas confined by a piston. We have discussed generalizations of this idea in the context of relativistic Dirac particles [26]. More recently, we proposed [27] a magnetically-driven QHEN, based on the combined effects of a parabolic confining potential, representing a semiconductor quantum dot, and an external magnetic field. In the single-particle picture, this configuration possesses an exact solution in terms of effective Landau levels, which constitute a discrete spectrum, and the effective confinement length is determined by the Landau radius. Here, the inter-level distance can be modulated by tuning the external magnetic field [27], with the effect of modulating the Landau radius. An interesting variation over this idea is the magneto-strain driven QHEN that, as we recently proposed, could be implemented on an ensemble of graphene flakes [28]. In this article, we will start with a brief introduction to the theoretical framework of quantum thermodynamics. We will present a microscopic formulation of the first law of thermodynamics, both from general out-of-equilibrium considerations, as well as in the more restricted quasi-static conditions. Finally, we will present two explicit realistic examples of the application of this formalism, in the analysis of magnetically-driven quantum engines [27,28].

2. General Theory

Along this section, we will present a brief introduction to the main definitions and physical quantities involved in the theory of quantum thermodynamics. These definitions will allow us to formulate the basic theory and assumptions required for the analysis of QHENs. First, we shall introduce the reduced density matrix operator characterizing the non-equilibrium states of an open quantum system, as well as the Markovian approximation leading to the master equation that governs its time evolution. From this description, we shall derive the non-equilibrium, microscopic version of the first law of thermodynamics for open quantum systems. At last, we shall focus on the quasi-static approximation to this law, which constitutes the theoretical basis for the two explicit examples of QHENs to be discussed in Section 3 and Section 4.

In general terms, a QHEN can be pictured as a quantum open system in contact with one or more thermal reservoirs, with which the system interacts and exchanges energy and information.

As depicted schematically in Figure 1, if the open system (S) and the reservoir (B) were isolated entities, they would be described by Hamiltonian operators ${\widehat{H}}_S$ and ${\widehat{H}}_B$, acting over states in Hilbert spaces ${\mathcal{H}}_S$ and ${\mathcal{H}}_B$, respectively. Due to the interaction ${\widehat{H}}_I$ that mixes the degrees of freedom of both subsystems S and B, the composed Hilbert space for the combination is ${\mathcal{H}}_S\otimes {\mathcal{H}}_B$, with the Hamiltonian:

$$\begin{array}{c}\hfill \widehat{H}={\widehat{H}}_S\otimes 1_B+1_S\otimes {\widehat{H}}_B+{\widehat{H}}_I.\end{array}$$

(1)

Let us define by $\widehat{\rho }$ the density matrix operator describing the mixed state of the combined system S and B. The reduced density matrix operator that characterizes the state of S, denoted by ${\widehat{\rho }}_S$, is obtained by tracing over the reservoir degrees of freedom,

$$\begin{array}{c}\hfill {\widehat{\rho }}_S={\mathrm{Tr}}_B\widehat{\rho }.\end{array}$$

(2)

The equation of motion for the reduced density matrix Equation (2) is then given by:

$$\begin{array}{c}\hfill \frac{d{\widehat{\rho }}_S\left(t\right)}{dt}=-\frac{i}{\hslash }{\mathrm{Tr}}_B\left[\widehat{H},\widehat{\rho }\left(t\right)\right].\end{array}$$

(3)

Under general physical assumptions, mainly concerning the short correlation times for observables in the reservoir B (to be discussed later in more detail), Equation (3) can be approximated by a Markovian master equation [1,2] of the Lindblad form. Mathematically, this involves the existence of a dynamical map from the space of density matrices of the system S onto itself. This map combines unitary evolution with a dissipative term described by a super-operator ${\mathcal{L}}_D$ and can be seen as a quantum mechanical version of the Liouville equation [1,2]:

$$\begin{array}{c}\hfill \frac{d{\widehat{\rho }}_S}{dt}=-\frac{i}{\hslash }[{\widehat{H}}_S^{{}^{\prime }},{\widehat{\rho }}_S]+{\mathcal{L}}_D{\widehat{\rho }}_S.\end{array}$$

(4)

Here, the super-operator ${\mathcal{L}}_D$ represents the dissipative contribution arising from the interaction between the open system S and the reservoir B. Notice that in general, ${\widehat{H}}_S^{{}^{\prime }}$ may include a renormalization of the Hamiltonian ${\widehat{H}}_S$ due to interaction between the system and the reservoir [2]. We shall come back to this point later in more detail.

The general form [1] of the super-operator ${\mathcal{L}}_D$ representing dissipation, within the assumptions leading to the Markovian master Equation (4), can be expressed by a set of Lindblad operators ${\widehat{V}}_j$, such that, for any observable $\widehat{O}$:

$$\begin{array}{c}\hfill {\mathcal{L}}_D\widehat{O}=\sum _j{\gamma }_j\left({\widehat{V}}_j\widehat{O}{\widehat{V}}_j^{†}-\frac{1}{2}\left\{{\widehat{V}}_j^{†}{\widehat{V}}_j,\widehat{O}\right\}\right).\end{array}$$

(5)

Let us now discuss the implications of the previous quantum-mechanical description of the open system dynamics in the context of non-equilibrium thermodynamics. The instantaneous ensemble-average energy for the system out-of-equilibrium is given by:

$$\begin{array}{c}\hfill E\left(t\right)={\mathrm{Tr}}_S\left({\widehat{\rho }}_S\left(t\right){\widehat{H}}_S\right),\end{array}$$

(6)

where the trace in this case is taken over the degrees of freedom of the open system S. The energy balance arising from this equation is:

$$\begin{array}{c}\hfill \frac{dE}{dt}={\mathrm{Tr}}_S\left(\frac{d{\widehat{\rho }}_S}{dt}{\widehat{H}}_S\right)+\langle \frac{\partial {\widehat{H}}_S}{\partial t}\rangle =\dot{Q}+\dot{W}.\end{array}$$

(7)

This is a non-equilibrium, quantum mechanical version of the first law of thermodynamics, where $\dot{Q}$ represents the heat dissipation rate exchanged between the system and the environment, whereas $\dot{W}$ is the power. These terms are defined by the expressions:

$$\begin{array}{ccc}\hfill \dot{Q}& =& {\mathrm{Tr}}_S\left(\widehat{H}{\mathcal{L}}_D{\rho }_S\right)={\mathrm{Tr}}_S\left({\mathcal{L}}_D{\widehat{H}}_S{\rho }_S\right)=\langle {\mathcal{L}}_D{\widehat{H}}_S\rangle ,\hfill \\ \hfill \dot{W}& =& {\mathrm{Tr}}_S\left({\widehat{\rho }}_S\frac{\partial {\widehat{H}}_S}{\partial t}\right)=\langle \frac{\partial {\widehat{H}}_S}{\partial t}\rangle .\hfill \end{array}$$

(8)

We shall now review the general assumptions required to obtain the Lindblad form of the Markovian master equation presented in Equation (4), starting from an explicit microscopic derivation. In our present discussion and along the examples that will be developed along this work, we shall assume that the reservoir B is a truly macroscopic one, as compared to the small size of the open system S. Under this assumption, it is safe to say that the relaxation time scale for the dynamics of the reservoir ${\tau }_B$ is much shorter that the corresponding one for the system S, i.e., ${\tau }_B\ll {\tau }_S$. Moreover, when the reservoir is considered as a macroscopic object, the eigenstates of ${\widehat{H}}_B$ form a continuum. Under these two physical assumptions, in the weak coupling limit between system and reservoir, the Born–Markov approximation [2] is valid,

$$\begin{array}{c}\hfill \widehat{\rho }\left(t\right)\sim {\widehat{\rho }}_S\left(t\right)\otimes {\widehat{\rho }}_B.\end{array}$$

(9)

The most general form of interaction Hamiltonian can be expressed as:

$$\begin{array}{c}\hfill {\widehat{H}}_I=\sum _{\alpha }{\widehat{A}}_{\alpha }\otimes {\widehat{B}}_{\alpha },\end{array}$$

(10)

with ${\widehat{A}}_{\alpha }$ and ${\widehat{B}}_{\alpha }$ two sets of Hermitian operators acting on the Hilbert spaces ${\mathcal{H}}_S$ and ${\mathcal{H}}_B$, respectively.

Let us assume the spectral decomposition for the open system Hamiltonian ${\widehat{H}}_S={\sum }_{ϵ}{\sum }_{i=1}^{g_{ϵ}}ϵ|{ϵ}^i\rangle \langle {ϵ}^i|$ and define the family of operators in frequency space [2]:

$$\begin{array}{ccc}\hfill {\widehat{A}}_{\alpha }\left(\omega \right)& =& \sum _{{ϵ}^{\prime }-ϵ=\omega }\widehat{P}\left(ϵ\right){\widehat{A}}_{\alpha }\widehat{P}\left({ϵ}^{\prime }\right),\hfill \\ \hfill {\widehat{B}}_{\alpha }\left(\omega \right)& =& \sum _{{ϵ}^{\prime }-ϵ=\omega }\widehat{P}\left(ϵ\right){\widehat{B}}_{\alpha }\widehat{P}\left({ϵ}^{\prime }\right),\hfill \end{array}$$

(11)

for $\widehat{P}\left(ϵ\right)={\sum }_{i=1}^{g_{ϵ}}|{ϵ}^i\rangle \langle {ϵ}^i|$ projection operators onto the subspace of eigenvectors of ${\widehat{H}}_S$ with eigenvalue ϵ and possible degeneracy $g_{ϵ}\ge 1$.

Under the assumption ${\tau }_B\ll {\tau }_S$ stated above, it is possible to separate two different time scales for the system and reservoir dynamics. Hence, a further approximation involves the averaging of the rapidly oscillating terms in the dynamics of the reservoir, also known as rotating wave approximation, where the environment operators average out to zero, $\langle {\widehat{B}}_{\alpha }\left(t\right)\rangle =\mathrm{Tr}\left({\widehat{B}}_{\alpha }\left(t\right){\widehat{\rho }}_B\right)=0$. Under these approximations, the dynamics [2] of the reduced density matrix operator, in the Schrödinger picture, is expressed by the Born–Markov master equation:

$$\begin{array}{c}\hfill \frac{d{\widehat{\rho }}_S}{dt}=-\frac{i}{\hslash }\left[{\widehat{H}}_S+{\widehat{H}}_{LS},{\widehat{\rho }}_S\left(t\right)\right]+{\mathcal{L}}_D{\widehat{\rho }}_S\left(t\right).\end{array}$$

(12)

Mathematically, this implies the existence of a one-parameter dynamical map, from the space of density matrices $V\left(t\right):\mathcal{S}\left({\mathcal{H}}_S\right)\to \mathcal{S}\left({\mathcal{H}}_S\right)$ of the system S onto itself, with $\widehat{\rho }\left(0\right)\to \widehat{\rho }\left(t\right)=V\left(t\right){\widehat{\rho }}_S\left(0\right)$ for $t>0$. This map possesses a semigroup property $V\left(t_1\right)V\left(t_2\right)=V(t_1+t_2)$. Its infinitesimal generator, as seen in Equation (12), combines unitary evolution with a dissipative term described by a super-operator ${\mathcal{L}}_D$ [1,2]. The terminology super-operator refers to the fact that the dynamical map $V\left(t\right)$ acts upon operators (reduced density matrix) instead of vectors in Hilbert space. In Equation (12), we have defined the Lamb-shift operator [2]:

$$\begin{array}{c}\hfill {\widehat{H}}_{LS}=\sum _{\omega }\sum _{\alpha \beta }{\sigma }_{\alpha \beta }\left(\omega \right){\widehat{A}}_{\alpha }^{†}\left(\omega \right){\widehat{A}}_{\beta }\left(\omega \right),\end{array}$$

(13)

whose physical meaning is a renormalization of the open system Hamiltonian due to its interaction with the rapidly fluctuating environment [2], by analogy with the Lamb shift [29] in atomic physics. On the other hand, the dissipation super-operator in Equation (12) is explicitly given by [2]:

$$\begin{array}{c}\hfill {\mathcal{L}}_D{\widehat{\rho }}_S=\sum _{\alpha ,\beta }\sum _{\omega }{\gamma }_{\alpha \beta }\left(\omega \right)\left({\widehat{A}}_{\beta }\left(\omega \right){\widehat{\rho }}_S\left(t\right){\widehat{A}}_{\alpha }^{†}\left(\omega \right)-\frac{1}{2}\left\{{\widehat{A}}_{\alpha }^{†}\left(\omega \right){\widehat{A}}_{\beta }\left(\omega \right),{\widehat{\rho }}_S\left(t\right)\right\}\right).\end{array}$$

(14)

The coefficients in the dissipator and Lamb-shift operator are given by the real and imaginary parts of the Fourier transform of the bath auto-correlation functions:

$$\begin{array}{c}\hfill {\Gamma }_{\alpha \beta }\left(\omega \right)={\int }_0^{\infty }ds{e}^{i\omega s}\langle {\widehat{B}}_{\alpha }^{†}\left(s\right){\widehat{B}}_{\beta }\left(0\right)\rangle =\frac{1}{2}{\gamma }_{\alpha \beta }\left(\omega \right)+i{\sigma }_{\alpha \beta }\left(\omega \right).\end{array}$$

(15)

Notice that Equation (15) can be expressed in the standard Lindblad form Equation (4) by diagonalizing the matrices ${\gamma }_{\alpha \beta }$.

2.1. Relaxation to Equilibrium

If the reservoir is a macroscopic system, it can be assumed to be in thermal equilibrium, with a Gibbsian canonical distribution at inverse temperature β,

$$\begin{array}{c}\hfill {\widehat{\rho }}_B=\frac{exp(-\beta {\widehat{H}}_B)}{{\mathrm{Tr}}_{B}exp(-\beta {\widehat{H}}_B)}.\end{array}$$

(16)

Under this assumption, a natural question to pose is if the small open system S would relax to a stationary distribution corresponding to its own thermal distribution at the same temperature as the reservoir, i.e., if:

$$\begin{array}{c}\hfill \underset{t\to \infty }{lim}{\widehat{\rho }}_S\left(t\right)={\widehat{\rho }}_{S,eq}=\frac{exp(-\beta {\widehat{H}}_S)}{{\mathrm{Tr}}_{S}exp(-\beta {\widehat{H}}_S)}.\end{array}$$

(17)

It turns out that for the approximate Markovian dynamics described by the Born–Markov master Equation (12), this is indeed the case. From the Kubo–Martin–Schwinger (KMS) relation for the correlators of the reservoir operators:

$$\begin{array}{c}\hfill \langle {\widehat{B}}_{\alpha }^{†}\left(t\right){\widehat{B}}_{\beta }\left(0\right)\rangle ={\mathrm{Tr}}_B\left({\widehat{\rho }}_B{\widehat{B}}_{\alpha }^{†}\left(t\right){\widehat{B}}_{\beta }\left(0\right)\right)=\langle {\widehat{B}}_{\beta }\left(0\right){\widehat{B}}_{\alpha }^{†}(t+i\beta )\rangle \end{array}$$

(18)

that follows directly from the thermal distribution of the reservoir Equation (16), one has the properties:

$$\begin{array}{ccc}\hfill {\gamma }_{\alpha \beta }(-\omega )& =& e^{-\omega \beta }{\gamma }_{\alpha \beta }\left(\omega \right),\hfill \\ \hfill {\widehat{\rho }}_{S,eq}{\widehat{A}}_{\alpha }\left(\omega \right)& =& e^{\beta \omega }{\widehat{A}}_{\alpha }\left(\omega \right){\widehat{\rho }}_{S,eq},\hfill \\ \hfill {\widehat{\rho }}_{S,eq}{\widehat{A}}_{\alpha }^{†}\left(\omega \right)& =& e^{-\beta \omega }{\widehat{A}}_{\alpha }^{†}\left(\omega \right){\widehat{\rho }}_{S,eq}.\hfill \end{array}$$

(19)

Using these relations, it is straightforward to prove that the dissipation term vanishes for the canonical Gibbsian distribution Equation (17):

$${\mathcal{L}}_D{\widehat{\rho }}_{S,eq}=0.$$

(20)

On the other hand, it is evident that the Lamb-shift operator Equation (13) commutes with the open system Hamiltonian ${\widehat{H}}_S$, which implies:

$$\begin{array}{c}\hfill \left[{\widehat{H}}_{LS},{\widehat{H}}_S\right]=0⟹\left[{\widehat{H}}_S+{\widehat{H}}_{LS},{\widehat{\rho }}_{S,eq}\right]=0.\end{array}$$

(21)

Therefore, it follows from Equations (20) and (21) that the thermal distribution Equation (17) is indeed a steady-state solution for the small system master Equation (12) in the Born–Markov approximation, given that the macroscopic reservoir B is in a thermal equilibrium state.

2.2. Quasi-Static Evolution

In the previous section, we have discussed sufficient conditions under which the dynamical evolution of the reduced density matrix operator ${\widehat{\rho }}_S\left(t\right)$, when in contact with a macroscopic reservoir in thermal equilibrium, will relax towards the thermal distribution ${\widehat{\rho }}_{S,eq}$. Assuming that the system is in such an equilibrium state, in analogy with the classical thermodynamics analysis, we can conceive of a sequence of quasi-static processes that drives the subsystem along a sequence of equilibrium states. This is a very special type of dynamics, where we shall assume that one or more physical parameters in the set $\left\{{\lambda }_j\right\}$ (such as geometrical dimensions or external fields), on which the Hamiltonian ${\widehat{H}}_S\left(\left\{{\lambda }_j\right\}\right)$ depends explicitly, can be varied at an arbitrary slow rate ${\dot{\lambda }}_j$. To be more precise, let us assume that $|n;\{{\lambda }_j\}\rangle$ constitutes the set of eigenvectors of ${\widehat{H}}_S$:

$${\widehat{H}}_S|n;\left\{{\lambda }_j\right\}\rangle =E_n\left(\left\{{\lambda }_j\right\}\right)|n;\left\{{\lambda }_j\right\}\rangle ,$$

(22)

where n represents a set of indexes (can be continuous, not necessarily discrete) that labels the spectrum of the Hamiltonian. The equilibrium density matrix operator is diagonal in the energy eigenbasis:

$${\widehat{\rho }}_{S,eq}=\sum _n{p}_n\left(\left\{{\lambda }_j\right\}\right)|n;\left\{{\lambda }_j\right\}\rangle \langle n;\left\{{\lambda }_j\right\}|.$$

(23)

Due to the normalization condition ${\mathrm{Tr}}_S{\widehat{\rho }}_{S,eq}=1$, we have:

$$\sum _n{p}_n\left(\left\{{\lambda }_j\right\}\right)=1.$$

(24)

In this representation, the von Neumann entropy [30] adopts a simple expression in terms of the probability coefficients:

$$S\left(\left\{{\lambda }_j\right\}\right)=-k_B{\mathrm{Tr}}_S\left({\widehat{\rho }}_{S}ln{\widehat{\rho }}_S\right)=-k_B\sum _n{p}_n\left(\left\{{\lambda }_j\right\}\right)ln\left(p_n\left(\left\{{\lambda }_j\right\}\right)\right).$$

(25)

The ensemble-average energy $E=\langle {\widehat{H}}_S\rangle$ of the open system S is given by:

$$E={\mathrm{Tr}}_S\left({\widehat{\rho }}_{S,eq}{\widehat{H}}_S\right)=\sum _n{p}_n\left(\left\{{\lambda }_j\right\}\right)E_n\left(\left\{{\lambda }_j\right\}\right).$$

(26)

The statistical ensemble just described can be submitted to an arbitrary quasi-static process, involving the modulation of one or more of the parameters $\left\{{\lambda }_j\right\}$ [26,27], and hence, the ensemble-average energy in Equation (26) changes accordingly:

$$\begin{array}{ccc}\hfill dE& =& {\mathrm{Tr}}_S\left({\widehat{H}}_S\delta {\widehat{\rho }}_{S,eq}\right)+{\mathrm{Tr}}_S\left({\widehat{\rho }}_{S,eq}\delta {\widehat{H}}_S\right)\hfill \\ \hfill & =& \sum _n\sum _j{E}_n\left(\left\{{\lambda }_j\right\}\right)\frac{\partial }{\partial {\lambda }_j}{p}_n\left(\left\{{\lambda }_j\right\}\right)\delta {\lambda }_j+\sum _n\sum _j{p}_n\left(\left\{{\lambda }_j\right\}\right)\frac{\partial }{\partial {\lambda }_j}{E}_n\left(\left\{{\lambda }_j\right\}\right)\delta {\lambda }_j\hfill \\ & =& \delta Q+\delta W.\hfill \end{array}$$

(27)

The first term in Equation (27) represents a process along which the spectrum remains unaltered, while the second one represents an isentropic process. These two terms are in correspondence with the macroscopic notions of work and heat, respectively. Therefore, Equation (27) represents a microscopic version of the first law of thermodynamics for the statistical ensemble of single-particle systems [26,27,28].

We notice that Equation (27) strongly parallels the non-equilibrium relations described along Equations (6)–(8). Indeed, the quasi-static process described above via Equation (27) can be considered as a very particular form of a dynamical process under two considerations. First, we assume that the process is performed arbitrarily close to equilibrium, such that the dynamics is uniquely determined by the rate of change of the set of parameters $\left\{{\dot{\lambda }}_j\right\}$. Concretely, in a given interval of time $\delta t$, we have $\delta {\lambda }_j={\dot{\lambda }}_j\delta t$. In the second place, for the system to remain arbitrarily close to equilibrium at any instant, we must assume that the rates are slow enough in order to satisfy $\delta {\lambda }_j/{\dot{\lambda }}_j\gg {\tau }_S\gg {\tau }_B$, with ${\tau }_S$ and ${\tau }_B$ the characteristic relaxation time scales for the open system S and the reservoir B, respectively. This is in practice an idealization, since such a process would take of course a nearly infinite amount of time to yield a finite change in any of the quantities $\delta W=\dot{W}\delta t$ and $\delta Q=\dot{Q}\delta t$.

In the following sections, we shall develop in detail two explicit examples of the application of this formalism, each of them constituting a proof of principle for the concrete realization of a quantum heat engine.

3. Magnetically-Driven Quantum Engine on a Quantum Dot Array

A system that is physically possible to construct in the laboratory and of which many examples exist indeed [31,32] is a two-dimensional array of nearly identical, cylindrical semiconductor quantum dots. For our purposes, this array need not be periodic, but a sufficient distance between the dots is necessary in order to ensure that no charge tunneling effects occur among them. Let us further assume that a constant gate potential is applied to each dot in the array, such that the local charge density can be controlled. In a fairly idealized limit, we shall further assume that each dot is charged with a single electron, and hence, the whole array of non-interacting dots constitutes an ensemble of identical replicas of a single-particle system. Under these assumptions, let us consider the action of an external, uniform magnetic field perpendicular to the plane of the array. The geometric confinement imposed by the structure of the dots, in combination with the imposed magnetic field, leads to the emergence of discrete Landau levels, with a characteristic confinement length-scale defined by an effective Landau radius. This characteristic length can be tuned by modulating the intensity of the external magnetic field, such that a sequence of processes can be devised in order to construct a QHEN, as discussed in [27]. In particular, along the next sections, we shall discuss two different alternatives, that we shall the call iso-energetic cycle and the Carnot cycle, respectively.

3.1. The Single-Particle Spectrum in a Cylindrical Quantum Dot under an External Magnetic Field

As discussed in [27], we shall conceive of a simplified model for a cylindrical semiconductor quantum dot, through the confining potential [27,31,32].

$$\begin{array}{c}\hfill V_{dot}(x,y)=\frac{m^{*}}{2}{\omega }_d^2\left(x^2+y^2\right).\end{array}$$

(28)

Here, $m^{*}$ is the effective mass of an electron in this system. For instance, in Gallium-Arsenide $m^{*}\sim 0.067m_0$ [31,32]. In addition to this “geometric” confinement, we shall impose a uniform magnetic field directed along the z-axis, $\mathbf{B}=\widehat{z}B$. For convenience, we choose the symmetric gauge $\mathbf{A}=\frac{B}{2}(-y,x,0)$, such that the system is represented by the single-particle Hamiltonian:

$$\widehat{H}=\frac{1}{2m^{*}}\left[{\left(p_x-\frac{eBy}{2}\right)}^2+{\left(p_y+\frac{eBx}{2}\right)}^2\right]+V_{dot}(x,y).$$

(29)

The corresponding single-particle eigenstates correspond to Landau levels [32] with the energy spectrum:

$$E_{n_{\rho },m}\left(B\right)=\hslash \Omega (2n_{\rho }+\mid m\mid +1)-m\frac{\hslash {\omega }_B}{2}.$$

(30)

Here, $n_{\rho }=0,1,2,...$ and $m=0,\pm 1,\pm 2,..$ are the radial and azimuthal quantum numbers, respectively. These levels are characterized by an effective frequency:

$$\Omega =\sqrt{{\omega }_d^2+\frac{{\omega }_B^2}{4}},$$

(31)

where ${\omega }_B=eB/m^{*}$ is the field-dependent cyclotron frequency. The eigenfunctions are expressed in terms of associated Laguerre polynomials [32]:

$$\begin{array}{ccc}\hfill \langle \mathbf{r}|{\psi }_{\mathbf{n}}\left(B\right)\rangle & =& \sqrt{\frac{2}{l_{e,B}^2}}\sqrt{\frac{n_{\rho }!}{\left(n_{\rho }+\left|m\right|\right)!}}{\left(\frac{\rho }{l_{e,B}}\right)}^{\left|m\right|}{e}^{-\frac{{\rho }^2}{2l_{e,B}^2}}{L}_{n_{\rho }}^{\left|m\right|}\left(\frac{{\rho }^2}{l_{e,B}^2}\right),\hfill \end{array}$$

(32)

where $l_{e,B}=\sqrt{\hslash /\left(m^{*}\Omega \right)}$ is the effective Landau radius that captures the confining effects of the geometric potential, as well as the external magnetic field, and $\mathbf{n}\equiv (n_{\rho },m)$.

The single-particle spectrum after Equation (30) is displayed in Figure 2, in units of the energy scale $\hslash {\omega }_d$, as a function of the external magnetic field, for the first fourteen eigenstates. From the figure, it is clear that, even at arbitrary large magnetic fields, the ground state $(n_{\rho }=0,m=0)$ is non-degenerate. The excited states, on the other hand, are non-degenerate except for a discrete set of crossovers as the field intensity increases. In particular, the excited state $(n_{\rho }=1,m=0)$ exhibits only a single crossover with the state $(n_{\rho }=0,m=3)$.

3.2. The Iso-Energetic Cycle

An elementary scheme for a quantum heat-engine, originally proposed by Bender et al. [4,5] in the context of a single Schrödinger particle and more recently extended by us to a relativistic Dirac particle [26], is the iso-energetic cycle . The cycle is composed of two iso-entropic and two iso-energetic trajectories. In particular, during the iso-energetic trajectories, the “working substance” must exchange energy with an energy reservoir [6,7], such as single-mode radiation in a cavity [6].

In [27], we proposed a series of quasi-static processes (see Figure 3), where the confinement length to be modified is the effective Landau radius:

$$\begin{array}{c}\hfill l_{e,B}=\sqrt{\hslash /(m^{*}\Omega )}={\left(l_d^{-4}+l_B^{-4}/4\right)}^{-1/4},\end{array}$$

(33)

with $l_d=\sqrt{\hslash /\left(m^{*}{\omega }_d\right)}$ representing the geometric confinement imposed by the cylindrical potential well. The effective Landau radius $l_{e,B}$ can be modified through $l_B$, by quasi-statically adjusting the intensity of the external magnetic field. Along these trajectories, according to the general theory presented in Section 2, the total change in the ensemble average energy of the system is given by:

$$\begin{array}{ccc}\hfill dE& =& \sum _{\mathbf{n}}{p}_{\mathbf{n}}\left(B\right)d{E}_{\mathbf{n}}\left(B\right)+\sum _{n}d{p}_{\mathbf{n}}\left(B\right)E_{\mathbf{n}}\left(B\right)=\delta W+\delta Q\hfill \end{array}$$

(34)

where we have introduced the two-valued index $\mathbf{n}\equiv (n_{\rho },m)$ to enumerate the eigenstates of the Hamiltonian, defined in Equation (30).

We shall analyze in detail each type of process. For an iso-entropic process defined by the constraint [33] $\left\{p_{\mathbf{n}}\left(B\right)\right\}=const.$, the quasi-static work performed in varying the external magnetic field B is related to the ensemble-average magnetization $M=-{\left(\partial E/\partial B\right)}_S$ (see Figure 4), $\delta W=-MdB$. Hence, if the magnetic field is varied from $B=B_{\gamma }$ to $B=B_{\delta }$, the total work performed by the system is given by the expression:

$$\begin{array}{ccc}\hfill W_{\gamma \to \delta }& =& {\int }_{B_{\gamma }}^{B_{\delta }}dB{\left(\frac{\partial E}{\partial B}\right)}_{\{p_{\mathbf{n}}\left(B_{\gamma }\right)=p_{\mathbf{n}}\left(B_{\delta }\right)\}=const.}=\sum _{\mathbf{n}}{p}_{\mathbf{n}}\left(B_{\gamma }\right)\left[E_{\mathbf{n}}\left(B_{\delta }\right)-E_{\mathbf{n}}\left(B_{\gamma }\right)\right].\hfill \end{array}$$

(35)

Here, we have chosen a sign convention, such that for an expansion process $l_{e,B_{\delta }}>l_{e,B_{\gamma }}$, the work performed by the system is negative [34], indicating that the ensemble-averaged energy is decreasing during expansion.

On the other hand, the parametric trajectory for an iso-energetic process is defined by the equation:

$$\begin{array}{c}\hfill dE=\delta W_{\gamma \to \delta }+\delta Q_{\gamma \to \delta }=0.\end{array}$$

(36)

The solution to this equation, for $B\in [B_{\delta },B_{\gamma }]$, is given by the path:

$$\begin{array}{c}\hfill \sum _{\mathbf{n}}{p}_{\mathbf{n}}\left(B\right)E_{\mathbf{n}}\left(B\right)=\sum _{\mathbf{n}}{p}_{\mathbf{n}}\left(B_{\gamma }\right)E_{\mathbf{n}}\left(B_{\gamma }\right),\end{array}$$

(37)

along with the normalization condition Equation (24). The integral of Equation (36) from $B_{\gamma }\to B_{\delta }$ yields:

$$\begin{array}{c}\hfill \Delta E=W_{\gamma \to \delta }+Q_{\gamma \to \delta }=0.\end{array}$$

(38)

Here, the first term $W_{\gamma \to \delta }$ corresponds to the magnetic work performed by the system when changing its induced magnetization due to the applied field, at constant total energy. The second term $Q_{\gamma \to \delta }=-W_{\gamma \to \delta }$ corresponds to the amount of energy exchanged by the system with the environment, in order to rearrange its internal level occupation. Therefore, for an iso-energetic process, the heat exchanged by the system with the environment is:

$$\begin{array}{c}\hfill Q_{\gamma \to \delta }=\sum _{\mathbf{n}}{\int }_{B_{\gamma }}^{B_{\delta }}{E}_{\mathbf{n}}\left(B\right)\frac{d{p}_{\mathbf{n}}\left(B\right)}{dB}dB.\end{array}$$

(39)

Evidently, Equation (37) combined with the normalization condition Equation (24) is not enough to uniquely define the coefficients $p_{\mathbf{n}}\left(B\right)$ along an iso-energetic trajectory. An exception is the case when the energy scale of all of the processes involved is such that only transitions between two adjacent levels are possible. Since the only driving force to induce transitions is a quasi-static variation of the magnetic field intensity, the azimuthal symmetry of the Hamiltonian Equation (29) is conserved at any moment along this process. Therefore, the transition probability amplitudes impose the selection rule ${\mathcal{T}}_{\mathbf{1}\to \mathbf{2}}\propto \int d\phi e^{i(m_1-m_2)\phi }={\delta }_{m_1,m_2}$, and hence, the azimuthal quantum numbers of initial and final states must be the same $m_1=m_2$, that is angular momentum $L_z$ is conserved for such a transition.

Let us focus on the particular case when the two states involved are the ground state and the first accessible excited state, respectively: $\mathbf{1}\equiv (n_{\rho }=0,m=0)$ and $\mathbf{2}\equiv (n_{\rho }=1,m=0)$ (see Figure 2 and Figure 3). It is clear that transitions between these two energy levels, when the quasi-static variation of the magnetic field intensity is the only driving force, are allowed by symmetry, since $m_1=m_2=0$ for both eigenstates. Moreover, by looking at Figure 2, it is clear that $(n_{\rho }=1,m=0)$ is the lowest excited state accessible from the ground state $(n_{\rho }=0,m=0)$ that respects the selection rule imposed by angular momentum conservation. From Figure 2, it is also evident that a crossover occurs between the excited states $(n_{\rho }=1,m=0)$ and $(n_{\rho }=0,m=3)$. However, transitions between the ground state $(n_{\rho }=0,m=0)$ and the excited state $(n_{\rho }=0,m=3)$ or between $(n_{\rho }=0,m=3)$ and $(n_{\rho }=1,m=0)$ are forbidden by angular momentum conservation, if the quasi-static variation of the magnetic field intensity is the only driving force, as discussed above. The next excited state that would respect the symmetry is $(n_{\rho }=2,m=0)$, but as clearly seen in Figure 2, this one is quite high in energy and does not cross at any point with $(0,0)$ nor $(1,0)$. Therefore, under these considerations, the ground state $(n_{\rho }=0,m=0)$ and the excited state $(n_{\rho }=1,m=0)$ constitute an effective two-level system.

In practice, single quantum dots can be prepared in the ground state $(0,0)$ by coupling them to a modified high-Q single-defect cavity [35] due to Purcell’s effect [36]. Once the system has been prepared in its ground state, the iso-energetic cycle can start by quasi-statically tuning the static magnetic field in the absence of external radiation sources.

The effective two-level system constituted by the states $(0,0)$ and $(1,0)$ as described before, along with the different trajectories involved in the cycle, is represented in Figure 3. Combining Equation (37) with the normalization condition Equation (24), the iso-energetic trajectories are described by:

$$\begin{array}{ccc}\hfill p_{\mathbf{1}}\left(B\right)& =& \frac{E_{\mathbf{2}}\left(B_1\right)-E_{\mathbf{2}}\left(B\right)}{E_{\mathbf{1}}\left(B\right)-E_{\mathbf{2}}\left(B\right)}+\frac{E_{\mathbf{1}}\left(B_1\right)-E_{\mathbf{2}}\left(B_1\right)}{E_{\mathbf{1}}\left(B\right)-E_{\mathbf{2}}\left(B\right)}{p}_{\mathbf{1}}\left(B_1\right),\hfill \end{array}$$

(40)

with $p_{\mathbf{2}}\left(B\right)=1-p_{\mathbf{1}}\left(B\right)$ after the normalization condition Equation (24). The heat exchanged by the system with the environment during the iso-energetic trajectory connecting the initial and final states $B_1\to B_2$ is given by the expression:

$$\begin{array}{ccc}\hfill Q_{1\to 2}& =& \left[E_{\mathbf{2}}\left(B_1\right)+\left(E_{\mathbf{1}}\left(B_1\right)-E_{\mathbf{2}}\left(B_1\right)\right)p_{\mathbf{1}}\left(B_1\right)\right]ln\left[\frac{E_{\mathbf{1}}\left(B_2\right)-E_{\mathbf{2}}\left(B_2\right)}{E_{\mathbf{1}}\left(B_1\right)-E_{\mathbf{2}}\left(B_1\right)}\right].\hfill \end{array}$$

(41)

For the effective two-level system described in Figure 3, we consider the cycle that starts in the ground state with $p_{\mathbf{1}}\left(B_1\right)=1$. Then, the system experiences an iso-energetic expansion from $l_{B_1}\to l_{B_2}>l_{B_1}$, followed by an iso-entropic expansion from $l_{B_2}\to l_{B_3}>l_{B_2}$. Then, it experiences an iso-energetic compression $l_{B_3}\to l_{B_4}<l_{B_3}$, to finally return to its initial ground state through an iso-entropic compression $l_{B_4}\to l_{B_1}$.

For definiteness, let us consider that after the iso-energetic process $l_{B_1}\to l_{B_2}$, the system ends completely localized in the excited state $\mathbf{n}=\mathbf{2}$. In this condition, we have:

$$\begin{array}{c}\hfill p_{\mathbf{1}}\left(B_2\right)=0,p_{\mathbf{2}}\left(B_2\right)=1.\end{array}$$

(42)

The conservation of total energy between the initial and final states connected through the iso-energetic process leads to the equation:

$$\begin{array}{c}\hfill p_{\mathbf{1}}\left(B_1\right)E_{\mathbf{1}}\left(B_1\right)=p_{\mathbf{2}}\left(B_2\right)E_{\mathbf{2}}\left(B_2\right),\end{array}$$

(43)

where $p_{\mathbf{1}}\left(B_1\right)=p_{\mathbf{2}}\left(B_2\right)=1$ for maximal expansion. Therefore, given the spectrum in Equation (30), Equation (43) implies that $l_{B_2}/l_{B_1}={\alpha }_1$, where ${\alpha }_1$ is determined by the condition:

$$\begin{array}{c}\hfill \hslash {\omega }_d\sqrt{1+N_{{\Phi }_1}^2}=3\hslash {\omega }_d\sqrt{1+\frac{N_{{\Phi }_1}^2}{{\alpha }_1^4}}.\end{array}$$

(44)

Here, we have defined $N_{{\Phi }_1}=l_d^2/\left(2l_{B_1}^2\right)={\Phi }_{B_1}/{\Phi }_0$ as the number of flux quanta ${\Phi }_0=h/\left(2e\right)$ piercing the area $\pi l_d^2$. Equation (44) possesses physically meaningful solutions:

$$\begin{array}{c}\hfill {\alpha }_1=\frac{\sqrt{3N_{{\Phi }_1}}}{{\left(N_{{\Phi }_1}^2-8\right)}^{1/4}}\end{array}$$

(45)

when $N_{{\Phi }_1}>2\sqrt{2}$. Therefore, the minimal initial value of the external field required to perform the cycle is $B_{1,min}=4\sqrt{2}\hslash /\left(e{l}_d^2\right)$. For instance, if one considers a typical size of a semiconductor quantum dot of $l_d=70$ nm [31], the minimal initial field is $B_{1,min}\sim 4.8$ T. The heat exchanged with the environment along this first stage of the cycle is calculated after Equation (39):

$$\begin{array}{c}\hfill Q_{1\to 2}=E_{\mathbf{1}}\left(B_1\right)ln\left[\frac{E_{\mathbf{2}}\left(B_1\right)-E_{\mathbf{1}}\left(B_1\right)}{E_{\mathbf{2}}\left({\alpha }_1^{-2}{B}_1\right)-E_{\mathbf{1}}\left({\alpha }_1^{-2}{B}_1\right)}\right],\end{array}$$

(46)

while the work performed is $W_{1\to 2}=-Q_{1\to 2}$.

The next process along the cycle is an iso-entropic expansion (see Figure 3), characterized by the condition $p_{\mathbf{2}}\left(B_2\right)=p_{\mathbf{2}}\left(B_3\right)=1$. We shall define the expansion parameter $\alpha \equiv l_{B_3}/l_{B_2}>1$. Notice that $\alpha >1$ can be arbitrarily chosen. Along this stage, we have $l_{B_2}={\alpha }_1{l}_{B_1}$.

The cycle continues with a maximal compression process from $l_{B_3}={\alpha }_1\alpha l_{B1}$ to $l_{B_4}={\alpha }_3{\alpha }_1\alpha l_{B1}$ under iso-energetic conditions (see Figure 3). The condition for energy conservation is in this case similar to Equation (43), implying $p_{\mathbf{2}}\left(B_3\right)=p_{\mathbf{1}}\left(B_4\right)=1$ and:

$$\begin{array}{c}\hfill \hslash {\omega }_d\sqrt{1+\frac{N_{{\Phi }_1}^2}{{\left({\alpha }_1{\alpha }_3\alpha \right)}^4}}=3\hslash {\omega }_d\sqrt{1+\frac{N_{{\Phi }_1}^2}{{\left({\alpha }_1\alpha \right)}^4}}.\end{array}$$

(47)

The solution to Equation (47) fixes the value for the compression coefficient ${\alpha }_3<1$:

$$\begin{array}{c}\hfill {\alpha }_3=\frac{N_{{\Phi }_1}^{1/2}}{{\left(8{\left({\alpha }_1\alpha \right)}^4+9N_{{\Phi }_1}^2\right)}^{1/4}}.\end{array}$$

(48)

The heat exchanged by the system with the environment along this process, applying Equation (39), is given by the expression:

$$\begin{array}{ccc}\hfill Q_{3\to 4}& =& E_{\mathbf{2}}\left(\frac{B_1}{{\left(\alpha {\alpha }_1\right)}^2}\right)ln\left[\frac{E_{\mathbf{1}}\left(\frac{B_1}{{\left(\alpha {\alpha }_1{\alpha }_3\right)}^2}\right)-E_{\mathbf{2}}\left(\frac{B_1}{{\left(\alpha {\alpha }_1{\alpha }_3\right)}^2}\right)}{E_{\mathbf{1}}\left(\frac{B_1}{{\left(\alpha {\alpha }_1\right)}^2}\right)-E_{\mathbf{2}}\left(\frac{B_1}{{\left(\alpha {\alpha }_1\right)}^2}\right)}\right],\hfill \end{array}$$

(49)

and the work performed is $W_{3\to 4}=-Q_{3\to 4}$. The last path along the cycle is an adiabatic process (see Figure 3), which returns the system to its initial ground state with $p_{\mathbf{1}}\left(B_4\right)=p_{\mathbf{1}}\left(B_1\right)=1$.

The efficiency of the cycle is defined by the ratio [27]:

$$\begin{array}{c}\hfill \eta \left(N_{{\Phi }_1},\alpha \right)=1-\left|\frac{Q_{3\to 4}}{Q_{1\to 2}}\right|=1-3\frac{{\Theta }_1\left(\alpha {\alpha }_1\right)}{{\Theta }_1\left(1\right)}\frac{ln\left[\frac{{\Theta }_1\left(\alpha {\alpha }_1{\alpha }_3\right)}{{\Theta }_1\left(\alpha {\alpha }_1\right)}\right]}{ln\left[\frac{{\Theta }_1\left(1\right)}{{\Theta }_1\left({\alpha }_1\right)}\right]},\end{array}$$

(50)

where we have defined ${\Theta }_1\left(\alpha \right)=\sqrt{1+N_{{\Phi }_1}^2/{\alpha }^4}$. The trend of the efficiency is shown in Figure 5 as a function of the expansion parameter α, for different values of the initial external field $B_1$ (expressed in terms of $N_{{\Phi }_1}$). For very large fields $N_{{\Phi }_1}\gg 1$, one has from Equations (45) and (48) that ${\alpha }_1=1/{\alpha }_3=\sqrt{3}$, and hence, the efficiency tends to the asymptotic limit [27]:

$$\begin{array}{c}\hfill \eta \to 1-1/{\alpha }^2,N_{{\Phi }_1}\gg 1.\end{array}$$

(51)

Remarkably, this asymptotic result has been obtained before for an iso-energetic cycle driven by a mechanical external force, both in the Schrödinger [4], as well as in the low-energy limit for the Dirac particle case [26]. This suggests that it may represent a universal maximal efficiency for any quantum mechanical engine based on the iso-energetic cycle construction.

3.3. The Quantum Carnot Cycle

A quantum mechanical version of the Carnot cycle, as applied to the statistical ensemble of quantum dots under consideration, is composed of four stages: two isothermal and two iso-entropic processes, as shown in Figure 6.

Along the first stage of the cycle, the system is brought into contact with a thermal reservoir at temperature $T_H$. Under isothermal conditions, the Landau radius is expanded from $l_{e,B_1}\to l_{e,B_2}$. Since thermal equilibrium with the reservoir is assumed along this process, the von Neumann entropy of the system achieves a maximum for the Boltzmann distribution [30,37]:

$$\begin{array}{c}\hfill p_{\mathbf{n}}(B,{\beta }_H)={\left[Z(B,{\beta }_H)\right]}^{-1}{e}^{-{\beta }_H{E}_{\mathbf{n}}\left(B\right)},\end{array}$$

(52)

with $\beta ={\left(k_{B}T\right)}^{-1}$, and the normalization factor is given by the partition function (see [27]):

$$Z(B,\beta )=\sum _{\mathbf{n}}{e}^{-\beta E_{\mathbf{n}}}=Z^{+}{Z}^{-}.$$

(53)

Here, we have defined:

$$Z^{\pm }=\frac{1}{2sinh\left(\frac{\hslash \beta {\omega }_{\pm }}{2}\right)},$$

(54)

with:

$${\omega }_{\pm }=\Omega \pm \frac{{\omega }_B}{2}.$$

(55)

From the general theory presented in Section 2, the heat absorbed by the system from the thermal reservoir is given by:

$$\begin{array}{ccc}\hfill Q_{1\to 2}& =& {\int }_{B_1}^{B_2}\sum _{\mathbf{n}}{E}_{\mathbf{n}}\left(B\right)\frac{d{p}_{\mathbf{n}}(B,{\beta }_H)}{dB}dB=E(B_2,{\beta }_H)-E(B_1,{\beta }_H)+{\beta }_H^{-1}ln\left[\frac{Z(B_2,{\beta }_H)}{Z(B_1,{\beta }_H)}\right].\hfill \end{array}$$

(56)

The final result follows from substituting the explicit expression for the partition function Equation (53) and the definition of the ensemble-averaged energy of the single-particle system $E=\langle \widehat{H}\rangle =-\partial lnZ/\partial \beta$,

$$\begin{array}{c}\hfill E(B,\beta )=\frac{\hslash {\omega }_{+}}{2}coth\left(\frac{\beta \hslash {\omega }_{+}}{2}\right)+\frac{\hslash {\omega }_{-}}{2}coth\left(\frac{\beta \hslash {\omega }_{-}}{2}\right).\end{array}$$

(57)

Similarly, during the third stage of the cycle (see Figure 6), the system is again brought into contact with a thermal reservoir, but at a lower temperature $T_C<T_H$. Therefore, the probability distribution of states in the ensemble is $p_{\mathbf{n}}(B,{\beta }_C)$, as defined in Equation (52), but with $T_C$ instead of $T_H$. The heat released to the reservoir during this stage is given by the expression:

$$\begin{array}{c}\hfill Q_{3\to 4}=E(B_4,{\beta }_C)-E(B_3,{\beta }_C)+{\beta }_C^{-1}ln\left[\frac{Z(B_4,{\beta }_C)}{Z(B_3,{\beta }_C)}\right].\end{array}$$

(58)

The second and fourth stages of the cycle constitute iso-entropic trajectories (see Figure 6). When substituting the Boltzmann distribution $p_{\mathbf{n}}(\beta ,B)={\left[Z(\beta ,B)\right]}^{-1}exp(-\beta E_{\mathbf{n}}\left(B\right))$ into the expression for the von Neumann entropy Equation (25), we obtain the relation:

$$\begin{array}{c}\hfill S/k_B=\beta E+lnZ(\beta ,B).\end{array}$$

(59)

The equation of state is obtained from Equation (57) as:

$$\begin{array}{ccc}\hfill M& =& -{\left(\frac{\partial E}{\partial B}\right)}_S=-{\mu }_B\frac{{\omega }_{+}}{\Omega }coth\left(\frac{\beta \hslash {\omega }_{+}}{2}\right)+{\mu }_B\frac{{\omega }_{-}}{\Omega }coth\left(\frac{\beta \hslash {\omega }_{-}}{2}\right),\hfill \end{array}$$

(60)

with M the ensemble-average magnetization as a function of the external magnetic field (see Figure 7) and ${\mu }_B=e\hslash /\left(2m^{*}\right)$ the Bohr magneton. In the last line, we made use of the explicit analytical expression Equation (53) for the partition function to calculate the derivative. The work performed during the second stage of the process is $W_{2\to 3}=E(B_2,{\beta }_H)-E(B_3,{\beta }_C)$. We are now in conditions to discuss the second and fourth stages of the Carnot cycle. These iso-entropic trajectories impose implicit conditions for the intensities of the magnetic field,

$$\begin{array}{ccc}\hfill \Delta S_{2\to 3}& =& S(B_3,{\beta }_C)-S(B_2,{\beta }_H)=0,\hfill \\ \hfill \Delta S_{4\to 1}& =& S(B_1,{\beta }_H)-S(B_4,{\beta }_C)=0.\hfill \end{array}$$

(61)

Expressing the entropies in terms of Equation (59), in combination with Equations (57) and (54), we obtain the condition:

$$\begin{array}{c}\hfill {\beta }_C\left[E(B_3,{\beta }_C)-E(B_4,{\beta }_C)\right]+ln\left[\frac{Z(B_3,{\beta }_C)}{Z(B_4,{\beta }_C)}\right]={\beta }_H\left[E(B_2,{\beta }_H)-E(B_1,{\beta }_H)\right]+ln\left[\frac{Z(B_2,{\beta }_H)}{Z(B_1,{\beta }_H)}\right].\end{array}$$

(62)

The fourth and final stage of the cycle also corresponds to an iso-entropic trajectory (see Figure 6) where $l_{B_4}\to l_{B_1}$, and the work performed by the system is given by $W_{4\to 1}=E(B_4,{\beta }_C)-E(B_1,{\beta }_H)$.

The efficiency of the quantum Carnot cycle is given by:

$$\begin{array}{ccc}\hfill {\eta }^C& =& 1-\frac{Q_{3\to 4}}{Q_{1\to 2}}=1-\frac{T_C}{T_H},\hfill \end{array}$$

(63)

where we have made use of Equations (56) and (58) to obtain the second equality. Remarkably, the efficiency is identical to the classical Carnot cycle. This result is in agreement with what we found in a recent work, where the efficiency for a mechanically-driven quantum heat engine based on a relativistic Dirac particle was studied [26].

4. A Magneto-Strain-Driven Quantum Engine on a Graphene Layer

Graphene, the ultimate slab of graphite [38], is an atomic monolayer of carbon atoms arranged in a honeycomb lattice that possesses remarkable mechanical, thermal and electronic properties [39,40,41,42,43,44]. From the electronic point of view, it is a semi-metal characterized by the presence of Dirac points in its band structure [39,40] that confer pseudo-relativistic properties to the charge carriers in this material. From the mechanical point of view [43], it is extremely flexible as seen by the presence of a phonon bending mode with quadratic dispersion near the Γ-point (see Figure 8), but at the same time extremely strong, with the highest Young modulus known to any material up to date [45]. The nearly ballistic propagation of phonons in graphene confer to this material a remarkably high thermal conductivity that can reach up to 4000 W/(m·K) [43,46].

The interplay between mechanics and electronics, determined by electron-phonon interactions [44] in graphene, induces interesting effects over its transport properties. In particular, the effect of strain-induced pseudo-magnetic fields on the electronic properties of graphene has been extensively discussed in the literature [39,40,41,48,49]. Moreover, by strain engineering, it is possible to generate nearly homogeneous pseudo-magnetic fields [40,41], thus allowing for the emergence of relativistic Landau levels, as confirmed experimentally [48]. In this chapter, we shall describe a conceptual model for a QHEN that is based on these singular properties.

4.1. The Single-Particle Spectrum

The components of the in-plane pseudo-vector potential ${\mathbf{A}}_S$ induced by mechanical strain are [40,41,49]:

$$\begin{array}{c}\hfill A_{S,1}=\frac{\beta }{2a}\left(u_{11}-u_{22}\right),A_{S,2}=\frac{\beta }{2a}\left(-2u_{12}\right).\end{array}$$

(64)

Here, the local displacement vector is defined as $\mathbf{u}=(u_1,u_2,z)$, with $u_i$ the in-plane components and z the out-of-plane component, while $\beta =\partial lnt/\partial lna$ is the relative change of the hopping parameter $\delta t/t$ with respect to the dilation of the lattice constant $\delta a/a$ [50]. The strain tensor $\left[u_{ij}\right]$ is defined by:

$$\begin{array}{c}\hfill u_{ij}=\frac{1}{2}\left({\partial }_i{u}_j+{\partial }_j{u}_i+{\partial }_{i}z{\partial }_{j}z\right).\end{array}$$

(65)

The pseudo-magnetic field generated by the strain-induced gauge potential is given by ${\mathbf{B}}_S={\widehat{e}}_3\left({\partial }_1{A}_{S,2}-{\partial }_2{A}_{S,1}\right)$.

It was recently shown [49,51,52] that the continuum expansion of the tight-binding Hamiltonian for strained graphene in the vicinity of a single Dirac point (valley), up to $O(u_{ij}^2)$, is given by the expression:

$$\begin{array}{ccc}\hfill \widehat{H}& =& -i\int d^{2}x{\widehat{\psi }}^{†}\left(x\right)\left[v_{ij}\left(x\right){\widehat{\sigma }}_i{\partial }_j+i{v}_F{\widehat{\sigma }}_i{A}_{S,i}+v_F{\widehat{\sigma }}_i{\Gamma }_i-\overline{\gamma }{B}_S{\widehat{\sigma }}_3\right]\widehat{\psi }\left(x\right).\hfill \end{array}$$

(66)

Here, $v_F\sim {10}^{6}m{s}^{-1}$ is the Fermi velocity for the undistorted graphene lattice and ${\widehat{\sigma }}_i$ the Pauli matrices. The presence of strain induces a local position dependence on the velocity [49,51], which then becomes a tensor:

$$\begin{array}{c}\hfill v_{ij}\left(x\right)=v_F\left({\delta }_{ij}-\frac{\beta }{4}\left(2u_{ij}+{\delta }_{ij}{u}_{kk}\right)+{\tilde{u}}_{ij}\right).\end{array}$$

(67)

Here, ${\tilde{u}}_{ij}=\frac{1}{2}\left({\partial }_i{u}_j+{\partial }_j{u}_i\right)$. On top of the gauge field ${\mathbf{A}}_S$, a “geometric” vector field ${\Gamma }_i$ appears, whose origin can be understood within the covariant picture [51] as the pseudo-spin connection for fermions propagating in a curved space [49]:

$$\begin{array}{c}\hfill {\Gamma }_i=\frac{1}{2v_F}{\partial }_j{v}_{ij}=-\frac{\beta }{4}\left({\partial }_j{u}_{ij}+\frac{1}{2}{\partial }_i{u}_{jj}\right)+\frac{1}{2}{\partial }_j{\tilde{u}}_{ij}.\end{array}$$

(68)

The last term in the Hamiltonian Equation (66) represents a pseudo-Zeeman coupling between the pseudo-spin degree of freedom and the strain pseudo-magnetic field ${\mathbf{B}}_S$ [52].

Let us consider the deformation field (see Figure 9):

$$\begin{array}{c}\hfill u_1=-2u_S{x}_1{x}_2,u_2=-u_S(x_1^2-x_2^2),z=0,\end{array}$$

(69)

where $u_S$ is a constant characterizing the magnitude of the in-plane displacement, in association through Equation (64) with the vector potential:

$$\begin{array}{c}\hfill {\mathbf{A}}_S=2\frac{\beta u_S}{a}\left(-x_2,x_1,0\right).\end{array}$$

(70)

This in turn implies a pseudo-magnetic field:

$$\begin{array}{c}\hfill {\mathbf{B}}_S=4\frac{\beta u_S}{a}{\widehat{e}}_3,\end{array}$$

(71)

which is constant in magnitude and points along the out-of-plane direction ${\widehat{e}}_3$, where in the last step, we used that for graphene $\beta =\partial lnt/\partial lna\sim 2$ [49,50]. Therefore, for the deformation field defined by Equation (69), the velocity is a constant diagonal tensor $v_{ij}=v_F{\delta }_{ij}$.

The Hamiltonian in Equation (66) describes the physics in the vicinity of a single Dirac point (valley). The inclusion of the pseudo-magnetic field terms breaks the valley degeneracy since, as opposed to a real magnetic field, the strain-induced gauge potential has opposite signs in the vicinity of each Dirac point [39]. Therefore, we define the Hamiltonian describing both valleys by the higher-dimensional structure:

$$\begin{array}{c}\hfill \widehat{H}=\int d^{2}x\left({\widehat{\Psi }}^{(+)†}\left(x\right),{\widehat{\Psi }}^{(-)†}\left(x\right)\right)\left[\begin{array}{cc}{\widehat{h}}^{+}& 0\\ 0& {\widehat{h}}^{-}\end{array}\right]\left(\begin{array}{c}{\widehat{\Psi }}^{(+)}\left(x\right)\\ {\widehat{\Psi }}^{(-)}\left(x\right)\end{array}\right).\end{array}$$

(72)

Here, $\xi =\pm$ represents each of the $K_{\xi }$ valleys, corresponding to the two inequivalent points $\xi \frac{4\pi }{3\sqrt{3}a}{\widehat{\mathbf{e}}}_1$ in the first Brillouin zone, respectively. By adopting the standard convention for the components of the spinor field at each sub-lattice $(A,B)$, we have defined the spinor:

$$\begin{array}{c}\hfill {\widehat{\Psi }}^{\left(\xi \right)}\left(x\right)={\left({\psi }_{A\uparrow }^{\left(\xi \right)},{\psi }_{B\uparrow }^{\left(\xi \right)},{\psi }_{A\downarrow }^{\left(\xi \right)},{\psi }_{B\downarrow }^{\left(\xi \right)}\right)}^T=\sum _{s=\uparrow ,\downarrow }{\chi }_s\otimes {\widehat{\Psi }}_s^{\left(\xi \right)}.\end{array}$$

(73)

Here, ${\widehat{\psi }}_s^{\left(\xi \right)}$ is the two-component pseudo-spinor at valley $K_{\xi }$ arising from the bipartite graphene lattice, whereas the electronic spinors ${\chi }_{\uparrow }={(1,0)}^T$ and ${\chi }_{\downarrow }={(0,1)}^T$ are eigenvectors of ${\widehat{\sigma }}_3$ with eigenvalues $s=\{\pm \}$.

The Hamiltonian operator defined in Equation (72) includes the effect of mechanical strain, via the pseudo-vector potential in Equation (70), whose curl is the pseudo-magnetic field in Equation (71). This effect on itself is sufficient to generate pseudo-relativistic Landau levels in the energy spectrum.

As we shall discuss in more detail in the next section, for the purpose of constructing a QHEN, the characteristic confinement length-scale given by the Landau radius must be tunable in a quasi-static way. In practice, it is difficult to adjust the strain applied to a nanoscale structure under realistic experimental conditions. Therefore, in addition to mechanical strain, we consider the effect of a uniform magnetic field $\mathbf{B}={\widehat{\mathbf{e}}}_{3}B$ that is easier to control experimentally. For convenience, we choose the gauge $\mathbf{A}=\frac{B}{2}\left(-x_2,x_1,0\right)$, which has the same form as the strain-induced pseudo-vector potential in Equation (70). For the combined magnetic field and strain, the effective Hamiltonian density at each valley $K_{\xi }$, including the electronic spin degree of freedom, becomes [28]:

$$\begin{array}{ccc}\hfill {\widehat{h}}^{\left(\xi \right)}& =& \xi v_F\left(\mathbf{1}\otimes \widehat{\mathit{\sigma }}\right)·\left[-i\nabla +\xi {\mathbf{A}}_S+\mathbf{A}\right]-\xi \overline{\gamma }{B}_S\left(\mathbf{1}\otimes {\widehat{\sigma }}_3\right)-\gamma B\left({\widehat{\sigma }}_3\otimes \mathbf{1}\right).\hfill \end{array}$$

(74)

The last term in the Hamiltonian Equation (74) is the Zeeman interaction, characterized by a coupling constant $\gamma =g{\mu }_B/2$, with $g\sim 1.8$ for graphene [53]. Notice that the Zeeman coupling only involves the real magnetic field $\mathbf{B}$, since the strain-induced pseudo-magnetic field ${\mathbf{B}}_S$ does not interact with the electronic spin.

The Hamiltonian Equation (72) determines a system of two decoupled Dirac equations in first quantization, one for each electronic spin component $s=\{\pm \}$,

$$\begin{array}{c}\hfill \left(\xi v_F\widehat{\mathit{\sigma }}·\left[-i\nabla +\xi {\mathbf{A}}_S+\mathbf{A}\right]-s\gamma B\mathbf{1}-\xi \overline{\gamma }{B}_S{\widehat{\sigma }}_3\right){\psi }_s^{\left(\xi \right)}\left(x\right)=E_s^{\xi }{\psi }_s^{\left(\xi \right)}\left(x\right).\end{array}$$

(75)

The pseudo-spinor eigenstates in Equation (75) are ${\psi }_{\lambda ,n,m,s}^{\left(\xi \right)}\left(x\right)\equiv \langle x|{\psi }_{\lambda ,n,m,s}^{\left(\xi \right)}\rangle$, with (see [28] for details):

$$\begin{array}{ccc}\hfill |{\psi }_{\lambda ,n,m,s}^{\left(\xi \right)}\rangle & =& \frac{1}{\sqrt{2}}\left(\begin{array}{c}{\alpha }_n|n-1,m\rangle \\ {\beta }_n|n,m\rangle \end{array}\right),n>0,\hfill \\ \hfill |{\psi }_{n=0,m,s}^{\left(\xi \right)}\rangle & =& \left(\begin{array}{c}0\\ |0,m\rangle \end{array}\right),n=0.\hfill \end{array}$$

(76)

Here, the coefficients ${\alpha }_n$ and ${\beta }_n$ are real constants [28]. The corresponding energy eigenvalues are [28]:

$$\begin{array}{c}\hfill E_{n,s}^{\xi }=\left\{\begin{array}{cc}\lambda \hslash {\Omega }_{\xi }\sqrt{n+{\Delta }_{\xi }^2}-s\gamma B,& n>0,\\ \xi \hslash {\Omega }_{\xi }{\Delta }_{\xi }-s\gamma B,& n=0.\end{array}\right.\end{array}$$

(77)

Here, $\lambda =\pm$ is the band index, while ${\Omega }_{\xi }=v_F\sqrt{2e|B_{\xi }|/\hslash }$ is the effective frequency, expressed in terms of the effective “total” magnetic field $B_{\xi }=B+\xi B_S$ that results from the combination of the strain-induced pseudo-magnetic field $B_S$ and the real magnetic field B at each valley $K_{\xi }$. We have also defined the dimensionless parameters ${\Delta }_{\xi }=\overline{\gamma }{B}_S{\hslash }^{-1}{\Omega }_{\xi }^{-1}$ at each valley. The two quantum numbers $(n,m)$ correspond to the quantization of the orbital Landau level $n\ge 0$, and the guiding-center (i.e., the center of the classical cyclotronic orbit) coordinate $m\ge 0$, respectively [50]. The energy levels described by Equation (77) are degenerate in the guiding-center quantum number m, with the same degeneracy factor $N_{n,s}^{\xi }=N_{\varphi }^{\xi }$ for each Landau level and each valley. On the other hand, $N_{\varphi }^{\xi }={\Phi }_{B_{\xi }}/{\varphi }_0$ is the number of magnetic flux quanta ${\varphi }_0=h/2e$ piercing the area $\mathcal{A}$ of the graphene flake [50], with ${\Phi }_{B_{\xi }}=B_{\xi }\mathcal{A}$ the “total” flux. As previously discussed, only the magnetic field B couples to the electronic spin $s=\{\pm \}$, as seen in the Zeeman term in Equation (77), which is proportional to $\gamma =g{\mu }_B/2$. Regarding the pseudo Zeeman term, when $B_S$ is expressed in Tesla, we estimate [40,48,52] $\overline{\gamma }=\frac{3\pi a^3}{4\beta }\frac{V^{\prime }}{{\varphi }_0}=9.788\times {10}^{-5}\mathrm{eV}{\mathrm{T}}^{-1}$, for $a=1.42\mathring{A}$ the carbon-carbon bond length, and $V^{\prime }=6{\mathrm{eV}\AA }^{-1}$ [52,54]. Remarkably, $\overline{\gamma }\sim 1.7{\mathsf{\mu }}_B$ is on the order of magnitude of the Bohr magneton.

4.2. The Quantum Engine Cycle

As the “working substance” for a QHE, let us consider a statistical ensemble of replicas of a single-particle system [26,27] consisting of a single electron in the conduction band ($\lambda =+$) of the graphene flake described by Equation (72). This can in principle be achieved by charging an otherwise neutral graphene flake with a positive gate potential. Each replica [26,27] may be in any of the different eigenstates of the Hamiltonian Equation (72). The single-particle system is then in a statistically mixed quantum state [30], described by the density matrix operator $\widehat{\rho }={\sum }_{n,m,s,\xi }{p}_{n,m,s,\xi }\left(B\right)|{\psi }_{n,m,s}^{\left(\xi \right)}\left(B\right)\rangle \langle {\psi }_{n,m,s}^{\left(\xi \right)}\left(B\right)|$, with $|{\psi }_{n,m,s}^{\left(\xi \right)}\left(B\right)\rangle$ a spinor eigenstate of the single-particle Hamiltonian Equation (72) for a given magnetic field intensity B and pseudo-magnetic field $B_S$. The indexes $(n,m,s,\xi )$ enumerate the eigenstates of the Hamiltonian defined by Equation (77), with Landau level n, at valley ξ and electronic spin component s, for $\lambda =+$. Here, the coefficient $0\le p_{n,m,s,\xi }\left(B\right)\le 1$ represents the probability for the system, within the statistical ensemble, to be in the particular state $|{\psi }_{n,m,s}^{\left(\xi \right)}\left(B\right)\rangle$. Therefore, the $\left\{p_{n,m,s,\xi }\left(B\right)\right\}$ satisfy the normalization condition:

$$\begin{array}{c}\hfill \mathrm{Tr}\widehat{\rho }=\sum _{n,m,s,\xi }{p}_{n,m,s,\xi }\left(B\right)=1.\end{array}$$

(78)

As discussed in Section 2, the von Neumann entropy reduces to the explicit expression:

$$\begin{array}{c}\hfill S\left(B\right)=-k_B\sum _{n,m,s,\xi }{p}_{n,m,s,\xi }\left(B\right)ln\left(p_{n,m,s,\xi }\left(B\right)\right).\end{array}$$

(79)

In our notation, we emphasize the explicit dependence of the energy eigenstates $\left\{\right|{\psi }_{n,m,s}^{\left(\xi \right)}\left(B\right)\rangle \}$, as well as the probability coefficients $\left\{p_{n,m,s,\xi }\left(B\right)\right\}$, on the intensity of the external magnetic field B.

The ensemble-average energy of the quantum single-particle system is:

$$\begin{array}{c}\hfill E=\mathrm{Tr}\left(\widehat{\rho }\widehat{H}\right)=\sum _{n,m,s,\xi }{p}_{n,m,s,\xi }\left(B\right)E_{n,s}^{\xi }\left(B\right)=\sum _{n,s,\xi }{p}_{n,s,\xi }\left(B\right)E_{n,s}^{\xi }\left(B\right),\end{array}$$

(80)

where we introduced the coefficients $p_{n,s,\xi }\equiv {\sum }_m{p}_{n,m,s,\xi }$ in order to take notational advantage of the degeneracy in the spectrum with respect to the quantum number m.

The statistical ensemble just described can be submitted to an arbitrary quasi-static process, either by modulating the magnetic field intensity or by exchanging energy with a reservoir. Along such a process, the ensemble-average energy will change accordingly [26,27], in a microscopic version of the first law of thermodynamics for the statistical ensemble of single-particle systems, as defined by Equation (27). Here, from Equation (27), the work along the process connecting states with magnetic fields $B_{\gamma }\to B_{\delta }$ is:

$$\begin{array}{ccc}\hfill W_{\gamma \to \delta }& =& {\int }_{B_{\gamma }}^{B_{\delta }}dB{\left(\frac{\partial E}{\partial B}\right)}_{\left\{p_{n,s,\xi }\left(B\right)\right\}=\mathrm{const}.}={\int }_{B_{\gamma }}^{B_{\delta }}dB{\left(\frac{\partial E}{\partial B}\right)}_S.\hfill \end{array}$$

(81)

On the other hand, the heat exchanged by the system with the environment while modifying its temperature from $T_{\gamma }\to T_{\delta }$ will be:

$$\begin{array}{ccc}\hfill Q_{\gamma \to \delta }& =& {\int }_{T_{\gamma }}^{T_{\delta }}dT{\left(\frac{\partial E}{\partial T}\right)}_B.\hfill \end{array}$$

(82)

For the statistical ensemble just defined, let us consider a cycle by devising a sequence of quasi-static trajectories as depicted in Figure 10. Initially, the single-particle system, while submitted to an external magnetic field of intensity $B_1$, is brought into thermal equilibrium with a macroscopic thermostat at temperature $T_1\equiv T_C$. The system then equilibrates to a Gibbsian ensemble at the temperature of the reservoir ${\beta }_1=\left(k_B{T}_1\right)$:

$$\begin{array}{c}\hfill p_{n,s,\xi }\left(B_1\right)={\left[Z(B_1,{\beta }_1)\right]}^{-1}{N}_{\varphi }^{\xi }{e}^{-{\beta }_1{E}_{n,s}^{\xi }\left(B_1\right)},\end{array}$$

(83)

where the normalization factor is defined by the partition function [28]:

$$\begin{array}{c}\hfill Z(B_1,{\beta }_1)=\sum _{n,s,\xi }{N}_{\varphi }^{\xi }{e}^{-\beta E_{n,s}^{\left(\xi \right)}\left(B_1\right)}\sim e^{-{\beta }_1\overline{\gamma }\xi B_S}+e^{-\beta \hslash {\Omega }_{\xi }\sqrt{1+{\Delta }_{\xi }^2}}.\end{array}$$

(84)

The ensemble-average energy for the Gibbsian distribution in Equation (83) is given by the expression:

$$\begin{array}{ccc}\hfill E({\beta }_1,B_1)& =& -{\left(\frac{\partial lnZ}{\partial \beta }\right)}_{B_1}\hfill \\ & =& -\gamma B_{1}tanh\left({\beta }_1\gamma B_1\right)+\frac{{\sum }_{\xi =\pm }{N}_{\varphi }^{\xi }\left(\xi \overline{\gamma }{B}_S{e}^{-\xi {\beta }_1\overline{\gamma }{B}_S}+\hslash {\Omega }_{\xi }\sqrt{1+{\Delta }_{\xi }^2}{e}^{-{\beta }_1\hslash {\Omega }_{\xi }\sqrt{1+{\Delta }_{\xi }^2}}\right)}{{\sum }_{\xi =\pm }{N}_{\varphi }^{\xi }\left(e^{-\xi {\beta }_1\overline{\gamma }{B}_S}+e^{-{\beta }_1\hslash {\Omega }_{\xi }\sqrt{1+{\Delta }_{\xi }^2}}\right)}\hfill \end{array}$$

(85)

It also follows from the definition Equation (79), along with Equation (83), that the entropy can be expressed in terms of the partition function [28]:

$$\begin{array}{ccc}\hfill S(\beta ,B)/k_B& =& \beta E(\beta ,B)+lnZ(B,\beta )\hfill \\ & =& \beta E(\beta ,B)+ln\left[2cosh\left(\beta \gamma B\right)\right]+ln\left[\sum _{\xi =\pm }{N}_{\varphi }^{\xi }\left(e^{-\xi \beta \overline{\gamma }{B}_S}+e^{-\beta \hslash {\Omega }_{\xi }\sqrt{1+{\Delta }_{\xi }^2}}\right)\right].\hfill \end{array}$$

(86)

The system performs work along the iso-entropic trajectory $1\to 2$, according to Equation (81), $W_{1\to 2}=E(T_2,B_2)-E(T_C,B_1)$, and along the iso-entropic trajectory $3\to 4$, $W_{3\to 4}=E(T_4,B_1)-E(T_H,B_2)$. A physical interpretation of the work performed by the engine is obtained by considering the statistical mechanical definition of the ensemble-average magnetization, that is $M=-{\left(\partial E/\partial B\right)}_S$. Therefore, we clearly have $W=-\int MdB$, a familiar expression from classical macroscopic thermodynamics.

Along the constant magnetic field trajectories $2\to 3$ and $4\to 1$, the system exchanges heat with the reservoirs. The heat absorbed by the system from the hot reservoir at $T_3=T_H$ is:

$$\begin{array}{c}\hfill Q_H={\int }_{T_2}^{T_H}dT{\left(\frac{\partial E}{\partial T}\right)}_{B_2}=E(T_H,B_2)-E(T_2,B_2).\end{array}$$

(87)

Similarly, the heat released by the system to the cold reservoir at $T_C$ is:

$$\begin{array}{c}\hfill Q_C={\int }_{T_4}^{T_C}dT{\left(\frac{\partial E}{\partial T}\right)}_{B_1}=E(T_C,B_1)-E(T_4,B_1).\end{array}$$

(88)

The efficiency of the engine is then given by the expression:

$$\begin{array}{c}\hfill \eta =\left|\frac{W_{1\to 2}+W_{3\to 4}}{Q_H}\right|=1-\left|\frac{Q_C}{Q_H}\right|.\end{array}$$

(89)

The intermediate temperatures $T_2$ and $T_4$ must be determined numerically from the condition that connects the initial and final states along each iso-entropic trajectory (see Figure 10):

$$\begin{array}{ccc}\hfill S(B_1,T_C)& =& S(B_2,T_2),\hfill \\ \hfill S(B_2,T_H)& =& S(B_1,T_4).\hfill \end{array}$$

(90)

For given values of the initial magnetic field $B_1$, the strain pseudo-magnetic field $B_S$ and the reservoir temperatures $T_C$ and $T_H$, the efficiency is a function of the magnetic field $B_2$. We choose to parametrize this dependency by defining the ratio [28]:

$$\begin{array}{c}\hfill r\left(B_2\right)=l_{B_1}/l_{B_2},\end{array}$$

(91)

where $l_B=\min \{l_{+},l_{-}\}$ is a characteristic confinement length for the semi-classical cyclotronic orbit, defined as the minimum Landau radius among the two inequivalent valleys. In Figure 11, we plot the result of our numerical calculation of the efficiency, as a function of the magnetic field expressed in terms of the ratio $r\left(B_2\right)$. In this particular example, we have chosen $B_1=4\mathrm{T}$, $B_S=20\mathrm{T}$ and the temperatures $T_C=30\mathrm{K}$ and $T_H=100\mathrm{K}$ at the cold and hot reservoirs, respectively. In Figure 11, we also compare the effect of the pseudo-Zeeman term, by calculating the efficiency when setting $\overline{\gamma }=1.7{\mu }_B$ and $\overline{\gamma }=0$, respectively. It is evident from the comparison of both curves that the pseudo-Zeeman effect produces a relative enhancement of the efficiency as compared to the case when this term is absent. We find that the numerical solution for the efficiency as a function of the compression ratio can be accurately represented by the parametric form:

$$\begin{array}{c}\hfill \eta =1-{\left[r\left(B_2\right)\right]}^{-\alpha }.\end{array}$$

(92)

Here, the exponent α depends on the temperatures $T_C$ and $T_H$, as well as on the strain field $B_S$. In particular, for the choice of parameters represented in Figure 11, we find that $\alpha =6.88$ for $\overline{\gamma }=1.7{\mathsf{\mu }}_B$, whereas $\alpha =4.2$ for $\overline{\gamma }=0$. Remarkably, this parametric expression is analogous to the well-known formula for the efficiency of the Otto cycle that works with a classical ideal gas, with r instead of the volumetric ratio that applies to the classical case. An even closer analogy between both cases can be put forward by noticing that $r\left(B_2\right)=l_{B_1}/l_{B_2}>1$ can be literally interpreted as a “compression ratio” between the effective Landau radii, which in practice defines a characteristic confinement length for the semi-classical cyclotronic orbit associated with each Landau level.

Numerical solutions for the efficiency only exist up to a maximum compression ratio, which in the example displayed in Figure 11 is $r^{max}=1.19$ for $\overline{\gamma }=1.7{\mu }_B$, whereas $r^{max}=1.33$ for $\overline{\gamma }=0$. At this point, the efficiency attains its maximum value, which exactly matches the Carnot efficiency for the same temperatures in the thermostats, i.e., $\eta \left(r^{max}\right)={\eta }_C=1-T_C/T_H=0.7$. More generally, using the parametric form Equation (92), one concludes that:

$$\begin{array}{c}\hfill r<{\left(\frac{T_H}{T_C}\right)}^{1/\alpha }\equiv r^{max}.\end{array}$$

(93)

5. Discussion

Along this work, we have presented a general theoretical formulation for quantum heat engines (QHEN). We have connected the more general, non-equilibrium thermodynamics aspects of the problem, with a specific formulation proposed by us based on the idealization of a quasi-static limit for the engine operation. Moreover, we have applied this quasi-static analysis to characterize the performance and to explicitly calculate the efficiency of two different QHENs recently proposed by us. The cyclic operation of those QHENs is driven by the modulation of an external magnetic field, as a control parameter for the characteristic confinement length-scale that determines the single-particle spectrum in both cases, defined as the Landau radius. In both examples, it was shown that the quasi-static efficiency of a hypothetical device constructed on this principle can be controlled by adjusting physical parameters, such as the magnetic field, strain and temperatures of the reservoirs that act as the heat source and sink, respectively. Remarkably, the Carnot efficiency imposes an upper limit to the theoretical efficiency obtained for both QHENs analyzed along this work, thus reflecting the robustness of classical thermodynamics.