1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2012, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Gallium Rubidium Oxygen Manganese Argon Carbon Silicon
36 #ifndef _nbnxn_pairlist_h
37 #define _nbnxn_pairlist_h
43 /* A buffer data structure of 64 bytes
44 * to be placed at the beginning and end of structs
45 * to avoid cache invalidation of the real contents
46 * of the struct by writes to neighboring memory.
50 } gmx_cache_protect_t
;
52 /* Abstract type for pair searching data */
53 typedef struct nbnxn_search
* nbnxn_search_t
;
55 /* Function that should return a pointer *ptr to memory
57 * Error handling should be done within this function.
59 typedef void gmx_nbat_alloc_t(void **ptr
,size_t nbytes
);
61 /* Function that should free the memory pointed to by *ptr.
62 * NULL should not be passed to this function.
64 typedef void gmx_nbat_free_t(void *ptr
);
67 int cj
; /* The j-cluster */
68 unsigned excl
; /* The exclusion (interaction) bits */
71 #define NBNXN_CI_SHIFT 127
72 #define NBNXN_CI_DO_LJ(subc) (1<<(7+3*(subc)))
73 #define NBNXN_CI_HALF_LJ(subc) (1<<(8+3*(subc)))
74 #define NBNXN_CI_DO_COUL(subc) (1<<(9+3*(subc)))
76 /* Simple pair-list i-unit */
78 int ci
; /* i-cluster */
79 int shift
; /* Shift vector index plus possible flags */
80 int cj_ind_start
; /* Start index into cj */
81 int cj_ind_end
; /* End index into cj */
84 /* Grouped pair-list i-unit */
86 int sci
; /* i-super-cluster */
87 int shift
; /* Shift vector index plus possible flags */
88 int cj4_ind_start
; /* Start index into cj4 */
89 int cj4_ind_end
; /* End index into cj4 */
93 unsigned imask
; /* The i-cluster interactions mask for 1 warp */
94 int excl_ind
; /* Index into the exclusion array for 1 warp */
98 int cj
[4]; /* The 4 j-clusters */
99 nbnxn_im_ei_t imei
[2]; /* The i-cluster mask data for 2 warps */
103 unsigned pair
[32]; /* Exclusion bits for one warp, *
104 * each unsigned has bit for 4*8 i clusters */
108 gmx_cache_protect_t cp0
;
110 gmx_nbat_alloc_t
*alloc
;
111 gmx_nbat_free_t
*free
;
113 gmx_bool bSimple
; /* Simple list has na_sc=na_s and uses cj *
114 * Complex list uses cj4 */
116 int na_ci
; /* The number of atoms per i-cluster */
117 int na_cj
; /* The number of atoms per j-cluster */
118 int na_sc
; /* The number of atoms per super cluster */
119 real rlist
; /* The radius for constructing the list */
120 int nci
; /* The number of i-clusters in the list */
121 nbnxn_ci_t
*ci
; /* The i-cluster list, size nci */
122 int ci_nalloc
; /* The allocation size of ci */
123 int nsci
; /* The number of i-super-clusters in the list */
124 nbnxn_sci_t
*sci
; /* The i-super-cluster list */
125 int sci_nalloc
; /* The allocation size of sci */
127 int ncj
; /* The number of j-clusters in the list */
128 nbnxn_cj_t
*cj
; /* The j-cluster list, size ncj */
129 int cj_nalloc
; /* The allocation size of cj */
131 int ncj4
; /* The total number of 4*j clusters */
132 nbnxn_cj4_t
*cj4
; /* The 4*j cluster list, size ncj4 */
133 int cj4_nalloc
; /* The allocation size of cj4 */
134 int nexcl
; /* The count for excl */
135 nbnxn_excl_t
*excl
; /* Atom interaction bits (non-exclusions) */
136 int excl_nalloc
; /* The allocation size for excl */
137 int nci_tot
; /* The total number of i clusters */
139 struct nbnxn_list_work
*work
;
141 gmx_cache_protect_t cp1
;
145 int nnbl
; /* number of lists */
146 nbnxn_pairlist_t
**nbl
; /* lists */
147 gmx_bool bCombined
; /* TRUE if lists get combined into one (the 1st) */
148 gmx_bool bSimple
; /* TRUE if the list of of type "simple"
149 (na_sc=na_s, no super-clusters used) */
150 int natpair_ljq
; /* Total number of atom pairs for LJ+Q kernel */
151 int natpair_lj
; /* Total number of atom pairs for LJ kernel */
152 int natpair_q
; /* Total number of atom pairs for Q kernel */
153 } nbnxn_pairlist_set_t
;
155 enum { nbatXYZ
, nbatXYZQ
, nbatX4
, nbatX8
};
158 real
*f
; /* f, size natoms*fstride */
159 real
*fshift
; /* Shift force array, size SHIFTS*DIM */
160 int nV
; /* The size of *Vvdw and *Vc */
161 real
*Vvdw
; /* Temporary Van der Waals group energy storage */
162 real
*Vc
; /* Temporary Coulomb group energy storage */
163 int nVS
; /* The size of *VSvdw and *VSc */
164 real
*VSvdw
; /* Temporary SIMD Van der Waals group energy storage */
165 real
*VSc
; /* Temporary SIMD Coulomb group energy storage */
166 } nbnxn_atomdata_output_t
;
168 /* LJ combination rules: geometric, Lorentz-Berthelot, none */
169 enum { ljcrGEOM
, ljcrLB
, ljcrNONE
, ljcrNR
};
172 gmx_nbat_alloc_t
*alloc
;
173 gmx_nbat_free_t
*free
;
174 int ntype
; /* The number of different atom types */
175 real
*nbfp
; /* Lennard-Jones 6*C6 and 12*C12 params, size ntype^2*2 */
176 int comb_rule
; /* Combination rule, see enum above */
177 real
*nbfp_comb
; /* LJ parameter per atom type, size ntype*2 */
178 real
*nbfp_s4
; /* As nbfp, but with stride 4, size ntype^2*4 */
179 int natoms
; /* Number of atoms */
180 int natoms_local
; /* Number of local atoms */
181 int *type
; /* Atom types */
182 real
*lj_comb
; /* LJ parameters per atom for combining for pairs */
183 int XFormat
; /* The format of x (and q), enum */
184 int FFormat
; /* The format of f, enum */
185 real
*q
; /* Charges, can be NULL if incorporated in x */
186 int na_c
; /* The number of atoms per cluster */
187 int nenergrp
; /* The number of energy groups */
188 int neg_2log
; /* Log2 of nenergrp */
189 int *energrp
; /* The energy groups per cluster, can be NULL */
190 gmx_bool bDynamicBox
; /* Do we need to update shift_vec every step? */
191 rvec
*shift_vec
; /* Shift vectors, copied from t_forcerec */
192 int xstride
; /* stride for a coordinate in x (usually 3 or 4) */
193 int fstride
; /* stride for a coordinate in f (usually 3 or 4) */
194 real
*x
; /* x and possibly q, size natoms*xstride */
195 int nout
; /* The number of force arrays */
196 nbnxn_atomdata_output_t
*out
; /* Output data structures */
197 int nalloc
; /* Allocation size of all arrays (for x/f *x/fstride) */