[gmx-users] Acetonitrile using CHARMM ff

[Sonia Milena Aguilera Segura]

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?

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:

; Run control
integrator               = sd       ; Langevin dynamics
tinit                    = 0
dt                       = 0.0005
nsteps                   = 2000000   ; 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                   = 40000000   ; 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

==================================================