Nonequilibrium fluctuations in quantum heat engines: Theory, example, and possible solid state experiments

We study the stochastic energetic exchanges in quantum heat engines. Due to microreversibility, these obey a fluctuation relation, called the heat engine fluctuation relation, which implies the Carnot bound: no machine can have an efficiency larger than Carnot's efficiency. The stochastic thermodynamics of a quantum heat engine (including the joint statistics of heat and work and the statistics of efficiency) is illustrated by means of an optimal two-qubit heat engine, where each qubit is coupled to a thermal bath and a two-qubit gate determines energy exchanges between the two qubits. We discuss possible solid state implementations with Cooper pair boxes and flux qubits, quantum gate operations, and fast calorimetric on-chip measurements of single stochastic events.


Introduction
The field of non equilibrium quantum thermodynamics has received a large impulse in the last two decades due to the discovery of a number of exact relations which characterise the response of physical (possibly small) systems, to external perturbations, namely applied mechanical forces or thermodynamic forces (e.g.temperature gradients, and chemical potential gradients).[1,2] Unlike traditional thermodynamics [3], which focusses on macroscopic quantities, fluctuation relations focus on their microscopic, fluctuating, counterparts.To exemplify this, consider the two fundamental objects of thermodynamic investigation, work and heat.A macroscopic thermal engine delivers a certain amount of work while withdrawing a corresponding amount of heat from a hot thermal reservoir.There can be variations in these amounts between different cycles, but typically these fluctuations are negligible.However as the machine size scales down, likewise will the work output and heat absorbed scale down.Accordingly, their fluctuations will become more and more relevant.It then becomes useful to investigate the stochastic properties of such fluctuating quantities.Fluctuation relations pose stringent constraints on the statistics of such fluctuating quantities like heat and work, due to the symmetries (in particular time-reversal symmetry) characterising the microscopic motions of atoms and molecules form which they originate.
Fluctuation relations have been reported for both classical and quantum systems [1,2,4,5,6].In fact identical fluctuation relations hold regardless of whether the same system is regarded as classical or quantum.Despite their formal identity, classical and quantum fluctuation relations are deeply different in the way they can be accessed experimentally.Concerning work, for example, while typically one can measure the fluctuating work applied to a classical nano system, e.g., a stretched RNA molecule, by continuously monitoring a displacement x and its conjugate force f (e.g.extension and tension in the molecule), and obtaining the work as W = − f dx, [7,8,9] this is typically impossible in a quantum system.In the quantum scenario, the situation is much complicated by the invasiveness of the measurement apparatus which can lead to a collapse of the wave function.The prescription accordingly is to measure the energy of the system twice (at the beginning and end of the forcing protocol), by means of two projective measurements and obtain the work as their difference [10,11,12,13].
This two-measurement scheme has proved however very challenging from the experimental point of view [14,15], so much that it has been realised only very recently [16].This occurrence has triggered the proposal of a number of alternative methods.One such method proposes to replace the two invasive projections with many less invasive measurements (POVM) carried on a smaller portion of the system [17].This method is particularly well suited for studies of transport induced by gradients of temperature and chemical potential [17].Some experiments already exist which can be explained in terms of these multiple measurements [18,19].They regard the full counting statistics of electrons transported through double quantum dot due to an applied chemical potential difference.We shall remark however that in those experiments all quantum coherences are suppressed.
Another very ingenious method, which is particularly well suited for obtaining the work statistics of a driven system, requires a special coupling of the driven system to an ancilla, e.g. a qubit, and replaces the two energy measurements with state tomography of the qubit at the sole final time [20,21].This method is a form of Ramsey interferometry and gives experimental access to the characteristic function of work, namely the Fourier transform of the probability density function of work.This has led to the first experimental measurement of quantum work statistics ever performed.It has been performed in a liquid-NMR set-up, and has reconstructed the work pdf of a driven two level system [22].A proposal for implementing the method with solid state quantum devices has been put forward in [23].The most promising aspect of this method is that it can be used not only to asses the work statistics of closed systems as in the performed experiment, but also of systems which stay in contact with a thermal bath [23].
Roncaglia et al. [24] have proposed to couple the system to a quantum pointer, e.g. a spin chain.The coupling is engineered so that that a single final projective measurement of the state of the pointer will contain information about the work performed on the system.Like the interferometric method, this method is best suited for the measurement of work.Its experimental realisation however appears extremely challenging.
In this work we focus on yet another method which has been discussed recently in [25,26] and is based on the calorimetric measurement of photon released and absorbed by thermal reservoirs.This quantum calorimeter is currently under development.The method is well suited for simultaneously measuring both heat and work in a driven quantum system which stays in contact with one or more baths.For this reason it is very promising for the experimental study of the stochastic energy exchange of quantum thermal machines.
Since the seminal work of [27], showing how the three level maser could be understood as a thermal machine, quantum thermal machines have been widely studied in the literature [28,29,30] and are still under vigorous investigation [31,32,33,34,35,36,37,38,39,40], see also the recent review [41] and references therein.However, while so far the focus was on the average value of heat and work, here we focus on their fluctuations as well.As recently reported [42], a special form of the fluctuation relation holds for quantum thermal machines.This form implies that no quantum thermal engine can over-perform the Carnot efficiency.This universal and exact result was anticipated long ago in [28] but only for those quantum mechanical open systems whose dynamics can be well approximated by a Markovian master equation in Lindblad form.
After revisiting the heat engine fluctuation relation, we shall introduce a model of thermal engine based on two-qubits each coupled to its own reservoir and subject to a unitary gate operation.We will identify the regimes when the engine works as heat engine, refrigerator, or heater (dud engine), and study its full stochastic characteristics, including the probability density function of its efficiency.The most intriguing features of the presented machine are (a) that at maximum power it can reach efficiency above the Curzon-Albhorn efficiency, and (b) that increasing the speed of its operation increases the power output without affecting its efficiency.
It is important to stress that the engine presented here can be implemented in a real solid state device and its stochastic energetic exchanges can be measured using the current and soon available technology.Below we discuss possible solid state implementations based on the calorimetric measurement scheme.

