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!
42 #include "gmx_omp_nthreads.h"
43 #include "nbnxn_kernel_ref.h"
44 #include "../nbnxn_consts.h"
45 #include "nbnxn_kernel_common.h"
47 /*! \brief Typedefs for declaring lookup tables of kernel functions.
49 typedef void (*p_nbk_func_ener)(const nbnxn_pairlist_t *nbl,
50 const nbnxn_atomdata_t *nbat,
51 const interaction_const_t *ic,
58 typedef void (*p_nbk_func_noener)(const nbnxn_pairlist_t *nbl,
59 const nbnxn_atomdata_t *nbat,
60 const interaction_const_t *ic,
65 /* Analytical reaction-field kernels */
68 /* Include the force+energy kernels */
70 #include "nbnxn_kernel_ref_outer.h"
73 /* Include the force+energygroups kernels */
76 #include "nbnxn_kernel_ref_outer.h"
80 /* Include the force only kernels */
81 #include "nbnxn_kernel_ref_outer.h"
86 /* Tabulated exclusion interaction electrostatics kernels */
89 /* Include the force+energy kernels */
91 #include "nbnxn_kernel_ref_outer.h"
94 /* Include the force+energygroups kernels */
97 #include "nbnxn_kernel_ref_outer.h"
101 /* Include the force only kernels */
102 #include "nbnxn_kernel_ref_outer.h"
104 /* Twin-range cut-off kernels */
105 #define VDW_CUTOFF_CHECK
107 /* Include the force+energy kernels */
108 #define CALC_ENERGIES
109 #include "nbnxn_kernel_ref_outer.h"
112 /* Include the force+energygroups kernels */
113 #define CALC_ENERGIES
114 #define ENERGY_GROUPS
115 #include "nbnxn_kernel_ref_outer.h"
119 /* Include the force only kernels */
120 #include "nbnxn_kernel_ref_outer.h"
122 #undef VDW_CUTOFF_CHECK
128 coultRF, coultTAB, coultTAB_TWIN, coultNR
131 p_nbk_func_ener p_nbk_c_ener[coultNR] =
133 nbnxn_kernel_ref_rf_ener,
134 nbnxn_kernel_ref_tab_ener,
135 nbnxn_kernel_ref_tab_twin_ener
138 p_nbk_func_ener p_nbk_c_energrp[coultNR] =
140 nbnxn_kernel_ref_rf_energrp,
141 nbnxn_kernel_ref_tab_energrp,
142 nbnxn_kernel_ref_tab_twin_energrp
145 p_nbk_func_noener p_nbk_c_noener[coultNR] =
147 nbnxn_kernel_ref_rf_noener,
148 nbnxn_kernel_ref_tab_noener,
149 nbnxn_kernel_ref_tab_twin_noener
153 nbnxn_kernel_ref(const nbnxn_pairlist_set_t *nbl_list,
154 const nbnxn_atomdata_t *nbat,
155 const interaction_const_t *ic,
164 nbnxn_pairlist_t **nbl;
168 nnbl = nbl_list->nnbl;
171 if (EEL_RF(ic->eeltype) || ic->eeltype == eelCUT)
177 if (ic->rcoulomb == ic->rvdw)
183 coult = coultTAB_TWIN;
187 #pragma omp parallel for schedule(static) num_threads(gmx_omp_nthreads_get(emntNonbonded))
188 for (nb = 0; nb < nnbl; nb++)
190 nbnxn_atomdata_output_t *out;
193 out = &nbat->out[nb];
195 if (clearF == enbvClearFYes)
197 clear_f(nbat, nb, out->f);
200 if ((force_flags & GMX_FORCE_VIRIAL) && nnbl == 1)
206 fshift_p = out->fshift;
208 if (clearF == enbvClearFYes)
210 clear_fshift(fshift_p);
214 if (!(force_flags & GMX_FORCE_ENERGY))
216 /* Don't calculate energies */
217 p_nbk_c_noener[coult](nbl[nb], nbat,
223 else if (out->nV == 1)
225 /* No energy groups */
229 p_nbk_c_ener[coult](nbl[nb], nbat,
239 /* Calculate energy group contributions */
242 for (i = 0; i < out->nV; i++)
246 for (i = 0; i < out->nV; i++)
251 p_nbk_c_energrp[coult](nbl[nb], nbat,
261 if (force_flags & GMX_FORCE_ENERGY)
263 reduce_energies_over_lists(nbat, nnbl, Vvdw, Vc);