[gmx-users] A question about simulating a box of ethane in GROMACS
Justin A. Lemkul
jalemkul at vt.edu
Tue Sep 6 16:45:39 CEST 2011
Phil (Yang) Song wrote:
> Hi, all
>
> I am learning GROMACS 4.5.1 recently and encountered a problem that
> puzzled me for weeks. I am wondering if someone can help me with this.
>
> Here is the problem: I was trying to perform a MD simulation for 216
> ethane molecules in a box. I have first generated a box of ethane with
> random position and orientation. Then, I wanted to run a minimization
> for the box of ethane to get rid of close contact of the randomly ethane
> molecules. However, I have found that some of the ethane molecules laid
> on top of each other after minimization and infinity forces were resulted.
>
Are they overlapping only after minimization, or does the starting configuration
contain clashes? Atoms should not be significantly displaced during EM, so I'd
be inclined to believe the starting structure is the problem.
> Since my mdp file was successfully used to generate other type of
> molecules such as acetylene and ethylene and they should be correct.
> Also, I have tried to use conjugate gradient instead of steepest decent
> for the minimization, but this effort came out to be in vain.
>
> Later, I thought the topology file should be the source of error. I have
> carefully checked the topology file for couple of times but did not find
> any obvious error. I have also tried to change the charges on each atom
> in the ethane to 0 and then run the minimization. Again, the ethane was
> attracted into each other. Hence, I think the vdw could be the problem.
> However, the parameters of the vdw interaction are all from the opls-aa
> files that come with gromacs and it should be fine.
>
Simplify the system. Rather than dealing with a full box of ethane, try
minimizing one molecule, then a few (like 2 or 4) in a box, then build up. If
these small systems are stable, then the topology and .mdp are fine (though some
of your cutoffs look strange to me). Then build your target system once you've
verified that all the pieces should be working together.
-Justin
> Eventually, I am puzzled...
>
> I am iterating the content of the content of mdp input file and itp file
> in the email and also attaching the input and output of my result as a
> tarball. Hope someone can help me to solve this problem.
>
> Thank you in advance!
>
> Best,
> Phil Song
>
> ================================================================================
> MDP file
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; preprocessing options
>
> ; title of the simulation
> title =
>
> ; include path
> include =
>
> ; cpp flag
> define =
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; energy minimization
> integrator = steep
>
> ; Start time and timestep in ps
> tinit = 0
> dt = 0.001
> nsteps = 500
>
> ; run continuation or reperforming
> init_step = 0
>
> emtol = 1000.0
> emstep = 0.01
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; periodic boundary conditions
> pbc = xyz
>
> periodic_molecules = no
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; electrostatics and vdw options
> ; electrostatics method
> coulombtype = PME
> rcoulomb = 1.3
>
> ; fourier setup for PME
> fourierspacing = 0.12
> fourier_nx = 0
> fourier_ny = 0
> fourier_nz = 0
> pme_order = 4
> ewald_rtol = 1e-5
> optimize_fft = yes
>
> ewald_geometry = 3d
>
> ; vdw cutoff
> rvdw-switch = 0
> rvdw = 1.31
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; neighbor list options
> ; neighbor list frequency
> nstlist = 10
>
> ; neighbor search algorithm
> ns_type = grid
>
> ; neighbor list cut-off
> rlist = 1.30
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; thermostat and barostat
> tcoupl = v-rescale
> tc_grps = System
> ref_t = 50.5
> tau_t = 0.5
>
> pcoupl = berendsen
> pcoupltype = isotropic
> nstpcouple = 10
> ref_p = 1.01325
> tau_p = 1.0
>
> ; Using compressibility of ethanol at 273.1 (approximately)
> ; from J. Phys. Chem., 1963, 67 (10), pp 2160.2164
>
> compressibility = 1.0e-5
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; velocty generation
> gen_vel = yes
> gen_temp = 50.5
> ; gen_seed = 173529
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; output control options
> ; output frequency for coords, velocities and forces
> nstxout = 5000
> nstvout = 5000
> nstfout = 5000
>
> ; output frequency for energy
> nstenergy = 5000
> ; energygrps =
>
> ; output frequency for log
> nstlog = 5000
>
> ; output options for trajectory
> nstxtcout = 5000
> xtc_precision = 100000
> ; xtc_grps =
>
>
> ================================================================================
> ITP file
>
> ; Ethane:
> ; Assign of the atom index
> ;
> ; H2 H6
> ; | |
> ; H3 - C1 - C5 - H7
> ; | |
> ; H4 H8
> ;
>
> [ moleculetype ]
> ; name nrexcl
> Ethane 3
>
> [ atoms ]
> ; nr type resnr residu atom cgnr charge mass
> 1 opls_135 1 ETH CE1 1 -0.18
> 2 opls_140 1 ETH HE2 1 0.06
> 3 opls_140 1 ETH HE3 1 0.06
> 4 opls_140 1 ETH HE4 1 0.06
> 5 opls_135 1 ETH CE5 2 -0.18
> 6 opls_140 1 ETH HE6 2 0.06
> 7 opls_140 1 ETH HE7 2 0.06
> 8 opls_140 1 ETH HE8 2 0.06
>
> [ bonds ]
> ; ai aj funct c0 c1
> 1 5 1
> 2 1 1
> 3 1 1
> 4 1 1
> 6 5 1
> 7 5 1
> 8 5 1
>
> [ pairs ]
> ; ai aj funct
> 2 6
> 2 7
> 2 8
> 3 6
> 3 7
> 3 8
> 4 6
> 4 7
> 4 8
>
> [ angles ]
> ; ai aj ak funct c0 c1
> ; H3-C-C
> 2 1 5 1
> 3 1 5 1
> 4 1 5 1
> ; H-C-H
> 2 1 3 1
> 2 1 4 1
> 3 1 4 1
> ; C-C-H3
> 6 5 1 1
> 7 5 1 1
> 8 5 1 1
> ; H-C-H
> 6 5 7 1
> 6 5 8 1
> 7 5 8 1
>
> [ dihedrals ]
> ; ai aj ak al funct c0 c1 c2
> c3 c4
> 2 1 5 6 3
> 3 1 5 6 3
> 4 1 5 6 3
> 2 1 5 7 3
> 3 1 5 7 3
> 4 1 5 7 3
> 2 1 5 8 3
> 3 1 5 8 3
> 4 1 5 8 3
>
>
>
> --
> Phil (Yang) Song
> Room 509, 590 Comm. Ave.
> Chem. Dept., Boston Univ.
> Boston, MA, 02215, USA
>
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
More information about the gromacs.org_gmx-users
mailing list