PHYSICAL REVIEW E 101, 042601 (2020)

Dynamics of a membrane coupled to an active fluid
Chia-Chun Liang,1,* Kento Yasuda,2 Shigeyuki Komura ,2 Kuo-An Wu,1 and Hsuan-Yi Chen

3,4,†

1

2

Department of Physics, National Tsing Hua University, Hsinchu 30013, Taiwan
Department of Chemistry, Graduate School of Science, Tokyo Metropolitan University, 192-0397 Tokyo, Japan3
Department of Physics, National Central University, Jhongli 32001, Taiwan
4
Institute of Physics, Academia Sinica, Taipei 11529, Taiwan(Received 13 December 2019; accepted 11 March 2020; published 1 April 2020)

The dynamics of a membrane coupled to an active fluid on top of a substrate is considered theoretically. It is
assumed that the director field of the active fluid has rotational symmetry in the membrane plane. This situation
is likely to be relevant for in vitro reconstructed actomyosin-membrane system. Different from a membrane
coupled to a polar active fluid, this model predicts that only when the viscosity of the fluid above the membrane issufficiently large, a contractile active fluid is able to slow down the relaxation of the membrane for perturbations
with wavelength comparable to the thickness of the active fluid. Hence, our model predicts a finite-wavelength
instability in the limit of strong contractility, which is different from a membrane coupled to a polar active
fluid. However, a membrane coupled to an extensile active fluid is always unstable against long-wavelength
perturbations due to active extensile stress enhanced membrane undulation.
DOI: 10.1103/PhysRevE.101.042601

I. INTRODUCTION

Investigating the dynamics of membranes is important bothin the advance of cell biophysics [1,2] and soft condensed
matter physics [3]. An interesting example that bridges theequilibrium soft matter, nonequilibrium active matter, and biophysics is the flickering phenomena of red blood cells. It was
found that the fluctuation power spectrum of the membrane
depends on the viscosity of the solvent [4], and the fluctuationdissipation relation is violated due to active processes in
the red blood cells [5,6]. Membranes that contain active
inclusions are also nonequilibrium soft matter. A membrane
containing active pumps have fluctuation power spectrum that
depends on the active force dipoles exerted by the pumps[7], and instabilities can occur when the membrane curvatureenhances the activities of the pumps [8]. Furthermore, it has
been predicted that finite size domain can be induced due tothe active conformational transitions of the inclusions [9,10].
In red blood cells, cell membranes are coupled to a relatively regular spectrin network [1]. In most biological cells,
however, the cell membranes are coupled to the actomyosin
cytoskeleton [2], in which actin filaments are densely crosslinked and myosin motors exert contractile stress [11]. The
spatial organization of actin filaments and the way that contractile stress is coupled to the membranes are different fromred blood cells and membranes containing active inclusions.
To study how membrane dynamics are influenced by theout-of-equilibrium cytoskeleton, we propose in this article aminimal model in which the membrane is in contact with apassive simple fluid on one side, whereas the other side is*
Present address: Department of Earth System Science, Universityof California, Irvine, California 92697, USA.
†hschen@phy.ncu.edu.tw

2470-0045/2020/101(4)/042601(10)filled with an active fluid on a solid substrate. The directionof contractile stress is dictated by the orientation of the actin
filaments which we assume to be distributed isotropically inthe plane that is parallel to the membrane. Thus, our model describes a fluid membrane on top of an active fluid which could
be reconstructed in vitro [12], or seen in vivo but away fromthe relatively highly polarized lamellipodium. The generalityof this model also allows it to describe the surface dynamics ofa fluid containing either contractile or extensile nonswimmingrod-shaped active particles on a solid substrate [13]. However,
it is important to notice that the system under consideration in
our model is slightly different from previous studies of the
dynamics for active fluid thin film [14,15]. In Refs. [14,15],
active particles with polar order were considered and the focus
was the coupling effect between membrane dynamics to thedirector fluctuations in the membrane plane.
The main result of our calculation shows that, for contractile active fluids, active stress generally helps to stabilize aperturbed membrane. However, when the passive fluid is sufficiently viscous, the contractility in the active fluid can slow
down perturbations with wavelengths comparable to the thickness of the active fluid. This means that at sufficiently strong
contractility, the membrane can undergo a finite-wavelength
instability. Similar to many active polar or nematic systems[16], the origin of this instability is the coupling between the
splay deformation of the director field and the active stressgenerated flow. However, the fact that the above instability
occurs at wavelengths close to the thickness of the active
fluid marks a significant difference from the long-wavelength
instability predicted for a membrane on top of a thin active
polar fluid film [14,15], and the active Fréedericksz transition[17]. For extensile active fluids, we further find that thesystem is always unstable in the long-wavelength limit due
to the same coupling between director and flow field. The
membrane dynamics presented in our model is valid from

042601-1©2020 American Physical Society

LIANG, YASUDA, KOMURA, WU, AND CHEN

PHYSICAL REVIEW E 101, 042601 (2020)

cell, in which the organization of actin filaments is polarized[18]. If we replace the membrane by an interface between theactive and the passive fluids, then this model can describe a
fluid containing nonswimming rod-shaped active particles that
follow the plane of the interface but without in-plane nematicorder.
When the membrane is deformed, there is an elastic energy
associated with the deformation. The elastic free energy of themembrane can be expressed as
1(1)dxdy [γ (∇h)2 + κ (∇ 2 h)2 ],F=
2FIG. 1. Schematics of the system. A simple fluid of viscosity η1
is above the membrane, below the membrane is an active viscous
fluid with viscosity η2 and director field n. The height of themembrane is l. The dashed lines represent the director field of theactive fluid. The angle between local director field and the xy plane isdenoted as θ. A solid surface is located at the xy plane. For simplicity
we consider perturbation from flat membrane state with a single
wave vector qx̂, and the system is translationally invariant in the y
direction.

