Minimal work protocols for inertial particles in non-harmonic traps

TL;DR

The paper addresses minimizing work for an underdamped Brownian particle during transitions between non-harmonic end-states. It develops a multiscale perturbation framework and a numerical method based on solving a cell problem for fields rho and sigma, enabling non-Gaussian boundary conditions. Key contributions include explicit near-overdamped forms of the optimal protocol, analysis of momentum and position cumulants showing symmetry-breaking due to inertia, and a demonstration that the overdamped entropy-production bound is tight in the adiabatic limit. This approach advances nanoscopic control by providing practical protocols and insights for inertia-dominated stochastic thermodynamics.

Abstract

Progress in miniaturized technology allows us to control physical systems at nanoscale with remarkable precision. Experimental advancements have sparked interest in control problems in stochastic thermodynamics, typically concerning a time-dependent potential applied to a nanoparticle to reach a target stationary state in a given time with minimal energy cost. We study this problem for a particle subject to thermal fluctuations in a regime that takes into account the effects of inertia, and, building on the results of a previous work, provide a numerical method to find optimal controls even for non-Gaussian initial and final conditions, corresponding to non-harmonic confinements. The control protocol and the time-dependent position distribution are qualitatively different from the corresponding overdamped limit: in particular, a symmetry of the boundary conditions, which is preserved in the absence of inertia, turns out to be broken in the underdamped regime. We also show that the momentum mean tends to a constant value along the trajectory, except close to the boundary, while the evolution of the position mean and of the second moments is highly non-trivial. Our results also support that the lower bound on the optimal entropy production computed from the overdamped case is tight in the adiabatic limit.

Figures (6)

• Figure 1: Optimal nucleation process. We compute the optimal evolution with boundary conditions \ref{['eq:nucleation']}, in the overdamped (orange) and underdamped (blue) case, in finite time. Panels (a)-(h) show the probability density \ref{['eq:fover']} in the overdamped dynamics, computed with the method outlined in Section \ref{['sec:numerical']}, and the position marginal distribution \ref{['eq:f']} in the underdamped dynamics. Panels (i)-(p) show the drift in the two cases, computed by Eqs. \ref{['eq:uover']} and \ref{['eq:control']}. Parameters of the dynamics: $\mathrm{t}_f=2$, $\varepsilon=0.2$, $\tau=\beta=1$, and $g=0.01$. For the numerical method we use $N=2 \times 10^5$ and $h=0.2$. The drifts are shown for area of mass greater than $10^{-4}$.
• Figure 2: Evolution of the cumulants. For the process of Fig. \ref{['fig:potentials']}, we compute the overdamped (orange) and the underdamped moments (dark blue). We also simulate the dynamics \ref{['eq:dyn']} under the predicted underdamped optimal control protocol \ref{['eq:control']} and compute the corresponding empirical cumulants, shown in light blue. Dynamics \ref{['eq:dyn']} is simulated through a Euler-Maruyama scheme with time-step $\Delta t=0.005$, and $5 \times 10^5$ independent trajectories are considered. Boundary conditions are applied on the drift for areas of density smaller than $10^{-5}$: on the left, we use the absolute value of the boundary value and its negative on the right.
• Figure 3: Dependence on $g$. Momentum mean (a) and variance (b) for different values of $g$ in the nucelation process, as functions of time. Assigned boundary conditions and all other numerical parameters and method used are the same as in Fig. \ref{['fig:potentials']}
• Figure 4: Solution in $(p,q)$ space. In each panel, we plot the joint distribution (center), and the marginal distribution of the position (left) $\tilde{f}_t(q)$ and momentum $\tilde{f}_t(p)$ (top). Different panels correspond to different times. The joint density is computed as a numerical average of weighted sample trajectories of \ref{['eq:dyn']}. For details of the method, see Section 3 of sanders2024numericalintegration. Trajectories are discretized by the Euler-Maruyama scheme with time step $\Delta t=0.005$ and $5\times 10^4$ independent realizations from $2.6 \times 10^3$ equally spaced points in $[-3,3]\times[-10,10]$. On the left inset, the perturbative prediction for marginal density of the position in the underdamped dynamics is shown in orange, and on the top inset, the Gaussian prediction for the marginal density of the momentum is shown by the dashed green line. This is estimated by a Gaussian distribution whose moments are the momentum mean and variance as shown in Fig. \ref{['fig:cumulants']}. Other parameters as in Fig. \ref{['fig:potentials']}.
• Figure 5: Lower bound on the cost. The total entropy production \ref{['model:ep']} (blue) of the nucleation process in the underdamped dynamics is shown for different values of $t_f$, as a function of $g$ (decreasing along the horizontal axis), and compared with the lower bound \ref{['eq:overdamped_lower_bound']} (orange). We also report, for comparison, the scaled $\mathcal{W}_2$-distance \ref{['eq:w2_dist']} (green dashed), computed from a numerical estimate of the squared Wasserstein-2 distance between the sampled histograms of initial and final distributions. At $g = 10^{-6}$, the difference between the entropy production cost and the lower bound is $1.4$ for $t_f=2$, $0.14$ for $t_f=5$, and $1.3 \times 10^{-3}$ for $t_f=50$. Assigned boundary conditions and all other numerical parameters are as in Fig. \ref{['fig:potentials']}.