\section{Introduction}
Quantum thermodynamics is concerned with how manifestly qfuantum phenomena, such as coherence and entanglement, affect thermal processes. A primary objective \cite{2016Goold143001,2016Vinjanampathy545} is understanding how these quantum phenomena can be harnessed to improve thermal processes, for example to increase work extraction from a heat engine. Conversely, what limitations does quantum mechanics place on thermal processes?
There has been research investigating the work extracted from a quantum measurement \cite{Erez2012,Kammerlander2016,Jacobs2012,Ding2018,Yi2017}, or measurements with feedback \citep{Funo2013,Maruyama2009,Elouard2017,Elouard2018}, known as a Maxwell's demon. Systems initialised in non-equilibrium states allow for more work extraction according to generalisations of the second law \cite{Hasegawa2010,Niedenzu2018}. Correlations between systems can also increase the extractable work \cite{Oppenheim2002,Perarnau-Llobet2015,Zurek2003,Jevtic2012,Rio2011}, and there is a thermodynamic cost to creating correlations in uncorrelated, thermal systems \cite{Esposito2010}.

Of particular interest is the effect of coherence in thermal processes. Although it contributes to the free energy, coherence cannot be extracted as work using thermal operations \cite{Lostaglio2015,Korzekwa2016}, but can be used as a catalyst \cite{Aberg2014} and thought of as a quantum resource \cite{Baumgratz2014,Lostaglio2015X,Marvian2016}.
Instead of considering coherence in the system itself Scully et al. \cite{2003scully862} studied a system which interacted with a coherent heat bath. Specifically they considered a photo-Carnot engine ---  a single mode cavity undergoing a Carnot cycle --- in which a heat bath is emulated by repeated weak interactions with identical three-level atoms with two degenerate ground states. Coherence between the two ground states can change the thermalisation temperature of the cavity, so a Carnot cycle can be implemented using a single bath and switching on coherence. There has been extensions to this setup, including adding cavity leakage \cite{Quan2006} and using entangled qubits \cite{Dillenschneider2009}, or more general multi-atom baths \cite{Turkpence2017,Dag2016,Dag2018,Hardal2015}.

In the quasistatic limit, the efficiency is maximal, but since this requires infinite cycle duration, the output power is zero and the engine is impractical; one must trade off engine efficiency for output power \cite{Curzon1975}. We investigate the effect on the power and efficiency of a quantum harmonic oscillator undergoing an Otto cycle with a coherent heat bath. Following Scully et al. we emulate a heat bath using repeated weak interactions with identical three-level systems with thermal populations and coherence between the two degenerate energy ground states. The interaction between the three-level atoms and the system is governed by a Jaynes-Cummings Hamiltonian, from which we derive a master equation for the state of the system.  This allows us to analyse the performance of the engine in a more practical regime; without requiring the duration of each stroke to be infinite.

The two isentropic strokes involve changing the harmonic frequency $\omega$. If this occurs rapidly then many off-diagonal coherences will be excited in the system. To correct for this effect we use shortcut-to-adiabaticity techniques (see \cite{ATOMTorrontegui2013} for a review). These techniques allow the state to undergo effectively adiabatic evolution in a finite time, and can be used to drive a thermal engine \cite{Campo2015,Beau2016,Abah2018,Abah2018PRE}. In \cite{EPLLutz2017}, the authors calculate a cost which must be paid to perform their shortcut-to-adiabaticity protocol, which increases with shorter stroke times. We adopt this cost as the reduction in work extracted during the Otto cycle.

This paper has the following structure. In section~\ref{sec:interaction} we describe the interaction between the bath and the system, showing the change in the thermalisation temperature of the system when there is coherence between the ground states of the atoms. In section~\ref{sec:Otto} we describe the four strokes of the Otto cycle, and in section~\ref{sec:finite} we describe some issues with running an engine in finite time. Section~\ref{sec:results} contains our results of the effect of coherence on the power and efficiency of the heat engine in the finite time regime; and in section~\ref{sec:optimise} we optimise the output power of the heat engine using coherence and by changing the relative duration of each stroke.

\section{Interaction with Bath}\label{sec:interaction}

The working fluid of the heat engine is a quantum harmonic oscillator with unit mass, described by the Hamiltonian
\begin{equation}
H_{S} = \frac{1}{2}\left(p^{2}+\omega^{2}x^{2}\right) = \hbar \omega \left( a^\dagger a + \frac{1}{2}\right).
\end{equation}
where $a^{\dagger}a$ is the usual photon number operator.
A sequence of $(\ell+1)$-level atoms interact for a short time with the system, and after many weak interactions this stream of atoms emulates a
thermal bath as it causes the system to thermalise. Here we leave the number of atomic energy levels general, but we will shortly assume that $\ell=2$. The atoms have $\ell$ degenerate ground states $\ket{g_j}$ and a single excited
state $\ket{e}$ with an energy level spacing tuned to the system frequency $\omega$,
\begin{equation}
H_{R} = \hbar \omega \ketbra{e}{e}.
\end{equation}
In this paper we consider loss to the environment to be negligible and focus on the interaction with the atoms as the dominant process.

The interaction between a single atom and the system is described by the familiar Hamiltonian in the rotating wave approximation
\begin{equation}\label{eq:Hi}
  H_{I} =  \hbar \Omega \left(\sum_{j=1}^\ell a^\dagger  \ketbra{g_j}{e} + a  \ketbra{e}{g_j}\right),
\end{equation}
where $\Omega$ is the interaction strength. Introducing \mbox{$J_-=\ketbra{G}{e}$} and $J_+=J_-^\dagger$, where \mbox{$\ket{G}=\sum_j  \ket{g_j}/\sqrt{\ell}$},
we can write the interaction Hamiltonian as
\mbox{$H_\mathrm{I} =  \hbar \sqrt{\ell}\Omega (a^\dagger  J_-+ a  J_+)$}. In the interaction picture with zero detuning, the total Hamiltonian reduces simply to the interaction Hamiltonian $H_{I}$. Thus the interaction with duration $\tau$ is governed by the unitary
\begin{align*}
  U &= e^{-i H_{\mathrm{I}} \tau/\hbar}=e^{-i\lambda\tau(a^\dagger  J_-+ a  J_+)} \\
  &\approx I - i\lambda\tau (a^\dagger  J_-+ a  J_+) -\frac{\lambda^2\tau^{2}}{2}(a^\dagger a J_-J_+ + aa^\dagger J_+J_-)
\end{align*}
where we have expanded to second order in $\lambda\tau=\sqrt{\ell}\Omega \tau$ which is assumed to be small.

Assuming the atom starts in the initial state $\rho_R$ at some time $t'$, each interaction induces a completely positive trace preserving (CPTP) map on the system-atom product state $\rho\otimes \rho_{R}$:
\begin{align}
	&\rho\left(t^{\prime}+\tau\right) = \mathcal{E}(\rho\left(t^{\prime}\right)) =  \mathrm{Tr}_R \{U (\rho\left(t^{\prime}\right)\otimes\rho_R) U^\dagger\} \label{eq:cptp-map}\\
	&\approx \rho\left(t^{\prime}\right) + \lambda^2 \tau^{2} (\bra{e}\rho_R\ket{e} \mathcal{D}[a^\dagger] +
	\bra{G}\rho_R\ket{G} \mathcal{D}[a])\rho\left(t^{\prime}\right), \nonumber
\end{align}
where $\mathcal{D}[b]\rho = b\rho b^\dagger -1/2 b^\dagger b \rho -1/2 \rho b^\dagger b$ is the Lindblad superoperator. In the derivation of \eqref{eq:cptp-map} we assumed that the atoms were block-diagonal in the energy basis, consequently
\begin{align}\label{eq:meanzero}
    \bra{e}\rho_{R}\ket{G}=0,
\end{align}
which removes the term of order $\lambda\tau$. It is straightforward to show that this map contains thermal states as fixed points.

To derive a master equation, we first note that since $\lambda^{2}\tau^{2}$ is assumed small, the map \eqref{eq:cptp-map} is equivalent to $N$ applications of a weaker imaginary map of time duration $\tau_N = \tau/N$:
\begin{multline*}
\rho\left(t^{\prime}+\tau\right) = \\
\left(I +\lambda^2 \tau\tau_N \left(\bra{e}\rho_R\ket{e} \mathcal{D}[a^\dagger] + \bra{G}\rho_R\ket{G} \mathcal{D}[a]\right)\right)^{N}\rho\left(t^{\prime}\right).
\end{multline*}
Considering a single step of this weaker map we can go to the continuum limit
\begin{equation}
\frac{d\rho(t^{\prime})}{dt^{\prime}} = \lim_{\tau_N\rightarrow 0} \frac{\rho\left(t^{\prime}+\tau_N\right)-\rho\left(t^{\prime}\right)}{\tau_N},
\end{equation}
which yields a dissipative Lindblad master equation
\begin{equation}
	\frac{d\rho(t)}{dt} = \bra{e}\rho_R\ket{e} \mathcal{D}[a^\dagger] \rho(t) +
	\bra{G}\rho_R\ket{G} \mathcal{D}[a] \rho(t), \label{eq:me}
\end{equation}
where for simplification we have rescaled time to the dimensionless parameter \mbox{$t=\lambda^2 \tau t^{\prime}$}.

If the initial state of the system is diagonal in photon number, the steady state of the master equation satisfies
\begin{align}
	\bra{n}\rho\ket{n} = \frac{(E/G)^n}{\sum_n (E/G)^n} = \frac{G-E}{G}  \left( \frac{E}{G} \right)^n ,
\end{align}
where the labeling $E=\bra{e}\rho_R\ket{e}$ and $G=\bra{G}\rho_R\ket{G}$ has been introduced to simplify the expressions.

If the atoms are prepared in a thermal state at inverse temperature $\beta_R$ then
\begin{align}
\rho_R = p_e \ketbra{e}{e} + \frac{1-p_e}{\ell}\sum_j \ketbra{g_j}{g_j},
\end{align}
where
\begin{equation}
p_e = \exp(-\hbar \omega \beta_R)/Z, \label{eq:pe}
\end{equation}
and the partition function is $Z = \ell+\exp(-\hbar \omega \beta_R)$. In this case $E=p_e$ and $G=(1-p_e)/\ell = 1/Z$, and the steady state of the master equation~\eqref{eq:me} is a thermal distribution at inverse temperature $\beta = \beta_R$,
\begin{align}
\rho \rightarrow \gamma_\beta = \frac{e^{-\beta \hbar a^\dagger a}}{\mathrm{Tr}(e^{-\beta \hbar a^\dagger a})}.
\end{align}
So the atoms, despite being a rather structured interaction, are emulating a thermal bath at a fixed temperature to which the system thermalises. It can be shown that if the system begins in a thermal state, then it remains in a thermal state under the evolution of this master equation, equilibrating to a particular temperature in the steady state.

In a thermodynamic analysis we are less interested in the explicit density matrix, than the internal energy. Changes in average energy are proportional to changes in average photon number, $\Delta \bar{H}_{S} = \hbar\omega \Delta \bar{n}$.
Master equation~(\ref{eq:me}) allows us to derive a time evolution equation for the mean photon number $\bar{n}=\mathrm{Tr}(a^\dagger a \rho)$,
\begin{align}
	\frac{d\bar{n}}{dt} = E(\bar{n}+1)-G\bar{n}. \label{eq:dndt}
\end{align}
This equation has the solution
\begin{align}
	\bar{n}(t) = \left( \bar{n}(0) -\frac{E}{G-E}\right)e^{-(G-E)t}+ \frac{E}{G-E}. \label{eq:nt}
\end{align}
The steady state limit exists provided $G > E$ and is
\begin{align}
\bar{n}_\mathrm{ss} = \lim_{t\rightarrow \infty} \bar{n}\left(t\right) = \frac{E}{G-E}.
\end{align}
If $G\leq E$ then the average photon number will increase unbounded and the system will not thermalise. While this can be achieved using coherence, the atoms cease to behave similarly to a thermal bath and so we exclude this case.\\
For simplicity, from here onwards we will use \mbox{$\Delta \equiv G-E$}.

In the steady state limit, the harmonic oscillator is in a thermal state, and the steady state average photon number is related to the inverse temperature as
\begin{equation}
\bar{n}_{\mathrm{ss}} = \frac{1}{e^{\beta\hbar\omega}-1},
\end{equation}
and therefore system will thermalise to the temperature
\begin{equation}
	\beta = \frac{1}{\hbar\omega}\ln\left(\frac{\bar{n}_{\mathrm{ss}}+1}{\bar{n}_{\mathrm{ss}}}\right) = \frac{1}{\hbar\omega}\ln\left(\frac{G}{E}\right).
\end{equation}
The temperature is  strictly monotonic in $\bar{n}_{\mathrm{ss}}$, so we can consider $\bar{n}_{\mathrm{ss}}=E/(G-E)$ as the effective temperature of the atoms. Hence in section~\ref{sec:results} when we compare different pairs of thermal baths with the same effective temperatures, we are holding $\bar{n}_{\mathrm{ss}}$ constant.

The value of $G=\bra{G}\rho_{R}\ket{G}$ (and hence $\bar{n}_{\mathrm{ss}}$ and $\beta$ too) is altered by coherence in the energy ground space of the atoms. In general the atoms are block-diagonal in energy
\begin{equation}
\rho_{R} = p_{e}\ketbra{e}{e} + (1-p_{e})\rho_{g},
\end{equation}
where $p_{e}$ is defined as it was previously \eqref{eq:pe}, and $\rho_{g}$ is the quantum state in the degenerate ground subspace. Hence \mbox{$0\leq G =\left(1-p_{e}\right)\bra{G}\rho_{g}\ket{G}\leq 1-p_{e}$}. Coherence in the energy ground space can be chosen to increase or decrease the thermalisation temperature $\beta$ of the system
\begin{equation}
\beta = \beta_{R} + \frac{1}{\hbar\omega}\ln\left(\ell \bra{G}\rho_{g}\ket{G}\right).
\end{equation}
In the case of no coherence then $\bra{G}\rho_{g}\ket{G} = \ell^{-1}$ and $\beta=\beta_{R}$. Coherence can be used to make the effective temperature arbitrarily hot; as $\bra{G}\rho_{g}\ket{G}$ approaches $\left(e^{-\beta_{R}\hbar\omega}\right)/\ell$ from above, the thermalisation temperature tends to infinity. However since $\bra{G}\rho_{g}\ket{G}\leq 1$, the effective temperature can only be lowered as cold as
\begin{equation}
\beta = \beta_{R} + \frac{1}{\hbar\omega}\ln(\ell).\label{eq:coollimit}
\end{equation}
A corollary is that when using maximal coherence to lower the effective temperature of a thermal bath, the final temperature is bounded from above
\begin{equation}
T \leq \frac{\hbar\omega}{\ln\left(\ell\right)}.
\end{equation}

