Multilayer Graphene as an Endoreversible Otto Engine

1. Introduction

There is rapidly growing interest in the development of quantum technologies, devices that take advantage of the unique properties of quantum systems to enhance their performance. This, in turn, has led to increased focus on the field of quantum thermodynamics [1,2,3,4]. Within the broad spectrum of topics that fall under the umbrella of quantum thermodynamics, significant attention is paid to the study of quantum heat engines - devices that extend the principles of classical heat engines to include working mediums made up of quantum systems [5,6,7,8,9].

Efficiency, defined as the ratio of the net work to the heat absorbed from the hot reservoir, is by far the most prominent metric of engine performance. To maximize engine efficiency the strokes of the cycle must be carried out quasistatically. However, truly quasistatic strokes require infinite time to implement, thus leading to vanishing power output. Practically useful metrics of heat engine performance must therefore account for cycles implemented in finite time. Endoreversible thermodynamics [10,11,12] provides a framework for introducing finite-time behavior by assuming that, while the working medium remains in a state of local equilibrium at all times during the cycle, the heating and cooling strokes occur quickly enough that the working medium never fully thermalizes with the hot and cold reservoirs. A prominent performance characteristic within endoreversible thermodynamics is the efficiency at maximum power (EMP) which corresponds to maximizing the power output with respect to the external control parameter and then determining the efficiency at that maximum power output. Endoreversible cycles have also been studied in the context of quantum heat engines. We draw particular attention to Ref. [13], where it was shown that the EMP of an endoreversible Otto cycle with a quantum harmonic oscillator as the working medium exceeds the Curzorn-Albhorn (CA) efficiency, the EMP achieved by the Otto cycle with a classical working medium.

When considering possible systems to serve as the working medium of a quantum heat engine, graphene stands out as intriguing candidate. Graphene’s optical, electronic, and mechanical properties have been extensively studied in recent years [14,15,16,17,18,19,20]. Furthermore, graphene is a prominent Dirac material, systems whose low energy excitations behave as relativistic massless fermions [21]. Thus the study of quantum heat engines with graphene as a working medium can give insight into the role of relativistic quantum features in engine performance [22,23]. In particular, the performance of a quasistatic Otto engine with twisted bilayer graphene was recently studied [24]. In this work it was found that the highest efficiency is reached when the twist angle corresponds to the magic angle of 0.96 degrees. These results show that for a heat engine with multilayer graphene as the working medium the number and configuration of the crystal sheets plays a significant role in the engine performance. Significant attention has also been given to graphene-based engines in the context of continuous, thermoelectric machines [25,26,27,28]. Notably, graphene has also been used in the construction of an experimental nanoscale cyclic heat engine [29].

As in Ref.[24] the calculation is quasi-static, assuming that the bilayer always reaches thermodynamic equilibrium at reservoir temperatures and only for the chosen temperatures ($T_{min}=30$ K and $T_{max}=100$ K) and minimum magnetic field ($B_{min}=5$ T). Therefore in this manuscript we analyze the finite time performance of an Otto cycle with multilayer graphene as the working medium using the framework of endoreversible thermodynamics. Additionally, we explore different temperatures and external magnetic fields to determine the efficiency reached at maximum power for each case, and compare these results with the ideal Courzon-Ahlborn efficiency.

In Section 2 we provide relevant background, including the analytical results of the energy spectrum for monolayer, bilayer, and trilayer graphene. In Section 3 we determine a closed form for the partition function and examine the equilibrium thermodynamic behavior of multilayer graphene. In Section 4 we introduce the endoreversible Otto cycle for multilayer graphene before presenting the results for the engine efficiency, power output, and EMP in Section 5.

2. Model

In multilayer graphene, the crystal sheets are placed on top of each other in different stacking configurations and are connected through weak van der Waals forces. The stacking configurations are determined by the orientation of the two triangular sublattices that make up the primary honeycomb lattice of a single sheet. For two stacked sheets, three possible orientations, A, B, and C, are possible, each corresponding to displacing one of the sublattice atoms along the edge of the honeycomb with respect to the neighboring sheet [30]. Subject to a perpendicular external magnetic field, these systems can be analyzed using a $\pi$-orbital continuum model. Such an analysis is described extensively in Ref. [30]. In our analysis we will focus on two particular stacking configurations. For bilayer graphene we consider Bernal stacking, also known as AB stacking. For the case of trilayer graphene we consider the rhombohedral configuration, also known as ABC stacking [31]. Significantly, an exact analytical result for the energy spectrum as a function of the external magnetic field can be found for these two cases.