thin film limit, where lubrication approximation holds, to bulk
limit where wavelength is small compared to the thicknessof active fluid. Extension of our model to polar active fluid
and systems with polymerization/depolymerization dynamics
will be interesting.
This article is organized as follows. In Sec. II we present
our model from the geometry of the system, the elasticity ofthe director field and membrane, to the equations of motion
and boundary conditions. Section III presents the membranerelaxation rates calculated from our model. Also presented
and discussed are phase diagrams showing the parameter
range where the system becomes unstable. In Section IV wecompare our model with previous studies on similar systems[14,15,17], and discuss possible experimental realizations.
Some details on the calculation in our model and the solution
in the thin film limit are reported in the Appendices.
II. MEMBRANE ON TOP OF AN ACTIVE VISCOUS FLUID

Consider an active fluid on a solid substrate, the upper surface of the active fluid is a membrane whereas the upper sideof the membrane is filled with a simple fluid, as represented
in Fig. 1. In the steady state, the membrane is flat at z = l, in
the presence of fluctuations the membrane position becomesz(x, y, t ) = l + h(x, y, t ).
We assume that the director field of the active fluid isparallel to the membrane due to surface anchoring effect butotherwise isotropic. This means that the director field is in
the isotropic phase and thickness of the active fluid is notvery large. When this system is in the steady state, the tensor
ni n j  vanishes for all i and j except nx nx  = ny ny  = 1/2.
Here ni is the ith component of the unit vector n indicating
the local director direction, and ni n j  is the average of ni n jover a small region large compared to molecular size. Thisgeometry mimics the situation for an artificial membrane on
top of a reconstructed actin-myosin network above a solid
surface. However, it is different from the lamellipodium of awhere γ is the membrane tension and κ is the bending rigidity.
Without loss of generality, we consider a perturbation on
the membrane with wave vector qx̂, and assume translational
invariance in the y direction. Then the membrane deformation
h(x, y, t ) can be written ash(x, y, t ) = h(q, t )eiqx + c.c.,(2)where c.c. is the complex conjugate of h(q, t )eiqx . Let us
denote the angle between the local director field of the active
fluid and the xy plane by θ . The symmetry of the system
indicates that the elastic free energy of the director field is     2 k⊥∂θ 2∂θ
kz ∂θ 2
Fd = dxdydz
+
+
,
2
∂x
∂y
2 ∂z
(3)
where k⊥ and kz are elastic constants analogous to the thesplay and bending constants of nematic liquid crystal [19].
Notice that the symmetry of our system does not allow thesystem to acquire elastic energy against twist deformation.
For simplicity we assume that the director field relaxes much
faster than other slow modes. Use θ ∼ eiqx and ∂θ /∂y = 0,minimize Fd with respect to θ , we obtain the equation for theangle θ in the active fluid in the limit of small θ askz

d 2θ− k⊥ q2 θ = 0.
dz2

(4)

Let the director field on the boundaries be parallel to themembrane and the solid surface, in the limit of small membrane deformation, θ (z = l, t ) = iqh(q, t ) and θ (z = 0, t ) =0. Using one-constant approximation k⊥ = kz , we obtain thefollowing director field in the active fluidθ (z, t ) = iqh(q, t )sinh(qz).
sinh(ql )(5)

The actomyosin network is driven away from equilibriumby several active processes. First, actin filaments grow by recruiting free actin monomers to the plus ends of the filaments,
and this is balanced by actin depolymerization in the minus
ends of the filaments. Second, myosin motors utilize energy of
ATP hydrolysis to exert forces to the actin filaments and surrounding solvent. They provide the cytoskeleton a contractile
stress that depends on local orientation of actin filaments andmyosin density. In this article, we focus on the effect of active
contractility, while actin polymerization/depolymerization isneglected. Moreover, the density of actin network in thesystem is treated as a constant for simplicity.

042601-2

DYNAMICS OF A MEMBRANE COUPLED TO AN ACTIVE …

The local active stress in the active fluid is given by [11]σiaj = χ Q̃i j ,(6)where χ > 0 (χ < 0) for a contractile [20] (extensile [21])active fluid, and Q̃i j is the nematic order parameter, i.e., thetraceless part of ni n j . Although there is also an isotropic
active stress, it does not affect the dynamics of the systembecause of the incompressibility condition imposed on theactive fluid. In the limit of small deformation, the relevant
components of Q̃i j areQ̃xz = Q̃zx = θ /2.

(7)

The dynamics of the membrane is coupled to the surrounding fluids. For z > l, there exists a passive fluid with viscosity
η1 . The dynamics of this passive fluid is described by the
Stokes equation with incompressibility condition
η1 ∇ 2 v − ∇ p + f1 δ(z − l ) = 0,(8)

∇ · v = 0,

(9)

PHYSICAL REVIEW E 101, 042601 (2020)

comes from the active stress, f2u and f2d are the forces from themembrane and the supporting substrate, respectively. Similar
to the passive fluid above the membrane, we assume that theactive fluid is incompressible.
In our model, the membrane is treated as a twodimensional surface with zero thickness and the following
conditions need to be satisfied at this surface. First, the flow
field needs to be continuous,v(z = l + ) = v(z = l − ).

Second, the shear stress is continuous across the membrane
because the membrane is a two-dimensional fluid,σxz (z = l + ) = σxz (z = l − ).

η2 ∇ 2 v − ∇ p + f a + f2u δ(z − l ) + f2d δ(z) = 0,(10)