It may seem paradoxical that the coherence in the ground state can raise the thermalisation temperature of the system arbitrarily; however, this is not a mysterious consequence of quantum coherence but simply due to the nature of the coupling between the atoms and field \eqref{eq:Hi}. By tuning the quantum state in the ground subspace $\rho_{g}$ we can prevent the system from releasing energy to the atoms and so each interaction can only increase the energy of the system.

For the rest of this paper we will specialise to the case in which $\ell=2$.

\section{The Otto cycle}\label{sec:Otto}

\begin{figure}[t]
    \centering
    \includegraphics[width=0.45\textwidth]{Otto_cycle.pdf}
    \caption{The Otto cycle shown for the variables $\bar{n}$ and $\omega$, and as $T$-$S$ diagram (inset). Strokes described in SEC.~\ref{sec:Otto}. Isotherms at temperatures $T_{h}$ (red) and $T_{c}$ (blue) are shown. The work extracted in a single cycle is the area enclosed by the rectangle, $W = \hbar\left(\omega_{h}-\omega_{c}\right)\left(\bar{n}_{h}-\bar{n}_{c}\right)$. The heat input to the system is the energy change of the system during the isochoric heating stroke, $Q_{h} = E_{1} = \hbar\omega_{h}\left(\bar{n}_{h}-\bar{n}_{c}\right)$.}
    \label{fig:ottocycle}
\end{figure}
In this section we examine the heat engine under the quasistatic evolution of the Otto  thermal cycle.
The Otto cycle (see FIG.~\ref{fig:ottocycle}) involves the following four strokes:
\begin{enumerate}
    \item  \emph{Isochoric heating}: the system interacts with the hot bath, eventually thermalising to temperature $T_h=\beta_{h}^{-1}$.
    \item  \emph{Isentropic expansion}: the system is isolated from the heat baths and the system frequency is decreased $\omega_{h}\rightarrow\omega_{c}$. Work is extracted from the system.
    \item  \emph{Isochoric cooling}: the system interacts with the cold bath, eventually thermalising to temperature $T_c=\beta_{c}^{-1}$.
    \item  \emph{Isentropic compression}: the system is isolated from the heat baths and the system frequency is increased $\omega_{c}\rightarrow\omega_{h}$. Work is performed on the system to close the cycle.
\end{enumerate}
The Otto cycle allows the easy quantification of heat and work for a cycle. Since the isochoric (constant volume) stokes involve no change of the system Hamiltonian, any change in energy is due to a heat flow between the system and thermal baths, $\Delta Q = \Delta E_{1}+\Delta E_{3}$. During the isentropic (constant entropy) strokes, the system is isolated, and so any change in energy is due to work done on the system, $-\Delta W = \Delta E_{2}+\Delta E_{4}$. Compare this to the Carnot cycle, where work is done on the system during every stroke. The internal energy is $U = \bar{H}_{S} = \hbar\omega \left( \bar{n} + \frac{1}{2} \right)$ and hence
\begin{equation}
    dU = \underbrace{\frac{dU}{d\omega}d\omega}_{-dW}+\underbrace{\frac{dU}{d\bar{n}}d\bar{n}}_{dQ},
\end{equation}
which is a statement of the first law for this system. During a quasi-static evolution the system remains in a thermal state. The entropy can be found by integrating
\begin{figure}[t]
    \centering
    \includegraphics[scale=0.4]{expansion_vs_compression.pdf}
    \caption{The energy cost of the isentropic expansion $V_{e}$ and compression $V_{c}$ strokes, using the shortcut-to-adiabaticity protocol. The energy cost is scaled by $1/\hbar\omega$, where $\omega = 10^{14}$Hz. Recall that $t_{\text{cycle}}$ is in terms of the dimensionless time parameter \eqref{eq:me}. The cost is large for short cycle times and approaches zero as they become longer in duration. In this plot we used a hot bath with temperature $\bar{n}_{\mathrm{ss}}=2$, and a cold bath with $\bar{n}_{\mathrm{ss}}=0.55$, both of which were incoherent. The system interacted with the hot bath for duration $t_{h}$ and with the cold bath for $t_{c}$; and we chose $t_{h}=t_{c}=t_{\text{cycle}}/4$. Expansion carries a larger cost than compression because it begins with the system at a hotter temperature, and $V\left(\tau\right)\propto \bar{n}\left(0\right)$ \eqref{eq:cost}. The compression stroke begins with the steady cycle initial photon number $\bar{n}_{sc}$ \eqref{eq:steadycycle}; the expansion stroke begins with initial photon number $\bar{n}_{h}$, which is found by evolving $\bar{n}_{sc}$ under \eqref{eq:nt} for time $t_{h}$. The integration in \eqref{eq:cost} was completed numerically.}
    \label{fig:expand_compress}
\end{figure}
\begin{equation}
    d S = \beta dQ = \ln\left(\frac{\bar{n}+1}{\bar{n}}\right) d\bar{n}
\end{equation}
from which we find
\begin{equation}
    S = \left(\bar{n}+1\right)\ln\left(\bar{n}+1\right)-\bar{n}\ln\left(\bar{n}\right).
\end{equation}
The efficiency of the Otto cycle is the ratio of the work extracted to the heat input from the hot bath
\begin{equation}
    \eta  \equiv\frac{W}{Q_{h}}=\frac{Q_{h}+Q_{c}}{Q_{h}}.
\end{equation}
We can easily read off from FIG.~\ref{fig:ottocycle} that the work done in a cycle is the area $W=\hbar\left(\omega_h-\omega_c\right)\left(\bar{n}_{h}-\bar{n}_{c}\right)$, and therefore
\begin{equation}
    \eta = \frac{\hbar\left(\omega_{h}-\omega_{c}\right)\left(\bar{n}_{h}-\bar{n}_{c}\right)}{\hbar\omega_{h}\left(\bar{n}_{h}-\bar{n}_{c}\right)}=1-\frac{\omega_{c}}{\omega_{h}}.
\end{equation}
We can relate this to the temperature (see FIG.~\ref{fig:ottocycle}), since during the work strokes, the product $\beta\omega$ is constant,
\begin{equation}
    \eta = 1-\frac{T_{2}}{T_{h}} = 1-\frac{T_{c}}{T_{1}}.
\end{equation}

\section{The Finite Time Regime}\label{sec:finite}

Previous works \cite{2003scully862,Quan2006,Dillenschneider2009,Turkpence2016} have considered a harmonic oscillator undergoing a Carnot or Otto cycle in the quasistatic limit. However the quasistatic limit is impractical as it produces zero output power; so we examine the harmonic oscillator undergoing a finite time Otto cycle.
In this section we discuss two issues in operating a heat engine for finite time. During the work strokes, we use shortcut-to-adiabaticity techniques to control the system energy populations; and during the heat strokes, we need to guarantee that the average energy returns to its original value at the end of each cycle.

\subsection{Work Strokes}

During the work strokes, the system evolves under the time-dependent Hamiltonian \cite{EPLLutz2017,CHEMDeff2010,PTPHusimi1953}
\begin{equation}
    H\left(t\right) = \frac{1}{2}\left(\omega\left(t\right)^{2}x^{2}+p^{2}\right).
