---
abstract: 'We study the elastic response of bilayer membranes with fixed projected area to both stretching and shape deformations. A surface tension is associated to each of these deformations. By using model amphiphilic membranes and computer simulations, we are able to observe both the types of deformation, and thus, both the surface tensions, related to each type of deformation, are measured for the same system. These surface tensions are found to assume different values in the same bilayer membrane: in particular they vanish for different values of the projected area. We introduce a simple theory which relates the two quantities and successfully apply it to the data obtained with computer simulations.'
author:
- Alberto Imparato
title: Surface tension in bilayer membranes with fixed projected area
---

Bilayer membranes are composed of amphiphilic molecules whose hydrophobic part is strongly insoluble in water. Such molecules tends to form interfaces with the solvent in order to reduce the interfacial energy [@Isr]: a bilayer is thus formed by two adjacent sheets of amphiphilic molecules separating two aqueous phases. The interest in such structures resides in the fact that they play an outstanding role in the organization of biological cells: bilayer membranes of lipids form the basic structure of a cell’s plasma membrane and of internal membranes, which surround the organelles in eukaryotic cells, such as the nucleus, mitochondria and the Golgi apparatus [@ABL; @L1].

In a fluid membrane fluctuating freely in a solvent, all internal degrees of freedom, related, e.g, to the hydrocarbon chain conformation, and to to the local molecular density, relax on a fast time scale. The membrane is then characterized by a zero shear modulus, i.e. any shear deformation induces a flow of the amphiphilic molecules within the membrane. Such a system is then subject to only two types of elastic deformations: bending deformations and stretching deformations. The first are deformations normal to the membrane plane which change the membrane shapes, while the latter are in-plane deformations which change the local area per molecule projected onto the membrane midplane.

The surface of a fluctuating membrane can thus be viewed as an interface with the solvent, and we will call [*fluctuation*]{} surface tension the surface tension of the bilayer membrane associated to shape fluctuations. The surface tension of an interface saturated by surfactant molecules is expected to vanish (see, e.g., [@Saf]), and thus the fluctuation surface tension of a fluid membrane fluctuating freely in a solvent is usually assumed to be zero. However in the case of a constrained bilayer membrane, such as those fixed on a frame, the geometrical constrain can lead to a non vanishing surface tension. Non vanishing fluctuation surface tension have been observed, e.g., in [@MM].

On a macroscopic length scale, a membrane fixed on a frame can be treated as an elastic sheet: if one of the sides of the frame is movable and one tries to change the frame area by pulling or pushing on this side, an elastic force will appear which will tend to restore the preferred value of the frame area. The parameter describing such elastic behaviour is the [*mechanical*]{} surface tension and can assume non negative values too [@GL1; @LE2].

Intuitively, one would expect that the two surface tensions associated with the two different membrane deformations, namely shape deformation and lateral stretching, are different quantities, see ref. [@DL] and discussion therein. However, in experiments, they are usually measured by observing the one or the other type of deformation, and so one cannot compare them experimentally for the same system. On the contrary, the two surface tensions can be measured using computer simulations of model membranes, and this has been done in two previous works [@FP0; @Ste]. But in one case, the two surface tensions have been found to be equal within the statistical errors [@FP0], while in the other case they have been found to be proportional with respect to each other [@Ste].

The aim of this paper is thus the systematical study and comparison between the fluctuation surface tension and the mechanical surface tension, measured by observing both the types of deformation in the same bilayer membrane. In the present paper we consider bilayer membrane with fixed projected area. The projected area correspond to the area of the frame in the picture drawn above. By using model bilayer membrane and computer simulations, we measure both the surface tension associated to shape deformation and that associated to stretching deformations. We show, for the first time, that in the same system these two quantities assume different values, and in particular they do not vanish for the same value of the projected area. The model used here for the computer simulations, have been largely used in previous works and has been quite successful in describing many dynamic and elastic properties of real bilayer membranes [@GL1; @GL2; @alb1; @alb2].

The paper is organized as follows: in section \[model\], we describe the model bilayer, together with the computer simulation techniques which we use in the present work. In section \[mis\_S\] we discuss the theory of the mechanical surface tension and the method adopted to measure it. We then show the results for such quantity obtained by computer simulations. In section \[mis\_s\], after reviewing the classical elasticity theory for the shape fluctuations of bilayer membranes, we show the results for the shape fluctuation surface tension, obtained by computer simulations. In section \[compare\_sec\], the two quantities are compared, and a simple theory which connects them is introduced and discussed. We discuss our results and conclude in section \[concl\].

In the following, we consider a bilayer membrane made up of $N$ amphiphilic molecules, at temperature $T$, spanning a square frame of area $A_p=L^2$, its effective area being $A$.

