Lorentz forces induce inhomogeneity and fluxes in active systems

TL;DR

This paper explores how space-dependent magnetic fields induce inhomogeneous density distributions and fluxes in active charged particles, revealing nonequilibrium steady states with circulating fluxes that differ from passive systems.

Contribution

It provides a theoretical framework and analytical results showing how Lorentz forces create inhomogeneity and fluxes in active matter, validated by simulations and with potential experimental implications.

Findings

Magnetic field gradients induce inhomogeneous density profiles.

Lorentz forces generate divergence-free bulk fluxes in active systems.

Reversible circulating fluxes observed in confined geometries.

Abstract

We consider the nonequilibrium dynamics of a charged active Brownian particle in the presence of a space dependent magnetic field. It has recently been shown that the Lorentz force induces a particle flux perpendicular to density gradients, thus preventing a diffusive description of the dynamics. Whereas a passive system will eventually relax to an equilibrium state, unaffected by the magnetic field, an active system subject to a spatially varying Lorentz force settles into a nonequilibrium steady state characterized by an inhomogeneous density and divergence-free bulk fluxes. A macroscopic flux of charged active particles is induced by the gradient of the magnetic field only and does not require additional symmetric breaking such as density or potential gradients. This stands in marked contrast to similar phenomena in condensed matter such as the classical Hall effect. In a confined geometry we observe circulating fluxes, which can be reversed by inverting the direction of the magnetic field. Our theoretical approach, based on coarse-graining of the Fokker-Planck equation, yields analytical results for the density, fluxes, and polarization in the steady state, all of which are validated by direct computer simulation. We demonstrate that passive tracer particles can be used to measure the essential effects of the Lorentz force on the active particle bath, and we discuss under which conditions the effects of the flux could be observed experimentally.

Equations

$${\mathbf{F}}_L=q\mathbf{v}\times \mathbf{B},$$

$$\frac{d\mathbf{r}(t)}{dt}=\mathbf{v}(t),$$

$$m\frac{d\mathbf{v}(t)}{dt}=-\gamma \Gamma (\mathbf{r}(t))\cdot \mathbf{v}(t)+f\mathbf{p}(t)+\sqrt{2\gamma T}\mathbf{\xi }(t),$$

$$\frac{d\mathbf{p}(t)}{dt}=\sqrt{2D_r}\mathbf{\eta }(t)\times \mathbf{p}(t),$$

$$\frac{\partial }{\partial t}Q(t)$$

$$\mathbf{J}(\mathbf{r},\mathbf{p},t)=-\frac{1}{\gamma }{\Gamma }^{-1}\cdot \left(T{\nabla }_{\mathbf{r}}-f\mathbf{p}\right)Q(t),$$

$${\mathbf{J}}_{rot}(\mathbf{r},\mathbf{p},t)=-D_r\mathcal{R}Q(t)$$

$${\Gamma }^{-1}(\mathbf{r})=1-\frac{\kappa (\mathbf{r})}{1+{\kappa }^2(\mathbf{r})}M+\frac{{\kappa }^2(\mathbf{r})}{1+{\kappa }^2(\mathbf{r})}{M}^2,$$

$$\kappa (\mathbf{r})=\frac{qB(\mathbf{r})}{\gamma }.$$

$$Q(\mathbf{r},\mathbf{p},t)=\rho +\mathbf{\sigma }\cdot \mathbf{p}+\mathbf{\tau }:\left(\mathbf{p}\mathbf{p}-1/3\right)+\Omega ,$$

$$\frac{\partial }{\partial t}\rho (\mathbf{r},t)=-\nabla \cdot \mathbf{j}(\mathbf{r},t),$$

$$\mathbf{j}(\mathbf{r},t)={\Gamma }^{-1}(\mathbf{r})\cdot \left[\frac{f}{\gamma }\bar{\mathbf{p}}(\mathbf{r},t)\rho (\mathbf{r},t)\right]-\frac{T}{\gamma }{\Gamma }^{-1}(\mathbf{r})\cdot \nabla \rho (\mathbf{r},t),$$

$$\bar{\mathbf{p}}(\mathbf{r},t)=-\frac{l}{3\rho (\mathbf{r},t)}\nabla \cdot \left[{\Gamma }^{-1}(\mathbf{r})\rho (\mathbf{r},t)\right],$$

$$l=\frac{f}{2D_r\gamma }$$

$$\rho (x)=N{\left[1+{\kappa }^2(x)\right]}^{\delta /2},$$

$$\delta =\frac{\Delta T}{T+\Delta T},$$

$${\bar{p}}_x(x)=\frac{(2-\delta )l}{3}\frac{\kappa (x){\kappa }^{\prime }(x)}{{\left[1+{\kappa }^2(x)\right]}^2},$$

$${\bar{p}}_y(x)=-\frac{l}{3}\frac{{\kappa }^{\prime }(x)}{{\left[1+{\kappa }^2(x)\right]}^2}\left[1-(1-\delta ){\kappa }^2(x)\right]$$

$$j_y(x)=\rho (x)\frac{f}{\gamma }{\bar{p}}_y(x).$$

$$⟨\Delta x(t{)}^2⟩=⟨\Delta y(t{)}^2⟩=⟨\Delta z(t{)}^2⟩={\left(\frac{ft}{\gamma }\right)}^2,$$

$$D_z=\frac{T+\Delta T}{\gamma }.$$

$$D_x(x)=D_y(x)=\frac{T+\Delta T}{\gamma }\frac{1}{1+\kappa (x{)}^2}.$$

$$D_x^{(eff)}=\frac{T+\Delta T}{\gamma }{\int }_0^{L}dx\frac{\rho (x)}{1+\kappa (x{)}^2},$$

$${\tau }_1=\frac{2{\gamma }^2{D}_x^{(eff)}}{f^2}.$$

$${\tau }_2=\frac{2{\gamma }^2{D}_z}{f^2}.$$

$${\tau }_3=\frac{L^2}{8D_x^{(eff)}}.$$

$$r_a=1+\frac{v_y^2{d}^4}{D_p^2{L}^2}.$$

$$v_y\approx \frac{4\pi }{3}\frac{v_0^2{\sigma }_b\lambda e|\nabla B|}{T}{R}^4.$$

$$\kappa =qB/\gamma \sim 0.1$$

$$\frac{\partial }{\partial t}P(t)=\left(L_{rev}+L_{irr}+L_r+L_a\right)P(t),$$

Full text

I Introduction

The Lorentz force modifies the trajectory of a charged particle without performing work on it. In general the Lorentz force $\boldsymbol{F}_{L}$ acts on a moving charge $q$ as

[TABLE]

where $\boldsymbol{v}$ denotes the velocity vector of the charge and $\boldsymbol{B}$ is the applied magnetic field. The Lorentz force gives rise to important phenomena in solid-state physics such as the classical Hall effect, in which electrical conductors develop a voltage difference transverse to the direction of the current if a magnetic field is applied perpendicular to the direction of the flow of charge carriers Hall (1879). Note that a macroscopic flux of charge carriers induced by Lorentz forces requires an explicitly broken symmetry that gives rise to a macroscopic velocity vector in addition to the symmetry breaking due to the magnetic field vector. Simply speaking the classical Hall effect leads to a deviation of an existing macroscopic flux. The main finding of this work is that in active systems driven by internal forces a macroscopic flux emerges from a flux-free system in the presence of a magnetic field gradient.

While the effect of the Lorentz force on the properties of materials has been thoroughly studied in the context of solid-state physics Altland and Simons (2010), much less is known about its influence on soft-matter systems, which are dominated by overdamped diffusive dynamics. Part of the reason is that Lorentz forces may appear to be too weak to affect slowly moving colloidal particles significantly. However, for Lorentz force to affect the motion of the particle, it is not its absolute value that is important but how its magnitude compares with the frictional force coming from the solvent. Since both the forces depend on velocity, the strength of Lorentz force relative to frictional force can be quantified and shown to be significant for highly charged mesoscopic particles in strong magnetic fields.

Before coming back to these points it is important to understand how the Lorentz force affects diffusive dynamics and how the diffusion equation can be extended to incorporate nonconservative velocity-dependent forces. It is known that as a consequence of the Lorentz force the corresponding Smoluchowski equation for the probability distribution acquires a tensorial diffusion coefficient Balakrishnan (2008), the components of which are determined by the applied magnetic field and the friction coefficient. Recently it has been shown that this tensor has an antisymmetric part, which gives rise to additional fluxes in the direction perpendicular to density gradients Chun et al. (2018); Vuijk et al. (2019). In a passive system, these fluxes affect the dynamics of a system as it relaxes from an initial nonequilibrium state to an equilibrium state. The equilibrium properties, however, are unaffected by the Lorentz force. A known consequence of the Lorentz force is a reduction of the diffusion coefficient in the plane perpendicular to the magnetic field, which can be understood by the deflection of the diffusing particle in the direction perpendicular to its trajectory according to the Eq. (1). This does not affect the equilibrium distribution of a system of passive particles.

A completely different situation arises when a driven system is subjected to an external magnetic field. A particularly interesting example of a driven system is active matter, which is ubiquitous in biology; examples include molecular motors inside cells Jülicher et al. (1997) and microscopic swimmers such as bacteria and algae Cates (2012). Besides biological active matter, there exists synthetic active matter, most notably Janus particles, which are ’two-faced’ particles with only one of them functionalized. When the functionalized half of the particle is activated, chemically or under illumination, the particle is propelled forward Jiang et al. (2010); Walther and Müller (2013); Buttinoni et al. (2012); Samin and van Roij (2015); Würger (2015). In contrast to externally driven systems, active matter has the hallmark feature that it is out of equilibrium without an externally broken symmetry. Much progress has been made in understanding the properties of active matter by using active Brownian particles (ABPs) as a model system, which captures the essential features of an internally driven system in a minimalistic fashion. ABPs violate time-reversal symmetry by consuming fuel to generate motion, often referred to as self-propulsion. Besides the application to biological systems, active matter serves as a paradigm to study the effect of broken time-reversal symmetry and nonequilibrium steady states in general Katz et al. (1983, 1984); Cates (2012); Fodor et al. (2016).

The absence of time-reversal symmetry can be exploited to obtain active fluxes by breaking the spatial symmetry of the system. For instance, sliding motion of active particles along an asymmetric wall Nikola et al. (2016); Katuri et al. (2018); Rodenburg et al. (2018); Reichhardt and Reichhardt (2017); Di Leonardo et al. (2009) can give rise to fluxes along a boundary. Similarly, in systems without walls or obstacles, fluxes can be induced by making the activity spatially inhomogeneous Stenhammar et al. (2016). Explicitly including alignment interactions, such as torque on an active particle, can also give rise to fluxes in a spatially inhomogeneous activity field as has been shown experimentally in Ref. Lozano et al. (2016). Whereas in these studies activity has been time independent, one can induce fluxes by making the swim speed space and time dependent, as has been shown in theoretical and numerical studies Geiseler et al. (2016, 2017a, 2017b); Merlitz et al. (2018); Zhu et al. (2018); Lozano and Bechinger (2019). In these systems the active particles ’surf’ on an activity wave resulting in active fluxes.

In this work, we show that the Lorentz force can induce macroscopic effects in a system of ABPs such as a magnetic-field dependent density distribution and flux-patterns. Note that, in contrast to externally driven system (for example, charge carriers in conductors under an applied voltage), ABPs do not have a preferential direction of motion. Thus, inhomogeneities are solely controlled by the spatial variation of the magnetic field. Whereas interaction between particles were necessary to induce bulk fluxes in a system with a space-dependent swim speed which is constant in time Stenhammar et al. (2016), this is not the case for the system studied here. Moreover, with Lorentz forces, the fluxes can be reversed by reversing the direction of the magnetic field.

For slowly moving particles, Lorentz forces can be extremely small. However, the effect of the Lorentz force on the motion of a particle is determined by ratio of the magnitudes of the Lorentz force and the frictional force. This ratio, which we refer to as the diffusive Hall constant, may become significant for highly charged particles in strong magnetic fields. We note that ABPs with magnetic interaction have been studied before; however, only by including a magnetic dipole moment without considering the Lorentz force. Vidal-Urquiza and Córdova-Figueroa (2017); Weber et al. (2011); Radtke and Schimansky-Geier (2012); Ao et al. (2014); Kline et al. (2005).

As an example of the effect of the Lorentz force on charged active particles we show the results of Brownian dynamics simulations where the ABPs are subjected to an inhomogeneous magnetic field in Fig. 1. The Gaussian shaped magnetic field is radially symmetric. The two signatures of the Lorentz force on ABPs are clearly visible: Circulating fluxes perpendicular to the gradient of the magnetic field and an inhomogeneous density distribution of the particles. The fluxes observed in ABPs under Lorentz force are reminiscent of the Corbino effect in electrical conductors Adams (1915). Whereas a radial flux of charges is necessary for the Corbino effect in conductors, this is not the case for ABPs which do not require explicit symmetry breaking in the system except the spatial inhomogeneity of the magnetic field.

On a microscopic level, fluxes in our system are related to the polarization of the active particles. Polar ordering can emerge in active systems without aligning interactions. This has been shown in previous studies, for instance, under gravity Enculescu and Stark (2011); Stark (2016) and inhomogeneous activity Stenhammar et al. (2016); Sharma and Brader (2017). As in these studies, the origin of polarization in our system is kinetic in nature. However, in contrast to previous studies, the particles in our system orient not only along the gradient of the applied field but also perpendicular to it. Related to this polar order are active fluxes in the same direction.

A common starting point of the theoretical description of ABPs is the Langevin equation. Most models of ABP ignore inertial effects and are based on overdamped equations of motion Sharma and Brader (2016, 2017); Cates and Tailleur (2013); Marini Bettolo Marconi and Maggi (2015); Vuijk et al. (2018); Merlitz et al. (2018). Due to the nonconservative nature of the Lorentz force, the overdamped equation of motion cannot be obtained by simply setting the mass of the particles to zero Lau and Lubensky (2007); Chun et al. (2018). Though there exists a nontrivial limiting procedure Hottovy et al. (2012, 2015); Volpe and Wehr (2016) that yields the small-mass limit of the (velocity) Langevin equation, the resulting overdamped equation is not suitable for determining velocity dependent variables such as flux and entropy Celani et al. (2012); Marino et al. (2016); Birrell (2018); Chun et al. (2018); Vuijk et al. (2019) despite the fact that it captures the position statistics accurately. In this work, we have avoided these problems by performing simulations using the (velocity) Langevin equation with a finite but small mass. Our theoretical approach is based on coarse-graining of the Fokker-Planck equation using a gradient expansion Schnitzer (1993); Cates and Tailleur (2013) , which yields an equation for the flux and density distribution as a function of time and the position variable alone. This equation also yields the space-dependent polarization in the system by identifying the corresponding contribution from the orientational degrees of freedom. We further show that one can improve the accuracy of the results by combining the coarse-graining approach with linear-response theory. The results of this method match those of Brownian dynamics simulations.

Our calculations lead to a straightforward estimation of the magnitude of this new kind of Hall effect in ABPs by comparing the amount of directed motion with the undirected motion of the active particles. We present a possible experimental realization of our findings: The motion of a large tracer particle, insensitive to Lorentz force, would become anisotropic due to the presence of the (active) flux patterns. Using reasonable estimates for the system parameters, our calculations suggest that this anisotropy could be measured experimentally.

II Active particles with Lorentz forces

We consider a single spherical, self-propelled, charged Brownian particle in a magnetic field $\boldsymbol{B}(\boldsymbol{r})$ with a spatially varying magnitude but constant direction. The active motion of the particle is modeled as a force in the direction of the orientation of the particle specified by an embedded unit vector $\boldsymbol{p}$ undergoing rotational diffusion. The state of the particle is determined by the position coordinate $\boldsymbol{r}(t)$, the velocity $\boldsymbol{v}(t)$ and the orientation $\boldsymbol{p}(t)$.

The dynamics are described by the following (Stratonovich) stochastic differential equations:

[TABLE]

[TABLE]

and

[TABLE]

where $m$ is the particle mass, $\gamma$ is the friction coefficient, $f$ is the self-propulsion force, $T$ is the temperature in units such that Boltzmann constant is equal to one and $D_{r}$ is the rotational-diffusion constant. The (dimensionless) matrix $\Gamma(\boldsymbol{r})=\boldsymbol{1}+\frac{qB(\boldsymbol{r})}{\gamma}M$, where $q$ is the charge of the particle and $M$ is a matrix such that $B(\boldsymbol{r})M\cdot\boldsymbol{v}=\boldsymbol{B}(\boldsymbol{r})\times\boldsymbol{v}$, and $B(\boldsymbol{r})\equiv|\boldsymbol{B}(\boldsymbol{r})|$. The stochastic vectors $\boldsymbol{\xi}$ and $\boldsymbol{\eta}$ are Gaussian distributed with zero mean and autocorrelation $\left<\boldsymbol{\xi}(t)\boldsymbol{\xi}^{T}(t^{\prime})\right>=\left<\boldsymbol{\eta}(t)\boldsymbol{\eta}^{T}(t^{\prime})\right>=\delta(t-t^{\prime})\boldsymbol{1}$. Note that there is no direct coupling between the orientation and the magnetic field.

We start with Eq. (3) with the inertia term because in presence of a Lorentz force, the overdamped equation of motion cannot be used to determine observables depending on the velocity process (such as flux), even though the position-statistics can be accurately described Hottovy et al. (2015); Chun et al. (2018); Vuijk et al. (2019). However, the Fokker-Planck equation (FPE) for the overdamped motion corresponding to Eqs. (2) and (3) can be obtained by a different method based on an expansion in powers of mass. This expansion removes the velocity degrees of freedom (see App. A) and yields the following FPE for the probability density $Q(t)\equiv Q(\boldsymbol{r},\boldsymbol{p},t)$:

[TABLE]

where

[TABLE]

is the translational flux,

[TABLE]

is the rotational flux, $\boldsymbol{\mathcal{R}}=\boldsymbol{p}\times\nabla_{\boldsymbol{p}}$ is the rotation operator, and

[TABLE]

with

[TABLE]

This parameter, which we refer to as the diffusive Hall constant, measures the strength of the Lorentz force relative to the frictional force in the system. When this parameter is comparable to unity, the Lorentz force affects the motion of the particle significantly. As we show later in Sec. V, the Lorentz force may become comparable to the frictional force for highly charged particles in strong magnetic fields.

Note that the matrix $\Gamma^{-1}$ in Eq. (8) is not symmetric, and therefore, $T\Gamma^{-1}/\gamma$ cannot be interpreted as a diffusion tensor Vuijk et al. (2019). Due to the antisymmetric terms in $\Gamma^{-1}$, there are additional fluxes in the system that are perpendicular to the density gradient.

Equation (5) can be coarse grained further by a gradient expansion to remove the orientational degrees of freedom. We consider the ansatz

[TABLE]

where the scalar $\rho$, the vector $\boldsymbol{\sigma}$ and the matrix $\boldsymbol{\tau}$ are function of $\boldsymbol{r}$ and $t$, and $\Omega$ contains all higher-order terms. Equations for the time evolution of the probability density $\rho$, the total orientation $\boldsymbol{\sigma}/3$ and higher-order tensors can be obtained by inserting Eq. (10) into Eq. (5) and integrating out the orientation degrees of freedom (see App. B). In the limit where gradients are small compared to the persistence length of the ABP Schnitzer (1993); Cates and Tailleur (2013), one obtains a FPE for the probability density $\rho(\boldsymbol{r},t)$ as a function of time and position alone:

[TABLE]

where

[TABLE]

and

[TABLE]

where

[TABLE]

is the persistence length of the ABP, and $\boldsymbol{\bar{p}}(\boldsymbol{r},t)$ is the polarization, which is defined as the average orientation per particle $\left<\boldsymbol{p}\delta\left(\boldsymbol{r}^{\prime}-\boldsymbol{r}\right)\right>/\rho(\boldsymbol{r},t)$ where the angle brackets indicate an average with respect to the probability density of $\boldsymbol{r}^{\prime}$, $\boldsymbol{p}$ and $\boldsymbol{v}$ at time $t$.

III Nonequilibrium Steady State

Here we focus on the steady-state bulk behavior of the system and therefore use periodic boundary conditions in all directions. Furthermore, the analysis is restricted to a magnetic field with a constant orientation, in this case the $z$ direction, and a magnitude that depends on a single coordinate, say $x$. Since the motion along the direction of the magnetic field is not affected by the Lorentz force, there is neither flux nor density variation along the $z$ direction in steady state. For the chosen magnetic field, the system is translationally invariant in the $y$ direction; therefore, the flux and the polarization in the $y$ direction are independent of the $y$ coordinate. This also implies that in steady state there can be no flux in the $x$ direction because of the periodicity and equivalence of the positive and negative $x$ directions. Moreover, there can be no local $x$ fluxes forming closed loops via the $y$ direction because of the translation invariance in the $y$ direction. The absence of fluxes in the $x$ direction has been corroborated by simulations.

By setting the $x$ component of the flux [Eq. (12)] to zero, we obtain the following analytical expression for density:

[TABLE]

where $N^{-1}=\int_{0}^{L}dx\left[1+\kappa^{2}(x)\right]^{\delta/2}$ and

[TABLE]

where $\Delta T=f^{2}/6D_{r}\gamma$. The polarization can be calculated by inserting Eq. (15) in Eq. (13), which yields

[TABLE]

for the $x$ component, and

[TABLE]

for the $y$ component, where $\kappa^{\prime}(x)=\frac{d\kappa(x)}{dx}$.

The flux has two contributions, one coming from the polarization and the other from the gradient in the density [see Eq. (12)]. In steady state, the two contributions to the flux in the $x$ direction are equal in magnitude, but in opposite direction, resulting in zero flux in this direction. In the $y$ direction there can be no density gradients because of the translation symmetry in this direction. This means that the contribution to the flux in the $y$ direction coming from a density gradient along this direction is zero; however, the contribution from the polarization is not zero, resulting in a flux in the $y$ direction. This flux can be calculated by using the solution for density [Eq. (15)] as input to the $y$ component of Eq. (12), yielding

[TABLE]

Note that the orientation and flux in the $y$ direction are perpendicular to the magnetic field gradient.

Before we test the predictions of the theory, a few remarks are in order. First, $\Delta T$ appearing in Eqs. (15)–(18) has units of energy (note that the Boltzmann constant is unity, so $T$ has units of energy as well) and has previously been interpreted as an effective increase in the temperature due to activity Solon et al. (2015); Farage et al. (2015). However, we do not interpret $\Delta T$ as increase in temperature because this system cannot be mapped to an equilibrium system at a higher temperature. Second, the expressions for the density, flux and polarization can be improved significantly by combining the gradient expansion with linear-response calculations of the polarization. Such a calculation results in a differential equation for the polarization, which can be solved numerically. This solution can be used to replace the linear-order contribution of the activity to the polarization calculated using the gradient expansion and using this polarization as input to the flux in Eq. (15). Details of the linear-response theory, its application to the polarization and the combined solution are shown, respectively, in Appendices C, D and E. Finally, it is interesting to note that $j_{y}(x)/\rho(x)$ is the macroscopic velocity of the ABP in the presence of the gradient field, while $f/\gamma$ gives the intrinsic microscopic velocity of the ABP. While the intrinsic velocity leads to random motion on large time scales, the field-induced motion persists (in our case in the $y$ direction) on any time scale. Thus the polarization [Eq. (18)] gives the fraction of the swim force that is converted into macroscopic directed motion.

Figure 2 shows the predictions of the theory for the density, flux and polarization for $\kappa(x)=4\sin(2\pi x/10)$ and $\gamma=1$ together with the results from Brownian dynamics (BD) simulations. Consistent with the prediction of equation (15), the density is high where the magnitude of the magnetic field is large (see Fig. 2 a). Qualitatively, the accumulation of particles in the high magnetic field region can be understood as follows. A magnetic field hinders the diffusive motion of Brownian particles because there is a force perpendicular to the velocity which results in circular motion. For passive Brownian particles this means that on time scales longer than the velocity autocorrelation time the effect of a uniform magnetic field is a decrease of the diffusion constant in the plane perpendicular to the magnetic field Balakrishnan (2008). For a spatially varying magnetic field, the diffusion coefficient of particles is small in the region where magnetic field is large. Whereas, for passive particles, the steady-state density distribution is independent of the diffusion coefficient, active particles accumulate in the region of low swim speed, which correspond to regions of high friction, or in this case high magnetic field Cates (2012); Sharma and Brader (2017).

Passive Brownian particles in a space dependent magnetic field have a uniform steady-state density (a Boltzmann distribution with a uniform potential). A system of ABPs with a homogeneous activity also has a uniform steady-state density. However, the combination of a space-dependent magnetic field and uniform activity results in an inhomogeneous steady-state density. An inhomogeneous density without a space dependent potential means that the density does not follow a Boltzmann distribution, a hallmark of nonequilibrium systems.

There is polarization in the plane perpendicular to the magnetic field; see Eqs. (17) (18) and Figs. 2 c and d. In steady state, the polarization and the gradient of the density in the $x$ direction are such that the fluxes resulting from them cancel each other. The polarization in the $y$ direction, however, cannot be canceled by a density gradient because when the magnetic field varies only along the $x$ direction, the system is translation invariant in the $y$ direction and therefore, the density, flux and polarization cannot depend on the $y$ coordinate. Because there is polarization in the $y$ direction, but no density gradient to counter the flux from the polarization, there is in steady state a flux in the $y$ direction [see Eq. (19) and Fig. 2 b]. Due to the periodicity of the field, there is flux in both the positive as the negative $y$ direction resulting in zero center-of-mass movement.

Note that there is no direct coupling between the magnetic field and the polarization [see Eqs. (3) and (4)], and the inhomogeneous density along the $x$ direction is not a result of the polarization. Both the polarization and the inhomogeneous density are kinetic in nature, and result from a filtering mechanism. The magnetic field hinders the motion of the ABP and a particle with an orientation parallel to the gradient of the magnetic field will get stuck. A particle with an orientation antiparallel to the gradient of the magnetic field can more easily move towards regions of lower magnetic field. These effects result in a force imbalance between regions of small and large magnetic fields. In steady state the force imbalance is canceled by gradient in the chemical potential due to the increase of density in the high magnetic region.

The presence of steady-state bulk fluxes is a key feature of the nonequilibrium steady state induced by Lorentz forces in a system of active particles. These fluxes are unique in the sense that the underlying mechanism is independent of inter-particle interactions Stenhammar et al. (2016), aligning torques Lozano et al. (2016) or sliding motion along asymmetric walls Nikola et al. (2016); Katuri et al. (2018); Rodenburg et al. (2018); Reichhardt and Reichhardt (2017); Di Leonardo et al. (2009). Moreover, the direction of fluxes is perpendicular to the gradient of the applied magnetic field. This particular property is what gives rise to the circulating steady-state fluxes shown in Fig. 1 in which the magnetic field is radially symmetric. The direction of the fluxes can be reversed by simply reversing the direction of the magnetic field.

In case of a finite system with walls, the bulk properties of the system are similar to that of an infinite system; however, there is also flux parallel to the wall resulting in closed-loop fluxes (see Fig. 3). A flux along walls has recently been investigated in the case of chiral active particle Caprini and Marconi (2019). The flux along the wall induced by Lorentz force can possibly be explained by the same mechanism.

IV Active diffusion under Lorentz force

Here we investigate the mean squared displacement (MSD) of an ABP in the three directions. The MSD in the $y$ direction is particularly interesting: for a magnetic field varying in the $x$ direction, the steady-state fluxes in the $y$ direction constitute layers of alternating flows (see inset to Fig. 4 and Fig. 2). How does the emergence of these layers of alternating flows influence the MSD in the $y$ direction?

Free ABPs have a swim speed of $f/\gamma$ Merlitz et al. (2017); therefore, on short time scales the motion of ABPs is ballistic and is the same in the three directions:

[TABLE]

where $\Delta x(t)$ is $x(t)-x(0)$, and similarly for $\Delta y(t)$ and $\Delta z(t)$. The diffusion parallel to the magnetic field orientation is unaffected by the field. This implies that the diffusion coefficient in the $z$ direction is that of a freely diffusing ABP, which is

[TABLE]

Since a magnetic field hinders the motion in the plane perpendicular to its direction, the diffusion in the $x$ and $y$ directions is reduced. For an inhomogeneous magnetic field, the spatially varying diffusion coefficient in the $x$ and $y$ directions can be obtained from the symmetric part of the tensor in Eq. (8) as

[TABLE]

On a length scale larger than the period of the magnetic field, the space dependence of the diffusion constant becomes irrelevant, and one can approximate the effective diffusion constant for MSD in the $x$ direction by averaging over space:

[TABLE]

where $L$ is the period of the magnetic field.

There are three time scales in the system. Under the assumption that the length scale of ballistic motion is much smaller than the period of the field, the shortest time scale, $\tau_{1}$, is the transition from ballistic to diffusive motion in the $x$ direction. This time scale can be obtained from $(f\tau_{1}/\gamma)^{2}=2D_{x}^{(eff)}\tau_{1}$, and is

[TABLE]

The second time scale, $\tau_{2}$, is the transition from ballistic to diffusive motion in the $z$ direction, which can be calculated from $(f\tau_{2}/\gamma)^{2}=2D_{z}\tau_{2}$ as

[TABLE]

The largest time scale, $\tau_{3}$, is the transition to diffusive motion in the $y$ direction. The diffusion in the $y$ direction is increased due to the flux lanes in the system. Because the net flux in the system is zero, the motion in the $y$ direction becomes diffusive on the same time scale as the time it takes to diffuse from on lane to another. This time scale can be calculated from $2D_{x}^{(eff)}\tau_{3}=(L/2)^{2}$ as

[TABLE]

The MSD in the $y$ direction becomes subballistic on the same time scale as the $x$ direction, namely $\tau_{1}$; however the diffusive regime starts later at $\tau_{3}$. This means that there is a super diffusive regime that persists over several decades in time. Figure  4 shows that when $t>\tau_{3}$ the MSD in the $y$ direction has the same effective diffusion coefficient as the $z$ direction. The magnetic field hinders diffusion in both $x$ and $y$ direction. Whereas the diffusion in $x$ direction is expectedly smaller than in the $z$ direction, the flux lanes increase the effective diffusion along $y$ and therefore the MSD in the $y$ direction ’catches up’ with $z$ direction.

V Passive diffusion of tracer particle

We showed in Sec III that the polarization [Eq. (18)] is the fraction of the swim force that is converted into macroscopic motion. Since polarization scales as $l\kappa^{\prime}\ll 1$, it is small in the limit where gradients are small compared to the persistence length of the ABP. To overcome this limitation, we consider a system in which the effect of the anisotropic motion of the ABPs on a tracer particle can accumulate.

Consider again the system of ABPs in steady state in a periodically varying magnetic field. As discussed above the ABPs will form a periodic array of lanes carrying an alternating macroscopic flux in $y$ direction, see Fig. 5. We now immerse an uncharged, passive tracer particle in the system and study its MSD. On can expect that the steady-state fluxes of the ABPs affect the diffusion in the $y$ direction of the tracer particle. Similar problems have been studied in different contexts, such as the diffusion of a passive particle in presence of random velocity fields Redner (1989, 1990); Bouchaud et al. (1990) and the motion of polymer chains in random layered flows Sommer and Blumen (1996).

Due to the (alternating) velocity field, one can expect that the MSD of the tracer particle along $y$ direction is enhanced. This implies that the motion of the particle is anisotropic. The enhancement of the diffusion along the $y$ direction can be estimated as follows: The particle diffuses freely with a diffusion rate $D_{p}$ between the velocity lanes. Assuming that inside the lane, the particle moves with the local velocity, the distance that the particle travels along the lane is $l_{f}=v_{y}d^{2}/2D_{p}$, where $v_{y}$ is the local velocity and $d$ is the width of the lane. The particle experiences these flow lanes every $\tau_{y}=L^{2}/8D_{p}$ seconds, where $L/2$ is the distance between two neighboring lanes. The increase in the diffusion rate along the $y$ direction can be estimated as $l_{y}^{2}/2\tau_{y}\sim v_{y}^{2}d^{4}/D_{p}L^{2}$. The anisotropy parameter, defined as the ratio of the MSDs in the $x$ and $y$ directions, is

[TABLE]

An experimental measurement of the anisotropy of the diffusion of a tracer particle, as described above, would provide an, albeit indirect, experimental verification of the active fluxes in the system. Since Lorentz force is relativistic in nature and ABPs are colloidal particles, one expects that $v_{y}$ is extremely small. As we show below, this is indeed the case for the known experimentally realizable ABPs. However, since the anisotropy parameter depends inversely on the diffusion coefficient of the tracer particle, it might still be observed experimentally for sufficiently large tracer particles on long time scales. The smaller the diffusion constant of the tracer particle the longer it will stay in one lane and thus accumulating the effect by drifting along with the velocity field. In the extreme case of a vanishing diffusion constant, $D_{p}$=0, the only observable displacement on a long time scale is the drift and thus $r_{a}\gg 1$.

The anisotropy parameter can be estimated as follows: The lane velocity can be obtained from Eq. (19) as $v_{y}=j_{y}/\rho=v_{0}\bar{p}_{y}$, where $v_{0}=f/\gamma$ is the self-propulsion speed and $\bar{p}_{y}$ is the $y$ component of the polarization vector inside a lane. It can be obtained from Eq. (18) as $\bar{p}_{y}\approx lq|\nabla B|/3\gamma$, where $l=v_{0}/2D_{r}$ is the persistence length and $|\nabla B|$ is the (maximum) gradient of the magnetic field. With this the lane velocity can be estimated as $v_{y}\sim v_{0}^{2}q|\nabla B|/(2D_{r}\gamma)$. To estimate the charge, we consider that the active particle is a Janus-shaped particle with half of its surface grafted with polymer chains, and the other half serves as the propulsion engine Jiang et al. (2010); Walther and Müller (2013); Buttinoni et al. (2012); Samin and van Roij (2015); Würger (2015). The total charge on such an ABP is $q=2\pi R^{2}\sigma_{b}\lambda e$, where $R$ is the radius of the ABP, $\sigma_{b}$ is the grafting density, $\lambda$ is the number of charges per polymer chain, and $e$ is the elementary charge. Using this together with the Stokes relations for the rotational diffusion coefficient $D_{r}$ and $\gamma$, one obtains

[TABLE]

In order to estimate the magnitude of lane velocity, we consider $T=298/k_{b}\,kgm^{2}s^{-2}$, where $k_{b}$ is the Boltzmann constant, a dense grafting density of $\sigma_{b}=1nm^{-2}$ and highly charged chains with $\lambda=100$. For the magnetic field we assume a period of $L\sim 1\,cm$. For the magnetic field gradient we use $10\,T/cm$, which can be obtained in magnetic matrices used for filtering fine and weakly magnetic particles Ge et al. (2017). For the ABPs we take the values $v_{0}\approx 3\,\mu ms^{-1}$ and $R\approx 2\,\mu m$ from Ref. Lozano et al. (2016). With these, we obtain the lane velocity $v_{y}\approx 2.5\times 10^{-3}\,\mu ms^{-1}$. The width of a lane can be approximated as $d=L/(2\pi)\sim 0.1\,cm$. Note that for the chosen parameters the diffusive Hall-constant is of the order,

[TABLE]

implying that Lorentz force becomes comparable to frictional force on the particle. Finally, assuming that the diffusion coefficient of the tracer particle is $D_{p}\sim 10^{-13}\,m^{2}/s$, a reasonable estimate for large tracer particles ($\sim\mu m$), we obtain from Eq. (27) $r_{a}\approx 2$. Of course long-time scales are necessary to detect the effect. However, the distance that a tracer particle drifts along a velocity lane in a day time turns out to be approximately $0.2\,mm$, which can be detected experimentally. These estimates suggest that Lorentz forces may be observable experimentally for charged, active particles. We note that it is not the magnitude of the magnetic field that is important for this Hall-like effect, but it is only the gradient which enters the calculation.

A few remarks are necessary with respect to a possible experimental setup. First, the total charge of the poly-electrolyte brush is compensated by counter-ions arising from the dissociation of the charged groups of polymers in water. In a salt-free solution the counter-ions will be strongly localized inside the brush and thus effectively cancel the effect of the Lorentz-force Pincus (1991). However, a superior aspect of the external magnetic field is that it is not screened by electric charges. Thus, at a high salt level in the solution the counter-ions are quasi-free Borisov et al. (1994); He et al. (2010). Furthermore, the exchange dynamics of them will be much faster than the Lorentz-induced drift. Therefore, brush-decorated active colloidal particles in a salted solution could be a promising experimental system.

VI Discussion and Conclusions

Lorentz force has the unique property that it depends on the velocity of the particle and is always perpendicular to it. Although this force generates particle currents, they are purely rotational and do no work on the system. As a consequence, the equilibrium properties of a Brownian system, for instance the steady-state density distribution, are independent of the applied magnetic field. The dynamics, however, are affected by Lorentz force: the Smoluchowski equation picks up a tensorial (diffusion) coefficient, which reflects the anisotropy of the particle’s motion Balakrishnan (2008). The diffusion rate perpendicular to the direction of the magnetic field decreases with increasing field whereas the rate along the field remains unaffected. In addition to this effect, Lorentz force gives rise to nondiffusive fluxes which are perpendicular to the density gradients Vuijk et al. (2019); Chun et al. (2018).

The effects caused by the Lorentz force, however, occur only in nonequilibrium and cease to exist when the distribution of charged particles reaches equilibrium. A system of ABPs, in contrast, is continuously driven out of equilibrium. Lorentz force couples to the nonequilibrium dynamics of an ABP via its self-propulsion. In this work, we showed that when subjected to an inhomogeneous magnetic field, a system of uniformly active ABPs settles into an unusual nonequilibrium steady state. This state is characterized by (a) inhomogeneously distributed density with the same spatial symmetry as that of the magnetic field and (b) macroscopic fluxes perpendicular to the gradient of the magnetic field.

The theoretical approach is based on coarse graining of the Smoluchowski equation and yields analytical results for the density distribution, active fluxes, and polarization in the steady state. We also showed that when combined with the linear-response theory, the theoretical predictions match the results from BD simulations.

That the density distribution in our system is inhomogeneous can be understood as follows. The swim force increases the ABP’s velocity in the direction of its orientation. The Lorentz force rotates the velocity vector. The effect of this rotation is a decrease in the speed of the particle in the direction of its orientation, which is similar as decreasing the effect of the swim force. From this it follows that a spatially varying magnetic field has a similar effect as ABP as subjecting them to an inhomgeneous activity field. However, the existence of nondiffusive fluxes perpendicular to the density gradient in the steady-state is a unique feature of this system and is a consequence of the antisymmetric coupling of the $x$ and $y$ directions due to the Lorentz force. In case of a finite system with walls, the bulk properties of the system are similar to that of an infinite system; however, there is also flux parallel to the wall resulting in closed-loop fluxes

The MSD of an ABP in the presence of an inhomogeneous magnetic field has nontrivial behavior related to the existence of flux lanes in the system. Whereas the motion along the direction of the field is unaffected by the field, it is hindered in the plane perpendicular to the field. The diffusion coefficient in the $z$ direction is, as expected, the same as that of a freely diffusing ABP. The diffusion in the $x$ direction is expectedly smaller than in the $z$ direction; however, the same does not hold for the $y$ direction. The diffusion in the $y$ direction is super diffusive at intermediate times. At long times, the motion along the $y$ direction becomes diffusive with the same rate as along the $z$ coordinate. The emergence of flux lanes in the $y$ direction gives rise to this nontrivial behavior of the MSD in the $y$ direction.

Our study naturally raises the question whether the effects induced by Lorentz force in a system of ABPs can be observed experimentally. A Lorentz force on a colloidal particle is extremely small due to the particle’s low velocity. However, for such a particle to be affected by a Lorentz force, it is not the absolute value of the Lorentz force that is important but how it compares with the frictional forces. Since both frictional and Lorentz force depend on the velocity, for Lorentz force to have any significant effect on the motion of the particle the diffusive Hall constant $\kappa=qB/\gamma$ must be of the order of unity. In Sec. V, we considered the possibility of experimental verification of this study where we have proposed to use the high total charge of a poly-electrolyte brush grafted to an ABP to obtain a value of $\kappa\simeq 0.1$. Direct experimental verification of the flux lanes might still be difficult. A possibility is to measure the effect of the flux lanes on the diffusive behavior of a passive tracer particle on which the effect of the flux lanes can accumulate. Using values for propulsion speed and particle size from experiments Lozano et al. (2016) and large but reasonable estimates for other parameters in the system, we estimated a drift of the tracer particle along the flux lines of the active system of the order of millimetre per day. This indicates that magnetic-field induced effect on a charged active system could be measured experimentally.

Finally, note some general aspects of the application of magnetic fields to soft-matter systems. First of all, magnetic fields are not screened by static charges such as electric field are. Strong electric field strength in particular combined with thin-film geometries are prone to short-cuts, irreversibly damaging the system. By contrast the application even of high magnetic fields is possible and is used in, for example, high-resolution NMR-experiments to medicine. In this work we have shown that Lorentz forces may be detectable in active soft-matter systems, which opens new possibilities to manipulate such systems. Interestingly, macroscopic effects are induced by the symmetry breaking of the magnetic field only and do not require a macroscopic flux of the particles which establishes a novel quality of the Hall-like effect, which in this case depends on the gradient of the field only.

Appendix A small-mass limit

In order to obtain the small-mass limit of the FPE corresponding to Eqs. (2)(3) and (4), one can use the method described in Chap. 10 of Ref. Risken (1989). This was done for a single passive particle in a homogeneous magnetic field in Ref. Chun et al. (2018). We start with the FPE corresponding to Eqs. (2), (3) and (4) for the time evolution of the probability density $P(t)\equiv P(\boldsymbol{r},\boldsymbol{v},\boldsymbol{p},t)$ is Gardiner (2009)

[TABLE]

where the time-evolution operator has been split up in a reversible part

[TABLE]

an irreversible part without the part coming from Eq. (4)

[TABLE]

the part representing the rotation of the orientation vector

[TABLE]

with the rotation operator $\boldsymbol{\mathcal{R}}=\boldsymbol{p}\times\nabla_{\boldsymbol{p}}$ Morse et al. (1953), and the active, nonequilibrium operator coming from the self-propulsion

[TABLE]

Equation (30) is equivalent to

[TABLE]

where

[TABLE]

and

[TABLE]

where $L$ can be any of the operators in Eq. (30), and

[TABLE]

is the solution to $L_{irr}R(\boldsymbol{v})=0$, normalized such that the integral over $\boldsymbol{v}$ is one. The transformed operators are

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

where $\boldsymbol{b}=\sqrt{T/m}\nabla_{\boldsymbol{v}}+\frac{1}{2}\sqrt{m/T}\boldsymbol{v}$ and $\boldsymbol{b}^{\dagger}=-\sqrt{T/m}\nabla_{\boldsymbol{v}}+\frac{1}{2}\sqrt{m/T}\boldsymbol{v}$.

The eigenfunctions of the operator $b^{\dagger}_{\alpha}b_{\alpha}$, where $\alpha$ is $x$,$y$ or $z$, are

[TABLE]

and

[TABLE]

where $H_{n}$ are Hermite polynomials. The operators $b_{\alpha}^{\dagger}$ and $b_{\alpha}$ are the raising and lowering operators of the eigenfunctions: $b_{\alpha}^{\dagger}\psi_{n}(v_{\alpha})=\sqrt{n+1}\psi_{n+1}(v_{\alpha})$ and $b_{\alpha}\psi_{n}(v_{\alpha})=\sqrt{n}\psi_{n-1}(v_{\alpha})$. Because the eigenfunctions are orthonormal,

[TABLE]

they can be used to expand $\bar{P}(t)$:

[TABLE]

where $c_{n_{x},n_{y},n_{z}}\equiv c_{n_{x},n_{y},n_{z}}(\boldsymbol{r},\boldsymbol{p},t)$.

The probability density for the position and orientation, $Q(t)\equiv Q(\boldsymbol{r},\boldsymbol{p},t)$, is given by the first expansion coefficient:

[TABLE]

The translational flux,

[TABLE]

can be evaluated by using Eqs. (36), (46) and $v_{\alpha}\psi_{0}(v_{\alpha})=\sqrt{\frac{T}{m}}\psi_{1}(v_{\alpha})$:

[TABLE]

Using this, together with the definition of the flux in orientation space

[TABLE]

one obtains the following equation for the time evolution of the probability density:

[TABLE]

Next, we need an expression for the $\boldsymbol{c}_{1}$ in terms of $Q(t)$. We are interested in the case where the magnetic field is oriented perpendicular to the gradient of the field, and take this to be the $z$ direction, so $\boldsymbol{B}(\boldsymbol{r})=B(\boldsymbol{r})\hat{\boldsymbol{z}}$. Equation (35) together with the orthonormality of the eigenfunctions yields an hierarchy of equations for the functions $c_{n_{x},n_{y},n_{z}}$ called a Brinkman hierarchy Brinkman (1956):

[TABLE]

where $\boldsymbol{D}=\sqrt{\frac{T}{m}}\nabla_{\boldsymbol{r}}$ and $\hat{\boldsymbol{D}}=\sqrt{\frac{T}{m}}\nabla_{\boldsymbol{r}}-\frac{1}{\sqrt{Tm}}f\boldsymbol{p}$.

The equation governing the time evolution of the first expansion coefficient is

[TABLE]

where $\boldsymbol{c}_{1}=(c_{1,0,0},c_{0,1,0},c_{0,0,1})^{T}$. This is of course the same as Eq. (53). The equation for $\boldsymbol{c}_{1}(t)$ is again a differential equation in $t$, which can be solved after a Laplace transformation:

