# Driving rapidly while remaining in control: classical shortcuts from
Hamiltonian to stochastic dynamics

David Guéry-Odelin1, Christopher Jarzynski2,3,4, Carlos A. Plata5, Antonio
Prados5 and Emmanuel Trizac6 1 Laboratoire Collisions, Agrégats, Réactivité,
IRSAMC, Université de Toulouse, CNRS, UPS, France
2 Department of Chemistry and Biochemistry, University of Maryland, College
Park, MD, USA
3 Institute for Physical Science and Technology, University of Maryland,
College Park, MD, USA
4 Department of Physics, University of Maryland, College Park, MD, USA
5 Física Teórica, Universidad de Sevilla, Apartado de Correos 1065, E-41080
Sevilla, Spain
6 Université Paris-Saclay, CNRS, LPTMS, 91405, Orsay, France

###### Abstract

Stochastic thermodynamics lays down a broad framework to revisit the venerable
concepts of heat, work and entropy production for individual stochastic
trajectories of mesoscopic systems. Remarkably, this approach, relying on
stochastic equations of motion, introduces time into the description of
thermodynamic processes—which opens the way to fine control them. As a result,
the field of finite-time thermodynamics of mesoscopic systems has blossomed.
In this article, after introducing a few concepts of control for isolated
mechanical systems evolving according to deterministic equations of motion, we
review the different strategies that have been developed to realize finite-
time state-to-state transformations in both over and underdamped regimes, by
the proper design of time-dependent control parameters/driving. The systems
under study are stochastic, epitomized by a Brownian object immersed in a
fluid; they are thus strongly coupled to their environment playing the role of
a reservoir. Interestingly, a few of those methods (inverse engineering,
counterdiabatic driving, fast-forward) are directly inspired by their
counterpart in quantum control. The review also analyzes the control through
reservoir engineering. Besides the reachability of a given target state from a
known initial state, the question of the optimal path is discussed. Optimality
is here defined with respect to a cost function, a subject intimately related
to the field of information thermodynamics and the question of speed limit.
Another natural extension discussed deals with the connection between
arbitrary states or non-equilibrium steady states. This field of control in
stochastic thermodynamics enjoys a wealth of applications, ranging from
optimal mesoscopic heat engines to population control in biological systems.

  * November 2022

##  1 Introduction

Thermodynamics originated in efforts to tame the motive power of fire [1].
Originally concerned with notions of heat and temperature, the field was
formalized during the 19th century into a set of universal principles that
govern the properties of macroscopic systems in thermal equilibrium, as well
as transformations between equilibrium states [2]. Beginning in the 1860s,
statistical physics revealed the microscopic roots of thermodynamics, greatly
enhancing its power to predict and explain systems’ material properties.
Reversible transformations have traditionally played a central role in the
foundations of thermodynamics. Such transformations occur in an idealized
adiabatic (infinitely slow [3]) limit, in which a system’s dynamical behavior
is essentially irrelevant. More recent developments, however, have focused on
nonequilibrium, finite-time processes, where dynamics become important. In
particular the growing field of stochastic thermodynamics [4, 5, 6] extends
the concepts of heat, work and entropy production to individual trajectories
of microscopic systems, evolving under stochastic equations of motion.

By introducing time into the description of thermodynamic processes, it
becomes possible to formulate new questions, and in particular to investigate
the control and optimization of finite-time thermodynamic transformations, see
Fig. 1. Recent years have seen a surge of activity in this area. In this
review, we survey finite-time control methods that have been developed within
the framework of stochastic thermodynamics. Many of these methods (though not
all!) have been inspired by developments in the control of quantum systems,
particularly in the field of shortcuts to adiabaticity [7]. For this reason,
we begin with a brief overview of quantum frameworks and features that are
especially relevant for control ideas in stochastic thermodynamics.

The goal of quantum shortcuts to adiabaticity is to steer a system to evolve
from an eigenstate $|n(0)\rangle$ of an initial Hamiltonian $\widehat{H}(0)$
to the corresponding eigenstate $|n(t_{\rm f})\rangle$ of a final Hamiltonian
$\widehat{H}(t_{\rm f})$. For infinitely slow driving, this is achieved
automatically by virtue of the quantum adiabatic theorem [8], which guarantees
that the system remains in the instantaneous eigenstate $|n(t)\rangle$ of
$\widehat{H}(t)$ at all times. For rapid driving, three broad approaches for
achieving the above-mentioned goal have emerged: inverse engineering [9, 10],
transitionless or counterdiabatic driving [11, 12, 13], and fast-forward
methods [14, 15].

In inverse engineering methods, instead of deducing a system’s evolution under
a given driving protocol (as is the customary approach in physics), one seeks
to engineer a driving protocol that produces the desired evolution. This is
accomplished by exploiting the equations of motion that govern the system’s
dynamics. In quantum mechanics, inverse engineering methods have been applied
to a variety of dynamical frameworks, including the Schrödinger equation, the
evolution operator, the dynamical invariants, and the density matrix formalism
[7].

In the counterdiabatic approach, for a given time-dependent reference
Hamiltonian $\widehat{H}(t)$, one seeks a Hamiltonian $\widehat{H}_{CD}(t)$
with the following property: if the system evolves unitarily under
$\widehat{H}+\widehat{H}_{CD}$ from an initial eigenstate $|n(0)\rangle$, then
throughout the process the system remains in the $n$’th eigenstate of
$\widehat{H}(t)$. In other words the system follows the adiabatic trajectory
$|n(t)\rangle$, even when the driving is rapid. Explicit expressions for
$\widehat{H}_{CD}$ are given by Eqs. (10) and Eq. (11lbd) below [11, 12, 13,
16]. Both results steer the system exactly along the adiabatic trajectory
$|n(t)\rangle$, hence the solution to the counterdiabatic problem is not
always unique.

In the fast-forward approach [14, 15], a potential
$\widehat{U}_{FF}(t)=U_{FF}(\widehat{\mathbf{x}},t)$ is designed such that, if
the system evolves under $\widehat{H}+\widehat{U}_{FF}$ from an initial
eigenstate $|n(0)\rangle$, then at $t=t_{\rm f}$ the system arrives at the
desired final state $|n(t_{\rm f})\rangle$. At intermediate times the state of
the system takes the form
$\psi(\mathbf{x},t)=e^{iS(\mathbf{x},t)/\hbar}\langle\mathbf{x}|n(t)\rangle$,
where $S(\mathbf{x},t)$ is real. Note that while $\psi(\mathbf{x},t)$ itself
is not an eigenstate of $\widehat{H}(t)$, its coordinate-space probability
distribution coincides with the eigenstate probability distribution:
$|\psi(\mathbf{x},t)|^{2}=|\langle\mathbf{x}|n(t)\rangle|^{2}$.

In the shortcuts described above, the goal is to make the system arrive
rapidly at a destination it would have reached naturally had the process been
carried out quasistatically. In the quantum case, the desired destination is
an energy eigenstate, but as we shall see the same goal can be reformulated
for classical systems governed by Hamilton’s equations (Sec. 2), and for
stochastic systems evolving under overdamped and underdamped Brownian dynamics
(Secs. 3.1 \- 3.3). In fact, not only is the goal the same—rapid evolution to
a quasistatic outcome—but there are close similarities between the various
quantum, classical and stochastic shortcuts designed to achieve this goal. A
number of strategies for constructing shortcuts can be unified within a
framework organized around the continuity equation [16]. In each case the
strategy involves identifying a velocity field $v(x,t)$, or else the
corresponding acceleration field $a(x,t)$, then using this field to construct
the counterdiabatic or fast-forward Hamiltonian or potential—see Eqs. (11lt),
(11lz), (11lbd) and (11lbr).

Figure 1: The White Rabbit is doing his best to overcome his unpunctuality,
which is not an easy task in Wonderland. In this review, we study a related
type of question, and first address the possibility of finding a path from A
to B (A and B are fixed). When several such paths may exist, a second level of
question deals with searching for the optimal one; for the White Rabbit, this
amounts to finding the fastest.

Although the term shortcuts to adiabaticity (STA) is widely used in the
context of quantum and Hamiltonian classical dynamics, for stochastic systems
other expressions such as engineered swift equilibration [17] and shortcuts to
isothermality [18] have been introduced. While these terms are descriptive
within particular contexts, they do not fully capture the broad scope of
methods that have been developed. In this review we will use the general
terminology, swift state-to-state transformations, and the acronym SST, to
embrace the entire catalogue of shortcutting—see A for a summary of acronyms
used throughout.

This review focuses primarily on systems evolving under stochastic dynamics,
but we begin in Sec. 2 by describing shortcuts developed for classical
mechanical systems governed by Hamiltonian dynamics. Mainly, we describe the
generalisation to classical mechanics of counterdiabatic and fast-forward
driving. These methods, both in their quantum and classical versions, are
based on the time manipulation of the system Hamiltonian. Also, we briefly
describe other shortcuts—for instance in the framework of the Boltzmann
equation. Section 3 is devoted to shortcuts for systems in contact with a
thermal bath and thus described by stochastic dynamics, codified for instance
by the Fokker-Planck equation—focusing on rapid transitions between
equilibrium states. Therein, we further generalize three main STA quantum
approaches (inverse engineering, counterdiabatic method, and fast-forward) to
the stochastic framework, by time manipulating the potential energy, i.e. by
applying a suitably chosen external force. In addition, we consider a new
possibility of tuning the evolution of these systems: engineering the thermal
environment.

Loosely speaking, it may argued that the finite-time driving associated with
SST, which involves the exploration of non-equilibrium states, requires more
resources than an infinitely slow adiabatic connection—for which the system is
at equilibrium for all times. To be concrete, for an isothermal process the
work performed on the system is minimized when the process is performed
infinitely slowly, and in this case the work is equal to the free energy
difference between the final and initial states. For SST, there emerges an
irreversible contribution to the work that depends on the path swept by the
system. It is therefore natural to explore how to design protocols that
optimally use the resources available by minimizing some “cost” function, e.g.
the irreversible work for the isothermal connection. Also, the system
parameters that become time-dependent in SST must often verify certain
constraints. The question of minimizing a certain cost function that is a
functional of the trajectory, with constrained parameters, is the central
problem of optimal control theory [19, 20]—a collection of tools that has been
employed in the applied mathematics and engineering literature for a long
time, but only recently in physics [21, 22]. In Sec. 4, we analyze how these
tools translate to the setting of stochastic thermodynamics, including a brief
discussion of information geometry ideas and the emergence of the so-called
classical speed limits. Then, in Sec. 5, we extend previously discussed
techniques to transitions between non-equilibrium states. Section 6 deals with
applications of stochastic shortcuts, including heat engines. Finally,
perspectives and conclusions are drawn in Sec. 7.

The stochastic dynamics considered in Secs. 3-6 account for the effects of
thermal noise in the surrounding bath. However, the resulting randomness
generally becomes negligible (that is, of relative order $N^{-1/2}$ where $N$
is the number of degrees of freedom of the system) in the macroscopic limit
$N\rightarrow\infty$. In this limit the effects of the bath can typically be
modeled by deterministic equations of motion. For this reason our review
generally focuses on systems with one or a few degrees of freedom, where the
role of noise can not be neglected.

##  2 Shortcuts for isolated classical systems

###  2.1 Background and Setup

As indicated above, the general goal of STA is to steer a system to arrive
quickly at a destination that it would have reached without external
assistance (so to speak) if the process were performed infinitely slowly. For
quantum shortcuts, this desired destination is an energy eigenstate. For
shortcuts involving isolated classical systems, the goal is most naturally
framed in terms of adiabatic invariants, therefore we begin this section with
a discussion of classical adiabatic invariants.

The study of these invariants traces back more than a century [23, 24] to the
problem of a simple pendulum whose length $\ell$ varies with time. If $\ell$
changes slowly, then so too do both the pendulum frequency $\nu$ and its
energy $E$. However, in the harmonic regime of small oscillations the ratio
$E/\nu$ remains fixed in the limit of infinitely slow variation of the
pendulum length.111 Note the correspondence between the classical adiabatic
invariant, $E/\nu$, and the quantum number $n$: since the $n$’th eigenenergy
of a harmonic oscillator is given by $E_{n}=h\nu\left(n+\frac{1}{2}\right)$,
the quantum adiabatic invariance of $n$ matches the classical adiabatic
invariance of $E/\nu$.

More generally, consider a classical particle of mass $m$ in one degree of
freedom, described by a Hamiltonian

$H(z,t)=\frac{p^{2}}{2m}+U(x,t),$ (1)

where $z=(x,p)$ denotes a point in phase space. For any fixed value of $t$, we
take $U(x,t)$ to be a confining potential, with closed (periodic) orbits in
phase space. We assume that $U(x,t)$ has a single minimum whose location may
depend on $t$;222For potentials with two or more local minima, the adiabatic
invariance of the action breaks down when phase space separatrices are crossed
[25] To exclude this complicating feature, we assume only a single minimum.;
that $U(x,t)$ varies with time only during the interval
$t_{\scriptsize\textnormal{i}}=0\leq t\leq t_{\scriptsize\textnormal{f}}$; and
that the time-dependence of the potential is turned on and off smoothly at
$t=0$ and $t=t_{\scriptsize\textnormal{f}}$, respectively—more precisely,
$U(x,t)$ is twice differentiable with respect to time. Under these conditions,
if we treat $t$ as fixed parameter then every trajectory evolving under
$H(z,t)$ is closed (i.e. periodic). If we instead let $t$ denote the running
time then the action

$I(E,t)=\oint_{E}dx\,p$ (2)

is an adiabatic invariant [26]. The right side indicates a clockwise line
integral around the energy shell $E$, that is the level set $H(z,t)=E$, which
forms a closed loop in phase space (see Fig. 2). We will use the notation
${\cal E}(t,I)$ to denote the energy shell of $H(z,t)$ whose action is $I$.

Imagine a classical trajectory $z(t)$ that evolves under Hamilton’s equations
of motion as $H(z,t)$ is varied extremely slowly over a long time interval
$0\leq t\leq t_{\scriptsize\textnormal{f}}$. Let $E(t)\equiv H(z(t),t)$ denote
the slowly evolving energy of this trajectory. In stating that the action $I$
is an adiabatic invariant, we mean that in the limit of infinitely slow
driving its value remains constant along this trajectory:

$I(E(t),t)=I(E(0),0)\equiv I_{\scriptsize\textnormal{i}}\quad,\quad 0\leq
t\leq t_{\scriptsize\textnormal{f}}$ (3)

even though, in general, $E(t)\neq E(0)$. For a harmonic oscillator with time-
dependent frequency $\nu(t)$, the action defined by Eq. (2) is equal to
$I(E,t)=E/\nu$, in agreement with the discussion of the simple pendulum,
above.

The integral in Eq. (2) gives the phase space volume 333By convention we use
the term “volume” rather than “area”, even though phase space is two-
dimensional. enclosed by the energy shell, hence we can equally well write

$I(E,t)=\int dz\,\theta\left[E-H(z,t)\right]$ (4)

where $\theta(\cdot)$ is the unit step function, and $\int dz=\int dx\int dp$
denotes integration over phase space. Thus the adiabatic invariance of the
action can be described as follows: if a trajectory is initially located on an
energy shell ${\cal E}(0,I_{\scriptsize\textnormal{i}})$ that encloses phase
space volume $I_{\scriptsize\textnormal{i}}$, then for any
$t\in[0,t_{\scriptsize\textnormal{f}}]$ the trajectory will be located on the
energy shell ${\cal E}(t,I_{\scriptsize\textnormal{i}})$, which encloses the
same amount of phase space.

Classical shortcuts to adiabaticity (SST) are concerned with the situation in
which the variation of the potential $H(z,t)$ is not slow. Consider again the
Hamiltonian given by Eq. (1), only now imagine that the interval
$[0,t_{\scriptsize\textnormal{f}}]$ over which $H$ varies with time is not
particularly long—indeed, it can be arbitrarily short, though finite. In this
situation, the action $I$ is generally not invariant:

$I(E(t),t)\neq I(E(0),0)\,.$ (5)

We will view the non-invariance of the action as a deficiency, to be corrected
using tools similar to those devised for quantum SST/ STA.

In what follows, the term slow driving refers to the adiabatic limit, while
fast driving denotes non-adiabatic time-dependence of $H(z,t)$.

