# Quantum Work Fluctuations in connection with Jarzynski Equality

###### Abstract

A result of great theoretical and experimental interest, Jarzynski equality predicts a free energy change of a system at inverse temperature from an ensemble average of non-equilibrium exponential work, i.e., . The number of experimental work values needed to reach a given accuracy of is determined by the variance of , denoted . We discover in this work that in both harmonic and an-harmonic Hamiltonian systems can systematically diverge in non-adiabatic work protocols, even when the adiabatic protocols do not suffer from such divergence. This divergence may be regarded as a type of dynamically induced phase transition in work fluctuations. For a quantum harmonic oscillator with time-dependent trapping frequency as a working example, any non-adiabatic work protocol is found to yield a diverging at sufficiently low temperatures, markedly different from the classical behavior. The divergence of indicates the too-far-from-equilibrium nature of a non-adiabatic work protocol and makes it compulsory to apply designed control fields to suppress the quantum work fluctuations in order to test Jarzynski equality.

Introduction. Work fluctuation theorems are a central topic in modern non-equilibrium statistical mechanics Hanggireview1 ; Hanggireview2 . One outstanding result is Jarzynski equality (JE) Jarzynski1 which makes use of an average of non-equilibrium exponential work, i.e., (where is the work and is the inverse temperature) to obtain important equilibrium information. It takes the form , where represents an ensemble average of all possible work values with the system initially prepared in a Gibbs state, is the free energy difference between initial and final systems at the same . Regarded as a recent breakthrough in non-equilibrium statistical mechanics, JE holds no matter how the work protocol is executed, fast or slow, adiabatic or highly non-adiabatic. The quantum version of JE takes precisely the same form Kurchan1 ; Tasaki1 , provided that the value of quantum work is interpreted as the energy difference obtained from two energy projective measurements. It is also convenient to define the so-called dissipated work . Then JE gives .

Early applications of JE focused on biomolecular systems, where it is a well-established method to compute Dellago2 . Proof-of-principle experiments were carried out Bustamante1 ; Kiang1 ; Ritort1 following an early proposal by Hummer and Szabo Szabo1 . Other experiments include those using classical oscillators Rabbiosi1 ; Bechinger1 . On the quantum side, the first experiment testing JE was first reported for a quantum spin Serra1 , relying on interferometric schemes to reconstruct work distributions Vedral1 ; Paternostro2 . While other settings have been proposed, e.g., circuit QED Campisi1 , a more recent experimental confirmation, following the proposal Lutz3 , uses a trapped ion subject to an external force Kihwan1 . These experiments are by no means straightforward. In both classical and quantum cases, the ensemble average of a highly nonlinear (in fact, exponential) function of work, as requested by JE, could become practically demanding in yielding a well-converged result for , as rare events can potentially dominate the average Jarzynski1 ; Jarzynski2 . This in fact has motivated a number of studies in the literature to investigate the errors in predicting based on JE Zuckerman1 ; Bustamante2 ; Zuckerman2 ; Dellago1 ; Ritort2 ; Zuckerman3 ; Yi1 ; Kirkpatrick1 . As suggested by the central limit theorem (CLT), the errors are related to the variance in the exponential work, i.e., . The larger is, the more realizations of we must collect to reach a given accuracy in . To obtain within an error of , the number of realizations needed was estimated as Jarzynski3 ; Dellago1 .

The quantity represents the ensemble average of exponential quantities and the Boltzmann exponential cutoff for high-energy states might not always be effective in suppressing the contributions from the high-energy tail. Indeed, in contrast to JE itself which is always well behaved, can diverge, suggesting higher-energy components from the initial Gibbs distribution may contribute even more to the variance. Surprisingly, prior to our work Deng1 , this possibility of divergence was rarely mentioned, with one exception Dellago1 treating an adiabatic work protocol. According to the principle of minimal work fluctuations (PMWF) Gaoyang1 ; Deng1 , once diverges for an adiabatic protocol, then it will suffer from analogous divergence under arbitrary work protocols sharing the same initial and final system Hamiltonians. The accuracy in estimating in such situations becomes problematic. One can no longer rely on the CLT to predict how the error scales with the number of experimental runs. A generalization of the CLT can prove useful Deng1 for a specific situation, but in general the error estimate under a diverging is unknown.