fia = χ ∂ j Q̃i j

(11)

where

(13)

Finally, the discontinuity of normal stress across membranebalances with the elastic force of the membrane,

andwhere v is the flow field, p is the pressure, and f1 is the force
exerted on the fluid by the membrane [3].
Between the membrane and the rigid surface, there exists
an active fluid with viscosity η2 . In this active fluid, actin
filaments are crosslinked thanks to the reversible bonds due tovarious linkers. This gives the network a mechanical responsethat is elastic on short time scale, but viscous on long timescale. We are interested in the dynamics on lengths greaterthan the typical molecular size in the active fluid (∼10 nm),and long time limit in which the cytoskeleton between themembrane and the rigid substrate can be modeled as a viscous
fluid under the influence of active forces. In this limit theactive fluid is described by the following equation:(12)−σzz (z = l + ) + σzz (z = l − ) = −

δF
.
δh

(14)

At the rigid surface z = 0, however, we assume the no-slip
condition for the flow field,v(z = 0) = 0.

(15)

III. LINEARIZED HYDRODYNAMICS CLOSE TO THE
STEADY STATE

With incompressibility condition, the velocity and pressurein the system can be expressed in terms of the membrane
conformation and the boundary forces f1 , f2u , and f2d . These
forces are further solved by imposing the boundary conditions
Eqs. (12)–(15). Time evolution of the membrane height is obtained from the kinematic boundary condition vz |z=l = ∂h/∂t.
Detailed intermediate steps of this calculation are presented inAppendix A.
The resulting equation for membrane dynamics is∂h(q, t )= −[λp (q) + λa (q)]h(q, t ),
∂t

(16)

where
λp (q) =(γ q + κq3 )[−1 + E + e4ql (1 + E ) − 2e2ql (E + 2ql + 2E q2 l 2 )]2η2 [(−1 + E )2 + e4ql (1 + E )2 − 2e2ql (−1 + E )(1 + 2q2 l 2 )](17)is the contribution from the restoring force provided by membrane elasticity. In the above, E = η1 /η2 specifies the viscosity
ratio of the passive fluid at z > l and the active fluid at z < l. However,χ
λa (q) =32η2 sinh(ql )[−q2 l 2 (−1 + E 2 ) + {cosh(ql ) + E sinh(ql )}2 ]×[(1 + 4q2 l 2 − 4E q3 l 3 ) sinh(ql ) − (E + 4ql + 4E q2 l 2 + q3 l 3 ) cosh(ql ) + sinh(3ql ) + E cosh(3ql )](18)

comes from the active stress associated with spatially nonuniform director field.
In the limit ql  1 the relaxation rate approaches
γ q + κq3,
2(η1 + η2 )χλa =.
8(η1 + η2 )λp =042601-3

(19)
(20)

LIANG, YASUDA, KOMURA, WU, AND CHEN

PHYSICAL REVIEW E 101, 042601 (2020)

FIG. 2. (a) The contribution of membrane elasticity to membrane relaxation rate for a tensionless membrane. In the limit of large ql it scales
as q3 l 3 . Curves for η1 /η2 = E = 0.75 (solid), 1.5 (long dashed), and 3.0 (short dashed) are shown. (b) Contribution of membrane elasticity
to membrane relaxation rate for a tensionless membrane in the long-wavelength limit. Solid curves: E = 0.75, dashed curves: E = 3.0. Thincurves represent the long-wavelength expansion presented in Eq. (22), thick curves represent the full solution presented in Eq. (17).

The passive contribution reduces to that for a membrane
separating two bulk fluids. The active contribution is independent of q, and it indicates that when ql  1 a contractility
assisted membrane relaxation is faster than a passive system.
However, for an extensile active fluid, active stress stretches
a curved membrane and prefers a curved membrane configuration, and hance the membrane relaxation is slower than apassive system [22]. Note that at sufficiently large q (for atensionless membrane this happens when q  (χ /κ )1/3 ), thepassive contribution dominates the membrane decay rate andthe membrane is basically passive.
In the limit ql  1, however, the relaxation rate isη1 4 4
γ q + κq3 q3 l 3−λp ≈q l + ··· ,
(21)
η2
3
2η2


η1
χ q2 l 2
1 − 3 ql + · · · .
(22)
λa ≈
12η2
η2
Here the leading contributions are independent of the viscosity of the passive fluid as expected from the lubrication
approximation [23], whereas the contributions from the passive fluid (η1 ) appear only in the next order in ql. It is worth
noticing that when ql is small, active stress dominates themembrane relaxation, λp < λa , and the system is stable (unstable) when the active fluid is contractile, χ > 0 (extensile,χ < 0). The mechanism for this long-wavelength instability
in extensile system is again due to active stretching induced
membrane undulation [24], which is also responsible for thesign of λa in Eq. (20).
The stability of the membrane in the long-wavelength
limit predicted by our model is different from a membrane
on top of a thin polar active fluid film. In the latter case,a long-wavelength instability due to the coupling between
the fluctuations of the director field in the membrane plane(the xy plane in our geometry) and membrane deformation
was predicted [14,15]. This instability can be stabilized by
the active force quadrupoles of the active fluid in the thinfilm limit [25]. In our model, such a coupling does not exist
because the system does not have any in-plane polar order. In

