[gmx-users] Trying to use user-defined potentials but getting out of range error

[Justin Lemkul]

On 2/9/14, 4:07 PM, Davide Mercadante wrote:
> Hello GMX users,
>
> I am trying to use user defined potentials and to create a table to give
> mdrun I am using a modified version of the C script given here as an
> example:
>
> http://www.gromacs.org/Documentation/How-tos/Tabulated_Potentials
>
> However, I have noticed that even the default script doesn't write out 0
> values for values of r lower than 0.04 as described by the loop and defined
> by the conditional ("if r<0.04 {.....}"). I am getting, by default the
> following output (only first lines are reported):
>
> #
>
> # Example LJ 6-12 Potential
>
> #
>
> 0.0000000000e+00   0.0000000000e+00 1.0000000000e+00   0.0000000000e+00
> 0.0000000000e+00   2.0000009537e+00 4.0000000000e-02
>
> 2.0000000000e-03   8.1920000000e-36 4.6985177795e-18   -1.2474001935e+291
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 4.0000000000e-03   6.7108864000e-32 3.0070513789e-16   -7.9833612381e+292
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 6.0000000000e-03   1.3060694016e-29 3.6700416182e-15   -1.2773377981e+294
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 8.0000000000e-03   5.4975581389e-28 3.8490257650e-14   -1.0218702385e+295
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 1.0000000000e-02   1.0000000000e-26 1.7532907184e-13   -4.0874809539e+295
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 1.2000000000e-02   1.0699320538e-25 4.6976532713e-13   -1.6349923816e+296
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 1.4000000000e-02   7.9371477325e-25 1.7418702204e-12   -3.2699847631e+296
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 1.6000000000e-02   4.5035996274e-24 2.4633764896e-12   -6.5399695263e+296
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 1.8000000000e-02   2.0822964866e-23 5.7024213310e-12   -1.3079939053e+297
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 2.0000000000e-02   8.1920000000e-23 1.1221060598e-11   -2.6159878105e+297
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 2.2000000000e-02   2.8281005788e-22 1.9389445473e-11   -5.2319756210e+297
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 2.4000000000e-02   8.7648833847e-22 3.0064980937e-11   -1.0463951242e+298
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 2.6000000000e-02   2.4811528732e-21 8.5036607586e-11   -2.0927902484e+298
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 2.8000000000e-02   6.5021114225e-21 1.1147969411e-10   -2.0927902484e+298
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 3.0000000000e-02   1.5943230000e-20 1.3695063438e-10   -4.1855804968e+298
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 3.2000000000e-02   3.6893488147e-20 3.1531219067e-10   -8.3711609936e+298
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 3.4000000000e-02   8.1138303246e-20 3.4759466729e-10   -8.3711609936e+298
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 3.6000000000e-02   1.7058172818e-19 7.2990993036e-10   -1.6742321987e+299
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 3.8000000000e-02   3.4449804052e-19 7.3785812839e-10   -1.6742321987e+299
> 0.0000000000e+00   -0.0000000000e+00 4.0000000000e-02
>
> 4.0000000000e-02   2.5000000000e+01 6.2500000000e+02   -2.4414062500e+08
> -3.6621093750e+10   5.9604644775e+16 1.7881393433e+19
>
> 4.2000000000e-02   2.3809523810e+01 5.6689342404e+02   -1.8218149332e+08
> -2.6025927617e+10   3.3190096508e+16 9.4828847166e+18
>
> 4.4000000000e-02   2.2727272727e+01 5.1652892562e+02   -1.3781101808e+08
> -1.8792411556e+10   1.8991876704e+16 5.1796027375e+18
>
>
> For values of r <0.04 the output doesn't seem to be 0 for some of the
> column..however, the script goes like this:
>
> if (r<0.04) {
>              fprintf(fout, "%12.10e   %12.10e %12.10e   %12.10e %12.10e
>    %12.10e %12.10e\n", r,0,0,0,0,0,0);
>          } else {
>              fprintf(fout, "%12.10e   %12.10e %12.10e   %12.10e %12.10e
>    %12.10e %12.10e\n", r,f,fprime,g,gprime,h,hprime);
>
>
>
> Can you please advice why the conditional doesn't seem to be working
> properly or I am overlooking something here?
>

I have fixed the code on the Gromacs site.  It was written a long time ago when 
0 would be interpreted as a floating-point number.  Newer compilers interpret it 
as an integer, generating junk if it even compiles.  I ran the corrected code 
and reproduced the provided table6-12.xvg exactly.

-Justin

-- 
==================================================

Justin A. Lemkul, Ph.D.
Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 601
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

==================================================