HBOND_MATRIX gives too many hydrogen bonds

450 views
Skip to first unread message

dan.k...@gmail.com

unread,
Jun 16, 2017, 2:49:40 AM6/16/17
to PLUMED users
Hello,

I'm am trying to use the number of backbone hydrogen bonds as a collective variable, but HBOND_MATRIX is giving me more than 500 hydrogen bonds when there should only be ~10. I read this thread: (https://groups.google.com/forum/#!msg/plumed-users/gZRPu_8HX7s/CV5MhJEZBAAJ) and looked at the example here (https://plumed.github.io/doc-v2.3/user-doc/html/_h_b_o_n_d__m_a_t_r_i_x.html) to get my input file as shown below:

# Activate metadynamics 

MOLINFO STRUCTURE=protein.pdb
WHOLEMOLECULES RESIDUES=all MOLTYPE=protein

hatoms: GROUP ATOMS=1,2,3,...
oatoms: GROUP ATOMS=4,5,6,...
hbond: HBOND_MATRIX ATOMS=oatoms HYDROGENS=hatoms SWITCH={RATIONAL R_0=3.2} HSWITCH={RATIONAL R_0=2.3} ASWITCH={RATIONAL R_0=0.167pi} SUM
rsums: ROWSUMS MATRIX=hbond MEAN
csums: COLUMNSUMS MATRIX=hbond MEAN

PRINT STRIDE=1 ARG=rsums.mean,csums.mean FILE=COLVAR


Does anyone have any suggestions as to what I'm doing wrong? I'm mainly confused about the SWITCH parameters since I thought the PLUMED default unit was nm and anything bigger than 1 seems like it would be too long for a hydrogen bond, but I'm also not really understanding how to properly process the contact matrix. 

Best,
Dan

Carlo Camilloni

unread,
Jun 16, 2017, 3:08:10 AM6/16/17
to plumed...@googlegroups.com
I think that R_0=3.2 is wrong here because standard unit for distance is nm. Try to use R_0=0.32

Carlo


--
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.
To post to this group, send email to plumed...@googlegroups.com.
Visit this group at https://groups.google.com/group/plumed-users.
To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/674d879a-f200-4474-9e39-cdaf49b1ef73%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Daniel Kozuch

unread,
Jun 16, 2017, 3:23:33 AM6/16/17
to plumed...@googlegroups.com
That was my initial reaction. I also changed 2.3 to 0.23, but these setting gave ~0 hydrogen bonds (actually something like 0.05, I assume because the switching is continuous). Somewhere in between 0.32 and 3.2 is where I get a reasonable number of hydrogen bonds, but I don't want to pick the cutoffs arbitrarily. I thought the problem might be that the example uses COLUMNSUMS with the option MEAN. Shouldn't this be some form of SUM (why would we want the mean of the column?).

Dan

On Fri, Jun 16, 2017 at 3:08 AM, Carlo Camilloni <carlo.c...@gmail.com> wrote:
I think that R_0=3.2 is wrong here because standard unit for distance is nm. Try to use R_0=0.32

Carlo
On 16 Jun 2017, at 08:49, dan.k...@gmail.com wrote:

Hello,

I'm am trying to use the number of backbone hydrogen bonds as a collective variable, but HBOND_MATRIX is giving me more than 500 hydrogen bonds when there should only be ~10. I read this thread: (https://groups.google.com/forum/#!msg/plumed-users/gZRPu_8HX7s/CV5MhJEZBAAJ) and looked at the example here (https://plumed.github.io/doc-v2.3/user-doc/html/_h_b_o_n_d__m_a_t_r_i_x.html) to get my input file as shown below:

# Activate metadynamics 

MOLINFO STRUCTURE=protein.pdb
WHOLEMOLECULES RESIDUES=all MOLTYPE=protein

hatoms: GROUP ATOMS=1,2,3,...
oatoms: GROUP ATOMS=4,5,6,...
hbond: HBOND_MATRIX ATOMS=oatoms HYDROGENS=hatoms SWITCH={RATIONAL R_0=3.2} HSWITCH={RATIONAL R_0=2.3} ASWITCH={RATIONAL R_0=0.167pi} SUM
rsums: ROWSUMS MATRIX=hbond MEAN
csums: COLUMNSUMS MATRIX=hbond MEAN