While the discussion above has focused on the evolution of a single trajectory
$z(t)$, the invariance of the action for slow driving and the breaking of that
invariance for rapid driving are conveniently visualized in terms of closed
loops evolving in phase space; see Fig. 2. Imagine, at $t=0$, a collection of
infinitely many initial conditions distributed over a single energy shell of
the initial Hamiltonian, with action $I_{\scriptsize\textnormal{i}}$. These
initial conditions define a loop ${\cal L}_{\scriptsize\textnormal{i}}$ that
coincides with the energy shell, ${\cal L}_{\scriptsize\textnormal{i}}={\cal
E}(0,I_{\scriptsize\textnormal{i}})$, as depicted schematically in Fig. 2(a).
From each of these initial conditions a trajectory $z(t)$ evolves under
$H(z,t)$. At any later time, $t>0$, a snapshot of these trajectories defines a
new closed loop ${\cal L}(t)$. Under slow driving, Eq. (3) implies that the
loop ${\cal L}(t)$ “clings” at all time to the instantaneous shell whose
action is $I_{\scriptsize\textnormal{i}}$, that is, ${\cal L}(t)={\cal
E}(t,I_{\scriptsize\textnormal{i}})$ for all
$t\in[0,t_{\scriptsize\textnormal{f}}]$. In particular, the initial loop
${\cal L}_{\scriptsize\textnormal{i}}={\cal L}(0)$ is mapped onto a final loop
${\cal L}_{\scriptsize\textnormal{f}}\equiv{\cal
L}({t_{\scriptsize\textnormal{f}}})={\cal
E}(t_{\scriptsize\textnormal{f}},I_{\scriptsize\textnormal{i}})$ that
coincides with an energy shell of the final Hamiltonian (see the gray loop in
Fig. 2(b)). In terms of the initial and final energies of the trajectory, we
have

$I(E_{\scriptsize\textnormal{i}},0)=I(E_{\scriptsize\textnormal{f}},t_{\scriptsize\textnormal{f}})=I_{\scriptsize\textnormal{i}}\,,$
(6)

where $E_{\scriptsize\textnormal{i}}=E(0)$ and
$E_{\scriptsize\textnormal{f}}=E(t_{\scriptsize\textnormal{f}})$. Under fast
driving, by contrast, the loop ${\cal L}(t)$ strays away from the energy shell
${\cal E}(t,I_{\scriptsize\textnormal{i}})$ as illustrated for
$t=t_{\scriptsize\textnormal{f}}$ in Fig. 2(b).

Figure 2: In (a), the gray loop depicts an energy shell ${\cal
E}(0,I_{\scriptsize\textnormal{i}})$ of $H(z,0)$, and the red loop ${\cal
L}_{\scriptsize\textnormal{i}}={\cal L}(0)$ depicts initial conditions for a
set of trajectories that subsequently evolve under $H(z,t)$, where $z=(x,p)$.
The two loops are identical at $t=0$, as shown. In (b), the gray loop depicts
the energy shell ${\cal
E}(t_{\scriptsize\textnormal{f}},I_{\scriptsize\textnormal{i}})$ of
$H(z,t_{\scriptsize\textnormal{f}})$, while the red loop ${\cal
L}_{\scriptsize\textnormal{f}}={\cal L}(t_{\scriptsize\textnormal{f}})$
depicts the final conditions for the set of trajectories. In the adiabatic
limit the two would coincide, ${\cal L}_{\scriptsize\textnormal{f}}={\cal
E}(t_{\scriptsize\textnormal{f}},I_{\scriptsize\textnormal{i}})$, as the
action $I$ becomes invariant. For any choice of
$t_{\scriptsize\textnormal{f}}$, the phase space volume enclosed by both loops
in (b) is equal to that enclosed by the loops in (a), by Liouville’s theorem,
even though the loop ${\cal L}_{\scriptsize\textnormal{f}}$ is not itself an
energy shell.

The goal of classical SST can be stated as follows: for a rapidly driven
Hamiltonian $H(z,t)$ and an action $I_{\scriptsize\textnormal{i}}$, devise a
strategy that evolves the loop ${\cal L}_{\scriptsize\textnormal{i}}={\cal
E}(0,I_{\scriptsize\textnormal{i}})$ to the loop ${\cal
L}_{\scriptsize\textnormal{f}}={\cal
E}(t_{\scriptsize\textnormal{f}},I_{\scriptsize\textnormal{i}})$, under
Hamiltonian evolution. This is analogous to the situation in quantum SST,
where the goal is to use unitary evolution to evolve a wavefunction from an
eigenstate of an initial Hamiltonian to the corresponding eigenstate of the
final Hamiltonian.

In the classical case, we attempt to meet this goal by designing an auxiliary
term $H_{\rm aux}(z,t)$ such that the desired evolution is generated by the
Hamiltonian

$H_{\rm SST}(z,t)=H(z,t)+H_{\rm aux}(z,t)\,.$ (7)

The auxiliary term effectively steers the evolving loop ${\cal L}(t)$ to the
desired target, namely the final energy shell ${\cal
E}(t_{\scriptsize\textnormal{f}},I_{\scriptsize\textnormal{i}})$. While the
dynamics are generated by $H_{\rm SST}(z,t)$, we emphasize that the energy
shells ${\cal E}(t,I)$ are always defined with respect to the original
Hamiltonian $H(z,t)$.

We now identify three different flavors of this problem, which differ from one
another in how ambitiously the above-mentioned goal is addressed.

  * •

In one version of the problem, $H_{\rm aux}$ is independent of the choice of
action $I_{\scriptsize\textnormal{i}}$, and the action remains invariant
throughout the entire process, i.e. Eq. (3) is satisfied at all times and for
any choice of $E_{\scriptsize\textnormal{i}}=E(0)$. This represents the
strongest version of the problem, and we will refer to it as global
counterdiabatic (GCD) driving.

  * •

In a somewhat more relaxed version, we continue to insist that the action be
invariant throughout the process (Eq. (3)), but now we allow $H_{\rm aux}$ to
depend on the choice of $I_{\scriptsize\textnormal{i}}$, equivalently on the
choice of $E_{\scriptsize\textnormal{i}}$; we will call this local
counterdiabatic (LCD) driving.

  * •

In the most relaxed version of the problem, we allow $H_{\rm aux}$ to depend
on the choice of $I_{\scriptsize\textnormal{i}}$ (as with LCD driving), and we
further allow the invariance of the action to be broken at intermediate times
$0<t<t_{\scriptsize\textnormal{f}}$, insisting only that
$I(E_{\scriptsize\textnormal{f}},t_{\scriptsize\textnormal{f}})=I(E_{\scriptsize\textnormal{i}},0)$
(Eq. (6)). We refer to this version as fast-forward (FF) driving.

The designations global and local reflect the distinction that in one case
(GCD) a single $H_{\rm aux}$ must succeed for all energy shells, i.e.
globally, while in the other (LCD) we are free to design $H_{\rm aux}$ based
on the energy shell under consideration, i.e. locally. The terms
counterdiabatic and fast-forward are taken from the literature on quantum
shortcuts—the former signifies the suppression of non-adiabatic excitations
[11], while the latter evokes a mechanism for rapidly arriving at a desired
final destination, regardless of the intermediate path taken to get there
[14].

All three versions of the classical SST problem—as defined above, for one
degree of freedom—have been solved, in the sense that in each version an
explicit recipe has been devised for constructing an auxiliary Hamiltonian
that achieves the desired steering of trajectories. As described in greater
detail below, for GCD driving the auxiliary Hamiltonian is generally a
complicated, non-linear function of both position ($x$) and momentum ($p$).
For LCD driving, the auxiliary Hamiltonian is simpler, taking the form $H_{\rm
aux}(z,t)=pv(x,t)$, with $v(x,t)$ given by Eq. (11ls) below. Finally, for FF
driving, the auxiliary Hamiltonian does not depend on momentum at all, in
other words the desired goal can be satisfied by adding a time dependent
potential energy function $U_{\rm FF}(x,t)$ to the original Hamiltonian.
Perhaps not surprisingly, as we relax our demands on the performance of the
auxiliary Hamiltonian, from global to local to fast-forward, the form of
$H_{\rm aux}$ becomes simpler, at least in its dependence on momentum: $H_{\rm
aux}(x,p,t)\rightarrow pv(x,t)\rightarrow U_{\rm FF}(x,t)$.

LCD and FF driving are closely related: the function $v(x,t)$ that appears in
the LCD auxiliary Hamiltonian (Eq. (11lt)) is also used to construct the fast-
forward potential $U_{\rm FF}$, as elaborated in Sec. 2.2. In general, GCD
driving is not closely related to the other two. However, for a particular
class of driving protocols that go by the name of scale-invariant driving
[27], GCD and LCD driving are identical, giving the same $H_{\rm aux}(z,t)$
(see Eq. (11lp) below). Scale-invariant protocols are those for which the
potential $U(x,t)$ in Eq. (1) has the form [28, 27]

$U(x,t)=\frac{1}{\sigma^{2}}U_{0}\left(\frac{x-\mu}{\sigma}\right)$ (8)

where $\sigma=\sigma(t)$ and $\mu=\mu(t)$ are functions of time. By varying
$\mu$ we translate the potential, while varying $\sigma$ stretches or squeezes
the potential along the $x$ axis, and rescales its magnitude, without
otherwise altering its profile. Potentials of this form give rise to
convenient scaling properties of both the classical dynamics and the quantum
energy eigenstates [27]. Power-law potentials $U(x,t)=\alpha\,[x/L(t)]^{c}$,
with $c=2,4,\cdots$, offer an illustrative example of scale-invariant driving
[29]. The harmonic oscillator with time-dependent stiffness ($c=2$) and the
particle-in-a-box, with time-dependent box length ($c\rightarrow\infty$)
belong to this class.

We note in passing that if $U(x,t)$ is scale-invariant (8) for some
$\sigma(t)$ and $\mu(t)$, then so is $U(x,t)+b/(x-\mu)^{2}$ for any real
constant $b$.444 This follows trivially, since
$b/y^{2}=(1/\sigma^{2})\,b/(y/\sigma)^{2}$, with $y=x-\mu$. Potentials of this
form will arise in Sec. 2.4.

###  2.2 Recipes for classical shortcuts to adiabaticity

We now describe how to construct $H_{\rm aux}(z,t)$ for the three situations
just outlined. For the derivations of these methods and further details, we
refer to the original papers, cited below.

In Secs. 2.2.1-2.2.3, we will assume that $U(x,t)$ is twice differentiable
with respect to time (see comments after Eq. (1)), hence both $U$ and
$\partial_{t}U$ are continuous functions of time. Since $\partial_{t}U=0$
outside the interval $[0,t_{\scriptsize\textnormal{f}}]$, this assumption
implies the boundary conditions

$\frac{\partial U}{\partial t}(x,0^{+})=\frac{\partial U}{\partial
t}(x,t_{\scriptsize\textnormal{f}}^{-})=0.$ (9)

In other words the time-dependence of $U$ is turned on and off smoothly rather
than abruptly. In Sec. 2.2.4 we will briefly consider the implications of
relaxing this assumption.

####  2.2.1 Global counterdiabatic driving

Classical GCD mirrors the quantum approach developed by Demirplak and Rice
[11] and Berry [13]. In the quantum case, the auxiliary or counterdiabatic
Hamiltonian is given as a sum over energy eigenstates of the original
Hamiltonian $\widehat{H}(t)$:

$\widehat{H}_{\rm GCD}(t)=i\hbar\sum_{m}\left(|\dot{m}\rangle\langle
m|-\langle m|\dot{m}\rangle|m\rangle\langle m|\right)\,.$ (10)

Here $|m\rangle=|m(t)\rangle$ denotes the $m$’th eigenstate of
$\widehat{H}(t)$, and $|\dot{m}\rangle=\partial_{t}|m(t)\rangle$ is its time
derivative.

We are interested in constructing the classical counterpart of
$\widehat{H}_{\rm GCD}(t)$. While Eq. (10) may not seem well-suited to this
end, the quantum operator defined (uniquely) by that equation can equivalently
be defined by [29]:

$\displaystyle\left[\widehat{H}_{\rm GCD},\widehat{H}\right]$ $\displaystyle=$
$\displaystyle i\hbar\left(\partial_{t}\widehat{H}-{\rm
diag}\,\partial_{t}\widehat{H}\right)$ (11a) $\displaystyle\langle
n|\widehat{H}_{\rm GCD}|n\rangle$ $\displaystyle=$ $\displaystyle
0\qquad\forall\,n$ (11b)

where ${\rm diag}\,\widehat{A}\equiv\sum_{m}|m\rangle\langle
m|\widehat{A}|m\rangle\langle m|$. Defining $\widehat{H}_{\rm GCD}$ in this
manner is convenient, as the correspondence between the quantum commutator
$[\widehat{A},\widehat{B}]$ and the classical Poisson bracket $\\{A,B\\}$,
together with the correspondence between quantum energy eigenstates and
classical energy shells, suggest a natural classical analogue:

$\displaystyle\left\\{H_{\rm GCD},H\right\\}$ $\displaystyle=$
$\displaystyle\partial_{t}H-\langle\partial_{t}H\rangle_{H}$ (11la)
$\displaystyle\langle H_{\rm GCD}\rangle_{E}$ $\displaystyle=$ $\displaystyle
0\qquad\forall\,E$ (11lb)

The notation

$\langle A\rangle_{E}=\frac{\int dz\,\delta\left(E-H\right)A(z)}{\int
dz\,\delta\left(E-H\right)}$ (11lm)

denotes a microcanonical average of an observable $A(z)$ over the energy shell
$E$ of $H(z,t)$. The left side of Eq. (11la) is evaluated at a phase point
$z$, and the notation $\langle\partial_{t}H\rangle_{H}$ on the right indicates
that the average is taken over the energy shell containing $z$; see Ref. [29]
for further details, as well as Ref. [30] where the same result is derived by
means of classical generating functions, and Ref. [31], where
$\widehat{H}_{\rm GCD}$ and $H_{\rm GCD}$ are framed as adiabatic gauge
potentials. Yet another interesting approach is taken in Ref. [32], where
$H_{\rm GCD}$ is constructed from the dispersionless Korteweg-de Vries
hierarchy, building on earlier work in the quantum context [33].

We have introduced Eqs. (11la) and (11lb) as the classical analogue of Eqs.
(11a) and (11b). It can be verified directly from classical analysis [29] that
the function $H_{\rm GCD}(z,t)$ defined by Eqs. (11la) and (11la) has the
counterdiabatic property we seek: along a trajectory $z_{G}(t)$ evolving under
$H+H_{\rm GCD}$, the action is preserved exactly:

$\frac{d}{dt}I\left[H(z_{G}(t),t),t\right]=0$ (11ln)

Hence the solution to a quantum problem, Eqs. (11a) and (11b), combined with
semiclassical reasoning, has yielded the exact solution to the analogous
classical problem, Eqs. (11la) and (11lb).

To translate Eqs. (11la) and (11lb) into an explicit function $H_{\rm
GCD}(z,t)$, note that Eq. (11la) implies the following relation, for any
points $z_{a}$ and $z_{b}$ on an energy shell $E$ of $H(z,t)$:

$H_{\rm GCD}(z_{b},t)-H_{\rm
GCD}(z_{a},t)=\int_{a}^{b}ds\,\left[\partial_{t}H(z(s),t)-\langle\partial_{t}H\rangle_{E}\right]\,,$
(11lo)

where $t$ is treated here as a parameter, and $z(s)$ is a trajectory of energy
$E$ that evolves with time $s$ under fixed $H(z,t)$, from $z(a)=z_{a}$ to
$z(b)=z_{b}$. This relation determines $H_{\rm GCD}(z,t)$ for all points on
the energy shell $E$, up to an additive constant whose value is in turn
determined by Eq. (11lb). This additive constant can depend on time, but it
has no dynamical relevance.

Equations (11la) and (11lb) thus provides an explicit recipe for constructing
$H_{\rm GCD}$. For scale-invariant driving (see Eq. (8)) this recipe leads to
the particularly simple expression [27]

$H_{\rm GCD}(z,t)\,=\,\frac{\dot{\sigma}}{\sigma}(x-\mu)p\,+\,\dot{\mu}\,p$
(11lp)

where the dots denote derivatives with respect to time. In the special cases
of the harmonic oscillator and the particle-in-a-box, Eq. (11lp) reproduces
results derived previously by other means [34, 35, 29].

