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,2014, 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.
40 #include "types/commrec.h"
41 #include "gromacs/utility/basedefinitions.h"
42 #include "gromacs/utility/smalloc.h"
60 typedef struct gmx_ga2la {
66 int start_space_search;
69 /* Clear all the entries in the ga2la list */
70 static void ga2la_clear(gmx_ga2la_t ga2la)
76 for (i = 0; i < ga2la->nalloc; i++)
78 ga2la->laa[i].cell = -1;
83 for (i = 0; i < ga2la->nalloc; i++)
85 ga2la->lal[i].ga = -1;
86 ga2la->lal[i].next = -1;
88 ga2la->start_space_search = ga2la->mod;
92 static gmx_ga2la_t ga2la_init(int nat_tot, int nat_loc)
98 /* There are two methods implemented for finding the local atom number
99 * belonging to a global atom number:
100 * 1) a simple, direct array
101 * 2) a list of linked lists indexed with the global number modulo mod.
102 * Memory requirements:
104 * 2) nat_loc*(2+1-2(1-e^-1/2))*4 ints
105 * where nat_loc is the number of atoms in the home + communicated zones.
106 * Method 1 is faster for low parallelization, 2 for high parallelization.
107 * We switch to method 2 when it uses less than half the memory method 1.
109 ga2la->bAll = (nat_tot < 1024 || 9*nat_loc >= nat_tot);
112 ga2la->nalloc = nat_tot;
113 snew(ga2la->laa, ga2la->nalloc);
117 /* Make the direct list twice as long as the number of local atoms.
118 * The fraction of entries in the list with:
119 * 0 size lists: e^-1/f
120 * >=1 size lists: 1 - e^-1/f
121 * where f is: the direct list length / #local atoms
122 * The fraction of atoms not in the direct list is: 1-f(1-e^-1/f).
124 ga2la->mod = 2*nat_loc;
125 ga2la->nalloc = over_alloc_dd(ga2la->mod);
126 snew(ga2la->lal, ga2la->nalloc);
134 /* Set the ga2la entry for global atom a_gl to local atom a_loc and cell. */
135 static void ga2la_set(gmx_ga2la_t ga2la, int a_gl, int a_loc, int cell)
137 int ind, ind_prev, i;
141 ga2la->laa[a_gl].la = a_loc;
142 ga2la->laa[a_gl].cell = cell;
147 ind = a_gl % ga2la->mod;
149 if (ga2la->lal[ind].ga >= 0)
151 /* Search the last entry in the linked list for this index */
153 while (ga2la->lal[ind_prev].next >= 0)
155 ind_prev = ga2la->lal[ind_prev].next;
157 /* Search for space in the array */
158 ind = ga2la->start_space_search;
159 while (ind < ga2la->nalloc && ga2la->lal[ind].ga >= 0)
163 /* If we are at the end of the list we need to increase the size */
164 if (ind == ga2la->nalloc)
166 ga2la->nalloc = over_alloc_dd(ind+1);
167 srenew(ga2la->lal, ga2la->nalloc);
168 for (i = ind; i < ga2la->nalloc; i++)
170 ga2la->lal[i].ga = -1;
171 ga2la->lal[i].next = -1;
174 ga2la->lal[ind_prev].next = ind;
176 ga2la->start_space_search = ind + 1;
178 ga2la->lal[ind].ga = a_gl;
179 ga2la->lal[ind].la = a_loc;
180 ga2la->lal[ind].cell = cell;
183 /* Delete the ga2la entry for global atom a_gl */
184 static void ga2la_del(gmx_ga2la_t ga2la, int a_gl)
190 ga2la->laa[a_gl].cell = -1;
196 ind = a_gl % ga2la->mod;
199 if (ga2la->lal[ind].ga == a_gl)
203 ga2la->lal[ind_prev].next = ga2la->lal[ind].next;
205 /* This index is a linked entry, so we free an entry.
206 * Check if we are creating the first empty space.
208 if (ind < ga2la->start_space_search)
210 ga2la->start_space_search = ind;
213 ga2la->lal[ind].ga = -1;
214 ga2la->lal[ind].cell = -1;
215 ga2la->lal[ind].next = -1;
220 ind = ga2la->lal[ind].next;
227 /* Change the local atom for present ga2la entry for global atom a_gl */
228 static void ga2la_change_la(gmx_ga2la_t ga2la, int a_gl, int a_loc)
234 ga2la->laa[a_gl].la = a_loc;
239 ind = a_gl % ga2la->mod;
242 if (ga2la->lal[ind].ga == a_gl)
244 ga2la->lal[ind].la = a_loc;
248 ind = ga2la->lal[ind].next;
255 /* Returns if the global atom a_gl available locally.
256 * Sets the local atom and cell,
257 * cell can be larger than the number of zones,
258 * in which case it indicates that it is more than one cell away
259 * in zone cell - #zones.
261 static gmx_bool ga2la_get(const gmx_ga2la_t ga2la, int a_gl, int *a_loc, int *cell)
267 *a_loc = ga2la->laa[a_gl].la;
268 *cell = ga2la->laa[a_gl].cell;
270 return (ga2la->laa[a_gl].cell >= 0);
273 ind = a_gl % ga2la->mod;
276 if (ga2la->lal[ind].ga == a_gl)
278 *a_loc = ga2la->lal[ind].la;
279 *cell = ga2la->lal[ind].cell;
283 ind = ga2la->lal[ind].next;
290 /* Returns if the global atom a_gl is a home atom.
291 * Sets the local atom.
293 static gmx_bool ga2la_get_home(const gmx_ga2la_t ga2la, int a_gl, int *a_loc)
299 *a_loc = ga2la->laa[a_gl].la;
301 return (ga2la->laa[a_gl].cell == 0);
304 ind = a_gl % ga2la->mod;
307 if (ga2la->lal[ind].ga == a_gl)
309 if (ga2la->lal[ind].cell == 0)
311 *a_loc = ga2la->lal[ind].la;
320 ind = ga2la->lal[ind].next;
327 /* Returns if the global atom a_gl is a home atom.
329 static gmx_bool ga2la_is_home(const gmx_ga2la_t ga2la, int a_gl)
335 return (ga2la->laa[a_gl].cell == 0);
338 ind = a_gl % ga2la->mod;
341 if (ga2la->lal[ind].ga == a_gl)
343 return (ga2la->lal[ind].cell == 0);
345 ind = ga2la->lal[ind].next;
356 #endif /* _gmx_ga2la_h */