PRINT STRIDE=1 ARG=rsums.mean,csums.mean FILE=COLVAR


Does anyone have any suggestions as to what I'm doing wrong? I'm mainly confused about the SWITCH parameters since I thought the PLUMED default unit was nm and anything bigger than 1 seems like it would be too long for a hydrogen bond, but I'm also not really understanding how to properly process the contact matrix. 

Best,
Dan


--
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users+unsubscribe@googlegroups.com.

--
You received this message because you are subscribed to a topic in the Google Groups "PLUMED users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/plumed-users/zoJ_wn6_qVY/unsubscribe.
To unsubscribe from this group and all its topics, send an email to plumed-users+unsubscribe@googlegroups.com.

To post to this group, send email to plumed...@googlegroups.com.
Visit this group at https://groups.google.com/group/plumed-users.

Gareth Tribello

unread,
Jun 16, 2017, 4:46:21 AM6/16/17
to plumed...@googlegroups.com
Hello

If you just want the total number of hydrogen bonds you might be better off doing just printing hbond.sum as this is the sum of the elements of the hbond matrix.  The values output by rowsums and columnsums are the average number of hydrogen bonds per atom in your system.

In terms of picking the cutoffs arbitrarily.  You have to look at what is the typical length of a hydrogen bond in your system and pick a cutoff to suit that.  To a certain extent this choice is always somewhat arbitrary.

Gareth



To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.

To post to this group, send email to plumed...@googlegroups.com.
Visit this group at https://groups.google.com/group/plumed-users.

Daniel Kozuch

unread,
Jun 16, 2017, 11:28:43 AM6/16/17
to plumed...@googlegroups.com
Thanks, that seems to be working.

Best,
Dan

On Fri, Jun 16, 2017 at 4:46 AM, Gareth Tribello <gareth....@gmail.com> wrote:
Hello

If you just want the total number of hydrogen bonds you might be better off doing just printing hbond.sum as this is the sum of the elements of the hbond matrix.  The values output by rowsums and columnsums are the average number of hydrogen bonds per atom in your system.

In terms of picking the cutoffs arbitrarily.  You have to look at what is the typical length of a hydrogen bond in your system and pick a cutoff to suit that.  To a certain extent this choice is always somewhat arbitrary.

Gareth


On 16 Jun 2017, at 08:22, Daniel Kozuch <dko...@princeton.edu> wrote:

That was my initial reaction. I also changed 2.3 to 0.23, but these setting gave ~0 hydrogen bonds (actually something like 0.05, I assume because the switching is continuous). Somewhere in between 0.32 and 3.2 is where I get a reasonable number of hydrogen bonds, but I don't want to pick the cutoffs arbitrarily. I thought the problem might be that the example uses COLUMNSUMS with the option MEAN. Shouldn't this be some form of SUM (why would we want the mean of the column?).

Dan

On Fri, Jun 16, 2017 at 3:08 AM, Carlo Camilloni <carlo.camilloni@gmail.com> wrote:
I think that R_0=3.2 is wrong here because standard unit for distance is nm. Try to use R_0=0.32

Carlo
On 16 Jun 2017, at 08:49, dan.k...@gmail.com wrote:

Hello,

I'm am trying to use the number of backbone hydrogen bonds as a collective variable, but HBOND_MATRIX is giving me more than 500 hydrogen bonds when there should only be ~10. I read this thread: (https://groups.google.com/forum/#!msg/plumed-users/gZRPu_8HX7s/CV5MhJEZBAAJ) and looked at the example here (https://plumed.github.io/doc-v2.3/user-doc/html/_h_b_o_n_d__m_a_t_r_i_x.html) to get my input file as shown below:

# Activate metadynamics 

MOLINFO STRUCTURE=protein.pdb
WHOLEMOLECULES RESIDUES=all MOLTYPE=protein

