Previous Thread
Next Thread
Print Thread
#38459 05/13/22 08:43 PM
Joined: May 2022
Posts: 2
T
txj Offline OP
Forum Member
OP Offline
Forum Member
T
Joined: May 2022
Posts: 2
I am using CHARMM to run simualtions of a pure water box with ~2000 TIP4 water. I found the force outputed in XYZ file is not correct. For example, here are the last several lines of the file:

OH -0.117556615442641E+02 -0.460737086741699E+01 -0.504966436489195E+01 -0.896962810587684E-01 -0.173079977399105E+00 -0.481165302110584E-01 0.000000000000000E+00 0.000000000000000E+00 0.000000000000000E+00
OM -0.118031743479324E+02 -0.451734875393761E+01 -0.502282548502928E+01 -0.871213924496163E-01 -0.144800493693594E+00 -0.138447897817385E+00 0.000000000000000E+00 0.000000000000000E+00 0.000000000000000E+00
H1 -0.113836050592983E+02 -0.372543913170627E+01 -0.505115624108677E+01 0.129661680997219E+00 -0.266193695140701E+00 -0.320327127050701E+00 0.000000000000000E+00 0.000000000000000E+00 0.000000000000000E+00
H2 -0.126565853167218E+02 -0.448726215730393E+01 -0.474942763475455E+01 -0.280393034970152E+00 0.234814024439068E+00 -0.781388701678074E+00 0.000000000000000E+00 0.000000000000000E+00 0.000000000000000E+00
OH -0.920377932247452E+01 -0.292489385283675E+01 0.171612180879107E+02 0.261133883253496E+00 0.641836712759974E-02 -0.224156577754866E+00 0.106645767748745E+01 -0.109439792523077E+01 -0.262052863952775E+01
OM -0.914211865733862E+01 -0.295232157244103E+01 0.172420107575801E+02 0.194001142856658E+00 -0.626893644332873E-01 -0.196464073114893E+00 -0.302238294894556E+01 0.950654595829137E+01 0.117020104788113E+02
H1 -0.829465994835501E+01 -0.266934800539263E+01 0.173175194750189E+02 0.280628971884755E+00 -0.510378296330748E+00 0.505309228793788E+00 -0.128182206452567E+01 -0.471861955630476E+01 -0.299333692575178E+01
H2 -0.942655088720442E+01 -0.348573895752148E+01 0.179042237666623E+02 -0.505618946457559E+00 -0.246026430692034E+00 -0.645375792911392E+00 0.240800395831427E+01 -0.208514960990031E+01 -0.348649064796158E+01


The forces (last three columns) on the second last water are all exact zeros, which is impossible. Actually, the force on them are always zeros throughout the XYZ file. Only the first ~60 water molecules seem to have correct force (I calculated based on leap-frog algorithm, as well as those OMs have zero force).

This is verly likely a bug of outputing force in parallel CHARMM MD:
1) I used 32 cores, and ~2000/32 ~= ~60 (number of molecules with correct forces).
2) When I used one core, everything is correct.


Some other points:

1) I found the direction of the force (in XYZ file) is reversed, i.e., it is actually the force acting on the system by the atom, rather than the force acting on the atom. This is even the case in COOR FORCE. Is this the default setup in CHARMM?
2) The coordinates of the XYZ file is shifted by one step, 1st step in XYZ stores coordinate of the 2nd step..., etc.

Last edited by txj; 05/13/22 09:05 PM.
txj #38460 05/14/22 06:01 PM
Joined: Sep 2003
Posts: 8,605
Likes: 23
rmv Online Content
Forum Member
Online Content
Forum Member
Joined: Sep 2003
Posts: 8,605
Likes: 23
The 'dynamc' info file mentions the CHARMM convention-- "The force on a particle is equal to the negative gradient of the potential energy of the particle."

The remaining issues appear to be bug in the placement of the write statement for the XYZ file within the dynamics control loop done for each integration step. To me, it also suggest the MXYZ output feature is older code that has not been used or rigorously tested for years. Depending on your goals, post-processing stored dynamics coordinates to get the forces might be a better option.


Rick Venable
computational chemist

txj #38462 05/19/22 06:22 PM
Joined: May 2022
Posts: 2
T
txj Offline OP
Forum Member
OP Offline
Forum Member
T
Joined: May 2022
Posts: 2
Thanks, Rick!

I fount that the results from COOR FORCE are a little different from XYZ's, but the resultant force and troque on one molecule are very close.


Moderated by  BRBrooks, lennart, rmv 

Link Copied to Clipboard
Powered by UBB.threads™ PHP Forum Software 7.7.5
(Release build 20201027)
Responsive Width:

PHP: 7.3.31-1~deb10u1 Page Time: 0.008s Queries: 20 (0.005s) Memory: 0.7390 MB (Peak: 0.7858 MB) Data Comp: Off Server Time: 2022-09-26 07:18:32 UTC
Valid HTML 5 and Valid CSS