Project

General

Profile

Bug #624

Error in nb_generic.c causes force miscalculation and/or explosion

Added by Ryan Hayes about 9 years ago. Updated almost 9 years ago.

Status:
Closed
Priority:
Normal
Assignee:
-
Category:
mdrun
Target version:
Affected version - extra info:
Affected version:
Difficulty:
uncategorized
Close

Description

Created an attachment (id=579)
gzipped folder of debugging files

src/gmxlib/nonbonded/nb_generic.c leaves out some forces during calculation leading to system explosion in some cases.

I've been modifying gromacs-4.0.7 to calculate physical quantities, so I've been working with 4.0.7 for consistency, but only one line of nb_generic.c has changed between 4.0.7 and 4.5, so this is likely still a problem in current code.

Substantial differences in the coulomb energy are observed between the optimized kernel and the nb_generic kernel. These differences are over two orders of magnitude are observed during the frist step in the output of md.log from mdrun when using vanilla coulomb cut-off forces. The errors arise because nb_generic only calculates the coulomb forces between the first atom of each charge group. In an explicit solvent bath, where the first atom of all charge groups is a negatively charged oxygen atom, this leads to system explosion after a few dozen steps. The work around is to put each atom in a different charge group, but this is unsatisfactory, especially if undocumented.

Reproducing the bug:
I'm attaching a gzipped folder with the necessary files.

To use nb_generic in my calculations, I set a flag near line 440 of nonbonded.c. Copy this file into src/gmxlib/nonbonded.c

To see a simple system (two water molecules) where the coulomb forces fail, rename spce_orig.itp to spce.itp, grompp mdin_amber_pbc.mdp with SOL_10.top and SOL_10.gro, (SOL_10.itp refers to this spce.itp) then run the result in mdrun. There should be a factor of -169.64 difference between the kernel result and this nb_generic result on Step 0 (in md.log). The discrepancy can be removed by using spce_new.itp where the molecule is split into three charge groups.

Following the same proceedure with the system SOL_8 (a bath with a protein sized hole in the middle) will result in a similar discrepancy (a factor of -360 on Step 0 in md.log), and an explosion after about 30 steps.

Other systems with varying effects are also in the folder, all the systems with only one atom per charge group (in this case protein only systems) give consistent results, but SOL_10 and SOL_8 are the best two systems to see the effect.

Ryan Hayes

nb_generic_files.tar.gz (1.49 MB) nb_generic_files.tar.gz gzipped folder of debugging files Ryan Hayes, 12/22/2010 11:57 AM

History

#1 Updated by Ryan Hayes about 9 years ago

My workaround was to split water atoms into three charge groups. This workaround turns out to be untenable because it causes an unacceptable increase in spontaneous system heating due to small errors introduced by using coulomb cutoff, so this bug is currently stopping progress on my project.

#2 Updated by Ryan Hayes about 9 years ago

I am sorry, this bug has been corrected in later versions (4.5.1) of gromacs. The problem seems to have stemmed from still using solvent optimizations that nb_generic didn't account for, and while nb_generic.c was not changed, functions calling it are smarter about not doing these optimizations when GMX_NB_GENERIC is set.

#3 Updated by Mark Abraham about 9 years ago

  • Assignee deleted (Erik Lindahl)

Ryan Hayes wrote:

I am sorry, this bug has been corrected in later versions (4.5.1) of gromacs. The problem seems to have stemmed from still using solvent optimizations that nb_generic didn't account for, and while nb_generic.c was not changed, functions calling it are smarter about not doing these optimizations when GMX_NB_GENERIC is set.

Yes, that sounds about right.

#4 Updated by Rossen Apostolov almost 9 years ago

  • Status changed from New to Closed

Closing since it's resolved in version 4.5.1.

Also available in: Atom PDF