The distribution of work performed on a NIS junction

The distribution of work performed on a NIS junction Jaime E. Santos Email: [email protected] Center for Nanostructured Graphene (CNG), Department of Micro and Nanotechnology, Technical University of Denmark, DK-2800 Kongens Lyngby, Denmark Centro de F´ısica and Departamento de F´ısica, Universidade do Minho, P-4710-057 Braga, Portugal Pedro Ribeiro Email: [email protected] Russian Quantum Center, Novaya Street 100 A, Skolkovo, Moscow Area, 143025 Russia CeFEMA, Instituto Superior T´ecnico, Universidade de Lisboa, Avenida Rovisco Pais, 1049-001 Lisboa, Portugal Stefan Kirchner Email: [email protected] Center for Correlated Matter, Zhejiang University, Hangzhou, Zhejiang 310058, China Max Planck Institute for the Physics of Complex Systems, N¨othnitzer Str. 38, D-01187 Dresden, Germany Max Planck Institute for Chemical Physics of Solids, N¨othnitzer Str. 40, D-01187 Dresden, Germany Abstract. We propose an experimental setup to measure the work performed in a normal-metal/insulator/superconducting (NIS) junction, subjected to a voltage change and in contact with a thermal bath. We compute the performed work and argue that the associated heat release can be measured experimentally. Our results are based on an equivalence between the dynamics of the NIS junction and that of an assembly of two-level systems subjected to a circularly polarised field, for which we can determine the work-characteristic function exactly. The average work dissipated by the NIS junction, as well as its fluctuations, are determined. From the work characteristic function, we also compute the work probability-distribution and show that it does not have a Gaussian character. Our results allow for a direct experimental test of the Crooks-Tasaki fluctuation relation. PACS numbers: 05.70.Ln,74.45.+c,74.50.+r arXiv:1506.02905v2 [cond-mat.stat-mech] 31 Jan 2016 The distribution of work performed on a NIS junction 2 1. Introduction Fluctuation relations are central to our present understanding of statistical mechanics. Their long and distinguished history goes back to at least the studies by Callen and Welton [1], Green [2] and Kubo [3], which were inspired by the works of Einstein on the Brownian movement [4] and of Johnson [5] and Nyquist [6] on noise in electrical circuits. The derivation by Jarzynski of a rather general non-equilibrium work relation [7], linking the average of the exponential of the work being performed on a system with the equilibrium free energy difference between initial and final equilibrium states of the system, is of particular interest. Subsequently, Crooks [8, 9] obtained the Jarzynski equality from a relation between the probability of a given amount of work being performed on a system and the probability that the system performs the same amount of work on its surroundings if the time-reversed protocol is undertaken. Such heightened interest in non-equilibrium work relations has not only been fuelled by novel theoretical advances in out-of-equilibrium dynamics, but also by experimental progress in preparing and probing non-equilibrium evolution (see [10, 11, 12, 13, 14] for reviews on the subject). In the present paper, we propose a relatively simple yet realistic experiment based on a proposal by Crooks [15], which would allow for a direct test of the Crooks-Tasaki fluctuation relation. The experimental system we consider consists of a normal-metal/insulator/superconducting (NIS) junction between a superconductor and a normal metal, which is initially short-circuited and is subjected to a given voltage protocol (see below). We establish that for the proposed protocol, the full work distribution-function has non-Gaussian character with a non-standard decay exponent. We also compute the first two moments of such distribution, which can be determined experimentally by measuring the average heat released and its variance. In order to obtain a closed expression for the characteristic function of the work distribution, we derive an equivalence between the dynamics of an NIS junction and the one of an assembly of two-level systems subjected to a circularly polarized field. Using this mapping, the work distribution is determined. For realistic parameters of experimentally available NIS junctions as e.g. those of Ref. [16], our findings show that at cryogenic temperatures, even the tiny amounts of heat released will result in a measurable volume change for a probe connected to the junction. The structure of the paper is as follows: Sect. 2 places the recent developments briefly outlined above into their proper historical context and reviews the theoretical tools necessary to study work-fluctuation relations at the quantum level. In Sect. 3, we introduce the Hamiltonian of the NIS junction and describe in detail the work protocol applied to the system. In Sect. 4, we use results that are derived in Appendix A and Appendix B to compute both the average work dissipated as well as its fluctuations and the full work distribution. Appendix C contains an extension of the results presented in Appendix B and Appendix D contains a discussion of a method of measurement of the tunnelling amplitude that describes the NIS junction. Section 5 contains the main The distribution of work performed on a NIS junction 3 result of our paper, where we compute the numerical value of the average heat released on a NIS junction such as the ones discussed in Ref. [16]. We leave to Appendix E the discussion of the experimental techniques that can be used to measure such a quantity. In particular, we show that the heat that is released can be measured by determining the volume change of a probe that absorbs it. We provide estimates of the relative volume change, which should be measurable using capacitive methods [17, 18]. Finally, in Sect. 6, we present our conclusions. For increased readability, most technical details have been relegated to the appendices. 2. Perspective and previous works In this section, we give a brief historical overview of the evolution of the subject of nonequilibrium fluctuation theorems. Besides the early contributions already mentioned above, a statement of a general fluctuation theorem involving the free energy of a gas of hard-spheres, interacting at large distances via the Lennard-Jones potential, can be found in the work of Zwanzig [19]. A subsequent important development was the work of Bochkov and Kuzovlev [20, 21], which also concerned the relation of probabilities between a given work protocol and its time reversed version, but which adopted a different perspective to that of Jarzynski and Crooks, namely the former authors considered the generalised force that performs work on the system as being external to the dynamical description of the system, rather than internal. Another significant contribution was the work of Evans and Searles (see [22] and references therein), which considers the relation between the probability for a dynamical trajectory characterised by a certain value of a dissipation function and its time-reversed counterpart. Regarding developments concerning the validity and extension of the Jarzynski-Crooks relations at the classical level, see the experimental works [23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33] and the theoretical ones [34, 35, 36, 37, 38, 39], the above list not being exhaustive. A quantum generalisation of the Crooks equality was obtained by Tasaki [40] and also by Kurchan [41] (for the case of cyclic protocols). A version of the Jarzynski equality valid for quantum systems was also derived by Yukawa [42], albeit treating work at a quantum level as an observable, a view that was shown not to be correct [43]. Such a view lead in the past to a debate concerning the validity of such fluctuation relations at the quantum level [44, 45, 46]. The extension of the Jarzynski equality to isolated systems, using the definition of the work distribution function of Tasaki and Kurchan, is due to Mukamel [47]. A quantum generalisation of the work of Bochkov and Kuzovlev referred above, as well as a clarification of the relation of their work to that of Jarzynski and Crooks, can be found in a paper by Campisi et al. [48]. Extensions of non-equilibrium work relations to isolated systems (microcanonical fluctuation theorems) are given in [49, 50], whereas work relations valid for arbitrary open quantum systems were obtained in [51]. The last reference shows the validity of such theorems for systems that can arbitrarily exchange heat with the bath while a work protocol is being applied to them (see also below). Other The distribution of work performed on a NIS junction 4 Figure 1. (a) NIS junction. (b) Schematic representation of the voltage protocol considered theoretical developments concerning quantum systems can be found in the references [52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72], where again the list is not exhaustive. At the experimental level, several proposals to demonstrate the validity of the Crooks-Tasaki relation in the quantum domain have been presented, including the use of a series of projective measurements [73], the measurement of this relation in the optical spectra of systems subjected to weak quenches [74], or the use of qubit interferometry [75, 76, 77]. The former as well as the latter proposal have been successfully implemented, see [78, 79]. An example of the experimental confirmation of a generalised version of the Jarzynski equality in the quantum domain was discussed in [71]. Finally, a discussion of possible solid state experiments performed on quantum heat engines is given in [72]. The above mentioned proposal by Crooks [15] relies on an indirect method of observation namely, the measurement of the heat released by a quantum system after a cyclic work protocol has been performed on the system. A specific set-up to measure the heat released by a circuit including a resistor and a Cooper-pair box, based on the temperature increase of a local probe, was already presented in [80], where the Cooperpair box acts as a single two-level system. The present proposal has the additional advantage that the heat released scales with the contact area between the normal-metal and the superconducting films. The work probability distribution, Pβ(W, τ), is defined, for a class of work protocols where the system is decoupled from the thermal bath during the application of the protocol (see below), of duration τ, as [40] Pβ(W, τ) = 1 Z(β)X m,n e−βE0 m| hn, τ |b U(τ, 0)|m, 0i |2δ(W−Eτ n+E0 m),(1) where |m, 0iare the eigenstates of the system’s Hamiltonian at the beginning of the protocol, c H(0) (which is given for the NIS junction by (5) at t= 0, see below), and E0 mare the corresponding eigenvalues, and where |n, τ iare the eigenstates of the system’s Hamiltonian at the end of the protocol, c H(τ) (which is given by the same equation at t=τ), with Eτ nbeing the corresponding eigenvalues, with β= 1/(kBT) and Z(β) = Tr(e−βb H(0)). The operator b U(τ, 0) is the time-evolution operator of the system The distribution of work performed on a NIS junction 5 during the protocol, i.e. in the interval [0, τ], which can be written as a time-ordered exponential of the full Hamiltonian of the system, given for a NIS junction by eq. (5). The work characteristic-function of the NIS junction is the Fourier transform of eq. (1), and is given by Gβ(υ, τ) = Tr b U†(τ, 0) eib H(τ)υb U(τ, 0) e−ib H(0)υb ρ(β),(2) where b ρ(β) = e−βb H(0)/Z(β) is the density matrix of the junction before the application of the protocol. The unitary evolution of the system in the interval [0, τ] reflects the work protocol’s restriction mentioned above, i.e. the NIS junction is disconnected from the thermal bath at t= 0 and remains adiabatically insulated while the protocol is being applied, being reconnected to the thermal bath at t=τand undergoing equilibration through the exchange of heat with the bath afterwards. It is possible to lift such a restriction on the protocol by explicitly considering the interaction of the system with the thermal bath [51], but we will not consider such an extension here, see however the discussion below. The knowledge of the work characteristic-function allows one to compute all work moments. If one performs the analytic continuation υ=iβ in the above definition and uses the cyclic character of the trace operation, one obtains Gβ(iβ, τ) = Zτ(β)/Z(β) = e−β∆F,(3) where Zτ(β) = Tr(e−βb H(τ)) is the partition function for a junction in equilibrium, described by the Hamiltonian c H(τ), ∆Fbeing the difference in free-energies between the final equilibrium state (after thermalisation at a time t > τ) and the initial equilibrium state at t= 0. Since in our case the Hamiltonian goes through a cycle, ∆F= 0. Expressing Gβ(iβ, τ) through its Fourier transform, we obtain he−βW i=Z+∞ −∞ dW e−βW Pβ(W, τ) = e−β∆F= 1 ,(4) which is a statement of the Jarzynski equality, which the work probability distribution satisfies by definition. Considering the relation between Gβ(υ, τ) and Gβ(υ+iβ, τ), using again the cyclic character of the trace, it is possible to derive the identity Pβ(−W, −τ) = e−βW Pβ(W, τ) where by Pβ(W, −τ) we denote the work probability distribution for the time-reversed protocol. This is a statement of the Crooks-Tasaki relation. Note that in the case of cyclic protocols, one also has ∆U= 0, where Udenotes the internal energy of the system. Therefore, from the first law of thermodynamics, ∆U=hWiβ−∆Q, where ∆Qdenotes the average heat released by the system on which the work protocol is applied, we obtain ∆Q=hWiβ, which shows that for cyclic protocols the measurement of the released heat allows for the determination of hWiβ, as stated. The distribution of work performed on a NIS junction 6 3. NIS junction model The NIS junction is composed of a normal-metal film, a thin insulating layer and a superconducting film. We consider that such a junction can be described by an Hamiltonian consisting of three terms: a Fermi-gas describing the normal-metal film, a BCS superconductor describing the superconducting film and a tunnelling Hamiltonian coupling the two [81, 82, 83, 84, 85, 86, 87], which gives rise to induced superconductivity in the normal-metal film through the proximity effect (see [88] and references therein). We use the Hamiltonian c H(t) = c Hn(t) + c Hs+c X. (5) The operator c Hn(t) = Pk,σ (ξn k−e φ(t)) b c† nkσb cnkσis a free-fermion Hamiltonian describing the normal-metal, with φ(t) being the voltage difference across the junction, whose value changes in time. The operator c Hs=X k,σ ξs kb c† skσb cskσ+X k∆kb c† sk↑b c† s−k↓+ ∆kb cs−k↓b csk↑,(6) is the Hamiltonian describing the BCS superconductor, where the pairing function ∆k is different from zero on a band of width 2¯hωDaround the Fermi-level, where ωDis the superconducting material’s Debye frequency. Finally, the third term of Eq. (5) is given by c X=X k,σ tσ(k)b c† skσb cnkσ+b c† n−k−σb cs−k−σ,(7) and represents the tunnelling process across the junction. We consider the overall tunnelling matrix element to be small with respect to the gap function ∆k. Such matrix element is taken to be invariant under time-reversal, i.e. tσ(k) = t−σ(−k) [85]. The fermions’ momentum is directed along the plane of the films and we assume it to be conserved in a tunnelling process. We will measure the kinetic energy of electrons in the normal-metal film and in the superconducting film with respect to their common Fermienergy µ. The superor subscripts nand srefer to operators or quantities pertaining to the normal-metal or to the superconductor, respectively. The magnitude of the tunnelling amplitude will depend on the width of the aluminium oxide layer and its value can be adequately measured for a given device using the method discussed in Appendix D. In this method, the power dissipated by the junction is measured at a constant applied voltage equal in value to the energy gap of the superconductor (divided by e), which is much larger than the value used in our protocol. We consider the NIS to be initially (t < 0) in equilibrium with an heat bath at temperature Tand that the voltage φacross the junction is zero. Such voltage is changed from 0 to a finite value Vmax and returned to zero within a finite time-interval starting at t= 0 and ending at t=τ. During such an interval the junction is decoupled from the heat bath. For simplicity, we consider that the voltage protocol is symmetric with respect to the mid-point of the interval t=τ/2. However, the calculation presented in Appendix A is valid for more general protocols. Also, φ(t) varies smoothly within a The distribution of work performed on a NIS junction 7 time-interval τ0τand is returned to zero at the end of the protocol within the same time frame (see figure 1 (b)). Physically, for small values of the tunnelling amplitude, the term given by eq. (7) leads to the opening of a gap between the energy bands of the normal metal, of magnitude ∆n∆, leaving the energy bands of the superconductor largely unaffected. This is the proximity effect [88]. In our protocol, we place ourselves within a particular adiabatic limit, i.e. we will consider that the switching-time τ0¯h/∆, where ∆ is the value of the gap-function within the superconductor. Thus, the voltage difference across the junction is held approximately constant within the microscopic time-scale of the superconductor and the change of voltage does not induce transitions which involve the superconducting bands. However, we also take τ¯h/∆nτ0¯h/∆. This implies that the dynamics is fully diabatic with respect to the separation between the two energy-bands in the normal metal and thus that the voltage protocol induces transitions between these bands, leading to work being performed on the system, see Appendix B. Using the Hamiltonian for an NIS junction and the voltage protocol specified, the work probability distribution can be obtained using Eq. (2), as will be shown in the next section. 4. Work dissipated due to the application of the voltage protocol to the junction In this section, we establish that the dynamics of the NIS junction is equivalent to that of an assembly of quantum spins in a time-dependent magnetic field. This entails an equality between the work probability distributions of both systems, from which it follows that one can compute the average heat released by the NIS junction, as well as its higher order fluctuations. Starting from Eq. (5), the explicit dependence of the Hamiltonian on the timedependent bias voltage can be eliminated by a gauge transformation at the expense of acquiring a time-dependent hopping term between the normal and superconducting films. Thus, the Hamiltonian becomes c H1(t) = c Hn(0) + c Hs+c X(t) (8) with c X(t) = X k,σ tσ(k)eieΦ(t)/¯hb c† skσb cnkσ+e−ieΦ(t)/¯hb c† n−k−σb cs−k−σ,(9) where Φ(t) = Rt 0du φ(u). Anticipating that, in general, tσ(k)∆, one can, in the specific adiabatic limit discussed above, perform a generalised SchriefferWolff transformation that takes into account the explicit time-dependence of c X(t). Furthermore, such a transformation can only be performed if eVmax ∆, see Appendix A. After such a transformation, the metal and superconducting films effectively decouple and only the transformed metallic excitations are subjected to the time-dependent bias. The distribution of work performed on a NIS junction 8 The Hamiltonian, c Hn(t), which describes the normal-metal film, explicitly reads c Hn(t) = X k,σ ˜ ξn kσb c† nkσb cnkσ +X k∆n ke−2ieΦ(t)/¯hb c† nk↑b c† n−k↓+ ∆n ke2ieΦ(t)/¯hb cn−k↓b cnk↑,(10) where ˜ ξn kσ=ξn k+|tσ(k)|2(ξn k+ξs k) (ξn k)2−(Es k)2,(11) ∆n k=−∆k(|t↑(k)|2+|t↓(k)|2) 2 [ (ξn k)2−(Es k)2],(12) are, respectively, the renormalized kinetic energy of electrons and the induced superconducting parameter, in the normal metal, due to the proximity effect, Es k= q(ξs k)2+|∆k|2being the energy of the superconducting excitations. This result is shown in detail in Appendix A. We thus conclude that for spin-independent tunnelling matrix elements, the Hamiltonian describing the normal-metal is equivalent, in each ((k,↑),(−k,↓)) subspace, to the Hamiltonian of a two-level system under the action of a circularly polarised field. Furthermore, the Hamiltonian that describes the superconductor does not contribute to the work-characteristic function, see Appendix A. One can obtain the work characteristic function for a NIS junction simply by considering the product of individual characteristic-functions for two-level systems subjected to a circularly polarised field, see Appendix B. The logarithm of the overall work characteristic-function, which is the generating function of the connected workmoments, is given by Wβ(υ, τ)≡ln Gβ(υ, τ) = X k ln n(1 −pk) +pk· cosh (β+ 2iυ)q|∆n k|2+ (˜ ξn k)2 cosh βq|∆n k|2+ (˜ ξn k)2       ,(13) with pk=|∆n k|2(eVmax)2 (|∆n k|2+(˜ ξn k)2)(|∆n k|2+(˜ ξn k−eVmax)2)sin2τ ¯hq|∆n k|2+ (˜ ξn k−eVmax)2. The detailed calculation can be found in Appendix B. The average work dissipated per atom in the normal-metal film, ¯w, is obtained from the derivative of Eq. (13), evaluated at υ= 0, ¯w=2(eVmax)2 Nat X k tanh βq|∆n k|2+ (˜ ξn k)2 ×|∆n k|2sin2τ ¯hq|∆n k|2+ (˜ ξn k−eVmax)2 q|∆n k|2+ (˜ ξn k)2(|∆n k|2+ (˜ ξn k−eVmax)2) ,(14) where Nat is the total number of atoms in the normal-metal film. From Eq. (11), one has that for weak-coupling between the two films, ˜ ξn k≈ξn k. Moreover, if, as in the The distribution of work performed on a NIS junction 9 system experimentally studied in [16], the dispersion relation is the same in both the normal-metal and the superconducting film (since the normal-metal film is composed of the same material as that of the superconducting one, but weakly doped with another metal), we have from Eq. (12) that ∆n k=|t(k)|2 ∆∗ k. Assuming that both the order parameter ∆kand |t(k)|2only depend on kthrough their dependence on ξn kand are only non-zero in a vicinity of width ¯hωDaround the Fermi level of the material composing the films, where ωDis the Debye frequency of the said material, one can transform the above summation into an integral and write the above equation as ¯w= 2(eVmax)2Z¯hωD −¯hωD dξ ρ(ξ) tanh βq|∆n(ξ)|2+ξ2 ×|∆n(ξ)|2sin2τ ¯hq|∆n(ξ)|2+ (ξ−eVmax)2 q|∆n(ξ)|2+ξ2(|∆n(ξ)|2+ (ξ−eVmax)2),(15) where we have introduced ρ(ξ) = 1 Nat Pkδ(ξ−ξn k), the density of states per atom (and per spin-species) in the normal metal film. Differentiating Eq. (13) twice with respect to υ, at υ= 0, one obtains the meansquare deviation and thus the relative deviation √δw2/w of the work performed by the NIS junction during the protocol. If one transforms such result into an integral using the same assumptions we have used to obtain Eq. (15), this integral can also be computed numerically. In a similar fashion, higher moments can be calculated. The full characterisation of the work fluctuations requires the determination of the work probability distribution, through the computation of the inverse Fourier transform of Gβ(υ, τ). Writing the probability distribution in terms of the logarithm of the characteristic function, as defined in Eq. (13), one has Pβ(W, τ) = Z+∞ −∞ dυ 2πe−iWυ +Wβ(υ,τ).(16) Since Wis an extensive quantity, i.e. it is proportional to Nat, the integral can be computed using a saddle-point approximation. Minimising the argument of the exponential in Eq. (16), one obtains the following saddle-point equation w=Z¯hωD −¯hωD dξ ρ(ξ)∂xS(ξ, x),(17) where w=W/Nat is the work per atom of the sample and S(ξ, x) is given by S(ξ, x) = ln     (1 −p(ξ)) + p(ξ)·cosh h(β+ 2x)q|∆n(ξ)|2+ξ2i cosh hβq|∆n(ξ)|2+ξ2i    ,(18) where p(ξ) = |∆n(ξ)|2(eVmax)2 (|∆n(ξ)|2+ξ2)(|∆n(ξ)|2+(ξ−eVmax)2)sin2τ ¯hq|∆n(ξ)|2+ (ξ−eVmax)2. Note that the solutions, x(w), of the saddle-point equation (17) correspond to the analytic continuation of the exponent of the integral given in Eq. (16) to imaginary values of υ=−ix(w). Numerically solving Eq. (17) for a set of values of w(with the same assumptions as in the two calculations previously performed) and substituting its solution in Eq. (16), one obtains the result plotted in figure 2 for the logarithm The distribution of work performed on a NIS junction 16 ×ei ¯h(eΦ(u)+ξn k(t−u)) ·((ξn k+ξs k)·cos(ωs k(t−u)) −iξs k(ξn k+ξs k) + |∆k|2 Es k·sin(ωs k(t−u))!,(A.18) and sσ k(t) = −i tσ(k)σ∆k (ξn k)2−(Es k)2·eie ¯hΦ(t)−tσ(k)σ∆k ¯h[(ξn k)2−(Es k)2]Zt 0du eφ(u) (A.19) ×ei ¯h(eΦ(u)+ξn k(t−u)) · cos(ωs k(t−u)) −iξn k Es k·sin(ωs k(t−u))!. We stated above that we are considering protocols in which the applied voltage difference φ(t) changes from 0 to Vmax in a time interval τ0i.e. φ(t) = Vmaxψ(t/τ0), where ψ(x) changes from 0 to 1 in an interval of length 1. The integrals which enter the second term of both (A.18) and (A.19) are all, after the change of variables u→u/τ0, of the form eVmaxτ0 ¯hZt/τ0 0du ψ(u)eieVmaxτ0 ¯hΨ(u)e−iτ0 ¯h(ξn k±Es k)u,(A.20) where Ψ(u) = Ru 0dv ψ(v). The term ψ(u)eieVmaxτ0 ¯hΨ(u)is a bounded function, whereas the exponential e−iτ0 ¯h(ξn k±Es k)uis at least equal to e−iτ0 ¯h|∆k|u, which in the adiabatic regime τ0¯h/∆ oscillates rapidly. Thus, provided that tτ0and eVmaxτ0/¯h≤1, then (A.20) is zero by the Riemann-Lebesgue lemma. From these conditions, we also obtain eVmax ∆, i.e. the tunnelling voltage applied should be well within the gap. Therefore, we conclude that in the adiabatic limit and for tτ0, (A.18) and (A.19) reduce to          rσ k(t)≈i tσ(k) (ξn k+ξs k) (ξn k)2−(Es k)2·eie ¯hΦ(t) sσ k(t)≈ − i tσ(k)σ∆k (ξn k)2−(Es k)2·eie ¯hΦ(t) .(A.21) From such a discussion and since φ(τ) = 0, it also follows that drσ k dt t=τ≈dsσ k dt t=τ≈0 and hence we obtain db S dt t=τ = 0. In such a case, we can write eq. (A.11) as Gβ(υ, τ) = Tr c V†(τ, 0) eib H(τ)υc V(τ, 0) e−ib H(0)υb ρ(β),(A.22) where c H(τ) is given by (A.9) (with tsubstituted by τin that equation). Substituting the operator b S(t) by its expression as given in (A.12), with rσ k(t) and sσ k(t) given by (A.21), in (A.9), we obtain that c H(t) = c Hs+c Hn(t), where the operators c Hsand c Hn(t), pertaining, respectively, to the superconductor and the normal-metal, are given by c Hs=X k,σ ˜ ξs kσb c† skσb cskσ+X k∆s kb c† sk↑b c† s−k↓+ ∆s kb cs−k↓b csk↑,(A.23) with ˜ ξs kσ=ξs k−λ2|tσ(k)|2(ξn k+ξs k) (ξn k)2−(Es k)2,(A.24) ∆s k= ∆k· 1−λ2 2·|t↑(k)|2+|t↓(k)|2 (ξn k)2−(Es k)2!,(A.25) The distribution of work performed on a NIS junction 17 being the renormalised kinetic energy and pairing function in the superconductor, and c Hn(t) = X k,σ ˜ ξn kσb c† nkσb cnkσ +X k∆n ke−2ieΦ(t)/¯hb c† nk↑b c† n−k↓+ ∆n ke2ieΦ(t)/¯hb cn−k↓b cnk↑,(A.26) with ˜ ξn kσ=ξn k+λ2|tσ(k)|2(ξn k+ξs k) (ξn k)2−(Es k)2,(A.27) ∆n k=−λ2 2·∆k(|t↑(k)|2+|t↓(k)|2) (ξn k)2−(Es k)2,(A.28) being the renormalised kinetic energy and induced pairing function in the normal-metal, due to the proximity effect. However, such an induced pairing function is multiplied by a time-dependent phase, due to the applied voltage across the junction. If one assumes that the tunnelling square amplitudes are independent of spin, i.e. |t↑(k)|2=|t↓(k)|2, the modified kinetic energies are also independent of the spin one can, introducing in each ((k,↑),(−k,↓)) subspace the operators, b σz k=b c† nk↑b cnk↑+b c† n−k↓b cn−k↓−1,(A.29) b σ+ k=b c† nk↑b c† n−k↓,(A.30) b σ− k=b cn−k↓b cnk↑,(A.31) whose algebra is isomorphic to the spin 1/2 algebra, write c Hn(t), up to a constant factor, as c Hn(t) = X k˜ ξn kb σz k+ ∆n ke−2ieΦ(t)/¯hb σ+ k+ ∆n ke2ieΦ(t)/¯hb σ− k.(A.32) We thus conclude that the Hamiltonian describing the normal-metal is equivalent, in each ((k,↑),(−k,↓)) subspace, to the Hamiltonian of a two-level system under the action of a circularly polarised field. Moreover, since c Hsis time-independent and commutes with c Hn(t), it does not contribute to (A.22). We can thus write this equation as Gβ(υ, τ) = Tr c U†(τ, 0) eib Hn(τ)υc U(τ, 0) e−ib Hn(0)υb ρn(β),(A.33) where c U(t, 0) = Texp −i ¯hZt 0du c Hn(u),(A.34) and b ρn(β) = e−βb Hn(0)/Zn(β), with Zn(β) = Tr(e−βb Hn(0)). Therefore, we have shown that there exists an equivalence, in what concerns the calculation of the work characteristic-function in the adiabatic limit, between the dynamics of the NIS junction and that of an assembly of independent two-level systems, subjected to a circularly polarised field. In Appendix B, we will consider the calculation of the work characteristic-function of the latter system. The distribution of work performed on a NIS junction 18 Appendix B. The work characteristic-function of a two-level system in contact with a thermal bath We now wish to compute (A.33) where c Hn(t) is given by (A.32). Since the different k modes are independent, we can restrict the calculation of the said function to that of a single mode. Furthermore, we will drop the klabel, as no confusion can arise. The derivation of Appendix A shows that one can reduce the dynamics of the system under the applied voltage protocol to the dynamics within the low lying energy states of the normal metal if one assumes that the transition time τ0¯h/∆, see eq. (A.32). We will now further assume that Φ(t)≈¯hωt/2, where ω= 2eVmax/¯h. Such an approximation, if it were exact, would actually imply that the voltage was being quenched. For it to be valid, one must have eVmaxτ0/¯h1, i.e. τ0¯h/∆n, as eVmax ≈∆n. Therefore, the time evolution of the system is only adiabatic with respect to the larger energy gap, being diabatic [102] with respect to the much smaller energy gap, i.e. one must have ¯h/∆τ0¯h/∆nτ. Therefore, for a clear separation between the two time-scales associated with the dynamics of the system to exist, one must have tσ(k)∆, as stated above. We will change somewhat the notation with regard to the previous section, so as to keep the result obtained as general as possible, rather than identifying it solely with the dynamics of the NIS junction. In the new notation, the initial Hamiltonian is given by c Hn(0) = h 2b σz+Γ 2b σx,(B.1) where in the NIS junction context h= 2˜ ξnand Γ = 2∆n, and where, without loss of generality, ∆ncan be chosen to be a real-number, as one can always fix the arbitrary superconducting phase ϕto be zero. The quantities hand Γ can be viewed as the components of a constant (pseudo) magnetic field applied to the two-level system. Since the system is initially in equilibrium with a thermal bath, the system’s partition function is given by Zn(β) = 2 cosh β√Γ2+h2/2. The time-dependent Hamiltonian c Hn(t) can be written as c Hn(t) = h 2b σz+Γ 2( cos(ωt)b σx+ sin(ωt)b σy).(B.2) The dynamics of the two-level system is such that the constant applied field is substituted by a circularly polarised field in [0, τ], i.e. by a Rabi dynamics during that interval. After the application of the protocol, the circularly polarised field is again replaced by a time-independent field and the system is described by the constant Hamiltonian c Hn(τ). In the process, the applied field has been rotated by an angle θ=ωτ around the zaxis. It is well known that for the Rabi dynamics, one can write the time-evolution operator c U(τ, 0) in (A.33) as c U(τ, 0) = c Rτc Uτ, with c Rτ=e−iωτbσz/2representing a The distribution of work performed on a NIS junction 19 pure rotation around the zaxis and c Uτ=e−ib H0τ/¯h, where c H0is a time-independent pseudo-Hamiltonian, which depends on ω, and is given by c H0=1 2(h−¯hω)b σz+Γ 2b σx.(B.3) Substituting this expression for c U(τ, 0) in (A.33) and noting that c R† τeib Hn(τ)uc Rτ= eib R† τb Hn(τ)b Rτu=eib Hn(0)u, since c R† τc Hn(τ)c Rτ=c Hn(0), one obtains for the characteristicfunction the result Gβ(υ, τ) = 1 Zn(β)Tr eib H0τ/¯heib Hn(0)υe−ib H0τ/¯he−ib Hn(0)υe−βb Hn(0).(B.4) Performing the trace over the complete set of states that diagonalises c Hn(0), we obtain for Gβ(υ, τ), the expression Gβ(υ, τ) = 1 Zn(β)X σ,σ0| hσ, b n|eib H0τ/¯h|σ0,b ni |2 ×ei(σ0−σ)υ√Γ2+h2/2e−βσ√Γ2+h2/2,(B.5) where the unit vector b nrefers to the direction of the applied field at t= 0, its components being given by b nx=Γ √Γ2+h2and b nz=h √Γ2+h2. Moreover, writing c H0=1 2qΓ2+ (h−¯hω)2(b n0·b σ), where b n0 x=Γ √Γ2+(h−¯hω)2and b n0 z=h−¯hω √Γ2+(h−¯hω)2, we can expand the exponential eib H0τ/¯has eib H0τ/¯h= cos  τqΓ2+ (h−¯hω)2 2¯h  1 +i(b n0·b n) sin  τqΓ2+ (h−¯hω)2 2¯h (b n·b σ) +isin  τqΓ2+ (h−¯hω)2 2¯h [b n×(b n0×b n)] ·b σ.(B.6) The first two terms of (B.6) are diagonal in the |σ, b nibasis, whereas the last term is only non-zero when evaluated between two states of the |σ, b nibasis with opposite spin. We thus obtain that the square moduli of the matrix elements that appear in eq. (B.5) are given by | hσ, b n|eib H0τ/¯h|σ0,b ni |2= cos2 τqΓ2+ (h−¯hω)2 2¯h  + (b n0·b n)2sin2 τqΓ2+ (h−¯hω)2 2¯h  δσ,σ0 +|b n×(b n0×b n)|2 ×sin2 τqΓ2+ (h−¯hω)2 2¯h δσ,−σ0.(B.7) The distribution of work performed on a NIS junction 20 Taking into account that |b n×(b n0×b n)|2= 1 −(b n0·b n)2and using the expressions for the components of b nand b n0given above in (B.7), we finally obtain the expression Gβ(υ, τ) =  cos2 τqΓ2+ (h−¯hω)2 2¯h  +(Γ2+h2−h¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)·sin2 τqΓ2+ (h−¯hω)2 2¯h   +Γ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)·sin2 τqΓ2+ (h−¯hω)2 2¯h  ×cosh h(β/2 + iυ)√Γ2+h2i cosh β√Γ2+h2/2.(B.8) It is trivial to check that Gβ(υ, τ), as given by (B.8), fulfils both the equality Gβ(0, τ) = 1 and Gβ(iβ, τ) = 1, which provides a check on the correctness of the result, since these equalities are built into the definition of Gβ(υ, τ) by construction, as pointed out above. Performing the inverse Fourier transform on Gβ(u, τ), we obtain for Pβ(W, τ) the result Pβ(W, τ) =  cos2 τqΓ2+ (h−¯hω)2 2¯h + (Γ2+h2−h¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h  δ(W) +Γ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h  ×eβW/2 2 cosh β√Γ2+h2/2 ·hδW−√Γ2+h2+δW+√Γ2+h2i.(B.9) The form of (B.9) is easily interpreted from the two-level structure of the problem, as either no work is performed on the system if the time-dependent field does not induce transitions between the levels (first term), or otherwise the time-dependent field induces a transition between the ground state and the excited state, involving an amount of work performed on the system equal to W=√Γ2+h2, or the inverse transition is induced involving a negative amount of work W=−√Γ2+h2being performed on the system (second term). Note that if one sets W→ −W,ω→ −ωand h→ −h, Γ → −Γ in (B.9), one can directly check that Pβ(−W, −τ) = e−βW Pβ(W, τ), i.e. the work distribution satisfies the Crooks-Tasaki relation. The average work performed on the system during the application of the protocol can be computed either by differentiating (B.8) with respect to υ, or directly from (B.9), The distribution of work performed on a NIS junction 21 and is given by hWiβ=Γ2(¯hω)2 √Γ2+h2(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h  ×tanh β√Γ2+h2/2.(B.10) Since the change in the free energy of the system is zero during the transformation, this quantity is equal to the energy dissipated by the system into the thermal bath during the equilibration process occurring after t > τ. One can easily check that, according to (B.10), hWiβis always larger or equal to zero in agreement with the second law of thermodynamics. It has a maximum at resonance, i.e. if ω= (Γ2+h2)/(¯hh) and if ωτ |Γ|/√Γ2+h2= (2n+1)π, with nbeing an integer. Note that for a single two-level system, one can, for each value of τ, choose ωsuch that hWiβ= 0, i.e. one takes ωsuch that τqΓ2+ (h−¯hω)2/¯h= 2nπ. The existence of such a minimum can be understood from that fact that for such choice of ω,eib H0τ/¯h= 1 , and thus no transitions between the levels are induced by the unitary transformation, see also (B.5). This behaviour regarding the dissipated work in a two-level system mirrors the corresponding behaviour of the Rabi formula for the transition probabilities of such a system under the influence of a circularly polarised field. The mean-square deviation of the work performed during the application of the protocol can be computed either by differentiating (B.8) twice with respect to υ, or directly from (B.9), and is given by h(δW )2iβ=Γ2(¯hω)2 (Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h ×  1−Γ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h  ×tanh β√Γ2+h2/2i.(B.11) In Sect. 4, we use (B.8) to compute the work dissipated by the NIS junction due to the application of the voltage protocol. Since it is also of interest, albeit not for the solution of the NIS problem, we leave to Appendix C the calculation of the work characteristic-function for an isolated two-level system due to the application of a circularly polarised field to such a system. Appendix C. The generating function of an isolated two-level system One can use the results obtained in Appendix B to compute the characteristic function or the work distribution function for an isolated system, described by a micro-canonical ensemble, which undergoes a transformation that is analogous to the one described in that section, i.e. the system is initially isolated and is coupled to the circularly polarised field at t= 0, being decoupled from it at t=τ. Such characteristic function is given The distribution of work performed on a NIS junction 22 by an inverse Laplace transform [49] that involves Gβ(υ, τ) and the partition function Z(β) = 2 cosh β√Γ2+h2/2 GE(υ, τ)ω0(E) = ZC dβ 2πi eβE Gβ(υ, τ)Z(β),(C.1) where Eis the energy of the system and ω0(E) = δE−√Γ2+h2/2+ δE+√Γ2+h2/2is the density of states of the system at t= 0. The contour C is the inverse Laplace transform contour from c−i∞and c+i∞, with cchosen such that all the singularities of the integrand are located to the left of c. In our case, c= 0. Note that these results cannot be generalised to an assembly of two-level systems, since the work characteristic-function Gβ(υ, τ) is in this case the product of characteristic-functions for the individual systems. The integral can be readily performed and after factoring the term ω0(E) out, we obtain for GE(υ, τ) the result GE(υ, τ) =  cos2 τqΓ2+ (h−¯hω)2 2¯h + (C.2) (Γ2+h2−h¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h  + Γ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h e−2iυE . From this expression, one can obtain, as above, by inverse Fourier transformation of GE(υ, τ), the work function distribution for an isolated two-level system. This is given by PE(W, τ) =  cos2 τqΓ2+ (h−¯hω)2 2¯h + (Γ2+h2−h¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h  δ(W) +Γ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2) ×sin2 τqΓ2+ (h−¯hω)2 2¯h δ(W+ 2E).(C.3) This result can be easily interpreted if one again notes that either no transitions occur between the two levels and in this case the work performed is zero (first term), or otherwise the work performed is −2Ewhere Eis the energy of the initial level (second term). Also, note that PE+W(−W, −τ) = PE(W, τ), which is the version of the CrooksTasaki relation appropriate for isolated systems in which the expressions for the microcanonical density of states corresponding to the initial and to the final Hamiltonian are identical, since the spectrum of the system does not change under application of the work protocol [49]. The distribution of work performed on a NIS junction 23 We can also compute the average work dissipated in the transformation, using the expression for GE(υ, τ) or that for PE(W, τ). This is given by hWiE=−2EΓ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h .(C.4) Again, one can minimise or maximise this quantity by an appropriate choice of the value of ω. This result may be relevant in the context of quantum phase-shift gates [103] where the transformation can be implemented at finite frequency without generation of heat. Finally, we can also compute the mean-square deviation of the work performed during the application of the protocol, using the expression for GE(υ, τ) or that for PE(W, τ). This quantity is given by h(δW )2iE=4E2Γ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h ·(C.5)  1−Γ2(¯hω)2 (Γ2+h2)(Γ2+ (h−¯hω)2)sin2 τqΓ2+ (h−¯hω)2 2¯h  . Appendix D. tunnelling matrix element extraction from experiment The power dissipated per atom in a NIS junction in the steady state, when submitted to a constant voltage, is given by [81, 82, 83, 84, 88] P=2πeV ¯hNat X kσ|tσ(k)|2(f(ξn k−eV )−f(ξn k)) h|uk|2δ(ξn k−Es k−eV ) +|vk|2δ(ξn k+Es k−eV )i,(D.1) where f() = 1/(eβ + 1) is the Fermi function and where            |uk|2=1 21 + ξs k Es k |vk|2=1 21−ξs k Es k,(D.2) are the squares of the occupation factors in the superconductor. Using the density of states per atom (and per spin) in the normal metal ρ(ξ) introduced above, we can write this equation as P=4πeV ¯hZ+∞ −∞ dξ |t(ξ)|2ρ(ξ) qξ2+ ∆2(ξ)(f(ξ−eV )−f(ξ))γ(ξ),(D.3) where the function γ(ξ) is given by γ(ξ) = (qξ2+ ∆2(ξ) + eV/2) δ(ξ−eV −qξ2+ ∆2(ξ)) + (qξ2+ ∆2(ξ)−eV/2) δ(ξ−eV +qξ2+ ∆2(ξ)) .(D.4) The distribution of work performed on a NIS junction 24 Assuming that the gap in the superconductor is constant within the Debye window of frequencies, one can perform the above integral, where only the second delta function in the definition of γ(ξ) gives a contribution. We obtain P=2π|t(ξ0)|2∆2ρ(ξ0) sinh(e|V|/2kBT) ¯h e|V|[ cosh(e|V|/2kBT) + cosh(∆2/2e|V|kBT) ] ,(D.5) where ξ0=eV/2−∆2/2eV . This implies that in first-order perturbation theory there is no power dissipated, and hence no current flowing for voltages well below the gap. Note, however, that eq. (15) corresponds to a second-order calculation, since ∆n∝ |t|2. Note that in the model described by eq. (5), the tunnelling process conserves momentum in the plane of the films, and hence the junction is not ohmic when both these films are in the normal state. However, one can measure the power dissipated when e|V|= ∆ in the superconducting state of the aluminium film. We obtain from (D.5) P=π ¯h|t(0)|2ρ(0)∆ tanh(∆/2kBT).(D.6) Therefore, a measurement of Pat e|V|= ∆, as a function of the temperature, allows the extraction of |t(0)|2. Note that ∆(T) is itself a function of the temperature, see [92]. Note that an expression of the tunnelling matrix element in terms of parameters of the system such as the thickness of the insulating layer (see comment on section 3) always requires a model of the potential barrier, see e.g. [104]. The above method avoids the need to resort to such modelling. Appendix E. Techniques for measurement of heat released in the NIS junction Several techniques can be used to measure the small heat release that was computed in section 5. The most straightforward technique would be the use of standard calorimetry [17, 18]. In such a case, if the conversion of the work performed on the film, into heat, occurs in a short time frame after the end of the work protocol, there is an increase ∆T in the temperature of the normal metal film given by ∆T=hWiβ ncv ,(E.1) where nis the number of moles of material contained in the probe, and cvthe material’s molar heat capacity. If we use hWiβ≈0.7 nJ/g, as computed in section 5, we obtain a variation of temperature ∆T= 1.28 ×10−5K, likely to be too small to be measured. A possible alternative technique would be the measurement of the released heat at constant temperature by coupling the system to an external probe that absorbs a quantity of heat, ∆Q=hWiβ, isothermally. This transformation would be analogous to the isothermal expansion phase of a Stirling engine [105]. If an infinitesimal amount of heat δQ is released in the NIS junction during the equilibration process, after application of the protocol (as stated in section 5, the equilibration time is of the order of the electron-phonon relaxation time τpe), we have from the second law of thermodynamics, The distribution of work performed on a NIS junction 25 that δQ/T =dSpr +dSHB, where dSpr is the variation of the entropy of the probe (which would be the elastic degrees of freedom of the normal metal film) and dSHB is the variation of entropy of the heat bath (surrounding helium bath). Since the temperature, volume and number of atoms of the bath do not change, dSHB = 0. Thus, we have δQ/T =dSpr =∂S ∂VT dV=∂P ∂T V dV=− ∂V ∂T |P ∂V ∂P |T dV=α kT dV,(E.2) where we have applied the Maxwell relation on going from the second to the third equality and the implicit function theorem on going from the third to the fourth equality in the equation above, and where κTis the isothermal compressibility of the probe and αits thermal expansion coefficient, as follows from equilibrium thermodynamics. The reasoning above is analogous to the one made when discussing the measurement of the magnetocaloric effect from magnetisation data in isothermal conditions, see Ref. [106]. Note that the stresses in the film and the pressure of the fluid that makes up the heat bath do not need to be equal and moreover, that the mechanical element that controls the expansion of the film would only move in time scales much larger than τpe, necessary for the equilibration of the elastic degrees of freedom of the film. Assuming that the thermodynamic quantities that enter eq. (E.2) are independent of the volume for the small variations measured, we obtain, integrating this equation ∆V=κT αT ∆Q . (E.3) Note that in this case the variation of volume is inversely proportional to the thermal expansion parameter, which implies that the smallness of this quantity in aluminium at low temperatures will actually amplify the effect to be measured. Therefore, the (average) relative variation of the volume of the probe is given by ∆V V=hWiβ γncvT,(E.4) where γ=αvm cvkTis the Gr¨uneisen parameter of the material that constitutes the probe and vmis its molar volume. For aluminium, γ≈1.7 [107] and cv≈1.5 mJ mol−1K−1 [18]. The measurement of the relative variation of the probe’s volume in isothermal conditions thus gives direct experimental access to the work produced during the prescribed protocol. A histogram of the distribution of data obtained for multiple realizations of the protocol will allow the reconstruction of the full work distribution. If we use hWiβ≈0.7 nJ/g, we obtain a relative deviation ∆V V≈7.5×10−6, on average, which should be measurable using, e.g. capacitive methods or SQUID dilatometers, as such a relative deviation is in the limit of precision of capacitive methods [17, 18], provided that other sources of heat dissipation can be minimized. References [1] Callen H B and Welton T A. Irreversibility and generalized noise. Phys. Rev., 83:34–40, 1951.