This work discovers the divergence of that is only present for non-adiabatic dynamics. Such a possibility is consistent with PMWF. This represents intriguing and more complex situations, where the “phase” boundary between a finite and a diverging is yet to be located case by case. Our working example below is mainly a quantum harmonic oscillator (QHO) with time-dependent trapping frequency because recent experiments to test quantum JE Kihwan1 and to construct quantum heat engines Singer1 were based on QHO (we consider an-harmonic Hamiltonians in the end). Work is done to the system by varying . For a QHO of fundamental and experimental importance, all quantum transition probabilities can be analytically obtained. To our great surprise, even in such a prototypic system with all dynamical aspects known for so many years, divergence of may be induced by non-adiabatic work protocols. The domain of divergence as a function of temperature and other dynamical parameters is found to possess complicated structures. In general, so long as the dynamics is non-adiabatic, quantum effects at lower temperatures tend to enlarge the domain of divergence in . Indeed, even when the classical is well behaved, its quantum counterpart is bound to diverge at sufficiently low temperatures. This itself constitutes a new aspect of quantum effects in non-equilibrium statistical mechanics. Our results indicate (i) that direct applicability of JE in free energy estimates can become much limited as a system of interest approaches the deep quantum regime, and (ii) that quantum non-adiabatic effects or “inner friction” italy1 alone can induce critical changes in work fluctuations. This work shall also motivate parallel studies on several generalized quantum fluctuation theorems genFT1 ; genFT2 ; genFT3 ; genFT4 ; genFT5 .

Work characteristic function and model system. Throughout we focus on work applied to an isolated system. Let be the time-dependent Hamiltonian subject to a work protocol. The work distribution function is assumed to be . The Fourier transformation of yields the so-called work characteristic function Talkner1

(1) |

For a work protocol starting at and ending at , can be expressed as a quantum correlation function , where is the initial Gibbs state and is the final Hamiltonian in the Heisenberg representation. While the value of at yields JE, the value of at determines . That is,

(2) |

This observation indicates that a divergence in is equivalent to a divergence in or in . The quantum correlation function on the imaginary axis is hence our central object. An operative expression of (1) is given by

(3) |

where represents the transition probability from th state of to th state of . Equations (2) and (3) indicate that, for a system with a finite-dimension Hilbert space such as a spin system, is obtained from a finite summation and will be always well behaved.

To lay a solid ground for our surprising findings, we aim to present our main results with computational and analytical calculations supporting each other. With this in mind, the parametric QHO system seems to be the best choice as a model with an infinite-dimensional Hilbert space. QHO belongs to the algebraic class of integrable systems, for which many solutions in the form of exact propagators can be found in the literature Nikonov ; Lohe . Systems in the same algebraic class share common dynamical features. For example, the Calogero-Sutherland model despite being an interacting system is dynamically analogous to QHO Jaramillo . In addition, in the history of quantum mechanics and statistical physics, QHO plays a pivotal role, providing general insights into quantum zero-point energy, low-temperature behavior of the heat capacity of a crystal, open quantum systems, and the relationship between classical mechanics and quantum mechanics, etc. In particular, in understanding quantum-classical correspondence, tuning the dimensionless ratio of the thermal energy over the energy level spacing is essential to make transitions between classical and quantum regimes.

Consider then a QHO with driving in its trapping frequency, with . Work is done as the trapping frequency changes from to . Thanks to Husimi, all are available as a function of the Husimi coefficient Husimi1 , which is defined as

High-temperature and low-temperature regimes. It is of interest to first examine the behavior of in the strictly adiabatic case, i.e., . Because the convergence criteria of the resulting geometric series in Eq. (3) becomes with , it is obvious that if , then and hence diverge. With PMWF, one further concludes that the divergence occurs in all non-adiabatic work protocols.

Next we investigate the behavior of under non-adiabatic work protocols, provided that adiabatic work protocols do not yield divergence in , i.e., in the regime of . To that end, we partially resort to a compact expression for Husimi1 . That is,

(4) |

Extra care is needed when one extends such a non-analytic result based on integration assuming real to the imaginary axis of . The Jarzynski equality is always recovered from (4) at Lutz1 , i.e., in the regime of . Equation (4) then takes us to the following compact expression for , . Encouraged by this, we cautiously use this compact expression to help us digest the possible critical boundary of divergence in

(5) |

The real denominator of in Eq. (5) approaching zero signifies a boundary separating finite values of from its diverging values. We stress however, all such “shortcut” conclusions are carefully checked against the original expression in Eq. (3), with in the sum series truncated at some large values of and (see Supplementary Material supp ).

For the high-temperature (classical) regime we have

(6) |

The boundary of divergence in is identified at

(7) |

Equation (7) for adiabatic cases still reproduces the known boundary. In generic non-adiabatic cases with , the divergence boundary is pushed to smaller values, i.e., , in agreement with the classical result (see Supplementary Material supp ). Later on we consider specific work protocols to further digest this result. One important feature is that this phase boundary at is no longer dependent upon temperature. Indeed, Eq. (6) shows that itself is temperature-independent in the high-temperature (classical) regime. Expanding the denominator of Eq. (6) around under a given , we find diverges as as approaches the critical value from below.

The behavior of is entirely different in the low-temperature (deep quantum) regime , where the compact expression for becomes

(8) |

For strictly adiabatic cases, i.e., , Eq. (8) is well behaved. However, for any work protocol with even slight non-adiabaticity, , one observes that the denominator in Eq. (8) hits zero or becomes imaginary at sufficiently low temperatures. This indicates that must diverge as temperature decreases, so long as the protocol is not strictly adiabatic. The boundary of divergence in , within the low-temperature regime, is located at , which does not explicitly depend on . Further, from Eq. (8) we find that as approaches , diverges as . To understand the divergence, one can actually focus on the contribution to made by the transitions from -th state of to the ground state of . After applying Sterling’s formula to for large , this contribution is found to scale as as above, because for , the more rare the initial state is (sampled from the initial Boltzmann distribution), the more contribution it makes to . In other words, quantum transitions associated with very negative work values, though exponentially suppressed by with , still not suppressed as sharply as in classical cases so as to lose the competition from the exponential increasing factor contained in at sufficiently low temperatures. This important insight has been confirmed by the behavior of classical and quantum deformed Jarzynzki equalities we just proposed dJE . . This predicts the same

All these predicted features are checked in numerics, with some examples shown in Fig. 1. There the results are obtained based on Eq. (3) directly as well as from Eq. (5). In the high-temperature regime (small ), all the plots vs (in dimensionless units) become flat, in agreement with Eq. (6). This also indicates that the classical does not diverge in these cases. For lower temperatures, however, all the plotted curves tend to blow up, except for the strict adiabatic case (bottom curve) whose approaches zero as temperature decreases, in agreement with Eq. (8). As seen from Fig. 1, local minima in as a function of might emerge, reflecting a competition between quantum fluctuations and thermal fluctuations. In addition, the divergent behavior of close to or (obtained from direct numerics and not shown) is also found to agree with the scaling laws of and .

Specific work protocols at intermediate temperatures. To further investigate the behavior of , we turn to a specific work protocol where the trapping frequency is suddenly quenched from to . In this case . In the high-temperature regime, Eq. (7) yields a critical , namely, diverges if . Fig. 2 depicts the numerically obtained domain of divergence in in terms of and , with panel (a) exactly showing this high-temperature behavior. Quantum effects however dramatically enlarge the domain of divergence (white area) in : From panel (b) to (d), temperature decreases and the domain of divergence in gradually invades the (classical) domain of finite in panel (a). For even lower temperatures than shown in Fig. 2, the domain of convergence (gray area) collapse into the line , suggesting any actual quench in will yield divergence in .

We have also examined a finite-time work protocol, where the parameter is chosen to be time-independent, with delCampo1 . Compared with the sudden quench case, the divergence domain is now fragmented into multiple domains with subtle phase boundaries. In addition, as the temperature decreases, the domain of divergence again grows. Similar disconnected regions of divergence are observed along the time of driving, where convergent domains are found around the adiabatic times delCampo1 ; Uzdin1 ; Rezek1 for which is close to 1. Detailed results are presented in Supplementary Material supp . Finally, it is worth mentioning the possibility of work protocols with , see for example adam1 ; Lutz2 . In these cases may even diverge under compression of the harmonic potential () at both high and low temperatures.