Refs. [14,15], a long-wavelength instability associated withsplay of the director field in the xz plane for a contractile
active fluid is also discussed, which is again different fromour prediction. To discuss this difference in details, we present
the lubrication approximation of our model in Appendix B by
keeping track of the order of magnitude of the sub-leading
terms. There we show that a consistent expansion at ql  1term by term leads to a stable (unstable) membrane in thelong-wavelength limit when the membrane is coupled to acontractile (extensile) active fluid.
To see the crossover from short-wavelength to longwavelength behavior, we plot in Figs. 2 and 3 the passive
and active contribution to the membrane relaxation rate asfunctions of ql for a tensionless membrane (γ = 0) withdifferent E . Figure 2 shows that the passive contribution
increases monotonically with q due to the bending free-energy
density κq4 /2. In the limit of large ql, this part scales as q3 l 3
as expected. In the limit of small ql, the long-wavelength
expansion Eq. (22) holds reasonably well as ql becomes
sufficiently small. Figure 3 shows that the active contribution
in the large ql limit also approaches the prediction provided
by Eq. (20). In the long-wavelength limit, the expansion
presented in Eq. (22) again holds better as ql becomes sufficiently small. For both passive and active contributions, thelong-wavelength expansion agrees better with the full solution
when E = η1 /η2 is smaller.
The active contribution also shows interesting features that
depend on the value of E , which can be clearly seen fromFig 3(a). For small E , the active contribution also increases
monotonically with q until it reaches the limiting magnitude provided by Eq. (20), indicating contractility-assisted
membrane relaxation at all wavelengths. However, as E becomes larger than Ec ≈ 1.5, contractility slows down relaxations for perturbations with ql ∼ 1. This indicates that eventhough contractility assists membrane relaxation in the shortand long-wavelength limits, there is a possible contractilityinduced membrane instability for an active fluid with sufficiently strong active stress when E > Ec . In contrast, fora membrane on top of an extensile active fluid, χ < 0, the042601-4

DYNAMICS OF A MEMBRANE COUPLED TO AN ACTIVE …

PHYSICAL REVIEW E 101, 042601 (2020)

FIG. 3. (a) Contribution of contractility to membrane relaxation rate. In the limit of large ql the relaxation rate approaches that provided
by the short-wavelength limit solution in Eq. (20). Solid curves: E = 0.75, dashed curves: E = 1.5, dotted curves: E = 3.0. Thick curves
represent the full solution presented in Eq. (18), thin lines represent the large ql limit presented in Eq. (20). Notice that when E > Ec ≈ 1.5
there is a region where λa becomes negative. (b) Contribution of contractility to membrane relaxation rate in the long wavelength limit. Thick
curves represent the full solution presented in Eq. (18), thin curves correspond to the results from the long-wavelength limit expansion Eq. (22).

system is always unstable in the long-wavelength limit regardless the magnitude of membrane tension, as can be seen from
Eq. (22).
The stability diagram for a membrane on top of an active
fluid can be plotted by introducing a dimensionless parameterχ l 3 /κ which characterizes the relative strength of the active
stress with respect to membrane elasticity. It can be seen from
Fig. 4(a) that when E < Ec , a system with a contractile active
fluid (χ > 0) is always stable. However, when the active fluid
is extensile (χ < 0), the membrane is always unstable in thelong-wavelength limit. Figure 4(b) shows that when E > Ec ,there is a finite-wavelength instability when the active fluid
is contractile (χ > 0), and the first unstable mode has ql oforder unity. When the active fluid is extensile (χ < 0), thesystem is always unstable in the short-wavelength limit, andthere is another finite-wavelength instability in the negative χregion when the activity is sufficiently strong. Two importantfeatures in E > Ec regime should be noticed. First, due tolimited range of ql that can be shown in the figure, it is noteasy to see that in the short-wavelength limit the system isalways stable because membrane elasticity is the dominating
drive of the relaxation dynamics of the system. Second, an
extensile system is stable in the wavelength regime wherea contractile system has a finite-wavelength instability. Thisis simply because λ changes sign when χ changes sign.
The complete phase diagram which summarizes the above
discussion of the stability of the system is shown in Fig. 5.
IV. DISCUSSION

The instability for a membrane on top of a contractile
active fluid is similar to the active Fréedericksz transition for
active nematics [17], because in both cases the instability isa result of the coupling between splay deformation and flow

FIG. 4. (a) For E < Ec (blue solid curve: E = 1.0, red dashed curve: E = 0.5, Ec ≈ 1.5), the system becomes unstable when the active
fluid is extensile and activity is sufficiently strong. Region below the curves is where the system becomes unstable. In the limit of small ql
the system is always unstable if the active fluid is extensile. (b) For E > Ec , (blue solid curve: E = 1.8, red dashed curve: E = 1.7), there isa finite-wavelength instability if the active fluid is contractile. If the active fluid is extensile, then in the small ql limit the system is always
unstable; with ql > 1 there is also an finite-wavelength instability.
042601-5

LIANG, YASUDA, KOMURA, WU, AND CHEN

PHYSICAL REVIEW E 101, 042601 (2020)

FIG. 5. Complete phase diagram for the system. The verticaldashed line indicates E = Ec ≈ 1.5. In regions S1 and S2, the systemis stable. In region LI1 and LI2, the system has a long-wavelength instability. In region F, there is a finite-wavelength instability. In region
LF, there is a long-wavelength instability and a finite-wavelength
instability.

