Gromacs
2019.3
|
#include "gmxpre.h"
#include "pme.h"
#include "config.h"
#include <cassert>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <algorithm>
#include <list>
#include "gromacs/domdec/domdec.h"
#include "gromacs/ewald/ewald-utils.h"
#include "gromacs/fft/parallel_3dfft.h"
#include "gromacs/fileio/pdbio.h"
#include "gromacs/gmxlib/network.h"
#include "gromacs/gmxlib/nrnb.h"
#include "gromacs/hardware/hw_info.h"
#include "gromacs/math/gmxcomplex.h"
#include "gromacs/math/invertmatrix.h"
#include "gromacs/math/units.h"
#include "gromacs/math/vec.h"
#include "gromacs/math/vectypes.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/md_enums.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/timing/cyclecounter.h"
#include "gromacs/timing/wallcycle.h"
#include "gromacs/timing/walltime_accounting.h"
#include "gromacs/topology/topology.h"
#include "gromacs/utility/basedefinitions.h"
#include "gromacs/utility/cstringutil.h"
#include "gromacs/utility/exceptions.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/gmxmpi.h"
#include "gromacs/utility/gmxomp.h"
#include "gromacs/utility/logger.h"
#include "gromacs/utility/real.h"
#include "gromacs/utility/smalloc.h"
#include "gromacs/utility/stringutil.h"
#include "gromacs/utility/unique_cptr.h"
#include "calculate-spline-moduli.h"
#include "pme-gather.h"
#include "pme-gpu-internal.h"
#include "pme-grid.h"
#include "pme-internal.h"
#include "pme-redistribute.h"
#include "pme-solve.h"
#include "pme-spline-work.h"
#include "pme-spread.h"
This file contains function definitions necessary for computing energies and forces for the PME long-ranged part (Coulomb and LJ).
Functions | |
static bool | addMessageIfNotSupported (const std::list< std::string > &errorReasons, std::string *error) |
Help build a descriptive message in error if there are errorReasons why PME on GPU is not supported. More... | |
bool | pme_gpu_supports_build (std::string *error) |
Checks whether the GROMACS build allows to run PME on GPU. TODO: this partly duplicates an internal PME assert function pme_gpu_check_restrictions(), except that works with a formed gmx_pme_t structure. Should that one go away/work with inputrec? More... | |
bool | pme_gpu_supports_hardware (const gmx_hw_info_t &hwinfo, std::string *error) |
Checks whether the detected (GPU) hardware allows to run PME on GPU. More... | |
bool | pme_gpu_supports_input (const t_inputrec &ir, const gmx_mtop_t &mtop, std::string *error) |
Checks whether the input system allows to run PME on GPU. TODO: this partly duplicates an internal PME assert function pme_gpu_check_restrictions(), except that works with a formed gmx_pme_t structure. Should that one go away/work with inputrec? More... | |
static bool | pme_gpu_check_restrictions (const gmx_pme_t *pme, std::string *error) |
Finds out if PME with given inputs is possible to run on GPU. This function is an internal final check, validating the whole PME structure on creation, but it still duplicates the preliminary checks from the above (externally exposed) pme_gpu_supports_input() - just in case. More... | |
PmeRunMode | pme_run_mode (const gmx_pme_t *pme) |
Returns the active PME codepath (CPU, GPU, mixed). More... | |
gmx::PinningPolicy | pme_get_pinning_policy () |
Return the pinning policy appropriate for this build configuration for relevant buffers used for PME task on this rank (e.g. running on a GPU). | |
static void | setup_coordinate_communication (pme_atomcomm_t *atc) |
Set up coordinate communication. | |
static int | mult_up (int n, int f) |
Round n up to the next multiple of f . | |
static double | estimate_pme_load_imbalance (struct gmx_pme_t *pme) |
Return estimate of the load imbalance from the PME grid not being a good match for the number of PME ranks. | |
static void | init_atomcomm (struct gmx_pme_t *pme, pme_atomcomm_t *atc, int dimind, gmx_bool bSpread) |
Initialize atom communication data structure. | |
static void | destroy_atomcomm (pme_atomcomm_t *atc) |
Destroy an atom communication data structure and its child structs. | |
static void | init_overlap_comm (pme_overlap_t *ol, int norder, int nnodes, int nodeid, int ndata, int commplainsize) |
Initialize data structure for communication. | |
int | minimalPmeGridSize (int pmeOrder) |
Return the smallest allowed PME grid size for pmeOrder . | |
bool | gmx_pme_check_restrictions (int pme_order, int nkx, int nky, int nkz, int numPmeDomainsAlongX, bool useThreads, bool errorsAreFatal) |
Check restrictions on pme_order and the PME grid nkx,nky,nkz. More... | |
static int | div_round_up (int enumerator, int denominator) |
Round enumerator . | |
gmx_pme_t * | gmx_pme_init (const t_commrec *cr, const NumPmeDomains &numPmeDomains, const t_inputrec *ir, int homenr, gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj, gmx_bool bReproducible, real ewaldcoeff_q, real ewaldcoeff_lj, int nthread, PmeRunMode runMode, PmeGpu *pmeGpu, const gmx_device_info_t *gpuInfo, PmeGpuProgramHandle pmeGpuProgram, const gmx::MDLogger &) |
Construct PME data. More... | |
void | gmx_pme_reinit (struct gmx_pme_t **pmedata, const t_commrec *cr, struct gmx_pme_t *pme_src, const t_inputrec *ir, const ivec grid_size, real ewaldcoeff_q, real ewaldcoeff_lj) |
As gmx_pme_init, but takes most settings, except the grid/Ewald coefficients, from pme_src. This is only called when the PME cut-off/grid size changes. | |
void | gmx_pme_calc_energy (struct gmx_pme_t *pme, int n, rvec *x, real *q, real *V) |
Calculate the PME grid energy V for n charges. More... | |
static void | calc_initial_lb_coeffs (struct gmx_pme_t *pme, const real *local_c6, const real *local_sigma) |
Calculate initial Lorentz-Berthelot coefficients for LJ-PME. | |
static void | calc_next_lb_coeffs (struct gmx_pme_t *pme, const real *local_sigma) |
Calculate next Lorentz-Berthelot coefficients for LJ-PME. | |
int | gmx_pme_do (struct gmx_pme_t *pme, int start, int homenr, rvec x[], rvec f[], real chargeA[], real chargeB[], real c6A[], real c6B[], real sigmaA[], real sigmaB[], matrix box, const t_commrec *cr, int maxshift_x, int maxshift_y, t_nrnb *nrnb, gmx_wallcycle *wcycle, matrix vir_q, matrix vir_lj, real *energy_q, real *energy_lj, real lambda_q, real lambda_lj, real *dvdlambda_q, real *dvdlambda_lj, int flags) |
Do a PME calculation on a CPU for the long range electrostatics and/or LJ. More... | |
void | gmx_pme_destroy (gmx_pme_t *pme) |
Destroys the PME data structure. | |
void | gmx_pme_reinit_atoms (const gmx_pme_t *pme, const int nAtoms, const real *charges) |
This function updates the local atom data on GPU after DD (charges, coordinates, etc.). TODO: it should update the PME CPU atom data as well. (currently PME CPU call gmx_pme_do() gets passed the input pointers for each computation). More... | |
Variables | |
const int | gmxCacheLineSize = 64 |
Number of bytes in a cache line. More... | |
|
static |
Help build a descriptive message in error
if there are errorReasons
why PME on GPU is not supported.
Calculate the PME grid energy V for n charges.
The potential (found in pme
) must have been found already with a call to gmx_pme_do() with at least GMX_PME_SPREAD and GMX_PME_SOLVE specified. Note that the charges are not spread on the grid in the pme struct. Currently does not work in parallel or with free energy.
bool gmx_pme_check_restrictions | ( | int | pme_order, |
int | nkx, | ||
int | nky, | ||
int | nkz, | ||
int | numPmeDomainsAlongX, | ||
bool | useThreads, | ||
bool | errorsAreFatal | ||
) |
Check restrictions on pme_order and the PME grid nkx,nky,nkz.
With errorsAreFatal=true, an exception or fatal error is generated on violation of restrictions. With errorsAreFatal=false, false is returned on violation of restrictions. When all restrictions are obeyed, true is returned. Argument useThreads tells if any MPI rank doing PME uses more than 1 threads. If at calling useThreads is unknown, pass true for conservative checking.
The PME GPU restrictions are checked separately during pme_gpu_init().
int gmx_pme_do | ( | struct gmx_pme_t * | pme, |
int | start, | ||
int | homenr, | ||
rvec | x[], | ||
rvec | f[], | ||
real | chargeA[], | ||
real | chargeB[], | ||
real | c6A[], | ||
real | c6B[], | ||
real | sigmaA[], | ||
real | sigmaB[], | ||
matrix | box, | ||
const t_commrec * | cr, | ||
int | maxshift_x, | ||
int | maxshift_y, | ||
t_nrnb * | nrnb, | ||
gmx_wallcycle * | wcycle, | ||
matrix | vir_q, | ||
matrix | vir_lj, | ||
real * | energy_q, | ||
real * | energy_lj, | ||
real | lambda_q, | ||
real | lambda_lj, | ||
real * | dvdlambda_q, | ||
real * | dvdlambda_lj, | ||
int | flags | ||
) |
Do a PME calculation on a CPU for the long range electrostatics and/or LJ.
The meaning of flags
is defined above, and determines which parts of the calculation are performed.
gmx_pme_t* gmx_pme_init | ( | const t_commrec * | cr, |
const NumPmeDomains & | numPmeDomains, | ||
const t_inputrec * | ir, | ||
int | homenr, | ||
gmx_bool | bFreeEnergy_q, | ||
gmx_bool | bFreeEnergy_lj, | ||
gmx_bool | bReproducible, | ||
real | ewaldcoeff_q, | ||
real | ewaldcoeff_lj, | ||
int | nthread, | ||
PmeRunMode | runMode, | ||
PmeGpu * | pmeGpu, | ||
const gmx_device_info_t * | gpuInfo, | ||
PmeGpuProgramHandle | pmeGpuProgram, | ||
const gmx::MDLogger & | mdlog | ||
) |
Construct PME data.
gmx::InconsistentInputError | if input grid sizes/PME order are inconsistent. |
GpuManager
that holds gmx_device_info_t
* and PmeGpuProgramHandle
and perhaps other related things whose lifetime can/should exceed that of a task (or perhaps task manager). See Redmine #2522. void gmx_pme_reinit_atoms | ( | const gmx_pme_t * | pme, |
int | nAtoms, | ||
const real * | charges | ||
) |
This function updates the local atom data on GPU after DD (charges, coordinates, etc.). TODO: it should update the PME CPU atom data as well. (currently PME CPU call gmx_pme_do() gets passed the input pointers for each computation).
[in] | pme | The PME structure. |
[in] | nAtoms | The number of particles. |
[in] | charges | The pointer to the array of particle charges. |
|
static |
Finds out if PME with given inputs is possible to run on GPU. This function is an internal final check, validating the whole PME structure on creation, but it still duplicates the preliminary checks from the above (externally exposed) pme_gpu_supports_input() - just in case.
[in] | pme | The PME structure. |
[out] | error | The error message if the input is not supported on GPU. |
bool pme_gpu_supports_build | ( | std::string * | error | ) |
Checks whether the GROMACS build allows to run PME on GPU. TODO: this partly duplicates an internal PME assert function pme_gpu_check_restrictions(), except that works with a formed gmx_pme_t structure. Should that one go away/work with inputrec?
[out] | error | If non-null, the error message when PME is not supported on GPU. |
bool pme_gpu_supports_hardware | ( | const gmx_hw_info_t & | hwinfo, |
std::string * | error | ||
) |
Checks whether the detected (GPU) hardware allows to run PME on GPU.
[in] | hwinfo | Information about the detected hardware |
[out] | error | If non-null, the error message when PME is not supported on GPU. |
bool pme_gpu_supports_input | ( | const t_inputrec & | ir, |
const gmx_mtop_t & | mtop, | ||
std::string * | error | ||
) |
Checks whether the input system allows to run PME on GPU. TODO: this partly duplicates an internal PME assert function pme_gpu_check_restrictions(), except that works with a formed gmx_pme_t structure. Should that one go away/work with inputrec?
[in] | ir | Input system. |
[in] | mtop | Complete system topology to check if an FE simulation perturbs charges. |
[out] | error | If non-null, the error message if the input is not supported on GPU. |
PmeRunMode pme_run_mode | ( | const gmx_pme_t * | pme | ) |
Returns the active PME codepath (CPU, GPU, mixed).
[in] | pme | The PME data structure. |
const int gmxCacheLineSize = 64 |
Number of bytes in a cache line.
Must also be a multiple of the SIMD and SIMD4 register size, to preserve alignment.