2.1. Monolayer Graphene

In monolayer graphene the application of a perpendicular magnetic field results in unevenly spaced Landau levels with an energy spectrum proportional to the root of the level quantum number n [32],

$$E_n=\pm \sqrt{2e\hslash v_f^{2}nB},n=0,1,2,...,$$

(1)

where B is the magnitude of the magnetic field, e is the electron charge, ℏ is Planck’s constant, and $v_f$ is the Fermi velocity ($\sim {10}^6$ m/s). Such an energy spectrum is characteristic of ultra-relativistic massless particles with the Fermi velocity playing the role of the speed of light. The positive energy branch corresponds to particle behavior and the negative energy branch to holes [30]. These energy levels are four times degenerate, including the zero energy state, where the factor of four arises from spin degeneracy and non-equivalent BZ points K and $K^{{}^{\prime }}$, known as valley degeneracy.

2.2. Bilayer Graphene: AB stacking

For a bilayer system, the most stable coupling corresponds to Bernal, or AB, stacking. This consists of displacing the A sublattice atoms of the upper layer so that they lie on top of the B sublattice atoms of the lower layer. Notably, the bilayer system has a quadratic dispersion relation, which gives rise to an interesting phenomenon. While the Dirac equation still models the dynamics of the low energy states, the quadratic dispersion relation indicates that the described charge carriers have mass. In this case, under a perpendicular external magnetic field, the low energy spectrum takes the form [33],

$$E_n=\pm \hslash {\omega }_c\sqrt{n(n-1)},n=0,1,2,..$$

(2)

where ${\omega }_c\equiv eB/m^{*}$ corresponds to the cyclotron frequency. This effective mass, $m^{*}$, is related to the Fermi velocity and the interlayer interaction parameter, $t_{\perp }$, by $m^{*}=t_{\perp }/2v_f^2$. This corresponds to a numerical value of $m^{*}\sim (0.039\pm 0.002)m_e$, where $m_e$ is the electron rest mass. Note that Eq. (2) has two zero energy levels, corresponding to $n=0$ and $n=1$.

2.3. Trilayer Graphene: ABC stacking

Trilayer graphene in the ABC configuration acts as a semiconductor with a gate tunable band gap. The energy spectrum has the form [30],

$$E_n=\pm \frac{{\left(2\hslash v_F^{2}eB\right)}^{3/2}}{t_{\perp }^2}\sqrt{n(n-1)(n-2)},n=0,1,2,...$$

(3)

Note that in the case of trilayer graphene the zero-energy state is 12 fold degenerate, while the other energy states remain fourfold degenerate just as in the bilayer and monolayer case.

3. Partition Function and Equilibrium Thermodynamics

Comparing the energy spectra presented in the previous section, we see a pattern emerging in how the energy scales with the magnetic field. Each energy is proportional to $B^{\mathcal{N}/2}$ where $\mathcal{N}$ is the number of layers. To illustrate this behavior, we plot the first ten positive and negative energy states as a function of the external field for monolayer, bilayer, and trilayer graphene in Fig. Figure 1. We also see the energy spectra follow a common structure in regard to the quantum number n, which takes the form,

$$f_{\mathcal{N}}\left(n\right)=\sqrt{{\displaystyle \prod _{k=0}^{\mathcal{N}-1}\left(n-k\right)}}.$$

(4)

Therefore, we can compactly write the energy spectrum for the multilayer system in the form [30],

$${\mathcal{E}}_{n,\mathcal{N}}={\theta }_{\mathcal{N}}{B}^{\frac{\mathcal{N}}{2}}{f}_{\mathcal{N}}\left(n\right),$$

(5)

where ${\theta }_{\mathcal{N}}\equiv {(2e\hslash v_f^2)}^{\mathcal{N}/2}{\left(t_{\perp }\right)}^{1-\mathcal{N}}$, is a constant that depends on the number of layers and the stacking structure of the system.

