Long-time persistence of hydrodynamic memory boosts microparticle transport

In a viscous fluid, the past motion of an accelerating particle is retained as an imprint on the vorticity field, which decays slowly as $t^{-3/2}$. At low Reynolds number, the Basset-Boussinesq-Oseen (BBO) equation correctly describes nonuniform particle motion, capturing hydrodynamic memory effects associated with this slow algebraic decay. Using the BBO equation, we numerically simulate driven single-particle transport to show that memory effects persist indefinitely under rather general driving conditions. In particular, when driving forces do not vary smoothly, hydrodynamic memory substantially lowers the effective transport friction. Remarkably, this enables coasting over a spatially uneven potential that otherwise traps particles modeled with pure Stokes drag. Our results provide direct physical insight into role of particle-fluid coupling in nonequilibrium microparticle transport.

Introduction.-In the pioneering work of Alder and Wainwright [1], molecular dynamics simulations of hardspheres liquid argon revealed a long-time tail in the velocity autocorrelation function, which decays as t −3/2 . Whereas standard Langevin theory predicts exponential decay over a Brownian relaxation time τ B , this algebraic decay is caused by viscous coupling between nonuniform particle motion and unsteady flow in the ambient fluid [1][2][3]. As a particle accelerates, the fluid it displaces also accelerates, causing a virtual mass force that augments the particle's apparent inertia. The particle also imparts momentum to the fluid as vorticity, which diffuses over a particle radius in a kinematic time τ ν proportional to the fluid density. A particle that translates over its own radius before the ambient fluid relaxes thus feels a force due to vorticity generated by its past motion [4]. This delayed self-interaction force-the Basset history forceunderlies the slow t −3/2 decay associated with hydrodynamic memory.
The Basset-Boussinesq-Oseen (BBO) equation [5][6][7] is a remarkably general description of nonuniform particle motion that captures the virtual mass and history force, both of which become significant when τ ν τ B [8]. Analytical solutions to the BBO equation, however, are generally infeasible; hence the importance of recent progress in addressing the considerable computational challenges of numerically integrating the history force [9,10]. Indeed, at low Reynolds number, it is tempting to neglect fluid inertia a priori, especially since memoryless models are ubiquitous, computationally accessible, and conceptually simple. The memoryless assumption, however, tacitly implies that τ ν τ B ; these timescales are not well separated in dense media (e.g., liquids), even when Re 1 [11][12][13]. This is indeed the case for the myriad mesoscopic phenomena at the micron scale and below, such as conformational changes in biomacromolecules [14,15], unsteady motion of microswimmers [16,17], colloidal flows [4,18], and particle dispersion and sedimentation in turbulence or environmental flows [19]. Yet, without a clear picture of the physical mechanism, a frugal modeler may still reason that memory effects can be safely ignored if relevant observation timescales, such as long-time steady-state behavior, are much longer than both τ ν and τ B .
In this letter, we use numerical simulations to unambiguously show that neglecting the history force can lead to qualitatively incorrect particle transport under general nonequilibrium conditions. Memory effects do not merely protract the transient relaxation to steady-state but persist indefinitely under a significant parameter regime, implying that the BBO description may be necessary regardless observation timescales. We begin by comparing particle dynamics with and without hydrodynamic memory under time-dependent square-wave driving. The history force, which models the transient retention of latent kinetic energy in a viscous fluid bath, is found to significantly reduce the effective friction for driving frequencies below τ −1 ν by smoothing the energy transfer from an abruptly varying force to a particle. We then show that this leads to enhanced transport through a tilted potential with gaps, whereas pure Stokes dissipation tends to suppress itinerant behavior. Additionally, our results demonstrate that BBO dynamics can be efficiently modeled using an extended phase space approach wherein the memory kernel is approximated by an exponential sum [20][21][22].
Problem overview.-Consider a microsphere [23] of radius R, mass m, and density ρ s in an incompressible fluid of density ρ and viscosity η, with no-slip boundary conditions. In the limit of low Reynolds number and Mach number, nonuniform motion is described by the BBO equation [5][6][7], The second term-the quasisteady component-is the familiar Stokes drag with friction coefficient ζ S = 6πηR, which defines the Brownian relaxation time, τ B = m e /ζ S . Unsteady flow introduces two additional forces: (1) the virtual mass force, where m a = 2 3 πR 3 ρ is the added mass and m e = m + m a is the effective mass; (2) the Basset history force, a convolution over past acceleration. Here, the kinematic time, τ ν = R 2 /ν, which depends on the kinematic viscosity of the fluid, ν = η/ρ, quantifies the time for vorticity to diffuse over a distance R.
It is instructive to consider Eq. (1) in dimensionless form. First, we introduce the following dimensional scales, which we use in the ensuing analyses (also cf. Arminski and Weinbaum [24]): where the dimensionless time ist = t/t c , dimensionless velocity isṽ = v/v c , and so on. Accordingly, length is scaled by L c = v th τ B , mass by M c = m e , force by F c = k B T /L c , and power by P c = E c /τ B . By introducing β = 9ρ 2ρs+ρ , one may show that the particle and fluid timescales are related by β = τ ν /τ B [24,25]. Redefining t as t,ṽ as v, and so on, for simplicity, one obtains the following dimensionless form: In the present analyses, we assume neutral buoyancypertinent to many processes in liquid water-so ρ s = ρ, β = 3, and the history term is of order unity by Eq. (3). In this case, the history force can be neglected only if F varies in such a way that a particle will always be moving at nearly its terminal speed andv ∼ 0.
There are a number of derivations of Eq. (1) [26-28], the common starting point being the unsteady Stokes equation [29], which must be used when τ ν /τ B is O(1) [13]. Generally, F may consist of deterministic timedependent and conservative components, as well as thermal fluctuations [30]; respectively, F(x, t) = f (t) − ∇U (x) + ξ(t). Due to the Basset memory kernel, the second fluctuation-dissipation theorem implies that ξ(t) is a colored noise process [31, 32]. In the limit of vanishing fluid inertia, ρ s ρ and τ ν τ B , the history force vanishes, ξ becomes a white noise process, and Eq. (1) reduces to a memoryless Markovian description.
However, the full non-Markovian dynamics of Eq. (1) obfuscates the physical picture considerably. To  bolster physical intuition, we focus exclusively on one-dimensional zero-temperature dynamics; finitetemperature dynamics will be treated in a follow-up study, where colored noise is naturally handled by the extended phase space approach used here [33], which is often called Markovian embedding [21,22,34]. For brevity, we refer to spherical particles obeying Eq. (1) as BBO beads; likewise, particles resisted only by Stokes drag (i.e., zero virtual mass and history forces, with m e → m) are called Stokes beads. Energy flow for a constant-force pulse.-We begin by studying BBO and Stokes beads acted upon by a finiteduration constant force. This simplest example builds intuition for the Basset history force, which can make the bead response quite complicated. Importantly, we plainly illustrate the non-dissipative nature of the history term, which can redistribute energy back to a BBO bead at later times, thereby smoothing its response to a changing external force.
Previously, Arminski and Weinbaum [24] studied the effect of pulse duration τ 0 and waveform (square and triangle pulses), finding that varying τ 0 altered the overall behavior more dramatically than varying the waveform did. They also showed that the maximum bead displacement ∆x max depends only on the total impulse J max = F dt and not the dynamics, i.e., ∆x max = J max ; we will use this important result later. Assuming beads start and end at rest, the history term only acts to alter the details of bead motion on intermediate timescales.
Here, we consider a constant-force square pulse, f (t) = f 0 for 0 ≤ t ≤ τ 0 and 0 otherwise, acting on a bead starting from rest, v(0) = 0. We expect the long-time displacement to be proportional to f 0 τ 0 = √ m e k B T . Figure 1 shows results for a "long" square pulse case of duration τ 0 = 10τ B , where f 0 = 1 10 F c and ∆x max = L c .
Due to the history force and, to a lesser extent, the virtual mass force, a BBO bead responds more gradually than a Stokes bead [ Fig. 1(b)], taking tens of τ B to reach 80% of its long-time displacement [ Fig. 1(a)]. On the other hand, since the force is constant during the driving phase (t < τ B ), the total input work and peak power is evidently lower for a BBO bead [ Fig. 1(c,d)]. Figure 1(c) details how input work is apportioned between the instantaneous kinetic energy, total heat loss (through Stokes dissipation), and fluid stress-energy lost through vorticity diffusion (through the history term). In particular, Stokes dissipation represents energy irreversibly lost to the bath as heat and is always positive. By contrast, the instantaneous power dissipation through the history force abruptly switches sign shortly into the relaxation phase [t > 10τ B , Fig. 1(d)], representing the transfer of stored latent energy in the fluid back to the BBO bead. This particular effect is less pronounced for short-and medium-duration pulses, though the delayed displacement at intermediate times is quite similar [35]. In any case, a BBO bead takes longer to displace than a Stokes bead but, since hydrodynamic memory improves momentum retention, less energy is required to achieve the same displacement.
Temporal square-wave forcing.-The single pulse results suggest that the Basset term improves the retention of momentum and energy as compared to Stokes drag alone, motivating an examination of beads driven by time-periodic forcing. However, to meaningfully compare differences between BBO and Stokes beads in transport, we should quantify the tradeoff between input energy, total displacement, and overall transport speed. To this end, we may define an effective friction coefficient ζ t by analogy with Stokes' law, defined using the following time-averaged quantities (indicated by overbar and subscript t): the effective velocity, and effective driving force, where W t is net input work, ∆x t = x(t) − x(0) is net displacement. Note that force balance in steady-state implicitly defines an effective drag via F t ∼ −F drag . Using Eqs. (4) to (6), it is easy to compute running estimates of the effective friction. We now consider beads driven by square waves with period τ = 10τ B , pulse duration τ 0 , and amplitude f 0 , holding fixed the impulse per cycle, an always-on force (τ 0 = τ ). Figure 2(a) shows that BBO and Stokes beads achieve around 90% and 100% of the maximum possible displacement by t = 250τ B , consistent with the expectation that long-time displacement depends only on the total impulse, ∆x max = 25L c . Moreover, since fixing J cyc determines the cycle-averaged momentum input rate, we expect BBO beads to attain comparable net velocities to Stokes beads after many forcing cycles; Fig. 2(b) confirms this intuition. For always-on forcing when the pulse duration is equal to the period, the effective driving force is identical for BBO and Stokes beads [black, Fig. 2(c)]. By contrast, increasing the abruptness of the pulses (shorter duration, larger magnitude) also increases the apparent effective driving force (i.e., increases the effective drag), but with a smaller effect on the BBO beads; however, the effective drag increases much more for a Stokes bead than a BBO bead as abruptness is increased [red/blue, Fig. 2(c)]. This behavior is expected because, as compared to a Stokes bead, a BBO bead has a lower instantaneous velocity during the forcing phase [cf. Fig. 1(c)], which reduces the total input work [cf. Fig. 1(b)]. It follows that a BBO bead's effective friction is comparable to a Stokes bead's when the driving force is always on but is lower when intermittently driven.
Tilted spatial potential with gaps.-We now consider beads in a tilted potential with spatial period λ and an overall tilt F 0 that acts as an average driving force. We use a simple, piecewise-linear construction that alternates between downward-sloped segments (i.e., constant, positive force) and gaps, flat segments of length d with zero force [Fig. 3, black/gray line]. This simple construction dovetails with the preceding square-wave analysis.
All beads are initialized with zero velocity at the top of a downhill segment [36]. A bead is therefore initially driven by a constant force of magnitude FIG. 3. Tilted periodic potential (black line) with wavelength λ, gaps of length d, and average slope F0 (red dashes). Here, an itinerant BBO bead progresses beyond the first gap, but a Stokes bead is trapped.
until reaching the first gap, by which time it has received a fixed input of energy, W = −∆U 0 , where ∆U 0 = −F 0 λ, and the bead begins coasting. It is clear a bead must coast past the first gap so as to achieve sustained transport; otherwise, it becomes localized and transport is halted. Whether a BBO or Stokes bead becomes itinerant depends on two factors dictated by its dynamics: the momentum carried by a bead upon entering the gap and the rate at which momentum is dissipated while coasting.
Once in sustained transport, an itinerant bead receives an average energy per unit displacement (over a wavelength λ) equal to the average potential slope F 0 , which in turn fixes the effective driving force. After an itinerant bead traverses many wavelengths Recall that time-periodic driving, by contrast, essentially fixes the effective velocity but not the effective driving force [cf. Fig. 2(b,c)]. Accordingly, we anticipate that the effective friction in a tilted potential is predominantly determined by the effective velocity, which in turn is dictated by the dynamics. . Varying the wavelength λ, slope F 0 , and gap size d over a wider parameter space demonstrates that this effect is general [37]. The point we wish to emphasize, however, is that it is possible construct a potential in which ordinary initial conditions lead to the localization of Stokes beads while BBO beads become itinerant.
Discussion.-In the tilted potential, a bead becomes itinerant when the momentum initially gained downhill is sufficient to coast past the first gap. The initial momentum input is given by the initial impulse, where the time t * to arrive at the first gap depends on the dynamics. In particular, t * BBO ≥ t * Stokes , because a Stokes bead travels at least as fast as its BBO counterpart [ Fig. 2 Since long-time displacement only depends on J * , BBO beads coast farther than Stokes beads; if BBO beads get trapped, so do Stokes beads. Hydrodynamic memory therefore always aids transport here by allowing larger gaps to be traversed.
The history term retards momentum decay, smoothing out the response to fluctuations in momentum/energy input [ Fig. 2(b)]. From a physical standpoint, the transient retention of momentum/energy in the fluid bath endows BBO beads with latent kinetic energy, manifesting upon rapid bead deceleration when the external force drops to zero [ Fig. 1(d)]; the remaining momentum imparted to the fluid is lost at domain boundaries through vorticity diffusion. Broadly speaking, the history force-by tempering sharp instantaneous velocity gains for abruptly changing driving forces-mitigates Stokes dissipation, which grows quadratically with velocity. Moreover, for driving periods larger than roughly τ ν [38], the history term accumulates latent bead energy over many cycles, causing the cycle-averaged work Our findings provide a foothold for understanding the compressible flow case [18,39,40], as well as full non-Markovian dynamics at finite temperature where hydro-dynamic memory drives long-time persistence of thermal fluctuations [41]. In view of Ref. [42], where it was shown that even tagged fluid particles closely follow BBO dynamics, the present results are germane to temperature and pressure control algorithms in molecular dynamics simulations, especially as extended phase space approaches show promise in improving the kinetic consistency of thermostats [43,44]. Indeed, such approaches can be brought to bear on the question of stochastic transport efficiency of linear molecular motors like kinesin [45][46][47] The total external force may be decomposed into conservative, time-dependent, and fluctuating thermal force components: In the present analyses, the thermal noise ξ is neglected in the zero-temperature limit. The physical parameters in Eq. (1) are summarized in Table I. There are two important timescales implied in Eq. (1), one corresponding to the bead and the other to the fluid (see Table I). The particle relaxation time of interest is τ B , the Brownian relaxation time. τ B characterizes the rate of momentum damping due to Stokes drag, which gives rise to exponential decay of velocity autocorrelation function. The competing fluid timescale is the kinematic time τ ν , which corresponds to the timescale of vorticity diffusion τ ν and characterizes the time to transport (shear) momentum in the fluid over a particle radius.

NUMERICAL METHOD
The BBO equation [Eq. (1)] is challenging to solve numerically due to complications arising from the integral and, in the finite-temperature case, time-correlations in thermal fluctuations that give rise to colored noise [14]: 1. The memory kernel t −1/2 is singular at t = t , which requires care when performing direct numerical integration of the convolution term.
2. Computing the convolution integral necessitates storage of a relatively large portion of the past velocity trajectory.
3. The colored noise necessitates the generation of correlated pseudo-random number sequences with the correct statistics.
Be that as it may, there is a powerful approach known as a Markovian embedding, which can be used to transform a history-dependent integrodifferential equation into a numerically tractable system of equations that are local in time. To wit, it has been known for some time that any non-Markovian process can be represented as a projection of a higher-dimensional Markov process [7,15,21]. Although we are concerned with the zero-temperature case in the present study, the Markovian embedding technique naturally address the above three challenges, including the generating of the correct colored noise. We leverage Markovian embedding to obtain an approximation to the non-Markovian GLE represented by Eq. (1) in terms of a finite set of auxiliary variables. There are a number of advantages to such an approach, as it: (1) avoids the computational (i.e. memory intensive) expense and complication of direct computation of the Basset history force (2) generates, on-the-fly, the statistically correct colored noise satisfying the fluctuation-dissipation theorem [12,16]; (3) allows the resulting (embedded) dynamical system to be propagated forward in time using standard, robust Langevin integrators, controlled to an arbitrary degree of precision; and (4) can achieve acceptable levels of statistical precision small using relatively few auxiliary variables (i.e. small embedding dimension).

Markovian embedding of BBO in the zero-temperature limit
We follow the approach taken by Goychuk [9], Siegle et al. [19] and, more generally, Baczewski and Bond [2], wherein the algebraic form of the memory kernel is approximated as a weighted sum of exponential terms. This is known as a Prony series approximation in constitutive modeling contexts [8,13]). This embedding leads to an extended phase space representation in terms of the usual position x and velocity v, as well as a finite set of N auxiliary variables s i , where i = 1, ..., N , each obeying a stochastic differential equation (SDE) with white Gaussian noise. Note that the numerical system studied here will be somewhat simplified, as we can neglect the noise terms in the zero-temperature limit.
Using the fractal scaling employed by Siegle et al. [19], we obtain an approximation to the full memory kernel (which includes both the Stokes drag-which were dropped, in an ad hoc fashion, in the equation of motion studied by Siegle et al. [19]-as well as the history term), where ν k = ν 0 /b k , b is a scaling parameter, ν 0 is a high-frequency cutoff, and C(b) is a fitting coefficient. The delta function in the first term comes directly from Stokes' law; the other delta function (in the sum) comes from the positive singularity in the memory kernel (when t 0 → 0 − ) and is necessary to ensure the integral of the kernel in the long-time limit vanishes [18,19]. 1 The slowly decaying power law, t −3/2 , can be approximated over N log 10 b − 2 decades through a suitable choice of auxiliary dimensions, N , and temporal spacing, b, of the exponential terms in Eq. (2) [10,19]. The primary restriction on the high-frequency cutoff is that the timestep be small enough to achieve numerical stability when an explicit integration method is used, i.e., ∆t ν −1 0 . In principle, decreasing b and increasing N allows the true memory kernel to be approximated to an arbitrary degree of precision. In most cases, however, one is likely to be interested in efficiency, whereby the memory kernel may be approximated by a small number of auxiliary variables that satisfactorily capture the relevant physical range of timescales.
The embedding specified via Eq. (2) transforms the linear, integro-differential BBO equation, Eq. (1), into a coupled system of D = N + 2 linear equations with N auxiliary variables of the exponential relaxation type. The auxiliary variables comprise an N -dimensional extended phase space that approximates both the convolution integral representing the history force:ẋ where N is the number of auxiliary variables used to approximate the power law in Eq. (1), k = 1, . . . , N , and D = N + 2 is the dimensionality of the embedded dynamics. For notational brevity, we have used γ B = ζ/m e , which is the Brownian collision frequency, and γ 0 = k γ k .

