Monitoring quantum Otto engines

Unlike classical systems, a measurement performed on a quantum system always alters its state. In this work, the impacts of two diagnostic schemes to determine the performance of quantum Otto heat engines are compared: In one scheme, the energy of the engine's working substance is measured after each stroke (repeated measurements), and in the other one, the energies after each stroke are recorded in one or two pointer states and measured only after the completion of a prescribed number of cycles (repeated contacts). A single pointer state suffices if one is only interested in either work or heat. For joint work and heat diagnostics, two pointers are needed. These schemes are applied to Otto engines, whose working substance consists of a two-level system. Depending on the engine protocol, the duration of a single cycle may be infinite or finite. Because in the repeated contact scheme, the number of measurements is drastically reduced compared to the repeated measurement scheme, the quantum coherence after and during the contact diagnostics is much better maintained than repeated measurements that destroy any coherence at the end of each stroke. We demonstrate that maximum power, reliability, and efficiency of the engine in the presence of repeated contacts typically outperform these figures of merit of repeated measurements. Due to the improved coherence persistence, heat engines with a finite cycle duration require a larger number of cycles to reach a periodically asymptotic state. Overall, our results document the importance of taking into account the particular nature of diagnostic tools for monitoring and testing purposes but also for feedback control, both in theory and experiment.


I. INTRODUCTION
Heat engines have not only been important devices from a technological perspective, but they also played a significant role in establishing and examining the laws of thermodynamics [1]. Remarkably, with the extension of thermodynamics from macroscopic to ever smaller systems during the past decades, the interest in the same paradigmatic model systems such as Carnot and Otto engines has revived, with efficiency and power as figures of merit. A particular aspect of small engines results from the generic randomness of their dynamics due to the coupling of the microscopically small working substance with heat baths at different temperatures [2]. Already shortly after the invention of the maser the similarity of its action principle with a heat engine was noted [3]. Apart from scattered early considerations [4][5][6], it took about forty years until this concept gained broader interest and the idea of quantum Carnot and Otto engines was further elaborated [7,8]. Since then basically all of the distinguishing features of quantum systems have been investigated in view of their potentially beneficial or possibly detrimental aspects for an engine, such as coherence [9][10][11][12][13], discreteness and spacing of energy levels [14][15][16], quantum correlations and entanglement [17,18], Fermi versus Bose statistics [19], and projective measurements as an energy source [20,21]. Experimental investigations of quantum heat engines were suggested in Ref. [22] and realized using trapped ions [23,24], optomechanical systems [25], nitrogen-vacancy centers [26], and NMR techniques [27].
So far, the monitoring of engines has found only relatively little attention, be it as a merely diagnostic tool or as a basic element of a feedback control device [28]. The functioning of four stroke Otto-like engines during a single cycle was tracked in Refs. [21,29,30] by projective energy measurements of the working substance immediately after the completion of each stroke. From these four measured values the energy changes during work strokes and thermalization strokes corresponding to work and heat, respectively, were characterized by their joint statistics. The reduction of the engine's working substance state with respect to the instantaneous energy basis specifying each measurement, in general alters the subsequent dynamics and hence may influence the resulting energy differences between subsequent measurements. This total suppression of coherence is specific for the two point projective measurement scheme of work [31][32][33] and has often been interpreted as a deficiency [34,35]. On the other hand, the traditional assignment of average work done on a thermally closed system by a control parameter change as the difference of the average energies calculated after and before the forcing [4,36] cannot be extended to the definition of a proper random variable [37][38][39][40].
Only if the interference of the system with a particular monitoring device is specified, e.g. in terms of two projective energy measurements in the two point measure-ment scheme or by two subsequent interactions of the system with the same pointer state of a measurement apparatus [38], or by any other generalized measurement scheme [41], work can be properly specified as a random variable at the price that its probability distribution depends on the chosen measurement procedure. This principle restriction still leaves the possibility to look for detection schemes that optimize particular features such as the accuracy or the amount of back-action on the system, to name just two possibly complementary criteria.
For the diagnostics of a reciprocating engine such as an Otto engine, the monitoring of the energy of the working substance several times during a single cycle is crucial in order to be able to infer the total performed work and the consumed heat during a single, and also over many cycles. With projective energy measurements any coherence possibly built up during a work stroke [11] or after a thermalization stroke [42] would be completely removed with an impact on the performance of the subsequent strokes. The coherencies of the bath have been identified as a potential boost to the performance of an engine [9], while those built up in the working substance have been found as potentially detrimental [8,11,43,44].
In the present paper we adapt a recently proposed method to determine the sum of an observable by means of repeated contacts [45] to find the total work and the heat supplied to an engine in a prescribed number of cycles. We compare the obtained results with those from repeated measurements and demonstrate that the repeated contact method is in general less invasive and also leads to a better performance of the engine.
We organize the paper as follows: In Sec. II we introduce the general quantum Otto engine and formulate the repeated contact and repeated measurement schemes. General expressions for the joint work and heat probability density functions (PDFs) are presented for both schemes. In Sec. III, a two-level system, i.e. a single qubit, is considered as the working substance for which general expressions for the work strokes are presented. For the heat strokes we consider perfect thermalization beyond weak system-bath coupling and alternatively finite time thermalization at weak system-bath coupling. We compare our analytic expressions of work and heat PDFs as well as the second moments resulting from the two measurement schemes for a single engine cycle. In Sec. IV numerical results for several cycles are presented. Finally, in Sec. V, we draw our conclusions.