Note that the energy spectra in Eqs. (1), (2), and (3) include both positive and negative energy solutions corresponding to electrons and holes, respectively. It is important to highlight that our approximation considers the energy spectra close to the Dirac points of the Brillouin zone, as additional energy solutions require extended approaches not yet fully solved. We focus on the electron excitations only, similar to established analysis for the case of a twisted bilayer graphene system [24]. This regime can be achieved by gating the monolayer graphene [34] and through moderate electron doping for the case of bilayer [35] and trilayer graphene samples [36]. This makes our restriction to the positive energy solutions experimentally available for the three graphene systems we study.

To accurately determine the partition function, we need to carefully consider the degeneracy of the energy levels, especially for the zero-energy state. The compact form of counting these degenerate states in the partition function is given by,

$$\mathcal{Z}=4(\mathcal{N}-1)+\sum _{n=0}^{\infty }4e^{-\beta {\mathcal{E}}_{n,\mathcal{N}}}.$$

(6)

The partition function for the energy spectrum given in Eq. (5) does not have a simple closed-form solution, except for the case of one layer. However, if we assume that the number of layers $\mathcal{N}$ is not very large compared with the number of states n, then the energy spectrum can be approximated as

$${\mathcal{E}}_{n,\mathcal{N}}\approx {\theta }_{\mathcal{N}}{B}^{\frac{\mathcal{N}}{2}}{n}^{\frac{\mathcal{N}}{2}}.$$

(7)

For large n we can approximate the partition function sum as an integral of the form,

$$\mathcal{Z}\approx 4(\mathcal{N}-1)+4{\int }_0^{\infty }dn{e}^{-\beta {\theta }_{\mathcal{N}}{B}^{\frac{\mathcal{N}}{2}}{n}^{\frac{\mathcal{N}}{2}}}.$$

(8)

Noting that,

$${\int }_0^{\infty }dx{e}^{-a{x}^{\frac{\mathcal{N}}{2}}}=a^{-\frac{2}{\mathcal{N}}}\Gamma \left(\frac{2+\mathcal{N}}{\mathcal{N}}\right),$$

(9)

we obtain a simple analytical form for the partition function,

$$\mathcal{Z}\left(T,B,\mathcal{N}\right)=4(\mathcal{N}-1)+4{\left(\frac{{\theta }_{\mathcal{N}}{B}^{\frac{\mathcal{N}}{2}}}{k_{B}T}\right)}^{-\frac{2}{\mathcal{N}}}\Gamma \left(\frac{2+\mathcal{N}}{\mathcal{N}}\right).$$

(10)

From the partition function, all relevant equilibrium thermodynamic properties can be determined as follows,

$$\mathcal{F}=-k_{B}Tln\mathcal{Z},\mathcal{S}=-{\left(\frac{\partial \mathcal{F}}{\partial T}\right)}_B,$$

(11)

$$\mathcal{U}=k_B{T}^2{\left(\frac{\partial ln\mathcal{Z}}{\partial T}\right)}_B,{\mathcal{C}}_B={\left(\frac{\partial \mathcal{U}}{\partial T}\right)}_B,$$

(12)

$$\mathcal{M}=-\left(\frac{\partial \mathcal{F}}{\partial B}\right),$$

(13)

where $\mathcal{F}$ is the free energy, $\mathcal{S}$ is the entropy, $\mathcal{U}$ is the internal energy, ${\mathcal{C}}_B$ is the heat capacity, and $\mathcal{M}$ is the magnetization. In Figs. Figure 2 and Figure 3 we compare the internal energy ($\mathcal{U}$) and entropy ($\mathcal{S}$), respectively, for monolayer, bilayer, and trilayer graphene. To ensure that our analytical approximation for the partition function is valid, we also plot the internal energy and entropy determined from numerical calculations of the partition function sum up to 50,000 terms.

4. The endoreversible Otto cycle

