On 23/05/2012 8:36 AM, Tom Dupree wrote:
Greetings all,
I now realise I asked too many questions with that one post.
I would like to check that mdrun -rerun can be used to do single point energy
calculations on the co-ordinate frames from my trajectory using a different
electrostatic model by changing the .mdp file used by grompp to build the .tpr
file?
Yes. You can observe differences from your original run, even with
identical .tpr and thus model, because the rerun does (for example)
neighbour searching for every frame, since it cannot assume atoms in
nearby frames share spatial locality. Obviously there is no way to know
if that frame would be present in the ensemble produced by the new .tpr,
but you can measure the energy if you want to.
Mark
All the best,
Tom
___________________________________________
Message: 2
Date: Mon, 21 May 2012 03:26:07 +0000
From: Tom Dupree<t.dup...@unsw.edu.au>
Subject: [gmx-users] LIE methodology check (mdrun -rerun cutoffs and
box vectors)
To: "gmx-users@gromacs.org"<gmx-users@gromacs.org>
Message-ID:
<0a0e67195a48dc4a96321968e2b384b024bf6...@infpwxm007.ad.unsw.edu.au>
Content-Type: text/plain; charset="windows-1252"
Greetings all,
Background:
I have a series of MD runs of a protein ligand complex where I have used PME,
they represent far more calculations than I have time to repeat in full. I
recently became aware that the coulombic interaction energy reported when PME
is used is not valid for use with g_lie (or any LIE method). I have read some
relevant LIE papers and it seems that reaction field is the electrostatic
method used.
http://pubs.acs.org/doi/abs/10.1021/jm7012198
Current methodology:
I found that using mdrun ?rerun and a new .tpr file I can change reported
energies. Are these energies corresponding to the settings in the new .mdp
file? Can I change it from PME to Generalised-Reaction-Field? Are there any
errors associated with changing it, as opposed to errors implicit in either
method? I am understanding ?rerun to only be calculated for the co-ordinate
frames that are outputted into the .xtc file irrespective of any settings in
the new .mdp file and/or old .mdp file for the energy.
As I only outputted .xtc files I don?t have velocities (or did I set that in
the .mdp file?) but I assume that errors in kinetic energy, temperature and
pressure are not important in the rerun as I am only concerned with interaction
energy?
Right assuming all the above is valid I would like to improve my accuracy as
much as possible.
My understanding of cut offs is that the bigger they are the more accurate and
the more computationally expensive your calculations are until they approach
half the shortest box vector where pbc issues cause errors. As I am only doing
reruns now I should use the largest possible cut offs that don?t hit the pbc
limit?
Is there an accuracy difference between the short range and long range
calculations for the LJ potential? (is there a reason to use long range at all?)
For the reruns can I ignore the problem of a water molecule seeing both sides
of my protein and use a cut off greater than the editconf -d setting (in my
case 1.6)? Or to be strictly accurate should I work out the distance to the
other side of the mirror protein from the binding site via the pbc boundary and
use half that?
Which brings me to my issue with the unit cell. I set -bt dodecahedron. But I
don't see it as one (vmd/pymol display of water atoms/unit cell). How do I
relate the vectors in the.gro file to the pbc arrangement?
------------------------------------
2584CL CL 7792 5.230 2.213 3.017 0.0075 0.3566 0.5809
2585CL CL 7793 1.448 2.710 1.489 0.5771 0.2054 0.0146
2586CL CL 7794 3.070 1.270 3.024 -0.1756 -0.0490 -0.1351
2587CL CL 7795 5.112 2.558 1.562 -0.1543 -0.0592 -0.0072
2588CL CL 7796 2.906 1.133 1.094 -0.2265 -0.3952 0.1190
4.80953 4.80953 3.40085 0.00000 0.00000 0.00000 0.00000
2.40476 2.40476
----------------------
In my trials I got the error
?ERROR 1 [file topol.top, line 50]:
ERROR: The cut-off length is longer than half the shortest box vector or
longer than the smallest box diagonal element. Increase the box size or
decrease rlist.?
So if I understand the .gro file correctly (now) these vectors are the long
diagonals of three rhombuses. Two are at right angles to each other and the
third connects two a lower edge on each of those. Trial and error lead me to
the maximum cut off that grompp will accept is 2.40476, which is d/2, if I read
this correctly. My question is, why don't my waters and/or unit cell appear to
be a dodecahedron when viewed in VMD/pymol? For that matter why don't the
waters and unit cell agree on a shape when viewed in VMD? (both are rhombic
prisims).
All the best,
Tom
--
gmx-users mailing list gmx-users@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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
