2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2021, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
35 /*! \libinternal \file
37 * \brief This file makes declarations used for building
38 * the reverse topology
41 * \ingroup module_domdec
44 #ifndef GMX_DOMDEC_REVERSETOPOLOGY_H
45 #define GMX_DOMDEC_REVERSETOPOLOGY_H
52 #include "gromacs/mdlib/vsite.h"
53 #include "gromacs/topology/idef.h"
54 #include "gromacs/utility/listoflists.h"
57 struct gmx_ffparams_t;
61 struct ReverseTopOptions;
65 class VirtualSitesHandler;
66 enum class DDBondedChecking : bool;
69 //! Options for linking atoms in make_reverse_ilist
70 enum class AtomLinkRule
72 FirstAtom, //!< Link all interactions to the first atom in the atom list
73 AllAtomsInBondeds //!< Link bonded interactions to all atoms involved, don't link vsites
76 struct MolblockIndices
84 struct reverse_ilist_t
86 std::vector<int> index; /* Index for each atom into il */
87 std::vector<int> il; /* ftype|type|a0|...|an|ftype|... */
88 int numAtomsInMolecule; /* The number of atoms in this molecule */
91 /*! \internal \brief Struct for thread local work data for local topology generation */
94 /*! \brief Constructor
96 * \param[in] ffparams The interaction parameters, the lifetime of the created object should not exceed the lifetime of the passed parameters
98 thread_work_t(const gmx_ffparams_t& ffparams) : idef(ffparams) {}
100 InteractionDefinitions idef; /**< Partial local topology */
101 std::unique_ptr<gmx::VsitePbc> vsitePbc = nullptr; /**< vsite PBC structure */
102 int numBondedInteractions = 0; /**< The number of bonded interactions observed */
103 gmx::ListOfLists<int> excl; /**< List of exclusions */
104 int excl_count = 0; /**< The total exclusion count for \p excl */
107 /*! \internal \brief Options for setting up gmx_reverse_top_t */
108 struct ReverseTopOptions
110 //! Constructor, constraints and settles are not including with a single argument
111 ReverseTopOptions(gmx::DDBondedChecking ddBondedChecking,
112 bool includeConstraints = false,
113 bool includeSettles = false) :
114 ddBondedChecking(ddBondedChecking),
115 includeConstraints(includeConstraints),
116 includeSettles(includeSettles)
120 //! \brief For which bonded interactions to check assignments
121 const gmx::DDBondedChecking ddBondedChecking;
122 //! \brief Whether constraints are stored in this reverse top
123 const bool includeConstraints;
124 //! \brief Whether settles are stored in this reverse top
125 const bool includeSettles;
128 //! \internal \brief Reverse topology class
129 class gmx_reverse_top_t
133 gmx_reverse_top_t(const gmx_mtop_t& mtop, bool useFreeEnergy, const ReverseTopOptions& reverseTopOptions);
135 ~gmx_reverse_top_t();
137 //! Gets the options that configured the construction
138 const ReverseTopOptions& options() const;
140 //! Gets the interaction list for the given molecule type
141 const reverse_ilist_t& interactionListForMoleculeType(int moleculeType) const;
143 //! Returns the total count of bonded interactions, used for checking partitioning
144 int expectedNumGlobalBondedInteractions() const;
146 //! Returns the molecule block indices
147 gmx::ArrayRef<const MolblockIndices> molblockIndices() const;
148 //! Returns whether the reverse topology describes intermolecular interactions
149 bool hasIntermolecularInteractions() const;
150 //! Gets the interaction list for any intermolecular interactions
151 const reverse_ilist_t& interactionListForIntermolecularInteractions() const;
152 //! Returns whether the reverse topology describes interatomic interactions
153 bool hasInterAtomicInteractions() const;
154 //! Returns whether there are interactions of type F_POSRES and/or F_FBPOSRES
155 bool hasPositionRestraints() const;
156 //! Returns the per-thread working structures for making the local topology
157 gmx::ArrayRef<thread_work_t> threadWorkObjects() const;
158 //! Returns whether the local topology interactions should be sorted
159 bool doSorting() const;
161 //! Private implementation definition
163 //! Private implementation declaration
164 std::unique_ptr<Impl> impl_;
167 /*! \brief Returns the number of atom entries for il in gmx_reverse_top_t */
168 int nral_rt(int ftype);
170 /*! \brief Return whether interactions of type \p ftype need to be assigned exactly once */
171 bool dd_check_ftype(int ftype, const ReverseTopOptions& rtOptions);
173 /*! \brief Return global topology molecule information for global atom index \p i_gl */
174 void global_atomnr_to_moltype_ind(gmx::ArrayRef<const MolblockIndices> molblockIndices,
181 /*! \brief Make the reverse ilist: a list of bonded interactions linked to atoms */
182 int make_reverse_ilist(const InteractionLists& ilist,
183 const t_atoms* atoms,
184 const ReverseTopOptions& rtOptions,
185 AtomLinkRule atomLinkRule,
186 reverse_ilist_t* ril_mt);
188 /*! \brief Generate and store the reverse topology */
189 void dd_make_reverse_top(FILE* fplog,
191 const gmx_mtop_t& mtop,
192 const gmx::VirtualSitesHandler* vsite,
193 const t_inputrec& inputrec,
194 gmx::DDBondedChecking ddBondedChecking);