The Heat Engine Fluctuation Relation (HEFR)
Consider a driven bi-partite system: with factorized initial condition, Without lack of generality we shall assume throughout this work β 1 ≤ β 2 , i.e., the first sub-system is assumed to be not colder than the second, at the initial time.Also we shall assume that at all times the Hamiltonian is time reversal symmetric [43].We further assume the compound system is thermally isolated and the driving is turned on at time t = 0 and turned off at time t = τ .At these two times simultaneous projective measurements of the energies of both sub-systems are performed, giving the results , where i = 1, 2 and E i k is the k-th eigenvalue of sub-system i.According to the quantum exchange fluctuation theorem [17,44,45], it is where ∆E i = E i m i − E i n i is the observed energy change in sub-system i, P (∆E 1 , ∆E 2 ) is the joint probability of observing ∆E 1 and ∆E 2 , and P (−∆E 1 , −∆E 2 ) is the joint probability of observing −∆E 1 and −∆E 2 when the reversed driving V (τ −t) is applied.The driving V (t) injects some amount of energy in the compound system: which is in fact the work performed by the external driving source to drive the system.Part of this energy, ∆E 1 goes into sub-system 1, and part of it, ∆E 2 , goes into subsystem 2. Using the above equation to make the change of variable ∆E 2 → W , we obtain a fluctuation relation for the joint probability of work W and ∆E 1 : Multiplying by P (−∆E 1 , −W ) and integrating in dW d∆E 1 , one obtains the integral form of the fluctuation relation . Scheme of a quantum thermal machine.An isolated system, (big green rectangle) is driven by an external time dependent field.The system is composed of two sub-systems (red and blue rectangles).Each subsystem is composed of a small quantum system (small circle) and a large system, namely a thermal reservoir (large circular section).The two small circles form the working substance, we shall call them working parts.The drive acts on the working substance, thus injecting work W in the whole system.A part of it, ∆E 1 , is delivered to subsystem 1, via the working part 1.The rest, ∆E 2 = W − ∆E 1 , is delivered to subsystem 2, via the working part 2. Each working part retains a part of the delivered energy ∆U i , and dumps the rest −Q i = ∆E i − ∆U i into its reservoir.These energetic exchanges are possible due to possibly time-dependent couplings between the two working parts, and between each working part and its reservoir (dashed lines).At the beginning of the driving each subsystem is at thermal equilibrium with a given temperature T i .
Using the Jensen's inequality as usual, one obtains from this where η C = 1 − β 1 /β 2 is Carnot's efficiency.The above equations hold regardless of size of the two subsystems, as long as the assumptions introduced are satisfied.In modelling a quantum thermal machine we shall consider each subsystem as composed of two parts, namely a heat reservoir and a small quantum system which constitutes a part of the working substance, see Fig. 1.We shall call the small quantum systems the working parts.The driving is applied on the working substance.The received work W , is shared between sub-system 1 and 2 as ∆E 1 and ∆E 2 .We allow for the possibility of a time dependence of the couplings between the reservoirs and the working parts, in which case we consider them as part of the time dependent part of V (t) of the Hamiltonian.This encompasses continuous mode thermal machines, where the couplings between the working parts and their respective reservoirs are constant in time and non-vanishing, and machines operating in discrete mode (via distinct strokes) where those couplings can be switched on and off during operation.The three level maser is an example of continuous mode engine while Carnot, Otto, Diesel engine etc. operate in discrete mode.
The average quantities ∆E 1 , ∆E 2 , W define the operation regime of the machine: When the thermal machine works as a heat engine, Eq. ( 7) gives This is the second law of thermodynamics as expressed for a heat engine.Our derivation proves its universality based on the time-reversal symmetric unitary dynamics of the whole system, and the initial bi-Gibbsian preparation.In a similar way, when the machine operates as refrigerator, one finds Before proceeding it is worth remarking that there is a freedom of arranging the position of the border between the two subsystems, i.e. to arrange the initial bi-Gibbsian equilibrium.In fact fluctuation relations for heat engines have been derived previously assuming the working substance is fully included in one of the two subsystem only, say subsystem 2 [42,46].In that case work W is delivered to subsystem 2, which retains a part ∆E 2 (shared between reservoir 2, −Q 2 , and working substance ∆U 2 ), an dumps the other other part ∆E 1 = −Q 1 directly into reservoir 1 as heat.That arrangement is particularly useful for a machine working as heat engine, because Eq. ( 7) would read − Q 1 / W ≤ η C as in standard thermodynamics books.
Here we adopt instead the scheme in Fig. 1 because we have in mind an implementation where the coupling of the working parts to the reservoirs is fixed and cannot be manipulated, while one can turn the interaction between the working parts on and off.By keeping this coupling off, it is then straightforward to prepare each working part in thermal equilibrium with its own bath.This corresponds to the scenario depicted in Fig. 1.With our arrangement the average energy ∆E 1 can be identified with the heat − Q 1 only when the energy ∆U 1 stored in the working part is null or negligible as compared to ∆E 1 and − Q 1 .This happens when the number of cycles is long and the working substance has a finite energy spectrum.Then ∆U 1 remains bounded while ∆E 1 and − Q 1 grow linearly in time.If the condition is met then Otherwise, if the condition is not met one can well have the ratio − Q 1 / W be larger than η C .This however does not have an impact on the second law of thermodynamics stating that a machine working in a cycle (implying ∆U 1 = ∆U 2 = 0) cannot have an efficiency larger than Carnot's efficiency.
It is important to stress that the choice of borders and appropriate associate thermodynamic quantities is the key to obtaining exact transient fluctuation relations, like Eq. ( 5), that is fluctuation relations that hold regardless of the time duration of the process under investigation [2,17,47,48].In the long time limit, steady-state fluctuation relations hold which are independent of the border choice.

