April 2016
EPL, 114 (2016) 28002
doi: 10.1209/0295-5075/114/28002

www.epljournal.org

Nano-domain formation in charged membranes:
Beyond the Debye-Hückel approximation
Ryuichi Okamoto1 , Naofumi Shimokawa2 and Shigeyuki Komura1
1

Department of Chemistry, Graduate School of Science and Engineering, Tokyo Metropolitan University
Tokyo 192-0397, Japan
2
School of Materials Science, Japan Advanced Institute of Science and Technology - Ishikawa 923-1292, Japan
received 4 March 2016; accepted in final form 27 April 2016
published online 12 May 201687.16.D- – Membranes, bilayers, and vesicles
82.45.Gj – Electrolytes
PACS 87.16.dt – Structure, static correlations, domains, and rafts
PACS
PACS

Abstract – We investigate the microphase separation in a membrane composed of charged lipids,
by taking into account explicitly the electrostatic potential and the ion densities in the surrounding solvent. While the overall (membrane and solvent) charge neutrality is assumed, themembrane can have a non-zero net charge. The static structure factor in the homogeneous stateis analytically obtained without using the Debye-Hückel approximation and is found to have apeak at an intermediate wave number. For a binary membrane composed of anionic and neutral
lipids, the characteristic wave number corresponds to a scale from several to tens of nanometers.
Our numerical calculation further predicts the existence of nano-domains in charged membranes.
c EPLA, 2016
Copyright ⃝Introduction. – Much attention has been paid tophase separations in artificial multi-component lipid membranes. In these systems, many degrees of freedom suchas lipid composition and membrane shape deformation
are coupled to each other, leading to complex phase behaviors [1]. In particular, long-lived small domains inmembranes may play important biological roles [2–4]. In
general, charge-induced microphase formation has been
intensively investigated in soft matter such as polyelectrolytes [5,6], electrolyte fluid mixtures [7–9], and charged
Langmuir monolayers [10,11]. Several authors have studied the microphase formation in membranes composed ofboth anionic and cationic lipids for which the net charge
within the membrane vanishes [12–14].
In the last decade, researchers have investigated phaseseparations in giant unilamellar vesicles (GUVs) composed of anionic and neutral lipids [15–19]. In theseexperiments, GUVs have a non-zero net charge. Froma theoretical viewpoint, Guttman and Andelman originally predicted a microphase separation in binary charged
membranes within the Debye-Hückel (DH) approximation [20]. However, no evidence for a charge-induced microphase formation on the scale of the optical resolution
was found in refs. [15–19]. Hence some theories for a
macrophase separation in such charged membranes havebeen developed [21–23]. Recently, Puﬀ et al. reported theformation of nanoscale domains, whose scale is smaller
than the optical resolution, with the addition of ganglioside GM1 [24]. In their mixtures, GM1 is anionic while theother components are all neutral. Given these experimental observations, it is necessary to study the competition
between the macrophase and microphase separations incharged membranes, as well as the characteristic lengthscale associated with the microphase separation.
The DH approximation is justified when bκ ≫ 1, whereb and κ are the Gouy-Chapman length and the Debye
wave number, respectively. Notice that b is inversely proportional to the surface charge density. In a strongly segregating charged membrane, each domain usually has a large
surface charge density, and thus the DH approximation is
no more valid [14]. Furthermore, for membranes that have
non-zero net charge, the DH approximation is inapplicable even to homogeneous nor to weakly segregated states.
For a binary membrane composed of anionic and neutral
lipids (as discussed later in more detail), we can estimateφb ∼ 1 Å with φ being the fraction of the anionic lipid,
whereas for a 1:1 electrolyte solution, κ is in the range
10−3 –10−1 Å−1 . Hence the DH condition bκ ≫ 1 is notsatisfied unless φ is very small, and it is imperative to go
beyond it.
In this letter, we investigate the microphase formation in binary membranes composed of charged lipids

28002-p1

