[gmx-users] A question about simulating a box of ethane in GROMACS
Phil (Yang) Song
ysong at bu.edu
Tue Sep 6 21:23:30 CEST 2011
Hi,
Thank everyone for your suggestions.
I have tried to put less molecule in a small box and find this works!
Eventually I have solve the problem by using a larger initial box size for
the random generated position and orientation. It seemed that the overall
interaction of OPLS-AA between ethanes is attractive rather than repulsive,
which leads to atomic clashes...
The original box size was calculated from experimental density but
apparently this does not work for the OPLS-AA. I had used larger initial box
size previously (20% larger than experimental calculated box size) but it
did not work out and then I gave up on this direction. Now, I have tried to
enlarge the box size by 40% and everything works!
Finally I should offer my gratitude to all of you. Thank you!
Best,
Phil
On Tue, Sep 6, 2011 at 10:59 AM, Mark Abraham <Mark.Abraham at anu.edu.au>wrote:
> On 7/09/2011 12:40 AM, 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.
>
>
> That's almost certainly inappropriate. Atomic clashes are inevitable, and
> EM will not fix severe problems.
>
> Get a single molecule in a small box, and use genconf to replicate it. Then
> use EM and equilibrate thoroughly to remove the residual ordering. Judicious
> use of NPT will fix your density to whatever you want.
>
> Mark
>
>
> 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.
>
> 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.
>
> 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
>
>
>
>
> --
> gmx-users mailing list gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
--
Phil (Yang) Song
Room 509, 590 Comm. Ave.
Chem. Dept., Boston Univ.
Boston, MA, 02215, USA
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20110906/a4c0e13d/attachment.html>
More information about the gromacs.org_gmx-users
mailing list