Optimal two-qubit engine
Our aim is to propose a minimal model of thermodynamic quantum engine which could be implemented and tested experimentally as a solid state quantum device.The simplest model one can think of is that of a single qubit coupled to two reservoirs at different temperatures.The qubit is driven by an external drive which changes its Hamiltonian in time for example by changing its energy spacing H qubit (t) = ω(t)σ z /2.If one has the further ability to couple and decouple the qubit from the two reservoirs one can implement a 4-stroke engine, e.g. a Otto cycle.This can be realised, e.g., by interfacing the qubit to the thermal reservoirs by means of band-pass filters, as proposed in Ref. [49].Here we focus instead on the case when the coupling to the reservoirs are fixed in time.In order to have a heat engine/refrigerator in continuous mode, a more complex working substance is necessary than a mere two level system.One needs a working substance that would be able to de-route the energy towards the wanted direction (from the hot bath to work source and cold bath for a heat engine; from the cold bath and the work source to the hot one for a refrigerator).For this reason we introduce a second qubit.Qubit one is in contact with the first bath and qubit 2 is contact with bath 2, as in Fig. 1.A time dependent coupling V (t) couples the two qubits for a time duration [0, τ ].The full Hamiltonian is: where H B,i , H int,i , i = 1, 2, are the i-th bath Hamiltonian and its interaction with qubit i, respectively, and is the the i-th qubit Hamiltonian.Here σ z i denotes the z Pauli sigma matrix of the i-th qubit.
To keep the discussion as simple and intuitive as possible we introduce a useful assumption, namely the coupling V (t) is turned on for a time period [0, τ ] that is very short compared to the relaxation time of each qubit to its own bath.The effect of the coupling V (t) can accordingly be modelled by a unitary operator U acting in the Hilbert space of the working substance, namely the two qubits.We shall call U the gate operation.The two qubits are initially each in thermal equilibrium with its own bath, i.e. their state is characterised by the density matrix with Z i = Tr e −β i H q,i = 2 cosh(β i ω i /2).The average work injected into the working substance by applying the unitary U is and the energy taken by each sub-system is: If after the application of the gate U each qubit is let interact with its respective reservoir for a sufficiently long time so as to reach the state of thermal equilibrium.During this thermalisation step they will give the heats − Q i = ∆E i to the baths.We are interested in the unitary that outputs the most work per cycle.Therefore we have searched for the unitary that maximises W .We have pursued this task by parametrising a 4×4 unitary by means of 15 angles as discussed in [50] and performing a maximisation over the corresponding 15 dimensional space.Numerics clearly indicates that maximum work output is achieved by means of the complex SWAP unitaries, reading in the {|+, + , |+, − , |−, + , |−, − } basis: With these U's we find In the following we fix the gate to be any complex swap gate in Eq. (15).Quite remarkably, the same unitaries also maximise the heat engine efficiency.