A. Mode of operation
For the sake of concreteness we consider an Otto engine with a quantum mechanical working substance, which later will be specified as a two-level system. The subsequent four strokes constituting a full cycle of the engine consist in compressing, heating, expanding, and cooling Each pointer is initialized in the same state σ. The measurements are taken at the beginnings and the end of the work strokes when the working substance is thermally isolated. The results of the measurements are symbolically presented at the right margin. The registered pointer states can be used to calculate heat and work. The lower panel illustrates the repeated contact scheme with two pointers, the upper one for the work and the lower one for the heat. Again both pointers are initialized in the same state σ. Here the actions of the contacts on the two pointers are accumulated with appropriate signs, see the Sec. II B for the details, and read out by a projective position measurement after a prescribed number of engine cycles is completed. Between subsequent contacts the states of the two pointers do not evolve. the working substance, see Fig. 1. The compression and expansion strokes are caused by changing a control parameter of the working substance in thermal isolation. Therefore, the energy change of the working substance during these strokes can be referred to as work. The parameter change in the compression stroke is designed such that the energy level distances of the working substance Hamiltonian increase. Hence, the time evolution of the density matrix caused by the compression is determined by a unitary operator U with ρ 2 = U ρ 1 U † , with ρ k = ρ(t k ), k = 1, 2, · · · , denoting the density matrices at the times indicated in Fig. 1, whereby times t k with k = 1, · · · , 4 specify instants within the first cycle, with k = 5, · · · , 8 within the second cycle, etc.. For the sake of simplicity, the expansion is assumed to be the time-reversed protocol of the compression. Therefore, the unitary operator giving the corresponding transformation can be written asŨ = K * U † K with K * being the adjoint of the antiunitary time-reversal operator K (K * K = KK * = 1), [46], such that ρ 4 =Ũ ρ 3Ũ † . The relation betweenŨ and U † holds as long as K * H(t)K = H(t) for all time t during the work strokes with H(t) being the time-dependent engine Hamiltonian. During the second and the fourth stroke the working substance is brought into contact with a hot and a cold heat bath, respectively, while its Hamiltonian remains unchanged. These two heat strokes are characterized by trace preserving, positive, linear maps Φ h and Φ c relating the respective density matrices of the working substance at the beginnings and the ends of these strokes according to ρ 3 = Φ h (ρ 2 ) and ρ 5 = Φ c (ρ 4 ). Later, we shall give explicit expressions specifying the unitary operator U ,Ũ and the maps Φ h and Φ c for the case of a two-level system as working substance.

B. Diagnostic means I: Repeated measurements
The central quantities characterizing the performance of an engine are the total work W performed on and the total heat Q supplied by the hot bath to the working substance during N cycles. Both quantities follow from the sequence of energies e = (e 1 , e 2 , · · · , e 4N ) of the working substance taken at the beginning and the end of each work stroke when the working substance is decoupled from either heat bath, see Fig. 1 yielding where µ 1+4l = µ 4+4l = 0, µ 2+4l = −1, µ 3+4l = 1 with l = 0, · · · , N − 1.
In the repeated measurement scheme the energies e k are measured individually by bringing the working substance for a short time into contact with the pointer of a measurement apparatus [47]. This contact causes a unitary transformation of the joint state of the pointer and the working substance given by where H k is the working substance Hamiltonian at the start of the kth stroke, i.e. H 1 = H 4 = H 5 = · · · = H c and H 2 = H 3 = · · · = H h , see also Fig. 1; furthermore, P k denotes the momentum operator conjugate to the pointer coordinate of the kth measurement apparatus. The parameter κ results as the product of the strength and the duration of the contact which is assumed to be so short that the free motion of the working substance appears as frozen. Above, and throughout this work we will set = k B = 1. Before the kth interaction, the working substance and the respective pointer are uncorrelated. The latter resides in a Gaussian pure state σ k with vanishing mean value and variance Σ 2 . Its position matrix element σ(x, y) = x|σ k |y is independent of k and is given by If, after the contact has taken place, the pointer position is measured projectively at a position x, the nonnormalized density matrix of the working substance becomes where e m k and P m k denote the eigenvalues and eigenprojectors of the Hamiltonian H k = m e m k P m k and ρ denotes the density matrix of the working substance before the contact. For the sake of simplicity and without restriction of generality, the pointer position is measured in units of energy such that κ = 1. The subsequent measurements of the energies e k during N cycles, giving rise to the pointer positions x N = (x 1 , · · · , x 4N ), lead to a mapping of the initial working substance density matrix ρ onto the non-normalized density matrix φ RM x N (ρ) given by k labelling the eigenstates of the Hamiltonian H k . The operator D m N , m N (ρ) describes the action of the alternating timeevolution operators characterizing the subsequent strokes and the interrupting projections due to the energy measurements on the initial density matrix ρ. They are recursively given by where the sequences m x (ρ) one obtains for the joint PDF p( x)to find the pointers at the values given by the components of the vector x the formal expression For Gaussian pointers as specified in Eq. (4) this PDF can be written as where ] denotes a Gaussian PDF with vanishing mean value and variance Σ 2 . The coefficient D m, m is defined as In order not to overburden the notation we omitted the index N in the vectors m ( ) and x.
Identifying the components x k as the measured energies, one obtains for the joint probability of work and heat, P RM (W, Q), the expression Going to the characteristic function G RM (u, v) = dW dQe iuW e ivQ P RM (W, Q) one can perform the 4Nfold x-integration to find with µ k as defined below Eq.
with δX = X − X (X = W, Q). Hence, the PDF is given by The more pronounced broadening of the Gaussian contributions in the marginal work PDF compared to the marginal heat PDF is caused by the twofold number of energy measurements required for the work compared to the heat.