\end{equation}
If the system frequency is changed very slowly, then by the adiabatic theorem, the state of the system is constant. Quickly changing the frequency will excite off-diagonal terms in the density matrix. Using \emph{shortcuts-to-adiabaticity} (STA) techniques \cite{ATOMTorrontegui2013}, the state of the system can be constant even when quickly changing the frequency (this has recently been demonstrated experimentally \cite{Deng2018,Diao2018}). This can be achieved by adding a \emph{counterdiabatic term} to the Hamiltonian, $H\left(t\right)\rightarrow H\left(t\right)+H_{\text{STA}}\left(t\right)$. For example, in \cite{DelCampo2013}, the author uses the following counterdiabatic term
\begin{equation}
    H_{\text{STA}}\left(t\right)=\frac{1}{2}\left(\frac{\ddot{\omega}}{2\omega}-\frac{3\dot{\omega}^{2}}{4\omega^{2}}\right)x^{2}.\label{eq:conterdiabatic}
\end{equation}
The author changes the system frequency according to
\begin{align}
    \omega\left(t\right) =& \;\omega_{i} +10\left(\omega_{f}-\omega_{i}\right)\left(\frac{t}{\tilde{\tau}}\right)^{3}-15\left(\omega_{f}-\omega_{i}\right)\left(\frac{t}{\tilde{\tau}}\right)^{4}\nonumber\\
    &+6\left(\omega_{f}-\omega_{i}\right)\left(\frac{t}{\tilde{\tau}}\right)^{5}.
\end{align}
This protocol will change the frequency from $\omega_{i}\rightarrow\omega_{f}$ in time $\tilde{\tau}$. The author shows that if the system begins in a thermal state, then under this counterdiabatic Hamiltonian \eqref{eq:conterdiabatic}, the energy populations will be constant. Such a protocol can be used to drive an thermal cycle \cite{Campo2015,Beau2016,Abah2018,Abah2018PRE}; we use this protocol as well.  The authors in \cite{EPLLutz2017} describe an energy cost associated with the addition of this counterdiabatic term, given by
\begin{align}
   V\left(\tilde{\tau}\right) &=  \frac{1}{\tilde{\tau}}\int_{0}^{\tilde{\tau}} \expect{H_{\text{STA}}\left(t\right)}dt\nonumber \\
   &=\frac{\hbar \bar{n}\left(0\right)}{\tilde{\tau}}\int_{0}^{\tilde{\tau}}\omega\left(\frac{\ddot{\omega}}{4\omega^{3}}-\frac{\dot{\omega}^{2}}{4\omega^{4}}\right)dt, \label{eq:cost}
\end{align}
where $\bar{n}\left(0\right)$ is the average photon number at $t=0$; which implies the isentropic expansion strokes will have a larger energy cost than the compression strokes (see FIG.~\ref{fig:expand_compress}; in all our plots, the integration was completed numerically). The energy cost is monotonically decreasing in the stroke duration. We note here that while we adopt the energy cost described in \cite{EPLLutz2017}, there exists a diversity of opinion about what constitutes the appropriate energy cost of running a shortcut-to-adiabaticity protocol \cite{Abah2018,Kosloff2017,Torrontegui2017,Zheng2016,Calzetta2018}.

The efficiency of the heat engine will therefore be diminished by two energy costs $V_{e}$ and $V_{c}$, one each for the isentropic expansion and compression strokes
\begin{equation}
    \eta = \frac{W-V_{e}-V_{c}}{Q_{h}}. \label{eq:efficiencywithcost}
\end{equation}

The authors of \cite{EPLLutz2017} consider the cost $V_{e}+V_{c}$ as an additional heat input into the system. We disagree with this. This energy cost is due to a unitary process when the system is isolated from other systems. Therefore we consider the energy cost as a reduction in the work extracted from the system.

Using this counterdiabatic term the energy populations of the system are constant, and therefore the average photon number $\bar{n}$ of the system is constant. We assume that the change in energy of the system during the isentropic expansion stroke $\Delta\bar{H}_{S}=\hbar\left(\omega_{h}-\omega_{c}\right)\bar{n}_{h}$ is entirely extracted as work. Hence, the net extracted work is proportional to the change in the frequency of the system $\omega_{h} - \omega_{c}$.

\subsection{Heat Strokes}

\begin{figure}[t]
    \centering
    \includegraphics[scale=0.38]{steady.pdf}
    \caption{The average photon number $\bar{n}$ as the system reaches a steady cycle. Regardless the initial average photon number $\bar{n}_{0}$, $\bar{n}$ eventually reaches a closed cycle \eqref{eq:steadylimit}. The photon number is constant during the work strokes, since the shortcut-to-adiabaticity protocol holds the energy populations constant. Each stroke takes a quarter of the total cycle time: $t_{\text{cycle}}/4$.}
    \label{fig:steady_cycle}
\end{figure}

During the heating (isochoric) strokes, the system frequency $\omega$ is constant, and the heat into the system is proportional to the change in the average photon number  $\bar{n}_{h}-\bar{n}_{c}$ of the system.

We assume the system begins in a thermal state. After a single cycle, the system has interacted with the hot bath for time $t_{h}$, and the cold bath for time $t_{c}$. If the systems begins with an arbitrary initial average photon number $\bar{n}_{0}$, it is not guaranteed that the cycle will be closed, since $\bar{n}$ may not return to its initial value. In this section we show that the system will always asymptotically approach a ``steady cycle", where the average photon number at the end of the cycle is the same as that at the beginning of the cycle.

Since $\bar{n}$ is constant during the work strokes, the average photon number after one cycle is given by two applications of \eqref{eq:nt},
\begin{align}
    \bar{n}_{1} = &\left(\left(\bar{n}_{0}-\frac{E_{h}}{\Delta_{h}}\right)e^{-\Delta_{h}t_{h}}+\frac{E_{h}}{\Delta_{h}}-\frac{E_{c}}{\Delta_{c}}\right) \nonumber\\
    &\times e^{-\Delta_{c}t_{c}} + \frac{E_{c}}{\Delta_{c}}.
    \label{eq:aftercycle}
\end{align}
Recall that we defined $\Delta \equiv G-E$, the difference in the transition rates of the master equation \eqref{eq:me}.
We have introduced subscripts $h$ and $c$ to denote the terms which arise due to interactions with the hot and cold baths respectively.
We want to know the initial photon number such that $\bar{n}_{1}=\bar{n}_{0}$. We denote it $\bar{n}_{sc}$,
\begin{align}
\bar{n}_{sc} &= \frac{\frac{E_{h}}{\Delta_{h}}e^{-\Delta_{c}t_{c}}\left(e^{-\Delta_{h}t_{h}}-1\right)+\frac{E_{c}}{\Delta_{c}}\left(e^{-\Delta_{c}t_{c}}-1\right)}{e^{-\Delta_{h}t_{h}}e^{-\Delta_{c}t_{c}}-1}. \nonumber \\
&= \frac{\frac{E_{h}}{\Delta_{h}}e^{\frac{-\Delta_{c}t_{c}}{2}}\sinh\left(\frac{\Delta_{h}t_{h}}{2}\right)+\frac{E_{c}}{\Delta_{c}}e^{\frac{\Delta_{h}t_{h}}{2}}\sinh\left(\frac{\Delta_{c}t_{c}}{2}\right)}{\sinh\left(\frac{\Delta_{h}t_{h}+\Delta_{c}t_{c}}{2}\right)}.
\label{eq:steadycycle}
\end{align}