Operation
The operation of the swap-machine is dictated by the relative signs of ∆E 1 , ∆E 2 , W .With β 1 ≤ β 2 the conditions for each mode of operation are: • HEAT ENGINE: The explanation of the above conditions is as follows.After the SWAP-gate operation is performed the two qubits are in the states where If , hence the cold qubit cools down and the hot qubit heats up ∆E 1 > 0, ∆E 2 < 0. Also, since ω 2 /ω 1 < β 1 /β 2 < 1, then W > 0. Hence we have the fridge operation.If , hence the hot qubit cools down and the and the cold qubit heats up ∆E 1 < 0, ∆E 2 > 0. In this case, depending on the relative size of ω 1 and ω 2 we will have either heat engine or heater.Let u then W < 0 and we have the heat engine.Otherwise dud engine.Fig. 2 shows a cartoon of the operation of the machine in the refrigerator mode.

Efficiency
For the heat engine operation, it is 0 For a fridge because for the fridge β 1 β 2 < ω 2 ω 1 < 1.Note that in order for the engine to function the two qubits must have different energy spacings ω i , otherwise the work intake (output), will be exactly null.Not also that the efficiency depends only on the ratio ω 2 /ω 1 and not on the temperatures β 1 , β 2 .This is a peculiar feature of the swap unitary.In the following we shall focus on heat engine operation.

