Ionic screening and dissociation are crucial for understanding chemical self-propulsion in polar solvents

Abstract

Polar solvents like water support the bulk dissociation of themselves and their solutes into ions, and the re-association of these ions into neutral molecules in a dynamic equilibrium, e.g., H₂O₂ ⇌ H⁺ + HO₂⁻. Using continuum theory, we study the influence of these association–dissociation reactions on the self-propulsion of colloids driven by surface chemical reactions (chemical swimmers). We find that association–dissociation reactions should have a strong influence on swimmers' behaviour, and therefore should be included in future modelling. In particular, such bulk reactions should permit charged swimmers to propel electrophoretically even if all species involved in the surface reactions are neutral. The bulk reactions also significantly modify the predicted speed of chemical swimmers propelled by ionic currents, by up to an order of magnitude. For swimmers whose surface reactions produce both anions and cations (ionic self-diffusiophoresis), the bulk reactions produce an additional reactive screening length, analogous to the Debye length in electrostatics. This in turn leads to an inverse relationship between swimmer radius and swimming speed, which could provide an alternative explanation for recent experimental observations on Pt-polystyrene Janus swimmers [S. Ebbens et al., Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2012, 85, 020401]. We also use our continuum theory to investigate the effect of the Debye screening length itself, going beyond the infinitely-thin-screening-length approximation used by previous analytical theories. We identify significant departures from this limiting behaviour for micron-sized swimmers under typical experimental conditions and find that the approximation fails entirely for nanoscale swimmers.

---

1 Introduction

The 20th century witnessed a revolution in condensed matter physics, due to the ready availability of well-characterised colloidal particles (1 nm to 10 µm in size). These particles are often viewed as ‘large atoms’: they are small enough to be subject to Brownian motion, and thus to all the machinery of equilibrium statistical physics, but large enough that their microscopic dynamics and interactions can be observed and tuned. Studying colloidal particles has led to fundamental breakthroughs. Most notably, the observation and subsequent understanding of Brownian motion in colloidal systems¹ led to acceptance of the molecular picture of matter.

Moving into the 21st century, physicists have started to recruit colloids to tackle systems that are intrinsically out-of-equilibrium, specifically where the components are themselves self-propelled. This is the field of ‘active matter’. A wide range of novel, self-propelled colloids^2–9 have been synthesised—see example sketches in Fig. 1 for two designs relevant to this work. Such self-propelled colloids are intrinsically out of equilibrium—they continuously transform chemical, thermal, or electromagnetic energy into directed motion—and recent work has focussed on using these systems to experimentally explore exciting non-equilibrium phenomena such as phase separation and collective motion.^10–13

Fig. 1 Cartoon of the two paradigmatic chemical swimmers discussed in the text. Both swimmers move at a few µm s⁻¹ in 10% H₂O₂ solution, powered by the decomposition of H₂O₂ on their surfaces. (a) Bimetallic (typically gold–platinum) rod,² of typical length 2 µm and width 300 nm. The accepted propulsion mechanism for these swimmers is via a H⁺ current, as shown. (b) Platinum–polymer (usually polystyrene) Janus sphere,³ of typical radius 1 µm.

In parallel with this research, much work has gone into understanding the experimental propulsion mechanisms at the level of surface chemical reactions.^{2–4,14–17} Working out how these tiny motors function is a fundamental problem in its own right. However, understanding the propulsion mechanism is also an essential first step in understanding the experimental collective behaviour. This is because unlike biological swimmers such as E. coli, where the biochemical reactions responsible for propulsion take place internally, synthetic swimmers are usually propelled by surface chemical reactions that directly modify the chemical, electrostatic, or temperature fields of their surroundings. These fields modify the propulsion speed of other swimmers, generating so-called ‘phoretic’ swimmer–swimmer interactions in addition to the hydrodynamic and contact interactions experienced by all swimmers.^18–21 As these phoretic interactions are directly coupled to the chemical reactions responsible for propulsion, knowing how artificial swimmers self-propel is essential for understanding their collective behaviour.

This bottom-up approach contrasts with the tactic employed by most theoretical modellers of active matter, which is to explore the phenomenology arising from minimal or effective models of swimmer–swimmer interactions, by which we mean models that have no explicit connection with microscopic parameters. Minimal models include the Vicsek model²² or the active Brownian model^23,24 which considers only contact forces, while more complex, but still effective, hydrodynamic^25–27 and phoretic models^18,28–31 are now also being explored. This research is valuable in unlocking generic non-equilibrium physics principles and is often able to reproduce experimental behaviour (phase separation, etc.) surprisingly well. However, it is often unclear whether such agreement is due to judicious tuning of parameters, or whether there is a real correspondence in generic physics so that details do not matter. In lieu of a general theory for non-equilibrium physics, only in-depth knowledge of the microscopic physics for specific swimmers is likely to allow the resolution of this issue. If a microscopically justified model predicts the same phenomenon as generic models, then details may indeed not matter. On the other hand, if a microscopic model fails, then new physics is indicated.

That said, current understanding of self-propulsion mechanisms is often very incomplete, and hence the necessary foundations for models which can reproduce multiparticle behaviour from a bottom-up perspective are lacking. In particular, most research so far has focussed on unravelling the surface chemistry of the swimmer itself.^{2,14–17,32–34} This is understandable and necessary, but it has meant that other aspects, such as the chemistry of the bulk solvent, have been neglected.

In this article, we show that taking into account the chemistry of common polar solvents such as water and hydrogen peroxide, significantly, and often qualitatively, modifies the predicted propulsion behaviour of almost all self-propelled synthetic swimmers. In this first study, we limit our discussion to the propulsion of single particles, because this is a necessary first step in understanding more complex behaviour. For most of the paper, we also do not go into the details of the surface chemistry, which are not well understood. This allows us to highlight the bulk solvent effects, at the expense of explicit predictions for the propulsion speed. However, in Section 7, in order to compare our results with experiments, we do examine the predicted propulsion speed obtained with suitably simple assumptions for the surface chemistry. One of the main conclusions of our paper is that more detailed experimental studies of both the surface and bulk chemistry are crucial for a detailed understanding of self-propulsion.

2 Chemical propulsion

2.1 Self-electrophoresis

We discuss here the most experimentally typical self-propelled colloids, which we term ‘chemical swimmers’. They are most easily defined by example. Fig. 1 shows two chemical swimmers, both powered by the catalytic decomposition of hydrogen peroxide on their surfaces. Because the colloid surface is anisotropic, this reaction produces chemical gradients which, via interaction with the particle surface, eventually lead to self-propulsion.

We say ‘eventually’ because the propulsion mechanism of these swimmers is somewhat involved. For the example given in Fig. 1a, H₂O₂ decomposition does not occur just by the simple chemical reaction

2H₂O₂ → 2H₂O + O₂, (R1)

but also occurs partially electrochemically, with two half reactions taking place preferentially on the Au or Pt surfaces^2,4 for exampleThese half reactions produce a proton gradient outside the colloid, which generates a local electric field. The colloid surface, like most surfaces in water, is charged, so this electric field causes electroosmotic flow over the colloid surface, leading to self-propulsion. This propulsion mechanism is called ‘self-electrophoresis’.⁴ The electric field also generates a proton current outside, which is balanced by an electric current inside the conductive swimmer.

A large body of experimental evidence confirms that self-electrophoresis is the appropriate propulsion mechanism for these bimetallic swimmers. For example, their propulsion speed scales inversely with salt concentration,^11,16 which is expected from a simple application of Ohm's law. Recent results^16,17 indicate that self-electrophoresis is also the appropriate propulsion mechanism for the type of colloid shown in Fig. 1b, which has a single metallic coating. This is at first surprising because there is no obvious mechanism for producing the ionic gradient needed for self-electrophoresis. However, geometrical differences between the equator and pole of the catalytic coating, such as thickness variation, may couple to the half-reaction rates in (R2) and so provide the necessary asymmetry.¹⁷ In this paper, we go further and show that these effects are not limited to swimmers that can support ionic currents. All swimmers in aqueous solution are likely to be self-electrophoretic to a major degree, whatever their surface reaction mechanism.

2.2 Surface chemistry

Before we discuss the effects of the bulk solution, we point out one general difficulty with self-electrophoresis that will also apply to other complex propulsion mechanisms. This is that the relevant surface reaction rates are extremely hard to measure. The overall reaction rate (R1) can be easily obtained by measuring reactant or product concentrations,^16,17 but for self-electrophoresis the important rate is the proton production rate (R2), and this is likely to make up only a tiny proportion of the overall reaction, with the remainder proceeding via(R1).² Measuring the rate of an individual reaction pathway like (R2) is challenging, and has not yet been done, to our knowledge, for any self-propelled particle.

This would not be a problem if we could predict these rates. However, surface catalysis is a notoriously sensitive phenomenon, and these surface reaction rates are likely to vary unpredictably with almost every parameter, e.g., pH, ionic strength, and surface roughness.³⁵ Reversing the argument, the only currently available method of estimating these reaction rates is from the particle propulsion velocity itself. That is, with a sufficiently accurate microscopic model, the surface reaction rates can be inferred from the propulsion speed.³² The catalytic chemistry of micro- and nanoparticles is of huge industrial importance, so this provides another major motivation for obtaining a detailed theoretical understanding of self-electrophoresis.

2.3 Bulk chemistry

At first glance, the chemistry of the bulk solution is much simpler than that of the surface. However, a polar solvent such as water presents two complications which have not yet been fully taken into account. The first of these is electrostatic screening. In so-called phoretic mechanisms, such as self-electrophoresis, fluid flow is generated in a layer around the particle surface. In self-electrophoresis, the thickness of this interaction layer is given by the electrostatic screening or Debye length κ⁻¹: outside this screening layer, the free charge density, which is responsible for fluid flow, decays rapidly to zero.

Analytical studies typically make use of a thin-screening approximation κa ≫ 1, where a is the swimmer radius, because this dramatically simplifies the calculation of propulsion speed, flow fields, etc.^16,17,33,36 However, this assumption is not generally valid: κ⁻¹ is of order 100 nm for an experimentally typical 3 M H₂O₂ solution,¹⁶ and active colloids typically range in size from 10 nm³⁷ to 10 µm.¹⁵ In this paper, we show that the thin-screening approximation can be dropped from analytical calculations, and that this dramatically reduces the predicted propulsion speeds, by up to several orders of magnitude for nanoscale swimmers. To the best of our knowledge, this result has not been shown even with numerics, such as the finite element method (FEM), for which the thin screening approximation is not employed.^14,34 Here, we use FEM calculations to verify the analytical results.

The second complication of aqueous and similar polar environments is that the solvent is not chemically inert. It is an ‘active fluid’ that can be driven out of chemical equilibrium by the reactions on the particle surface. This consideration has been appreciated for biological fluids such as the cytoplasm,²⁰ where biomolecules are continuously synthesized and broken down, but it is also true for simple fluids like water which permit the ionic dissociation of both themselves and any polar solutes, e.g., H₂O ⇌ H⁺ + OH⁻ and H₂O₂ ⇌ H⁺ + HO₂⁻. The implications of these reactions for self-electrophoresis are the main focus of this paper. The most striking implication is that a gradient of a neutral molecule like H₂O₂ will result in ionic gradients, here of H⁺ and HO₂⁻ ions, which will themselves produce electric fields. This means that a surface reaction with only uncharged species like H₂O₂ or H₂O can itself generate self-electrophoretic propulsion.

It is worth highlighting that these bulk reactions should also qualitatively modify interparticle interactions. To demonstrate this, we describe an effect called ‘reactive screening’,²⁰ which will underlie much of our later discussion. We illustrate this effect with a simple 1D model, see Fig. 2. Let an uncharged molecule of diffusivity D be produced uniformly at a plane surface z = 0 and consumed in the bulk (z > 0) with rate γ. The steady-state concentration profile c(z) then obeys

where diffusion, on the left, is balanced by consumption, on the right. The solution to eqn (1) is an exponentially screened concentration profile, c = c₀ exp(−qz), with c₀ the concentration at the surface, and the reactive screening length q⁻¹ = (D/γ)^1/2. This uncharged model has been applied to the diffusiophoresis of small particles inside a biological cell, where the relevant bulk reactions are the breakdown of biomolecules in the cytoplasm.²⁰ With self-electrophoresis, as we shall see, reactive screening can also exponentially screen the electrostatic potential, effectively turning off the long-range electrostatic interactions, which would otherwise be inevitable. Hence, this reactive screening is a qualitatively new effect of bulk reactions, which cannot be ignored a priori even at the level of phenomenological models of collective motion.

Fig. 2 Schematic representation of the simple 1D molecular screening model introduced in the text. A wall releases molecules (red triangles) which diffuse with D (black arrow) and are consumed in the bulk with a rate γ (blue symbol). This leads to an exponential decay of the concentration, as indicated using the dashed red line and continuum-level red gradient.

3 Overview of main results