[TABLE]

where the tilde indicates a Laplace transformation from the variable $t$ to $s$, $D_{1,2}$ is a matrix of which the entries are linearly proportional to $D_{x}$, $D_{y}$ or $D_{z}$, $\boldsymbol{c}_{2}(t)$ is a vector of which the elements are all $c_{n_{x},n_{y},n_{z}}$ such that $n_{x}+n_{y}+n_{z}=2$, and

[TABLE]

Note that $\widetilde{\Gamma}_{m}^{(1)}(s)=\Gamma(\boldsymbol{r})+\mathcal{O}(m)$. The Laplace transformation of the equation for $\boldsymbol{c}_{2}(t)$ is

[TABLE]

where $\widetilde{\Gamma}_{m}^{(2)}$ is a matrix of which the entries are at least of order $m^{0}$, $D_{2,3}$ is a matrix whose elements are proportional to $D_{x}$, $D_{y}$ or $D_{z}$, $\hat{D}_{2,1}$ is a matrix whose elements are proportional to $\hat{D}_{x}$, $\hat{D}_{y}$ or $\hat{D}_{z}$, and $\boldsymbol{c}_{3}(t)$ is a vector of which the elements are all $c_{n_{x},n_{y},n_{z}}$ such that $n_{x}+n_{y}+n_{z}=3$.

By only retaining the zeroth-order contribution in $m$ to Eq. (55), one obtains an equation that is independent of $v$ and holds for $sm/\gamma<<1$ (that is, for times longer than the velocity autocorrelation time). The vector $\boldsymbol{D}$ is of order $m^{-1/2}$, so only the terms of order $m^{1/2}$ of $\boldsymbol{c}_{1}$ are needed, which means that one does not need the matrices $D_{1,2}$, $D_{2,3}$ and $\hat{D}_{2,1}$. Equation (A) can be used to replace the $\boldsymbol{c}_{2}$ in Eq. (A), which, after inverting the Laplace transformation and only retaining terms of order $m^{1/2}$ or lower, becomes

[TABLE]

which gives for the translational flux

[TABLE]

where

[TABLE]

The result for the translational flux together with Eqs. (49) and (55) gives following Fokker-Planck equation:

[TABLE]

Appendix B Gradient expansion

The goal of the gradient expansion is to integrate out the orientational degrees of freedom of the probability density $Q(\boldsymbol{r},\boldsymbol{p},t)$ and obtain a drift-diffusion equation for the probability density $\rho(\boldsymbol{r},t)$ of the position variable Tailleur and Cates (2008); Cates and Tailleur (2013). This expansion is based on a decomposition of the probability density $Q(\boldsymbol{r},\boldsymbol{p},t)$ in spherical harmonics:

[TABLE]

where the scalar $\rho$, the vector $\boldsymbol{\sigma}$ and the matrix $\boldsymbol{\tau}$ are function of $\boldsymbol{r}$ and $t$. These functions are linear combinations of, respectively, the zeroth, first and second spherical harmonic components of $Q$, and $\Omega$ is the projection onto the third and higher-order spherical harmonics. Using this decomposition in Eqs. (52), (53) and (60) and integrating over the orientational degrees of freedom gives

[TABLE]

where

[TABLE]

is the flux in position space. An equation for $\boldsymbol{\sigma}$ can be obtained by multiplying Eq. (53) by $\boldsymbol{p}$ and then integrating over the orientational degrees of freedom. This gives

[TABLE]

where $a$, $b$ and $c$ are the indices of the vectors and matrices, and $\Gamma^{-1}_{ab}$ denotes the $a$ $b$ component of the matrix $\Gamma^{-1}$. Similarly, an equation for $\boldsymbol{\tau}$ can be obtained by multiplying Eq. (53) by $\boldsymbol{p}\boldsymbol{p}-\boldsymbol{1}/3$ and then integrating over the orientational degrees of freedom:

[TABLE]

where the $\partial_{c}\Upsilon_{cab}$ comes from the projection on to the third and higher-order harmonics. Equations (64) (B) (B) are exact, but not closed since $\Upsilon$ is not specified.

We now consider the limit in which the gradients in the system are small in comparison to the persistence length of the ABP Schnitzer (1993); Cates and Tailleur (2013). The relaxation time of $\rho$ is of order $\sim(\nabla)^{-1}$, whereas all the other harmonics relax in times of order $\sim 1$. This implies that in the limit of small gradients, $\rho$ is the slowest mode. We thus assume that the time derivative of $\boldsymbol{\sigma}$ and terms with both a $\boldsymbol{\sigma}$ and a $\nabla$ can be neglected. This yields

[TABLE]

which is related to the polarization by $\boldsymbol{\bar{p}}=\frac{1}{3}\boldsymbol{\sigma}/\rho$. With this expansion the flux becomes

[TABLE]

In a system in steady state with periodic boundary conditions and a magnetic field of the form $B(x)\hat{e}_{z}$, the flux in the x direction must be zero (see Sec. III). This condition together with Eq. (69) yields the steady-state density:

[TABLE]

where the proportionality constant is determined by the normalization of the density, $\delta=\frac{\Delta T}{T+\Delta T}$ and $\Delta T=\frac{f^{2}}{6D_{r}\gamma}$. Using this density to calculate the $y$ component of the steady-state flux gives

[TABLE]

The steady-state polarization in the $x$ and $y$ directions can be determined from Eq. (68) and the density:

[TABLE]

[TABLE]

Appendix C Linear-response theory

Linear-response theory for ABPs has been developed in Ref. Sharma and Brader (2016) to study the average swim speed, and was later applied to a system of interacting ABPs with inhomogeneous activity to study the density profile Sharma and Brader (2017) and flux Merlitz et al. (2018). Here we adapt the linear-response theory for ABPs and apply it to an charged ABP in a space-dependent magnetic field.

[TABLE]

where

[TABLE]

is the operator corresponding to the particle’s diffusive behavior, $\mathcal{L}_{r}=L_{r}$ is the rotational time-evolution operator and

[TABLE]

is the nonequilibrium time-evolution operator. Next, we split the probability density $Q(t)$ in an equilibrium ($Q_{eq}$) and a nonequilibrium ($\delta Q(t)$) part, where the equilibrium part is the steady-state solution to Eq. (74) with $f=0$, which is constant in space. The equation for the nonequilibrium part is

[TABLE]

because $\partial_{t}Q_{eq}=\left(\mathcal{L}_{d}+\mathcal{L}_{r}\right)Q_{eq}=0$.

We are interested in nonequilibrium steady state averages, and therefore assume that the system started out at $t=-\infty$ in the equilibrium state. The previous equation can be integrated from $t=-\infty$ to $t=0$ by treating $\mathcal{L}_{a}Q_{eq}$ as an inhomogeneity:

[TABLE]

where $\delta Q_{ss}=\delta Q(t=\infty)$ is the steady-state solution of the nonequilibrium part of the probability density. This equation is exact and can be used to calculate the steady-state nonequilibrium average of a generic function $g=g(\boldsymbol{r},\boldsymbol{p})$:

[TABLE]

where $Q_{ss}(\boldsymbol{r},\boldsymbol{p})=Q_{eq}+\delta Q_{ss}(\boldsymbol{r},\boldsymbol{p})$ is the nonequilibrium steady-state probability density, the super script $\boldsymbol{r},\boldsymbol{p}$ denotes the variables that are averaged over, the subscript $eq$ indicates that the average is taken with respect to the equilibrium distribution, and a $\dagger$ indicates the adjoint of an operator. The adjoint operators acting on a function $g$ are

[TABLE]

[TABLE]

and

[TABLE]

where the superscript $T$ indicates a transpose.

