1 /* 2 * This file is part of the GROMACS molecular simulation package. 3 * 4 * Copyright (c) 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 37 * Implements a force calculator based on GROMACS data structures. 38 * 39 * Intended for internal use inside the ForceCalculator. 40 * 41 * \author Victor Holanda <victor.holanda@cscs.ch> 42 * \author Joe Jordan <ejjordan@kth.se> 43 * \author Prashanth Kanduri <kanduri@cscs.ch> 44 * \author Sebastian Keller <keller@cscs.ch> 45 */ 46 47 #ifndef NBLIB_GMXCALCULATOR_H 48 #define NBLIB_GMXCALCULATOR_H 49 50 #include <memory> 51 52 #include "nblib/vector.h" 53 54 struct nonbonded_verlet_t; 55 struct t_forcerec; 56 struct t_nrnb; 57 struct interaction_const_t; 58 struct gmx_enerdata_t; 59 60 namespace gmx 61 { 62 template<typename T> 63 class ArrayRef; 64 class StepWorkload; 65 } // namespace gmx 66 67 namespace nblib 68 { 69 class Box; 70 class NbvSetupUtil; 71 class SimulationState; 72 struct NBKernelOptions; 73 74 /*! \brief GROMACS non-bonded force calculation backend 75 * 76 * This class encapsulates the various GROMACS data structures and their interplay 77 * from the NBLIB user. The class is a private member of the ForceCalculator and 78 * is not intended for the public interface. 79 * 80 * Handles the task of storing the simulation problem description using the internal 81 * representation used within GROMACS. It currently supports short range non-bonded 82 * interactions (PP) on a single node CPU. 83 * 84 */ 85 86 class GmxForceCalculator final 87 { 88 public: 89 GmxForceCalculator(); 90 91 ~GmxForceCalculator(); 92 93 //! Compute forces and return 94 void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput, gmx::ArrayRef<gmx::RVec> forceOutput); 95 96 //! Puts particles on a grid based on bounds specified by the box (for every NS step) 97 void setParticlesOnGrid(gmx::ArrayRef<const int> particleInfoAllVdw, 98 gmx::ArrayRef<const gmx::RVec> coordinates, 99 const Box& box); 100 101 private: 102 //! Friend to allow setting up private members in this class 103 friend class NbvSetupUtil; 104 105 //! Non-Bonded Verlet object for force calculation 106 std::unique_ptr<nonbonded_verlet_t> nbv_; 107 108 //! Only nbfp and shift_vec are used 109 std::unique_ptr<t_forcerec> forcerec_; 110 111 //! Parameters for various interactions in the system 112 std::unique_ptr<interaction_const_t> interactionConst_; 113 114 //! Tasks to perform in an MD Step 115 std::unique_ptr<gmx::StepWorkload> stepWork_; 116 117 //! Energies of different interaction types; currently only needed as an argument for dispatchNonbondedKernel 118 std::unique_ptr<gmx_enerdata_t> enerd_; 119 120 //! Non-bonded flop counter; currently only needed as an argument for dispatchNonbondedKernel 121 std::unique_ptr<t_nrnb> nrnb_; 122 123 //! Legacy matrix for box 124 matrix box_{ { 0 } }; 125 }; 126 127 } // namespace nblib 128 129 #endif // NBLIB_GMXCALCULATOR_H 130