C. Diagnostic means II: Repeated contacts
Instead of registering each of the 4N energy values in a separate pointer state one may follow the idea of Ref. [48] and subsequently transfer the actual energy value at the beginning of each stroke with the appropriate sign according to the Eqs. (1) and (2) to two pointer states, one for work and the other one for heat. The transfer of information from the working substance to the pointer states is mediated in a similar way as for repeated measurements by a sufficiently short interaction causing a unitary transformation V k of the form acting on the joint Hilbert space of the working substance and the two pointer states. Here P W and P Q denote the momentum operators that are canonically conjugate to the position operators of the work and heat pointers, respectively. In contrast to the repeated measurement scheme, the pointers are only read out by a projective position measurements after the completion of N cycles. Again, the pointers are assumed to be gauged in units of energy. Between two subsequent contacts the working substance does not interact with the two pointers but their states remain correlated as a result of the previous contacts. To simplify the further analysis we assume that the pointers have no own dynamics, i.e. they change their states only due to the interactions as specified by the unitary operators defined in Eq. (20). The subsequent strokes and contacts, resulting after the completion of N cycles in pointer state positions at specific values W and Q, lead to a non-normalised reduced density matrix of the working substance, which is given by where ρ denotes the initial density matrix of the working substance and individual work and heat outcomes are denoted as with µ k as defined below Eq. (2). The joint PDF P RC (W, Q) of finding these work and heat values is given by the trace of φ RC W,Q (ρ) yielding As for repeated measurements it is a linear superposition of Gaussian PDFs with weights D m, m , however, with the following differences: 1. for repeated contacts the contributing Gaussians for work and heat factorize while for repeated measurements they are correlated; 2. the widths of the contributing Gaussians are independent of N given by Σ whereas those for repeated measurements are proportional to √ N Σ; For repeated contacts with a non-selective heat measurement the marginal PDF P RC,2P (W ) defined as differs by the extra factor of exp[−(Q m − Q m ) 2 /(8Σ 2 )] from the PDF resulting from a mere work measurement performed with a single pointer state. In the latter case one obtains For repeated contacts the same type of contextuality, i.e. the dependence of the observed results on the presence of a non-selective measurement, also holds for the marginal heat PDF. In both cases it is more pronounced for rather imprecise registration, i.e. for pointer states with relatively large variances Σ 2 . For repeated measurements the joint work and heat PDF [see Eq. (17)] results according to Eq. (12) in a classical way from the joint PDF p( x) of all individual energies. Therefore repeated measurements do not present this particular contextuality feature.
Once the work and heat PDFs are known for a specified diagnostic method, the performance of the engine can be characterized by its efficiency η and reliability R which are defined as where · denotes the average of the indicated quantity with respect to the PDF of the considered diagnostic method. Higher order efficiencies as introduced in Ref. [49] can be obtained from higher order moments of work and heat. Also a fluctuating efficiency, W/Q, can in principle be determined [50], but will not be considered here.

A. Work strokes
The general formulation outlined in Sec. II can be applied to any working substance undergoing adiabatic or non-adiabatic protocols. For a numerically feasible but yet representative example, we will use in this and subsequent sections a two-level system as the working substance and allow for non-adiabatic work strokes, i.e. strokes comprising shifts of the instantaneous energy eigenstates and transitions between them during the unitary evolution. Without specifying the time-dependent Hamiltonian in detail, the system Hamiltonians H c and H h before and after the first work stroke, respectively, are written as where |− u and |+ u denote the ground and the excited states, respectively, of the Hamiltonian H u with u = h, c.
The most general form of the unitary evolution, during the work stroke, then reads AboveŨ ≡ C * U † C is the evolution operator for the timereversed protocol with the antiunitary operator K chosen to be the complex conjugation operator C. Here α ∈ [0, 1] is the transition probability between the ground and the excited states and ϕ is the corresponding phase (both parameters depend on the stroke protocol and in particular on its duration). The limiting case α = 0 refers to an adiabatic process while α = 1 corresponds to a swap. For example, a Landau-Zener protocol is described by the time-dependent Hamiltonian where the velocity of the drive, v = 4 2 h − 2 c /T 1 , is chosen such that the eigenvalues of the Hamiltonian change from ± c at t = 0 to ± h at t = T 1 /2, the time corresponding to the duration of a work stroke. The parameters specifying the unitary operator U are given in Ref. [51] by where Γ(z) denotes the Gamma function. For other protocols characterized by a single parameter see Ref. [52] and for more general cases Ref. [53,54].
In the first instance we shall not refer to a particular work protocol and use the general form of a unitary time evolution U as specified by Eq. (30). Before we do so, we consider the heating and cooling strokes of the two-level working substance that are achieved by a coupling to the heat reservoirs at temperatures T h and T c , respectively. Depending on the duration and the strength of the respective contacts perfect or imperfect thermalization of the two-level system can be reached.

B. Perfect thermalization
Provided that the contact between the working substance and the heat bath is sufficiently long, the finally reached reduced state of the two-level system is always the same, dependent on the temperature of the heat bath and the strength of the interaction but independent of its state ρ at the beginning of the contact. Hence, the operation Φ u , u = h, c describing a perfect thermalization stroke is given by the projection onto the final state ρ β and can be written as Only for weak coupling ρ u , u = h, c coincides with the Gibbs state ρ u ∝ exp[−β u H u ] defined by the Hamiltonian of the two-level system at the inverse temperature β u . While these states lack coherence the latter may be found in generalized Gibbs states resulting from the partial trace of the total system's Gibbs state with respect to the bath [40,42]. For a perturbative treatment in the coupling strength up to second order see Append. A. In the asymptotic limit of an initially sharp pointer state with Σ = 0 the joint work and heat PDFs charac-terizing a single cycle are identical for repeated measurements and repeated contacts (two-pointer scheme), i.e., P RM 1 (W, Q) = P RC 1 (W, Q). The first two moments can be expressed as with Here, d c and d h denote the populations of the exited states of the cold and the hot density matrices, respectively, see the Eqs. (A4) and (A5) in Append. A. Even when the Gaussian pointer has a finite width Σ = 0, the differences between repeated measurements and contacts in the average values of work and heat are of the order exp[− 2 c(h) /(2Σ 2 )] which could be extremely small for a precise measurement apparatus, i.e. if c(h) < Σ. On the other hand, the second moments differ by which is dominated by the 3Σ 2 term arising due to the difference in the number of measurements. Accordingly, also the reliability [see Eq. (28)] differs for repeated contacts and repeated measurements. Overall, however, even though the two measurement schemes are quite different they do not show major significant differences for a perfectly thermalizing quantum Otto engine.

