Re: equilibrate the bulk water system

From: Hasan haska (hasanhaska_at_yahoo.com)
Date: Sat Mar 14 2015 - 14:38:03 CDT

Dear Maxim Belkin,
Thank you for your valuable explanation. I used a spring constant of 450 for OH bond and a spring constant of 55.00 for HOH angle in my parameter file as you advised. And the problem was solved. But should I also use “rigidBonds water”  or  “rigidBonds all” definition in conf file ? 
My definitions for hydrogen and oxygen molecules are 0S1, HS2, HS3 and types are OS and HS ? I think that these are not the same in general charmm prm file. Can ıt be known by namd that I use water molecules in my simulation system ?
And I also found another namd prm file for SPC/E water model which is in below. And I realized that this prm file is different form my prm file. Could you please help me about choosing the suitable SPC/E prm values for my bulk water simulation ?
*****************My prm File:   SPC/E
BONDS
!V(bond) = Kb(b - b0)**2!Kb: kcal/mole/A**2!b0: A!atom type Kb          b0
OS HS         450.0 1.000

ANGLES
!V(angle) = Ktheta(Theta - Theta0)**2!V(Urey-Bradley) = Kub(S - S0)**2!Ktheta: kcal/mole/rad**2!Theta0: degrees!Kub: kcal/mole/A**2 (Urey-Bradley)!S0: A!atom types     Ktheta    Theta0   Kub     S0

HS OS HS      55.0    109.47

 NONBONDED 
!V(Lennard-Jones) = Eps,i,j[(Rmin,i,j/ri,j)**12 - 2(Rmin,i,j/ri,j)**6]!epsilon: kcal/mole, Eps,i,j = sqrt(eps,i * eps,j)!Rmin/2: A, Rmin,i,j = Rmin/2,i + Rmin/2,j!atom  ignored    epsilon      Rmin/2   ignored   eps,1-4       Rmin/2,1-4

OS 0.     -0.1553943     1.7766 0.    -0.1553943    1.7766HS 0. 0.0000       0.0 0.     0.0000       0.0
END*************

**************Another PRM File : SPC/E
BONDS
!V(bond) = Kb(b - b0)**2!Kb: kcal/mole/A**2!b0: A!atom type Kb          b0
OS HS   600.000     1.000

ANGLES
!V(angle) = Ktheta(Theta - Theta0)**2!V(Urey-Bradley) = Kub(S - S0)**2!Ktheta: kcal/mole/rad**2!Theta0: degrees!Kub: kcal/mole/A**2 (Urey-Bradley)!S0: A!atom types     Ktheta    Theta0   Kub     S0

HS OS HS   124.00    109.50

 NONBONDED 
!V(Lennard-Jones) = Eps,i,j[(Rmin,i,j/ri,j)**12 - 2(Rmin,i,j/ri,j)**6]!epsilon: kcal/mole, Eps,i,j = sqrt(eps,i * eps,j)!Rmin/2: A, Rmin,i,j = Rmin/2,i + Rmin/2,j!atom  ignored    epsilon      Rmin/2   ignored   eps,1-4       Rmin/2,1-4

HS     0.000000  -0.000000     0.000000OS     0.000000  -0.155298     1.776857

END
***********
 

     On Friday, March 13, 2015 7:20 PM, Maxim Belkin <mbelkin_at_ks.uiuc.edu> wrote:
   

 Hasan,
The problem is with your parameter file:
1. Bond:  OS HS 0.0 1.0Which means: spring constant 0.0, distance 1.0
2. Angle:  HS OS HS 0.0 109.47Which means: spring constant 0.0, angle 109.47
3. You don't use rigid bonds (i.e., in config file you don't specify: rigidBonds water)
During minimization bonds in water molecules collapse, angles become 180, and you get this error message in the end. So, you get exactly what you ask for given the parameter file you feed in.
Try using a spring constant of 450 for OH bond, a spring constant of 55.00 for HOH angle.
OS HS 450.0 1.0HS OS HS 55.0 109.47
Another, minor note: edit your parameter file in a different text editor...
Maxim

On Mar 13, 2015, at 6:57 AM, Hasan haska <hasanhaska_at_REMOVE_yahoo.com> wrote:
Dear Norman Geist,

I increased the min. step form 100000 to 150000  as you adviced and run the conf file in below again. But I got the same error "Reason: FATAL ERROR: Stray PME grid charges detected!" . What should I do to solve this error ? Could you please suggest a right procedure to equilibriate the bulk water simulation system ?

Thaks in advance. 

#--- This is a test namd configuration file
#############################################################
## 
#############################################################
 
paratypeCharmm  on
parameters    water.PRM
 
structure              water.psf  #
coordinates            water_.pdb  #
outputName             water_1out  #
 
set temperature    300
  
temperature $temperature
                  
firsttimestep      0

 
# Integrator Parameters
timestep            1.0
nonbondedFreq       2
fullElectFrequency  4    
stepspercycle       16  
 
 
# Force-Field Parameters
exclude         scaled1-4 
1-4scaling      1.0
cutoff          18.
switching       on
switchdist      15.            
pairlistdist    20.
 
 
# Constant Temperature Control
langevin            on   ;# do langevin dynamics
langevinDamping     1    ;# damping coefficient (gamma) of 5/ps
langevinTemp        $temperature
langevinHydrogen    on    ;#  couple langevin bath to hydrogens
 
 
# Constant Pressure Control (variable volume)
 
