[gmx-developers] tabulated potential energy=nan for r=0 nm and charged atoms

[Berk Hess]

That depends on the kernel you are using.
If you are using tabulated Coulomb only, no LJ, this is kernel number 300.
In which directory you have to look depends on your architecture.
It might be sse assembly kernels.
In that case you might want to use the generic kernel instead.

Berk

On 05/09/2011 05:30 PM, chris.neale at utoronto.ca wrote:
> Could you please point me to the file and line number so that I can 
> change this to a simple sqrt? If it impacts the timing too much, then 
> I may just leave r=r^2 and remake my table to provide values based on 
> r^2.
>
> Thank you very much for your assistance,
> Chris.
>
> Quoting Berk Hess <hess at cbr.su.se>:
>
>> Hi,
>>
>> In the innerloop r is calculated as rsquared*invsqrt(rsquared).
>>
>> Berk
>>
>> On 05/09/2011 04:54 PM, chris.neale at utoronto.ca wrote:
>>> Dear Developers:
>>>
>>> I find that even when tabulated potentials specify f(x)=0 and  
>>> -f'(x)=0 for x=0, I still get a coulomb(SR) energy=nan when two  
>>> charges occupy the exact same spot.
>>>
>>> Both the manual 4.5.4 page 151 equation 6.23 and the gromacs_nb.pdf 
>>>  tabulated potentials tutorial (  
>>> http://www.gromacs.org/Documentation/How-tos/Tabulated_Potentials ) 
>>>  indicate that setting f(x)=0 and -f'(x)=0 for x=0 should lead to 
>>> an  energy of zero in this case.
>>>
>>> Is this a bug or a case of WYSIWYG? Can anybody offer any hints on  
>>> where I might start searching through the source code with printf  
>>> statements to see where the nan develops? I'm guessing there is  
>>> some if(r>0) in the code that doesn't technically need to be there.
>>>
>>> More information and a test system is available in my users-list  
>>> post for which I received no responses:
>>> http://lists.gromacs.org/pipermail/gmx-users/2011-May/061140.html
>>>
>>> Thank you very much,
>>> Chris.
>>>
>>
>> -- 
>> gmx-developers mailing list
>> gmx-developers at gromacs.org
>> http://lists.gromacs.org/mailman/listinfo/gmx-developers
>> Please don't post (un)subscribe requests to the list. Use the www
>> interface or send it to gmx-developers-request at gromacs.org.
>
>
>