Euler-Maruyama integration
The equations of motion were solved using forward-Euler integration, which was sufficient for the purposes of the present study. The BBO code is implemented in the Cython language [5], which is available open source (https: \bitbucket.org/sseyler/glsimulator/src/bbo/) and will be refined in the near future to enable more general use.
The forward-Euler discretization of Eq. (4) is where ∆t is the timestep, and a superscript n denotes evaluation at the n th timestep, and so on. We have implemented more sophisticated integration schemes, including a stochastic symplectic scheme that might be suitable for more sophisticated numerical experiments [2,20]. A dimensionless timestep of ∆t = 10 −5 was used and found to be more than adequate in terms of accuracy. We used N = 12 auxiliary variables for the single-pulse simulations [ Fig. 1 in the main text], and up to N = 21 for longer simulations with time-periodic forcing or simulations in the tilted potential [ Fig. 2-4 in the main text]. A highfrequency cutoff ν 0 = 10 4 provided sufficient accuracy in all simulations. Following Siegle et al. [19], we set the scaling parameter b = 5 and fitting coefficient C(b) = 1.78167, which were sufficient to generate essentially indistinguishable trajectories from the analytical solutions (in the single-pulse case).
Finally, consider the following substitutions in the equation for velocity, and, for the auxiliary equations, These substitutions, which are helpful for implementing the algorithm, enable the conversion of Eqs. (5) to (7) to the final simple form:

