Gromacs
2016.6
|
#include "gmxpre.h"
#include "runner.h"
#include "config.h"
#include <assert.h>
#include <signal.h>
#include <stdlib.h>
#include <string.h>
#include <algorithm>
#include "gromacs/commandline/filenm.h"
#include "gromacs/domdec/domdec.h"
#include "gromacs/domdec/domdec_struct.h"
#include "gromacs/essentialdynamics/edsam.h"
#include "gromacs/ewald/pme.h"
#include "gromacs/fileio/checkpoint.h"
#include "gromacs/fileio/oenv.h"
#include "gromacs/fileio/tpxio.h"
#include "gromacs/gmxlib/md_logging.h"
#include "gromacs/gmxlib/network.h"
#include "gromacs/gpu_utils/gpu_utils.h"
#include "gromacs/hardware/cpuinfo.h"
#include "gromacs/hardware/detecthardware.h"
#include "gromacs/listed-forces/disre.h"
#include "gromacs/listed-forces/orires.h"
#include "gromacs/math/calculate-ewald-splitting-coefficient.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/calc_verletbuf.h"
#include "gromacs/mdlib/constr.h"
#include "gromacs/mdlib/force.h"
#include "gromacs/mdlib/forcerec.h"
#include "gromacs/mdlib/gmx_omp_nthreads.h"
#include "gromacs/mdlib/integrator.h"
#include "gromacs/mdlib/main.h"
#include "gromacs/mdlib/md_support.h"
#include "gromacs/mdlib/mdatoms.h"
#include "gromacs/mdlib/mdrun.h"
#include "gromacs/mdlib/minimize.h"
#include "gromacs/mdlib/nbnxn_search.h"
#include "gromacs/mdlib/qmmm.h"
#include "gromacs/mdlib/sighandler.h"
#include "gromacs/mdlib/sim_util.h"
#include "gromacs/mdlib/tpi.h"
#include "gromacs/mdrunutility/threadaffinity.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/md_enums.h"
#include "gromacs/mdtypes/state.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/pulling/pull.h"
#include "gromacs/pulling/pull_rotation.h"
#include "gromacs/timing/wallcycle.h"
#include "gromacs/topology/mtop_util.h"
#include "gromacs/trajectory/trajectoryframe.h"
#include "gromacs/utility/cstringutil.h"
#include "gromacs/utility/exceptions.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/gmxmpi.h"
#include "gromacs/utility/pleasecite.h"
#include "gromacs/utility/smalloc.h"
#include "deform.h"
#include "md.h"
#include "membed.h"
#include "repl_ex.h"
#include "resource-division.h"
Implements the MD runner routine calling all integrators.
Macros | |
#define | NNSTL sizeof(nstlist_try)/sizeof(nstlist_try[0]) |
Number of elements in the neighborsearch list trials. | |
Functions | |
static void | increase_nstlist (FILE *fp, t_commrec *cr, t_inputrec *ir, int nstlist_cmdline, const gmx_mtop_t *mtop, matrix box, gmx_bool bGPU, const gmx::CpuInfo &cpuinfo) |
Try to increase nstlist when using the Verlet cut-off scheme. | |
static void | prepare_verlet_scheme (FILE *fplog, t_commrec *cr, t_inputrec *ir, int nstlist_cmdline, const gmx_mtop_t *mtop, matrix box, gmx_bool bUseGPU, const gmx::CpuInfo &cpuinfo) |
Initialize variables for Verlet scheme simulation. | |
static void | override_nsteps_cmdline (FILE *fplog, gmx_int64_t nsteps_cmdline, t_inputrec *ir, const t_commrec *cr) |
Override the nslist value in inputrec. More... | |
static integrator_t * | gmx::my_integrator (unsigned int ei) |
Return the correct integrator function. | |
int | gmx::mdrunner (gmx_hw_opt_t *hw_opt, FILE *fplog, struct t_commrec *cr, int nfile, const t_filenm fnm[], const gmx_output_env_t *oenv, gmx_bool bVerbose, int nstglobalcomm, ivec ddxyz, int dd_rank_order, int npme, real rdd, real rconstr, const char *dddlb_opt, real dlb_scale, const char *ddcsx, const char *ddcsy, const char *ddcsz, const char *nbpu_opt, int nstlist_cmdline, gmx_int64_t nsteps_cmdline, int nstepout, int resetstep, int nmultisim, int repl_ex_nst, int repl_ex_nex, int repl_ex_seed, real pforce, real cpt_period, real max_hours, int imdport, unsigned long Flags) |
Driver routine, that calls the different methods. More... | |
Variables | |
gmx_int64_t | deform_init_init_step_tpx |
First step used in pressure scaling. | |
matrix | deform_init_box_tpx |
Initial box for pressure scaling. | |
tMPI_Thread_mutex_t | deform_init_box_mutex = TMPI_THREAD_MUTEX_INITIALIZER |
MPI variable for use in pressure scaling. | |
static const int | nbnxnReferenceNstlist = 10 |
Cost of non-bonded kernels. More... | |
const int | nstlist_try [] = { 20, 25, 40 } |
The values to try when switching. | |
static const float | nbnxn_cpu_listfac_ok = 1.05 |
Max OK performance ratio beween force calc and neighbor searching. | |
static const float | nbnxn_cpu_listfac_max = 1.09 |
Too high performance ratio beween force calc and neighbor searching. | |
static const float | nbnxn_knl_listfac_ok = 1.22 |
Max OK performance ratio beween force calc and neighbor searching. | |
static const float | nbnxn_knl_listfac_max = 1.3 |
Too high performance ratio beween force calc and neighbor searching. | |
static const float | nbnxn_gpu_listfac_ok = 1.20 |
Max OK performance ratio beween force calc and neighbor searching. | |
static const float | nbnxn_gpu_listfac_max = 1.30 |
Too high performance ratio beween force calc and neighbor searching. | |
|
static |
Override the nslist value in inputrec.
with value passed on the command line (if any)
|
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).