The theory of self-electrophoretic propulsion is mathematically involved, even without the introduction of additional bulk reactions, so we will use this section to sketch out our main results in advance. This will necessarily skim over or simplify many relevant details that will be addressed fully in later sections.

3.1 Overall framework

Our main mathematical result is that the self-electrophoretic propulsion speed of an arbitrary, uniformly charged, spherical swimmer can be written, if a suitable linearization is applied, in the form

U = U^SM(j^s, c^salt, σ,…)F(κa)B(qa,…), (2)

where U^SM is the ‘standard model’ propulsion speed assuming the thin screening limit without bulk reactions, explored for example in ref. 33 and 38. U^SM depends on, among other parameters, the surface reaction rates j^s, the salt concentration c^salt, and the surface charge density σ. Here we introduce two new factors, F and B, to account for realistic electrostatic screening and bulk reactions, respectively. These factors depend on the dimensionless parameters κa and qa, with κ and q the inverse electrostatic and inverse reactive screening lengths, and a the swimmer radius.

3.2 Electrostatic screening

The new dimensionless factor F(κa) is exactly analogous to the well-known function f(κa) (Henry's function³⁹) that controls the speed of a particle undergoing electrophoresis in an external field via U_ext = μ_EE^∞, with μ_E = ζεf(κa)/η the electrophoretic mobility, ε the dielectric constant, ζ the particle's surface potential, η the solution viscosity, and E^∞ the external electric field.^39,40 Both f and F are plotted in Fig. 3. For κa ≪ 1, i.e., for small particles or low salt concentration, F decreases rapidly, scaling as (κa)³. This is different from external electrophoresis because of the different geometries of the driving fields—a uniform field for external electrophoresis compared with a dipole for self-electrophoresis.

Fig. 3 (a) Schematic showing the difference in boundary conditions between external electrophoresis (upper) and self electrophoresis (lower). Thick, coloured arrows show the direction of motion of positively charged particles. (b) Henry's function f(κa) which determines the mobility of a particle in an external electric field (—), and F(κa), the equivalent function for self-electrophoresis ().

The implication of this a³ scaling is that, other things being equal, nanoswimmers should swim much slower than microswimmers. Experimentally, however, nanoswimmers are found to swim faster than equivalent microswimmers.³⁷ From this we conclude that other things are not equal: either the surface reaction rates are much larger for nanoswimmers, or the standard self-electrophoresis theory does not apply for these small swimmers. If this issue can be resolved, which we do not attempt here, it will likely also give insight into the related phenomenon of directed motion in nanoscale biological enzymes.^41,42

3.3 Bulk reactions

The effects of ionic dissociation depend upon the nature of the surface reaction responsible for propulsion. A common feature is the importance of the reactive screening length q⁻¹ which controls the propulsion behaviour through the parameter qa. We can understand why qa is the relevant parameter as follows: for qa ≪ 1, the swimmer is smaller than the reactive screening length, so any molecules produced at the swimmer surface will diffuse away or return to the swimmer surface before they have time to react. In this ‘reactionless limit’, the swimmer will behave as though there are no bulk reactions, which is the usual, tacit assumption. For swimmers larger than the reactive screening length, qa ≫ 1, we are in a ‘reactive limit’ where the bulk ionic reactions dominate the behaviour. For typical experimental conditions, e.g., 3 mol L⁻¹ H₂O₂, we find a reactive screening length q⁻¹ ≈ 70 nm, which is in the centre of the experimental range of swimmer radii.^15,37 Both the reactionless and reactive limits, and the intermediate regime (qa ≈ 1), should therefore be experimentally relevant.

We now explain the effect of bulk reactions on specific types of swimmer. The overall surface reaction we focus on is the H₂O₂ decomposition reaction (R1). As we have seen, this overall reaction can occur through several different pathways. We therefore define three model swimmers, shown in Fig. 4 with surface reactions that are representative of these different pathways. A real swimmer might exhibit any or all of these.

Fig. 4 Schematic of the effect of bulk ionic reactions on the propulsion of three model swimmers. The upper panel shows the system without and the lower panel with bulk reactions. Coloured arrows indicate fluxes of three chemical species H₂O₂ (purple), H⁺ (red) and HO₂⁻ (blue). The thickness of the arrows corresponds very roughly to the relative intensity of the fluxes. White arrows denote the direction of particle propulsion (x = no propulsion). Arrow length indicates relative speed. Dashed semicircles show the approximate extent of the reactive screening length q⁻¹.

The upper panels (a–c) of Fig. 4 show these model swimmers without bulk reactions. In (a), there is a single surface flux of neutral H₂O₂ molecules. This models the purely neutral decomposition of H₂O₂ in reaction (R1). Here we first make three general points: first, O₂ and H₂O are not included here for modelling simplicity. This is justified because O₂ does not dissociate, and H₂O dissociates much less than H₂O₂, so neither species should contribute strongly to self-electrophoresis. Second, we assume that the chemical fluxes are dipolar rather than monopolar, i.e., they both leave and enter the particle surface. For a uniformly charged particle, which we assume for simplicity, only the dipolar component of the flux contributes to the propulsion speed,³³ so our choice of a dipolar surface-flux profile does not affect our results, and simplifies the argument. Third, we consider only self-electrophoresis here, and ignore ‘neutral self-diffusiophoresis’, which is propulsion generated by a direct non-electrostatic interaction between a neutral species, such as O₂ and the swimmer surface,^3,43 and which is typically much weaker than self-electrophoresis.^16,34,44 Hence the model swimmer in (a) does not move, because, without bulk reactions, a surface reaction involving only uncharged species cannot generate electric fields, and therefore cannot produce self-electrophoresis.

In (b) a surface proton flux generates an electric field viareaction (R2). We assume the particle is positively charged, and it then swims in the direction indicated by the white arrow. This corresponds to the standard self-electrophoresis model, e.g., for Au–Pt swimmers.² The electric field has a dipolar form, like the proton flux.

In (c), we have a third mechanism, with equal fluxes of H⁺ and HO₂⁻ ions. There is no net electrical current for this swimmer since there are equal positive and negative fluxes. However, self-propulsion still occurs. This is because the two ions diffuse at different rates (H⁺ faster than HO₂⁻), and this creates a so-called diffusion potential, which acts to prevent net charge separation. The diffusion potential leads to an associated (self-generated) electric field, which then produces motion via electrophoresis in the usual way (white arrow). This propulsion mechanism is called ‘ionic diffusiophoresis’,⁴⁵ and is typically used to model swimmers composed of solid salts, which generate propulsion through dissolution of the swimmer itself, e.g., AgCl(s) → Ag⁺(aq) + Cl⁻(aq).⁴⁶ We include this model here because ionic diffusiophoresis may contribute to the propulsion of Pt-Janus swimmers, for example via

with subsequent recombination of H⁺ and OH⁻ in the bulk. However, note that reaction (R3) is not the reaction shown in Fig. 4c, where HO₂⁻ is used instead of OH⁻ for modelling simplicity.

The lower half of Fig. 4(d–f) shows the effect on each of these swimmers of a single ionic reaction occurring in the bulk, aqueous phase

H₂O₂ ⇌ H⁺ + HO₂⁻. (R4)

As we mentioned before, this reaction will only begin to have a significant effect when we are in the reactive, qa > 1 regime. In Fig. 4, the reactive screening length q⁻¹ is indicated by the dashed line: for these particular swimmers, qa ≈ 3. The white arrows show the qualitative effect of this reaction on the propulsion speed, which is different for each of the swimmers. For (a→d), the reaction generates propulsion, for (b→e) the reaction increases the propulsion speed, and for (c→f) the propulsion speed falls. We now briefly explain the reason for these effects.

In the absence of a swimmer, reaction (R4) is in a state of dynamic equilibrium. If a swimmer consumes or produces molecules on either side of this equilibrium, then this will push the reaction out of equilibrium, and the system will respond so as to reduce the effect of that perturbation: this is Le Chatelier's principle. Thus, in Fig. 4d, the H₂O₂ flux injected from the particle surface is to the left of the equilibrium, so the H₂O₂ partially dissociates into ions, producing ionic fluxes. In Fig. 4e the proton flux is to the right of equilibrium, so there is net ionic recombination in the bulk to give an H₂O₂ flux (this small H₂O₂ flux is not shown because it does not significantly contribute to self-electrophoresis) and an HO₂⁻ flux in the opposite direction to the original proton flux. In Fig. 4f, the proton and HO₂⁻ fluxes are both to the right of equilibrium, so these both recombine with counterions in the bulk to give an H₂O₂ flux instead of the ionic fluxes.

For (a→d), the new ionic fluxes produce a diffusion potential, which generates motion. Hence a swimmer without any electrochemical reactions on its surface can still exhibit self-electrophoretic propulsion. Crucially, it will also display the experimental behaviour that would be expected of a self-electrophoretic swimmer, e.g., propulsion speed scaling inversely with salt concentration (via the U^SM factor in eqn (2)). This means that the kind of ionic behaviour observed in ref. 16, 17 and 32 does not a priori require an electrochemical surface reaction. In practice, however, we find that, because of the weak dissociation of H₂O₂, the simple non-electrochemical surface reaction mechanism in Fig. 4a→d cannot account for the magnitude of the experimentally observed propulsion in, e.g., ref. 16: genuine self-electrophoretic propulsion is still required.

In (b→e), the important point is that chemical reactions conserve charge, and this also implies the conservation of electrical current. There is a net electrical current in (b), and because this current is conserved it will have the same magnitude with or without bulk reactions. It is only the identity of the current-carrying ions which changes: in this case, the current becomes partially carried by HO₂⁻ ions travelling in the opposite direction, see Fig. 5. As we discuss later, the propulsion speed scales inversely with the diffusivity of the current-carrying ion. In this case HO₂⁻ diffuses approximately 10 times slower than H⁺, and this is why the speed increases. In fact, in the appropriate environment of high pH (= high HO₂⁻ concentration), the predicted speed increases ten-fold because the current becomes entirely carried by HO₂⁻ ions.

Fig. 5 The effect of bulk reactions on ionic currents. Regions with excess H⁺ ions become depleted in HO₂⁻ ions and vice versa. This means that an initial H⁺ current flowing from excess to depleted regions (left) is partially replaced by a HO₂⁻ current flowing in the opposite direction (right).

In (c→f), on the other hand, there is no net electric current to be conserved and both anions and cations react freely with their counterions in the bulk. Hence, far from the swimmer, the ionic gradients become vanishingly small, with a resultant drop in propulsion speed compared to the case without bulk reactions. In detail, the presence of ions in the bulk, due to the surface reactions, generates a diffusion potential (similar to the situation in Fig. 4d). However, since the ions in Fig. 4f can recombine through bulk reactions, the further one is from the swimmer surface, the fewer ions generated by the surface reaction remain to induce the diffusion potential. This shows up as an exponentially screened potential, with screening length q⁻¹. This also affects the swimming speed, because we find that the magnitude of the diffusion potential scales with the thickness of the screening layer, leading to a scaling of U ∝ 1/(qa) for large swimmers. This lowers the predicted propulsion speed by up to a factor of approximately 100 for the largest ≈10 µm radius swimmers.

3.4 Comparison with experiments

As previously discussed, measurements of relevant surface reaction rates are not currently available. This makes it difficult for us (or others) to predict propulsion speeds that can immediately be compared with experiments. Nevertheless, we will present some speed calculations with the simplest assumption of fixed surface reaction rates in Section 7. These comparisons indicate that the U ∝ 1/a scaling observed with Pt-polymer Janus particles¹⁵ might come from an ionic diffusiophoretic mechanism like that discussed above, see Fig. 4c→f. In addition, our results indicate that the speed of nanoscale swimmers is too high to be explained by self-electrophoresis with fixed surface reaction rates, see Section 3.2.

3.5 Summary

In brief, we find that ionic reactions and electrostatic screening should have very significant and system-dependent effects on the propulsion of a wide range of chemical swimmers. These effects include increasing or decreasing the predicted speed by several orders of magnitude, as well as the qualitatively new behaviour of reactive screening. In addition, we find that even swimmers with no ionic surface reaction can behave as though they are self-electrophoretic. The remainder of this paper provides a detailed account of the theory that gives rise to these results, and compares them to experiments as far as current data allow.

4 Theoretical model

In this and the following two sections, we present a quantitative model of self-electrophoresis. Here in Section 4, we will lay out the general theoretical model and detail how this will be applied to the specific H₂O₂ reaction system described above. In order to obtain analytical results we also linearize our theory. We will then apply this model to obtain explicit results, first for a system with only surface reactions (Section 5), and then with bulk reactions (Section 6).

4.1 General model

The standard theoretical approach to self-electrophoresis involves coupling the chemical fluxes arising from reactions on the particle surface to bulk differential equations (Nernst–Planck, Poisson, and Navier–Stokes).^14,34 This treatment generally ignores bulk chemical reactions by assuming that each chemical species is conserved. We adopt the standard approach, but include bulk reactions by coupling chemical fluxes to local reaction rates. We solve this model numerically using COMSOL. Separately, and in common with previous work,³³ we also linearize the model to obtain an analytical approximation. Unlike in previous work, the analytical solution does not require the assumption of a thin electrostatic screening layer.

We consider a spherical swimmer of radius a and uniform surface charge density σ, see Fig. 6. The electrostatic boundary condition of such a particle is

with ϕ the electrostatic potential field and ε the dielectric constant of the fluid (the dielectric constant of the particle is assumed to be zero). Here (s) and n̂ indicate evaluation at, and the normal out of, the particle surface. We have chosen a uniform, dielectric boundary condition for simplicity. However, in Appendix A.5 we show formally that, with an appropriate choice of surface potential, an equipotential (conducting) surface gives the same swimming speed as a dielectric. We do not deal with mixed dielectric/conducting particles here, but this should not qualitatively affect the basic physics of the self-propulsion.

Fig. 6 Diagram of a model swimmer, highlighting the distinction between bulk and surface parameters. In the bulk, we have an electrostatic potential field ϕ, chemical concentration fields c_l, a pressure field p, and fluid velocity u. Far from the particle, these fields approach uniform values (superscript ∞). On the particle surface, the uniform surface charge density σ and the nonuniform molecular fluxes out of the surface j^s_l set boundary conditions for the bulk potential and concentration fields, respectively. The particle, of radius a, is axisymmetric around the z axis, so the surface fluxes are parameterized by the polar angle θ. The swimming velocity, which we calculate, is U = Uẑ.

Propulsion is generated by reactions on the swimmer surface. These reactions produce and consume N different chemical species, labelled l = 1…N. The surface production (or, if negative, consumption) rate per unit area of each species is j^s_l(θ), and is a function only of θ, the polar angle with respect to the symmetry axis ẑ, see Fig. 6. The surface reaction rates can be equated to the bulk flux j_l of each species out of the particle surface, giving the boundary condition

n̂·j_l(s) = j^s_l(s). (4)

These bulk fluxes obey the classical Nernst–Planck equation³⁴with z_l, D_l, and c_l, respectively the valence, diffusivity, and concentration field of each chemical species, u the fluid flow field, e the fundamental charge, k_B Boltzmann's constant, and T temperature. Physically, eqn (5) expresses the bulk fluxes as linear sums of advective, diffusive, and conductive terms respectively. Eqn (5) is the standard flux expression used in studies of self-electrophoresis. Its main simplification is the neglect of cross-coupling terms between the molecular fluxes, and this is valid as long as we are in the dilute limit with relatively small ionic gradients,⁴⁴ which is true here.

Without bulk reactions, conservation of chemical species would require that the bulk fluxes are incompressible vector fields, i.e., ∇·j_l = 0. This is the standard approach, see ref. 14 and 34, and Section 5 here. Bulk chemical reactions can be incorporated by writing instead

∇·j_l = R_l(c₁…c_N), (6)

where R_l is the local rate at which each chemical is produced (if negative, consumed), in chemical reactions. In general, R_l depends on the local concentration of all chemical species involved in reactions with species l. Note that chemical reactions are charge-conserving, i.e.,everywhere, and combining this condition with eqn (6) implies the conservation of electrical current

∇·i = 0. (8)

where the electrical current .

Infinitely far from the particle, the chemical concentrations are labelled c^∞_l, and are determined by equilibrium equations and charge neutrality. The other boundary conditions at infinity are j^∞_l = 0, ϕ^∞ = 0, u^∞ = 0 (in the lab frame), and p^∞ = p_atm, where p is the hydrostatic pressure field and p_atm is the atmospheric pressure, whose absolute value does not affect the calculations.

The electrostatic potential ϕ is determined by the Poisson equation

ε∇²ϕ = −ρ_e, (9)

with charge density . The interaction of the electrostatic potential and the unbalanced charge density (ρ_e ≠ 0) generates a force density

f = −ρ_e∇ϕ, (10)

and this drives fluid flow via the Stokes equations for low-inertia, incompressible flow

η∇²u = ∇p − f, (11)

∇·u = 0. (12)

Finally, the swimmer is not held in place, so fluid flow around it will cause it to move with some propulsion velocity U. This propulsion velocity is determined by the condition that there is no net force acting on the total system of swimmer plus fluid out to infinity.⁴⁷ This force-free condition is simply a reflection of the fact that all the forces are internal to this total system—there are no long-range, external forces like gravity. The force-free condition can be translated into an expression for U by using the Lorentz reciprocal theorem, which is a restatement of the Stokes equations in integral form. This gives a closed-form expression for the propulsion velocity⁴⁸where r is the distance from the particle centre, r̂ and [small straight theta, Greek, circumflex] are unit vectors in the r and θ directions, and the scalar speed U is defined by U = Uẑ. The volume integral is over the region outside the particle.

4.2 Numerical solution

We solve the full non-linear model numerically using FEM implemented in COMSOL. To do this, we make several modifications to the above equation system. In particular, we define a new force density f_FEM to replace f in eqn (10). The two quantities are related byThis redefinition does not influence the result of the calculation, but limits spurious flows related to numerical artefacts in the electrostatics.⁴⁹ Our calculations are performed on an axisymmetric spherical domain of radius L = 10a + 25κ⁻¹, which we verified to be sufficient to eliminate most finite-size effects in our speed calculations. This frame co-moves with the colloid. We impose no-slip at the colloid surface, and on the edge of the domain we employ the same boundary conditions as the theory has at infinity. For the fluid velocity we impose a no-stress condition on the edge of the domain

[η(∇u + (∇u)*) − p[Doublestruck I]]·n̂ = 0, (15)

with n̂ the normal to the domain, * denoting transposition, and [Doublestruck I] the 3D identity matrix. This is equivalent to imposing a force-free condition on the swimmer–fluid system.^44,50

Our technique is to first obtain approximate numerical solutions for the electrostatic and concentration fields in the absence of advection, so neglecting the first term in eqn (5). This approach is justified because experimental swimmers generally have low Péclet numbers, i.e., molecular diffusion D dominates over advection. The Péclet number is defined as Pe = Ua/D, with U ≈ 10 µm s⁻¹, a ≈ 1 µm, and D ≈ 10⁻⁹ m² s⁻¹ typical for experiments on microswimmers, leading to Pe ≈ 0.01.

The flow field is then computed self-consistently on the domain by employing the force density, eqn (14), following from the concentration and potential fields. The speed of the swimmer is determined by taking the average of the fluid velocity on the edge of the domain: U = −〈u〉_r=L, where U is in the lab frame and u in the co-moving frame. We subsequently verified the low-Pe approximation by solving the fully coupled equations (with advection) directly in a limited number of cases, which gave agreement to within a few per cent. See Appendix C for full details of the numerical calculations.

4.3 Analytical solution

We also linearize the model to provide an analytical solution. To do this we assume that the fields ϕ and c_l have only small deviations from their values in the uncharged, unreactive state where ϕ = 0 and c_l = c^∞_l everywhere (for ϕ, this assumption corresponds to the usual Debye–Hückel approximation, ϕ ≪ k_BT/e). We then expand the model to linear order in the small dimensionless parameters ψ = ϕe/(k_BT) and x_l = (c_l − c^∞_l)/c^∞_l. Applying this linearization to eqn (5) gives

j_l = −c^∞_lD_l[∇x_l + z_l∇ψ], (16)

where the advection term has been dropped entirely because u scales quadratically with the small parameters (eqn (10) contains a product of ρ_e and ϕ, which are both small). We must also Taylor expand the production rates R_l to linear order, i.e.,where [scr O, script letter O] (·) means ‘of order ·’, and the elementsare components of a matrix which we can call the linear reaction matrix. Its meaning will become clearer when we consider specific reactions. From eqn (6) and (16), we have, to linear orderThis set of N equations, together with the Poisson equation, which we rewrite asmakes up a system of N + 1 linear differential equations in N + 1 fields (x_l, l ∈ 1…N and ψ).

This system of equations is soluble in a spherical geometry by standard spectral methods, and the electrostatic potential field so obtained can then be used to calculate the propulsion speed by evaluating the integral in eqn (13). In doing this, we make the further usual assumption of a relatively small driving field.³⁹ That is, if we define ϕ = ϕ^eq + ϕ^sr where ϕ^eq is the electrostatic potential in the absence of surface reactions and ϕ^sr is the additional potential generated by these reactions then ϕ^sr ≪ ϕ^eq. As a result, the surface reaction rates j^s_l, which only come into ϕ^sr, contribute linearly to the final velocity. The algebra required to solve eqn (19) and (20) is significant, so we go through this explicitly in Appendix A.

4.4 Specific H₂O₂ reaction model

The chemical reaction system we consider is the simplified version of the H₂O₂ reaction system described in Section 3. On the particle surface, H₂O₂ decomposes into O₂ and H₂O. For simplicity, however, we ignore both products of this reaction: O₂ because it is electrically neutral and does not dissociate, and H₂O because it dissociates much less than H₂O₂—the respective equilibrium constants⁵¹ are K_eq,H2O = 1.0 × 10⁻¹⁴ mol L⁻¹ (pH = 7) and K_eq,H2O2 = 2.5 × 10⁻¹² mol L⁻¹. In the bulk, we ignore any slow decomposition of H₂O₂viareaction (R1), and the only bulk reaction we consider is the ionic dissociation reaction (R4), which we rewrite hereWe therefore have only three chemically active species, with associated subscripts in brackets: H₂O₂ (∘), H⁺ (+), and HO₂⁻ (−). Protonation reactions like (R4′) are normally extremely rapid, with kinetics controlled by the diffusion and collision of the ions,⁵² and with simple first order rate expressions

k_forward = k_disc_∘,

k_reverse = k_asc₊c₋, (21)

where we estimate the association rate constant to be k_as = 4.9 × 10¹⁰ mol⁻¹ L s⁻¹ using the Smoluchowski–Debye theory for diffusion-limited reactions, see Appendix B. This theoretical value agrees closely with experimentally measured rates for similar reactions,⁵²e.g., H⁺ + HCO₃⁻ ⇌ H₂CO₃ in water has k_as = 5 × 10¹⁰ mol⁻¹ L s⁻¹ (ref. 53).

The dissociation rate constant k_dis = 0.12 s⁻¹ is then determined from the equilibrium constant K_eq = k_dis/k_as = 2.5 × 10⁻¹² mol L⁻¹ (ref. 51). Far from the particle, the system is in equilibrium, so we have

c^∞₊c^∞₋ = K_eqc^∞_∘, (22)

The production rates are R₊ = R₋ = k_forward − k_reverse and R_∘ = k_reverse − k_forward, and linearizing using eqn (18) gives the linear reaction matrixwhere the order of rows and columns is ∘, +, −.

The three reactive species have diffusivities D_∘ = 1.7 × 10⁻⁹ m² s⁻¹ (ref. 54), D₊ = 9.3 × 10⁻⁹ m² s⁻¹ (ref. 55), and D₋ = 0.9 × 10⁻⁹ m² s⁻¹ (ref. 56). We also have two unreactive ions, which we take to be Na⁺ and Cl⁻ with diffusivities D_Na⁺ = 1.3 × 10⁻⁹ m² s⁻¹ and D_Cl⁻ = 2.0 × 10⁻⁹ m² s⁻¹ (ref. 57). Because these ions are not involved in chemical reactions at the surface or in the bulk, their concentration fields are in equilibrium with the electrostatic potential. The implication is that the diffusivity of these ions does not contribute to the propulsion speed in the linear regime. We show this mathematically in Appendix A.1.

The chemical concentrations at infinity are determined by the chemical equilibrium, eqn (22) and by charge balance

These two equations (eqn (22) and (24)) connect five concentrations, so we can set three concentrations freely. In practice, we choose instead to set the H₂O₂ concentration, the total ionic strength, and the pH. The reaction scheme presented here is the simplest possible that gives the necessary freedom: bulk ionic dissociation reactions require at least three reactive species, and the two non-reactive ions are necessary to allow the ionic strength to be modified independent of other parameters.

For the variable parameters, our base set, used unless specified otherwise, is 1 mmol L⁻¹ salt, i.e., , a = 500 nm, and c^∞_∘ = 3 mol L⁻¹. For these parameters, κ⁻¹ = 10 nm, and c^∞₊ = c^∞₋ = 3 × 10⁻⁶ mol L⁻¹. These values were chosen because micron-sized particles and H₂O₂ concentrations of order 3 mol L⁻¹ are experimentally typical,¹⁶ while the 1 mmol L⁻¹ baseline salt concentration allows us to scan a wide range of the important parameter κa for realistically sized particles.

Meanwhile, the surface reactions are specified by surface fluxes j^s_l, l ∈ {∘, +, −}, of the three active species. We consider the three model swimmers shown in Fig. 4, referred to as: S_∘, the nominally neutral swimmer; S₊, powered by a proton current; and S₌, powered by ionic diffusiophoresis. As mentioned above, only the dipolar part of the fluxes, that is the 1st Legendre component, contributes to the propulsion speed of uniformly charged swimmers,³³ so we include only this term by setting j^s_l ≡ j^s_l,1 cos θ for each surface flux where j^s_l,1 is a constant coefficient. For S_∘, only j^s_∘,1 is finite; for S₊, only j^s_+,1 is finite; and for S₌, j^s_+,1 and j^s_−,1 are equal and finite, with j^s_∘,1 = 0.

Our model makes a number of simplifications. This includes those chemical simplifications already discussed, as well as the neglect of potential contaminants such as CO₂, which also undergo ionic dissociation. The main purpose of this paper is to illustrate the physical principles behind the effect of bulk reactions on self-propulsion. These physical principles will also apply to a more complex and realistic H₂O₂ reaction system, as well as to other chemical systems.^58,59 We also neglect any dependence of the surface parameters on environmental conditions, so, for example, we take σ = constant, independent of pH, salt concentration etc. This does not imply that surface parameters are independent of the environment; it is just that detailed knowledge of this dependence is currently lacking. The ‘pure’ effect of bulk reactions which we capture will occur in addition to any such interdependence.

5 Electrophoresis without bulk reactions

Before discussing the effect of bulk ionic reactions, it is important to set out the basic theory for propulsion by self-electrophoresis without such reactions. This theory has been set out multiple times before for the limit of vanishing electrostatic screening length, κa ≫ 1, ref. 14, 33, 34, 36, and 60. Here we extend the theory to include the effect of a finite κa.

For comparison, we first write down the standard results for electrophoresis in an external, linear field,⁴⁰Fig. 3a (top). We align the external field along the z axis and consider the velocity U_extẑ of a uniformly charged spherical particle with radius a and small surface (ζ) potential⁴⁰

in an externally imposed electric field E = −∇ϕ. Far from the particle, E is a constant linear field E = E^∞ẑ. The particle is suspended in an aqueous solution of a monovalent salt, e.g., NaCl. In a weak field, particle velocity is proportional to electric field strength, U_ext = μ_EE^∞, with μ_E called the electrophoretic mobility. The standard expression for μ_E for small ζ is^39,40where f(κa) is Henry's function,³⁹ which accounts for electrostatic screening and depends only on κa, the ratio between particle radius a and the electrostatic screening length κ⁻¹. The function f is plotted in Fig. 3b: it has constant limits of f(∞) = 1, (the Debye or Smoluchowski limit) corresponding to high salt concentration or large particles, and f(0) = 2/3, (the Hückel limit) corresponding to small particles or non-polar solvents. Eqn (26), typically in either the high or low κa limit, is the expression commonly used to compute colloidal ζ potentials from mobility measurements in, e.g., commercial Zetasizers.

In self-electrophoresis, the independent parameters are the surface reaction rates, and therefore the ionic fluxes, rather than the electric field. To facilitate understanding, we translate the expression for external electrophoresis into these terms. We write down expressions for the inverse electrostatic screening length

the ionic conductivityand the concentration-averaged diffusivitytogether with Ohm's lawwhich relates the electric field to the ionic current density at infinity i^∞ẑ and which we can rewrite asCombining eqn (31) with eqn (26) we then haveNote that for electrophoresis in an external field, the particle speed is inversely proportional to the concentration-averaged diffusivity D̅.

We now compare eqn (32) with the analogous expression for the most well-studied self-electrophoretic swimmer, a proton-powered bimetallic swimmer.^14,36 Consider a spherical swimmer of radius a, with surface charge density σ, a surface proton flux j^s₊(θ) = j^s_+,1 cos θ and no bulk reactions: we call this model S^NBR₊. In this case eqn (19) and (20) can be easily solved to yield, after some algebra

where the first and second terms are ϕ^eq and ϕ^sr, the potentials generated, respectively, by the surface charge and the surface reactions. The ⋯ indicate additional, electrostatically screened terms that are necessary to match the electrostatic boundary conditions, but which make no contribution to the propulsion, see Appendix A.5. The propulsion speed is obtained by evaluating eqn (13) with eqn (33) to givewhere F is, like f, a function of κa only. The full form of F is given in Appendix A.3.

Eqn (34) corresponds closely to eqn (32), the particle velocity with external electrophoresis, and we compare these expressions factor by factor:

(I) The relevant current density i^∞ becomes −j^s_+,1e/3 because of the exclusive dependence of the propulsion speed on the first Legendre component of the flux³³ discussed in Section 4.4.

(II) The relevant diffusivity D̅ becomes D₊ because, for self-electrophoresis in steady state, the ionic current can only be carried by the active ion involved in reactions at the particle surface, in this case H⁺. There can be no net flux of the other ions, or they would build up at the particle surface. In fact, the other ions are in local equilibrium with the electrostatic potential ϕ, i.e., etc., from standard Debye–Hückel theory,⁴⁰ and the swimmer behaviour therefore cannot depend on their dynamic properties at all. The appropriate version of Ohm's law for the self-electrophoretic swimmer is thus not eqn (31), but instead⁶¹

which depends only on the mobility of the active ion, and the electrical current. Eqn (33) and the propulsion speed have the same dependencies. This difference between external- and self-electrophoresis has been confirmed in numerical calculations^60,61 and the inverse scaling of self-propulsion speed with the diffusivity of the active ion will be crucial for understanding the effects of bulk ionic reactions in Section 6.

(III) We have chosen to parameterize our model in terms of σ rather than ζ because this is the most natural choice from a microscopic point of view. Much of the charge on the surface, both of conducting and dielectric particles, is due to surface-absorbed groups, leaving σ fixed as other parameters, such as κ, vary. This has been demonstrated experimentally for dielectric particles.⁶² Nevertheless, the experimental evidence indicates that self-electrophoretic propulsion speeds scale with κ⁻² (ref. 11, 16 and 32) which is consistent with a fixed ζ, not a fixed σ, though to our knowledge there is no microscopic justification for this. Since we are most interested in bulk effects, we do not insist upon a particular surface parameterization, and eqn (25) can be used to translate our results into a parameterization where ζ is fixed, which gives at small radius, a speed scaling as a²/κ rather than a³ as in eqn (36) below. Note that this point is distinct from the choice between conducting and dielectric boundary conditions on the particle surface, which is discussed in Section 4.1.

(IV) We have replaced f(κa) with an equivalent expression for self-electrophoresis, F(κa), shown in Fig. 3b. In the thin-screening limit, F(κa → ∞) = 1, and eqn (34) then agrees with previous self-electrophoresis results in the thin-screening limit,³³ except that ref. 33 incorrectly assumes that the propulsion is controlled by the total ionic diffusivity D̄, (see point II above). In the opposite limit, F(κa → 0) = (κa/2)³, so that for small κa the propulsion speed scales with a³

The reason that F(0) → 0, while f(0) is finite, is the different geometry of the driving currents, as illustrated in Fig. 3a. For self-electrophoresis, the driving potential is a local, dipolar field which decays over a length of order a, see eqn (33), whereas in external electrophoresis the driving potential is infinite in extent. Therefore, in self-electrophoresis, additional factors of a in the propulsion speed are to be expected.†

Several of the features of eqn (34) have been verified experimentally for bimetallic swimmers^2,4 like the Au–Pt swimmers in Fig. 1a, which explains the wide-acceptance of the self-electrophoretic model for this system. As discussed above, this equation has also been found to be applicable to single-catalyst swimmers such as Pt-polystyrene Janus particles,^11,16,17 suggesting that these swimmers are also powered by proton currents.^16,17

The additional screening parameter F(κa), is more problematic. It predicts that the speed of a swimmer will drop off sharply as κa decreases, Fig. 3. This drop-off is significant for surprisingly large κa: F(10) ≈ 0.5, while F(1) < 0.1. This shows that the thin-screening limit, which is commonly employed,^{16,17,21,33,63} is not justified even for the common situation of a 1 µm radius swimmer in 3 mol L⁻¹ H₂O₂, where κa ≈ 10 (ref. 16). However, to the best of our knowledge, there is no experimental evidence for this drop off. In fact, a small number of experiments show a larger speed for nanoswimmers^37,64 than is typical for microswimmers.^3,65 We discuss this experimental comparison in more detail in Section 7.

6 Electrophoresis with bulk reactions

We now examine the effect of bulk reactions on the propulsion of model swimmers, in particular of reaction (R4′), H₂O₂ ⇌ H⁺ + HO₂⁻, on the three model swimmers depicted again for convenience in Fig. 7a. In Section 6.1 we write down the general form of expressions for the swimming speed when bulk reactions are included, before focussing on the effect of two experimentally relevant parameters—swimmer radius and H⁺ concentration—on the swimming speed, in Sections 6.2–6.5. In Section 7 we will compare our predictions with experimental observations.

Fig. 7 (a) Recap of the model swimmers S_∘, S₊, and S₌ and the effect of bulk reactions, from Fig. 4. (b and c) Magnitude of the dimensionless bulk reaction factors |B| for type S_∘ (, insets), S₊ (), and S₌ (●) propulsion, from analytical theory with (solid curves) and without (broken curves) bulk reactions; and FEM simulations (symbols). × indicates the base parameter set defined in text. For (b) the particle radius and (c) the proton concentration c^∞₊, at fixed κ.

6.1 General form of the solutions

Mathematically, the bulk reactions make it impractical to solve even the linearized problem by hand. Instead, we solve the system of equations, eqn (19) and (20) symbolically in MATLAB, see Appendix A. The final solution is very similar to the reactionless solution, with an extra bulk reaction factor B(qa,…), so for each of our model swimmers, we can writein the form of eqn (2). Here, † indicates a particular swimmer type, i.e., † ∈ {∘, +, =} and j^s_† is the appropriate surface flux density for that swimmer. We define j^s_† = j^s_∘,1 for the S_∘ swimmers, and j^s_† = j^s_+,1 for the S₊ and S₌ swimmers. The use of D₊ in the denominator of eqn (37) is an arbitrary definition. Under this definition, the bulk reaction factors in the absence of bulk reactions have the constant values B^NBR_∘ = 0, B^NBR₊ = −1 and B^NBR₌ = (d₊ − d₋)/(d₊d₋) = 9.3. Here, d_l is a rescaled diffusivity, d_l = D_l/D₊. By definition, d₊ = 1, but we retain d₊ for symmetry of notation.

Note that the expression in square brackets in eqn (37) is identical to eqn (34). This emphasizes that all the propulsion mechanisms, S_∘, S₊, and S₌, are really forms of self-electrophoresis, and display all the responses to, e.g., salt-concentration, particle radius, and surface charge, which standard self-electrophoresis models predict. The inclusion of bulk reactions just adds a new layer of phenomena on top of this behaviour.

6.2 Influence of bulk reactions on swimming

We study the bulk reactions by varying two common experimental parameters: particle radius, and proton concentration, i.e., pH, Fig. 7b and c, respectively. Including bulk reactions (solid curves = analytic; solid symbols = FEM numerics) introduces a range of effects compared to the case with no bulk reactions (broken, horizontal lines), with qualitatively different behaviour for the three swimmer models.

Examining Fig. 7b first, the bulk reactions permit propulsion of the neutral swimmer S_∘ (inset), and B_∘ increases with radius, saturating for large radii. However, the magnitude of B_∘ always remains smaller than that of the other swimmers by a factor of order 10⁻⁶. In practice, this is typically partially compensated for by the much larger flux of the neutral species. For the proton-current-driven swimmer S₊, B₊ shows plateaux at both large and small radius, with the large-radius plateau approximately twice as high. For the ionic-diffusiophoretic swimmer S₌, B₌ scales inversely with radius for large radius, but has a plateau at small radius.

Meanwhile, varying the proton concentration c^∞₊, as in Fig. 7c, produces a peak in B_∘ and B₌, and decreases the overall value of B₌ by at least a factor of 5 compared to without reactions. For S₊ there are again two plateaux, at high and low c^∞₊, with the low c^∞₊ plateau now a factor of approximately 10 higher than the other. The main control parameter for all these effects is qa, and there is a qualitative change of behaviour for all three swimmers at qa ≈ 1: the vertical lines on Fig. 7b are for qa = 1. In Table 1, we write down the bulk parameters for each of the model swimmers in the limits qa ≪ 1 and qa ≫ 1. The full analytical expressions, which are lengthy, are provided in Appendix A.4, but the basic physics can be understood from the limiting behaviour. For the table, we have also assumed weak ionic dissociation, i.e., K_eq ≪ c^∞_∘ which is valid here, and thin electrostatic screening, κa ≫ 1. These assumptions also apply to the analytical expressions given in the rest of this section. Table 1 matches Fig. 7 in all but one respect, which is the scaling of B_∘ at qa ≪ 1, and this difference occurs because the assumption κa ≫ 1 does not hold for small a in Fig. 7b. The parameters α and d* will be defined below. For Fig. 7c, qa ≈ 7 or larger, so we will assume that this figure is always in the qa ≫ 1 limit.

Table 1 The bulk mobility factors predicted in the thin-screening κa ≫ 1, κa ≫ qa, and low dissociation c^∞₊,c^∞₋ ≪ c^∞_∘ limits, for low, qa ≪ 1 and high qa ≫ 1 reaction rates. In both limits the prefactor should be multiplied by the relevant expression in the right-hand columns. The full expressions are given in Appendix A.4

	Prefactor ×	qa ≪ 1	qa ≫ 1
B ∘			1
B +	1
B =		1

---

To understand the results shown in Fig. 7 and Table 1, we will examine the bulk reactions in terms of three physical principles: reactive screening, the composition of the electrical current, and the dissociation of the neutral flux. Though we focus on these underlying principles, which are crucial for understanding the effect of bulk chemical reactions on any swimmer, this structure also allows us to discuss the three model swimmers in a logical order: S₌, S₊ and S_∘.

6.3 Reactive screening (model S₌)

If an ion is released from the particle surface, it will react and come into local equilibrium with the surrounding solution. The characteristic distance over which this approach to equilibrium occurs can be called a ‘reactive screening length’ q⁻¹. As for the simple model discussed in Section 2.3, the reactive screening length is a balance between molecular diffusion and the reaction rate. However, the expression for q is more complex than in the simple model. We findMathematically, q corresponds to one of the eigenvalues of the linear system of equations, eqn (19) and (20), see Appendix A.2. For our base parameter set, we obtain a screening length q⁻¹ = 74 nm.

Just as for the simple model, reactive screening gives an exponential decay of chemical concentrations with distance from the particle surface. The inclusion of charged species means that the electrostatic potential can now also become screened. However, we observe this reactive electrostatic screening only for S₌ swimmers, where the two ions released from the surface both react with oppositely charged ions in the bulk solution, causing an exponential decay in the resulting diffusion potential, see Fig. 8a.

Fig. 8 Normalized surface-reaction-generated potential ϕ^sr for (a) S₌ and (c) S₊ swimmers, with (right) and without (left) bulk reactions. In each case, the potentials are normalized by the largest value of |ϕ^sr| without bulk reactions. (b) Normalized radial decay of ϕ^sr along X–X′ for (b) S₌ and (d) S₊. Solid curves are with bulk reactions, dashed curves without.

For S₊ swimmers, no such screening is observed, see Fig. 8b. This is because the electrical current is conserved, so cannot be screened, and hence the associated electrical field also retains its unscreened dipolar form. Similarly, the H₂O₂ concentration field c_∘ remains unscreened because this field is approximately conserved in the weak-dissociation limit. This results in an unscreened electrostatic potential field for S_∘ swimmers (not shown).

For S₌ swimmers, reactive screening also explains the 1/(qa) scaling of B₌ at high qa, as we show with a simple scaling argument: from eqn (4) and (5), we expect a fixed ratio between the surface reaction rates and the concentration gradients normal to the surface. For example, at r = a

independent of other parameters. For qa ≫ 1, the concentration decays exponentially away from the surface

c₊ ∝ exp(−qr). (40)

Differentiating this equation with respect to r gives ∂c₊/∂r ≈ −qc₊ and comparing this with eqn (39) yieldsSince the diffusion potential is proportional to the ionic concentrations (c₊, c₋), and the propulsion speed U₌ is proportional to the diffusion potential, we have U₌ ∝ q⁻¹. On the other hand, without bulk reactions the only relevant length scale is a, so a similar argument gives U^NBR₌ ∝ a. Therefore, one obtains B₌ ∝ U₌/U^NBR₌ = 1/(qa). Physically, for qa ≫ 1, the concentration flux only has the small screening length q⁻¹ over which to set up a diffusion potential, whereas without bulk reactions a length of order a is available.

This 1/(qa) scaling immediately explains the B₌ ∝ 1/a scaling in Fig. 7b. We can also understand the peak in B₌ in Fig. 7c by noting that the screening length q⁻¹ vanishes both for high c^∞₊ and for high c^∞₋ in eqn (38). Since c^∞₋ scales inversely with c^∞₊ due to the ionic equilibrium of eqn (22), this means that q⁻¹ vanishes at either end of the c^∞₊ scale. As B₌ ∝ q⁻¹, it too vanishes at either extreme and is peaked for intermediate c^∞₊. Physically, at either end of the c^∞₊ scale, the high concentration of ions screens electric fields, preventing the formation of a diffusion potential.

6.4 Composition of the electrical current (model S₊)

The total electrical current in the bulk is a conserved quantity, and is therefore not screened. However, the individual ionic fluxes making up that current are not conserved, and the bulk reactions modify the identity of the current-carrying ions. As discussed in Section 5, this is important because the swimming speed scales inversely with the diffusivity of the current-carrying ion or ions. In our system, an initially pure proton current will be partially replaced by HO₂⁻ ions travelling in the opposite direction, as illustrated in Fig. 5. In the reactive, qa ≫ 1 limit we can calculate the composition of this electrical current relatively simply. From this, we will obtain the propulsion speed of the S₊ swimmer.

In the qa ≫ 1 limit, at any point outside the thin reactive screening layer, the ions released from the surface will have had time to come into equilibrium with each other. This is equivalent to requiring that the chemical production rates vanish, i.e., R_l = 0. In the linear approximation, see eqn (17), this means

In other words, the deviations in concentration x_l of each of the reactive chemical species are coupled by the reaction matrix given in eqn (23). This concentration coupling also implies a coupling of the chemical fluxes, in the same way that charge conservation implies the conservation of electrical current. Consider the linearized flux equation, eqn (16). Multiplying both sides by k_ml/(D_lc^∞_l) and summing over l yieldsThen, from eqn (42), with l and m exchanged, the first term on the right vanishes, while charge conservation in reactions implies (see eqn (7)), so the second term vanishes too. Hence, the general flux coupling equation isFor our specific system, substituting the expression for from eqn (23) into eqn (44) then givesThe physical meaning of eqn (45) is that each of the molecular fluxes has a characteristic scale set by D_lc^∞_l, and that, with this scaling, the relationship between the currents is set by the stoichiometry of the bulk reactions. Eqn (45) can be rearranged to give each of the ionic fluxes j_± in terms of the conserved quantities i and j_∘where D* is the concentration-averaged diffusivity of the active ions‡and the dimensionless factor α, which specifies the equilibrium decomposition of a neutral current into ionic currents, isThe meaning of the first term in eqn (46), which is relevant for S₊ swimmers, is that the electric current is carried by a fixed proportion of H⁺ ions travelling in one (positive) direction, and a counter current of HO₂⁻ ions in the opposite (negative) direction. In the second term, which is relevant for S_∘ swimmers, the neutral flux j_∘ continuously dissociates into H⁺ and HO₂⁻ ions, producing small, equal fluxes of these ions, which travel with the neutral flux.

If we are also outside the electrostatic screening length, which is the case in the κa ≫ 1 limit, we also have a zero-charge-density condition, which reads

Just as above, but now multiplying eqn (16) by z_l/D_l, we can derive a direct relationship between the electric field and the chemical fluxes⁶¹which, combined with eqn (46), yields a version of Ohm's law for the reactive limitComparing the first term of this equation with eqn (35) for self-electrophoresis without bulk reactions, we see that they are identical apart from the switch from D₊ to D*.

For an S₊ swimmer, j_∘ = 0, and for high qa all of the electric field outside a thin screening layer will be determined by the first term of eqn (51). From this we can understand the 1/D* factor which appears in B₊, see Table 1. Just as without bulk reactions (see eqn (34)) the propulsion speed is inversely proportional to the diffusivity of the current-carrying ion. However, the current is now made up of two ions, with a total effective diffusivity D*. This explains the 2× speed increase in Fig. 7b: D₊ at low a is replaced by D* at high a, and for c^∞₊ = c^∞₋, which is the case in Fig. 7b, D* = (D₊ + D₋)/2 ≈ D₊/2.

To understand the effect of varying c^∞₊, Fig. 7c, we examine the form of D* in eqn (47). At high c^∞₊, D* ≈ D₊, while at low c^∞₊ (=high c^∞₋), D* ≈ D₋. Physically, this is again simply a result of the relative number of each ionic species: if there is an overwhelming number of protons in solution, then the ionic current must be carried predominately by protons. This explains the factor of D₊/D₋ ≈ 10 speed difference between the two plateaux for B₊ in Fig. 7c.

6.5 Dissociation of the neutral flux (model S_∘)

To understand the dissociation of a purely neutral flux, we examine the parameter α in eqn (48). The form of α can be explained by the fact that in the absence of a net electrical current, e.g., for S_∘ swimmers, the ionic currents are constrained by j₊ = j₋. This means that the total ionic flux will be limited by whichever ion has the lower value of D_lc^∞_l, as this ion will contribute most to the flux balance in eqn (45). Hence the parameter α, like q⁻¹, vanishes at the extreme ends of the c^∞₊ scale: at low c^∞₊ it is limited by the low proton concentration, and at high c^∞₊, by the low HO₂⁻ concentration.

The dissociation of the neutral flux generates a diffusion potential. Hence, the prefactor for B_∘ in Table 1 is made up of two factors: α and (d₊ − d₋)/(d₊d₋), the latter of which controls the diffusion potential just as for B₌. The peak in B_∘ as a function of c^∞₊ then follows directly from the behaviour of α.

Interestingly, both S_∘ and S₌ show peaks in speed at intermediate c^∞₊, but for two different reasons. For S₌, the reason is that the reassociation of ions is slowest at intermediate concentrations. For S_∘, the reason is that the least conductive fraction of the solution limits the total carrying capacity, and this effect is strongest at either extreme in pH.

7 Comparison with experiments

We now compare our theoretical predictions with experimental results, in so far as this is possible at present. We stress here again the lack of understanding of the surface chemistry, and in particular of the effect of experimental parameters on the ionic surface reaction rates. We have not attempted to predict these reaction rates, so we cannot immediately test our theoretical predictions. In this section alone, we will make the simple and typical^3,15,16 assumption that the surface properties, i.e., surface reaction rates and surface charge densities, vary only with fuel concentration and are otherwise constant, so that U ∝ κ⁻³F(κa)B(qa) for all swimmers (eqn (37)). This allows us to make some suggestive comparisons with experiments. We note that more complex reaction rate dependencies based on electrochemical modelling of the surface have been proposed previously.^14,17,36 However, to the best of our knowledge, these models do not have independent experimental justification, or experimental validation from speed experiments beyond that achieved by the assumption of uniform surface properties.

Independent of bulk reactions and swimmer type, we predict a speed scaling with a³ for small particles. In particular, for a proton powered swimmer U^NBR₊, the predicted speed dependence is as shown in Fig. 9a (solid curve), due mostly to the new electrostatic screening parameter F(κa)—the bulk reactions do not significantly modify the form of this curve. Here, we have matched the solution parameters to those of the typical microswimmer experiments of ref. 65, choosing the (constant) surface parameters to give a speed of 7 µm s⁻¹ (blue circle) for a 1 µm radius bimetallic sphere, as in that paper (see Appendix B for the surface parameters used). We also plot nanoswimmer data (black triangle) from ref. 37, where a swimmer of radius 10 nm had U = 650 µm s⁻¹. These two experiments used similar concentrations of H₂O₂, but differed in salt concentration. Thus, to predict the nanoswimmer results, we also plot a theoretical curve (dashed), with the salt concentration modified to match those of the nanoswimmer experiment,³⁷ but keeping the surface parameters the same as in the microswimmer experiments.⁶⁵

Fig. 9 Comparison between theory and experiments. (a) The predicted speed of a swimmer powered by a proton current, in the presence of chemical reactions, with parameters chosen to match typical measurements on microparticles (red solid, theory; blue circle, experiment⁶⁵); and, with the same surface parameters, but bulk solution parameters chosen to match experiments on a nanoswimmer (red dashed, theory; black triangle, experiment³⁷). (b) The speed of an S₌ swimmer plotted against experimental data on Janus-Pt microswimmers.¹⁵ The experimental error bars are smaller than the data points.

If the assumption of constant surface properties holds true, then this dashed curve should agree with the experimental value for nanoswimmers. Instead, there is a clear disparity amounting to several orders of magnitude. This is not the result of our linear approximation: we find a good match between analytics and numerics up to values of σ and j^s₊ higher than those used in plotting Fig. 9a (Appendix C). The disparity is also not significantly reduced if we assume uniform ζ rather than uniform σ, see Section 5. The discrepancy could be explained in at least two ways. It may be that self-electrophoresis is not the correct propulsion mechanism for bimetallic nano-swimmers. It has recently been found that nanometre scale biological enzymes also exhibit self-propulsion,^41,42 and a range of mechanisms has been proposed for this propulsion,⁶⁶ some of which might also apply to bimetallic nano-swimmers. Alternatively, it may be that the assumption of constant surface properties is inappropriate. That is, the proton current density could be much higher for these nano-swimmers than for micro-swimmers. Whatever the explanation, this discrepancy highlights the need for more systematic studies of identical or comparable swimmers over wide parameter ranges, as in ref. 15–17, and for independent measurements of the relevant ionic reaction rates.

For Pt-polystyrene Janus particles, such systematic studies do exist.¹⁵ These show a U ∝ a⁻¹ scaling for 0.2 µm < a < 5 µm (this scaling has also been observed over a narrower range for some bimetallic swimmers⁶⁵). Self-electrophoresis S₊, reaction (R2), is currently the preferred mechanism for Pt-polystyrene Janus swimmers,^16,17 but comparison of this 1/a scaling with Fig. 7b suggests self-ionic diffusiophoresis S₌ as an alternative mechanism, corresponding to reaction (R3). This is plausible: ion release without net electrical currents, which would correspond to reaction (R3), has previously been observed for H₂O₂ decomposition on Pt.⁶⁷ This mechanism would also avoid the conceptual difficulty of producing a net ionic current in single-catalyst systems.^16,17 However, when we plot the experimental data from ref. 15 against our theoretical predictions for S₌ propulsion, Fig. 9b—which is again scaled to match the experimental data for 1 µm radius swimmers, see Appendix B—we see that the fit fails at small a, again due to the F(κa) parameter. It is possible that evaluation of the complete H₂O₂–H₂O reaction system would provide a better fit, but this goes beyond the scope of this work.

Note that the 1/a scaling has previously been explained by postulating that the overall surface reaction rate j^s_∘ is limited by diffusion,¹⁵ and therefore scales as 1/a just from geometrical arguments. However, the diffusion-limit implies a large flux density j^s_∘ ≈ D_∘c^∞_∘/a, which for a 1 µm radius swimmer in 3 M H₂O₂, as in ref. 15 requires j^s_∘ ≈ 3 × 10²⁴ m⁻² s⁻¹. So far, only much smaller rates, j^s_∘ ≈ 10²² m⁻² s⁻¹, have been measured, both by us¹⁶ and by the authors¹⁷ of ref. 15. Therefore, these swimmers do not appear to be in the diffusion limited regime, so this explanation for the 1/a scaling cannot hold.

Next, we have previously calculated the values of the uncharged flux j^s_∘ and the charge density σ for Pt-coated Janus swimmers.¹⁶ We estimated that the propulsion speed of such swimmers was too high to be explained by a purely uncharged reaction like (R1).¹⁶ This estimate did not allow for bulk ionic reactions. However, including these reactions, we calculate in Appendix C that such a mechanism could still only account for ≈5% of the observed speed of these swimmers. Hence, a model with just surface reaction (R1), even with bulk dissociation, cannot explain the propulsion of H₂O₂-powered swimmers, so that such swimmers probably still require more complex ionic surface reaction schemes like (R2) and (R3). Nevertheless, purely neutral-surface-reaction mechanisms could still be relevant for swimmers powered by more dissociative fuels, such as hydrazine.⁵⁹

Turning to the effect of pH, there have been two suggestive studies,^16,59 but no systematic investigation. First, we found that NaOH reduced the swimming speed of Pt-polystyrene Janus swimmers, but that this effect was much weaker than the speed reduction due to NaCl.¹⁶ This is consistent with our prediction that increasing pH at fixed Debye length should raise the swimming speed for any of the 3 swimming mechanisms discussed. Raising the pH corresponds to moving left from the × symbol in Fig. 7c; for all swimming mechanisms the value of B increases in this direction.

Second, the silica–iridium swimmers of ref. 59 show a clear spike in speed as a function of fuel (hydrazine) concentration similar in form to the peaks in Fig. 7c. This spike could be due to modulation of pH by the reaction product ammonia, which would imply that either neutral self-diffusiophoresis or self-ionic diffusiophoresis dominates this swimmer's self-propulsion. In both these experiments, however, variation of the reaction rates with pH could also explain the results,⁶⁸ so further systematic study is necessary.

Finally, we note that the bulk ionic dissociation rates are also not well known, see Section 4.4. Our estimate of the association constant k_as is a diffusion limited rate, and is therefore an upper limit, though for a proton exchange reaction, this is likely an accurate estimate.⁵² Nevertheless, it might be argued that the effects of bulk reactions will not be observed in practice for lower reaction rates, as we require qa ⪆ 1. However, since q⁻¹ ∝ k_as^−1/2, see eqn (38), for a much lower value of k_as, say 1000 times lower, we still find q⁻¹ ∼ 2 µm, well within the experimental microswimmer size range. Further, the value of q⁻¹ becomes arbitrarily small in both high and low pH solutions, see eqn (38).

8 Conclusion

In this article, we have theoretically explored the influence of common polar solvents such as water and hydrogen peroxide on the propulsion behaviour of chemically-propelled, synthetic microswimmers. We have focussed on two unavoidable properties of such polar solvents—electrostatic screening and ionic dissociation—and calculated their effect on the swimming speed of a wide range of microswimmers propelled by chemical reactions on their surfaces. These effects have not been studied systematically before; nevertheless, they are highly significant, and including these effects can modify predicted swimming speeds by several orders of magnitude.

By ionic dissociation, we mean the breaking up of neutral molecules, including water, into charged species, for example, H₂O₂ ⇌ H⁺ + HO₂⁻. One of our main predictions is that this kind of ionic dissociation reaction allows even microswimmers whose surface chemistry does not involve any ions, e.g., swimmers propelled by the simple decomposition of hydrogen peroxide, 2H₂O₂ → 2H₂O + O₂, to generate ionic gradients and thereby electric fields. The implication of this is that all microswimmers in water should experience some degree of self-electrophoresis, i.e., propulsion via self-generated electric fields. This is significant because self-electrophoresis is much more efficient than other putative propulsion mechanisms and is likely to dominate over them. Put simply: our results imply that all swimmers in aqueous solution are likely to be self-electrophoretic to a large degree.

The second major prediction of our work is that for some chemically-propelled swimmers, ionic dissociation reactions may result in a kind of exponential ‘reactive-screening’.²⁰ Electrical and chemical concentration fields generated by surface reactions on microswimmers are usually taken to decay slowly into the bulk solution, that is, as a power-law with distance. Ionic dissociation can instead produce a short-ranged exponential decay of these fields, just as in electrostatic screening. This is significant since these chemical and electrical fields are implicated in inter-swimmer interactions and collective behaviour, and the interaction range will play a crucial role in this behaviour.

Our third prediction relates to electrostatic screening itself. Most analytical treatment of microswimmers has focussed on the thin-screening limit, where the electrostatic screening length is much smaller than the swimmer size. For very small swimmers, this limit does not apply, and we find that this massively reduces the predicted swimming speed. This is important because experiments on nanoscale swimmers³⁷ show that these in fact swim faster than microswimmers, in apparent contradiction to our predictions. This opens up the exciting possibility that nanoscale swimmers move by entirely novel mechanisms compared to their microscopic counterparts.

Finally, the general conclusion that we draw from our results is that more experimental work is required to understand self-propulsion mechanisms. The effect of ionic dissociation in particular depends crucially on the type of surface reactions which are responsible for propulsion—and the details of these reactions remain almost universally unknown. What is most urgently required in this regard is the independent measurement of surface reaction rates, which is challenging, and has so far only been achieved in the simplest of cases. However, recent results with electroosmotic pumps⁶⁹ suggest that such measurements will not long remain beyond our reach. We particularly hope that our theoretical results will lead to renewed efforts in this direction.

Looking ahead, our results suggest that a deeper understanding of self-propulsion will lead to greater insights into swimmer–swimmer interactions and collective effects. This is particularly relevant to synthetic swimmers, as their propulsion is closely coupled to their interactions through self-generated electrostatic, chemical, and hydrodynamic flow fields. We have shown here that reactive screening can qualitatively change the electrostatic interactions between swimmers. A detailed follow-up study will look explicitly at such interactions. Further theoretical work will focus on applying our calculations to fully realized experimental systems, e.g., mixed metal-dielectric swimmers.

Appendix

A Calculation of the analytical solution

In this appendix, we explicitly calculate the propulsion speed for a general swimmer. This calculation is based on the linearized model described in Section 4, but now with a more concise notation, A.1. From this linearized model, we determine first the electrostatic potential fields, A.2, then the propulsion speed, A.3, as set out in Section 4. In A.4, we apply this general calculation to determine the speed of the model swimmers presented in the main text. Finally, in A.5, we demonstrate the equivalence of uniform charge and uniform potential boundary conditions for the calculation of the propulsion speed.

A.1 The linearized model

We begin with the linearized model described in Section 4.3. For notational convenience, we define a composite dimensionless parameter y_l by combining the linearized potential ψ and concentration x_l fieldsWith this notation, the linear system of equations, eqn (19) and (20) is given byEqn (53) represents a system of N + 1 linear equations. However, several of the species, typically inactive ions such as Na⁺ or Cl⁻, may not be involved in any bulk or surface reactions, and we will now show that these inactive species can be eliminated. We specify that the first N′ indices (N′ < N) correspond to the reactive species. For the remaining, unreactive species, all the bulk reaction coefficients, k_lm are zero, and there is no surface flux, so eqn (53) can only be satisfied if

y_l = −z_ly₀, l > N′. (54)

This is the linear approximation to the Boltzmann distribution, which one expects, since these unreactive species should be in equilibrium with the electric field. Using eqn (54), these ions can be eliminated from the remaining N′ + 1 parts of eqn (53) to yieldwhere κ is the inverse Debye screening lengthwith the Bjerrum length, l_B = e²/(4πεk_BT), and where χ_l is a dimensionless ionic concentration

χ_l = 4πl_Bκ⁻²z_l²c^∞_l. (57)

Eliminating the inactive ions makes it clear that the motion of the swimmer cannot depend on the diffusivity of these ions, and is only affected by them through the value of κ and through charge balance.

Finally, linearizing the boundary conditions in eqn (3) and (4) gives

A.2 Calculation of the electrostatic potential

Eqn (55) has the form of a matrix equation with components corresponding to the chemical concentrations and the electrostatic potential, so it is convenient to introduce some additional matrix notation. The bold font is reserved for real-space vectors, such as the fluid velocity u, while vectors in this concentration-potential space will be underlined. A general vector t̲ will have N′ + 1 components labelled t_l, while a matrix will have (N′ + 1) × (N′ + 1) components labelled T_lp. A point in the concentration-potential space is specified by the vector y̲, with components y_l, as defined in eqn (52). Using this notation, we can rewrite eqn (55) aswhich can be solved by finding the N′ + 1 eigenvectors of the matrix , with eigenvalues μ_p. These eigenvectors define a new basis, in which is diagonal. Defining w̲ as the representation of y̲ in this basis, we havewhere the matrix is diagonal, with components G_lp = δ_lpg_p, where δ_lp is the Kronecker delta and is an inverse screening length. For the model system described in the main text, the unique values of g_p are κ, 0 and q, as given in eqn (38). Here, for clarity, we use the index p to refer to the screening lengths, and the indices l or m to refer to the concentrations and potentials, even where these are dummy indices.

Eqn (60) is a series of N′ + 1 independent Helmholtz equations, and the full solution to this equation is just a vector of individual solutions to the Helmholtz equation. In spherical polar coordinates, these solutions have the form^20,70

with w_p,n an as yet undetermined surface coefficient, P_n the Legendre polynomial of order n,⁷¹θ the polar angle, andWe refer to the Legendre components by the subscript n throughout, and where we have multiple subscripts, the Legendre subscript shall be preceded by a comma. Transforming back into the original coordinate frame linearly combines the solutions in eqn (61), so that the final form for the electrostatic potential iswith analogous expressions for each concentration field. Here, ϕ_p,n are surface coefficients which we will now determine.

Transformation back into the original coordinate system is achieved with a transformation matrix

where each element K_lp of is equal to the lth component (in the original coordinate system) of the pth eigenvector. Applying this transformation to eqn (61) givesThe boundary conditions specified in eqn (58) can also be rearranged into a matrix equationwhere b̲ is a vector specifying each of the boundary fluxes or charge density. We define the harmonic components b̲_n of b̲ by , with analogous expressions defining . The solution to the boundary conditions is found by inverting eqn (66) to yieldwherein which the diagonal matrix has elementsInserting the boundary conditions into eqn (65) then giveswhere the surface coefficients areIn particular, this yields for the surface coefficients of the electrostatic potential (for which the index l = 0)Eqn (72), together with eqn (63) completely determines the electrostatic potential field.

We can also determine the equilibrium and non-equilibrium components of the potential. Writing ϕ = ϕ^eq + ϕ^sr, where ϕ^eq = ϕ({j^s_l} → 0) is the equilibrium potential distribution without any surface chemical reactions (here, {j^s_l} is the complete set of surface fluxes), and ϕ^sr = ϕ(σ → 0) is the additional potential generated by the surface reactions, we have

with surface coefficientsHere b̲^eq = b̲({j^s_l} → 0) is the vector specifying the boundary conditions for a charged but unreactive particle, and b̲^sr = b̲(σ → 0) specifies the boundary conditions for an uncharged but reactive particle.

A.3 Calculation of the propulsion speed

Having determined the electrostatic potential, we calculate the fluid flow by making use of the Lorentz reciprocal theorem.³⁸ This allows one to transform the Stokes equation, eqn (11), from a 3D partial differential equation into an integral equation on the 2D domain boundary (the swimmer surface). Using this approach, a general formula for the propulsion velocity U of a non-slip sphere generated by an axisymmetric distribution of force density f has been derived⁴⁸where the volume integral is over the region outside the sphere, and the scalar speed U used in the main text is defined by U = Uẑ. Here f = −ρ_e∇ϕ from eqn (10).

For a uniformly charged sphere, the equilibrium potential distribution is

Making the usual assumption of a small driving field, i.e., ϕ^sr ≪ ϕ^eq then giveswhere the ϕ^sr_p,1 are to be read out from eqn (76) andwhich is the self-electrophoretic equivalent of the Henry function for electrophoresis in an external field.³⁹ We have verified that eqn (79) is also obtained by solving the 3D Stokes equations directly, following Henry's methods.^39,40,70 We also write down a single-argument form of the self-electrophoretic function F(x) = F(x,0), which is useful when considering swimmers without bulk reactions, for which y = qa = 0This is the function discussed in Section 6.

A.4 Propulsion speed for the model swimmers

We now write down the propulsion speed for the 3 model swimmers discussed in the main text. These expressions were determined by solving eqn (79) symbolically in MATLAB, and making the further assumption of weak ionic dissociation. We find that the bulk reaction factors, as defined in eqn (37) arewhere Θ(κa,qa) depends on the relationship between the 3 length scales a, κ⁻¹ and q⁻¹ and isWith the limits κ ≫ q, and either qa ≫ 1 or qa ≪ 1, we obtain the expressions given in Table 1.

A.5 A note on electrostatic boundary conditions

In this section, we show that a particle with fixed, uniform surface charge σ has the same propulsion velocity as an equivalent particle with fixed, uniform surface potential ζ, as long asTo do this, we first need to show that modifying the electrostatic boundary conditions of the particle has only a limited effect on the fields of concentration and potential; namely, that such modifications can only introduce electrostatic fields corresponding to the equilibrium Debye–Hückel solutions around passive colloids.

We take a swimmer, in a given chemical environment, and apply to it three sets of boundary conditions. Boundary conditions (1) and (2), with corresponding solutions y̲⁽¹⁾ and y̲⁽²⁾, have equal chemical flux boundary conditions (equal surface reaction rates), but have arbitrary, different electrostatic boundary conditions. Boundary condition (3) consists of a no flux condition on all species (no surface reactions), and the electrostatic boundary condition

y⁽³⁾₀(s) = y⁽²⁾₀(s) − y⁽¹⁾₀(s). (85)

Since there are no fluxes through this particle's surface, each chemical species is in equilibrium, and the solution to this boundary condition is just the equilibrium, Debye–Hückel solutionwhere the equilibrium potential field ψ must satisfy both the electrostatic boundary condition, eqn (85), and the Debye–Hückel equation

∇²ψ = κ²ψ. (87)

One can then show by direct substitution of eqn (86) into eqn (55), that the solutions to the three boundary problems are related by y̲⁽²⁾ − y̲⁽¹⁾ = y̲⁽³⁾. In particular, ϕ⁽²⁾ − ϕ⁽¹⁾ = ψ, which implies, from eqn (87)

∇²[ϕ⁽²⁾ − ϕ⁽¹⁾] = κ²[ϕ⁽²⁾ − ϕ⁽¹⁾]. (88)

Hence the difference ϕ⁽²⁾ − ϕ⁽¹⁾ between the electric potential fields of particles (1) and (2) corresponds to an equilibrium Debye–Hückel solution around a passive colloid.

As before, we make the assumption of a small driving field, ϕ^sr ≪ ϕ^eq, where ϕ = ϕ^eq + ϕ^sr. Now, consider two particles (1′) and (2′), with equal surface reactions, but where (1′) has uniform surface charge density σ, and (2′) has uniform surface potential ζ, with σ and ζ satisfying eqn (84). In this case, the two equilibrium fields are equal, i.e., ϕ^(1′)eq = ϕ^(2′)eq, and are given by eqn (78). Subtracting this equality from eqn (88) yields, for the remaining, non-equilibrium part of the potential

∇²[ϕ^(2′)sr − ϕ^(1′)sr] = κ²[ϕ^(2′)sr − ϕ^(1′)sr]. (89)

In other words, the difference in the reaction-generated electrostatic potential field between (1′) and (2′) is an equilibrium, Debye–Hückel type field, which has an inverse screening length g_p = κ. Since there is a (κ − g_p) factor in eqn (79), we see that such a field can have no effect on the propulsion speed. This proves the assertion that, to linear order, a particle with fixed, uniform surface charge σ will have the same propulsion velocity as an equivalent particle with fixed, uniform surface potential ζ, as long as eqn (84) is satisfied.

In fact, one can make a more general statement, which we will not prove. For any two particles (1′) and (2′), with equal arbitrary shape, surface reactions, and equilibrium (possibly non-uniform) fields ϕ^eq, not only the propulsion speed but the entire flow field will be the same. A physical justification for this conclusion is that if the interaction between one equilibrium field (ϕ^eq) and another (the difference field between (1′) and (2′)) could generate fluid flow, then this would constitute a perpetual-motion machine. Analogous conclusions have also been drawn for electrophoresis in an external field.⁷²

B Experimental parameters

B.1 The ionic association constant

An important parameter in our calculations is the ionic reaction association constant k_as in (R4′), H⁺ + HO₂⁻ ⇌ H₂O₂ in water, see Section 4.4. We were unable to find a value for this constant in the literature. However, reactions involving the transfer of a proton or a hydroxyl ion are normally sufficiently fast to be diffusion limited.⁵² It has been shown,⁷³ that the diffusion-limited rate constant between two species, A and B, with diffusivities D_A, D_B, and valences z_A, z_B, which react at a short distance r_AB is⁷⁴

k_as = [4π(D_A + D_B)r_AB]f(z_Az_B,r_AB). (90)

Here, f(z_Az_B,r_AB) is a modifier for charged speciesFor reactions between oppositely charged species, over a typical reaction distance in water of r_AB = 0.2 nm,⁷⁴f(−1,r_AB) = 3.59. For the reaction between H⁺ (species A) and HO₂⁻ (species B), this yields, using the diffusivities quoted in the main text, k_as = 4.9 × 10¹⁰ mol⁻¹ L s⁻¹.

B.2 Comparison with experiments

For the comparison of self-electrophoretic micro- and nano- swimmers in Fig. 9, we take experimental parameters from ref. 65 (microswimmers) and ref. 37 (nanoswimmers). For microswimmers, we use a = 1 µm, c^∞_∘ = 1.5 mol L⁻¹ and no added salt. For nanoswimmers, we use a = 15 nm, c^∞_∘ = 1.5 mol L⁻¹, pH 7 and 5 mmol L⁻¹ NaCl (this has the same Debye length as 1 mmol L⁻¹ trisodium citrate, which was used in practice³⁷). For both swimmers, we take σ = 10⁻²e nm⁻² (ref. 21) and j^s_+,1 = 1.66 × 10⁻⁷ mol m⁻² s⁻¹, which is chosen to match the microswimmer speed in ref. 65.

For comparison between S₌ and polystyrene-Pt Janus particles, we take experimental parameters from ref. 15, which are c^∞_∘ = 3 mol L⁻¹ with no added salt. We used σ = 10⁻²e nm⁻² again and j^s_+,1 = 6.42 × 10⁻⁶ mol m⁻² s⁻¹, which is chosen to match the experimental propulsion speed at a = 1 µm s⁻¹.

C Finite element method calculations

In this section we give additional details for the numerical FEM calculations discussed in Section 4.2. FEM calculations are performed using the COMSOL 5.1 Multiphysics Modelling package.

We employed the following strategies to accelerate the calculations and obtain high quality results. (i) The solutions were obtained in a 2D cylindrically symmetric geometry. (ii) We ignored the advective coupling term in eqn (5). This allowed us to split the problem into electrostatic plus hydrodynamic parts, as for the linear theory, and thus solve the uncoupled equations more efficiently. This approach is justified, since the Péclet number (Pe) ≲ 10⁻² for typical experimental systems. We also verified this directly, by including the advective coupling term in a subset of the data points, finding good agreement. (iii) We created a physics-specific mesh, see Fig. 10, on which we solved the system. Quadrilateral elements were used out to a distance of 3κ⁻¹ from the colloid surface. These elements grow exponentially in size with increasing distance, whilst maintaining a constant number along the tangential direction. The remainder of the domain was meshed with triangular elements which grow larger with distance from the colloid. This approach is necessary to ensure convergence of the model. (iv) The following polynomial orders were employed for the test functions: electrostatics (3), diffusion (5) and hydrodynamics (2 + 3). These higher orders proved necessary to reduce spurious flow (see also ref. 44). (v) Finite-size scaling was employed to check for artefacts arising from the finite extent of the simulation domain, we found that for L = 10a + 25κ⁻¹ the effects on the speed of the particle were negligible. (vi) Mesh refinement was used for several simulations to determine the dependence of our result on the element size. (vii) We also varied the tolerance on the residual for a few cases to verify that our solutions had sufficiently converged.

Fig. 10 The mesh on which the FEM calculations are performed. This particular mesh was generated for radius a = 0.5 µm and a salt concentration of 10⁻⁵ mol L⁻¹, but illustrates the generic features of all the meshes. The rotational symmetry of the simulation domain is exploited to calculate on a quasi-2D domain: the symmetry axis is indicated by the dashed red line. The domain typically has a radius L = 10a + 25κ⁻¹ in size. This domain is subdivided into two pieces on which triangular and quadrilateral elements are used. In a range of 3κ⁻¹ around the colloid the domain consists of quadrilaterals, which grow in size geometrically, see the zoom-in (blue box). Beyond this range the elements are triangular and are allowed to grow out linearly to best fit the domain boundary and reduce the overall number of elements.

To verify the analytic results, we first performed calculations with sufficiently low values of the surface charge density and flux to remain in the linear regime. These j^s and σ are given in Table 2. Fig. 7 in the main text shows that there is excellent correspondence between the theory and FEM calculations in this regime. Different fluxes were used for the different propulsion models because the low efficiency of type S_∘ and S₌ propulsion mean that numerical errors become significant more quickly as the flux density is reduced for these models.

Table 2 The charge densities and the first Legendre components of the surface flux densities used in Fig. 7 in the main text and Fig. 12 here. The flux densities have units mol m⁻² s⁻¹, and the charge densities have units e nm⁻². The final column gives the product of σ and the relevant non-zero flux density, with units e mol m⁻² nm⁻² s⁻¹

Fig.	Type	j ^s∘,1	j ^s+,1	j ^s−,1	σ	σj ^s
7	S ∘	3 × 10^−1	0	0	10^−4	3 × 10^−5
	S +	0	3 × 10^−7	0	10^−4	3 × 10^−11
	S =	0	3 × 10^−5	3 × 10^−5	10^−4	3 × 10^−9
12	S ∘	1.5 × 10^−2	0	0	10^−2	1.5 × 10^−4
	S +	0	1.5 × 10^−5	0	10^−2	1.5 × 10^−7
	S =	0	1.5 × 10^−5	1.5 × 10^−5	10^−2	1.5 × 10^−7

---

In addition, the FEM calculations and the linearized theory produce essentially identical electrostatic potential fields. Fig. 11a illustrates this for type S₊ electrophoresis. Note that we had to use a much smaller computational domain than we typically use (L = 6a here, rather than L = 10a + 25κ⁻¹), in order to show details in Fig. 11a. This means that the deviation from the theory, which stems from the ϕ = 0 boundary condition on the edge of the domain (Fig. 11b), occurs closer to the particle than in our regular calculations. However, the potential and flow-fields decay sufficiently rapidly that this does not affect the potential near the particle, or the propulsion speed beyond a few per cent.

Fig. 11 (a) Comparison of normalized surface-reaction-generated potential fields ϕ^sr for (left) FEM and (right) linearized theory, for S₊ self-electrophoresis, with j^s_+,1 = 3 × 10⁻⁷ mol m⁻² s⁻¹, and with other conditions as in the base parameter set (see main text). The radius of the simulation domain L = 3 µm = 6a here. (b) Normalized radial decay of ϕ^sr for linearized theory () and FEM calculations () along X–X′ in (a).

We can also use the FEM to go beyond the linear approximation. We defer to future work a systematic investigation of the non-linear behaviour, and here focus on the propulsion speed for selected experimentally relevant values of the surface charge density and chemical fluxes. These values are taken from measurements on the Pt-polystyrene Janus swimmers in ref. 21, and are listed in Table 2. The neutral flux density j^s_∘,1 is that which would be produced by a Janus particle which uniformly consumes H₂O₂ on one hemisphere at a rate Γ = 8 × 10¹⁰ molecules per second per particle. This rate was measured for a = 1 µm radius particles in 3 mol L⁻¹ H₂O₂.¹⁶ The surface charge density is taken from the electrophoretic mobility measurements made on the same particles in ref. 16. The ionic fluxes are unknown, but we arbitrarily set j^s_±,1 = 10⁻³j^s_∘,1, so that S₊ electrophoresis gives a speed of order 100 µm s⁻¹, which is somewhat larger than typical experimental values, ≈10 µm s⁻¹ for Au–Pt spherical microswimmers.^10,65 Hence, our results should overestimate the non-linear behaviour of the propulsion speed. Note that though the ionic flux densities for the experimentally realistic case are sometimes lower than those for the linear case, the product of charge density and surface flux is always greater in the experimentally realistic case, Table 2.

Fig. 12a and b, both with 1 mmol L⁻¹ NaCl, correspond to Fig. 7 in the main text. We see that the analytical theory continues to match the FEM calculations well even for these realistic values of the flux and charge densities. However, many experiments are performed with no added salt, and as shown in Fig. 12c, the agreement worsens as the salt concentration falls. This is to be expected, since it is low salt that generates a high-ζ, large-screening-length regime where linear approximations break down.⁷⁵ In fact, with 0 mmol L⁻¹ NaCl, the dimensionless zeta-potential ζe/(k_BT) = 5.6 for these particles, well beyond the Debye–Hückel regime of ζe/(k_BT) ≪ 1. Nevertheless, for all propulsion types, the agreement remains semi-quantitative between simulations and theory over the whole radius range for 0 mmol L⁻¹ NaCl, Fig. 12d.

Fig. 12 Propulsion speed for realistic parameters (given in this text) for type S₀ (, inset), S₊ (), and S₌ (●) propulsion, from analytical theory (solid curves) and FEM simulations (symbols). For (a) [H⁺] at fixed κ equivalent to 1 mmol L⁻¹ NaCl, (b) particle radius with 1 mmol L⁻¹ NaCl (c) NaCl concentration, (d) particle radius with 0 mmol L⁻¹ NaCl. In (d), the black arrow indicates the experimental point from ref. 16 referred to in the text.

From Fig. 12d, we obtain a speed of 0.5 µm s⁻¹ for type S_∘ electrophoresis with particles of radius a = 1 µm, no salt, and 3 M H₂O₂ (the black arrow indicates the relevant data point). As stated in the main text, this predicted speed can account for at most 5% of the experimentally measured propulsion speed of 15–20 µm s⁻¹ obtained for Pt-polystyrene Janus particles under the same conditions.¹⁶

Acknowledgements

This work was funded by UK EPSRC grant EP/J007404/1 (WP); ERC Advanced Grant ERC-2013-AdG 340877-PHYSAP (AB and WP); NWO Rubicon Grant #680501210 and Marie Skłodowska-Curie Intra European Fellowship (G.A. No. 654916) within Horizon 2020 (JdG); and DFG SPP 1726 ‘Microswimmers—from single particle motion to collective behaviour’ (JdG and CH; HO1108/24-1). We thank Patrick Kreissl for checking many of the calculations; and Mike Cates, Michael Kuron, Davide Marenduzzo, Alexander Morozov, Mihail Popescu, Georg Rempfer, Joakim Stenhammar, and Teun Vissers for useful comments on the manuscript.

References

1. A. Einstein, Ann. der Physik, 1905, 17, 549.
2. 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, J. Am. Chem. Soc., 2004, 126, 13424.
3. J. R. Howse, R. A. L. Jones, A. J. Ryan, T. Gough, R. Vafabakhsh and R. Golestanian, Phys. Rev. Lett., 2007, 99, 048102.
4. Y. Wang, R. M. Hernandez, D. J. Bartlett, J. M. Bingham, T. R. Kline, A. Sen and T. E. Mallouk, Langmuir, 2006, 22, 10451.
5. A. A. Solovev, Y. Mei, E. Bermúdez Ureña, G. Huang and O. G. Schmidt, Small, 2009, 5, 1688.
6. H.-R. Jiang, N. Yoshinaga and M. Sano, Phys. Rev. Lett., 2010, 105, 268302.
7. I. Buttinoni, G. Volpe, F. Kümmel, G. Volpe and C. Bechinger, J. Phys.: Condens. Matter, 2012, 24, 284129.
8. S. Ahmed, W. Wang, L. Mair, R. Fraleigh, S. Li, L. Castro, M. Hoyos, T. Huang and T. Mallouk, Langmuir, 2013, 29, 16113.
9. S. J. Ebbens and J. R. Howse, Soft Matter, 2010, 6, 726–738.
10. I. Theurkauff, C. Cottin-Bizonne, J. Palacci, C. Ybert and L. Bocquet, Phys. Rev. Lett., 2012, 108, 268303.
11. J. Palacci, S. Sacanna, A. P. Steinberg, D. J. Pine and P. M. Chaikin, Science, 2013, 339, 936.
12. I. Buttinoni, J. Bialké, F. Kümmel, H. Löwen, C. Bechinger and T. Speck, Phys. Rev. Lett., 2013, 110, 238301.
13. F. Ginot, I. Theurkauff, D. Levis, C. Ybert, L. Bocquet, L. Berthier and C. Cottin-Bizonne, Phys. Rev. X, 2015, 5, 011004.
14. J. L. Moran and J. D. Posner, J. Fluid Mech., 2011, 680, 31.
15. S. Ebbens, M.-H. Tu, J. R. Howse and R. Golestanian, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2012, 85, 020401.
16. A. T. Brown and W. C. K. Poon, Soft Matter, 2014, 10, 4016.
17. S. Ebbens, D. A. Gregory, G. Dunderdale, J. R. Howse, Y. Ibrahim, T. B. Liverpool and R. Golestanian, Europhys. Lett., 2014, 106, 58003.
18. R. Soto and R. Golestanian, Phys. Rev. Lett., 2014, 112, 068301.
19. W. E. Uspal, M. N. Popescu, S. Dietrich and M. Tasinkevych, Soft Matter, 2015, 11, 434.
20. E. J. Banigan and J. F. Marko, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2016, 93, 012611.
21. A. Brown, I. Vladescu, A. Dawson, T. Vissers, J. Schwarz-Linek, J. Lintuvuori and W. Poon, Soft Matter, 2016, 12, 131.
22. T. Vicsek, A. Czirók, E. Ben-Jacob, I. Cohen and O. Shochet, Phys. Rev. Lett., 1995, 75, 1226.
23. J. Stenhammar, A. Tiribocchi, R. J. Allen, D. Marenduzzo and M. E. Cates, Phys. Rev. Lett., 2013, 111, 145702.
24. X. Zheng, B. ten Hagen, A. Kaiser, M. Wu, H. Cui, Z. Silber-Li and H. Löwen, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2013, 88, 032304.
25. R. Matas-Navarro, R. Golestanian, T. Liverpool and S. Fielding, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2014, 90, 032304.
26. A. Zöttl and H. Stark, Phys. Rev. Lett., 2014, 112, 118101.
27. X. Yang and M. Marchetti, Phys. Rev. Lett., 2015, 115, 258101.
28. S. Thakur and R. Kapral, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2012, 85, 026121.
29. O. Pohl and H. Stark, Phys. Rev. Lett., 2014, 112, 238303.
30. T. Bickel, G. Zecua and A. Würger, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2014, 89, 050303.
31. B. Liebchen, D. Marenduzzo, I. Pagonabarraga and M. E. Cates, Phys. Rev. Lett., 2015, 115, 258301.
32. W. F. Paxton, P. T. Baker, T. R. Kline, Y. Wang, T. E. Mallouk and A. Sen, J. Am. Chem. Soc., 2006, 128, 14881.
33. R. Golestanian, T. B. Liverpool and A. Ajdari, New J. Phys., 2007, 9, 126.
34. J. Moran, P. Wheat and J. Posner, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2010, 81, 065302.
35. E. Gileadi, Physical Electrochemistry: Fundamentals, Techniques and Applications, Wiley-VCH Weinheim, Germany, 2011.
36. B. Sabass and U. Seifert, J. Chem. Phys., 2012, 136, 214507.
37. T.-C. Lee, M. Alarcón-Correa, C. Miksch, K. Hahn, J. G. Gibbs and P. Fischer, Nano Lett., 2014, 14, 2407.
38. B. Sabass and U. Seifert, J. Chem. Phys., 2012, 136, 064508.
39. D. Henry, Proc. R. Soc. London, Ser. A, 1931, 133, 106.
40. S. Kim and S. J. Karrila, Microhydrodynamics: principles and selected applications, Courier Dover Publications, 2013.
41. H. S. Muddana, S. Sengupta, T. E. Mallouk, A. Sen and P. J. Butler, J. Am. Chem. Soc., 2010, 132, 2110.
42. S. Sengupta, K. K. Dey, H. S. Muddana, T. Tabouillot, M. E. Ibele, P. J. Butler and A. Sen, J. Am. Chem. Soc., 2013, 135, 1406.
43. R. Golestanian, T. B. Liverpool and A. Ajdari, Phys. Rev. Lett., 2005, 94, 220801.
44. J. de Graaf, G. Rempfer and C. Holm, IEEE Trans. Nanobiosci., 2015, 14, 272.
45. J. L. Anderson, Ann. Rev. Fluid Mech., 1989, 21, 61.
46. W. Duan, M. Ibele, R. Liu and A. Sen, Eur. Phys. J. E: Soft Matter Biol. Phys., 2012, 35, 1.
47. H. A. Stone and A. D. T. Samuel, Phys. Rev. Lett., 1996, 77, 4102.
48. M. Teubner, J. Chem. Phys., 1982, 76, 5564.
49. G. Rempfer, G. Davies, C. Holm and J. de Graaf, J. Chem. Phys., 2016, 145, 044901.
50. P. Kreissl, C. Holm and J. de Graaf, J. Chem. Phys., 2016, 144, 204902.
51. A. J. Everett and G. J. Minkoff, Trans. Faraday Soc., 1953, 49, 410.
52. E. Caldin, Fast reactions in solution, Blackwell Scientific Publications, Oxford, 1964.
53. Y. Pocker and D. W. Bjorkquist, J. Am. Chem. Soc., 1977, 99, 6537.
54. D. M. H. Kern, J. Am. Chem. Soc., 1954, 76, 4208.
55. CRC handbook of chemistry and physics, ed. W. M. Haynes, CRC press, Boca Raton, USA, 93rd edn, 2013.
56. F. van den Brink, W. Visscher and E. Barendrecht, J. Electroanal. Chem., 1984, 172, 301.
57. E. Samson, J. Marchand and K. Snyder, Mater. Struct., 2003, 36, 156.
58. M. Ibele, T. E. Mallouk and A. Sen, Angew. Chem., Int. Ed., 2009, 48, 3308.
59. W. Gao, A. Pei, R. Dong and J. Wang, J. Am. Chem. Soc., 2014, 136, 2276.
60. J. L. Moran and J. D. Posner, Phys. Fluids, 2014, 26, 042001.
61. M. Esplandiu, A. Farniya and D. Reguera, J. Chem. Phys., 2016, 144, 124702.
62. B. R. Midmore, G. V. Pratt and T. M. Herrington, J. Colloid Interface Sci., 1996, 184, 170.
63. S. Das, A. Garg, A. I. Campbell, J. Howse, A. Sen, D. Velegol, R. Golestanian and S. J. Ebbens, Nat. Commun., 2015, 6, 8999.
64. M. Xuan, J. Shao, X. Lin, L. Dai and Q. He, ChemPhysChem, 2014, 15, 2255.
65. P. M. Wheat, N. A. Marine, J. L. Moran and J. D. Posner, Langmuir, 2010, 26, 13052.
66. R. Golestanian, Phys. Rev. Lett., 2015, 115, 108102.
67. Y. Ono, T. Matsumura, N. Kitajima and S. Fukuzumi, J. Phys. Chem., 1977, 81, 1307.
68. D. W. McKee, J. Catal., 1969, 14, 355.
69. A. Afshar Farniya, M. J. Esplandiu, D. Reguera and A. Bachtold, Phys. Rev. Lett., 2013, 111, 168301.
70. J. Y. Kim and B. J. Yoon, J. Colloid Interface Sci., 2003, 262, 101.
71. K. F. Riley, M. P. Hobson and S. J. Bence, Mathematical methods for physics and engineering: a comprehensive guide, Cambridge University Press, 2006.
72. R. W. O'Brien and L. R. White, J. Chem. Soc., Faraday Trans. 2, 1978, 74, 1607.
73. P. Debye, Trans. Electrochem. Soc., 1942, 82, 265.
74. Chemical Kinetics, ed. K. J. Laidler, Harper Collins, New York, USA, 3rd edn, 1987.
75. R. J. Hunter, Foundations of colloid science, Oxford University Press, 2001.

---

Footnotes

† The precise a³ factor in eqn (36) can be obtained from a scaling argument. We do not include this, as it is involved, and barely more informative than this qualitative explanation.

‡ Note that D* ≠ D̅ (eqn (29)), as D̅ includes the diffusivity of the inactive ions.

---