If the interaction time with the cold bath $t_{c}$ is large, then $\bar{n}_{sc}$ is the thermalisation average photon number of the cold bath.
\begin{equation}
\lim_{t_{c}\rightarrow\infty} \bar{n}_{sc} = \frac{E_{c}}{\Delta_{c}}.
\end{equation}
So if the system begins with $\bar{n}_{0}=n_{sc}$, then after one cycle, the average photon number will return to $\bar{n}_{0}$. However even if $\bar{n}_{0}\ne \bar{n}_{sc}$, after multiple cycles the the average photon number at the beginning of a cycle will tend to $\bar{n}_{sc}$. After $\delta$ cycles (see Appendix~\ref{sec:Steady_cycle} for details)
\begin{equation}
    \bar{n}_{\delta} = e^{-\delta\Delta_{h}t_{h}}e^{-\delta\Delta_{c}t_{c}}\left(\bar{n}_{0}-\bar{n}_{sc}\right)+\bar{n}_{sc},
\end{equation}
so, as seen in FIG.~\ref{fig:steady_cycle},
\begin{equation}
    \lim_{\delta\rightarrow\infty} \bar{n}_{\delta} = \bar{n}_{sc},\label{eq:steadylimit}
\end{equation}
and the system always reaches a `steady cycle', where the average photon number returns to its initial value after each cycle. This has been previously studied in \cite{SCFeldmann1996,Kosloff2017}. We assume for the rest of this paper that the system has reached a steady cycle.

\section{Performance}\label{sec:results}

\subsection{The effect of coherence}
The work extracted in a single cycle is $\hbar\left(\omega_{h}-\omega_{c}\right)\left(\bar{n}_{h}-\bar{n}_{c}\right)$. We can use \eqref{eq:nt} to write
\begin{equation}
\bar{n}_{h}-\bar{n}_{c} = \left(\bar{n}_{c}-\frac{E_{h}}{\Delta_{h}}\right)\left(e^{-\Delta_{h} t_{h}}-1\right),
\end{equation}
where the average photon number begins from $\bar{n}\left(0\right)=\bar{n}_{c}$ and evolves to $\bar{n}\left(t_{h}\right)=\bar{n}_{h}$.
In a steady cycle, $\bar{n}_{c}$ is given by \eqref{eq:steadycycle}, and so
\begin{align}
\bar{n}_{h}-\bar{n}_{c} = \left(\frac{E_{h}}{\Delta_{h}}-\frac{E_{c}}{\Delta_{c}}\right)\frac{2\sinh\left(\frac{\Delta_{c}t_{c}}{2}\right)\sinh\left(\frac{\Delta_{h}t_{h}}{2}\right)}{\sinh\left(\frac{\Delta_{c}t_{c}+\Delta_{h}t_{h}}{2}\right)}. \label{eq:nh-nc}
\end{align}
Keeping the temperature of the baths constant, \mbox{$\bar{n}_{h}-\bar{n}_{c}$} is monotonically increasing in both $\Delta_{h}$ and $\Delta_{c}$ (see Appendix~\ref{sec:AppWork}, for the calculation).

We have seen that coherence in the energy ground space of the incoming atoms can increase or decrease the effective temperature of atomic thermal bath (see section~\ref{sec:interaction}).
To investigate the effect of coherence on the power of the engine, we study three pairs (hot and cold) of heat baths, all at the same effective inverse temperatures, $\beta_{h}$ and $\beta_{c}$. The first is a pair of incoherent (\textsc{I}) baths, with thermalisation photon numbers
\begin{equation}
\begin{matrix}
\dfrac{E_{h}^{\textsc{i}}}{\Delta_{h}^{\textsc{i}}} &= \dfrac{1}{e^{\beta_{h}\hbar\omega}-1}, \\[12pt]
\dfrac{E_{c}^{\textsc{i}}}{\Delta_{c}^{\textsc{i}}} &= \dfrac{1}{e^{\beta_{c}\hbar\omega}-1}.
\end{matrix}
\end{equation}
\begin{figure}[b]
\includegraphics[scale=0.38]{Bath_diagram.pdf}
\begin{tabular}{l | c | c | c | c }
Incoherent Baths (\textsc{I}) & $E_{h}$ & $\Delta_{h}$ & $E_{c}$ & $\Delta_{c}$ \\ \hline
Coherent Hot Bath (\textsc{CH}) & $E_{c}$ & $\Delta_{h}\frac{E_{h}}{E_{c}}$ & $E_{c}$ & $\Delta_{c}$ \\ \hline
Coherent Cold Bath (\textsc{CC}) & $E_{h}$ & $\Delta_{h}$ & $E_{h}$ & $\Delta_{c}\frac{E_{c}}{E_{h}}$
\end{tabular}
\caption{Diagram and table of three different pairs of baths, all with the same effective temperatures. Temperature according to thermal populations in hard line, effective temperature in broken lines. The coherent hot (\textsc{CH}) baths involve using coherence to create a hot bath from two cold baths; the coherent cold (\textsc{CC}) baths use coherence to create a cold bath from two hot baths. Without the use of coherence, these latter two bath systems would not be usable for running a thermal engine.}
\label{fig:bath_table}
\end{figure}
We compare this with two other pairs of baths with the same effective temperatures. The first pair (coherent hot or \textsc{CH}) have cold thermal populations consistent with inverse temperature $\beta_{c}$, however there is coherence in the energy ground space of one of the baths such that the effective inverse temperature is that of the incoherent hot bath $\beta_{h}$. Explicitly
\begin{equation}
E_{h}^{\textsc{ch}}=E_{c}^{\textsc{ch}}=E_{c}^{\textsc{i}},\quad \Delta_{c}^{\textsc{ch}} = \Delta_{c}^{\textsc{i}},\quad \Delta_{h}^{\textsc{ch}} = \Delta_{h}^{\textsc{i}}\frac{E_{c}^{\textsc{i}}}{E_{h}^{\textsc{i}}}.
\end{equation}
The second pair of thermal baths (coherent cold or \textsc{CC}) is the reverse: hot  thermal populations at inverse temperature $\beta_{h}$, but with sufficient coherence in the ground space of one of the bath such that it has an effective inverse temperature of the incoherent cold bath $\beta_{c}$. That is
 \begin{equation}
E_{h}^{\textsc{cc}}=E_{c}^{\textsc{cc}}=E_{h}^{\textsc{i}},\quad \Delta_{h}^{\textsc{cc}} = \Delta_{h}^{\textsc{i}},\quad \Delta_{c}^{\textsc{cc}} = \Delta_{c}^{\textsc{i}}\frac{E_{h}^{\textsc{i}}}{E_{c}^{\textsc{i}}}.
\end{equation}
These three different pairs of baths are summarised in FIG.~\ref{fig:bath_table}.
Without the addition of coherence, these latter two pairs of baths would be unusable; the addition of coherence creates an effective temperature difference, and allows them to be used to extract work, as seen in \cite{2003scully862}. In the first case the coherence is being used to turn a cold bath into a hot bath, and in the second case coherence is being used to change a hot bath into a cold bath. We note here that
\begin{equation}
\Delta_{h}^{\textsc{ch}} < \Delta_{h}^{\textsc{i}}\text{ and } \Delta_{c}^{\textsc{cc}} > \Delta_{c}^{\textsc{i}}.
\label{eq:newbaths}
\end{equation}
From \eqref{eq:nt}, $\Delta=G-E$ effects not only the effective temperature of the bath but also the rate at which the system thermalises. Faster thermalisation allows more work to be extracted per cycle. Creating the coherent hot bath requires a smaller $\Delta_{h}$ as compared with the incoherent bath pair, and less work will be extracted in a cycle. Conversely, the creation of the coherent cold bath requires a larger $\Delta_{c}$ than the incoherent baths, and more work will be extracted per cycle.

\subsection{Cost}\label{sec:cost}

\begin{figure}[t]
\includegraphics[scale=0.4]{Cost_with_inset.pdf}
\caption{The cost of the shortcut to adiabaticity protocol $V_{e}+V_{c}$ (see \eqref{eq:cost} and FIG.~\ref{fig:expand_compress}) for the three pairs of thermal baths. The total energy cost is scaled by $1/\hbar\omega$, where $\omega= 10^{14}$Hz. Recall that $t_{\text{cycle}}$ is in terms of the dimensionless time parameter \eqref{eq:me}. The three pairs of thermal baths are described in FIG.~\ref{fig:bath_table}. Both the coherent cold \textsc{CC} and the coherent hot \textsc{CH} systems have a reduced cost derived from their coherence (see Appendix~\ref{sec:AppCost} for details). The cost is large for $t_{\text{cycle}}\rightarrow 0$. As $t_{\text{cycle}}\rightarrow\infty$, the cost $V_{e}+V_{c}\rightarrow 0$. Here we used a hot bath with temperature $\bar{n}_{\mathrm{ss}}=2$, and a cold bath with $\bar{n}_{\mathrm{ss}}=0.55$.}\label{fig:cost}
\end{figure}

The cost of the shortcut-to-adiabaticity protocol \eqref{eq:cost} is proportional to the average photon number $\bar{n}$ of the system at the beginning of the protocol. This will be different for the different bath systems, since this depends not only on the temperature of the baths but also on the thermalisation rates $\Delta_{h}$ and $\Delta_{c}$. The cost of the compression stroke $V_{c}$ is proportional to the average photon number $\bar{n}_{c}$. Since we are assuming the system has reached a steady cycle, this is given by $\bar{n}_{sc}$ \eqref{eq:steadycycle}. If we keep the temperature of the baths constant but change the absorption rate $\Delta_{h}$ and $\Delta_{c}$, we can see that $\bar{n}_{sc}$ is monotonically increasing in $\Delta_{h}$, but monotonically decreasing in $\Delta_{c}$. This is also true for the expansion stroke cost $V_{e}$, during which the average photon number is $\bar{n}_{h}$ (see Appendix~\ref{sec:AppCost} for details).\\
This implies, from \eqref{eq:newbaths} that the use of coherence either to increase the temperature of a hot bath, or to decrease the temperature of a cold bath, both reduce the cost of the shortcut-to-adiabaticity protocol. As can be seen in FIG.~\ref{fig:cost}, the \textsc{CC} and \textsc{CH} bath systems using coherence have a smaller energy cost than the incoherent (\textsc{I}) bath system.

\subsection{Efficiency}\label{sec:efficiency}

\begin{figure}[t]
    \centering
    \includegraphics[scale=0.4]{efficiency.pdf}
    \caption{Efficiency curves \eqref{eq:efficiency} of the different baths. Recall that $t_{\text{cycle}}$ is in terms of the dimensionless time parameter \eqref{eq:me}. The three pairs of thermal baths are described in FIG.~\ref{fig:bath_table}. The coherent cold \textsc{CC} always has the highest efficiency \eqref{eq:effCC}. The coherent hot \textsc{CH} has the worst \eqref{eq:effCH}. Before $t^{\prime}$, the work extracted is less than the energy cost of the shortcut-to-adiabaticity protocol $W<V_{e}+V_{c}$. With longer cycle times, the cost of implementing the shortcut-to-adiabaticity protocol becomes negligible, and the efficiency approaches its maximum. Here we used a hot bath with temperature $\bar{n}_{\mathrm{ss}}=2$, and a cold bath with $\bar{n}_{\mathrm{ss}}=0.55$.}
    \label{fig:efficiency}
\end{figure}

The efficiency \eqref{eq:efficiencywithcost} is reduced from maximum because of the cost of implementing the shortcut-to-adiabaticity protocol.
\begin{equation}
\eta = \frac{W-V_{e}-V_{c}}{Q_{h}} = 1-\frac{\omega_{c}}{\omega_{h}}-\frac{V_{e}+V_{c}}{\hbar\omega_{h}\left(\bar{n}_{h}-\bar{n}_{c}\right)}.\label{eq:efficiency}
\end{equation}
Since the cost tends toward zero with increasing cycle time $t_{\text{cycle}}$, engines in which the expansion and compression strokes have a longer duration will be more efficient (as shown in FIG.~\ref{fig:efficiency}). There is a time $t^{\prime}$, before which, the energy cost of the shortcut-to-adiabaticity protocol will exceed the work extracted. This will result in a negative efficiency, rendering the engine too expensive to run.

\begin{figure}[t]
\includegraphics[scale=0.4]{power.pdf}
\caption{The power curves \eqref{eq:power} for the three pairs of heat baths. The power is scaled by $1/\hbar\omega$, where $\omega= 10^{14}$Hz. Recall that $t_{\text{cycle}}$ is in terms of the dimensionless time parameter \eqref{eq:me}. The three pairs of thermal baths are described in FIG.~\ref{fig:bath_table}. The \textsc{CC} baths always have a higher output power than the incoherent baths \eqref{eq:powerCC}. The \textsc{CH} baths do not necessarily have a lesser output, but they do for these temperatures (see Appendix~\ref{sec:AppPower} for details). As $t_{h},t_{c}\rightarrow\infty$, the difference $\bar{n}_{h}-\bar{n}_{c}\rightarrow\frac{E_{h}}{\Delta_{h}}-\frac{E_{c}}{\Delta_{c}}$, $W$ approaches a constant value and $P\rightarrow 0$. Before $t^{\prime}$, the engine is too costly to run as the cost of using shortcut-to-adiabaticity techniques is larger than the work extracted $W<V_{e}+V_{c}$. Here we used a hot bath with temperature $\bar{n}_{\mathrm{ss}}=2$, and a cold bath with $\bar{n}_{\mathrm{ss}}=0.55$.}
\label{fig:power}
\end{figure}

For the baths that we have hitherto examined, coherence has reduced the cost of the shortcut-to-adiabaticity protocol. However the reduction in the efficiency is also inversely dependent upon the heat input to the system by the hot bath $Q_{h}$, which is proportional to $\bar{n}_{h}-\bar{n}_{c}$ \eqref{eq:nh-nc}. We saw previously that $\bar{n}_{h}-\bar{n}_{c}$ is monotonically increasing in both $\Delta_{h}$ and $\Delta_{c}$ (see Appendix~\ref{sec:AppWork}). In the case of the $\textsc{CC}$ thermal baths the thermalisation rate $\Delta^{\textsc{cc}}_{c}>\Delta^{\textsc{i}}_{c}$, so the heat input from the hot bath $Q_{h}^{\textsc{cc}}$, is greater than the heat input using incoherent thermal baths $Q_{h}^{\textsc{i}}$,
\begin{equation}
Q_{h}^{\textsc{cc}} > Q_{h}^{\textsc{i}}.
\end{equation}
\begin{figure*}[t]
\subfigure[]{\includegraphics[scale=0.29]{maxpowerP.pdf}\label{fig:maxp}}
\subfigure[]{\includegraphics[scale=0.29]{maxpowerQ.pdf}}
\subfigure[]{\includegraphics[scale=0.29]{maxpowerR.pdf}}
\caption{The power curves of the engines can be increased by changing the relative times of the four strokes of the Otto cycle. This will change the amount of heat transfer to the system in each cycle, and alter the energy cost of using shortcut-to-adiabaticity techniques. The power is scaled by $1/\hbar\omega$, where $\omega= 10^{14}$Hz. The three pairs of thermal baths are described in FIG.~\ref{fig:bath_table}. The cycle time $t_{\text{cycle}}$ is parametrised in \eqref{eq:Ptcycle} and \eqref{eq:paramqr}. Here we numerically find the optimal values of $p,q,r$ for the \text{CC} baths. These are three cross-sections of the function Max Power $(p,q,r)$, one for each variable \eqref{eq:maxpqr}. In (a), $q=r=0.5$, and the peak occurs at $p\approx 0.58$ for the \textsc{CC} baths. In (b), $p=0.58$, $r=0.5$ and the peak occurs at $q\approx 0.69$. In (c), $p=0.58, q=0.69$ and the peak occurs at $r\approx 0.59$. See FIG.~\ref{fig:powerMAX} for the increase in the output power of the heat engine using these values. Here we used a hot bath with temperature $\bar{n}_{\mathrm{ss}}=2$, and a cold bath with $\bar{n}_{\mathrm{ss}}=0.55$.}
\label{fig:maxpower}
\end{figure*}
We can immediately deduce an improvement in the efficiency using the $\textsc{CC}$ thermal baths. Since the $\textsc{CC}$ thermal baths have a smaller shortcut-to-adiabaticity energy cost $V_{e}+V_{c}$ and a larger heat input $Q_{h}^{\textsc{cc}}$, then from \eqref{eq:efficiency} we must have
\begin{equation}
\eta^{\textsc{cc}} > \eta^{\textsc{i}}.\label{eq:effCC}
\end{equation}

For the $\textsc{CH}$ thermal baths, the coherence provides a reduction in the cost of the shortcut-to-adiabaticity protocol, however since $\Delta_{h}^{\textsc{ch}}<\Delta_{h}^{\textsc{i}}$, there is a reduction in the heat input to the system
\begin{equation}
Q_{h}^{\textsc{ch}} < Q_{h}^{\textsc{i}}.
\end{equation}
It can be shown (see Appendix~\ref{sec:AppEff} for details) that the reduction in heat into the system outweighs the benefit of coherence in reducing the cost of the shortcut-to-adiabaticity protocol, and we have (see FIG.~\ref{fig:efficiency})
\begin{equation}
\eta^{\textsc{ch}} < \eta^{\textsc{i}}.\label{eq:effCH}
\end{equation}

\subsection{Power}

The power generated from a single thermal cycle (of duration $t_{\text{cycle}}$) is the work per unit time
\begin{equation}
P = \frac{W-V_{e}-V_{c}}{t_{\text{cycle}}}.
\end{equation}
The work extracted is $W=\hbar\left(\omega_{h}-\omega_{c}\right)\left(\bar{n}_{h}-\bar{n}_{c}\right)$. Using \eqref{eq:nh-nc}, we can rewrite the power as
\begin{align}
P = &\;\frac{\hbar\left(\omega_{h}-\omega_{c}\right)}{t_{\text{cycle}}}\left(\frac{E_{h}}{\Delta_{h}}-\frac{E_{c}}{\Delta_{c}}\right)\frac{2\sinh\left(\frac{\Delta_{c}t_{c}}{2}\right)\sinh\left(\frac{\Delta_{h}t_{h}}{2}\right)}{\sinh\left(\frac{\Delta_{c}t_{c}+\Delta_{h}t_{h}}{2}\right)}\nonumber\\
&-\frac{V_{e}+V_{c}}{t_{\text{cycle}}}. \label{eq:power}
\end{align}
FIG.~\ref{fig:power} compares $P\left(t_{\text{cycle}}\right)$ for the different thermal baths. As we saw in the previous section, before time $t^{\prime}$ the cost $V_{e}+V_{c}$ is greater than the work extracted, and the engine becomes too costly to run. As $t_{h},t_{c}\rightarrow \infty$, the work approaches a constant value
\begin{equation}
\lim_{t_{h},t_{c}\rightarrow\infty} W = \hbar\left(\omega_{h}-\omega_{c}\right)\left(\frac{E_{h}}{\Delta_{h}}-\frac{E_{c}}{\Delta_{c}}\right),
\end{equation}
and $P\rightarrow 0$. Since $P$ is the product of transcendental functions, we cannot find an analytic expression for the cycle time $t_{\text{cycle}}$ at which $P$ is maximal. \\
To reiterate, the work $W$ is proportional to \mbox{$\bar{n}_{h}-\bar{n}_{c}$} which is monotonic in $\Delta_{h}$ and $\Delta_{c}$ (see Appendix~\ref{sec:AppWork}). Since $\Delta_{c}^{\textsc{cc}}>\Delta_{c}^{\textsc{I}}$, the work extracted from the \textsc{CC} baths in a cycle is greater than the incoherent (\textsc{I}) baths.
The \textsc{CC} baths extract more work, and have a smaller shortcut-to-adiabaticity energy cost, so we can immediately deduce an increase in power (FIG.~\ref{fig:power})
\begin{equation}
P^{\textsc{cc}}>P^{\textsc{i}}.\label{eq:powerCC}
\end{equation}
The \textsc{CH} baths also have a smaller energy cost, however $\Delta_{h}^\textsc{ch}<\Delta_{h}^{\textsc{i}}$, which implies that less work is extracted in a cycle from the \textsc{CH} baths than the incoherent (\textsc{I}) baths. The power of the \textsc{CH} baths need not be less than the incoherent baths (\textsc{I}), but in our simulations we find it typically is (see Appendix~\ref{sec:AppPower} for details). In FIG.~\ref{fig:power}, the decrease in work extracted in a single cycle is large enough such that there is a net decrease in the power $P^{\textsc{ch}}$ of the \textsc{CH} baths, i.e. $P^{\textsc{ch}}<P^{\textsc{i}}$.

\section{Optimising the Power}\label{sec:optimise}
\subsection{Using coherence}

Can we further improve the output power?
So far we have been looking at two extremes: using coherence to generate a hot bath (\textsc{CH}) and a cold bath (\textsc{CC}), with the other bath in each pair having no coherence (refer to FIG.~\ref{fig:bath_table}). But what about in-between these two extremes? We could compare these baths with others in which both hot and cold baths are created using coherence. We could look at the continuum in between these two extremes, parametrised by $\pi$:
\begin{align}
E_{h}^{\pi} = E_{c}^{\pi} &= \pi E_{h}^{\textsc{i}} + (1-\pi) E_{c}^{\textsc{i}}, \nonumber \\
\Delta_{h}^{\pi} &= \left(\pi E_{h}^{\textsc{i}} + (1-\pi) E_{c}^{\textsc{i}}\right)\frac{\Delta_{h}^{\textsc{i}}}{E_{h}^{\textsc{i}}}, \\
\Delta_{c}^{\pi} &= \left(\pi E_{h}^{\textsc{i}} + (1-\pi) E_{c}^{\textsc{i}}\right)\frac{\Delta_{c}^{\textsc{i}}}{E_{c}^{\textsc{i}}}, \nonumber
\end{align}
where $0\leq \pi\leq 1$. The limiting cases are those studied in the previous section: $\pi=0$ corresponds to the \textsc{CH} baths and $\pi=1$ corresponds to the \textsc{CC} baths. Along this continuum, the maximal work per cycle is achieved at $\pi=1$. In the previous section we mentioned numerous times (shown in Appendix~\ref{sec:AppWork}) that the work extracted in a single cycle is monotonically increasing with $\Delta_{h}$ and $\Delta_{c}$. Since $E_{h}^\textsc{i}>E_{c}^{\textsc{i}}$, we see that $\Delta_{h}^{\pi}$ and $\Delta_{c}^{\pi}$ are maximised when $\pi=1$. Using coherence to decrease the temperature of a bath increases the work extracted per cycle, and using coherence to increase the temperature of a bath decreases the work extracted per cycle. Thus it is expected to be that the optimal bath combination would be to use coherence to cool the cold bath, and to use no coherence in heating the hot bath.

This reasoning suggests the way to maximise the work extracted per cycle using coherence is to create both the hot and cold baths from hotter baths, using as much coherence as possible to decrease their effective temperatures. However as we saw in section~\ref{sec:interaction}, if a system is cooled maximally, then the effective temperature $T_{\text{eff}}\leq \hbar\omega/\ln\left(2\right)$. Such a bath can be created using coherence to cool a thermal bath with temperature
\begin{equation}
T_{R}=\frac{\hbar\omega}{\beta_{\text{eff}}\hbar\omega-\ln\left(2\right)}.
\end{equation}
A bath with temperature $T_{\text{eff}}\geq \hbar\omega/\ln\left(2\right)$ can be created by using coherence to cool a bath of infinite temperature.

\subsection{Changing stroke times}

Hitherto, we have been assuming that for a given cycle time $t_{\text{cycle}}$, each stroke had equal duration $t_{\text{cycle}}/4$. The system interacts with the heat baths for a total time $t_{Q}=t_{h}+t_{c}$, and work is done on the system for a time $t_{W} = t_{\text{cycle}}-t_{Q}$. We can break down $t_{W}$ further, as $t_{W} = t_{W_{e}}+t_{W_{c}}$, where $t_{W_{e}}$ is the duration of the isentropic expansion, and $t_{W_{c}}$ is the duration of the isentropic compression.

We can improve the output power of the heat engine by altering all the relative times of the four cycle strokes. We parametrise this using three variables.

First, we can alter the relative durations of the heat and work strokes within a single cycle. We parametrise the total cycle time $t_{\text{cycle}}$ as follows
\begin{equation}
\begin{matrix}
t_{Q} = p \times t_{\text{cycle}} \\
t_{W} = (1-p) \times t_{\text{cycle}}
\end{matrix}.
\label{eq:Ptcycle}
\end{equation}

Increasing the proportion $p$ of $t_{\text{cycle}}$ in which the system interacts with the heat bath causes more heat to flow in and out of the system, allowing more work to be extracted per cycle, at the price of a larger energy cost in running the shortcut-to-adiabaticity protocol.
We want to find the value of $p$ in which the peak power is maximised.
We cannot find an analytic expression for an optimal value for $p$, so we proceed numerically.

We can further parametrise $t_{Q}$ and $t_{W}$,
\begin{equation}
\begin{matrix}
t_{h} = q\times t_{Q} & t_{W_{e}} = r\times t_{W} \\
t_{c} = (1-q)\times t_{Q} & t_{W_{c}} = (1-r)\times t_{W}
\end{matrix}.\label{eq:paramqr}
\end{equation}
Increasing the proportion $q$ of $t_{Q}$ which the system interacts with the hot bath allows more heat to enter the system, increasing the work that can be extracted. However less heat is transferred to the cold bath, so more work is required during the compression stroke.
The costs of expansion ($V_{e}$) and compression ($V_{c}$) differ due to the different initial average photon numbers (see FIG.~\ref{fig:expand_compress}). Increasing $r$ reduces the cost of the isentropic expansion strokes; however, it increases the cost of the isentropic compression strokes.

\begin{figure}[t]
\includegraphics[scale=0.4]{powerMAX.pdf}
\caption{The output power can be improved by altering the relative times of each stroke of the Otto cycle. See FIG.~\ref{fig:maxpower} for the numerical search for the optimal values of $p,q,r$ for the \textsc{CC} baths. Power with $p=0.58$, $q=0.69$ and $r=0.59$ in thick lines, showing noticeable enhancement. Dotted lines show $p=q=r=0.5$ (where the duration of all strokes is $t_{\text{cycle}}/4$, i.e. FIG.~\ref{fig:power}). The power is scaled by $1/\hbar\omega$, where $\omega= 10^{14}$Hz. Recall that $t_{\text{cycle}}$ is in terms of the dimensionless time parameter \eqref{eq:me}. The three pairs of thermal baths are described in FIG.~\ref{fig:bath_table}. Here we used a hot bath with temperature $\bar{n}_{\mathrm{ss}}=2$, and a cold bath with $\bar{n}_{\mathrm{ss}}=0.55$.}
\label{fig:powerMAX}
\end{figure}
We want to find the values of $p,q,r$ such that the peak power is maximised. That is, we want to maximise the following function
\begin{equation}
\text{Max Power}\left(p,q,r\right) = \max_{t_{\text{cycle}}} P\left(t_{\text{cycle}},p,q,r\right).
\label{eq:maxpqr}
\end{equation}
In FIG.~\ref{fig:maxpower}, three cross sections of this function are shown, one for each variable. We see the values which maximise the power are different for each of the baths. We numerically found values of $p,q,r$ which maximised the Max Power for the \textsc{CC} baths with arbitrarily chosen temperatures. In this case we found a bias towards increasing the duration of the heat strokes over the work strokes, and in particular a higher proportion of the interaction time $t_{Q}$ to be with the hot bath. The power is also improved by slightly slowing the expansion strokes.
In FIG.~\ref{fig:powerMAX} the improvement in the power curve with these values is shown.

\section{Conclusion}

Coherence in the energy ground space of the atoms used to emulate a thermal bath does not change the average energy of the atoms, but changes the interaction of the atoms with the system; specifically the rate $G$ at which the system loses photons to the incoming atoms. This has two effects. First, it changes the thermalisation temperature of the system. So coherence can be used to create a temperature difference in two baths which previously were in equilibrium. Secondly, it changes the rate of thermalisation $\Delta$. These effects are inversely correlated, in that decreasing the thermalisation temperature causes an increase in the rate of thermalisation, and vice versa. In this paper we have shown how the coherence can be used to increase (or decrease) the power and efficiency of a heat engine. We used coherence to create a temperature difference between two bath systems in equilibrium; creating both a cold bath, and a hot bath. In both cases, the coherence reduced the cost of implementing the shortcut-to-adiabaticity protocol. Using coherence to create the cold bath improved the power of the engine, because it increased the thermalisation rate $\Delta$, allowing more heat transfer into the system in a cycle, and more work to be extracted. The efficiency was also improved, because the increased heat into the system reduced the effect of the cost of the shortcut-to-adiabaticity protocol. Conversely, creating a hot bath requires the thermalisation rate $\Delta$ to be reduced, which can in turn reduce the output power of the engine, and always reduces the efficiency. This demonstrates a more nuanced and subtle perspective of the benefit of coherence in the finite-time regime, it is a quantity which requires careful optimisation.
The output power can be further improved by using the maximum amount of coherence in cooling the heat baths, also by changing the relative durations of the four strokes in the thermal cycle.

This work can be extended by considering a quantum trajectory unravelling of the master equation \eqref{eq:me}. This would involve continuously monitoring the photon number of the system, and using that information to reduce the duration of the heat strokes. Of course the memory which stores the measurement results would need to be erased to close the cycle.
Considering atoms with a larger number of degenerate ground states, might allow more efficient energy transfer between the bath and system.
Furthermore, correlated states between the atoms would allow the study of more general quantum resources, extending beyond Markovian interactions.

\begin{acknowledgements}
T. Guff would like to thank J. D. Cresser for his many useful conversations.
This research was funded in part by the Australian Research Council Centre of Excellence for Engineered Quantum Systems (Project number CE110001013, CE170100009).
\end{acknowledgements}