1 /*
2 * This file is part of the GROMACS molecular simulation package.
3 *
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2012,2014,2015,2016,2018 by the GROMACS development team.
7 * Copyright (c) 2019,2020, by the GROMACS development team, led by
8 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
9 * and including many others, as listed in the AUTHORS file in the
10 * top-level source directory and at http://www.gromacs.org.
11 *
12 * GROMACS is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU Lesser General Public License
14 * as published by the Free Software Foundation; either version 2.1
15 * of the License, or (at your option) any later version.
16 *
17 * GROMACS is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 * Lesser General Public License for more details.
21 *
22 * You should have received a copy of the GNU Lesser General Public
23 * License along with GROMACS; if not, see
24 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
25 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 *
27 * If you want to redistribute modifications to GROMACS, please
28 * consider that scientific software is very special. Version
29 * control is crucial - bugs must be traceable. We will be happy to
30 * consider code for inclusion in the official distribution, but
31 * derived work must not be called official GROMACS. Details are found
32 * in the README & COPYING files - if they are missing, get the
33 * official version at http://www.gromacs.org.
34 *
35 * To help us fund GROMACS development, we humbly ask that you cite
36 * the research papers on the package. Check out http://www.gromacs.org.
37 */
38 #ifndef GMX_TOPOLOGY_IFUNC_H
39 #define GMX_TOPOLOGY_IFUNC_H
40
41 #include "gromacs/math/vectypes.h"
42
43 struct t_fcdata;
44 struct t_graph;
45 union t_iparams;
46 struct t_mdatoms;
47 struct t_pbc;
48
49 /* TODO: Remove this typedef when t_ilist is removed */
50 typedef int t_iatom;
51
52 /* Real vector type with an additional, unused 4th element.
53 * This type is used to allow aligned 4-wide SIMD loads and stores.
54 */
55 typedef real rvec4[4];
56
57 /*
58 * The function type t_ifunc() calculates one interaction, using iatoms[]
59 * and iparams. Within the function the number of atoms to be used is
60 * known. Within the function only the atomid part of the iatoms[] array
61 * is supplied, not the type field (see also t_ilist). The function
62 * returns the potential energy. If pbc==NULL the coordinates in x are
63 * assumed to be such that no calculation of PBC is necessary,
64 * If pbc!=NULL a full PBC calculation is performed.
65 * If g!=NULL it is used for determining the shift forces.
66 * With domain decomposition ddgatindex can be used for getting global
67 * atom numbers for warnings and error messages.
68 * ddgatindex is NULL when domain decomposition is not used.
69 */
70
71 constexpr unsigned int IF_NULL = 0;
72 constexpr unsigned int IF_BOND = 1 << 0;
73 constexpr unsigned int IF_VSITE = 1 << 1;
74 constexpr unsigned int IF_CONSTRAINT = 1 << 2;
75 constexpr unsigned int IF_CHEMBOND = 1 << 3;
76 constexpr unsigned int IF_BTYPE = 1 << 4;
77 constexpr unsigned int IF_ATYPE = 1 << 5;
78 constexpr unsigned int IF_DIHEDRAL = 1 << 6;
79 constexpr unsigned int IF_PAIR = 1 << 7;
80 constexpr unsigned int IF_TABULATED = 1 << 8;
81 constexpr unsigned int IF_LIMZERO = 1 << 9;
82 /* These flags tell to some of the routines what can be done with this
83 * item in the list.
84 * With IF_BOND a bonded interaction will be calculated.
85 * With IF_BTYPE grompp can convert the bond to a Morse potential.
86 * With IF_BTYPE or IF_ATYPE the bond/angle can be converted to
87 * a constraint or used for vsite parameter determination by grompp.
88 * IF_LIMZERO indicates that for a bonded interaction the potential
89 * does goes to zero for large distances, thus if such an interaction
90 * it not assigned to any node by the domain decompostion, the simulation
91 * still continue, if mdrun has been told so.
92 */
93
94 struct t_interaction_function // NOLINT (clang-analyzer-optin.performance.Padding)
95 {
96 const char* name; /* the name of this function */
97 const char* longname; /* The name for printing etc. */
98 int nratoms; /* nr of atoms needed for this function */
99 int nrfpA, nrfpB; /* number of parameters for this function. */
100 /* this corresponds to the number of params in */
101 /* iparams struct! (see idef.h) */
102 /* A and B are for normal and free energy components respectively. */
103 unsigned int flags; /* Flags (see above) */
104 };
105
106 #define NRFPA(ftype) (interaction_function[(ftype)].nrfpA)
107 #define NRFPB(ftype) (interaction_function[(ftype)].nrfpB)
108 #define NRFP(ftype) (NRFPA(ftype) + NRFPB(ftype))
109 #define NRAL(ftype) (interaction_function[(ftype)].nratoms)
110
111 #define IS_CHEMBOND(ftype) \
112 (interaction_function[(ftype)].nratoms == 2 && (interaction_function[(ftype)].flags & IF_CHEMBOND))
113 /* IS_CHEMBOND tells if function type ftype represents a chemical bond */
114
115 /* IS_ANGLE tells if a function type ftype represents an angle
116 * Per Larsson, 2007-11-06
117 */
118 #define IS_ANGLE(ftype) \
119 (interaction_function[(ftype)].nratoms == 3 && (interaction_function[(ftype)].flags & IF_ATYPE))
120 #define IS_VSITE(ftype) (interaction_function[(ftype)].flags & IF_VSITE)
121
122 #define IS_TABULATED(ftype) (interaction_function[(ftype)].flags & IF_TABULATED)
123
124 /* this MUST correspond to the
125 t_interaction_function[F_NRE] in gmxlib/ifunc.cpp */
126 enum
127 {
128 F_BONDS,
129 F_G96BONDS,
130 F_MORSE,
131 F_CUBICBONDS,
132 F_CONNBONDS,
133 F_HARMONIC,
134 F_FENEBONDS,
135 F_TABBONDS,
136 F_TABBONDSNC,
137 F_RESTRBONDS,
138 F_ANGLES,
139 F_G96ANGLES,
140 F_RESTRANGLES,
141 F_LINEAR_ANGLES,
142 F_CROSS_BOND_BONDS,
143 F_CROSS_BOND_ANGLES,
144 F_UREY_BRADLEY,
145 F_QUARTIC_ANGLES,
146 F_TABANGLES,
147 F_PDIHS,
148 F_RBDIHS,
149 F_RESTRDIHS,
150 F_CBTDIHS,
151 F_FOURDIHS,
152 F_IDIHS,
153 F_PIDIHS,
154 F_TABDIHS,
155 F_CMAP,
156 F_GB12_NOLONGERUSED,
157 F_GB13_NOLONGERUSED,
158 F_GB14_NOLONGERUSED,
159 F_GBPOL_NOLONGERUSED,
160 F_NPSOLVATION_NOLONGERUSED,
161 F_LJ14,
162 F_COUL14,
163 F_LJC14_Q,
164 F_LJC_PAIRS_NB,
165 F_LJ,
166 F_BHAM,
167 F_LJ_LR_NOLONGERUSED,
168 F_BHAM_LR_NOLONGERUSED,
169 F_DISPCORR,
170 F_COUL_SR,
171 F_COUL_LR_NOLONGERUSED,
172 F_RF_EXCL,
173 F_COUL_RECIP,
174 F_LJ_RECIP,
175 F_DPD,
176 F_POLARIZATION,
177 F_WATER_POL,
178 F_THOLE_POL,
179 F_ANHARM_POL,
180 F_POSRES,
181 F_FBPOSRES,
182 F_DISRES,
183 F_DISRESVIOL,
184 F_ORIRES,
185 F_ORIRESDEV,
186 F_ANGRES,
187 F_ANGRESZ,
188 F_DIHRES,
189 F_DIHRESVIOL,
190 F_CONSTR,
191 F_CONSTRNC,
192 F_SETTLE,
193 F_VSITE1,
194 F_VSITE2,
195 F_VSITE2FD,
196 F_VSITE3,
197 F_VSITE3FD,
198 F_VSITE3FAD,
199 F_VSITE3OUT,
200 F_VSITE4FD,
201 F_VSITE4FDN,
202 F_VSITEN,
203 F_COM_PULL,
204 F_DENSITYFITTING,
205 F_EQM,
206 F_EPOT,
207 F_EKIN,
208 F_ETOT,
209 F_ECONSERVED,
210 F_TEMP,
211 F_VTEMP_NOLONGERUSED,
212 F_PDISPCORR,
213 F_PRES,
214 F_DVDL_CONSTR,
215 F_DVDL,
216 F_DKDL,
217 F_DVDL_COUL,
218 F_DVDL_VDW,
219 F_DVDL_BONDED,
220 F_DVDL_RESTRAINT,
221 F_DVDL_TEMPERATURE, /* not calculated for now, but should just be the energy (NVT) or enthalpy (NPT), or 0 (NVE) */
222 F_NRE /* This number is for the total number of energies */
223 };
224
IS_RESTRAINT_TYPE(int ifunc)225 static inline bool IS_RESTRAINT_TYPE(int ifunc)
226 {
227 return ifunc == F_POSRES || ifunc == F_FBPOSRES || ifunc == F_DISRES || ifunc == F_RESTRBONDS
228 || ifunc == F_DISRESVIOL || ifunc == F_ORIRES || ifunc == F_ORIRESDEV
229 || ifunc == F_ANGRES || ifunc == F_ANGRESZ || ifunc == F_DIHRES;
230 }
231
232 /* Maximum allowed number of atoms, parameters and terms in interaction_function.
233 * Check kernel/toppush.c when you change these numbers.
234 */
235 constexpr int MAXATOMLIST = 6;
236 constexpr int MAXFORCEPARAM = 12;
237 constexpr int NR_RBDIHS = 6;
238 constexpr int NR_CBTDIHS = 6;
239 constexpr int NR_FOURDIHS = 4;
240
241 extern const t_interaction_function interaction_function[F_NRE];
242 /* initialised interaction functions descriptor */
243
244 #endif
245