Efficiency at maximum power
Given the two temperatures T 1 , T 2 the maximal efficiency, i.e., Carnot's efficiency is reached when ω 2 /ω 1 → β 1 /β 2 .In this regime however, the work tends to zero, see Eq. (18).It is interesting that here the power at Carnot efficiency is zero, as with standard strokes engines, but not because of slow operation.
On the other hand, given the two temperatures T 1 , T 2 one can find the value of ω 1 and ω 2 for which the power output, − W , is maximum.This can be achieved by maximising the work output in Eq. ( 18).The maximum depends indeed only on the ratio Ω = ω 2 /ω 1 .This can be best seen by setting ω 1 as the unit of energy, so that ω 1 = 1, and all energies are measured as multiples of ω 1 .With these units We denote the value of Ω for which − W is maximum at given β 1 , β 2 as Ω * (β 1 , β 2 ).The corresponding efficiency, namely the efficiency at maximum power is: For example, the value of Ω * is Ω * = 0.83 for k B T 1 = 3/2, k B T 2 = 1 (in units of ω 1 as explained above).The corresponding efficiency at maximum power is η * ≃ 0.17.
The figure shows that the the maximum power efficiency can be both larger and smaller than the Curzon-Albhorn efficiency.However for sufficiently low β 1 (hotter hot reservoir), η * < η CA , while, for sufficiently high β 1 , (colder hot reservoir), η * > η CA , that is at very low temperature the efficiency at maximum power is above the Curzon-Albhorn efficiency.
We have performed an analysis of the maximum power efficiency η * for β 1 ≃ β 2 , i.e., in the low η C limit.In accordance to linear response theory we expect the linear coefficient of the expansion to match the value 1/2 [52].Since our engine is not a thermoelectric engine (work is provided by time-dependent pulses, rather than by a DC electric potential difference), and does not have the left-right symmetry (in order for it to output some work the energy spacings of the two qubits, ω 1 and ω 2 , should be different), we do not expect that the value 1/8 for the quadratic coefficient, predicted in those cases Ref. [53], to be obeyed.The results of the low η C analysis, reported in Fig. 4, corroborate these expectations.The figure presents plots of η * /η C for various values of β 2 as function of η C .The plots clearly show that that is, the linear coefficient 1/2 is obeyed whole the quadratic coefficient is a function f (β 2 ) whose value may differ from 1/8.

Modelling: Quantum jumps
The above analysis based on the simplified assumption of unitary gate followed by thermalisation, allowed us to make predictions about the average work and heats that go in the two reservoirs.It does not suffice however for the full stochastic characterisation of the engine.In order to achieve that we need to model the dynamics of the thermalization.We assume then that the effect of the thermal environment on each qubit can be modelled by means of a master equation of Lindblad form. where and σ i = σ x i + iσ y i is the annihilation operator for the qubit i, and σ i , its adjoint, is the creation operator.
To obtain the statistics of energy exchanges with the bath during the thermalisation step, we proceed to un-ravel the master equation [54], as proposed in [55] and [26].This results in a stochastic differential equation in the Hilbert space of each qubit: The deterministic part is given by: while the stochastic Poisson increments have the ensemble expectations The stochastic equations can be solved by means of the Monte Carlo Wave Function (MCWV) method [56].In the present case of an undriven single qubit they result in a dichotomic Poisson process governed by the two rates: depending on whether the qubit is in the down state |− or up state |+ .Note that these rates are detailed balanced: Ref. [55] has studied the fluctuation relations for such quantum trajectories but for systems being in contact with a single thermal reservoir.Here our working substance, the two qubits, is in contact with two distinct reservoires.The analysis performed in [55] can however be extended to multiple reservoirs.It results in the following fluctuation relation for the probability of a given quantum trajectory γ: where γ is the time reverse of γ.We remark that in the case of multiple reservoirs γ is not only specified by the temporal evolution of the state of the central system (the working substance in our case), call it χ t , but also by the succession i n indicating which bath (labelled by i) was responsible for each of the N jumps (labelled by n), γ = ({χ t }, {i n }) Accordingly the time reversed trajectory results by requiring that the temporal evolution of the central system state is inverted and if the n-th jump of the forward trajectory γ was caused by the i-th bath, so was the last n-th jump of the backward trajectory γ.That is γ t = ({χ T −t }, {i N −n }).In our case the trajectory γ has two components γ t = (γ 1,t , γ 2,t ), each specifying the temporal evolution of the state of each qubit.No extra indexes are necessary because all jumps in γ 1,t are caused by reservoir 1 and all jumps in γ 2,t are caused by reservoir 2. Accordingly , means the heat ceded to the i-th reservoir during the realisation of γ.Specifically )ds.Obviously in our case Q i is a functional of γ i only.In Eq. ( 41) a, b denote the initial and final state of the trajectory γ, i.e. γ 0 = a, γ T = b, and p a,b are the respective probability that these states are observed.With our choice (12) Multiplying by ) and performing a path integral over all trajectories γ one recovers Eq. ( 3).Accordingly all subsequent relations in Sec. 2 are obeyed within our quantum jump modelling.