C. Imperfect thermalization
Perfect thermalization requires long times during which the working substance is coupled to one of the two heat baths, leading to long cycle periods and correspondingly small power output of the engine. The dynamics of the relaxation process in general is very complicated. A substantial simplification occurs only at week coupling in the Markovian limit in which the dynamics is governed by a Lindblad equatioṅ , and the Lindblad jump operators L ± u = |± u ∓ u |. The Hamiltonians H u are given by Eq. (29). Even though the stationary solutions of the master equations for the hot and the cold heat baths are given by the Gibbs states, corresponding to the Hamiltonians H u at the respective temperatures β u and hence do not possess any coherence, the solutions of the master equations at any finite time in general do contain coherences. This leads to differences between repeated measurements and repeated contacts even in the limit of infinite precision (Σ → 0). For a single cycle, starting at some initial state c.] one finds the following expressions for the first two moments of work and heat for a small pointer state variance including leading order O(Σ 2 ) corrections: where θ = h τ h with τ h being the duration of the hot bath stroke. A dependence of the moments on the initial coherence parameter q is only seen in higher order Σ 2 corrections. From the expressions of these moments the according variances and the covariance of work and heat can be obtained. Instead of these bulky expressions we present the variances and covariances for the two monitoring schemes in Fig. 2a, b as functions of θ. The difference between repeated measurements and repeated contacts is governed by damped oscillations stemming from the e −γθ cos (2(θ + ϕ))-term in the deviation of the first two work moments, see Eqs. (40), (42) and Fig. 2c, d. Yet, the maximal amplitude of these oscillations is rather small rendering the (co)variances for repeated measurements and repeated contacts almost identical. The moments of heat can be obtained by setting c = 0 and multiplying appropriate global signs (negative for odd order moments and positive for even ones) in the expressions for the moments of work. As seen from the above equations, this also implies that the average heat for two-pointer repeated contacts and repeated measurements are the same, i.e., Q = Q RM = Q RC . The variances of the heat are not displayed. They behave in a similar way as the those of the work as presented in Fig. 2 with solely a constant difference of 3Σ 2 between repeated measurements and repeated contacts and therefore are not extra displayed. Note that the engine may starts at an arbitrary initial state specified by the parameters d and q. Here and in the following examples, however, we will initiate the monitored sequence of N cycles with the invariant state ρ * of the map F, i.e. where describes the dynamics of a complete cycle in the absence of any monitoring. Furthermore we note that the master equation (38), which governs the thermalization strokes of the engine, leads to decoupled equations of motion for the populations and coherences. Hence, the thermalization maps Φ u resulting as the solutions of the master equation (38) obey the conditions with swap-operators L ± u = |± u ∓ u |, as defined below Eq. (38), and P k u = |k u k u | with k = ±, denoting the projection operators onto the eigenstates of the Hamiltonian H u . Due to this decoupling, the work PDFs for oneand two-pointer repeated contacts become identical and we denote the average value as W RC (see Append. B for more details). In the limit θ → ∞, Σ = 0, and ρ = ρ c , which leads to perfect thermalization, Eqs. (39)-(44) match Eqs. (35)- (37). For a highly non-adiabatic work stroke (α large) the cosine terms in the expressions for the first two moments play a significant role. Then the sign of cos(2(θ + ϕ)) determines which one of the two monitoring schemes leads to a larger work on average for a single cycle. This affects also the efficiency of the engine because the average heat is the same for either scheme.