The Otto cycle consists of four strokes, illustrated graphically in Fig. Figure 4 using an entropy (S) - magnetic field (B) diagram. The first stroke ($\mathrm{A}\to \mathrm{B}$) is an isentropic compression in which the external field is varied from $B_1$ to $B_2$ while the working medium is isolated from the thermal reservoirs. During this stroke an amount of work, $W_{\mathrm{comp}}$, must be supplied to compress the working medium. The second stroke ($\mathrm{B}\to \mathrm{C}$) is an isochoric heating stroke in which the working medium draws an amount of heat, $Q_{\mathrm{in}}$, from the hot reservoir while the external field is held constant. The third stroke ($\mathrm{C}\to \mathrm{D}$) is an isentropic expansion where the working medium is again disconnected from the thermal reservoirs and the external field is varied from $B_2$ back to $B_1$. During this stroke an amount of work, $W_{\exp }$, is extracted from the expansion of the working medium. The final stroke ($\mathrm{D}\to \mathrm{A}$) is an isochoric cooling stroke in which the working medium expels and amount of heat, $Q_{\mathrm{out}}$ to the cold reservoir while the external field is held constant. Note that the work parameter (B) plays the role of an inverse volume, increasing during the compression stroke ($\mathrm{A}\to \mathrm{B}$) and decreasing during the expansion stroke ($\mathrm{C}\to \mathrm{D}$).

Characteristic of the framework of endoreversibility, we will assume the working medium remains in a state of local equilibrium throughout the cycle, but, due to finite-time thermalization strokes, never achieves global equilibrium with the reservoirs. The thermodynamic equation of state for the internal energy of working medium at each corner of the cycle must thus be expressed in terms of the corresponding temperature, $T_A$, $T_B$, $T_C$, or $T_D$, and the external field strength, $B_1$ or $B_2$. Ultimately, we want to determine expressions for the engine performance figures of merit solely in terms of the experimentally controllable parameters, namely the temperatures of the thermal reservoirs, $T_l$ and $T_h$, the magnetic field strengths $B_1$ and $B_2$, and the durations of the heating and cooling strokes, ${\tau }_h$ and ${\tau }_l$. In order to do so, we must model the thermal conduction during the isochoric strokes and apply the constraint that the entropy remains constant during the isentropic strokes. For this endoreversible analysis we will follow the procedure established in Ref. [13].

During the isentropic compression stroke ($\mathrm{A}\to \mathrm{B}$) the working medium is decoupled from the thermal reservoirs. As such, all change in the working medium’s internal energy can be associated with work,

$$W_{\mathrm{comp}}=U_{\mathrm{B}}(T_{\mathrm{B}},B_2)-U_{\mathrm{A}}(T_{\mathrm{A}},B_1).$$

(14)

During the isochoric heating stroke ($\mathrm{B}\to \mathrm{C}$), the external field is held constant. Thus the difference in internal energy can be associated entirely with heat,

$$Q_{\mathrm{in}}=U_{\mathrm{C}}(T_{\mathrm{C}},B_2)-U_{\mathrm{B}}(T_{\mathrm{B}},B_2).$$

(15)

As mentioned above, unlike in the quasistatic case, $T_{\mathrm{C}}\ne T_h$ since the working medium does not fully thermalize with the hot reservoir. As the heating stroke is now carried out in finite time, we must determine how the temperature of the working medium changes during the duration of the stroke. The temperatures $T_{\mathrm{B}}$ and $T_{\mathrm{C}}$, corresponding to the temperature of the working medium at the beginning and ending of the heating stroke, respectively, must satisfy the conditions,

$$T\left(0\right)=T_{\mathrm{B}},T\left({\tau }_h\right)=T_{\mathrm{C}}\mathrm{and}{T}_{\mathrm{B}}<T_{\mathrm{C}}\le T_h,$$

(16)

where ${\tau }_h$ is the duration of the heating stroke. Consistent with the assumptions of endoreversibility, we model thermal conduction in the working medium using Fourier’s law. In this case, the temperature change from $T_{\mathrm{B}}$ to $T_{\mathrm{C}}$ can be found by applying Newton’s law of cooling,

$$\frac{dT}{dt}=-{\alpha }_h\left(T\left(t\right)-T_h\right),$$

(17)

where ${\alpha }_h$ is a constant that depends on the thermal conductivity and heat capacity of the working medium. Solving Eq. (17) yields,

$$T_{\mathrm{C}}-T_h=(T_{\mathrm{B}}-T_h)e^{-{\alpha }_h{\tau }_h}.$$

(18)

Just as in the compression stroke, the work extracted during the isentropic expansion stroke ($\mathrm{C}\to \mathrm{D}$) is found from,

$$W_{\exp }=U_{\mathrm{D}}(T_{\mathrm{D}},B_1)-U_{\mathrm{C}}(T_{\mathrm{C}},B_2).$$

(19)

During the isochoric cooling stroke ($\mathrm{D}\to \mathrm{A}$) the heat exchanged with the cold reservoir is given by,