DECREASING THE DURATION OF A SINGLE SQUARE PULSE
In the main text, Fig. 1 shows the results for a single pulse where the pulse duration was an order of magnitude longer than the Brownian relaxation time, i.e., 10τ B . This result corresponds to the long-pulse results presented in Fig. 3 of Arminski and Weinbaum [1]. When the pulse is long, the sign-change of the Basset history force is rather abrupt once the external force is removed and the bead enters the relaxation phase. We also show results for an abrupt pulse with a short duration of 1 10 τ B [ Fig. 1], which also corresponds to the short-pulse result presented in Fig.  2 of Arminski and Weinbaum [1]. For a very abrupt pulse whose duration is an order of magnitude shorter than τ B , the history force clearly dominates the Stokes drag. Thus, the energy imparted by short-timescale forces on a BBO bead will be dissipated primarily through the action of vorticity diffusion, which carries the energy away from the bead and toward the fluid boundaries [ Fig. 1(b,d)]. That the Basset history force dominates Stokes drag on timescales shorter than the Brownian relaxation time agrees with the conclusions reached by Kheifets et al. [11], where the history force is seen to be dominant on the timescale of thermal velocity fluctuations. Conversely, when a force is applied to a BBO bead gradually, dissipation of the input energy is shared more equally between the Stokes drag and history terms.
For completeness, we also present results for the case of a pulse of "moderate" abruptness, having a duration exactly equal to τ B [Fig. 2]. In this case, the Basset history force does in fact appear to change sign between 1 and 1.5 τ B , albeit only slightly [ Fig. 2(b,d)]. The shift is much less dramatic than for when the pulse is longer than a Brownian relaxation time (i.e., 10τ B ). We conclude that applying a longer duration pulse therefore has the effect of storing energy in the vorticity field of the imaginary fluid bath surrouding the particle. This energy, which is really the energy associated with the stress field of the fluid, can be thought of as latent bead kinetic energy, since some of it will be returned to the bead in the future (in the form of momentum). Our simulations suggest that a square pulse of duration ∼ 3τ B = τ ν causes the Basset term to maximize the return of latent energy to the bead during the relaxation phase. Our simulations suggest that the value of the pulse duration where latent energy return becomes important is around τ ν = 3τ B .

SQUARE-WAVE DRIVING: VARYING DRIVING PERIOD AND PULSE DURATION
Using a longer-period forcing frequency, we see the same overall behavior as in Fig. 2 of the main text, where BBO beads exhibit a lower friction under time-periodic driving conditions [ Fig. 3]. Again, the difference between the Stokes and BBO beads' effective friction is clearly established after only several forcing cycles.
Decreasing the forcing period to exactly one kinematic time (τ = 3τ B = τ ν ), however, qualitatively changes the behavior [ Fig. 4]. When the force is always on (i.e., pulse duration is equal to the period, τ 0 = τ = τ ν ), the effective friction for Stokes beads rapidly approaches a minimum value of 1 [dotted black line, Fig. 4(d)]. BBO beads, on the other hand, slowly approach this value, which is indicative of the algebraic decay of the memory kernel [solid black line, Fig. 4(d)]; the effective friction is about 1.5 times the value for a Stokes bead after several tens of Brownian relxation times, around 1.25 times the value after around 100τ B , and within 10% at 500τ B . However, more interesting behavior results when pulses are applied for a brief portion of a driving period. In particular, while BBO beads initially have a larger effective friction than Stokes beads, there are crossover points around roughly 50τ B and 150τ B ; after 500τ B , the effective friction is about 15% lower for short pulses [τ 0 = 1 10 τ , blue, Fig. 4(d)] and 10% lower with very short pulses [τ 0 = 1 100 τ , red, Fig. 4(d)]. For shorter-duration pulses, the situation is more interesting. When the pulse duration is τ 0 = 1 10 τ = 1 10 τ , the effective friction for BBO beads starts off above that of Stokes beads [blue lines, Fig. 4(d)]. However, after around 60τ B , the BBO beads actually have a lower effective friction than Stokes beads. This value continues decreasing up to and beyond the duration of the simulation (i.e., 500τ B ), at which point the effective friction is around 30% lower for BBO beads than Stokes beads; Stokes beads attain their minimum effective friction after several tens of τ B . The situation is similar in the case of very short pulses, τ 0 = 1 100 τ = 1 100 τ , where Stokes beads have a lower effective friction than BBO beads up to around 100τ B 30τ ν , after which the effective friction of BBO beads drops below that of Stokes beads [red lines, Fig. 4(d)]. After 500τ B , BBO beads have about a 15% lower effective friction than Stokes beads. When the period is equal to the kinematic time, hydrodynamic memory still appears to produce a measurable advantage that manifest in a lower effective friction for BBO beads versus Stokes beads.
Once the forcing period is decreased well below a kinematic time, τ = 1 10 τ ν = 3τ B [Fig. 5], the differences between BBO and Stokes beads are much more muted. In particular, the effective friction for all beads is seen to settle to ζ ≈ 1 [ Fig. 5(d)] rather quickly. The friction of BBO beads versus Stokes beads is only substantially larger for the first 100-200τ B , after which values are with 30% of the minimum friction value achieved by Stokes beads in the always-on case [i.e., when τ 0 = τ , dotted black lines, Fig. 5(d)]. The reason for the differences in the high-frequency case can be rationalized by directly comparing the effective driving forces and velocities. Whereas the effective driving forces are comparable between BBO and Stokes beads subjected to identical forcing protocols in time [ Fig. 5(c); values are within roughly 20% of one another after only ∼ 10τ B ], BBO beads have a measurably lower effective velocity than Stokes beads [ Fig. 5(b); values differ by over 40% around 10τ B and are almost within 20% after 50τ B ]. However, if the BBO beads driven under this forcing protocol for a longer time (i.e., more cycles), it is expected that the effective friction achieved in steady state will eventually be the same as that achieved by Stokes beads. Overall, the differences between BBO and Stokes beads are nevertheless small as compared to the differences that arise for driving forces with a much longer period [as in Fig. 2 in the main text and Fig. 3 above]. Therefore, if Stokes beads have any transport advantage at all (by not having the Basset term present), it is for the scenario in which a high frequency force with period τ τ ν is applied for a total length of time somewhat shorter than 100τ B -before the effective friction of a BBO bead has had time to relax. Indeed, in the long-time limit, the effective friction for BBO beads would converge to the value achieved by Stokes beads. However, the overall differences between BBO and Stokes beads for high frequency forcing is relatively small as compared to when driving frequencies are below τ −1 ν ; at high frequencies, the effective friction is less sensitive to the pulse duration (i.e., abruptness) and is relatively close to the theoretical minimum value for always-on forcing.
To recapitulate, when the forcing frequency is below τ −1 ν , then there is a measurable drop in the effective friction when memory is included. When the frequency is right at τ −1 ν , there is a crossover after around 100 Brownian relxation times, where the BBO bead effective friction drops below that of the Stokes bead (for a given forcing function). In the case where the frequency is ten times above τ −1 ν , the values for the friction all appear to approach a long-time value of ζ = ζ S , with the effective friction of the BBO beads taking slightly longer than that for the Stokes beads, though the differences are relatively minor. It therefore seems that, for non-smoothly-varying driving forces, BBO beads are somewhat more difficult to transport for high-frequency driving (i.e., when BBO beads do not have enough time to reach a relatively constant steady-state velocity), but substantially easier to transport for low-frequency driving.

