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 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
9 * Copyright (c) 2001-2009, The GROMACS Development Team
11 * Gromacs is a library for molecular simulation and trajectory analysis,
12 * written by Erik Lindahl, David van der Spoel, Berk Hess, and others - for
13 * a full list of developers and information, check out http://www.gromacs.org
15 * This program is free software; you can redistribute it and/or modify it under
16 * the terms of the GNU Lesser General Public License as published by the Free
17 * Software Foundation; either version 2 of the License, or (at your option) any
19 * As a special exception, you may use this file as part of a free software
20 * library without restriction. Specifically, if other files instantiate
21 * templates or use macros or inline functions from this file, or you compile
22 * this file and link it with other files to produce an executable, this
23 * file does not by itself cause the resulting executable to be covered by
24 * the GNU Lesser General Public License.
26 * In plain-speak: do not worry about classes/macros/templates either - only
27 * changes to the library have to be LGPL, not an application linking with it.
29 * To help fund GROMACS development, we humbly ask that you cite
30 * the papers people have written on it - you can find them on the website!
33 #define UNROLLI NBNXN_CPU_CLUSTER_I_SIZE
34 #define UNROLLJ NBNXN_CPU_CLUSTER_I_SIZE
36 /* We could use nbat->xstride and nbat->fstride, but macros might be faster */
39 /* Local i-atom buffer strides */
44 /* All functionality defines are set here, except for:
45 * CALC_ENERGIES, ENERGY_GROUPS which are defined before.
46 * CHECK_EXCLS, which is set just before including the inner loop contents.
49 /* We always calculate shift forces, because it's cheap anyhow */
50 #define CALC_SHIFTFORCES
53 #define NBK_FUNC_NAME(x,y) x##_rf_##y
56 #ifndef VDW_CUTOFF_CHECK
57 #define NBK_FUNC_NAME(x,y) x##_tab_##y
59 #define NBK_FUNC_NAME(x,y) x##_tab_twin_##y
65 NBK_FUNC_NAME(nbnxn_kernel_ref,noener)
68 NBK_FUNC_NAME(nbnxn_kernel_ref,ener)
70 NBK_FUNC_NAME(nbnxn_kernel_ref,energrp)
74 (const nbnxn_pairlist_t *nbl,
75 const nbnxn_atomdata_t *nbat,
76 const interaction_const_t *ic,
79 #ifdef CALC_SHIFTFORCES
90 const nbnxn_ci_t *nbln;
91 const nbnxn_cj_t *l_cj;
98 #ifdef VDW_CUTOFF_CHECK
106 gmx_bool do_LJ,half_LJ,do_coul;
107 int cjind0,cjind1,cjind;
110 real xi[UNROLLI*XI_STRIDE];
111 real fi[UNROLLI*FI_STRIDE];
115 #ifndef ENERGY_GROUPS
120 int egp_sh_i[UNROLLI];
137 const real *tab_coul_FDV0;
139 const real *tab_coul_F;
140 const real *tab_coul_V;
151 sh_invrc6 = ic->sh_invrc6;
162 tabscale = ic->tabq_scale;
164 halfsp = 0.5/ic->tabq_scale;
168 tab_coul_FDV0 = ic->tabq_coul_FDV0;
170 tab_coul_F = ic->tabq_coul_F;
171 tab_coul_V = ic->tabq_coul_V;
176 egp_mask = (1<<nbat->neg_2log) - 1;
180 rcut2 = ic->rcoulomb*ic->rcoulomb;
181 #ifdef VDW_CUTOFF_CHECK
182 rvdw2 = ic->rvdw*ic->rvdw;
185 ntype2 = nbat->ntype*2;
190 shiftvec = shift_vec[0];
196 for(n=0; n<nbl->nci; n++)
202 ish = (nbln->shift & NBNXN_CI_SHIFT);
203 /* x, f and fshift are assumed to be stored with stride 3 */
205 cjind0 = nbln->cj_ind_start;
206 cjind1 = nbln->cj_ind_end;
207 /* Currently only works super-cells equal to sub-cells */
209 ci_sh = (ish == CENTRAL ? ci : -1);
211 /* We have 5 LJ/C combinations, but use only three inner loops,
212 * as the other combinations are unlikely and/or not much faster:
213 * inner half-LJ + C for half-LJ + C / no-LJ + C
214 * inner LJ + C for full-LJ + C
215 * inner LJ for full-LJ + no-C / half-LJ + no-C
217 do_LJ = (nbln->shift & NBNXN_CI_DO_LJ(0));
218 do_coul = (nbln->shift & NBNXN_CI_DO_COUL(0));
219 half_LJ = ((nbln->shift & NBNXN_CI_HALF_LJ(0)) || !do_LJ) && do_coul;
222 #ifndef ENERGY_GROUPS
226 for(i=0; i<UNROLLI; i++)
228 egp_sh_i[i] = ((nbat->energrp[ci]>>(i*nbat->neg_2log)) & egp_mask)*nbat->nenergrp;
233 for(i=0; i<UNROLLI; i++)
237 xi[i*XI_STRIDE+d] = x[(ci*UNROLLI+i)*X_STRIDE+d] + shiftvec[ishf+d];
238 fi[i*FI_STRIDE+d] = 0;
248 Vc_sub_self = 0.5*c_rf;
252 Vc_sub_self = 0.5*tab_coul_V[0];
254 Vc_sub_self = 0.5*tab_coul_FDV0[2];
259 for(i=0; i<UNROLLI; i++)
261 qi[i] = facel*q[ci*UNROLLI+i];
264 if (l_cj[nbln->cj_ind_start].cj == ci_sh)
267 Vc[egp_sh_i[i]+((nbat->energrp[ci]>>(i*nbat->neg_2log)) & egp_mask)]
271 -= qi[i]*q[ci*UNROLLI+i]*Vc_sub_self;
278 while (cjind < cjind1 && nbl->cj[cjind].excl != 0xffff)
285 #include "nbnxn_kernel_ref_inner.h"
289 /* cppcheck-suppress duplicateBranch */
293 #include "nbnxn_kernel_ref_inner.h"
298 #include "nbnxn_kernel_ref_inner.h"
304 for(; (cjind<cjind1); cjind++)
310 #include "nbnxn_kernel_ref_inner.h"
314 /* cppcheck-suppress duplicateBranch */
318 #include "nbnxn_kernel_ref_inner.h"
323 #include "nbnxn_kernel_ref_inner.h"
326 ninner += cjind1 - cjind0;
328 /* Add accumulated i-forces to the force array */
329 for(i=0; i<UNROLLI; i++)
333 f[(ci*UNROLLI+i)*F_STRIDE+d] += fi[i*FI_STRIDE+d];
336 #ifdef CALC_SHIFTFORCES
339 /* Add i forces to shifted force list */
340 for(i=0; i<UNROLLI; i++)
344 fshift[ishf+d] += fi[i*FI_STRIDE+d];
351 #ifndef ENERGY_GROUPS
359 printf("atom pairs %d\n",npair);
363 #undef CALC_SHIFTFORCES