#include "gmxpre.h"
#include "pairlist_tuning.h"
#include <cassert>
#include <cmath>
#include <cstdlib>
#include <algorithm>
#include <string>
#include "gromacs/domdec/domdec.h"
#include "gromacs/hardware/cpuinfo.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/calc_verletbuf.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/interaction_const.h"
#include "gromacs/mdtypes/multipletimestepping.h"
#include "gromacs/mdtypes/state.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/topology/topology.h"
#include "gromacs/utility/cstringutil.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/logger.h"
#include "gromacs/utility/strconvert.h"
#include "gromacs/utility/stringutil.h"
#include "nbnxm_geometry.h"
#include "pairlistsets.h"
Implements functions for tuning adjustable parameters for the nbnxn non-bonded search and interaction kernels.
- Author
- Berk Hess hess@.nosp@m.kth..nosp@m.se
|
static bool | supportsDynamicPairlistGenerationInterval (const t_inputrec &ir) |
| Returns if we can (heuristically) change nstlist and rlist. More...
|
|
void | increaseNstlist (FILE *fp, t_commrec *cr, t_inputrec *ir, int nstlist_cmdline, const gmx_mtop_t *mtop, const matrix box, bool useOrEmulateGpuForNonbondeds, const gmx::CpuInfo &cpuinfo) |
| Try to increase nstlist when using the Verlet cut-off scheme. More...
|
|
static void | setDynamicPairlistPruningParameters (const t_inputrec &inputrec, const gmx_mtop_t &mtop, const matrix box, const bool useGpuList, const VerletbufListSetup &listSetup, const bool userSetNstlistPrune, const interaction_const_t &interactionConst, PairlistParams *listParams) |
| Set the dynamic pairlist pruning parameters in ic . More...
|
|
static std::string | formatListSetup (const std::string &listName, int nstList, int nstListForSpacing, real rList, real interactionCutoff) |
| Returns a string describing the setup of a single pair-list. More...
|
|
void | setupDynamicPairlistPruning (const gmx::MDLogger &mdlog, const t_inputrec &inputrec, const gmx_mtop_t &mtop, matrix box, const interaction_const_t &interactionConst, PairlistParams *listParams) |
| Set up the dynamic pairlist pruning. More...
|
|
static std::string formatListSetup |
( |
const std::string & |
listName, |
|
|
int |
nstList, |
|
|
int |
nstListForSpacing, |
|
|
real |
rList, |
|
|
real |
interactionCutoff |
|
) |
| |
|
static |
Returns a string describing the setup of a single pair-list.
- Parameters
-
[in] | listName | Short name of the list, can be "" |
[in] | nstList | The list update interval in steps |
[in] | nstListForSpacing | Update interval for setting the number characters for printing nstList |
[in] | rList | List cut-off radius |
[in] | interactionCutoff | The interaction cut-off, use for printing the list buffer size |
void increaseNstlist |
( |
FILE * |
fplog, |
|
|
t_commrec * |
cr, |
|
|
t_inputrec * |
ir, |
|
|
int |
nstlistOnCmdline, |
|
|
const gmx_mtop_t * |
mtop, |
|
|
const matrix |
box, |
|
|
bool |
useOrEmulateGpuForNonbondeds, |
|
|
const gmx::CpuInfo & |
cpuinfo |
|
) |
| |
Try to increase nstlist when using the Verlet cut-off scheme.
- Parameters
-
[in,out] | fplog | Log file |
[in] | cr | The communication record |
[in] | ir | The input parameter record |
[in] | nstlistOnCmdline | The value of nstlist provided on the command line |
[in] | mtop | The global topology |
[in] | box | The unit cell |
[in] | useOrEmulateGpuForNonbondeds | Tells if we are using a GPU for non-bondeds |
[in] | cpuinfo | Information about the CPU(s) |
static void setDynamicPairlistPruningParameters |
( |
const t_inputrec & |
inputrec, |
|
|
const gmx_mtop_t & |
mtop, |
|
|
const matrix |
box, |
|
|
const bool |
useGpuList, |
|
|
const VerletbufListSetup & |
listSetup, |
|
|
const bool |
userSetNstlistPrune, |
|
|
const interaction_const_t & |
interactionConst, |
|
|
PairlistParams * |
listParams |
|
) |
| |
|
static |
Set the dynamic pairlist pruning parameters in ic
.
- Parameters
-
[in] | inputrec | The input parameter record |
[in] | mtop | The global topology |
[in] | box | The unit cell |
[in] | useGpuList | Tells if we are using a GPU type pairlist |
[in] | listSetup | The nbnxn pair list setup |
[in] | userSetNstlistPrune | The user set ic->nstlistPrune (using an env.var.) |
[in] | interactionConst | The nonbonded interactions constants |
[in,out] | listParams | The list setup parameters |
void setupDynamicPairlistPruning |
( |
const gmx::MDLogger & |
mdlog, |
|
|
const t_inputrec & |
inputrec, |
|
|
const gmx_mtop_t & |
mtop, |
|
|
matrix |
box, |
|
|
const interaction_const_t & |
interactionConst, |
|
|
PairlistParams * |
listParams |
|
) |
| |
Set up the dynamic pairlist pruning.
- Parameters
-
[in,out] | mdlog | MD logger |
[in] | inputrec | The input parameter record |
[in] | mtop | The global topology |
[in] | box | The unit cell |
[in] | interactionConst | The nonbonded interactions constants |
[in,out] | listParams | The list setup parameters |
static bool supportsDynamicPairlistGenerationInterval |
( |
const t_inputrec & |
ir | ) |
|
|
static |
Returns if we can (heuristically) change nstlist and rlist.
- Parameters
-
[in] | ir | The input parameter record |
const int c_nbnxnDynamicListPruningMinLifetime = 4 |
|
static |
The minimum nstlist for dynamic pair list pruning.
In most cases going lower than 4 will lead to a too high pruning cost. This value should be a multiple of c_nbnxnGpuRollingListPruningInterval
const int c_nbnxnGpuRollingListPruningInterval = 2 |
|
static |
The interval in steps at which we perform dynamic, rolling pruning on a GPU.
Ideally we should auto-tune this value. Not considering overheads, 1 would be the ideal value. But 2 seems a reasonable compromise that reduces GPU kernel launch overheads and also avoids inefficiency on large GPUs when pruning small lists. Because with domain decomposition we alternate local/non-local pruning at even/odd steps, which gives a period of 2, this value currenly needs to be 2, which is indirectly asserted when the GPU pruning is dispatched during the force evaluation.
const int nbnxnReferenceNstlist = 10 |
|
static |
Cost of non-bonded kernels.
We determine the extra cost of the non-bonded kernels compared to a reference nstlist value of 10 (which is the default in grompp).