Ryuichi Okamoto et al.
for general ionic strength. We assume the overall (membrane and solvent) charge neutrality, while the membrane can have a non-zero net charge. Solving the full
non-linear Poisson-Boltzmann equation (PBE), we discussi) the scale of the microphase structures, and ii) the conditions for the microphase formation. Our theory predicts
a microphase endpoint (MEP) [25] in the compositiontemperature plane at which the macrophase spinodal line
meets an end of the microphase spinodal line. This point
cannot be obtained within the DH approximation. We
also find that the characteristic length scale of the microphase separation is in the range from several to tens ofnanometers except in the vicinity of the MEP.
Free energy and PBE. – As shown in fig. 1, we consider a flat fluid membrane composed of A- and B-lipid
molecules having electric charges eZJ (J = A, B), where eis the elementary charge and ZJ is the valence number. In
water, these lipid molecules form a bilayer structure. Herethe hydrophobic tails face each other and the hydrophilic
head groups are in contact with water. Any interactions
between diﬀerent monolayers are neglected although there
are situations in which inter-monolayer coupling plays arole [22,23]. Then we are allowed to consider only a twodimensional (2D) monolayer located at z = 0 which is incontact with the solvent occupying the region of z > 0 ina three-dimensional (3D) space. We use the abbreviations
x = (x, y) and X = (x, y, z).
We assume that the molecular area a2 occupied by an
A-lipid is the same as that occupied by a B-lipid. The
area fraction of A-lipid is denoted by φ(x), and that of
B-lipid is 1 − φ(x) under the areal incompressibility condition. For simplicity, we assume that the solvent is a 1:1
dilute electrolyte solution. The cation and anion number
densities are denoted by n+ (X) and n− (X), respectively.
In
neutrality is assumed;! 2our work, the overall! charge
3
Z
(1
−
φ)]
+
d
X
(n
d x [ZA φ ++ − n− ) = 0, where the! 2B
integrationd x is taken over the range −∞ < x, y < ∞,!
while d3 X is for the range z > 0 and −∞ < x, y < ∞.
The total free energy functional F = Fm + Fb is given by
the sum of the membrane contribution, Fm (in the absenceof the electrostatics contribution), and the bulk contribution, Fb . The former is given by
"
$
#
Fm
c
(1)
= d2 x f (φ) + (∇∥ φ)2 ,T
2where T is the temperature (we have set the Boltzmann
constant kB to unity), f (φ) is the scaled free energy (without the electrostatic contribution) per unit area for a homogeneous state, and the second term is the standard
squared-gradient form with a positive dimensionless coeﬃcient c. Note that ∇∥ = (∂x , ∂y ) is the 2D gradient
operator.
The charged lipids in the membrane and the ions inthe solvent generate a gradient of the electrostatic potential Ψ(X) satisfying the Poisson equation, ε∇2 Ψ =−4πe(n+ − n− ), where ε is the dielectric constant of

1:1 electrolyte solution
dielectric constantA-lipid

B-lipid

Fig. 1: Schematic illustration of a cross section of a charged
membrane in an electrolyte solution of the dielectric constant ε.
The electric charges of A-lipid and B-lipid are eZA and eZB ,where ZJ (J = A, B) is the valence number and e is the elementary charge.

water, and ∇ = (∂x , ∂y , ∂z ) is the 3D gradient operator.
The areal charge density of the membrane is given by
σ(φ) = [ZB + (ZA − ZB )φ]/a2 in units of e. Then theboundary condition at z = 0 is ε∂z Ψ|z=0 = −4πeσ(φ).
The bulk part of the free energy Fb consists of the entropyof ions and the electrostatic energy,
%
'
"
&
Fb
ε
3
3
2
= d X
(∇Ψ) , (2)ni {ln(ni λi ) − 1} +
T
8πT
i=±
where λi is the thermal de Broglie length of the i-th ion
species.
In equilibrium, the ion densities n1 and n2 approach toa common constant value n∞ as z → 0. The electrostatic
potential Ψ also becomes constant far from the membrane (without loss of generality, we can assume Ψ → 0
as z → !∞). (We then introduce
the grand potential!
Ω = F − d2 x i ni µi∞ − h d3 X (φ − φ0 ), where µi∞ =T ln(n∞ λ3i ) is the ion chemical potential far from the membrane and φ0 is the average composition. In the mean-field
theory, the equilibrium state is obtained by minimizing Ωwith respect to both ni and φ. With the aid of the relation δ(ε|∇Ψ|2 ) = 2ε∇ · (Ψ∇δΨ) + 8πeΨ(δn+ − δn− ), we
derive the nonlinear PBE as the minimization condition
with respect to ni [26]∇2 ψ = κ2 sinh ψ,(3)where ψ = eΨ/T is the dimensionless potential, andκ = [8πe2 n∞ /(εT )]1/2 is the Debye wave number.
Minimization with respect φ yields
h
∂f
∂σ
≡
− c∇2∥ φ + ψ= const.
T
∂φ
∂φ

(4)

Fluctuations about the homogeneous state. – At
high temperatures, the translational entropy of the lipid
molecules dominates inter-molecular interactions, leading
to a homogeneous state, φ = φ0 . In this state, all thevariables are constant in the lateral xy-direction, whileni and ψ depend on z because of the electric charges ofthe lipid molecules. Therefore ψ obeys a one-dimensional

28002-p2

Nano-domain formation in charged membranes: Beyond the Debye-Hückel approximation(1D) PBE ∂z2 ψ = κ2 sinh ψ with the boundary condition∂z ψ|z=0 = −4πℓσ0 , where σ0 = σ(φ0 ) and ℓ = e2 /εTis the Bjerrum length. It is convenient to introduce a
dimensionless number η = κ/(2πℓσ0 ) which can be either
positive or negative. Note here that the Gouy-Chapman
length is given by b = 1/(2πℓ|σ0 |) = |η|/κ.
The 1D PBE has a well-known exact solution [26],ψ0 (z) = 2 ln

1 + Γe−κz.
1 − Γe−κz(5)

In the above, the dimensionless number Γ is theroot of )Γ2 + 2ηΓ − 1 = 0 and is given byΓ = −η ± η 2 + 1 (η ≷ 0). The ion densities in equilibrium are expressed as n±0 (z) = n∞ e∓ψ0 (z) . In the DH
condition |η| ≫ 1 which corresponds to a high salt and/or
small surface charge condition, we have Γ ≃ 1/(2η) andobtain ψ0 ≃ 2e−κz /η.
In order to see the fluctuations about the homogeneous state, we superimpose the variations, φ0 → φ0 +
δφ(x) and ni0 (z) → ni0 (z) + δni (X). We then examine the free energy deviation ∆F up to the bilinear orderin the variations, δφ and δni . Since we are interested in
the fluctuations of δφ, we further minimize ∆F with respect to δni . Here we introduce the! in-plane Fourier transform of a function g(x) as gk = d2 x e−ik·x g(x), wherek = (kx , ky ) is the 2D wave vector. In the Fourier space,∆F is written as*+"
1
d2 k ∂ 2 f (φ0 )∆F
2
2=+ ck + θ0 Pk (0) |δφk |2 ,
T
2
(2π)2
∂φ20
(6)

Pk (z) = R(z̄)P0 (z) into eq. (7) and obtain
R′′ + (ln P02 )′ R′ − ϵ̂k̄ 2 R = 0,(9)where the prime denotes the derivative with respect toz̄, and we have introduced a “book keeping parameter” ϵ̂which will be set equal to unity at the end. We expand R
in powers of ϵ̂ such that R = R0 + ϵ̂R1 + · · · . Then up to
the first order in ϵ̂, we haveR0′′ + (ln P02 )′ R0′ = 0,R1′′ + (ln P02 )′ R1′ = k̄ 2 R0 ,

(10)
(11)

and we may obtain the perturbation solution. However, itturns out that a secular term appears in ϵ̂R1 , and hence
the perturbation solution is only locally valid in the vicinity of the boundary, but not uniformly valid in the entire
region z > 0.
In order to cure the breakdown of such a perturbation
calculation, apparently diﬀerent but almost equivalent
ways have been developed [27–30]. Among these we use
the renormalization group (RG) method in refs. [29,30].
Given the solution G(z̄0 ) at any point z̄0 > 0, we first seek
the solution of eqs. (10) and (11) such that Pk = RP0 → 0
as z̄ → ∞. Some calculation yields
R(z̄; z̄0 ) =*03+2ϵ̂k̄ 2
Γ2 1 −2z̄
G 1−e− e−2z̄0(z̄ − z̄0 ) ++ O(ϵ̂2 ).
2
2
(12)

Here the term proportional to ϵ̂(z̄ − z̄0 ) is a secular termthat becomes large as z̄ − z̄0 is increased, leading to thebreakdown of the regular perturbation calculation. To obtain a uniformly valid solution, we impose in eq. (12) theRG equation, ∂R/∂ z̄0 |z̄=z̄0 = 0. This yields the diﬀerenwith the boundary condition ∂!z Pk (z)|z=0 = −4πℓ. One tial equation for G as
can easily find that ψ0 + θ0 d2 k/(2π)2 eik·x Pk (z) δφk2∂Gϵ̂k̄ 2 1 2 −2z̄0
is a solution of the PBE under the boundary condition
Γ e+− 1 G = 0.
(13)∂z ψ(z)|z=0 = −4πℓ[σ0 + θ0 δφ(x)] when the surface charge∂ z̄02
heterogeneity θ0 δφ(x) is suﬃciently small.
In the DH condition |η| ≫ 1 , we can set cosh ψ0 ≃ 1 in Using the solution G(z̄0 ), we obtain the improved solutionas R(z̄; z̄0 = z̄) = G(z̄), which is uniformly valid up toeq. (7) and obtain
the order of ϵ̂ [30]. Geometrically, the improved solution
. )
/
4πℓ
exp −z k 2 + κ2 .
(8) G(z̄) is the envelope of the family of curves {R(z̄; z̄0 )}z̄0
Pk (z) ≃ √parametrized by z̄0 . It is tangent at each point z̄ = z̄0 to
k 2 + κ2a member of the family R(z̄; z̄0 ) that is locally valid in the
The above expression was obtained in ref. [20] and used in vicinity of z̄ [29].
0simulations [13]. Notice that θ02 Pk (0) in eq. (6) does not
Solving eq. (13) under the boundary condition
depend on φ0 within the DH condition.
∂Pk /∂ z̄|z̄=0 = −2/(ησ0 ), we finally obtain Pk (z) as
Perturbation solution of nonlinear PBE. – In
* 20
3+
2
P0 (z)k̄Γ2 1 −2z̄
order to go beyond the DH approximation and discuss
Pk (z) ≃
e
exp −
−1
z̄ +
, (14)
general values of η, we may seek the solution of eq. (7)22
1 + γ k̄ 2
perturbatively in powers of k 2 . This is because the solution is P0 (z) = ∂ψ0 /∂σ0 for k = 0. Introducing the di- where γ = η 2 Γ/(η + Γ) ≥ 0. For the DH condition |η| ≫mensionless quantities z̄ = κz and k̄ = k/κ, we substitute 1, eq. (14) reduces to eq. (8) within the approximation
where k = |k| and θ0 = ∂σ(φ0 )/∂φ0 = (ZA − ZB )/a2 . In
the above, Pk (0) = Pk (z = 0) where Pk (z) satisfies, 2∂z − k 2 − κ2 cosh ψ0 (z) Pk (z) = 0,(7)28002-p3

Ryuichi Okamoto et al.
(k 2 + κ2 )1/2 ≃ κ(1 + k̄ 2 /2). For a general value of η,
it can be shown that eq. (14) provides a good approximation for any z > 0 if k̄ ≪ 1. Even if k̄ ≪ 1 is notsatisfied, it is still a good approximation for suﬃciently
small z satisfying k̄ 2 (ηΓ + Γ2 z̄)2 ≪ 1. This region, where
eq. (14) is accurate for not so small k̄, indeed exists if(k̄ηΓ)2 ≪ 1. For a low salt and/or large surface charge
condition, |η| ≪ 1, the region is given by k 2 (b + z)2 ≪ 1
if (kb)2 ≪ 1. Equation (14) is enough for our purpose because we only need Pk (0) to calculate the compositional
structure factor (see eq. (6)).
Conditions for microphase separation. – Using thedefinition p = P0 (z = 0) = 2/[σ0 (Γ + η)], we substitutePk (0) ≃ p/(1 + γ k̄ 2 ) into eq. (6) and obtain the structure
factor S(k) = ⟨|δφk |2 ⟩. Here both p and γ are positive.
The positivity of p means that the electrostatic interaction
tends to prevent the instability towards a macrophase separation [23]. In addition, the positivity of γ implies thepossibility of the microphase formation when the temperature is decreased. If c < θ02 pγ/κ2 is satisfied, S(k) takes
a maximum value at an intermediate wave number
k∗ =

κ
γ 1/2

%4

θ02 pγ
cκ2

51/2

'1/2

−1

.

(15)

When the temperature is decreased to a certain value,S(k ∗ ) diverges and the modes δφk with |k| = k ∗ become
unstable. This leads to the microphase separation characterized by a typical wave length 2π/k ∗ .
For the sake of further discussion, we assume thatthe free energy density f (φ) is given by the BraggWilliams form, f (φ) = [φ ln φ + (1 − φ) ln(1 − φ) +χφ(1 − φ)]/a2 , where χ is dimensionless and is roughly
proportional to 1/T . With this choice, the instability
condition towards the microphase separation is written as
χ > χ∗ with
a2 cγ ∗ 4χ∗ = χs −(k ) .
(16)2κ2Here χs = [{φ0 (1 − φ0 )}−1 + a2 θ02 p]/2 is the spinodal forthe macrophase separation. Equation (16) defines the microphase endpoint (MEP) [25], χE = χ∗ = χs , at which
an end of the macrophase spinodal meets the microphasespinodal. The MEP can be a Lifshitz point when it is also
the critical point of the macrophase separation. We notethat the critical point has been
√ discussed, and it√ previously
is located at (φ0 , χ)c = ((3− 3)/2, 2+ 3) ≈ (0.634, 3.73)
in the low salt limit [21–23].
Membranes composed of anionic and neutral
lipids. – Following the previous experiments [15–19,24],
we hereafter set ZA = −1 and ZB = 0 (note that cationiclipids are not so common in biomembranes). For typicallipid membranes in water, we set ℓ = 7 Å and a = 8 Å [31,32]. Then we have φ0 b = a2 /(2πℓ) = 1.46 Å. The Debye
wave number κ is usually in the range 10−3 –10−1 Å−1 for