IV. NUMERICAL RESULTS
For a working substance with a d-dimensional Hilbert space, the number of possible realizations of ( m, m ) increases exponentially as d 8N , and so do the coefficients D m, m . Hence, the computational resources (CPU time and RAM) required to track the individual coefficients become inaccessible even for d = 2 already at a mod- In both of the plots, the differences, which are almost imperceptible in the panels a and b, can be divided into two parts. The first part is the constant shift proportional to Σ 2 in Eqs. (42) and (44). The remaining nontrivial differences (NTD), which are plotted in panels c (for δW 2 ) and d (for δW δQ ), result from the terms proportional to cos(2(θ + ϕ)) in the Eqs. (40) and (42). The engine parameters are as follows: The initial half energy gap c = 1, the final half energy gap h = 3.7, and the Gaussian width of the measurement is Σ = 0. erately large number of cycles. To avoid this numerical impossibility, we device a scheme to reduce the computational complexity to N 2 by grouping terms according to their work values. The scheme (see Append. C for details) is restricted to a two-level working substance undergoing thermalization strokes that are governed by decoupled population and coherence dynamics as specified by Eq. (47). Under these conditions the method is exact and allows us to obtain the work distributions for a large number of cycles such that the asymptote can be reached.
FIG. 3. Panels a-b display the probability density functions (PDFs) of work for adiabatic work strokes (α = 0) and different measurement schemes (a) as well as for different nonadiabatic work strokes (b). Solid lines represent two-pointer repeated contacts PDFs, and the dashed line presents the PDF for repeated measurements. One-pointer results are not drawn since they are exactly the same as two-pointer scheme for α = 0 and have a negligible difference for α = 0. Panel c displays the joint PDF of work and heat calculated for repeated measurements, and d is the joint PDF calculated for two-pointer repeated contacts. The PDFs are computed for perfect thermalization with thermal baths at inverse temperatures βc = 0.25, β h = 0.025, both characterized by the the same Ohmic spectal density with Lorentz-Drude cut-off frequency ωD = 0.2, system-operator S = σz + σx, and damping coefficient γ = 0.5. This results in density matrices with an excited-state probability dc = 0.37759 and d h = 0.45388, and coherence qu = +u|ρu|−u given by q h = −1.9205 × 10 −6 and qc = 5.0813 × 10 −5 (see Append. A for details). The phase during the non-adiabatic evolution is set to ϕ = 0 in all cases, and the joint PDF is shown for α = 0.1. The engine is started in the generalized Gibbs state of the cold bath with parameters dc and qc. All other parameters are the same as in Fig. 2. A. One cycle, perfect thermalization Figure 3 presents the joint and marginal work and heat PDFs of an engine with perfect thermalization for a single cycle (N = 1). Panel a displays the work PDFs for repeated measurements and repeated contacts for adiabatic work strokes, i.e., for α = 0. In this case, there is no difference between the marginal two-pointer work PDF given by Eq. (25) and the one-pointer work PDF (26) as explained in Append. D. In accordance with the Eq. (26) the PDF consists of three peaks located at W = 0 and W = ±2( c − h ) because transitions between different energy levels are impossible for adiabatic work strokes. These peaks are broader for repeated measurements compared to repeated contacts. This feature prevails for nonadiabatic work strokes. In panel b we present the results for moderately and strongly non-adiabatic work strokes. In order not to clutter the plot we present only the work PDFs resulting from repeated contacts. With decreasing adiabaticity (increasing α) six further peaks appear in accordance with Eq. (26) accompanied by a decrease of the weights of the three adiabatic peaks. At the same time the weights of those peaks at positive work values may increase to such an extent that the average work becomes positive. That means that the engine only consumes energy rather than to perform work. As for classical engines any non-adiabaticity is detrimental with respect to the efficiency but must be accepted to achieve finite power. The same features as described for repeated contacts also hold for the repeated measurement PDFs with the only difference that their peaks are broader.
Figures 3c and d present the joint work and heat PDFs for repeated measurements and contacts, respectively. In the case of repeated measurements, each individual Gaussian contribution contains a negative covariance indicating an anti-correlation of the according work and heat contributions. In contrast, for repeated contacts the individual Gaussian contributions are isotropic. In both cases the main contribution to the covariance of work and heat results from the linear superposition of the individual Gaussian distributions.

B. N cycles, imperfect thermalization
Numerical results for adiabatic and non-adiabatic work strokes and imperfect thermalization are presented in Fig. 4. The panels a and b display the work PDFs for repeated measurements and repeated contacts for two different numbers of cycles N . For the smaller number of cycles, N = 5, still several peaks are visible in the repeated measurements work PDF. Because their widths grow as √ N this multiply peaked structure disappears with increasing number of cycles to finally approach a Gaussian shape, which in the present case is virtually indistinguishable from the exact PDF for N = 30. In contrast, the fine-structure of the repeated contact work PDF remains visible also for large numbers of cycles due to the fact that the widths of the individual Gaussian contributions are solely determined by the initial pointer state and hence are N independent.
In Fig. 4c the average work per cycle is displayed as it results as a function of the number of cycles both for repeated measurements (red closed circles) and repeated contacts (blue crosses). The initial state of the working substance is chosen as the asymptotic state of an unobserved engine reached after sufficiently many cycles. Only for a single cycle (N = 1) the average work for repeated measurements is larger than that for repeated contacts. For any larger number of cycles the repeated contact scenario outperforms the repeated measurements both with respect to the average work (see panel c) and the reliability (see panel d). As a function of the number of cycles, the work behaves like a biased random walk in that its first and second moments asymptotically grow in proportion to the number of cycles yielding a constant work per cycle and a reliability that increases as the square root of the number of cycles.
For repeated measurements a map can be constructed that assigns the density matrix of the working substance after a cycle to any initial density matrix. The map possesses a unique invariant state describing the asymptotic regime of a large number of cycles. This state can be used to calculate the asymptotic average work per cycle, see Append. B. Even though such a map of the reduced working substance density matrix does not exist for repeated contacts, in the special cases of the Lindblad equation (38) that dynamically decouples populations and coherences, the asymptotic value of the average work can be determined from the asymptotic density matrix in the absence of any contacts or measurements, see Append. B. This unobserved asymptotic density matrix of the working substance can again be obtained as the invariant state of a quantum map in the absence of contacts as considered in Ref. [55]. The resulting average works are presented in panel c as dashed black lines.

