# Self-avoidant memory effects on enhanced diffusion
in a stochastic model of environmentally responsive swimming droplets

Katherine Daftari  Katherine A. Newhall<EMAIL_ADDRESS>University of
North Carolina at Chapel Hill

(January 2021)

###### Abstract

Enhanced diffusion is an emergent property of many experimental microswimmer
systems that usually arises from a combination of ballistic motion with random
reorientations. A subset of these systems, autophoretic droplet swimmers that
move as a result of Marangoni stresses, have additionally been shown to
respond to local, self-produced chemical gradients that can mediate self-
avoidance or self-attraction. Via this mechanism, we present a mathematical
model constructed to encode experimentally observed self-avoidant memory and
numerically study the effect of this particular memory on the enhanced
diffusion of such swimming droplets. To disentangle the enhanced diffusion due
to the random reorientations from the enhanced diffusion due to the self-
avoidant memory, we compare to the widely-used active Brownian model.
Paradoxically, we find that the enhanced diffusion is substantially suppressed
by the self-avoidant memory relative to that predicted by only an equivalent
reorientation persistence timescale in the active Brownian model. We attribute
this to transient self-caging that we propose is novel for self-avoidant
systems. Additionally, we further explore the model parameter space by
computing emergent parameters that capture the velocity and reorientation
persistence, thus finding a finite parameter domain in which enhanced
diffusion is observable.

##  I Introduction

Active particles are a class of nonliving nonequilibrium systems that derive
their motility from environmental energy consumption that is transformed into
self-propulsion. There is great variety among motility-inducing forces,
including chemical forces [1, 2, 3, 4], photoelectric forces [5, 6, 7, 8], and
autophoresis [9, 10, 11, 12, 13]. Such active systems are designed to mimic
self-propulsion seen in microscale living systems, including the run-and-
tumble behavior exhibited by bacteria [14, 15], chemotactic responses to
environmental stimuli [16, 17, 18], gravitactic responses [19, 20], and
photoelectric responses [21, 22, 23]. For a comprehensive review of micro-
scale active systems and current research developments, see Refs. [24, 25,
26].

A hallmark feature of active systems is a ballistic movement, or “swimming”
that when interrupted by random and frequent directional changes gives rise to
enhanced diffusion [6, 8, 15, 27]. The biological advantage of enhanced
diffusion is greater exploration of an area in a shorter period of time when
compared to passive diffusion. Consequently, biological effective diffusion
and other single-particle emergent behaviors such as micro-scale transport,
bacterial motion, and cell migration patterns and their biomimetic
applications are research areas of great interest [24].

In parallel, complete understanding of these phenomena via mathematical
modeling provides design inspiration and permits cost-effective testing of
novel systems; the most common model for active particles is the active
Brownian particle (ABP) model. ABP combines directed motion resulting from a
velocity dependent on the amount of available energy or “fuel” with a
rotational diffusion dependent on a defined persistence timescale, resulting
in enhanced diffusion at time scales longer than the correlation time of the
rotational diffusion [27]. This model of competing ballistic and diffusive
motion accurately predicts the enhanced diffusion of many experimental
systems, such as those found in Refs. [2, 3, 4, 6, 28, 26, 27].

We are interested in the additional effect of spatio-temporal memory observed
in slowly-dissolving autophoretic droplets [10, 11, 12, 13], in addition to
the persistence memory seen in ABP. As these autophoretic droplets interact
with the surfactant suspension, the particular physics induces a self-avoidant
memory response. Above a critical surfactant concentration, the leaking solute
from the droplets creates oily micelles. Local heterogeneities in the
surfactant concentration induce Marangoni stresses that cause the droplets to
spontaneously swim in a random direction. As they move, each particle
continually emits a diffusing wake of surfactant, thus changing the local
concentration profile. Via response to the field, such changes mediate
repulsive inter-particle and intra-particle interactions; the particles veer
to avoid other particles and places where the other particles, or itself, have
recently been. Due to the slow dissipation of the wake, the surfactant
concentration field contains the slowly-diffusing history of the particles’
past locations.

Despite being too large for the effects of thermal noise to be visible, the
ballistic motion of the autophoretic experimental droplets is still punctuated
by randomized directional changes. Such changes reflect a transition between a
dipolar (swimming) and a quadrupolar (stopped) hydrodynamic mode and the
distribution of re-orientation events increases with Péclet number, droplet
size, and the viscosity of the surrounding suspension [13]. While this run-
and-tumble-like behavior produces an enhanced diffusion that is consistent
with the ABP model [12], we seek an understanding of the additional self-
avoidant memory effect at play, particularly on the enhanced diffusion.

We propose a model with a tunable memory response (which we distinguish from
directional persistence) while ignoring hydrodynamic effects for a particle’s
modification of and response to its environment by considering such effects as
built-in, rather than emergent, features. In this model, the particle is a
mobile source of diffusing surfactant that descends its self-produced
concentration gradient, resulting in a sustained “swimming” state. To recover
the stochasticity produced by local fluctuations in the surfactant gradient,
we introduce additive Brownian noise. The interaction of the noise and self-
produced concentration gradient result in directional reorientations mimicking
those produced by the experimental particles. This results in enhanced
diffusion that intuitively one might expect the encoded self-avoidant memory
to amplify as the particle evades its own past locations. However, we find the
opposite: a suppression of enhanced diffusion over that predicted by an ABP
with the same velocity and orientational persistence. We find evidence of
transient self-caging as a possible explanation for this behavior.
Furthermore, we surmise that this contrast from the motivating experimental
systems is likely due to the origin of the noise: in the experiment, local
hydrodynamic effects produce the random reorientations, but in our model, the
additive Brownian noise perturbs the memory-induced swimming direction.

In this paper, we begin in Sec. II by presenting the mathematical details of
our model for self-avoidant swimming droplets. We investigate the memory
effects of these model swimmers at long time scales by comparing the mean
square displacement (MSD) to that of ABPs with the same velocity and
orientational persistence in Sec. III. To make these comparisons, we
analytically derived an expression for the velocity in our model and
numerically compute its orientational persistence timescale. We find that the
equivalent ABP overestimates the enhanced diffusion of the model self-avoidant
droplets, which we attribute to an unexpected side-effect of self-avoidant
memory: transient self-caging. In Sec. IV we further investigate the parameter
space of the model, finding that with fixed noise strength, there is a limited
regime of self-avoidant-memory strength within which enhanced diffusion is
observable; the zero-memory limit of our model is not ABP. We conclude the
paper in Sec. V.

##  II A Model for self-avoidant Memory

Motivated by the experimental system described previously, we propose a
coupled model of a diffusion partial differential equation (PDE) for the
surfactant concentration $c(\mathbf{x},t)$ and a stochastic differential
equation (SDE) for the particle’s location $\mathbf{X}(t)$. These equations
are

$\begin{split}\partial_{t}c(\mathbf{x},t)&=D\Delta c(\mathbf{x},t)+\alpha
DR^{2}\delta_{R}(\mathbf{x}-\mathbf{X}(t)),\end{split}$ (1a) for
$\mathbf{x}\in\Omega\subseteq\mathbb{R}^{2},t\geq 0$, and
$\begin{split}d\mathbf{X}(t)&=-\beta
R\left(\int_{\Omega}\delta_{R}(\mathbf{x}-\mathbf{X}(t))\nabla_{\mathbf{x}}c(\mathbf{x},t)d\mathbf{x}\right)dt\\\
&\quad\quad\quad\quad\quad\quad+\sqrt{\sigma}d\mathbf{W}(t),\end{split}$ (1b)

with prescribed initial conditions ($c(\mathbf{x},0)=0$ and
$\mathbf{X}(t)=\mathbf{0}$ unless otherwise noted) and boundary conditions
(reflecting boundary conditions on $\partial\Omega$ unless otherwise noted).