2

(a)

(b)
1.2

1.5
1

0.8

0.5

0.4

0

10 -2

10 -1

1

0

10 -2

10 -1

1

Fig. 2: (Colour online) The scaled characteristic wave numberk∗ a as a function of the scaled Debye wave number κa for (a)c = 0.15 and (b) 0.4. The filled circles are the microphase
endpoints (MEPs) at which k∗ a vanishes, and correspondto κE .

1:1 electrolyte solutions. The parameter η < 0 is then estimated to be φ0 |η| ∼ 10−3 –10−1 . Hence the DH condition|η| ≫ 1 is not satisfied even for a high salt solution, except
for a very small fraction φ0 of the charged lipid. In neutral
membranes, the line tension τ ∼ T (c/a2 )1/2 between thecoexisting phases has been measured to be several pN [33].
Then we may estimate that c is in the range 0.1–1.
In the present case, we can rewrite eq. (15) as
∗

k =

κ
γ 1/2

%0

|Γ|(φ∗ /φ0 )3(Γ + η)2

31/2

'1/2

−1

(17)

with φ∗ = [a2 /(2cπ 2 ℓ2 )]1/3 . Here the factor |Γ|/(Γ + η)2
is a monotonic decreasing function of |η| and is less thanunity. Therefore, k ∗ cannot exist for any κ when φ0 > φ∗ .
In fig. 2, the maximum wave number k ∗ is plotted as afunction of κ, where we set (a) c = 0.15 (corresponding toφ∗ = 0.76) and (b) c = 0.4 (φ∗ = 0.4). For κb ≪ 1, we
may set η ≃ 0, Γ ≃ −1 and γ ≃ η 2 in eq. (17) to obtaink ∗ b ≃ [(φ∗ /φ0 )3/2 − 1]1/2 , which is independent of κ. As
κ is increased, k ∗ decreases and eventually vanishes atκ = κE (φ0 ) (marked with filled circles in fig. 2) when the
MEP condition (φ∗ /φ0 )3 = (Γ + η)2 /|Γ| holds. It should
be noted that the predicted value of 2π/k ∗ corresponds tothe scale in the range from several to tens of nanometers
unless κ is very close to κE .
In figs. 3(a) and (b), we plot χ∗ as a function of κ forc = 0.15 and 0.4, respectively. All the curves for different φ0 values exhibit non-monotonic dependence on κ
(though it is not apparent for φ0 = 0.74 in (a) and for
φ0 = 0.53 in (b)). When κ is varied for fixed φ0 values,
one can show that χ∗ takes a maximum value when κsatisfies c = 2b2 θ02 Γ3 /σ0 if the inequality c < 2b2 θ02 /|σ0 |holds. In these plots, the MEP is located at χE = χ∗ (κE )(marked with filled circles). For κ < κE , the inequality χs > χ∗ always holds for each curve, indicating thatthe macroscopic (k = 0) mode is more stable than themodes with the characteristic non-zero wave number k ∗(see eq. (16)).