The response of the system to the activity up to linear-order in $f$ is obtained by removing from the full time-evolution operator the nonequilibrium operator $\mathcal{L}_{a}^{\dagger}$:

[TABLE]

This Green-Kubo relation defines the active transport coefficient $\alpha$ corresponding to an functions $g$ by $\alpha=\lim_{f\to 0}\left(<g>_{lr}-<g>_{eq}\right)/f$. Equation (C) shows that only quantities that are odd functions of $\boldsymbol{p}$ have a nonzero linear response because the angular average in equilibrium any function of an odd power of $\boldsymbol{p}$ yields zero by symmetry.

Appendix D Polarization from linear-response theory

Here we calculate the polarization using the general linear-response theory developed in the previous appendix. The steady-state polarization is defined as the average orientation per particle:

[TABLE]

where $\rho(\boldsymbol{r}^{\prime})$ is the average steady-state density, which is the same as the steady-state probability density of particles as a function of the position alone. Because the density operator $\hat{\rho}(\boldsymbol{r}^{\prime})=\delta(\boldsymbol{r}^{\prime}-\boldsymbol{r})$ is independent of the orientation, the density has no linear response to the activity; therefore, when one only considers contributions up to linear order in the swim force to the polarization, the position-dependent density can be replaced by the bulk density $\rho_{b}$.

The orientation up to linear order in $f$ can be calculated using Eq. (C) to evaluate the average in Eq. (85):

[TABLE]

where the superscript $(1)$ indicates that it is correct up to first order in the self-propulsion force $f$. Because in equilibrium the position is not correlated with the orientation, one can integrate out the orientational degrees of freedom in the average. The orientational integral in the equilibrium average in Eq. (D) is

[TABLE]

where the $\boldsymbol{p}(t)$ is the solution of the stochastic process resulting from Eq. (4) with initial orientation $\boldsymbol{p}(0)$, and the angle brackets indicate an average over realizations of the noise. This is the autocorrelation function for Brownian motion on a unit sphere and can be calculated by expanding the corresponding FPE in spherical harmonics Farage et al. (2015) as

[TABLE]

With the result for the autocorrelation of $\boldsymbol{p}(t)$, Eq. (D) becomes

[TABLE]

where average is taken with respect to the equilibrium probability density for the position variable. By adding an other integral, one can take the gradient of $\Gamma^{-1}$ out of the average:

[TABLE]

which can be written as

[TABLE]

where the response function $\boldsymbol{\chi}\left(\boldsymbol{r}^{\prime},\boldsymbol{r}^{\prime\prime},t\right)$ is

[TABLE]

and

[TABLE]

is similar to a Van Hove function, but with a space-dependent diffusion constant Van Hove (1954); Hansen and McDonald (1990). In the case of a space-independent diffusion tensor, Eq. (93) would become a Van Hove function after integrating out one of the coordinates.

By working out what $\nabla_{\boldsymbol{r}}\cdot\Gamma^{-1}(\boldsymbol{r})$ is want finds that there is no orientation along the direction of the magnetic field. Working out the derivative gives

[TABLE]

Without loss of generality one can consider the case where the magnetic field is oriented along the $z$ direction. In this case the matrices $M$ and $M^{2}$ are

[TABLE]

so clearly there is no orientation along the direction of the magnetic field.

The operator $\mathcal{L}_{d}^{\dagger}$ acting on $\boldsymbol{r}$ in Eq. (D) has the same effect as acting with the adjoint on $\boldsymbol{r}^{\prime}$ (see, for example, Chap. 4.2 of Ref. Risken (1989)); therefore

[TABLE]

where

[TABLE]

and

[TABLE]

was used. In the linear-response calculation of the orientation $\rho(\boldsymbol{r}^{\prime})$ can be replaced with $\rho_{b}$.

Equation (D) can be simplified to

[TABLE]

and after the $\boldsymbol{r}^{\prime\prime}$ integral this becomes

[TABLE]

Multiplying both sides of this equation by $2D_{r}-\mathcal{L}_{d}(\boldsymbol{r})$ results in

[TABLE]

which is equivalent to the following differential equation:

[TABLE]

This equation is the main result of this appendix. Solutions can be obtained numerically, or by further approximations.

For large $D_{r}$ and a magnetic field with small gradients, the contribution from the operator $\mathcal{L}_{d}$ to the left hand side of Eq. (102) is negligible compared to the contribution coming from the part with $D_{r}$; therefore, one can approximate, in this case, the solution to this equation by

[TABLE]

Note that this is the same Eqs. (72) and (73) without the terms of cubic order in $f$.

Even though the density does not have a linear response to the swim speed, the linear response calculations for the orientation can be used to calculate the density and fluxes Sharma and Brader (2017). First, one integrates out the orientational degrees of freedom from Eq. (53):

[TABLE]

with

[TABLE]

where $\boldsymbol{\bar{p}}(\boldsymbol{r},t)\equiv\frac{1}{\rho(\boldsymbol{r},t)}\int d\boldsymbol{p}Q(\boldsymbol{r},\boldsymbol{p},t)$ is the polarization.

When the magnetic field is oriented in the $z$ direction and depends only on the $x$ coordinate, equating the steady-state flux in the $x$ direction to zero gives the following density:

[TABLE]

where the proportionality constant is determined by normalization, and the superscript $(2)$ indicates that it is correct up to second order in the swim force (there is no linear-order response in the density). With this density the steady-state flux in the $y$ direction can be calculated from Eq. (105):

[TABLE]

where $\rho_{b}$ is the bulk density.

Appendix E Improved solution

Only using the linear order contributions in $f$ to the orientation gives poor agreement with simulations for the density; see Fig. 6. The terms of cubic order in $f$ contribute to the diffusion in the system, and by ignoring them there has to be a larger density gradient to compensate for the flux resulting from the orientation. This problem can be overcome by combining the first order result from linear-response theory with the cubic terms from the gradient expansion.

Using the gradient expansion (App. B) we found that the flux is

[TABLE]

where $\boldsymbol{\bar{p}}=\boldsymbol{\bar{p}}(\boldsymbol{r},t)$ is the average polarization at position $\boldsymbol{r}$, which is

[TABLE]

The second term on the right hand side is higher order in the swim force because a gradient in the density is at least second order. The first term on the right hand side of this equation is proportional to the swim force and is equal to the linear-response solution in the limit of small gradients. However, Eq. (102) gives the linear-order response of the polarization without assuming that the gradients are small. We therefore replace the linear order term with the numerical solution to the linear-response equation for the polarization [Eq. (102)] and rewrite the previous equation as

[TABLE]

When the magnetic field is oriented in the $z$ direction and varies in the $x$ direction, one can use this expression for the polarization together with the expression for the flux [Eq. (108)] to obtain the steady-state density by equating the flux in the $x$ direction to zero:

[TABLE]

where $N$ is a normalization constant such that $\int_{0}^{L}dx\rho(x)=1$ and $L$ the size of the periodic box. Note that this expression is the same as the density calculated only using linear-response theory [Eq. (106)], but with $T$ replaced by $T+\Delta T$. The flux in the $y$ direction can be calculated from Eq. (108) together with the result from the density as

[TABLE]

And finally, the density can be used in Eq. (110) to obtain the polarization in the $x$ and $y$ directions in terms of the linear-response polarization as

[TABLE]

and

[TABLE]

Appendix F simulation details and numerics

The Brownian dynamics simulations were done by using the Euler algorithm to discretize the equations of motion in time [Eqs. (2)-(4)]. After each time increment of the orientation vector [Eq. (4)], the vector is rescaled such that the length is always unity. For all simulations a time step of $\Delta t=0.5\ 10^{-6}$ was used. The flux in the system was obtained by sampling the velocity of the ABP. The linear-response equation for the orientation [Eq. (102)] was integrated using the boundary value problem solver from the Scipy Python library Jones et al. (01).