Suppression of work fluctuations by optimal control: An approach based on Jarzynski’s equality

Suppression of work fluctuations by optimal control: An approach based on Jarzynski’s equality

Gaoyang Xiao and Jiangbin Gong Department of Physics, National University of Singapore, Singapore 117542
Department of Physics and Centre for Computational Science and Engineering, National University of Singapore, Singapore 117542
July 15, 2019

Understanding and manipulating work fluctuations in microscale and nanoscale systems are of both fundamental and practical interest. For example, aspects of work fluctuations will be an important factor in designing nanoscale heat engines. In this work, an optimal control approach directly exploiting Jarzynski’s equality is proposed to effectively suppress the fluctuations in the work statistics, for systems (initially at thermal equilibrium) subject to a work protocol but isolated from a bath during the protocol. The control strategy is to minimize the deviations of individual values of from their ensemble average given by , where is the work, is the inverse temperature, and is the free energy difference between two equilibrium states. It is further shown that even when the system Hamiltonian is not fully known, it is still possible to suppress work fluctuations through a feedback loop, by refining the control target function on the fly through Jarzynski’s equality itself. Numerical experiments are based on linear and nonlinear parametric oscillators. Optimal control results for linear parametric oscillators are also benchmarked with early results based on accelerated adiabatic processes.

I Introduction

Aspects of thermodynamics at the nanoscale have attracted great attention in recent years. In a small system with few number of particles, thermal fluctuations and quantum fluctuations in work, heat, and other quantities can be comparable to their ensemble mean values. Understanding and manipulating these fluctuations hence become an interesting and important topic. From a fundamental point of view, the connections between nonequilibrium statistical fluctuations with equilibrium properties, such as those established through fluctuation theorems (e.g., Jarzynski equality Jarzynski.97.PRL (); Jarzynski.97.PRE () and crooks theorem crooks ()), have laid a solid and fruitful foundation for nanoscale thermodynamics Campisi.11.RMP (). It is interesting to note that the existence of such type of fluctuation theorems does not rule out the possibility to control statistical fluctuations. That is, under appropriate conditions, it is still possible to manipulate nonequilibrium processes and alter the statistical fluctuations in thermodynamic quantities, without violating any fundamental law in thermodynamics or modifying the fluctuation theorems themselves. In the case of Jarzynski’s equality, i.e., Jarzynski.97.PRL (); Jarzynski.97.PRE (), where is work done on a system initially at thermal equilibrium with inverse temperature and is the free energy difference between the final and initial equilibrium states of the same , a change in fluctuations of can have an impact on how fast the statistics of may converge towards . This is the first motivation of this study.

Advances in our knowledge of thermodynamics of small systems have stimulated studies of efficient energy devices (classical or quantum) at microscale and nanoscale, sometimes consisting of few particles or few degrees of freedom. Of particular interest here is that several designs of microscale or nanoscale heat engines have been proposed theoretically Abah.12.PRL (); Bergenfeldt.14.PRL (); Zhang.14.PRL (). For this type of applications where heat-to-work conversion efficiency and the power of work output are apparently crucial, the stability or reliability of work output also becomes an important performance indicator Deng.13.PRE (); dario (). This is because the cycle-to-cycle fluctuations in the work output are an unavoidable characteristic of nanoscale heat engines. Given two heat engines with the same mean work output per cycle, it seems desirable to prefer the one with less fluctuations and hence more uniform output. This understanding further motivates us to ask the following question: how to systematically suppress work fluctuations in a given protocol (such as in one step of a heat engine cycle that does not involve any heat exchange with a bath)?

It has been realized that work fluctuations in some controlled processes can be relatively smaller than those in bare uncontrolled processes. In particular, work fluctuations in an adiabatic process (here “adiabatic” means very slow, as compared with the system’s own natural time scale) are relatively small Deng.13.PRE (). However, an adiabatic process in the standard sense is simply too slow, and as such the power of work output based on adiabatic processes vanishes. For this reason it is natural to consider accelerated adiabatic processes to suppress work fluctuations and increase the power of work output at the same time Deng.13.PRE (); campoarxiv (). Certainly, accelerated adiabatic processes just represent a special type of controlled dynamics, and extensions based on alternative constructions of a control field should be possible. This paper is precisely one extension of Ref. Deng.13.PRE () in order to suppress work fluctuations in a wider variety of systems. The extension is based on an optimal control framework Shi.90.JCP (); Peirce.88.PRA (); Shi.91.CPC (); Sundermann.99.JCP (); Ohtsuki.04.JCP () that is far more versatile and far more flexible. We shall also use known results based on accelerated adiabatic processes to benchmark our optimal control results.

To see the necessity and advantages of an optimal-control-based extension, let us first make brief comments on early studies of accelerated adiabatic processes Demirplak.08.JCP (); Demirplak.05.JPCB (); Demirplak.03.JPCA (); Berry.09.JPAMT (); Masuda.10.PRSLSA (); Chen.10.PRL (); Chen.10.PRLa (); Choi.11.PRA (); Schaff.11.NJP (); Bason.12.NP (); Zhang.13.PRL (); Jarzynski.13.PRA (); Campo.11.PRA (). For the sake of discussions later on, we divide accelerated adiabatic processes known to date into two types (though these two types can even be regarded as being equivalent upon a transformation  Deffner.14.PRX () and there is no clear distinction in the literature). In the first type, an additional control Hamiltonian is introduced to drive a system (within a short time scale), such that the evolution of the system, either classical or quantum, still follows the adiabatic evolution of the original bare system. Here we call this type of control as fast-forward adiabatic driving (FFAD) Masuda.10.PRSLSA (). The extra control Hamiltonian in realizing FFAD may be found for very simple systems but in general, its analytical form is not available and and it is also challenging to precisely implement it experimentally. The second type of accelerated adiabatic processes, which we call “shortcuts to adiabaticity” (STA) mugareview (), relies on special time dependence of system parameters. In this sense, the required form of the control Hamiltonian is already contained in the original bare Hamiltonian. Within a finite time, the system under STA control can reach the same final states as those reached in conventional adiabatic processes. However, during an STA process the system may make nonadiabatic transitions Tu.14.PRE (); campoarxiv (). The realization of STA has been established in several simple systems such as two-level systems or linear parametric oscillators. As Ref. Deffner.14.PRX () showed, STA in known examples can be understood as a result of a canonical transformation of FFAD (or vice versa). Though Campo and Jarzynski extended STA to more complicated systems such as many-body systems Campo.13.PRL (), a general discussion by them also indicates that this method can only apply to systems with the so-called “scale-invariant driving” Deffner.14.PRX (). For these reasons both FFAD and STA require full knowledge of the system Hamiltonian, and this requirement presents a limitation when we attempt to use accelerated adiabatic processes to suppress work fluctuations in general situations.

In this work we restrict ourselves to classical systems (but our optimal control approach can be easily extended to quantum systems). As such only work fluctuations across a thermal ensemble are to be suppressed. Our classical results should be also useful in guiding possible manipulations of work fluctuations in quantum systems, especially in cases with a relatively high temperature. Specifically, we develop a method to design control fields to suppress work fluctuations based on the well-known optimal control theory (OCT). The peculiarities of our OCT here lie in two aspects. First, we need to handle a thermal ensemble in order to consider work fluctuations. Thus our control is ensemble-based optimal control. Second, we need to design a useful control target function in order to reach our goal. As shown below, our control target function directly exploits Jarzynski’s equality, and the constructed target function is biased against deviations of individual values of from their ensemble average, i.e., . We show that our optimal control method built in this manner can be applied to a wide variety of systems, including highly nonlinear systems. When it is applied to previously known simple systems, it has essentially the same quantitative performance as that based on accelerated adiabatic processes, but now the control fields have many possible forms. Furthermore, our OCT method can be extended to implement a feedback mechanism, so as to suppress work fluctuations in systems with unknown system parameters. This is done by guessing the control target function iteratively through Jarzynski’s equality again.

This paper is arranged as follow. In Sec. II, we apply OCT to the topic of suppression of classical work fluctuations, with necessary details. In Sec. III, we test our approach using a simple system and compare our results with those obtained based on accelerated adiabatic processes. In Sec. IV, we consider a nonlinear oscillator and show that our OCT approach can still operate well. In Sec. V, we consider a nonlinear oscillator with some unknown system parameters. We show that, in the spirit of feedback optimal control, it is possible to refine the control target function on the fly, and then in the end it is still possible to find a control field to suppress the work fluctuations. The last section of this paper gives a brief summary. An appendix regarding numerical implementation of our OCT is also given.

Ii Theory of Optimal Control of work fluctuations: use of Jarzynski’s equality to construct a control target function

Let us consider a general time-dependent classical system with Hamiltonian , where represents phase space coordinates (we assume one-dimensional problems), and is a time-dependent parameter due to a work protocol starting from to . For example, this can represent the frequency of a parametric oscillator, the length of a pendulum, or the width of a box, etc. Such a classical system is initially prepared at thermal equilibrium, with the phase space probability density given by


where is the partition function associated with . To suppress work fluctuations in this protocol defined by , an extra control field is applied. Then the total Hamiltonian of the system becomes


