From: Mark Abraham Date: Wed, 27 Aug 2014 12:15:24 +0000 (+0200) Subject: Move pair-interaction code into listed-forces module X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?p=alexxy%2Fgromacs.git;a=commitdiff_plain;h=4488887e17c2ad98c5885a8181fc460158e1476a Move pair-interaction code into listed-forces module Did C++ conversion, added Doxygen, renamed functions to reflect that they are now local to a source file, and put them into an anonymous namespace. Made glatnr() local to the listed-forces module. Converted IS_LISTED_LJ_C to static function and made it local to bonded module. Removed temporary cycle-warning suppression. Change-Id: If3a3045021fe7a3f1d67d14f2294551fc41b5422 --- diff --git a/docs/doxygen/cycle-suppressions.txt b/docs/doxygen/cycle-suppressions.txt index c49d9ed7ea..b87e5469ca 100644 --- a/docs/doxygen/cycle-suppressions.txt +++ b/docs/doxygen/cycle-suppressions.txt @@ -10,4 +10,3 @@ topology -> fileio topology -> legacyheaders pbcutil -> fileio pbcutil -> legacyheaders -listed-forces -> gmxlib diff --git a/src/gromacs/listed-forces/bonded.cpp b/src/gromacs/listed-forces/bonded.cpp index 6afb559bba..2415f501a6 100644 --- a/src/gromacs/listed-forces/bonded.cpp +++ b/src/gromacs/listed-forces/bonded.cpp @@ -55,7 +55,6 @@ #include -#include "gromacs/gmxlib/nonbonded/pairs.h" #include "gromacs/legacyheaders/disre.h" #include "gromacs/legacyheaders/force.h" #include "gromacs/legacyheaders/macros.h" @@ -75,6 +74,8 @@ #include "gromacs/utility/fatalerror.h" #include "gromacs/utility/smalloc.h" +#include "listed-internal.h" +#include "pairs.h" #include "restcbt.h" /*! \brief Mysterious CMAP coefficient matrix */ @@ -98,25 +99,6 @@ const int cmap_coeff_matrix[] = { }; -/* TODO This function should go and live in nonbonded.c where it is - really needed. Here, it only supports giving a fatal error message - with FENE_bonds */ -int glatnr(int *global_atom_index, int i) -{ - int atnr; - - if (global_atom_index == NULL) - { - atnr = i + 1; - } - else - { - atnr = global_atom_index[i] + 1; - } - - return atnr; -} - /*! \brief Compute dx = xi - xj, modulo PBC if non-NULL * * \todo This kind of code appears in many places. Consolidate it */ @@ -4308,6 +4290,15 @@ real tab_dihs(int nbonds, //! \endcond +/*! \brief Return true if ftype is an explicit pair-listed LJ or + * COULOMB interaction type: bonded LJ (usually 1-4), or special + * listed non-bonded for FEP. */ +static bool +isPairInteraction(int ftype) +{ + return ((ftype) >= F_LJ14 && (ftype) <= F_LJC_PAIRS_NB); +} + gmx_bool ftype_is_bonded_potential(int ftype) { @@ -4512,7 +4503,7 @@ static real calc_one_bond(int thread, nb0 = idef->il_thread_division[ftype*(idef->nthreads+1)+thread]; nbn = idef->il_thread_division[ftype*(idef->nthreads+1)+thread+1] - nb0; - if (!IS_LISTED_LJ_C(ftype)) + if (!isPairInteraction(ftype)) { if (ftype == F_CMAP) { @@ -4575,8 +4566,8 @@ static real calc_one_bond(int thread, } else { - v = do_nonbonded_listed(ftype, nbn, iatoms+nb0, idef->iparams, x, f, fshift, - pbc, g, lambda, dvdl, md, fr, grpp, global_atom_index); + v = do_pairs(ftype, nbn, iatoms+nb0, idef->iparams, x, f, fshift, + pbc, g, lambda, dvdl, md, fr, grpp, global_atom_index); } if (thread == 0) diff --git a/src/gromacs/listed-forces/bonded.h b/src/gromacs/listed-forces/bonded.h index b3acfd3513..1b6d87f21e 100644 --- a/src/gromacs/listed-forces/bonded.h +++ b/src/gromacs/listed-forces/bonded.h @@ -72,13 +72,6 @@ extern "C" { struct t_graph; struct t_pbc; -/*! \brief Returns the global topology atom number belonging to local atom index i. - * This function is intended for writing ascii output - * and returns atom numbers starting at 1. - * When global_atom_index=NULL returns i+1. - */ -int glatnr(int *global_atom_index, int i); - /*! \brief Return whether this is an interaction that actually * calculates a potential and works on multiple atoms (not e.g. a * connection or a position restraint). diff --git a/src/gromacs/listed-forces/listed-internal.cpp b/src/gromacs/listed-forces/listed-internal.cpp new file mode 100644 index 0000000000..fd69378866 --- /dev/null +++ b/src/gromacs/listed-forces/listed-internal.cpp @@ -0,0 +1,62 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 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. + * + * 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. + * + * 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. + * + * 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. + * + * 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. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * + * \brief This file defines functions needed internally by the module. + * + * \author Mark Abraham + * \ingroup module_listed-forces + */ +#include "gmxpre.h" + +#include "listed-internal.h" + +#include + +int glatnr(int *global_atom_index, int i) +{ + int atnr; + + if (global_atom_index == NULL) + { + atnr = i + 1; + } + else + { + atnr = global_atom_index[i] + 1; + } + + return atnr; +} diff --git a/src/gromacs/listed-forces/listed-internal.h b/src/gromacs/listed-forces/listed-internal.h new file mode 100644 index 0000000000..7e2982bc1d --- /dev/null +++ b/src/gromacs/listed-forces/listed-internal.h @@ -0,0 +1,55 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 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. + * + * 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. + * + * 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. + * + * 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. + * + * 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. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * + * \brief This file contains declarations for functions needed + * internally by the module. + * + * \author Mark Abraham + * \ingroup module_listed-forces + */ +#ifndef GMX_LISTED_FORCES_LISTED_INTERNAL_H +#define GMX_LISTED_FORCES_LISTED_INTERNAL_H + +/*! \brief Returns the global topology atom number belonging to local + * atom index i. + * + * This function is intended for writing ascii output and returns atom + * numbers starting at 1. When global_atom_index=NULL returns i+1. + */ +int +glatnr(int *global_atom_index, int i); + +#endif diff --git a/src/gromacs/gmxlib/nonbonded/pairs.c b/src/gromacs/listed-forces/pairs.cpp similarity index 82% rename from src/gromacs/gmxlib/nonbonded/pairs.c rename to src/gromacs/listed-forces/pairs.cpp index d24d280cf3..1cb683f7b4 100644 --- a/src/gromacs/gmxlib/nonbonded/pairs.c +++ b/src/gromacs/listed-forces/pairs.cpp @@ -1,9 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 1991-2000, University of Groningen, The Netherlands. - * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014, by the GROMACS development team, led by + * Copyright (c) 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. @@ -34,14 +32,22 @@ * To help us fund GROMACS development, we humbly ask that you cite * the research papers on the package. Check out http://www.gromacs.org. */ +/*! \internal \file + * + * \brief This file defines functions for "pair" interactions + * (i.e. listed non-bonded interactions, e.g. 1-4 interactions) + * + * \author Mark Abraham + * + * \ingroup module_listed-forces + */ #include "gmxpre.h" #include "pairs.h" -#include +#include #include "gromacs/legacyheaders/types/group.h" -#include "gromacs/listed-forces/bonded.h" #include "gromacs/math/vec.h" #include "gromacs/pbcutil/ishift.h" #include "gromacs/pbcutil/mshift.h" @@ -49,8 +55,14 @@ #include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/fatalerror.h" -static void -nb_listed_warning_rlimit(const rvec *x, int ai, int aj, int * global_atom_index, real r, real rlimit) +#include "listed-internal.h" + +namespace +{ + +/*! \brief Issue a warning if a listed interaction is beyond a table limit */ +void +warning_rlimit(const rvec *x, int ai, int aj, int * global_atom_index, real r, real rlimit) { gmx_warning("Listed nonbonded interaction between particles %d and %d\n" "at distance %.3f which is larger than the table limit %.3f nm.\n\n" @@ -73,15 +85,10 @@ nb_listed_warning_rlimit(const rvec *x, int ai, int aj, int * global_atom_index, } } - - -/* This might logically belong better in the nb_generic.c module, but it is only - * used in do_nonbonded_listed(), and we want it to be inlined there to avoid an - * extra functional call for every single pair listed in the topology. - */ -static real -nb_evaluate_single(real r2, real tabscale, real *vftab, - real qq, real c6, real c12, real *velec, real *vvdw) +/*! \brief Compute the energy and force for a single pair interaction */ +real +evaluate_single(real r2, real tabscale, real *vftab, + real qq, real c6, real c12, real *velec, real *vvdw) { real rinv, r, rtab, eps, eps2, Y, F, Geps, Heps2, Fp, VVe, FFe, VVd, FFd, VVr, FFr, fscal; int ntab; @@ -90,7 +97,7 @@ nb_evaluate_single(real r2, real tabscale, real *vftab, rinv = gmx_invsqrt(r2); r = r2*rinv; rtab = r*tabscale; - ntab = rtab; + ntab = static_cast(rtab); eps = rtab-ntab; eps2 = eps*eps; ntab = 12*ntab; @@ -127,24 +134,27 @@ nb_evaluate_single(real r2, real tabscale, real *vftab, return fscal; } -static real -nb_free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, real alpha_vdw, - real tabscale, real *vftab, - real qqA, real c6A, real c12A, real qqB, real c6B, real c12B, - real LFC[2], real LFV[2], real DLF[2], - real lfac_coul[2], real lfac_vdw[2], real dlfac_coul[2], real dlfac_vdw[2], - real sigma6_def, real sigma6_min, real sigma2_def, real sigma2_min, - real *velectot, real *vvdwtot, real *dvdl) +/*! \brief Compute the energy and force for a single pair interaction under FEP */ +real +free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, + real alpha_vdw, real tabscale, real *vftab, + real qqA, real c6A, real c12A, real qqB, + real c6B, real c12B, real LFC[2], real LFV[2], real DLF[2], + real lfac_coul[2], real lfac_vdw[2], real dlfac_coul[2], + real dlfac_vdw[2], real sigma6_def, real sigma6_min, + real sigma2_def, real sigma2_min, + real *velectot, real *vvdwtot, real *dvdl) { - real r, rp, rpm2, rtab, eps, eps2, Y, F, Geps, Heps2, Fp, VV, FF, fscal; - real qq[2], c6[2], c12[2], sigma6[2], sigma2[2], sigma_pow[2], sigma_powm2[2]; + real rp, rpm2, rtab, eps, eps2, Y, F, Geps, Heps2, Fp, VV, FF, fscal; + real qq[2], c6[2], c12[2], sigma6[2], sigma2[2], sigma_pow[2]; real alpha_coul_eff, alpha_vdw_eff, dvdl_coul, dvdl_vdw; real rpinv, r_coul, r_vdw, velecsum, vvdwsum; real fscal_vdw[2], fscal_elec[2]; real velec[2], vvdw[2]; int i, ntab; - const real half = 0.5; + const real minusOne = -1.0; + const real oneThird = 1.0 / 3.0; const real one = 1.0; const real two = 2.0; const real six = 6.0; @@ -172,7 +182,7 @@ nb_free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, real a } else { - rp = pow(r2, half*sc_r_power); /* not currently supported as input, but can handle it */ + rp = std::pow(r2, half * sc_r_power); /* not currently supported as input, but can handle it */ rpm2 = rp/r2; } @@ -185,7 +195,7 @@ nb_free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, real a * Correct for this by multiplying with (1/12.0)/(1/6.0)=6.0/12.0=0.5. */ sigma6[i] = half*c12[i]/c6[i]; - sigma2[i] = pow(half*c12[i]/c6[i], 1.0/3.0); + sigma2[i] = std::pow(half*c12[i]/c6[i], oneThird); /* should be able to get rid of this ^^^ internal pow call eventually. Will require agreement on what data to store externally. Can't be fixed without larger scale changes, so not 5.0 */ if (sigma6[i] < sigma6_min) /* for disappearing coul and vdw with soft core at the same time */ @@ -202,19 +212,16 @@ nb_free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, real a if (sc_r_power == six) { sigma_pow[i] = sigma6[i]; - sigma_powm2[i] = sigma6[i]/sigma2[i]; } else if (sc_r_power == fourtyeight) { sigma_pow[i] = sigma6[i]*sigma6[i]; /* sigma^12 */ sigma_pow[i] = sigma_pow[i]*sigma_pow[i]; /* sigma^24 */ sigma_pow[i] = sigma_pow[i]*sigma_pow[i]; /* sigma^48 */ - sigma_powm2[i] = sigma_pow[i]/sigma2[i]; } else - { /* not really supported as input, but in here for testing the general case*/ - sigma_pow[i] = pow(sigma2[i], sc_r_power/2); - sigma_powm2[i] = sigma_pow[i]/(sigma2[i]); + { /* not really supported as input, but in here for testing the general case*/ + sigma_pow[i] = std::pow(sigma2[i], sc_r_power/2); } } @@ -243,11 +250,11 @@ nb_free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, real a { /* Coulomb */ rpinv = one/(alpha_coul_eff*lfac_coul[i]*sigma_pow[i]+rp); - r_coul = pow(rpinv, -one/sc_r_power); + r_coul = std::pow(rpinv, minusOne / sc_r_power); /* Electrostatics table lookup data */ rtab = r_coul*tabscale; - ntab = rtab; + ntab = static_cast(rtab); eps = rtab-ntab; eps2 = eps*eps; ntab = 12*ntab; @@ -264,10 +271,10 @@ nb_free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, real a /* Vdw */ rpinv = one/(alpha_vdw_eff*lfac_vdw[i]*sigma_pow[i]+rp); - r_vdw = pow(rpinv, -one/sc_r_power); + r_vdw = std::pow(rpinv, minusOne / sc_r_power); /* Vdw table lookup data */ rtab = r_vdw*tabscale; - ntab = rtab; + ntab = static_cast(rtab); eps = rtab-ntab; eps2 = eps*eps; ntab = 12*ntab; @@ -322,23 +329,24 @@ nb_free_energy_evaluate_single(real r2, real sc_r_power, real alpha_coul, real a return fscal; } +} // namespace + real -do_nonbonded_listed(int ftype, int nbonds, - const t_iatom iatoms[], const t_iparams iparams[], - const rvec x[], rvec f[], rvec fshift[], - const t_pbc *pbc, const t_graph *g, - real *lambda, real *dvdl, - const t_mdatoms *md, - const t_forcerec *fr, gmx_grppairener_t *grppener, - int *global_atom_index) +do_pairs(int ftype, int nbonds, + const t_iatom iatoms[], const t_iparams iparams[], + const rvec x[], rvec f[], rvec fshift[], + const struct t_pbc *pbc, const struct t_graph *g, + real *lambda, real *dvdl, + const t_mdatoms *md, + const t_forcerec *fr, gmx_grppairener_t *grppener, + int *global_atom_index) { - int ielec, ivdw; real qq, c6, c12; rvec dx; ivec dt; - int i, j, itype, ai, aj, gid; + int i, itype, ai, aj, gid; int fshift_index; - real r2, rinv; + real r2; real fscal, velec, vvdw; real * energygrp_elec; real * energygrp_vdw; @@ -347,7 +355,7 @@ do_nonbonded_listed(int ftype, int nbonds, gmx_bool bFreeEnergy; real LFC[2], LFV[2], DLF[2], lfac_coul[2], lfac_vdw[2], dlfac_coul[2], dlfac_vdw[2]; real qqB, c6B, c12B, sigma2_def, sigma2_min; - + const real oneThird = 1.0 / 3.0; switch (ftype) { @@ -380,8 +388,8 @@ do_nonbonded_listed(int ftype, int nbonds, DLF[1] = 1; /* precalculate */ - sigma2_def = pow(fr->sc_sigma6_def, 1.0/3.0); - sigma2_min = pow(fr->sc_sigma6_min, 1.0/3.0); + sigma2_def = pow(fr->sc_sigma6_def, oneThird); + sigma2_min = pow(fr->sc_sigma6_min, oneThird); for (i = 0; i < 2; i++) { @@ -459,7 +467,7 @@ do_nonbonded_listed(int ftype, int nbonds, * disadvantage is that the warning is printed twice */ if (warned_rlimit == FALSE) { - nb_listed_warning_rlimit(x, ai, aj, global_atom_index, sqrt(r2), fr->tab14.r); + warning_rlimit(x, ai, aj, global_atom_index, sqrt(r2), fr->tab14.r); warned_rlimit = TRUE; } continue; @@ -472,15 +480,15 @@ do_nonbonded_listed(int ftype, int nbonds, c6B = iparams[itype].lj14.c6B*6.0; c12B = iparams[itype].lj14.c12B*12.0; - fscal = nb_free_energy_evaluate_single(r2, fr->sc_r_power, fr->sc_alphacoul, fr->sc_alphavdw, - fr->tab14.scale, fr->tab14.data, qq, c6, c12, qqB, c6B, c12B, - LFC, LFV, DLF, lfac_coul, lfac_vdw, dlfac_coul, dlfac_vdw, - fr->sc_sigma6_def, fr->sc_sigma6_min, sigma2_def, sigma2_min, &velec, &vvdw, dvdl); + fscal = free_energy_evaluate_single(r2, fr->sc_r_power, fr->sc_alphacoul, fr->sc_alphavdw, + fr->tab14.scale, fr->tab14.data, qq, c6, c12, qqB, c6B, c12B, + LFC, LFV, DLF, lfac_coul, lfac_vdw, dlfac_coul, dlfac_vdw, + fr->sc_sigma6_def, fr->sc_sigma6_min, sigma2_def, sigma2_min, &velec, &vvdw, dvdl); } else { /* Evaluate tabulated interaction without free energy */ - fscal = nb_evaluate_single(r2, fr->tab14.scale, fr->tab14.data, qq, c6, c12, &velec, &vvdw); + fscal = evaluate_single(r2, fr->tab14.scale, fr->tab14.data, qq, c6, c12, &velec, &vvdw); } energygrp_elec[gid] += velec; diff --git a/src/gromacs/gmxlib/nonbonded/pairs.h b/src/gromacs/listed-forces/pairs.h similarity index 70% rename from src/gromacs/gmxlib/nonbonded/pairs.h rename to src/gromacs/listed-forces/pairs.h index 4582b08aef..c895f81b85 100644 --- a/src/gromacs/gmxlib/nonbonded/pairs.h +++ b/src/gromacs/listed-forces/pairs.h @@ -1,9 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 1991-2000, University of Groningen, The Netherlands. - * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014, by the GROMACS development team, led by + * Copyright (c) 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. @@ -34,9 +32,17 @@ * To help us fund GROMACS development, we humbly ask that you cite * the research papers on the package. Check out http://www.gromacs.org. */ - -#ifndef _pairs_h -#define _pairs_h +/*! \internal \file + * + * \brief This file declares functions for "pair" interactions + * (i.e. listed non-bonded interactions, e.g. 1-4 interactions) + * + * \author Mark Abraham + * + * \ingroup module_listed-forces + */ +#ifndef GMX_LISTED_FORCES_PAIRS_H +#define GMX_LISTED_FORCES_PAIRS_H #include "gromacs/legacyheaders/types/forcerec.h" #include "gromacs/legacyheaders/types/mdatom.h" @@ -44,28 +50,19 @@ #include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/real.h" -#ifdef __cplusplus -extern "C" { -#endif -#if 0 -} /* fixes auto-indentation problems */ -#endif - struct t_graph; struct t_pbc; -/* Calculate VdW/charge listed pair interactions (usually 1-4 interactions). +/*! \brief Calculate VdW/charge listed pair interactions (usually 1-4 + * interactions). + * * global_atom_index is only passed for printing error messages. */ real -do_nonbonded_listed(int ftype, int nbonds, const t_iatom iatoms[], const t_iparams iparams[], - const rvec x[], rvec f[], rvec fshift[], - const struct t_pbc *pbc, const struct t_graph *g, - real *lambda, real *dvdl, const t_mdatoms *md, const t_forcerec *fr, - gmx_grppairener_t *grppener, int *global_atom_index); - -#ifdef __cplusplus -} -#endif +do_pairs(int ftype, int nbonds, const t_iatom iatoms[], const t_iparams iparams[], + const rvec x[], rvec f[], rvec fshift[], + const struct t_pbc *pbc, const struct t_graph *g, + real *lambda, real *dvdl, const t_mdatoms *md, const t_forcerec *fr, + gmx_grppairener_t *grppener, int *global_atom_index); #endif diff --git a/src/gromacs/topology/idef.h b/src/gromacs/topology/idef.h index aaac95175f..e243e453f1 100644 --- a/src/gromacs/topology/idef.h +++ b/src/gromacs/topology/idef.h @@ -152,12 +152,6 @@ enum { #define IS_RESTRAINT_TYPE(ifunc) (((ifunc == F_POSRES) || (ifunc == F_FBPOSRES) || (ifunc == F_DISRES) || (ifunc == F_RESTRBONDS) || (ifunc == F_DISRESVIOL) || (ifunc == F_ORIRES) || (ifunc == F_ORIRESDEV) || (ifunc == F_ANGRES) || (ifunc == F_ANGRESZ) || (ifunc == F_DIHRES))) -/* A macro for checking if ftype is an explicit pair-listed LJ or COULOMB - * interaction type: - * bonded LJ (usually 1-4), or special listed non-bonded for FEP. - */ -#define IS_LISTED_LJ_C(ftype) ((ftype) >= F_LJ14 && (ftype) <= F_LJC_PAIRS_NB) - typedef union t_iparams { /* Some parameters have A and B values for free energy calculations.