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.
37 * Declares input generator class for CP2K QMMM
39 * \author Dmitry Morozov <dmitry.morozov@jyu.fi>
40 * \ingroup module_applied_forces
42 #ifndef GMX_APPLIED_FORCES_QMMMINPUTGENERATOR_H
43 #define GMX_APPLIED_FORCES_QMMMINPUTGENERATOR_H
47 #include "gromacs/pbcutil/pbc.h"
48 #include "gromacs/utility/arrayref.h"
50 #include "qmmmtypes.h"
56 * Class that takes QMMMParameters, Coordinates, Point charges, Box dimensions, pbcType.
57 * Generates QM/MM sample input parameters and pdb-style coordinates for CP2K.
58 * Input are generated as std::string objects which can be stored in tpr KVT
59 * and/or flushed into the files.
61 class QMMMInputGenerator
64 /*! \brief Construct QMMMInputGenerator from its parameters
66 * \param[in] parameters Structure with QMMM parameters
67 * \param[in] pbcType Periodic boundary conditions
68 * \param[in] box Matrix with full box of the system
69 * \param[in] q Point charges of each atom in the system
70 * \param[in] x Coordinates of each atom in the system
72 QMMMInputGenerator(const QMMMParameters& parameters,
75 ArrayRef<const real> q,
76 ArrayRef<const RVec> x);
78 /*! \brief Generates sample CP2K input file
81 std::string generateCP2KInput() const;
83 /*! \brief Generates PDB file suitable for usage with CP2K.
84 * In that PDB file Point Charges of MM atoms are provided with Extended Beta field
86 std::string generateCP2KPdb() const;
88 //! \brief Returns computed QM box dimensions
89 const matrix& qmBox() const;
91 //! \brief Returns computed translation vector in order to center QM atoms inside QM box
92 const RVec& qmTrans() const;
95 //! \brief Check if atom belongs to the global index of qmAtoms_
96 bool isQMAtom(index globalAtomIndex) const;
98 /*!\brief Calculates dimensions and center of the QM box.
99 * Also evaluates translation for the system in order to center QM atoms inside QM box
101 * \param[in] scale Factor of how much QM box would be bigger than the radius of QM system
102 * \param[in] minNorm Minimum norm of the QM box vector
104 void computeQMBox(real scale, real minNorm);
106 //! \brief Generates &GLOBAL section of CP2K Input
107 static std::string generateGlobalSection();
109 //! \brief Generates &DFT section of CP2K Input
110 std::string generateDFTSection() const;
112 //! \brief Generates &QMMM section of CP2K Input
113 std::string generateQMMMSection() const;
115 //! \brief Generates &MM section of CP2K Input
116 static std::string generateMMSection();
118 //! \brief Generates &SUBSYS section of CP2K Input
119 std::string generateSubsysSection() const;
121 //! QMMM Parameters structure
122 const QMMMParameters& parameters_;
123 //! Simulation PbcType
129 //! PBC-aware center of QM subsystem
131 //! Translation that shifts qmCenter_ to the center of qmBox_
133 //! Set containing indexes of all QM atoms
134 std::set<index> qmAtoms_;
135 //! Atoms point charges
136 ArrayRef<const real> q_;
137 //! Atoms coordinates
138 ArrayRef<const RVec> x_;
140 //! Scale of the generated QM box with respect to the QM-subsystem size
141 static constexpr real sc_qmBoxScale = 1.5;
142 //! Minimum length of the generated QM box vectors in nm
143 static constexpr real sc_qmBoxMinLength = 1.0;
146 /*! \brief Transforms vector a such as distance from it to the plane defined by vectors b and c
147 * will be h minimum length will be milL and maximum length maxL
149 * \param[in] a Vector which should be scaled
150 * \param[in] b First vector that forms the plane
151 * \param[in] c Second vector that forms the plane
152 * \param[in] h Distance from the end of a to the plane of (b,c)
153 * \param[in] minNorm Minimum norm of vector
154 * \param[in] maxNorm Maximum norm of vector
156 RVec computeQMBoxVec(const RVec& a, const RVec& b, const RVec& c, real h, real minNorm, real maxNorm);
160 #endif // GMX_APPLIED_FORCES_QMMMINPUTGENERATOR_H