diff --git a/sPuReMD/src/ffield.c b/sPuReMD/src/ffield.c
index cab5dcacbe11b8798d19ca6f61f33f7efedf6590..e4a7d144fe59d355caf579145c437be70d8f14dc 100644
--- a/sPuReMD/src/ffield.c
+++ b/sPuReMD/src/ffield.c
@@ -886,6 +886,17 @@ void Read_Force_Field( const char * const ffield_file,
         Tokenize( s, &tmp, MAX_TOKEN_LEN );
         l = sstrtol( tmp[0], __FILE__, __LINE__ );
 
+        for ( i = 0; i < reax->max_num_atom_types; ++i )
+        {
+            for ( j = 0; j < reax->max_num_atom_types; ++j )
+            {
+                for ( k = 0; k < reax->max_num_atom_types; ++k )
+                {
+                    reax->hbp[i][j][k].is_valid = FALSE;
+                }
+            }
+        }
+
         for ( i = 0; i < l; i++ )
         {
             if ( fgets( s, MAX_LINE, fp ) == NULL )
@@ -914,6 +925,7 @@ void Read_Force_Field( const char * const ffield_file,
                 val = sstrtod( tmp[6], __FILE__, __LINE__ );
                 reax->hbp[j][k][m].p_hb3 = val;
 
+                reax->hbp[j][k][m].is_valid = TRUE;
             }
         }
 
diff --git a/sPuReMD/src/hydrogen_bonds.c b/sPuReMD/src/hydrogen_bonds.c
index 14d30b2489cd11dceb1fce0838eb7242914477ae..45392a52c9307325afe8a452db78931670ab7212 100644
--- a/sPuReMD/src/hydrogen_bonds.c
+++ b/sPuReMD/src/hydrogen_bonds.c
@@ -51,7 +51,7 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
         int i, j, k, pi, pk, itr, top;
         int type_i, type_j, type_k;
         int start_j, end_j, hb_start_j, hb_end_j;
-        int hblist[MAX_BONDS];
+        int *hblist, hblist_size;
         real r_ij, r_jk, theta, cos_theta, sin_xhz4, cos_xhz1, sin_theta2;
         real e_hb, exp_hb2, exp_hb3, CEhb1, CEhb2, CEhb3;
         rvec dcos_theta_di, dcos_theta_dj, dcos_theta_dk;
@@ -70,6 +70,8 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
         int tid = omp_get_thread_num( );
 #endif
 
+        hblist = NULL;
+        hblist_size = 0;
         bonds = lists[BONDS];
         bond_list = bonds->bond_list;
         hbonds = lists[HBONDS];
@@ -110,6 +112,13 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
                 }
 
                 top = 0;
+                if ( Num_Entries( j, bonds ) > hblist_size )
+                {
+                    hblist_size = Num_Entries( j, bonds );
+                    hblist = srealloc( hblist, sizeof(int) * hblist_size,
+                            __FILE__, __LINE__ );
+                }
+
                 for ( pi = start_j; pi < end_j; ++pi )
                 {
                     pbond_ij = &bond_list[pi];
@@ -149,12 +158,13 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
                         pbond_ij = &bond_list[pi];
                         i = pbond_ij->nbr;
 
-                        if ( i != k )
+                        if ( i != k
+                                && system->reax_param.hbp[type_i][type_j][type_k].is_valid == TRUE )
                         {
                             bo_ij = &pbond_ij->bo_data;
                             type_i = system->atoms[i].type;
                             r_ij = pbond_ij->d;
-                            hbp = &system->reax_param.hbp[ type_i ][ type_j ][ type_k ];
+                            hbp = &system->reax_param.hbp[type_i][type_j][type_k];
 #if defined(_OPENMP)
                             f_i = &workspace->f_local[tid * system->N + i];
 #else
@@ -320,6 +330,11 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
                 }
             }
         }
+
+        if ( hblist != NULL )
+        {
+            sfree( hblist, __FILE__, __LINE__ );
+        }
     }
 
     data->E_HB += e_hb_total;
diff --git a/sPuReMD/src/reax_types.h b/sPuReMD/src/reax_types.h
index 61da0a8186b33c02aefb441aa7b19ad4a8ec1917..09b61506a5eebba76620e9eb5ad5a917584fd26b 100644
--- a/sPuReMD/src/reax_types.h
+++ b/sPuReMD/src/reax_types.h
@@ -701,7 +701,7 @@ struct three_body_header
 /* hydrogen-bond parameters */
 struct hbond_parameters
 {
-    /**/
+    /* ideal H-bond distance */
     real r0_hb;
     /**/
     real p_hb1;
@@ -709,6 +709,8 @@ struct hbond_parameters
     real p_hb2;
     /**/
     real p_hb3;
+    /* TRUE if parameters are set for this triplet of atom types, FALSE otherwise */
+    int is_valid;
 };