-/* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
+/*
+ * This file is part of the GROMACS molecular simulation package.
*
+ * Copyright (c) 2012,2013,2014, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
*
- * This source code is part of
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
*
- * G R O M A C S
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
*
- * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2009, The GROMACS Development Team
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
*
- * Gromacs is a library for molecular simulation and trajectory analysis,
- * written by Erik Lindahl, David van der Spoel, Berk Hess, and others - for
- * a full list of developers and information, check out http://www.gromacs.org
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
*
- * This program is free software; you can redistribute it and/or modify it under
- * the terms of the GNU Lesser General Public License as published by the Free
- * Software Foundation; either version 2 of the License, or (at your option) any
- * later version.
- * As a special exception, you may use this file as part of a free software
- * library without restriction. Specifically, if other files instantiate
- * templates or use macros or inline functions from this file, or you compile
- * this file and link it with other files to produce an executable, this
- * file does not by itself cause the resulting executable to be covered by
- * the GNU Lesser General Public License.
- *
- * In plain-speak: do not worry about classes/macros/templates either - only
- * changes to the library have to be LGPL, not an application linking with it.
- *
- * To help fund GROMACS development, we humbly ask that you cite
- * the papers people have written on it - you can find them on the website!
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
*/
#define UNROLLI NBNXN_CPU_CLUSTER_I_SIZE
#define CALC_SHIFTFORCES
#ifdef CALC_COUL_RF
-#define NBK_FUNC_NAME(x,y) x##_rf_##y
+#define NBK_FUNC_NAME2(ljt, feg) nbnxn_kernel ## _ElecRF ## ljt ## feg ## _ref
#endif
#ifdef CALC_COUL_TAB
#ifndef VDW_CUTOFF_CHECK
-#define NBK_FUNC_NAME(x,y) x##_tab_##y
+#define NBK_FUNC_NAME2(ljt, feg) nbnxn_kernel ## _ElecQSTab ## ljt ## feg ## _ref
+#else
+#define NBK_FUNC_NAME2(ljt, feg) nbnxn_kernel ## _ElecQSTabTwinCut ## ljt ## feg ## _ref
+#endif
+#endif
+
+#if defined LJ_CUT && !defined LJ_EWALD
+#define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJ, feg)
+#elif defined LJ_FORCE_SWITCH
+#define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJFsw, feg)
+#elif defined LJ_POT_SWITCH
+#define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJPsw, feg)
+#elif defined LJ_EWALD
+#ifdef LJ_EWALD_COMB_GEOM
+#define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJEwCombGeom, feg)
#else
-#define NBK_FUNC_NAME(x,y) x##_tab_twin_##y
+#define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJEwCombLB, feg)
#endif
+#else
+#error "No VdW type defined"
#endif
static void
#ifndef CALC_ENERGIES
-NBK_FUNC_NAME(nbnxn_kernel_ref,noener)
+NBK_FUNC_NAME(_F)
#else
#ifndef ENERGY_GROUPS
-NBK_FUNC_NAME(nbnxn_kernel_ref,ener)
+NBK_FUNC_NAME(_VF)
#else
-NBK_FUNC_NAME(nbnxn_kernel_ref,energrp)
+NBK_FUNC_NAME(_VgrpF)
#endif
#endif
#undef NBK_FUNC_NAME
- (const nbnxn_pairlist_t *nbl,
- const nbnxn_atomdata_t *nbat,
- const interaction_const_t *ic,
- rvec *shift_vec,
- real *f
+#undef NBK_FUNC_NAME2
+(const nbnxn_pairlist_t *nbl,
+ const nbnxn_atomdata_t *nbat,
+ const interaction_const_t *ic,
+ rvec *shift_vec,
+ real *f
#ifdef CALC_SHIFTFORCES
- ,
- real *fshift
+ ,
+ real *fshift
#endif
#ifdef CALC_ENERGIES
- ,
- real *Vvdw,
- real *Vc
+ ,
+ real *Vvdw,
+ real *Vc
#endif
- )
+)
{
const nbnxn_ci_t *nbln;
const nbnxn_cj_t *l_cj;
const real *shiftvec;
const real *x;
const real *nbfp;
- real rcut2;
+ real rcut2;
#ifdef VDW_CUTOFF_CHECK
- real rvdw2;
+ real rvdw2;
#endif
- int ntype2;
- real facel;
- real *nbfp_i;
- int n,ci,ci_sh;
- int ish,ishf;
- gmx_bool half_LJ,do_coul;
- int cjind0,cjind1,cjind;
- int ip,jp;
+ int ntype2;
+ real facel;
+ real *nbfp_i;
+ int n, ci, ci_sh;
+ int ish, ishf;
+ gmx_bool do_LJ, half_LJ, do_coul, do_self;
+ int cjind0, cjind1, cjind;
+ int ip, jp;
- real xi[UNROLLI*XI_STRIDE];
- real fi[UNROLLI*FI_STRIDE];
- real qi[UNROLLI];
+ real xi[UNROLLI*XI_STRIDE];
+ real fi[UNROLLI*FI_STRIDE];
+ real qi[UNROLLI];
#ifdef CALC_ENERGIES
#ifndef ENERGY_GROUPS
- real Vvdw_ci,Vc_ci;
+ real Vvdw_ci, Vc_ci;
#else
int egp_mask;
int egp_sh_i[UNROLLI];
#endif
- real sh_invrc6;
+#endif
+#ifdef LJ_POT_SWITCH
+ real swV3, swV4, swV5;
+ real swF2, swF3, swF4;
+#endif
+#ifdef LJ_EWALD
+ real lje_coeff2, lje_coeff6_6, lje_vc;
+ const real *ljc;
#endif
#ifdef CALC_COUL_RF
real k_rf2;
#ifdef CALC_ENERGIES
- real k_rf,c_rf;
+ real k_rf, c_rf;
#endif
#endif
#ifdef CALC_COUL_TAB
int ninner;
#ifdef COUNT_PAIRS
- int npair=0;
+ int npair = 0;
#endif
-#ifdef CALC_ENERGIES
- sh_invrc6 = ic->sh_invrc6;
+#ifdef LJ_POT_SWITCH
+ swV3 = ic->vdw_switch.c3;
+ swV4 = ic->vdw_switch.c4;
+ swV5 = ic->vdw_switch.c5;
+ swF2 = 3*ic->vdw_switch.c3;
+ swF3 = 4*ic->vdw_switch.c4;
+ swF4 = 5*ic->vdw_switch.c5;
+#endif
+
+#ifdef LJ_EWALD
+ lje_coeff2 = ic->ewaldcoeff_lj*ic->ewaldcoeff_lj;
+ lje_coeff6_6 = lje_coeff2*lje_coeff2*lje_coeff2/6.0;
+ lje_vc = ic->sh_lj_ewald;
+
+ ljc = nbat->nbfp_comb;
#endif
#ifdef CALC_COUL_RF
l_cj = nbl->cj;
ninner = 0;
- for(n=0; n<nbl->nci; n++)
+ for (n = 0; n < nbl->nci; n++)
{
- int i,d;
+ int i, d;
nbln = &nbl->ci[n];
ci = nbln->ci;
ci_sh = (ish == CENTRAL ? ci : -1);
- half_LJ = (nbln->shift & NBNXN_CI_HALF_LJ(0));
+ /* We have 5 LJ/C combinations, but use only three inner loops,
+ * as the other combinations are unlikely and/or not much faster:
+ * inner half-LJ + C for half-LJ + C / no-LJ + C
+ * inner LJ + C for full-LJ + C
+ * inner LJ for full-LJ + no-C / half-LJ + no-C
+ */
+ do_LJ = (nbln->shift & NBNXN_CI_DO_LJ(0));
do_coul = (nbln->shift & NBNXN_CI_DO_COUL(0));
+ half_LJ = ((nbln->shift & NBNXN_CI_HALF_LJ(0)) || !do_LJ) && do_coul;
+#ifdef LJ_EWALD
+ do_self = TRUE;
+#else
+ do_self = do_coul;
+#endif
#ifdef CALC_ENERGIES
#ifndef ENERGY_GROUPS
Vvdw_ci = 0;
Vc_ci = 0;
#else
- for(i=0; i<UNROLLI; i++)
+ for (i = 0; i < UNROLLI; i++)
{
egp_sh_i[i] = ((nbat->energrp[ci]>>(i*nbat->neg_2log)) & egp_mask)*nbat->nenergrp;
}
#endif
#endif
- for(i=0; i<UNROLLI; i++)
+ for (i = 0; i < UNROLLI; i++)
{
- for(d=0; d<DIM; d++)
+ for (d = 0; d < DIM; d++)
{
xi[i*XI_STRIDE+d] = x[(ci*UNROLLI+i)*X_STRIDE+d] + shiftvec[ishf+d];
fi[i*FI_STRIDE+d] = 0;
}
+
+ qi[i] = facel*q[ci*UNROLLI+i];
}
- /* With half_LJ we currently always calculate Coulomb interactions */
- if (do_coul || half_LJ)
- {
#ifdef CALC_ENERGIES
+ if (do_self)
+ {
real Vc_sub_self;
#ifdef CALC_COUL_RF
#else
Vc_sub_self = 0.5*tab_coul_FDV0[2];
#endif
-#endif
#endif
- for(i=0; i<UNROLLI; i++)
+ if (l_cj[nbln->cj_ind_start].cj == ci_sh)
{
- qi[i] = facel*q[ci*UNROLLI+i];
-
-#ifdef CALC_ENERGIES
- if (l_cj[nbln->cj_ind_start].cj == ci_sh)
+ for (i = 0; i < UNROLLI; i++)
{
+ int egp_ind;
#ifdef ENERGY_GROUPS
- Vc[egp_sh_i[i]+((nbat->energrp[ci]>>(i*nbat->neg_2log)) & egp_mask)]
+ egp_ind = egp_sh_i[i] + ((nbat->energrp[ci]>>(i*nbat->neg_2log)) & egp_mask);
#else
- Vc[0]
+ egp_ind = 0;
#endif
- -= qi[i]*q[ci*UNROLLI+i]*Vc_sub_self;
- }
+ /* Coulomb self interaction */
+ Vc[egp_ind] -= qi[i]*q[ci*UNROLLI+i]*Vc_sub_self;
+
+#ifdef LJ_EWALD
+ /* LJ Ewald self interaction */
+ Vvdw[egp_ind] += 0.5*nbat->nbfp[nbat->type[ci*UNROLLI+i]*(nbat->ntype + 1)*2]/6*lje_coeff6_6;
#endif
+ }
}
}
+#endif /* CALC_ENERGIES */
cjind = cjind0;
while (cjind < cjind1 && nbl->cj[cjind].excl != 0xffff)
{
#define CALC_COULOMB
#define HALF_LJ
-#include "nbnxn_kernel_ref_inner.h"
+#include "gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_ref_inner.h"
#undef HALF_LJ
#undef CALC_COULOMB
}
- /* cppcheck-suppress duplicateBranch */
else if (do_coul)
{
#define CALC_COULOMB
-#include "nbnxn_kernel_ref_inner.h"
+#include "gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_ref_inner.h"
#undef CALC_COULOMB
}
else
{
-#include "nbnxn_kernel_ref_inner.h"
+#include "gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_ref_inner.h"
}
#undef CHECK_EXCLS
cjind++;
}
- for(; (cjind<cjind1); cjind++)
+ for (; (cjind < cjind1); cjind++)
{
if (half_LJ)
{
#define CALC_COULOMB
#define HALF_LJ
-#include "nbnxn_kernel_ref_inner.h"
+#include "gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_ref_inner.h"
#undef HALF_LJ
#undef CALC_COULOMB
}
- /* cppcheck-suppress duplicateBranch */
else if (do_coul)
{
#define CALC_COULOMB
-#include "nbnxn_kernel_ref_inner.h"
+#include "gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_ref_inner.h"
#undef CALC_COULOMB
}
else
{
-#include "nbnxn_kernel_ref_inner.h"
+#include "gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_ref_inner.h"
}
}
ninner += cjind1 - cjind0;
/* Add accumulated i-forces to the force array */
- for(i=0; i<UNROLLI; i++)
+ for (i = 0; i < UNROLLI; i++)
{
- for(d=0; d<DIM; d++)
+ for (d = 0; d < DIM; d++)
{
f[(ci*UNROLLI+i)*F_STRIDE+d] += fi[i*FI_STRIDE+d];
}
if (fshift != NULL)
{
/* Add i forces to shifted force list */
- for(i=0; i<UNROLLI; i++)
+ for (i = 0; i < UNROLLI; i++)
{
- for(d=0; d<DIM; d++)
+ for (d = 0; d < DIM; d++)
{
fshift[ishf+d] += fi[i*FI_STRIDE+d];
}
*Vc += Vc_ci;
#endif
#endif
- }
+ }
#ifdef COUNT_PAIRS
- printf("atom pairs %d\n",npair);
+ printf("atom pairs %d\n", npair);
#endif
}
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff