Asymptotics of work distributions in a stochastically driven system

TL;DR

This paper extends the EN theory to stochastically driven overdamped Langevin systems, providing exact asymptotic work distribution forms at arbitrary times and validating results with numerical simulations.

Contribution

The paper develops an extension of EN theory to stochastic driving systems, deriving new analytic solutions for work distributions in such models.

Findings

Derived asymptotic forms of work distributions at arbitrary times.

Validated analytic solutions with numerical simulations.

Compared steady-state dissipated work behavior with previous results.

Abstract

We determine the asymptotic forms of work distributions at arbitrary times $T$, in a class of driven stochastic systems using a theory developed by Engel and Nickelsen (EN theory) (arXiv:1102.4505v1 [cond-mat.stat-mech]), which is based on the contraction principle of large deviation theory. In this paper, we extend the theory, previously applied in the context of deterministically driven systems, to a model in which the driving is stochastic. The models we study are described by overdamped Langevin equations and the work distributions in the path integral form, are characterised by having quadratic actions. We first illustrate EN theory, for a deterministically driven system - the breathing parabola model, and show that within its framework, the Crooks flucutation theorem manifests itself as a reflection symmetry property of a certain characteristic polynomial function. We then extend our analysis to a stochastically driven system, studied in ( arXiv:1212.0704v2 [cond-mat.stat-mech], arXiv:1402.5777v1 [cond-mat.stat-mech]) using a moment-generating-function method, for both equilibrium and non - equilibrium steady state initial distributions. In both cases we obtain new analytic solutions for the asymptotic forms of (dissipated) work distributions at arbitrary $T$. For dissipated work in the steady state, we compare the large $T$ asymptotic behaviour of our solution to that already obtained in ( arXiv:1402.5777v1 [cond-mat.stat-mech]). In all cases, special emphasis is placed on the computation of the pre-exponential factor and the results show excellent agreement with the numerical simulations. Our solutions are exact in the low noise limit.

Tables

Table 1. Table 1: Asymptotic behaviour of work distributions. List of asymptotic forms of P ​ ( W ) 𝑃 𝑊 P(W) known so far adapted from Holubec . A , B , C 𝐴 𝐵 𝐶 A,B,C are constants that depend upon the explicit form of the driving protocol and the duration of the protocol T 𝑇 T .

Potential U(x(t), $\lambda \left(t\right)$)	Tail behaviour
$\frac{\kappa }{2}{\left[x\left(t\right)-\lambda \left(t\right)\right]}^2$	$A{e}^{-{\left(BW-C\right)}^2}$
$\frac{1}{2}\lambda \left(t\right)x{\left(t\right)}^2$	$A{\left|W\right|}^{-1/2}{e}^{-B\left|W\right|}$
$-g\text{log }\left|x\right|+\frac{1}{2}\lambda \left(t\right)x{\left(t\right)}^2$	$A{\left|W\right|}^{-\left(1-\beta g\right)/2}{e}^{-B\left|W\right|}$
$\lambda \left(t\right)\left|x\left(t\right)\right|$	$A{e}^{-B\left|W\right|}$

Table 2. Table 2: Asymptotic Forms of P ​ ( W , T ) 𝑃 𝑊 𝑇 P(W,\;T) . The exact asymptotic form of P ​ ( W ) 𝑃 𝑊 P(W) can be determined for any value of T 𝑇 T . Table 2 gives the exactly computed form for a few T 𝑇 T values, which may be compared with experiments / numerical simulations.

T	$P\left(W^{+}\right)$	$P\left(W^{-}\right)$
0.3	$\frac{0.49}{\sqrt{\left|W\right|}}{e}^{-1.04\left|W\right|}$	$\frac{0.49}{\sqrt{\left|W\right|}}{e}^{-2.04\left|W\right|}$
0.5	$\frac{0.48}{\sqrt{\left|W\right|}}{e}^{-0.819\left|W\right|}$	$\frac{0.48}{\sqrt{\left|W\right|}}{e}^{-1.819\left|W\right|}$
0.7	$\frac{0.48}{\sqrt{\left|W\right|}}{e}^{-0.72\left|W\right|}$	$\frac{0.48}{\sqrt{\left|W\right|}}{e}^{-1.72\left|W\right|}$
1	$\frac{0.52}{\sqrt{\left|W\right|}}{e}^{-0.65\left|W\right|}$	$\frac{0.52}{\sqrt{\left|W\right|}}{e}^{-1.65\left|W\right|}$

Table 3. Table 3: Asymptotic Forms for P ​ ( W d , T ) 𝑃 subscript 𝑊 𝑑 𝑇 P(W_{d},\;T) . The exact asymptotic forms of P ​ ( W d ) 𝑃 subscript 𝑊 𝑑 P(W_{d}) for different values of T 𝑇 T , obtained using Eq. ( 78 ).

T	$P\left(W_d^{+}\right)$	$P\left(W_d^{-}\right)$
1	$\frac{1.13}{\sqrt{\left|W_d\right|}}{e}^{-1.08\left|W_d\right|}$	$\frac{1.13}{\sqrt{\left|W_d\right|}}{e}^{-2.08\left|W_d\right|}$
2	$\frac{3.30}{\sqrt{\left|W_d\right|}}{e}^{-0.88\left|W_d\right|}$	$\frac{3.30}{\sqrt{\left|W_d\right|}}{e}^{-1.88\left|W_d\right|}$
3	$\frac{7.44}{\sqrt{\left|W_d\right|}}{e}^{-0.78\left|W_d\right|}$	$\frac{7.44}{\sqrt{\left|W_d\right|}}{e}^{-1.78\left|W_d\right|}$
4	$\frac{10.44}{\sqrt{\left|W_d\right|}}{e}^{-0.72\left|W_d\right|}$	$\frac{10.44}{\sqrt{\left|W_d\right|}}{e}^{-1.72\left|W_d\right|}$

Equations

$$P(W)\sim \exp (-\beta \tilde{S}),$$

$$V(x)=\frac{1}{2}\lambda (t)x(t{)}^2,$$

$$\dot{x}(t)=-V^{\prime }(x)+\sqrt{\frac{2}{\beta }}\eta (t),$$

$$W\left[x(\cdot )\right]\equiv {\int }_0^{T}dt\frac{\partial V}{\partial \lambda }\dot{\lambda }=\frac{1}{2}{\int }_0^{T}dt\dot{\lambda }(t)x(t{)}^2.$$

$$P(W)=\frac{\text{N}}{Z_0}\int d{x}_0\int d{x}_T\int \frac{dq}{4\pi /\beta }{\int }_{x(0)=x_0}^{x(T)=x_T}\mathcal{D}[x]e^{-\beta S[x,q]},$$

$$S[x,q]=V_0(x_0)+{\int }_0^{T}dt\left(\frac{1}{4}[\dot{x}+V^{\prime }(x){]}^2+\frac{iq}{4}\dot{\lambda }(t)x(t{)}^2\right)-\frac{iq}{2}W\equiv S_W-\frac{iq}{2}W.$$

$$S[x,q]=\tilde{S}+S_{lin}+S_{quad},$$

$$\tilde{S}$$

$$\begin{array}{rl}{S}_{lin}& =\frac{\partial {\tilde{V}}_0}{\partial x_0}{y}_0-\frac{1}{2}[{\dot{\tilde{x}}}_0+{\lambda }_0{\tilde{x}}_0]y_0+\frac{1}{2}\int dt\text{A}\tilde{x}y+\frac{1}{2}[{\dot{\tilde{x}}}_T+{\lambda }_T{\tilde{x}}_T]y_T\\ & +\frac{ir}{2}(\int dt\frac{\partial \tilde{V}}{\partial \lambda }\dot{\lambda }-W),\end{array}$$

$$\displaystyle\begin{split}S_{quad}&=\frac{1}{2}\frac{\partial^{2}\tilde{V}_{0}}{\partial x_{0}^{2}}y_{0}^{2}-\frac{1}{4}(\;\lambda_{0}y_{0}+\dot{y}_{0}\;)y_{0}+\frac{1}{4}\int dt\;y\;\textbf{A}\;y+\frac{1}{4}(\;\lambda_{T}y_{T}+\dot{y}_{T}\;)y_{T}\\[5.0pt] &+\frac{ir}{2}\int dt\;\bigg{(}\;2\;\dot{\lambda}\;\tilde{x}(t)\;y(t)\;\bigg{)},\end{split}$$

$$\text{A}=-\frac{d^2}{d{t}^2}-\left(\left(1-iq\right)\dot{\lambda }-{\lambda }^2\right).$$

$$P(W)\sim e^{-\beta \tilde{S}},$$

$$\text{A}\tilde{x}=-\ddot{\tilde{x}}(t)-\left(\left(1-i\tilde{q}\right)\dot{\lambda }-{\lambda }^2\right)\tilde{x}=0,$$

$${\dot{\tilde{x}}}_0={\lambda }_0{\tilde{x}}_0,{\dot{\tilde{x}}}_T=-{\lambda }_T{\tilde{x}}_T.$$

$$W=W[\tilde{x}(\cdot )].$$

$$\tilde{S}=S[\tilde{x},\tilde{q}]=-\frac{i{\tilde{q}}^{\ast }}{2}W.$$

$$P(W)\sim e^{\beta \frac{i{\tilde{q}}^{\ast }}{2}W}.$$

$$\text{I}:=\int d{y}_0\int d{y}_T\int \frac{dr}{4\pi /\beta }{\int }_{y_0}^{y_T}Dy(\cdot )e^{-\beta S_{quad}}.$$

$$y(\cdot )={\Sigma }_n{c}_n{\varphi }_n(t)$$

$$\text{I}=\frac{\text{J}}{\sqrt{\frac{4\pi }{\beta }}}\frac{1}{\sqrt{d_0^2\det {\text{A}}_{i\tilde{q}=i{\tilde{q}}^{\ast }}^{\prime }}}.$$

$$d_0\equiv \frac{1}{||\tilde{x}||}{\int }_0^{T}dt\left(\dot{\lambda }(t){\tilde{x}}^2(t)\right),$$

$$\text{J}=\frac{Z_0}{\text{N}}\times \sqrt{\det {\text{A}}_{iq=0}}.$$

$$P(W)=2\times \sqrt{\frac{\beta }{4\pi d_0^2}}\times \sqrt{\frac{\det {\text{A}}_{iq=0}}{\det {\text{A}}_{iq=i{\tilde{q}}^{\ast }}^{\prime }}}\times e^{\beta \frac{i{\tilde{q}}^{\ast }}{2}W}(1+O(1/\beta )).$$

$$P_F(W)=e^{\beta W-\beta \Delta F}{P}_R(W).$$

$$\lambda (t)=\frac{1}{2-t},t=0\text{ to }1,$$

$$\lambda (t)=\frac{1}{1+t},t=0\text{ to }1.$$

$$\displaystyle M\;\left[\begin{array}[]{c}x_{0}\\ \dot{x}_{0}\end{array}\right]$$

$$M$$

$$F(k=1-iq)=\text{Det}\left[M+N{H}_k(T)\right],$$

$$H_{k}(t)=\left[\begin{array}[]{cc}x_{1}(t)&x_{2}(t)\\ \dot{x}_{1}(t)&\dot{x}_{2}(t)\end{array}\right].$$

$$k\equiv 1-iq.$$

Full text

Key words:

Large deviation theory, Fluctuation theorems, Functional determinants, Stochastic thermodynamics.

1 Introduction

Stochastic thermodynamics extends the definition of thermodynamic quantities such as entropy, heat and work, to the level of stochastic trajectories, and has become an important area of research in non-equilibrium statistical mechanics of small systems Seifert . Non-equilibrium fluctuations are particularly relevant for systems with a small number of degrees of freedom, where there are large fluctuations around the average and considerable contributions from rare events. Hence a major advancement in the field was due to the discovery of fluctuation theorems (FTs) which extend equilibrium fluctuation-dissipation relations to far-from-equilibrium regimes. Depending on the initial conditions of the system under consideration, FTs put constraints on various thermodynamic distributions and demonstrate that the second law holds statistically in systems with stochastic dynamics; a positive entropy production is exponentially more likely than its negative counterpart. Among the various versions of the FT, the notable ones are the Crooks Fluctuation Theorem (CFT) Crooks , and the Jarzynski Equality (JE) Jarzynski1 ; Jarzynski2 which relate non-equilibrium quantities such as the thermodynamic (Jarzynski) work to equilibrium quantities such as the free energy difference. Advancements in technology have allowed experimental verifications of these fluctuation theorems in a variety of systems (see references in Seifert ).

Relatedly, there has been significant interest in calculating work distributions for different models. A well studied example is that of a colloidal particle in a harmonic optical trap, where either the mean position of the trap Vanzone1 ; Vanzone2 ; van3 ; van4 or the stiffness speck ; Nickelsen2011 ; bretpjar ; bretpprl ; Park is externally modulated. In the literature these models are referred to as the sliding parabola or the breathing parabola respectively. For the sliding parabola, if the driving is deterministic, the work distribution is known to be a Gaussian Vanzone1 . For the breathing parabola, the solution given in speck is formally exact for an arbitrary driving protocol. In loghar , an exact calculation of the work distribution has been carried out for the Brownian particle in a logarithmic-harmonic potential. In Saha2015 , exact work statistics have been obtained for a Brownian particle in the presence of non-conservative forces such as torques.

It is however hard to find an exact expression for the full work distribution $P(W)$ except in the few cases mentioned above. This problem is resolved to a certain extent by using techniques based on large deviation theory Touchette . One of the recent developments in this regard is a theory developed in Nickelsen2011 by Engel and Nickelsen (EN Theory), which is used to compute the tail forms of work distributions including the pre-exponential factor at arbitrary times $T$. To derive the asymptotic probability for a certain rare work value, the probability of the most likely trajectory that gives rise to this specific work value is considered. The constraint to rare work values is effectively equivalent to a low temperature limit. The advantage of this method over other similar large deviation techniques sabha1 ; Sabhapandit is that the results for the tails are also valid for very small time durations of the process which are typical for experimental situations. EN theory also has the additional benefit that it reduces the problem of computing the asymptotic form of the work distributions to solving a system of ordinary differential equations. One can therefore use available BVP solvers bvp to obtain the full asymptotic form Nickelsen2011 ; Holubec . In addition to the problems studied in Nickelsen2011 , a few more systems have been studies using EN theory Nickelsen2 ; Ryabov . In Holubec , Holubec et al put forward a functional-form conjecture to classify the exact asymptotic form of $P(W)$ based on the form of the driving used. In addition, they also obtain the analytic solution to the asymptotic form of $P(W)$ in an absolute value potential (V-potential). Table 1 is adapted from Holubec , where the so far known results for the asymptotic forms of work distributions for various driving protocols are listed. In all the above cases the driving protocol $\lambda(t)$ considered, is a deterministic function of time.

In this paper, we will first revisit a model with a deterministic driving protocol, the breathing parabola Nickelsen2011 , to familiarize the reader with the methods of EN theory. Then, in this simple model, by considering a particular forward and reverse protocol, we will illustrate the fluctuation theorem for the Jarzynski work. We will show that, within the framework of the EN theory, the fluctuation theorem manifests itself as a symmetry property of a certain characteristic polynomial, which also determines the exact moment generating function at finite times. We then extend this analysis to a stochastically driven system (i.e. $\lambda(t)$ is stochastic). Stochastic driving protocols have been looked at previously in sabha1 ; Sabhapandit ; Varley and recently in unc . In sabha1 , $\lambda(t)$ is taken to be a Gaussian random noise, and the large-time asymptotic form of the distribution of work done by the stochastic force $\lambda(t)$ (the product of the stochastic force $\lambda$ times the displacement $dx$) has been computed, and fluctuation theorems analysed. In Sabhapandit , $\lambda(t)$ is considered instead to be the Ornstein Uhlenbeck process. This is also the stochastic protocol that we study in this paper. The model describes the dynamics of a colloidal particle in a harmonic potential whose mean position is stochastically modulated; we call this model the stochastic sliding parabola (SSP) model. The model studied in unc in the context of the recently discussed finite time thermodynamics uncertainty relation, is equivalent to a discrete version of the SSP. In (Sabhapandit, ), it was claimed that the work fluctuation theorem was valid only in certain regions of the parameter space. Subsequently it was shown in Varley , that when the dissipated work, $W_{d}$ ( defined as the difference between the Jarzynski work, that differs from the definition in sabha1 ; Sabhapandit by a boundary term, and the equilibrium free energy difference $\Delta F$) is considered, the Crooks fluctuation theorem is satisfied in all regions of the parameter space. In all these cases, the authors used a moment generating function method, and the $P(W)$ obtained is valid only in the large $T$ limit. In this paper, we study the SSP for both equilibrium and non-equilibrium initial conditions using EN theory. We show that when the initial points are sampled from an equilibrium distribution, the Jarzynski work ($W$) satisfies a transient fluctuation theorem. To our knowledge, this has not been noted before in the context of this model. We compute the exact forms of the tails of the PDF for different time durations, and compare them with numerical simulations. For the special case when $\lambda(t=0)$ is unconstrained, we are able to obtain the closed asymptotic form of $P(W)$ as a function of $T$. For non-equilibrium steady state initial conditions, we obtain the asymptotic form of $P(W_{d})$ in a similar manner and compare its large-time behaviour with the asymptotic form obtained using the results in Varley . Our comparison shows that, for very large times, the leading-order asymptotic form of both methods match but they disagree in the sub-leading pre-exponential behaviour. We validate all our results using numerical simulations. In all cases, special emphasis is put on the computation of the pre-exponential factor, which involves calculating a fluctuation determinant feynman ; schulman that is a functional determinant ( of a matrix differential operator Kirsten ; functional in the case of the SSP). We do this by using a technique developed in Kirsten , which is based on the spectral -$\zeta$ functions of Sturm-Liouville type operators.

The paper is organized as follows. In Section 2, we introduce various techniques and methods used in this article and illustrate EN theory in the context of the breathing parabola problem. Particularly, in Section 2.2, we discuss the FT for a specific choice of forward and reverse protocols and obtain the asymptotic forms of $P(W)$ in each case. In Section 3 we extend this analysis to the SSP model. Exact asymptotic forms of PDFs of the Jarzynski work ($W$) and the dissipated work ($W_{d}$) are obtained for equilibrium and non-equilibrium steady state initial conditions in Sections 3.1 and 3.2 respectively. In Section 3.2.2, we compare the asymptotic forms of $P(W_{d})$ with the results from Varley . In Section 4 we present our conclusions. Various technical details of the paper including the computation of the pre-exponential factor, are presented in Appendices A - E.

2 Basic methods

For systems modelled using overdamped Langevin equations, EN theory may be used to compute the asymptotic form of the work probability distribution analytically Nickelsen2011 ; engel . The theory was initially developed for equilibrium initial conditions, and later generalized to non-equilibrium initial conditions as well Ryabov . The method can be summarized as follows:

Let $W[x(\cdot)]$ be any functional of the stochastic trajectory $x(\cdot)$. Write down $P(W[x(\cdot)]=W)$ as a path integral, with the corresponding action $S[x(\cdot)]$, by constraining the trajectories to have a work value $W$.

• 2.

The first order approximation of $P(W)$ for rare $W$ is obtained as

[TABLE]

where $\tilde{S}$ is the action $S$ evaluated along the optimal trajectory that minimizes the action, found by solving the corresponding Euler Lagrange equations together with natural boundary conditions.

• 3.

An improved estimate is then obtained by including the pre-exponential factor, which also takes into account the contributions from the trajectories lying close to the optimal trajectory. This is done by expanding $S$ to second order in variations around the optimal trajectory and by performing Gaussian integrations over the variations.

As is usual for Gaussian integrals, the computation of the pre-exponential factor requires the evaluation of the ratio of determinants of certain differential operators (functional determinants). In this article we use the generalization of a method developed in Kirsten , based on the spectral -$\zeta$ function of Sturm-Liouville operators to determine this ratio. In section 2.1, we illustrate EN theory using the breathing parabola model. This model has been studied by Engel and Nickelsen Nickelsen2011 and the full asymptotic form of the work distribution including the pre-exponential factor has been computed. This model serves as a good starting point to the discussions that follow in the paper since the techniques we will use later are the generalizations of the techniques presented here. We will stick to the notations used in Nickelsen2011 unless required otherwise.

2.1 Illustrative example: The breathing parabola

The breathing parabola potential is given by

[TABLE]

where the stiffness of the trap $\lambda(t)$ is the driving protocol which is varied deterministically from a value $\lambda(t=0)=\lambda_{0}$ to $\lambda(t=T)=\lambda_{T}$ during each realization of the process. The motion of a colloidal particle in this potential can be described by the overdamped Langevin equation,

[TABLE]

where $V^{\prime}(x)=dV/dx$ and $\beta$ is the inverse temperature defined as $1/(k_{B}T)$. $\eta(t)$ is a Gaussian white noise, with $\langle\eta(t)\rangle=0$, and $\langle\eta(t)\;\eta(s)\rangle=\delta(t-s)$. The Jarzynski work done along a trajectory $x(\cdot)$ of this system for a time interval $\left[0,T\right]$ is defined sekimoto1 ; sekimoto2 as,

[TABLE]

$W[x(\cdot)]$ being a functional of a stochastic trajectory, is a random variable by itself. Using the Onsager Machlup formalism Onsager1 ; Onsager2 , the probability of $W[x(\cdot)]$ taking a value $W$ can be written down as a path integral Nickelsen2011 ,

[TABLE]

where the augmented action $S$ is given by,

[TABLE]

N is the normalization constant corresponding to mid point discretization in the functional integral and $Z_{0}$ is the initial equilibrium partition function. EN theory uses the contraction principle of large deviation theory Touchette to calculate the asymptotic behaviour of $P(W)$ for large $W$. This is implemented by first evaluating the integrals in Eq. (6) using the saddle point approximation, and finding the optimal trajectory $(\tilde{x}(\cdot),\tilde{q})$ which minimizes $S$. To find the optimal trajectory, $S$ is studied near the vicinity of a trajectory $\tilde{x}(t)$ and a value $\tilde{q}$ of $q$ by writing $x(t)=\tilde{x}(t)+y(t)$ and $q=\tilde{q}+r$ and by expanding $S$ to second order333 Notice that $S$ will not have expansion terms of order $>2$ in $y(\cdot)$ and $r$. This is because of the quadratic form of the action Eq. (6). in $y(\cdot)$ and $r$ as,

[TABLE]

where

[TABLE]

where A is a second order Sturm-Liouville type differential operator,

[TABLE]

2.1.1 Leading-order behaviour

The leading-order asymptotic behaviour of $P(W)$ is obtained by approximating,

[TABLE]

where $\tilde{S}$ is the action evaluated along the optimal trajectory $\tilde{x}(\cdot)$ and $\tilde{q}$, obtained by demanding that $S_{lin}$ vanishes for an arbitrary variation $\left(y(\cdot),r\right)$. This results in the Euler-Lagrange equations,

[TABLE]

with boundary conditions,

[TABLE]

The above equations constitute a Sturm-Liouville eigenvalue problem with the parameter $i\tilde{q}$. Therefore there are infinitely many values $i\tilde{q}^{(n)}$ for which a non trivial solution exists, all of which are saddle points of the action $S$. It can be verified that each $i\tilde{q}^{(n)}$ value corresponds to two solutions, $\pm\tilde{x}^{n}(t)$. It is possible to show that all saddle points except the one corresponding to the smallest $i\tilde{q}(\equiv i\tilde{q}^{*})$ value are unstable because they fail the Hessian (second derivative) test444Since the action is quadratic, the Hessian operator in this case, is the same operator as the operator for the optimal trajectory itself. This means that $\tilde{x}^{(n)}$ is an eigenfunction of the Hessian operator with eigenvalue [math]. According to the Courant nodal theorem courant , $\tilde{x}^{(n)}$ has $n$ nodes. Taken together, these two points imply that the Hessian operator corresponding to $i\tilde{q}^{(n)}$ will have $(n-1)$ eigenfunctions with negative eigenvalues. Therefore the saddle points $x^{(n)}$ for $n>1$, are unstable and do not contribute to the asymptote of $P(W)$.. Therefore they do not contribute to the asymptote of $P(W)$. The term proportional to $r$ in Eq. (8b) gives us the constraint equation,

[TABLE]

Using Eq. (11), (12) and Eq. (13), it can be shown that along the optimal trajectory,

[TABLE]

Therefore using Eq. (10), the leading-order asymptotic behaviour of $P(W)$ is determined to be,

[TABLE]

2.1.2 The pre-exponential factor

The next step is to improve this estimate by also taking into account the fluctuations around the optimal trajectory. This is done by retaining $S_{quad}$ (Eq. (8c)) in the exponent of Eq. (5) by computing the Gaussian integrals in

[TABLE]

As already mentioned, since the action is quadratic, the fluctuation governing operator A that appear in $S_{quad}$ is the same as the operator for the optimal trajectory. The operator A therefore has a zero-mode which is nothing but the optimal trajectory itself. Writing integration variable $y$ as a series expansion in terms of the normalized eigenfunctions $\phi_{n}(t)$ of A as,

[TABLE]

and then performing the Gaussian integrals over the expansion parameters $c_{n}$ and $r$, it can be shown that the zero-mode of A gets omitted naturally and does not cause any problems to the integral in Eq. (16). This gives us a compact expression for the pre-exponential factor,

[TABLE]

The notation $\det\textbf{A}^{\prime}_{\;i\tilde{q}=i\tilde{q}^{*}\;}$ stands for a determinant omitting the zero-mode. The factor $d_{0}$ is defined as,

[TABLE]

where $\tilde{x}(t)$ is the zero-mode and $||\tilde{x}||$ stands for its norm. J is a factor stemming from the Jacobian of the transformation of integration variables. The value of J for the breathing parabola problem and also a class of similar potentials has been determined in Nickelsen2011 , and can be shown to be equal to,

[TABLE]

A derivation of Eq. (18) and (19) can be found in Nickelsen2011 . Now using Eqns. (20), (18), (15) and (5) we finally obtain,

[TABLE]

The factor 2 corresponds to the two equipotent saddle points $(i\tilde{q}^{*},\pm\tilde{x}(t))$. Computing the pre-exponential factor now reduces to evaluating $d_{0}$ and the ratio of functional determinants appearing in the expression above. As we will see, computing $d_{0}$ is rather straightforward. The evaluation of the determinant ratio is however more involved and is carried out using a technique developed in (Kirsten, ). Notice also that in Eq. (21), the protocol $\lambda(t)$ is not specified. In the next Section, we apply the method to a specific forward and the corresponding reverse protocol and exactly compute the asymptotic form of $P(W)$ including the pre-exponential factor, for both cases. We also analyse the FT within the framework of EN theory.

2.2 Forward and reverse protocols: Fluctuation relation

In this Section, we will compute the exact asymptotic form of the work distribution in the breathing parabola problem, for a specific choice of forward (F) and reverse (R) protocols. The work distributions are then known to satisfy the Crooks fluctuation theorem:

[TABLE]

We consider the following forward protocol,

[TABLE]

and the corresponding reverse protocol,

[TABLE]

This particular reverse protocol was studied as an illustrative example in Nickelsen2011 . In order to find the leading-order behaviour of $P(W)$ in each case, one has to obtain the smallest $iq$ values for which the ELE Eqns. (11), (12) have a non trivial solution. Here we use a formalism used in Kirsten . As we will see, this method will turn out to be useful for the computation of the pre-exponential factor as well. First we write the boundary conditions as two matrix equations (for notational simplicity, we will use $x(t)$ instead of $\tilde{x}(t)$.),

[TABLE]

In this case, M and N are matrices,

[TABLE]

It can then be shown that the $iq$ values for which the equations (11), (12) have a non trivial solution, can be obtained as the roots of a function (the characteristic polynomial),

[TABLE]

where $H_{k}(t)$ is the matrix of fundamental solutions of ELE Eqns. (11),(12) defined as Kirsten ; Kirsten1 ,

[TABLE]

We will also make a particular choice of $H_{k}(t)$, namely that $H_{k}(0)=\textbf{I}_{2}$. For convenience, we have defined a variable,

[TABLE]

The $i\tilde{q}^{*}$ value may then be found by looking at the smallest roots of the characteristic polynomial $F(k)$. In Fig.  1 we plot $F(k)$ vs. $k$ for both the forward and reverse protocols for a time period $T=1$.

Notice the interesting symmetry of the characteristic polynomials namely that, $F(k)$ for the forward protocol is the mirror image of $F(k)$ of the corresponding reverse protocol. As we will show in Appendix D, this symmetry is a consequence of the FT Crooks ; Park ; linDiff . As a consequence of this form of $F(k)$, we have

[TABLE]

For $T=1$, we get $k_{F}^{*}=3.66=-k_{R}^{*}$ and therefore $P_{F/R}(W)$ has the leading-order asymptotic form ($\beta=1$),

[TABLE]

Now we move on to the computation of the pre-exponential factor. Obtaining the factor $d_{0}^{2}$ is straightforward, using Eq. (13) in Eq. (19) we see that,

[TABLE]

where $||x||$ is the norm of the zero-mode, and is defined as the usual inner product,

[TABLE]

The $*$ in the equation above stands for complex conjugation. The computation of the second factor, which is the ratio of the functional determinants, is rather involved, and may be evaluated using a technique that is developed based on the spectral $\zeta$ function of Sturm-Liouville type operators Kirsten . The power of the method is that it enables one to compute the determinant ratio in terms of the zero-mode itself. In terms of suitably normalized zero-mode solutions, the determinant ratio become

[TABLE]

The subscript $N$ denotes that a particular normalization is chosen. Explicit details of the calculation including the choice of normalization will be discussed in Appendix B and C.

Using Eq. (21), Eq. (40) and Eq. (42), the exact asymptotic form of $P(W)$ with the pre-exponential factor can now be determined for both forward and reverse protocols. Doing explicit computations for $T=1$, (see Appendix C) we get,

[TABLE]

The asymptotic forms obtained above are consistent with the Crooks fluctuation relation (Eq. (22)), with $\beta\Delta F=\frac{1}{2}\text{ ln }\frac{k_{f}}{k_{i}}=\frac{1}{2}\text{ ln }\frac{1}{2}$ for the reverse process speck ; Park . We have also calculated $P(W)$ numerically, and the results are in good agreement with the theoretical predictions.

2.3 Remarks

The main points in the calculation of the asymptotic form of $P(W)$ for the breathing parabola hold in general for PDFs with a class555where all the terms in the action $S_{W}$ are of degree $2$ in variables $x(t)$ and $\dot{x}(t)$ ( For example, the sliding parabola with deterministic driving Vanzone1 also has a quadratic action. But the action also contains terms of degree less than 2 in $x(t)$ and $\dot{x}(t)$. As a result of this, the the operator appearing in $S_{quad}$ is different from the one that appear in $S_{lin}$. In this case, EN theory gives the full work distribution which is a Gaussian for all values of $W$.). of quadratic augmented action $S_{W}$ (Eq. (6)), and may be summarised as follows.

• •

$P(W)$ has, to leading-order, the functional form,

[TABLE]