For non-scale-invariant driving, it is difficult to obtain closed-form
expressions for $H_{\rm GCD}$, and one must resort to solving Eqs. (11la) and
(11lb) or (11lo) numerically. An exception is the case of the “tilted piston”,
which involves a particle inside a one-dimensional box with hard walls
separated by a distance $L$, with a potential inside the box that is linear in
$x$, with slope $\alpha$. Either $L$ or $\alpha$, or both, can be made time-
dependent. Even for this relatively modest extension of the particle-in-a-box,
the exact expression for $H_{\rm GCD}$ is complicated and non-linear in $p$
[36].

Finally, although Eqs. (11la) (11lb) uniquely defines $H_{\rm GCD}$, Eq.
(11ln) remains satisfied for a trajectory $z_{G}(t)$ evolving under $H+H_{\rm
GCD}+f(H)$, for any differentiable function $f(\cdot)$. As a result, to
construct a globally counterdiabatic auxiliary Hamiltonian we need only
satisfy Eq. (11la) and not necessarily Eq. (11lb) – imposing the latter (i.e.
setting $f=0$) amounts to a kind of gauge choice. Analogous comments apply in
the quantum case to Eqs. (11a) and (11b).

####  2.2.2 Local counterdiabatic driving

In the case of local counterdiabatic driving, we select a value of action,
$I_{\scriptsize\textnormal{i}}$. For any
$t\in[0,t_{\scriptsize\textnormal{f}}]$, ${\cal
E}(t,I_{\scriptsize\textnormal{i}})$ is the energy shell of $H(z,t)$ with the
same action $I_{\scriptsize\textnormal{i}}$; we refer to this shell as the
adiabatic energy shell for our choice of $I_{\scriptsize\textnormal{i}}$. This
energy shell represents the desired evolution we wish to generate: if a
trajectory begins on the energy shell ${\cal
E}(0,I_{\scriptsize\textnormal{i}})$ at $t=0$, we want to guarantee that it
will be found on the shell ${\cal E}(t,I_{\scriptsize\textnormal{i}})$ at all
$t\in[0,t_{\scriptsize\textnormal{f}}]$. Unlike with GCD, this goal can be
accomplished with an auxiliary Hamiltonian that is linear in momentum, $p$.

The LCD auxiliary Hamiltonian involves a function $v(x,t)$ that is constructed
as follows. At time $t$, the shell ${\cal E}(t,I_{\scriptsize\textnormal{i}})$
forms a loop in phase space, with left and right turning points, $x_{L}(t)$
and $x_{R}(t)$. Let

$\bar{p}_{\pm}(x,t)=\pm\left[2m({\cal
E}(t,I_{\scriptsize\textnormal{i}})-U(x,t))\right]^{1/2}$ (11lq)

denote the upper and lower branches of this loop, and let

${\cal S}(x,t)=2\int_{x_{L}(t)}^{x}dx^{\prime}\,\bar{p}_{+}(x^{\prime},t)$
(11lr)

denote the volume of phase space enclosed by ${\cal
E}(t,I_{\scriptsize\textnormal{i}})$ between the left turning point and a
vertical line located at position $x\in[x_{L},x_{R}]$. (To avoid confusion, we
stress that ${\cal S}$ is not Hamilton’s principal function, which appears in
the Hamilton-Jacobi equation [26].) Since ${\cal S}$ increases monotonically
with $x$ we can invert it, writing $x=x({\cal S},t)$, with ${\cal
S}\in[0,I_{\scriptsize\textnormal{i}}]$. We then define a velocity field

$v(x,t)=\frac{\partial}{\partial t}x({\cal S},t)=-\frac{\partial{\cal
S}/\partial t}{\partial{\cal S}/\partial x}$ (11ls)

using the cyclic identity of partial derivatives. From Eq. (9) it follows that
$v(x,0)=v(x,t_{\scriptsize\textnormal{f}})=0$. The LCD auxiliary Hamiltonian
is then given by

$H_{\rm LCD}(x,p,t)=pv(x,t).$ (11lt)

Adding $H_{\rm LCD}$ to the original Hamiltonian $H$ (as per Eq. (7)) gives us

$H_{\rm SST}(z,t)=\frac{p^{2}}{2m}+U(x,t)+pv(x,t)$ (11lu)

which generates the equations of motion

$\dot{x}=\frac{p}{m}+v(x,t)\quad,\quad\dot{p}=-\frac{\partial U}{\partial
x}-p\frac{\partial v}{\partial x}(x,t).$ (11lv)

Under these equations, the action is conserved exactly for any trajectory
$z_{L}(t)$ (where $L$ is short for $LCD$) launched from the initial adiabatic
energy shell [37]. That is,

$z_{L}(t)\in{\cal E}(t,I_{\scriptsize\textnormal{i}})$ (11lw)

for all $t\in[0,t_{\scriptsize\textnormal{f}}]$. For such a trajectory, the
first terms appearing on the right in Eq. (11lv) generate motion along the
instantaneous energy shell, while the second terms (involving $v$) force the
trajectory to remain attached to the evolving energy shell.

For scale-invariant driving, Eq. (11ls) gives

$v(x,t)=\frac{\dot{\sigma}}{\sigma}(x-\mu)+\dot{\mu}$ (11lx)

As mentioned shortly before Eq. (8), in this special case the GCD (Eqs. (11la)
and (11lb)) and LCD (Eq. (11lt)) prescriptions yield the same auxiliary
Hamiltonian, Eq. (11lp).

####  2.2.3 Fast-forward driving

Having defined $v(x,t)$ above, it is straightforward to construct the FF
auxiliary potential $U_{\rm FF}(x,t)$. We first introduce an acceleration
field

$a(x,t)=\frac{\partial v}{\partial x}v+\frac{\partial v}{\partial
t}=\frac{\partial^{2}}{\partial t^{2}}x({\cal S},t).$ (11ly)

$U_{\rm FF}$ is then defined, up to an arbitrary function of time, by

$-\frac{\partial U_{\rm FF}}{\partial x}=ma(x,t).$ (11lz)

Adding $U_{\rm FF}$ to $H$ gives

$H_{\rm SST}(z,t)=\frac{p^{2}}{2m}+U(x,t)+U_{\rm FF}(x,t)$ (11laa)

which generates the equations of motion

$\dot{x}=\frac{p}{m}\quad,\quad\dot{p}=-\frac{\partial U}{\partial
x}+ma(x,t).$ (11lab)

If two trajectories are launched from identical initial conditions on the
energy shell ${\cal E}(0,I_{\scriptsize\textnormal{i}})$, and one of them,
$z_{L}(t)$, evolves under Eq. (11lv), while the other, $z_{F}(t)$, evolves
under Eq. (11lab) (where $F$ is short for $FF$), then the two are related by
[37]

$x_{F}(t)=x_{L}(t)\quad,\quad p_{F}(t)=p_{L}(t)+mv(x_{L}(t),t)\,.$ (11lac)

The two trajectories reunite at $t=t_{\scriptsize\textnormal{f}}$, as
$v(x,t_{\scriptsize\textnormal{f}})=0$. Thus the fast-forward trajectory
$z_{F}(t)$ starts on the energy shell ${\cal
E}(0,I_{\scriptsize\textnormal{i}})$ at $t=0$, then strays from ${\cal
E}(t,I_{\scriptsize\textnormal{i}})$ at intermediate times, but ultimately
arrives at the final adiabatic energy shell:

$z_{F}(t_{\scriptsize\textnormal{f}})=z_{L}(t_{\scriptsize\textnormal{f}})\in{\cal
E}(t_{\scriptsize\textnormal{f}},I_{\scriptsize\textnormal{i}})$ (11lad)

In Sec. 2.5 of Ref. [31], the close relationship between LCD and FF driving is
described in terms of canonical gauge transformations that map Hamiltonians of
the form given by Eq. (11lu) into those given by Eq. (11laa). A similar
approach was used in Ref. [27], and analogous unitary transformations were
earlier introduced in the quantum context in Ref. [38].

For scale-invariant driving, Eqs. (11lx) - (11lz) lead to the following fast-
forward potential [27]:

$U_{\rm
FF}(x,t)=-\frac{m}{2}\frac{\ddot{\sigma}}{\sigma}(x-\mu)^{2}-m\ddot{\mu}x.$
(11lae)

For the special cases of the harmonic oscillator and particle-in-a-box,
equivalent results were obtained using an inverse engineering approach, in
Refs. [39, 40].

Let us focus briefly on the harmonic oscillator, for which the original
Hamiltonian is given by

$H(z,t)=\frac{p^{2}}{2m}+\frac{m}{2}\omega^{2}(t)x^{2},$ (11laf)

with
$\dot{\omega}(0)=\ddot{\omega}(0)=\dot{\omega}(t_{\scriptsize\textnormal{f}})=\ddot{\omega}(t_{\scriptsize\textnormal{f}})=0$.
The potential $U(x,t)$ in Eq. (11laf) can be cast into scale-invariant form
(see Eq. (8)) by setting $U_{0}(x)=mx^{2}/2$, $\sigma=\omega^{-1/2}$ and
$\mu=0$. Using Eq. (11lae), we can then combine $U$ and $U_{\rm FF}$ into a
single quadratic potential (see Eq. (11laa)):

$H_{\rm SST}(z,t)=\frac{p^{2}}{2m}+\frac{m}{2}\Omega^{2}(t)x^{2}$ (11lag)

with

$\Omega^{2}=\omega^{2}-\frac{3}{4}\frac{\dot{\omega}^{2}}{\omega^{2}}+\frac{1}{2}\frac{\ddot{\omega}}{\omega}\,.$
(11lah)

Rewriting the right side of Eq. (11lah) in terms of $\sigma$ rather than
$\omega$, and re-ordering terms, we get

$\ddot{\sigma}+\Omega^{2}\sigma=\sigma^{-3}\,.$ (11lai)

This is the Ermakov equation, which was used in Ref. [39] as the starting
point for designing a SST protocol for the harmonic oscillator. When $\Omega$
is constant, this equation admits a simple solution, since $\sigma^{2}$
follows a harmonic oscillator equation, see Eq. (11laz) below.

####  2.2.4 Boundary conditions in time

To this point, we have assumed that the time-dependence of $U$ is turned on
and off smoothly; see Eq. (9). For GCD and LCD driving this assumption is
unnecessary: as long as $U$ is once-differentiable with respect to time, the
results described in Secs. 2.2.1 and 2.2.2 remain valid. A discontinuity in
$\partial_{t}U$ at $t=0$ or $t=t_{\scriptsize\textnormal{f}}$ merely implies
that $H_{\rm aux}$ is turned on or off abruptly.

With FF driving the situation is subtler. If $\partial_{t}U$ is discontinuous
at $t=0$, then so is the velocity field $v(x,t)$, hence the acceleration field
$a(x,t)$ (Eq. (11ly)) is undefined. In this situation, the classical fast-
forward method of Sec. 2.2.3 can be salvaged if we apply the following
impulsive potential at $t=0$:

$U_{imp}(x,t)=-m\,\delta(t)\int^{x}dx^{\prime}\,v(x^{\prime},0^{+}).$ (11laj)

(The choice of the lower limit of integration is arbitrary.) This potential
causes a trajectory with initial phase space conditions $(x,p)$ at $t=0^{-}$
to jump suddenly to $(x,p+mv(x,0^{+}))$ at $t=0^{+}$. If $\partial_{t}U$ is
discontinuous at $t=t_{\scriptsize\textnormal{f}}$, then another impulse is
required:

$U_{imp}(x,t)=+m\,\delta(t-t_{\scriptsize\textnormal{f}})\int^{x}dx^{\prime}\,v(x^{\prime},t_{\scriptsize\textnormal{f}}^{-}).$
(11lak)

Once these impulses are included, the FF driving works as described in Sec.
2.2.3. In particular all points located on the initial energy shell ${\cal
E}(0,I_{\scriptsize\textnormal{i}})$ at $t=0^{-}$ evolve to points on the
final energy shell ${\cal
E}(t_{\scriptsize\textnormal{f}},I_{\scriptsize\textnormal{i}})$ at
$t=t_{\scriptsize\textnormal{f}}^{+}$.555If $\partial_{t}U$ is discontinuous
at intermediate times $t_{1},t_{2},\cdots$, then similar impulses at those
times can “rescue” the fast-forward method.

For an illustration of the effects of such impulses in the context of the
particle-in-a-box under fast-forward driving, see Sec. III.A of Ref. [27].

####  2.2.5 Beyond one degree of freedom

It is natural to ask whether the classical GCD, LCD and FF recipes discussed
above can be generalized to systems with $n>1$ degrees of freedom. This
question has largely remained unexamined in the literature.

When $n>1$, the definition of the adiabatic invariant itself depends on the
dynamics generated by the system Hamiltonian $H(\mathbf{x},\mathbf{p},t)$. If
the dynamics are integrable at all values of $t$, then the phase space
coordinates can be written in terms of action-angle variables
$(\mathbf{I},\mathbf{w})$, and the action variables $(I_{1},\cdots I_{n})$ are
adiabatic invariants [26]. It is plausible that the methods described above
could be applied to construct a separate auxiliary Hamiltonian $H_{\rm
aux}^{k}$ for each action-angle pair $(I_{k},w_{k})$, such that evolution
under $H+\sum_{k}H_{\rm aux}^{k}$ would preserve the value of each invariant.

If the dynamics are ergodic over the energy shell at all values of $t$, then
the phase space volume

$\Omega(E,t)=\int d^{n}\mathbf{x}\int
d^{n}\mathbf{p}\,\,\theta\left[E-H(\mathbf{x},\mathbf{p},t)\right]$ (11lal)

is the sole adiabatic invariant [41, 42, 43, 44, 45, 46, 47, 48]. Note that
the phase space volume is directly related to the Gibbs definition of entropy
in the microcanonical ensemble. In the special case when the ergodic
Hamiltonians $H(\mathbf{x},\mathbf{p},t)$, with
$t\in[0,t_{\scriptsize\textnormal{f}}]$, constitute a canonical
family,666$H(\mathbf{x},\mathbf{p},t)$ is a canonical family if $H(t_{1})$ can
be mapped to $H(t_{2})$ by a canonical transformation, for any
$t_{1},t_{2}\in[0,t_{\scriptsize\textnormal{f}}]$ [49]. the $n$-dimensional
analogue of Eqs. (11la) and (11lb) is known to have a solution [50],
suggesting that $H_{\rm GCD}$ can be constructed for canonical families of
ergodic Hamiltonians. However, it is not clear how to extend this approach to
the typical situation in which $H(\mathbf{x},\mathbf{p},t)$ is not a canonical
family.

Finally, integrable and ergodic Hamiltonians are themselves somewhat special:
a generic classical Hamiltonian in $n>1$ degrees of freedom has a mixed phase
space, with some trajectories evolving regularly and others irregularly
(chaotically) [51]. Adiabatic invariants for such systems are not as
unambiguously defined as for integrable or ergodic Hamiltonians, hence the
very notion of what constitutes a shortcut to adiabaticity becomes murky.

###  2.3 Classical shortcuts and microcanonical ensembles

The discussion in Sec. 2.2 focused on the evolution of individual trajectories
(Eqs. (11ln), (11lw), (11lad)). It is also useful to analyze this problem at
the statistical level.

In this section, we consider a microcanonical ensemble of initial conditions
on the energy shell ${\cal E}(0,I_{\scriptsize\textnormal{i}})$. This ensemble
can be understood as a probability density that is distributed over the loop
${\cal L}_{\scriptsize\textnormal{i}}$ shown in Fig. 2(a). Specifically, if we
use the angle $w\in[0,2\pi)$ of action-angle variables $(I,w)$ [26] to label
points around the loop, then the microcanonical distribution is uniform in
$w$. We emphasize that the microcanonical distribution represents a particular
choice of initial conditions for the ensemble; it does not arise from coupling
with a thermal bath. We now describe what happens to this ensemble, as
trajectories evolve from these initial conditions under global
counterdiabatic, local counterdiabatic, and fast-forward dynamics.

In the GCD case, the ensemble of trajectories clings to the adiabatic energy
shell ${\cal E}(t,I_{\scriptsize\textnormal{i}})$ during the process, as the
action is preserved (Eq. (11ln)). Moreover the ensemble remains microcanonical
at all times: for any $t\in[0,t_{\scriptsize\textnormal{f}}]$, a snapshot of
the ensemble of trajectories would show them to be distributed uniformly, with
respect to the angle variable $w$, on the adiabatic energy shell.

