[gmx-users] Segmentation fault while using external water model
Irem Altan
irem.altan at duke.edu
Fri Feb 19 02:57:49 CET 2016
Hi,
I’m not entirely sure. How can I check? Below is my .mdp file (with misleading/old comments):
title = OPLS Lysozyme NVT equilibration
define = -DPOSRES ; position restrain the protein
; Run parameters
integrator = md ; leap-frog integrator
nsteps = 50000000 ; 2 * 50000 = 10 ns
dt = 0.002 ; 2 fs
; Output control
nstxout = 5000 ; save coordinates every 10 ps
nstvout = 5000 ; save velocities every 10 ps
nstxtcout = 5000
nstenergy = 5000 ; save energies every 10 ps
nstlog = 5000 ; update log file every 10 ps
; Bond parameters
continuation = no ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
lincs_iter = 1 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
; Neighborsearching
cutoff-scheme = Verlet ;NEED TO RUN ON GPU
ns_type = grid ; search neighboring grid cells
nstlist = 5 ; 10 fs
rlist = 1.0 ; short-range neighborlist cutoff (in nm)
rcoulomb = 1.0 ; short-range electrostatic cutoff (in nm)
rvdw = 1.0 ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
; Temperature coupling is on
tcoupl = V-rescale ; modified Berendsen thermostat
tc-grps = Protein Non-Protein ; two coupling groups - more accurate
tau_t = 0.1 0.1 ; time constant, in ps
ref_t = 277 277 ; reference temperature, one for each group, in K
; Pressure coupling is off
pcoupl = no ; no pressure coupling in NVT
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes ; assign velocities from Maxwell distribution
gen_temp = 277 ; temperature for Maxwell distribution
gen_seed = -1 ; generate a random seed
Best,
Irem
> On Feb 18, 2016, at 8:21 PM, Mark Abraham <mark.j.abraham at gmail.com> wrote:
>
> Hi,
>
> Are these bonds being converted to constraints, like the original message
> suggested?
>
> Mark
>
> On Fri, 19 Feb 2016 02:14 Irem Altan <irem.altan at duke.edu> wrote:
>
>> Hi,
>>
>> I am still having trouble using TIP4P2005. When I start running my MD
>> simulation, it “runs” for a minute or so (although it does not start the
>> simulation itself), and then crashes with a segmentation fault.
>>
>> I have (modified) copy of the amber99sb forcefields in my working
>> directory, and I have made the following changes:
>>
>> - added a tip4p2005.itp file with the correct parameters
>> — named OW as OW_tip4p2005 and HW as HW_tip4p2005
>> - added the LJ parameters to ffnonbonded.itp:
>> OW_tip4p2005 8 16.00 0.0000 A 3.15890e-01 7.74908e-01
>> HW_tip4p2005 1 1.008 0.0000 A 0.00000e+00 0.00000e+00
>> - added the new atom types to atomtypes.atp:
>> HW_tip4p2005 1.00080 ; tip4p/2005 HW
>> OW_tip4p2005 16.00000 ; tip4p/2005 OW
>> - added the following line to watermodels.dat
>> tip4p2005 TIP4P/2005 tip4p/2005
>>
>> My system consists of a protein in its unit cell, unstripped of its
>> crystal waters. Below is how my topology file looks like (after solvation,
>> adding ions, and energy minimization):
>>
>> ; Include forcefield parameters
>> #include "./amber99sb.ff/forcefield.itp"
>>
>> ; Include chain topologies
>> #include "topol_Protein_chain_A.itp"
>> #include "topol_Protein_chain_B.itp"
>> #include "topol_Other_chain_A2.itp"
>> #include "topol_Other_chain_B2.itp"
>>
>> ; Include water topology
>> #include "./amber99sb.ff/tip4p2005.itp"
>>
>> #ifdef POSRES_WATER
>> ; Position restraint for each water oxygen
>> [ position_restraints ]
>> ; i funct fcx fcy fcz
>> 1 1 1000 1000 1000
>> #endif
>>
>> ; Include topology for ions
>> #include "./amber99sb.ff/ions.itp"
>> #include "./tip4p2005cw.itp"
>>
>> [ system ]
>> ; Name
>> Protein in water
>>
>> [ molecules ]
>> ; Compound #mols
>> Protein_chain_A 1
>> Protein_chain_B 1
>> Other_chain_A2 1
>> Other_chain_B2 1
>> HO4 96
>> HO4 62
>> SOL 5763
>> NA 6
>> CL 10
>>
>> where tip4p2005cw.itp is the same file, only with [settles] removed, and
>> molname changed to HO4 for labeling purposes.
>>
>> What could be the problem? (PS, see below for the contents of
>> tip4p2005.itp and tip4p2005cw.itp)
>>
>> Best,
>> Irem
>>
>>
>> =======================================================================================================
>> tip4p2005.itp (taken from here and slightly modified:
>> https://github.com/wesbarnett/trappeua/blob/master/trappeua.ff/tip4p2005.itp
>> )
>> ————————
>>
>> TIP4P/2005 Water Model
>>
>> ; J. L. F. Abascal and C. Vega, J. Chem. Phys 123 234505 (2005)
>> ; http://dx.doi.org/10.1063/1.2121687
>>
>>
>> [ moleculetype ]
>> ; molname nrexcl
>> SOL 2
>>
>> [ atoms ]
>> ; nr type resnr residu atom cgnr charge mass
>> 1 OW_tip4p2005 1 SOL OW 1 0.0000 15.994
>> 2 HW_tip4p2005 1 SOL HW1 1 0.5564 1.008
>> 3 HW_tip4p2005 1 SOL HW2 1 0.5564 1.008
>> 4 MW 1 SOL MW 1 -1.1128 0.0
>>
>> #ifndef FLEXIBLE
>> [ settles ]
>> ; i funct doh dhh
>> 1 1 0.09572 0.15139
>> #else
>> [ bonds ]
>> ; i j funct length force.c.
>> 1 2 1 0.09572 502416.0
>> 1 3 1 0.09572 502416.0
>>
>> [ angles ]
>> ; i j k funct angle force.c.
>> 2 1 3 1 104.52 628.02
>> #endif
>>
>> [ exclusions ]
>> 1 2 3 4
>> 2 1 3 4
>> 3 1 2 4
>> 4 1 2 3
>>
>> ; The position of the dummy is computed as follows:
>> ;
>> ; O
>> ;
>> ; D
>> ;
>> ; H H
>> ;
>> ; a = b = (1/2) * (distance(OD) / [ cos (angle(DOH)) * distance(OH) ] )
>> ; = (1/2) * (0.01546 nm / [ cos (52.26 deg) * 0.09572 nm ] )
>> ; = 0.13193828
>> ; Dummy pos x4 = x1 + a*(x2-x1) + b*(x3-x1)
>>
>> [ virtual_sites3 ]
>> ; Vsite from funct a b
>> 4 1 2 3 1 0.13193828 0.13193828
>>
>>
>> ==============
>> tip4p2005cw.itp
>> ————————
>>
>> ; TIP4P/2005 Water Model
>>
>> ; J. L. F. Abascal and C. Vega, J. Chem. Phys 123 234505 (2005)
>> ; http://dx.doi.org/10.1063/1.2121687
>>
>>
>> [ moleculetype ]
>> ; molname nrexcl
>> HO4 2
>>
>> [ atoms ]
>> ; nr type resnr residu atom cgnr charge mass
>> 1 OW_tip4p2005 1 SOL OW 1 0.0000 15.994
>> 2 HW_tip4p2005 1 SOL HW1 1 0.5564 1.008
>> 3 HW_tip4p2005 1 SOL HW2 1 0.5564 1.008
>> 4 MW 1 SOL MW 1 -1.1128 0.0
>>
>> ;#ifndef FLEXIBLE
>> ;[ settles ]
>> ; i funct doh dhh
>> ; 1 1 0.09572 0.15139
>> ;#else
>> [ bonds ]
>> ; i j funct length force.c.
>> 1 2 1 0.09572 502416.0
>> 1 3 1 0.09572 502416.0
>>
>> [ angles ]
>> ; i j k funct angle force.c.
>> 2 1 3 1 104.52 628.02
>> ;#endif
>>
>> [ exclusions ]
>> 1 2 3 4
>> 2 1 3 4
>> 3 1 2 4
>> 4 1 2 3
>>
>> ; The position of the dummy is computed as follows:
>> ;
>> ; O
>> ;
>> ; D
>> ;
>> ; H H
>> ;
>> ; a = b = (1/2) * (distance(OD) / [ cos (angle(DOH)) * distance(OH) ] )
>> ; = (1/2) * (0.01546 nm / [ cos (52.26 deg) * 0.09572 nm ] )
>> ; = 0.13193828
>> ; Dummy pos x4 = x1 + a*(x2-x1) + b*(x3-x1)
>>
>> [ virtual_sites3 ]
>> ; Vsite from funct a b
>> 4 1 2 3 1 0.13193828 0.13193828
>>
>> --
>> Gromacs Users mailing list
>>
>> * Please search the archive at
>> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
>> posting!
>>
>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>
>> * For (un)subscribe requests visit
>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>> send a mail to gmx-users-request at gromacs.org.
> --
> Gromacs Users mailing list
>
> * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-request at gromacs.org.
More information about the gromacs.org_gmx-users
mailing list