where is the control Hamiltonian and is assumed to be the time-dependent amplitude of a control field. The evolution of obeys the following Hamilton’s equation of motion: and .

Next we consider a certain physical quantity evaluated at , which is written as . Consistent with this notation and for convenience, we also define and from now on. The minimization of the ensemble-averaged value of will be important in our OCT. The thermal ensemble average of is given by


where represents the entire phase space. We call defined above a target function since it will be a quantity we hope to minimize. For a control problem, typically a cost function is also needed to reflect a cost-related constraint. Here we define the cost function using the amplitude of the control field, i.e.,


where is a weightage factor that may depend on time. The larger is, the more constraints posted on the control field due to cost considerations. The overall target function can then be defined as . That is, the problem becomes to minimize under the general dynamical constraints reflected by Hamilton’s equations of motion. To that end we introduce two Lagrange multiplier as functions of , denoted by and . We then minimize instead, with


Let be the variation in due to , an arbitrary variation in , . Then the variation in due to is found to be


To minimize we let . Since the variation is arbitrary, one has the following relations:


The list of relations in Eq. (7) can be numerically solved by an iteration procedure and some details are presented in Appendix.

After outlining the general steps in OCT, we turn to our main objective, which is to minimize the statistical work fluctuations if a work protocol is applied to a thermal ensemble. For the above-defined protocol of during which there is no interaction with a heat bath, the inclusive work as a function of is given by Jarzynski.97.PRL (); Jarzynski.07.CRP (); Campisi.11.RMP ():


where is simply the phase space coordinates at the end of the protocol. The ensemble average of is


In order to suppress work fluctuations across the thermal ensemble, one may now choose an appropriate target function to minimize via OCT. A naive and simple choice of the target function can be just the variance squared of . That is, one may choose the ensemble average as the target function. However, it turns out that, different control fields can yield different values of , and as a result one does not know beforehand. Thus such type of intuitive choice of the control target function is not adopted here. At the point, Jarzynski’s equality comes as a rescue. Note first that


where refers to the free energy difference associated with the total Hamiltonian . So long as the initial and final values of are zero (this is easily realized in the next section) , then will be the same as , the free energy difference associated with and . This being the case, no matter what the to-be-found control field is, the ensemble-averaged value of is fixed, which is given by . We are thus motivated to design a control target function as follows:


The above-defined form of the control target function is intriguing, because minimization of this function is then to directly suppress the derivations of possible individual values of from their ensemble-average value predicted by Jarzynski’s equality. In this regard, our OCT framework directly exploits Jarzynski’s equality. It is also expected that the found optimal control field can remove those rare trajectories with rare values of , which slow down the convergence of simulation results towards Jarzynski’s equality.

Iii Optimal Control of Work Fluctuations in Linear Parametric Oscillators

As a benchmark step, in this section we consider parametric (linear) oscillator systems. In particular, for such systems, the control Hamiltonian to realize accelerated adiabatic processes can be found easily, for both scenarios of fast-forward adiabatic driving (FFAD) and shortcuts to adiabaticity (STA). We can hence compare the performance of our optimal control fields with those found in FFAD and STA.

The Hamiltonian of a parametric oscillator with a time-dependent frequency, all in dimensionless units, is described by


For FFAD, the extra control field is Deng.13.PRE (); Jarzynski.13.PRA ()


For STA, the control field is found to be Campo.13.PRL ()


As an example, we choose a frequency protocol Campo.13.PRL ()


such that the above two control fields are indeed zero at or . Besides, note that


So during the entire protocol if we choose . In our calculations we choose and in dimensionless units. Then, it is easy to show that


That is, along an arbitrary trajectory, the work value is always positive. We further set to be as small as (as compared with ), such that the process will be highly nonadiabatic were there no control fields.

To benchmark our OCT approach outlined above with FFAD and STA, we introduce two kinds of control fields, i.e., and (with to be found numerically), in parallel with the above two control fields and to realize FFAD and STA. The weightage factor in Eq. (4) is set to be small, since here we are not much concerned with the cost of the control field. Nevertheless, the time-dependence of can be designed to further alter the profile of the field amplitude. Here we propose to use , where is a small constant. Such time-dependence of makes diverge at and , and as a consequence the numerically found will become zero automatically at these two boundary times due to the cost function in the OCT.

With the form of the control field and the weightage factor in OCT both specified, numerical iterations based on Eq. (7) can then yield explicit solutions of . The results are shown in Fig. 1 for both -type and -type optimal control. Their respective time dependence is also compared with the corresponding results to realize FFAD and STA. Note that in all these cases the control field amplitude is zero in the beginning or at the end. It is also seen that found from our optimal control is different from those in FFAD or STA, though the difference for the shown computational example is not yet dramatic.

Figure 1: (color online) Time dependence of optimal control field for a parametric linear oscillator whose frequency time dependence is given in Eq. (15), with , , and . All the plotted quantities here and in other figures are in scaled and hence dimensionless units. The inverse temperature is set to be . (a) The control field amplitude of -type optimal control as compared with that based on fast-forward adiabatic driving (FFAD). (b) The control field amplitude of -type optimal control as compared with that in the shortcuts to adiabaticity (STA) approach. Note that the field amplitudes in (b) are much higher than those in (a).

Next we attempt to quantitatively characterize the performance in suppressing the work fluctuations, for our OCT approach along with FFAD and STA. In particular, we randomly sample initial phase space points (through a standard importance sampling prcedure) according to the initial thermal probability distribution and then evolve them under the total Hamiltonian . Individual values of are denoted , and the fluctuations in and in are characterized by


where the total number of simulation trajectories is chosen to be . The variance in both and in itself under different control schemes are all presented in Table 1. A few interesting observations are in order. First, within expected statistical error due to a finite , the bare system, FFAD, STA, OCT with -type field, and OCT with -type field all yield the same . Second, though the found time dependence from OCT is different from FFAD or STA, the variances in and in obtained from OCT with -type field are all the same as those obtained in FFAD and STA. This confirms that our OCT framework is doing an excellent job in suppressing work fluctuations, to the same degree as accelerated adiabatic processes can reach. This also hints that in the parametric oscillator example here, accelerated adiabatic processes already reach an optimized suppression of work fluctuations. Third, even the mean work from -type OCT agrees with those obtained from FFAD and STA. Fourth, the two variances shown in the last row of Table 1, which is for OCT with -type control, are however slightly above those in other cases. This is because for -type control, as also shown in Fig. 1, the required amplitude of the control field is very high, so a small weightage factor chosen for the cost function can still cause a minor difference. This is also manifested in the last row of Table 1, which is again relatively higher than those obtained in OCT with -type field, in FFAD, or in STA.