For LCD driving, the ensemble of trajectories also clings to the adiabatic
energy shell (Eq. (11lw)) only now the distribution does not necessarily
remain microcanonical. In particular, a snapshot at
$t=t_{\scriptsize\textnormal{f}}$ would reveal a collection of final
conditions that are distributed non-uniformly with respect to $w$ – see Fig.
3c of Ref. [37] for an illustration.

The result that GCD driving preserves the microcanonical distribution while
LCD driving (in general) does not, can be traced back to the fact that in the
former case the auxiliary Hamiltonian is the same for all
$I_{\scriptsize\textnormal{i}}$, whereas in the latter case $H_{\rm aux}$
generally depends on $I_{\scriptsize\textnormal{i}}$, as discussed in Ref.
[37], Appendix D.

For FF driving, the trajectories depart from the adiabatic energy shell at
intermediate times but return to that shell at the final time. From Eq.
(11lad) it follows that the final conditions are generally distributed non-
uniformly (i.e. non-microcanonically), as they coincide with the final
conditions achieved under LCD driving.

At intermediate times the LCD and FF phase space distributions differ.
However, the projections of these distributions onto the $x$-axis are
identical, by Eq. (11lac). An analogous situation holds in the quantum case:
the LCD and FF wavefunctions differ for
$t\in(0,t_{\scriptsize\textnormal{f}})$, but their $x$-space probability
distributions coincide, $|\psi_{FF}(x,t)|^{2}=|\psi_{LCD}(x,t)|^{2}=|\langle
x|n(t)\rangle|^{2}$.777 This conclusion follows by comparing Eqs.(6) and (22)
of Ref. [16], which give the exact time-dependent wavefunctions for LCD and FF
dynamics, respectively.

Under scale-invariant driving, GCD and LCD trajectories are identical, as
mentioned earlier. In this special situation, all three flavors of shortcuts
(GCD, LCD, FF) map an initial microcanonical distribution to a final
microcanonical distribution.

###  2.4 Shortcuts in classical kinetic theory

We now discuss a situation in which classical shortcuts arise in the context
of kinetic theory.

The discussion so far has focused on a single particle in a one-dimensional
potential, Eq. (1). The situation becomes more complicated when dealing with
$N>1$ mutually interacting particles. Under appropriate conditions, however,
SST-inspired tools have found interesting applications to many-body,
interacting systems.

Ref. [52] considers a dilute gas of identical particles of mass $m$, governed
by the Boltzmann equation

$\frac{\partial f}{\partial
t}+\mathbf{v}\cdot\nabla_{\mathbf{r}}f+\frac{1}{m}\mathbf{F}(\mathbf{r},t)\cdot\nabla_{\mathbf{v}}f=I_{\rm
coll}[\mathbf{v}|f,f]$ (11lam)

where $f(\mathbf{r},\mathbf{v},t)$ is the single-particle density,
$\mathbf{r}$ and $\mathbf{v}$ denote three-dimensional position and velocity,
$\mathbf{F}=-\nabla_{\mathbf{r}}U(\mathbf{r},t)$ is a conservative force, and
$I_{\rm coll}$ is a term that models the effect of two-body collisions. When
$I_{\rm coll}=0$ we get the collisionless equation

$\partial_{t}f+\mathbf{v}\cdot\nabla_{\mathbf{r}}f+\frac{\mathbf{F}}{m}\cdot\nabla_{\mathbf{v}}f=0\,,$
(11lan)

describing a gas of non-interacting particles. If the potential $U$ is time-
independent, then under Eq. (11lam) the gas generically evolves to a canonical
distribution, by Boltmann’s H theorem [53, 54]. However, Boltzmann realized
already in the 1870’s that when $U\propto r^{2}$ (with $r=|\mathbf{r}|$),
oscillatory “breathing mode” solutions of Eq. (11lam) also exist.

It is known that $I_{\rm coll}=0$ when $f$ has the form

$f(\mathbf{r},\mathbf{v},t)=\exp\left(-\alpha-\eta
v^{2}-\boldsymbol{\gamma}\cdot\mathbf{v}\right)$ (11lao)

where $\alpha$, $\eta$ and $\boldsymbol{\gamma}$ are arbitrary functions of
position and time [52]. In other words the distribution $f$ given by Eq.
(11lao) belongs to the kernel of $I$. Thus any solution of the collisionless
equation (Eq. (11lan)) that has the form given by Eq. (11lao), is also a
solution of the Boltzmann equation (Eq. (11lam)). Building on this insight,
the authors of Ref. [52] discovered novel exact solutions of the Boltzmann
equation. In particular, they constructed time-dependent potentials of the
form

$U(\mathbf{r},t)=\frac{m}{2}\omega^{2}(t)r^{2}+\frac{b}{r^{2}}\,,$ (11lap)

with appropriately engineered $\omega(t)$, that rapidly steer the dilute gas
from an initial canonical distribution at temperature
$T_{\scriptsize\textnormal{i}}$, to a final canonical distribution at
temperature $T_{\scriptsize\textnormal{f}}$, under the assumption that the
evolution of the gas is accurately described by the Boltzmann equation. Here,
$b\geq 0$ should be time independent.

We now extend these results to include scale-invariant potentials as well as
time-periodic driving. Consider the time-dependent, spherically symmetric
potential

$\overline{U}(\mathbf{r},t)=\frac{1}{\sigma^{2}}U_{0}\left(\frac{r}{\sigma}\right)+\frac{b}{r^{2}}\,,$
(11laq)

where $\sigma(t)>0$ is a twice-differentiable function of time, and $b\geq 0$
is a constant. Next, let $\beta(0)>0$ denote an inverse temperature, and
define

$\beta(t)=\frac{\sigma^{2}(t)}{\sigma^{2}(0)}\beta(0)\quad,\quad\bar{\mathbf{v}}(\mathbf{r},\mathbf{v},t)=\mathbf{v}-\frac{\dot{\sigma}(t)}{\sigma(t)}\mathbf{r}$
(11lar)

and

$U(\mathbf{r},t)=\overline{U}(\mathbf{r},t)-\frac{m}{2}\frac{\ddot{\sigma}(t)}{\sigma(t)}r^{2}$
(11las)

(compare with Eqs. (11lx), (11lae)). We claim that the time-dependent
distribution

$f(\mathbf{r},\mathbf{v},t)=\frac{1}{Z_{0}}e^{-\beta(t)[m\bar{v}^{2}/2+\overline{U}(\mathbf{r},t)]}$
(11lat)

is an exact solution of the Boltzmann equation (Eq. (11lam)), when the force
$\mathbf{F}(\mathbf{r},t)$ is obtained from the potential given by Eq.
(11las).888The time-in-dependence of the normalization constant $Z_{0}$
follows from the form of $\overline{U}(\mathbf{r},t)$, as can be shown by
inspection.

To establish this claim, note that Eq. (11lat) has the form given by Eq.
(11lao), therefore $I_{\rm coll}=0$. Also, it follows from direct substitution
that Eq. (11lat) satisfies the collisionless equation, Eq. (11lan). Thus Eq.
(11lat) solves Eq. (11lam) exactly.

We can extend these results further by replacing Eq. (11laq) with

$\overline{U}(\mathbf{r},t)=\frac{1}{\sigma^{2}}\sum_{j=1}^{3}U_{0j}\left(\frac{x_{j}}{\sigma}\right)+\,U_{h}(x_{1},x_{2},x_{3})$
(11lau)

where the $x_{j}$’s are the Cartesian components of $\mathbf{r}$, the
$U_{0j}$’s are three (generally unrelated) potential functions, and $U_{h}$ is
a homogeneous function of degree -2:

$U_{h}(x_{1},x_{2},x_{3})\,=\,\frac{1}{(x_{1}x_{2}x_{3})^{2/3}}\,\psi\left(\frac{x_{1}}{x_{2}},\frac{x_{1}}{x_{3}}\right),$
(11lav)

where $\psi$ is an arbitrary function. This form for $U_{h}$ includes a wide
range of potentials, including as special cases $b/r^{2}$ (see Eq. (11laq))
and $a_{1}/x_{1}^{2}+a_{2}/x_{2}^{2}+a_{3}/x_{3}^{2}$. The potential $U_{h}$
should bear no explicit time dependence, unlike $\overline{U}$, which depends
explicitly on time through $\sigma$. If we again define
$U(\mathbf{r},t)=\overline{U}-m\ddot{\sigma}r^{2}/2\sigma$ (see Eq. (11las)),
then $f(\mathbf{r},\mathbf{v},t)$ given by Eq. (11lat) remains an exact
solution of the Boltzmann equation (11lam). This result can be established as
in the previous paragraph.

If we set $U_{0}(x)=U_{0j}(x)=\kappa x^{2}/2$, for some fixed $\kappa>0$ and
for $j=1,2,3$, and $U_{h}=b/r^{2}$ for some $b\geq 0$, then Eqs. (11laq) and
(11lau) describe the same potential, and $U(\mathbf{r},t)$ has the form given
by Eq. (11lap). If $U_{0j}(x)=\kappa_{j}x^{2}/2$ for different fixed values
$\kappa_{1},\kappa_{2},\kappa_{3}>0$, and $U_{h}=0$, then
$\overline{U}(\mathbf{r},t)$ is an anisotropic harmonic potential of the form
considered in Ref. [55], where exact solutions of the Boltzmann equation were
derived for such potentials. Note that, in general, even if the three
$U_{0j}$’s are identical functions, the resulting potential
$\overline{U}(\mathbf{r},t)$ (Eq. (11lau)) is not spherically symmetric.

We can use these results to design a protocol for driving the gas from a
canonical distribution at initial inverse temperature $\beta(0)^{-1}$ to a
canonical distribution at final inverse temperature
$\beta(t_{\scriptsize\textnormal{f}})^{-1}$, in an arbitrary, finite time
$t_{\scriptsize\textnormal{f}}$. To do this, we simply choose $\sigma(t)$ to
satisfy the boundary conditions $\dot{\sigma}=\ddot{\sigma}=0$ at $t=0$ and
$t=t_{\scriptsize\textnormal{f}}$, and

$\sigma(t_{\scriptsize\textnormal{f}})/\sigma(0)=\sqrt{\beta(t_{\scriptsize\textnormal{f}})/\beta(0)}\,.$
(11law)

This protocol can be used with either of the potential forms given by Eqs.
(11las) or (11lau). In the isotropic and anisotropic harmonic cases discussed
in the previous paragraph, these protocols reduce to the ones obtained in
Refs. [52] and [55], respectively.

Alternatively, we can design a driven breathing mode of the Boltzmann equation
by choosing the protocol

$\sigma(t)=\bar{\sigma}+\Delta\sigma\cos^{2}(\omega t)\quad,\quad
0<\Delta\sigma<\bar{\sigma}$ (11lax)

again using either Eq. (11las) or Eq. (11lau). The exact solution given by Eq.
(11lat) then oscillates periodically in time, with frequency $2\omega$.

Finally, let us consider how the above results relate to the undriven
breathing modes discovered by Boltzmann [52]. In Eq. (11laq), take
$U_{0}(x)=mx^{2}/2$ and $b=0$ so that $\overline{U}(\mathbf{r},t)$ becomes a
spherically symmetric harmonic oscillator with a time-dependent stiffness
$\sigma^{-4}(t)$. By Eq. (11las) this choice leads to the total potential

$U(\mathbf{r},t)=\frac{m}{2}\left(\frac{1}{\sigma^{4}}-\frac{\ddot{\sigma}}{\sigma}\right)r^{2}\equiv\frac{m}{2}\Omega^{2}(t)r^{2}\,.$
(11lay)

Note that the relationship between $\sigma$ and $\Omega$ here is described by
the Ermakov equation (11lai). If we now take $\Omega(t)$ to be constant rather
than time-dependent, and we solve for $\sigma(t)$ (see Eq. (11laz) below),
then Eq. (11lat) becomes a time-periodic solution of the Boltzmann equation
for a dilute gas in a fixed harmonic potential
$U(\mathbf{r},t)=m\Omega^{2}r^{2}/2$. This class of solutions coincides with
Boltzmann’s breathing modes. A related question pertains to the type of static
confining potential for which a breathing solution can exist. It is addressed
in C.

When $\Omega$ is constant, the general solution of the Ermakov equation
(11lai) is given by

$\sigma(t)=\left[\tau-\sqrt{\tau^{2}-\frac{1}{\Omega^{2}}}\cos(2\Omega
t+\phi)\right]^{1/2}$ (11laz)

where $\phi$ and $\tau\geq 1/\Omega$ are constants. Substituting this result
into Eq. (11lar) for $\beta(t)$ reveals that the effective inverse temperature
$\beta(t)$ oscillates harmonically with frequency $2\Omega$. Moreover, taking
the first three derivatives of $\beta(t)$ with respect to time, we
straightforwardly obtain

$\frac{d^{3}}{dt^{3}}{\beta}(t)+4\Omega^{2}\dot{\beta}(t)=0$ (11lba)

which is equivalent (for static $\Omega$) to Eq. (9b) of Ref. [52].

###  2.5 Relations to quantum shortcuts to adiabaticity

The classical shortcuts described above have quantum counterparts as we now
briefly discuss.

We have already seen this correspondence in the global counterdiabatic case:
just as the classical term $H_{\rm GCD}(z,t)$ defined by Eqs. (11la) and
(11lb) generates the desired shortcut for any choice of initial energy shell
(Eq. (11ln)), so too its quantum counterpart $\widehat{H}_{\rm GCD}(t)$ given
by Eqs. (11a) and (11b) generates the desired shortcut for any energy
eigenstate [11, 13].

For local counterdiabatic driving, the auxiliary term $H_{\rm LCD}=pv(x,t)$ is
designed for a specific choice of initial energy shell, as discussed in Sec.
2.2.2. In the quantum case, let us instead choose an initial energy eigenstate
$|n(0)\rangle$ and construct the cumulative distribution

${\cal
F}(x,t)=\int_{-\infty}^{x}dx^{\prime}\,\left|\phi_{n}(x^{\prime},t)\right|^{2}\quad,\quad\phi_{n}(x,t)\equiv\langle
x|n(t)\rangle.$ (11lbb)

Now define a velocity field $v(x,t)$ in a manner analogous to Eq. (11ls), but
with ${\cal F}(x,t)$ playing the role of ${\cal S}(x,t)$:

$v(x,t)=-\frac{\partial_{t}{\cal F}}{\partial_{x}{\cal F}}\,.$ (11lbc)

Using this field we construct a quantum auxiliary term

$\widehat{H}_{\rm
LCD}(t)=\frac{\widehat{p}\widehat{v}+\widehat{v}\widehat{p}}{2}\quad,\quad\widehat{v}\equiv
v(\widehat{x},t)$ (11lbd)

(compare with Eq. (11lt)). If the system begins in the state $|n(0)\rangle$
and then evolves under $\widehat{H}(t)+\widehat{H}_{\rm LCD}(t)$, it will
remain in the instantaneous eigenstate $|n(t)\rangle$ (up to an overall time-
dependent phase) for all $t\in[0,t_{\scriptsize\textnormal{f}}]$ [16].

Using the velocity field $v(x,t)$ given by Eq. (11lbc), we next construct an
acceleration field $a(x,t)=v\partial_{x}v+\partial_{t}v$ and a corresponding
fast-forward potential $U_{\rm FF}(x,t)$ via $-\partial_{x}U_{\rm FF}=ma$ (see
Eqs. (11ly) and (11lz)). If the system evolves under
$\widehat{H}(t)+\widehat{U}_{\rm FF}(x,t)$, from an initial state
$|n(0)\rangle$, then at the final time $t=t_{\scriptsize\textnormal{f}}$ it
will arrive at the eigenstate $|n(t_{\scriptsize\textnormal{f}})\rangle$ (up
to a phase), though at intermediate times it generally will not be in the
state $|n(t)\rangle$ [16]. The potential $U_{\rm FF}(x,t)$ obtained in this
manner – using ${\cal F}$ rather than ${\cal S}$ to construct $v$ – is
equivalent to the fast-forward potential originally derived by Masuda and
Nakamura [14] and further studied in Refs. [15, 56, 57].

It is interesting to note that both quantum and classical auxiliary terms are
constructed using a velocity field determined from a cumulative function,
either ${\cal F}(x,t)$ or ${\cal S}(x,t)$. As we shall see later, this pattern
applies as well to stochastic shortcuts to adiabaticity, where the cumulative
function ${\cal F}(x,t)$ is given in terms of a canonical probability
distribution.