C. Power of a Landau-Zener engine
So far we have characterized the two work strokes by the Hamiltonians of the expanded and the compressed working substances [see Eqs. (32) and (33)] and unitary operators specifying the changes of its state, however, without specifying a particular forcing protocol. As a consequence, the duration of the work strokes is not specified. In general, determining the unitary time evolution operator for a particular forcing protocol of given duration presents a nontrivial task [52,54]. In order to also have full control over the duration of an engine cycle we consider a forcing protocol with a linear driving as specified in Eq. (32). The unitary time evolution for a work stroke of the duration T 1 /2 is then determined by the parameters as given in the Eq. (33). The duration of the thermalization strokes, τ h and τ c are assumed to be fixed by a the single parameter θ = c τ c = h τ h , resulting in the total thermalization time Therefore, a single cycle takes the time T 1 + T 2 . The power of the engine is given by where W denotes the average work performed on the working substance per cycle. In Fig. 5 the power and the efficiency for repeated contacts and repeated measurements as well as for different numbers of cycles are presented as functions of the duration of the work and thermalization strokes. The variance of the pointer states is chosen relatively small such that coherences that may build up during the nonadiabatic work strokes and also partially during an imperfect thermalization, are almost completely suppressed by repeated measurements while they can partly survive in the presence of repeated contacts. Because the generation of coherences is further suppressed in the limit of adiabatic work processes as well as for perfect thermalization, one finds substantial differences between repeated measurements and repeated contacts for small values of T 1 and small to intermediate values of T 2 . As already observed in Sec. IV A, the total average work performed on the working substance in a single cycle becomes positive for perfect thermalization and work strokes with a sufficiently large parameter α and hence for a sufficiently short duration of the work strokes. This feature of a dud engine with formally negative power and efficiency persists for imperfect thermalization and also for any numbers of cycles at sufficiently small values of the work stroke duration T 1 indicated as white regions in  The dependence of the maximal power P * = max T1,T2 P on the number of cycles for repeated measurements and contacts is displayed in Fig. 6a. This maximum is chosen from the set of all engines with variable stroke times T 1 and T 2 and with all other parameters fixed. In both cases, the power approaches an asymptotic value at large numbers of cycles. While the asymptotic limit is reached from above already after a few cycles for repeated measurements, it is approached considerably slower from below for repeated contacts. The asymptotic values are indicated as broken horizontal lines. In both cases, these values can be determined from the invariant maps that propagate the density matrix cycle by cycle. For the repeated contact scenario, the map describes the dynamics in the absence of any monitoring, while for the repeated measurement scenario it describes the singlecycle engine dynamics in the presence of non-selective energy measurements at the times indicated in Fig. 1. For a more detailed description we refer to Append. B.
For the considered parameter values of the engine, the position where the maximum power is assumed in the T 1 -T 2 plane turns out to be almost independent of the number of cycles. Therefore, the convergence speed of the map as it is iterated, mainly determines the convergence of the power with increasing number of cycles as it can be seen in Fig. 6a. For sufficiently many cycles the convergence speed is ruled by the second largest eigenvalue of the respective map, while its largest eigenvalue is given by Λ 1 = 1. It is non-degenerate for all engines with a positive thermalization time T 2 > 0. The dependence of the second largest eigenvalue Λ 2 as a function of T 1 and T 2 is exemplified in the panels Fig. 6b,c. In both scenarios the second eigenvalue converges for large thermalization times to zero, however much faster for repeated measurements than for repeated contacts in agreement with the faster convergence of the maximal power towards its asymptotic limit, i.e., Λ RM 2 < Λ RC 2 . Outside the regime of short work strokes, the second eigenvalue is almost independent of T 1 as it is primarily determined by the duration and strength of the thermalization strokes.

V. CONCLUSIONS
In this paper, we compared two different monitoring schemes, both designed to record the statistics of work performed on and heat supplied to a quantum Otto engine over a specified number of engine cycles. Based on these statistics, we investigated the performance of such engines running subject to different modes of operation, such as different work strokes ranging from adiabatic to sudden compression and expansion of the working substance, and different types of heating and cooling strokes. The method was formulated for general working substances and exemplified for a two-level system acting as working substance.
Apart from the joint and marginal work and heat PDFs, we studied performance measures such as efficiency, power and reliability that are all determined by the first two moments of work and heat. For an engine that performs N cycles, the two monitoring schemes are based on contacts of the working substance with one, two, or 4N pointers, depending on the particular scheme.
These engine-pointer contacts are supposed to be so short that the dynamics of the working substance is frozen during these moments but effective enough to induce a change of the pointer state depending on the actual energy of the system.
For the repeated measurement scheme each contact is immediately followed by a projective measurement of the pointer state whose value corresponds to the energy of the working substance at the instant of the contact, hence representing a generalized energy measurement according to the von Neumann scheme. The precision of the measurement can be controlled by the position variance Σ 2 of the initial pointer states which we supposed to be Gaussian and identical for all employed pointers. From the joint statistics of all measured energy values, the joint and marginal PDFs of work and heat can be found. These PDFs result as linear superpositions of Gaussian contributions with centers at work and heat values given by proper sums and differences of the 4N energy values but also fictitious values resulting from a double set of energies. The latter though are exponentially suppressed provided the variance Σ 2 is sufficiently small. All contributing Gaussians are of the same shape with a variance growing proportionally to the number of cycles leading to a smeared-out, almost Gaussian PDF for sufficiently large numbers of cycles.
For repeated contacts there are just two pointers one being associated to work and the other one to heat. These pointers are "kept on" during the whole monitored period of cycles. At each contact the pointer state is modified according to the present energetic state of the working substance, left unchanged during the strokes and read out by projective position measurements only after the completion of N cycles. The energetic states at the subsequent contacts are inscribed with alternating signs in the work pointer. For the heat pointer only those contacts before and after a heat stroke are registered with appropriate signs. It is interesting to note that the marginal distributions of work obtained from the joint work and heat PDF in general differ from the work PDF of a repeated contact scheme with only the work pointer but without the heat pointer. The same holds for the heat PDFs obtained either as the marginal of the joint work and heat PDF or as the outcome of a single heat pointer repeated contact scheme. This aspect of quantum contextuality is absent in the repeated measurement scheme. It indicates a higher sensitivity of the repeated contact scheme to the quantumness of a process as compared to the repeated measurement scheme.
The higher sensitivity of the repeated contact scheme is closely related to the lesser back-action of the sole contacts as compared to the individual read-out of the pointer-state position for repeated measurements leading to different suppression factors making long-lived coherences possible. In the example of an Otto engine with a two-level system working substance and imperfect thermalization strokes this manifests itself in an average work per cycle that first increases with the number of cycles and finally saturates at a constant value for repeated contacts. In contrast, for repeated measurements the average work per cycle decreases towards an asymptotic value that is smaller than that for repeated contacts. Accordingly, the power but also the efficiency as well as the reliability of the work, defined as the ratio of the average work and its root-mean-square deviation, are larger for repeated contacts than for repeated measurements provided that the engine protocol allows coherences being still present after a cycle. Due to the less strong back-action imposed by repeated contacts in comparison to repeated measurements these coherences may subsist leading to a better performance of a repeated-contactsmonitored engine. We expect that this behaviour is not restricted to the special case of a two-level working substance but that it prevails in general.
In a recent work, Watanabe et. al. [56] studied the amount of energy a quantum Otto engine may deposit into an external system over several engine cycles. The comparison of the amount of out-coupled energy into a "flywheel" with the energy produced by the free engine according to a repeated contact diagnostic scheme may provide further insight into the effectiveness of charging and discharging strategies of quantum batteries.
The generalized Gibbs state can be obtained up to second order in system-bath interaction using canonical perturbation theory [42]. In this Appendix we omit the subscripts c or h used in the main text for the cold and hot reservoir for notational simplicity. The generalized Gibbs state can be interpreted as the reduced density matrix of an extended system governed by the total Hamiltonian H tot = H S + H B + V in a canonical state at the inverse temperature β. Here, H S and H B , denote the Hamiltonians of the proper system and the bath, respectively and V denotes the interaction. We assume that the bath is composed of independent Harmonic oscillators of frequency ω n and the system-bath interaction takes the form V = −S ⊗ n c n x n , with S being a system operator and x n being the position operator of the nth mode of the bath coupling to the system via strength c n . Performing a perturbative expansion upto second-order in V we obtain, Above Z S = Tr[exp[−βH S ]] is the system partition function of the bare system and the operator D encodes information about finite system-bath coupling, i.e., The operatorS x = e xH S Se −xH S is the imaginary time free evolution of the system operator and the force-force correlator C(τ ) = n,n c n c n Tr[e −βH Bx n (τ )x n ]/Tr[e −βH B ] with x n (τ ) = e iτ H B x n e −iτ H B . The force-force correlator can be expressed in terms of the spectral density of the bath, as For our perfectly thermalizing two-level Otto engine discussed in Sec. III we use the system operator S = σ z +σ x such that the density matrix with Z S = e β + e −β . Throughout this work we have used an Ohmic spectral density with Lorentz-Drude cutoff as J(ω) = γω/(1 + ω 2 /ω 2 D ) for all our numerical simulations.