bare system
optimal control of -type
optimal control of of -type
Table 1: The performance of suppressing work fluctuations in the absence or presence of several different control fields, mainly characterized by the variance in work and the variance in , using trajectories. The system is a parametric linear oscillator whose frequency time dependence is given in Eq. (15, with , , and (duration of the protocol). The inverse temperature is set to be . Note that the obtained values of are all around the theoretical value theoretically obtained from Jarzyski’s equality.

Having successfully benchmarked our OCT approach, we now shed light on the flexibility of our optimal control approach. First of all, we can introduce different time dependence to the weightage factor when accounting for the cost of the control field so as to obtain different , the time dependence of a control field. For example, instead of , we have considered


for our -type optimal control. The obtained is presented in Fig. 2 in comparison with our previous result. Interestingly, although the time-dependence of the control field varies significantly with changes in the cost function weightage factor , all the three cases shown in Fig. 2 yield the same variance of , i.e., . This is a clear demonstration that there are many possible solutions in suppressing work fluctuations to a certain level.

So far we have chosen , which is more or less to simulate an instantaneous limit for the bare system. That is, the frequency of the parameter oscillator is changed rapidly as compared with the system’s own time scale. As shown above, in such a parameter regime our optimal control can suppress work fluctuations to the same degree as that achieved in accelerated adiabatic processes. To further check the usefulness of OCT, we now consider a slower protocol in which and . In this case, during the protocol the system’s own bare Hamiltonian will be important in the time evolution. Our optimal control fields are found to perform also very well in this regime. In particular, without a control field decreases to . This is because nonadiabatic effects and hence work fluctuations are weaker for a slower protocol. Interestingly, with FFAD or optimal control, we still find . One example of found from our -type optimal control is presented in Fig. 3 in comparison with the field for FFAD. As seen from Fig. 3, for a slow protocol here the required control fields to suppress fluctuations are much weaker than that presented in Fig. 1

Figure 2: (color online) Time dependence of -type optimal control field, choosing different weightage function (indicated on the panel) in the cost function defined in Eq. (4). The system considered here and other system parameters are the same as in Fig. 1. The three choices of lead to three different control fields, but all of them yield from trajectories.

Figure 3: (color online) Time dependence of a -type optimal control field, as compared with that in FFAD, for a parametric linear oscillator. All system parameters are the same as in Fig. 1, except that , , and .

Iv Optimal control of work fluctuations in nonlinear oscillators

We are now ready to apply our OCT approach to systems of nonlinear oscillators under certain protocols. Explicit solutions to realize FFAD and STA cannot be found for general nonlinear systems (excluding scale-invariant systems), but our OCT framework equally applies. In particular, let us consider the following system Hamiltonian with a time-dependent nonlinear term:


where is still the time-varying system parameter to model the same protocol defined in Eq. (15). In this system, the nonlinear term vanishes at . This choice is just for computational convenience when we sample the initial states from a thermal ensemble ensemble (which is still Gaussian). Below only -type optimal control are presented (note that, as observed earlier, the required field strength for -type control is much larger). In addition, for the sake of comparison, we also examine the parallel performance of the previous FFAD field obtained in the absence of the nonlinear term. To stress that the previous FFAD will no longer strictly give rise to accelerated adiabatic processes here due to the nonlinear term in the Hamiltonian, we call such a control approximate FFAD.

As the first computational example, we set and , , , and then numerically compute the fluctuations of . For approximate FFAD, we have ; whereas in our optimal control we have . Since these results are almost identical with our previous result in the absence of a nonlinear term, there must be a reason. The reason is simple. When the temperature is not high enough, the system is largely confined to a small neighborhood of the lowest energy state and hence the nonlinear term only plays a negligible role. This is consistent with a recent experiment Rohringer.13.a1 (), which showed that approximate FFAD can perform well in the presence of some degree of nonlinearity. For the same value of and , we have also considered a slow protocol with and . In this case both optimal control and approximate FFAD still yield .

Figure 4: (color online) Time dependence of a -type optimal control field obtained for a parametric nonlinear oscillator defined in Eq. (20), as compared with that of an approximate FFAD field obtained in the absence of the quartic term in the system Hamiltonian. The difference between the solid and dashed lines indicates the impact of the nonlinearity.

To enhance the nonlinear effect, we next consider a case with a much higher temperature, i.e., (still with , , and ). The time-dependence of the found optimal control field is presented in Fig. 5, which is seen to be very different from that of approximate FFAD. To understand the details better, we also count the number of trajectories that may give a negative work output. As shown previously, for the protocol here the work would be always positive if the nonlinear term were not introduced. So the presence of negative work values in the bare system does indicate the presence of nonlinearity. Detailed computational results are shown in Table 2. It is seen that the performance of optimal control is much better than that achieved by approximate FFAD. While the optimal control field has essentially suppressed almost all negative work values (but one out of one million), approximate FFAD increases probabilities of negative work. The reason why our optimal control approach is so effective in removing negative work values is quite simple. For cases with a positive , the value of with a negative would be too drastically larger than its “target” value , so it will be rejected by our optimization algorithm as extreme or rare values of . Note also that the variance in work under the optimal control field is less than one quarter of that for the bare system, and less than half of that obtained with approximate FFAD. The mean work under the optimal control field is also significantly smaller than that in the bare system or in the FFAD case. Because is an interesting quantity called the dissipated work, a decrease in indicates less dissipated work. We also show in Fig. 5 the work probability distribution . There it is seen that the optimal control field most effectively suppresses very large or very small work values, in addition to an almost complete removal of negative work values.

 Process  Probability of negative work
bare system
FFAD control
optimal control of -type
Table 2: The performance of work fluctuation suppression in the absence or presence of several different control fields, mainly characterized by the variance in work and the variance in , using trajectories. The system is a parametric nonlinear oscillator defined in Eq. (20), and the time dependence of is still given in Eq. (15), with , , and . The nonlinear parameter and the inverse temperature is set to be to enhance anharmonic effects.

Figure 5: (color online) Probability density of work distribution, in the absence of presence of different types of control fields. The results are obtained along with those presented in Table 2, with all the system parameters of a parametric nonlinear oscillator the same as that described in Table 2. Note that the optimal control field has suppressed a long tail of the work distribution and has also almost completely suppressed negative work values.

V Optimal Control of Work Fluctuations in Systems with Unknown Parameters

The flexibility in optimal control and early studies of feedback control motivated us to ask whether OCT can be used to suppress work fluctuations in those systems that have unknown system parameters. At the first thought, since with unknown system parameters, it is impossible to calculate and thus our optimal control approach seems not applicable. However, we show in this section that it is possible to refine the control target function on the go, by iteratively guessing the target function through Jarzynski’s equality. This possibility should be good news for the application of Jarzynski’s equality itself, as we can now indeed predict from non-equilibrium work values, with suppressed work fluctuations and hence better performance in reaching the correct based on a finite number of trajectories.

Figure 6: Variance of obtained from optimal control applied to a parametric nonlinear oscillator described in Eq. (20), as a function of incorrectly preassumed values of the nonlinear parameter . The actual nonlinear parameter . is still given by Eq. (15). Other system parameters are , , and . It is seen that despite wrong values of are used in searching for optimal control fields, the obtained variance of does not change significantly.

Let us first examine how our OCT framework depends on the knowledge of the bare system Hamiltonian. The first two relations in Eq. (7) require information of the initial thermal distribution and the target function. In other words, the full knowledge of the Hamiltonian in the beginning and at the end of the protocol is needed there. The third and fourth relations in Eq. (7) describe continuous time evolution under the total Hamiltonian including the control field. However, if the protocol is fast, then the main component of the total Hamiltonian can be the control field, so the bare system Hamiltonian may not play an important role there.

To have some quantitative ideas, we consider again the anharmonic Hamiltonian in Eq. (20), i.e., , with and . For this designed case, the quartic term vanishes at and , and hence the target function in our OCT is still fully known. Next we set , but in our construction and implementation of OCT we do not use this piece of knowledge. Instead, we use some wrong values of when computationally search for the OCT field. The knowledge of is used for checking only, after we have obtained the OCT field and start to look into the actual work fluctuations with the control field thus obtained.

As shown in Fig. 6, with used in our optimal control algorithm varying from to , the fluctuations in is not changing significantly. That is, for a wide range of incorrectly assumed values of , the associated optimal control field can still effectively suppress work fluctuations, with ranging from to . Additional numerical investigations of the evolving trajectories further indicate that the total time duration is too short for to play a role, thus confirming our qualitative insights above.

This computational example is enlightening, but we assume there that the (bare) system Hamiltonian is fully known at two boundary times and . Since in such situations can be exactly calculated from the bare system Hamiltonian, we are still one step away from predicting using the work values in a nonequilirium protocol. So our final question is the following, if some parameter in the (bare) system Hamiltonian is indeed unknown to us, how to construct the optimal control target function, suppress the work fluctuations, and eventually predict ?

Borrowing the idea from feedback optimal control theory Judson.92.PRL (), we aim to propose a useful computational feedback procedure as illustrated in Fig. 7. It consists of the following steps:

Figure 7: Procedure of executing a feedback loop to construct OCT fields iteratively, in order to apply OCT to systems with unknown system parameters. The feedback is to refine the control target function based on Jarzynski’s equality, by use of Jarzynski’s equality itself to guess the free energy difference.
  • A certain small number of initial states are first sampled according to system’s thermal distribution, and then evolved in accord with a guessed control field.

  • Based on the previous step, a rough estimate of is obtained to yield , which is then used to yield/update the control target function.

  • To optimize the control field based on Eq. (7), one may neglect the effect of the bare system Hamiltonian or preassume some wrong parameter values because the evolution is mainly dictated by the control field during a rapid work protocol.

  • The control field is then updated by the output from an optimal control algorithm and then all the previous steps are repeated until some convergence threshold is met.

bare system
Feedback OCT 
Feedback OCT 
Feedback OCT 
Table 3: The performance of work fluctuation suppression in the absence or presence of several different control fields, mainly characterized by the variance in work and the variance in , using trajectories. The system is a parametric nonlinear oscillator defined in Eq. (21), whose quartic term does not vanish at the end of the protocol. The coefficient of the quartic term is never used when searching for the control field via feedback OCT. The control field is obtained by five iterations of the feedback loop illustrated in Fig. 7. Other system parameters are , , and , is given by Eq. (15).

To demonstrate our strategy we consider the following system


Note that in this case, the quartic term is not vanishing at , with assumed to be “unknown” when we seek the optimal control field. That is, we will not use this value to construct our control field. Instead, we use some pre-assumed wrong values of when solving Eq. (7). We still set the quartic term to be zero at for the convenience in initial state sampling. In the first iteration of our numerical experiment, we use a null control field to start with. We update the target function and the optimal control field based on trajectories only. Remarkably, the optimal control field can already converge well after only five iterations of the above four steps. In particular, we compare in Table 3 the results from trajectories, in terms of , , , and . In obtaining the numerical results we have used three different preassumed values of (all are much different from the real value). As seen from the third column of Table 3, for all three cases, the fluctuations in are effectively well suppressed (as compared with the bare case without a control field). The ensemble-average in all the three cases are also close to the true theoretical value 0.3272 (obtained by numerically computing the partition function of with ). Interestingly, the value of for the bare system case (second column, second row of Table 3) is still slightly away from this theoretical value. This suggests that in the bare system case a high-quality convergence towards Jarzynski’s equality has not been achieved with trajectories. Thus, the presence of a control field suppressing work fluctuations is seen to have, albeit slightly, enhanced the convergence of the simulation towards Jarzynski’s equality. To see this more clearly, we show in Fig. 8 how the numerically obtained average value of gradually converges towards the theoretical value , as the number of trajectories increases to . The horizontal line in Fig. 8 represents the true theoretical value. It is seen that all the three cases of feedback control have converged to the true theoretical value with about trajectories, but the bare system case still has a non-negligible shift from the theoretical result.

Returning to Table 3, from the last two columns it is also seen that the mean work (hence also the dissipated work) as well as the variance of is suppressed by the optimal control field by about one order of magnitude. This significant control over the work output and its fluctuations is achieved even though part of the system parameters is unknown to us! Interestingly, this does not mean that the variance of (see second column of Table 3) will be also significantly reduced by the control field. Qualitatively, note that a higher temperature induces a wider initial probability distribution and hence larger thermal fluctuations in work values, but on the other hand, when calculating the larger work values are still scaled down by the inverse temperature . As such, from the explicit computational example here it is learned that one should not underestimate the implications of a seemingly “small” suppression in the variance of .

Figure 8: (color online) The convergence of towards the theoretical value (indicated by the horizontal line) versus the number of trajectories used in the simulations. The system and all the system parameters are all specified in Table 3. It is seen that the three cases with optimal control fields obtained from feedback mechanism can yield better results, despite that the actual value of the nonlinear parameter is not needed in constructing the control fields.

Vi Concluding Remarks

In this paper we have proposed an optimal control approach aiming at the suppression of work fluctuations associated with a protocol applied to a thermal ensemble. Our approach is original in the sense that we directly exploits Jarzynski’s equality in defining our control target function. Indeed, the control target function is constructed to bias against the deviation of individual values of from their ensemble average , where is the free energy difference. This approach is shown to be very effective. In the case of a parametric oscillator, the performance of our optimal control approach is simply as good as previous methods based on accelerated adiabatic processes. More importantly, our optimal control approach can equally apply to rather arbitrary nonlinear systems. We also note that a recent study Lu.14.PRA () considered a somewhat related optimal control approach for nonlinear systems, but assuming weak anharmonic potential to validate a perturbative treatment. Our approach can however be applied to systems with strong nonlinear effects, as we have shown through nonlinear oscillator systems with a large quartic term. One might also think that since is needed in the construction of a target function, we may have to know the full information of the system Hamiltonian. Through simple numerical experiments, we showed that this intuition is not true, because a feedback mechanism can help us to refine the control target function on the fly. As such, we claim that the suppression of work fluctuations via a control field is possible even when we do not have full knowledge of the system.

We believe that this study is just a motivator and a starting point along the general issue of controlling and manipulating nonequilibrium thermodynamic properties. Many follow up studies are on the way. For example, it is of immediate interest to extend this study to quantum systems in order to account for quantum effects and quantum fluctuations. So far we have assumed sufficiently rapid work protocols during which system-bath interaction is neglected. In the future it is necessary to extend our approach here to open systems. To that direction we note an early study considering classical overdamped motion under a special type of optimization Schmiedl.07.PRL ().

Vii Ackowledgements

We are grateful to Peter Hanggi, Dario Poletti, and Manas Mukherjee for interesting discussions.

Viii Appendix: The Iteration Procedure of OCT

In this section we discuss numerical iteration details in seeking optimal control fields from Eq. (7). The procedure is rather standard in OCT Shi.90.JCP (). First, we divide the relations in Eq. (7) together with Hamilton’s equations of motion into four parts.


Next, in numerical calculations, we can consider a sufficiently large but finite phase space area, by considering the initial thermal distribution in Eq. (1). For each initial phase space point sampled, one integrates Hamilton’s equation from to , using a guessed field time dependence . Then the final values of the Lagrange multipliers and can be computed from Eq. (22. With this, Eq. (23) can be integrated backwards. Finally, with , , and all solved numerically, the field time dependence for the next iteration can be obtained by referring to Eq. (24). This iterative procedure is repeated until certain convergence criteria are satisfied. In our implementation of this procedure we simply set the starting to be zero. At the end of each iteration, we update the field time dependence for the next iteration as follows:


Basically, the iteration equation means that a difference between final state and target state is calculated after each iteration, and then a modification of the control field is made.


  • (1) C. Jarzynski, Phys. Rev. Lett. 78, 2690 (1997).
  • (2) C. Jarzynski, Phys. Rev. E 56, 5018 (1997).
  • (3) G. E. Crooks, Phys. Rev. E60, 2721 (1999).
  • (4) M. Campisi, P. Hänggi, and P. Talkner, Rev. Mod. Phys. 83, 771 (2011).
  • (5) O. Abah, J. Roßnagel, G. Jacob, S. Deffner, F. Schmidt-Kaler, K. Singer, and E. Lutz, Phys. Rev. Lett. 109, 203006 (2012).
  • (6) C. Bergenfeldt, P. Samuelsson, B. Sothmann, C. Flindt, and M. B¨¹ttiker, Phys. Rev. Lett. 112, 076803 (2014).
  • (7) K. Zhang, F. Bariani, and P. Meystre, Phys. Rev. Lett. 112, 150602 (2014).
  • (8) J. Deng, Q. Wang, Z. Liu, P. Hänggi, and J. B. Gong, Phys. Rev. E 88, 062122 (2013).
  • (9) Y. Zheng and D. Poletti, Phys. Rev. E 90, 012145 (2014).
  • (10) A. del Campo, J. Goold, and M. Paternostro, arXiv:1305.3223.
  • (11) S. Shi and H. Rabitz, J. Chem. Phys. 92, 364 (1990).
  • (12) A. P. Peirce, M. A. Dahleh, and H. Rabitz, Phys. Rev. A 37, 4950 (1988).
  • (13) S. Shi and H. Rabitz, Comput. Phys. Commun. 63, 71 (1991).
  • (14) K. Sundermann and R. de Vivie-Riedle, J. Chem. Phys. 110, 1896 (1999).
  • (15) Y. Ohtsuki, G. Turinici, and H. Rabitz, J. Chem. Phys. 120, 5509 (2004).
  • (16) M. Demirplak and S. A. Rice, J. Chem. Phys. 129, 154111 (2008).
  • (17) M. Demirplak and S. A. Rice, J. Phys. Chem. B 109, 6838 (2005).
  • (18) M. Demirplak and S. A. Rice, J. Phys. Chem. A 107, 9937 (2003).
  • (19) M. V. Berry, J. Phys. A: Math. Theor. 42, 365303 (2009).
  • (20) S. Masuda and K. Nakamura, Proc. R. Soc. London Ser. A 466, 1135 (2010).
  • (21) X. Chen, I. Lizuain, A. Ruschhaupt, D. Guéry-Odelin, and J. G. Muga, Phys. Rev. Lett. 105, 123003 (2010).
  • (22) X. Chen, A. Ruschhaupt, S. Schmidt, A. del Campo, D. Guéry-Odelin, and J. G. Muga, Phys. Rev. Lett. 104, 063002 (2010).
  • (23) S. Choi, R. Onofrio, and B. Sundaram, Phys. Rev. A 84, 051601 (2011).
  • (24) J. F. Schaff, P. Capuzzi, G. Labeyrie, and P. Vignolo, New J. Phys. 13, 113017 (2011).
  • (25) M. G. Bason, M. Viteau, N. Malossi, P. Huillery, E. Arimondo, D. Ciampini, R. Fazio, V. Giovannetti, R. Mannella, and O. Morsch, Nat. Phys. 8, 147 (2012).
  • (26) J. Zhang, J. H. Shim, I. Niemeyer, T. Taniguchi, T. Teraji, H. Abe, S. Onoda, T. Yamamoto, T. Ohshima, J. Isoya, and D. Suter, Phys. Rev. Lett. 110, 240501 (2013).
  • (27) C. Jarzynski, Phys. Rev. A 88, 040101 (2013).
  • (28) A. del Campo, Phys. Rev. A 84, 031606 (2011).
  • (29) S. Deffner, C. Jarzynski, and A. del Campo, Phys. Rev. X 4, 021013 (2014).
  • (30) E. Torrontegui, S. Ibanez, S. Martí nez-Garaot, M. Modugno, A. del Campo, D. Guéry-Odelin, A. Ruschhaupt, X. Chen, and J. G. Muga, Adv. At. Mol. Opt. Phys. 62, 117-169 (2013).
  • (31) Z. C. Tu, Phys. Rev. E 89, 052148 (2014).
  • (32) A. del Campo, Phys. Rev. Lett. 111, 100502 (2013).
  • (33) C. Jarzynski, C. R. Phys. 8, 495 (2007).
  • (34) W. Rohringer, D. Fischer, F. Steiner, I. E. Mazets, J. Schmiedmayer, and M. Trupke, arXiv 1312.5948 (2013).
  • (35) R. S. Judson and H. Rabitz, Phys. Rev. Lett. 68, 1500 (1992).
  • (36) X. J. Lu, X. Chen, J. Alonso, and J. G. Muga, Phys. Rev. A 89, 023627 (2014).
  • (37) T. Schmiedl and U. Seifert, Phys. Rev. Lett. 98, 108301 (2007).
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description