The discussion in the previous three paragraphs requires a rather strong
caveat. For excited states $n>0$, the field $v(x,t)$ given by Eq. (11lbc)
generically diverges at the nodes of the eigenstate, i.e. where
$\phi_{n}(x,t)=0$, leading to ill-behaved auxiliary terms $\widehat{H}_{\rm
LCD}$ and $\widehat{U}_{\rm FF}$. Thus in general the quantum LCD and FF
approaches described above are limited to ground states $n=0$, which have no
nodes; scale-invariant systems are an exception to this statement [27].

Ref. [58] develops a semiclassical fast-forward approach that avoids the
problem posed by the nodes of $\phi_{n}(x)$. In this approach, the fields $v$
and $a$ are constructed using ${\cal S}$ (Eq. (11ls)) rather than ${\cal F}$
(Eq. (11lbc)). The resulting auxiliary term $\widehat{U}_{\rm FF}$ is well-
behaved, but no longer guides the wavefunction exactly to the desired final
state $|n(t_{\scriptsize\textnormal{f}})\rangle$. Instead it provides an
approximate shortcut that is expected to work well in the semiclassical regime
of large $n$. Numerical simulations support this expectation [58].

##  3 Shortcuts for classical systems in contact with a thermal bath

In the previous section, we have obtained finite-time protocols for driving an
isolated system, with the same final state as with an infinitely long driving.
In the remainder, our interest goes to systems that are (strongly) coupled to
an environment; the systems we have in mind are epitomized by a Brownian
object such as a colloid in a fluid. The latter plays the role of a thermal
bath [59], which introduces fluctuations. Therefore, the system necessitates a
stochastic/probabilistic description, a new feature compared to the treatment
in Sec. 2. This goes with a point-of-view-change, from single trajectories to
probability distributions. From an experimental perspective, two routes can be
proposed, where the question, formulated in a probabilistic fashion, becomes
meaningful. a) One may be interested in repetitions of an experiment involving
a Brownian object, such as a macromolecule, or a nano device; statistics is
then gathered by automatizing the protocol, as e.g. in [60] or [61]; b) A
single experiment can allow for measuring a distribution function, provided it
involves a collection of Brownian objects, manipulated simultaneously. When
they do not interact, or when interactions are weak such as in a low density
colloidal system, routes a) and b) are equivalent.

In a nutshell, what we now aim at suitably driving with an external force is a
probability density, which shall be denoted $\rho(x,t)$: What should the drive
be to evolve the system from an initial density
$\rho(x,t_{\scriptsize\textnormal{i}})$ at time
$t_{\scriptsize\textnormal{i}}$ to a target distribution
$\rho(x,t_{\scriptsize\textnormal{f}})$ at time
$t_{\scriptsize\textnormal{f}}$? Such a question of fast driving a system from
a given equilibrium state to another appeared under different names in the
literature: Engineered Swift Equilibration [17],111Also variants such as
Engineered Swift Relaxation for the connection of non-equilibrium states [62],
see Sec. 5. or Shortcuts to Isothermality [18].222The latter work is akin to
the counterdiabatic method presented in Sec. 3.2. It has been investigated
experimentally for small systems (e.g. colloids [17] and AFM tips [63]) in
contact with a thermostat. There, one searches for the proper variation of the
control parameters to ensure a transfer from the initial state at time
$t_{\scriptsize\textnormal{i}}$ to the desired final state at time
$t_{\scriptsize\textnormal{f}}$.

###  3.1 Inverse Engineering

A method that often proves efficient in practice is the commonly-called
inverse engineering technique. This approach has been successfully applied to
ordinary differential equations in classical and quantum physics [7]. If we
denote by $\mathbf{X}$ the dynamical variables (state vector with $n$
components), the control problem is encapsulated in a set of coupled
differential equations

$\dot{\mathbf{X}}={\mathbf{f}}({\mathbf{X}},{\mathbf{\lambda}}(t)).$ (11lbe)

There is no general statement for the reachability of the desired state under
the driving provided by ${\mathbf{\lambda}}(t)$ (control vector with $r$
components) except when ${\mathbf{f}}$ is a linear function with time-
invariant coefficients, for which the Kalman rank condition applies [64,
19].333For a linear system of the type
${\mathbf{f}}({\mathbf{X}},{\mathbf{\lambda}(t)})=A{\mathbf{X}}(t)+B{\mathbf{\lambda}}(t)$,
where $A$ and $B$ are time-independent matrices of size $n\times n$ and
$n\times r$, respectively, we define the $n\times nr$ controllability matrix
$C=[B\;AB\;...\;A^{n-1}B]$. The system is controllable if $C$ has $n$ linearly
independent columns, i.e. rank$(C)=n$ [64, 19]. When the desired state belongs
to the reachable set of solutions, there exists a wide variety of possible
protocols to reach the final state.

The inverse engineering method provides a convenient way to find out a proper
driving. It amounts to impose the evolution of the dynamical variable
${\mathbf{X}}$ and to infer from Eq. (11lbe), the expression for
${\mathbf{\lambda}}(t)$. However, this inverse use of the differential
equation that governs the dynamics is not always easy to handle. The
mathematical property which allows for such an inverse use of the dynamical
(including nonlinear) equations is known as the flatness property, and can be
considered as an extension of the Kalman’s controllability criterion [65].
This strategy has been successfully used to transport a particle in a moving
harmonic potential, both in classical and quantum physics [66, 67, 68], or to
shuttle the particle, i.e. to set a given velocity to the particle—see [69]
and references therein.

For a system in contact with a thermal bath, the dynamical variables obey
stochastic differential equations (Langevin-like) that cannot be directly
written as a set of continuous equations such as that of Eqs. (11lbe). For
Brownian motion, this difficulty has been circumvented in two steps. First, as
emphasized above, the discussion is made not on the individual trajectories
but on the density distribution that obeys the Smoluchowski, or Fokker-Planck
(FP), equation. (See B for a crash recapitulation of the essential aspects of
the Langevin and frameworks for stochastic processes.) Second, an ansatz
depending on a set of a few effective dynamical variables, such as the moments
of the PDF, is usually proposed to get a finite set of equations in the form
of Eqs. (11lbe). The coupled equations on the moments can be alternatively
derived from the Langevin equation.

Consider the overdamped motion of a bead of micron size immersed in a fluid
and trapped by a general confining potential $U(x,t)$. The density $\rho(x,t)$
obeys the overdamped —or Smoluchowski—equation

$\gamma\partial_{t}\rho(x,t)=\partial_{x}\left[\partial_{x}U(x,t)\rho(x,t)\right]+\beta^{-1}\partial^{2}_{x}\rho(x,t),$
(11lbf)

where $\gamma$ is the friction coefficient and $\beta=(k_{B}T)^{-1}$ refers to
the inverse of the temperature, $k_{B}$ being Boltzmann’s constant—see B for a
detailed discussion of the stochastic dynamics framework. Particular attention
has been paid to the harmonic potential, which models paradigmatic systems
such as optical tweezers. Equation (11lbf) then becomes

$\gamma\partial_{t}\rho(x,t)=\partial_{x}\left[\kappa
x\rho(x,t)\right]+\beta^{-1}\partial^{2}_{x}\rho(x,t).$ (11lbg)

In the harmonic problem, the control (time-dependent) parameters are a priori
the temperature $T$ (see Sec. 3.4) and the stiffness $\kappa$ of the
potential—typically, the friction coefficient is assumed to be constant, see B
for details. Here, the statistics of the bead’s position remains Gaussian if
it is initially, with a standard deviation $\sigma(t)$ that plays the role of
an effective dynamical variable and obeys

$\dot{\sigma}=-\frac{\kappa}{\gamma}\sigma+\frac{k_{B}T}{\gamma}\frac{1}{\sigma}.$
(11lbh)

For both time-independent temperature and stiffness, the system approaches the
canonical equilibrium distribution, for which the variance of the position is

$\sigma_{eq}^{2}=\frac{k_{B}T}{\kappa},$ (11lbi)

i.e. the equilibrium equation of state.

Let us assume for instance that the temperature is constant and we vary the
stiffness of the trap. The initial state corresponds to thermal equilibrium
state with the initial stiffness
$\kappa_{\scriptsize\textnormal{i}}=\kappa(t_{\scriptsize\textnormal{i}})$.
The objective is to reach the equilibrium state with the desired final value
of $\kappa_{\scriptsize\textnormal{f}}=\kappa(t_{\scriptsize\textnormal{f}})$
in a chosen amount of time $t_{\scriptsize\textnormal{f}}$. Those boundaries
conditions define the values
$\sigma_{\scriptsize\textnormal{i}}=(k_{B}T/\kappa_{\scriptsize\textnormal{i}})^{1/2}$
and
$\sigma_{\scriptsize\textnormal{f}}=(k_{B}T/\kappa_{\scriptsize\textnormal{f}})^{1/2}$
at initial and final time. The inverse engineering technique involves the
choice of an interpolation function $\sigma(t)$ between those two values, and
to subsequently infer the time-dependent stiffness $\kappa(t)$ to be applied,
directly from Eq. (11lbh), since both $\sigma$ and $\dot{\sigma}$ are then
known. For the very same problem in the underdamped regime, one shall solve
the Kramers equation for the phase space distribution. In this latter case,
the effective dynamical system boils down to a set of 3 coupled linear
equations for the time evolution of the 3 moments $\langle x^{2}\rangle$,
$\langle v^{2}\rangle$, and $\langle xv\rangle$. As a result, the strategy to
extract $\kappa(t)$ is slightly more involved [70].

The inverse engineering method applied here is specific to the manipulation of
Gaussian states. As already commented, a remarkable property is that such an
initial condition, under time-dependent harmonic forcing, remains Gaussian at
all times, and thus preserves its shape.444This can be viewed as a consequence
of the conservation of the Gaussian character when summing Gaussian variables,
and dwells on the fact that the solution to Langevin equation is then linear
in the noise history. Yet, while it is difficult to compute analytically the
potential required to connect two arbitrary non-Gaussian states, it turns out
that if both have the same shape, a simple solution can be found, therefore
generalizing the aforementioned Gaussian result. We here impose the following
shape for the PDF $\rho(x,t)$:

$\rho(x,t)=\frac{1}{\sigma(t)Z_{\scriptsize\textnormal{i}}}\exp\left[-\beta
U_{\scriptsize\textnormal{i}}\left(\frac{x-\mu(t)}{\sigma(t)}\right)\right],$
(11lbj)

with $\mu(t)$ and $\sigma(t)$ continuous real functions, such that
$\mu(t_{\scriptsize\textnormal{i}})=0$ and
$\sigma(t_{\scriptsize\textnormal{i}})=1$, which accounts for shifting and re-
scaling the space dependence. Note that, although we used similar ideas in the
scale-invariant protocols introduced in Eq. (8)—justifying the same notation
for the parameters, processes preserving the shape are different from those
scale-invariant protocols.555Specifically, one can be mapped onto the other if
$U_{\scriptsize\textnormal{i}}(x)$ is proportional to a power of $x$. For
example, if we consider the harmonic potential
$U_{\scriptsize\textnormal{i}}(x)=k_{\scriptsize\textnormal{i}}x^{2}/2$, the
potential $U$ remains harmonic for both cases, shape-preserving and scale-
invariant, but with a different stiffness,
$k_{\scriptsize\textnormal{i}}/\sigma^{2}$ and
$k_{\scriptsize\textnormal{i}}/\sigma^{4}$ respectively. For general
$U_{\scriptsize\textnormal{i}}$, the map is not guaranteed. The quantity
$Z_{\scriptsize\textnormal{i}}=\int_{\cal D}dx\,e^{-\beta
U_{\scriptsize\textnormal{i}}(x)}$ is the partition function that guarantees
the correct normalization of the distribution. Such distributions only connect
states that belong to the same family of potentials, whatever this family is.
The final potential then reads [71]

$U_{\scriptsize\textnormal{f}}(x)=U_{\scriptsize\textnormal{i}}\left[\frac{x-\mu(t_{\scriptsize\textnormal{f}})}{\sigma(t_{\scriptsize\textnormal{f}})}\right].$
(11lbk)

The time-dependent driving potential required to ensure shape preservation is
found by introducing the ansatz (11lbj) into the equation. One finds [71]

$\displaystyle U(x,t)=$ $\displaystyle
U_{\scriptsize\textnormal{i}}\left[\frac{x-\mu(t)}{\sigma(t)}\right]$ (11lbl)
$\displaystyle-\gamma\frac{2\dot{\mu}(t)\sigma(t)[x-\mu(t)]+[x-\mu(t)]^{2}\dot{\sigma}(t)}{2\sigma(t)}.$

Such a driving potential involves two different contributions: the shape-
bearing potential itself and an additional harmonic potential, whose stiffness
and center are determined by certain combinations of the shift and scaling
functions $\mu(t)$ and $\sigma(t)$. The extra time-dependent harmonic
potential is nothing but the counterdiabatic term, to be discussed further
below. To enforce the smoothness of the potential, one can add the following
extra conditions on the parameters at initial and final time:
$\dot{\mu}(t_{\scriptsize\textnormal{i}})=\dot{\mu}(t_{\scriptsize\textnormal{f}})=0$
and
$\dot{\sigma}(t_{\scriptsize\textnormal{i}})=\dot{\sigma}(t_{\scriptsize\textnormal{f}})=0$.

Inverse engineering techniques have also been employed to address the
underdamped situation [70]. The authors worked out eligible conservative,
velocity-independent, drivings $U(x,t)$. However, the problem becomes involved
and limitations appear, presumably inherent to the functional forms chosen for
constructing explicit solutions.

Related in spirit to inverse engineering are stochastic methods that generate
Brownian paths conditioned to start and end at prescribed (ensemble of) points
[72]. The conditioning, that endows these paths with precise statistical
properties, can be of various types: a bridge, meaning a path that starts at
some $x_{\scriptsize\textnormal{i}}$ at $t=t_{\scriptsize\textnormal{i}}$, and
ends at a given $x_{\scriptsize\textnormal{f}}$ at
$t=t_{\scriptsize\textnormal{f}}$; an excursion, meaning a bridge with
$x_{\scriptsize\textnormal{f}}=x_{\scriptsize\textnormal{i}}$ that is
furthermore constrained to lie at all times to the right of
$x_{\scriptsize\textnormal{i}}$, etc. These methods can be used to generate
computationally the constrained paths in an efficient manner; a naive variant
would amount to pruning an ensemble of unconstrained paths, keeping only those
trajectories that fulfil the imposed constraints, a highly inefficient way of
proceeding. For overdamped Langevin dynamics, it was shown on general grounds
that irrespective of $x_{\scriptsize\textnormal{i}}$, addition of an external
harmonic force centered at $x_{\scriptsize\textnormal{f}}$ and with stiffness
$\gamma/(t_{\scriptsize\textnormal{f}}-t)$, which thus diverges for $t\to
t_{\scriptsize\textnormal{f}}$, generates a bona fide bridge [72]. More
generally, as might be anticipated from the previous example, the gist of the
approach is to add a time-dependent entropic potential of the form
$-2k_{B}T\log Q(x,t)$; $Q$, encoding the constraints considered (such as
remaining in the allowed half-space for an excursion), is the probability
density to be at $x_{\scriptsize\textnormal{f}}$ at time
$t_{\scriptsize\textnormal{f}}$, having started at point $x$ at time
$t<t_{\scriptsize\textnormal{f}}$ [72, 73]:

$Q(x,t)\equiv
P(x_{\scriptsize\textnormal{f}},t_{\scriptsize\textnormal{f}}|x,t).$ (11lbm)

While this result can be obtained by inverse engineering, it should be
stressed that it does feature an important difference with short-cutting ideas
discussed here, in the sense that there does not exist an infinitely slow
process that is being accelerated in some form. Note that $Q$ introduced above
fulfills the backwards equation, and is as such intimately related to first
passage problems [74].

###  3.2 Counterdiabatic method

In this section, we extend the counterdiabatic method (see Secs. 1 and 2) to
systems in contact with a thermal bath. As above, we model the evolution of
such systems with overdamped Langevin dynamics at the single-trajectory level,
and with the FP equation at the ensemble level, and we restrict ourselves to
systems with a single degree of freedom, $x$. Because momentum is ignored in
the overdamped limit, we will use the generic notation $U(x)$ rather than
$H(x,p)$ to denote the system’s energy function. This potential, used to drive
the system, can be viewed as a Hamiltonian.