Appendix B: Reduced maps for measured engines
In this section, we present the methods to determine the map for the average work performed in the (N + 1)th cycle based on the knowledge of the average work performed during N cycles. In case of repeated measurements, the map fully captures the distribution of work performed at any (N + 1)th cycle, whereas for repeated contacts we introduce a fictitious work distribution whose average is equivalent to the asymptotic average work per cycle. The method is explained for work performed in the last cycle but can likewise be applied for the heat absorbed in the last cycle.

Repeated measurements
First, we consider the case of repeated measurements. As described in Sec. II B, in this case after every stroke a new Gaussian pointer is adopted, brought into contact with the engine, and measured. We begin with the non-normalised density matrix for (N + 1) cycles, using Eq. (6), Utilising the relation in Eq. (7), splitting the product above into one for N cycles and the other for (N + 1)th cycle, and using the non-normalized density matrix for N cycles [Eq. (6)] we can rewrite the above equation as, The work performed in the (N + 1)th cycle is given by and its corresponding PDF is with p( x N +1 ) = Tr[φ RM x N +1 (ρ)] according to Eq. (9). Not-ing that the integration variables within the delta function depend only on the (N + 1)th cycle pointer positions x k , the integral can be split as, All information about the previous N cycles is contained in the rightmost integral above. Hence, using Eq. (B2) this integral can be expressed as, The simplification on the r.h.s. implies that as long as we know the map S m,m and the state it acts on, i.e., d 4N x φ RM x N (ρ), we can evaluate the work distribution P RM (W (N +1) ). Thus, the task reduces to finding the map F RM , such that, Using Eq. (B2) with N → N − 1 and integrating both sides we can easily obtain the map, Note the difference betweenm ( ) and m ( ) is that N is replaced by (N + 1) in the concatenation. The map F RM is Markovian since it gives us the state at the N th cycle given the state at N −1 [see Eq. (B6)]. Overall, in case of repeated measurements since the map F RM is known it is straightforward to obtain the work PDF for the work performed in the (N + 1)th cycle using Eqs. (B5)-(B7).

