2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2016,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.
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.
38 * \brief This file contains internal function definitions for performing the PME calculations on GPU.
39 * These are not meant to be exposed outside of the PME GPU code.
40 * As of now, their bodies are still in the common pme_gpu.cpp files.
42 * \author Aleksei Iupinov <a.yupinov@gmail.com>
43 * \ingroup module_ewald
46 #ifndef GMX_EWALD_PME_GPU_INTERNAL_H
47 #define GMX_EWALD_PME_GPU_INTERNAL_H
49 #include "gromacs/fft/fft.h" // for the gmx_fft_direction enum
50 #include "gromacs/gpu_utils/devicebuffer_datatype.h"
51 #include "gromacs/gpu_utils/gpu_macros.h" // for the GPU_FUNC_ macros
52 #include "gromacs/utility/arrayref.h"
54 #include "pme_gpu_types_host.h"
55 #include "pme_output.h"
57 class GpuEventSynchronizer;
58 struct DeviceInformation;
61 struct gmx_pme_t; // only used in pme_gpu_reinit
64 enum class PmeForceOutputHandling;
68 struct PmeGpuSettings;
76 //! Type of spline data
77 enum class PmeSplineDataType
80 Derivatives, // dtheta
81 }; // TODO move this into new and shiny pme.h (pme-types.h?)
83 //! PME grid dimension ordering (from major to minor)
84 enum class GridOrdering
90 /*! \libinternal \brief
91 * Returns the size of the block size requirement
93 * The GPU version of PME requires that the coordinates array have a
94 * size divisible by the returned number.
96 * \returns Number of atoms in a single GPU atom data chunk, which
97 * determines a minimum divisior of the size of the memory allocated.
99 int pme_gpu_get_atom_data_block_size();
101 /*! \libinternal \brief
102 * Returns the number of atoms per chunk in the atom spline theta/dtheta data layout.
104 * \param[in] pmeGpu The PME GPU structure.
105 * \returns Number of atoms in a single GPU atom spline data chunk.
107 int pme_gpu_get_atoms_per_warp(const PmeGpu* pmeGpu);
109 /*! \libinternal \brief
110 * Synchronizes the current computation, waiting for the GPU kernels/transfers to finish.
112 * \param[in] pmeGpu The PME GPU structure.
114 GPU_FUNC_QUALIFIER void pme_gpu_synchronize(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu)) GPU_FUNC_TERM;
116 /*! \libinternal \brief
117 * Allocates the fixed size energy and virial buffer both on GPU and CPU.
119 * \param[in,out] pmeGpu The PME GPU structure.
121 void pme_gpu_alloc_energy_virial(PmeGpu* pmeGpu);
123 /*! \libinternal \brief
124 * Frees the energy and virial memory both on GPU and CPU.
126 * \param[in] pmeGpu The PME GPU structure.
128 void pme_gpu_free_energy_virial(PmeGpu* pmeGpu);
130 /*! \libinternal \brief
131 * Clears the energy and virial memory on GPU with 0.
132 * Should be called at the end of PME computation which returned energy/virial.
134 * \param[in] pmeGpu The PME GPU structure.
136 void pme_gpu_clear_energy_virial(const PmeGpu* pmeGpu);
138 /*! \libinternal \brief
139 * Reallocates and copies the pre-computed B-spline values to the GPU.
141 * \param[in,out] pmeGpu The PME GPU structure.
143 void pme_gpu_realloc_and_copy_bspline_values(PmeGpu* pmeGpu);
145 /*! \libinternal \brief
146 * Frees the pre-computed B-spline values on the GPU (and the transfer CPU buffers).
148 * \param[in] pmeGpu The PME GPU structure.
150 void pme_gpu_free_bspline_values(const PmeGpu* pmeGpu);
152 /*! \libinternal \brief
153 * Reallocates the GPU buffer for the PME forces.
155 * \param[in] pmeGpu The PME GPU structure.
157 void pme_gpu_realloc_forces(PmeGpu* pmeGpu);
159 /*! \libinternal \brief
160 * Frees the GPU buffer for the PME forces.
162 * \param[in] pmeGpu The PME GPU structure.
164 void pme_gpu_free_forces(const PmeGpu* pmeGpu);
166 /*! \libinternal \brief
167 * Copies the forces from the CPU buffer to the GPU (to reduce them with the PME GPU gathered
168 * forces). To be called e.g. after the bonded calculations.
170 * \param[in] pmeGpu The PME GPU structure.
172 void pme_gpu_copy_input_forces(PmeGpu* pmeGpu);
174 /*! \libinternal \brief
175 * Copies the forces from the GPU to the CPU buffer. To be called after the gathering stage.
177 * \param[in] pmeGpu The PME GPU structure.
179 void pme_gpu_copy_output_forces(PmeGpu* pmeGpu);
181 /*! \libinternal \brief
182 * Checks whether work in the PME GPU stream has completed.
184 * \param[in] pmeGpu The PME GPU structure.
186 * \returns True if work in the PME stream has completed.
188 bool pme_gpu_stream_query(const PmeGpu* pmeGpu);
190 /*! \libinternal \brief
191 * Reallocates the buffer on the GPU and copies the charges/coefficients from the CPU buffer.
192 * Clears the padded part if needed.
194 * \param[in] pmeGpu The PME GPU structure.
195 * \param[in] h_coefficients The input atom charges/coefficients.
197 * Does not need to be done for every PME computation, only whenever the local charges change.
198 * (So, in the beginning of the run, or on DD step).
200 void pme_gpu_realloc_and_copy_input_coefficients(PmeGpu* pmeGpu, const float* h_coefficients);
202 /*! \libinternal \brief
203 * Frees the charges/coefficients on the GPU.
205 * \param[in] pmeGpu The PME GPU structure.
207 void pme_gpu_free_coefficients(const PmeGpu* pmeGpu);
209 /*! \libinternal \brief
210 * Reallocates the buffers on the GPU and the host for the atoms spline data.
212 * \param[in,out] pmeGpu The PME GPU structure.
214 void pme_gpu_realloc_spline_data(PmeGpu* pmeGpu);
216 /*! \libinternal \brief
217 * Frees the buffers on the GPU for the atoms spline data.
219 * \param[in] pmeGpu The PME GPU structure.
221 void pme_gpu_free_spline_data(const PmeGpu* pmeGpu);
223 /*! \libinternal \brief
224 * Reallocates the buffers on the GPU and the host for the particle gridline indices.
226 * \param[in,out] pmeGpu The PME GPU structure.
228 void pme_gpu_realloc_grid_indices(PmeGpu* pmeGpu);
230 /*! \libinternal \brief
231 * Frees the buffer on the GPU for the particle gridline indices.
233 * \param[in] pmeGpu The PME GPU structure.
235 void pme_gpu_free_grid_indices(const PmeGpu* pmeGpu);
237 /*! \libinternal \brief
238 * Reallocates the real space grid and the complex reciprocal grid (if needed) on the GPU.
240 * \param[in] pmeGpu The PME GPU structure.
242 void pme_gpu_realloc_grids(PmeGpu* pmeGpu);
244 /*! \libinternal \brief
245 * Frees the real space grid and the complex reciprocal grid (if needed) on the GPU.
247 * \param[in] pmeGpu The PME GPU structure.
249 void pme_gpu_free_grids(const PmeGpu* pmeGpu);
251 /*! \libinternal \brief
252 * Clears the real space grid on the GPU.
253 * Should be called at the end of each computation.
255 * \param[in] pmeGpu The PME GPU structure.
257 void pme_gpu_clear_grids(const PmeGpu* pmeGpu);
259 /*! \libinternal \brief
260 * Reallocates and copies the pre-computed fractional coordinates' shifts to the GPU.
262 * \param[in] pmeGpu The PME GPU structure.
264 void pme_gpu_realloc_and_copy_fract_shifts(PmeGpu* pmeGpu);
266 /*! \libinternal \brief
267 * Frees the pre-computed fractional coordinates' shifts on the GPU.
269 * \param[in] pmeGpu The PME GPU structure.
271 void pme_gpu_free_fract_shifts(const PmeGpu* pmeGpu);
273 /*! \libinternal \brief
274 * Copies the input real-space grid from the host to the GPU.
276 * \param[in] pmeGpu The PME GPU structure.
277 * \param[in] h_grid The host-side grid buffer.
279 void pme_gpu_copy_input_gather_grid(const PmeGpu* pmeGpu, float* h_grid);
281 /*! \libinternal \brief
282 * Copies the output real-space grid from the GPU to the host.
284 * \param[in] pmeGpu The PME GPU structure.
285 * \param[out] h_grid The host-side grid buffer.
287 void pme_gpu_copy_output_spread_grid(const PmeGpu* pmeGpu, float* h_grid);
289 /*! \libinternal \brief
290 * Copies the spread output spline data and gridline indices from the GPU to the host.
292 * \param[in] pmeGpu The PME GPU structure.
294 void pme_gpu_copy_output_spread_atom_data(const PmeGpu* pmeGpu);
296 /*! \libinternal \brief
297 * Copies the gather input spline data and gridline indices from the host to the GPU.
299 * \param[in] pmeGpu The PME GPU structure.
301 void pme_gpu_copy_input_gather_atom_data(const PmeGpu* pmeGpu);
303 /*! \libinternal \brief
304 * Waits for the grid copying to the host-side buffer after spreading to finish.
306 * \param[in] pmeGpu The PME GPU structure.
308 void pme_gpu_sync_spread_grid(const PmeGpu* pmeGpu);
310 /*! \libinternal \brief
311 * Does the one-time GPU-framework specific PME initialization.
312 * For CUDA, the PME stream is created with the highest priority.
314 * \param[in] pmeGpu The PME GPU structure.
316 void pme_gpu_init_internal(PmeGpu* pmeGpu);
318 /*! \libinternal \brief
319 * Initializes the CUDA FFT structures.
321 * \param[in] pmeGpu The PME GPU structure.
323 void pme_gpu_reinit_3dfft(const PmeGpu* pmeGpu);
325 /*! \libinternal \brief
326 * Destroys the CUDA FFT structures.
328 * \param[in] pmeGpu The PME GPU structure.
330 void pme_gpu_destroy_3dfft(const PmeGpu* pmeGpu);
332 /* The PME stages themselves */
334 /*! \libinternal \brief
335 * A GPU spline computation and charge spreading function.
337 * \param[in] pmeGpu The PME GPU structure.
338 * \param[in] xReadyOnDevice Event synchronizer indicating that the coordinates are ready in the device memory;
339 * can be nullptr when invoked on a separate PME rank or from PME tests.
340 * \param[in] gridIndex Index of the PME grid - unused, assumed to be 0.
341 * \param[out] h_grid The host-side grid buffer (used only if the result of the spread is expected on the host,
342 * e.g. testing or host-side FFT)
343 * \param[in] computeSplines Should the computation of spline parameters and gridline indices be performed.
344 * \param[in] spreadCharges Should the charges/coefficients be spread on the grid.
346 GPU_FUNC_QUALIFIER void pme_gpu_spread(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
347 GpuEventSynchronizer* GPU_FUNC_ARGUMENT(xReadyOnDevice),
348 int GPU_FUNC_ARGUMENT(gridIndex),
349 real* GPU_FUNC_ARGUMENT(h_grid),
350 bool GPU_FUNC_ARGUMENT(computeSplines),
351 bool GPU_FUNC_ARGUMENT(spreadCharges)) GPU_FUNC_TERM;
353 /*! \libinternal \brief
354 * 3D FFT R2C/C2R routine.
356 * \param[in] pmeGpu The PME GPU structure.
357 * \param[in] direction Transform direction (real-to-complex or complex-to-real)
358 * \param[in] gridIndex Index of the PME grid - unused, assumed to be 0.
360 void pme_gpu_3dfft(const PmeGpu* pmeGpu, enum gmx_fft_direction direction, int gridIndex);
362 /*! \libinternal \brief
363 * A GPU Fourier space solving function.
365 * \param[in] pmeGpu The PME GPU structure.
366 * \param[in,out] h_grid The host-side input and output Fourier grid buffer (used only with testing or host-side FFT)
367 * \param[in] gridOrdering Specifies the dimenion ordering of the complex grid. TODO: store this information?
368 * \param[in] computeEnergyAndVirial Tells if the energy and virial computation should be performed.
370 GPU_FUNC_QUALIFIER void pme_gpu_solve(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
371 t_complex* GPU_FUNC_ARGUMENT(h_grid),
372 GridOrdering GPU_FUNC_ARGUMENT(gridOrdering),
373 bool GPU_FUNC_ARGUMENT(computeEnergyAndVirial)) GPU_FUNC_TERM;
375 /*! \libinternal \brief
376 * A GPU force gathering function.
378 * \param[in] pmeGpu The PME GPU structure.
379 * reductions. \param[in] h_grid The host-side grid buffer (used only in testing mode)
381 GPU_FUNC_QUALIFIER void pme_gpu_gather(PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
382 const float* GPU_FUNC_ARGUMENT(h_grid)) GPU_FUNC_TERM;
384 /*! \brief Sets the device pointer to coordinate data
385 * \param[in] pmeGpu The PME GPU structure.
386 * \param[in] d_x Pointer to coordinate data
388 GPU_FUNC_QUALIFIER void pme_gpu_set_kernelparam_coordinates(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
389 DeviceBuffer<gmx::RVec> GPU_FUNC_ARGUMENT(d_x)) GPU_FUNC_TERM;
391 /*! \brief Return pointer to device copy of force data.
392 * \param[in] pmeGpu The PME GPU structure.
393 * \returns Pointer to force data
395 GPU_FUNC_QUALIFIER void* pme_gpu_get_kernelparam_forces(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu))
396 GPU_FUNC_TERM_WITH_RETURN(nullptr);
398 /*! \brief Return pointer to GPU stream.
399 * \param[in] pmeGpu The PME GPU structure.
400 * \returns Pointer to stream object.
402 GPU_FUNC_QUALIFIER const DeviceStream* pme_gpu_get_stream(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu))
403 GPU_FUNC_TERM_WITH_RETURN(nullptr);
405 /*! \brief Return pointer to the sync object triggered after the PME force calculation completion
406 * \param[in] pmeGpu The PME GPU structure.
407 * \returns Pointer to sync object
409 GPU_FUNC_QUALIFIER GpuEventSynchronizer* pme_gpu_get_forces_ready_synchronizer(
410 const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu)) GPU_FUNC_TERM_WITH_RETURN(nullptr);
412 /*! \libinternal \brief
413 * Returns the PME GPU settings
415 * \param[in] pmeGpu The PME GPU structure.
416 * \returns The settings for PME on GPU
418 inline const PmeGpuSettings& pme_gpu_settings(const PmeGpu* pmeGpu)
420 return pmeGpu->settings;
423 /*! \libinternal \brief
424 * Returns the PME GPU staging object
426 * \param[in] pmeGpu The PME GPU structure.
427 * \returns The staging object for PME on GPU
429 inline const PmeGpuStaging& pme_gpu_staging(const PmeGpu* pmeGpu)
431 return pmeGpu->staging;
434 /*! \libinternal \brief
435 * Sets whether the PME module is running in testing mode
437 * \param[in] pmeGpu The PME GPU structure.
438 * \param[in] testing Whether testing mode is on.
440 inline void pme_gpu_set_testing(PmeGpu* pmeGpu, bool testing)
444 pmeGpu->settings.copyAllOutputs = testing;
445 pmeGpu->settings.transferKind = testing ? GpuApiCallBehavior::Sync : GpuApiCallBehavior::Async;
449 /* A block of C++ functions that live in pme_gpu_internal.cpp */
451 /*! \libinternal \brief
452 * Returns the energy and virial GPU outputs, useful for testing.
454 * It is the caller's responsibility to be aware of whether the GPU
455 * handled the solve stage.
457 * \param[in] pme The PME structure.
458 * \param[out] output Pointer to output where energy and virial should be stored.
460 GPU_FUNC_QUALIFIER void pme_gpu_getEnergyAndVirial(const gmx_pme_t& GPU_FUNC_ARGUMENT(pme),
461 PmeOutput* GPU_FUNC_ARGUMENT(output)) GPU_FUNC_TERM;
463 /*! \libinternal \brief
464 * Returns the GPU outputs (forces, energy and virial)
466 * \param[in] pme The PME structure.
467 * \param[in] computeEnergyAndVirial Whether the energy and virial are being computed
468 * \returns The output object.
470 GPU_FUNC_QUALIFIER PmeOutput pme_gpu_getOutput(const gmx_pme_t& GPU_FUNC_ARGUMENT(pme),
471 bool GPU_FUNC_ARGUMENT(computeEnergyAndVirial))
472 GPU_FUNC_TERM_WITH_RETURN(PmeOutput{});
474 /*! \libinternal \brief
475 * Updates the unit cell parameters. Does not check if update is necessary - that is done in pme_gpu_prepare_computation().
477 * \param[in] pmeGpu The PME GPU structure.
478 * \param[in] box The unit cell box.
480 GPU_FUNC_QUALIFIER void pme_gpu_update_input_box(PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
481 const matrix GPU_FUNC_ARGUMENT(box)) GPU_FUNC_TERM;
483 /*! \libinternal \brief
484 * Finishes the PME GPU computation, waiting for the output forces and/or energy/virial to be copied to the host.
485 * If forces were computed, they will have arrived at the external host buffer provided to gather.
486 * If virial/energy were computed, they will have arrived into the internal staging buffer
487 * (even though that should have already happened before even launching the gather).
488 * Finally, cudaEvent_t based GPU timers get updated if enabled. They also need stream synchronization for correctness.
489 * Additionally, device-side buffers are cleared asynchronously for the next computation.
491 * \param[in] pmeGpu The PME GPU structure.
493 void pme_gpu_finish_computation(const PmeGpu* pmeGpu);
495 //! A binary enum for spline data layout transformation
496 enum class PmeLayoutTransform
502 /*! \libinternal \brief
503 * Rearranges the atom spline data between the GPU and host layouts.
504 * Only used for test purposes so far, likely to be horribly slow.
506 * \param[in] pmeGpu The PME GPU structure.
507 * \param[out] atc The PME CPU atom data structure (with a single-threaded layout).
508 * \param[in] type The spline data type (values or derivatives).
509 * \param[in] dimIndex Dimension index.
510 * \param[in] transform Layout transform type
512 GPU_FUNC_QUALIFIER void pme_gpu_transform_spline_atom_data(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
513 const PmeAtomComm* GPU_FUNC_ARGUMENT(atc),
514 PmeSplineDataType GPU_FUNC_ARGUMENT(type),
515 int GPU_FUNC_ARGUMENT(dimIndex),
516 PmeLayoutTransform GPU_FUNC_ARGUMENT(transform)) GPU_FUNC_TERM;
518 /*! \libinternal \brief
519 * Gets a unique index to an element in a spline parameter buffer (theta/dtheta),
520 * which is laid out for GPU spread/gather kernels. The index is wrt the execution block,
521 * in range(0, atomsPerBlock * order * DIM).
522 * This is a wrapper, only used in unit tests.
523 * \param[in] order PME order
524 * \param[in] splineIndex Spline contribution index (from 0 to \p order - 1)
525 * \param[in] dimIndex Dimension index (from 0 to 2)
526 * \param[in] atomIndex Atom index wrt the block.
527 * \param[in] atomsPerWarp Number of atoms processed by a warp.
529 * \returns Index into theta or dtheta array using GPU layout.
531 int getSplineParamFullIndex(int order, int splineIndex, int dimIndex, int atomIndex, int atomsPerWarp);
533 /*! \libinternal \brief
534 * Get the normal/padded grid dimensions of the real-space PME grid on GPU. Only used in tests.
536 * \param[in] pmeGpu The PME GPU structure.
537 * \param[out] gridSize Pointer to the grid dimensions to fill in.
538 * \param[out] paddedGridSize Pointer to the padded grid dimensions to fill in.
540 GPU_FUNC_QUALIFIER void pme_gpu_get_real_grid_sizes(const PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
541 gmx::IVec* GPU_FUNC_ARGUMENT(gridSize),
542 gmx::IVec* GPU_FUNC_ARGUMENT(paddedGridSize)) GPU_FUNC_TERM;
544 /*! \libinternal \brief
545 * (Re-)initializes the PME GPU data at the beginning of the run or on DLB.
547 * \param[in,out] pme The PME structure.
548 * \param[in] deviceInfo The GPU device information structure.
549 * \param[in] pmeGpuProgram The PME GPU program data
550 * \throws gmx::NotImplementedError if this generally valid PME structure is not valid for GPU runs.
552 GPU_FUNC_QUALIFIER void pme_gpu_reinit(gmx_pme_t* GPU_FUNC_ARGUMENT(pme),
553 const DeviceInformation* GPU_FUNC_ARGUMENT(deviceInfo),
554 const PmeGpuProgram* GPU_FUNC_ARGUMENT(pmeGpuProgram)) GPU_FUNC_TERM;
556 /*! \libinternal \brief
557 * Destroys the PME GPU data at the end of the run.
559 * \param[in] pmeGpu The PME GPU structure.
561 GPU_FUNC_QUALIFIER void pme_gpu_destroy(PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu)) GPU_FUNC_TERM;
563 /*! \libinternal \brief
564 * Reallocates the local atoms data (charges, coordinates, etc.). Copies the charges to the GPU.
566 * \param[in] pmeGpu The PME GPU structure.
567 * \param[in] nAtoms The number of particles.
568 * \param[in] charges The pointer to the host-side array of particle charges.
570 * This is a function that should only be called in the beginning of the run and on domain
571 * decomposition. Should be called before the pme_gpu_set_io_ranges.
573 GPU_FUNC_QUALIFIER void pme_gpu_reinit_atoms(PmeGpu* GPU_FUNC_ARGUMENT(pmeGpu),
574 int GPU_FUNC_ARGUMENT(nAtoms),
575 const real* GPU_FUNC_ARGUMENT(charges)) GPU_FUNC_TERM;
577 /*! \brief \libinternal
578 * The PME GPU reinitialization function that is called both at the end of any PME computation and on any load balancing.
580 * This clears the device-side working buffers in preparation for new computation.
582 * \param[in] pmeGpu The PME GPU structure.
584 void pme_gpu_reinit_computation(const PmeGpu* pmeGpu);
587 * Blocks until PME GPU tasks are completed, and gets the output forces and virial/energy
588 * (if they were to be computed).
590 * \param[in] pme The PME data structure.
591 * \param[in] computeEnergyAndVirial Tells if the energy and virial computation should be performed.
592 * \param[out] wcycle The wallclock counter.
593 * \return The output forces, energy and virial
595 GPU_FUNC_QUALIFIER PmeOutput pme_gpu_wait_finish_task(gmx_pme_t* GPU_FUNC_ARGUMENT(pme),
596 bool GPU_FUNC_ARGUMENT(computeEnergyAndVirial),
597 gmx_wallcycle* GPU_FUNC_ARGUMENT(wcycle))
598 GPU_FUNC_TERM_WITH_RETURN(PmeOutput{});