On 11/29/12 8:18 AM, Ariana Torres Knoop wrote:
Dear GROMACS users,
We are trying to determine the electric fields at the C and N atoms in a
protein. To do so, we thought of performing a one step energy minimization with
the Lennard-Jones as well as the bonded interactions turned off, in order to
get only the coulomb forces.
For that, according to
www.gromacs.org/Documentation/How-tos/Tabulated_Potentials , we made a table
(table.xvg) in which the Lennard-Jones potential was set to zero and we also
made three tables (table_b1.xvg, table_a5.xvg and table_d5.xvg) according to
section 4.2.13 of the Manual (Tabulated interaction functions), in which the
f(x)=0 for the whole range of distances, angles, dihedrals present in the
simulation box. The numbers 1, 5, 2 correspond to the function number referred
to in the topology file. For example for the bonds interaction, in the topology
file we have:
[ bonds ]
; ai aj funct c0 c1 c2 c3
1 2 1
1 3 1
1 4 1
1 5 1
5 6 1
.
.
.
so we named the table "table_b1.xvg".
We changed the run.mdp file as follows:
vdw-type = User
coulombtype = PME-User
And then we try to run it using:
mdrun -v -deffnm run.tpr
When checking the energy contribution using g_energy, Lennard-Jones
contribution was indeed zero, but
for the bonded interactions this is not the case:
Energy Average Err.Est. RMSD Tot-Drift
-------------------------------------------------------------------------------
U-B 18180.1 -- 276.416 -540.638 (kJ/mol)
Proper Dih. 11732.3 -- 27.6829 -54.0396 (kJ/mol)
Improper Dih. 1027.26 -- 47.1598 -91.7097 (kJ/mol)
LJ-14 0 -- 0 0 (kJ/mol)
LJ (SR) 0 -- 0 0 (kJ/mol)
We though that maybe the tables for the bonded interactions were not taken, so
we tried:
mdrun -v -tableb table_b1.xvg table_a5.xvg table_d2.xvg -deffnm run.tpr
But it also did not work.
We are not sure whether the problem is in the table definition, the way of
running the program, the run type, the name of the tables, or if there should
be some other modifications in the *mdp file.
We would be very grateful if anyone could help us out here.
I think the problem comes from the fact that the function types you are
specifying in the topology correspond to actual potentials and thus the tables
are being ignored. Per Table 5.5, it seems that any value for the table is
acceptable (n >= 0) but using numbers corresponding to actual function types may
not work. The specified format also requires that both the table number and
force constant (k) be specified in the topology, but the [bonds] section you
show above has no value of k and the function type corresponds to a simple
harmonic bond. There is no "bond" energy in the output, which suggests you are
using constraints - is that the case?
-Justin
--
========================================
Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
--
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