Repeated contacts
In case of repeated contacts, one cannot assign a value of work to each individual cycle as we did with W (j) for repeated measurements by using the classically stored information about the outcomes of the sequence of energy measurements for N cycles. To circumvent this problem we restrict our consideration to engines with a working substance consisting of a two-level system and moreover being subject to thermalization strokes obeying the decoupled population and coherences dynamics specified by the Eq. (47). In this case we are able to introduce a fictitious work PDF that can be computed recursively and that determines the average work per cycle in the limit of large cycle numbers when the engine has reached a steady mode of operation.
In order to proceed, we begin with the work PDFs defined in Eqs. (25) and (26), In case of dynamical decoupling between populations and coherences for a two-level working substance, [Eq. (47) Thus, for the models considered in Secs. IV B-IV C the one-and two-pointer work PDFs always coincide. Moreover, using the above simplifications one can easily identify the non-normalized reduced density matrix for the working substance as, Here we also omitted the only remaining con- to the exponential factor vanish in the final trace operation applied to obtain the work PDFs from the non-normalized reduced density matrices. In addition, we made use of the fact that the innermost part of the telescopic action of the map D m, m on ρ is given by P m1 1 ρP m 1 1 , see Eq. (7). Splitting the density matrix ρ into its diagonal part ρ d = m=1,2 P m 1 ρP m 1 and the remaining coherence contribution ρ c = ρ − ρ d one obtains exp −(e m1 1 − e Note that the map ρ →ρ preserves the trace and the positivity.
Next we introduce the fictitious work distributioñ which exactly resembles Eq. (B5) with, The average work from the fictitious distribution, represents the difference of the average work performed by the engine in the (N + 1)th and N th cycle, i.e., W (N +1) = W N +1 − W N ( this relation follows from the definitions of the individual quantities as explained below). As the monitored engine reaches its asymptotic state, W (∞) represents the average work performed by the engine per cycle, similarly to the repeated measurements case, because the monitoring process on average has the same effect on the engine state in the N th and (N + 1)th cycle after sufficiently many cycles.
The average work performed in (N +1) and in N cycles can be written as The third line holds because the non-selective map F RC that is defined as preserves the trace so that where we shifted the integrand by the amount of W m, m and used the normalization of the Gaussian distributions. On the other hand, .
Integrating x first, we get Furthermore, using a fictitious non-normalised density matrix, similar to Eq. (B1), with the fictitious pointer statẽ where F RC is defined in Eq. (B17). Comparing the above equation with Eq. (B7) we notice the missing exponential terms in the map. This is simply because the exponential terms in case of repeated contacts get simplified, see Eq. (B12), and they only affect the initial state of the engine. In other words, in case of repeated measurements to obtain the state of the engine after N cycles we apply the map F RM N times on the initial density matrix ρ whereas in case of repeated contacts the map F RC is applied N times on the fictitious initial density matrixρ. As trace and positivity preserving maps, F RC and F RM possess at least one invariant solution. The uniqueness of the invariant solution is guaranteed by the dissipative nature of the thermalization strokes. As a conseqence, the invariant densities can be iteratively reached as ρ * = lim N →∞ (F X ) N (ρ), where X = RC, RM . Accordingly, the spectra of these maps contain one non-degenerate eigenvalue Λ 1 = 1 and the absolute values of all other eigenvales are smaller than one.
In Fig. 6b-c, the second-largest eigenvalues Λ 2 of the maps F RC(RM) are displayed as functions of the durations T 1 and T 2 . When Λ 2 is close to 1, the engine needs a long time to converge to the asymptotic cycle and vice versa.
2i )] which is zero due to Eq. (47). Hence, for case (i) only when Eq. (B9) is satisfied we have a contribution to the work distributions.
For case (ii), we divide into three sub-cases: (iia) m 2i+1 = m 2i+1 for which A = Tr[P This exhausts all possible cases and overall only when Eq. (B9) is satisfied we get a contribution to the work distribution.

Appendix C: Reduction in computational costs
In this section, we device a computational scheme to reduce the complexity from exponential (d 8N with d being the dimension of working substance Hilbert space and N being the number of cycles) to quadratic (N 2 ) as proposed in the main text. Our proof below focuses on the work PDFs given by Eqs. (18) and (26), Above we have reintroduced the subscript N indicating the cycle number which was suppressed in the main text.
The scheme for the heat PDFs can be derived in exactly the same way. The above equations can be rewritten as, (C4) In case of repeated contacts, the general expression reduces to the above due to the dynamical decoupling between populations and coherence [see Append. B 2 Eq. (B12)]. The reduction in computational resources occurs because we only need to compute the coefficients for a fixed W N that scale as N 2 rather than the exponential scaling of d 8N .
In order to see this reduction in computational resources, let us apply the map S m,m , defined in Eq. Here we derive the scaling reduction for the case of repeated contacts and the proof for repeated measurements can be obtained in the same spirit. The exponential terms appearing in Eq. (C4) have been absorbed in the effective initial density matrixρ [see Eq. (B13)] since they only depend on the first contact (m 1 and m 1 ), hence, only alters the initial state. It is important to note here that our proof works only if the exponentials appearing in Eq. (C1) can be expressed as a product of exponentials implying that such a dramatic reduction is always possible for repeated measurements and is impossible if the populations and coherence do not dynamically decouple in the case of repeated contacts. Using Eq. (7) (C7) Therefore, the operator D RC W N +1 (ρ) is obtained recursively from D RC W N (ρ) in terms of a sum over m, m and W N . The required number of computations does not scale with N (for m ( ) we have 2 4 terms for d = 2) and the sum over W N scales as N 2 . Thus, the computations performed, due to the grouping with fixed work values, reduces from an exponential scaling of 2 8N (d = 2) to a quadratic N 2 .
Appendix D: Equivalence between one-and two-pointer work PDFs for a perfectly thermalizing adiabatic engine In this section, we show that for a perfectly thermalizing adiabatic engine diagnosed via repeated contacts the one-and two-pointer work PDFs are equal independent of the dimensionality d of the Hilbert space of the working substance. In order to achieve our objective, we begin with the one-and two-pointer work PDFs as defined in Eqs. (25) and (26) Comparing the above equations, it is easy to see that the PDFs are equal when either Q m = Q m or when D m, m = 0.
We restrict ourselves to a single cycle (N = 1) perfectly thermalizing engine, such that the elements of m ( ) = (m For case (iiia), since the evolution during the work strokes is adiabatic, i.e., if the system starts in an eigenstate of the initial Hamiltonian, it will end in the corre-sponding eigenstate of the final Hamiltonian, either Overall, with the three cases above, all possible combinations of m and m are excluded that would lead to a difference between the one-and the two-pointer work PDFs. Moreover, our proof does not depend on the working substance of the engine and hence holds for all thermalizing adiabatic quantum Otto engines diagnosed via repeated contacts.