Bug #177
position and velocity set to nan
Description
it is an customized FF with section
[ dihedraltypes ]
; i l func q0 cq
;added by Dongsheng
OM OM 0 0.000 1 2
OA OA 0 0.000 1 2
OM OA 0 0.000 1 2
in ffpolymerbon.itp
If there is entries like following in topol.top
[ dihedrals ]
1 2 3 4 1
It will result in nan in position and velocity and cause error
-------------------------------------------------------
Program mdrun_s, VERSION 3.3.1
Source code file: nsgrid.c, line: 220
Fatal error:
Number of grid cells is zero. Probably the system and box collapsed.
-------------------------------------------------------
I have verified in bondfree.c that coordinates are the same after calling pdihs.
If there is no dihedral defined in topol.top. The simulation is alright.
Is there any other place which changes the coordinate?
History
#1 Updated by Yang Ye about 12 years ago
Created an attachment (id=257)
test case
test file including force field and run files.
type "grompp" and "mdrun" to run the case.
#2 Updated by David van der Spoel about 12 years ago
This is because your dihedrals are exactly 0, that is all 4 atoms are co-linear. With more "physical" input this wouldn't happen.
#3 Updated by David van der Spoel about 12 years ago
OK, I have made a fix where the dihedral gives zero force if it is completely linear. That means the simulation will not crash but the forces may be discontinuous close to a linear molecule. In addition the linear angles will also give you a discontinuous force close to 180, basically because we disallow 180 degrees angles altogether. This is not a pretty solution but it's the best I can do.
The dihedral angle is defined as the angle between two planes, however we cannot make these planes since the vectors making up the planes are on a line, hence the planes are undefined.
Erik, shall I commit this fix?
#4 Updated by Erik Lindahl about 12 years ago
Hi,
That sounds like the only reasonably alternative, so please do!
#5 Updated by David van der Spoel about 12 years ago
Created an attachment (id=265)
gmx/include/vec.h
New vec.h
#6 Updated by David van der Spoel about 12 years ago
Created an attachment (id=266)
gmx/src/gmxlib/bondfree.c
#7 Updated by David van der Spoel about 12 years ago
With the two file attached the problem is solved, but please keep in mind my previous comments.
#8 Updated by Yang Ye about 12 years ago
Thanks. I am aware for the effects of the fix.