Equation (1a) is a diffusion equation with diffusivity $D$ and a source term
at the particle’s current location. We note that the inclusion of diffusion on
the source term differs from similar models [33, 34, 35]. Motivation for this
decision and the resulting effects are discussed later in this section. We
introduce a “size” $R$ to the particle using the radially-symmetric mollified
delta function, $\delta_{R}(\mathbf{x}-\mathbf{X}(t))=\frac{1}{2\pi
R^{2}}e^{-\frac{|\mathbf{x}-\mathbf{X}(t)|^{2}}{2R^{2}}}$. Our particle
produces a Gaussian emission pattern with rate $\alpha D$ and remains fixed in
size, thus we ignore depletion effects. (For the treatment of the particle as
a point source with the Dirac delta function, see App. A.2; interestingly the
particle does not swim.) The time evolution of the concentration field holds
the temporally-decaying memory of the particle’s spatial history. (For an
elegant characterization of the effects of the evolving concentration field as
experienced by the particle, refer to [29].)

Equation (1b) is a modified Langevin equation for a Brownian particle in a
force field in the strong friction limit. We again mimic the size of the
particle by convolving the gradient of the concentration field with the
mollified delta function. We define the particle’s response strength to the
concentration gradient to be $\beta R$. As is conventional for an overdamped
Langevin equation, $\mathbf{W}$ is a two-dimensional Weiner process scaled by
the noise strength $\sqrt{\sigma}$. Recall, this noise is to reproduce the
stochasticity introduced by local fluctuations in the surfactant gradient. The
experimental particles are athermal.

The stated model in Eq. (1) articulates the explicit relationship between the
evolving concentration field and the particle trajectories. As the particle
moves, its emissions induce changes in the local concentration field and it
leaves behind diffusing physical evidence of its trajectory. Thus, the
historical information or memory of the particle’s past locations is contained
within the current state of the evolving concentration field. The memory
encoded in the concentration field allows each particle to “remember” where it
has been (hotspots in the concentration field) and avoid its past trajectory
with response strength decreasing as the time lag increases. Via integration
of the whole gradient field at each time point, the particle becomes spatially
omniscient as it moves with dependence on the affecting forces from every
spatial point on the domain. In time, the particles are pseudo-omniscient as
their ability to “see” into the past through interaction with the
concentration gradient diminishes exponentially in time. This unique behavior
abolishes time-reversal symmetry although the coupled configuration in Eq.(1)
is Markovian since there is no explicit dependence on the trajectory’s past
steps.

To limit the number of parameters under investigation, we nondimensionalize
Eq. (1). We choose $R$ as a natural length scale and non-dimensionalize $c$
without any scaling for simplicity. Temporarily, we leave time scale $T$
arbitrary. Under the scalings $\mathbf{y}=\frac{\mathbf{x}}{R}$,
$\mathbf{Y}=\frac{\mathbf{X}}{R}$, $\tau=\frac{t}{T}$,and
$\mathbf{B}=\frac{\mathbf{W}}{\sqrt{T}}$, we arrive at

$\begin{split}\partial_{t}c(\mathbf{y},t)&=\mu\Delta
c(\mathbf{y},t)+\mu\phi\exp\left[-\frac{|\mathbf{y}-\mathbf{Y}(t)|^{2}}{2}\right],\end{split}$
(2a) for $\mathbf{y}\in\Omega,\quad t\geq 0$, and
$\begin{split}d\mathbf{Y}(t)&=-\nu\left(\int_{\Omega}\exp\left[-\frac{|\mathbf{y}-\mathbf{Y}(t)|^{2}}{2}\right]\nabla_{\mathbf{y}}c(\mathbf{y},t)d\mathbf{y}\right)dt\\\
&\quad\quad\quad\quad\quad\quad+\sqrt{\epsilon}d\mathbf{B}(t)\end{split}$ (2b)

where $c$, $t$ and $\Omega$ are re-used for their non-dimensional versions for
convenience. We have mapped the dimensional parameters as follows:
$D\rightarrow\mu=\frac{DT}{R^{2}}$, $\alpha\rightarrow\phi=\frac{\alpha
R^{2}}{2\pi}$, $\beta\rightarrow\nu=\frac{\beta T}{2\pi R}$, and
$\sigma\rightarrow\epsilon=\frac{\sigma T}{R^{2}}$.

We note that a typical time scale for the diffusion equation is
$T=\frac{L^{2}}{D}$. Although traditional, this choice would prevent us from
seeing directly the effects of changing $D$, which encodes the memory
timescale. Increasing diffusivity would contract time such that the past-
history wake of the particle would adjust to decay at the same rate. Thus, to
observe the effects of this memory-encoding diffusivity, we choose to fix the
stochastic diffusivity $\epsilon=0.75^{2}$, thereby choosing
$T=0.75^{2}\frac{R^{2}}{\sigma}$, and focus this study on the deterministic
parameters $\mu$ and $\nu$ relative to the noise.

We can simplify the system by taking the Fourier transform and solving the PDE
(2a) on an infinite domain, $\Omega=\mathbb{R}^{2}$, explicitly. Incorporating
this solution into the SDE (2b) we arrive at the mathematically equivalent
system for the particle in an infinite domain

$\begin{split}d\mathbf{Y}&=\frac{\pi}{2}\mu\nu\phi\int_{0}^{t}\exp\left[-\frac{|\mathbf{Y}(t)-\mathbf{Y}(s)|^{2}}{4(1+\mu(t-s))}\right]\\\
&\cdot(1+\mu(t-s))^{-2}(\mathbf{Y}(t)-\mathbf{Y}(s))dsdt+\sqrt{\epsilon}d\mathbf{B};\end{split}$
(3)

see App. A.1 for details. This non-Markovian SDE explicitly reveals the
dependence on all the particle’s previous locations via integration in time;
the exponential kernel decays in both time and space, revealing $\mu^{-1}$ to
be the self-avoidant memory timescale.

This formulation is convenient for simulation since it does not require
boundary conditions and solving the PDE both removes computations requiring
the Laplacian and transforms the SDE integral from $\mathbb{R}^{2}$ to $t$. We
integrate in time with the Euler-Maruyama method while using Simpson’s rule to
integrate the memory kernel at each step. This algorithm is a first order
method.

The limiting behavior of these two equivalent systems, Eq. (2) and (3),
foreshadows their distinction from the active Brownian model since it reveals
that removing the distinguishing feature of memory by taking $D\to\infty$ will
not reduce our model to ABP. The parameter $D$ was added to the source term in
Eq. (1a) to achieve balance between the rate at which the oil diffuses and the
rate at which the oil is expelled in this limit. (If instead the source term
remained constant relative to $D$, then it would effectively vanish in the
limit of $D\to\infty$.) The nondimensional parameter $\mu$ therefore appears
on the source term in Eq. (2a), and $\partial_{t}c(\mathbf{y},t)\to\infty$ as
$\mu\to\infty$. To leading order, the concentration field satisfies the
Poisson equation

$\Delta
c(\mathbf{y})=-\phi\exp\left[-\frac{|\mathbf{y}-\mathbf{Y}(t)|^{2}}{2}\right]\qquad\mu\to\infty.$
(4)

The concentration field is now memory-less since it instantaneously
equilibrates as the particle moves. On an infinite domain, the solution to Eq.
(4) will be radially symmetric around the particle’s location, and therefore
the integral in Eq. (2b) will always be zero. As a result, particles
experience motility solely from thermal fluctuation, namely simple Brownian
motion.

Also noteworthy is the “full memory” limit of Eq. (2a) which is
$\partial_{t}c\to 0$ as $\mu\to 0$. The concentration field remains fixed at
its initial conditions as the source term and the diffusion term vanish in
this limit. The particle experiences thermal fluctuations while responding to
the concentration gradient of the fixed concentration landscape, thereby
statistically preferring concentration minima. The steady state (if one
exists) would be almost solely determined by the initial topography of
$c(\mathbf{y})$ and the relative size of $\epsilon$. In fact, the entire
coupled system in Eq. (2) reduces to simple memoryless Brownian motion

$d\mathbf{Y}=\sqrt{\epsilon}d\mathbf{B}$ (5)

in this limit, as the source term which encodes the memory vanishes. This is
consistent with Eq. (3) which also reduces to simple Brownian motion in the
same limit $\mu\to 0$ when it is assumed that the initial concentration field
is constant. Therefore we focus our study on intermediate range of $\mu$ where
the effects of the noise, the swimming, and the effects of memory are all
observable. The limiting behaviors of our model as described above can all be
traced back to the addition of a second mu on the source term. As stated
previously, inclusion of a diffusive scaling on the source term was required
to ensure that the source term remained in the limit as $D\rightarrow\infty$,
in the hope that the model would revert to active Brownian motion with no
self-avoidant memory. In our results, we discuss the role of this diffusive
scaling in generating previously unseen behaviors in this class of models.

##  III Comparative Analysis

Numerically simulated trajectories of the coupled model given by Eq. (2) are
shown in Fig. 1. These trajectories illustrate the main features of active
matter: a swimming velocity with a slowly diffusing direction. Increasing
$\nu$, and therefore the response to the concentration gradient causes the
particle to swim faster, shown in Fig. 1a, while increasing $\mu$, and
therefore shortening the timescale of the diffusion (decreasing the memory),
has a secondary effect on the velocity, but also causes the particles to turn
faster, shown in Fig. 1b. An increase in turning frequency was also seen
experimentally in [12] as surfactant concentration was increased, which
prompted a transition from ballistic motion to diffusion. (See Figure 2 in
[12].)

We seek to look beyond the combined effects of swimming and random directional
changes in producing enhanced diffusion and understand the additional effects
of self-avoidant memory. Specifically, we compare our model to ABP given by
the equations

$dX=V\cos(\theta(t))dt+\sqrt{\epsilon}dW_{x}$ (6a)
$dY=V\sin(\theta(t))dt+\sqrt{\epsilon}dW_{y}$ (6b)
$d\theta=\frac{1}{\sqrt{\tau}}dW_{\theta}$ (6c)

where $\sqrt{\epsilon}$ is the strength of the additive noise in each spatial
component (consistent with the model in Eq. (3)), $V$ is the swimming
velocity, and $\tau$ is the persistence timescale of the rotational diffusion
[12, 27, 30]. These latter two parameters do not explicitly appear in our
model; we will compute them and compare the MSD of the two models to
understand the effects of self-avoidant memory on enhanced diffusion.

In Sec. III.1 we present an analytic equation that is numerically solved for
the velocity of the swimming solution to Eq. (3). This velocity is consistent
with the intermediate ballistic regime of the MSD, computed numerically for
Eq. (3) and given by

$\mathbb{E}[\mathbf{X}(t)^{2}]=4V^{2}\tau^{2}\left[2\left(e^{-\frac{t}{2\tau}}-1\right)+\frac{t}{\tau}\right]+2\epsilon
t$ (7)

for Eq. (6). In Sec. III.2 we determine $\tau$ by numerically computing the
orientation correlation function but find that the memory induced from
modifying the environment causes a reduced effective diffusion as compared to
ABP with identical angular persistence. Given that in the experimental
systems, the angular persistence is a byproduct of local hydrodynamic effects
while in our model the angular reorientation is a globally-derived effect of
the past spatial history contained in the concentration field, it is perhaps
not surprising that active Brownian motion accurately predicts the enhanced
diffusion of the experimental system, while here it does not.

Figure 1: Numerical simulations of Eq. (2) were carried out using a forward-
time centered-space finite difference scheme for the PDE and the Euler-
Maruyama method for integrating the SDE where the trapezoid rule was used to
compute the integral therein. We confine $\mathbf{X}(t)$ to a box
$B=\\{x\in(-5,5),y\in(-5,5)\\}$ with insulating boundary conditions such that
$\frac{\partial c(\mathbf{x},t)}{\partial
t}=0\quad\forall\mathbf{x}\in\\{\partial B\\}$. This has the effect of
reflecting the particle back into the domain when it reaches the boundary. The
initial condition is $c(\mathbf{x},0)=0$. Note that these trajectories are
visually indistinguishable from active Brownian motion. (a) For $\phi=1$,
$\mu=5$ and a noise level of $\sqrt{\epsilon}=0.75$, we see the dominant
effect of $\nu$ which is to increase the velocity. (b) For $\phi=1$, $\nu=7$,
and a noise level of $\sqrt{\epsilon}=0.75$, we see the effects of $\mu$ which
primarily increases the turning frequency and has a secondary effect on
velocity.

###  III.1 Intermediate Time Scale Behaviors: Ballistic Motion

Unlike active Brownian models, the proposed model has a non-explicit intrinsic
velocity; directed motion at this velocity may become observable at
intermediate time scales under appropriate conditions for $\nu$ and $\mu$. To
find an analytic form for the velocity, we seek a deterministic constant
velocity (“steady state”) solution to the combined model Eq. (3). Without loss
of generality, suppose $\mathbf{Y}(t)=\langle Vt,0\rangle$; thus $V$ must
solve

$\begin{split}\frac{d\mathbf{Y}}{dt}=V=&\frac{\pi}{2}\mu\nu\phi\int_{-\infty}^{t}\exp\left[-\frac{|Vt-
Vs|^{2}}{4(1+\mu(t-s))}\right]\\\ &\cdot(1+\mu(t-s))^{-2}(Vt-
Vs)ds.\end{split}$ (8)

Under the transformation $z=\mu(t-s)$, the constant velocity $V$ therefore
satisfies

$1=\frac{\pi}{2}\frac{\nu}{\mu}\phi\int_{0}^{\infty}\frac{z}{(1+z)^{2}}\exp\left[-\left(\frac{V}{\mu}\right)^{2}\frac{z^{2}}{4(1+z)}\right]dz.$
(9)

For each value of $\frac{\nu}{\mu}$, we solve for the value of $\frac{V}{\mu}$
that makes the above integral equal to 1 numerically with the `ptimize.r`ot
function in the Python

    cipy  library.
    Under the change of variable

$x=\frac{2}{\pi}\arctan(z)$, we map the domain $(0,\infty)$ to $(0,1)$ for
ease of numerical integration. The resulting monotonically increasing
dependence of $\frac{V}{\mu}$ on $\frac{\nu}{\mu}$ is plotted as the solid
black line in Fig. 2b.

We can directly compare this theory to ABP on the timescale at which ballistic
motion is dominant. It is evident from Fig. 2a that the ballistic portion of
the simulated MSD aligns with the computed velocity from Eq. (9). At such
small times, the MSD of ABP asymptotically reduces to

$\mathbb{E}[\mathbf{X}(t)^{2}]\sim V^{2}t^{2}+2\epsilon t$ (10)

as $t\rightarrow 0$ (see App. A.3 for details.) Fitting $V$ from the ballistic
portion of the MSD of our particles is also in good agreement with the theory,
as shown in Fig. 2b.

We point out that the existence of an observable ballistic regime in the MSDs
from our model requires a sufficient swimming velocity $V$ to dominate the
additive noise. In the ABP model, this can be guaranteed by changing the
stated parameter $V$, whereas in our model, there must be consideration for
the parameters $\mu$ and $\nu$ due to the explicit functional relationship
$V=f(\mu,\nu)$ given by Eq. (9). To see this functional relationship more
clearly, the contours of constant velocity are plotted in Fig. 3a and the
contours of constant $\nu$ are plotted in Fig. 3c. These figures agree with
the limits from Sec. II in that $V=f(\mu,\nu)\rightarrow 0$ when taking either
$\mu\to 0$ or $\mu\to\infty$ with fixed $\nu$, and the model system Eq. (2) or
Eq. (3) reduces the particle motion to simple Brownian motion. Furthermore,
taking $V\to 0$ in Eq. (8) results in a divergent integral; for the integral
to converge, either $\mu$ or $\nu$ in the prefactor must also go to zero. The
result is no transition to swimming at a small finite values of these
parameters. Similarly, the integral also diverges as $\mu\to\infty$. Figure 3c
most clearly shows the relevant intermediate values of $\mu$ for which a
significant velocity exists and ABP-like motion with a ballistic regime is
expected for the model system.

Figure 3a most clearly shows that for any given velocity, there exists a
minimal $(\mu,\nu)$ pair. If we interpret this in the context of the coupled
model given by Eq. (2), it suggests the existence of an optimal response
strength and diffusivity pairing which act on the particle to produce directed
motion at a specific speed. Moving off of this minimum illustrates the
parameter couplings which must balance to keep the particle moving at a given
speed. For example, decreasing memory (increasing $\mu$) allows the particle’s
trail to diffuse faster which weakens local gradients, and thus requires that
the response strength to the weakened gradient be increased (increasing
$\nu$).

Figure 2: (color online) (a) Empirically computed MSDs of the model in
comparison to benchmark pure ballistic motion. Three distinct regions of the
MSD are visible: classical diffusion, directed motion in alignment with the
benchmark ballistic lines, and enhanced diffusion where the MSD pulls away
from the ballistic motion. Legend colors are consistent with other panels
throughout the rest of the figure. Legend colors correspond to all panels. (b)
Collapsed model velocity theory from solving Eq. (9) in comparison to fitted
values extracted from the ballistic regime of the MSD of the form
$V^{2}t^{2}$. (c) Model OCFs (noisy solid lines) with fitted $\tau$ (solid
smooth curves) with comparison to $\tau$ fitted from the enhanced diffusion
regime of the MSD consistent with ABP Eq. (11) (dashed smooth curves). Values
of $\tau$ for each corresponding color are given in the legend. (d) Both
fitted $\tau$ values from ABP MSD (dashed) and OCF of our model (solid) as
they vary over $\mu$. (e) Model MSD (solid) and ABP MSD Eq. 7 (dashed) under
the same $\tau$ fitted to the model OCF and theoretical $V$. (f) Signed
difference between model and ABP MSDs in panel (e).  Figure 3: (color online)
Visualizations of the four dimensional parameter space of $V$ from solving Eq.
(9), $\tau$ from fitting the model OCF, along with model parameters $\nu$ and
$\mu$. (a) and (b) depict constant velocity $V$ contours, while in (c) and (d)
constant $\nu$ contours are shown. Panels (b) and (d) show sub-regions of
panels (a) and (c), respectively, with added individual points depicting
values of $\tau$ given by the color bar. Note that $\tau$ values larger than
the experiment length (greater than 30) are colored black.

###  III.2 Long Time Scale Behaviors: Enhanced Diffusion

Figure 2a shows a departure of the MSD from ballistic motion at longer time
scales. For ABP, this departure happens at time scales $t\gg\tau$ for which
the MSD (7) is asymptotic to

$\mathbb{E}[\mathbf{X}(t)^{2}]\sim(4V^{2}\tau+2\epsilon)t$ (11)

as $t\rightarrow\infty$. Particle reorientations that decorrelate with
timescale $\tau$ enhance the diffusion term $2\epsilon t$ with the term
$4V^{2}\tau$.

To estimate $\tau$ from trajectories given by our model we numerically compute
the normalized orientation correlation function (OCF)

$C(\Delta t)=\left<\cos\left(\frac{\mathbf{v}(t)\cdot\mathbf{v}(t+\Delta
t)}{|\mathbf{v}(t)||\mathbf{v}(t+\Delta t)|}\right)\right>_{t}$ (12)

where $\mathbf{v}(t)=\mathbf{Y}(t)-\mathbf{Y}(t-\Delta t)$ gives the
directional displacement between times $t$ and $t-\Delta t$ (see App. A.4 for
details). The portion of this function displaying exponential decay is fit by
a single exponential given by $C(\Delta t)=e^{-\frac{\Delta t}{\tau}}$ as is
consistent with ABP [12, 30]. (Note that this correlation $C(\Delta t)\to 0$
as $\Delta t\to 0$ due to the additive uncorrelated fluctuations in the model;
OCF values near 1 are indicative of directed motion.) These fits are shown by
the smooth solid lines in Fig. 2c and the resulting values of $\tau$ as a
function of $\mu$ in Fig. 2d.

While the model OCF is well fit by an exponential decay, the long time
asymptotics of the MSD given by Eq. 11 in Fig. 2e (dashed lines) shows that
ABP substantially overestimates the enhanced diffusion of our model (solid
lines). This overestimation is larger for the two larger values of $\mu$ that
correspond to weaker self-avoidant memory, as shown in Fig. 2f. Alternatively,
using $V$ from Eq. (9), $\tau$ is determined by fitting the long time MSD to
Eq. (11). These values of $\tau$, plotted as the dashed lines in Fig. 2d,
substantially underestimate the decorrelation time scale of our model, also
shown by the corresponding dashed lines of exponential decay in Fig. 2c.
Although the form of exponential decay of the orientational persistence is
consistent with ABP and quantifiable by $\tau$, it alone is not enough to
predict the enhanced diffusion of our model. There are additional effects of
self-avoidant memory beyond the persistence memory, which is the only memory
present in ABP.

At constant velocity, the effect of increasing self-avoidant memory
(decreasing $\mu$) is seen in Fig. 2a. The MSDs with smaller $\mu$ in both
cases depart from the ballistic regime earlier, and thus exhibit less enhanced
diffusion. This corresponds to Fig. 2d where for smaller $\mu$ the OCF decays
more rapidly as measured by a smaller value of $\tau$. This is further
illustrated in Fig. 3b, showing a decrease of $\tau$ with decreasing $\mu$
along contours of constant velocity. For fixed $\mu$, $\tau$ increases with
decreasing $\nu$ (although velocity decreases). Thus we see that one effect of
self-avoidant memory as it is present in our model is to decrease
orientational persistence: swimmers with high memory experience weak
orientational persistence and vice versa.

A more exotic effect of self-avoidant memory is shown by the trajectories in
Fig. 4 and provides a plausible explanation for the surprising fact that ABP
overestimates the enhanced diffusion of the model. To avoid crossing their own
self history, paths turn back on themselves and continue turning inward,
becoming caged for a while before enough diffusion has occurred for them to
leave this self-created trap. This transient self-trapping perhaps explains
the reduced enhanced diffusion as compared to ABP with equivalent orientation
persistence time. Self-trapping has been studied for autophoretic systems like
that which we model, although it has only been found in chemoattractant
experimental systems [31, 32] and model systems [33, 34, 35] having self-
attracting memory.

Figure 4: With $\mu=0.01$ and $V=6$, four sample paths are shown which
illustrate the caging of enhanced diffusion experienced due to high memory.

##  IV Limiting Behavior

Since the relevant experimental systems are well described by ABP and we can
explicitly tune memory in our model, we anticipated that we could find a
parameter regime (low memory, high $\mu$) in which the enhanced diffusion
produced by our model was well described by ABP. However, as discussed at the
end of Sec. II, the limiting behavior of systems Eq. 2 and Eq. 3 as
$\mu\to\infty$ is simple Brownian motion, indicating that enhanced diffusion
with low self-avoidant memory may not be possible. We revisit this limit with
further simulations in light of the emergent parameters $V$ and $\tau$,
considering both $\nu\to\infty$ with $V$ fixed as well as dynamic $V$ with
$\nu$ fixed. Additionally, we investigate the high memory limit ($\mu\to 0$),
and find it consistent with Eq. (5) describing classical Brownian motion with
(unfortunately) no further memory effects to investigate.

From Figs. 3a and b, we can consider the limit most likely to produce enhanced
diffusion consistent with ABP: removing the memory via the limit
$\mu\rightarrow\infty$ while keeping the particle at constant velocity by
fixing $V$. Visually we observe that as $\mu\rightarrow\infty$, the contours
of $V$ become flatter, reproducing the behavior seen in Fig. 1 which shows
that ballistic motion is sensitive to changes in the gradient response $\nu$.
Moreover, remaining on one $V$ contour requires $\nu\to\infty$ much slower
than $\mu\to\infty$. To investigate the enhanced diffusion in this limit, we
look at $\tau$ in Fig. 3b. Following a $V$ contour as $\mu\to\infty$ results
in an increase in $\tau$ corresponding to longer orientational persistence (or
less change in direction).

As a result of decreasing the self-avoidant memory timescale ($\mu\to\infty$
while maintaining a constant velocity $V$), we find that both the past history
of the trajectory and the Brownian noise become less important in influencing
the future location of the trajectory. Furthermore, with the addition of an
increased gradient response by taking $\nu\to\infty$ (as required to keep the
particle at constant $V$), the deterministic gradient response force in Eq.
(2) dominates the Brownian noise, and consecutive steps become more
correlated. This increases the persistence time $\tau$, and the trajectories
approach purely ballistic motion with no enhanced diffusion at observable
finite times.

Figs. 3a and b also allow for considering infinite memory ($\mu\rightarrow
0$), again with constant velocity. In Fig. 3a, we see that the gradient
response required (given by the size of $\nu$) to keep the particle moving at
constant velocity $V$ rapidly blows up to $\infty$. This largely unsurprising
as the prefactor on the deterministic term in Eq. (3) contains the product
$\mu\nu$; taking $\mu\to 0$ while keeping this integral response factor
relatively constant would necessitate $\nu\to\infty$. Fig. 3b shows a
corresponding decrease in $\tau$, limiting towards pure diffusion. Returning
to Eq. (2a), as $\mu\to 0$ both the diffusion and the source term go to zero,
thus the concentration field would remain fixed in time. If this initial
concentration field was constant, then the particle would have no gradient to
respond to and therefore only undergo pure Brownian motion in this infinite
memory regime, corresponding to $\tau=0$. This suggests that rather than
trying to start at finite $\mu$ and witness the effects of self-avoidant
memory fade as $\mu\to\infty$, as this model was set up to do, future work
should remove $\mu$ from the source term in Eq. (2a) and start at $\mu=0$ to
witness the effects of self-avoidant memory fade as $\mu$ increases away from
zero.

The limiting ballistic motion when taking both $\mu$ and $\nu$ to infinity is
in contrast to the limiting Brownian motion behavior of Eq. (3) as
$\mu\to\infty$ while keeping $\nu$ fixed. By following contours of $\nu$ in
Fig. 3c, we see that the velocity first increases with $\mu$ and then
decreases, approaching zero velocity as $\mu\to\infty$, which is consistent
with the trend shown in Fig. 1b. It is interesting to observe in Fig. 3d that
$\tau$ appears to be relatively static along the contours of $\nu$. Note that
these values of $\tau$ were mainly computed at points to the left of the
maximum velocity of the fixed $\nu$ contours. Figure 1b indicated that $\tau$
decreases with increasing $\mu$ and fixed $\nu$. When we also consider a
decrease velocity, the trajectories can be assumed to approach Brownian
motion.

In summary, increasing $\mu$ to decrease the effects of memory either results
in increasing $\tau$ (by fixing $V$) and therefore creating straighter
trajectories that do not display enhanced diffusion in the MSD over the
timescale of the simulation, or in decreasing $V$ to zero (by fixing $\nu$)
which results in a purely diffusive MSD. The memory is responsible for both
the ballistic motion measured by velocity $V$ and the effective rotational
diffusion measured by orientational persistence time $\tau$, so naturally
follows that these effects are both lost with increasing $\mu$. If the
concentration field diffuses infinitely fast by taking $\mu\rightarrow\infty$
with $\nu$ fixed, we lose deterministic motion since the gradient of the
concentration field is always zero at the particle’s center and radially
isotropic around the particle, thus the net force acting on the particle is
always zero. The effective angular diffusion is lost when taking
$\mu\rightarrow\infty$ with $V$ fixed because this requires large $\nu$ such
that the immediate deterministic forces overwhelm the noise and any past
history, and so reduce the MSD to almost exclusively ballistic motion. Thus,
incorporation of self-avoidant memory is not simply an addendum to the active
Brownian model that can be removed without consequence; by its complex
interactions with the enhanced diffusion we see that it makes for a
categorically unique model.

##  V Conclusions

We have analyzed the self-avoidant memory effects of a model coupling an
active swimmer and an environmental chemical field. Like the experimental
system it was inspired by, it can exhibit ABP-like behavior with the MSD
having both a ballistic and a long-time enhanced diffusion regime. With an
analytical formula for the velocity, $V$, we faithfully reproduced the
ballistic regime. The enhanced diffusion in our model is a result of both
angular persistence and the self-avoidant memory, whereas ABP only includes
orientational persistence. We found that numerically computing the orientation
decorrelation (or persistence) time, $\tau$, enhanced diffusion predicted by
ABP overestimates the enhanced diffusion in our model. Thus, our proposed
model did not faithfully capture the dynamics of the experimental system at
long time scales in the same way that ABP did. (We hypothesize this is due to
the origin of the noise, which is hydrodynamic and not thermal in the
experimental system.) Instead, we discovered that the self-avoidant memory in
our model led to transient self-trapping that suppressed the enhanced
diffusion. This self-trapping has, to date, been suggested to occur only in
self-attracting experimental systems [31, 32] and computational models [33,
34, 35]. Determining if the addition of the diffusive scaling to the source
term or if the modeling of hydrodynamic effects as thermal noise produced this
new behavior will be left to future work.

Through these investigations, we kept the noise parameter $\epsilon$ fixed,
while changing the gradient response parameter $\nu$ and the diffusion $\mu$
to find that both latter parameters control the implicit parameters $V$ and
$\tau$. With only two control parameters, we were unable to independently tune
each timescale of behavior: the velocity $V$, the memory timescale $\mu^{-1}$,
and the angular persistence timescale $\tau$. Taking $\mu\to\infty$ to remove
memory effects, we either arrived at simple Brownian motion by fixing $\nu$ or
purely ballistic motion by fixing $V$ and allowing $\nu\to\infty$; the memory
is responsible for both the ballistic motion and the effective rotational
diffusion. Taking $\mu\to 0$, we again arrive at simple Brownian motion,
having removed all self-avoidant memory with our choice of scaling the source
term in the concentration field by $\mu$. We thereby identified an
intermediate regime of $\mu$ for which enhanced diffusion is present on a
finite timescale, but at a lower magnitude than expected for ABP with
equivalent angular persistence. This regime will be used in future work to
study self-avoidant memory effects in many-particle simulations, investigating
motility induced phase separation and associated dynamic pattern formation,
which is commonly observed in active systems with particles that are repulsive
to one another.

###### Acknowledgements.

KAN thanks Eric Vanden-Eijnden for initial discussions concerning the model.

##  Appendix A Supplementary Material

###  A.1 Solving the Diffusion Equation to Combine and Nondimensionalize the
Coupled System

Consider the dimensional, 2D system given in Eq. (1). Taking the Fourier
Transform of Eq. (1a) we arrive at the ODE

$\hat{c}_{t}+D(\mathbf{k}^{2})\hat{c}=\frac{\alpha
DR^{2}}{2\pi}\exp\left[-\frac{R^{2}}{2}\mathbf{k}^{2}+i(\mathbf{k}\cdot\mathbf{X}(t))\right].$
(13)

We compute the integrating factor of Eq. 13 which is $I=e^{\int
d\overline{k}^{2}dt}=e^{Dt\overline{k}^{2}}$. From this Eq. 13 can be
rewritten as

$\frac{d}{dt}\left(e^{Dt\mathbf{k}^{2}}\hat{c}\right)=\frac{\alpha
DR^{2}}{2\pi}e^{Dt\mathbf{k}^{2}}e^{\left[-\frac{R^{2}}{2}\mathbf{k}^{2}+i(\mathbf{k}\cdot\mathbf{X}(t))\right]}.$
(14)

Integrating both sides of Eq. 14 gives the solution to Eq. 13

$\hat{c}=\frac{\alpha
DR^{2}}{2\pi}\int_{0}^{t}e^{-D(t-s)\mathbf{k}^{2}}e^{\left[-\frac{R^{2}}{2}\mathbf{k}^{2}+i(\mathbf{k}\cdot\mathbf{X}(s))\right]}ds.$
(15)

Taking the inverse Fourier Transform of Eq. 15 yields the solution,
$c(\mathbf{x},t)$, to Eq. (1a), which is

$c=\frac{\alpha DR^{2}}{2\pi}\int_{0}^{t}(R^{2}+2D(t-s))^{-1}\cdot
e^{-\frac{|\mathbf{x}-\mathbf{X}(s)|^{2}}{2(R^{2}+2D(t-s))}}ds.$ (16)

We can incorporate the solution to Eq. (1a), which is Eq. 16, into Eq. (1b) by
taking the gradient, $\nabla c$, which is

$\nabla c=-\frac{\alpha
DR^{2}}{2\pi}\int_{0}^{t}(R^{2}+2D(t-s))^{-2}(\mathbf{x}-\mathbf{X}(s))\cdot
e^{-\frac{|\mathbf{x}-\mathbf{X}(s)|^{2}}{2(R^{2}+2D(t-s))}}ds.$ (17)

The SDE path evolution equation Eq. (1b) then becomes

$d\mathbf{X}=\frac{\alpha D\beta
R}{(2\pi)^{2}}\left[\int_{0}^{t}(R^{2}+2D(t-s))^{-2}\left(\int_{\mathbb{R}^{2}}(\mathbf{x}-\mathbf{X}(s))\cdot
e^{-\frac{|\mathbf{x}-\mathbf{X}(s)|^{2}}{2(R^{2}+2D(t-s))}-\frac{|\mathbf{x}-\mathbf{X}(t)|^{2}}{2R^{2}}}d\mathbf{x}\right)ds\right]dt+\sqrt{\sigma}d\mathbf{W}.$
(18)

Evaluation of the spatial integral over $\mathbb{R}^{2}$ reduces Eq. 18 to

$d\mathbf{X}=\frac{\alpha D\beta
R^{3}}{2^{3}\pi}\left[\int_{0}^{t}\left((\mathbf{X}(t)-\mathbf{X}(s))e^{-\frac{|\mathbf{X}(t)-\mathbf{X}(s)|^{2}}{4(R^{2}+D(t-s))}}(R^{2}+D(t-s))^{-2}\right)ds\right]dt+\sqrt{\sigma}d\mathbf{W}.$
(19)

By nondimensionalizing under the scalings $\mathbf{y}=\frac{\mathbf{x}}{R}$,
$\mathbf{Y}=\frac{\mathbf{X}}{R}$, $\tau=\frac{t}{T}$,
$E=\frac{c}{\overline{c}}$ with $\overline{c}=1$ and
$\mathbf{B}=\frac{\mathbf{W}}{\sqrt{T}}$, Eq. 18 becomes

$Rd\mathbf{Y}=\frac{\alpha D\beta
R^{3}}{2^{3}\pi}\left[\int_{0}^{\tau}\left(e^{-\frac{|R\mathbf{Y}(\tau)-R\mathbf{Y}(\zeta)|^{2}}{4(R^{2}+DT(\tau-\zeta))}}(R\mathbf{Y}(\tau)-R\mathbf{Y}(\zeta))(R^{2}+DT(\tau-\zeta))^{-2}\right)Td\zeta\right]Td\tau+\sqrt{\sigma
T}d\mathbf{B}.$ (20)

The SDE path evolution given by Eq. 20 then simplifies to

$d\mathbf{Y}=\frac{\alpha D\beta
R^{3}T^{2}}{2^{3}\pi}\left[\int_{0}^{\tau}\left(e^{-\frac{|R\mathbf{Y}(\tau)-R\mathbf{Y}(\zeta)|^{2}}{4(R^{2}+DT(\tau-\zeta))}}(\mathbf{Y}(\tau)-\mathbf{Y}(\zeta))(R^{2}+DT(\tau-\zeta))^{-2}\right)d\zeta\right]d\tau+\frac{\sqrt{\sigma
T}}{R}d\mathbf{B}$ (21)

where we also introduce the nondimensional dummy time variable
$\zeta=\frac{s}{T}.$ Since random variables are by definition dimensionless, a
random variable $Z(t)$ is not a function of $t$; it is a draw from the
distribution of $Z$ at time $t$. We want to then be able to draw from the same
random variable at nondimensional time $\tau$, thus we map
$\mathbf{X}(t)\rightarrow\mathbf{U}(\tau)$. Similarly, we want to evaluate the
now-dimensionless integrand over its associated dimensionless time interval,
so we set the upper limit of integration to $\tau$.

Incorporating the nondimensional parameters
$D\rightarrow\mu=\frac{DT}{R^{2}}$, $\alpha\rightarrow\phi=\frac{\alpha
R^{2}}{2\pi}$, $\beta\rightarrow\nu=\frac{\beta T}{2\pi R}$, and
$\sigma\rightarrow\epsilon=\frac{\sigma T}{R^{2}}$ and exchanging $s$ for
$\zeta$ and $t$ for $\tau$ for notational convenience we have the
nondimensional SDE path evolution equation

$d\mathbf{Y}=\frac{\pi}{2}\mu\nu\phi\left[\int_{0}^{t}\left(e^{-\frac{|\mathbf{Y}(t)-\mathbf{Y}(s)|^{2}}{4(1+\mu(t-s))}}(\mathbf{Y}(t)-\mathbf{Y}(s))(1+\mu(t-s))^{-2}\right)ds\right]dt+\sqrt{\epsilon}d\mathbf{B}$
(22)

in agreement with Eq. 3.

###  A.2 Computation of the Velocity Integral Formulation Using a Dirac Delta
Function

To assess the case in which the particle is considered a point source, we
substitute the mollified delta function,
$\delta_{R}(\mathbf{x}-\mathbf{X}(t))=\frac{1}{2\pi
R^{2}}e^{-\frac{|\mathbf{x}-\mathbf{X}(t)|^{2}}{2R^{2}}}$, in Eq. 2 for a
Dirac delta function,

$\frac{\partial c}{\partial t}=D\Delta c+\alpha
D\delta^{2}(\mathbf{x}-\mathbf{X}(t))$ (23a) $d\mathbf{X}(t)=-\beta
R\left(\int_{\Omega}\delta^{2}(\mathbf{x}-\mathbf{X}(t))\nabla
cd\mathbf{x}\right)dt+\sqrt{\sigma}d\mathbf{W}.$ (23b)

Here, $\delta^{2}(\mathbf{x}-\mathbf{X}(t))$ is a 2-dimensional Dirac delta
function centered at $\mathbf{X}(t)$. The $R^{2}$ in the source term of the
original PDE given by Eq. 1a is no longer necessary. Accordingly, the units of
$\alpha$ are $[\alpha]=c$ and the units of $\beta$ remain
$[\beta]=\frac{L}{cT}.$

Nondimensionalizing Eq. 23 with the scalings
$\mathbf{y}=\frac{\mathbf{x}}{R}$, $\mathbf{Y}=\frac{\mathbf{X}}{R}$,
$\tau=\frac{t}{T}$, and $\mathbf{B}=\frac{\mathbf{W}}{\sqrt{T}}$ and where
$\mu=\frac{DT}{R^{2}}$, $\phi=\frac{\alpha}{2\pi}$, $\nu=\frac{\beta
T}{R2\pi}$ and $\epsilon=\frac{\sigma T}{R^{2}}$ we arrive at the new system

$\frac{\partial c}{\partial t}=\mu\Delta
c+2\pi\mu\phi\delta^{2}(\mathbf{y}-\mathbf{Y}(t))$ (24a)
$d\mathbf{Y}(t)=-2\pi\nu\left(\int_{\Omega}\delta^{2}(\mathbf{y}-\mathbf{Y}(t))\nabla
cd\mathbf{y}\right)dt+\sqrt{\epsilon}d\mathbf{B}$ (24b)

where $c$, $t$ and $\Omega$ are re-used for their non-dimensional versions for
convenience.

As in the case with the sized particle, we take the Fourier Transform of the
Eq. 24a to arrive at the ODE

$\hat{c}_{t}+\mu\mathbf{k}^{2}\hat{c}=\mu\phi
e^{i\mathbf{k}\cdot\mathbf{Y}(t)}.$ (25)

We compute the integrating factor of Eq. 25 which is
$I=e^{\int\mu\mathbf{k}^{2}dt}=e^{\mu t\mathbf{k}^{2}}$. From this, Eq. 25 can
be rewritten as

$\frac{d}{dt}\left(\hat{c}\cdot e^{\mu t\mathbf{k}^{2}}\right)=\mu\phi
e^{i\mathbf{k}\cdot\mathbf{Y}(t)}\cdot e^{\mu t\mathbf{k}^{2}}.$ (26)

Integrating both sides of Eq. 26 gives

$\hat{c}=\mu\phi\int_{0}^{t}e^{-\mu(t-s)\mathbf{k}^{2}+i\mathbf{k}\cdot\mathbf{Y}(s)}ds.$
(27)

We take the inverse Fourier Transform of Eq. 27 to find the solution to Eq.
24a, which is

$c=\mu\phi\int_{0}^{t}\left(2\mu(t-s)\right)^{-1}e^{-\frac{|\mathbf{y}-\mathbf{Y}(s)|^{2}}{4(\mu(t-s))}}ds.$
(28)

We can incorporate the solution to Eq. 24a into Eq. 24b by computing the
gradient $\nabla c$ of Eq. 28, which is

$\nabla
c=-\mu\phi\int_{0}^{t}\frac{(\mathbf{y}-\mathbf{Y}(s))}{(2\mu(t-s))^{2}}\exp\left[{-\frac{|\mathbf{y}-\mathbf{Y}(s)|^{2}}{4(\mu(t-s))}}\right]ds.$
(29)

Eq. 24b then becomes

$d\mathbf{Y}(t)=\nu\mu\phi\frac{\pi}{2}\int_{0}^{t}(\mu(t-s))^{-2}\int_{\Omega}\delta^{2}(\mathbf{y}-\mathbf{Y}(t))(\mathbf{y}-\mathbf{Y}(s))\exp\left[{-\frac{|\mathbf{y}-\mathbf{Y}(s)|^{2}}{4(\mu(t-s))}}\right]d\mathbf{y}dsdt+\sqrt{\epsilon}d\mathbf{B}.$
(30)

Evaluation of the spatial integral over $\mathbb{R}^{2}$ reduces (30) to

$d\mathbf{Y}(t)=\nu\mu\phi\frac{\pi}{2}\int_{0}^{t}(\mu(t-s))^{-2}(\mathbf{Y}(t)-\mathbf{Y}(s))\exp\left[{-\frac{|\mathbf{Y}(t)-\mathbf{Y}(s)|^{2}}{4(\mu(t-s))}}\right]dsdt+\sqrt{\epsilon}d\mathbf{B}.$
(31)

Now, suppose that: $\mathbf{Y}(t)=\langle Vt,0\rangle$. This simplifies Eq. 31
to

$\frac{d\mathbf{Y}(t)}{dt}=V=\nu\mu\phi\frac{\pi}{2}\int_{-\infty}^{t}(\mu(t-s))^{-2}(Vt-
Vs)\exp\left[{-\frac{|Vt-Vs)|^{2}}{4(\mu(t-s))}}\right]ds.$ (32)

By making the variables given by $z=\mu(t-s)$ and $ds=-\frac{1}{\mu}dz$, we
see that Eq. 32 is considerably reduced to

$V=\nu\mu\phi\frac{\pi}{2}\int_{0}^{\infty}z^{-1}\frac{V}{\mu}\exp\left[-\left(\frac{V}{\mu}\right)^{2}\frac{z}{4}\right]dz.$
(33)

After some further simplification we arrive at the following expression,

$1=\frac{\pi}{2}\frac{\nu}{\mu}\phi\int_{0}^{\infty}\frac{1}{z}\exp\left[-\left(\frac{V}{\mu}\right)^{2}\frac{z}{4}\right]\frac{1}{\mu}dz.$
(34)

This integral equation is not pointwise convergent for finite $V$ and thus
indicates that when we consider the particle to be a point source with the
self-avoidant memory that we have defined, the particle does not swim.

###  A.3 Computing the Small Time Asymptotics of the Active Brownian MSD

Recall the MSD given for an active Brownian model with translational noise and
rotational diffusion given in Eq (6):

$\mathbb{E}[\mathbf{X}(t)^{2}]=4V^{2}\tau^{2}\left[2\left(e^{-\frac{t}{2\tau}}-1\right)+\frac{t}{\tau}\right]+2\epsilon
t$ (35)

Starting from the MSD in Eq. (7) for the active Brownian model, we rewrite the
exponential as an infinite series to arrive at

$\mathbb{E}[\mathbf{X}(t)^{2}]=4V^{2}\tau^{2}\left[2\left(\sum_{n=0}^{\infty}\left(\frac{t}{2\tau}\right)^{n}-1\right)+\frac{t}{\tau}\right]+2\epsilon
t.$ (36)

This is asymptotic to

$\mathbb{E}[\mathbf{X}(t)^{2}]\approx
4V^{2}\tau^{2}\left[2\left(\left(1-\frac{t}{2\tau}+\frac{t^{2}}{8\tau^{2}}\right)-1\right)+\frac{t}{\tau}\right]+2\epsilon
t$ (37)

as $t\to 0$ by just retaining a few leading order terms.

In the small time scale regime where $t\ll\tau$, we see that

$\mathbb{E}[\mathbf{X}(t)^{2}]\approx V^{2}t^{2}+2\epsilon t$ (38)

in agreement with Eq. (10). This expression is dominated by the diffusion-
generated term $2\epsilon t$ at the smallest time scales (where $t\ll t^{2}$)
and dominated by the directed motion term $V^{2}t^{2}$ when $t^{2}$ becomes
sufficiently larger than $t$. In the large time scale regime where $t\gg\tau
t$, we see that

$e^{-\frac{t}{2\tau}}\rightarrow 0$

and therefore

$\mathbb{E}[\mathbf{X}(t)^{2}]\approx(4V^{2}\tau^{2}+2\epsilon)t$

as given by Eq. (11). This expression contains the amount of enhanced
diffusion, $4V^{2}\tau^{2}$.

###  A.4 Computing MSD and OCF from Position Time Series Generated by the
Model

Absent a closed form expression for the mean square displacement of our model,
we compute the empirical MSD from the position time series of length $N+1$
given by $\mathbf{X}(t)$: $\\{\mathbf{X}(0),\ldots,\mathbf{X}(N)\\}$. To avoid
introducing any correlations into the increment averages, we use non-
overlapping increments. To achieve statistical accuracy, then, we simulate
many trajectories and average over multiple particles. We denote the integer
lag time as $\Delta L$ indicating the displacement traveled by the particle
between observations $j$ and $j+\Delta L$ and given by $\mathbf{X}(j+\Delta
L)-\mathbf{X}(j)$. The total number of non-overlapping increments of length
$\Delta L$ in a time series of length $N+1$ is
$k=\lfloor\left(\frac{N+1}{\Delta L}\right)\rfloor$. (In the event that the
index lag length $\Delta L$ does not evenly divide the number of increments
$N+1$, we remove the extra data from the beginning of the time.) Thus, the
empirical formula for the mean square displacement over the lag time $\Delta
L$ of a single particle is given by

$\Delta L^{2}=\frac{1}{k-1}\sum_{i=1}^{k}(\mathbf{X}(N-(i-1)\cdot\Delta
L)-\mathbf{X}(N-i\cdot\Delta L))^{2}.$ (39)

As shown in Fig. 5, successive increases in $\Delta L$ result in a sampling
process which coarse grains the position time series.

Using the same partitioning process described above and shown in Fig. 5 we can
compute the non-overlapping displacements and find the cosine between
consecutive pairs. The resulting time average of these computed cosines gives
the orientation correlation function, for which the formula is given in Eq.
(12).

Figure 5: Position time series of a sample trajectory with coarse grained lag
times of 10 increments. Panel (a) includes entire trajectory and panel (b) is
the inset identified with the dashed square.

## References

  * [1] W. F. Paxton, K. C. Kistler, C. C. Olmeda, A. Sen, S. K. St. Angelo, Y. Cao, T. E. Mallouk, P. E. Lammert, and V. H. Crespi.  Catalytic nanomotors: Autonomous movement of striped nanorods.  Journal of the American Chemical Society, 126(41):13424–13431, September 2004.
  * [2] H. Ke, S. Ye, R. L. Carroll, and K. Showalter.  Motion analysis of self-propelled Pt-silica particles in hydrogen peroxide solutions.  Journal of Physical Chemistry A, 114(17):5462–5467, April 2010\.
  * [3] I. Theurkauff, C. Cottin-Bizonne, J. Palacci, C. Ybert, and L. Bocquet.  Dynamic clustering in active colloidal suspensions with chemical signaling.  Physical Review Letters, 108(26):268303, June 2012.
  * [4] A.-Y. Jeea, Y.-K. Choa, S. Granicka, and T. Tlustya.  Catalytic enzymes are active matter.  PNAS, 115(46):E10812–E10821, November 2018.
  * [5] A. Sen, M. Ibele, Y. Hong, and D. Velegol.  Chemo and phototactic nano/microbots.  Faraday Discussions, 143(0):15–27, 2009.
  * [6] I. Buttinoni, G. Volpe, F. Kümmel, G. Volpe, and C. Bechinger.  Active brownian motion tunable by light.  Journal of Physics: Condensed Matter, 24(28):284129, Jun 2012.
  * [7] J. Palacci, S. Sacanna, A. P. Steinberg, D. J. Pine, and P. M. Chaikin.  Living crystals of light-activated colloidal surfers.  Science, 339(6122):936–940, February 2013.
  * [8] J. Palacci, S. Sacanna, S.-H. Kim, G.-R. Yi, D. J. Pine, and P. M. Chaikin.  Light-activated self-propelled colloids.  Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 372(2029), November 2014.
  * [9] S. Thutupalli, R. Seemann, and S. Herminghaus.  Swarming behavior of simple model squirmers.  New Journal of Physics, 13(7):073021, July 2011.
  * [10] P. G. Moerman, H. W. Moyses, E. B. van der Wee, D. G. Grier, A. van Blaaderen, W. K. Kegel, J. Groenewold, and J. Brujic.  Solute-mediated interactions between active droplets.  Phys. Rev. E, 96:032607, Sep 2017.
  * [11] C. Jin, C. Kruger, and C. C. Maass.  Chemotaxis and autochemotaxis of self-propelling droplet swimmers.  Proceecings of the National Academy of Sciences, 114(20):5089–5094, May 2017.
  * [12] A. Izzet, P. G. Moerman, P. Gross, J. Groenewold, A. D. Hollingsworth, J. Bibette, and J. Brujic.  Tunable persistent random walk in swimming droplets.  Physical Review X, 10(2):021035, May 2020.
  * [13] B. V. Hokmabad, M. Jalaal R. Dey, M. Almukambetova D. Mohanty, K. A. Baldwin, D. Lohse, and C. C. Maass.  Emergence of bimodal motility in active droplets.  arXiv:2005.12721v2, 2020.
  * [14] H. C. Berg and D. A. Brown.  Chemotaxis in escherichia coli analysed by three-dimensional tracking.  Nature, 239:500–504, october 1972.
  * [15] Patteson, A. E., Gopinath, A., Goulian, M., Arratia, and P. E.  Running and tumbling with e. coli in polymeric solutions.  Scientific Reports, 5(1):15761, 2015.
  * [16] J. Adler.  Chemotaxis in bacteria.  Science, 153(3737):708–716, August 1966.
  * [17] D. D. Thomas and A. P. Peterson.  Chemotactic auto-aggregation in the water mould achlya.  Microbiology, 136(5):847–853, May 1990.
  * [18] Z. Alirezaeizanjani, R. Großmann, V. Pfeifer, M. Hintsche, and C. Beta.  Chemotaxis strategies of bacteria with multiple run modes.  Science Advances, 6(22), May 2020.
  * [19] D. P. Häder.  Polarotaxis, gravitaxis and vertical phototaxis in the green flagellate, euglena gracilis.  Archives of Microbiology, 147(2):179–183, 1987.
  * [20] W.-L. Chen, H. Ko, H.-S. Chuang, H. H. Bau, and D. Raizen.  Caenorhabditis elegans exhibits positive gravitaxis.  bioRxiv, 2019.
  * [21] R. L. Stavis and R. Hirschberg.  Phototaxis in chlamydomonas reinhardtii.  Journal of Cell Biology, 59:367–377, 1973.
  * [22] E. Hildebrand and N. Dencher.  Two photosystems controlling behavioural responses of halobacterium halobium.  Nature, 257(5521):46–48, 1975.
  * [23] Y. Yang, V. Lam, M. Adomako, R. Simkovsky, A. Jakob, N. C. Rockwell, S. E. Cohen, A. Taton, J. Wang, J. C. Lagarias, A. Wilde, D. R. Nobles, J. J. Brand, and S. S. Golden.  Phototaxis in a wild isolate of the cyanobacterium synechococcus elongatus.  PNAS, 115(52):E12378–E12387, December 2018.
  * [24] J. Zhang, E. Luijten, B. A. Grzybowski, Bartosz A., , and S. Granick.  Active colloids with collective mobility status and research opportunities.  Chem. Soc. Rev, 46(18):5551–5569, 2017.
  * [25] S. J. Ebbens.  Active colloids: Progress and challenges towards realising autonomous applications.  Current Opinion in Colloid Interface Science, 21:14–23, 2016\.
  * [26] S. J. Ebbens and J. R. Howse.  In pursuit of propulsion at the nanoscale.  Soft Matter, 6(4):726–738, 2010.
  * [27] J. R. Howse, R. A. L. Jones, A. J. Ryan, T. Gough, R. Vafabakhsh, and R. Golestanian.  Self-motile colloidal particles: From directed propulsion to random walk.  Physical Review Letters, 99(4):048102–, 07 2007.
  * [28] D. P. Singh, U. Choudhury, P. Fischer, and A. G. Mark.  Non‐equilibrium assembly of light‐activated colloidal mixtures.  Advanced Materials, 29(32):1701328, June 2017.
  * [29] B. M. Schulz, S. Trimper, and M. Schulz.  Feedback-controlled diffusion: self-trapping to true self-avoiding walks.  Physics Letters A, 339:224–231, January 2005.
  * [30] H. Löwen.  Inertial effects of self-propelled particles: From active brownian to active langevin motion.  The Journal of Chemical Physics, 152(4):040901, January 2020.
  * [31] B. Liebchen and H. Löwen.  Synthetic chemotaxis and collective behavior in active matter synthetic chemotaxis and collective behavior in active matter.  Accounts of Chemical Research, 51(12):2982–2990, October 2018.
  * [32] Y. Tsori and P.-G. de Gennes.  Self-trapping of a single bacterium in its own chemoattractant.  Europhysics Letters (EPL), 66(4):599–602, 2004.
  * [33] A. Sengupta, S. v. Teeffelen, and H. Löwen.  Dynamics of a microorganism moving by chemotaxis in its own secretion.  Physical Review E, 80(3):031122–, 09 2009.
  * [34] R. Grima.  Strong-Coupling Dynamics of a Multicellular Chemotactic System.  Physical Review Letters, 95(4):128103, 09 2005.
  * [35] W. T. Kranz, A. Gelimson, K. Zhao, G. C. L. Wong, and R. Golestanian.  Effective Dynamics of Microorganisms That Interact with Their Own Trail.  Physical Review Letters, 117(6):038101, 08 2016.