field. However, in active Fréedericksz transition, the upperboundary does not deform and the maximum deformation of
the director field happens not at the upper boundary of theactive fluid but in the middle. When the transition occurs, thesystem develops a spontaneous flow that is independent ofx. Furthermore, for this transition to occur, the director fieldneeds to have a finite relaxation rate.
However, our model focuses on the limit of fast director
relaxation. With a membrane that can deform on the top, in
our model the maximum deformation of the director fieldoccurs on the upper boundary of the active fluid. Thus, theinstability in our model manifests itself as a membrane undulation instability, and it strongly depends on the viscosityof the passive fluid. By comparing Ref. [17] with our model,
we expect that for a system with a sufficiently viscous passive
fluid (E > Ec ) and relatively slow director dynamics, active
Fréedericksz transition should happen first; but for a systemwith E > Ec and fast director relaxation, one should firstobserve membrane undulation instability. For a system witha passive fluid that is less viscous (E < Ec ), active Fréedericksz transition may be observed, but not finite-wavelength
membrane undulation instability.
Previous studies of thin active fluid film applied lubrication approximation to obtain the coupled dynamics ofmembrane fluctuations and the deformation of polar order inthe membrane plane [14,15]. Recently, it has been pointed
out in Ref. [25] that the next order active force multipole can restore the stability of the system considered inRefs. [14,15]. Although these results are important progress
in two-dimensional active nematics, their main focus was the
dynamics of the director in the plane parallel to the membrane.
In this aspect, the focus of these studies is not related to ourmodel.
It is interesting to point out that in both Refs. [14,15], thepossible long-wavelength instability for director splay in thexz plane are also discussed. This is closely related to ourmodel, but their prediction is different from our result forql  1. To find out the origin of this difference, we present
in Appendix B a small-ql expansion and show that, for amembrane on top of a contractile active fluid, the system isstable in the long-wavelength limit if one consistently takes all
terms of lowest order in ql into account. The finite-wavelength
instability predicted by our model appears only when higher
order terms in small ql expansion are included. These terms
depend on the viscosity of the passive fluid and are clearly
beyond the regime where lubrication approximation holds.
This further shows that our general theory is a good approachto study the membrane dynamics coupled to an active fluid,
ranging from ql  1 to ql  1 limits.
In typical in vitro experiments, the viscosity of the passive fluid is at best as large as the viscosity of the active
counterpart. Our model predicts that in these systems themembrane relaxation is faster (slower) than a passive system
if the active fluid is contractile (extensile). To observe acontractility-induced finite-wavelength instability, one may
add passive polymers to the passive fluid to increase E . For
E ≈ 1.7–2.0, finite-wavelength instability can be observed forχ l 3 /κ  10–30.
The parameter χ has dimension of energy per unit volume,
since the contractility comes from the force dipoles of themolecular motors in the active fluid, the magnitude of χ isof the order of the typical magnitude of single force dipole
times the density of the motors. For a typical biomembraneκ ≈ 20–40 kB T [1], and the magnitude of a typical force
dipole generated by a molecular motor is around 10 kB T , thedimensionless parameter χ l 3 /κ has magnitude around c0 l 3
where c0 is the average number density of molecular motors.
Since we typically have c0−1/3  l, it is easy for the instability
to be observed.
We assume the director field to be isotropic in the membrane and substrate planes. Therefore, the active fluid is in an
isotropic phase in the bulk, and the orientational order is due
to surface anchoring. In real experiments, such spatial director
distribution can be achieved when the thickness of active fluid
is not very large compared to the nematic correlation length(roughly, a few times the “mesh size” of the network formedby the active particles). For example, in vitro actyomyosin
systems have mesh size at least on the order of ∼50 nm, andthe director field can be prepared to be parallel to the substrate
for active fluid with thickness up to ∼1 μm. For nonmotile
rod-shaped bacteria, the thickness of the active fluid can beas large as ∼10 μm. This means that, first, it is appropriate
to apply hydrodynamic theory to describe our model systems.
Furthermore, Figs. 2(a) and 3(a) show that, even when ql ∼ 5,
the large ql expressions Eqs. (19) and (20) already describe
the membrane dynamics well. Hence, our analysis, even thelarge ql limit, is experimentally accessible.
In summary, we have studied the dynamics of a membrane
above a finite-thickness active fluid beyond lubrication approximation. Our model predicts that, in the short-wavelength
limit (ql  1, q  (χ /κ )1/3 ), the membrane relaxation isdominated by the passive driving of membrane elasticity. For
longer wavelengths, the system is affected by active stress.
For extensile active fluid, the membrane is always unstablein the long-wavelength limit. For contractile active fluid, themembrane can have a finite-wavelength instability when thepassive fluid above the membrane is sufficiently viscous andthe strength of active stress exceeds a threshold that depends
on E . The long-wavelength instability predicted for extensile
active fluid and finite-wavelength instability predicted for the042601-6

DYNAMICS OF A MEMBRANE COUPLED TO AN ACTIVE …

PHYSICAL REVIEW E 101, 042601 (2020)

contractile active fluid are both within typical experimental
parameter range.
In this study we considered a relatively simple active fluidin which the density of active particles is uniform and the
director field is isotropic in the membrane plane. It will bevery interesting to extend current model, for example, by including the dynamics of director field, one can obtain a phase
diagram which includes the membrane undulation instabilities
and active Fréedericksz transition. By considering situations
where concentration variations in the system becomes important, especially when polymerization/depolymerization [18](or birth/death) of the active particles contribute significantly
to the dynamics of the system, one can explore the physics
of cytoskeleton (or biofilm). Finally, systems with different
membrane-director coupling, for example, the situation when
the director field is perpendicular to the membrane [26], andthe active fluid has finite thickness, is also interesting andworth studying in the future.

Science (JSPS). S.K. acknowledges support by a Grant-inAid for Scientific Research (C) (Grants No. 18K03567 and
No. 19K03765) from the JSPS. H.-Y.C. is supported by the
Ministry of Science and Technology, Taiwan (Grants No.
MOST 107-2112-M-008-024- and No. MOST 108-2112-M008-016) and K.-A.W. acknowledges the support of Ministryof Science and Technology, Taiwan (Grants No. MOST 1052112-M-007-031-MY3 and No. 108-2112-M-007-024-). H.Y.C. and K.-A.W. acknowledge the support from National
Center for Theoretical Sciences, Taiwan. The authors also
thank an anonymous referee for the helpful suggestions on
Eq. (3).

