Measurement of work and heat in the classical and quantum regimes
MMeasurement of work and heat in the classical and quantum regimes
P. Solinas,
1, 2
M. Amico, ∗ and N. Zangh`ı
1, 2 Dipartimento di Fisica, Universit`a di Genova, via Dodecaneso 33, I-16146, Genova, Italy INFN - Sezione di Genova, via Dodecaneso 33, I-16146, Genova, Italy The Graduate School and University Center, The City University of New York, New York, NY 10016, USA (Dated: February 3, 2021)Despite the increasing interest, the research field which studies the concepts of work and heat at quantumlevel has suffered from two main drawbacks: first, the difficulty to properly define and measure the work, heatand internal energy variation in a quantum system and, second, the lack of experiments. Here, we report for thefirst time the full experimental characterization of the dissipated heat, work and internal energy variation in atwo-level quantum system interacting with an engineered environment. We use the IBMQ quantum computerto implement the driven system’s dynamics in a dissipative environment. The experimental data allow us toconstruct quasi-probability distribution functions from which we recover the correct averages of work, heatand internal energy variation in the dissipative processes. Interestingly, by increasing the environment couplingstrength, we observe a reduction of the pure quantum features of the energy exchange processes that we interpretas the emergence of the classical limit. This makes the present approach a privileged tool to study, understandand exploit quantum effects in energy exchanges.
The oddities of quantum mechanics such as particle entan-glement, superposition of states, interference between evolu-tion paths and so on, have proven to be a valuable asset toenvision quantum devices able to outperform the correspond-ing classical ones. From metrology [1], to the detection ofgravitational waves [2] passing through quantum computationand communication [3], the use of quantum effects has givena decisive impulse towards new discoveries. With the adventof quantum technologies this trend will be reinforced allowingfor mass production of quantum-based devices. In this direc-tion, studying the energy exchange processes of a quantumsystem with an external drive and an environment could haveimportant implication for the future development of quantumtechnologies. This is the aim of a relatively new area of re-search referred to as quantum thermodynamics.Although numerous results have been obtained, there is stillno a clear consensus about how to determine what is the workand the heat in a quantum system [4–9]. In a closed quantumsystem, the work done on the system is equal to the variation of the internal energy of the system. However, the need forinformation at different times makes it impossible to envisiona proper measurement protocol for these quantities [10].This limitation has delayed experimental verification, es-pecially in case of a quantum system undergoing dissipativedynamics. Indeed, apart from a single experiment with closedquantum systems [11] and interesting proposals for open ones[12–18], a full and convincing measurement of the the dissi-pated heat is still missing. In this article, we fill this gap byshowing that it is possible to obtain the correct and expectedaverage values of work, heat and internal energy variation,and important information about the underlying quantum pro-cesses.To avoid any confusion and interpretative pitfalls, wechoose a more practical approach by focusing on simple andprecise questions. Given a quantum system controlled by anexternal time-dependent field, how much energy does the sys-tem absorb during the evolution? How much work does the external field do on the quantum system? And what is theheat dissipated by the quantum system?To answer to these questions, we implement the detectionscheme proposed in Refs. [8, 9, 19] on a IBMQ device [20]. Inparticular, we study a two-level quantum system, i.e., a qubit,subject to an external driving field and interacting with an en-gineered environment. The quantum detector and the engi-neered environment are represented by three additional qubits.The advantage of using an engineered environment is that wecan tune the system-environment coupling strength and ex-plore different dissipative regimes.Our physical observable is the phase of the detector qubitthat is measured with standard techniques [20, 21] from whichwe recover the information about the average work, heat andinternal energy variation while preserving the full quantumfeatures of the evolution [8, 9, 19, 21].Furthermore, from the measured detector phase we are ableto construct a quasi-characteristic generating function and aquasi-probability density function (QPDF) for these physicalobservables. The QPDF reproduces the statistics of the two-measurement protocol (TMP) [4] when the system is initiallyin an eigenstate of the Hamiltonian and keeps much more in-formation about the evolution of the system. In a direct anal-ogy with the Wigner function [22], the negative regions of thederived QPDF are associated to the violation of the Leggett-Garg inequalities and are the signature of a pure quantum phe-nomenon [19, 23–26]. The disappearance of these regions inpresence of strong dissipation can be seen as a proof of theemergence of the classical limit in energy exchange processesinduced by the presence of an environment.We start considering a two level system (denoted by S )that evolves under unitary evolution U S = U z U x with U z = exp ( − i β σ z ) and U x = exp ( − i ασ x ) where σ i ( i = x , y , z ) arethe usual Pauli operators. The system is initially in the state | ψ (cid:105) = cos θ | (cid:105) S + sin θ e i φ | (cid:105) S where | (cid:105) S and | (cid:105) S are theeigenstates of σ z .This evolution is generated by the time-dependent Hamilto- a r X i v : . [ qu a n t - ph ] F e b SDE E UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 SDE E UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 SDE E UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 FIG. 1. Three schemes for the detection of the variation of internalenergy ∆ U , work W and dissipated heat Q . The letter S, D, E andE denote the system, the detector, the environment qubit 1 and 2, re-spectively. The unitary dynamics are represented with red and light-blue filled squares corresponding to U x and U z , respectively. Thedissipative gates are represented with circles connecting the S qubitwith the E and E qubits. They corresponds to R x (dark red) and R z (dark blue). The system-detector coupling operations are representedwith rectangles spanning the S-D qubit space. They are U − χ , x and U χ , x represented with dashed red lines from up-left to down-right andup-right to down-left, respectively. The operators U − χ , z and U χ , z arerepresented with dashed blue lines from up-left to down-right andup-right to down-left, respectively. To simplify the presentation thegates needed to initialize the system and the detector are not shown(see [21] for details). nian H S = εσ x / ≤ t < t and H S = εσ z / t < t ≤ T with appropriate t and T . The fact that the Hamiltonianchanges in time assures that the external field does work onthe system [8].The detector is represented by an additional two-level sys-tem (denoted by D ). Its Hamiltonian is H D = ω Σ z / Σ i with i = x , y , z are the Pauli operators actingon the detector) and it is time-independent. The detector isinitialized in an equal superposition of eigenstates of H D , i.e., ( | (cid:105) D + | (cid:105) D ) / √ H SD = f ( χ , t ) H S ( t ) ⊗ H D al-lows us to store information about the system energy intothe accumulated phase of the detector [9, 19]. The func-tion f ( χ , t ) in H SD determines the time at which the system-detector coupling is active and its coupling strength χ . Ifthe system-detector coupling occurs on time scales muchsmaller than all of the other time-scales, we can assume that f ( χ , t i ) = χ / ε δ ( t − t i ) which generates the transformation U χ , t i = exp { i χε H S ( t i ) ⊗ H D } .We model the environment as an amplitude-damping chan-nel [3, 27] in which the environment can only induce relax- ation from the excited to the ground state. We also assumethat the unitary evolution occurs on shorter time scales withrespect to the relaxation time scale. Hence, we can imaginethe system evolution as a sequence of unitary dynamics fol-lowed by relaxation.To access to the whole spectrum of dissipative evolution,from unitary to strongly dissipative, we use an engineered en-vironment that can be simulated with two additional two-levelsystems [21]. On a NISQ device this translates into two ad-ditional qubits, allowing us to describe the relaxation of thesystem in the σ z and σ x channels with two operators denotedwith R z and R x , respectively (see figure 1).The building block of the measurement schemes is the op-erator sequence U χ , t i U U − χ , t j where U is a unitary operatoracting on the system. As discussed in Ref. [9], this allowsus to have information about the system energy changes in thetime interval t i − t j . Here, it is of fundamental importance tonote that work and heat are associated to different kinds of dy-namical evolution [8]. The work is associated to the changesin the Hamiltonian while the heat is associated to the changeof the system state when the Hamiltonian is constant.Since the system Hamiltonian changes only at time t , tomeasure the work we couple the system and the detectorshortly before and after t (Fig. 1).In the present experiment, the (dissipative) dynamics of thesystem is given by the operator U diss = R z U z R x U x . The evo-lution generated by the latter operator can be divided in twoparts, R x U x and R z U z , where the system Hamiltonian can beconsidered constant. Before and after the time t , when theHamiltonian is constant, the change in the system energy canbe associated to the action of the environment and, thus, inter-preted as dissipated heat. Analogously, we follow the schemein Fig. 1 to store information about the dissipated heat in the σ x and σ z system basis, respectively. Note that with this cou-pling scheme, we obtain information about the heat, i.e., thedissipated energy, supplied to the environment [21]. If we areinterested in the variation of the internal energy of the system,we need to couple the system and the detector at the beginningand at the end of the evolution only (see Fig. 1).The physical observable measured by the experiments is thephase accumulated in the detector using the different schemes.For a given system-detector coupling strength χ , it reads [9,19, 28] G χ , F = D (cid:104) | ρ D ( T ) | (cid:105) DD (cid:104) | ρ D | (cid:105) D = Tr S , E (cid:104) U χ , F ρ S U † − χ , F (cid:105) . (1)where F = ∆ U , W , Q , ρ ( T ) , ρ S ( T ) = Tr D , E [ ρ ( T )] and ρ D ( T ) = Tr S , E [ ρ ( T )] are the total, the final system and detec-tor density operators, respectively; U ± χ , F represents the full(system, detector and environments) evolution [9, 21].By changing the system-detector coupling strength χ , weobtain the quasi-characteristic function (QCGF) that is relatedto quasi-moments of the observable F [9, 19, 24, 28]. TheFourier Transform of the QCGF allows us to obtain a quasi-probability distribution function: P ( F ) = (cid:82) d χ G χ , F e i χ F . By - - - - - - - - - p=0 p=1p=0.5 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3) P U P W P Q E /✏ (4) h U ih W ih Q i p (5)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3) P U P W P Q E /✏ (4) h U ih W ih Q i p (5)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3) P U P W P Q E /✏ (4) h U ih W ih Q i p (5)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3) P U P W P Q E /✏ (4) h U ih W ih Q i p (5)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3) P U P W P Q E /✏ (4) h U ih W ih Q i p (5)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3) P U P W P Q E /✏ (4) h U ih W ih Q i p (5)1 FIG. 2. The quasi-probability distribution functions for the variation of internal energy, work and dissipated heat (y-axis) with respect to theenergy ( x -axis), normalized to the energy gap of the system E / ε . The parameter p represents the strength of the dissipation from p = p = θ = . φ = . α = β = . using the prefix ”quasi”, we stress that the above quantitiesare not obtained by direct physical measurements but derivedfrom the analysis of the detector phase measurements. In-deed, P ( F ) presents negative probability regions [9, 29] thatare related to the violation of the Leggett-Garg inequalities[9, 23–26, 29] and can therefore be seen as a signature of thequantum properties of the work.The main results of the experiments on the IBMQ deviceare shown in Fig. 2 (see [21] for the details about the imple-mentation) for different the dissipative parameters p . Whenno relaxation is present, i.e., p =
0, no heat is dissipated andthe heat distributions is peaked around Q =
0. The internalenergy and work QPDFs are identical and show the classi-cal energy peaks at E / ε = ± ,
0. However, there are alsoquantum energy peaks at E / ε = ± / p = . W and ∆ U have some notable differences but the features describedabove persist. However, in this case some heat is being dissi-pated by the system as pointed out by the presence of a peakat energy Q = E / ε =
1. Notice that, while the internal energyis bounded between ±
1, the dissipated heat is not. Because ofthe chosen dynamics, the system can dissipate energy at twotimes (formally when we apply R x and R y ) corresponding to amaximum energy exchange of two energy quanta.The strong dissipation case p = W and ∆ U distribution are now both positively defined and noquantum energy exchange at E / ε = ± / Q distribution is alwaysclassical, i.e., positively defined and with no half-quantumenergy exchange. This is a feature we expect from a largeMarkovian environment that is always in equilibrium and, inthis sense, a classical environment.These results and trends are confirmed by the behavior ofthe average values shown in Fig. 3 as a function of the dis- - - UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3) P U P W P Q E /✏ (4) h U ih W ih Q i p (5)1 h U i /✏ h W i /✏ h Q i /✏ h U i /" h W i /" h Q i /"p (6)+ ¯ ¯ Im[ G , F ] (7)2 h U i /✏ h W i /✏ h Q i /✏ h U i /" h W i /" h Q i /"p (6)+ ¯ ¯ Im[ G , F ] (7)2 h U i /✏ h W i /✏ h Q i /✏ h U i /" h W i /" h Q i /"p (6)+ ¯ ¯ Im[ G , F ] (7)2 FIG. 3. The average values of the variation of internal energy, workand dissipated heat normalized to ε as a function of the dissipativeparameter p . The light blue dots represent the experimental data.The dark blue curve is the theoretical prediction of the discussed ap-proach. The red curve is the prediction of the TMP obtained by sim-ulations with the IBMQ computers. The error in these simulations iswithin the curve thickness. As can be seen, for strong dissipation thethree curves converge demonstrating that the quantum features aredestroyed and the classical limit is reached. sipative parameter p . The experimental points (in light blue)are in good agreement with the expected theoretical prediction(blue curve). For no or weak dissipation ( p = p = . p = (cid:104) ∆ U (cid:105) + (cid:104) Q (cid:105) − (cid:104) W (cid:105) =
0. In particular, we find from the experimentaldata that the energy conservation law is satisfied within a theexperimental error [21].The present approach has several advantages with respectone presented in Refs [4, 30]. First, it allows us to obtainthe dissipated heat by acting on the system degrees of free-dom. This is an important difference with respect to the theo-retical proposal to measure the variation of the energy of theenvironment, which is practically unrealizable since it wouldrequire an insurmountable number of measurements. Otherviable schemes have been designed in order to measure thedissipated energy quanta [12–16]. However, these can be usedonly in specific physical systems while our scheme is system-independent and, thus, can be used with any quantum system.More importantly, the information about the classical TMPdistributions is included in the more general QPDFs that,at the same time, contains much more valuable information about the quantumness of the process. In Wigner’s spirit [22]these can be see as quantum correction to an underlying clas-sical process. The Wigner function has become an invalu-able tool to understand quantum phenomena [31, 32]. Analo-gously, the presented approach could shed light on the energyexchanges allowing us to identify their quantum features. Thisis the first step toward understanding if and when quantum en-ergy exchanges can be more efficient that classical ones andexploit their quantum advantage.
ACKNOWLEDGEMENTS
The authors would like to thank E. De Vito, S. Gasparinettiand A. Toigo for fruitful discussions. PS and NZ acknowledgefinancial support from INFN. The views expressed are thoseof the author and do not reflect the official policy or positionof Q-CTRL. ∗ Mirko Amico is currently employed at Q-CTRL inc.[1] V. Giovannetti, S. Lloyd, and L. Maccone, Nature Photonics ,222 (2011).[2] J. Aasi, J. Abadie, B. P. Abbott, R. Abbott, T. D. Abbott, M. R.Abernathy, C. Adams, T. Adams, P. Addesso, R. X. Adhikari,et al., Nature Photonics , 613 (2013).[3] M. A. Nielsen and I. L. Chuang, Quantum computation andquantum information (Cambridge University Press, Cambridge,UK, 2010).[4] M. Campisi, P. H¨anggi, and P. Talkner, Rev. Mod. Phys. , 771(2011).[5] P. Talkner, E. Lutz, and P. H¨anggi, Phys. Rev. E , 050102(2007).[6] R. Dorner, S. R. Clark, L. Heaney, R. Fazio, J. Goold, and V. Ve-dral, Phys. Rev. Lett. , 230601 (2013).[7] L. Mazzola, G. De Chiara, and M. Paternostro, Phys. Rev. Lett. , 230602 (2013).[8] P. Solinas, D. V. Averin, and J. P. Pekola, Phys. Rev. B ,060508 (2013).[9] P. Solinas and S. Gasparinetti, Phys. Rev. E , 042150 (2015).[10] M. Perarnau-Llobet, E. B¨aumer, K. V. Hovhannisyan, M. Hu-ber, and A. Acin, Phys. Rev. Lett. , 070601 (2017).[11] T. B. Batalh˜ao, A. M. Souza, L. Mazzola, R. Auccaise, R. S.Sarthour, I. S. Oliveira, J. Goold, G. De Chiara, M. Paternostro,and R. M. Serra, Phys. Rev. Lett. , 140601 (2014).[12] J. Pekola, P. Solinas, A. Shnirman, and D. Averin, New Journalof Physics , 115006 (2013).[13] S. Gasparinetti, K. L. Viisanen, O.-P. Saira, T. Faivre, M. Arzeo,M. Meschke, and J. P. Pekola, Phys. Rev. Applied , 014007(2015).[14] K. L. Viisanen, S. Suomela, S. Gasparinetti, O.-P. Saira,J. Ankerhold, and J. P. Pekola, New Journal of Physics ,055014 (2015).[15] O.-P. Saira, M. Zgirski, K. L. Viisanen, D. S. Golubev, and J. P.Pekola, Phys. Rev. Applied , 024005 (2016).[16] V. Vesterinen, O.-P. Saira, I. R¨ais¨anen, M. M¨ott¨onen,L. Gr¨onberg, J. Pekola, and J. Hassel, Superconductor Scienceand Technology , 085001 (2017). [17] B. Karimi, D. Nikoli´c, T. Tuukkanen, J. T. Peltonen, W. Belzig,and J. P. Pekola, Phys. Rev. Applied , 054001 (2020).[18] B. Karimi, F. Brange, P. Samuelsson, and J. P. Pekola, NatureCommunications , 367 (2020).[19] P. Solinas and S. Gasparinetti, Phys. Rev. A , 052103 (2016).[20] Ibm quantum experience https://quantum-computing.ibm.com/ .[21] P. Solinas, M. Amico, and N. Zangh`ı,
Supplementary informa-tion (2020).[22] E. Wigner, Phys. Rev. , 749 (1932).[23] A. J. Leggett and A. Garg, Phys. Rev. Lett. , 857 (1985).[24] A. A. Clerk, Phys. Rev. A , 043824 (2011).[25] M. Lostaglio, Phys. Rev. Lett. , 040602 (2018).[26] A. Levy and M. Lostaglio, PRX Quantum , 010309 (2020). [27] G. Garc´ıa-P´erez, M. A. C. Rossi, and S. Maniscalco, npj Quan-tum Information , 1 (2020).[28] P. Solinas, H. J. D. Miller, and J. Anders, Phys. Rev. A ,052115 (2017).[29] P. P. Hofer and A. A. Clerk, Phys. Rev. Lett. , 013603(2016).[30] S. Gasparinetti, P. Solinas, A. Braggio, and M. Sassetti, New J.Phys. , 115001 (2014).[31] P. Bertet, A. Auffeves, P. Maioli, S. Osnaghi, T. Meunier,M. Brune, J. M. Raimond, and S. Haroche, Phys. Rev. Lett. ,200402 (2002).[32] S. Del´eglise, I. Dotsenko, C. Sayrin, J. Bernu, M. Brune, J.-M.Raimond, and S. Haroche, Nature , 510 (2008). Supplemental Information
OPERATOR DECOMPOSITION IN IBMQ FUNDAMENTAL GATES
To implement the protocol described in the main text on the IBM quantum computer we need the Hadamard gate, the phasegate u ( θ ) and the controlled-NOT gate [S1]. The first two read in the {| (cid:105) , | (cid:105)} basis H = √ (cid:18) − (cid:19) (S1)and u ( θ ) = (cid:18) e i θ (cid:19) . (S2)The initial state of the system | (cid:105) S can be transformed in | ψ (cid:105) = cos θ | (cid:105) S + sin θ e i φ | (cid:105) S by the application of the operators U in = u (cid:16) ϕ + π (cid:17) H u ( θ ) H . (S3)The dynamical operator acting on the system is U S ( α , β ) = exp {− i β σ z } exp {− i ασ x } which can be obtained (a part from airrelevant exp { i ( β + α ) } factor) from the gate sequence U S ( α , β ) = u ( β ) H u ( α ) H . In the experiments and simulations, wehave set θ = . φ = . α = β = . U χ = e i χ e − i χ e − i χ
00 0 0 e i χ (S4)in the {| (cid:105) , | (cid:105) , | (cid:105) , | (cid:105)} basis.We denote with CX i , j the controlled-NOT gates that has the i -th qubit as the control and j has the target qubit [S1]. If thecoupled qubits are the system ( i =
1) and the detector one ( j = , gate in the same basis readsCX , = . (S5)When H S = εσ z , U χ = exp { i χσ z ⊗ σ z } and it can be implemented as U χ , z = CX , [ ⊗ u ( − χ )] CX , (S6)(a part a exp { i χ } factor).We also need another system-detector coupling operator when H S = εσ x ; that is U χ = exp { i χσ x ⊗ σ z } . The logical gatesequence reads U χ , x = H CX , [ ⊗ u ( − χ )] CX , H . (S7)Finally, we need to measure the phase accumulated by the detector due to the interaction with the system. The real andimaginary part of the off-diagonal matrix elements are measured in separate experiments. For the real part, this is done byapplying the Hadamard operator followed by a measurement. For the imaginary part, applying the operator exp {− i πσ x / } followed by a measurement. The second operator can be implemented by the use of the IBMQ u ( π / , − π / ) gate that, in the {| (cid:105) , | (cid:105)} basis, is [S1] u ( a , − a ) = √ (cid:18) − e − ia e ia (cid:19) . (S8) ENGINEERED ENVIRONMENT
We make the working hypothesis that the main source of decoherence is the relaxation of the system qubit from the excited tothe ground state. To describe this type of relaxation we use an engineered amplitude-damping channel [S2, S3]. This allows usto study the work and heat distribution at different level of environmental coupling.The system-environment qubit register is composed by a system qubit q s and an environment qubit q env so that the full quantumstate describing the register is | q s , q env (cid:105) . The transformation we implement is | (cid:105) → | (cid:105)| (cid:105) → (cid:112) − p | (cid:105) − √ p | (cid:105)| (cid:105) → (cid:112) − p | (cid:105) + √ p | (cid:105)| (cid:105) → | (cid:105) (S9)and it describes the the exchange of an energy quantum between the system and the environment. In the cold environment limit,i.e., k B T (cid:28) ε where T is the bath temperature, the system cannot be excited by the environment but can only relax to the groundstate. In this case, assuming that the environment is in the ground state, the above transformation is | (cid:105) → | (cid:105)| (cid:105) → (cid:112) − p | (cid:105) + √ p | (cid:105) . (S10)If the system is initially in the state | ψ s (cid:105) = α | (cid:105) + β e i ϕ | (cid:105) , the total initial state evolves into | ψ in (cid:105) = α | (cid:105) + β e i ϕ | (cid:105) → | ψ fin (cid:105) = α | (cid:105) + β e i ϕ ( (cid:112) − p | (cid:105) + √ p | (cid:105) ) (S11)which corresponds in terms of the system density matrix reduced dynamics (in the {| (cid:105) , | (cid:105)} basis) to ρ S , in = (cid:18) α αβ e − i ϕ αβ e i ϕ β (cid:19) → ρ S , fin = (cid:18) α + p β αβ e − i ϕ √ − p αβ e i ϕ √ − p ( − p ) β (cid:19) . (S12)As we can see, for p = p = | (cid:105) state.The same results can be obtained with the the operator sum-representation [S2, S3]. The operator elements in the system {| (cid:105) , | (cid:105)} basis are M = (cid:18) √ − p (cid:19) M = (cid:18) √ p (cid:19) . (S13)The final system density operator is ρ S , fin = M ρ S , in M + M ρ S , in M . By direct calculation we obtain the result in Eq. (S12). Implementation in Qiskit
For the implementation in Qiskit, we use the C u , i , j ( θ ) gate that induces a controlled rotation of an angle θ on the j -th qubit ifthe i -th qubit is active [S1]. The rotation on the target qubit generates | (cid:105) → cos θ | (cid:105) + sin θ | (cid:105) and | (cid:105) → cos θ | (cid:105) − sin θ | (cid:105) .The C u , i , j ( θ ) acting on a generic | i , j (cid:105) state produces | (cid:105) → | (cid:105)| (cid:105) → | (cid:105)| (cid:105) → cos θ | (cid:105) + sin θ | (cid:105)| (cid:105) → − sin θ | (cid:105) + cos θ | (cid:105) . (S14)In our case, the working space is composed by the system qubit q s and by one environmental qubit q env so that the full statereads | q s , q env (cid:105) . The logical gate sequence needed isRelaxation q s , q env = CX q env , q s C u , q s , q env ( θ ) CX q env , q s (S15)(Notice the inversion in the controlled operations, i.e., for the CX q env , q s gates the control is on the environment qubit and the targetis the system qubit while for the C u operation the control is on the system qubit while the rotation is applied to the environmentqubit). The total transformation on a | q s , q env (cid:105) state is | (cid:105) → | (cid:105)| (cid:105) → − sin θ | (cid:105) + cos θ | (cid:105)| (cid:105) → cos θ | (cid:105) + sin θ | (cid:105)| (cid:105) → | (cid:105) . (S16)This is reduced to the amplitude-damping channel [Eq. (S9)] when θ = [ √ − p ] and to Eq. (S10) for a cold environment.Therefore, the gate sequence in Eq. (S15) describes the relaxation of the system qubit in the in the {| (cid:105) , | (cid:105)} basis , i.e., when H S = εσ z .The relaxation in the {| + (cid:105) , |−(cid:105)} basis has the same structure with a Hadamard rotation on the system qubit, i.e.,RelaxationX q s , q env = H q s CX q env , q s C u , q s , q env ( θ ) CX q env , q s H q s . (S17) FULL DYNAMICS IMPLEMENTATION
The operation sequences needed to build the quasi-probability distribution functions of the work, heat and internal energyvariation are discussed in details in Refs. [S4, S5] and depicted in Figure 1 of the main text.However, for the sake of clarity, in Fig. S1 we explicitly show the logical gate sequence used to simulate the driven anddissipative dynamics on the IBMQ device. Here, the logical gate sequence in Eq. (S3) used to initialize the system and the onesused to measure the real and imaginary part of the detector phase are not shown.The U x and U z (shown as red and light-blue filled squares in the main text) are the external drive dynamics discussed above.The rectangular boxes represent the system-detector coupling operation U ± χ , x and U ± χ , z . In the main text, they are associatedto boxes with i ) dashed red lines from up-left to down-right ( U − χ , x ), ii ) dashed red lines from up-right to down-left ( U χ , x ), iii ) dashed blue lines from up-left to down-right ( U − χ , z ) and iv ) dashed blue up-right to down-left U χ , z . Their implementations interms of logical gate sequenced are given in Eqs. (S7) and (S6), respectively.The dissipative evolution is denoted by R x and R z which are obtained with the logical gates sequence (S17) and (S15). In themain text they were shown as a dark red ( R x ) and a dark blue ( R z ) circle.Depending on the physical observable, we fix the value of the system-detector coupling parameter χ , implement the cor-responding logical gate sequence and perform a measure of the detector phase (not shown in the scheme). By repeating thisscheme changing χ , we are able to construct the quasi-characteristic generating function and through a Fourier transform thequasi-probability density function. S D E E U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 S D E E UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 S D E E UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 U ,x U ,x U ,z U ,z U x U z R x R z (1)1 UWQ (1) U ,x U ,x U ,z U ,z U x U z R x R z (2)SDE E (3)1 FIG. S1. Logical gate sequence for the measurement of the average values of work W , heat Q and internal energy variation ∆ U . The gatesequenced needed to initialize the system and to measure the detector phase are not shown. DATA ANALYSIS AND PRESENTATION
The experimental data are obtained from the IBM quantum computer ibmq vigo device [S1]. For the physical observable ∆ U , W and Q , we have implemented the corresponding gate sequence presented in Fig. 1 of the main text.We have taken the system-detector coupling parameter 0 ≤ χ ≤ χ max =
100 with a discretization of ∆ χ = .
1. For each valueof χ k = k ∆ χ (with 0 ≤ k ≤ G χ , F as in Eq. (1) of the main text and the implementation discussed above. The noise level in the experiments on the IBMquantum computers can change depending on different protocol implementations or can increase in time (the quantum computersneed periodic recalibrations). For these reasons, whenever high noise level was observed, we have rerun the experiments (upto 24000 experiments for every χ k ) and taken the average G χ , F . We have used the function symmetries, Re [ G χ ] = Re [ G − χ ] andIm [ G − χ ] = − Im [ G χ ] , to extend the function in the interval − χ max ≤ χ ≤ χ max . From G χ , F we have calculated the correspondingquasi-probability distribution functions as P ( F ) = (cid:82) χ max − χ max d χ G χ , F e i χ F .The experiments on IBM quantum computer are flanked by numerical simulations. These are obtained simulating the systemdynamics with the same parameter and precision, i.e., the same χ max and ∆ χ .Finally, the TMP distribution results are obtained by running a numerical simulation of the TMP protocol on the ibmq vigodevice [S1]. Since in this case we perform a direct measurement of the probability distribution, we do not need a detector qubitand we used a three-qubit register: one qubit for the system and two for the environment.To simulate the initial measurement, we prepare the system in one of the two two eigenstates of H S ( ) , i.e., ( | (cid:105) ± | (cid:105) ) / √ U x , z and relaxation operators R x , z as in the QPDF scheme. At the end of thesimulation, all three qubits are measured. The knowledge of the final measured state | q s , q env , q env (cid:105) (in the {| (cid:105) , | (cid:105)} basis)allows us to determine both the variation of internal energy and the dissipated heat. In this way, we can calculate the variationof internal energy ∆ U as the difference between the energy of the initial state and the final measured state. We can also calculatethe dissipated heat. Since at the beginning we have q env = q env =
0, the normalized dissipated heat is Q / ε = q env + q env . Thework can now be calculated as a derived quantity , equal to W = ∆ U + Q . The results obtained are weighted with the probabilityto start from the corresponding eigenstate when the measurement is performed on | ψ (cid:105) = cos θ | (cid:105) S + sin θ e i φ | (cid:105) S .To compare the experimental results for the QPDFs with the ones obtained from the simulation of the TMP, we focus only onthe classical energy peaks (integer quanta exchange). These correspond to the case in which the system has no initial coherencesand the energy processes becomes classical and, thus, reproduced by the TMP [S4, S5].The relevant quantities to be compared are the probability to obtain a process with energy E . For the present protocol, it h U i /✏ h W i /✏ h Q i /✏p (6)+ ¯ (7)2 h U i /✏ h W i /✏ h Q i /✏p (6)+ ¯ (7)2 h U i /✏ h W i /✏ h Q i /✏p (6)+ ¯ (7)2 h U i /✏ h W i /✏ h Q i /✏p (6)+ ¯ (7)2 h U i /✏ h W i /✏ h Q i /✏p (6)+ ¯ ¯ (7)2 h U i /✏ h W i /✏ h Q i /✏p (6)+ ¯ ¯ (7)2 h U i /✏ h W i /✏ h Q i /✏p (6)+ ¯ ¯ Im[ G , F ] (7)2 FIG. S2. Geometric calculation of the error. The imaginary part of G χ , F is linearly approximated around χ = χ = . δ determined by the statistical uncertainty in the determination of the ground and excited state populations. The error determines a change inthe line slope and, consequently, a change in the average values of the physical observables. corresponds to the integral (cid:82) E + δ E − δ d FP ( F ) , that is the probability to exchange an energy between E − δ and E + δ . On thecontrary, for the TMP we obtain a probability mass functions so the probability to obtain E is obtained evaluating the distributionat the point. These statistical distributions are intrinsically different, i.e., the experimental data are continuous (probabilitydistribution) functions while the TMP are discrete functions, so to present them together, we need to renormalize the experimentaldistributions.The oscillations in the experimental functions are due to the limitation of the experimental data. As χ max increases, thedistributions are more and more peaked around the discrete energy exchange ( E = { , ± , ± } ). Accordingly, we assume thatall the information about the physical processes is in these discrete points. Then, for every probability distribution, we adoptedthe following procedure:1. We calculate the values of the distribution at the peaks and sum them to obtain a norm: N = ∑ E P ( E ) .2. We then rescale the original (quasi-)probability distribution function with respect to this norm: P (cid:48) ( F ) = P ( F ) / N .The final effect of this procedure is to preserve the relative strength of the peaks give a direct visual comparison between theTMP distribution peaks and the classical peaks of the quasi-probability distributions. ERROR CALCULATION
As discussed above, in the experiments we obtain the real and imaginary part of the phase accumulated by the detector duringthe evolution by measuring the population of the ground p and the excited state p of the detector qubit. These are related to G χ , F in Eq. ( ) of the main text by the expressions Re [ G χ , F ] = p − p = p − H gate, and Im [ G χ , F ] = p − p = p − {− i πσ x / } gate(we recall that the initial state of the detector has D (cid:104) | ρ D | (cid:105) D = / p and p are affected by errors. Since it is not possible to characterize the error induced by thedynamics, i.e., the sequence of logical gates, we consider only the statistical error due to the final measurement.The interesting physical quantities are the average values of ∆ U , W or Q therefore we only discuss the experimental erroron these. The propagation of the experimental error to the QPDF and then the average values is not straightforward because itrequires the definition of a ”distance” between the two QPDFs (the measured one and the one with the included errors). Here,we use a more direct and simple way to calculate the error on the average quantities of interest.The G χ , F is built so that the first derivative (calculated in χ =
0) reproduces the average values of the physical observables F [S4–S6], i.e., (cid:104) F (cid:105) = − i d G χ , F d χ (cid:12)(cid:12)(cid:12) χ = . (S18)Note that, given the symmetries of G χ , F , only the imaginary part contributes to the first moment.As shown in [S4–S6], the accumulated phase on the detector is exp { i χ ( ε i − ε j ) } where ε i − ε j is the normalized energyvariation during the process. In our case, ε i − ε j = , ± ∆ U and W and ε i − ε j = , ± , ± Q . Therefore, the imaginarypart of G χ , F behaves as sin [ χ ( ε i − ε j )] and can be linearly approximated for χ (cid:28) / ( ε i − ε j ) ≤ /
2. In the linear approximation(see in Figure S2), (cid:104) F (cid:105) is the slope of the line G + (cid:104) F (cid:105) χ . In this geometric interpretation, the error on the populations p and p determine an error δ in Im [ G χ , F ] resulting, eventually, in a change of slope and the average value of the physical observable(see Figure S2).More quantitatively, the change in the slope due to δ at point χ = ¯ χ can be determined by geometric means and reads ∆ (cid:104) F (cid:105) = (cid:12)(cid:12)(cid:12) δ ¯ χ (cid:12)(cid:12)(cid:12) . (S19)For every value χ , we performed N = G χ , F . Since the experiments areindependent, the error on the determination of the populations p i can be estimated to be 1 / √ N (a more precise description ofthe error does not change significantly the results).The smallest positive value of χ at which the G χ , F was estimated is ¯ χ = .
1. We have repeated the experiments 3 N times for G χ , ∆ U and p =
1, 2 N times for G χ , ∆ U and p = . G χ , Q and p = . ,
1. For all the other experiments we have performed N experiments. These data with δ = / √ N allow us to calculate the errors for the physical averages presented in the main text.The theoretical prediction of the averages have been calculated with numerical simulations. We have plotted the averagescalculated directly as the derivative of the QCGFs. The theoretical calculation through the QPDFs give the same results.Finally, the averages of ∆ U , W and Q for the TMP, are obtained by numerical simulation with the IBMQ computers. Thestatistical error in this case is 1 / √ N TMP where N TMP are the number of simulated experiments performed. In this case, the errorbar are so small that are not visible in the plot.
NOTATION OF THE HEAT SIGN
The coupling sequence U χ , t i U U − χ , t j where U ± χ , t = exp { i χ H S ( t ) ⊗ H D } , allows us to determine the variation of the internalenergy of the system, i.e., ∆ U = ε f − ε i [S4, S5]. This is the energy supplied to the system ; if the system absorbs energy ∆ U > ∆ U <
0. This can be seen directly by calculating the CGF and the average value as discussed in Ref. [S4].On the contrary the transformation U − χ , t i U diss U χ , t j measures the energy supplied by the system to the environment. In otherwords we are ”measuring” q = ε i − ε f . Therefore, if q = ε i − ε f > q = ε i − ε f < q is the energy supplied to the system as heat.As a consequence, the energy conservation law takes the form (cid:104) ∆ U (cid:105) + (cid:104) Q (cid:105) − (cid:104) W (cid:105) = ∗ Mirko Amico is currently employed at Q-CTRL inc.[S1]
Ibm quantum experience https://quantum-computing.ibm.com/ .[S2] M. A. Nielsen and I. L. Chuang,
Quantum computation and quantum information (Cambridge University Press, Cambridge, UK, 2010).[S3] G. Garc´ıa-P´erez, M. A. C. Rossi, and S. Maniscalco, npj Quantum Information , 1 (2020).[S4] P. Solinas and S. Gasparinetti, Phys. Rev. E , 042150 (2015).[S5] P. Solinas and S. Gasparinetti, Phys. Rev. A , 052103 (2016).[S6] P. Solinas, H. J. D. Miller, and J. Anders, Phys. Rev. A96