Gromacs  2020.4
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Macros | Functions
#include "gmxpre.h"
#include "domdec_setup.h"
#include <cassert>
#include <cinttypes>
#include <cmath>
#include <cstdio>
#include "gromacs/domdec/domdec.h"
#include "gromacs/domdec/domdec_struct.h"
#include "gromacs/domdec/options.h"
#include "gromacs/ewald/pme.h"
#include "gromacs/gmxlib/network.h"
#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/perf_est.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/md_enums.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/topology/topology.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/logger.h"
#include "gromacs/utility/stringutil.h"
#include "box.h"
#include "domdec_internal.h"
#include "utility.h"
+ Include dependency graph for domdec_setup.cpp:

Description

This file defines functions used by the domdec module in its initial setup phase.

Author
Berk Hess hess@.nosp@m.kth..nosp@m.se

Macros

#define DD_GRID_MARGIN_PRES_SCALE   1.05
 Margin for setting up the DD grid.
 

Functions

static void factorize (int n, std::vector< int > *fac, std::vector< int > *mfac)
 Factorize n. More...
 
static int largest_divisor (int n)
 Find largest divisor of n smaller than n.
 
static int lcd (int n1, int n2)
 Compute largest common divisor of n1 and n2.
 
static gmx_bool fits_pme_ratio (int nrank_tot, int nrank_pme, float ratio)
 Returns TRUE when there are enough PME ranks for the ratio.
 
static gmx_bool fits_pp_pme_perf (int ntot, int npme, float ratio)
 Returns TRUE when npme out of ntot ranks doing PME is expected to give reasonable performance.
 
static int guess_npme (const gmx::MDLogger &mdlog, const gmx_mtop_t &mtop, const t_inputrec &ir, const matrix box, int nrank_tot)
 Make a guess for the number of PME ranks to use.
 
static int div_up (int n, int f)
 Return n divided by f rounded up to the next integer.
 
real comm_box_frac (const gmx::IVec &dd_nc, real cutoff, const gmx_ddbox_t &ddbox)
 Returns the volume fraction of the system that is communicated.
 
static gmx_bool inhomogeneous_z (const t_inputrec &ir)
 Return whether the DD inhomogeneous in the z direction.
 
static float comm_pme_cost_vol (int npme, int a, int b, int c)
 Estimate cost of PME FFT communication. More...
 
static float comm_cost_est (real limit, real cutoff, const matrix box, const gmx_ddbox_t &ddbox, int natoms, const t_inputrec &ir, float pbcdxr, int npme_tot, const gmx::IVec &nc)
 Estimate cost of communication for a possible domain decomposition.
 
static void assign_factors (const real limit, const bool request1D, const real cutoff, const matrix box, const gmx_ddbox_t &ddbox, int natoms, const t_inputrec &ir, float pbcdxr, int npme, int ndiv, const int *div, const int *mdiv, gmx::IVec *irTryPtr, gmx::IVec *opt)
 Assign penalty factors to possible domain decompositions, based on the estimated communication costs.
 
static gmx::IVec optimizeDDCells (const gmx::MDLogger &mdlog, const int numRanksRequested, const int numPmeOnlyRanks, const real cellSizeLimit, const bool request1DAnd1Pulse, const gmx_mtop_t &mtop, const matrix box, const gmx_ddbox_t &ddbox, const t_inputrec &ir, const DDSystemInfo &systemInfo)
 Determine the optimal distribution of DD cells for the simulation system and number of MPI ranks. More...
 
real getDDGridSetupCellSizeLimit (const gmx::MDLogger &mdlog, const bool request1DAnd1Pulse, const bool bDynLoadBal, const real dlb_scale, const t_inputrec &ir, const real systemInfoCellSizeLimit)
 Return the minimum cell size (in nm) required for DD.
 
void checkForValidRankCountRequests (const int numRanksRequested, const bool usingPme, const int numPmeRanksRequested, const bool checkForLargePrimeFactors)
 Checks that requests for PP and PME ranks honor basic expectations. More...
 
static int getNumPmeOnlyRanksToUse (const gmx::MDLogger &mdlog, const DomdecOptions &options, const gmx_mtop_t &mtop, const t_inputrec &ir, const matrix box, const int numRanksRequested)
 Return the number of PME-only ranks used by the simulation. More...
 
static int set_dd_dim (const gmx::IVec &numDDCells, const DDSettings &ddSettings, ivec *dims)
 Sets the order of the DD dimensions, returns the number of DD dimensions.
 
DDGridSetup getDDGridSetup (const gmx::MDLogger &mdlog, const t_commrec *cr, const int numRanksRequested, const DomdecOptions &options, const DDSettings &ddSettings, const DDSystemInfo &systemInfo, const real cellSizeLimit, const gmx_mtop_t &mtop, const t_inputrec &ir, const matrix box, gmx::ArrayRef< const gmx::RVec > xGlobal, gmx_ddbox_t *ddbox)
 Determines the DD grid setup. More...
 

Function Documentation

void checkForValidRankCountRequests ( int  numRanksRequested,
bool  usingPme,
int  numPmeRanksRequested,
bool  checkForLargePrimeFactors 
)

Checks that requests for PP and PME ranks honor basic expectations.

Issues a fatal error if there are more PME ranks than PP, or if the count of PP ranks has a prime factor that is too large to be likely to have good performance.

static float comm_pme_cost_vol ( int  npme,
int  a,
int  b,
int  c 
)
static

Estimate cost of PME FFT communication.

This only takes the communication into account and not imbalance in the calculation. But the imbalance in communication and calculation are similar and therefore these formulas also prefer load balance in the FFT and pme_solve calculation.

static void factorize ( int  n,
std::vector< int > *  fac,
std::vector< int > *  mfac 
)
static

Factorize n.

Parameters
[in]nValue to factorize
[out]facVector of factors (to be allocated in this function)
[out]mfacVector with the number of times each factor repeats in the factorization (to be allocated in this function)
DDGridSetup getDDGridSetup ( const gmx::MDLogger mdlog,
const t_commrec *  cr,
int  numRanksRequested,
const gmx::DomdecOptions options,
const DDSettings &  ddSettings,
const DDSystemInfo &  systemInfo,
real  cellSizeLimit,
const gmx_mtop_t &  mtop,
const t_inputrec &  ir,
const matrix  box,
gmx::ArrayRef< const gmx::RVec xGlobal,
gmx_ddbox_t *  ddbox 
)

Determines the DD grid setup.

Either implements settings required by the user, or otherwise chooses estimated optimal number of separate PME ranks and DD grid cell setup, DD cell size limits, and the initial ddbox.

static int getNumPmeOnlyRanksToUse ( const gmx::MDLogger mdlog,
const DomdecOptions options,
const gmx_mtop_t &  mtop,
const t_inputrec &  ir,
const matrix  box,
const int  numRanksRequested 
)
static

Return the number of PME-only ranks used by the simulation.

If the user did not choose a number, then decide for them.

static gmx::IVec optimizeDDCells ( const gmx::MDLogger mdlog,
const int  numRanksRequested,
const int  numPmeOnlyRanks,
const real  cellSizeLimit,
const bool  request1DAnd1Pulse,
const gmx_mtop_t &  mtop,
const matrix  box,
const gmx_ddbox_t &  ddbox,
const t_inputrec &  ir,
const DDSystemInfo &  systemInfo 
)
static

Determine the optimal distribution of DD cells for the simulation system and number of MPI ranks.

Returns
The optimal grid cell choice. The latter will contain all zeros if no valid cell choice exists.