1 /* 2 * This file is part of the GROMACS molecular simulation package. 3 * 4 * Copyright (c) 2018,2019,2020, 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. 8 * 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. 13 * 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. 18 * 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. 23 * 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. 31 * 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. 34 */ 35 /*! \libinternal \file 36 * \brief Declares interface to SHAKE code. 37 * 38 * \author David van der Spoel <david.vanderspoel@icm.uu.se> 39 * \author Berk Hess <hess@kth.se> 40 * \author Mark Abraham <mark.j.abraham@gmail.com> 41 * \ingroup module_mdlib 42 * \inlibraryapi 43 */ 44 45 #ifndef GMX_MDLIB_SHAKE_H 46 #define GMX_MDLIB_SHAKE_H 47 48 #include "gromacs/math/vec.h" 49 #include "gromacs/topology/block.h" 50 #include "gromacs/utility/real.h" 51 52 struct InteractionList; 53 class InteractionDefinitions; 54 struct t_inputrec; 55 struct t_nrnb; 56 struct t_pbc; 57 58 namespace gmx 59 { 60 template<typename T> 61 class ArrayRef; 62 63 enum class ConstraintVariable : int; 64 65 /*! \libinternal 66 * \brief Working data for the SHAKE algorithm 67 */ 68 struct shakedata 69 { 70 //! Returns the number of SHAKE blocks */ numShakeBlocksshakedata71 int numShakeBlocks() const { return sblock.size() - 1; } 72 73 //! The reference constraint vectors 74 std::vector<RVec> rij; 75 //! The reduced mass of the two atoms in each constraint times 0.5 76 std::vector<real> half_of_reduced_mass; 77 //! Multiplicative tolerance on the difference in the square of the constrained distance 78 std::vector<real> distance_squared_tolerance; 79 //! The reference constraint distances squared 80 std::vector<real> constraint_distance_squared; 81 /* SOR stuff */ 82 //! SOR delta 83 real delta = 0.1; 84 //! SOR omega 85 real omega = 1.0; 86 //! SOR gamma 87 real gamma = 1000000; 88 //! The SHAKE blocks, block i contains constraints sblock[i]/3 to sblock[i+1]/3 */ 89 std::vector<int> sblock = { 0 }; 90 /*! \brief Scaled Lagrange multiplier for each constraint. 91 * 92 * Value is -2 * eta from p. 336 of the paper, divided by the 93 * constraint distance. */ 94 std::vector<real> scaled_lagrange_multiplier; 95 }; 96 97 //! Make SHAKE blocks when not using DD. 98 void make_shake_sblock_serial(shakedata* shaked, InteractionDefinitions* idef, int numAtoms); 99 100 //! Make SHAKE blocks when using DD. 101 void make_shake_sblock_dd(shakedata* shaked, const InteractionList& ilcon); 102 103 /*! \brief Shake all the atoms blockwise. It is assumed that all the constraints 104 * in the idef->shakes field are sorted, to ascending block nr. The 105 * sblock array points into the idef->shakes.iatoms field, with block 0 106 * starting 107 * at sblock[0] and running to ( < ) sblock[1], block n running from 108 * sblock[n] to sblock[n+1]. Array sblock should be large enough. 109 * Return TRUE when OK, FALSE when shake-error 110 */ 111 bool constrain_shake(FILE* log, /* Log file */ 112 shakedata* shaked, /* Total number of atoms */ 113 const real invmass[], /* Atomic masses */ 114 const InteractionDefinitions& idef, /* The interaction def */ 115 const t_inputrec& ir, /* Input record */ 116 ArrayRef<const RVec> x_s, /* Coords before update */ 117 ArrayRef<RVec> xprime, /* Output coords when constraining x */ 118 ArrayRef<RVec> vprime, /* Output coords when constraining v */ 119 const t_pbc* pbc, /* PBC information */ 120 t_nrnb* nrnb, /* Performance measure */ 121 real lambda, /* FEP lambda */ 122 real* dvdlambda, /* FEP force */ 123 real invdt, /* 1/delta_t */ 124 ArrayRef<RVec> v, /* Also constrain v if not empty */ 125 bool bCalcVir, /* Calculate r x m delta_r */ 126 tensor vir_r_m_dr, /* sum r x m delta_r */ 127 bool bDumpOnError, /* Dump debugging stuff on error*/ 128 ConstraintVariable econq); /* which type of constraint is occurring */ 129 130 /*! \brief Regular iterative shake */ 131 void cshake(const int iatom[], 132 int ncon, 133 int* nnit, 134 int maxnit, 135 ArrayRef<const real> dist2, 136 ArrayRef<RVec> xp, 137 const t_pbc* pbc, 138 ArrayRef<const RVec> rij, 139 ArrayRef<const real> half_of_reduced_mass, 140 real omega, 141 const real invmass[], 142 ArrayRef<const real> distance_squared_tolerance, 143 ArrayRef<real> scaled_lagrange_multiplier, 144 int* nerror); 145 146 } // namespace gmx 147 148 #endif 149