where $i\tilde{q}^{*}$ is the smallest root of some characteristic polynomial $F(k\equiv 1-iq)$. The Crook’s fluctuation theorem is seen to manifest itself as a reflection-symmetry property of this characteristic polynomial $F$ around $k=0$. This symmetry is no surprise when we realize that the exact moment generating function of work in this case is given by,

[TABLE]

The $i\tilde{q}^{*}$ value which determines the asymptotics correspond to the singularity of this moment generating function lying close to zero. We will derive Eq. (45) in Appendix D.

• •

Including the pre-exponential factor, $P(W)$ takes the form,

[TABLE]

where all the factors may be computed in terms of the fundamental solutions of the Euler Lagrange equations. This form for the tail is exact in the low noise limit $\beta\to\infty$, since there are no higher order expansion terms in (7) that we are neglecting.

In the following Sections, we will apply the techniques developed here to a stochastically driven system.

3 The stochastic sliding parabola

Consider the dynamics of a colloidal particle in a harmonic trap where the mean position of the trap is externally modulated Vanzone1 ; Vanzone2 ; engel ; Sabhapandit ; Varley ; sabha1 . Such potentials go by the name sliding parabola, and have the general form :

[TABLE]

where $x(t)$ is the position variable and $\lambda(t)$ is the externally modulated mean position. We have set the stiffness of the trap to $1$. The dynamics of the colloidal particle in this potential can be described by the Langevin equation,

[TABLE]

where $\eta(t)$ is a thermal noise and $D$ is the diffusion coefficient. $\eta$ is assumed to be Gaussian with $\langle\eta(t)\rangle=0$, and $\langle\eta(t)\;\eta(s)\rangle=\delta(t-s)$. One of the natural time scales in the system is given by the relaxation time in the harmonic trap ${{\tau}_{r}}=\gamma/\kappa$ where $\gamma$ is the friction coefficient and $\kappa$ is the stiffness of the trap. In the case of a deterministic driving protocol, the exact statistics of the Jarzynski work done on the colloidal particle is known Vanzone1 . The work distribution is a Gaussian and satisfies the transient fluctuation theorem. The sliding parabola with a deterministic driving was also looked at in Nickelsen2011 ; engel , as a test example for the EN Theory, and in this case the method gives the full probability distribution (not just the large-$W$ form). Eq. (48) has also been studied both experimentally expt and analytically Sabhapandit ; sabha1 ; Varley when $\lambda(t)$ is a stochastic driving protocol. One of the cases studied is when $\lambda(t)$ is the Ornstein-Uhlenbeck process given by,

[TABLE]

$\xi(t)$ is again assumed to be a Gaussian noise with $\langle\xi\rangle=0$ and $\langle\xi(t)\xi(s)\rangle=\delta(t-s)$. The noise $\xi$ is usually athermal in origin with a diffusion coefficient $A$ as given in Eq. (49). $\tau_{0}$ gives the second natural time scale in the system in terms of the relaxation time of $\lambda$ correlations. The two noises are assumed to not have cross correlations, i.e. $\langle\eta(t)\;\xi(s)\rangle=0$. Hereafter, we will refer to the coupled equations (48) and (49) as the Stochastic Sliding Parabola (SSP).

In the remaining Sections of this paper, we study the SSP model for both equilibrium and non-equilibrium steady state initial conditions using EN theory. For equilibrium initial conditions, the dissipation function that satisfies a fluctuation theorem of the SSP can be identified with the Jarzynski work666This result can be derived from the ratio of the net probabilities of the forward trajectory and the corresponding time-reversed trajectory as discussed in Chernyk .. The form of the dissipated work in the steady state has been obtained in Varley . For both situations, the exact form of work distributions at arbitrary times is not known. Here we show that the discussions in Section 2.3 can be applied to this system, and we can hence compute the exact asymptotic form of both the transient and steady state work distributions including the pre-exponential factor, at arbitrary times $T$. For steady state initial conditions, we compare our results with Varley , in the appropriate limits. Without loss of generality, for the calculations that follow, we set $A=D=k_{B}\mathcal{T}$ and $\tau_{0}=\tau_{r}=1$.

3.1 Equilibrium initial condition. : Transient fluctuations

In the first case that we look at, we compute the asymptotic form of the distribution of the Jarzynski work done on the colloidal particle ($W$) by the stochastic force Eq. (49) starting from an initial equilibrium distribution given by,

[TABLE]

The particle is assumed be in thermal equilibrium initially for a fixed value of $\lambda_{0}$ and the partition function $Z_{0}$ is computed accordingly. The Jarzynski work done on the colloidal particle along each trajectory is defined in the same way as before :

[TABLE]

In terms of the joint probability density functional of trajectories $\left\{x\left(\cdot\right),\lambda\left(\cdot\right)\right\}_{0}^{T}$, the probability density function of work can be written down as,

[TABLE]

with the action

[TABLE]

The normalization constant for this case is pathintegral ,

[TABLE]

In order to find the large $W$ asymptotic behaviour of $P(W)$, we will adopt the methods discussed in Section 2.1. Here that means, we need to identify the optimal choice of both $\tilde{x}(t)$ and $\tilde{\lambda}(t)$ that minimizes the action $S$ for a given value of $W$. We follow the same procedure as before and put $x(t)=\tilde{x}(t)+y(t)$, $\lambda(t)=\tilde{\lambda}(t)+z(t)$ and $q=\tilde{q}+r$ and expand $S$ to second order in $y(\cdot),z(\cdot)$ and $r$.

[TABLE]

Notice again that since $S$ is quadratic, the expansion will not contain terms of order greater than quadratic in $y(\cdot)$, $z(\cdot)$ and $r$.

3.1.1 Leading-order form of $P(W)$

As in the case of the Breathing Parabola problem, The leading-order form of $P(W)$ can be computed in terms of the optimal trajectory $(\tilde{x}(t),\tilde{\lambda}(t))$ as,

[TABLE]

where $(\tilde{x}(t),\tilde{\lambda}(t))$ solve the following Euler Lagrange equations,

[TABLE]

together with Robin-type boundary conditions,

[TABLE]

and the constraint equation,

[TABLE]

Notice that (as we have seen in section 2.1.1,) a non trivial solution to the ELEs exists only for some specific values of $iq\equiv i\tilde{q}$, and only the smallest $i\tilde{q}$ value $(\equiv i\tilde{q}^{*})$ is relevant. Based on the discussion in Section 2.3, we can conclude that, to leading order,

[TABLE]

The $i\tilde{q}^{*}$ value may be found by looking at the roots of the function,

[TABLE]

corresponding to this problem (see Appendix C). In Fig. 3 we plot $F(k)$ vs $k$ for $T=1$.

It may be verified that $F(k)$ is a symmetric function under the transformations $k\rightarrow\;-k$ for any value of $T$. This indicates that the work fluctuations satisfy the fluctuation theorem

[TABLE]

This form of the Crooks fluctuation theorem is a consequence of the fact that $\Delta F=0$ for this particular choice of initial conditions. The smallest roots of $F(k)$ are found to be $k^{*}=\pm\;2.300$. Solving Eq. (61) with Eq. (62) for $k=\pm\;2.3$, and then using the constraint equation (63), it can be verified that $k=+\;2.3$ corresponds to the positive tail of $P(W)$ and $k=-\;2.3$ corresponds to the negative tail of $P(W)$. Hence using Eq.(64), the leading-order asymptotic forms of $P(W)$ may be written down as (For $\beta=1$),

[TABLE]

Notice that the asymptotic forms are again consistent with the fluctuation theorem (Eq. (66)).

3.1.2 $P(W)$ including the pre-exponential factor

As in Section 2.1.2, the asymptotic estimate for $P(W)$ obtained above can be improved by also taking into account the contributions coming from $S_{quad}$ in Eq. (55). Since the action in Eq. (53) is quadratic, it can be shown that the fluctuations around the optimal trajectory are again governed by the same operator A (Eq. (61)) which determines the optimal trajectory. $P(W)$ including the pre-exponential factor therefore takes the form:

[TABLE]

