2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2010, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
59 typedef struct gmx_ga2la {
65 int start_space_search;
68 /* Clear all the entries in the ga2la list */
69 static void ga2la_clear(gmx_ga2la_t ga2la)
75 for (i = 0; i < ga2la->nalloc; i++)
77 ga2la->laa[i].cell = -1;
82 for (i = 0; i < ga2la->nalloc; i++)
84 ga2la->lal[i].ga = -1;
85 ga2la->lal[i].next = -1;
87 ga2la->start_space_search = ga2la->mod;
91 static gmx_ga2la_t ga2la_init(int nat_tot, int nat_loc)
97 /* There are two methods implemented for finding the local atom number
98 * belonging to a global atom number:
99 * 1) a simple, direct array
100 * 2) a list of linked lists indexed with the global number modulo mod.
101 * Memory requirements:
103 * 2) nat_loc*(2+1-2(1-e^-1/2))*4 ints
104 * where nat_loc is the number of atoms in the home + communicated zones.
105 * Method 1 is faster for low parallelization, 2 for high parallelization.
106 * We switch to method 2 when it uses less than half the memory method 1.
108 ga2la->bAll = (nat_tot < 1024 || 9*nat_loc >= nat_tot);
111 ga2la->nalloc = nat_tot;
112 snew(ga2la->laa, ga2la->nalloc);
116 /* Make the direct list twice as long as the number of local atoms.
117 * The fraction of entries in the list with:
118 * 0 size lists: e^-1/f
119 * >=1 size lists: 1 - e^-1/f
120 * where f is: the direct list length / #local atoms
121 * The fraction of atoms not in the direct list is: 1-f(1-e^-1/f).
123 ga2la->mod = 2*nat_loc;
124 ga2la->nalloc = over_alloc_dd(ga2la->mod);
125 snew(ga2la->lal, ga2la->nalloc);
133 /* Set the ga2la entry for global atom a_gl to local atom a_loc and cell. */
134 static void ga2la_set(gmx_ga2la_t ga2la, int a_gl, int a_loc, int cell)
136 int ind, ind_prev, i;
140 ga2la->laa[a_gl].la = a_loc;
141 ga2la->laa[a_gl].cell = cell;
146 ind = a_gl % ga2la->mod;
148 if (ga2la->lal[ind].ga >= 0)
150 /* Search the last entry in the linked list for this index */
152 while (ga2la->lal[ind_prev].next >= 0)
154 ind_prev = ga2la->lal[ind_prev].next;
156 /* Search for space in the array */
157 ind = ga2la->start_space_search;
158 while (ind < ga2la->nalloc && ga2la->lal[ind].ga >= 0)
162 /* If we are at the end of the list we need to increase the size */
163 if (ind == ga2la->nalloc)
165 ga2la->nalloc = over_alloc_dd(ind+1);
166 srenew(ga2la->lal, ga2la->nalloc);
167 for (i = ind; i < ga2la->nalloc; i++)
169 ga2la->lal[i].ga = -1;
170 ga2la->lal[i].next = -1;
173 ga2la->lal[ind_prev].next = ind;
175 ga2la->start_space_search = ind + 1;
177 ga2la->lal[ind].ga = a_gl;
178 ga2la->lal[ind].la = a_loc;
179 ga2la->lal[ind].cell = cell;
182 /* Delete the ga2la entry for global atom a_gl */
183 static void ga2la_del(gmx_ga2la_t ga2la, int a_gl)
189 ga2la->laa[a_gl].cell = -1;
195 ind = a_gl % ga2la->mod;
198 if (ga2la->lal[ind].ga == a_gl)
202 ga2la->lal[ind_prev].next = ga2la->lal[ind].next;
204 /* This index is a linked entry, so we free an entry.
205 * Check if we are creating the first empty space.
207 if (ind < ga2la->start_space_search)
209 ga2la->start_space_search = ind;
212 ga2la->lal[ind].ga = -1;
213 ga2la->lal[ind].cell = -1;
214 ga2la->lal[ind].next = -1;
219 ind = ga2la->lal[ind].next;
226 /* Change the local atom for present ga2la entry for global atom a_gl */
227 static void ga2la_change_la(gmx_ga2la_t ga2la, int a_gl, int a_loc)
233 ga2la->laa[a_gl].la = a_loc;
238 ind = a_gl % ga2la->mod;
241 if (ga2la->lal[ind].ga == a_gl)
243 ga2la->lal[ind].la = a_loc;
247 ind = ga2la->lal[ind].next;
254 /* Returns if the global atom a_gl available locally.
255 * Sets the local atom and cell,
256 * cell can be larger than the number of zones,
257 * in which case it indicates that it is more than one cell away
258 * in zone cell - #zones.
260 static gmx_bool ga2la_get(const gmx_ga2la_t ga2la, int a_gl, int *a_loc, int *cell)
266 *a_loc = ga2la->laa[a_gl].la;
267 *cell = ga2la->laa[a_gl].cell;
269 return (ga2la->laa[a_gl].cell >= 0);
272 ind = a_gl % ga2la->mod;
275 if (ga2la->lal[ind].ga == a_gl)
277 *a_loc = ga2la->lal[ind].la;
278 *cell = ga2la->lal[ind].cell;
282 ind = ga2la->lal[ind].next;
289 /* Returns if the global atom a_gl is a home atom.
290 * Sets the local atom.
292 static gmx_bool ga2la_get_home(const gmx_ga2la_t ga2la, int a_gl, int *a_loc)
298 *a_loc = ga2la->laa[a_gl].la;
300 return (ga2la->laa[a_gl].cell == 0);
303 ind = a_gl % ga2la->mod;
306 if (ga2la->lal[ind].ga == a_gl)
308 if (ga2la->lal[ind].cell == 0)
310 *a_loc = ga2la->lal[ind].la;
319 ind = ga2la->lal[ind].next;
326 /* Returns if the global atom a_gl is a home atom.
328 static gmx_bool ga2la_is_home(const gmx_ga2la_t ga2la, int a_gl)
334 return (ga2la->laa[a_gl].cell == 0);
337 ind = a_gl % ga2la->mod;
340 if (ga2la->lal[ind].ga == a_gl)
342 return (ga2la->lal[ind].cell == 0);
344 ind = ga2la->lal[ind].next;
355 #endif /* _gmx_ga2la_h */