[gmx-users] Tabulated potentials make newbies crazy

[ms]

Dear Mark,

First of all, thanks for your patient and thorough reply!

Mark Abraham ha scritto:
> Yeah, that's an unfortunate fact of life. The manual can't be organized
> so that everybody finds all the information they want in one location
> with the detail they need right now. Using tabulated potentials is an
> advanced feature and if the user needs to search the manual diligently
> to find the information, that's probably good for them in the long run.
> There's a good table of contents, index and some cross references, and
> text-searching is also available...

Yes, that's what I am doing now but it seems I am not good enough at it
:) .Don't get me wrong, the manual is truly, amazingly, incredibly good
-but probably a wiki page on the subject could help people tie the
strings together more easily. When I'll get all of this correct I wish
to draft one so that you can correct it.

>> Just to try, therefore, I removed all occurences of energy groups in a
>> .ndx files, except for system and C-alpha (I have no solvent). And there
>> it complains that atoms belong to different energy groups, but belong to
>> the same charge group.
>>
>> I understand therefore that energy groups and charge groups must be
>> identically defined.
> 
> No. You should read in the manual about groups - sections 3.3 and 8.1.1

I've read them -read both of them *now* again, actually- and still I see
no mention of charge groups in the paragraph you linked -apart perhaps
for this sentence in 3.3: " This is done separately for Lennard-Jones
and Coulomb terms." -They are described in 3.4 ,and also in 4.6 for
example, but I find no information on how they're differently described
with respect to energy groups. It tells me they're described in the
topology, but nothing more. If I look at my topology, I find no charge
group definition.

The problem is that since I have a single molecule now, and the single
molecule must be neutral, so it must be all a single charge group
("Therefore we have to keep groups of atoms with total charge 0
together. These groups are called charge groups.", 4.6.2).

> Get the energy group stuff working before you worry about providing user
> tables. Divide-and-conquer makes troubleshooting much easier.

Right, thanks...

>> Anyway, just for the sake of learning how to get gmx use a tabulated
>> potential, I apply the tabulated potential to all the peptide chain,
>> (and I rename the table table_protein_protein.xvg). At this point, it
>> stops asking me for a "table.xvg" file. I realize that it probably wants
>> to know how to define other non bonded interaction. Since there is only
>> the protein, I rename the table_protein_protein.xvg as table.xvg and
>> retry.
> 
> Group names are probably case sensitive, so your filename must match
> your group names... the above doesn't look right.

The group name was in lower case, that is why I called it in lower case
(for C-alpha I used upper case as in the group name). It recognized it
somehow. But thanks for pointing it, I'll retry.

>> And now I cry because mdrun doesn't start, telling me:
>>
>>    Fatal error:
>>    Out of range force value 2.44081e+41 in file 'table.xvg'
>>
>> The "out of range force value" error does not appear at all on google. I
>> understand it's "kinda high" - but it's a very steep repulsion and so it
>> is naturally like that at distances close to 0.
> 
> That looks like a bug with mdrun trying to read from a non-existent
> file. It is always useful to provide your exact command lines, rather
> than leave us to infer from text.

Sorry. This is from the script I use to launch the job on the cluster
(it is not, however, a parallel job):

mdrun -v -s $STARTDIR/fullmd_xab10_01.tpr -o $STARTDIR/xab10_traj01.trr
-x $STARTDIR/xab10_traj01.xtc -c $STARTDIR/xab10_final01.gro -e
$STARTDIR/xab10_ener01.edr 2>&1 > $STARTDIR/output_20091127_xab01.txt

Now however, thanks to your suggestion, I checked my table and I found
that I defined a derivative (the 2.44081e+41 value indeed) for r=0.

Adjusting that (putting it to 0, since r=0 makes little sense), it
complains for the next, huge value, at r=0.002. The values it complains
about are always real ones in the table, so it finds it and reads it. I
can now regenerate it however such that it stays constant when energy is
larger, say, than 50.000. This is probably the least problem, but some
advice on it, just for the sake of completeness, could be nice.

>> Now, my questions are:
>> - What is the accepted range of values in tables?
> 
> I don't think this is the problem.

It is the least problem probably, given my confusion on energy-charge
groups, but it seems it is too...

>> - How do I define a steep repulsion potential correctly?
> 
> It's terse, but manual 6.7 seems to have the necessary information.

6.7 is one of the references I am obviously using, but it gives only
general (even if essential!!) information, nothing speficic on "good" or
"bad" potential shapes/values. But probably that's the least problem :)

>> - Is there a way to define a table for specific atoms in a chain,
>> leaving the other atoms using the standard VdW potential? I guess it
>> requires splitting the chain in two energy groups, but I don't
>> understand if it a)makes sense/creates problems b)how to deal with the
>> fact that the chain should be a single charge group.
> 
> Use energygrp-table like you've been trying. See 7.3.12

That is also something I keep under my pillow, but the charge/energy
group confusion I must resolve, before going on...

Thanks a lot and sorry for my stubborness -learning MD almost by
yourself is not the best of worlds, but that's how my laboratory seems
to work...

cheers,

m.