useGroupPressure      no ;# needed for 2fs steps
useFlexibleCell       yes  ;# no for water box, yes for membrane
useConstantArea       yes  ;# no for water box, yes for membrane
 
langevinPiston        on
langevinPistonTarget  1.01325 ;#  in bar -> 1 atm
langevinPistonPeriod  200.
langevinPistonDecay   50.
langevinPistonTemp    $temperature
 
 
#--- PBC
 
cellBasisVector1  55.0   0.0    0.0
cellBasisVector2   0.0  55.0    0.0
cellBasisVector3   0.0   0.0  100.0
 
 
#--- PME 
PME                 yes
PMEGridSpacing      1.0
 
 
#--- Output & Restart
 
binaryoutput    no
 
binaryrestart   yes
 
restartname     water_1restart      #
DCDfile         water_1out.dcd      #
 
restartfreq        5000      #
dcdfreq            5000
xstFreq            5000
outputEnergies     5000
outputPressure     5000
outputtiming       5000
 
minimize 150000reinitvels          $temperature
run 10000000
# 10 ns 2002 record

On Wednesday, March 11, 2015 5:19 PM, Norman Geist <norman.geist_at_uni-greifswald.de> wrote:

>From the Error message one can conclude instability of your system. The error message points out, that the grid charges for PME have changed more than expected between timesteps, which results from really fast atoms. The most likely reason are bad contacts in the initial structure and too short minimization. So increase the number of minimization steps to approx. 10000 steps. Also check the TOTAL energy of the system for convergence during the minimization run (VMDs Namd Plot plugin is useful for that).  For the latter NPT run, you are interested in: useflexiblecell yesuseconstantarea yes along with the piston pressure to have only Z changing. Norman Geist. From: owner-namd-l_at_ks.uiuc.edu [mailto:owner-namd-l_at_ks.uiuc.edu] On Behalf Of Hasan haska
Sent: Tuesday, March 10, 2015 10:52 PM
To: NAMD List
Subject: namd-l: equilibrate the bulk water system Dear NAMD users, I want to minimize and equilibrate a bulk water system that includes 7048 water molecules (SPC/E water model) in a rectangular box with the d=0.996 g/cm3 experimental density at 300K. Therefore I prepared a namd conf file below and started to run it. But after minimization was completed I got the error, which is in below. My aim is to run minimization and NPT simulation with keeping the dimension of the unit cell in the x-y plane constant while allowing fluctuations along the z axis in order to get the real cell sizes ( 55 X 55 X 70 Angs ) according to the density. Then run a NVT simulation with these cell size. How can I solve this problem ?  Could you please share your opinion and experience ?   Thanks.  http://www.filedropper.com/water_1 My simulation files are in the link above.  *****************************Error:------------- Processor 0 Exiting: Called CmiAbort ------------Reason: FATAL ERROR: Stray PME grid charges detected! FATAL ERROR: See http://www.ks.uiuc.edu/Research/namd/bugreport.html Charm++ fatal error:FATAL ERROR: Stray PME grid charges detected! **************************   #--- This is a test namd configuration file ############################################################### ############################################################# paratypeCharmm  on parameters    water.PRM  structure              water.psf  # coordinates            water_.pdb  # outputName             water_1out  #  set temperature    300   temperature $temperature                   firsttimestep      0   # Integrator Parameters timestep            1.0 nonbondedFreq       2 fullElectFrequency  4     stepspercycle       16     # Force-Field Parameters exclude         scaled1-4  1-4scaling      1.0 cutoff          18. switching       on switchdist      15.             pairlistdist    20.   # Constant Temperature Control langevin            on   ;# do langevin dynamics langevinDamping     1    ;# damping coefficient (gamma) of 5/ps langevinTemp        $temperature langevinHydrogen    on    ;#  couple langevin bath to hydrogens   # Constant Pressure Control (variable volume)  useGroupPressure      no ;# needed for 2fs steps useFlexibleCell       yes  ;# no for water box, yes for membrane useConstantArea       yes  ;# no for water box, yes for membrane  langevinPiston        on langevinPistonTarget  1.01325 ;#  in bar -> 1 atm langevinPistonPeriod  200. langevinPistonDecay   50. langevinPistonTemp    $temperature   #--- PBC  cellBasisVector1  55.0   0.0    0.0 cellBasisVector2   0.0  55.0    0.0 cellBasisVector3   0.0   0.0  100.0   #--- PME  PME                 yes PMEGridSpacing      1.0   #--- Output & Restart  binaryoutput    no  binaryrestart   yes  restartname     water_1restart      # DCDfile         water_1out.dcd      #  restartfreq        5000      # dcdfreq            5000 xstFreq            5000 outputEnergies     5000 outputPressure     5000 outputtiming       5000  minimize 10000 reinitvels          $temperature run 10000000 # 10 ns 2002 record

  

This archive was generated by hypermail 2.1.6 : Tue Dec 27 2016 - 23:20:59 CST