Fix PME run mode without PME
[alexxy/gromacs.git] / src / gromacs / taskassignment / decidegpuusage.h
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 2017,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 functionality for deciding whether tasks will run on GPUs.
37  *
38  * \author Mark Abraham <mark.j.abraham@gmail.com>
39  * \ingroup module_taskassignment
40  * \inlibraryapi
41  */
42
43 #ifndef GMX_TASKASSIGNMENT_DECIDEGPUUSAGE_H
44 #define GMX_TASKASSIGNMENT_DECIDEGPUUSAGE_H
45
46 #include <vector>
47
48 struct gmx_hw_info_t;
49 struct gmx_mtop_t;
50 struct t_inputrec;
51 enum class PmeRunMode;
52
53 namespace gmx
54 {
55
56 class MDLogger;
57
58 //! Record where a compute task is targetted.
59 enum class TaskTarget : int
60 {
61     Auto,
62     Cpu,
63     Gpu
64 };
65
66 //! Help pass GPU-emulation parameters with type safety.
67 enum class EmulateGpuNonbonded : bool
68 {
69     //! Do not emulate GPUs.
70     No,
71     //! Do emulate GPUs.
72     Yes
73 };
74
75 class MDAtoms;
76
77 /*! \brief Decide whether this thread-MPI simulation will run
78  * nonbonded tasks on GPUs.
79  *
80  * The number of GPU tasks and devices influences both the choice of
81  * the number of ranks, and checks upon any such choice made by the
82  * user. So we need to consider this before any automated choice of
83  * the number of thread-MPI ranks.
84  *
85  * \param[in]  nonbondedTarget             The user's choice for mdrun -nb for where to assign short-ranged nonbonded interaction tasks.
86  * \param[in]  gpuIdsToUse                 The compatible GPUs that the user permitted us to use.
87  * \param[in]  userGpuTaskAssignment       The user-specified assignment of GPU tasks to device IDs.
88  * \param[in]  emulateGpuNonbonded         Whether we will emulate GPU calculation of nonbonded interactions.
89  * \param[in]  buildSupportsNonbondedOnGpu Whether GROMACS was built with GPU support.
90  * \param[in]  nonbondedOnGpuIsUseful      Whether computing nonbonded interactions on a GPU is useful for this calculation.
91  * \param[in]  numRanksPerSimulation       The number of ranks in each simulation.
92  *
93  * \returns    Whether the simulation will run nonbonded tasks on GPUs.
94  *
95  * \throws     std::bad_alloc          If out of memory
96  *             InconsistentInputError  If the user requirements are inconsistent. */
97 bool decideWhetherToUseGpusForNonbondedWithThreadMpi(TaskTarget              nonbondedTarget,
98                                                      const std::vector<int>& gpuIdsToUse,
99                                                      const std::vector<int>& userGpuTaskAssignment,
100                                                      EmulateGpuNonbonded     emulateGpuNonbonded,
101                                                      bool buildSupportsNonbondedOnGpu,
102                                                      bool nonbondedOnGpuIsUseful,
103                                                      int  numRanksPerSimulation);
104
105 /*! \brief Decide whether this thread-MPI simulation will run
106  * PME tasks on GPUs.
107  *
108  * The number of GPU tasks and devices influences both the choice of
109  * the number of ranks, and checks upon any such choice made by the
110  * user. So we need to consider this before any automated choice of
111  * the number of thread-MPI ranks.
112  *
113  * \param[in]  useGpuForNonbonded        Whether GPUs will be used for nonbonded interactions.
114  * \param[in]  pmeTarget                 The user's choice for mdrun -pme for where to assign
115  *                                       long-ranged PME nonbonded interaction tasks.
116  * \param[in]  gpuIdsToUse               The compatible GPUs that the user permitted us to use.
117  * \param[in]  userGpuTaskAssignment     The user-specified assignment of GPU tasks to device IDs.
118  * \param[in]  hardwareInfo              Hardware information
119  * \param[in]  inputrec                  The user input
120  * \param[in]  mtop                      Global system topology
121  * \param[in]  numRanksPerSimulation     The number of ranks in each simulation.
122  * \param[in]  numPmeRanksPerSimulation  The number of PME ranks in each simulation.
123  *
124  * \returns    Whether the simulation will run PME tasks on GPUs.
125  *
126  * \throws     std::bad_alloc          If out of memory
127  *             InconsistentInputError  If the user requirements are inconsistent. */
128 bool decideWhetherToUseGpusForPmeWithThreadMpi(bool                    useGpuForNonbonded,
129                                                TaskTarget              pmeTarget,
130                                                const std::vector<int>& gpuIdsToUse,
131                                                const std::vector<int>& userGpuTaskAssignment,
132                                                const gmx_hw_info_t&    hardwareInfo,
133                                                const t_inputrec&       inputrec,
134                                                const gmx_mtop_t&       mtop,
135                                                int                     numRanksPerSimulation,
136                                                int                     numPmeRanksPerSimulation);
137
138 /*! \brief Decide whether the simulation will try to run nonbonded
139  * tasks on GPUs.
140  *
141  * The final decision cannot be made until after the duty of the rank
142  * is known. But we need to know if nonbonded will run on GPUs for
143  * setting up DD (particularly rlist) and determining duty. If the
144  * user requires GPUs for the tasks of that duty, then it will be an
145  * error when none are found.
146  *
147  * With thread-MPI, calls have been made to
148  * decideWhetherToUseGpusForNonbondedWithThreadMpi() and
149  * decideWhetherToUseGpusForPmeWithThreadMpi() to help determine
150  * the number of ranks and run some checks, but the final
151  * decision is made in this routine, along with many more
152  * consistency checks.
153  *
154  * \param[in]  nonbondedTarget             The user's choice for mdrun -nb for where to assign short-ranged nonbonded interaction tasks.
155  * \param[in]  userGpuTaskAssignment       The user-specified assignment of GPU tasks to device IDs.
156  * \param[in]  emulateGpuNonbonded         Whether we will emulate GPU calculation of nonbonded interactions.
157  * \param[in]  buildSupportsNonbondedOnGpu Whether GROMACS was build with GPU support.
158  * \param[in]  nonbondedOnGpuIsUseful      Whether computing nonbonded interactions on a GPU is useful for this calculation.
159  * \param[in]  gpusWereDetected            Whether compatible GPUs were detected on any node.
160  *
161  * \returns    Whether the simulation will run nonbonded and PME tasks, respectively, on GPUs.
162  *
163  * \throws     std::bad_alloc          If out of memory
164  *             InconsistentInputError  If the user requirements are inconsistent. */
165 bool decideWhetherToUseGpusForNonbonded(TaskTarget              nonbondedTarget,
166                                         const std::vector<int>& userGpuTaskAssignment,
167                                         EmulateGpuNonbonded     emulateGpuNonbonded,
168                                         bool                    buildSupportsNonbondedOnGpu,
169                                         bool                    nonbondedOnGpuIsUseful,
170                                         bool                    gpusWereDetected);
171
172 /*! \brief Decide whether the simulation will try to run tasks of
173  * different types on GPUs.
174  *
175  * The final decision cannot be made until after the duty of the rank
176  * is known. But we need to know if nonbonded will run on GPUs for
177  * setting up DD (particularly rlist) and determining duty. If the
178  * user requires GPUs for the tasks of that duty, then it will be an
179  * error when none are found.
180  *
181  * With thread-MPI, calls have been made to
182  * decideWhetherToUseGpusForNonbondedWithThreadMpi() and
183  * decideWhetherToUseGpusForPmeWithThreadMpi() to help determine
184  * the number of ranks and run some checks, but the final
185  * decision is made in this routine, along with many more
186  * consistency checks.
187  *
188  * \param[in]  useGpuForNonbonded        Whether GPUs will be used for nonbonded interactions.
189  * \param[in]  pmeTarget                 The user's choice for mdrun -pme for where to assign long-ranged PME nonbonded interaction tasks.
190  * \param[in]  userGpuTaskAssignment     The user-specified assignment of GPU tasks to device IDs.
191  * \param[in]  hardwareInfo              Hardware information
192  * \param[in]  inputrec                  The user input
193  * \param[in]  mtop                      Global system topology
194  * \param[in]  numRanksPerSimulation     The number of ranks in each simulation.
195  * \param[in]  numPmeRanksPerSimulation  The number of PME ranks in each simulation.
196  * \param[in]  gpusWereDetected          Whether compatible GPUs were detected on any node.
197  *
198  * \returns    Whether the simulation will run nonbonded and PME tasks, respectively, on GPUs.
199  *
200  * \throws     std::bad_alloc          If out of memory
201  *             InconsistentInputError  If the user requirements are inconsistent. */
202 bool decideWhetherToUseGpusForPme(bool                    useGpuForNonbonded,
203                                   TaskTarget              pmeTarget,
204                                   const std::vector<int>& userGpuTaskAssignment,
205                                   const gmx_hw_info_t&    hardwareInfo,
206                                   const t_inputrec&       inputrec,
207                                   const gmx_mtop_t&       mtop,
208                                   int                     numRanksPerSimulation,
209                                   int                     numPmeRanksPerSimulation,
210                                   bool                    gpusWereDetected);
211
212 /*! \brief Determine PME run mode.
213  *
214  * Given the PME task assignment in \p useGpuForPme and the user-provided
215  * FFT task target in \p pmeFftTarget, returns a PME run mode for the
216  * current run. It also checks the compatibility of the two.
217  *
218  * \note Aborts the run upon incompatible values of \p useGpuForPme and \p pmeFftTarget.
219  *
220  * \param[in]  useGpuForPme              PME task assignment, true if PME task is mapped to the GPU.
221  * \param[in]  pmeFftTarget              The user's choice for -pmefft for where to assign the FFT
222  * work of the PME task. \param[in]  inputrec                  The user input record
223  * */
224 PmeRunMode determinePmeRunMode(bool useGpuForPme, const TaskTarget& pmeFftTarget, const t_inputrec& inputrec);
225
226 /*! \brief Decide whether the simulation will try to run bonded tasks on GPUs.
227  *
228  * \param[in]  useGpuForNonbonded        Whether GPUs will be used for nonbonded interactions.
229  * \param[in]  useGpuForPme              Whether GPUs will be used for PME interactions.
230  * \param[in]  bondedTarget              The user's choice for mdrun -bonded for where to assign tasks.
231  * \param[in]  canUseGpuForBonded        Whether the bonded interactions can run on a GPU
232  * \param[in]  usingLJPme                Whether Vdw interactions use LJ-PME.
233  * \param[in]  usingElecPmeOrEwald       Whether a PME or Ewald type method is used for electrostatics.
234  * \param[in]  numPmeRanksPerSimulation  The number of PME ranks in each simulation, can be -1 for auto.
235  * \param[in]  gpusWereDetected          Whether compatible GPUs were detected on any node.
236  *
237  * \returns    Whether the simulation will run bondeded tasks on GPUs.
238  *
239  * \throws     std::bad_alloc          If out of memory
240  *             InconsistentInputError  If the user requirements are inconsistent. */
241 bool decideWhetherToUseGpusForBonded(bool       useGpuForNonbonded,
242                                      bool       useGpuForPme,
243                                      TaskTarget bondedTarget,
244                                      bool       canUseGpuForBonded,
245                                      bool       usingLJPme,
246                                      bool       usingElecPmeOrEwald,
247                                      int        numPmeRanksPerSimulation,
248                                      bool       gpusWereDetected);
249
250 /*! \brief Decide whether to use GPU for update.
251  *
252  * \param[in]  forceGpuUpdateDefault        If update should run on GPU by default.
253  * \param[in]  isDomainDecomposition        Whether there more than one domain.
254  * \param[in]  useUpdateGroups              If the constraints can be split across domains.
255  * \param[in]  pmeRunMode                   PME running mode: CPU, GPU or mixed.
256  * \param[in]  havePmeOnlyRank              If there is a PME-only rank in the simulation.
257  * \param[in]  useGpuForNonbonded           Whether GPUs will be used for nonbonded interactions.
258  * \param[in]  updateTarget                 User choice for running simulation on GPU.
259  * \param[in]  gpusWereDetected             Whether compatible GPUs were detected on any node.
260  * \param[in]  inputrec                     The user input.
261  * \param[in]  mtop                         The global topology.
262  * \param[in]  useEssentialDynamics         If essential dynamics is active.
263  * \param[in]  doOrientationRestraints      If orientation restraints are enabled.
264  * \param[in]  useReplicaExchange           If this is a REMD simulation.
265  * \param[in]  doRerun                      It this is a rerun.
266  * \param[in]  mdlog                        MD logger.
267  *
268  * \returns    Whether complete simulation can be run on GPU.
269  * \throws     std::bad_alloc            If out of memory
270  *             InconsistentInputError    If the user requirements are inconsistent.
271  */
272 bool decideWhetherToUseGpuForUpdate(bool                 forceGpuUpdateDefault,
273                                     bool                 isDomainDecomposition,
274                                     bool                 useUpdateGroups,
275                                     PmeRunMode           pmeRunMode,
276                                     bool                 havePmeOnlyRank,
277                                     bool                 useGpuForNonbonded,
278                                     TaskTarget           updateTarget,
279                                     bool                 gpusWereDetected,
280                                     const t_inputrec&    inputrec,
281                                     const gmx_mtop_t&    mtop,
282                                     bool                 useEssentialDynamics,
283                                     bool                 doOrientationRestraints,
284                                     bool                 useReplicaExchange,
285                                     bool                 doRerun,
286                                     const gmx::MDLogger& mdlog);
287
288
289 } // namespace gmx
290
291 #endif