28002-p4

Nano-domain formation in charged membranes: Beyond the Debye-Hückel approximation
4.3

4

(a)

(b)

4.2
3.9

80 (A)

80 (B)

60

60

40

40

20

20

4.1
4

3.8
3.9

0
3.7
10 -2

10 -1

1

3.8 -2
10

10 -1

4.2

4.2

(c)

(D)

4

(C)

(A)

4

3.7
0.6

40

60

0

80

0

20

40

60

80

0

20

40

60

80

0.6

0.7 0.8

0.9

80 (C)

80 (D)

60

60

40

40

20

20

(d)
4.1

3.8

3.6
0.5

20

3.9

(B)

3.8

0

1

0.7

0.8 0.4

0
0.5

0.6

0

20

40

60

80

0

0.7

0.2

Fig. 3: (Colour online) Stability diagrams in the (κ, χ)-plane((a) and (b)) and in the (φ0 , χ)-plane ((c) and (d)). The parameters are c = 0.15 in (a) and (c), and c = 0.4 in (b) and(d). The spinodal lines of the microphase separation are plotted with solid lines, while those of the macrophase separation
are shown by the dotted lines. The MEP for each curve ismarked with a filled circle. The four crosses in (c) correspondto the numerical simulations in fig. 4.

In figs. 3(c) and (d), we plot χ∗ as a function of φ0
for c = 0.15 and 0.4, respectively. For κa = 0.01 andfor φ0 " 0.01, we are allowed to set Γ ≃ −1 and η ≃ 0
in eq. (17). In this case, the MEP is almost located at(φ0 , χ)E = (φ∗ , χ∗ (φ∗ )). For a larger value of κ, φE exhibits a downward shift. Because k ∗ does not dependon φ0 within the DH approximation, it cannot predict
any MEP on the (φ0 , χ)-plane (see also below eq. (8)).
Indeed, within the DH approximation, the condition for
the existence of k ∗ is given by c < 2πℓθ02 /κ3 that is
independent of φ0 . With our choice of ℓ and a, this inequality becomes κa < 1.76c−1/3 and is always satisfied
for the parameter values used in figs. 2 and 3. In general,the DH approximation overestimates the possibility of the
microphase separation.
Numerical simulation. – Finally, we numerically integrate the equilibrium conditions, eqs. (3) and (4). We
set ZA = −1 and ZB = 0 as before. We prepare in
the xy-plane a membrane of size (85.3a)2 that is in contact with an electrolyte solution at z = 0. The size Lzof the solvent container is set to 16.25a. At z = Lz ,we impose the boundary condition, ψ(Lz ) = 0, which
is justified when κLz ≫ 1. A periodic boundary condition is employed in the lateral xy-directions. We thensolve simultaneously the fictitious dynamic equations,∂φ/∂t = −h + ⟨h⟩ and ∂ψ/∂t = ∇2 ψ − κ2 sinh ψ, where

0.3

0.4

0.5

Fig. 4: (Colour online) Equilibrium composition profiles φ(x).
The membrane size is measured in units of a. The values ofχ and φ0 are (A) (φ0 , χ) = (0.55, 3.72), (B) (0.6, 3.67), (C)(0.68, 3.9) and (D) (0.73, 4.05) while c = 0.15 and κa = 1.0
are fixed. These four cases are marked with crosses (A)–(D) in
fig. 3(c). Corresponding to these cases, we have (A) (k∗ a, χ∗ ) =(1.90, 3.65), (B) (1.65, 3.64) and (C) (0.219, 3.72), whereas in(D) the characteristic wave number k∗ does not exist and themacrophase separation occurs at χs = 3.87.

!
!
⟨· · · ⟩ = d2 x (· · · )/ d2 x denotes the areal average in themembrane. Although the dynamics itself has no physicalmeaning, the equilibrium pattern can be eﬃciently obtained as a stationary state.
In fig. 4, we present the numerically obtained equilibrium profiles of φ(x). We examine the phase separation
by varying the composition φ0 and the χ values as marked
with the four crosses (A)–(D) in fig. 3(c). As expected,
we clearly see microphase separations in (A)–(C) and amacrophase separation in (D). In (A) and (B), hexagonal
and stripe patterns are obtained, respectively. In thesetwo cases, the parameters are chosen so that the modes
with k = k ∗ are unstable wheareas the macroscopic mode(k = 0) is linearly stable. The characteristic lengths areapproximately (A) 2.6 nm and (B) 3.0 nm within our parameters. In fig. 4(C), where both macroscopic (k = 0)
and intermediate (k = k ∗ ) modes are unstable, there is no
periodic pattern but ring-like aggregates are formed. For
the parameter values in (C), a macrophase separated stateis also a metastable state where the composition profile issimilar to that in (D). We have found that the value ofthe grand potential Ω for the microphase is smaller thanthat of the macrophase. This suggests that the transition between the micro and the macro phases is discontinuous when the parameters (φ0 , χ and κ) are varied. In
ref. [12], the authors have discussed the strong segregation

28002-p5

Ryuichi Okamoto et al.
limit of membranes with no net charge. Assuming the
DH approximation, they also found a discontinuous transition between the micro and the macro phases as κ isvaried.

from the Ministry of Education, Culture, Sports, Science,
and Technology of Japan (MEXT). SK acknowledges
support from the Grant-in-Aid for Scientific Researchon Innovative Areas “Fluctuation and Structure” (Grant
No. 25103010) from the MEXT, the Grant-in-Aid for
Scientific Research (C) (Grant No. 15K05250) from theJSPS, and the JSPS Core-to-Core Program “International
Research Network for Non-equilibrium Dynamics of Soft
Matter ”.