ACKNOWLEDGMENTS

H.-Y.C. thanks J.-F. Joanny for interesting discussions.
C.-C.L. thanks the hospitality and the financial support of
TMU from NTHU-TMU cotutorial program. K.Y. acknowledges support by a Grant-in-Aid for JSPS Fellows (Grant
No. 18J21231) from the Japan Society for the Promotion ofp(q, z, t ) =APPENDIX A: DERIVATION OF THE MEMBRANEEVOLUTION EQUATION

To obtain the membrane evolution Eqs. (16), (17), and(18), we first use the incompressibility condition to obtain thepressure. This is done by taking the Fourier transform of thedivergence of the momentum Eqs. (8) and (10). The Fouriertransform of the velocity field can be solved by substituting
the pressure back to the Fourier transform of the momentum
equations.
Next we take the inverse Fourier transform in the z direction to obtain the pressure and flow field in terms of theboundary forces and active force. The resulting expression forthe pressure in the region z > l ise−q(z−l )[−i f1x (q, t ) + f1z (q, t )],2(A1)

and that in the region z < l isp(q, z, t ) = −χ qh[sinh(qz) + e−ql sinh(q(l − z)) − e−qz qz − eqz q(l − z)]
4 sinh(ql )

e−q(l−z)e−qzud− i f2x− i f2x(q, t ) − f2zu (q, t ) +(q, t ) + f2zd (q, t ) .
2
2
The flow field in the region z > l is+e−q(z−l )[iq(l − z) f1z (q, t ) + (1 + q(l − z)) f1x (q, t )],4η1 q
e−q(z−l )[iq(l − z) f1x (q, t ) + (1 − q(l − z)) f1z (q, t )],vz (q, z, t ) =4η1 q

vx (q, z, t ) =

(A2)

(A3)
(A4)

and the flow field in the region z < l isvx (q, z, t ) =1 iχ qh(q, t ){(−1 − qz + q2 z2 ) sinh(qz) + (qz − q2 z2 ) cosh(qz)4η2 q 4 sinh(ql )− [eql (2q(l − z) − q2 (l − z)2 ) + e−ql (1 − q(l − z))] sinh(q(l − z))− [−eql (2q(l − z) − q2 (l − z)2 ) + e−ql q(l − z)] cosh(q(l − z))}e−qze−q(l−z)udiq(l − z) f2zu (q, t ) + (1 − q(l − z) f2x−iqz f2zd (q, t ) + (1 − qz) f2xu(q, t ) +d (q, t ) ,
4η2 q
4η2 q


χ qh(q, t )1{(1 + qz + q2 z2 ) sinh(qz) − (qz + q2 z2 ) cosh(qz)vz (q, z, t ) = −4η2 q 4 sinh(ql )+(A5)− [eql q2 (l − z)2 ) + e−ql (1 + q(l − z))] sinh(q(l − z)) + [eql q2 (l − z)2 + e−ql q(l − z)] cosh(q(l − z))}+e−qze−q(l−z)udiq(l − z) f2x− iqz f2x
(q, t ) + (1 + q(l − z) f2zu (q, t ) +(q, t ) + (1 + qz) f2zd (q, t ) .
4η2 q
4η2 q
042601-7(A6)

LIANG, YASUDA, KOMURA, WU, AND CHEN

PHYSICAL REVIEW E 101, 042601 (2020)

The boundary forces are solved by imposing the boundary conditions. From the continuity of flow field at the membrane, wehave
η2 f1x (q, t ) = η1iχ qh(q, t )[(−1 − ql + q2 l 2 ) sinh(ql ) + (ql − q2 l 2 ) cosh(ql )]4 sinh(ql )

u
d
+ f2x
(q, t ) + e−ql (1 − ql ) f2xd (q, t ) − iql f2zd (q, t )η2 f1z (q, t ) = −η1

,

(A7)

χ qh(q, t )[(1 + ql + q2 l 2 ) sinh(ql ) − (ql + q2 l 2 ) cosh(ql )]4 sinh(ql )d+ f2zu (q, t ) + e−ql (1 + ql ) f2zd (q, t ) − iql f2x(q, t ).

(A8)

From the no slip condition at z = 0, we havedf2x(q, t ) =iχ qh(q, t ){((1 − ql )e−ql + (2ql − q2 l 2 )eql ) sinh(ql ) + (qle−ql + (−2ql + q2 l 2 )eql ) cosh(ql )}4 sinh(ql )u(q, t )) ,− e−ql iql f2zu (q, t ) + (1 − ql ) f2x

f2zd (q, t ) =(A9)χ qh(q, t ) 2 2 ql[(q l e + qle−ql ) cosh(ql ) − (q2 l 2 eql + (1 + ql )e−ql ) sinh(ql )]4 sinh(ql )u(q, t ) .
− e−ql (1 + ql ) f2zu (q, t ) + iql f2x(A10)

The continuity of shear stress across the membrane leads to
f1x (q, t ) = −1 iχ qh(q, t )[(1 + ql − q2 l 2 ) sinh(ql ) − (ql − q2 l 2 ) cosh(ql )]4 4 sinh(q, t )

u
d
+ f2x
(q, t ) + e−ql iql f2zd (q, t ) + (−1 + ql ) f2x(q, t ),

(A11)

and the generalized Laplace condition for the normal stress across the membrane gives−

1 χ qh(q, t )δF=[(1 + ql + q2 l 2 ) sinh(ql ) − (ql + q2 l 2 ) cosh(ql )]δh(−q, t )2 4 sinh(ql )f2zu (q, t )f1z (q, t )1 + ql d−ql iql d+++ef (q, t ) −f2z (q, t ) .
2
2
2 2x
2

After solving Eqs. (A7), (A8), (A9), (A10), (A11), and (A12),
the kinematic condition ∂h/∂t = vz |z=l gives the evolution
Eq. (16) for membrane height.

By using ∂z2 ∼ 1/l 2  q2 ∼ ∂x2 , the momentum equation for
the active fluid can be approximated as
η2 ∂z2 vx − ∂x p +APPENDIX B: APPROXIMATE SOLUTION
OF MEMBRANE DYNAMICS AT ql  1

z < l.

(B1)

In this limit, the angle of the director field isθ (q, z, t ) = iqh(q, t )zsinh(qz)≈ iqh(q, t ) .
sinh(ql )liχ qh+ O(q2 l 2 ) = 0,
2l

(B3)

and

In this Appendix, we present the calculation for membrane dynamics in the limit ql  1. To be consistent and tocompare to the ql  1 limit provided by Eq. (22), we markthe magnitude of the sub-leading terms in each step of ourcalculation, and only keep the leading terms. This also helps
us to compare our result to the calculation based on lubrication
approximation [14,15].
Let the typical magnitude of vx in the active fluid be
U . Then the incompressibility condition gives the typicalmagnitude of vz in the active fluid asvz ∼ qlU,

(A12)

(B2)

−∂z p −

χ q2 hz + O(q2 l 2 ) = 0.
2l

(B4)

Here O(q2 l 2 ) means that terms that are dropped are of orderq2 l 2 or smaller compared to the leading terms.
At z = 0, the no slip condition applies. The continuity ofshear stress at z = l leads toη2 ∂z vx −iχ qh
z + O(E ql ) = 0,
2l

(B5)

where O(E ql ) means that the term that is dropped comesfrom shear stress of the passive fluid, it is of order ql compare to the leading terms. This is because the passive fluid
above the membrane satisfies ∂z vx ∼ qvx , and the continuity of flow field at z = x means that η1 (∂z vx + ∂x vz )|z=l + ∼E ql η2 (∂z vx + ∂x vz )|z=l− . Similar analysis for the normal

042601-8

DYNAMICS OF A MEMBRANE COUPLED TO AN ACTIVE …

PHYSICAL REVIEW E 101, 042601 (2020)

stress across the membrane leads to∂F= −p + O(q2 l 2 ).
(B6)
−
∂h
From the momentum equation, the typical magnitude of pressure in the active fluid satisfies qp ∼ η2U/l 2 , leading to p ∼η2U/(ql 2 ) in the active fluid. From the momentum equationin the passive fluid, the pressure in the passive fluid is of orderq2 l 2 compared to the pressure in the active fluid. The viscousstress in the active fluid is η2 ∂z vz ∼ η2U q. From incompressibility and continuity of flow field across the membrane, thetypical magnitude of flow in the passive fluid is vz ∼ vx ∼ U .
This gives us η1 ∂z vz ∼ η1 qU . Therefore, the viscous stress inboth active and passive fluids, like the pressure in the passive
fluid, is of order q2 l 2 compared to the pressure in the active
fluid.
To obtain the evolution equation for the membrane height,
we first integrate the z-component of the momentum equationfrom z to l and use the boundary condition for normal stress.
After that the pressure in the active fluid is solved as

two terms. Second, the “leading dropped term” is O(E ql ),which comes from the viscous stress of the passive fluid thathas been neglected in the shear stress continuity boundary
condition Eq. (B5).
Finally, using the kinematic boundary condition ∂t h =vz |z=l and incompressibility condition, we obtain the evolution equation of the membrane height in the ql  1 limit as l∂h(q, t )= vz (q, l, t ) =∂z vz (q, z, t )dz∂t0 lvx (q, z, t )dz
= −iq
0


3 3
q2 l 21
3 q l(γ q+κq )+χ
+O(E ql ) h(q, t ).
=−
η2
3
12
(B9)

Note that, first, the dominant contribution from the active
stress to vx comes from the active stress in the shear stresscontinuity condition Eq. (B5) and iχ qh/(2l ) term in the
momentum Eq. (B3). The contribution from active stress termin the pressure in Eq. (B7) is of order q2 l 2 compared to theseThis gives the leading terms in Eq. (22). The leading correction is of order ql compared to the leading term, and it comesfrom the shear stress of the passive fluid.
It is interesting to compare our approach with the lubrication approximation presented in Refs. [14,15]. AlthoughRefs. [14,15] considered active fluid with polar order andfluctuating active particle density, when the in-plane (xy)fluctuation of the director field and density fluctuations are
neglected, their model is basically identical to our model.
However, in Refs. [14,15], they neglected the contribution ofactive stress to the shear stress at z = l. This makes the leading
active term in the resulting membrane equation different from
Eq. (B8). Furthermore, according to our calculation, termsof order ql or smaller compared to the leading terms cancome from the passive fluid at z > l, or vz in the active fluid.
These are neglected in lubrication approximation. Therefore,
the subleading terms in the membrane dynamic equation ofRefs. [14,15] are also different from Eq. (22). Because of
these differences, they predicted a long-wavelength instabilityof the system against director splay in the xz plane, but oursystematic thin film calculation predicts that the system isstable in the long-wavelength limit.