Before ending this section, we also mention our numerical studies of work fluctuations in quantum anharmonic oscillators and other systems such as particle in a box. There we gain similar insights. In particular, we consider a time-dependent quartic potential with the total Hamiltonian , where is the time of driving, as well as a time-independent quartic potential, with the total Hamiltonian . There quantum effects are also found to induce the divergence of even if it does not diverge in the classical (high-temperature) domain. Fig. 3 depicts a few such computational examples. In all these cases, the effect of a quartic potential is to bring the onset of divergence at higher temperatures as compared with the harmonic cases.

Conclusions. JE was hoped to take advantage of non-equilibrium work protocols to estimate the change of free energy in nanoscale and mesoscopic systems. The divergence in presents a hurdle to a direct application of JE. Note that if the divergence is solely induced by non-adiabatic quantum effects, then one promising solution is to use shortcuts to adiabaticity (STA) Deng2 ; delCampo1 ; Paternostro1 ; Paternostro1 ; Muga1 , to go around this divergence and yet still realizing speedy work protocols (but with some price Deffner1 ; Ueda1 ). Indeed, one can even use as a minimization target to design control fields Gaoyang1 ; gaoyang2 . Because quantum effects are seen to enlarge the domain of divergence drastically in the low temperature regime, using JE for free energy estimates in the deep quantum regime will not be fruitful in the absence of a designed control field. Nevertheless, the divergence in offers a new angle to study quantum work fluctuations Hanggireview2 ; Quan1 ; Rahav1 and to characterize the too-far-from-equilibrium nature of a work protocol. One fascinating question is to study how the divergence of may be reflected in the behavior of in the real- domain. Finally, one might wonder if the recent experiment testing quantum JE using a trapped ion Kihwan1 already suffered from the divergence issue exposed here. It turns out that this belongs to a fortunate case without any divergence in (see Supplementary Material supp ). Thus, analyzing possible divergence in exponential work fluctuations will be crucial in guiding the design of future quantum experiments testing JE.

###### Acknowledgements.

Acknowledgements: This work is funded by Singapore MOE Academic Research Fund Tier-2 project (Project No. MOE2014-T2-2- 119 with WBS No. R-144-000-350-112). J.G. also acknowledges encouraging discussions with many colleagues.Appendix

As explained in the main text, all divergences in are traced back to divergences in

(9) |

where are the state-to-state transition probabilities and is the Gibbs distribution of the initial state. All our numeric calculations involve a truncation in the series (9).

## Appendix A Transition Probabilities

The transition probabilities between instant eigenstates at and for the quantum harmonic oscillator (QHO) with arbitrary driving of the trapping frequency are given by Husimi1

(10a) | |||||

(10b) |

where, , and is the so called Husimi coefficient. The selection rules prevent transitions between energy levels with different parity. The definition of is Husimi1

(11) |

where and are the two solutions of the equation of motion of the corresponding classical oscillator, with initial conditions: and . Any non-adiabaticity is captured by deviations of from its adiabatic value, . The used in the main text is assumed to be .

## Appendix B Work fluctuations for the classical harmonic oscillator

For classical systems the work distribution is given by

(12) |

where is the initial classical Hamiltonian and represents the value of the final Hamiltonian for the trajectory emanating from .
For the classical harmonic oscillator (CHO) this integral results in two different expressions, associated to the cases: and , where is the Husimi coefficient defined in (11) and stands for the Husimi coefficient in sudden quench cases.

We first consider the regime , which will be the case if, for example, we have a monotonically increasing or decreasing from to . This gives rise to positive definite work or negative definite work. For compression of the harmonic trap the work distribution reads Deng2

(13) |

where is the compression ratio, is the step function and is the modified Bessel function of the first kind. This work distribution coincides with that of the QHO in the semiclassical regime vanZon1 . The asymptotic behavior determines the convergence of

(14) |

One can then use the approximation: , to get the tail

(15) |

Evaluation of requires us to integrate along from to . A divergence of will take place when the exponent in (15) is an increasing function of . Fortunately, the divergence does not show up for compression of the trap. However, for expansion of the trap the work is always negative, the work function reading

(16) |

Similarly, the tail is

(17) |

Since the term is an integration along from to , the possibility of a negative coefficient, i.e.,

(18) |

