ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/mdRipple/mdripple.tex
(Generate patch)

Comparing trunk/mdRipple/mdripple.tex (file contents):
Revision 3174 by xsun, Fri Jul 13 22:01:52 2007 UTC vs.
Revision 3265 by xsun, Fri Oct 19 18:40:38 2007 UTC

# Line 1 | Line 1
1   %\documentclass[aps,pre,twocolumn,amssymb,showpacs,floatfix]{revtex4}
2 < \documentclass[aps,pre,preprint,amssymb,showpacs]{revtex4}
2 > %\documentclass[aps,pre,preprint,amssymb]{revtex4}
3 > \documentclass[12pt]{article}
4 > \usepackage{times}
5 > \usepackage{mathptm}
6 > \usepackage{tabularx}
7 > \usepackage{setspace}
8 > \usepackage{amsmath}
9 > \usepackage{amssymb}
10   \usepackage{graphicx}
11 + \usepackage[ref]{overcite}
12 + \pagestyle{plain}
13 + \pagenumbering{arabic}
14 + \oddsidemargin 0.0cm \evensidemargin 0.0cm
15 + \topmargin -21pt \headsep 10pt
16 + \textheight 9.0in \textwidth 6.5in
17 + \brokenpenalty=10000
18 + \renewcommand{\baselinestretch}{1.2}
19 + \renewcommand\citemid{\ } % no comma in optional reference note
20  
21   \begin{document}
22 < \renewcommand{\thefootnote}{\fnsymbol{footnote}}
23 < \renewcommand{\theequation}{\arabic{section}.\arabic{equation}}
22 > %\renewcommand{\thefootnote}{\fnsymbol{footnote}}
23 > %\renewcommand{\theequation}{\arabic{section}.\arabic{equation}}
24  
25 < %\bibliographystyle{aps}
25 > \bibliographystyle{achemso}
26  
27 < \title{Dipolar Ordering of the Ripple Phase in Lipid Membranes}
28 < \author{Xiuquan Sun and J. Daniel Gezelter}
29 < \email[E-mail:]{gezelter@nd.edu}
30 < \affiliation{Department of Chemistry and Biochemistry,\\
31 < University of Notre Dame, \\
27 > \title{Dipolar ordering in the ripple phases of molecular-scale models
28 > of lipid membranes}
29 > \author{Xiuquan Sun and J. Daniel Gezelter \\
30 > Department of Chemistry and Biochemistry,\\
31 > University of Notre Dame, \\
32   Notre Dame, Indiana 46556}
33  
34 + %\email[E-mail:]{gezelter@nd.edu}
35 +
36   \date{\today}
37  
38 < \begin{abstract}
38 > \maketitle
39  
40 + \begin{abstract}
41 + Symmetric and asymmetric ripple phases have been observed to form in
42 + molecular dynamics simulations of a simple molecular-scale lipid
43 + model. The lipid model consists of an dipolar head group and an
44 + ellipsoidal tail.  Within the limits of this model, an explanation for
45 + generalized membrane curvature is a simple mismatch in the size of the
46 + heads with the width of the molecular bodies.  The persistence of a
47 + {\it bilayer} structure requires strong attractive forces between the
48 + head groups.  One feature of this model is that an energetically
49 + favorable orientational ordering of the dipoles can be achieved by
50 + out-of-plane membrane corrugation.  The corrugation of the surface
51 + stabilizes the long range orientational ordering for the dipoles in the
52 + head groups which then adopt a bulk anti-ferroelectric state. We
53 + observe a common feature of the corrugated dipolar membranes: the wave
54 + vectors for the surface ripples are always found to be perpendicular
55 + to the dipole director axis.  
56   \end{abstract}
57  
58 < \pacs{}
59 < \maketitle
58 > %\maketitle
59 > \newpage
60  
61   \section{Introduction}
62   \label{sec:Int}
63 + Fully hydrated lipids will aggregate spontaneously to form bilayers
64 + which exhibit a variety of phases depending on their temperatures and
65 + compositions. Among these phases, a periodic rippled phase
66 + ($P_{\beta'}$) appears as an intermediate phase between the gel
67 + ($L_\beta$) and fluid ($L_{\alpha}$) phases for relatively pure
68 + phosphatidylcholine (PC) bilayers.  The ripple phase has attracted
69 + substantial experimental interest over the past 30 years. Most
70 + structural information of the ripple phase has been obtained by the
71 + X-ray diffraction~\cite{Sun96,Katsaras00} and freeze-fracture electron
72 + microscopy (FFEM).~\cite{Copeland80,Meyer96} Recently, Kaasgaard {\it
73 + et al.} used atomic force microscopy (AFM) to observe ripple phase
74 + morphology in bilayers supported on mica.~\cite{Kaasgaard03} The
75 + experimental results provide strong support for a 2-dimensional
76 + hexagonal packing lattice of the lipid molecules within the ripple
77 + phase.  This is a notable change from the observed lipid packing
78 + within the gel phase.~\cite{Cevc87} The X-ray diffraction work by
79 + Katsaras {\it et al.} showed that a rich phase diagram exhibiting both
80 + {\it asymmetric} and {\it symmetric} ripples is possible for lecithin
81 + bilayers.\cite{Katsaras00}
82  
83 < As one of the most important components in the formation of the
84 < biomembrane, lipid molecules attracted numerous studies in the past
85 < several decades. Due to their amphiphilic structure, when dispersed in
86 < water, lipids can self-assemble to construct a bilayer structure. The
87 < phase behavior of lipid membrane is well understood. The gel-fluid
88 < phase transition is known as main phase transition. However, there is
89 < an intermediate phase between gel and fluid phase for some lipid (like
90 < phosphatidycholine (PC)) membranes. This intermediate phase
91 < distinguish itself from other phases by its corrugated membrane
92 < surface, therefore this phase is termed ``ripple'' ($P_{\beta '}$)
93 < phase. The phase transition between gel-fluid and ripple phase is
94 < called pretransition. Since the pretransition usually occurs in room
95 < temperature, there might be some important biofuntions carried by the
96 < ripple phase for the living organism.
83 > A number of theoretical models have been presented to explain the
84 > formation of the ripple phase. Marder {\it et al.} used a
85 > curvature-dependent Landau-de~Gennes free-energy functional to predict
86 > a rippled phase.~\cite{Marder84} This model and other related continuum
87 > models predict higher fluidity in convex regions and that concave
88 > portions of the membrane correspond to more solid-like regions.
89 > Carlson and Sethna used a packing-competition model (in which head
90 > groups and chains have competing packing energetics) to predict the
91 > formation of a ripple-like phase.  Their model predicted that the
92 > high-curvature portions have lower-chain packing and correspond to
93 > more fluid-like regions.  Goldstein and Leibler used a mean-field
94 > approach with a planar model for {\em inter-lamellar} interactions to
95 > predict rippling in multilamellar phases.~\cite{Goldstein88} McCullough
96 > and Scott proposed that the {\em anisotropy of the nearest-neighbor
97 > interactions} coupled to hydrophobic constraining forces which
98 > restrict height differences between nearest neighbors is the origin of
99 > the ripple phase.~\cite{McCullough90} Lubensky and MacKintosh
100 > introduced a Landau theory for tilt order and curvature of a single
101 > membrane and concluded that {\em coupling of molecular tilt to membrane
102 > curvature} is responsible for the production of
103 > ripples.~\cite{Lubensky93} Misbah, Duplat and Houchmandzadeh proposed
104 > that {\em inter-layer dipolar interactions} can lead to ripple
105 > instabilities.~\cite{Misbah98} Heimburg presented a {\em coexistence
106 > model} for ripple formation in which he postulates that fluid-phase
107 > line defects cause sharp curvature between relatively flat gel-phase
108 > regions.~\cite{Heimburg00} Kubica has suggested that a lattice model of
109 > polar head groups could be valuable in trying to understand bilayer
110 > phase formation.~\cite{Kubica02} Bannerjee used Monte Carlo simulations
111 > of lamellar stacks of hexagonal lattices to show that large head groups
112 > and molecular tilt with respect to the membrane normal vector can
113 > cause bulk rippling.~\cite{Bannerjee02}
114  
115 < The ripple phase is observed experimentally by x-ray diffraction
116 < ~\cite{Sun96,Katsaras00}, freeze-fracture electron microscopy
117 < (FFEM)~\cite{Copeland80,Meyer96}, and atomic force microscopy (AFM)
118 < recently~\cite{Kaasgaard03}. The experimental studies suggest two
119 < kinds of ripple structures: asymmetric (sawtooth like) and symmetric
120 < (sinusoidal like) ripple phases. Substantial number of theoretical
121 < explaination applied on the formation of the ripple
122 < phases~\cite{Marder84,Carlson87,Goldstein88,McCullough90,Lubensky93,Misbah98,Heimburg00,Kubica02,Bannerjee02}.
123 < In contrast, few molecular modelling have been done due to the large
124 < size of the resulting structures and the time required for the phases
125 < of interest to develop. One of the interesting molecular simulations
126 < was carried out by De Vries and Marrink {\it et
127 < al.}~\cite{deVries05}. According to their dynamic simulation results,
128 < the ripple consists of two domains, one is gel bilayer, and in the
129 < other domain, the upper and lower leaves of the bilayer are fully
130 < interdigitated. The mechanism of the formation of the ripple phase in
131 < their work suggests the theory that the packing competition between
132 < head group and tail of lipid molecules is the driving force for the
133 < formation of the ripple phase~\cite{Carlson87}. Recently, the ripple
134 < phase is also studied by using monte carlo simulation~\cite{Lenz07},
135 < the ripple structure is similar to the results of Marrink except that
136 < the connection of the upper and lower leaves of the bilayer is an
137 < interdigitated line instead of the fully interdigitated
138 < domain. Furthermore, the symmetric ripple phase was also observed in
69 < their work. They claimed the mismatch between the size of the head
70 < group and tail of the lipid molecules is the driving force for the
71 < formation of the ripple phase.
115 > In contrast, few large-scale molecular modeling studies have been
116 > done due to the large size of the resulting structures and the time
117 > required for the phases of interest to develop.  With all-atom (and
118 > even unified-atom) simulations, only one period of the ripple can be
119 > observed and only for time scales in the range of 10-100 ns.  One of
120 > the most interesting molecular simulations was carried out by de~Vries
121 > {\it et al.}~\cite{deVries05}. According to their simulation results,
122 > the ripple consists of two domains, one resembling the gel bilayer,
123 > while in the other, the two leaves of the bilayer are fully
124 > interdigitated.  The mechanism for the formation of the ripple phase
125 > suggested by their work is a packing competition between the head
126 > groups and the tails of the lipid molecules.\cite{Carlson87} Recently,
127 > the ripple phase has also been studied by Lenz and Schmid using Monte
128 > Carlo simulations.\cite{Lenz07} Their structures are similar to the De
129 > Vries {\it et al.} structures except that the connection between the
130 > two leaves of the bilayer is a narrow interdigitated line instead of
131 > the fully interdigitated domain.  The symmetric ripple phase was also
132 > observed by Lenz {\it et al.}, and their work supports other claims
133 > that the mismatch between the size of the head group and tail of the
134 > lipid molecules is the driving force for the formation of the ripple
135 > phase. Ayton and Voth have found significant undulations in
136 > zero-surface-tension states of membranes simulated via dissipative
137 > particle dynamics, but their results are consistent with purely
138 > thermal undulations.~\cite{Ayton02}
139  
140 < Although the organizations of the tails of lipid molecules are
141 < addressed by these molecular simulations, the ordering of the head
142 < group in ripple phase is still not settlement. We developed a simple
143 < ``web of dipoles'' spin lattice model which provides some physical
144 < insight in our previous studies~\cite{Sun2007}, we found the dipoles
78 < on head groups of the lipid molecules are ordered in an
79 < antiferroelectric state. The similiar phenomenon is also observed by
80 < Tsonchev {\it et al.} when they studied the formation of the
81 < nanotube\cite{Tsonchev04}.
140 > Although the organization of the tails of lipid molecules are
141 > addressed by these molecular simulations and the packing competition
142 > between head groups and tails is strongly implicated as the primary
143 > driving force for ripple formation, questions about the ordering of
144 > the head groups in ripple phase have not been settled.
145  
146 < In this paper, we made a more realistic coarse-grained lipid model to
147 < understand the primary driving force for membrane corrugation and to
148 < elucidate the organization of the anisotropic interacting head group
149 < via molecular dynamics simulation. We will talk about our model and
150 < methodology in section \ref{sec:method}, and details of the simulation
151 < in section \ref{sec:experiment}. The results are shown in section
152 < \ref{sec:results}. At last, we will discuss the results in section
146 > In a recent paper, we presented a simple ``web of dipoles'' spin
147 > lattice model which provides some physical insight into relationship
148 > between dipolar ordering and membrane buckling.\cite{Sun2007} We found
149 > that dipolar elastic membranes can spontaneously buckle, forming
150 > ripple-like topologies.  The driving force for the buckling of dipolar
151 > elastic membranes is the anti-ferroelectric ordering of the dipoles.
152 > This was evident in the ordering of the dipole director axis
153 > perpendicular to the wave vector of the surface ripples.  A similar
154 > phenomenon has also been observed by Tsonchev {\it et al.} in their
155 > work on the spontaneous formation of dipolar peptide chains into
156 > curved nano-structures.\cite{Tsonchev04,Tsonchev04II}
157 >
158 > In this paper, we construct a somewhat more realistic molecular-scale
159 > lipid model than our previous ``web of dipoles'' and use molecular
160 > dynamics simulations to elucidate the role of the head group dipoles
161 > in the formation and morphology of the ripple phase.  We describe our
162 > model and computational methodology in section \ref{sec:method}.
163 > Details on the simulations are presented in section
164 > \ref{sec:experiment}, with results following in section
165 > \ref{sec:results}.  A final discussion of the role of dipolar heads in
166 > the ripple formation can be found in section
167   \ref{sec:discussion}.
168  
169 < \section{Methodology and Model}
169 > \section{Computational Model}
170   \label{sec:method}
171  
95 Our idea for developing a simple and reasonable lipid model to study
96 the ripple phase of lipid bilayers is based on two facts: one is that
97 the most essential feature of lipid molecules is their amphiphilic
98 structure with polar head groups and non-polar tails. Another fact is
99 that dominant numbers of lipid molecules are very rigid in ripple
100 phase which allows the details of the lipid molecules neglectable. The
101 lipid model is shown in Figure \ref{fig:lipidMM}. Figure
102 \ref{fig:lipidMM}a shows the molecular strucure of a DPPC molecule. The
103 hydrophilic character of the head group is the effect of the strong
104 dipole composed by a positive charge sitting on the nitrogen and a
105 negative charge on the phosphate. The hydrophobic tail consists of
106 fatty acid chains. In our model shown in Figure \ref{fig:lipidMM}b,
107 lipid molecules are represented by rigid bodies made of one head
108 sphere with a point dipole sitting on it and one ellipsoid tail, the
109 direction of the dipole is fixed to be perpendicular to the tail. The
110 breadth and length of tail are $\sigma_0$, $3\sigma_0$. The diameter
111 of heads varies from $1.20\sigma_0$ to $1.41\sigma_0$.  The model of
112 the solvent in our simulations is inspired by the idea of ``DPD''
113 water. Every four water molecules are reprsented by one sphere.
114
172   \begin{figure}[htb]
173   \centering
174 < \includegraphics[width=\linewidth]{lipidMM}
175 < \caption{The molecular structure of a DPPC molecule and the
176 < coars-grained model for PC molecules.\label{fig:lipidMM}}
174 > \includegraphics[width=4in]{lipidModels}
175 > \caption{Three different representations of DPPC lipid molecules,
176 > including the chemical structure, an atomistic model, and the
177 > head-body ellipsoidal coarse-grained model used in this
178 > work.\label{fig:lipidModels}}
179   \end{figure}
180  
181 < Spheres interact each other with Lennard-Jones potential
182 < \begin{eqnarray*}
183 < V_{ij} = 4\epsilon \left[\left(\frac{\sigma_0}{r_{ij}}\right)^{12} -
184 < \left(\frac{\sigma_0}{r_{ij}}\right)^6\right]
185 < \end{eqnarray*}
186 < here, $\sigma_0$ is diameter of the spherical particle. $r_{ij}$ is
187 < the distance between two spheres. $\epsilon$ is the well depth.
188 < Dipoles interact each other with typical dipole potential
189 < \begin{eqnarray*}
190 < V_{ij} = \frac{|\mu|^2}{4 \pi \epsilon_0 r_{ij}^3}
191 < \left[ \hat{\bf u}_i \cdot \hat{\bf u}_j - 3 (\hat{\bf u}_i \cdot
192 < \hat{\bf r}_{ij})(\hat{\bf u}_j \cdot \hat{\bf r}_{ij}) \right]
193 < \end{eqnarray*}
194 < In this potential, $\mathbf{\hat u}_i$ is the unit vector pointing
195 < along dipole $i$ and $\mathbf{\hat r}_{ij}$ is the unit vector
196 < pointing along the inter-dipole vector $\mathbf{r}_{ij}$. Identical
197 < ellipsoids interact each other with Gay-Berne potential.
198 < \begin{eqnarray*}
181 > Our simple molecular-scale lipid model for studying the ripple phase
182 > is based on two facts: one is that the most essential feature of lipid
183 > molecules is their amphiphilic structure with polar head groups and
184 > non-polar tails. Another fact is that the majority of lipid molecules
185 > in the ripple phase are relatively rigid (i.e. gel-like) which makes
186 > some fraction of the details of the chain dynamics negligible.  Figure
187 > \ref{fig:lipidModels} shows the molecular structure of a DPPC
188 > molecule, as well as atomistic and molecular-scale representations of
189 > a DPPC molecule.  The hydrophilic character of the head group is
190 > largely due to the separation of charge between the nitrogen and
191 > phosphate groups.  The zwitterionic nature of the PC headgroups leads
192 > to abnormally large dipole moments (as high as 20.6 D), and this
193 > strongly polar head group interacts strongly with the solvating water
194 > layers immediately surrounding the membrane.  The hydrophobic tail
195 > consists of fatty acid chains.  In our molecular scale model, lipid
196 > molecules have been reduced to these essential features; the fatty
197 > acid chains are represented by an ellipsoid with a dipolar ball
198 > perched on one end to represent the effects of the charge-separated
199 > head group.  In real PC lipids, the direction of the dipole is
200 > nearly perpendicular to the tail, so we have fixed the direction of
201 > the point dipole rigidly in this orientation.  
202 >
203 > The ellipsoidal portions of the model interact via the Gay-Berne
204 > potential which has seen widespread use in the liquid crystal
205 > community.  Ayton and Voth have also used Gay-Berne ellipsoids for
206 > modeling large length-scale properties of lipid
207 > bilayers.\cite{Ayton01} In its original form, the Gay-Berne potential
208 > was a single site model for the interactions of rigid ellipsoidal
209 > molecules.\cite{Gay81} It can be thought of as a modification of the
210 > Gaussian overlap model originally described by Berne and
211 > Pechukas.\cite{Berne72} The potential is constructed in the familiar
212 > form of the Lennard-Jones function using orientation-dependent
213 > $\sigma$ and $\epsilon$ parameters,
214 > \begin{equation*}
215   V_{ij}({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j}, {\mathbf{\hat
216   r}_{ij}}) = 4\epsilon ({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j},
217   {\mathbf{\hat r}_{ij}})\left[\left(\frac{\sigma_0}{r_{ij}-\sigma({\mathbf{\hat u}_i},
218   {\mathbf{\hat u}_j}, {\mathbf{\hat r}_{ij}})+\sigma_0}\right)^{12}
219   -\left(\frac{\sigma_0}{r_{ij}-\sigma({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j},
220   {\mathbf{\hat r}_{ij}})+\sigma_0}\right)^6\right]
221 < \end{eqnarray*}
222 < where $\sigma_0 = \sqrt 2\sigma_s$, and orietation-dependent range
223 < parameter is given by
224 < \begin{eqnarray*}
225 < \sigma({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j}, {\mathbf{\hat r}_{ij}}) =
226 < {\sigma_{0}} \left(1-\frac{1}{2}\chi\left[\frac{({\mathbf{\hat r}_{ij}}
227 < \cdot {\mathbf{\hat u}_i} + {\mathbf{\hat r}_{ij}} \cdot {\mathbf{\hat
228 < u}_j})^2}{1+\chi({\mathbf{\hat u}_i} \cdot {\mathbf{\hat u}_j})} +
229 < \frac{({\mathbf{\hat r}_{ij}} \cdot {\mathbf{\hat u}_i} - {\mathbf{\hat r}_{ij}}
230 < \cdot {\mathbf{\hat u}_j})^2}{1-\chi({\mathbf{\hat u}_i} \cdot
231 < {\mathbf{\hat u}_j})} \right] \right)^{-\frac{1}{2}}
232 < \end{eqnarray*}
233 < and the strength anisotropy function is,
234 < \begin{eqnarray*}
235 < \epsilon ({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j}, {\mathbf{\hat r}_{ij}}) =
236 < {\epsilon^\nu}({\mathbf{\hat u}_i}, {\mathbf{\hat
237 < u}_j}){\epsilon'^\mu}({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j},
238 < {\mathbf{\hat r}_{ij}})
164 < \end{eqnarray*}
165 < with $\nu$ and $\mu$ being adjustable exponent, and
166 < $\epsilon({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j})$,
167 < $\epsilon'({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j}, {\mathbf{\hat
168 < r}_{ij}})$ defined as
169 < \begin{eqnarray*}
170 < \epsilon({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j}) =
171 < \epsilon_0\left[1-\chi^2({\mathbf{\hat u}_i} \cdot {\mathbf{\hat
172 < u}_j})^2\right]^{-\frac{1}{2}}
221 > \label{eq:gb}
222 > \end{equation*}
223 >
224 > The range $(\sigma({\bf \hat{u}}_{i},{\bf \hat{u}}_{j},{\bf
225 > \hat{r}}_{ij}))$, and strength $(\epsilon({\bf \hat{u}}_{i},{\bf
226 > \hat{u}}_{j},{\bf \hat{r}}_{ij}))$ parameters
227 > are dependent on the relative orientations of the two molecules (${\bf
228 > \hat{u}}_{i},{\bf \hat{u}}_{j}$) as well as the direction of the
229 > intermolecular separation (${\bf \hat{r}}_{ij}$).  $\sigma$ and
230 > $\sigma_0$ are also governed by shape mixing and anisotropy variables,
231 > \begin {eqnarray*}
232 > \sigma_0 & = & \sqrt{d_i^2 + d_j^2} \\ \\
233 > \chi & = & \left[ \frac{\left( l_i^2 - d_i^2 \right) \left(l_j^2 -
234 > d_j^2 \right)}{\left( l_j^2 + d_i^2 \right) \left(l_i^2 +
235 > d_j^2 \right)}\right]^{1/2} \\ \\
236 > \alpha^2 & = & \left[ \frac{\left( l_i^2 - d_i^2 \right) \left(l_j^2 +
237 > d_i^2 \right)}{\left( l_j^2 - d_j^2 \right) \left(l_i^2 +
238 > d_j^2 \right)}\right]^{1/2},
239   \end{eqnarray*}
240 < \begin{eqnarray*}
241 < \epsilon'({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j}, {\mathbf{\hat r}_{ij}}) =
242 < 1-\frac{1}{2}\chi'\left[\frac{({\mathbf{\hat r}_{ij}} \cdot {\mathbf{\hat
243 < u}_i} + {\mathbf{\hat r}_{ij}} \cdot {\mathbf{\hat
244 < u}_j})^2}{1+\chi'({\mathbf{\hat u}_i} \cdot {\mathbf{\hat u}_j})} +
245 < \frac{({\mathbf{\hat r}_{ij}} \cdot {\mathbf{\hat u}_i} - {\mathbf{\hat r}_{ij}}
246 < \cdot {\mathbf{\hat u}_j})^2}{1-\chi'({\mathbf{\hat u}_i} \cdot
247 < {\mathbf{\hat u}_j})} \right]
240 > where $l$ and $d$ describe the length and width of each uniaxial
241 > ellipsoid.  These shape anisotropy parameters can then be used to
242 > calculate the range function,
243 > \begin{equation*}
244 > \sigma({\bf \hat{u}}_{i},{\bf \hat{u}}_{j},{\bf \hat{r}}_{ij}) = \sigma_{0}
245 > \left[ 1- \left\{ \frac{ \chi \alpha^2 ({\bf \hat{u}}_i \cdot {\bf
246 > \hat{r}}_{ij} ) + \chi \alpha^{-2} ({\bf \hat{u}}_j \cdot {\bf
247 > \hat{r}}_{ij} ) - 2 \chi^2 ({\bf \hat{u}}_i \cdot {\bf
248 > \hat{r}}_{ij} )({\bf \hat{u}}_j \cdot {\bf
249 > \hat{r}}_{ij} ) ({\bf \hat{u}}_i \cdot {\bf \hat{u}}_j)}{1 - \chi^2
250 > \left({\bf \hat{u}}_i \cdot {\bf \hat{u}}_j\right)^2} \right\}
251 > \right]^{-1/2}
252 > \end{equation*}
253 >
254 > Gay-Berne ellipsoids also have an energy scaling parameter,
255 > $\epsilon^s$, which describes the well depth for two identical
256 > ellipsoids in a {\it side-by-side} configuration.  Additionally, a well
257 > depth aspect ratio, $\epsilon^r = \epsilon^e / \epsilon^s$, describes
258 > the ratio between the well depths in the {\it end-to-end} and
259 > side-by-side configurations.  As in the range parameter, a set of
260 > mixing and anisotropy variables can be used to describe the well
261 > depths for dissimilar particles,
262 > \begin {eqnarray*}
263 > \epsilon_0 & = & \sqrt{\epsilon^s_i  * \epsilon^s_j} \\ \\
264 > \epsilon^r & = & \sqrt{\epsilon^r_i  * \epsilon^r_j} \\ \\
265 > \chi' & = & \frac{1 - (\epsilon^r)^{1/\mu}}{1 + (\epsilon^r)^{1/\mu}}
266 > \\ \\
267 > \alpha'^2 & = & \left[1 + (\epsilon^r)^{1/\mu}\right]^{-1}
268   \end{eqnarray*}
269 < the diameter dependent parameter $\chi$ is given by
270 < \begin{eqnarray*}
271 < \chi = \frac{({\sigma_s}^2 -
272 < {\sigma_e}^2)}{({\sigma_s}^2 + {\sigma_e}^2)}
273 < \end{eqnarray*}
274 < and the well depth dependent parameter $\chi'$ is given by
275 < \begin{eqnarray*}
276 < \chi' = \frac{({\epsilon_s}^{\frac{1}{\mu}} -
277 < {\epsilon_e}^{\frac{1}{\mu}})}{({\epsilon_s}^{\frac{1}{\mu}} +
278 < {\epsilon_e}^{\frac{1}{\mu}})}
279 < \end{eqnarray*}
280 < $\sigma_s$ is the side-by-side width, and $\sigma_e$ is the end-to-end
281 < length. $\epsilon_s$ is the side-by-side well depth, and $\epsilon_e$
282 < is the end-to-end well depth. For the interaction between
283 < nonequivalent uniaxial ellipsoids (in this case, between spheres and
284 < ellipsoids), the range parameter is generalized as\cite{Cleaver96}
285 < \begin{eqnarray*}
286 < \sigma({\mathbf{\hat u}_i}, {\mathbf{\hat u}_j}, {\mathbf{\hat r}_{ij}}) =
287 < {\sigma_{0}} \left(1-\frac{1}{2}\chi\left[\frac{(\alpha{\mathbf{\hat r}_{ij}}
288 < \cdot {\mathbf{\hat u}_i} + \alpha^{-1}{\mathbf{\hat r}_{ij}} \cdot {\mathbf{\hat
289 < u}_j})^2}{1+\chi({\mathbf{\hat u}_i} \cdot {\mathbf{\hat u}_j})} +
290 < \frac{(\alpha{\mathbf{\hat r}_{ij}} \cdot {\mathbf{\hat u}_i} - \alpha^{-1}{\mathbf{\hat r}_{ij}}
291 < \cdot {\mathbf{\hat u}_j})^2}{1-\chi({\mathbf{\hat u}_i} \cdot
292 < {\mathbf{\hat u}_j})} \right] \right)^{-\frac{1}{2}}
293 < \end{eqnarray*}
208 < where $\alpha$ is given by
209 < \begin{eqnarray*}
210 < \alpha^2 =
211 < \left[\frac{\left({\sigma_e}_i^2-{\sigma_s}_i^2\right)\left({\sigma_e}_j^2+{\sigma_s}_i^2\right)}{\left({\sigma_e}_j^2-{\sigma_s}_j^2\right)\left({\sigma_e}_i^2+{\sigma_s}_j^2\right)}
212 < \right]^{\frac{1}{2}}
213 < \end{eqnarray*}
214 < the strength parameter is adjusted by the suggestion of
215 < \cite{Cleaver96}. All potentials are truncated at $25$ \AA\ and
216 < shifted at $22$ \AA.
269 > The form of the strength function is somewhat complicated,
270 > \begin {eqnarray*}
271 > \epsilon({\bf \hat{u}}_{i}, {\bf \hat{u}}_{j},{\bf \hat{r}}_{ij}) & = &
272 > \epsilon_{0}  \epsilon_{1}^{\nu}({\bf \hat{u}}_{i}.{\bf \hat{u}}_{j})
273 > \epsilon_{2}^{\mu}({\bf \hat{u}}_{i},{\bf \hat{u}}_{j},{\bf
274 > \hat{r}}_{ij}) \\ \\
275 > \epsilon_{1}({\bf \hat{u}}_{i},{\bf \hat{u}}_{j}) & = &
276 > \left[1-\chi^{2}({\bf \hat{u}}_{i}.{\bf
277 > \hat{u}}_{j})^{2}\right]^{-1/2} \\ \\
278 > \epsilon_{2}({\bf \hat{u}}_{i},{\bf \hat{u}}_{j},{\bf \hat{r}}_{ij}) &
279 > = &
280 > 1 - \left\{ \frac{ \chi' \alpha'^2 ({\bf \hat{u}}_i \cdot {\bf
281 > \hat{r}}_{ij} ) + \chi' \alpha'^{-2} ({\bf \hat{u}}_j \cdot {\bf
282 > \hat{r}}_{ij} ) - 2 \chi'^2 ({\bf \hat{u}}_i \cdot {\bf
283 > \hat{r}}_{ij} )({\bf \hat{u}}_j \cdot {\bf
284 > \hat{r}}_{ij} ) ({\bf \hat{u}}_i \cdot {\bf \hat{u}}_j)}{1 - \chi'^2
285 > \left({\bf \hat{u}}_i \cdot {\bf \hat{u}}_j\right)^2} \right\},
286 > \end {eqnarray*}
287 > although many of the quantities and derivatives are identical with
288 > those obtained for the range parameter. Ref. \citen{Luckhurst90}
289 > has a particularly good explanation of the choice of the Gay-Berne
290 > parameters $\mu$ and $\nu$ for modeling liquid crystal molecules. An
291 > excellent overview of the computational methods that can be used to
292 > efficiently compute forces and torques for this potential can be found
293 > in Ref. \citen{Golubkov06}
294  
295 < \section{Experiment}
295 > The choices of parameters we have used in this study correspond to a
296 > shape anisotropy of 3 for the chain portion of the molecule.  In
297 > principle, this could be varied to allow for modeling of longer or
298 > shorter chain lipid molecules. For these prolate ellipsoids, we have:
299 > \begin{equation}
300 > \begin{array}{rcl}
301 > d & < & l \\
302 > \epsilon^{r} & < & 1
303 > \end{array}
304 > \end{equation}
305 > A sketch of the various structural elements of our molecular-scale
306 > lipid / solvent model is shown in figure \ref{fig:lipidModel}.  The
307 > actual parameters used in our simulations are given in table
308 > \ref{tab:parameters}.
309 >
310 > \begin{figure}[htb]
311 > \centering
312 > \includegraphics[width=4in]{2lipidModel}
313 > \caption{The parameters defining the behavior of the lipid
314 > models. $l / d$ is the ratio of the head group to body diameter.
315 > Molecular bodies had a fixed aspect ratio of 3.0.  The solvent model
316 > was a simplified 4-water bead ($\sigma_w \approx d$) that has been
317 > used in other coarse-grained (DPD) simulations.  The dipolar strength
318 > (and the temperature and pressure) were the only other parameters that
319 > were varied systematically.\label{fig:lipidModel}}
320 > \end{figure}
321 >
322 > To take into account the permanent dipolar interactions of the
323 > zwitterionic head groups, we have placed fixed dipole moments $\mu_{i}$ at
324 > one end of the Gay-Berne particles.  The dipoles are oriented at an
325 > angle $\theta = \pi / 2$ relative to the major axis.  These dipoles
326 > are protected by a head ``bead'' with a range parameter ($\sigma_h$) which we have
327 > varied between $1.20 d$ and $1.41 d$.  The head groups interact with
328 > each other using a combination of Lennard-Jones,
329 > \begin{equation}
330 > V_{ij}(r_{ij}) = 4\epsilon_h \left[\left(\frac{\sigma_h}{r_{ij}}\right)^{12} -
331 > \left(\frac{\sigma_h}{r_{ij}}\right)^6\right],
332 > \end{equation}
333 > and dipole-dipole,
334 > \begin{equation}
335 > V_{ij}({\bf \hat{u}}_{i},{\bf \hat{u}}_{j},{\bf
336 > \hat{r}}_{ij})) = \frac{|\mu|^2}{4 \pi \epsilon_0 r_{ij}^3}
337 > \left[ \hat{\bf u}_i \cdot \hat{\bf u}_j - 3 (\hat{\bf u}_i \cdot
338 > \hat{\bf r}_{ij})(\hat{\bf u}_j \cdot \hat{\bf r}_{ij}) \right]
339 > \end{equation}
340 > potentials.  
341 > In these potentials, $\mathbf{\hat u}_i$ is the unit vector pointing
342 > along dipole $i$ and $\mathbf{\hat r}_{ij}$ is the unit vector
343 > pointing along the inter-dipole vector $\mathbf{r}_{ij}$.  
344 >
345 > For the interaction between nonequivalent uniaxial ellipsoids (in this
346 > case, between spheres and ellipsoids), the spheres are treated as
347 > ellipsoids with an aspect ratio of 1 ($d = l$) and with an well depth
348 > ratio ($\epsilon^r$) of 1 ($\epsilon^e = \epsilon^s$).  The form of
349 > the Gay-Berne potential we are using was generalized by Cleaver {\it
350 > et al.} and is appropriate for dissimilar uniaxial
351 > ellipsoids.\cite{Cleaver96}
352 >
353 > The solvent model in our simulations is identical to one used by
354 > Marrink {\it et al.}  in their dissipative particle dynamics (DPD)
355 > simulation of lipid bilayers.\cite{Marrink04} This solvent bead is a
356 > single site that represents four water molecules (m = 72 amu) and has
357 > comparable density and diffusive behavior to liquid water.  However,
358 > since there are no electrostatic sites on these beads, this solvent
359 > model cannot replicate the dielectric properties of water.
360 >
361 > \begin{table*}
362 > \begin{minipage}{\linewidth}
363 > \begin{center}
364 > \caption{Potential parameters used for molecular-scale coarse-grained
365 > lipid simulations}
366 > \begin{tabular}{llccc}
367 > \hline
368 >  & &  Head & Chain & Solvent \\
369 > \hline
370 > $d$ (\AA) & & varied & 4.6  & 4.7 \\
371 > $l$ (\AA) & & $= d$ & 13.8 & 4.7 \\
372 > $\epsilon^s$ (kcal/mol) & & 0.185 & 1.0 & 0.8 \\
373 > $\epsilon_r$ (well-depth aspect ratio)& & 1 & 0.2 &  1 \\
374 > $m$ (amu) & & 196 & 760 & 72.06 \\
375 > $\overleftrightarrow{\mathsf I}$ (amu \AA$^2$) & & & & \\
376 > \multicolumn{2}c{$I_{xx}$} & 1125 & 45000 & N/A \\
377 > \multicolumn{2}c{$I_{yy}$} & 1125 & 45000 & N/A \\
378 > \multicolumn{2}c{$I_{zz}$} &  0 &    9000 & N/A \\
379 > $\mu$ (Debye) & & varied & 0 & 0 \\
380 > \end{tabular}
381 > \label{tab:parameters}
382 > \end{center}
383 > \end{minipage}
384 > \end{table*}
385 >
386 > \section{Experimental Methodology}
387   \label{sec:experiment}
388  
389 < To make the simulations less expensive and to observe long-time
390 < behavior of the lipid membranes, all simulations were started from two
391 < separate monolayers in the vaccum with $x-y$ anisotropic pressure
389 > The parameters that were systematically varied in this study were the
390 > size of the head group ($\sigma_h$), the strength of the dipole moment
391 > ($\mu$), and the temperature of the system.  Values for $\sigma_h$
392 > ranged from 5.5 \AA\ to 6.5 \AA\ .  If the width of the tails is taken
393 > to be the unit of length, these head groups correspond to a range from
394 > $1.2 d$ to $1.41 d$.  Since the solvent beads are nearly identical in
395 > diameter to the tail ellipsoids, all distances that follow will be
396 > measured relative to this unit of distance.  Because the solvent we
397 > are using is non-polar and has a dielectric constant of 1, values for
398 > $\mu$ are sampled from a range that is somewhat smaller than the 20.6
399 > Debye dipole moment of the PC head groups.
400 >
401 > To create unbiased bilayers, all simulations were started from two
402 > perfectly flat monolayers separated by a 26 \AA\ gap between the
403 > molecular bodies of the upper and lower leaves.  The separated
404 > monolayers were evolved in a vacuum with $x-y$ anisotropic pressure
405   coupling. The length of $z$ axis of the simulations was fixed and a
406   constant surface tension was applied to enable real fluctuations of
407 < the bilayer. Periodic boundaries were used. There were $480-720$ lipid
408 < molecules in the simulations depending on the size of the head
409 < beads. All the simulations were equlibrated for $100$ ns at $300$
410 < K. The resulting structures were solvated in water ($6$ DPD
411 < water/lipid molecule). These configurations were relaxed for another
412 < $30$ ns relaxation. All simulations with water were carried out at
413 < constant pressure ($P=1$ atm) by $3$D anisotropic coupling, and
414 < constant surface tension ($\gamma=0.015$). Given the absence of fast
415 < degrees of freedom in this model, a timestep of $50$ fs was
416 < utilized. Simulations were performed by using OOPSE
417 < package\cite{Meineke05}.
418 <
419 < \section{Results and Analysis}
407 > the bilayer. Periodic boundary conditions were used, and $480-720$
408 > lipid molecules were present in the simulations, depending on the size
409 > of the head beads.  In all cases, the two monolayers spontaneously
410 > collapsed into bilayer structures within 100 ps. Following this
411 > collapse, all systems were equilibrated for $100$ ns at $300$ K.
412 >
413 > The resulting bilayer structures were then solvated at a ratio of $6$
414 > solvent beads (24 water molecules) per lipid. These configurations
415 > were then equilibrated for another $30$ ns. All simulations utilizing
416 > the solvent were carried out at constant pressure ($P=1$ atm) with
417 > $3$D anisotropic coupling, and constant surface tension
418 > ($\gamma=0.015$ N/m). Given the absence of fast degrees of freedom in
419 > this model, a time step of $50$ fs was utilized with excellent energy
420 > conservation.  Data collection for structural properties of the
421 > bilayers was carried out during a final 5 ns run following the solvent
422 > equilibration.  All simulations were performed using the OOPSE
423 > molecular modeling program.\cite{Meineke05}
424 >
425 > A switching function was applied to all potentials to smoothly turn
426 > off the interactions between a range of $22$ and $25$ \AA.
427 >
428 > \section{Results}
429   \label{sec:results}
430  
431 < Snapshots in Figure \ref{fig:phaseCartoon} show that the membrane is
432 < more corrugated increasing size of the head groups. The surface is
433 < nearly flat when $\sigma_h=1.20\sigma_0$. With
434 < $\sigma_h=1.28\sigma_0$, although the surface is still flat, the
435 < bilayer starts to splay inward; the upper leaf of the bilayer is
436 < connected to the lower leaf with an interdigitated line defect. Two
437 < periodicities with $100$ \AA\ width were observed in the
438 < simulation. This structure is very similiar to the structure observed
439 < by de Vries and Lenz {\it et al.}. The same basic structure is also
440 < observed when $\sigma_h=1.41\sigma_0$, but the wavelength of the
441 < surface corrugations depends sensitively on the size of the ``head''
442 < beads. From the undulation spectrum, the corrugation is clearly
443 < non-thermal.
431 > The membranes in our simulations exhibit a number of interesting
432 > bilayer phases.  The surface topology of these phases depends most
433 > sensitively on the ratio of the size of the head groups to the width
434 > of the molecular bodies.  With heads only slightly larger than the
435 > bodies ($\sigma_h=1.20 d$) the membrane exhibits a flat bilayer.
436 >
437 > Increasing the head / body size ratio increases the local membrane
438 > curvature around each of the lipids.  With $\sigma_h=1.28 d$, the
439 > surface is still essentially flat, but the bilayer starts to exhibit
440 > signs of instability.  We have observed occasional defects where a
441 > line of lipid molecules on one leaf of the bilayer will dip down to
442 > interdigitate with the other leaf.  This gives each of the two bilayer
443 > leaves some local convexity near the line defect.  These structures,
444 > once developed in a simulation, are very stable and are spaced
445 > approximately 100 \AA\ away from each other.
446 >
447 > With larger heads ($\sigma_h = 1.35 d$) the membrane curvature
448 > resolves into a ``symmetric'' ripple phase.  Each leaf of the bilayer
449 > is broken into several convex, hemicylinderical sections, and opposite
450 > leaves are fitted together much like roof tiles.  There is no
451 > interdigitation between the upper and lower leaves of the bilayer.
452 >
453 > For the largest head / body ratios studied ($\sigma_h = 1.41 d$) the
454 > local curvature is substantially larger, and the resulting bilayer
455 > structure resolves into an asymmetric ripple phase.  This structure is
456 > very similar to the structures observed by both de~Vries {\it et al.}
457 > and Lenz {\it et al.}.  For a given ripple wave vector, there are two
458 > possible asymmetric ripples, which is not the case for the symmetric
459 > phase observed when $\sigma_h = 1.35 d$.
460 >
461   \begin{figure}[htb]
462   \centering
463 < \includegraphics[width=\linewidth]{phaseCartoon}
464 < \caption{A sketch to discribe the structure of the phases observed in
465 < our simulations.\label{fig:phaseCartoon}}
463 > \includegraphics[width=4in]{phaseCartoon}
464 > \caption{The role of the ratio between the head group size and the
465 > width of the molecular bodies is to increase the local membrane
466 > curvature.  With strong attractive interactions between the head
467 > groups, this local curvature can be maintained in bilayer structures
468 > through surface corrugation.  Shown above are three phases observed in
469 > these simulations.  With $\sigma_h = 1.20 d$, the bilayer maintains a
470 > flat topology.  For larger heads ($\sigma_h = 1.35 d$) the local
471 > curvature resolves into a symmetrically rippled phase with little or
472 > no interdigitation between the upper and lower leaves of the membrane.
473 > The largest heads studied ($\sigma_h = 1.41 d$) resolve into an
474 > asymmetric rippled phases with interdigitation between the two
475 > leaves.\label{fig:phaseCartoon}}
476   \end{figure}
477  
478 < When $\sigma_h=1.35\sigma_0$, we observed another corrugated surface
479 < morphology. This structure is different from the asymmetric rippled
480 < surface; there is no interdigitation between the upper and lower
481 < leaves of the bilayer. Each leaf of the bilayer is broken into several
482 < hemicylinderical sections, and opposite leaves are fitted together
483 < much like roof tiles. Unlike the surface in which the upper
484 < hemicylinder is always interdigitated on the leading or trailing edge
485 < of lower hemicylinder, the symmetric ripple has no prefered direction.
486 < The corresponding cartoons are shown in Figure
487 < \ref{fig:phaseCartoon} for elucidation of the detailed structures of
488 < different phases. Figure \ref{fig:phaseCartoon} (a) is the flat phase,
489 < (b) is the asymmetric ripple phase corresponding to the lipid
490 < organization when $\sigma_h=1.28\sigma_0$ and $\sigma_h=1.41\sigma_0$,
491 < and (c) is the symmetric ripple phase observed when
492 < $\sigma_h=1.35\sigma_0$. In symmetric ripple phase, the bilayer is
493 < continuous everywhere on the whole membrane, however, in asymmetric
494 < ripple phase, the bilayer is intermittent domains connected by thin
495 < interdigitated monolayer which consists of upper and lower leaves of
496 < the bilayer.
478 > Sample structures for the flat ($\sigma_h = 1.20 d$), symmetric
479 > ($\sigma_h = 1.35 d$, and asymmetric ($\sigma_h = 1.41 d$) ripple
480 > phases are shown in Figure \ref{fig:phaseCartoon}.  
481 >
482 > It is reasonable to ask how well the parameters we used can produce
483 > bilayer properties that match experimentally known values for real
484 > lipid bilayers.  Using a value of $l = 13.8$ \AA for the ellipsoidal
485 > tails and the fixed ellipsoidal aspect ratio of 3, our values for the
486 > area per lipid ($A$) and inter-layer spacing ($D_{HH}$) depend
487 > entirely on the size of the head bead relative to the molecular body.
488 > These values are tabulated in table \ref{tab:property}.  Kucera {\it
489 > et al.}  have measured values for the head group spacings for a number
490 > of PC lipid bilayers that range from 30.8 \AA (DLPC) to 37.8 (DPPC).
491 > They have also measured values for the area per lipid that range from
492 > 60.6
493 > \AA$^2$ (DMPC) to 64.2 \AA$^2$
494 > (DPPC).\cite{NorbertKucerka04012005,NorbertKucerka06012006} Only the
495 > largest of the head groups we modeled ($\sigma_h = 1.41 d$) produces
496 > bilayers (specifically the area per lipid) that resemble real PC
497 > bilayers.  The smaller head beads we used are perhaps better models
498 > for PE head groups.
499 >
500   \begin{table*}
501   \begin{minipage}{\linewidth}
502   \begin{center}
503 < \caption{}
504 < \begin{tabular}{lccc}
503 > \caption{Phase, bilayer spacing, area per lipid, ripple wavelength
504 > and amplitude observed as a function of the ratio between the head
505 > beads and the diameters of the tails.  Ripple wavelengths and
506 > amplitudes are normalized to the diameter of the tail ellipsoids.}
507 > \begin{tabular}{lccccc}
508   \hline
509 < $\sigma_h/\sigma_0$ & type of phase & $\lambda/\sigma_0$ & $A/\sigma_0$\\
509 > $\sigma_h / d$ & type of phase & bilayer spacing (\AA) & area per
510 > lipid (\AA$^2$) & $\lambda / d$ & $A / d$\\
511   \hline
512 < 1.20 & flat & N/A & N/A \\
513 < 1.28 & asymmetric flat & 21.7 & N/A \\
514 < 1.35 & symmetric ripple & 17.2 & 2.2 \\
515 < 1.41 & asymmetric ripple & 15.4 & 1.5 \\
512 > 1.20 & flat & 33.4 & 49.6 & N/A & N/A \\
513 > 1.28 & flat & 33.7 & 54.7 & N/A & N/A \\
514 > 1.35 & symmetric ripple & 42.9 & 51.7 & 17.2 & 2.2 \\
515 > 1.41 & asymmetric ripple & 37.1 & 63.1 & 15.4 & 1.5 \\
516   \end{tabular}
517   \label{tab:property}
518   \end{center}
519   \end{minipage}
520   \end{table*}
521  
522 < The membrane structures and the reduced wavelength $\lambda/\sigma_0$,
523 < reduced amplitude $A/\sigma_0$ of the ripples are summerized in Table
524 < \ref{tab:property}. The wavelength range is $15~21$ and the amplitude
525 < is $1.5$ for asymmetric ripple and $2.2$ for symmetric ripple. These
526 < values are consistent to the experimental results. Note, the
527 < amplitudes are underestimated without the melted tails in our
528 < simulations.
522 > The membrane structures and the reduced wavelength $\lambda / d$,
523 > reduced amplitude $A / d$ of the ripples are summarized in Table
524 > \ref{tab:property}. The wavelength range is $15 - 17$ molecular bodies
525 > and the amplitude is $1.5$ molecular bodies for asymmetric ripple and
526 > $2.2$ for symmetric ripple. These values are reasonably consistent
527 > with experimental measurements.\cite{Sun96,Katsaras00,Kaasgaard03}
528 > Note, that given the lack of structural freedom in the tails of our
529 > model lipids, the amplitudes observed from these simulations are
530 > likely to underestimate of the true amplitudes.
531  
532 < The $P_2$ order paramters (for molecular bodies and head group
533 < dipoles) have been calculated to clarify the ordering in these phases
534 < quantatively. $P_2=1$ means a perfectly ordered structure, and $P_2=0$
535 < implies orientational randomization. Figure \ref{fig:rP2} shows the
536 < $P_2$ order paramter of the dipoles on head group rising with
537 < increasing head group size. When the heads of the lipid molecules are
538 < small, the membrane is flat. The dipolar ordering is essentially
539 < frustrated on orientational ordering in this circumstance. Another
540 < reason is that the lipids can move independently in each monolayer, it
541 < is not nessasory for the direction of dipoles on one leaf is
542 < consistant to another layer, which makes total order parameter is
543 < relatively low. With increasing head group size, the surface is
318 < corrugated, and dipoles do not move as freely on the
319 < surface. Therefore, the translational freedom of lipids in one layer
320 < is dependent upon the position of lipids in another layer, as a
321 < result, the symmetry of the dipoles on head group in one layer is tied
322 < to the symmetry in the other layer. Furthermore, as the membrane
323 < deforms from two to three dimensions due to the corrugation, the
324 < symmetry of the ordering for the dipoles embedded on each leaf is
325 < broken. The dipoles then self-assemble in a head-tail configuration,
326 < and the order parameter increases dramaticaly. However, the total
327 < polarization of the system is still close to zero. This is strong
328 < evidence that the corrugated structure is an antiferroelectric
329 < state. The orientation of the dipolar is always perpendicular to the
330 < ripple wave vector. These results are consistent with our previous
331 < study on dipolar membranes.
532 > \begin{figure}[htb]
533 > \centering
534 > \includegraphics[width=4in]{topDown}
535 > \caption{Top views of the flat (upper), symmetric ripple (middle),
536 > and asymmetric ripple (lower) phases.  Note that the head-group
537 > dipoles have formed head-to-tail chains in all three of these phases,
538 > but in the two rippled phases, the dipolar chains are all aligned {\it
539 > perpendicular} to the direction of the ripple.  Note that the flat
540 > membrane has multiple vortex defects in the dipolar ordering, and the
541 > ordering on the lower leaf of the bilayer can be in an entirely
542 > different direction from the upper leaf.\label{fig:topView}}
543 > \end{figure}
544  
545 < The ordering of the tails is essentially opposite to the ordering of
546 < the dipoles on head group. The $P_2$ order parameter decreases with
547 < increasing head size. This indicates the surface is more curved with
548 < larger head groups. When the surface is flat, all tails are pointing
549 < in the same direction; in this case, all tails are parallel to the
550 < normal of the surface,(making this structure remindcent of the
551 < $L_{\beta}$ phase. Increasing the size of the heads, results in
545 > The principal method for observing orientational ordering in dipolar
546 > or liquid crystalline systems is the $P_2$ order parameter (defined
547 > as $1.5 \times \lambda_{max}$, where $\lambda_{max}$ is the largest
548 > eigenvalue of the matrix,
549 > \begin{equation}
550 > {\mathsf{S}} = \frac{1}{N} \sum_i \left(
551 > \begin{array}{ccc}
552 >        u^{x}_i u^{x}_i-\frac{1}{3} & u^{x}_i u^{y}_i & u^{x}_i u^{z}_i \\
553 >        u^{y}_i u^{x}_i  & u^{y}_i u^{y}_i -\frac{1}{3} & u^{y}_i u^{z}_i \\
554 >        u^{z}_i u^{x}_i & u^{z}_i u^{y}_i  & u^{z}_i u^{z}_i -\frac{1}{3}
555 > \end{array} \right).
556 > \label{eq:opmatrix}
557 > \end{equation}
558 > Here $u^{\alpha}_i$ is the $\alpha=x,y,z$ component of the unit vector
559 > for molecule $i$.  (Here, $\hat{\bf u}_i$ can refer either to the
560 > principal axis of the molecular body or to the dipole on the head
561 > group of the molecule.)  $P_2$ will be $1.0$ for a perfectly-ordered
562 > system and near $0$ for a randomized system.  Note that this order
563 > parameter is {\em not} equal to the polarization of the system.  For
564 > example, the polarization of a perfect anti-ferroelectric arrangement
565 > of point dipoles is $0$, but $P_2$ for the same system is $1$.  The
566 > eigenvector of $\mathsf{S}$ corresponding to the largest eigenvalue is
567 > familiar as the director axis, which can be used to determine a
568 > privileged axis for an orientationally-ordered system.  Since the
569 > molecular bodies are perpendicular to the head group dipoles, it is
570 > possible for the director axes for the molecular bodies and the head
571 > groups to be completely decoupled from each other.
572 >
573 > Figure \ref{fig:topView} shows snapshots of bird's-eye views of the
574 > flat ($\sigma_h = 1.20 d$) and rippled ($\sigma_h = 1.35, 1.41 d$)
575 > bilayers.  The directions of the dipoles on the head groups are
576 > represented with two colored half spheres: blue (phosphate) and yellow
577 > (amino).  For flat bilayers, the system exhibits signs of
578 > orientational frustration; some disorder in the dipolar head-to-tail
579 > chains is evident with kinks visible at the edges between differently
580 > ordered domains.  The lipids can also move independently of lipids in
581 > the opposing leaf, so the ordering of the dipoles on one leaf is not
582 > necessarily consistent with the ordering on the other.  These two
583 > factors keep the total dipolar order parameter relatively low for the
584 > flat phases.
585 >
586 > With increasing head group size, the surface becomes corrugated, and
587 > the dipoles cannot move as freely on the surface. Therefore, the
588 > translational freedom of lipids in one layer is dependent upon the
589 > position of the lipids in the other layer.  As a result, the ordering of
590 > the dipoles on head groups in one leaf is correlated with the ordering
591 > in the other leaf.  Furthermore, as the membrane deforms due to the
592 > corrugation, the symmetry of the allowed dipolar ordering on each leaf
593 > is broken. The dipoles then self-assemble in a head-to-tail
594 > configuration, and the dipolar order parameter increases dramatically.
595 > However, the total polarization of the system is still close to zero.
596 > This is strong evidence that the corrugated structure is an
597 > anti-ferroelectric state.  It is also notable that the head-to-tail
598 > arrangement of the dipoles is always observed in a direction
599 > perpendicular to the wave vector for the surface corrugation.  This is
600 > a similar finding to what we observed in our earlier work on the
601 > elastic dipolar membranes.\cite{Sun2007}
602 >
603 > The $P_2$ order parameters (for both the molecular bodies and the head
604 > group dipoles) have been calculated to quantify the ordering in these
605 > phases.  Figure \ref{fig:rP2} shows that the $P_2$ order parameter for
606 > the head-group dipoles increases with increasing head group size. When
607 > the heads of the lipid molecules are small, the membrane is nearly
608 > flat. Since the in-plane packing is essentially a close packing of the
609 > head groups, the head dipoles exhibit frustration in their
610 > orientational ordering.
611 >
612 > The ordering trends for the tails are essentially opposite to the
613 > ordering of the head group dipoles. The tail $P_2$ order parameter
614 > {\it decreases} with increasing head size. This indicates that the
615 > surface is more curved with larger head / tail size ratios. When the
616 > surface is flat, all tails are pointing in the same direction (normal
617 > to the bilayer surface).  This simplified model appears to be
618 > exhibiting a smectic A fluid phase, similar to the real $L_{\beta}$
619 > phase.  We have not observed a smectic C gel phase ($L_{\beta'}$) for
620 > this model system.  Increasing the size of the heads results in
621   rapidly decreasing $P_2$ ordering for the molecular bodies.
622 +
623   \begin{figure}[htb]
624   \centering
625   \includegraphics[width=\linewidth]{rP2}
626 < \caption{The $P_2$ order parameter as a funtion of the ratio of
627 < $\sigma_h$ to $\sigma_0$.\label{fig:rP2}}
626 > \caption{The $P_2$ order parameters for head groups (circles) and
627 > molecular bodies (squares) as a function of the ratio of head group
628 > size ($\sigma_h$) to the width of the molecular bodies ($d$). \label{fig:rP2}}
629   \end{figure}
630  
631 < We studied the effects of the interactions between head groups on the
632 < structure of lipid bilayer by changing the strength of the dipole.
633 < Figure \ref{fig:sP2} shows how the $P_2$ order parameter changes with
634 < increasing strength of the dipole. Generally the dipoles on the head
635 < group are more ordered by increase in the strength of the interaction
636 < between heads and are more disordered by decreasing the interaction
637 < stength. When the interaction between the heads is weak enough, the
638 < bilayer structure does not persist; all lipid molecules are solvated
639 < directly in the water. The critial value of the strength of the dipole
640 < depends on the head size. The perfectly flat surface melts at $5$
641 < debye, the asymmetric rippled surfaces melt at $8$ debye, the
642 < symmetric rippled surfaces melt at $10$ debye. The ordering of the
643 < tails is the same as the ordering of the dipoles except for the flat
644 < phase. Since the surface is already perfect flat, the order parameter
645 < does not change much until the strength of the dipole is $15$
646 < debye. However, the order parameter decreases quickly when the
647 < strength of the dipole is further increased. The head groups of the
648 < lipid molecules are brought closer by stronger interactions between
649 < them. For a flat surface, a large amount of free volume between the
650 < head groups is available, but when the head groups are brought closer,
651 < the tails will splay outward, forming an inverse micelle. For rippled
652 < surfaces, there is less free volume available between the head
653 < groups. Therefore there is little effect on the structure of the
654 < membrane due to increasing dipolar strength. Unlike other systems that
655 < melt directly when the interaction is weak enough, for
656 < $\sigma_h=1.41\sigma_0$, part of the membrane melts into itself
657 < first. The upper leaf of the bilayer becomes totally interdigitated
658 < with the lower leaf. This is different behavior than what is exhibited
659 < with the interdigitated lines in the rippled phase where only one
660 < interdigitated line connects the two leaves of bilayer.
631 > In addition to varying the size of the head groups, we studied the
632 > effects of the interactions between head groups on the structure of
633 > lipid bilayer by changing the strength of the dipoles.  Figure
634 > \ref{fig:sP2} shows how the $P_2$ order parameter changes with
635 > increasing strength of the dipole.  Generally, the dipoles on the head
636 > groups become more ordered as the strength of the interaction between
637 > heads is increased and become more disordered by decreasing the
638 > interaction strength.  When the interaction between the heads becomes
639 > too weak, the bilayer structure does not persist; all lipid molecules
640 > become dispersed in the solvent (which is non-polar in this
641 > molecular-scale model).  The critical value of the strength of the
642 > dipole depends on the size of the head groups.  The perfectly flat
643 > surface becomes unstable below $5$ Debye, while the  rippled
644 > surfaces melt at $8$ Debye (asymmetric) or $10$ Debye (symmetric).
645 >
646 > The ordering of the tails mirrors the ordering of the dipoles {\it
647 > except for the flat phase}. Since the surface is nearly flat in this
648 > phase, the order parameters are only weakly dependent on dipolar
649 > strength until it reaches $15$ Debye.  Once it reaches this value, the
650 > head group interactions are strong enough to pull the head groups
651 > close to each other and distort the bilayer structure. For a flat
652 > surface, a substantial amount of free volume between the head groups
653 > is normally available.  When the head groups are brought closer by
654 > dipolar interactions, the tails are forced to splay outward, first forming
655 > curved bilayers, and then inverted micelles.
656 >
657 > When $\sigma_h=1.28 d$, the $P_2$ order parameter decreases slightly
658 > when the strength of the dipole is increased above $16$ Debye. For
659 > rippled bilayers, there is less free volume available between the head
660 > groups. Therefore increasing dipolar strength weakly influences the
661 > structure of the membrane.  However, the increase in the body $P_2$
662 > order parameters implies that the membranes are being slightly
663 > flattened due to the effects of increasing head-group attraction.
664 >
665 > A very interesting behavior takes place when the head groups are very
666 > large relative to the molecular bodies ($\sigma_h = 1.41 d$) and the
667 > dipolar strength is relatively weak ($\mu < 9$ Debye). In this case,
668 > the two leaves of the bilayer become totally interdigitated with each
669 > other in large patches of the membrane.   With higher dipolar
670 > strength, the interdigitation is limited to single lines that run
671 > through the bilayer in a direction perpendicular to the ripple wave
672 > vector.
673 >
674   \begin{figure}[htb]
675   \centering
676   \includegraphics[width=\linewidth]{sP2}
677 < \caption{The $P_2$ order parameter as a funtion of the strength of the
678 < dipole.\label{fig:sP2}}
677 > \caption{The $P_2$ order parameters for head group dipoles (a) and
678 > molecular bodies (b) as a function of the strength of the dipoles.
679 > These order parameters are shown for four values of the head group /
680 > molecular width ratio ($\sigma_h / d$). \label{fig:sP2}}
681   \end{figure}
682  
683 < Figure \ref{fig:tP2} shows the dependence of the order parameter on
684 < temperature. The behavior of the $P_2$ order paramter is
685 < straightforward. Systems are more ordered at low temperature, and more
686 < disordered at high temperatures. When the temperature is high enough,
687 < the membranes are instable. Since our model lacks the detailed
688 < information on lipid tails, we can not simulate the fluid phase with
689 < melted fatty acid chains. Moreover, the formation of the tilted
690 < $L_{\beta'}$ phase also depends on the organization of fatty groups on
691 < tails.
683 > Figure \ref{fig:tP2} shows the dependence of the order parameters on
684 > temperature.  As expected, systems are more ordered at low
685 > temperatures, and more disordered at high temperatures.  All of the
686 > bilayers we studied can become unstable if the temperature becomes
687 > high enough.  The only interesting feature of the temperature
688 > dependence is in the flat surfaces ($\sigma_h=1.20 d$ and
689 > $\sigma_h=1.28 d$).  Here, when the temperature is increased above
690 > $310$K, there is enough jostling of the head groups to allow the
691 > dipolar frustration to resolve into more ordered states.  This results
692 > in a slight increase in the $P_2$ order parameter above this
693 > temperature.
694 >
695 > For the rippled surfaces ($\sigma_h=1.35 d$ and $\sigma_h=1.41 d$),
696 > there is a slightly increased orientational ordering in the molecular
697 > bodies above $290$K.  Since our model lacks the detailed information
698 > about the behavior of the lipid tails, this is the closest the model
699 > can come to depicting the ripple ($P_{\beta'}$) to fluid
700 > ($L_{\alpha}$) phase transition.  What we are observing is a
701 > flattening of the rippled structures made possible by thermal
702 > expansion of the tightly-packed head groups.  The lack of detailed
703 > chain configurations also makes it impossible for this model to depict
704 > the ripple to gel ($L_{\beta'}$) phase transition.
705 >
706   \begin{figure}[htb]
707   \centering
708   \includegraphics[width=\linewidth]{tP2}
709 < \caption{The $P_2$ order parameter as a funtion of
710 < temperature.\label{fig:tP2}}
709 > \caption{The $P_2$ order parameters for head group dipoles (a) and
710 > molecular bodies (b) as a function of temperature.
711 > These order parameters are shown for four values of the head group /
712 > molecular width ratio ($\sigma_h / d$).\label{fig:tP2}}
713   \end{figure}
714  
715 + Fig. \ref{fig:phaseDiagram} shows a phase diagram for the model as a
716 + function of the head group / molecular width ratio ($\sigma_h / d$)
717 + and the strength of the head group dipole moment ($\mu$).  Note that
718 + the specific form of the bilayer phase is governed almost entirely by
719 + the head group / molecular width ratio, while the strength of the
720 + dipolar interactions between the head groups governs the stability of
721 + the bilayer phase.  Weaker dipoles result in unstable bilayer phases,
722 + while extremely strong dipoles can shift the equilibrium to an
723 + inverted micelle phase when the head groups are small.   Temperature
724 + has little effect on the actual bilayer phase observed, although higher
725 + temperatures can cause the unstable region to grow into the higher
726 + dipole region of this diagram.
727 +
728 + \begin{figure}[htb]
729 + \centering
730 + \includegraphics[width=\linewidth]{phaseDiagram}
731 + \caption{Phase diagram for the simple molecular model as a function
732 + of the head group / molecular width ratio ($\sigma_h / d$) and the
733 + strength of the head group dipole moment
734 + ($\mu$).\label{fig:phaseDiagram}}
735 + \end{figure}
736 +
737 + We have computed translational diffusion coefficients for lipid
738 + molecules from the mean square displacement,
739 + \begin{eqnarray}
740 + \langle {|\left({\bf r}_{i}(t) - {\bt r}_{i}(0) \right)|}^2 \rangle \\ \\
741 + & = & 6Dt
742 + \end{eqnarray}
743 + of the lipid bodies. The values of the translational diffusion
744 + coefficient for different head-to-tail size ratio are shown in table
745 + \ref{tab:relaxation}.
746 +
747 + We have also computed orientational diffusion constants for the head
748 + groups from the relaxation of the second-order Legendre polynomial
749 + correlation function,
750 + \begin{eqnarray}
751 + C_{\ell}(t) & = & \langle P_{\ell}\left({\bf \mu}_{i}(t) \cdot {\bf
752 + \mu}_{i}(0) \right) \rangle  \\ \\
753 + & \approx & e^{-\ell(\ell + 1) \theta t},
754 + \end{eqnarray}
755 + of the head group dipoles.  In this last line, we have used a simple
756 + ``Debye''-like model for the relaxation of the correlation function,
757 + specifically in the case when $\ell = 2$.   The computed orientational
758 + diffusion constants are given in table \ref{tab:relaxation}.  The
759 + notable feature we observe is that the orientational diffusion
760 + constant for the head group exhibits an order of magnitude decrease
761 + upon entering the rippled phase.  Our orientational correlation times
762 + are substantially in excess of those provided by...
763 +
764 +
765 + \begin{table*}
766 + \begin{minipage}{\linewidth}
767 + \begin{center}
768 + \caption{Rotational diffusion constants for the head groups
769 + ($\theta_h$) and molecular bodies ($\theta_b$) as well as the
770 + translational diffusion coefficients for the molecule as a function of
771 + the head-to-body width ratio.  The orientational mobility of the head
772 + groups experiences an {\it order of magnitude decrease} upon entering
773 + the rippled phase, which suggests that the rippling is tied to a
774 + freezing out of head group orientational freedom.  Uncertainties in
775 + the last digit are indicated by the values in parentheses.}
776 + \begin{tabular}{lccc}
777 + \hline
778 + $\sigma_h / d$ & $\theta_h (\mu s^{-1})$ & $\theta_b (1/fs)$ & $D (
779 + \times 10^{-11} m^2 s^{-1} \\
780 + \hline
781 + 1.20 & $0.206(1) $ & $0.0175(5) $ & $0.43(1)$ \\
782 + 1.28 & $0.179(2) $ & $0.055(2)  $ & $5.91(3)$ \\
783 + 1.35 & $0.025(1) $ & $0.195(3)  $ & $3.42(1)$ \\
784 + 1.41 & $0.023(1) $ & $0.024(3)  $ & $7.16(1)$ \\
785 + \end{tabular}
786 + \label{tab:relaxation}
787 + \end{center}
788 + \end{minipage}
789 + \end{table*}
790 +
791   \section{Discussion}
792   \label{sec:discussion}
793  
794 + Symmetric and asymmetric ripple phases have been observed to form in
795 + our molecular dynamics simulations of a simple molecular-scale lipid
796 + model. The lipid model consists of an dipolar head group and an
797 + ellipsoidal tail.  Within the limits of this model, an explanation for
798 + generalized membrane curvature is a simple mismatch in the size of the
799 + heads with the width of the molecular bodies.  With heads
800 + substantially larger than the bodies of the molecule, this curvature
801 + should be convex nearly everywhere, a requirement which could be
802 + resolved either with micellar or cylindrical phases.
803 +
804 + The persistence of a {\it bilayer} structure therefore requires either
805 + strong attractive forces between the head groups or exclusionary
806 + forces from the solvent phase.  To have a persistent bilayer structure
807 + with the added requirement of convex membrane curvature appears to
808 + result in corrugated structures like the ones pictured in
809 + Fig. \ref{fig:phaseCartoon}.  In each of the sections of these
810 + corrugated phases, the local curvature near a most of the head groups
811 + is convex.  These structures are held together by the extremely strong
812 + and directional interactions between the head groups.
813 +
814 + Dipolar head groups are key for the maintaining the bilayer structures
815 + exhibited by this model.  The dipoles are likely to form head-to-tail
816 + configurations even in flat configurations, but the temperatures are
817 + high enough that vortex defects become prevalent in the flat phase.
818 + The flat phase we observed therefore appears to be substantially above
819 + the Kosterlitz-Thouless transition temperature for a planar system of
820 + dipoles with this set of parameters.  For this reason, it would be
821 + interesting to observe the thermal behavior of the flat phase at
822 + substantially lower temperatures.
823 +
824 + One feature of this model is that an energetically favorable
825 + orientational ordering of the dipoles can be achieved by forming
826 + ripples.  The corrugation of the surface breaks the symmetry of the
827 + plane, making vortex defects somewhat more expensive, and stabilizing
828 + the long range orientational ordering for the dipoles in the head
829 + groups.  Most of the rows of the head-to-tail dipoles are parallel to
830 + each other and the system adopts a bulk anti-ferroelectric state.  We
831 + believe that this is the first time the organization of the head
832 + groups in ripple phases has been addressed.
833 +
834 + Although the size-mismatch between the heads and molecular bodies
835 + appears to be the primary driving force for surface convexity, the
836 + persistence of the bilayer through the use of rippled structures is a
837 + function of the strong, attractive interactions between the heads.
838 + One important prediction we can make using the results from this
839 + simple model is that if the dipole-dipole interaction is the leading
840 + contributor to the head group attractions, the wave vectors for the
841 + ripples should always be found {\it perpendicular} to the dipole
842 + director axis.  This echoes the prediction we made earlier for simple
843 + elastic dipolar membranes, and may suggest experimental designs which
844 + will test whether this is really the case in the phosphatidylcholine
845 + $P_{\beta'}$ phases.  The dipole director axis should also be easily
846 + computable for the all-atom and coarse-grained simulations that have
847 + been published in the literature.\cite{deVries05}
848 +
849 + Although our model is simple, it exhibits some rich and unexpected
850 + behaviors.  It would clearly be a closer approximation to reality if
851 + we allowed bending motions between the dipoles and the molecular
852 + bodies, and if we replaced the rigid ellipsoids with ball-and-chain
853 + tails.  However, the advantages of this simple model (large system
854 + sizes, 50 fs time steps) allow us to rapidly explore the phase diagram
855 + for a wide range of parameters.  Our explanation of this rippling
856 + phenomenon will help us design more accurate molecular models for
857 + corrugated membranes and experiments to test whether or not
858 + dipole-dipole interactions exert an influence on membrane rippling.
859 + \newpage
860   \bibliography{mdripple}
861   \end{document}

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines