2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2014, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
40 #include "gmx_fatal.h"
46 /* This include file implements the simplest hash table possible.
47 * It is limited to integer keys and integer values.
48 * The purpose is highest efficiency and lowest memory usage possible.
50 * The type definition is placed in types/commrec.h, as it is used there:
51 * typedef struct gmx_hash *gmx_hash_t
60 typedef struct gmx_hash {
67 int start_space_search;
70 /* Clear all the entries in the hash table */
71 static void gmx_hash_clear(gmx_hash_t hash)
75 for (i = 0; i < hash->nalloc; i++)
77 hash->hash[i].key = -1;
78 hash->hash[i].next = -1;
80 hash->start_space_search = hash->mod;
85 static void gmx_hash_realloc(gmx_hash_t hash, int nkey_used_estimate)
87 /* Memory requirements:
88 * nkey_used_est*(2+1-2(1-e^-1/2))*3 ints
89 * where nkey_used_est is the local number of keys used.
91 * Make the direct list twice as long as the number of local keys.
92 * The fraction of entries in the list with:
94 * >=1 size lists: 1 - e^-f
95 * where f is: the #keys / mod
96 * The fraction of keys not in the direct list is: 1-1/f(1-e^-f).
97 * The optimal table size is roughly double the number of keys.
99 /* Make the hash table a power of 2 and at least double the number of keys */
101 while (2*nkey_used_estimate > hash->mod)
105 hash->mask = hash->mod - 1;
106 hash->nalloc = over_alloc_dd(hash->mod);
107 srenew(hash->hash, hash->nalloc);
111 fprintf(debug, "Hash table mod %d nalloc %d\n", hash->mod, hash->nalloc);
115 /* Clear all the entries in the hash table.
116 * With the current number of keys check if the table size is still good,
117 * if not optimize it with the currenr number of keys.
119 static void gmx_hash_clear_and_optimize(gmx_hash_t hash)
121 /* Resize the hash table when the occupation is < 1/4 or > 2/3 */
122 if (hash->nkey > 0 &&
123 (4*hash->nkey < hash->mod || 3*hash->nkey > 2*hash->mod))
127 fprintf(debug, "Hash table size %d #key %d: resizing\n",
128 hash->mod, hash->nkey);
130 gmx_hash_realloc(hash, hash->nkey);
133 gmx_hash_clear(hash);
136 static gmx_hash_t gmx_hash_init(int nkey_used_estimate)
143 gmx_hash_realloc(hash, nkey_used_estimate);
145 gmx_hash_clear(hash);
150 /* Set the hash entry for global atom a_gl to local atom a_loc and cell. */
151 static void gmx_hash_set(gmx_hash_t hash, int key, int value)
153 int ind, ind_prev, i;
155 ind = key & hash->mask;
157 if (hash->hash[ind].key >= 0)
159 /* Search the last entry in the linked list for this index */
161 while (hash->hash[ind_prev].next >= 0)
163 ind_prev = hash->hash[ind_prev].next;
165 /* Search for space in the array */
166 ind = hash->start_space_search;
167 while (ind < hash->nalloc && hash->hash[ind].key >= 0)
171 /* If we are at the end of the list we need to increase the size */
172 if (ind == hash->nalloc)
174 hash->nalloc = over_alloc_dd(ind+1);
175 srenew(hash->hash, hash->nalloc);
176 for (i = ind; i < hash->nalloc; i++)
178 hash->hash[i].key = -1;
179 hash->hash[i].next = -1;
182 hash->hash[ind_prev].next = ind;
184 hash->start_space_search = ind + 1;
186 hash->hash[ind].key = key;
187 hash->hash[ind].val = value;
192 /* Delete the hash entry for key */
193 static void gmx_hash_del(gmx_hash_t hash, int key)
198 ind = key & hash->mask;
201 if (hash->hash[ind].key == key)
205 hash->hash[ind_prev].next = hash->hash[ind].next;
207 /* This index is a linked entry, so we free an entry.
208 * Check if we are creating the first empty space.
210 if (ind < hash->start_space_search)
212 hash->start_space_search = ind;
215 hash->hash[ind].key = -1;
216 hash->hash[ind].val = -1;
217 hash->hash[ind].next = -1;
224 ind = hash->hash[ind].next;
231 /* Change the value for present hash entry for key */
232 static void gmx_hash_change_value(gmx_hash_t hash, int key, int value)
236 ind = key & hash->mask;
239 if (hash->hash[ind].key == key)
241 hash->hash[ind].val = value;
245 ind = hash->hash[ind].next;
252 /* Change the hash value if already set, otherwise set the hash value */
253 static void gmx_hash_change_or_set(gmx_hash_t hash, int key, int value)
257 ind = key & hash->mask;
260 if (hash->hash[ind].key == key)
262 hash->hash[ind].val = value;
266 ind = hash->hash[ind].next;
270 gmx_hash_set(hash, key, value);
275 /* Returns if the key is present, if the key is present *value is set */
276 static gmx_bool gmx_hash_get(const gmx_hash_t hash, int key, int *value)
280 ind = key & hash->mask;
283 if (hash->hash[ind].key == key)
285 *value = hash->hash[ind].val;
289 ind = hash->hash[ind].next;
296 /* Returns the value or -1 if the key is not present */
297 static int gmx_hash_get_minone(const gmx_hash_t hash, int key)
301 ind = key & hash->mask;
304 if (hash->hash[ind].key == key)
306 return hash->hash[ind].val;
308 ind = hash->hash[ind].next;
319 #endif /* _gmx_hash_h */