In the isolated quantum and classical cases described in Secs. 1 and 2,
counterdiabatic driving aims to preserve an adiabatic invariant under rapid
driving. In the present context the role of the adiabatic invariant is played
by the functional form of the probability distribution function (PDF),
$\rho(x,t)$. Specifically, when the potential $U$ is driven very slowly, its
PDF evolves through a continuous sequence of equilibrium states
$\rho_{eq}(x,t)$,

$\displaystyle\rho_{eq}(x,t)=e^{\beta[F(t)-U(x,t)]},$ $\displaystyle
F(t)=-\beta^{-1}\ln\left[\int dx\,e^{-\beta U(x,t)}\right].$ (11lbn)

In the counterdiabatic method we seek to construct a term $U_{CD}(x,t)$ such
that under the full driving potential $U_{\rm SST}=U+U_{CD}$, the system
evolves through the same equilibrium states $\rho_{eq}(x,t)$, even when the
time dependence of $U(x,t)$ is not slow.

Note that, although the fine details may be different, the idea behind the
counterdiabatic method could be cast under the umbrella of inverse engineering
techniques. That is, we want to find an auxiliary $H_{CD}$ that enforces the
conservation of the adiabatic invariant. The separation into categories of the
techniques analyzed in this review has to be understood as a choice made more
by pedagogical purposes than by a rigorous categorization.

Consider an overdamped Brownian particle whose statistical state $\rho(x,t)$
obeys the FP equation (11lbf)—see B for details—with the potential $U_{\rm
SST}(x,t)=U(x,t)+U_{CD}(x,t)$. The instantaneous equilibrium distribution
associated with $U$, $\rho_{eq}$ given by Eq. (11lbn), plays the role of
adiabatic invariant. In other words, we aim to steer the system so that the
equilibrium PDF of the unperturbed potential $U(x,t)$ is maintained throughout
the process.

When $U(x,t)$ varies at finite rate, $U_{CD}\neq 0$ is needed to preserve the
prescribed evolution $\rho_{eq}$. Specifically, the counterdiabatic term
escorting the adiabatic evolution is given by [18]

$-\partial_{x}U_{CD}(x,t)=-\gamma\frac{\int^{x}dx^{\prime}\partial_{t}\rho_{eq}(x^{\prime},t)}{\rho_{eq}(x,t)}.$
(11lbo)

This expression can be derived from the cumulative distribution

$\mathcal{F}(x,t)=\int_{-\infty}^{x}dx^{\prime}\rho(x^{\prime},t)$ (11lbp)

and its velocity field

$v(x,t)=\partial_{t}x(\mathcal{F},t)=-\frac{\partial_{t}\mathcal{F}(x,t)}{\partial_{x}\mathcal{F}(x,t)},$
(11lbq)

where $x(\mathcal{F},t)$ is the function obtained by inverting
$\mathcal{F}(x,t)$ for fixed $t$. Specifically, using the velocity field, we
obtain

$-\partial_{x}U_{CD}(x,t)=\gamma v(x,t).$ (11lbr)

This approach is particularly appealing since it has been shown to be useful
as a unified procedure to derive counterdiabatic terms in quantum, classical
and thermal systems [16]—as in Eqs. (11ls) and (11lt).

The counterdiabatic force given in Eq. (11lbo) allows to recover the inverse
engineering results derived for shape preserving potentials in Sec. 3.1.
Indeed, injecting relation (11lbj) into (11lbo) (with $\rho$ playing the role
of $\rho_{eq}$) yields, rather unexpectedly, a counterdiabatic force that does
not depend on $U$ and thus holds irrespective of the functional form chosen
for the PDF [71], provided this form is conserved (and is thus a preserved
“adiabatic invariant”):

$\partial_{x}U_{CD}(x,t)\,=\,-\gamma\,\frac{\dot{\mu}(t)\sigma(t)+[x-\mu(t)]\dot{\sigma}(t)}{\sigma(t)},$
(11lbs)

which is consistent with the potential given in Eq. (11lbl), remembering
$U_{\rm SST}=U+U_{CD}$.

The counterdiabatic method has been also employed in the underdamped situation
[75, 18, 76]. At variance with the conservative drivings $U(x,t)$ derived with
the inverse engineering technique [70], the counterdiabatic potential contains
a term depending on the momentum of the particle. This entails serious
difficulties for the experimental implementation of such a counterdiabatic
driving. Counterdiabatic protocols have also been devised for systems
described by master equations, i.e. Markov processes with discrete states—see
Sec. 6.2. Finally, we present in F a derivation of the work fluctuation
relation [77] from counterdiabatic type of arguments.

###  3.3 Fast-forward

Now, let us consider another strategy for the swift connection of equilibrium
states: the so-called fast-forward procedure. In the quantum case [78, 14],
fast-forward refers to a protocol that makes it possible to reach a desired
final state, independently of the path swept to do so. For the isolated
classical systems analyzed in Sec. 2 of this review, fast-forward has been
employed in a similar sense. Therein, the fast-forward protocol led the system
to the target state, with the same value of the adiabatic invariant as the
initial state, by adding a velocity-independent potential $U_{FF}(x,t)$. The
price of such a procedure was the adiabatic invariant no longer preserved at
intermediate times. Derivation of fast-forward protocols are unified in both
quantum and classical mechanics by using the acceleration flow field [16], as
discussed in Sec. 2.2.3.

Here, we consider the extension of the fast-forward idea to the context of
systems with stochastic dynamics, again described by the overdamped FP
equation (11lbf). In contrast to the counterdiabatic driving just described,
there is no underlying adiabatic transformation over which the shortcut is
built. The idea is the following: one considers a certain reference
process—not necessarily slow—that connects two given states and then searches
for a tailor-made external potential that accelerates this reference process
[71]. In this way, the “frames” of the “movie” are fixed, given by the
reference process, but are played at a higher rate in the fast-forward
protocol. As explained below, the same idea allows for reproducing the frames
at a lower rate (slow-forward) or even play the movie backwards, at a higher
or a lower rate (fast-backward or slow-backward, respectively)—somehow
generalizing the shortcuts described in this review.

Going into specifics, and following Ref. [71], we consider a reference
solution $\rho_{r}(x,t)$ of Eq. (11lbf) under a reference potential
$U_{r}(x,t)$. Then, we introduce a time distortion $\Lambda(t)$ of the
reference, i.e., $\rho(x,t)=\rho_{r}(x,\Lambda(t))$, and look for the
potential $U(x,t)$ required to drive the evolution following the given
prescription $\rho_{r}(x,\Lambda(t))$. It is important to remark that, at
variance with the classical case, the connecting path is fixed and given by
$\rho_{r}(x,\Lambda(t))$.

The FP equation (11lbf) is rewritten as a continuity equation

$\displaystyle\partial_{t}\rho=-\partial_{x}\left(\rho v\right),$ (11lbta)
$\displaystyle
v(x,t)=-\gamma^{-1}\left[\partial_{x}U(x,t)+\beta^{-1}\partial_{x}\ln\rho(x,t)\right]$
(11lbtb)

where $v(x,t)$ is a velocity field.666This velocity field is exactly the same
defined in the counterdiabatic driving, since spatial integration of the FP
equation, written in the form (11lbta), from $-\infty$ to $x$ leads precisely
to Eq. (11lbq). The derivation of the fast-forward protocol is based on the
relation between the velocity fields in the reference and the manipulated
process,

$v(x,t)=\dot{\Lambda}(t)v_{r}(x,\Lambda(t)).$ (11lbtbu)

Hence, the driving potential can be solved. For the sake of clarity, we
display the solution for the auxiliary potential,

$\displaystyle-\partial_{x}U_{FF}(x,t)=$ $\displaystyle[1-\dot{\Lambda}(t)]$
(11lbtbv)
$\displaystyle\times[\partial_{x}U_{r}(x,\Lambda(t))+\beta^{-1}\partial_{x}\ln\rho_{r}(x,\Lambda(t))].$

The total driving potential is $U(x,t)=U_{r}(x,\Lambda(t))+U_{FF}(x,t)$. The
solution is given in terms of the reference process and the time map
$\Lambda(t)$. As expected, when the reference and manipulated dynamics
coincide, i.e., $\Lambda(t)=t$, one simply has
$\partial_{x}U(x,t)=\partial_{x}U_{r}(x,t)$ or, equivalently,
$\partial_{x}U_{FF}(x,t)=0$.

This type of driving not only allows for acceleration ($\dot{\Lambda}>1$) but
also for deceleration ($0<\dot{\Lambda}<1$) and even for the inversion of
time’s arrow ($\dot{\Lambda}<0$, meaning that the reference dynamics can be
“played backwards”). Combining simple reference processes, it is possible to
build up an operational welding protocol that connects arbitrary states [71].
Specifically, one can produce a welding connection between an initial state
$\rho_{\scriptsize\textnormal{i}}(x)$ and the target state
$\rho_{\scriptsize\textnormal{f}}(x)$ through an intermediate state
$\rho_{\scriptsize\textnormal{int}}(x)$, distorting two consecutive reference
processes. This construction relies on the acceleration of a first reference
relaxation process from $\rho_{\scriptsize\textnormal{i}}(x)$ to
$\rho_{\scriptsize\textnormal{int}}(x)$; and accelerating and reversing a
second relaxation process from $\rho_{\scriptsize\textnormal{f}}(x)$ to
$\rho_{\scriptsize\textnormal{int}}(x)$—see Fig. 3. If one seeks an operating
time $t_{\scriptsize\textnormal{f}}$ (with a starting time
$t_{\scriptsize\textnormal{i}}=0$), one can assign a time duration
$t_{\scriptsize\textnormal{f}}/2$ for each of the two steps, but other choices
are possible.

Figure 3: Sketch of the welding strategy to connect two arbitrary
distributions $\rho_{\scriptsize\textnormal{i}}$ and
$\rho_{\scriptsize\textnormal{f}}$. The connection is made in two steps. In
the first step, a reference process from $\rho_{\scriptsize\textnormal{i}}$ to
$\rho_{\scriptsize\textnormal{int}}$ is submitted to a fast-forward evolution.
In the second step, a reference process starting at
$\rho_{\scriptsize\textnormal{f}}$ and finishing at
$\rho_{\scriptsize\textnormal{int}}$ is time-reversed and accelerated (fast
backward). In doing so, one achieves the desired connection in a chosen time.

Note that the counterdiabatic method could be understood as a limit of the
fast-forward protocol presented above. Let us consider a reference process of
duration $t_{r}$ and take the limit $t_{r}\to\infty$, so that the reference
process becomes quasistatic. Hence, the limit reference process would sweep
equilibrium states. In this way, the resulting limit of the fast-forward
protocol would converge to the counteradiabatic method.

###  3.4 Engineering the thermal environment

In previous developments, we addressed Brownian objects in some environment at
thermal equilibrium, meaning the temperature entering the Langevin and FP
equations is fixed. For colloidal systems, the environment is usually water,
and it may be difficult to impose a chosen time dependence for its
temperature, especially if a massive heating is sought [79]. Yet, it is
possible to “fool” the colloidal beads, subjecting them to a random forcing
that will emulate an effective temperature exceeding several thousand kelvins
[79, 59]. The method is quite robust, and essentially requires that the
forcing frequency be large compared to the bead inverse relaxation time. It is
then possible to finely control the time dependence of the effective
temperature, by playing on the amplitude of the forcing, which opens new means
for driving the system [70], and paves the way towards the more general goal
of reservoir engineering.

Following this idea, a micrometric silica sphere has been driven in [70] by
the joint monitoring of a harmonic trap stiffness ($\kappa$), and the point of
zero force ($x_{0}$). In other words, the confining potential is of the form

$U(x,t)\,=\,\frac{1}{2}\kappa(t)\left(x-x_{0}(t)\right)^{2},$ (11lbtbw)

and the idea is to impose the proper time dependence jointly on both $\kappa$
and $x_{0}$. Compared to the more usual situation where $x_{0}$ is fixed, a
new contribution $\kappa(t)x_{0}(t)$ arises in the force balance. It is
important that a) $x_{0}$ remain small compared to the bead size, in order not
to affect the effective stiffness and b) that the correlation time of the
signal $x_{0}(t)$ be small compared to the protocol duration (itself by
construction smaller than the intrinsic relaxation time). Then, the forcing of
$x_{0}$ results in an effective heat bath for the colloidal degree of freedom:
this forcing has a time-dependent amplitude and, for practical purposes, can
be viewed as delta-correlated in time. In Ref. [70], the control of $x_{0}$
was achieved with an acousto-optic deflector, and the bath engineering made it
possible to quickly deconfine a colloidal state. Should one be able to play
only on the stiffness $\kappa(t)$, equivalent transformations would require
transiently expulsive forces with $\kappa<0$, which represent an experimental
challenge [80, 81]. A limitation of the approach is that the extra stochastic
forcing applied results in enhanced Brownian fluctuations, and in an effective
temperature increase. Other techniques would have to be applied when it comes
to cooling the center-of-mass motion of trapped beads, such as feedback-based
approaches [82].

Recently, the effective heating of optically trapped object allowed to devise
finite-time adiabatic processes [83]. Here, we stress that “adiabatic” is
understood in its usual thermodynamics meaning, of heat-exchange free [3]. A
Brownian object is inherently fluctuating, and strongly coupled to its
environment. Consider the compression at fixed temperature $T$ of a colloidal
bead, where for instance the stiffness $\kappa$ of a harmonic potential is
increased. When the bead has relaxed, the internal energy difference $\Delta
U$ vanishes on average between the initial and final states, which have the
same temperature. The first principle of thermodynamics [6] states that, since
the bead received work from the confining force, heat flew towards the bath,
on average [4]. If one seeks a vanishing heat exchange on average, it is
mandatory that the environment temperature increase. In the quasi-static
limit, this increase has to be proportional to $\sqrt{\kappa}$. This can be
seen as a consequence of Laplace relation between temperature $T$ and volume
$V$ for the adiabatic reversible transformation of an ideal gas:
$TV^{2/3}=\textnormal{const}$ for a monoatomic gas. Here, the confinement
length is $\sigma_{eq}$, so that the role of the volume is played by
$\sigma_{eq}^{3}$. Since $\sigma^{2}_{eq}=k_{B}T/\kappa$, this yields a
Laplace condition $T\sigma_{eq}^{2}=\textnormal{const}$ or, equivalently,
$T^{2}/\kappa=\textnormal{const}$. This can be viewed as the statement that
$n\Lambda^{3}=\textnormal{const}$, where $n$ is the typical density and
$\Lambda$ is De Broglie wavelength, which guarantees that volume in phase
space is conserved. This volume is computed from the (cubed) product of the
typical length in real space, $\sigma$, times the typical velocity, scaling
like $\sqrt{T}$. This discussion also illustrates that the quasi-static
criterion $T^{2}/\kappa=\textnormal{const}$ is space-dimension independent
[84, 85].

For finite-time adiabatic processes, not only does the mean heat released to
the thermal environment vanish between the initial and final states of the
transformation, but it also does at any time in between. The operating time of
these irreversible adiabats can be optimized by jointly controlling the
potential and the temperature. The condition of zero heat involves the kinetic
contribution to the energy: the only assumption being, consistently with the
overdamped description, that the velocity degree of freedom is always at
equilibrium with the time-dependent value of the temperature. Some general
results emerge [83], like (i) forbidden regions, i.e., final states that
cannot be reached adiabatically, and (ii) a speed limit, the existence of a
minimum, in general nonvanishing, time $t_{\scriptsize\textnormal{f}}^{*}$ for
the adiabatic connection. For the specific case of a harmonic confining
potential, it was shown that $T\sigma^{2}$ is nondecreasing over the adiabats.
This implies that
$T_{\scriptsize\textnormal{f}}/T_{\scriptsize\textnormal{i}}\geq\sqrt{\kappa_{\scriptsize\textnormal{f}}/\kappa_{\scriptsize\textnormal{i}}}$,
taking into account that the system is at equilibrium at the initial and final
times—see Eq. (11lbi). It is only in the quasistatic limit that $T\sigma^{2}$
remains constant and, moreover, the equilibrium equation of state (11lbi)
holds for all times, which leads to recover Bo and Celani’s result of constant
$T^{2}/\kappa$ [84]. Advantage was taken of these finite-time adiabatic
transformations to construct an irreversible Carnot engine featuring
interesting efficiency properties [86], see also Sec. 6.1. In Fig. 4, the
fastest possible adiabatic connection is illustrated for a $20\%$ compression,
$\sigma_{\scriptsize\textnormal{f}}=0.8\sigma_{i}$.

