[gmx-users] Re: surface tension
André Farias de Moura
moura at ufscar.br
Thu Mar 17 16:19:19 CET 2011
Dear Elisabeth,
I'm going to make some general remarks on what I would consider a
reasonable set of choices for such simulations:
(1) 3x3 nm in the xy plane seems perfectly fine to me if your model systems
is made up of some well-behaved molecular liquid, alkanes for instance.
(2) GROMACS manual describes the Ewald geometry option "3dc" for such
slab model systems if long range corrections are to be applied
(3) whether or not the NVT ensemble is adequate for such simulations
is arguable,
but definitely applying a plain isotropic NpT coupling introduces
errors. take a look
at the size of your system using g_energy, probably it shrank in all
three directions
due to the vacuum in the z-direction. again: either try NVT or use a
surface tension
coupling (setting compressibility to zero in the z direction to keep
the box length
constant in that direction)
(4) surface tension should be a positive number, meaning that creating a surface
requires some amount of energy.
best regards,
André
On Wed, Mar 16, 2011 at 8:31 PM, Elisabeth <katesedate at gmail.com> wrote:
> Dear Dr. Andre F. de Moura,
>
> Thanks for your invaluable help. I folloed your instructions and made a 3 nm
> box of alkane (125 molecules) which gives the actual density of material,
> and increased the size in Z direction. i. e 3X3X6 then did a 8 ns NPT with
> parrinello-rahman algorithm. I dont know why I am getting negative surface
> tension. However, the gamma I get from formula Lz=6nm n=2, gamma: 6/2(Pzz -
> 1/2(Pxx + Pyy)) is almost identical to "#Surf*SurfTen divided by 2" values.
>
> 0- My question is, the system I have now 3x3x6 is in no way similar to
> bilayers ( where the structure of system shows a planar surface parrale to
> XY. Is this because I am studying a pure system? making it different from
> bilayer study?
>
> 1- negative results are because I did not perform NVT? I am interested in
> NPT to see the effect of pressure on surface tension. does NPT work with
> surface tenio studies in gromcas?
>
> 2- In your previous message you pointed that if the empty region is large
> enough "periodicity no longer affects the system during a typical MD
> simulation, except
> for the other two directions". My question is given the 3 nm box size, the
> empty region in X and Y are also larger than cut-off ! Do I need to work
> with a 1 nm box size ( because of using 1.1 nm cut offs) instead of 3nm to
> make sure only in Z direction BC is not working?
>
> 3- I am not sure if ~ -3 or ~ -1.7 is the actual surface tension value I
> should consider. Please see the averages for different time periods. In
> other words I dont know over how many steps I need to read the average.
>
> 4- You suggested to turn off long range interactions, Does that mean PME is
> not needed? Please see the mdp file below.
>
>
> Thank you so much for you time.
> Best regards,
>
> ====================================================================================
>
> 8 ns (4 000 000 steps, 0.002 dt) NPT on a box of 3 nm replicated by nbox 112
> in Z:
>
>
>
> Statistics over 3500001 steps [ 1000.0001 thru 8000.0005 ps ], 3 data sets
>
> All averages are exact over 3500001 steps
>
>
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
>
> -------------------------------------------------------------------------------
>
> Pressure (bar) 20.2098 375.159 375.158 -0.00039324
> -2.75268
>
> Density (SI) 616.449 9.13759 9.13759 -3.21055e-06
> -0.0224739
>
> #Surf*SurfTen -1.69698 5074.59 5074.08 -0.0354603
> -248.222
>
>
>
> Statistics over 3000001 steps [ 2000.0001 thru 8000.0005 ps ], 2 data sets
>
> All averages are exact over 3000001 steps
>
>
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
>
> -------------------------------------------------------------------------------
>
> Pressure (bar) 20.2281 375.189 375.18 -0.00145544
> -8.73263
>
> #Surf*SurfTen -1.18619 5062.76 5062.55 -0.0265989
> -159.593
>
>
>
> Statistics over 2500001 steps [ 3000.0002 thru 8000.0005 ps ], 2 data sets
>
> All averages are exact over 2500001 steps
>
>
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
>
> -------------------------------------------------------------------------------
>
> Pressure (bar) 20.1216 374.81 374.795 -0.00236499
> -11.825
>
> #Surf*SurfTen -3.98186 5055.47 5055.43 -0.013803
> -69.0151
>
>
>
> From formula: -1.89 bar nm
>
>
>
>
>
> Statistics over 2000001 steps [ 4000.0002 thru 8000.0005 ps ], 2 data sets
>
> All averages are exact over 2000001 steps
>
>
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
>
> -------------------------------------------------------------------------------
>
> Pressure (bar) 20.0417 376.748 376.725 -0.00359389
> -14.3756
>
> #Surf*SurfTen -6.07893 5074.92 5074.89 -0.0160507
> -64.2028
>
>
>
> Statistics over 1500001 steps [ 5000.0000 thru 8000.0005 ps ], 2 data sets
>
> All averages are exact over 1500001 steps
>
>
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
>
> -------------------------------------------------------------------------------
>
> Pressure (bar) 20.0348 377.054 377.046 -0.00277188
> -8.31566
>
> #Surf*SurfTen -6.4344 5085.67 5085.63 -0.0216084
> -64.8252
>
>
>
> Gamma from Formula: -3.12 bar.nm
>
>
>
> Statistics over 1000001 steps [ 6000.0005 thru 8000.0005 ps ], 2 data sets
>
> All averages are exact over 1000001 steps
>
>
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
>
> -------------------------------------------------------------------------------
>
> Pressure (bar) 20.0461 377.419 377.398 -0.00690407
> -13.8082
>
> #Surf*SurfTen -6.27837 5096.53 5096.51 0.0228292
> 45.6584
>
>
>
> From formula: -3.075 bar.nm
>
>
>
>
>
> Statistics over 500001 steps [ 7000.0005 thru 8000.0005 ps ], 2 data sets
>
> All averages are exact over 500001 steps
>
>
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
>
> -------------------------------------------------------------------------------
>
> Pressure (bar) 20.0929 377.928 377.925 -0.00557213
> -5.57214
>
> #Surf*SurfTen -3.82994 5097.87 5097.75 0.119265
> 119.266
>
>
>
> ====================================================
>
> mdp file:
>
> integrator = md
>
> dt = 0.002
>
> nsteps = 4000000
>
> nstenergy = 100
>
> nstxout = 100
>
> nstvout = 0
>
>
>
>
>
> nstlist = 10
>
> ns_type = grid
>
>
>
> coulombtype = PME
>
> vdw-type = Shift
>
> rcoulomb-switch = 0
>
> rvdw-switch = 0.9
>
>
>
> rlist = 1.1
>
> rcoulomb = 1.1
>
> rvdw = 1.0
>
> fourierspacing = 0.12
>
> fourier_nx = 0
>
> fourier_ny = 0
>
> fourier_nz = 0
>
> pme_order = 4
>
> ewald_rtol = 1e-5
>
> optimize_fft = yes
>
>
>
> ; Temperature coupling
>
> Tcoupl = v-rescale
>
> tc-grps = System
>
> tau_t = 0.1
>
> ref_t = 300
>
>
>
> ; Pressure coupling
>
> Pcoupl = parrinello-rahman
>
> Pcoupltype = isotropic
>
> tau_p = 0.5
>
> compressibility = 4.5e-5
>
> ref_p = 20
>
>
>
> gen_vel = yes
>
> gen_temp = 300.0
>
> gen_seed = 173529
>
>
>
> constraints = all-bonds
>
> constraint-algorithm = lincs
>
> ==================================
>
> 2011/3/15 André Farias de Moura <moura at ufscar.br>
>>
>> Dear Elisabeth,
>>
>> PBC are still there when you increase the box length in one direction, but
>> that
>> increase creates an empty region between the periodic images. provided the
>> empty region is large enough (larger than cutoff values, for
>> instance), periodicity
>> no longer affects the system during a typical MD simulation, except
>> for the other
>> two directions.
>>
>> as regards long range corrections, either do not use them at all or choose
>> a two
>> dimensional correction (see the manual).
>>
>> and try some google search on the simulation of Langmuir films. for
>> instance, I
>> wrote this paper a while ago:
>>
>> Molecular Dynamics Simulation of a Perylene-Derivative Langmuir Film
>> André F. de Moura and Milan Trsic
>> J. Phys. Chem. B, 2005, 109 (9), pp 4032–4041
>> DOI: 10.1021/jp0452711
>>
>> best
>>
>> Andre
>>
>>
>> On Tue, Mar 15, 2011 at 1:52 AM, Elisabeth <katesedate at gmail.com> wrote:
>> > Dear Andre,
>> >
>> > Thanks for the helpful information.
>> >
>> > I need to do some text reading to understand the periodic BC effect you
>> > are
>> > talking about. I dont see why increasing length in z direction does not
>> > lead
>> > to periodic BC in z and only for x, y ? does that mean the thickness of
>> > layer would be the Z dimension? then how much increase in one direction
>> > is
>> > reasonable? (If I have a 2 nm box).
>> >
>> > Also Can you please introduce some text book?
>> >
>> > Thank you,
>> > Best regards,
>> >
>> > 2011/3/15 André Farias de Moura <moura at ufscar.br>
>> >>
>> >> Dear Elisabeth,
>> >>
>> >> actually, it is the other way around, you need increase the box length
>> >> in
>> >> one direction, thus keeping periodic boundary conditions in the other
>> >> two
>> >> directions while a (infinitely periodic) surface is created. and notice
>> >> that
>> >> using genconf with -nbox 3 3 1 will increase your system but will not
>> >> make
>> >> it a surface system, unless you increase the box length in one
>> >> direction.
>> >>
>> >> as regards the size, the larger the model system, the smaller should be
>> >> the fluctuations. but mind that you should increase the size by a few
>> >> order of magnitude for any noticeable decrease on the (huge) RMSD
>> >> values you are getting. if you want to check the results for
>> >> convergence
>> >> maybe you could try either a block averaging or a running average
>> >> (grace
>> >> can do that for you).
>> >>
>> >> best regards
>> >>
>> >> Andre
>> >>
>> >> On Mon, Mar 14, 2011 at 7:57 PM, Elisabeth <katesedate at gmail.com>
>> >> wrote:
>> >> > Hello,
>> >> >
>> >> > Thank you for your answer.
>> >> >
>> >> > 1- If I am right I have to increase the length in two directions
>> >> > rather
>> >> > than
>> >> > one, to create a plane parallel to XY for example?
>> >> >
>> >> > 2- Can you please give me an idea on how many molecules I need to
>> >> > have
>> >> > in
>> >> > the box and also what should be the thickness of layer? I have now
>> >> > 3nm X
>> >> > 9 X
>> >> > 9 dimensions. That is thickness of 3nm. What I did was replicating a
>> >> > 3nm
>> >> > box
>> >> > using genconf -nbox 3 3 1. I dont know what is the correct way of
>> >> > creating a
>> >> > layer for surface tension calculation.
>> >> >
>> >> > I appreciate any comments about number of molecules, box dimensions
>> >> > for
>> >> > such
>> >> > a study.
>> >> >
>> >> > 3- my last question is how can I make sure surface tension reported
>> >> > by
>> >> > g_energy is the equilibrated one. RMSD is very big compared to surf.
>> >> > ten. !
>> >> >
>> >> > Thanks for your time.
>> >> > Elisabeth
>> >> >
>> >> > ******************************************
>> >> >
>> >> > if you are interested in the surface tension of a pure liquid, which
>> >> > I
>> >> > assume is
>> >> > true from your message, then you need to create at least one surface,
>> >> > since
>> >> > periodic boundary conditions make the model system infinite, i.e.,
>> >> > without a
>> >> > surface whatsoever.
>> >> >
>> >> > the easiest way to make that happen is to increase the length of the
>> >> > box
>> >> > in
>> >> > one direction, say the z direction. that way you will end up with a
>> >> > system
>> >> > that
>> >> > resemble a (thin) liquid film with vacuum below and above, meaning
>> >> > that
>> >> > you
>> >> > now have two surfaces. run a regular simulation (NVT) e use g_energy
>> >> > to
>> >> > get
>> >> > the surface tension.
>> >> >
>> >> > btw: as any other pressure related property, fluctuations are huge.
>> >> >
>> >> > best
>> >> >
>> >> > Andre
>> >> >
>> >> > On Wed, Mar 9, 2011 at 12:25 PM, Elisabeth <katesedate at gmail.com>
>> >> > wrote:
>> >> >> Dear gmx users,
>> >> >>
>> >> >> Since I am new to surface tension topic I need to ask very trivial
>> >> >> questions. Please help me out with these simple questions.
>> >> >>
>> >> >> As a starting point I am going to calculate surface tension of a
>> >> >> pure
>> >> >> alkane
>> >> >> in a cubic box and compare with experimental values.
>> >> >>
>> >> >> 1- g_energy is giving #Surf*SurfTen by default. On the other hand
>> >> >> surface
>> >> >> tension can be obtained by gamma = (Pzz - (Pxx+Pyy)/2) / Lz. i.e
>> >> >> Pres-XX-(bar), Pres-YY(bar), Pres--(bar)
>> >> >>
>> >> >> Can anyone tell me what the difference between these two is?
>> >> >>
>> >> >> 2- In pressure coupling settings there is surface_tension option
>> >> >> which
>> >> >> I
>> >> >> guess is applicable where surface tension needs to be kept fixed. If
>> >> >> one
>> >> >> want to calculate surface tension I dont think this option make
>> >> >> sense.
>> >> >> Am
>> >> >> I
>> >> >> right?
>> >> >>
>> >> >> 3- I am using the following setting: I calculate the average for a
>> >> >> 2ns
>> >> >> run
>> >> >> and different start times as shown below. Although T, P and other
>> >> >> quantities
>> >> >> are equilibrated after 200ps, surface tension is not giving a
>> >> >> constant
>> >> >> value. Is that because I am not using berenden P coupling? (As
>> >> >> mentioned
>> >> >> in
>> >> >> the manual surface tension works with berendsen)
>> >> >>
>> >> >> Pcoupl = Parrinello-Rahman
>> >> >> Pcoupltype = isotropic
>> >> >> tau_p = 1 1
>> >> >> compressibility = 4.5e-5
>> >> >> ref_p = 40
>> >> >>
>> >> >>
>> >> >> time period for which average is calculated Average RMSD
>> >> >> Fluct. Drift Tot-Drift
>> >> >>
>> >> >>
>> >> >>
>> >> >> -------------------------------------------------------------------------------
>> >> >> 1-2000 ps run: #Surf*SurfTen 6.43844 3588.74
>> >> >> 3588.35 0.091406 182.721
>> >> >> 500-2000 ps #Surf*SurfTen 12.8518 3605.72
>> >> >> 3605.26 0.132126 198.189
>> >> >> 1000-2000ps #Surf*SurfTen 18.8821 3610.97
>> >> >> 3610.8 0.11819 118.191
>> >> >> 1500-2000ps #Surf*SurfTen 23.0072 3585.51
>> >> >> 3584.93 -0.444037 -222.019
>> >> >>
>> >> >>
>> >> >>
>> >> >> 4- Assuming I am getting surface tension for a cubic box, to compare
>> >> >> this
>> >> >> with reported values in literature I need to divide by 6 (no. of
>> >> >> surfaces)?
>> >> >>
>> >> >> 5- Does box six affect the results? (mine is 3.3 nm ).
>> >> >>
>> >> >>
>> >> >> Thank you,
>> >> >>
>> >
>> >
>
>
More information about the gromacs.org_gmx-users
mailing list