CHAPTER FIVE

Dynamics of Heterogeneity
in Fluid Membranes
Shigeyuki Komura*1, Sanoop Ramachandran{, Kazuhiko Seki{,
Masayuki Imai}*Department of Chemistry, Graduate School of Science and Engineering, Tokyo Metropolitan University,
Tokyo, Japan{
Physique des Polymères, Université Libre de Bruxelles, Brussels, Belgium{
National Institute of Advanced Industrial Science and Technology (AIST), Tsukuba, Japan}
Department of Physics, Faculty of Science, Tohoku University, Sendai, Japan
1
Corresponding author: e-mail address: komura@tmu.ac.jp

Contents
1. Introduction
2. Membrane Hydrodynamics
3. Dynamics of Concentration Fluctuations
3.1 Time-dependent Ginzburg–Landau model
3.2 Effective diffusion coefficient
3.3 Membrane as a 2D microemulsion
3.4 Related works
4. Phase Separation Dynamics
4.1 Domain coarsening
4.2 Model and simulation technique
4.3 Domain growth dynamics
4.4 Correlated diffusion
4.5 Related works
5. Conclusion and Outlook
Acknowledgments
References

130
133
136
136
138
143
145
146
146
147
148
151
154
159
161
161

Abstract
In this chapter, we investigate the dynamics of heterogeneity observed in multicomponent fluid membranes, particularly focusing on hydrodynamic effects due tothe membrane and solvent. Two situations are discussed separately: above and below
the miscibility transition temperature. In the former case, we calculate the wave number
dependence of the effective diffusion coefficient by changing the temperature and/or
the thickness of the bulk fluid. We also obtain the effective diffusion coefficient when
multicomponent membranes are regarded as two-dimensional microemulsions. For
the case below the transition temperature, we study the domain growth exponent
in a binary membrane using a particle-based simulation method. With the addition

Advances in Planar Lipid Bilayers and Liposomes, Volume 16
ISSN 1554-4516
http://dx.doi.org/10.1016/B978-0-12-396534-9.00005-2# 2012 Elsevier Inc.

All rights reserved.

129

130

Shigeyuki Komura et al.

of bulk solvent, a change in the growth exponent from two- to three-dimensional nature
is observed. Along with the measurement of correlated diffusion, we conclude thatthe phase separation takes place through the Brownian coagulation process in our
simulation. We shall review some other works which are related to the present subject.

1. INTRODUCTION
Biological membranes typically contain various components such as
lipid mixtures, sterols, and proteins that are indispensable to cell functions[1,2]. Rather than being uniformly distributed in the membrane, there is
growing evidence that some cellular components are incorporated indomains arising from lateral lipid segregation in membranes. In 1997,
Simons and Ikonen proposed a hypothesis which suggested that some lipids
organize themselves into submicron-sized domains termed “lipid rafts” [3].
It was postulated that lipid rafts serve as platforms for proteins which
attribute certain functionality to each domain.
In 2006, the definition of lipid rafts was proposed as follows [4]: “Membrane rafts are small (10–200 nm), heterogeneous, highly dynamic, steroland sphingolipid-enriched domains that compartmentalize cellular processes.
Small rafts can sometimes be stabilized to form larger platforms throughprotein-protein and protein-lipid interactions.” In spite of such a definition,
the existence of lipid rafts is still debatable [5]. This is because lipid rafts have
not yet been directly observed in vivo. One of the key issues that has been
repeatedly asked is concerned with the lateral size and lifetime of lipid domains
or clusters. The reported value ranges from 20 to 200 nm for the size, and from10 2 to 103 s for the lifetime [6]. Recently, a high-resolution observation of thelipid dynamics in living cells is achieved using stimulated emission depletion
fluorescent microscopy [7]. They reported that sphingomyeline and membraneproteins are transiently (10–20 ms) trapped in about 20 nm diameter areas.
Stimulated by the lipid raft hypothesis, studies on artificial model membranes
have been developed intensively in the last decade. Especially, the phase behavior of ternary giant vesicles composed of saturated lipids with high chain-melting
temperature, unsaturated lipids with low chain-melting temperature, and cholesterol has been investigated for various mixtures [8,9]. It is known that such aternary membrane is homogeneously mixed in the high-temperature region. By
decreasing the temperature, the membrane undergoes a phase separation
between the coexisting liquid-ordered (Lo) and liquid-disordered (Ld)phases. As a typical example, we show in Fig. 5.1 the phase diagram of aternary vesicle consisting of DOPC (dioleoylphosphatidylcholine), PSM

131

Dynamics of Heterogeneity in Fluid Membranes

CholesterolA

47/23/30

B

3/57/40

One liquid

D

C
s
uid

iq

ol
Tw

68/22/10
DOPCdan
lid id
o
s iqu
l

F

E

30/60//10

27/53/20
PSM

8/72/20

Figure 5.1 Phases observed by fluorescence microscopy of giant unilamellar vesicles containing mixture of DOPC, PSM, and cholesterol at 298 K. The dark liquid
phase is rich in PSM and cholesterol, while the bright liquid phase is rich in DOPC.
Adapted from Ref. [10].

(N-palmitoyl-D-sphingomyelin) and cholesterol together with several
fluorescence micrographs [10]. In these pictures, dark liquid domains arerich in PSM and cholesterol, while bright liquid ones are rich in DOPC.
The minority liquid phase forms circular domains which undergo lateral
Brownian motion on the vesicle. In order to determine the ternary phase
diagrams and to identify the domain morphologies, a large amount ofresearch has been conducted using various experimental techniques [11].
On the other hand, behaviors of multicomponent membranes above thetransition temperature have also gained much attention. Performing deuterium NMR experiments, Veatch et al. made a notable attempt to investigate
critical fluctuations in lipid mixtures consisting of DOPC, DPPC
(dipalmitoylphosphatidylcholine), and cholesterol [12]. A broadening of
NMR resonances in the vicinity of the critical point was attributed to the compositional fluctuations on the scale less than 50 nm in the membrane. A morequantitative analysis of critical fluctuations using fluorescence microscopy wasaddressed by Honerkamp-Smith et al. for ternary mixtures of DPPC, diPhyPC(diphytanoylphosphatidylcholine), and cholesterol [13]. Typical microscope
pictures of concentration fluctuations are shown in Fig. 5.2 for higher temperatures. From the measurement of the critical exponents, the authors concluded
that the critical behavior in ternary membranes belongs to the universality class132

Shigeyuki Komura et al.

33.6 °C

32.8 °C

32.6 °C

32.4 °C

32.2 °C

32.0 °C

31.75 °C

31.0 °C

28.0 °C

Figure 5.2 Giant vesicle consisting of DPPC, diPhyPC, and cholesterol passing through acritical temperature at Tc  32.5  C. For temperatures above Tc, the vesicle exhibits concentration fluctuations. For temperatures below Tc, the vesicle undergoes a macroscopic phase separation. The scale bar corresponds to 20 mm. Adapted from Ref. [13].

of the two-dimensional (2D) Ising model [14]. Furthermore, by analyzing thecorrelation length and the line tension, it was shown that giant plasma membrane vesicles extracted from that of living rat basophil leukemia cells also exhibit a critical behavior [15]. These results allow us to speculate that lateralheterogeneity present in real cell membranes at physiological conditions could
correspond to critical fluctuations.
As mentioned in the first paragraph, it is important to note that lipid rafts
are highly dynamical objects. Being stimulated by various experiments, weshall discuss in this chapter the dynamics of heterogeneity in multicomponent membranes. It will be stressed that the hydrodynamic interaction mediated not only by the fluid membrane itself but also by the bulk
solvent plays an essential role in the dynamics of multicomponent membranes at large scales. We discuss the cases of above and below the miscibility
transition temperature separately.
In Section 2, we first provide a general framework of the membrane hydrodynamics. In Section 3, we investigate the dynamics of critical concentration fluctuations above the transition temperature. Based on the
Ginzburg–Landau approach with full hydrodynamics, we calculate the wave

133

Dynamics of Heterogeneity in Fluid Membranes

number dependence of the effective diffusion coefficient. We shall also consider the situation when the multicomponent membranes form microemulsions [16]. In Section 4, we study the domain growth dynamics below
the transition temperature using dissipative particle dynamics (DPD) simulations. We show that the presence of a bulk fluid will alter the domain growth
exponent from that of 2D to 3D, indicating a significant role played by themembrane–solvent coupling. In order to elucidate the underlying physicalmechanism of this effect, we look into the diffusion properties in the membrane by measuring two-particle correlated diffusion. We show that quasi2D phase separation proceeds by the Brownian coagulation (BC) mechanism
which reflects the 3D nature of the bulk solvent. We will also compare ourresults with the related works in the literatures.

2. MEMBRANE HYDRODYNAMICS
In this section, we first establish the governing equations for the fluid
membrane and its surrounding environment. Our aim is to derive the membrane mobility tensors which will be used in the later sections. More detailsof the calculation are given in Refs. [17,18].
As shown in Fig. 5.3, we assume that the membrane is an infinite planar
sheet of liquid, and its out-of-plane fluctuations are totally neglected, which
is justified for typical bending rigidities of bilayers. The liquid membrane isembedded in a bulk fluid such as water or solvent that is bounded by hardwalls. Such a situation is worth considering because biological membranesinteract strongly with other cells, substrates, or even the underlying cytoskeleton which can affect the structural and transport properties of the
z

Wall

h+
h+s
0
h–s

h

–h–

Wall

Figure 5.3 Schematic picture showing a planar liquid membrane having 2D viscosity 
located at z ¼ 0. It is sandwiched by a solvent of 3D viscosity 
s . Two impenetrable walls
are located at z ¼  h bounding the solvent.

134

Shigeyuki Komura et al.

membrane. Let v(r) be the 2D velocity of the membrane fluid and the 2D
vector r ¼ (x, y) represents a point in the plane of the membrane. We first
assume the membrane to be incompressiblerv ¼ 0,

ð5:1Þ

where r is a 2D differential operator. We work in the low-Reynolds number regime of the membrane hydrodynamics so that the inertial effects can beneglected. This allows us to use the 2D Stokes equation given by [19]
r2 v  rp þ f s þ F ¼ 0,

ð5:2Þwhere  is the 2D membrane viscosity, p(r) the 2D in-plane pressure, fs(r)the force exerted on the membrane by the surrounding fluid (“s” stands forthe solvent), and F(r) is any other force acting on the membrane which weshall discuss in later sections.
As presented in Fig. 5.3, the membrane is fixed in the xy-plane at z ¼ 0.
The upper (z > 0) and the lower (z < 0) fluid regions are denoted by “þ”
and “,” respectively. The velocities and pressures in these regions are written as v(r, z) and p(r, z), respectively. Since the 3D viscosity of the upper
and the lower solvent can be different, we denote them as s , respectively.
Consider the situation in which impenetrable walls are located at z ¼  h,
where hþ and h can be different in general. Similar to the liquid membrane,
the solvent in both regions are taken to be incompressible~  ¼ 0,
rv

ð5:3Þ

~ represents a 3D differential operator. We also neglect the solvent
where rinertia, and hence, it obeys the 3D Stokes equations~2  ~ s r v  rp ¼ 0:

ð5:4ÞThe presence of the surrounding solvent is important because it exerts force
on the fluid membrane. This force, indicated as fs in Eq. (5.2), is given by the
projection of ðsþ  s Þz¼0 ^ez on the xy-plane. Here ^ez is the unit vector
along the z-axis, and s are the stress tensors due to the solventh  T i~  þ rv
~
rv
:
ð5:5Þ
s ¼ p I þ s
In the above, I is the identity tensor and the superscript “T” indicates thetranspose.
Using the stick boundary conditions at z ¼ 0 and z ¼  h, we solve thehydrodynamic equations (5.3) and (5.4) to obtain fs. Then we calculate themembrane velocity from Eq. (5.2) asDynamics of Heterogeneity in Fluid Membranes

v½k ¼ G½kF½k,

135

ð5:6Þ

where v[k] and others are the Fourier components defined such as byðdkvðrÞ ¼v½kexpðikrÞ,
ð5:7Þ
ð2pÞ2


with k ¼ kx ,ky . After some calculations, one can show that the mobility
tensor G[k] in Fourier space is given by [17,18]


1
ka kb

 dab  2 , ð5:8Þ
Gab ½k ¼ 2
þkk þ k þs cothðkh Þ þ s cothðkh Þwith a, b ¼ x, y and k ¼ jkj. For simplicity, we consider the case when thetwo walls are located at equal distances from the membrane, that is,hþ ¼ h ¼ h. Then the above mobility tensor becomes


1
ka kb
ð5:9Þ
Gab ½k ¼ 2
dab  2 ,k
½k þ nk cothðkhÞ
where n 1 ¼ /2s with s ¼ (þs þ s )/2. An almost equivalent expression
to Eq. (5.9) has been derived for Langmuir monolayers in which there is
only one wall or a substrate [20,21]. In the following, we will employ
Eq. (5.9) as the general membrane mobility tensor.
We now discuss the two limiting situations of Eq. (5.9). Saffman andDelbrück (SD) investigated the case when the two walls are located infinitely
away from the membrane, which is called as the free membrane case [22–24].
Taking the limit of kh  1 in Eq. (5.9), the mobility tensor becomes [25,26]


1
ka kb
free
Gab ½k ¼
ð5:10Þ
dab  2 :kðk2 þ nkÞ

The quantity n 1 is called as the SD hydrodynamic screening length.
The real space expression of this mobility tensor is obtained by the Fouriertransform of Eq. (5.10) [21,25,26]
"
#
1
2
H
ð
nr
Þ
Y
ð
nr
Þ
1
1
free
Gab
ðrÞ ¼
H0 ðnr Þ  Y0 ðnr Þ þ 2 2 
þ
dab
4
pn r
nr
nr
"
#
1
4
2H1 ðnr Þ 2Y1 ðnr Þra rb
þ
 2 2þ

 H0 ðnr Þ þ Y0 ðnr Þ 2 ,
r
4
pn r
nr
nr
ð5:11Þ
where r ¼ jrj. In the above, Hn ðzÞ are Struve functions and Yn(z) are Bessel
functions of the second kind.

136

Shigeyuki Komura et al.

Evans and Sackmann (ES) considered the opposite kh 1 limit, forwhich the membrane is confined between the two walls although the solvent is still left [27]. In this case, Eq. (5.9) takes the following form:


1
ka kb
con
ð5:12Þ
dab  2 :Gab ½k ¼
kðk2 þ k2 Þ
pﬃﬃﬃﬃﬃﬃﬃﬃﬃ
In the above, k1 ¼ n1 h is called as the ES hydrodynamic screening
length, and it is the geometric mean of n 1 and h [28]. The above mobility
tensor was used in a phenomenological membrane hydrodynamic modelwith momentum decay [26,29–31]. The real space representation of theabove mobility tensor becomes"#1
K1 ðkr Þ
1
con
Gab ðrÞ ¼
K0 ðkr Þ þ
 2 2 dab
2p
kr
kr
"
#
ð5:13Þ
1
2K1 ðkr Þ
2 ra rbþ
K0 ðkr Þ þ 22 2 ,
2p
kr
kr r
where Kn(z) are the modified Bessel functions of the second kind. In
Section 3, we shall mainly use the general mobility tensor Eq. (5.9), whereas
either Eq. (5.10) or (5.12) is used in Section 4.

3. DYNAMICS OF CONCENTRATION FLUCTUATIONS
3.1. Time-dependent Ginzburg–Landau model
In this section, we use the idea of critical dynamics to calculate the effective
diffusion coefficient in multicomponent lipid membranes. Based on
Ginzburg–Landau approach with full hydrodynamics, we calculate the decay rateof the concentration fluctuations occurring in membranes [32,33]. We deal withthe general case where the membrane is surrounded by a bulk solvent and
two walls as depicted in Fig. 5.3. We also study the situation when themulticomponent membranes form 2D microemulsions [16]. This interesting
viewpoint is motivated by a recent work which predicts the reduction of theline tension in membranes containing saturated, unsaturated, and hybrid lipids(one tail saturated and the other unsaturated) [34–37]. We shall explore theconcentration fluctuations in 2D microemulsion with full hydrodynamics.
Consider a two-component fluid membrane composed of lipid A andlipid B whose local area fractions are denoted by fA ðrÞ and fB ðrÞ, respectively. Since the relation fA ðrÞ þ fB ðrÞ ¼ 1 holds, we introduce a new

137

Dynamics of Heterogeneity in Fluid Membranes

variable defined by cðrÞ ¼ fA ðrÞ  fB ðrÞ. Then the simplest form of thefree-energy functional F fcg describing the fluctuation around the homogeneous state is
ð h
i
a
c
F fcg ¼ dr c2 þ ðrcÞ2 ,ð5:14Þ
22
where a > 0 is proportional to the temperature difference from the criticaltemperature Tc and c > 0 is related to the line tension.
The time evolution of concentration in the presence of hydrodynamic
flow is given by the time-dependent Ginzburg–Landau equation for aconserved order parameter [38]@cdFþ rðvcÞ ¼ Lr2 ,@tdc

ð5:15Þ

where L is the kinetic coefficient. In the membrane hydrodynamic equation (5.2), we need to incorporate the thermodynamic force due to the concentration fluctuations. Hence we have
F ¼ crdF:
dc

ð5:16Þ

We implicitly assume that the relaxation of the velocity v is much fasterthan that of concentration c. The membrane velocity can be formally solved
as follows using the appropriate 2D mobility tensor Gab ðr, r0 Þ derived in theprevious sectionðdF
va ðr, tÞ ¼ dr0 Gab ðr, r0 Þ r0b c:
ð5:17Þ
dcðr0 ÞSince our interest is in the concentration fluctuations around the homogeneous where the bar indicates the spatialstate, we define dcðr, tÞ ¼ cðr, tÞ  c,average. The free-energy functional expanded in powers of dc becomes
ð h
i
a
c
ð5:18Þ
F fdcg ¼ dr ðdcÞ2 þ ðrdcÞ2 :2
2
Substituting Eq. (5.17) into Eq. (5.15), we get@dcðr,t Þ
dF
¼ Lr2
@t
dðdcÞ
ð


 dr0 ðra dcðrÞÞGab ðr,r0 Þ r0 b dcðr0 Þ

dF:dðdcðr0 ÞÞð5:19Þ

138

Shigeyuki Komura et al.

We now consider the dynamics of the time-correlation function defined bySðr, tÞ ¼ hdcðr1 , tÞdcðr2 , 0Þi,

ð5:20Þ

where r ¼ r2  r1 . Within the factorization approximation [39], the spatial
Fourier transform of Sðr, tÞ defined byð
S½k,t  ¼ drSðr, tÞ expðikrÞ
ð5:21Þ
satisfies the following equation@S½k, t 
¼  Gð1Þ ½k þ Gð2Þ ½k S½k, t :@t

ð5:22Þ

The first term, Gð1Þ ½k, denotes the van Hove part of the relaxation rate
given byGð1Þ ½k ¼ LkB Tk2 w1 ½k:

ð5:23Þ

Here, the static correlation function is defined bykB T,w½k ¼ hdc½kdc½ki ¼  2
c k þ x2

ð5:24Þ

where x ¼ (c/a)1/2 is the correlation length, kB the Boltzmann constant, andT the temperature.
As for the second term in Eq. (5.22), Gð2Þ ½k denotes the hydrodynamic
part of the decay rate given by
ð
1
dq
ð2Þ
G ½k ¼
ka Gab ½qkb w½k þ q:ð5:25Þ
w½k ð2pÞ2
When we use Eq. (5.9) for the mobility tensor Gab, the hydrodynamic partof the decay rate is expressed with an integral as
ð
kB T
dq
w½q
k2 q2  ðkqÞ2
ð2Þ
G ½k ¼:
w½k ð2pÞ2 jk  qj2 þ njk  qj cothðjk  qjhÞ jk  qj2
ð5:26Þ

3.2. Effective diffusion coefficient
We now introduce an effective diffusion coefficient (due only to the hydrodynamic part) D[k] defined by

139

Dynamics of Heterogeneity in Fluid Membranes

Gð2Þ ½k ¼ k2 D½k:ð5:27Þ

In order to deal with dimensionless quantities, we rescale all the lengths by
the SD hydrodynamic screening length n 1 ¼ /2s such that K ¼ k/n,
Q ¼ q/n, X ¼ xn, and H ¼ hn. Then D[k] can be rewritten as
ð
ð 2p
 1
kB T 
D½K;X, H  ¼ 2 1 þ K 2 X 2
dQ
dy
4p 
0
0
ð5:28Þ
Q3 sin2 ypﬃﬃﬃﬃ  ,
ð1 þ Q2 X 2 Þ G2 þ G3=2 coth GH
with G ¼ K2 þ Q2  2KQ cos y. Since this integral cannot be performed
analytically, we evaluate it via a numerical method. We explore the dependencies of D on the variable K, and the parameters X and H. Notice thatthe dimensionless correlation length X also measures the proximity of thetemperature with respect to the critical temperature Tc.
In Fig. 5.4, we plot the diffusion coefficient D (scaled by kBT/4p) as afunction of dimensionless wave number K for a different solvent thickness H
while the correlation length is fixed to X ¼ 1 (i.e., fixed temperature). In thelimit of K 1, D is almost a constant. The calculated D starts to increase101

4phD/kBT

100

H = 100

H=1

10–1

H = 0.01
X=1

10–2
10–3

10–2

10–1

100
K

101

102

103

Figure 5.4 Scaled effective diffusion coefficient D as a function of K for H ¼ 0.01, 1, and100 when X ¼ 1. The solid lines are from the analytical expression given in Eq. (5.29)
obtained in the limit of small H.

140

Shigeyuki Komura et al.

around K  1, and a logarithmic behavior (extracted via numerical fitting) isseen for K  1. In this plot, we see that D becomes smaller for smaller valuesof H. Figure 5.5 shows the diffusion coefficient D as a function of wave
number K for different X (i.e., different temperatures) while the solvent
height is fixed to H ¼ 1. Here D is nearly constant for K 1/X and follows
an S-shaped curve with increasing K. Finally, a logarithmic dependence isobserved for large K. We note that this logarithmic behavior for K  1/X isin contrast to that of 3D critical fluids given by the Kawasaki function which
increases linearly with the wave number [40].
In Fig. 5.6, we explore the effect of the correlation length X on D fordifferent values of H when K ¼ 10 3. The quantity X is a measure of an effective size of the correlated region formed transiently in the membrane dueto thermal fluctuations. When X 1, the diffusion coefficient D decreases
only logarithmically, which is typical for a pure 2D system [22–24]. The
proximity to the walls results in a loss of momentum from the membrane[29]. This leads to a rapid suppression of the velocity field within themembrane. Consequently, the values of D are smaller for smaller H. The
flattening of the curves for X  1 is due to the dominance of the X2terms in the numerator and denominator in Eq. (5.28). In Fig. 5.7,

101

100

X = 0.01

4phD/kBT

X=1
10–1

10–2

10–3
X = 100
10–4
10–3

10–2

10–1

H=1
100
K

101

102

103

Figure 5.5 Scaled effective diffusion coefficient D as a function of K for X ¼ 0.01, 1, and100 when H ¼ 1.

141

Dynamics of Heterogeneity in Fluid Membranes

101
100
10–1

H = 100

4phD/kBT

10–2
10–3
10–4

H=1

10–5
H = 0.01

10–6
10–7

K = 0.001

10–8 – 3
10

10 – 2

10 – 1

10 0

10 1

10 2

10 3

10 4

X

Figure 5.6 Scaled effective diffusion coefficient D as a function of X for H ¼ 0.01, 1, and100 when K ¼ 10 3. The solid lines are from the analytical expression given in Eq. (5.29)
obtained in the limit of small H.

102
101

X = 0.01

100

4phD/kBT

10–1
X=1

10–2
10–3
10–4

X = 100

10–5
10–6
10

K = 0.001

–7

10–3

10–2

10–1

100

101

102

103

H

Figure 5.7 Scaled effective diffusion coefficient D as a function of H for X ¼ 0.01, 1, and100 when K ¼ 10 3.

142

Shigeyuki Komura et al.

we plot D as a function of H for different values of X when K ¼ 10 3. Herethe values of D are smaller for larger X, as it should be. In general, there is amonotonic increase of D followed by a saturation to a constant value forlarger H. It is interesting to note that the crossover occurs when H  X holds.
For confined membranes (small h), we showed in Section 2 that the generalmobility tensor Eq. (5.9) reduces to Gcon
ab given by Eq. (5.12). In this case, one
can obtain an analytical expression for the effective diffusion coefficient [32]. In
terms of the dimensionless quantities K, X, and H, it is written as
2
0
1
2 2
kB T 1 þ K X 4
X
D½K;X, H  ¼
ln@pﬃﬃﬃﬃﬃA
4p 2K 2 X 2
H
0
1


H
X
A
ð5:29Þ
þ 2 1 þ K 2 X 2 ln@pﬃﬃﬃﬃﬃ
XH ð1 þ K 2 X 2 Þ!
Kþ4 þ K2 þ Kþ2 OHO
þ 2 ln,O  K2  1
2X
where
O¼

qﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
ðK 2 X 2 þ X 2 =H  1Þ2 þ 4K 2 X 2 ,ð5:30Þ

pﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
K 2 X 2  X 2 =H :

ð5:31Þ

andK ¼

Equation (5.29) is plotted using solid lines in Fig. 5.4 for H ¼ 0.01 and 1 with
X ¼ 1. For H ¼ 0.01, the analytical and numerical data coincide, giving
credence to accuracy of the numerical solutions. It is seen that even for
H ¼ 1, the agreement is still acceptable. For H ¼ 100, however, a significant
deviation is observed (not shown), which is expected as this limit is beyond
the valid range of Eq. (5.29).
The solid lines in Fig. 5.6 also represent the analytical result of Eq. (5.29).
It is seen that the analytical and the numerical data points almost coincide for
H ¼ 0.01 and 1. The agreement between the numerical result and the analytical expression is beyond the expected range of H 1 and reaches up to
H  1, as pointed out by Stone and Ajdari [28]. Hence Eq. (5.29) is useful inanalyzing the experimental data in many situations.
When the critical temperature Tc is approached from above, the correlation length diverges according to x  jT  Tc jn , where n is the critical

Dynamics of Heterogeneity in Fluid Membranes

143

exponent. When Tc is approached from below, on the other hand, the orderparameter given by the difference in lipid compositions vanishes as
dc  ðTc  T Þb . From the experiments on model multicomponent vesicles, ¼ 1=8 and n ¼ 1,
these critical exponents were found to have values close to b
respectively [13]. Furthermore, the experiment on giant plasma membrane
vesicles measured the critical exponent g ¼ 7=4 which characterizes the criticalbehavior of the osmotic compressibility [15]. These static exponents seem tocoincide with the exact results of the 2D Ising model [14]. The description
presented in this section uses a mean-field approach, and therefore, the ¼ 1=2, n ¼ 1=2, and g ¼ 1, respectively. The
corresponding exponents are b
discrepancies between these values are still under discussion.

3.3. Membrane as a 2D microemulsion
The role of surfactant molecules in 3D microemulsions is to reduce the surface tension at the interface between oil and water. In an analogy to 3D
microemulsions, hybrid lipids (one chain unsaturated and the other saturated)
act as lineactant molecules which stabilize finite-sized domains in 2D. In otherwords, hybrid lipids play a similar role to surfactant molecules at the interface
between Lo and Ld domains. It should be also noticed that hybrid lipids form amajor percentage of all naturally existing lipids [41]. Based on a simple model ofhybrid lipids, Brewster et al. showed that finite-sized domains can be formed inequilibrium [34,35]. A subsequent model predicted that domains are evenmore stabilized in a system of saturated/hybrid/cholesterol lipid membranes[36,37]. Being motivated by this idea, we calculate the decay rate ofconcentration fluctuations when the free energy of the multicomponent
membrane has the form of a 2D microemulsion.
The free-energy functional for a microemulsion includes a higher orderderivative term and is expressed in terms of dc as [16]
ð h
2 i
a
c
g
F ME fdcg ¼ dr ðdcÞ2 þ ðrdcÞ2 þ r2 dc ,
ð5:32Þ
2
2
2
with a, g > 0 and c < 0. The negative value of c creates 2D interfaces, while
the term with positive g is a stabilizing term. This form of the free energy has
been used previously to study coupled modulated bilayers [42,43]. As in theprevious section, the decay rate of the correlation function can be split into
two parts. First, the van Hove part now becomesð1Þ

GME ½k ¼ LkB Tk2 w1
ME ½k,ð5:33Þ

144

Shigeyuki Komura et al.

where the static correlation function wME ½k is [44]wME ½k ¼

kB T:
4
gk þ ck2 þ a

ð5:34Þ

By defining
c
,
2g
 2
a
c
4
,
s ¼ 
g
2g
k20 ¼ 

ð5:35Þ
ð5:36Þ

we can write the static correlation function askB Ti:wME ½k ¼ h
2
g ðk2  k20 Þ þ s4

ð5:37Þ

On plotting wME as a function of k, a peak appears at k ¼ k0 followed by
a 1/k4 decay. The width of the peak is given by s. A lamellar phase appears
when s ¼ 0. Notice that c ¼ 0 is called the Lifshitz point at which the peakoccurs for k ¼ 0. Using the form of Eq. (5.37), we can write Eq. (5.33) as
h
i
2
ð1Þ
GME ½k ¼ Lgk2 k2  k20 þ s4 :
ð5:38Þ
Similar to the previous section, we next write the hydrodynamic part ofthe decay rate in terms of the effective diffusion coefficient DME ½k asð2Þ

GME ½k ¼ k2 DME ½k:

ð5:39Þ

Using the mobility tensor given by Eq. (5.9), we can write DME as
DME ½K;K0 ,S, H  ¼
ð 2p

ð1
dQ
0

i

kB T h 2
2 2
4

K
þ
S
K
0
4p2 

Q3 sin2 yi
dy hpﬃﬃﬃﬃ  ,20
ðQ2  K02 Þ þ S4 G2 þ G3=2 coth GH

ð5:40Þ

where K ¼ q/n, K0 ¼ q0/n, Q ¼ q/n, S ¼ s/n, H ¼ hn, and G ¼ K2 þ Q2  2
KQ cos y.
In Fig. 5.8, we plot DME as a function of K for different values of H when
K0 ¼ S ¼ 1 are fixed. When K 1, D shows a constant value. We also
observe the 2D characteristic of logarithmic behavior of D for K  1.

145

Dynamics of Heterogeneity in Fluid Membranes

102

101

4phDME/kBT

H = 100
100
H=1
10–1

10–2
H = 0.01
10–3

10–3

10–2

10–1

K 0= S = 1
100
K

101

102

103

Figure 5.8 Scaled effective diffusion coefficient DME as a function of K for H ¼ 0.01, 1,
and 100 when K0 ¼ S ¼ 1.

For K 1, the effect of the outer environment is felt with the suppression of
the diffusion coefficient with smaller H. The curves almost overlap whenK  1, indicating the negligible effect of the outer environment at large
wave numbers. An interesting feature of DME is the dip occurring at
K  K0 which does not exist for binary critical fluids. This can be attributed
to the peak at k ¼ k0 in wME ½k [39,44]. We also note that, for 3D
microemulsions, the effective diffusion coefficient varies linearly with thewave number when it is large enough [39], which is in contrast to thepresent 2D microemulsions.

3.4. Related works
In addition to the present theory, there are some other theoretical works
on concentration fluctuations in multicomponent membranes. Using
renormalization group techniques, Tserkovnyak and Nelson calculated protein diffusion in a multicomponent membrane close to a rigid substrate [45].
They pointed out that, in the vicinity of the critical point, the effective
protein diffusion coefficient acquires a power-law behavior. Inaura andFujitani first discussed concentration fluctuations in free membranes surrounded by a 3D solvent [46]. They used the mobility tensor Gfree
ab in
Eq. (5.10) and calculated numerically the effective diffusion coefficient.

146

Shigeyuki Komura et al.

Hence their calculation corresponds to the special case of the present theory.
Haataja also discussed critical dynamics in multicomponent lipid membranes
and showed that the effective diffusion coefficient exhibits a crossover from alogarithmic behavior to an algebraic dependence (inversely proportional tothe correlation length) for larger length scales [47]. However, since an approximate empirical relation for the diffusion coefficient of a moving object
was employed [48], his theory should be distinguished from ours. Concentration fluctuations in membranes in the absence of hydrodynamic effect wasrecently discussed by McConnell [49].

4. PHASE SEPARATION DYNAMICS
4.1. Domain coarsening
Although biomembranes can be regarded as 2D viscous fluids, they are not isolated pure 2D systems since lipids are coupled to the adjacent solvent. Hence it isof great interest to investigate the phase separation dynamics in such a quasi-2D
liquid membrane in the presence of hydrodynamic interaction. (We use theword “quasi-2D” whenever the membrane is coupled to the bulk fluid.) To
address this problem, we consider a 2D binary viscous membrane in contact
with a bulk solvent. We employ a simple model in which the membrane is confined to a plane with the bulk fluid particles added above and below. In ourmodel using DPD simulation technique, the exchange of momentum between
the membrane and the bulk solvent is naturally taken into account. We particularly focus on the effect of bulk solvent on the quasi-2D phase separation.
Before explaining our simulation, let us briefly review here the generalknowledge about phase separation of binary fluids following a quench [50].
The dynamic scaling hypothesis assumes that there exists a scaling regime
characterized by the average domain size R that grows with time t as R ta
with an universal exponent a. For 3D off-critical binary fluids, there is an
initial growth by the BC process [51], followed by the Lifshitz–Slyozov(LS) evaporation–condensation process [52]; both mechanisms show a
growth exponent a ¼ 1/3. For critical mixtures, there is an intermediatea¼ 1 regime owing to hydrodynamic flow effect [53]. This is followed bya late time inertial regime of a ¼ 2/3 [54]. The scenario is slightly different
for pure 2D systems [55]. For an off-critical mixture, it was predicted thatafter the initial formation of domains, they grow by the BC mechanism witha different exponent a ¼ 1/2 (as will be explained later), followed by a crossover to the LS mechanism which gives a ¼ 1/3 even in 2D. For critical mixtures, on the other hand, the initial quench produces an interconnectedDynamics of Heterogeneity in Fluid Membranes

147

structure which coarsens and then breaks up due to the interface diffusion
with an exponent a ¼ 1/2. After the breakup processes, coarsening takes
place through BC that is again characterized by the a ¼ 1/2 scaling [51].
These predictions were confirmed by molecular dynamics simulations in
2D [56]. The exponent a ¼ 1/2 was also observed in 2D lattice-Boltzmann
simulations in the presence of thermal noise for a critical mixture [57].

4.2. Model and simulation technique
We use a structureless model of the 2D fluid membrane within the DPD
framework [58,59]. As shown in Fig. 5.9, the 2D membrane is represented
by a single layer of particles confined to a plane. In order to study phaseseparation, we introduce two species of particles, A and B. The bulk fluidwhich we call as “solvent” (S) is also represented by single particles of same
size as that of the membrane particles. All particles have the same mass m.
In DPD, the interaction between any two particles, within a range r0, is
linearly repulsive. The pairwise interaction leads to full momentum conservation, which in turn brings out the correct fluid hydrodynamics. The force
on a particle i is given by i
dvi X h C  RmþFð5:41Þ
¼
Fij rij þ FD
r
,v
ij ij
ij
ij rij ,dt
j6¼i
where rij ¼ ri  rj and vij ¼ vi  vj . Of the three types of forces acting ontherij ,particles, the conservative force on particle i due to j is FC
ij ¼ aij o rij ^Figure 5.9 Image of the fluid membrane with the bulk fluid called solvent. The yellow(A) and red (B) particles represent the two components constituting the membrane,
while blue ones (S) represent the solvent. For clarity, only a fraction of the solvent particles are shown.

148

Shigeyuki Komura et al.

where aij is an interaction strength and ^rij ¼ rij =rij with rij ¼ jr
 ij j. The second
D2
type of force is the dissipative force Fij ¼ Gij o rij ^rij vij ^rij , where Gij isthe dissipative strength
for the pair (i, j). The last is the random force1=2  
z
FR
¼
s
ð
Dt
Þ
o
r
^
r
ij
ij ij ij , where sij is the amplitude of the random noiseij
for the pair (i, j), and zij is a random variable with zero mean and unit variance which is uncorrelated for different pairs of particles and different time
steps. The dissipative and random forces act as a thermostat, provided thefluctuation–dissipation theorem s2ij ¼ 2GijkBT is satisfied. The weight factor
is chosen as o(rij) ¼ 1  rij/r0 up to the cutoff radius r0 and zero thereafter.
The particle trajectories are obtained by solving Eq. (5.41) using the velocity-Verlet integrator. In the simulation, r0 and m set the scales for length andmass, respectively, while kBT sets the energy scale. The time is measured
in units of t ¼ (mr20/kBT)1/2. The numerical value of the amplitude of
the random force is assumed to be the same for all pairs such that sij ¼ 3.0[(kBT)3m/r20]1/4, and the fluid density is set as r ¼ 3.0. We set kBT ¼ 1 andthe integration time step is chosen to be Dt ¼ 0.01t.
The membrane is constructed by placing particles in the xy-plane in the
middle of the simulation box (see Fig. 5.9). Owing to the structureless representation of the constituent particles, we apply an external potential so as
to maintain the membrane integrity. This is done by fixing the z-coordinatesof all the membrane particles. The work involves the systematic variation of
the height of the simulation box starting from the pure 2D case. In the absenceof solvent, we work with a 2D-box of dimensions Lx Ly ¼ 80 80 with
19,200 particles constituting the membrane. For the quasi-2D studies, we
add solvent particles S above and below the membrane and increase the heightof the box as Lz ¼ 5, 20 and 40. For all the cases, there are 19,200 membrane
particles. The largest box size (Lz ¼ 40) has 748,800 solvent particles. The box
with height Lz ¼ 40 is found to be sufficiently large enough to prevent thefinite-size effect which affects the membrane–solvent interaction. The systemis then subject to periodic boundary conditions in all the three directions. For
phase separation simulations, we introduce two species of membrane particles
A and B. The interaction parameter between various particles are given by
aAA ¼ aBB ¼ aSS ¼ aAS ¼ aBS ¼ 25 and aAB ¼ 50. In order to do a quench, themembrane is first equilibrated with a single component, following which afraction of the particles are instantaneously changed to the B type.

4.3. Domain growth dynamics
First, we describe the results of the phase separation dynamics. The snapshots
for A:B composition set to 70:30 (off-critical mixture) are shown in Fig. 5.10
for both pure 2D case (A) and quasi-2D case with Lz ¼ 40 (B). Qualitatively,

Dynamics of Heterogeneity in Fluid Membranes

149

Figure 5.10 The snapshots of a 70:30 (off-critical) mixture undergoing phase separation
at t ¼ 0, 150, and 1000 (top to bottom) for a pure 2D (A) and quasi-2D system with Lz ¼ 40
(B). The same sequences of a 50:50 (critical) mixture for a pure 2D (C) and quasi-2D system (D).

it is seen that the domains for the quasi-2D case are smaller in size whencompared at the same time step. We also monitor the average domain sizeR(t) which can be obtained from the total interface length L(t) between thetwo components. This is because R(t) and L(t) are related by L(t) ¼ 2pN(t)R(t), where N(t) is the number of domains. The area occupied by the
B-component is given by A ¼ pN ðt ÞR2 ðtÞ which is a conserved quantity.
Then we haveRðtÞ ¼ 2A=L ðtÞ:

ð5:42Þ

When the domain size grows as R ta, one has L t a and N t 2a. The
domain size R(t) for 70:30 mixture is shown in Fig. 5.11. In this plot, anaverage over 10 independent trials has been taken. It can be seen that thepure 2D case has a growth exponent a ¼ 1/2. Upon the addition of solvent,
we observe that the exponent shifts to a lower value of a ¼ 1/3. This exponent is reminiscent of the phase separation dynamics of an off-critical mixture in 3D. By systematically increasing the amount of solvent in the system
by changing the height Lz, we can see a clear deviation from the pure 2D

150

Shigeyuki Komura et al.

1

0.8

log10 Ra = 1/2
0.6

0.4

a = 1/3

0.2
70:30
0

1

1.5

2

2.5

3

log10 tFigure 5.11 The average domain size R as a function of time t for a 70:30 off-criticalmixture. The upper curve is the pure 2D case showing an a ¼ 1/2 scaling, and the lower
curve is the quasi-2D case when Lz ¼ 40 showing a distinct a ¼ 1/3 scaling.
1.5

a = 1/2

log10 R1

a = 1/3
0.5

50:50
0

1

1.5

2
log10 t2.5

Figure 5.12 The average domain size R as a function of time t for a 50:50 critical mixture.
The upper curve is the pure 2D case showing an a ¼ 1/2 scaling, and the lower curve isthe quasi-2D case when Lz ¼ 40 showing a distinct a ¼ 1/3 scaling.

behavior (not shown). There is no further change if Lz is increased beyond
40. A larger system size Lx Ly ¼ 200 200 also produced the same scaling
for the pure 2D case, which demonstrates that finite-size effects are small.
In Fig. 5.12, we show the result for a component ratio of 50:50 (criticalmixture), whereas the corresponding patterns are given in Fig. 5.10C and D
for 2D and quasi-2D cases, respectively. In this case, the growth exponent

Dynamics of Heterogeneity in Fluid Membranes

151

for the pure 2D case is less obvious owing to rapid coarsening of thedomains. However, by simulating a bigger system 200 200 with the same
areal density, an a ¼ 1/2 exponent is indeed obtained. Similar to the offcritical case, the growth of the domains is slowed down by the additionof solvent and the exponent is reduced to a ¼ 1/3. These results indicate thatsolvent is responsible for slowing down the growth dynamics.
The observed exponent a ¼ 1/2 in pure 2D systems can be explained interms of the BC mechanism [55]. From dimensional analysis, the 2D diffusion
coefficient of the domain is given by D2 kBT/, where  is the membrane
2D viscosity. Using the relation


kB T
2
R
D2 t
t,
ð5:43Þ

we find R t1/2. For 3D systems, on the other hand, the diffusion coefficientof the droplet is inversely proportional to its size, D3 1/R, a well-known
Stokes–Einstein relation. Hence the BC mechanism in 3D gives rise to an
exponent a ¼ 1/3. (In general, the exponent is a ¼ 1/d, where d is the spacedimension.) The change in the exponent from a ¼ 1/2 to 1/3 due to theaddition of solvent implies the crossover from 2D to 3D behaviors of
the phase separation dynamics even though the lateral coarsening takes place
only within the 2D geometry. Notice that the LS mechanism shows an
exponent of a ¼ 1/3 in both 2D and 3D. We thus conclude that oursimulations are still in the early time of the coarsening dynamics.

4.4. Correlated diffusion
In order to justify our argument, it is necessary to examine the size dependence of the domain diffusion coefficient in quasi-2D systems. This can becalculated by tracking the mean-squared displacement of domains of various
radii. The equivalent information can be more efficiently obtained by
calculating the two-particle longitudinal coupling diffusion coefficient ina single component membrane rather than in a binary system.
Consider a pair of particles separated by a 2D vector r, undergoing diffusion in the fluid membrane. The two-particle mean-squared displacement
is given by [25,26]
D
E
Drak Drbl ¼ 2Dklab ðrÞt,ð5:44Þ
where Drka is the displacement of the particle k(¼1, 2) along the axis a(¼x, y),and Dkl
ab is the diffusion tensor giving self-diffusion when k ¼ l and the152

Shigeyuki Komura et al.

coupling between them when k ¼
6 l. The x-axis is defined along the line connecting a pair of particles 1 and 2, that is, r ¼ r x^. Hence, we have D12
xy ¼ 0 bysymmetry. The longitudinal coupling diffusion coefficient, DL ðr Þ ¼ D12
xÞ,xx ðr^
gives the coupled diffusion along the line of centers of the particles.
We first describe the analytical expression of DL(r) for the free membranecase, which can be essentially obtained from Eq. (5.11) (see Ref. [17] fordetails). Using the Einstein relation for over-damped dynamics, we obtain
Dfree
L ðr Þ ¼
At short distances r
becomeskB T
2
pH1 ðnr Þ pY1 ðnr Þ
 2 2þ

:
4p
nr
nr
nr

ð5:45Þ

n 1, the asymptotic form of the above expression
 
2
1
ln
gþ ,
4p
nr
2

kB T
Dfree
L ðr Þ 

ð5:46Þ

where g ¼ 0.5772 . . . is Euler’s constant. At large inter-particle separationsr  n 1, on the other hand, Eq. (5.45) reduces to
Dfree
L ðr Þ 

kB T
kB T
,
¼
2pnr 4ps r

ð5:47Þ

showing the asymptotic 1/r decay which reflects the 3D nature of this limit.
Notice that Eq. (5.47) depends only on the solvent viscosity s but not on themembrane viscosity  any more.
In Fig. 5.13, we plot the measured longitudinal coupling diffusion coefficient DL as a function of 2D distance r. In these simulations, we have worked
with only single component membranes with the same system sizes and number of particles as those used for the phase separation simulations. We have also
taken average over 20 independent trials. In the pure 2D case without any
solvent, DL shows a logarithmic dependence on r. This is consistent with
Eq. (5.46) obtained when the coupling between the membrane and solvent
is very weak so that the membrane can be regarded almost as a pure 2D system.
Using Eq. (5.46) as an approximate expression, we get from the fitting as
kBT/4p  0.89 10 2 and n 1  20. In an ideal case, the SD screening
length should diverge due to the absence of solvent. The obtained finite value
for n 1 is roughly set by the half of the system size in the simulation.
When we add solvent (Lz ¼ 40), the DL is decreased and no longer
behaves logarithmically. In this case, we use the full expression Eq. (5.45)
for the fitting and obtained kBT/4p  1.35 10 2 and n 1  1. In the above

153

Dynamics of Heterogeneity in Fluid Membranes

10–1

DL

10–2

10–3

1/r
10–4

1

10

100

r

Figure 5.13 Longitudinal coupling diffusion DL as a function of particle separation r. The
upper circles are data for the pure 2D case. The lower squares correspond to the casewith solvent when Lz ¼ 40. The upper solid line is the fit by Eq. (5.46), and the lower solidline is the fit by Eq. (5.45). The dashed line shows the 1/r dependence.

two fits, we have neglected the first two points as they lie outside the range ofvalidity, r  1, of Eq. (5.45). Since n 1  1 when the solvent is present, thedata shown in Fig. 5.13 are in the crossover region, r ≳ n 1, showing an approach toward the asymptotic 1/r dependence as in Eq. (5.47). Hence we conclude that the solvent brings in the 3D hydrodynamic property to the diffusion
in membranes. This is the reason for the exponent a ¼ 1/3 in the phase separation dynamics and justifies that it is mainly driven by BC mechanism.
In our simulations, the membrane and the solvent have very similar viscosities. This sets the SD length scale to be of the order of unity (in units ofparticle size), which is consistent with the value n 1  1 obtained from thefitting. As explained above, the fit also provides the 2D membrane viscosity
as   6, and hence we obtain as s  3. In real biomembranes sandwiched bywater, the value of the SD length is much larger than the lipid size and is in
the order of submicron scale [23]. Hence the 3D nature of hydrodynamics
should be observed for large enough domains [60].
For confined membranes, on the other hand, the appropriate mobility
tensor is Gcon
ab given by Eq. (5.12). In this case, the longitudinal coupling
diffusion coefficient can be obtained as [17]DconL ðr Þ ¼

kB T 1K1 ðkr Þ
,2p k2 r 2
kr

ð5:48Þ

154

Shigeyuki Komura et al.

where K1(z) is modified Bessel function of the second kind (see also Eq. 13).
At short distances r k 1, we have
 
kB T
2
1
con
DL ð r Þ lngþ ,ð5:49Þ
4p
kr
2
which is almost identical to Eq. (5.46) except v is replaced now by k. At long
distances r  k 1, on the other hand, we get
Dcon
L ðr Þ 

kB T
kB Th
¼
,
2
2
2pk r4ps r 2

ð5:50Þ

which exhibits a 1/r2 dependence. This is in contrast to Eq. (5.47). Repeating the similar scaling argument, we predict that, in the presence of walls, thedomain growth exponent should be a ¼ 1/4 within the BC mechanism. A
similar simulation with two confining hard walls is required as a next step. In
biological systems, the above situation with solid walls can be relevant because the cell membranes are strongly anchored to the underlying cytoskeleton or are tightly adhered to other cells.

4.5. Related works
We mention here several related works on the phase separation dynamics
below the transition temperature. Veatch and Keller observed the kineticsof domain growth on ternary vesicles as presented in Fig. 5.14 [9]. For a

A
Domainripening
0s

1s

4s

20 s

30 s

50 s

85 s

95 s

0s

1s

2s

3s

5s

12 s

14 s

64 s

0s

7s

8s

9s

0s

7s

12 s

16 s

B
Spinodaldecomposition

C
Viscousfingering

Figure 5.14 Phase separation dynamics on ternary giant vesicles. (A) Brownian coagulation observed for 1:1 DOPC/DPPC þ 25% Chol., (B) spinodal decomposition observedfor 1:1 DOPC/DPPC þ 35% Chol., (C) viscous fingering observed for 1:9 DOPC/DPPC þ 25% Chol. (left series) and 1:1 DOPC/DMPC þ 25% Chol. (right series). Adaptedfrom Ref. [9].

Dynamics of Heterogeneity in Fluid Membranes

155

vesicle with an off-critical composition, dark circular domains grow bycolliding and coalescing with each other rather than through the LS
evaporation–condensation process (Fig. 5.14A). When the composition isnearly critical, spinodal decomposition takes place when the temperature
is decreased through transition point (Fig. 5.14B). For a highly asymmetriccomposition as in Fig. 5.14C, striped domains are observed when thetemperature is raised; possibly viscous fingering.
The domain growth exponent a was first measured by Saeki et al. for ternary vesicles composed of DOPC, DPPC, and cholesterol [61]. They indeed
observed a power-law behavior and reported a rather small value a ¼ 0.15.
The reason for this slow dynamics was partially attributed to the curvatureof the domain portion. Subsequently, Yanagisawa and some of the present
authors found that the domain coarsening processes are classified into two
types, that is, normal coarsening and trapped coarsening [62]. In the former case,the domains having flat circular shape grow through BC process, although thegrowth exponent was found to be a ¼ 2/3. For the trapped coarsening, on theother hand, the domain growth is suppressed at a certain domain size because
the repulsive inter-domain interactions obstruct the coalescence of domains.
A two-color imaging technique of the trapped domains revealed that thisrepulsive interactions are induced by the budding of domains. By observing
the motion of domains, they also found that each domain is attracted toward
the largest one following the flow around it. This is a strong indication of
the hydrodynamic interactions acting between the domains.
As for computer simulations, Laradji and Kumar performed a DPD study
on phase separation dynamics in both two-component vesicles and open
membranes using a coarse-grained model for the membrane lipids[63–65]. In their model, the self-assembly of a bilayer in the presence ofsolvent is naturally taken into account. For off-critical case, they obtained
an exponent a  0.3 which is in good agreement with our simulation
result, a  1/3. Similar to our observation, they also reported that the BC
process takes place rather than the LS evaporation–condensation process.
However, when the excess area of a vesicle is large enough, domains are
no more flat and reshape into caps. In such a case, Eq. (5.42) is not valid
and it is more convenient to measure the interface length L(t). According
to their simulation, it scales as L(t) t 0.4, indicating that the interfacial
length decreases faster as compared to the case of small excess area. For
membranes with critical composition, the growth exponent was found to
be a  0.5 for both vesicles and open membranes. This value is incontrast to our result a  1/3 as shown in Fig. 5.12. In Fig. 5.15, we156

Shigeyuki Komura et al.

Figure 5.15 The snapshots for a 50:50 critical mixture undergoing phase separation foropen membranes with different projected areas (I and II). Adapted from Ref. [64].

present a sequence of spinodal decomposition for open membranes obtained
from their simulation. One of the differences between their simulation andours is that the membrane is not allowed to exhibit out-of-plane fluctuations in
our case. In their subsequent work, the effect of asymmetry in the bilayer lipid
composition was studied [66]. Notice that flip-flop events are extremely rare in
their simulation. This asymmetry sets a spontaneous curvature for domains
which are capped. Interestingly, they found that the coupling between the
spontaneous curvature and composition effectively leads to micro-phaseseparation. The observed scaling behavior was L(t) t 0.13, indicative of a

Dynamics of Heterogeneity in Fluid Membranes

157

logarithmic growth. Obviously, such a bilayer nature of fluid membranes is nottaken into account in our simulation.
Camley et al. [67,68] and Fan et al. [69] performed numerical simulationsof 2D version of model H using the mobility tensor of Eq. (5.10) for freemembranes. The basic formulation of this model is given by that
described in Section 3.1 except that the free-energy functional in
Eq. (5.14) is replaced withða1c
F fcg ¼ dr c2 þ c4 þ ðrcÞ2 ,ð5:51Þ
2
42
where a < 0 is now negative in order to study the domain growth below thetransition temperature. In the dynamical equations, Camley et al. included
Gaussian white thermal noise terms which are distributed with variances set
by the fluctuation–dissipation theorem [67,68]. By using reasonable
parameter values, they were almost successful in reproducing theexperimentally observed pattern evolution given in Fig. 5.14B. For thedomain growth law, it was pointed out that the scaling behavior R ta
would hold only when R n 1 or R  n 1 because it requires the factthat R is the only relevant length scale.
For off-critical mixtures, there are at least two coarsening mechanism inthis model: BC process and LS evaporation–condensation process, as mentioned before. The LS mechanism is essentially driven by the line tension l
between coexisting phases, and the corresponding growth law is given by [50]RðLlt Þ1=3 ,ð5:52Þ

where L is the kinetic coefficient given in Eq. (5.15). Notice that this exponent a ¼ 1/3 does not depend on the dimensionality of the system. (In3D, for example, the line tension l is replaced with a surface tension.)On the other hand, the BC mechanism leads to a ¼ 1/2 for R n 1 anda ¼ 1/3 for R  n 1, as we have discussed in detail. All these different scaling
regimes have been confirmed by their numerical simulation [68], in agreement with our DPD simulation.
For critical mixtures, however, the results are rather complex.
Figure 5.16 shows the time evolution of the phase-separated patterns fora mixture of critical composition [69]. Comparing Fig. 5.16E and F, for example, we notice that the morphologies and sizes of the elongated domains
have changed considerably, whereas more isotropic ones remain almost thesame. As a result, isolated circular domains coarsen slower than elongated158

Shigeyuki Komura et al.

Figure 5.16 Time sequences of concentration maps during spinodal decomposition
for different combinations of membrane and solvent viscosities (A–C and D–F). Time
increases from left to right. Adapted from Ref. [69].

ones. From this observation together with quantitative analysis, they concluded that dynamical scaling breaks down since a single length scale cannot
represent the domain distribution for different times [68,69]. Such a
violation of the scaling behavior was reported before for pure 2D binary
fluids [70,71]. Nevertheless, Camley et al. also found a region of apparent
scaling where a ¼ 1/2 holds [68]. Although this may not be a universal
feature, the power-law behavior seems to appear due to the competition
between thermal and hydrodynamic effects. Furthermore, even a value
close to a ¼ 1/3 was observed when the kinetic coefficient L is large[68]. This is what they call the “Cahn-Hilliard” regime in which bulk
diffusion described by L is dominant. As a whole, the domain coarsening
in critical mixtures seems to be far from universal and may dependent on
the membrane viscosity.
In some other models, membranes undergoing simultaneously a phaseseparation and permanent exchange of lipids with the surrounding medium
were considered [72,73]. For example, Foret proposed the following
dynamical equation in the absence of flow [74]


@c
dF
 ,
¼ Lr2
g cc@tdc

ð5:53Þ

Dynamics of Heterogeneity in Fluid Membranes

159

 is the average value of c, and F iswhere g is the typical exchange rate, cgiven by Eq. (5.51). This equation should be contrasted with Eq. (5.15). It isinteresting to note that the above equation also describes the dynamics ofblock copolymer systems which undergo micro-phase separation [75]. It
is well known for block copolymers that the domain growth will eventuallystop, and the domain size reaches a stationary value. Hence Eq. (5.53) gives
rise to finite-sized raft domains in equilibrium. We point out here that 2D
microemulsion discussed in Section 3.3 also exhibits a micro-phase separation below the transition temperature. Gómez et al. proposed a similar
model in which only cholesterol is recycled [76,77]. As a generalizationof these models, Fan et al. discussed the influence of nonequilibrium lipid
transport, membrane compartmentalization, and membrane proteins on
the phase separation dynamics [78–80]. Recently, micro-phase separationin nonequilibrium membrane was discussed in Ref. [81]. It should be
noted, however, that the hydrodynamic effect is not taken into accountin these models.

5. CONCLUSION AND OUTLOOKIn this chapter, we have discussed the dynamics of heterogeneity inmulticomponent fluid membranes particularly focusing on hydrodynamic
effects. For membranes above the transition temperature, we haveinvestigated the kinetics of critical concentration fluctuations in twocomponent mixtures. Based on the time-dependent Ginzburg–Landau
approach with full hydrodynamics, the wave number dependence of theeffective diffusion coefficient was calculated as a function of the temperature and/or the thickness of the bulk fluid. For membranes below thetransition temperature, we have discussed the effects of an embedding bulk
solvent on the phase separation dynamics using DPD simulations. We haveshown that the presence of a bulk fluid alters the domain growth exponent
from that of 2D to 3D, indicating the significant role played by themembrane–solvent coupling. We further demonstrated that quasi-2D
phase separation proceeds by the BC mechanism which reflects the 3D
nature of the bulk solvent.
Even though the models presented in this chapter capture the essential
physics, they are simplistic in several respects. The bending stiffness of typicalmembranes is of the order of 10kBT which may be not always large enoughto neglect the out-of-plane displacements of the membrane. The dynamics

160

Shigeyuki Komura et al.

of out-of-plane fluctuations of a homogeneous membrane was previously
studied by Levine and MacKintosh [82]. It is worthwhile to comment thatthe presence of a substrate or the second membrane would suppress out-ofplane membrane fluctuations [83,84]. In our treatment, we have also
neglected the effects of membrane curvature which can be significant
when the radius of curvature becomes close to the hydrodynamic
screening length. Henle et al. considered the diffusion of a point object
on a spherically closed membrane [85,86]. The extension of this work tofinite-sized objects is an interesting but difficult proposition.
In some theories, viscoelasticity of an infinitely large flat membrane has
been considered [87,88]. Although it turned out that membranes are
purely viscous in the latest report [89], their experimental techniqueusing particle tracking microrheology provides us with a new clue toinvestigate the dynamical responses of lipid bilayers coupled to thesurrounding environments. Recently, viscoelasticity of phospholipid
Langmuir monolayers in a liquid-condensed phase was measured using
active microrheology technique [90]. From our point of view, however,the viscoelasticity of the surrounding media is more important for largescale dynamics in membranes. This is indeed the case according to therecent work by Granek who calculated the frequency-dependenttransverse fluctuations of a membrane surrounded by viscoelastic media[91]. Currently, we are investigating the lateral dynamics in a purely
viscous lipid membrane surrounded by viscoelastic media whose viscosity
is a frequency-dependent one s[o] [92]. For example, the frequencydependent membrane mobility tensor is obtained as


1
ka kb
Gab ½k, o ¼ 2
ð5:54Þ
dab  2 ,kk þ 2s ½ok0 cothðk0 hÞ þ ior
where


iors 1=2:k ¼k 1þs ½ok20

ð5:55Þ

In the above, r and rs are the densities of the membrane and solvent, respectively. In the limit of o ! 0, these expressions reduce to Eq. (5.9). Using
Eq. (5.54), we can show that the mean square displacement of a disk embedded in the membrane exhibits an anomalous diffusion [92]. An useful
relation which connects the mean square displacement and the solvent
modulus is obtained.

Dynamics of Heterogeneity in Fluid Membranes

161

ACKNOWLEDGMENTS
We thank D. Andelman, H. Diamant, Y. Fujitani, G. Gompper, T. Hamada, Y. Hirose,
T. Kato, P. B. S. Kumar, S. L. Keller, C.-Y. D. Lu, N. Oppenheimer, and Y. Sakuma
for useful discussions. This work was supported by KAKENHI (Grant-in-Aid for
Scientific Research) on Priority Area “Soft matter physics” and Grant No. 21540420
from the MEXT of Japan. This work was also supported by the JSPS Core-to-Core
Program “International research network for nonequilibrium dynamics of soft matter.”

REFERENCES[1] B. Alberts, A. Johnson, P. Walter, J. Lewis, M. Raff, Molecular Biology of the Cell,
Garland Science, New York, 2008.
[2] R. Lipowsky, E. Sackmann, Structure and Dynamics of Membranes, Elsevier, Amsterdam,
1995.
[3] K. Simons, E. Ikonen, Functional rafts in cell membranes, Nature 387 (1997) 569–572.
[4] L.J. Pike, Rafts defined: a report on the Keystone symposium on lipid rafts and cell function, J. Lipid Res. 47 (2006) 1597–1598.
[5] M. Leslie, Do lipid rafts exist? Science 334 (2011) 1046–1047.
[6] J. Fan, M. Sammalkorpi, M. Haataja, Formation and regulation of lipid microdomains incell membranes: theory, modeling and speculation, FEBS Lett. 584 (2010) 1678–1684.
[7] C. Eggeling, C. Ringemann, R. Medda, G. Schwarzmann, K. Sandhoff, S. Polyakova,
V.N. Belov, B. Hein, C. von Middendorff, A. Schönle, S.W. Hell, Direct observationof the nanoscale dynamics of membrane lipids in a living cell, Nature 457 (2009)
1159–1162.
[8] S.L. Veatch, S.L. Keller, Organization in lipid membranes containing cholesterol, Phys.
Rev. Lett. 89 (2002) 268101.
[9] S.L. Veatch, S.L. Keller, Separation of liquid phases in giant vesicles of ternary mixturesof phospholipids and cholesterol, Biophys. J. 85 (2003) 3074–3083.
[10] S.L. Veatch, S.L. Keller, Miscibility phase diagrams of giant vesicles containing
sphingomyelin, Phys. Rev. Lett. 94 (2005) 148101.
[11] S.L. Veatch, S.L. Keller, Seeing spots: complex phase behavior in simple membranes,
Biochim. Biophys. Acta 1746 (2005) 172–185.
[12] S.L. Veatch, O. Soubias, S.L. Keller, K. Gawrisch, Critical fluctuations in domain forming lipid mixtures, Proc. Natl. Acad. Sci. U.S.A. 104 (2007) 17650–17655.
[13] A.R. Honerkamp-Smith, P. Cicuta, M. Collins, S.L. Veatch, M. Schick, M. den Nijs,
S.L. Keller, Line tensions, correlation lengths and critical exponents in lipid membranes
near critical points, Biophys. J. 95 (2008) 236–246.
[14] L. Onsager, Crystal statistics: a two-dimensional model with an order-disorder transition,
Phys. Rev. 65 (1944) 117–149.
[15] S.L. Veatch, P. Cicuta, P. Sengupta, A. Honerkamp-Smith, D. Holowka, B. Baird,
Critical fluctuations in plasma membrane vesicles, ACS Chem. Biol. 3 (2008) 287–293.
[16] G. Gompper, M. Schick, Self-assembling Amphiphilic Systems, Academic, New York,
1994.
[17] S. Ramachandran, S. Komura, K. Seki, G. Gompper, Dynamics of a polymer chain confined in a membrane, Eur. Phys. J. E 34 (2011) 46.
[18] K. Seki, S. Ramachandran, S. Komura, Diffusion coefficient of an inclusion in a liquid
membrane supported by a solvent of arbitrary thickness, Phys. Rev. E 84 (2011) 021905.
[19] L.D. Landau, E.M. Lifshitz, Fluid Mechanics, Pergamon Press, Oxford, 1987.
[20] D.K. Lubensky, R.E. Goldstein, Hydrodynamics of monolayer domains at the air-waterinterface, Phys. Fluids 8 (1996) 843–854.

162

Shigeyuki Komura et al.

[21] Th.M. Fischer, The drag on needles moving in a Langmuir monolayer, J. Fluid Mech.
498 (2004) 123–137.
[22] P.G. Saffman, M. Delbrück, Brownian motion in biological membranes, Proc. Natl.
Acad. Sci. U.S.A. 72 (1975) 3111–3113.
[23] P.G. Saffman, Brownian motion of thin sheets of viscous fluid, J. Fluid Mech. 73 (1976)
593–602.
[24] B.D. Hughes, B.A. Pailthorpe, L.R. White, The translational and rotational drag on acylinder moving in a membrane, J. Fluid Mech. 110 (1981) 349–372.
[25] N. Oppenheimer, H. Diamant, Correlated diffusion of membrane proteins and theireffect on membrane viscosity, Biophys. J. 96 (2009) 3041–3049.
[26] N. Oppenheimer, H. Diamant, Correlated dynamics of inclusions in a supported membrane, Phys. Rev. E 82 (2010) 041912.
[27] E. Evans, E. Sackmann, Translational and rotational drag coefficients for a disk moving
in a liquid membrane associated with a rigid substrate, J. Fluid Mech. 194 (1988)
553–561.
[28] H. Stone, A. Ajdari, Hydrodynamics of particles embedded in a flat surfactant layer
overlying a subphase of finite depth, J. Fluid Mech. 369 (1998) 151–173.
[29] K. Seki, S. Komura, Brownian dynamics in a thin sheet with momentum decay, Phys.
Rev. E 47 (1993) 2377–2383.
[30] S. Komura, K. Seki, Diffusion constant of a polymer chain in biomembranes, J. Phys. II
5 (1995) 5–9.
[31] S. Ramachandran, S. Komura, M. Imai, K. Seki, Drag coefficient of a liquid domain in a
two-dimensional membrane, Eur. Phys. J. E 31 (2010) 303–310.
[32] K. Seki, S. Komura, M. Imai, Concentration fluctuations in binary fluid membranes,
J. Phys. Condens. Matter 19 (2007) 072101.
[33] S. Ramachandran, S. Komura, K. Seki, M. Imai, Hydrodynamic effects on concentration fluctuations in multicomponent membranes, Soft Matter 7 (2011) 1524–1531.
[34] R. Brewster, P.A. Pincus, S.A. Safran, Hybrid lipids as a biological surface-active component, Biophys. J. 97 (2009) 1087–1094.
[35] R. Brewster, S.A. Safran, Line active hybrid lipids determine domain size in phaseseparation of saturated and unsaturated lipids, Biophys. J. 98 (2010) L21–L23.
[36] T. Yamamoto, R. Brewster, S.A. Safran, Chain ordering of hybrid lipids can stabilize
domains in saturated/hybrid/cholesterol lipid membranes, Europhys. Lett. 91 (2010)
28002.
[37] T. Yamamoto, S.A. Safran, Line tension between domains in multicomponent membranes is sensitive to degree of unsaturation of hybrid lipids, Soft Matter 7 (2011)
7021–7033.
[38] P.M. Chaikin, T.C. Lubensky, Principles of Condensed Matter Physics, Cambridge
University Press, Cambridge, 1995.
[39] M. Nonomura, T. Ohta, Decay rate of concentration fluctuations in microemulsions,J. Chem. Phys. 110 (1999) 7516–7523.
[40] K. Kawasaki, Kinetic equations and time correlation functions of critical fluctuations,Ann. Phys. 61 (1970) 1–56.
[41] G. van Meer, D. Voelker, G. Feigenson, Membrane lipids: where they are and how they
behave, Nat. Rev. Mol. Cell Biol. 9 (2008) 112–124.
[42] Y. Hirose, S. Komura, D. Andelman, Coupled modulated bilayers: a phenomenologicalmodel, Chemphyschem 10 (2009) 2839–2846.
[43] Y. Hirose, S. Komura, D. Andelman, Concentration fluctuations and phase transitions
in coupled modulated bilayers, submitted.
[44] M. Teubner, R. Strey, Origin of the scattering peak in microemulsions, J. Chem. Phys.
87 (1987) 3195–3200.

Dynamics of Heterogeneity in Fluid Membranes

163[45] Y. Tserkovnyak, D.R. Nelson, Conditions for extreme sensitivity of protein diffusion inmembranes to cell environments, Proc. Natl. Acad. Sci. U.S.A. 103 (2006) 15002–15007.
[46] K. Inaura, Y. Fujitani, Concentration fluctuation in a two-component fluid membranesurrounded with three-dimensional fluids, J. Phys. Soc. Jpn. 77 (2008) 114603.
[47] M. Haataja, Critical dynamics in multicomponent lipid membranes, Phys. Rev. E 80
(2009) 020902.
[48] E.P. Petrov, P. Schwille, Translational diffusion in lipid membranes beyond the
Saffman-Delbrück approximation, Biophys. J. 94 (2008) L41–L43.
[49] H. McConnell, Composition fluctuations, correlated response, and protein solvation inmembranes, J. Chem. Phys. 132 (2010) 205102.
[50] A.J. Bray, Theory of phase-ordering kinetics, Adv. Phys. 51 (2002) 481–587.
[51] K. Binder, D. Stauffer, Theory for the slowing down of the relaxation and spinodal decomposition of binary mixtures, Phys. Rev. Lett. 33 (1974) 1006–1009.
[52] E.M. Lifshitz, L.P. Pitaevskii, Physical Kinetics, Pergamon Press, Oxford, 1981.
[53] E.D. Siggia, Late stages of spinodal decomposition in binary mixtures, Phys. Rev. A 20
(1979) 595–605.
[54] H. Furukawa, Effect of inertia on droplet growth in a fluid, Phys. Rev. A 31 (1995)
1103–1108.
[55] M.S. Miguel, M. Grant, J.D. Gunton, Phase separation in two-dimensional binary
fluids, Phys. Rev. A 31 (1985) 1001–1005.
[56] P. Ossadnik, M.F. Gyure, H.E. Stanley, S.C. Glotzer, Molecular dynamics simulation of
spinodal decomposition in a two-dimensional binary fluid mixture, Phys. Rev. Lett. 72(1994) 2498.
[57] G. Gonnella, E. Orlandini, J.M. Yeomans, Phase separation in two-dimensional fluids:
the role of noise, Phys. Rev. E 59 (1999) R4741–R4744.
[58] R.D. Groot, P.B. Warren, Dissipative particle dynamics: bridging the gap between
atomistic and mesoscopic simulation, J. Chem. Phys. 107 (1997) 4423–4435.
[59] P. Espanol, P. Warren, Statistical mechanics of dissipative particle dynamics, Europhys.
Lett. 30 (1995) 191–196.
[60] P. Cicuta, S.L. Keller, S.L. Veatch, Diffusion of liquid domains in lipid bilayer membranes, J. Phys. Chem. B 111 (2007) 3328–3331.
[61] D. Saeki, T. Hamada, K. Yoshikawa, Domain-growth kinetics in a cell-sized liposome,J. Phys. Soc. Jpn. 75 (2006) 013602.
[62] M. Yanagisawa, M. Imai, T. Masui, S. Komura, T. Ohta, Growth dynamics of domains
in ternary fluid vesicles, Biophys. J. 92 (2007) 115–125.
[63] M. Laradji, P.B.S. Kumar, Dynamics of domain growth in self-assembled fluid vesicles,
Phys. Rev. Lett. 93 (2001) 198105.
[64] M. Laradji, P.B.S. Kumar, Domain growth, budding and fission in phase-separating
self-assembled fluid bilayers, J. Chem. Phys. 123 (2005) 224902.
[65] S. Ramachandran, M. Laradji, P.B.S. Kumar, Lateral organization of lipids in
multi-component liposomes, J. Phys. Soc. Jpn. 78 (2009) 041006.
[66] M. Laradji, P.B.S. Kumar, Anomalously slow domain growth in fluid membranes withasymmetric transbilayer lipid distribution, Phys. Rev. E 73 (2006) 040901(R).
[67] B.A. Camley, F.L.H. Brown, Dynamic simulations of multicomponent lipid
membranes over long length and time scales, Phys. Rev. Lett. 105 (2010) 148102.
[68] B.A. Camley, F.L.H. Brown, Dynamic scaling in phase separation kinetics forquasi-two-dimensional membranes, J. Chem. Phys. 135 (2011) 225106.
[69] J. Fan, T. Han, M. Haataja, Hydrodynamic effects on spinodal decomposition kinetics
in planar lipid bilayer membranes, J. Chem. Phys. 133 (2010) 235101.
[70] A.J. Wagner, J.M. Yeomans, Breakdown of scale invariance in the coarsening of phaseseparating binary fluids, Phys. Rev. Lett. 80 (1998) 1429.

164

Shigeyuki Komura et al.

[71] A.J. Wagner, M.E. Cates, Phase ordering of two-dimensional symmetric binary fluids: adroplet scaling state, Europhys. Lett. 56 (2001) 556.
[72] M.S. Turner, P. Sens, N.D. Socci, Nonequilibrium raftlike membrane domains undercontinuous recycling, Phys. Rev. Lett. 95 (2005) 168301.
[73] L. Foret, Aggregation on a membrane of particles undergoing active exchange with areservoir, Eur. Phys. J. E 35 (2012) 12.
[74] L. Foret, A simple mechanism of raft formation in two-component fluid membranes,
Europhys. Lett. 71 (2005) 508.
[75] I.W. Hamley, The Physics of Block Copolymers, Oxford University Press, Oxford, 1998.
[76] J. Gómez, F. Sagués, R. Reigada, Actively maintained lipid nanodomains inbiomembranes, Phys. Rev. E 77 (2008) 021907.
[77] J. Gómez, F. Sagués, R. Reigada, Nonequilibrium patterns in phase-separating ternary
membranes, Phys. Rev. E 80 (2009) 011920.
[78] J. Fan, M. Sammalkorpi, M. Haataja, Domain formation in the plasma membrane: rolesof nonequilibrium lipid transport and membrane proteins, Phys. Rev. Lett. 100 (2008)
178102.
[79] J. Fan, M. Sammalkorpi, M. Haataja, Lipid microdomains: structural correlations,fluctuations, and formation mechanisms, Phys. Rev. Lett. 104 (2010) 118101.
[80] J. Fan, M. Sammalkorpi, M. Haataja, Influence of nonequilibrium lipid transport, membrane compartmentalization, and membrane proteins on the lateral organization of theplasma membrane, Phys. Rev. E 81 (2010) 011908.
[81] P. Sens, M.S. Turner, Microphase separation in nonequilibrium biomembranes, Phys.
Rev. Lett. 106 (2011) 238101.
[82] A.J. Levine, F.C. MacKintosh, Dynamics of viscoelastic membranes, Phys. Rev. E 66
(2002) 061606.
[83] S. Sankararaman, G.I. Menon, P.B.S. Kumar, Two-component fluid membranes near
repulsive walls: linearized hydrodynamics of equilibrium and nonequilibrium states,
Phys. Rev. E 66 (2002) 031914.
[84] N. Gov, A.G. Zilman, S. Safran, Hydrodynamics of confined geometries, Phys. Rev. E
70 (2004) 011104.
[85] M.L. Henle, R. McGotry, A.B. Schofield, A.D. Dinsmore, A.J. Levine, The effect ofcurvature and topology on membrane hydrodynamics, Europhys. Lett. 84 (2008)
48001.
[86] M.L. Henle, A.J. Levine, Hydrodynamics in curved membranes: the effect of geometry
on particulate mobility, Phys. Rev. E 81 (2010) 011905.
[87] B.A. Camley, F.L.H. Brwon, Beyond the creeping viscous flow limit for lipid bilayer
membranes: theory of single-particle microrheology, domain flicker spectroscopy, andlong-time tails, Phys. Rev. E 84 (2011) 021904.
[88] T. Han, M. Haataja, Comprehensive analysis of compositional interface fluctuations inplanar lipid bilayer membranes, Phys. Rev. E 84 (2011) 051903.
[89] C.W. Harland, M.J. Bradley, R. Parthasarathy, Phospholipid bilayers are viscoelastic,
Proc. Natl. Acad. Sci. U.S.A. 107 (2010) 19146 108 (2011) 14705.
[90] S.Q. Choi, S. Steltenkamp, J.A. Zasadzinski, T.M. Squires, Active microrheology andsimultaneous visualization of sheared phospholipid monolayers, Nat. Commun. 2(2011) 312.
[91] R. Granek, Membrane surrounded by viscoelastic continuous media: anomalous diffusion and linear response to force, Soft Matter 7 (2011) 5281.
[92] S. Komura, S. Ramachandran, K. Seki, EPL 97 (2012) 68007.