allows divergence. In particular, for adiabatic expansion of the trap () the criteria of divergence is .

Now consider the regime where the nonadiabaticity parameter is even larger than that associated with a sudden quench, i.e., . This is possible if we introduce multiple quenches to the trapping frequency . In this case, the integral (12) for compression or expansion becomes

(19) |

where is the Macdonald function or modified Bessel function of the third kind. This result is consistent with that reported in Lutz1 ; vanZon1 using a different approach starting from the characteristic function of the quantum harmonic oscillator (QHO) and taking a semiclassical limit. Similarly, with the approximation , one obtains the tail

(20) |

Because now can take both positive and negative values, the tail must go to zero at both in order to have a finite integral over . The condition for divergence after integration over is

(21a) | |||

(21b) |

Briefly,

(22) |

Combining these results and after a few straightforward steps, we finally reach a compact condition for divergence

(23) |

Classical characteristic function. First we consider the case . From the integral formula

(24) |

one obtains from (14) the compact expression

(25) |

With some algebra is easy to show that this result coincides with the high-temperature limit of the quantum characteristic function at , depicted in Eq. (7) of the main text. Next we consider the case . From the integral

(26) |

and Eq. (19) it is easy to see that takes the same form as , but in the domain .

## Appendix C Work protocols

The existence of additional domains of divergence in away of adiabaticity is allowed by the PMWF. We study such behavior for two characteristic work protocols: (i) Sudden quench, , and (ii) a protocol with independent of the time of evolution, covering the regime for different values of the time of driving . On the other hand, protocols involving multiple quenches Lutz2 ; adam1 report greater than sudden quench, . Although not treated here, it is worth mentioning that such protocols with can lead to divergences in in the high-temperature regime not only for expansion protocols but also for compression.

We briefly recall that

(27) |

where is the characteristic function of work. For the quantum harmonic oscillator

(28) |

where are the corresponding transition probabilities. Our numerical approach is based on truncation of the quantum numbers and in the latter expression. Using the generating function method one arrives at the compact expression (to be derived in the next section)

(29) |

In particular,

(30) |

For the high-temperature (classical) regime , we have

(31) |

while for the low-temperature (deep quantum) regime ,

(32) |

We now discuss the specific work protocols.

### c.1 Sudden quench

For the case of sudden quench, , the domain of divergence in for the high-temperature regime (classical) is , as predicted by Eq. (31). As we depart from the high-temperature limit, an extra domain of divergence emerges from the regime of , reaching smaller values with lower temperatures, see panels (b) to (e) in Fig. 4. This new domain divergence originates from quantum effects as they are absent in the limit of . We can further digest the divergence by looking for the zero and imaginary values of the denominator in the characteristic function in Eq. (30), i.e.,

(33) |

To gain useful insights we consider the limit . Using , Eq. (33) reduces to

(34) |

Noting the obvious inequality

(35) |

the existence of a divergence requires

(36) |

This divergence is not present in the high-temperature (classical) regime, being associated to . Panels (c) to (e) in Fig. 4 show such boundary extending to , but always within the compression sector, . As an example, the inequality (36) predicts an onset of divergence around for the temperature , in the limit . This is in agreement with Fig. 4(d).

### c.2 A protocol with constant .

The second protocol allows us to study the impact of a varying protocol duration , between sudden quench () and the slow driving limit (). The explicit protocol is specified by delCampo1

(37) |

Thanks to a time-independent , one can readily obtain the corresponding Husimi coefficient

(38) |

where , and . Depicted in Fig. 5 is the behavior of for a time of driving and temperatures .

Compared with the sudden quench case, the divergent domains (see Fig. 5(b) and Fig. 5(d)) are now fragmented into multiple domains with subtle phase boundaries. This fragmentation can be partially traced back to the non-monotonic behavior of as a function of and (see dashed lines in Fig. 5(a) and Fig. 5(c)). Indeed, for , then , where .

In addition, as the temperature decreases, the domain of divergence tends to grow (the white area in Fig. 5(c) is greater than in Fig. 5(b)). In the low temperature limit, the domain for having finite (gray area) shrinks to almost zero area, leaving only zero-measure boundaries separating different domains of divergence. Such zero-measure cases with finite arise from finite-time adiabatic dynamics whose , which is possible because for