Joshua Ballanco wrote:
On Feb 1, 2009, at 8:22 PM, Mark Abraham wrote:
Joshua Ballanco wrote:
On Feb 1, 2009, at 6:48 PM, Mark Abraham wrote:
Joshua Ballanco wrote:
Hi,
I'm attempting to model a system involving a small DNA 3-mer.
Without any explicit constraints, the helix begins to come apart
around 0.75 ns to 1 ns into the simulation.
Presumably you have a 3-mer of helices, of which at least one comes
apart. Does a single helix in water survive? (Giving a better
description of your simulation system would be a good idea!)
Apologies for not being more descriptive... It is a single strand of
DNA containing 3 A-T base pairs. The system also contains a single
Arginine residue. Simulating it in water leads to the single DNA
strand gradually coming apart. With the DNA and water alone, the
helix stays together much longer, but still eventually comes apart.
OK, so your model physics for DNA is intrinsically broken. Where did
you get it?
The coordinates are from 3DNA. I'm using the terms from the G53a6 force
field for DADE and DTHY. As for the H-bond physics, I've thus far been
unable to find a good suggestion for how to handle these explicitly
using Gromacs. Do you have a recommendation as to where I should be
looking? (None of the primary literature I've looked through thus far
seems concerned with MD of such short DNA fragments).
No, I've no idea since I don't simulate DNA.
So I'm now attempting to add restraints for the base-pair H-bonds,
but I'm having trouble. It seems like no matter what I try, my
system reliably explodes within the first 1 ns. My constraints look
like this:
[ distance_restraints ]
; ai aj type index type’ low up1 up2 fac
18 136 1 0 2 0.0 2.0 2.1 1.0
14 134 1 0 2 0.0 2.0 2.1 1.0
43 114 1 0 2 0.0 2.0 2.1 1.0
39 112 1 0 2 0.0 2.0 2.1 1.0
68 92 1 0 2 0.0 2.0 2.1 1.0
64 90 1 0 2 0.0 2.0 2.1 1.0
I've tried pre-equilibrating for up to 100 ps, but even that
doesn't prevent the system from eventually exploding.
Your .mdp settings for distance restraints may also be relevant here
- not least in setting the existence and magnitude of these restraints.
As I understand, the only relevant lines are:
constraints = all-bonds
integrator = md
disre = simple
disre-fc and others are also relevant. See manual chapter 7.
Thanks for the pointer. I had overlooked most of the options there,
since I'm not actually doing anything related to NMR. (That'll teach me
to read more carefully!) Unfortunately, playing around with this,
disre-tau, disre-weighting, and the weighting factors for each bond have
not, so far, avoided the explosion.
OK, that's no longer surprising - distance restraints will not usefully
fix a broken model physics.
For PME I was using:
coulombtype = PME
rlist = 0.55
rcoulomb = 0.55
rvdw = 0.55
fourierspacing = 0.1375
I agree with Justin that these are very weird for normal usage.
Thanks for pointing that out. I'm relatively new with Gromacs, and
hastily reduced these values to fix the relatively small box my system
fits in. I doubled the short box dimension (triclinic; was --> 2.0, 2.1,
1.1 now --> 2.0, 2.1, 2.0) and increased the radii to the (as far as I
can tell) more recommended values:
coulombtype = PME
rlist = 0.9
rcoulomb = 0.9
rvdw = 0.9
fourierspacing = 0.12
Well, that's more like it. Values for these parameters are actually
intrinsic to the forcefield parametrization process, and one should vary
them only with caution. This algorithmic constraint sets a minimum size
for the simulation, of course.
When using PBC, just fitting your system into a box doesn't address the
real issue. In a real solution this 3-mer would be close to infinite
dilution, which can't be modeled without a serious chunk of solvent
around it. This consideration dwarfs the algorithmic one I refer to above.
Unfortunately, even with all of these changes, I'm still getting an
explosion (and my simulation is quite a bit slower).
Anyone can get quick random numbers - you don't even need a simulation
package :-P There's no substitute for background literature reading,
doing tutorials, and experimenting with preparedness for failure. :-)
Thanks again for the pointers. I'm going to try running everything with
ffamber to see if it does a better job with the DNA (without the added
restraints). I presume the port validated with Gromacs 3.3.1 is still
good for 4.0?
Don't know - check out the documentation about the forcefield port -
search the web.
Mark
_______________________________________________
gmx-users mailing list [email protected]
http://www.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to [email protected].
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