hatoms: GROUP ATOMS=1,2,3,...
oatoms: GROUP ATOMS=4,5,6,...
hbond: HBOND_MATRIX ATOMS=oatoms HYDROGENS=hatoms SWITCH={RATIONAL R_0=3.2} HSWITCH={RATIONAL R_0=2.3} ASWITCH={RATIONAL R_0=0.167pi} SUM
rsums: ROWSUMS MATRIX=hbond MEAN
csums: COLUMNSUMS MATRIX=hbond MEAN

PRINT STRIDE=1 ARG=rsums.mean,csums.mean FILE=COLVAR


Does anyone have any suggestions as to what I'm doing wrong? I'm mainly confused about the SWITCH parameters since I thought the PLUMED default unit was nm and anything bigger than 1 seems like it would be too long for a hydrogen bond, but I'm also not really understanding how to properly process the contact matrix. 

Best,
Dan


-- 
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users+unsubscribe@googlegroups.com.

To post to this group, send email to plumed-users@googlegroups.com.

-- 
You received this message because you are subscribed to a topic in the Google Groups "PLUMED users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/plumed-users/zoJ_wn6_qVY/unsubscribe.
To unsubscribe from this group and all its topics, send an email to plumed-users+unsubscribe@googlegroups.com.
To post to this group, send email to plumed-users@googlegroups.com.

-- 
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users+unsubscribe@googlegroups.com.
To post to this group, send email to plumed-users@googlegroups.com.

--
You received this message because you are subscribed to a topic in the Google Groups "PLUMED users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/plumed-users/zoJ_wn6_qVY/unsubscribe.
To unsubscribe from this group and all its topics, send an email to plumed-users+unsubscribe@googlegroups.com.
To post to this group, send email to plumed...@googlegroups.com.
Visit this group at https://groups.google.com/group/plumed-users.
To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/F8AA7171-9818-4E26-BFAB-4FA2482E2BA5%40gmail.com.

Lakshman Verma

unread,
Feb 10, 2021, 2:22:04 AM2/10/21
to PLUMED users
Hi,
I am trying to calculate hydrogen bond in my system as a function of distance from a surface. I followed previous threads https://groups.google.com/g/plumed-users/c/zoJ_wn6_qVY/m/gjFa-OHJAAAJ .This solved problem of getting too many h-bonds. However, the number still doesn't make sense. Each water can have max 4 hbonds in a frame (2 donated and 2 accepted), however, the column sum and row sum both give more than 4 hbonds per atom for a test bulk water droplet. I am not able to understand why this is happening. 

My hbond parameters look like following:
mat: HBOND_MATRIX ATOMS=1-3000:3 HYDROGENS=2-3000:3,3-3000:3 SWITCH={RATIONAL R_0=0.350} HSWITCH={RATIONAL R_0=0.230} ASWITCH={RATIONAL R_0=0.167pi} SUM
rsums: ROWSUMS MATRIX=mat HIGHEST
csums: COLUMNSUMS MATRIX=mat HIGHEST
DUMPMULTICOLVAR DATA=rsums FILE=donors.xyz
DUMPMULTICOLVAR DATA=csums FILE=acceptors.xyz
PRINT ARG=* FILE=colvar

A sample output of acceptors.xyz
1000
 7.000000 7.000000 7.000000
X 2.520000 2.928000 2.848000 6.008639
X 2.422000 3.794000 3.493000 4.220039
X 4.033000 4.558000 2.363000 4.164962
X 4.182000 4.676000 3.248000 4.601371
X 2.041000 3.148000 2.921000 4.462403
X 4.626000 2.586000 3.895000 5.487307
X 3.275000 4.226000 2.344000 3.814607
X 4.197000 4.230000 5.448000 3.302398
X 3.845000 3.840000 5.278000 3.802542
X 3.938000 2.382000 3.106000 5.198168
X 4.444000 3.623000 2.889000 4.401865

Also, I noticed that the sum of the number of hbonds dumped over all dumped atoms in column sum is the same as that in the row sum for a frame. Is it because the Hbond matrix would be symmetrical? 
Also, is there a way to dump the actual number of hbonds and not the mean for each water with their coordinates, as in the sample output? 