Figure 4: Fastest control and evolution for an adiabatic (in the sense of zero
average heat) $20\%$ compression of a harmonically trapped particle. Time
$\tau=t/t_{\scriptsize\textnormal{f}}^{*}$ has been made dimensionless with
the shortest possible duration of the process
$t_{\scriptsize\textnormal{f}}^{*}$, which is reached for a linear evolution
of $\sigma_{x}^{2}$. The fastest connection requires the stiffness to be
discontinuous at the initial and final times. The example shown corresponds to
$T_{\scriptsize\textnormal{f}}/T_{\scriptsize\textnormal{i}}=16$,
$\kappa_{\scriptsize\textnormal{f}}/\kappa_{\scriptsize\textnormal{i}}=25$ and
entails an acceleration of a factor around $5.7$ with respect to the
relaxation time scale
$t_{\scriptsize\textnormal{rel}}=\gamma/\kappa_{\scriptsize\textnormal{f}}$.
See Sec. 4.2 and Ref. [83] for further details.

##  4 Optimal control theory

Hitherto, this review has focused on techniques that allow for connecting two
given states. Naturally, having established that going from A to B in a finite
time is feasible raises the question of what the best path is. Answering this
kind of question is the main goal of optimal control theory [19, 20], which
combines well with inverse engineering problems. However, the choice of the
best pathway depends on the quantity to be optimized (time, some kind of cost
function…) that should thus be defined in the first place. (See for instance
Refs. [68, 62] for specific examples, ranging from ultracold atoms to granular
systems.)

###  4.1 Minimization of the mean work

####  4.1.1 Harmonic connections

Although optimal connection problems in the context of finite-time
thermodynamics date back to the 1970’s and 1980’s [87, 88, 89], the first
solution of an optimal connection problem in the “modern” context of
stochastic thermodynamics and SST is due to Schmiedl and Seifert [90]. They
considered an overdamped Brownian particle submitted to harmonic trapping,
where either the position or the stiffness of the trap is controlled. We
recall that experiments with colloidal particles are usually well described by
the overdamped FP, or Smoluchowski, equation (11lbf) [59]—see also B.

In Ref. [90], the optimal control needed to minimize the mean work for an
isothermal process

$\displaystyle\left\langle W\right\rangle$
$\displaystyle=\int_{0}^{t_{\scriptsize\textnormal{f}}}\\!\\!dt\int_{-\infty}^{+\infty}\\!\\!\\!dx\,\partial_{t}U(x,t)\,\rho(x,t)$
(11lbtbx)
$\displaystyle=\int_{0}^{t_{\scriptsize\textnormal{f}}}\\!\\!dt\int_{-\infty}^{+\infty}\\!\\!\\!dx\,\dot{\boldsymbol{\lambda}}(t)\cdot\partial_{\boldsymbol{\lambda}}U(x,\boldsymbol{\lambda}(t))\,\rho(x,t).$

was derived. Note that we are assuming that the potential $U$ depends on time
through some externally controlled parameters $\boldsymbol{\lambda}$.777See B
for a brief account of the definitions of work and heat in the context of
stochastic thermodynamics. Dealing with the variation of the work, an Euler-
Lagrange equation for the control parameter was obtained and analytically
solved. The assumption of harmonic potential, and therefore Gaussian states,
has much to do with the fact that the problem is analytically solvable. The
dynamics of the system, in principle codified in the FP equation, can be
simplified to an ordinary differential equation for the only relevant moment
of the distribution—e.g. its standard deviation, see Eq. (11lbh).

It is natural to try to transpose optimal protocols to actual experiments, for
example the compression or decompression of a harmonically trapped Brownian
particle—by controlling the stiffness of the trap. Still, the experimental
implementation of the optimal protocols presents some difficulties.
Specifically, the negative values of the stiffness needed for
decompression—for short enough connecting times 888From Eq. (11lbh), it is
seen that for a fast decompression, where $\dot{\sigma}$ is strongly positive
in some time window, then $\kappa\simeq-\gamma\dot{\sigma}/\sigma$, and is
thus negative in the same time window.—are experimentally challenging. A step
forward to solve this issue has been made by employing an optical feedback
trap [80]. Although the connection considered therein is not optimal, it is
neatly shown that it is possible to decompress the Brownian particle in a
finite-time with a potential that becomes repulsive—i.e., with negative
stiffness—inside a certain time window. Also, the optimal control—the
stiffness of the trap—possesses finite discontinuities at the initial and
final times in the overdamped limit. These discontinuities have the same
formal origin as in the classic problems of finite-time thermodynamics [88],
the linearity of the Lagrangian in its highest derivative [91].

Discontinuities in the “control functions” are better rationalized in the
context of Pontryagin’s maximum principle of optimal control theory than
within the framework of variational calculus [19, 20]. In optimal control
theory, the control function only has to be piecewise continuous, and thus
discontinuities in the control like those appearing in the stiffness of the
harmonic trap are treated in a natural and mathematical rigorous way, see e.g.
Ref. [92]. More general potentials, beyond the harmonic case, represent a
challenge because solving the FP equation cannot be mapped onto solving an
ordinary differential equation. However, numerical minimization of the mean
work has been carried out, which shows that the predicted discontinuities are
robust features of the optimal control [93]. In the underdamped case, the
discontinuities of the control become harsher, they do not involve finite
jumps but delta peaks [94].

Optimal harmonic connections considering the limitation stemming from bounded
stiffness have also been investigated in the overdamped case [92].
Specifically, the stiffness has been assumed to be bounded between 0 and a
maximum value, $0<\kappa(t)<\kappa_{\max}$. Pontryagin’s principle provides
the adequate framework to solve such constrained optimal problem. The time
evolution of the stiffness turns out to be built by two pieces. In the first
one, the equations for the protocol are similar to those coming out from the
unconstrained problem [90] whereas, in the second piece, the control is kept
fixed and equal to one of its limiting values. These two pieces smoothly
match, in the sense that the dynamical variable—the variance of the position
of the Brownian particle—is continuous and has continuous time derivative. As
a consequence of the bounds, and depending on the target value of the
stiffness and the desired connection time, the target state may become
inaccessible. When the connection is possible, the minimum work is greater
than the minimum one found for the unconstrained case, the difference between
them becomes large in some situations. Similarly to the unconstrained case,
the control develops finite jumps at the initial and final times.

####  4.1.2 Beyond the harmonic case

After the problem of optimally connecting Gaussian states—in the sense of
minimizing the mean work—our interest goes to the analytical derivation of the
optimal connection for an arbitrary nonlinear potential—still in the
isothermal case. This challenging problem has been first addressed in a series
of related papers [95, 96, 97, 98, 99, 100]. Starting from the Langevin
description, instead of the equivalent FP equation for the PDF, they show that
the above question has quite a general answer in the overdamped limit.
Building on long-established relations between dynamical systems and
stochastic control theory [101, 102], the minimization of the average heat
released to the reservoir can be mapped onto an optimal mass transport
problem, ruled by the Burgers equation, and explicitly solved for several
physical situations [95, 96]. Since optimizing average heat and work is the
same problem,111Note that the first principle states that $\left\langle
Q\right\rangle+\left\langle W\right\rangle=\Delta U$. As the initial and final
states of the SST are fixed, $\Delta U$ is also fixed and minimizing
$\left\langle W\right\rangle$ thus entails minimizing $-\left\langle
Q\right\rangle$, i.e. the average heat released to the thermal bath. for the
particular case of a harmonic trap the optimal protocol that minimizes the
mean work in Refs. [90, 103] is recovered. Once more, the optimal control
presents discontinuities at the initial and final times.222The authors
addressed the regularization of these discontinuities by introducing a penalty
for the current acceleration [97]. A different “surgery” procedure, which also
introduces boundary layers, has been considered to avoid these jumps [92]. It
is worth stressing that the minimization of heat is shown to be equivalent to
the minimization of entropy production, providing a refined version of the
second law of thermodynamics [96]. In general, this line of research evidences
that differential geometry concepts can also be useful to investigate
optimization problems in stochastic thermodynamics [99, 104].

The generalization of the above results to Markov jump processes, governed by
a master rather than a FP equation, has also been carried out [98]. In the
continuum limit, the results converge to those previously described. It has
been shown that optimal protocols—in the same sense of minimizing entropy
production—for systems with discrete states may involve non-conservative
forces [105]. Here, non-conservative means that the logarithm of the ratio of
forward and backward rates cannot be written as a difference of state
functions. Also, the relevance of information geometry concepts like the
Wasserstein distance has been unveiled in this context [106]. Interestingly,
the Wasserstein distance is closely related to entropy production on quite
general grounds, also for continuous states [107, 108].

Furthermore, the possibility of extending the above results to the
underdamped, Langevin-Kramers, case has been investigated [100]. Therein, the
emergence of singularities and also of momentum dependence in the optimal
driving potential makes the situation less clear-cut than that found in the
overdamped situation. Recently, the minimization of the work for the
particular case of a counterdiabatic connection in the underdamped regime has
been investigated [76], but the driving is once more velocity-dependent. In
fact, the derivation of the optimal conservative driving potential $U(x,t)$
that provides the minimum irreversible work in the underdamped situation is,
to the best of our knowledge, an open question. An advance in this direction
has been accomplished by working out a momentum-independent protocol that
approximates the momentum-dependent optimal protocol for counterdiabatic
driving [109]. Therein, the momentum-dependent terms of the auxiliary
Hamiltonian are removed by combining a variational method with a gauge
transformation.

In the overdamped limit, the general problem for the optimization of the mean
work performed during the connection between arbitrary states has also been
worked in the FP framework [110, 111]. Specifically, the starting point is the
FP equation for the probability density $\rho(x,t)$, written as a continuity
equation, Eq. (11lbta). Still, the main role is played by the cumulative
distribution introduced in Eq. (11lbp). Making use of the method of
characteristics, the general results for the optimal connecting potential
obtained from the Langevin equation [95] are recovered. Here, we give the main
results for deriving the minimum work and the associated optimal protocol—for
a more detailed derivation thereof, see D.

For a quasi-static process in which the system remains at equilibrium for all
times, the average work $\left\langle W\right\rangle$ equals the free energy
difference $\Delta F$ between the final and initial states. For a finite-time
process, the second principle implies $\left\langle W\right\rangle>\Delta F$
and the irreversible (or excess) work is defined as

$W_{\scriptsize\textnormal{irr}}\equiv\left\langle W\right\rangle-\Delta F\geq
0.$ (11lbtby)

Starting from Eq. (11lbtbx) for the mean work, repeated use of integration by
parts and the FP equation leads to

$\displaystyle
W_{\scriptsize\textnormal{irr}}=\int_{0}^{t_{\scriptsize\textnormal{f}}}dt\,P_{\scriptsize\textnormal{irr}}(t),$
(11lbtbza) $\displaystyle
P_{\scriptsize\textnormal{irr}}(t)=\gamma\int_{-\infty}^{+\infty}\\!\\!\\!dx\,v^{2}(x,t)\rho(x,t)\geq
0,$ (11lbtbzb)

i.e. $P_{\scriptsize\textnormal{irr}}$ stands for the irreversible power in
the considered finite-time process and $v$ is given in (11lbtb). Note that,
consistently with our discussion above, $W_{\scriptsize\textnormal{irr}}$ (or
$P_{\scriptsize\textnormal{irr}}$) vanishes for a reversible process only:
$v(x,t)$ then identically equals zero—and thus the PDF has the equilibrium
shape $\rho(x,t)\propto e^{-\beta U(x,t)}$—for all $(x,t)$.

The problem of minimizing the average work $\left\langle W\right\rangle$ is
then equivalent to finding the fields $\rho(x,t)$ and $v(x,t)$ such that
$W_{\scriptsize\textnormal{irr}}$ becomes minimum, while verifying the FP
equation (11lbta). This can be done by employing the method of Lagrange
multipliers, introducing an auxiliary field $\psi(x,t)$ and seeking the
unconstrained minimum of

$J[\rho,v,\psi]\equiv\int_{0}^{t_{\scriptsize\textnormal{f}}}\\!dt\\!\int_{-\infty}^{+\infty}\\!\\!\\!dx\,L,\;L=v^{2}\rho+\psi\left[\partial_{t}\rho+\partial_{x}(\rho
v)\right].$ (11lbtbzca)

The Euler-Lagrange equations for this problem read

$\partial_{t}\psi+v\partial_{x}\psi=v^{2},\quad\partial_{x}\psi=2v,$
(11lbtbzcb)

plus Eq. (11lbta). Combining them, one gets the closed equation
$\partial_{t}\psi+\frac{1}{4}\left(\partial_{x}\psi\right)^{2}=0$, which is
nothing but the Burgers equation for the auxiliary field $\psi(x,t)$ [95].
Equivalently, one may write the Burgers equation in terms of the velocity
field $v(x,t)$,

$\partial_{t}v+v\partial_{x}v=0.$ (11lbtbzcc)

In order to find the optimal profiles for $\rho(x,t)$ and $v(x,t)$—note that
the latter provides us with the driving potential $U(x,t)$, the FP equation
(11lbta) and the Burgers equation (11lbtbzcc) must be simultaneously solved.
The solution is

$\displaystyle v=\varphi(x-vt),$ (11lbtbzcda)
$\displaystyle\rho(x,t)=\frac{\rho_{\scriptsize\textnormal{i}}(x-vt)}{1+t\varphi^{\prime}(x-vt)},$
(11lbtbzcdb)

where

$\rho_{\scriptsize\textnormal{i}}(x)\equiv\rho(x,t=0)$ (11lbtbzcdce)

is the initial distribution and $\varphi$ is—for the time being—an arbitrary
function.

The function $\varphi(x)$ is determined by the system reaching the target
distribution at the final time $t=t_{\scriptsize\textnormal{f}}$, i.e.
$\rho(x,t_{\scriptsize\textnormal{f}})=\rho_{\scriptsize\textnormal{f}}(x)$.
This condition is easier to implement by employing the cumulative distribution
introduced in Eq. (11lbp). Specifically, the initial and final cumulative
distributions
$\mathcal{F}_{\scriptsize\textnormal{i},\scriptsize\textnormal{f}}$,
corresponding to the initial and final probability distributions
$\rho_{\scriptsize\textnormal{i},\scriptsize\textnormal{f}}$, are needed.
After defining

$\Xi(x)\equiv\mathcal{F}_{\scriptsize\textnormal{i}}^{-1}\left(\mathcal{F}_{\scriptsize\textnormal{f}}(x)\right),$
(11lbtbzcdcf)

one gets

$\varphi(x)=\frac{\Xi^{-1}(x)-x}{t_{\scriptsize\textnormal{f}}},$
(11lbtbzcdcg)

where $A^{-1}(x)$ stands for the inverse function of $A$, i.e.
$A^{-1}(A(x))=x$. Note that finding $\varphi$ makes it possible to obtain the
driving potential, making use of Eq. (11lbtbzcda) and the definition of the
velocity field $v(x,t)$, Eq. (11lbtb),

$\partial_{x}U=-\gamma v-k_{B}T\partial_{x}\ln\rho.$ (11lbtbzcdch)

The irreversible power over the optimal protocol is shown to be

$\displaystyle
P_{\scriptsize\textnormal{irr}}^{*}\,=\,\frac{\gamma}{t_{\scriptsize\textnormal{f}}^{2}}\int_{-\infty}^{+\infty}\\!\\!\\!dx\,\rho_{\scriptsize\textnormal{i}}(x)\left[\Xi^{-1}(x)-x\right]^{2}.$
(11lbtbzcdci)

Note that, as emphasized by our notation, the optimal irreversible power does
not depend on time; it is a constant proportional to
$t_{\scriptsize\textnormal{f}}^{-2}$. The irreversible work immediately
follows,

$\displaystyle
W_{\scriptsize\textnormal{irr}}^{*}=t_{\scriptsize\textnormal{f}}P_{\scriptsize\textnormal{irr}}^{*}=\frac{\gamma}{t_{\scriptsize\textnormal{f}}}\int_{-\infty}^{+\infty}\\!\\!\\!dx\,\rho_{\scriptsize\textnormal{i}}(x)\left[\Xi^{-1}(x)-x\right]^{2},$
(11lbtbzcdcj)

which is then proportional to $t_{\scriptsize\textnormal{f}}^{-1}$. Quite
expectedly, it vanishes in the limit $t_{\scriptsize\textnormal{f}}\to\infty$,
where the optimal process tends to be quasi-static.

