On 7/10/17 10:52 AM, Sonia Milena Aguilera Segura wrote:
Dear Justin,
Thank you for the answer. I changed the two parameters suggested in the mdp
file and I ran again a minimization, 200 ps NVT, 200 ps NPT, and 1 ns MD for
the two cases of the previous mail, and now I am getting densities around 771
g/m3 which is slightly underestimated, but close to what other authors have
obtained (774 others and 777 experimental). Also, I got slightly higher values
for dielectric constants and diffusivities, also closer to another simulation
with CHARMM. The energies also changed, but I guess that was expected. It looks
like reproducing the dielectric constant with the current parameters is not
possible. Is there anything that can be changed in order to get a better
description? In terms of simulation, what is the dielectric constant depending
of?
Structure, charge distribution, etc. This may just be a limitation of an
additive force field model. We typically see neat liquid properties better
reproduced with polarizable models.
Moreover, I observed that this time I got lower values for P during the NPT
equilibration, but still is too far from 1 bar. I really don't understand why
for the NVT simulation I get a T around 298, but as soon as I turn on the
pcoupl, the T rises to 300-301 K and the P gets average values of 7 and 4 bar
(vs 8 and 14 for the previous simulations). Then at the end of the 1-ns MD the
temperature remains around 301 and the P is -1 and 2.7 bar. Considering the
parameters I am using, is there anything I can change to make the P coupling
better? I am running a 3 nm box with 308 molecules. This is the full mdp file:
http://www.gromacs.org/Documentation/Terminology/Pressure
Your box is very small and will be subject to large fluctuations.
-Justin
; Run control
integrator = sd ; Langevin dynamics
tinit= 0
dt = 0.0005
nsteps = 200 ; 1 ns
nstcomm = 100
;energygrps = non-Water
; Neighborsearching and short-range nonbonded interactions
cutoff-scheme= verlet
nstlist = 20
ns_type = grid
pbc = xyz
rlist= 1.2
; Electrostatics
coulombtype = PME
rcoulomb = 1.2
; van der Waals
vdwtype = cutoff
vdw-modifier = force-switch
rvdw-switch = 1.0
rvdw = 1.2
; Apply long range dispersion corrections for Energy and Pressure
DispCorr = no
; Spacing for the PME/PPPM FFT grid
fourierspacing = 0.12
; EWALD/PME/PPPM parameters
pme_order= 6
ewald_rtol = 1e-06
epsilon_surface = 0
; Temperature coupling
; tcoupl is implicitly handled by the sd integrator
tc_grps = system
tau_t= 1.0
ref_t= 298.15
; Pressure coupling is on for NPT
Pcoupl = Parrinello-Rahman
tau_p= 1.0
compressibility = 4.5e-05
ref_p= 1.0
; Do not generate velocities
gen_vel = no
; options for bonds
constraints = none ; we only have C-H bonds here
; Type of constraint algorithm
constraint-algorithm = lincs
; Constrain the starting configuration
; since we are continuing from NPT
continuation = yes
; Highest order in the expansion of the constraint coupling matrix
lincs-order = 12
Thank you very much,
Sonia Aguilera
PhD student
ENSCM
; Run control
integrator = sd ; Langevin dynamics
tinit= 0
dt = 0.0005
nsteps = 4000 ; 20 ns
nstcomm = 100
; Neighborsearching and short-range nonbonded interactions
cutoff-scheme= verlet
nstlist = 20
ns_type = grid
pbc = xyz
rlist= 1.2
; Electrostatics
coulombtype = PME
rcoulomb = 1.2
; van der Waals
vdwtype = cutoff
vdw-modifier = potential-switch
rvdw-switch = 1.0
rvdw = 1.2
; Apply long range dispersion corrections for Energy and Pressure
DispCorr = EnerPres
CHARMM uses a force switch, and only apply dispersion correction in the case of
lipid monolayers.
http://www.gromacs.org/Documentation/Terminology/Force_Fields/CHARMM
-Justin
==
--
==
Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow
Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201
jalem...@outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul
