*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2018, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2018 by the GROMACS development team.
+ * Copyright (c) 2019,2020,2021, 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.
#ifndef GMX_TOPOLOGY_IDEF_H
#define GMX_TOPOLOGY_IDEF_H
-#include <stdio.h>
+#include <cstdio>
+
+#include <array>
+#include <vector>
#include "gromacs/math/vectypes.h"
-#include "gromacs/utility/basedefinitions.h"
+#include "gromacs/topology/ifunc.h"
#include "gromacs/utility/real.h"
-#ifdef __cplusplus
-extern "C" {
-#endif
-
-/* check kernel/toppush.c when you change these numbers */
-#define MAXATOMLIST 6
-#define MAXFORCEPARAM 12
-#define NR_RBDIHS 6
-#define NR_CBTDIHS 6
-#define NR_FOURDIHS 4
-
-typedef int t_iatom;
-
-/* this MUST correspond to the
- t_interaction_function[F_NRE] in gmxlib/ifunc.c */
-enum {
- F_BONDS,
- F_G96BONDS,
- F_MORSE,
- F_CUBICBONDS,
- F_CONNBONDS,
- F_HARMONIC,
- F_FENEBONDS,
- F_TABBONDS,
- F_TABBONDSNC,
- F_RESTRBONDS,
- F_ANGLES,
- F_G96ANGLES,
- F_RESTRANGLES,
- F_LINEAR_ANGLES,
- F_CROSS_BOND_BONDS,
- F_CROSS_BOND_ANGLES,
- F_UREY_BRADLEY,
- F_QUARTIC_ANGLES,
- F_TABANGLES,
- F_PDIHS,
- F_RBDIHS,
- F_RESTRDIHS,
- F_CBTDIHS,
- F_FOURDIHS,
- F_IDIHS,
- F_PIDIHS,
- F_TABDIHS,
- F_CMAP,
- F_GB12_NOLONGERUSED,
- F_GB13_NOLONGERUSED,
- F_GB14_NOLONGERUSED,
- F_GBPOL_NOLONGERUSED,
- F_NPSOLVATION_NOLONGERUSED,
- F_LJ14,
- F_COUL14,
- F_LJC14_Q,
- F_LJC_PAIRS_NB,
- F_LJ,
- F_BHAM,
- F_LJ_LR_NOLONGERUSED,
- F_BHAM_LR_NOLONGERUSED,
- F_DISPCORR,
- F_COUL_SR,
- F_COUL_LR_NOLONGERUSED,
- F_RF_EXCL,
- F_COUL_RECIP,
- F_LJ_RECIP,
- F_DPD,
- F_POLARIZATION,
- F_WATER_POL,
- F_THOLE_POL,
- F_ANHARM_POL,
- F_POSRES,
- F_FBPOSRES,
- F_DISRES,
- F_DISRESVIOL,
- F_ORIRES,
- F_ORIRESDEV,
- F_ANGRES,
- F_ANGRESZ,
- F_DIHRES,
- F_DIHRESVIOL,
- F_CONSTR,
- F_CONSTRNC,
- F_SETTLE,
- F_VSITE2,
- F_VSITE3,
- F_VSITE3FD,
- F_VSITE3FAD,
- F_VSITE3OUT,
- F_VSITE4FD,
- F_VSITE4FDN,
- F_VSITEN,
- F_COM_PULL,
- F_EQM,
- F_EPOT,
- F_EKIN,
- F_ETOT,
- F_ECONSERVED,
- F_TEMP,
- F_VTEMP_NOLONGERUSED,
- F_PDISPCORR,
- F_PRES,
- F_DVDL_CONSTR,
- F_DVDL,
- F_DKDL,
- F_DVDL_COUL,
- F_DVDL_VDW,
- F_DVDL_BONDED,
- F_DVDL_RESTRAINT,
- F_DVDL_TEMPERATURE, /* not calculated for now, but should just be the energy (NVT) or enthalpy (NPT), or 0 (NVE) */
- F_NRE /* This number is for the total number of energies */
-};
-
-#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)))
+struct gmx_ffparams_t;
typedef union t_iparams
{
* The harmonic type is used for all harmonic potentials:
* bonds, angles and improper dihedrals
*/
- struct {
+ struct
+ {
real a, b, c;
} bham;
- struct {
+ struct
+ {
real rA, krA, rB, krB;
} harmonic;
- struct {
+ struct
+ {
real klinA, aA, klinB, aB;
} linangle;
- struct {
+ struct
+ {
real lowA, up1A, up2A, kA, lowB, up1B, up2B, kB;
} restraint;
/* No free energy supported for cubic bonds, FENE, WPOL or cross terms */
- struct {
+ struct
+ {
real b0, kb, kcub;
} cubic;
- struct {
+ struct
+ {
real bm, kb;
} fene;
- struct {
+ struct
+ {
real r1e, r2e, krr;
} cross_bb;
- struct {
+ struct
+ {
real r1e, r2e, r3e, krt;
} cross_ba;
- struct {
+ struct
+ {
real thetaA, kthetaA, r13A, kUBA, thetaB, kthetaB, r13B, kUBB;
} u_b;
- struct {
+ struct
+ {
real theta, c[5];
} qangle;
- struct {
+ struct
+ {
real alpha;
} polarize;
- struct {
+ struct
+ {
real alpha, drcut, khyp;
} anharm_polarize;
- struct {
+ struct
+ {
real al_x, al_y, al_z, rOH, rHH, rOD;
} wpol;
- struct {
- real a, alpha1, alpha2, rfac;
+ struct
+ {
+ real a, alpha1, alpha2;
} thole;
- struct {
+ struct
+ {
real c6, c12;
} lj;
- struct {
+ struct
+ {
real c6A, c12A, c6B, c12B;
} lj14;
- struct {
+ struct
+ {
real fqq, qi, qj, c6, c12;
} ljc14;
- struct {
+ struct
+ {
real qi, qj, c6, c12;
} ljcnb;
/* Proper dihedrals can not have different multiplicity when
* doing free energy calculations, because the potential would not
* be periodic anymore.
*/
- struct {
- real phiA, cpA; int mult; real phiB, cpB;
+ struct
+ {
+ real phiA, cpA;
+ int mult;
+ real phiB, cpB;
} pdihs;
- struct {
+ struct
+ {
real dA, dB;
} constr;
/* Settle can not be used for Free energy calculations of water bond geometry.
* Use shake (or lincs) instead if you have to change the water bonds.
*/
- struct {
+ struct
+ {
real doh, dhh;
} settle;
- struct {
+ struct
+ {
real b0A, cbA, betaA, b0B, cbB, betaB;
} morse;
- struct {
+ struct
+ {
real pos0A[DIM], fcA[DIM], pos0B[DIM], fcB[DIM];
} posres;
- struct {
- real pos0[DIM], r, k; int geom;
+ struct
+ {
+ real pos0[DIM], r, k;
+ int geom;
} fbposres;
- struct {
+ struct
+ {
real rbcA[NR_RBDIHS], rbcB[NR_RBDIHS];
} rbdihs;
- struct {
+ struct
+ {
real cbtcA[NR_CBTDIHS], cbtcB[NR_CBTDIHS];
} cbtdihs;
- struct {
+ struct
+ {
real a, b, c, d, e, f;
} vsite;
- struct {
- int n; real a;
+ struct
+ {
+ int n;
+ real a;
} vsiten;
/* NOTE: npair is only set after reading the tpx file */
- struct {
- real low, up1, up2, kfac; int type, label, npair;
+ struct
+ {
+ real low, up1, up2, kfac;
+ int type, label, npair;
} disres;
- struct {
+ struct
+ {
real phiA, dphiA, kfacA, phiB, dphiB, kfacB;
} dihres;
- struct {
- int ex, power, label; real c, obs, kfac;
+ struct
+ {
+ int ex, power, label;
+ real c, obs, kfac;
} orires;
- struct {
- int table; real kA; real kB;
+ struct
+ {
+ int table;
+ real kA;
+ real kB;
} tab;
- struct {
+ struct
+ {
int cmapA, cmapB;
} cmap;
- struct {
+ struct
+ {
real buf[MAXFORCEPARAM];
- } generic; /* Conversion */
+ } generic; /* Conversion */
} t_iparams;
typedef int t_functype;
-/*
- * The nonperturbed/perturbed interactions are now separated (sorted) in the
- * ilist, such that the first 0..(nr_nonperturbed-1) ones are exactly that, and
- * the remaining ones from nr_nonperturbed..(nr-1) are perturbed bonded
- * interactions.
+/* List of listed interactions, see description further down.
+ *
+ * TODO: Consider storing the function type as well.
+ * TODO: Consider providing per interaction access.
*/
-typedef struct t_ilist
+struct InteractionList
{
+ /* Returns the total number of elements in iatoms */
+ int size() const { return static_cast<int>(iatoms.size()); }
+
+ /* Returns whether the list is empty */
+ bool empty() const { return iatoms.empty(); }
+
+ /* Adds one interaction to the list */
+ template<std::size_t numAtoms>
+ void push_back(const int parameterType, const std::array<int, numAtoms>& atoms)
+ {
+ const std::size_t oldSize = iatoms.size();
+ iatoms.resize(iatoms.size() + 1 + numAtoms);
+ iatoms[oldSize] = parameterType;
+ for (std::size_t i = 0; i < numAtoms; i++)
+ {
+ iatoms[oldSize + 1 + i] = atoms[i];
+ }
+ }
+
+ /* Adds one interaction to the list */
+ void push_back(const int parameterType, const int numAtoms, const int* atoms)
+ {
+ const std::size_t oldSize = iatoms.size();
+ iatoms.resize(iatoms.size() + 1 + numAtoms);
+ iatoms[oldSize] = parameterType;
+ for (int i = 0; i < numAtoms; i++)
+ {
+ iatoms[oldSize + 1 + i] = atoms[i];
+ }
+ }
+
+ /* Appends \p ilist at the back of the list */
+ void append(const InteractionList& ilist)
+ {
+ iatoms.insert(iatoms.end(), ilist.iatoms.begin(), ilist.iatoms.end());
+ }
+
+ /* Clears the list */
+ void clear() { iatoms.clear(); }
+
+ /* List of interactions, see explanation further down */
+ std::vector<int> iatoms;
+};
+
+/* List of interaction lists, one list for each interaction type
+ *
+ * TODO: Consider only including entries in use instead of all F_NRE
+ */
+using InteractionLists = std::array<InteractionList, F_NRE>;
+
+/* Deprecated list of listed interactions */
+struct t_ilist
+{
+ /* Returns the total number of elements in iatoms */
+ int size() const { return nr; }
+
+ /* Returns whether the list is empty */
+ bool empty() const { return nr == 0; }
+
int nr;
- int nr_nonperturbed;
- t_iatom *iatoms;
+ t_iatom* iatoms;
int nalloc;
-} t_ilist;
+};
+
+/* TODO: Remove t_ilist and remove templating on list type in mshift.cpp */
/*
- * The struct t_ilist defines a list of atoms with their interactions.
+ * The structs InteractionList and t_ilist defines a list of atoms with their interactions.
* General field description:
* int nr
* the size (nr elements) of the interactions array (iatoms[]).
* identifier is an index in a params[] and functype[] array.
*/
-typedef struct
+/*! \brief Type for returning a list of InteractionList references
+ *
+ * TODO: Remove when the function type is made part of InteractionList
+ */
+struct InteractionListHandle
{
- real *cmap; /* Has length 4*grid_spacing*grid_spacing, */
+ const int functionType; //!< The function type
+ const std::vector<int>& iatoms; //!< Reference to interaction list
+};
+
+/*! \brief Returns a list of all non-empty InteractionList entries with any of the interaction flags in \p flags set
+ *
+ * \param[in] ilists Set of interaction lists
+ * \param[in] flags Bit mask with one or more IF_... bits set
+ */
+static inline std::vector<InteractionListHandle> extractILists(const InteractionLists& ilists, int flags)
+{
+ std::vector<InteractionListHandle> handles;
+ for (size_t ftype = 0; ftype < ilists.size(); ftype++)
+ {
+ if ((interaction_function[ftype].flags & flags) && !ilists[ftype].empty())
+ {
+ handles.push_back({ static_cast<int>(ftype), ilists[ftype].iatoms });
+ }
+ }
+ return handles;
+}
+
+/*! \brief Returns the stride for the iatoms array in \p ilistHandle
+ *
+ * \param[in] ilistHandle The ilist to return the stride for
+ */
+static inline int ilistStride(const InteractionListHandle& ilistHandle)
+{
+ return 1 + NRAL(ilistHandle.functionType);
+}
+
+struct gmx_cmapdata_t
+{
+ std::vector<real> cmap; /* Has length 4*grid_spacing*grid_spacing, */
/* there are 4 entries for each cmap type (V,dVdx,dVdy,d2dVdxdy) */
-} gmx_cmapdata_t;
+};
-typedef struct gmx_cmap_t
+struct gmx_cmap_t
{
- int ngrid; /* Number of allocated cmap (cmapdata_t ) grids */
- int grid_spacing; /* Grid spacing */
- gmx_cmapdata_t *cmapdata; /* Pointer to grid with actual, pre-interpolated data */
-} gmx_cmap_t;
+ int grid_spacing = 0; /* Grid spacing */
+ std::vector<gmx_cmapdata_t> cmapdata; /* Lists of grids with actual, pre-interpolated data */
+};
-typedef struct gmx_ffparams_t
+enum
{
- int ntypes;
- int atnr;
- t_functype *functype;
- t_iparams *iparams;
- double reppow; /* The repulsion power for VdW: C12*r^-reppow */
- real fudgeQQ; /* The scaling factor for Coulomb 1-4: f*q1*q2 */
- gmx_cmap_t cmap_grid; /* The dihedral correction maps */
-} gmx_ffparams_t;
-
-enum {
- ilsortUNKNOWN, ilsortNO_FE, ilsortFE_UNSORTED, ilsortFE_SORTED
+ ilsortUNKNOWN,
+ ilsortNO_FE,
+ ilsortFE_SORTED
+};
+
+/* Struct with list of interaction parameters and lists of interactions
+ *
+ * TODO: Convert to a proper class with private data members so we can
+ * ensure that the free-energy sorting and sorting setting is consistent.
+ */
+class InteractionDefinitions
+{
+public:
+ /* Constructor
+ *
+ * \param[in] ffparams The interaction parameters, the lifetime of the created object should not exceed the lifetime of the passed parameters
+ */
+ InteractionDefinitions(const gmx_ffparams_t& ffparams);
+
+ // Clears data not read in from ffparams
+ void clear();
+
+ // The interaction parameters
+ const std::vector<t_iparams>& iparams;
+ // The function type per type
+ const std::vector<int>& functype;
+ // Position restraint interaction parameters
+ std::vector<t_iparams> iparams_posres;
+ // Flat-bottomed position restraint parameters
+ std::vector<t_iparams> iparams_fbposres;
+ // The list of interactions for each type. Note that some, such as LJ and COUL will have 0 entries.
+ std::array<InteractionList, F_NRE> il;
+ /* The number of non-perturbed interactions at the start of each entry in il */
+ std::array<int, F_NRE> numNonperturbedInteractions;
+ // The sorting state of interaction in il
+ int ilsort = ilsortUNKNOWN;
+ // The dihedral correction maps
+ gmx_cmap_t cmap_grid;
};
-typedef struct t_idef
+/* Deprecated interation definitions, used in t_topology */
+struct t_idef
{
int ntypes;
int atnr;
- t_functype *functype;
- t_iparams *iparams;
+ t_functype* functype;
+ t_iparams* iparams;
real fudgeQQ;
- gmx_cmap_t cmap_grid;
- t_iparams *iparams_posres, *iparams_fbposres;
- int iparams_posres_nalloc, iparams_fbposres_nalloc;
+ t_iparams * iparams_posres, *iparams_fbposres;
- t_ilist il[F_NRE];
- int ilsort;
- int nthreads;
- int *il_thread_division;
- int il_thread_division_nalloc;
-} t_idef;
+ t_ilist il[F_NRE];
+ int ilsort;
+};
/*
* The struct t_idef defines all the interactions for the complete
* such as LJ and COUL will have 0 entries.
* int ilsort
* The state of the sorting of il, values are provided above.
- * int nthreads
- * The number of threads used to set il_thread_division.
- * int *il_thread_division
- * The division of the normal bonded interactions of threads.
- * il_thread_division[ftype*(nthreads+1)+t] contains an index
- * into il[ftype].iatoms; thread th operates on t=th to t=th+1.
- * int il_thread_division_nalloc
- * The allocated size of il_thread_division,
- * should be at least F_NRE*(nthreads+1).
*/
-void pr_iparams(FILE *fp, t_functype ftype, const t_iparams *iparams);
-void pr_ilist(FILE *fp, int indent, const char *title,
- const t_functype *functype, const t_ilist *ilist,
- gmx_bool bShowNumbers,
- gmx_bool bShowParameters, const t_iparams *iparams);
-void pr_ffparams(FILE *fp, int indent, const char *title,
- const gmx_ffparams_t *ffparams, gmx_bool bShowNumbers);
-void pr_idef(FILE *fp, int indent, const char *title, const t_idef *idef,
- gmx_bool bShowNumbers, gmx_bool bShowParameters);
-
-#ifdef __cplusplus
-}
-#endif
+namespace gmx
+{
+class TextWriter;
+} // namespace gmx
+
+void printInteractionParameters(gmx::TextWriter* writer, t_functype ftype, const t_iparams& iparams);
+void pr_iparams(FILE* fp, t_functype ftype, const t_iparams& iparams);
+void pr_ilist(FILE* fp,
+ int indent,
+ const char* title,
+ const t_functype* functype,
+ const InteractionList& ilist,
+ bool bShowNumbers,
+ bool bShowParameters,
+ const t_iparams* iparams);
+void pr_idef(FILE* fp, int indent, const char* title, const t_idef* idef, bool bShowNumbers, bool bShowParameters);
+
+/*! \brief
+ * Properly initialize idef struct.
+ *
+ * \param[in] idef Pointer to idef struct to initialize.
+ */
+void init_idef(t_idef* idef);
+
+/*! \brief
+ * Properly clean up idef struct.
+ *
+ * \param[in] idef Pointer to idef struct to clean up.
+ */
+void done_idef(t_idef* idef);
#endif
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff