[gmx-users] Negative pressure in NVT simulation

Neda Rafiee nerafiee at ipm.ir
Sun Jul 15 20:33:05 CEST 2018


Dear Gromacs users,
I have a system part of which is frozen in three dimensions, so I cannot do an NPT. To add solvent to my system, I used solvate but since I am using a user-defined vdwradii.dat to inhibit water enter my solid slab (my slab has a thickness about 3 nm), after some steps I understood that density of water is less than 1000 and so I calculated available volume that water molecules can be added to my box by subtravting the volume of my slab and then I used "gmx insert-molecules" command to add as much water molecules as required to get a density of 1000 Kg/m^3. Now, I am doing NVT and I get correct density for water but I see the pressure in my system is a very large negative number (minus ten to power 24). I want to discuss with you if it is a real problem in my system or not! I am using these mdp options:

constraints         =  none
constraint_algorithm = LINCS
lincs_order         =  4
lincs_iter          =  1
lincs_warnangle     =  30
integrator          =  md
emstep              =  0.0001
emtol               =  50.0
;
; Center of mass removal ------------------------------------------
comm_mode           =  none; linear    ; Change for no pbc
comm_grps           =  system; 
nstcomm             =  10        ; freq of center of mass removal
;
; Time steps ------------------------------------------------------
dt                  =  0.001     ;1 fs
nsteps              =  1000000     ;1 ns
;
; Output control ---------------------------------------------------
nstxout             =  10000       ; freq of coord output to trr
nstvout             =  10000       ; freq of velocity output to trr
nstfout             =  0           ; freq of force output to trr
nstlog              =  10000       ; freq of output to log file
nstenergy           =  10000       ; default of -1 sets equal to nstlist
nstxtcout           =  10000       ; freq of coord output to xtc
;
; Neighbor list and pbc's ------------------------------------------
cutoff-scheme       =  group 
nstlist             =  10        ; default is 10 -- try more frequent
ns_type             =  grid     ; simple (for no pbc) or grid
;
pbc                 =  xyz      ; CHANGE
periodic_molecules  =  no
;
; Interaction parameters -------------------------------------------
rlist               =  1
rlistlong           =  1
rcoulomb            =  1
;
coulombtype         =  PME
epsilon_r           =  1
;
vdwtype             =  user; cut-off
rvdw                =  1
DispCorr            = no ;  no or EnerPres
;
energygrps          = LYS LEU VAL PHE1 PHE2 ALA GLU ASP NaB ClB OXB OSil OGem OMB SIB SiSil SiGem HSil HGem OW HW
energygrp_table = NaB OXB NaB OSil NaB OGem NaB OMB OMB OW OW SIB OW SiSil OW SiGem OW HSil OW HGem HW OMB HW OXB OW OSil OW OGem HW OSil HW OGem OXB OW

freezegrps          = SIB OXB OMB SiSil OSil SiGem OGem
freezedim           = Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y
;
;  temperature coupling is on ----------------------------------
Tcoupl              =  v-rescale  ;or v-rescale berendsen nose-hoover
tau_t               =  0.3
tc-grps             =  System
ref_t               =  300.0
; Pressure coupling is not on ----------------------------------
Pcoupl              =  no ; no berendsen Parrinello-Rahman
tau_p               =  4
pcoupltype          =  anisotropic ;   isotropic anisotropic
compressibility     =  10.e-5  10.e-5  10.e-5  0.0  0.0  0.0
ref_p               =  1.0 1.0 1.0 0.0 0.0 0.0
; distance restraint force constant
;disre               = simple
;disre_fc            = 10000

; Generate velocites -------------------------------------------
gen_vel             =  yes
gen_temp            =  300.0
gen_seed            =  173529

[ahmadabadi_i.ch.sharif at compute-0-1 1-1_max_ROG_over_Sil]$ cat  run.mdp 
constraints         =  none
constraint_algorithm = LINCS
lincs_order         =  4
lincs_iter          =  1
lincs_warnangle     =  30
integrator          =  md
emstep              =  0.0001
emtol               =  50.0
;
; Center of mass removal ------------------------------------------
comm_mode           =  none; linear    ; Change for no pbc
comm_grps           =  system; 
nstcomm             =  10        ; freq of center of mass removal
;
; Time steps ------------------------------------------------------
dt                  =  0.001     ;1 fs
nsteps              =  1000000     ;1 ns
;
; Output control ---------------------------------------------------
nstxout             =  10000       ; freq of coord output to trr
nstvout             =  10000       ; freq of velocity output to trr
nstfout             =  0           ; freq of force output to trr
nstlog              =  10000       ; freq of output to log file
nstenergy           =  10000       ; default of -1 sets equal to nstlist
nstxtcout           =  10000       ; freq of coord output to xtc
;
; Neighbor list and pbc's ------------------------------------------
cutoff-scheme       =  group 
nstlist             =  10        ; default is 10 -- try more frequent
ns_type             =  grid     ; simple (for no pbc) or grid
;
pbc                 =  xyz      ; CHANGE
periodic_molecules  =  no
;
; Interaction parameters -------------------------------------------
rlist               =  1
rlistlong           =  1
rcoulomb            =  1
;
coulombtype         =  PME
epsilon_r           =  1
;
vdwtype             =  user; cut-off
rvdw                =  1
DispCorr            = no ;  no or EnerPres
;
energygrps          = LYS LEU VAL PHE1 PHE2 ALA GLU ASP NaB ClB OXB OSil OGem OMB SIB SiSil SiGem HSil HGem OW HW
energygrp_table = NaB OXB NaB OSil NaB OGem NaB OMB OMB OW OW SIB OW SiSil OW SiGem OW HSil OW HGem HW OMB HW OXB OW OSil OW OGem HW OSil HW OGem OXB OW

freezegrps          = SIB OXB OMB SiSil OSil SiGem OGem
freezedim           = Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y
;
;  temperature coupling is on ----------------------------------
Tcoupl              =  v-rescale  ;or v-rescale berendsen nose-hoover
tau_t               =  0.3
tc-grps             =  System
ref_t               =  300.0
; Pressure coupling is not on ----------------------------------
Pcoupl              =  no ; no berendsen Parrinello-Rahman
tau_p               =  4
pcoupltype          =  anisotropic ;   isotropic anisotropic
compressibility     =  10.e-5  10.e-5  10.e-5  0.0  0.0  0.0
ref_p               =  1.0 1.0 1.0 0.0 0.0 0.0
; distance restraint force constant
;disre               = simple
;disre_fc            = 10000

; Generate velocites -------------------------------------------
gen_vel             =  yes
gen_temp            =  300.0
gen_seed            =  173529
------------------

Thanks in advance.
Neda


More information about the gromacs.org_gmx-users mailing list