Conclusion. – In this letter, we have investigated themicrophase formation in charged membranes, where themembrane can have non-zero net charge. Without assuming the DH condition, we obtained the solution eq. (14)of the nonlinear PBE when the charge heterogeneity inthe membrane is small. Using the solution, we have discussed the microphase separation in a binary membrane REFERENCES
composed of anionic and neutral lipids, for which the DH[1] Komura S. and Andelman D., Adv. Colloid. Interface
approximation is not justified except for a very small fracSci., 208 (2014) 34.
tion of the charged lipid. Our theory reveals that the[2] Simons K. and Ikonen E., Nature, 387 (1997) 569.
characteristic wave number k ∗ , at which S(k) takes a max[3] Simons K. and Gerl M. J., Nat. Rev. Mol. Cell Biol.,
imum value, corresponds to the scale in the range from11 (2010) 688.
[4] Arumugam S. and Bassereau P., Essays Biochim., 57
several to tens of nanometers except in the vicinity of the(2015) 109.
MEP. This explains why microphase separated structures[5]
Borue V. Y. and Erukhimovich I. Y., Macromolecules,
have not been observed by optical microscopy measure21(1988) 3240.
ments. We further predict a charge-induced MEP in the[6] Joanny J. F. and Leibler L., J. Phys., 51 (1990) 545.
composition-temperature (φ0 , χ)-plane, which cannot be[7] Sadakane K. et al., Phys. Rev. Lett., 103 (2009) 167803.
obtained within the DH approximation. The numerical[8] Onuki A. et al., Bull. Chem. Soc. Jpn., 84 (2011) 569.
simulation also shows that our model exhibits both the[9] Onuki A. et al., J. Phys.: Condens. Matter, 23 (2011)
microphase and macrophase separations depending on the284113.
composition and/or the temperature.
[10] Iimura K. et al., Langmuir, 18 (2002) 10183.
We make further remarks. i) For small scales corre- [11] Caruso B. et al., J. Phys. Chem. B, 118 (2014) 519.
sponding to our predicted values of k ∗ , it would be prefer- [12] Solis F. J. et al., J. Chem. Phys., 122 (2005) 054905.
able to use a microscopic density functional free energy [13] Velichko Y. S. and Olvera de la Cruz M., J. Chem.
Phys., 124 (2006) 214705.
for Fm rather than the mean field f (φ) and the gradient[14]
Naydenov
A. et al., Langmuir, 23 (2007) 12016.
expansion form in eq. (1). Nevertheless, we believe that[15] Shimokawa N. et al., Chem. Phys. Lett., 496 (2010) 59.
the present approximation provides reliable predictions of[16] Vequi-Suplicy C. C. et al., Biochim. Biophys. Acta
the nano-domain formation. ii) The DH approximation
Biomembr., 1798 (2010) 1338.
is valid for a very small charge density σ0 in the homo- [17] Blosser M. C. et al., Biophys. J., 104 (2013) 2629.
geneous state and/or in the weak segregation regime, as [18] Pataraia S. et al., Biochim. Biophys. Acta Biomembr.,discussed in this letter. However, even if σ0 is very small,1838 (2014) 2036.
it is not justified in the strong segregation regime where [19] Himeno H. et al., Soft Matter, 10 (2014) 7959.
the charge density in each domain becomes large. With- [20] Guttman G. and Andelman D., J. Phys. II, 3 (1993)
1411.
out assuming the DH condition, Naydenov et al. discussed
the strong segregation regime (as well as the weak segre- [21] May S. et al., Phys. Rev. Lett., 89 (2002) 268102.
gation regime) of a membrane that has no net charge, [22] Baciu C. L. and May S., J. Phys.: Condens. Matter, 16(2004) S2455.
i.e., σ0 = 0 [14]. In the weak segregation limit where the[23] Shimokawa N. et al., Phys. Rev. E, 84 (2011) 031919.
DH approximation is justified for membranes with σ0 = 0,[24] Puff N. et al., Biochim. Biophys. Acta Biomembr., 1838
their result coincides with eq. (8). iii) The solution eq. (14)(2014) 2105.
has further applications, such as the electrostatic contribu- [25] Olmsted P. D. and Hamley I. W., Europhys. Lett., 45
tion to the bending rigidity [20], and the charge regulation(1999) 83.
eﬀect [26] on a surface which has a spatially heterogeneous [26] Israelachivili J. N., Intermolecular and Surface Forces,
ionizable group distribution.
3rd edition (Academic Press, London) 2011.
∗∗∗
RO and SK thank C. Watanabe for informative
discussion. NS acknowledges support from the Grant-inAid for Young Scientist (B) (Grant No. 26800222) fromthe Japan Society for the Promotion of Science (JSPS)and the Grant-in-Aid for Scientific Research on Innovative Areas “Molecular Robotics” (Grant No. 15H00806)[27] Bender C. M. and Orszag S. A., Advanced Mathematial Methods for Scientists and Engineers (Springer,
New York) 1999.
[28] Chen L. Y. et al., Phys. Rev. E, 54 (1996) 376.
[29] Kunihiro T., Prog. Theor. Phys., 94 (1995) 503.
[30] Kunihiro T., Prog. Theor. Phys., 97 (1997) 179.
[31] Kučerka N. et al., Biochim. Biophys. Acta Biomembr.,
1808 (2011) 2761.
[32] Petrache H. I. et al., Biophys. J., 88 (2004) 1574.
[33] Tian A. et al., Phys. Rev. Lett., 98 (2007) 208102.

28002-p6