$$Q_{\mathrm{out}}=U_{\mathrm{A}}(T_{\mathrm{A}},B_1)-U_{\mathrm{D}}(T_{\mathrm{D}},B_1),$$

(20)

where, in analogy to the heating stroke, $T_{\mathrm{A}}$ and $T_{\mathrm{D}}$ satisfy the conditions,

$$T\left(0\right)=T_{\mathrm{D}}\mathrm{and}T\left({\tau }_l\right)=T_{\mathrm{A}}\mathrm{with}{T}_{\mathrm{D}}>T_{\mathrm{A}}\ge T_l.$$

(21)

We again apply Fourier’s law and Newton’s law of cooling to model the temperature change during the stroke,

$$\frac{dT}{dt}=-{\alpha }_l\left(T\left(t\right)-T_l\right),$$

(22)

which after solving yields,

$$T_{\mathrm{A}}-T_l=\left(T_{\mathrm{D}}-T_l\right)e^{-{\alpha }_l{\tau }_l}.$$

(23)

With expressions for the work done and heat exchanged during each stroke of the cycle we can now determine the cycle efficiency,

$$\eta =-\frac{W_{\mathrm{comp}}+W_{\exp }}{Q_{\mathrm{in}}},$$

(24)

and power output,

$$P=-\frac{W_{\mathrm{comp}}+W_{\exp }}{\gamma ({\tau }_h+{\tau }_l)}.$$

(25)

Note that $\gamma$ is a multiplicative factor that implicitly incorporates the duration of the isentropic strokes [13].

By definition, the entropy remains constant during the isentropic strokes. We can use this fact to obtain a relationship between the initial and final temperatures and magnetic field strengths during the isentropic strokes. Using $dS(T,B)=0$ we obtain the following first order differential equation,

$$\frac{dB}{dT}=-\frac{{\left(\frac{\partial S}{\partial T}\right)}_B}{{\left(\frac{\partial S}{\partial B}\right)}_T}.$$

(26)

Taking the partial derivatives of the entropy found from Eq. (11) we arrive at,

$$\frac{dB}{dT}=\frac{2B}{\mathcal{N}T}.$$

(27)

Solving Eq. (27) for the compression stroke we find,

$$\frac{T_{\mathrm{A}}}{T_{\mathrm{B}}}={\left(\frac{B_1}{B_2}\right)}^{\frac{\mathcal{N}}{2}}.$$

(28)

Similarly, solving Eq. (27) for the expansion stroke we have,

$$\frac{T_{\mathrm{C}}}{T_{\mathrm{D}}}={\left(\frac{B_2}{B_1}\right)}^{\frac{\mathcal{N}}{2}}.$$

(29)

This relationship between the temperature, external field, and number of layers can be seen graphically in Fig. Figure 5, where we have plotted curves of constant entropy as a function of the temperature and external field for monolayer, bilayer, and trilayer graphene.

5. Results

5.1. Efficiency

We are now in a position where we can determine all of our characterizations of engine performance in terms of experimentally controllable parameters. First, combining Eq. (24) with Eqs. (15), (19), (14), (28) and (29) we arrive at a simple expression for the engine efficiency,

$$\eta =1-{\left(\frac{B_1}{B_2}\right)}^{\frac{\mathcal{N}}{2}}.$$

(30)

We note that this expression is strikingly similar to the classical expression for the Otto efficiency, with the layer number $\mathcal{N}$ playing the role of the ratio of heat capacities.

5.2. Power

Similarly, by combining Eq. (25) with (19), (14), (28) and (29) we arrive at an analytical expression for the power output,

$$\begin{array}{cc}\hfill P=& \frac{2\left(1-{\kappa }^{-\mathcal{N}/2}\right)}{\mathcal{N}\gamma \left({\tau }_l+{\tau }_h\right)}[\frac{\Sigma }{B_2\kappa (\mathcal{N}-1){\Sigma }^{-2/\mathcal{N}}+\Gamma \left(\frac{\mathcal{N}+2}{\mathcal{N}}\right)}\hfill \\ \hfill & -\frac{\Lambda }{B_2\kappa (\mathcal{N}-1){\Lambda }^{-2/\mathcal{N}}+\Gamma \left(\frac{\mathcal{N}+2}{\mathcal{N}}\right)}]\Gamma \left(\frac{\mathcal{N}+2}{\mathcal{N}}\right)\hfill \end{array}$$