Model bilayers and computer simulations {#model}
=======================================

Coarse grained model of amphiphilic bilayers
--------------------------------------------

We adopt the same coarse grained model introduced and used in refs.[@GL1; @GL2; @alb1; @alb2], which turned out to be an effective model to study several properties of bilayer membranes such as surface tension, bending rigidity and diffusion characteristics. The amphiphilic molecules are represented as linear chains of beads, a single bead representing the molecule hydrophilic head (H), or several $\mathrm{CH}_2/ \mathrm{CH}_3$ groups of the amphiphile hydrocarbon chain (C) see fig. \[mol\]. The water molecule is also represented as a single bead (W). The hydrophobic interaction of C with W and H particles is modelled by soft core potential $$U{_{SC}}({r_{ij}})=4 {\epsilon_{\alpha \beta}}{\left({\frac{ \ell'}{{r_{ij}}}}\right)}^9 \, ,
\label{usc}$$ while attractive interactions W-W, W-H, H-H, C-C, are modelled by a Lennard-Jones potential $$U{_{LJ}}({r_{ij}})=4 {\epsilon_{\alpha \beta}}{\left[{ {\left({\frac{\ell}{{r_{ij}}}}\right)}^{12}-{\left({\frac{\ell}{{r_{ij}}}}\right)}^{6}}\right]} \, ,
\label{ulj}$$ where $\alpha/\beta\, \epsilon \, \{ \text W, \text H, \text C \} $. Adjacent beads along a single model molecule interact via the harmonic potential $$U_2(r_{i,i+1})=k_2{\left({r_{i,i+1}-\ell}\right)}^2\, ,
\label{u2}$$ where $i$ and $i+1$ indicate two successive particles along the chain. The effects of hydrocarbon chain stiffness is modelled by letting all the particles within a single amphiphile molecule interact via the three-body bending potential $$U_3({\mathbf{r}}_{i-1,i},{\mathbf{r}}_{i,i+1})=k_3{\left({1-\frac{{\mathbf{r}}_{i-1,i}\cdot {\mathbf{r}}_{i,i+1}}{r_{i-1,i}\, r_{i,i+1}}}\right)}=k_3(1-\cos \phi_i) \, ,
\label{u3}$$ where ${\mathbf{r}}_{i,i+1}={\mathbf{r}}_{i+1}-{\mathbf{r}}_i$, see fig. \[mol\].

Simulation method and model parameters
--------------------------------------

In the present work we combine both Monte-Carlo (MC) and Molecular Dynamics (MD) simulation methods. Our simulations are carried out for a cuboidal boxes of constant volume and periodic boundary conditions. The MC algorithm is used to make the system relax towards configurations of minimal potential energy, the output of $500\times N$ MC steps are used as starting configurations for the MD part of the simulations.

For the parameters characterizing the above potentials we choose the values used in refs. [@GL1; @GL2; @alb1; @alb2]. For the interaction ranges of the LJ potentials we take $\ell=1/3$ nm which is of the same order as the LJ lengths of interactions for pairs of ${\text{CH}_2}$ groups, ${\text{CH}_3}$ groups or water molecules, as discussed in [@EB]. The characteristic length of the repulsive SC potential (\[u3\]), is taken to be $\ell'=1.05 \ell$ as in [@GL1]: with this choice the hard-core repulsion of the SC potential is of the same order of the repulsive part of the LJ potential. We take $\epsilon=2/N_{A_{V}}$ kJ ($N_{A_{V}}$ is the Avogadro number): this value is bigger than the LJ energies for pairs of ${\text{CH}_2}$ and/or water molecules reported in [@EB], but it takes into account that in our model one C particle corresponds to three or four ${\text{CH}_2}$ groups [@GL1]. The strength of the harmonic bond potential (\[u2\]) and of the three-body bending potentials (\[u3\]), are taken to be $k_2=5000 \epsilon/\ell^2$ and $k_3=2 \epsilon$ [@GL1].

Since we run MD simulations, two additional parameters are involved in the present model: the masses of the different beads, which enter the equations of motions, and the time step $\delta t$, used to discretize such equations. For simplicity’s sake, all the beads are taken to have the same mass $m=0.036/N_{A_{V}}$ kg, which is approximately the mean value of the mass of a water molecule and the mass of four ${\text{CH}_2}$ groups. Using this value of $m$ one obtains the characteristic time scale $\tau=\sqrt{m \ell^2 /\epsilon}=1.4$ ps. In the MD simulations, the integration time step $\delta t$ is taken to be $\delta t=\tau/2000=0.7$ fs. The [*leap-frog*]{} algorithm at constant temperature [@AT] with $k_B T=1.35 \epsilon$ (corresponding to $T=325$ K) is used for the integration of the equations of motions. In the following, if not differently specified, all the quantities will be expressed in units of the tree basic parameters $\ell$, $\epsilon$ and $m$.

In the present paper we consider bilayers with three different values of the number of amphiphilic molecules $N=512$, $N=768$ and $N=1152$. For each value of $N$, the simulations are run with different values of the projected area $A_p$, which corresponds to the area of the simulation box side parallel to the bilayer. In changing the projected area, we keep constant the number of water particles $N_W$, as well as the overall simulation box volume $V$, i.e., we keep the overall system density constant. The values of $N$, $N_W$ and $V$ used in the simulations described in the following are reported in table \[tab1\].

[|c|c|c|]{} $N$& $N_W$& $V(\ell^3)$\
512 & 3200 & 8640\
768 & 4800& 12960\
1152 & 7200 & 19440\

Stress tensor and mechanical surface tension {#mis_S}
============================================

On a macroscopic length scale, a membrane fixed on a frame can be viewed as an elastic sheet: the system will be characterized by some equilibrium value of the frame area where the force exerted on the frame sizes vanishes, and small changes in the frame area cause restoring forces to appear. If $F$ is the free energy of a framed bilayer, the surface tension $\Sigma$ conjugated to its projected area $A_p$ is defined as $$\Sigma={\frac{\partial F}{\partial A_p}}.
\label{defSigma}$$

In the present work, the bilayer surface tension is measured using a method first introduced by Schofield and Henderson [@SH], and then extended by Goetz and Lipowsky [@GL1]. A fluctuating membranes can be considered isotropic and homogeneous along any direction parallel to the membrane plane, if the amplitude of the fluctuations is not too high. Thus, the system stress tensor will be a function of the coordinate $z$ perpendicular to the bilayer plane. The surface tension is related to the system stress tensor via $$\Sigma=\int {{\mathrm d}}z {\left[{\Sigma_T(z) -\Sigma_N(z)}\right]},
\label{defS}$$ where $\Sigma_T(z)$ and $\Sigma_N(z)$ are the components of the stress tensor perpendicular and normal to the bilayer surface, respectively. The integral can be extended to infinity because the stress tensor in isotropic in the bulk water and so $\Sigma_T(z)=\Sigma_N(z)$. By noting that the stress tensor is defined as the negative of the pressure tensor, one sees that eq. (\[defS\]) is equivalent to the original expression of the surface tension of a planar liquid-vapor interface, as formulated by Kirkwood and Buff [@KB]. The macroscopic stress tensor $\mathbf \Sigma$ can be expressed in terms of the microscopic stress $\mathbf s$ which depends on the momenta and on the positions of the system particles within a small volume. The microscopic and the macroscopic stress tensor are related via $${\mathbf \Sigma}={\left\langle  \mathbf s \right\rangle},
\label{medias}$$ where the brackets denote ensemble average. The microscopic stress tensor ${\mathbf s}$ consists of two contributions one arising from the kinetic energy and the other from the interaction energy of the system particles within a small volume ${\mathbf s}={\mathbf s}_K+{\mathbf s}_{in}$. The kinetic contribution ${\mathbf s}_K$ to the macroscopic stress tensor, as given by eq. (\[medias\]), can be neglected, since it is isotropic, and must vanish on average. Thus, one is left with the interaction part of the microscopic stress tensor ${\mathbf s}_{in}$ only, and eq. (\[medias\]) becomes $${\mathbf \Sigma}={\left\langle  \mathbf s \right\rangle}_{in},
\label{mediasin}$$ Furthermore a planar bilayer membrane in solvent can be considered homogeneous along any direction parallel to the the bilayer plane, and so is the stress tensor. The system can thus be divided in thin slices parallel to the bilayer plane, and the stress tensor can be averaged over each of these slices. Thus the microscopic stress tensor is given by [@GL1; @SH; @LE2] $${\mathbf s}_{in}(z)=f(z_i,z_j,z) \frac{1}{\Delta V}\sum_{i>j} {\mathbf F}_{ij}\otimes {\mathbf r}_{ij}, 
\label{defsz}$$ where the sum is extended to all the system particle pairs, $\Delta V$ is the volume of the slice, the function $f(z_i,z_j,z)$ determines the actual contribution of the pair $i,j$ to the stress tensor in the slice of coordinate $z$, and the symbol $\otimes$ denotes the tensorial product. If both the particles are inside the current slice we take $f(z_i,z_j,z)=1$ . Let $\Delta z$ be the thickness of each slice, if both the particles are external to the current slice and on opposite sides, we take $f(z_i,z_j,z)=\Delta z/|z_i-z_j|$, where, using periodic boundary conditions, the shortest distance $|z_i-z_j|$ between them is considered. If just one of the two particles lies within the current slice, we take $f(z_i,z_j,z)=\Delta z/(2|z_i-z_j|)$. In all the other cases, the contribution of the particles $i,j$ to the stress tensor associated to the slice of coordinate $z$ vanishes, and thus we take $f(z_i,z_j,z)=0$. Inspection of eq. (\[defsz\]) suggests that its rhs can be viewed as a local density of the macroscopic virial tensor ${\mathbf W}=\sum_{i>j} {\mathbf F}_{ij}\otimes {\mathbf r}_{ij}$, and thus the function $f(z_i,z_j,z)$ determines the fraction of virial tensor density to be associated to the slice of coordinate $z$.

Measurements of the microscopic stress tensor
---------------------------------------------

Using the model amphiphilic bilayer described in section \[model\], the microscopic stress tensor $\mathbf s_{in}(z)$, as given by eq. (\[defsz\]), was measured for the three values of the number of amphiphilic molecules $N$ here considered, $N=512$, $N=768$, and $N=1152$. For each value of $N$, $\mathbf s_{in}(z)$ was measured for different values of the bilayer projected area $A_p$ keeping constant the total volume $V$, see table \[tab1\]. The simulation boxes have been divided into 140 $z$-slices, which corresponds to $\Delta z\sim 0.1\ell$, for the simulation box dimensions here used. The off-diagonal elements of the system stress tensor are found to vanish as expected, as well as the average of its kinetic part (data not shown). For each value of $A_p$, the mechanical surface tension $\Sigma$ is obtained using eq. (\[defS\]), by performing a discrete integration. For each value of $N$ and $A_p$ we run 5 simulations of $2\times 10^5$ MD steps, and a run of $6\times 10^6$ MC steps was interposed after each run of $2\times 10^5$ MD steps. The MC steps are inserted to allow the system to reach regions of its phase space not easily accessible with a single MD trajectory.

It is worth noting that the code for the stress tensor measuring is extremely time-consuming: a single run of $2\times 10^5$ MD steps, for a system with $N=1152$ amphiphiles, needs about one month of machine time, on a computer equipped with a 1 GHz processor. The relative fluctuations of the surface tension with respect to its average value are very large, as found in other work [@GL1]. Following ref. [@GL1], in order to reduce the effect of short time fluctuations, the total simulation time of each run is divided in blocks of 5000 MD steps, and the surface tension is subaveraged over each of these blocks. This average value over the 5000 MD steps is then used as the sample value for the current block.

The surface tension obtained from simulations is plotted in fig. \[compare\_S\], as a function of the projected area per amphiphile $a_p=A_p/(N/2)$, for the three values of $N$ here considered.

According to the classical elasticity theory, the parameter describing the effective compressibility of an elastic sheet to a change in its projected area is the area compressibility, $K_A$ defined as $$K_A=A_p{\frac{\partial \Sigma}{\partial A_p}} .$$ By integrating this last equation, one obtains $$\Sigma=K_A\ln{\left({\frac{A_p}{A^\dagger_p}}\right)}\simeq K_A\frac{A_p-A^\dagger_p}{A^\dagger_p},
\label{SigKA}$$ where the last equality holds for values of $A_p$ close to $A^\dagger_p$, which is the area at which $\Sigma$ vanishes. Upon integration of eq. (\[SigKA\]), we obtain the classical expression for the free energy of a stretched (or compressed) elastic sheet as a a function of $A_p$ around the equilibrium value $A^\dagger_p$ $$F=\frac{K_A}{2 A^\dagger_p} {\left({A_p-A^\dagger_p}\right)}^2,$$ where the reference free energy at $A=A^\dagger_p$ has been chosen to be equal to zero. Inspection of figure \[compare\_S\] indicates that the surface tension $\Sigma$ obtained by simulations follows the behaviour predicted by eq. (\[SigKA\]): this quantity is linear within a range of values around the equilibrium area $A^\dagger_p$. We find that for $A_p\ll A^\dagger_p$, ($a_p\simeq2.065$) the surface tension $\Sigma$ is no longer a monotonous function of $a_p$: this is probably due to the fact that, for such small projected area per molecule, the system exhibits buckling, and the membrane cannot be considered flat on average. In this case, the basic hypothesis that the membrane is isotropic and homogeneous along the the $z$-axis is not fulfilled. However, such hypothesis is required for measuring $\Sigma$ with the method discussed in this section, and thus the results obtained for very small values of $a_p$ must be inaccurate. Thus in the following we will consider values of $a_p$ such that $\Sigma$ is a monotonic function of $a_p$, i.e. $a_p>2.065$. The values of the equilibrium area, of the area compressibility and of the slope of $\Sigma$ as a function of $a_p$, for the three system sizes here considered, are listed in table \[tabS\]. In the same table, we also report the results for a smaller system which was considered in a previous work [@GL1].

[|c|c|c|c|]{} $N$& $a^\dagger_p\, (\ell^2)$& $K_A\, (\epsilon/\ell^2)$& $\partial \Sigma/ \partial a_p\, (\epsilon/\ell^3)$\
128 & $2.15\pm0.02$ & $11.8\pm1.5$& $5.5\pm 0.7$\
512 & $2.12\pm0.01$ & $11.7\pm1.5$& $5.5 \pm 0.7$\
768 & $ 2.12 \pm 0.01$ & $9.4 \pm 1.1$ & $4.5 \pm 0.5$\
1152 & $2.12\pm0.01$& $8.8\pm0.8$& $4.1 \pm 0.4$\

The results listed in table \[tabS\] suggests that the equilibrium area $A^\dagger_p$ does not depend on the system size, but is an intrinsic properties of the model amphiphilic molecule used here. However, the area compressibility decreases as a function of the system size, as already found in [@MM]: this decrease is due to the fact that a larger number of oscillation modes are introduced in the system as its projected area is increased. Thus, the bilayer becomes easier to compress or to stretch if its total projected area $A_p$ increases while the projected area per amphiphile $a_p$ is kept constant. Note that the values of the area compressibility found for the present model range between 263 and 360 $mJ/m^2$, see table \[tabS\]. Such values are slightly greater than those observed for real lipid bilayers, that are in the range 140-240 $mJ/m^2$ [@Ma], but are very similar to those values found in other works on computer simulations of amphiphilic bilayer with atomic resolution [@MM; @LE2]

Elasticity Hamiltonian and the fluctuation surface tension. {#mis_s}
===========================================================

In the Monge representation, the shape of the bilayer is described by the height function $h(x,y)$ which measures the distance of its midsurface from the reference $(x,y)$ plane. The classical elasticity Hamiltonian of fluctuating membranes, for small fluctuations, reads [@Can; @Helf] $$\begin{aligned}
H&=&\sigma A_p +\int_{A_p} {{\mathrm d}}x {{\mathrm d}}y{\left[{ \frac 1 2  \sigma (\nabla h)^2 + 
 \frac 1 2 \kappa (\nabla^2 h)^2}\right]}\nonumber \\
  &=&\sigma A_p + \frac{1} {2 L^2} \sum _{{\mathbf{q}}}{\left[{\sigma q^2+ \kappa q^4}\right]} |{\tilde h_{\mathbf{q}}}|^2 \, ,
\label{fapp}\end{aligned}$$ where ${\tilde h_{\mathbf{q}}}$ define the Fourier coefficients of $h(x,y)$ $${\tilde h_{\mathbf{q}}}=\frac  {1} {A_p}\int_{A_p} {{\mathrm d}}x {{\mathrm d}}y\, \exp{\left[{\mathrm{i} {\mathbf{r}}\cdot {\mathbf{q}}}\right]} h(x,y).$$ The two parameters appearing in eq. (\[fapp\]) are the [*bending rigidity*]{} $\kappa$, describing the resistance of the system to bending, and the surface tension $\sigma$, which takes into account the contribution of the bilayer total area $A$ to the system energy. In fact, eq. (\[fapp\]) can be rewritten as $$H=\sigma A+ \int {{\mathrm d}}x {{\mathrm d}}y \frac 1 2 \kappa (\nabla^2 h)^2.$$ The membrane fluctuation spectrum, defined as $S(q)\equiv|{\tilde h_{\mathbf{q}}}|^2$, depends only on $q=|{\mathbf{q}}|$ and exhibits the functional form [@LB] $$S(q)=\frac{k_B T}{A_p {\left({\sigma q^2 +\kappa q^4}\right)}}.
\label{fluspe}$$ Equation (\[fluspe\]) has been obtained under the hypothesis that the bilayer can be modelled as a geometrical surface, neglecting its thickness. Thus eq. (\[fluspe\]) holds only for wavenumbers smaller than an upper bond ${{q^*}}$, which is usually taken to be ${{q^*}}\simeq 2 \pi /d_0$, where $d_0$ is the bilayer thickness. We find that the model bilayer thickness is $d_0\simeq 6 \ell$, independently of the system size. For $q>q^*$, the shape of the bilayer is characterized by the local protrusions of the amphiphiles, and thus the fluctuation spectrum is dominated by a local surface tension term, see refs [@Ste; @GL2; @alb2]. However these sub-optical modes are not accessible in experiments, and thus, in the present work, we will focus on the small wavenumber, and on the effective values of $\kappa$ and $\sigma$ that can be obtained by eq. (\[fluspe\]). By sampling the height field $h(x,y)$ during the simulations, one can obtain the mean fluctuation spectrum (\[fluspe\]). By fitting $S(q)$ to the rhs of eq. (\[fluspe\]), for $q\lesssim q^*$, the values of the bending rigidity $\kappa$ and of the surface tension $\sigma$ can thus be estimated. The fluctuation spectrum shape depends on the projected area: and thus both $\kappa$ and $\sigma$ depend on this parameter. By changing the value of $A_p$ one can achieve the tensionless ($\sigma=0$) state [@GL2; @alb2]. The surface tension has been measured for the three system sizes here considered $N=521, \, 768,\, 1152$, for different projected area, and keeping constant the overall volume, see table \[tab1\]. For each value of the projected area, we run $10\times 10^5$ MD steps: a run of $6\times 10^6$ MC steps has been interposed after each run of $10^5$ steps. The fluctuation spectrum has been sampled during the MD simulations, every $5000$ MD steps. In fig.\[sq2\] the fluctuation spectrum of the larger system here considered is plotted as a function of the wavenumber $q$ for different value of the projected area per amphiphile $a_p$: the surface tension $\sigma$ vanishes for a given values of $a_p=a^*_p$, which depends on the microscopic details of the system (on the model parameters in the case of simulations), while for $a_p>a_p^*$ ($a_p<a_p^*$) it is greater (smaller) than zero. This is clearly shown in figure \[compare\_s\], where the surface tension $\sigma$ is plotted as a function of $a_p$ for the three values of $N$ here considered.

By linear fit of the data in figure \[compare\_s\], we obtain the tensionless area $a^*_p$ and the slope of the surface tension $\sigma$ as a function of $a_p$, the results are listed in table \[tabsig\]. As in the case of the mechanical surface tension $\Sigma$, the slope of the curve $\sigma(a_p)$ decreases as $N$ is increased. Comparison of tables \[tabS\] and \[tabsig\] indicates that the slopes of $\Sigma$ and $\sigma$ as functions of $a_p$ are different for a given system size.

[|c|c|c|]{} $N$& $a^\dagger_p\, (\ell^2)$& $\partial \sigma/ \partial a_p\, (\epsilon/\ell^3)$\
512 & $2.085\pm0.005$ & $5.4 \pm 0.7$\
768 & $ 2.09 \pm 0.01$ & $4.0 \pm 0.4$\
1152 & $2.09\pm0.01$& $3.6 \pm 0.4$\

Comparison between the two surface tension $\sigma$ and $\Sigma$ {#compare_sec}
================================================================

Comparison of figure \[compare\_s\] with figure \[compare\_S\] shows clearly that the projected area per amphiphile $a_p^*\simeq2.09$ where the surface tension $\sigma$ vanishes is different from the projected area per amphiphile $a_p^\dagger\simeq2.12$, at which the mechanical tension $\Sigma$ vanishes. Since these two quantities have been measured independently, for [*each*]{} value of the projected area, this last result suggests that one is dealing with two different quantities. In addition, the two vanishing areas $a_p^*$ and $a_p^\dagger$ do not exhibit any dependence on the system size $N$, and thus one can argue that they are two independent intrinsic properties of the amphiphilic molecule.

Comparison of tables \[tabS\] and \[tabsig\], and inspection of figure \[compare\_sS\], clearly indicate that the curves $\Sigma(a_p)$ and $\sigma(a_p)$ exhibit different slopes. Furthermore, if one plots on the same figure $\sigma$ and $\Sigma$ as functions of $a_p$, in the range where $\Sigma$ is linear with respect to $a_p$, the two data sets appear clearly shifted the one respect to the other, see figure \[compare\_sS\], where these two quantities have been plotted for the $N=1152$ case.

In order to relate the two surface tension $\sigma$ and $\Sigma$, we need first to define the statistical ensemble which characterizes the system considered here, i.e., the minimal set of thermodynamical variables which fully describe the system state. The classical argument that the system area readjusts to its preferred value (see, e.g., [@Saf]) cannot be invoked here because of the constrain represented by the frame area (simulation box area). The effective bilayer area $A$ can be evaluated using a discrete version of the formula $$A=\int_{A_p} {{\mathrm d}}x{{\mathrm d}}y \sqrt{1+{\left({\nabla h}\right)}^2},
\label{defA}$$ where the field $h(x,y)$ is sampled during the MD simulations. The effective area per amphiphile $a$ is plotted in fig. \[effa\] as a function of $a_p$. Inspection of this figure clearly indicates that the effective area $a$ is not independent of the projected area $a_p$, but is rather a function of $a_p$. Thus, the thermodynamical ensemble characterizing the framed bilayer here considered is the $(T,A_p=L^2, N)$ ensemble, where $L$ is the size of the simulation box face parallel to the bilayer.

The free energy in the $(T,A_p=L^2,N)$ ensemble for a bilayer membrane system is $$\begin{aligned}
F&=&-k_B T\ln Z=-k_B T\ln{\left({\int {\mathcal D} h\,  e^{-\beta H(\{h\})} }\right)}\nonumber\\
 &\simeq& -k_B T\ln{\left({ \prod_{i=0}^{N-1}\int \frac{ d h_i}{\lambda} e^{-\beta H(\{h_i\})} }\right)}\nonumber \\
       &=&-k_B T\ln {\left({ \prod_{{\mathbf{q}}}\int  \frac{d h_{{\mathbf{q}}}}{\sqrt N a_p\lambda} e^{-\beta H(\{h_q\})} }\right)}\nonumber \\
    &=&\sigma A_p-k_B T\ln {\left[{ \prod_{{\mathbf{q}}}\int  \frac{d h_{{\mathbf{q}}}}{\sqrt N a_p\lambda}
 e^{-  \frac \beta{2 L^2} \sum_{{\mathbf{q}}}(\sigma q^2+ \kappa q^4)  |h_{\mathbf{q}}|^2 } }\right]}\nonumber \\
    &=&\sigma A_p-k_B T \ln {\left[{ \prod_{{\mathbf{q}}} \frac{2\pi}{\lambda^2 a_p \beta (\sigma q^2+ \kappa q^4)} }\right]}^{\frac 1 2}\nonumber \\
    &=&\sigma A_p+{\frac{1}{2}}k_B T\sum_{\mathbf{q}}\ln {\left[{\frac{\lambda^2 a_p (\sigma q^2+ \kappa q^4)}{2\pi k_B T} }\right]}\, .
\label{free}\end{aligned}$$ The parameter $\lambda$ has the dimension of a length and is the analogous of the De Broglie wave-length in the ideal gas partition function. Without introducing this parameter, the partition function would be dimensionful. The equation of state for the intrinsic area is $$\begin{aligned}
{\left\langle  A \right\rangle}&=&{\frac{\partial F(T,A_p)}{\partial \sigma}} =A_p+\frac 1 2 k_B T \sum_{\mathbf{q}}\frac 1 {\sigma+\kappa q^2},
\label{area}
$$ which is a well known result, see, e.g, [@Rou]. Let us now introduce the quantity $N'$, which represents the number of membrane patches which fluctuate independently ($N'\le N/2$): the summations in equations (\[free\]) and (\[area\]) run over these $N'$ wavemodes, and thus we have $\sum_{{\mathbf{q}}}=N'$. Note that, in the large $A_p=L^2$ limit, $N'$ can be estimated as follows $$\begin{aligned}
N'&=& \sum_{{\mathbf{q}}}\simeq {\left({\frac{L}{2 \pi}}\right)}^2 \int {{\mathrm d}}q^2 =\frac{L^2}{2 \pi} \int^{{q^*}}_{{q_{min}}}q {{\mathrm d}}q\nonumber \\
&=& \frac{L^2}{4 \pi} {\left({{{q^*}}^2-{{q_{min}}}^2}\right)} \simeq \frac{A_p}{4 \pi}, 
\label{evaln1}\end{aligned}$$ where we have taken ${{q^*}}\simeq 1$ (see discussion is section \[mis\_s\]), and ${{q^*}}\gg{{q_{min}}}\simeq 2 \pi/L\sim 0$, in the large system size limit. Taking into account that $q^2=(2  \pi)^2 (n_x^2+n_y^2)/L^2$, and since the surface tension $\Sigma$ is the thermodynamic conjugate of $A_p=L^2$, as defined by eq. (\[defSigma\]), we have $$\begin{aligned}
\Sigma(T,A_p)&=&{\frac{\partial F}{\partial A_p}}=\sigma -{\frac{1}{2}}\frac{k_B T}{A_p} \sum_{{\mathbf{q}}} \frac{\kappa q^4}{\kappa q^4+\sigma q^2}\\
      &=&\sigma +{\frac{1}{2}}\frac{k_B T}{A_p} \sum_{{\mathbf{q}}} \frac{\sigma }{\kappa q^2+\sigma } -\frac{k_B T}{2} \frac{N'}{A_p}\label{Sigma2}\\
      &=& \frac{{\left\langle  A \right\rangle}}{A_p} \sigma -\frac{k_B T}{2} \frac{N'}{A_p}\label{Sigma4}\, ,\end{aligned}$$ where we substituted eq. (\[area\]) into eq. (\[Sigma2\]). Equation (\[Sigma4\]) relates thus the frame surface tension $\Sigma$ with two of the variables which define the actual ensemble $T,A_p$, and with the fluctuation surface tension $\sigma$. Comparison of tables \[tabS\] and \[tabsig\] indicates that, for a given system size, the slope of $\Sigma(a_p)$ is larger than the slope of $\sigma(a_p)$: since $A>A_p$ by definition of intrinsic area, the first term on the rhs of eq. (\[Sigma4\]) accounts for the difference of the slope between the two curves. Note that the second term on the rhs of eq. (\[Sigma4\]) takes into account the effect of the entropic elasticity on the system surface tension. This term is proportional to the system temperature, and inversely proportional to the system effective area per amphiphile $a_p'=A_p/N'$: the analogy with the case of polymers can be immediately drawn. It is worth noting that a similar results was obtained by Farago and Pincus in ref. [@FP], where a continuous expression for the free energy (\[free\]) was used, and where the surface tension $\sigma$ was taken to depend explicitly on the effective area $A$. A first rapid check of equation (\[Sigma4\]) can be done by noting that it predicts that, if the surface tension $\sigma$ vanishes, the frame surface tension $\Sigma$ is negative. This is confirmed by inspection of figures \[compare\_S\] and \[compare\_s\].

In equation (\[Sigma4\]) the only adjustable parameter is $N'$, while $\sigma$ can be estimated as described in section \[mis\_s\], and $A$ can be sampled during the MD simulations using eq. (\[defA\]). One can thus compare the values for $\Sigma$ predicted by eq.(\[Sigma4\]) with those directly measured as described in section \[mis\_S\]. In the following, the optimal value of the parameter $N'$ will be determined by fitting eq. (\[Sigma4\]) to the measured values plotted in fig. \[compare\_S\], for the three values of $N$ here considered. With this fit we find $N'=155$ for $N=512$, $N'=181$ for $N=768$ and $N'=213$ for $N=1152$. Using this values for $N'$, we plot the measured and the calculated values for $\Sigma$, as a function of $a_p$, for $N=512$ fig. \[simul1\], for $N=768$ fig. \[simulm\], and for $N=1152$ fig. \[simul2\]. Inspection of these figures indicates a good agreement between the values of $\Sigma$ predicted by eq.(\[Sigma4\]) and those calculated as described in section \[mis\_S\].

We now consider the scaling behaviour of the quantity $N'$ as a function of the system size. In fig. \[n1\], $N'$ is plotted as a function of the tensionless projected area $A^*_p$ (projected area at which $\sigma=0$), and as a function of the number of molecules $N$. Inspection of this figure, suggests that $N'$ is a linear function of $A^*_p$, and the slope of such function is in good agreement with that predicted by eq. (\[evaln1\]). It is worth noting that, as discussed is section \[mis\_s\], the projected area per molecule $a^*_p=A^*_p/(N/2)$ is independent of the system size, and thus we have $N'\sim A^*_p/(4 \pi)\sim N$, as one would expect in the large $N$ (or large $A_p$) limit.

Note that the intercept of the line plotted in fig. \[n1\] is non-zero, while eq. (\[evaln1\]) predicts a vanishing intercept in the large projected area limit. However the argument used to obtain eq. (\[evaln1\]) is no longer valid for small value of $A_p$ ($N$), and thus in this limit such an equation is incorrect.

Conclusion {#concl}
==========

In the present paper we have measured both the surface tension $\sigma$ appearing in the Hamiltonian which governs the bilayer shape fluctuations (\[fapp\]), and the surface tension $\Sigma$ which governs the elastic response of the system to a change in its projected area $A_p$. Using computer simulations, we have measured independently these two quantities for bilayer with different projected area. Our results indicate that the two surface tensions have different values for a given value of $A_p$. In particular the two projected area per amphiphile $a^*_p$ and $a^\dagger_p$ are different, where $a^*_p$ is the projected area per amphiphile at which the surface tension $\sigma$ vanishes, while $a^\dagger_p$ is the projected area per amphiphile where the surface tension $\Sigma$ vanishes. Furthermore, the two surface tensions are found to exhibit different slopes as functions of the projected area per molecule $a_p$, for a given system size. Using a simple thermodynamic argument, we manage to relate the two quantities $\sigma$ and $\Sigma$, and the relation between them which we found, eq. (\[Sigma4\]), nicely fits the data that we obtain from simulations. This is the most important result of the present work: it indicates that eq. (\[Sigma4\]) succeeds to capture the basic relation between the two quantities. Note that in a previous work [@Ste], where computer simulations of model membrane were considered, the two surface tensions were found to be proportional the one respect to the other, but no argument was introduced to explain such an effect.

The difference between the two quantities, does not appear to be due to some size effect, since the vanishing areas $a^*_p$ and $a^\dagger_p$ do not depend on the system size. The scaling behaviour of the effective wavemode number $N'$ as a function of the system size is a rather delicate issue. In the present work we find that $N'$ scales linearly with the number of molecule $N$ (the tensionless area $A^*_p$), but in a restricted range of values of $N$ ($A_p$). We were limited in the choice of the values of $N$ by the considerable amount of computation time required by the simulations. For larger values of $N$ one would expect that $N'$ grows either as $N$ or more slowly than $N$. If this were the case, since $A_p\sim N$, the second term of the rhs of eq. (\[Sigma4\]) would vanish in the large $N$ limit. In this case eq. (\[Sigma4\]) would read $\Sigma\simeq A/A_p \cdot \sigma$. We plan to consider larger system size in a future work. In conclusion, the results contained in this paper puts in evidence that in a bilayer membrane, the surface tension $\sigma$ appearing in the elasticity Hamiltonian (\[fapp\]) and the mechanical surface tension $\Sigma$ are two distinct thermodynamic quantities, which have to be measured independently for a full characterization of the system elastic properties, although for a constrained system, the geometrical constraints impose a relation between them.

It would be interesting to measure these two quantities in a real system, e.g, a lipid vesicle manipulated with a micropipette. On the basis of the results presented in this paper, one would expect that the geometrical constrains imposed by the micropipette would lead to a difference between the two surface tensions.

The author is indebted to R. Lipowsky for introducing him to this topic, for many interesting discussions and for a critical reading of this manuscript. I am grateful to J.B. Fournier for the long and interesting discussions on the two tension issue. I also thank L. Peliti for many discussions, and for his encouragements. This work was partially supported by MIUR-PRIN 2004. Finally, I thank the CRdC AMRA for the use of its computational resources.

[99]{} J. Israelachvili, Intermolecular and Surface Forces, (Academic Press, London, 1992). B. Alberts, A. Johnson, J. Lewis, M. Raff, K. Roberts, and P. Walter, [*Molecular Biology of the Cell*]{}, (Garland, New York, 2002). R. Lipowsky, Vesicles and biomembranes in [*Encyclopedia of Applied Physics*]{} Vol. 23, pages 199-222, (VCH Publishers, Weinheim and New York 1998). Safran S. A., Statistical Thermodynamics of Surfaces, Interfaces and Membranes (Addison-Wesley, Reading, Mass. 1994). S.J. Marrink and A.E. Mark, [*J. Phys. Chem. B*]{} [**105**]{}:6122-6127 (2001). R. Goetz and R. Lipowsky, [*J. Chem. Phys.*]{} [**108**]{}:7397-7409 (1998). E. Lindhal and O. Edholm, [*J. Chem. Phys.*]{} [**113**]{}:3882 (20 00). F. David and S. Leibler, [*J. Phys. France*]{} [**1**]{}:959-976 (1991). O. Farago and P. Pincus, [*J. Chem. Phys.*]{} [**120**]{}: 2934–2950 (2004). J. Stecki, [*J. Chem. Phys.*]{} [**120**]{}: 3508–3516 (2004). R. Goetz, G. Gompper and R. Lipowsky, [*Phys. Rev. Lett.*]{} [**82**]{}:221-224 (1999). A. Imparato, J.C. Shillcock and R. Lipowsky, [*Eur. Phys. J. E* ]{} [**11**]{}:21-28 (2003). A. Imparato, J.C. Shillcock and R. Lipowsky, [*Europhys. Lett.*]{} , [**69**]{}: 650-656 (2005). E. Egberts and H.J.C. Berendsen [*J. Chem. Phys.* ]{} [**89**]{}:3718-3732 (1988). M.P. Allen and D.J. Tildesley, [*Computer Simulations of Liquids*]{} (Oxford, Oxford, 1987). J.G. Kirkwood and F.P. Buff, [*J. Chem. Phys.*]{}, [**17**]{}:338 (1 949). P. Schofield , J.R. Henderson, Proc. R. Soc Lond. A MAT [**379**]{} ( 1776): 231-246 1982. D. Marsh , CRC handbook of lipid bilayers, Boca Raton, CRC Press, (1990). P.B. Canham, [*J. Theor. Biol*]{} [**26**]{}:61-81 (1970). W. Helfrich, [*Z. Naturforsch.*]{}, [**28c**]{}:693 (1973). J.F. Lennon and F. Brochrard, [*J. Phys. (Paris)*]{} [**36**]{}, 1035 (1975). D. Roux, [*Physica A*]{} [**172**]{}, 242–257 (1991). O.Farago and P. Pincus, [*Eur. Phys. J. E*]{}, [**11**]{}:399 (2003).

\[.9\][$\phi_i$]{} \[.9\][$i-1$]{} \[.9\][$i$]{} \[.9\][$i+1$]{} ![Cartoon of the model amphiphile molecule used in this paper.[]{data-label="mol"}](molecules.eps "fig:"){width="3cm"}

\[bl\]\[bl\]\[.8\][$N=1152$]{} \[bl\]\[bl\]\[.8\][$N=768$]{} \[bl\]\[bl\]\[.8\][$N=512$]{} \[br\]\[br\]\[1.\]\[-90\][$\Sigma$]{} \[ct\]\[ct\]\[1.\][$a_p$]{} ![Plot of $\Sigma$ as a function of the projected area per amphiphile $a_p$ for the three values of $N$ here considered, $N=512$, $N=768$ and $N=1152$.[]{data-label="compare_S"}](confronto_Sigma.eps "fig:"){width="12cm"}

\[bl\]\[bl\]\[1.\]\[-90\][$S$]{} \[br\]\[br\]\[1.\][$q$]{} \[tr\]\[tr\]\[.8\][$a_p=2.066$]{} \[tr\]\[tr\]\[.8\][$a_p=2.095$]{} \[tr\]\[tr\]\[.8\][$a_p=2.127$]{} \[tr\]\[tr\]\[.8\][$a_p=2.157$]{} \[tr\]\[tr\]\[.8\][$a_p=2.203$]{} ![Fluctuation spectrum $S$ as a function of the wavenumber $q$, for different values of the projected area per amphiphile $a_p$, for bilayers with $N=1152$. The tensionless state ($\sigma=0$) corresponds to the projected area per amphiphile $a_p^*=2.095$. The dashed line is obtained by fitting the fluctuation spectrum of the tensionless system to eq. (\[fluspe\]), with $\sigma=0$, for $q\lesssim 1$.[]{data-label="sq2"}](sq.eps "fig:"){width="12cm"}

\[bl\]\[bl\]\[.8\][$N=1152$]{} \[bl\]\[bl\]\[.8\][$N=768$]{} \[bl\]\[bl\]\[.8\][$N=512$]{} \[br\]\[br\]\[1.\]\[-90\][$\sigma$]{} \[ct\]\[ct\]\[1.\][$a_p$]{} ![Plot of $\sigma$ as a function of the projected area per amphiphile $a_p$ for the three values of $N$ here considered.[]{data-label="compare_s"}](confronto_sigma.eps "fig:"){width="12cm"}

\[ct\]\[ct\]\[1.\][$a_p$]{} \[ct\]\[ct\]\[.8\][$\sigma$]{} \[ct\]\[ct\]\[.8\][$\Sigma$]{} ![Plot of the fluctuation surface tension $\sigma$ and of the mechanical surface tension $\Sigma$ as functions of the projected area per amphiphile $a_p$, for the $N=1152$ system. The lines are linear fits to the two sets of data. The two sets of data appear to be clearly shifted and tilted with respect to each other.[]{data-label="compare_sS"}](confr_sS_1152.eps "fig:"){width="12cm"}

\[ct\]\[ct\]\[1.\][$a_p$]{} \[bl\]\[bl\]\[1.\]\[-90\][$a$]{} \[br\]\[br\]\[.8\][$N=1152$]{} \[br\]\[br\]\[.8\][$N=768$]{} \[br\]\[br\]\[.8\][$N=512$]{} ![Plot of the measured effective area per amphiphile $a$, as a function of the projected area per amphiphile $a_p$ (simulation box area), for the three values of $N$ here considered. []{data-label="effa"}](a_ap.eps "fig:"){width="12cm"}

\[br\]\[br\]\[1.\][$a_p$ ]{} \[br\]\[br\]\[.8\][$\Sigma$ ]{} \[br\]\[br\]\[.8\][$\Sigma(T,A_p,\sigma)$]{} ![Plot of the measured surface tension $\Sigma$ and of the estimate of the same quantity, as given by eq. (\[Sigma4\]), as functions of $a_p$, for a system with $N=512$ amphiphiles.[]{data-label="simul1"}](confr_2_sigma "fig:"){width="12cm"}

\[br\]\[br\]\[1.\][$a_p$ ]{} \[br\]\[br\]\[.8\][$\Sigma$ ]{} \[br\]\[br\]\[.8\][$\Sigma(T,A_p,\sigma)$]{} ![Plot of the measured surface tension $\Sigma$ and of the estimate of the same quantity, as given by eq. (\[Sigma4\]), as functions of $a_p$, for a system with $N=768$ amphiphiles.[]{data-label="simulm"}](confr_2_sigma_768.eps "fig:"){width="12cm"}

\[br\]\[br\]\[1.\][$a_p$ ]{} \[br\]\[br\]\[.8\][$\Sigma$ ]{} \[br\]\[br\]\[.8\][$\Sigma(T,A_p,\sigma)$]{} ![Plot of the measured surface tension $\Sigma$ and of the estimate of the same quantity, as given by eq. (\[Sigma4\]), as functions of $a_p$, for a system with $N=1152$ amphiphiles.[]{data-label="simul2"}](confr_2_sigma_1152 "fig:"){width="12cm"}

\[ct\]\[ct\]\[1.\][$N$]{} \[br\]\[br\]\[1.\]\[-90\][$N'$]{} \[ct\]\[ct\]\[1.\][$A^*_p$]{} ![Plot of the effective fluctuation modes $N'$ as a function of the vanishing tension projected area $A^*_p$ ($\sigma=0$) and as a function of the number of molecules $N$ (upper $x$-axis). The line is a linear fit of the data, whose slope, as a function of $A^*_p$, is $0.086\pm 0.003\simeq 1/(4 \pi)$, as predicted by eq. (\[evaln1\]). []{data-label="n1"}](n.eps "fig:"){width="12cm"}
