[gmx-users] RE: LIE methodology check
Mark.Abraham at anu.edu.au
Wed May 23 03:34:50 CEST 2012
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.
> All the best,
> Message: 2
> Date: Mon, 21 May 2012 03:26:07 +0000
> From: Tom Dupree<t.dupree at unsw.edu.au>
> Subject: [gmx-users] LIE methodology check (mdrun -rerun cutoffs and
> box vectors)
> To: "gmx-users at gromacs.org"<gmx-users at gromacs.org>
> <0A0E67195A48DC4A96321968E2B384B024BF646D at INFPWXM007.ad.unsw.edu.au>
> Content-Type: text/plain; charset="windows-1252"
> Greetings all,
> 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.
> 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,
More information about the gromacs.org_gmx-users