Stochastic thermodynamics of the SWAP engine
We operate the machine in the following manner.At time t = 0 we pick up a state randomly from the initial bi-Gibbsian distribution, Eq. ( 12).We apply the complex SWAP gate, Eq. ( 15), and generate the stochastic dynamics of each qubit using the MCWF method until time τ 2 , when we apply the complex SWAP again, and let evolve stochastically until time 2τ 2 , and so on for a total duration T = Nτ 2 .Our assumption is that the swap gate is much faster than the stochastic evolution time: τ ≪ τ 2 .Figure 4 shows a sketch of the resulting quantum trajectories of the two qubits, along with the energetic exchanges the various jumps signal.
The first important observation from Fig. 4 is that any time the energy ∆E 1 is given to subsystem 1, accordingly the energy ∆E 2 = −(ω 2 /ω 1 )∆E 1 is taken from subsystem 2. This implies that all trajectories have the same efficiency η = W/∆E 1 = (∆E 1 + ∆E 2 )/∆E 1 = 1 − ω 2 /ω 1 = η.In other words there are no efficiency fluctuations.This is because the gate swaps the eigenstates of the double qubit without creating superpositions thereof.Hence each swap pulse k deterministically and univocally results in well defined values of ∆E k 1,2 depending on the state of each qubit before its application.A generic unitary will typically create a superposition of the eigenstates, which can collapse either in the up state or down state of each qubit with according probability.The value of ∆E k 1,2 would be accordingly not uniquely defined by the state before a generic gate.
The constraint ∆E 2 = −(ω 2 /ω 1 )∆E 1 allows to express the heat engine fluctuation relation (5) as a relation for a single variable, say W . Since W is an integer multiple N W of ω = ω 1 − ω 2 , the fluctuation relation can be conveniently expressed as a relation for the probability P (N W ) that N W of work quanta are given off by the work source.We obtain then Figure 6 shows P (N W ) and the corresponding logarithmic ratio log P (N W )/P (−N W ) for one simulation of our engine.In an experimental realisation the probability P (N W ), can be constructed by recording the number and sign of the swaps occurred during each of many realisations in just one of the two qubits.
In Fig. 7, left panel, we report a plot of the joint probability distribution of heat and work P (Q 1 , W ). Note how W/(ω 1 − ω 2 ) differs from Q 1 /ω 1 at most by one unit.This is because ∆E 1 /ω 1 differs from Q 1 /ω 1 at most by one unit, i.e. one quantum of energy stored in qubit one as ∆U 1 .Because Q 1 is not exactly equal to ∆E 1 , the heat-efficiency η Q = −W/Q 1 has some fluctuations, in contrast to the ∆E-efficiency η = −W/∆E 1 .The statistics of η Q corresponding to the plot in Fig. 7 is reported in the right panel of Fig. 7.
Note the very pronounced peak at η.Note also a second peak η C .We observe that there is a finite probability that η Q is infinite.Because of the peak at infinity the quantity η Q is not well defined.As the number of cycles increases the spot in Fig 5 drifts and diffuses in the diagonal direction, but not in the transverse direction.A  ).The operation time of the machine is fixed and equal to T op = 30τ relax (τ relax is defined here as the longest among the thermal relaxation times, i.e., τ relax = γ −1 max[e β1ω1 − 1, e β2ω2 − 1]).At N = 10 cycles, the engine has plenty of time to relax to equilibrium, because each pulse is followed by a rest time of 3τ relax .By increasing the pulse frequency one can greatly enhance the work output.Solid line, ω 1 = 1, ω 2 = 0, 7, corresponding to efficiency η = 0.3.Dashed line: ω 1 = 1, ω 2 ≃ 0.83 corresponding to efficiency at max power η * =≃ 0.17.
consequence of this is that the peak at η in the efficiency probability increases while all other peaks decay.That is for large operation time the probability of η Q coincides with the probability of η as expected.‡

Increasing the power
As discussed above the swap heat engine works at the efficiency η = 1−ω 2 /ω 1 regardless of the power output.This is a great advantage over traditional engines because increasing the power has no cost in terms of reducing the efficiency for our engine.The power of our engine can be increased simply by increasing the swap-pulses frequency.This is illustrated in Fig. 8. Figure 8 suggests that the power output saturates at a regime value as the pulse frequency increases.The saturation value gives the maximum power for the given T 1 , T 2 , ω 1 , ω 2 .It is important to stress that for too frequent swap pulses, namely when their temporal separation τ 2 is of the same order as the temporal duration τ of the swap pulse, our simplifying assumption (namely that dynamics can be modelled separately by a unitary followed by the stochastic relaxation) does not hold any more.In a real experiment the power output is expected to decay in the range of highly frequent pulses.‡ According to large deviation theory, all peaks but the most likely fall with an exponential rate, the largest of which is for the Carnot efficiency [57].

