Merge release-2018 into master
[alexxy/gromacs.git] / src / gromacs / topology / idef.h
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) 2013,2014,2015,2016,2018, by the GROMACS development team, led by
7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8  * and including many others, as listed in the AUTHORS file in the
9  * top-level source directory and at http://www.gromacs.org.
10  *
11  * GROMACS is free software; you can redistribute it and/or
12  * modify it under the terms of the GNU Lesser General Public License
13  * as published by the Free Software Foundation; either version 2.1
14  * of the License, or (at your option) any later version.
15  *
16  * GROMACS is distributed in the hope that it will be useful,
17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19  * Lesser General Public License for more details.
20  *
21  * You should have received a copy of the GNU Lesser General Public
22  * License along with GROMACS; if not, see
23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
25  *
26  * If you want to redistribute modifications to GROMACS, please
27  * consider that scientific software is very special. Version
28  * control is crucial - bugs must be traceable. We will be happy to
29  * consider code for inclusion in the official distribution, but
30  * derived work must not be called official GROMACS. Details are found
31  * in the README & COPYING files - if they are missing, get the
32  * official version at http://www.gromacs.org.
33  *
34  * To help us fund GROMACS development, we humbly ask that you cite
35  * the research papers on the package. Check out http://www.gromacs.org.
36  */
37 #ifndef GMX_TOPOLOGY_IDEF_H
38 #define GMX_TOPOLOGY_IDEF_H
39
40 #include <stdio.h>
41
42 #include "gromacs/math/vectypes.h"
43 #include "gromacs/utility/basedefinitions.h"
44 #include "gromacs/utility/real.h"
45
46 /* check kernel/toppush.c when you change these numbers */
47 #define MAXATOMLIST 6
48 #define MAXFORCEPARAM   12
49 #define NR_RBDIHS   6
50 #define NR_CBTDIHS   6
51 #define NR_FOURDIHS     4
52
53 typedef int t_iatom;
54
55 /* this MUST correspond to the
56    t_interaction_function[F_NRE] in gmxlib/ifunc.c */
57 enum {
58     F_BONDS,
59     F_G96BONDS,
60     F_MORSE,
61     F_CUBICBONDS,
62     F_CONNBONDS,
63     F_HARMONIC,
64     F_FENEBONDS,
65     F_TABBONDS,
66     F_TABBONDSNC,
67     F_RESTRBONDS,
68     F_ANGLES,
69     F_G96ANGLES,
70     F_RESTRANGLES,
71     F_LINEAR_ANGLES,
72     F_CROSS_BOND_BONDS,
73     F_CROSS_BOND_ANGLES,
74     F_UREY_BRADLEY,
75     F_QUARTIC_ANGLES,
76     F_TABANGLES,
77     F_PDIHS,
78     F_RBDIHS,
79     F_RESTRDIHS,
80     F_CBTDIHS,
81     F_FOURDIHS,
82     F_IDIHS,
83     F_PIDIHS,
84     F_TABDIHS,
85     F_CMAP,
86     F_GB12_NOLONGERUSED,
87     F_GB13_NOLONGERUSED,
88     F_GB14_NOLONGERUSED,
89     F_GBPOL_NOLONGERUSED,
90     F_NPSOLVATION_NOLONGERUSED,
91     F_LJ14,
92     F_COUL14,
93     F_LJC14_Q,
94     F_LJC_PAIRS_NB,
95     F_LJ,
96     F_BHAM,
97     F_LJ_LR_NOLONGERUSED,
98     F_BHAM_LR_NOLONGERUSED,
99     F_DISPCORR,
100     F_COUL_SR,
101     F_COUL_LR_NOLONGERUSED,
102     F_RF_EXCL,
103     F_COUL_RECIP,
104     F_LJ_RECIP,
105     F_DPD,
106     F_POLARIZATION,
107     F_WATER_POL,
108     F_THOLE_POL,
109     F_ANHARM_POL,
110     F_POSRES,
111     F_FBPOSRES,
112     F_DISRES,
113     F_DISRESVIOL,
114     F_ORIRES,
115     F_ORIRESDEV,
116     F_ANGRES,
117     F_ANGRESZ,
118     F_DIHRES,
119     F_DIHRESVIOL,
120     F_CONSTR,
121     F_CONSTRNC,
122     F_SETTLE,
123     F_VSITE2,
124     F_VSITE3,
125     F_VSITE3FD,
126     F_VSITE3FAD,
127     F_VSITE3OUT,
128     F_VSITE4FD,
129     F_VSITE4FDN,
130     F_VSITEN,
131     F_COM_PULL,
132     F_EQM,
133     F_EPOT,
134     F_EKIN,
135     F_ETOT,
136     F_ECONSERVED,
137     F_TEMP,
138     F_VTEMP_NOLONGERUSED,
139     F_PDISPCORR,
140     F_PRES,
141     F_DVDL_CONSTR,
142     F_DVDL,
143     F_DKDL,
144     F_DVDL_COUL,
145     F_DVDL_VDW,
146     F_DVDL_BONDED,
147     F_DVDL_RESTRAINT,
148     F_DVDL_TEMPERATURE, /* not calculated for now, but should just be the energy (NVT) or enthalpy (NPT), or 0 (NVE) */
149     F_NRE               /* This number is for the total number of energies      */
150 };
151
152 #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)))
153
154 typedef union t_iparams
155 {
156     /* Some parameters have A and B values for free energy calculations.
157      * The B values are not used for regular simulations of course.
158      * Free Energy for nonbondeds can be computed by changing the atom type.
159      * The harmonic type is used for all harmonic potentials:
160      * bonds, angles and improper dihedrals
161      */
162     struct {
163         real a, b, c;
164     } bham;
165     struct {
166         real rA, krA, rB, krB;
167     } harmonic;
168     struct {
169         real klinA, aA, klinB, aB;
170     } linangle;
171     struct {
172         real lowA, up1A, up2A, kA, lowB, up1B, up2B, kB;
173     } restraint;
174     /* No free energy supported for cubic bonds, FENE, WPOL or cross terms */
175     struct {
176         real b0, kb, kcub;
177     } cubic;
178     struct {
179         real bm, kb;
180     } fene;
181     struct {
182         real r1e, r2e, krr;
183     } cross_bb;
184     struct {
185         real r1e, r2e, r3e, krt;
186     } cross_ba;
187     struct {
188         real thetaA, kthetaA, r13A, kUBA, thetaB, kthetaB, r13B, kUBB;
189     } u_b;
190     struct {
191         real theta, c[5];
192     } qangle;
193     struct {
194         real alpha;
195     } polarize;
196     struct {
197         real alpha, drcut, khyp;
198     } anharm_polarize;
199     struct {
200         real al_x, al_y, al_z, rOH, rHH, rOD;
201     } wpol;
202     struct {
203         real a, alpha1, alpha2, rfac;
204     } thole;
205     struct {
206         real c6, c12;
207     } lj;
208     struct {
209         real c6A, c12A, c6B, c12B;
210     } lj14;
211     struct {
212         real fqq, qi, qj, c6, c12;
213     } ljc14;
214     struct {
215         real qi, qj, c6, c12;
216     } ljcnb;
217     /* Proper dihedrals can not have different multiplicity when
218      * doing free energy calculations, because the potential would not
219      * be periodic anymore.
220      */
221     struct {
222         real phiA, cpA; int mult; real phiB, cpB;
223     } pdihs;
224     struct {
225         real dA, dB;
226     } constr;
227     /* Settle can not be used for Free energy calculations of water bond geometry.
228      * Use shake (or lincs) instead if you have to change the water bonds.
229      */
230     struct {
231         real doh, dhh;
232     } settle;
233     struct {
234         real b0A, cbA, betaA, b0B, cbB, betaB;
235     } morse;
236     struct {
237         real pos0A[DIM], fcA[DIM], pos0B[DIM], fcB[DIM];
238     } posres;
239     struct {
240         real pos0[DIM], r, k; int geom;
241     } fbposres;
242     struct {
243         real rbcA[NR_RBDIHS], rbcB[NR_RBDIHS];
244     } rbdihs;
245     struct {
246         real cbtcA[NR_CBTDIHS], cbtcB[NR_CBTDIHS];
247     } cbtdihs;
248     struct {
249         real a, b, c, d, e, f;
250     } vsite;
251     struct {
252         int  n; real a;
253     } vsiten;
254     /* NOTE: npair is only set after reading the tpx file */
255     struct {
256         real low, up1, up2, kfac; int type, label, npair;
257     } disres;
258     struct {
259         real phiA, dphiA, kfacA, phiB, dphiB, kfacB;
260     } dihres;
261     struct {
262         int  ex, power, label; real c, obs, kfac;
263     } orires;
264     struct {
265         int  table; real kA; real kB;
266     } tab;
267     struct {
268         int cmapA, cmapB;
269     } cmap;
270     struct {
271         real buf[MAXFORCEPARAM];
272     } generic;                                               /* Conversion */
273 } t_iparams;
274
275 typedef int t_functype;
276
277 /*
278  * The nonperturbed/perturbed interactions are now separated (sorted) in the
279  * ilist, such that the first 0..(nr_nonperturbed-1) ones are exactly that, and
280  * the remaining ones from nr_nonperturbed..(nr-1) are perturbed bonded
281  * interactions.
282  */
283 typedef struct t_ilist
284 {
285     int      nr;
286     int      nr_nonperturbed;
287     t_iatom *iatoms;
288     int      nalloc;
289 } t_ilist;
290
291 /*
292  * The struct t_ilist defines a list of atoms with their interactions.
293  * General field description:
294  *   int nr
295  *      the size (nr elements) of the interactions array (iatoms[]).
296  *   t_iatom *iatoms
297  *  specifies which atoms are involved in an interaction of a certain
298  *       type. The layout of this array is as follows:
299  *
300  *        +-----+---+---+---+-----+---+---+-----+---+---+---+-----+---+---+...
301  *        |type1|at1|at2|at3|type2|at1|at2|type1|at1|at2|at3|type3|at1|at2|
302  *        +-----+---+---+---+-----+---+---+-----+---+---+---+-----+---+---+...
303  *
304  *  So for interaction type type1 3 atoms are needed, and for type2 and
305  *      type3 only 2. The type identifier is used to select the function to
306  *      calculate the interaction and its actual parameters. This type
307  *      identifier is an index in a params[] and functype[] array.
308  */
309
310 typedef struct
311 {
312     real *cmap; /* Has length 4*grid_spacing*grid_spacing, */
313     /* there are 4 entries for each cmap type (V,dVdx,dVdy,d2dVdxdy) */
314 } gmx_cmapdata_t;
315
316 typedef struct gmx_cmap_t
317 {
318     int             ngrid;        /* Number of allocated cmap (cmapdata_t ) grids */
319     int             grid_spacing; /* Grid spacing */
320     gmx_cmapdata_t *cmapdata;     /* Pointer to grid with actual, pre-interpolated data */
321 } gmx_cmap_t;
322
323
324 typedef struct gmx_ffparams_t
325 {
326     int         ntypes;
327     int         atnr;
328     t_functype *functype;
329     t_iparams  *iparams;
330     double      reppow;    /* The repulsion power for VdW: C12*r^-reppow   */
331     real        fudgeQQ;   /* The scaling factor for Coulomb 1-4: f*q1*q2  */
332     gmx_cmap_t  cmap_grid; /* The dihedral correction maps                 */
333 } gmx_ffparams_t;
334
335 enum {
336     ilsortUNKNOWN, ilsortNO_FE, ilsortFE_UNSORTED, ilsortFE_SORTED
337 };
338
339 typedef struct t_idef
340 {
341     int         ntypes;
342     int         atnr;
343     t_functype *functype;
344     t_iparams  *iparams;
345     real        fudgeQQ;
346     gmx_cmap_t  cmap_grid;
347     t_iparams  *iparams_posres, *iparams_fbposres;
348     int         iparams_posres_nalloc, iparams_fbposres_nalloc;
349
350     t_ilist     il[F_NRE];
351     int         ilsort;
352 } t_idef;
353
354 /*
355  * The struct t_idef defines all the interactions for the complete
356  * simulation. The structure is setup in such a way that the multinode
357  * version of the program  can use it as easy as the single node version.
358  * General field description:
359  *   int ntypes
360  *      defines the number of elements in functype[] and param[].
361  *   int nodeid
362  *      the node id (if parallel machines)
363  *   int atnr
364  *      the number of atomtypes
365  *   t_functype *functype
366  *      array of length ntypes, defines for every force type what type of
367  *      function to use. Every "bond" with the same function but different
368  *      force parameters is a different force type. The type identifier in the
369  *      forceatoms[] array is an index in this array.
370  *   t_iparams *iparams
371  *      array of length ntypes, defines the parameters for every interaction
372  *      type. The type identifier in the actual interaction list
373  *      (ilist[ftype].iatoms[]) is an index in this array.
374  *   gmx_cmap_t cmap_grid
375  *      the grid for the dihedral pair correction maps.
376  *   t_iparams *iparams_posres, *iparams_fbposres
377  *      defines the parameters for position restraints only.
378  *      Position restraints are the only interactions that have different
379  *      parameters (reference positions) for different molecules
380  *      of the same type. ilist[F_POSRES].iatoms[] is an index in this array.
381  *   t_ilist il[F_NRE]
382  *      The list of interactions for each type. Note that some,
383  *      such as LJ and COUL will have 0 entries.
384  *   int ilsort
385  *      The state of the sorting of il, values are provided above.
386  */
387
388 void pr_iparams(FILE *fp, t_functype ftype, const t_iparams *iparams);
389 void pr_ilist(FILE *fp, int indent, const char *title,
390               const t_functype *functype, const t_ilist *ilist,
391               gmx_bool bShowNumbers,
392               gmx_bool bShowParameters, const t_iparams *iparams);
393 void pr_ffparams(FILE *fp, int indent, const char *title,
394                  const gmx_ffparams_t *ffparams, gmx_bool bShowNumbers);
395 void pr_idef(FILE *fp, int indent, const char *title, const t_idef *idef,
396              gmx_bool bShowNumbers, gmx_bool bShowParameters);
397
398 #endif