In general, the optimal potential stemming from the optimal velocity field, by
employing Eq. (11lbtb), has discontinuities at both the initial and final
times. For the harmonic case, this approach leads to the optimal connections
that can be worked out by using simpler methods, as described in Sec. 4.1.2.
For the general non-harmonic case, finding the optimal driving potential in an
explicit closed form can only be done in a few examples [110, 111]. The main
difficulty stems from the expression for the mapping $\Xi(x)$ defined in Eq.
(11lbtbzcdcf): only in simple cases is it possible to calculate explicitly the
cumulative distributions $\mathcal{F}_{\scriptsize\textnormal{i}}$ and
$\mathcal{F}_{\scriptsize\textnormal{f}}$, and to invert
$\mathcal{F}_{\scriptsize\textnormal{i}}$ and compose it with
$\mathcal{F}_{\scriptsize\textnormal{f}}$.

The above difficulties limit the usefulness of the exact optimal protocol for
practical implementations. Non-optimal driving potentials but with values of
$W_{\scriptsize\textnormal{irr}}$ close to the optimal one, which can be
expressed in closed form and do not present discontinuities, can be derived
with the welding procedure introduced in Ref. [71]—see also Sec. 3.3. Therein,
the authors look for the (sub)optimal connection belonging to fast-forward
protocols that leads to minimal work. Remarkably, such a protocol conserves
the property of the global optimum of delivering work at constant power.

Linear response theory [112] is a standard framework for understanding
nonequilibrium fluctuations. In fact, the problem above—i.e. optimization of
the mean work performed during the connection between arbitrary states—has
also been addressed using linear response theory [104, 113]. Sivak and Crooks
obtained the optimal protocol minimizing the irreversible work in the linear
regime, making use of information geometry concepts—see Sec. 4.3 for further
details. Once more, this approach leads to an optimal work with constant power
[114]. However, the assumption of linear response prevents this analysis from
capturing the discontinuities of the optimal potential at the initial and
final times. Bonança and Deffner deepened in the linear response approach,
thoroughly discussing its range of validity and comparing its predictions with
exact results [113].

###  4.2 Optimization of other figures of merit

In the previous subsection, we have gone over the literature related to the
optimization of average work. Nonetheless, the optimal approach in the context
of SST is not limited to the minimization of the average work. Herein, we
present different studies where the optimization of other relevant quantities
has been examined.

The minimization of the statistical error of the free energy—an estimator of
the difference between the free energy obtained after averaging over a certain
number of individual measurements and the real free energy change—has been
addressed [115]. Studying this quantity makes it possible to estimate the
number of experiments needed to attain a certain accuracy when performing free
energy measurements. Analytical results are not available even for simple
(harmonic) cases, but numerical optimization provides step-like protocols with
a significant reduction of the statistical error. Interestingly, this problem
is intimately related to the reverse process [116, 117].

Minimization of the average work carries no insight on work fluctuations.
Hence, interest aroused in the analysis of alternative figure of merits,
combining information of both average and fluctuations. Solon and Horowitz
studied the minimization of an objective function, which is a linear
combination of the work average $\langle W\rangle$ and the work standard
deviation $\sigma_{W}$, specifically $\mathcal{J}_{\alpha}=\alpha\langle
W\rangle+(1-\alpha)\sigma_{W}$, [118]. By varying the coefficient $\alpha$
from $0$ to $1$, the weight of work fluctuations (mean work) is reduced
(increased) in the search of the optimum protocol. The notion of Pareto-
optimal solutions is applied to classify all possible optimal protocols. Due
to the mathematical complexity of the optimization problem, this is carried
out numerically by tuning the value of the control in a finite set of times. A
first-order phase transition is found when illustrating the optimization
procedure above in a quantum dot.

The physical approach to information and memory storage [119, 120, 121] may
also give rise to optimization problems in the general framework described in
this review. Landauer’s principle states the minimum heat cost exhausted to
erase one random bit stored in a memory device [122, 123]. Such a bound refers
to a quasistatic process. Finite-time processes have been considered as well,
posing new relevant problems, not only from the theoretical point of view but
also from the technological one, where fast computational operations with
memory devices are necessary [124, 125, 126, 127, 128, 129, 130].

All previous instances belong in the optimization of energetic observables.
Remarkably, optimal problems involving time have also been investigated. The
minimization of time related observables is closely linked to information
geometry concepts and the so-called speed limits, which are discussed in Sec.
4.3. Below, we briefly report some results that have been obtained outside the
general framework of information geometry ideas.

The optimal static external potential required to minimize the escape time
$\tau$, of a Brownian particle confined in a box has been investigated [131,
132]. The escape time is defined here as the mean first passage time to the
end of the box, and it is considered that the external potential does not
introduce any bias between the starting and final points. Rather surprisingly,
in the overdamped regime, the escape time can approach zero arbitrarily close
[131], which requires divergent and strongly “squeezing” potentials [131,
132]. However, when some constraints are considered (e.g., the maximal
potential difference is below a certain threshold $\Delta U$) an expression,
which reminds Heisenberg’s time-energy uncertainty principle, is obtained
$\tau_{opt}\Delta U=\textnormal{const}$. Related to this, the optimization
problem for static external potential minimizing the first passage time to a
certain target, distributed according to a certain symmetric probability
distribution with respect the initial position, has been also worked out
[133]. On a different note, in the context of stochastic resetting, optimal
first passage time is a hot topic [134, 135, 136, 61, 137]. Usually, in those
studies, the external potential is fixed and the resetting rate is the object
that plays the role of the external control.

Another relevant time optimization problem is the minimization of the
connection time between the initial and target states. For the SST between
equilibrium states, the minimum connection time is zero for the unconstrained
problem—similarly to the situation described above for the escape time of the
Brownian particle. A different situation arises when the connection problem
has additional restrictions, for example when the protocol has to be adiabatic
in the thermodynamic sense of zero average heat. Therein, the adiabatic
constraint together with the second principle gives rise to the emergence of a
speed limit, i.e. the emergence of a minimum time for the adiabatic connection
[83]. In general, the instantaneous adiabat does not exist, and moreover there
appear forbidden regions beyond the quasi-static curve
$T^{2}/\kappa=\textnormal{const}$ at which the minimum connection time
diverges—as discussed in Sec. 3.4. Another problem in which a non-vanishing
minimum connection time emerges is the thermal bath engineering of
harmonically confined Brownian particles [138]. Therein, the restriction comes
from the bath temperature being bounded—in a certain interval, from a
practical point of view, and to non-negative values, from a fundamental point
of view. The optimal protocol, i.e. the brachistochrone, is of bang-bang type
and comprises as many bangs as the dimension of the system, with the
temperature alternating between its maximum and minimum value. As a
consequence, the minimum connection time increases with the dimension of the
system $d$, even for spherically symmetric confinement.

###  4.3 Information thermodynamics and speed limits

In the previous sections, we have discussed the optimization of the work and
other figures of merit over SST protocols—including the minimization of the
connection time for some specific situations. Here, we focus on the emergence
of the so-called classical speed limits, which are closely related to
information thermodynamics concepts, and their relevance in the context of
SST.

The acceleration of the connection entailed by SST protocols comes at a price:
for example, we have already discussed that there appears a non-vanishing
irreversible contribution to the average work, which only vanishes for
infinite connection time, for the SST connecting equilibrium states at
constant temperature. In fact, the minimum irreversible work, as given by Eq.
(11lbtbzcdcj), is proportional to $t_{\scriptsize\textnormal{f}}^{-1}$ and
then blows up for very short connecting times. Therefore, a natural question
arises, whether or not there exists a speed limit for SST, i.e. a minimum
value for the connection time $t_{\scriptsize\textnormal{f}}$.

The existence of a speed limit in the context of quantum mechanics dates back
to the 1940s [139]. It is related to the time-energy uncertainty relation; for
review of the subject, see Ref. [140]. The simplest situation is that of the
time evolution of a pure state in a conservative system. In this case, there
appear two inequalities for the time $t_{\scriptsize\textnormal{f}}$ necessary
to evolve from an initial state $\ket{\psi(0)}$ to an orthogonal target state
$\ket{\psi(t_{\scriptsize\textnormal{f}})}$,
$\braket{\psi(0)}{\psi(t_{\scriptsize\textnormal{f}})}=0$. First, the
Mandelstamm-Tamm bound [139], $t_{\scriptsize\textnormal{f}}\Delta
H\geq\pi\hbar/2$, and, second, the Margolus-Levitin bound [141],
$t_{\scriptsize\textnormal{f}}\braket{H}\geq\pi\hbar/2$. In these expressions,
$\braket{H}$ and $\Delta H$ are the (time-independent) expectation value and
standard deviation of the energy, respectively. In fact, it has been proven
that the combination of the Mandelstamm-Tamm and Margolus-Levitin inequalities
gives the tightest bound for $t_{\scriptsize\textnormal{f}}$ [142], i.e.

$t_{\scriptsize\textnormal{f}}\geq\max\left(\frac{\pi\hbar}{2\Delta
H},\frac{\pi\hbar}{2\braket{H}}\right).$ (11lbtbzcdck)

This inequality expresses the existence of a natural time scale, which cannot
be beaten for the evolution of a quantum system with time-independent
Hamiltonian $H$.

In statistical mechanics, only more recently the possible existence of
inequalities resembling the quantum time-energy uncertainty relation—and the
possible associated emergence of a speed limit—has been investigated [104,
143, 144, 145, 146, 147, 148, 149, 150, 83, 151, 62, 107, 152]. One of the
first derivations of an inequality for the connection time is probably that of
Ref. [104], within the linear response regime. For a review of these ideas,
see Ref. [153]. Let us consider a general system with control functions
$\boldsymbol{\lambda(t)}$, and a protocol that drives the system from an
initial state with
$\boldsymbol{\lambda}(t=0)=\boldsymbol{\lambda}_{\scriptsize\textnormal{i}}$
to the target state with
$\boldsymbol{\lambda}(t_{\scriptsize\textnormal{f}})=\boldsymbol{\lambda}_{\scriptsize\textnormal{f}}$;
the initial and final states must be close enough in order to use linear
response theory. The irreversible power can be written as a bilinear function
of the time derivative of the control functions $\boldsymbol{\lambda}(t)$,
with coefficients provided by the—positive definite—friction tensor
$\boldsymbol{\zeta}(\boldsymbol{\lambda})$,

$P_{\scriptsize\textnormal{irr}}(t)=\frac{d\boldsymbol{\lambda^{T}(t)}}{dt}\cdot\boldsymbol{\zeta}(\boldsymbol{\lambda(t)})\cdot\frac{d\boldsymbol{\lambda(t)}}{dt},$
(11lbtbzcdcl)

where the superindex $T$ stands for transpose. This expression naturally
introduces a metric in the problem and, in fact, a statistical length
$\mathcal{L}_{\scriptsize\textnormal{lin}}$ 333The subindex “lin” stresses
that this definition of statistical length is valid in linear response, see
below for a more general definition based on the Fisher information. is
defined as

$\mathcal{L}_{\scriptsize\textnormal{lin}}=\int_{0}^{t_{\scriptsize\textnormal{f}}}dt\sqrt{P_{\scriptsize\textnormal{irr}}(t)}.$
(11lbtbzcdcm)

This length somehow measures the distance swept by the system in parameter
space along the SST path from the initial to the target state. The connection
time $t_{\scriptsize\textnormal{f}}$ and the irreversible work
$W_{\scriptsize\textnormal{irr}}$ are shown to verify the inequality

$t_{\scriptsize\textnormal{f}}\,W_{\scriptsize\textnormal{irr}}\geq\mathcal{L}_{\scriptsize\textnormal{lin}}^{2}.$
(11lbtbzcdcn)

The equality only holds when $P_{\scriptsize\textnormal{irr}}$ is constant, as
it was for the protocol that minimizes the irreversible work—but note that Eq.
(11lbtbzcdci) holds for arbitrary initial and final equilibrium states, not
necessarily close, and thus it is not restricted to linear response.444A
perturbative approach to solve the FP equation has been introduced with the
aim of extending the range of validity of Sivak and Crook’s results [154].

Speed limits have been discussed in a broader context, beyond the linear
response regime [145, 147, 146, 155, 156, 148, 150, 83, 151, 149, 148, 62,
107, 152]. Here, we focus on processes that involve a net transformation of
states; for time-periodic or stationary processes, there are specific
inequalities showing that the entropy production rate bounds the rate at which
physical processes can be carried out [150]—which has been experimentally
checked [152]. On the one hand, the first results were derived under the
assumption of Markovian dynamics [145, 147] and are thus restricted to system
with dynamics described by master, for discrete variables, and FP (or
Langevin) equations, for continuous variables. On the other hand, approaches
based on information geometry concepts hold for general dynamics [149, 148],
not necessarily Markovian. Key to the latter results is the Fisher information
$I(t)$,

$I(t)\equiv\left\langle\left(\partial_{t}\ln P\right)^{2}\right\rangle=\int
d\mathbf{x}\frac{\left(\partial_{t}\rho(\mathbf{x},t)\right)^{2}}{\rho(\mathbf{x},t)},$
(11lbtbzcdco)

which is the curvature of the Kullback-Leibler divergence,

$\int\\!d\mathbf{x}\rho(\mathbf{x},t+dt)\ln\\!\left[\\!\frac{\rho(\mathbf{x},t+dt)}{\rho(\mathbf{x},t)}\right]\\!=\frac{1}{2}(dt)^{2}I(t)+O(dt)^{3}.$
(11lbtbzcdcp)

The Fisher information can be connected with entropy production and the so-
called thermodynamic uncertainty relations [146, 145, 155, 148], and it is
directly linked with the thermodynamic length $\mathcal{L}$, first introduced
in the 1980s [157, 158, 159]. For equilibrium systems, this relation was
addressed in a pioneering work by Crooks [114], which showed that

$\mathcal{L}=\int_{0}^{t_{\scriptsize\textnormal{f}}}dt\,\sqrt{I(t)}.$
(11lbtbzcdcq)

Note that, in general,
$\mathcal{L}\neq\mathcal{L}_{\scriptsize\textnormal{lin}}$.

Also relevant to our discussion is the divergence—also called the
thermodynamic cost [148]—of the path,

$\mathcal{C}=\frac{1}{2}\int_{0}^{t_{\scriptsize\textnormal{f}}}dt\,I(t).$
(11lbtbzcdcr)

The Cauchy-Schwarz inequality makes it possible to establish the following
lower bound for the connection time,

$2t_{\scriptsize\textnormal{f}}\,\mathcal{C}\geq\mathcal{L}^{2}.$
(11lbtbzcdcs)

Interestingly, Eq. (11lbtbzcdcs) is implicitly written in Ref. [114] [see Eq.
(9) therein], although it has not been explicitly stated as establishing a
speed limit for finite-time SST until recently—see Ref. [145] for the case of
Markovian dynamics and Ref. [148] for arbitrary dynamics. At variance with the
quantum case, it must be remarked that the tightness of the bound in Eq.
(11lbtbzcdcs) has not been proven, to the best of our knowledge. Indeed, Eq.
(11lbtbzcdcs) is valid for an arbitrary dynamics—see Ref. [62] for an analysis
thereof for granular fluids described at the kinetic level, where the PDF
obeys the non-linear (inelastic) Boltzmann equation—but there may exist a
larger bound for the connection time.

Equations (11lbtbzcdcn) and (11lbtbzcdcs) provide us with two inequalities
that must hold in SST. In principle, one could argue that Eq. (11lbtbzcdcn)
was derived under the framework of linear response but it is a direct
consequence of the definition of irreversible work and the Cauchy-Schwarz
inequality—it is linking $P_{\scriptsize\textnormal{irr}}(t)$ with a
Riemannian metric that linear response ensures 555Very recently, the geometric
ideas introduced in Ref. [104] have been extended beyond linear response
[109]. Therein, it has been shown that, for connecting equilibrium states, the
protocol that minimizes the energetic cost corresponds to the geodesic path..
A comparison between the predictions of both inequalities for the harmonic
case is given in E.

##  5 Beyond equilibrium: connecting NESS or arbitrary states

In the statistical mechanics context, SST have been mainly employed to connect
equilibrium states [90, 95, 96, 97, 98, 99, 100, 79, 17, 63, 59, 18, 70, 160,
80, 83, 81, 71, 92, 93, 94]. The analysis of accelerated connection between
NESSs has been initiated in [161, 62]. The delay in the engineering of SST
connections between NESSs stems from some difficulties that are inherent to
the initial and final states being non-equilibrium, as explained below.

For the sake of concreteness, we revisit the already addressed paradigmatic