The form of the Jacobian that is required to derive Eq. (68) is obtained in Appendix A. As in the previous case, $\det\textbf{A}^{\prime}$ in Eq. (68) is the determinant of the operator A omitting the zero-mode. $d_{0}$ is given by,

[TABLE]

where ||\left[\begin{array}[]{c}x\\ \lambda\end{array}\right]|| is the norm of the zero-mode. The factor 2 in Eq. (68) again accounts for the two equipotent saddle points, $(\pm x(t),\pm\lambda(t),i\tilde{q}^{*}_{\pm})$. Notice that one significant difference from the case of the breathing parabola is that the functional operators appearing in the determinant ratio in Eq. (68) are $2D$ functional operators. We hence generalise the method of functional determinants Kirsten that we used for the $1D$ Sturm-Liouville operator in the previous case, for this $2D$ case, to compute the ratio of the determinants appearing in Eq. (68). The explicit details are given in Appendices B and C. For the case $T=1$ one can hence compute,

[TABLE]

and therefore using Eq. (68), the improved estimate to the positive and negative tails are,

[TABLE]

In Table 2 we give the exact asymptotic forms of $P(W)$ for different values of $T$.

We have also numerically integrated the Langevin equations to get an estimate for $P(W)$. For a step size of $\Delta t=10^{-3}$, averaging over $10^{6}$ realizations gives a reasonably good agreement with our theoretical predictions. The results are plotted in Figure 4.

The results agree for large values of $W$, except for the very far tails where there are not sufficiently many sample points. In Figure 5 we present the numerical result, verifying the fluctuation theorem.

3.1.3 $\lambda_{0}$ arbitrary.

The case when we leave $\lambda_{0}$ unconstrained, (equivalent to sampling $\lambda_{0}$ from a Gaussian distribution with very large variance), is a special case in which one can obtain a closed asymptotic form of $P(W)$ as a function of $T$ - the time period of duration of the protocol. Following similar calculations as for the previous case (Appendix C), but with modified boundary conditions (62), we find for the positive tail,

[TABLE]

and for the negative tail,

[TABLE]

For $T\rightarrow\infty$ we find that the time-dependent part of the the pre-exponential factor converges to a value $\sqrt{\frac{7}{2}}$. Notice that the positive tail of the probability distribution decays as a power-law in $W$, and therefore the mean and higher moments do not exist. This behaviour of the tails can be attributed to the large fluctuations in the system which contribute to positive work values. However, for any finite variance of the initial Gaussian distribution of $\lambda_{0}$, the work distributions can be shown to have tails of the form (71), and well defined moments.

3.2 The steady state fluctuations.

In this section we study the steady state work fluctuations in the SSP. The functional form of the dissipated work ($W_{d}$) in the steady state of the SSP was identified in Varley . It was then shown to satisfy the FT by computing the corresponding moment generating function in the large $T$ limit. Here we look at the exact asymptotic form of $P(W_{d})$ using EN theory. We will then make a comparison with results from Varley in the appropriate limits.

3.2.1 Asymptotic form of $P(W_{d})$ for large $W_{d}$.

The dissipated work in the steady state in the time interval $[0,T]$ for the SSP is given by Varley ,

[TABLE]

The form of the dissipation function may also be identified from the formalism presented in Chernyk . Sampling the initial points from the stationary probability distribution (from Sabhapandit ),

[TABLE]

one can write down the probability distribution for dissipated work in the steady state as,

[TABLE]

with the augmented action

[TABLE]

Notice that the action in Eq. (77) is again quadratic. Therefore one can infer that in the asymptotic regime, $P(W_{d})$ must again be of the form (Section 2.3):

[TABLE]

It can be verified that the functional operator that determines the optimal trajectory, as well as the fluctuations around the optimal trajectory, is again A, as in the previous Sections 2.1.2 and 3.1.2. The difference comes in the boundary conditions; matrices $M$ and $N$ get modified accordingly.

As before, the leading-order behaviour of $P(W_{d})$ can be found from the smallest roots of the corresponding characteristic polynomial $F(k)$ ( Eq. (65) ). In Fig. 6 is a plot of the characteristic polynomial for $T=1$.

Notice that $F(k)$ is again symmetric about $k=0$, and this indicates that the dissipated work given by Eq. (74) indeed satisfies the fluctuation relation

[TABLE]

The smallest roots of $F(k)$, are found to be $k^{*}=\pm 3.16$ and therefore the leading-order asymptotic form is given by (for $\beta=1$):

[TABLE]

where the $\pm$ denotes the positive and negative tails respectively. As done earlier, one can also improve this estimate by including the pre-exponential factor. The steps involved are identical to the previous cases considered. The quadratic nature of the augmented action in Eq. (77) leads to a sub-leading power-law behaviour $\sim|W_{d}|^{-\frac{1}{2}}$ , and a numerical factor which is completely determined in terms of the zero-mode of the operator A, and depends only on the time duration $T$ of the protocol. Explicit calculations may be carried out in the same way as in the previous case (see Appendix C). For $T=1$ we find:

[TABLE]

In Table 3 we present the exact asymptotic forms of $P(W_{d})$ for different values of $T$

3.2.2 Comparison with results in Varley .

Verley et al Varley , obtain the exact form of the generating function of the probability distribution of the dissipation function for large values of the time duration $T$ of the protocol, and show that the Crooks fluctuation theorem is satisfied in this limit Varley . In order to compare the two methods, we have first inverted the generating function from Varley (details in Appendix E). The exact form of $P(W_{d},T)$ for $T\gg 1$ is obtained as :

[TABLE]

In order to compare this result with our calculations, we do an asymptotic expansion of Eq. (82) for large $W_{d}$. To leading-order we find that,

[TABLE]

We compare this with the leading-order form in Eq. (78) by checking how $k^{*}(T)$ behaves as $T$ becomes large. Knowing the exact form $F(k)$ for any value of $T$, this behaviour may be readily found. As we show in Fig. 7, we find,

[TABLE]

and this leads to the asymptotic forms given in Eq. (83). Therefore in the large-$T$ limit, the leading-order behaviour predicted by both methods agree.

Now we look at the sub-leading pre-exponential behaviour predicted by both methods. The asymptotic expansion of the pre-exponential factor of $P(W_{d},T)$ in Eq. (82) gives a sub-leading pre-exponential behaviour $\sim|W_{d}|^{-5/2}$, which is a much faster decay than predicted by Eq. (78) which suggests a sub-leading behaviour $\sim|W_{d}|^{-1/2}$ for any value of $T$. In Figure 8, we compare the results from both methods with our simulations of the Langevin dynamics.

The asymptotic forms obtained in Table 3 are good fits to the tails for all values of $T$. As $T$ gets larger, Eq. (82) becomes a good fit to the numerical data, and the two methods agree to a good extent at the tails (as expected from the same leading-order behaviour). In order to see the difference coming from the disagreement in the pre-exponential factors of both methods, extensive Langevin simulations will be needed.

4 Conclusion.

In this paper we have determined the exact asymptotic form of the work distribution, including the pre-exponential factor, in a class of stochastically driven systems, using a theory developed by Engel and Nickelsen (EN theory) Nickelsen2011 . In cases where the exact finite time work distribution is not known, EN theory can be applied to obtain good analytical approximations for the tails of the distributions, which are generically hard to observe in experiments or numerical simulations. This can then be combined with data from experimentally / numerically viable regimes to construct the full probability distribution. The extension of EN theory to stochastically driven systems has not been done previously and the analytic solutions that we obtain here are new.

EN theory involves writing an augmented action which carries all the information about initial conditions as well as the functional - $W[x(\cdot)]$ whose probability distribution is to be computed. The asymptotic form of work distributions are then computed by using the saddle point approximation, which formally corresponds to the small noise limit. For a class of quadratic (augmented) actions, we have shown that the smallest roots of a certain characteristic polynomial function $F(k)$ determines both the leading-order asymptotic form as well as the pre-exponential factor. The Crooks fluctuation theorem is then shown to manifest itself as the reflection symmetry property of this function. These features can be explained using the relation (proved in Appendix D),

