1. Introduction
A semi-permeable membrane allows the passage of solvent through its pores, but blocks the solute molecules. When the solute concentration differs on the two sides of such a membrane, osmosis takes place: solvent seeps through the membrane, flowing from the side of low solute concentration to the side of high solute concentration. In a linearised non-equilibrium description (De Groot & Mazur Reference De Groot and Mazur2013), the volumetric solvent flux is proportional to the jump in the osmotic pressure, offset by the mechanical pressure (Kedem & Katchalsky Reference Kedem and Katchalsky1958). The proportionality coefficient, known as the hydraulic coefficient of the membrane, is given by the ratio of the hydraulic conductivity
$k/\eta$
(where
$k$
is the Darcy permeability of the membrane, and
$\eta$
is the solvent viscosity) to the membrane thickness
$d$
. For small solute concentrations, the van’t Hoff expression gives the osmotic pressure as the product of
$RT$
with the solute molar concentration, where
$R$
is the universal gas constant, and
$T$
is the absolute temperature.
Certain lipid vesicles are made of semi-permeable membranes that bound a liquid core. When placed under an imposed solute concentration gradient, say of magnitude
$\alpha$
, the non-uniform jump of osmotic pressure results in an inward trans-membrane solvent flux in the solute-depleted side (‘front’) of the vesicle, and an outward flux in its solute-enriched side (‘back’). In the absence of external forces, a recoil-like mechanism necessitates that the freely suspended vesicle drifts down the gradient, a phenomenon known as ‘osmophoresis’ (Gordon Reference Gordon1981; Pope Reference Pope1982). With
$a$
as the vesicle size, typical jumps of the osmotic pressure scale as
$RT\alpha a$
. It follows that the characteristic velocity scale is
\begin{equation} \mathcal U = \frac {kRT\alpha a}{\eta {d}}. \end{equation}
Osmophoresis is relevant to cellular migration (Li et al. Reference Li, Boileau, Kung and Adler1988) and has been proposed as a possible mechanism for vesicle trafficking (Lipchinsky Reference Lipchinsky2015). In modelling, it is a standard practice to assume a spherical membrane (radius
$a$
) that, while permeable to the solvent, does not allow slippage. The fluid velocity therefore satisfies a no-slip condition, with the standard no-penetration condition replaced by a Darcy-like condition. In a linearised non-equilibrium continuum description, the problem governing the solute transport is quasi-steady in a co-moving reference frame. That problem was analysed by Anderson (Reference Anderson1983). When solute advection is neglected altogether, Anderson found that the vesicle drifts down the gradient with velocity
$\mathcal U/2$
.
In studying the motion of DMPC (dimyristoyl phosphatidycholine) vesicles in sucrose solutions, Nardi et al. (Reference Nardi, Bruinsma and Sackmann1999) found drift velocities more than three orders of magnitude larger than
$\mathcal U$
. They used rigid latex particles, of size comparable to that of the vesicles, as a ‘control group’. These impermeable particles did not exhibit any directed motion, thus excluding the possibility that other mechanisms animated by the imposed concentration gradient (in particular, buoyancy-driven convection) are responsible for the observed anomaly. In that experiment, the chamber was mounted on the stage of an inverted microscope, with the observed vesicles moving over the glass cover slip. Nardi et al. (Reference Nardi, Bruinsma and Sackmann1999) took great care to make sure that the vesicles did not adhere to the cover slip, but did not consider the possibility that the proximity to a boundary may be the very reason for the observed discrepancy.
In particular, Nardi et al. (Reference Nardi, Bruinsma and Sackmann1999) were apparently unaware of numerical observations of velocity amplification by Keh & Yang (Reference Keh and Yang1993). These authors used eigenfunction expansions in bi-spherical coordinates to solve the linearised problem of vesicle osmophoresis in the presence of an impermeable boundary, with the solute concentration gradient applied parallel to the boundary. In Stokes flow problems, the introduction of a no-slip boundary typically results in the diminution of the Stokes mobility (Goldman et al. Reference Goldman, Cox and Brenner1967) – an obvious consequence of the excess dissipation associated with the enhanced velocity gradients. Surprisingly, Keh & Yang (Reference Keh and Yang1993) numerically observed that the osmophoretic mobility actually increases monotonically with decreasing distance between the vesicle and the wall. Later, Chen & Keh (Reference Chen and Keh2003) revisited that problem using a collocation scheme, which is more suitable for small vesicle–wall separations.
The anomalously large velocities observed by the numerical computations are reminiscent of amplified electrophoretic mobility near a boundary (Keh & Chen Reference Keh and Chen1988). In that physiochemical problem, the amplification was rationalised by the intensification of the electric field in the narrow gap. Since electrophoresis is animated by electro-osmotic slip at the colloid boundary, with the velocity proportional to the electric field, such intensification results in a strong localised propulsion mechanism that overcomes the diminished Stokes mobility. The asymptotic analysis of that problem was carried out by Yariv & Brenner (Reference Yariv and Brenner2003). The prediction of mobility enhancement has been verified in two separate experiments where near-contact conditions were maintained (Xuan et al. Reference Xuan, Ye and Li2005; Liang et al. Reference Liang, Zhao and Yang2015).
When solute advection is negligible, there is a perfect analogy between the electric potential in the electrophoretic problem and the solute concentration in the osmophoretic problem: in particular, both scalar fields satisfy a no-flux condition on the particle (colloid or vesicle, respectively) boundary and on the wall (Schnitzer & Yariv Reference Schnitzer and Yariv2012). The fundamental difference between the two problems is the underlying flow mechanism. In the electrophoretic problem, flow is driven by a slip mechanism, with a velocity proportional to the gradient of the scalar field at a rigid boundary (Helmholtz–Smoluchowski slip); in the osmophoretic problem, it is driven by a trans-membrane seepage mechanism, with the velocity proportional to the scalar field itself.
We intend to investigate the osmophoretic problem using the linearised non-equilibrium description of Anderson (Reference Anderson1983). In that framework, the symmetry properties of Stokes flows (Jeffrey Reference Jeffrey and Torry1996; Leal Reference Leal2007) imply that the vesicle translates parallel to the wall and rotates about an axis that is perpendicular to the imposed gradient and parallel to the wall. Under that constraint, the particle maintains its distance from the wall. The fixed distance parametrises the vesicle–wall geometry.
Prior to formulating the problem, it is useful to provide typical values for the membrane parameters appearing in (1.1). Following Anderson (Reference Anderson1983), these are estimated as
\begin{equation} k\approx 10^{-16}\,\mbox{cm}^2, \quad {d} \approx 10^{-6}\, \mbox{cm}, \end{equation}
where the estimate for
$k$
assumes a membrane with pores of radius
$10$
Å that comprise
$10\, \%$
of the membrane volume. This is probably an upper limit; indeed, realistic data of lipid membranes yield significantly lower values (Bloom et al. Reference Bloom, Evans and Mouritsen1991; Mathai et al. Reference Mathai, Tristram-Nagle, Nagle and Zeidel2008).
2. Physical problem and governing equations
A spherical vesicle (radius
$a$
) bounded by a semi-permeable membrane (Darcy permeability
$k$
, thickness
$d$
) is freely suspended in a liquid solution (viscosity
$\eta$
) adjacent to an impermeable rigid wall, the distance between the vesicle centroid and the wall being
$a(1+\delta )$
. A solute concentration gradient, of magnitude
$\alpha$
, is externally imposed parallel to the wall. Our goal is the calculation of the rectilinear velocity
$\mathcal V$
of the vesicle centroid and the vesicle angular velocity
$\unicode {x1D6FA}$
; see figure 1(a).
Figure 1. (a) Schematic of the geometry in a laboratory reference frame. (b) Dimensionless geometry in a co-moving reference frame, with the thin arrows indicating the trans-membrane seepage.
We employ the simplest linearised non-equilibrium model, neglecting solute advection. Thus solute transport is purely diffusive, satisfying a no-flux condition at the membrane. Moreover, we neglect the solvent flux due to trans-membrane differences of out-of-equilibrium (‘dynamic’) pressure: with osmotic pressure differences being of order
$\eta \,\mathcal U {d}/ k$
and the dynamic pressure variations scaling as
$\eta \,\mathcal U / a$
, this assumption is justified for
$k / a{d}\ll 1$
. For the estimates (1.2),
$k / a{d}$
is approximately
$10^{-6}$
for a typical 10 micron vesicle. Finally, with small solute concentration (as tacitly assumed in using the van’t Hoff expression), the volumetric solvent flux effectively coincides with the liquid velocity.
In what follows, we normalise all length variables by
$a$
. We use a co-moving reference frame that translates with the centroid
$\mathcal C$
but does not rotate. In that reference frame, we employ the
$(x,y,z)$
Cartesian coordinates, with corresponding unit vectors denoted by
$(\hat {\boldsymbol{\imath }},\hat {\boldsymbol{\jmath }},\hat {\boldsymbol{{k}}})$
. The
$x$
-coordinate points down the gradient, the wall coincides with
$z=0$
, and the liquid domain is within the half-space
$z\gt 0$
; see figure 1(b). The membrane surface
$\mathcal M$
is accordingly given by the spherical surface
\begin{equation} x^2+y^2+(z-1-\delta )^2=1. \end{equation}
The solute concentration outside the vesicle, measured relative to the ambient concentration at
$\mathcal C$
and normalised by
$\alpha a$
, is denoted by
$c$
. Using linearised non-equilibrium transport with negligible solute advection, it is governed by: (i) Laplace’s equation,
\begin{equation} \nabla ^2c=0\quad \text{outside} \ \mathcal M ; \end{equation}
(ii) the no-flux condition at the membrane,
\begin{equation} \frac {\partial c}{\partial n} =0 \quad \text{at} \ \mathcal M, \end{equation}
with
$\partial /\partial n=\hat {\boldsymbol{n}}\boldsymbol{\cdot }\boldsymbol{\nabla }$
, where
$\hat {\boldsymbol{n}}$
is the outward-pointing unit normal to
$\mathcal M$
; (iii) the no-flux condition at the wall,
\begin{equation} \frac {\partial c}{\partial z} = 0 \quad \text{at} \ z=0; \end{equation}
and (iv) the imposed gradient condition,
\begin{equation} c \sim -x \quad \text{as} \ |\boldsymbol{x}|\to \infty . \end{equation}
In analysing the flow problem, we normalise velocity variables by
$\mathcal U$
, as defined by (1.1), and stress variables by
$\eta \, \mathcal U / a$
. The velocity field
$\boldsymbol{u}$
and pressure field
$p$
are governed by: (i) the Stokes equations,
\begin{equation} \boldsymbol{\nabla }\boldsymbol{\cdot } \boldsymbol{u} = 0, \quad \boldsymbol{\nabla } p = {\nabla }^2\boldsymbol{u}; \end{equation}
(ii) the no-slip condition on the membrane,
\begin{equation} \boldsymbol{u} = \hat {\boldsymbol{n}} c + \frac {a\unicode {x1D6FA}}{\mathcal U}\hat {\boldsymbol{\jmath }}\times \hat {\boldsymbol{n}} \quad \text{on}\ \mathcal M, \end{equation}
where the first term represents osmotic seepage, in which we have employed the van’t Hoff expression for osmotic pressure, and the second term represents rigid-body rotation; (iii) the no-slip condition at the wall,
\begin{equation} \boldsymbol{u} = -\hat {\boldsymbol{\imath }}\mathcal{V}/{\mathcal U} \quad \text{at} \ z=0; \end{equation}
and (iv) the far-field condtion,
\begin{equation} \lim _{|\boldsymbol{x}|\to \infty }\boldsymbol{u} = -\hat {\boldsymbol{\imath }}\mathcal{V}/{\mathcal U}. \end{equation}
The flow problem, animated by the osmotic term in (2.7), is ‘closed’ by the requirement that the vesicle is force- and torque-free. These conditions serve to determine
$\mathcal V$
and
$\unicode {x1D6FA}$
.
We note that the solute concentration inside the membrane does not appear in the above formulation. With the membrane impermeable to solute and with advection neglected, the interior concentration is simply retained at the value it had upon introducing the vesicle into the solution; that value affects neither the transport problem nor the flow problem.
The solute transport problem is decoupled from the flow, and can be solved independently. It may superficially appear that
$c$
is defined by the solute transport problem up to an additive constant, but such a freedom is incompatible with the dependence of (2.7) upon
$c$
. In fact, it is that very dependence that serves to resolve the indeterminacy: since the membrane core is incompressible,
$\oint _{\mathcal M} {\textrm{d}} A \, \hat {\boldsymbol{n}}\boldsymbol{\cdot }\boldsymbol{u}=0$
. Substituting (2.7) yields the constraint
\begin{equation} \oint _{\mathcal M} {\textrm{d}} A \, c = 0 \end{equation}
that serves to uniquely determine
$c$
.
3. Simplifications
3.1. Decomposition
For convenience, we decompose the flow problem into two subproblems, both satisfying (2.6). The first represents a stationary vesicle, where the boundary condition at the membrane accounts for the osmotic term in (2.7),
\begin{equation} \boldsymbol{u}^{{\textrm{I}}} = \hat {\boldsymbol{n}} c \quad \text{on} \ \mathcal M, \end{equation}
while the counterparts of (2.8) and (2.9) are homogeneous:
\begin{align} \boldsymbol{u}^{{\textrm{I}}} = \boldsymbol{0} \quad \text{at} \ z=0, \qquad \lim _{|\boldsymbol{x}|\to \infty }\boldsymbol{u}^{{\textrm{I}}} = \boldsymbol{0} . \end{align}
Due to the underlying linearity of the governing equations and the symmetry properties of the Stokes equations, the associated hydrodynamic force must act in the
$x$
-direction while the associated hydrodynamic torque about
$\mathcal C$
must act the
$y$
-direction. We write the force (normalised by
$\eta a \,\mathcal U$
) as
$f \hat {\boldsymbol{\imath }}$
, and the torque (normalised by
$\eta a^2 \,\mathcal U$
) as
$g \hat {\boldsymbol{\jmath }}$
.
In the second subproblem, the boundary condition at the membrane accounts for the rigid-body rotation in (2.7),
\begin{equation} \boldsymbol{u}^{{\textrm{II}}} = \frac {a\unicode {x1D6FA}}{\mathcal U}\hat {\boldsymbol{\jmath }}\times \hat {\boldsymbol{n}} \quad \text{on} \ \mathcal M, \end{equation}
while the conditions on the wall and at large distances are given by (2.8) and (2.9). In the laboratory frame, the fluid velocity vanishes on the wall and at large distances, while on
$\mathcal M$
it is given by
\begin{equation} \frac {\mathcal V}{\mathcal U}\hat {\boldsymbol{\imath }} + \frac {a\unicode {x1D6FA}}{\mathcal U}\hat {\boldsymbol{\jmath }}\times \hat {\boldsymbol{n}}. \end{equation}
The boundary velocity (3.4) represents a superposition of rigid-body translation and rotation. Given the linearity of the flow problem, it results in a force (normalised by
$\eta a \,\mathcal U$
) in the
$x$
-direction of the form
$-f^{\textit{tr}}\mathcal{V}/\mathcal U - f^{\textit{rot}} a\unicode {x1D6FA}/\mathcal U$
, and a torque (normalised by
$\eta a^2 \,\mathcal U$
) in the
$y$
-direction of the form
$-g^{\textit{rot}} a\unicode {x1D6FA}/\mathcal U - g^{\textit{tr}} \mathcal{V}/\mathcal U$
. The symmetry properties of Stokes flow (Happel & Brenner Reference Happel and Brenner1965) necessitate that the coupling coefficients
$f^{\textit{rot}}$
and
$g^{\textit{tr}}$
coincide, =
$h$
, say. The force and torque coefficients
$f^{\textit{tr}}$
and
$g^{\textit{rot}}$
, as well as the coupling coefficient
$h$
, may be considered as known functions of
$\delta$
.
Applying the conditions of a force- and torque-free particle to the combined velocity field
$\boldsymbol{u}^{{\textrm{I}}}+\boldsymbol{u}^{{\textrm{II}}}$
gives
\begin{equation} {\unicode{x1D64D}} \left (\begin{matrix} \mathcal{V}/\mathcal U\\ {a\unicode {x1D6FA}}/{\mathcal U} \end{matrix}\right ) = \left (\begin{matrix} f \\ g \end{matrix}\right ), \end{equation}
wherein
\begin{equation} {\unicode{x1D64D}} = \left (\begin{matrix} f^{\textit{tr}} & h \\ h & g^{\textit{rot}} \end{matrix}\right ) \end{equation}
is the grand resistance matrix (Kim & Karrila Reference Kim and Karrila2005). Inversion of (3.5) gives
\begin{equation} \left (\begin{matrix} \mathcal{V}/\mathcal U\\ {a\unicode {x1D6FA}}/{\mathcal U} \end{matrix}\right ) = {\unicode{x1D64D}}^{-1} \left (\begin{matrix} f \\ g \end{matrix}\right ). \end{equation}
At this stage, the calculation of
$\mathcal V$
and
$\unicode {x1D6FA}$
requires: (i) the solution of the concentration problem; (ii) the solution of the first flow subproblem, animated by (3.1); and (iii) the evaluation of the resulting loads,
$f$
and
$g$
. The rigid-body velocities are then obtained from (3.7).
3.2. Azimuthal symmetry
In what follows, we make use of cylindrical
$(\rho, \psi, z)$
coordinates, with
$\psi =0$
in the
$x$
-direction. The associated unit vectors are denoted by
$(\hat {\boldsymbol{e}}_\rho, \hat {\boldsymbol{e}}_\psi, \hat {\boldsymbol{{k}}})$
. The membrane surface (2.1) becomes
\begin{equation} \rho ^2+(z-1-\delta )^2=1. \end{equation}
Since
$x=\rho \cos \psi$
, the linearity of the governing equations in conjunction with (2.5) implies the following symmetry of the concentration,
\begin{equation} c(\rho, \psi, z) = \grave {c}(\rho, z)\cos \psi, \end{equation}
which trivially satisfies the consistency condition (2.10). Writing the velocity
$\boldsymbol{u}^{{\textrm{I}}}$
as
$\hat {\boldsymbol{e}}_\rho u + \hat {\boldsymbol{e}}_\psi v +\hat {\boldsymbol{{k}}}w$
, we anticipate the factorised form (O’Neill Reference O’Neill1964)
\begin{align} u(\rho, \psi, z) = \grave {u}(\rho, z)\cos \psi, \quad v(\rho, \psi, z) = \grave {v}(\rho, z)\sin \psi, \quad w(\rho, \psi, z) = \grave {w}(\rho, z)\cos \psi, \end{align}
and then also
\begin{equation} p^{{\textrm{I}}}(\rho, \psi, z) = \grave {p}(\rho, z)\cos \psi . \end{equation}
Laplace’s equation (2.2) becomes
\begin{equation} \frac {\partial ^2\grave {c}}{\partial z^2} + {\mathscr{L}}_1\grave {c} = 0, \end{equation}
wherein
\begin{equation} {\mathscr{L}}_m = \frac {\partial ^2}{\partial \rho ^2} + \frac {1}{\rho }\frac {\partial }{\partial \rho }-\frac {m^2}{\rho ^2}. \end{equation}
The continuity equation (2.6a ) becomes
\begin{equation} \frac {\partial \grave {u}}{\partial \rho } + \frac {\grave {u}+\grave {v}}{\rho } + \frac {\partial \grave {w}}{\partial z}=0, \end{equation}
while the Stokes equation (2.6b ) reads, in the radial, azimuthal and axial directions, respectively,
\begin{align} \frac {\partial \grave {p}}{\partial \rho } = \frac {\partial ^2\grave {u}}{\partial z^2} + {\mathscr{L}}_0\grave {u} - 2\frac {\grave {u}+\grave {v}}{\rho ^2}, \quad -\frac {\grave {p}}{\rho } = \frac {\partial ^2\grave {v}}{\partial z^2} + {\mathscr{L}}_0\grave {v} - 2\frac {\grave {u}+\grave {v}}{\rho ^2}, \quad \frac {\partial \grave {p}}{\partial z} = \frac {\partial ^2\grave {w}}{\partial z^2}+{\mathscr{L}}_1\grave {w}. \end{align}
4. Near-contact limit
As discussed in the Introduction, the osmophoretic problem was solved by Anderson (Reference Anderson1983) for a vesicle in free solution, corresponding to
$\delta =\infty$
in the present formulation. In that limit,
$\mathcal{V}/\mathcal{U}=1/2$
. Motivated by the numerical observations of Chen & Keh (Reference Chen and Keh2003), we here address the opposite limit,
\begin{equation} \delta \ll 1. \end{equation}
As a preliminary step to the subsequent analysis, we consider the geometry of the narrow gap separating
$\mathcal M$
from the wall. It follows from (3.8) that the ‘lower’ portion of
$\mathcal M$
is given by the equation
\begin{equation} z = 1+\delta -\sqrt {1-\rho ^2}, \end{equation}
which in the gap region, where
$\rho \ll 1$
, becomes
\begin{equation} z= \delta + \rho ^2/2 + \cdots . \end{equation}
Since
$z$
scales as
$\delta$
in the narrow gap, it follows that
$\rho$
scales there as
$\delta ^{1/2}$
. The gap region is accordingly handled using the stretched variables
\begin{equation} R = \rho /\delta ^{1/2}, \quad Z = z/\delta . \end{equation}
In terms of these, the surface (4.3) becomes
\begin{equation} Z = H(R;\delta ), \end{equation}
wherein
\begin{equation} H(R;\delta ) = H_0(R) + O(\delta ), \end{equation}
in which
\begin{equation} H_0(R) = 1 + \frac {R^2}{2}. \end{equation}
The outward (non-unit) vector normal to
$\mathcal M$
,
\begin{equation} -\hat {\boldsymbol{{k}}} + \delta ^{1/2} \hat {\boldsymbol{e}}_\rho \left [ \frac {{\textrm{d}} H_0}{{\textrm{d}} R} + O(\delta ) \right ], \end{equation}
is primarily in the negative
$z$
-direction.
4.1. Transport problem
Consider first the transport problem governing
$\grave {c}$
in the gap. It consists of: (i) the partial differential equation (cf. (3.12))
\begin{equation} \frac {\partial ^{2}\grave {c}}{\partial Z^{2}}+\delta \tilde {\mathscr{L}} \grave {c} = 0, \end{equation}
wherein (cf. (3.13))
\begin{equation} \tilde {\mathscr{L}}= \frac {\partial ^2}{\partial R^2} + \frac {1}{R}\frac {\partial }{\partial R}-\frac {1}{R^2}; \end{equation}
(ii) the boundary condition at the wall (cf. (2.4))
\begin{equation} \frac {\partial \grave {c}}{\partial Z}=0\quad \textrm{at}\ Z=0; \end{equation}
(iii) the boundary condition on
$\mathcal M$
, obtained from (2.3) and (4.8),
\begin{equation} \frac {\partial \grave {c}}{\partial Z}=\delta \left ( \frac {{\textrm{d}} H_0}{{\textrm{d}} R}+\cdots \right ) \frac {\partial \grave {c}}{\partial R}\quad \textrm{at}\ Z=H(R;\delta ); \end{equation}
and (iv) the requirement of regularity at
$R=0$
.
Since the imposed gradient (2.5) does not affect the gap problem, the problem posed by (4.9)–(4.12) is homogeneous. As such, it defines
$\grave {c}$
up to an arbitrary pre-factor, say
$\mu (\delta )$
. We therefore write
\begin{equation} \grave {c} = \mu (\delta )\, C, \end{equation}
where
$C$
is
${\textrm {ord}}(1)$
as
$\delta \to 0$
. At this stage, then, the goal is the calculation of
$C$
at leading order, with the pre-factor
$\mu (\delta )$
remaining to be determined. This task was carried out by Solomentsev et al. (Reference Solomentsev, Velegol and Anderson1997), who investigated the analogous problem of heat conduction between two spheres. Here, we summarise the key steps for completeness (see also Yariv & Brenner Reference Yariv and Brenner2003).
We employ the asymptotic expansion
\begin{equation} C(R,Z;\delta ) \sim C_0(R,Z) + \delta\, C_1(R,Z) + \cdots . \end{equation}
It readily follows from (4.9)–(4.12) at
${\textrm {ord}}(1)$
that
$C_0$
is independent of
$Z$
, say
$C_0(R)$
. To determine that function, we need a solvability condition. Consider then (4.9) and (4.11)–(4.12) at
${\textrm {ord}}(\delta )$
:
\begin{gather} \frac {\partial ^{2}C_1}{\partial Z^{2}} =-\tilde {\mathscr{L}} C_0 \quad \textrm{for}\ 0 \lt Z \lt H_0, \end{gather}
\begin{gather} \frac {\partial C_1}{\partial Z} =0\quad \textrm{at}\ Z=0, \end{gather}
\begin{gather} \frac {\partial C_1}{\partial Z} =\frac {{\textrm{d}} H_0}{{\textrm{d}} R}\frac {{\textrm{d}} C_1}{{\textrm{d}} R}\quad \textrm{at}\ Z=H_0. \end{gather}
By integrating (4.15a
) from
$Z=0$
to
$Z=H_0$
, and making use of (4.15b
) and (4.15c
), we obtain the second-order ordinary differential equation
\begin{equation} H_0\left (\frac {{\textrm{d}} ^2C_0}{{\textrm{d}} R^2} + \frac {1}{R}\frac {{\textrm{d}} C_0}{{\textrm{d}} R}-\frac {C_0}{R^2}\right ) +R\frac {{\textrm{d}} C_0}{{\textrm{d}} R}=0. \end{equation}
One of the solutions of (4.16) is
\begin{equation} C_0=-{D}R\,\left [ H_0(R)\right ] ^{\beta -1} {}_2F_1\left ( 1-\beta, 1-\beta, 2,\frac {R^{2}}{2H_0(R)}\right ), \end{equation}
wherein
$\beta =1/\sqrt {2}$
,
${}_2F_1$
is the hypergeometric function, and
$D$
is an integration constant (which, just like
$\mu (\delta )$
, remains undetermined at this stage). The other solution is
${\textrm {ord}}(R^{-1})$
as
$R\to 0$
, and is rejected by the requirement of regularity.
For future reference, we note that
\begin{equation} C_0 \sim -D(R +\cdots ) \quad \text{for} \ R\ll 1 \end{equation}
and
\begin{equation} C_0 \sim -{D}\,\frac {2^{1-\beta }\,\Gamma \left ( \sqrt {2}\right ) }{\left [ \Gamma \left ( 1+\beta \right ) \right ]^{2}}R^{\sqrt {2}-1} \quad \text{for} \ R\gg 1. \end{equation}
4.2. Flow problem: scaling arguments
In the gap region, (3.1) gives (recall (4.8))
\begin{align} u = \delta ^{1/2} \frac {{\textrm{d}} H_0}{{\textrm{d}} R} c \left [ 1 + O(\delta ) \right ], \quad v = 0, \quad w = - c \left [ 1 + O(\delta ) \right ] \quad \text{at} \ Z=H(R;\delta ). \end{align}
Using (3.9) and (3.10), these conditions become
\begin{align} \grave {u} = \delta ^{1/2} \frac {{\textrm{d}} H_0}{{\textrm{d}} R} \grave {c} \left [ 1 + O(\delta ) \right ], \quad \grave {v} = 0, \quad \grave {w} = - \grave {c} \left [ 1 + O(\delta ) \right ] \quad \text{at} \ Z=H(R;\delta ). \end{align}
It follows from (4.13) and (4.21c
) that
$\grave {w}$
scales as
$\mu (\delta )$
in the gap.
Consider now the flow problem. In the lubrication limit, the right-hand sides of the Stokes equations (3.15) are dominated by the
$\partial ^2/\partial z^2$
term. We see from (3.14) that both
$\grave {u}$
and
$\grave {w}$
scale as the product of
$\delta ^{-1/2}$
and
$\grave {w}$
, namely as
$\mu (\delta )\,\delta ^{-1/2}$
. Thus both (3.15a
) and (3.15b
) imply that
$\grave p$
scales as the product of
$\delta ^{-3/2}$
and these two fields, namely as
$\mu (\delta )\,\delta ^{-2}$
. (Equation (3.15c
) implies that the leading-order pressure is independent of
$z$
.)
The above scalings suggest writing
\begin{align} \grave {u} = \mu (\delta )\,\delta ^{-1/2} U, \quad \grave {v} = \mu (\delta )\,\delta ^{-1/2} V, \quad \grave {w} = \mu (\delta )\, W, \quad \grave p = \mu (\delta )\,\delta ^{-2} P. \end{align}
With the pertinent area of
$\mathcal M$
scaling as
$\delta$
(recall (4.3) et seq.), the associated contributions to the force
$f$
and torque
$g$
scale as
$\mu (\delta )\,\delta ^{-1/2}$
. Anticipating
\begin{equation} \mu (\delta )\,\delta ^{-1/2} \gg 1, \end{equation}
these ‘inner’ contributions dominate the corresponding contributions from the ‘outer’ region outside the gap, on the scale of the vesicle, which are presumably
${\textrm {ord}}(1)$
. We therefore write
\begin{align} f = \mu (\delta )\,\delta ^{-1/2} F, \quad g = \mu (\delta )\,\delta ^{-1/2} G. \end{align}
4.3. Flow problem: formulation
The variables
$U$
,
$V$
,
$W$
and
$P$
are functions of
$R$
and
$Z$
, as well as
$\delta$
. For each of these rescaled variables, we employ an expansion of the form (4.14). Our interest is in the leading-order fields.
Substituting (4.22) into (3.14) and (3.15) yields
\begin{equation} \frac {\partial U_0}{\partial R} + \frac {U_0+V_0}{R}+\frac {\partial W_0}{\partial Z}=0 \end{equation}
and
\begin{align} \frac {\partial P_0}{\partial R} =\frac {\partial ^2U_0}{\partial Z^2}, \quad -\frac {P_0}{R} = \frac {\partial ^2V_0}{\partial Z^2}, \quad \frac {\partial P_0}{\partial Z} = 0. \end{align}
The boundary conditions at the wall follow from (3.2a ) as
\begin{align} U_0=0, \quad V_0=0, \quad W_0=0 \quad \text{at} \ Z=0. \end{align}
The boundary conditions at the membrane follow from (4.21) and (4.22) as
\begin{align} U_0=0, \quad V_0=0, \quad W_0=-C_0(R) \quad \text{at} \ Z=H_0(R). \end{align}
The associated force is
\begin{equation} F_0 = -\unicode{x03C0} \int _0^\infty \left (RP_0 + \frac {\partial U_0}{\partial Z}-\frac {\partial V_0}{\partial Z}\right )_{Z=H_0} R\, {\textrm{d}} R, \end{equation}
while the torque about
$\mathcal C$
is
\begin{equation} G_0 = \unicode{x03C0} \int _0^\infty \left (\frac {\partial U_0}{\partial Z}-\frac {\partial V_0}{\partial Z}\right )_{Z=H_0} R\, {\textrm{d}} R. \end{equation}
4.4. Flow problem: lubrication analysis
From (4.26c
), we conclude that
$P_0$
is a function of
$R$
alone, say
$P_0(R)$
. The solution of (4.26a
), (4.27a
) and (4.28a
) is
\begin{equation} U_0 = \frac {1}{2}\frac {{\textrm{d}} P_0}{{\textrm{d}} R}Z[Z-H_0(R)]. \end{equation}
Similarly, the solution of (4.26b ), (4.27b ) and (4.28b ) is
\begin{equation} V_0 = -\frac {P_0}{2R}Z[Z-H_0(R)]. \end{equation}
Substitution of (4.31)–(4.32) into (4.25), followed by integration from
$Z=0$
to
$Z=H_0$
, yields, upon using (4.27c
) and (4.28c
),
\begin{equation} R^2H_0^3 \frac {{\textrm{d}} ^2P_0}{{\textrm{d}} R^2} + (3R^3H_0^2 + RH_0^3) \frac {{\textrm{d}} P_0}{{\textrm{d}} R} - H_0^3P_0 = -12 R^2 C_0(R), \end{equation}
which is an appropriate Reynolds equation. Also, substitution of (4.31)–(4.32) into (4.29) yields
\begin{equation} F_0 = -\frac {\unicode{x03C0} }{2}\int _0^\infty \left [(2R^2+H_0)P_0(R) +RH_0\frac {{\textrm{d}} P_0}{{\textrm{d}} R}\right ] \, {\textrm{d}} R. \end{equation}
For the integral to exist, it is required that
\begin{equation} P_0\ll R^{-1} \quad \text{for} \ R\ll 1 \end{equation}
and
\begin{equation} P_0\ll R^{-3} \quad \text{for} \ R\gg 1. \end{equation}
Integration by parts using (4.35)–(4.36) yields
\begin{equation} F_0 = -\frac {\unicode{x03C0} }{2}\int _0^\infty R^2\, P_0(R)\, {\textrm{d}} R. \end{equation}
A similar manipulation of (4.29) reveals that
\begin{equation} G_0 = F_0 . \end{equation}
4.5. Integrating the Reynolds equation
We now need to solve (4.33), with
$C_0$
given by (4.17). We immediately encounter three obstacles. The first is that we do not know at this stage the value of
$D$
, which appears as a multiplicative factor in (4.17). This obstacle is readily dealt with, since the differential equation (4.33) is linear. The second obstacle is technical: given the form (4.17), we have not been able to solve (4.33) in closed form. The Reynolds equation must therefore be integrated numerically. The third obstacle is fundamental: there are no auxiliary conditions for the second-order equation (4.33).
In the absence of such conditions, we need to rely upon the asymptotic behaviour of
$P_0$
for both small and large
$R$
. Using (4.18), we see that the particular integral of (4.33) satisfies
\begin{equation} \sim \frac {3}{2}DR^3 + \cdots \quad \text{for} \ R\ll 1. \end{equation}
Similarly, using (4.19) we see that the particular integral also satisfies
\begin{equation} \sim -24{D}(\sqrt {2}-1)\,\frac {2^{1-\beta }\,\Gamma \left ( \sqrt {2}\right ) }{\left [ \Gamma \left ( 1+\beta \right ) \right ] ^{2}}R^{\sqrt {2}-5} \quad \text{for} \ R\gg 1. \end{equation}
For small
$R$
, the behaviour of the complementary solutions is obtained from a local analysis about the regular singular point
$R=0$
. One complementary solution is
$P_{{c}}^{(1)} = R-3R^3/8 + \cdots$
. The other solution,
$P_{{c}}^{(2)} =1/R + \cdots$
, is rejected by (4.35). At large
$R$
, the behaviour of the complementary solutions is obtained from a local analysis about the regular singular point
$R=\infty$
. One complementary solution,
$R^{-\sqrt {3}-10} + \cdots$
, is subdominant to (4.40); the other,
$R^{-\sqrt {3}+10} + \cdots$
, is rejected by (4.36).
We integrated (4.33) numerically with
$D$
set to unity. As a substitute for auxiliary conditions, we used the small-
$R$
behaviour of
$P_0$
to construct effective initial conditions at a small value of
$R$
. This procedure was carried out twice. In the first instance, we employed (4.39), with
$D$
set to unity; in the second instance, where the right-hand side of (4.33) was set to zero, we employed
$P_{{c}}^{(1)}$
. The resulting particular integral and complementary solution are denoted by
$P_{{p}}$
and
$P_{{c}}$
, respectively. The permissible solution is
$P_{{p}}+EP_{{c}}$
, where the constant
$E$
is chosen such that
$P_0$
behaves like (4.40) at large
$R$
(again with
$D$
set to unity). This choice effectively eliminates the complementary mode that diverges as
$R^{-\sqrt {3}+10}$
at large
$R$
. Note that two terms have been required in the small-
$R$
expansion of
$P_{{c}}^{(1)}$
to be consistent with the concurrent use of (4.39). Due to the problem linearity, the requisite solution
$P_0$
is provided by
$D(P_{{p}}+EP_{{c}})$
.
Substitution of the resulting pressure distribution into (4.37)–(4.38) gives
\begin{equation} \frac {F_0}{\unicode{x03C0} D} =\frac {G_0}{\unicode{x03C0} D} =6.7004\ldots . \end{equation}
4.6. Reciprocal theorem
By making use of a methodology due to Brenner (Reference Brenner1964), it is possible to obtain
$F_0$
and
$G_0$
without actually solving the flow in the stationary particle subproblem. The starting point is the Lorentz reciprocal theorem of Stokes flow (Happel & Brenner Reference Happel and Brenner1965; Masoud & Stone Reference Masoud and Stone2019), involving two solutions of the Stokes equations. One is chosen as
$\boldsymbol{u}^{{\textrm{I}}}$
. The other is chosen as either the ‘translational’ flow field, associated with the rectilinear motion of a rigid particle with unit velocity in the
$x$
-direction, or the ‘rotational’ flow field, associated with the angular motion of a rigid particle about
$\mathcal C$
with unit velocity in the
$y$
-direction.
We prefer the ‘direct’ approach, where the flow field is calculated in the lubrication approximation; we believe that it is more illuminating. Moreover, use of the reciprocal theorem requires some familiarity with the details of the near-contact solutions of the translational (O’Neill & Stewartson Reference O’Neill and Stewartson1967) and rotational (Cooley & O’Neill Reference Cooley and O’Neill1968) problems. We therefore just present the final results of the alternative calculation using the reciprocal method:
\begin{equation} F_0 = -\unicode{x03C0} \int _0^\infty R\, P_0^{\textit{tr}}(R)\, C_0(R)\, {\textrm{d}} R, \quad G_0 = -\unicode{x03C0} \int _0^\infty R\, P_0^{\textit{rot}}(R)\, C_0(R) \, {\textrm{d}} R. \end{equation}
Here,
$P_0^{\textit{tr}}$
is the leading-order (factorised) pressure field appearing in the gap scale solution of the translational problem, while
$P_0^{\textit{rot}}$
is the corresponding field appearing in the gap scale solution of the rotational problem. The field
$P_0^{\textit{tr}}$
was calculated by O’Neill & Stewartson (Reference O’Neill and Stewartson1967) as
\begin{equation} P_0^{\textit{tr}} = \frac {6R}{5H_0^2}. \end{equation}
Later, Cooley & O’Neill (Reference Cooley and O’Neill1968) found that, by what seems like a happy coincidence,
$P_0^{\textit{rot}} = P_0^{\textit{tr}}$
. This provides an independent confirmation of (4.38). Substitution of (4.17) and (4.43) into (4.42) reproduces (4.41).
4.7.
Pre-factor scaling
$\mu (\delta )$
Up to this point, we managed to proceed without knowledge of the coefficient
$D$
and the pre-factor
$\mu (\delta )$
. To obtain these quantities, it is necessary to consider the outer transport problem, on the scale of the vesicle. That problem was solved by Yariv & Brenner (Reference Yariv and Brenner2003) in a different physical scenario using tangent-spheres coordinates, with the solution given in the form of a Hankel transform. That solution becomes singular at the origin, where the sphere appears to be touching the wall. Asymptotic matching using (4.13) and (4.19) gives
\begin{equation} \mu (\delta ) = \delta ^{1/\sqrt {2}-1/2} \end{equation}
and
\begin{equation} {D}=1.067437\ldots . \end{equation}
Similar results have been obtained by Solomentsev et al. (Reference Solomentsev, Velegol and Anderson1997) by comparing the gap scale solution with an exact solution in the form of an eigenfunction expansion. Note that (4.44) provides an a posteriori justification of (4.23).
Figure 2. Osmophoretic force
$f$
or osmophoretic torque
$g$
on a stationary vesicle. The solid line is the asymptotic prediction (4.24a
), using (4.41) and (4.45). The symbols show numerical results, obtained from table 1 of Chen & Keh (Reference Chen and Keh2003). The circles indicate
$f$
values; the squares indicate
$g$
values.
In figure 2, we compare our asymptotic approximation for the osmophoretic loads with the results predicted for
$f$
and
$g$
by the numerical computations of Chen & Keh (Reference Chen and Keh2003) in the absence of solute advection. The solid line is the leading-order approximation to (4.24), obtained using (4.41) and (4.44)–(4.45). That approximation is a power law, which appears on the log-log scale as a straight line of slope
$1/\sqrt {2}-1$
. The symbols indicate the numerical results, obtained from table 1 of Chen & Keh (Reference Chen and Keh2003). (Only results for
$\delta \lt 1$
are shown here.) Note that the coefficient
$g$
appearing in Chen & Keh (Reference Chen and Keh2003) is scaled differently.
5. Vesicle velocities
The rigid-body velocities of a freely suspended vesicle are obtained from (3.7), representing the force- and torque-free conditions. Substituting (4.24) and (4.38) into (3.7) gives
\begin{equation} \left (\begin{matrix} \mathcal{V}/\mathcal U\\ {a\unicode {x1D6FA}}/{\mathcal U} \end{matrix}\right ) \sim \mu (\delta )\,\delta ^{-1/2} F_0 {\unicode{x1D64D}}^{-1} \left (\begin{matrix} 1 \\ 1 \end{matrix}\right ). \end{equation}
The resistance coefficients in the limit
$\delta \ll 1$
were calculated by Goldman et al. (Reference Goldman, Cox and Brenner1967) using a lubrication analysis, where an ingenious argument eliminates the need for matching with the outer region (Jeffrey Reference Jeffrey and Torry1996):
\begin{equation} f^{\textit{tr}}/6\unicode{x03C0} \sim \frac {8}{15}\ln \frac {1}{\delta }, \quad h/6\unicode{x03C0} \sim -\frac {2}{15}\ln \frac {1}{\delta }, \quad g^{\textit{rot}}/6\unicode{x03C0} \sim \frac {8}{15}\ln \frac {1}{\delta }. \end{equation}
That is (recall (3.6)),
\begin{equation} {\unicode{x1D64D}} \sim \frac {4\unicode{x03C0} }{5}\ln \frac {1}{\delta } \left (\begin{matrix} 4 & -1 \\ -1 & 4 \end{matrix}\right ), \end{equation}
whose inverse is
\begin{equation} {\unicode{x1D64D}}^{-1} \sim \frac {1}{12\unicode{x03C0} \ln ({1}/{\delta })} \left (\begin{matrix} 4 & 1 \\ 1 & 4 \end{matrix}\right ). \end{equation}
Substituting (4.44) and (5.4) into (5.1) yields the requisite velocities
\begin{equation} \mathcal{V}/\mathcal U \sim \frac {5F_0 \delta ^{1/\sqrt {2}-1}}{12\unicode{x03C0} \ln ({1}/{\delta })} \quad \text{with} \quad \lim _{\delta \to 0}\frac {a\unicode {x1D6FA}}{\mathcal{V}}=1. \end{equation}
The ratio of angular velocity to rectilinear velocity effectively corresponds to pure rolling.
It is of interest to contrast these results with those pertaining to the classical problem of sedimentation of a rigid sphere due to an external force, say of magnitude
$\mathcal{F}$
, parallel to a rigid wall (Goldman et al. Reference Goldman, Cox and Brenner1967). In that problem, the appropriate velocity scale
$\mathcal U$
is simply
$\mathcal{F}/\eta a$
. Thus the counterpart of (5.1) is
\begin{equation} \left (\begin{matrix} \mathcal{V}/\mathcal U\\ {a\unicode {x1D6FA}}/{\mathcal U} \end{matrix}\right ) \sim {\unicode{x1D64D}}^{-1} \left (\begin{matrix} 1 \\ 0 \end{matrix}\right ). \end{equation}
Substituting (5.4) then gives
\begin{equation} \mathcal{V}/\mathcal U \sim \frac {1}{3\unicode{x03C0} \ln ({1}/{\delta })} \quad \text{with} \quad \lim _{\delta \to 0}\frac {a\unicode {x1D6FA}}{\mathcal{V}}=\frac {1}{4}. \end{equation}
6. Concluding remarks
Consider a 10 micron vesicle in aqueous solution (
$\eta = 10^{-2}\, \mbox{g} \, \mbox{cm}^{-1}\,\mbox{s}^{-1}$
) at room temperature (
$RT= 2.5 \times 10^{10}\,\mbox{g} \, \mbox{cm}^2 \,\mbox{s}^{-2} \,\mbox{mol}^{-1}$
). When exposed to a typical concentration gradient (
$\alpha \approx 10^{-3}\,\mbox{mol}\, \mbox{cm}^{-4}$
), we find from (1.1) and the estimates (1.2) that
\begin{equation} \mathcal U\approx 2.5 \times 10^{-4}\, \mbox{cm}\,\mbox{s}^{-1}. \end{equation}
Given the gross over-estimate in the value of
$k$
, osmophoretic velocities are rather moderate compared to those attained under similar concentration gradients due to other physicochemical mechanisms, such as diffusiophoresis (Anderson et al. Reference Anderson, Lowell and Prieve1982). The experimental observations of orders-of-magnitude discrepancy between theory and experiment (Nardi et al. Reference Nardi, Bruinsma and Sackmann1999) in conjunction with the numerical observations of mobility amplification (Chen & Keh Reference Chen and Keh2003) are therefore of practical interest.
In this paper, we employed the simplest linearised non-equilibrium description to analyse vesicle motion in close proximity to an impermeable wall. Our key results are the osmophoretic loads on a stationary vesicle, and in particular their anomalous power-law dependence upon the particle–wall separation. Upon making use of the well-known mobility relations of a rigid particle near a wall, we obtained an approximation for the rectilinear and angular velocities of the vesicle. At the leading order considered, these rigid-body velocities effectively satisfy the condition of pure rolling.
Admittedly, the derived approximation (5.5) for the rigid-body velocities is somewhat crude because of the ‘logarithmic error’ incurred by the use of (5.2). It is not difficult to improve (5.5) by incorporating the
${\textrm {ord}}(1)$
resistance corrections to the logarithmic terms (5.2). These corrections, obtained originally by Goldman et al. (Reference Goldman, Cox and Brenner1967) using comparison with exact solutions of the resistance problem (Dean & O’Neill Reference Dean and O’Neill1963; O’Neill Reference O’Neill1964), were later obtained using matched asymptotic expansions (O’Neill & Stewartson Reference O’Neill and Stewartson1967; Cooley & O’Neill Reference Cooley and O’Neill1968). In principle, then, we could achieve an ‘algebraic accuracy’ in the velocity approximations, with the error terms being smaller than some positive power of
$\delta$
. Since we have already achieved algebraic accuracy in the calculation of the osmophoretic loads, this may appear to be a natural course of action.
We have elected to avoid that refined calculation for two reasons. The first is the desire to maintain the focus on the fundamental analysis with the simplest possible model. The second reason is more pragmatic. While formally subdominant, the
${\textrm {ord}}(1)$
outer contributions to the osmophoretic loads from stresses on the vesicle scale could be numerically significant, given the mild singularity in (4.24). (With the rather marginal algebraic accuracy, we found the agreement with numerical computations of the exact problem, as illustrated in figure 2, pleasantly surprising.) The calculation of these outer contributions is a daunting task that we prefer to avoid. Consequently, we saw no point in introducing improved resistance expressions.
While our key result (5.5) does indicate a velocity blow-up as
$\delta \to 0$
, the associated divergence is rather mild. Decreasing the separation distance
$a\delta$
from 1 micron to 10 angstroms for a 10 micron vesicle, we find from (4.41) and (4.45) that the velocity pre-factor (5.5) merely increases from approximately 2.4 to approximately 4.5. This effect is too weak to explain the large (three orders of magnitude) discrepancy between the reported velocity measurements of Nardi et al. (Reference Nardi, Bruinsma and Sackmann1999) and the characteristic scale (1.1). It is important to emphasise, however, that the dramatic anomaly observed by Nardi et al. (Reference Nardi, Bruinsma and Sackmann1999) has not been reproduced by later experiments. Indeed, the controlled experiments of Gu et al. (Reference Gu, Tran, Lee, Zhang and Bishop2023) did not display any velocity amplification. Due to the scarcity of experiments, it is difficult to judge from the reported results whether the difference between the two experimenatal studies has to do with proximity to the wall.
Assuming an existing wall effect, it is desirable to improve the theoretical model. A natural refinement of the present model may involve the incorporation of solute advection. Employing a typical diffusivity in aqueous systems,
$D\approx 10^{-5}\,\mbox{cm}^{2}\,\mbox{s}^{-1}$
, and using the velocity scale (6.1) obtained for a 10 micron vesicle, the Péclet number
$a\,\mathcal{U}/D$
is approximately
$2.5\times 10^{-2}$
(Anderson Reference Anderson1983). It was observed by Anderson (Reference Anderson1983), however, that solute advection associated with the equilibrium solute concentration (both inside and outside the vesicle) may play a dominant role even for small Péclet numbers (‘a peculiar effect of convection’). Since the incorporation of the associated terms does not disturb the problem linearity, it may be possible to solve the associated near-contact problem. When these terms were incorporated into the numerical model of Chen & Keh (Reference Chen and Keh2003), they resulted (in certain cases – see table 2 in that paper) in the prediction of larger (albeit not by orders of magnitude) velocities than those predicted by the purely diffusive transport problem.
Other directions of interest involve a comparable two-dimensional analysis (Yariv Reference Yariv2016a ), where conformal-mapping methods may be able to handle the above-mentioned advection effects for arbitrary vesicle–wall separations. Finally, the motion of vesicles has received recent interest in the active matter community due to the possibility of self-generated concentration gradients (Peng et al. Reference Peng, Zhou and Brady2022). This suggests natural future extensions of the present contribution, following existing analyses of diffusiophoretic self-propulsion in close proximity to boundaries (Yariv Reference Yariv2016b ; Brandão Reference Brandão2024).
Declaration of interests.
The authors report no conflict of interest.
Open access