(31)

where we have defined,

$$\begin{array}{cc}\hfill & \Sigma \equiv \frac{e^{{\alpha }_h{\tau }_h}\left(e^{{\alpha }_l{\tau }_l}-1\right)k_B{T}_l+\left(e^{{\alpha }_h{\tau }_h}-1\right)k_B{T}_h{\kappa }^{\mathcal{N}/2}}{e^{{\alpha }_l{\tau }_l+{\alpha }_h{\tau }_h}-1},\hfill \\ & \Lambda \equiv \frac{\left(e^{{\alpha }_l{\tau }_l}-1\right)k_B{T}_l+e^{{\alpha }_l{\tau }_l}\left(e^{{\alpha }_h{\tau }_h}-1\right)k_B{T}_h{\kappa }^{\mathcal{N}/2}}{e^{{\alpha }_l{\tau }_l+{\alpha }_h{\tau }_h}-1}.\hfill \end{array}$$

(32)

Examining Eq. (31), we see that the power will vanish under the condition that $\Sigma =\Lambda$. This occurs under three conditions. The first is that ${\kappa }^{\mathcal{N}/2}\to T_l/T_h$. We see from Eq. 30 that this corresponds to limit of Carnot efficiency, under which we would expect the power to vanish. The second and third conditions are when $exp\left({\alpha }_h{\tau }_h\right)\to 1$ and $exp\left({\alpha }_l{\tau }_l\right)\to 1$, respectively. These conditions correspond to the limits of instantaneous thermalization strokes or vanishing thermal conductivity, both of which would prevent heat transfer and thus result in zero power. We also note that the power vanishes in the quasistatic limit of ${\tau }_l+{\tau }_h\to \infty$. From Eqs. (18) and (23) we see that this limit yields $T_3=T_h$ and $T_1=T_l$, which in turn maximizes the efficiency. This is a demonstration of the well-established trade-off between efficiency and power.

The efficiency and power are plotted as a function of the compression ratio, $\kappa$, in Fig. Figure 6 for the monolayer, bilayer, and trilayer systems. We see that the monolayer system has the lowest efficiency, but highest power output while the opposite is true for the trilayer system. The efficiency and power of the bilayer system falls between the monolayer and trilayer results.

5.3. Efficiency at maximum power

Due to the inherent trade-off between efficiency and power mentioned above, efficiency alone does not provide the most practically useful metric of engine performance. Instead, this role is played by the efficiency at maximum power. In this case the EMP is found by maximizing the power output with respect to the compression ratio, $\kappa$, and then determining the efficiency corresponding to this value of $\kappa$. For a classical Otto cycle, the EMP is given by the Curzon-Ahlborn efficiency,

$${\eta }_{\mathrm{CA}}=1-\sqrt{\frac{T_l}{T_h}}$$

(33)

Due to the complicated expression for power in Eq. (31) we maximize the power numerically. The EMP as a function of the ratio of bath temperatures is shown in Fig. Figure 7. We see that for the monolayer case, the EMP is identical to the Curzon-Ahlborn efficiency. This result can be confirmed analytically by taking the derivative of Eq. (31) and confirming that it vanishes for $\mathcal{N}=1$ and $\kappa =T_l/T_h$.

For the bilayer and trilayer systems, however, we see that the EMP exceeds the CA efficiency. The EMP is largest for the bilayer system, decreasing slightly in the trilayer case. This trend continues, with the EMP of larger layer numbers converging back towards the CA efficiency. However, it is important to note that if we increase the number of graphene layers significantly beyond the trilayer case, the assumptions made in determining the closed form of the partition function begin to break down.

It has been previously shown for both classical and quantum working mediums that, within the regime of linear response, EMP is bounded by the CA efficiency [37,38,39,40,41,42]. To achieve higher EMP requires going beyond the linear regime or by breaking time-reversal symmetry [43,44,45,46,47,48,49,50,51]. For a cyclic engine, the regime of linear response occurs near the equilibrium limit $T_l\approx T_h$.

To probe the behavior of the EMP for a multilayer graphene working medium in and around the linear response regime we define $T_l\equiv T$ and $T_h\equiv ϵT$. In Fig. Figure 8 we plot the EMP for the monolayer, bilayer and trilayer working mediums in comparison to the CA efficiency for $ϵ=1.1$, 2, and 10. As expected, for the monolayer system we see that at all examined values of $ϵ$ the EMP is identical to the CA efficiency. For the bilayer and trilayer working mediums we see that at $ϵ=1.1$, close to the equilibrium limit, the EMP is identical to CA, consistent with the results in the works mentioned above. As we move away from the equilibrium limit by increasing $ϵ$ we see that, at low bath temperatures, the EMP exceeds CA but that as the temperature increases the EMP converges back to the CA efficiency. As observed in Fig. Figure 7, in the low temperature regime the bilayer working medium EMP exceeds the CA efficiency by a greater amount than the trilayer working medium. However, as temperature increases the bilayer EMP converges to CA faster than the trilayer EMP.

From these results we see that two conditions must be met for the EMP to exceed CA. First, the difference in bath temperatures must be sufficiently far from the equilibrium limit of $T_l\approx T_h$. Second, the cycle must be operating in the low-temperature, quantum regime which for the multilayer graphene working medium is determined by the condition ${\theta }_{\mathcal{N}}{B}_2/k_B{T}_l\gg 1$. This second condition is consistent with the results shown in Ref. [13] for harmonic working mediums.

5.4. Engine vs refrigerator

For any arbitrary choice of parameters it is not guaranteed that the Otto cycle will function as an engine. In general, there are four possible types of thermal machines, corresponding to all possible combinations of directions heat and work flow consistent with the first and second laws of thermodynamics. An engine corresponds to positive work output, along with heat flow from the hot bath into the working medium, and from the working medium into the cold bath. A refrigerator corresponds to negative work output, along with heat flow from the cold bath into the working medium and from the working medium into the hot bath. A heater corresponds to negative work output, along with heat flow from the working medium into both baths. Finally, an accelerator corresponds to negative work output along with heat flow from the hot bath into the working medium and from the working medium into the cold bath.

By examining the signs of Eqs. (15), (19), (20), and (14) across the parameter space we can determine the regions where cycle will function as each type of thermal machine. Note, both the heater and accelerator are fundamentally nonequilibrium devices and thus we would not expect to find regions of parameter space corresponding to these devices under the assumption of endoreversible behavior.

In Fig. Figure 9 we show the regions where the cycle functions as either an engine or a refrigerator as a function of the hot bath temperature and compression ratio. We see that over the same region of parameter space the layer number has a significant impact on the apportionment between engine and refrigerator. In the monolayer case we see that majority of the examined region corresponds to the engine regime, while in the trilayer case the opposite is true, with a larger portion of the explored space corresponding to the refrigerator regime. The origins of this behavior can be understood from the plot of the power in Fig. Figure 6. As the layer number decreases we see that the power reduced power output leads to a smaller region where the work is positive and thus a reduced engine regime.

6. Conclusions

We have examined the performance of an endoreversible Otto cycle with a working medium of a multilayer graphene system. We have found that all examined performance metrics, including the efficiency, power, EMP and parameter regions under which the cycle functions as an engine or refrigerator all depend significantly on the number of layers. Most notably we have found that the EMP for bilayer and trilayer graphene working mediums exceeds the Curzon-Ahlborn efficiency. Conversely, we have demonstrated that the EMP of a monolayer graphene working medium is identical to the CA efficiency. As the energy spectrum of monolayer graphene can be mapped to that of the relativistic Dirac oscillator, this result is consistent with previous work examining the performance of endoreversible Otto engines with relativistic oscillators as the working medium [23]. We have also found two conditions necessary for the graphene EMP to exceed the CA efficiency, namely the cycle must be sufficiently far from the equilibrium limit and must be operating in the low temperature regime corresponding to ${\theta }_{\mathcal{N}}{B}_2/k_B{T}_l\gg 1$.

Experimental implementation of the multilayer graphene engine requires precise control of the number of layers as well as a tunable magnetic field of sufficient strength to induce Landau quantization. Fine control over multilayer structures has been demonstrated by folding monolayer graphene nanoribbons [52] and precise multilayer thickness measurements can be accomplished with electron spectroscopy [53]. Strong external magnetic fields can be generated either by direct application [54] or through strain induced psuedo-magnetic fields [55,56]. Future work involves examining the role of holes in the thermodynamics within the temperature range where the Dirac model for graphene remains valid, by utilizing Fermi-Dirac statistics.