Bug #788
Hydrogen bond output from CHARMM is wrong
Description
In version 4.5.4, running g_hbond on a trajectory using CHARMM27 results in incorrect output. g_hbond only detects the terminal amino group as an H-bond donor and no others. I have found that the reason is due to the fact that in CHARMM27, the amide proton is named "HN" for all residues. If one changes the names of these protons in the topology to simply "H" then g_hbond works correctly. I guess the fix would either be to expand the names that g_hbond accepts when identifying H-bond donors, or to change the way the atoms are named in the CHARMM27 residues (.rtp file).
Associated revisions
History
#1 Updated by Justin Lemkul over 9 years ago
- Assignee deleted (
David van der Spoel)
#2 Updated by Mark Abraham over 9 years ago
- Assignee set to Mark Abraham
#3 Updated by Mark Abraham over 9 years ago
- Status changed from New to Feedback wanted
I couldn't reproduce this with my own CHARMM27 trajectory+tpr. 8 donors and 15 acceptors were found for ACE-ARG-ALA-ALA-GLU-NAC, which sounds about right.
Moreover, the test in search_donors()
in src/tools/gmx_hbond.c
looks for the first letter of various atom names, and so "H" and "HN" both satisfy.
#4 Updated by Justin Lemkul over 9 years ago
- File correct.tpr correct.tpr added
- File wrong.tpr wrong.tpr added
- File frame.gro frame.gro added
It finds some, but not all. I'm attaching a test case to see if you can replicate it. With "correct.tpr," the H names are changed to "H" instead of "HN." For "wrong.tpr," the original "HN" name is retained.
g_hbond -s wrong.tpr -f frame.gro ... (select MainChain+H for both groups) Average number of hbonds per timeframe 4.000 out of 16.5 possible
Whereas:
g_hbond -s correct.tpr -f frame.gro (again, MainChain+H) Average number of hbonds per timeframe 29.000 out of 280.5 possible
#5 Updated by Mark Abraham over 9 years ago
This is a problem with how the default index groups are constructed, not with g_hbond. See make_ndx on frame.gro
#6 Updated by Mark Abraham over 9 years ago
- Status changed from Feedback wanted to 3
I fixed the construction of the default groups to be sensitive to this CHARMM-specific detail. This resolves the issue with Justin's input.
#7 Updated by Rossen Apostolov over 9 years ago
- Status changed from 3 to Closed
Adjusted construction of default index groups
The CHARMM convention of naming the H atom bonded to backbone N as
"HN" led to erroneous construction the MainChain+H and (complementary)
SideChain groups. Protein, Protein-H and SideChain-H were fortuitously
correct.
IssueID #788
Change-Id: Ib8525dba61e4f072ae5074a9f6b523ac13892833