TILTED SPATIAL POTENTIAL WITH GAPS: VARYING WAVELENGTH AND STEEPNESS
In order to test the robustness of our results, we performed long simulations t f = 2500τ B for both BBO and Stokes beads in the tilted spatial potential across five different values of the wavelength λ and five values for the overall tilt (i.e., overall slope) F 0 of the tilted potential, for a total of 25 wavelength-tilt combinations. Then, for each pair of values λ and F 0 , we then tested the sensitivity of bead transport to the gap width by varying d from 1 20 λ to λ in increments of 1 20 λ, simulating one bead for each value of the gap width. In total, we performed 25 × 20 = 500 simulations for each type of bead, with 1000 simulations in total.  Figure 6 and Figure 7 summarize the results for Stokes and BBO beads, respectively. A value of 1 means that, for the given combination of wavelength and tilt, 100% of the beads were able to clear the gaps and become itinerant. Equivalently, this means that the maximum gap width that was cleared by a bead was 95% of the wavelength of the potential. Similarly, a value of 0.1, for instance, means that only 10% of the beads (2 out of 20) acheived sustained transport; one of the beads encountered no gap (constant slope throughout the potential) and one of the beads cleared a gap that was 5% of the wavelength, i.e., a gap width d = 1 20 λ. Since the gap width is set as a fraction of the wavelength, as the wavelength of the tilted potential is increased, the gap size d increases proportionally. As a result, beads must clear a larger gap for a given gap fraction, so a larger amplitude is required to generate itinerant behavior for a given gap size. For the smallest wavelengths and steepest slopes, all beads, whether BBO beads or Stokes beads, were able to escape [top left, Fig. 6 and Fig. 7]. As the wavelength is increased, holding the potential slope constant, Stokes beads become localized relatively easily. For instance, when the potential slope is large 4 5 F 0 and the wavelength is λ, Stokes beads can only cross the gap when the gap width is less than 40% of a full wavelength, whereas BBO beads can cross a gap that is 90% of the wavelength; when the wavelength is 2λ, the maximum crossable gap width for Stokes beads is less than 20% of the wavelength, but 75% for BBO beads [top row, Fig. 6 and Fig. 7].