[TABLE]

for the exact moment generating function (MGF). For a colloidal particle in a harmonic trap, where the mean position of the trap is modulated according to the Ornstein-Uhlenbeck process, we have shown that the (dissipated) work distributions have the asymptotic behaviour,

[TABLE]

in both the transient and the steady state. Here the constants, $C_{1}$ and $C_{2}$ are fixed by the time duration of the driving and the noise coefficients, and are explicitly determined for all cases. The asymptotic form given by Eq. (86) can be shown to be universal for quadratic augmented actions. A rigorous discussion of this point and also the relation of the asymptotic form to the singularities of the MGF may be found in linDiff ; Park . In linDiff , the authors have shown that for the same class of systems for which the above discussions apply, determination of various PDFs and MGFs reduce to finding solutions of certain nonlinear differential equations (NLDEs), which in many cases need to be solved numerically. The method of functional determinants simplifies this problem to instead determining the solutions of ELEs, which are linear ordinary differential equations. In loghar , for the special case of a Brownian particle in a logarithmic harmonic potential, the authors obtained the asymptotic form of work distribution in terms of the solution of a Riccati differential equation. We are not aware of any other analytic methods for performing the exact finite-time computation of the asymptotic form including the pre-factor, in Langevin systems. Although we have restricted ourselves to the computation of asymptotics of work distributions in this work, Eq. (85) contains more information, such as the complete moment hierarchy. These aspects will be discussed in a future publication.

We believe that the methods that we discuss here has potential applications in the context of finite time thermodynamics of stochastic systems. It should be interesting to look into more applications of this theory, particularly in other stochastic potentials typical for experimental situations. For example, the theory can be developed to include stochastic driving governed by discrete time stochastic processes such as the one studied in unc , by using an appropriate path integral representation Doi ; Doii ; peliti . However, identifying universality classes for asymptotics of work distributions still remains an open problem.

Acknowledgement

We would like to thank Daniel Nickelsen for very helpful discussions, comments and a critical reading of an earlier version of this manuscript. We would also like to thank Viktor Holubec and Dominik Lips for pointing out an error in the references, and Gatien Verley for helpful comments on reference Varley .

Author contribution

Both the authors have contributed equally to this manuscript.

Appendix A The Jacobian

In this Appendix, we obtain the exact form of the Jacobian of the transformations that is required in deriving Eq. (68) in Section 3.1. We generalise the derivation of Engel and Nickelsen in Nickelsen2011 , to the SSP studied in Section 3 described by the Langevin equations,

[TABLE]

The propagator of the corresponding Fokker-Planck equation may be written down as,

[TABLE]

From the normalization condition

[TABLE]

where,

[TABLE]

we have,

[TABLE]

After several partial integrations in the RHS of Eq. (91) we obtain,

[TABLE]

The above Gaussian integral is straightforward to compute and gives,

[TABLE]

where J is the Jacobian to be determined. Note that $\textbf{A}=\textbf{A}_{iq=0}$ as given by Eqns. (61) and (62). Therefore we have,

[TABLE]

This is the result used in deriving Eq. (68). In a similar manner one can show that all the other Jacobians appearing in the main text, have the above form.

Appendix B Functional determinants

The necessity of computing functional determinants of certain differential operators arises in many different situations. As we have seen in the main text, computing the leading-order contribution to path integrals is one of them. In many cases it is not an absolute functional determinant that is required, but a ratio where one or both of the operators can in principle have a zero-mode. Profound mathematical techniques have been developed to compute functional determinants (or the ratio of determinants) even in situations when the operators have zero-modes. Here we apply the contour integration method suggested in Kirsten to the SSP problem discussed in Section 3. As we have already seen, the operators that appear in the ratio of determinants in Eq. (68) are $2\times 2$ matrix differential operators. In Kirsten Kirsten et al, have discussed the possible generalization of their techniques to such matrix differential operators as well. Recently in functional , Falco et al have also looked at a similar generalization. In contrast to these previous studies, the operator A that we study here has differential operator entries in the off diagonal terms as well. However, the methods discussed in Kirsten can also be generalized to this situation. In the case of the SSP, the matrix differential operator that we need to find the determinant of, is defined by the following problem:

[TABLE]

together with Robin-type boundary conditions:

[TABLE]

The form of the matrices $M$ and $N$ can be deduced from Eq. (62). Using the results from Kirsten , one can then write down the determinant ratio as,

[TABLE]

Here $H_{k}$ is the matrix of fundamental solutions of the homogeneous equation,

[TABLE]

defined as,

[TABLE]

The function B appear due to the presence of the zero-mode and need to be determined using the self adjointness property of the differential operator A in each case. The normalized zero-mode, $\textbf{u}_{N}(t)$ is defined as,

[TABLE]

where the constants are determined by,

[TABLE]

The inner product is the usual one, given by

[TABLE]

In case of the SSP discussed in Section 3.1, we find,

[TABLE]

In the next Section, we show how this theory can be used to compute the functional determinants that appear in the main text.

Appendix C Explicit Computations

Here we provide the explicit calculations using EN theory for two of the cases considered in the main text, the breathing parabola in Section 2.1 and the SSP in Section 3.1. For the breathing parabola problem, in Nickelsen2011 , for a specific choice of (reverse) protocol, EN theory was used to compute the exact asymptotic form of $P(W)$. Here we present the calculations for a particular choice of forward protocol, and give only the final solution for the corresponding reverse protocol. The solutions we obtain for the SSP problem in Section C.2 are new, and are generalizations of the calculations in Section C.1.

C.1 Breathing Parabola: $P_{F/R}(W)$

For the breathing parabola, we have considered the specific forward protocol,

[TABLE]

In terms of the shifted variable $k=1-iq$, the corresponding ELE reads (for simplicity we will use $x(t)$ instead of $\tilde{x}(t)$ everywhere.),

[TABLE]

Two independent solutions of this 2nd order differential equations are given by,

[TABLE]

Together with the boundary conditions,

[TABLE]

where M and N are matrices,

[TABLE]

The above system constitutes a second order Sturm-Liouville eigenvalue problem in $k$. A non trivial solution exists only for some specific values of $k$, which are given by the roots of the corresponding characteristic polynomial,

[TABLE]

As we have seen in the main text, the asymptotic behaviour of $P(W)$ is determined by the smallest value of $k(\equiv k^{*})$ for which $F(k)=0$. The value of $k^{*}$ may be obtained numerically. For the case $T=1$ we find $k^{*}=3.67$ (see Figure 1 in the main text). The leading-order asymptotic form of $P(W)$ is therefore,

[TABLE]

In order to improve this estimate one has to compute the pre-exponential factor as well. The computation of the factor $d_{0}^{2}$ is rather straightforward. From (19) we see that,

[TABLE]

where $x(t)$ is the zero-mode. For $k=-3.67$, we find that the zero-mode is

[TABLE]

Here $C_{1}$ is the undetermined constant in the solution, which is to be fixed using the constraint equation (13). Using the above form of $x(t)$ and the constraint equation (13), we find,

[TABLE]

The other factor which appears in the calculation of the pre-exponential factor is the determinant ratio of the two functional differential operators,

[TABLE]

Let us first compute the factor appearing in the numerator of the RHS. Using Eq. (124), the matrix of normalized fundamental solutions ($H(0)=\textbf{I}_{2}$) when $iq=0$ may be found as,

[TABLE]

Therefore,

[TABLE]

Notice that, this is also the limiting value of $F(k)$ as $k\rightarrow 1$ in Figure 1. Let us now look at the term in the denominator of the RHS of Eq. (140). The appropriately normalized solutions $x^{N}(t)$ can be identified using (124) and Eq. (119) adapted to this problem. We obtain,

[TABLE]

Using the methods discussed in Kirsten we find that for this problem, $\textbf{B}=-\frac{1}{x_{N}(T)}$. Together with this, one finds

[TABLE]