Solid state implementation and measurement scheme
For our implementation proposal we follow the scheme presented in Ref. [25], with the necessary modifications and extensions.Ref. [25] presents an experimental scheme where a single Cooper-Pair-Box (CPB), namely a qubit, is coupled to a resistor at some temperature T .The resistor comprises an electronic system coupled to a phononic one.When a photon is emitted (absorbed) into the resistor, the fast electronic system responds by heating up (cooling down) abruptly, and then relaxing to the thermal equilibrium set by the phononic substrate.A nano calorimeter can then be used to monitor the temperature of the electronic system, in order to detect absorbed/emitted photons.As reported in [58], sufficienty fast and sensitive calorimeters for this pourpose are currently under development and should be soon available.
In order to realise the SWAP engine, two such CPB + resistor systems should be realised on the same microchip, which does not seem to pose any particular difficulty.Each resistor is then monitored by an on chip calorimeter of the type in Ref. [58].At variance with the set-up proposed in Ref. [25] here the two CPBs have fixed energy gaps, hence they exchange photons of well defined energy ω i .This simplifies the measurement, because each calorimeter needs not measure the energy of the absorbed/emitted photon, but should just detect that a photon has been absorbed/emitted.The gate operation can be implemented by coupling the two CPB using two tunnel junctions connected in parallel as described in [59].This allows for the implementation of the iSWAP gate, namely the complex SWAP gate Fig. 9, top panel, shows the scheme of the implementation.An alternative implementation uses flux qubits operating at the optimal point, which have much longer coherence and relaxation times as compared to CPBs [60,61].The switchable coupling is realised by means of a third qubit sandwiched between the two qubits as demonstrated in Ref. [62].iSWAP gate can be realised by means of microwave driving for a targeted time duration [62].As demonstrated in Ref. [62], the minimum time for a universal gate is, in that set-up, about 22ns, while decoherence and relaxation times are at least 0.2 µs.This is in agreement with our assumption of fast gate operation, as compared to thermal relaxation.Fig. 9, bottom panel, shows this alternative implementation.
By means of the calorimetric measurement one can experimentally access the quantum trajectories of the type shown in Fig. 4.This is achieved in the following way.The calorimeters can only detect the heat quanta ∆Q i k ceded to the resistors i = 1, 2. If two consecutive emissions (absorptions) are observed to occur in the same bath i, it means that meanwhile a quantum of energy ∆E i k = +(−)ω i has been given Bottom panel: Two flux qubits are coupled via a third flux qubit as from Ref. [62].Both setups can be used to realise the iSWAP gate, Eq. ( 44).The two qubits exchange photons each with a different resistor kept at a given temperature.Two on-chip fast calorimeters detect single photon emission/absorption in each resistor.
to (taken from) the i th qubit by the work source.Summing up all the ∆E i k one obtains the total energy given to each subsystem, ∆E 1 , ∆E 2 and the work W = ∆E 1 + ∆E 2 .Having Q 1 , Q 2 , ∆E 1 , ∆E 2 , W one can address the full statistics of energetic exchanges of the engine, and accordingly can check the validity of the fluctuation relations (3,5).Note that the measurement apparatus can also be employed to check the coincidence of swap-induced jumps in the two qubits, thus quantifying the goodness of the swap operation.
The employment of flux qubit for implementing on chip coolers have been discussed also in Ref. [49].In the work of Ref. [49] the working substance is a single flux qubit which is alternatively coupled and decoupled from the two baths.This is attained by embedding the two bath-resistors each in a LCR circuit, acting as band-pass filters centred at different frequencies ω 1 , ω 2 .As the qubit level spacing is switched between these two values the qubit interacts primarily with one resistor or the other so as to realise a Otto cycle, where interactions with the cold and hot bath occur in alternation and are separated by slow, adiabatic drives.This realises the same average heat and work exchanges as the present engine, hence same efficiency, with the difference that the present engine works in continuous mode.Heat exchanges with the two baths occur here simultaneously, and no adiabatic drive is employed.This makes it more promising in regard to the delivered power.

Conclusions
Based on a previous work [42], we have here presented a detailed discussion of fluctuation relations for heat and work in quantum heat engines.These fluctuations are illustrated by means of an optimal two-qubit engine working in continuous mode.We studied its full stochastic energetic exchanges including the statistics of its efficiency.At the average level this engine realises the same thermodynamics as the single qubit Otto engine of Ref. [62] but is expected to deliver a higher power due to its continuous mode of operation (no adiabatic sweeps needed).We have presented possible implementations which employ Cooper pair boxes and flux qubits as working substances, two-qubit quantum gates, and on-chip fast calorimetry for the detection of single exchanged energy quanta.The proposed experiment would allow for the first fully stochastic characterisation of a quantum heat engine.

Figure 2 .
Figure 2. Illustration of the functioning of the SWAP engine in the refrigerator mode.a) Each qubit is in thermal equilibrium with its respective bath.b) The instantaneous SWAP gate is applied resulting in the injected work W .The hot qubit becomes hotter while the cold qubit becomes colder.c) Qubit 1 cedes heat to the hot bath.Qubit 2 withdraws heat from the cold bath.The initial equilibrium a) is re-established.In heat engine mode the swap cools the hot qubit and heats the cold qubit, while outputting work.The sign of the heat and work arrows gets inverted accordingly.

Figure 3 .
Figure 3. Efficiency at maximum power as function of β 1 for various values of β 1 .

Figure 3
Figure3shows the maximum power efficiency η * (β 1 , β 2 ) as a function of β 2 for various fixed values of β 1 .The figure also report the corresponding Carnot efficiency and the Curzon-Albhorn efficiency[51]

Figure 5 .
Figure 5.Typical quantum trajectory of the working substance of the SWAP engine.Top trajectory is for qubit 1. Bottom trajectory is for qubit 2. Green vertical lines indicate times when the instantaneous SWAP gate is applied.Transitions occurring at these time are due to the work done by the external work source.Transitions occurring between the SWAP pulses signal heat exchanges with the heat reservoirs.

Figure 6 .
Figure 6.Left panel: Probability P (N W ) of work quanta N W given off by the work source.Right Panel: Corresponding logarithmic fluctuation ratio ln[P (N W )/P (−N W )]. Dots: numerics.Solid line: theoretical line (β 1 ω 1 − β 2 ω 2 )N W . Discerpancy at large N W is ascribed to bad corresponding statistics.The histogram P (N W ) was constructed from a sample of 10 6 trajectories.Here k B T 1 = 1.5,kB T 2 = 1, ω 1 = 1, ω 2 = 5/6, corresponding to heat engine operation.The time between swaps 2 ≃ 0.65 was about half the relaxation time τ relax and N = 100 swap gates were applied.

Figure 8 .
Figure 8. Efficiency and work output as a function of number of applied pulses for a fixed operation time.Here k B T 1 = 3/2, k B T 2 = 1 (ini units of ω 1 ).The operation time of the machine is fixed and equal to T op = 30τ relax (τ relax is defined here as the longest among the thermal relaxation times, i.e., τ relax = γ −1 max[e β1ω1 − 1, e β2ω2 − 1]).At N = 10 cycles, the engine has plenty of time to relax to equilibrium, because each pulse is followed by a rest time of 3τ relax .By increasing the pulse frequency one can greatly enhance the work output.Solid line, ω 1 = 1, ω 2 = 0, 7, corresponding to efficiency η = 0.3.Dashed line: ω 1 = 1, ω 2 ≃ 0.83 corresponding to efficiency at max power η * =≃ 0.17.

Figure 9 .
Figure 9. Scheme of two possible experimental set-up.Top Panel: Two Cooper Pair Boxes are coupled by means of two Josephson junctions in parallel as from Ref. [59].Bottom panel: Two flux qubits are coupled via a third flux qubit as from Ref. [62].Both setups can be used to realise the iSWAP gate, Eq. (44).The two qubits exchange photons each with a different resistor kept at a given temperature.Two on-chip fast calorimeters detect single photon emission/absorption in each resistor.