Add C++ version of t_ilist
[alexxy/gromacs.git] / src / gromacs / topology / topology.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) 2011,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_TOPOLOGY_H
38 #define GMX_TOPOLOGY_TOPOLOGY_H
39
40 #include <cstdio>
41
42 #include <vector>
43
44 #include "gromacs/math/vectypes.h"
45 #include "gromacs/topology/atoms.h"
46 #include "gromacs/topology/block.h"
47 #include "gromacs/topology/idef.h"
48 #include "gromacs/topology/symtab.h"
49 #include "gromacs/utility/unique_cptr.h"
50
51 enum
52 {
53     egcTC,    egcENER,   egcACC, egcFREEZE,
54     egcUser1, egcUser2,  egcVCM, egcCompressedX,
55     egcORFIT, egcQMMM,
56     egcNR
57 };
58 /* Names corresponding to groups */
59 extern const char *gtypes[egcNR+1];
60
61 /*! \brief Molecules type data: atoms, interactions and exclusions */
62 struct gmx_moltype_t
63 {
64     gmx_moltype_t();
65
66     ~gmx_moltype_t();
67
68     /*! \brief Deleted copy assignment operator to avoid (not) freeing pointers */
69     gmx_moltype_t &operator=(const gmx_moltype_t &) = delete;
70
71     /*! \brief Default copy constructor */
72     gmx_moltype_t(const gmx_moltype_t &) = default;
73
74     char          **name;         /**< Name of the molecule type            */
75     t_atoms         atoms;        /**< The atoms in this molecule           */
76     InteractionList ilist[F_NRE]; /**< Interaction list with local indices  */
77     t_block         cgs;          /**< The charge groups                    */
78     t_blocka        excls;        /**< The exclusions                       */
79
80     /* TODO: Change ilist to InteractionLists */
81 };
82
83 /*! \brief Block of molecules of the same type, used in gmx_mtop_t */
84 struct gmx_molblock_t
85 {
86     int                    type = -1; /**< The molecule type index in mtop.moltype  */
87     int                    nmol = 0;  /**< The number of molecules in this block    */
88     std::vector<gmx::RVec> posres_xA; /**< Position restraint coordinates for top A */
89     std::vector<gmx::RVec> posres_xB; /**< Position restraint coordinates for top B */
90 };
91
92 /*! \brief Indices for a gmx_molblock_t, derived from other gmx_mtop_t contents */
93 struct MoleculeBlockIndices
94 {
95     int     numAtomsPerMolecule; /**< Number of atoms in a molecule in the block */
96     int     globalAtomStart;     /**< Global atom index of the first atom in the block */
97     int     globalAtomEnd;       /**< Global atom index + 1 of the last atom in the block */
98     int     globalResidueStart;  /**< Global residue index of the first residue in the block */
99     int     residueNumberStart;  /**< Residue numbers start from this value if the number of residues per molecule is <= maxres_renum */
100     int     moleculeIndexStart;  /**< Global molecule indexing starts from this value */
101 };
102
103 typedef struct gmx_groups_t
104 {
105     t_grps            grps[egcNR];  /* Groups of things                     */
106     int               ngrpname;     /* Number of groupnames                 */
107     char           ***grpname;      /* Names of the groups                  */
108     int               ngrpnr[egcNR];
109     unsigned char    *grpnr[egcNR]; /* Group numbers or NULL                */
110 } gmx_groups_t;
111
112 /*! \brief
113  * Returns group number of an input group for a given atom.
114  *
115  * Returns the group \p type for \p atom in \p group, or 0 if the
116  * entries for all atoms in the group are 0 and the pointer is thus null.
117  *
118  * \param[in] group Group to check.
119  * \param[in] type  Type of group to check.
120  * \param[in] atom  Atom to check if it has an entry.
121  */
122 int getGroupType (const gmx_groups_t *group, int type, int atom);
123
124 /* The global, complete system topology struct, based on molecule types.
125  * This structure should contain no data that is O(natoms) in memory.
126  *
127  * TODO: Find a solution for ensuring that the derived data is in sync
128  *       with the primary data, possibly by converting to a class.
129  */
130 struct gmx_mtop_t //NOLINT(clang-analyzer-optin.performance.Padding)
131 {
132     gmx_mtop_t();
133
134     ~gmx_mtop_t();
135
136     char                            **name; /* Name of the topology                 */
137     gmx_ffparams_t                    ffparams;
138     std::vector<gmx_moltype_t>        moltype;
139     std::vector<gmx_molblock_t>       molblock;
140     gmx_bool                          bIntermolecularInteractions; /* Are there intermolecular
141                                                                     * interactions?            */
142     std::unique_ptr<InteractionLists> intermolecular_ilist;        /* List of intermolecular interactions
143                                                                     * using system wide atom indices,
144                                                                     * either NULL or size F_NRE           */
145     int              natoms;
146     int              maxres_renum;                                 /* Parameter for residue numbering      */
147     int              maxresnr;                                     /* The maximum residue number in moltype */
148     t_atomtypes      atomtypes;                                    /* Atomtype properties                  */
149     gmx_groups_t     groups;                                       /* Groups of atoms for different purposes */
150     t_symtab         symtab;                                       /* The symbol table                     */
151     bool             haveMoleculeIndices;                          /* Tells whether we have valid molecule indices */
152
153     /* Derived data */
154     std::vector<MoleculeBlockIndices> moleculeBlockIndices;  /* Indices for each molblock entry for fast lookup of atom properties */
155 };
156
157 /* The mdrun node-local topology struct, completely written out */
158 typedef struct gmx_localtop_t
159 {
160     t_idef        idef;         /* The interaction function definition  */
161     t_atomtypes   atomtypes;    /* Atomtype properties                  */
162     t_block       cgs;          /* The charge groups                    */
163     t_blocka      excls;        /* The exclusions                       */
164 } gmx_localtop_t;
165
166 /* The old topology struct, completely written out, used in analysis tools */
167 typedef struct t_topology
168 {
169     char          **name;                        /* Name of the topology                 */
170     t_idef          idef;                        /* The interaction function definition  */
171     t_atoms         atoms;                       /* The atoms                            */
172     t_atomtypes     atomtypes;                   /* Atomtype properties                  */
173     t_block         cgs;                         /* The charge groups                    */
174     t_block         mols;                        /* The molecules                        */
175     gmx_bool        bIntermolecularInteractions; /* Inter.mol. int. ?   */
176     t_blocka        excls;                       /* The exclusions                       */
177     t_symtab        symtab;                      /* The symbol table                     */
178 } t_topology;
179
180 void init_mtop(gmx_mtop_t *mtop);
181 void init_top(t_topology *top);
182 void done_gmx_groups_t(gmx_groups_t *g);
183 void done_top(t_topology *top);
184 // Frees both t_topology and gmx_mtop_t when the former has been created from
185 // the latter.
186 void done_top_mtop(t_topology *top, gmx_mtop_t *mtop);
187 void done_localtop(gmx_localtop_t *top);
188 void done_and_sfree_localtop(gmx_localtop_t *top);
189
190 bool gmx_mtop_has_masses(const gmx_mtop_t *mtop);
191 bool gmx_mtop_has_charges(const gmx_mtop_t *mtop);
192 bool gmx_mtop_has_perturbed_charges(const gmx_mtop_t &mtop);
193 bool gmx_mtop_has_atomtypes(const gmx_mtop_t *mtop);
194 bool gmx_mtop_has_pdbinfo(const gmx_mtop_t *mtop);
195
196 void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop,
197              gmx_bool bShowNumbers, gmx_bool bShowParameters);
198 void pr_top(FILE *fp, int indent, const char *title, const t_topology *top,
199             gmx_bool bShowNumbers, gmx_bool bShowParameters);
200
201 void cmp_top(FILE *fp, const t_topology *t1, const t_topology *t2, real ftol, real abstol);
202 void cmp_groups(FILE *fp, const gmx_groups_t *g0, const gmx_groups_t *g1,
203                 int natoms0, int natoms1);
204
205 //! Deleter for gmx_localtop_t, needed until it has a proper destructor.
206 using ExpandedTopologyPtr = gmx::unique_cptr<gmx_localtop_t, done_and_sfree_localtop>;
207
208 #endif