Thank you,
Lakshman

On Friday, June 16, 2017 at 10:28:43 AM UTC-5 dko...@princeton.edu wrote:
Thanks, that seems to be working.

Best,
Dan

On Fri, Jun 16, 2017 at 4:46 AM, Gareth Tribello <gareth....@gmail.com> wrote:
Hello

If you just want the total number of hydrogen bonds you might be better off doing just printing hbond.sum as this is the sum of the elements of the hbond matrix.  The values output by rowsums and columnsums are the average number of hydrogen bonds per atom in your system.

In terms of picking the cutoffs arbitrarily.  You have to look at what is the typical length of a hydrogen bond in your system and pick a cutoff to suit that.  To a certain extent this choice is always somewhat arbitrary.

Gareth


On 16 Jun 2017, at 08:22, Daniel Kozuch <dko...@princeton.edu> wrote:

That was my initial reaction. I also changed 2.3 to 0.23, but these setting gave ~0 hydrogen bonds (actually something like 0.05, I assume because the switching is continuous). Somewhere in between 0.32 and 3.2 is where I get a reasonable number of hydrogen bonds, but I don't want to pick the cutoffs arbitrarily. I thought the problem might be that the example uses COLUMNSUMS with the option MEAN. Shouldn't this be some form of SUM (why would we want the mean of the column?).

Dan

On Fri, Jun 16, 2017 at 3:08 AM, Carlo Camilloni <carlo.c...@gmail.com> wrote:
I think that R_0=3.2 is wrong here because standard unit for distance is nm. Try to use R_0=0.32

Carlo
On 16 Jun 2017, at 08:49, dan.k...@gmail.com wrote:

Hello,

I'm am trying to use the number of backbone hydrogen bonds as a collective variable, but HBOND_MATRIX is giving me more than 500 hydrogen bonds when there should only be ~10. I read this thread: (https://groups.google.com/forum/#!msg/plumed-users/gZRPu_8HX7s/CV5MhJEZBAAJ) and looked at the example here (https://plumed.github.io/doc-v2.3/user-doc/html/_h_b_o_n_d__m_a_t_r_i_x.html) to get my input file as shown below:

# Activate metadynamics 

MOLINFO STRUCTURE=protein.pdb
WHOLEMOLECULES RESIDUES=all MOLTYPE=protein

hatoms: GROUP ATOMS=1,2,3,...
oatoms: GROUP ATOMS=4,5,6,...
hbond: HBOND_MATRIX ATOMS=oatoms HYDROGENS=hatoms SWITCH={RATIONAL R_0=3.2} HSWITCH={RATIONAL R_0=2.3} ASWITCH={RATIONAL R_0=0.167pi} SUM
rsums: ROWSUMS MATRIX=hbond MEAN
csums: COLUMNSUMS MATRIX=hbond MEAN

PRINT STRIDE=1 ARG=rsums.mean,csums.mean FILE=COLVAR


Does anyone have any suggestions as to what I'm doing wrong? I'm mainly confused about the SWITCH parameters since I thought the PLUMED default unit was nm and anything bigger than 1 seems like it would be too long for a hydrogen bond, but I'm also not really understanding how to properly process the contact matrix. 

Best,
Dan


-- 
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.

To post to this group, send email to plumed...@googlegroups.com.

-- 
You received this message because you are subscribed to a topic in the Google Groups "PLUMED users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/plumed-users/zoJ_wn6_qVY/unsubscribe.
To unsubscribe from this group and all its topics, send an email to plumed-users...@googlegroups.com.
To post to this group, send email to plumed...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

-- 
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.
To post to this group, send email to plumed...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to a topic in the Google Groups "PLUMED users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/plumed-users/zoJ_wn6_qVY/unsubscribe.
To unsubscribe from this group and all its topics, send an email to plumed-users...@googlegroups.com.

To post to this group, send email to plumed...@googlegroups.com.
Visit this group at https://groups.google.com/group/plumed-users.
Reply all
Reply to author
Forward
0 new messages