Putting the factors together in Eq. (21), we finally get,

[TABLE]

Similarly for the reverse protocol,

[TABLE]

it can be shown that,

[TABLE]

We compare these results with numerical simulations, and as we show in the main text, Figure 2, the prediction for the tail region is in excellent agreement with the theoretical predictions.

C.2 The Stochastic Sliding Parabola

In this Appendix, we do explicit computations to obtain the asymptotic form of the probability distribution Eq. (71), found in Section 3.1 for equilibrium initial conditions and $T=1$ (explicit calculations for the other cases discussed in section 3.1.3 and section 3.2 can be carried out in a similar manner). For simplicity we will use the notation $x(t)$ and $\lambda(t)$ instead of $\tilde{x}(t)$ and $\tilde{\lambda}(t)$.

First we note that when $|k|\neq 1$, following some algebra, the system of Euler-Lagrange equations for $\left(x,\lambda\right)$ given by Eq. (61) in the main text can be reduced to a fourth order ordinary differential equation for one of the variables (for example, $\lambda$) as,

[TABLE]

In terms of the solution $\lambda(t)$, $x(t)$ is then given by,

[TABLE]

Eq. (148) has four independent solutions given by,

[TABLE]

A general solution for a specific optimal trajectory ($k=k^{*}$) can always be written as a linear combination of these four independent solutions, where the coefficients are fixed by the boundary conditions and the constraint equation. As we discussed in Section 3.1.1, in order to compute the leading-order behaviour of $P(W)$, we only require the $k^{*}$ values and not the explicit solution. In order to find $k^{*}$, we look at the smallest roots of the function:

[TABLE]

$M$ and $N$ corresponding to the boundary conditions in (62) can be written down as,

[TABLE]

$H_{k}(t)$ has the form as given in Eq. (117). From Figure 3 in the main text, we see that the relevant $k^{*}$ values are given by $k^{*}=\pm 2.3$. Solving the ELEs (61) along with the boundary conditions (62), for $k=2.3$ yields,

[TABLE]

If we evaluate the work done along this optimal trajectory, we get,

[TABLE]

The work done is positive; this indicates that this pair of trajectories correspond to the positive tail of P(W). Similarly, solving the ELEs (61),(62) for $k=-2.3$ gives,

[TABLE]

The work done along this trajectory become,

[TABLE]

This value is negative, and therefore $k^{*}=-2.3$ corresponds to the negative tail of $P(W)$. With this we find that to leading-order, the positive and negative tails of P(W) have the functional form,

[TABLE]

respectively.

In order to improve this estimate, we next include the pre-exponential factor. The first factor which goes into the pre-exponential is $d_{0}^{2}$ defined in Eq. (69). This can be computed relatively easily as in the case of the breathing parabola. For both $k=\pm\;2.3$ we find using Eq. (161) and (164),

[TABLE]

(The superscript $\pm$ is used to denote the solutions for positive or negative tails). The next factor to be computed is the square root of the ratio of two functional determinants, for which we will use the formula,

[TABLE]

where,

[TABLE]

is the appropriately normalized solution to the ELEs, which have to be found using Eq. (119). For this particular case, We find that,

[TABLE]

Similarly

[TABLE]

Using the self-adjointness property of A, one can again compute,

[TABLE]

Also using Figure 3 to compute $F(1)$, we find,

[TABLE]

Therefore the full pre-exponential factor is,

[TABLE]

and the asymptotic form of $P(W)$, including the pre-exponential factor becomes,

[TABLE]

In a similar manner, the exact asymptotic forms discussed in Section 3.1.3 and Section 3.2 can be computed for any value of $T$. In particular, for the case that we discussed in Section 3.1.3, the exact asymptotic form can be obtained as a function of $T$.

Appendix D Origin of the symmetry of $F(k)$

In this section we will show the relation between $F(k)$ and the exact moment generating function of dissipated work, which explains the reflection symmetry of $F(k)$. First, using the path integral representation, an exact relation for the moment generating function can be written down as,

[TABLE]

In all the cases we have considered, the augmented action $S[\;x,\;\lambda,\;q\;]$ is quadratic, therefore by doing several partial integrations, it can be shown that it reduces to

[TABLE]

where the kernel $\textbf{A}_{k}$ is defined by the same operator that determines the optimal trajectory (Eq. (11), (61)) with the same boundary terms. Therefore the integral in Eq. (177) is a standard Gaussian integral which can be computed as,

[TABLE]

This determinant ratio can then be computed using the techniques developed in Kirsten . In terms of the function $F(k)$, we find,

[TABLE]

Due to Crooks fluctuation theorem Crooks , the moment generating function of dissipated work ($G$) must satisfy the relation,

[TABLE]

By writing $iq$ as $1-k$ and using Eq. (183), it can be immediately verified that for the above relation to hold, the function $F$ must satisfy,

[TABLE]

Hence the symmetry property of $F(k)$ is a consequence of the fluctuation theorem.

Appendix E Comparison with the results in Varley

In this Appendix, we will invert the generating function of $P(W_{d})$ obtained in Varley using the methods discussed in Sabhapandit and sabha1 . The solution that we obtain here will be used for comparison with the asymptotic form of $P(W_{d})$ calculated using the EN theory in Section 3.2.

Verley et al, in Varley , have shown that for very large $T$, the probability generating function of the dissipated work has the form,

[TABLE]

The notation $W_{T}$ stands for the dissipated work $W_{d}$ over a time duration $T$ of the driving. For the SSP considered in Section3.2, the functions $\phi$ and $g$ are given by Varley ,

[TABLE]

The probability density function can be obtained from the moment generating function $Z(\mu,t)$ by taking the inverse Fourier (two-sided Laplace) transform:

[TABLE]

where the integration is done along the imaginary axis in the complex-$\mu$ plane sabha1 . Using the large T form of $Z(\mu,T)$ given by Eq. (186) and (187) we write,

[TABLE]

where

[TABLE]

The large-$T$ form of $P(W_{T})$ can be obtained from Eq. (189) by using the method of steepest descent (for completeness, we reproduce the method and discussion from sabha1 here ). The saddle point $\mu^{*}$ is obtained from the solution of the condition $f_{w}^{\prime}(\mu^{*})=0$ as

[TABLE]

From the above expression, one finds that $\mu^{*}(w\rightarrow^{+}_{-}\infty)\rightarrow\mu_{\pm}$, where

[TABLE]

Therefore $\mu^{*}\in\left(\mu_{-},\mu_{+}\right)$. It is useful to notice that in terms of $\mu_{\pm}$,

[TABLE]

On the real axis, outside the interval $\left[\mu_{-},\mu_{+}\right]$, $\nu(\mu)$ is therefore imaginary. However, in order for the the integral in the definition of $Z$ (Eq. (186)) to converge, $Z(\mu,T)$ must be real for real values of $\mu$. For this reason, it is only within the range $\mu_{-}<\mu<\mu_{+}$ ( for which $\nu(\mu)$ is real and analytic), that the analytic continuation of $Z(\mu,T)$ to real $\mu$ is allowed. We hence expect the saddle point to also lie between these values. As we have already seen in Eq. (192), this is indeed the case. Since for $\mu_{-}<\mu<\mu_{+}$, $g(\mu)$ is analytic (the denominator is positive for all $\mu$ in the specified range), it can be neglected in the saddle point calculation as a sub-leading contribution.

The saddle point calculation relates $\phi(\mu)$ to the large deviation function $h_{s}(w)$ by the Legendre transform,

[TABLE]

We also see that,

[TABLE]

This means that $f_{w}(\mu)$ has a minimum at $\mu^{*}$ along real $\mu$. Now since $g(\mu)$ is analytic, the usual saddle point approximation method Sabhapandit gives,

[TABLE]

Using, Eq. (194), (191), (190) and also the relation $w=W_{d}/T$ in Eq. (196), we finally get

[TABLE]

We use this form of $P(W_{d})$ in Section 3.2.2 to compare with numerical results as well as the analytic forms obtained using the EN theory.