[1] R. Lipowsky and E. Sackmann, The Structure and Dynamics ofMembranes (Elsevier, Amsterdam, 1995).
[2] H. Lodish et al., Molecular Cell Biology, 8th ed. (W. H.
Freeman, New York, 2016).
[3] U. Seifert, Adv. Phys. 46, 13 (1997).
[4] S. Tuvia, A. Almagor, A. Bitler, S. Levin, R. Korenstein,
and S. Yedgar, Cell membrane fluctuation are regulated
by medium macroviscosity: Evidence for a metabolic
driving force, Proc. Natl. Acad. Sci. USA 94, 5045
(1997)[5] H. Turlier, D. A. Fedosov, B. Audoly, T. Auth, N. S. Gov,
C. Sykes, J.-F. Joanny, G. Gompper, and T. Betz, Equilibrium
physics breakdown reveals the active nature of red blood cell
flickering, Nat. Phys. 12, 513 (2016).
[6] K. Yasuda, S. Komura, and R. Okamoto, Dynamics of a membrane interacting with an active wall, Phys. Rev. E 93, 052407(2016).

[7] J.-B. Manneville, P. Bassereau, S. Ramaswamy, and J. Prost,
Active membrane fluctuations studied by micropipette aspiration, Phys. Rev. E 64, 021908 (2001).
[8] S. Ramaswamy, J. Toner, and J. Prost, Nonequilibrium Fluctuations, Traveling Waves, and Instabilities in Active Membranes,
Phys. Rev. Lett. 84, 3494 (2000).
[9] H.-Y. Chen, The Internal States of Active Inclusions and the
Dynamics of an Active Membrane, Phys. Rev. Lett. 92, 168101(2004).
[10] H.-Y. Chen and A. S. Mikhailov, Dynamics of biomembranes
with active multiple-state inclusions, Phys. Rev. E 81, 031901(2010).
[11] F. Jülicher, K. Kruse, J. Prost, and J.-F. Joanny, Active behaviorof the cytoskeleton, Phys. Rep. 449, 3 (2007).
[12] F.-C. Tsai, B. Stuhrmann, and G. H. Koendrink, Encapsulationof active cytoskeletal protein networks in cell-sized liposomes,
Langmuir 27, 10061 (2011).

p(q, z, t ) =∂Fχ q2 (l 2 − z2 )+h(q, t ) + O(q2 l 2 ).
∂h(−q, t )
4
l
(B7)

Next, we substitute pressure into the x-component of the
momentum equation, integrate over z two times, use thecontinuity of shear stress at z = l and no slip condition atz = 0. Then vx in the active fluid becomes



∂F
χ qh(q, t )
z2
1
iq
−i
lz −
vx = −
η2∂h(q, t )2l2iχ qh(q, t )z + O(E ql ) .
(B8)
+
2

042601-9

LIANG, YASUDA, KOMURA, WU, AND CHEN

PHYSICAL REVIEW E 101, 042601 (2020)[13] W. F. Paxton, K. C. Kistler, C. C. Olmeda, A. Sen, S. K. S.
Angelo, Y. Cao, T. E. Mallouk, P. E. Lammert, and V. H.
Crespi, Catalytic nanomotors: Autonomous movement ofstriped nanorods, J. Am. Chem. Soc. 126, 13424 (2004).
[14] S. Sankararaman and S. Ramaswamy, Instabilities and Waves
in Thin Films of Living Fluids, Phys. Rev. Lett. 102, 118107(2009).
[15] N. Sarkar and A. Basu, Instabilities and diffusion in a hydrodynamic model of a fluid membrane coupled to a thin active fluidlayer, Eur. Phys. J. E 35, 115 (2012).
[16] R. A. Simha and S. Ramaswamy, Hydrodynamic Fluctuations
and Instabilities in Ordered Suspensions of Self-Propelled Particles, Phys. Rev. Lett. 89, 058101 (2002).
[17] R. Voituriez, J.-F. Joanny, and J. Prost, Spontaneous flow transition in active polar gels, Europhys. Lett. 70, 404 (2005).
[18] S. Marbach, A. L. Godeau, D. Riveline, J.-F. Joanny, and J.
Prost, Theoretical study of actin layers attachment and separation, Eur. Phys. J. E 38, 122 (2015).
[19] P. G. de Gennes and J. Prost, The Physics of Liquid Crystals,
2nd ed. (Clarendon Press, Oxford, 1993).
[20] For the molecular origin of contractility in actyomyosin system,
see, for example, D. Bray, Cell Movements, 2nd ed. (Talyor andFrancis, New York, 2001).

[21] For the experimental measurement of the active extensile force
dipole of swimming Escherichia coli, see, for example, S.
Chattopadhyay, R. Moldovan, C. Yeung, and X. L. Wu, Swimming efficiency of bacterium Escherichia coli, Proc. Natl. Acad.
Sci. USA 103, 13712 (2006).
[22] This effect is implicitly presented, for example, in Eq. (6a)of Ref. [24]. In Ref. [24], the dynamics of a membrane ina bulk active solvent was studied, and the active contribution to the membrane dynamics is qualitatively (but not completely) the same as our ql  1 limit, because their model,
different from our model, has a finite director relaxationrate.
[23] G. K. Batchelor, An Introduction to Fluid Dynamics (Cambridge
University Press, Cambridge, 1967).
[24] A. Maitra, P. Srivastava, M. Rao, and S. Ramaswamy, Activating Membranes, Phys. Rev. Lett. 112, 258101 (2014).
[25] A. Maitra, P. Srivastava, M. C. Marchetti, J. S. Lintuvuori,
S. Ramaswamy, and M. Lenz, A nonequilibrium force canstabilize 2D active nematics, Proc. Natl. Acad. Sci. USA 115,6934 (2018).
[26] N. Sarkar and A. Basu, Generic instabilities in a fluid membrane
coupled to a thin layer of ordered active polar fluid, Eur. Phys.
J. E 36, 86 (2013).

042601-10