From 6f1f895fcbdf44d4dfe01ac237ebe45fc4130b98 Mon Sep 17 00:00:00 2001 From: Mark Abraham Date: Wed, 28 Feb 2018 16:20:03 +0100 Subject: [PATCH] Simplified handling of simulation runners Argments to runners are now passed via struct members. That struct is expected to change over the medium term, as most of its contents should either be in a container of modules, or data with such a module. As such, it does not make sense to do a lot of coding to manage its contents in a sound way. The practical effects of using this struct are almost the same as the old approach of passing a forest of arguments - some variables are declared in a scope, initialized with values, then used once. The form is chosen so that * the do_md() etc. functions do not need a large number of useless textual changes now * future changes to names and types of the former parameters do not need to be edited in matching way in multiple places of code and their doxygen * we don't need to mark things gmx_unused anywhere Useless unchecked return values have been removed. Some excessive includes in integrator.h were removed, which generates a few minor fixes elsewhere. Refs #1793 Change-Id: I678598175192c9c68113fdd79fcee17f8e5c504e --- src/gromacs/mdrun/CMakeLists.txt | 1 + src/gromacs/mdrun/{minimize.h => integrator.cpp} | 69 ++++++--- src/gromacs/mdrun/integrator.h | 159 +++++++++++++------- src/gromacs/mdrun/md.cpp | 53 ++----- src/gromacs/mdrun/minimize.cpp | 178 ++--------------------- src/gromacs/mdrun/minimize.h | 12 -- src/gromacs/mdrun/runner.cpp | 71 +++------ src/gromacs/mdrun/tpi.cpp | 46 +----- 8 files changed, 206 insertions(+), 383 deletions(-) copy src/gromacs/mdrun/{minimize.h => integrator.cpp} (55%) diff --git a/src/gromacs/mdrun/CMakeLists.txt b/src/gromacs/mdrun/CMakeLists.txt index 0654e1cc75..3b3b9b3f15 100644 --- a/src/gromacs/mdrun/CMakeLists.txt +++ b/src/gromacs/mdrun/CMakeLists.txt @@ -33,6 +33,7 @@ # the research papers on the package. Check out http://www.gromacs.org. gmx_add_libgromacs_sources( + integrator.cpp md.cpp minimize.cpp runner.cpp diff --git a/src/gromacs/mdrun/minimize.h b/src/gromacs/mdrun/integrator.cpp similarity index 55% copy from src/gromacs/mdrun/minimize.h copy to src/gromacs/mdrun/integrator.cpp index 00ed45ffd7..f5158bcae7 100644 --- a/src/gromacs/mdrun/minimize.h +++ b/src/gromacs/mdrun/integrator.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2015,2018, by the GROMACS development team, led by + * Copyright (c) 2018, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -32,31 +32,64 @@ * To help us fund GROMACS development, we humbly ask that you cite * the research papers on the package. Check out http://www.gromacs.org. */ -/*! \brief Declares the integrators for energy minimization and NMA +/*! \internal + * \brief Defines the dispatch function for the .mdp integrator field. * * \author David van der Spoel + * \author Mark Abraham * \ingroup module_mdrun */ -#ifndef GMX_MDRUN_MINIMIZE_H -#define GMX_MDRUN_MINIMIZE_H +#include "gmxpre.h" #include "integrator.h" +#include "gromacs/mdtypes/md_enums.h" +#include "gromacs/utility/exceptions.h" + namespace gmx { -//! Steepest descents energy minimization -integrator_t do_steep; - -//! Conjugate gradient energy minimization -integrator_t do_cg; - -//! Conjugate gradient energy minimization using the L-BFGS algorithm -integrator_t do_lbfgs; - -//! Normal mode analysis -integrator_t do_nm; - -} // namespace gmx +//! \brief Run the correct integrator function. +void Integrator::run(unsigned int ei) +{ + switch (ei) + { + case eiMD: + case eiBD: + case eiSD1: + case eiVV: + case eiVVAK: + if (!EI_DYNAMICS(ei)) + { + GMX_THROW(APIError("do_md integrator would be called for a non-dynamical integrator")); + } + do_md(); + break; + case eiSteep: + do_steep(); + break; + case eiCG: + do_cg(); + break; + case eiNM: + do_nm(); + break; + case eiLBFGS: + do_lbfgs(); + break; + case eiTPI: + case eiTPIC: + if (!EI_TPI(ei)) + { + GMX_THROW(APIError("do_tpi integrator would be called for a non-TPI integrator")); + } + do_tpi(); + break; + case eiSD2_REMOVED: + GMX_THROW(NotImplementedError("SD2 integrator has been removed")); + default: + GMX_THROW(APIError("Non existing integrator selected")); + } +} -#endif // GMX_MDRUN_MINIMIZE_H +} // namespace diff --git a/src/gromacs/mdrun/integrator.h b/src/gromacs/mdrun/integrator.h index d5be821410..04ff880cbf 100644 --- a/src/gromacs/mdrun/integrator.h +++ b/src/gromacs/mdrun/integrator.h @@ -32,9 +32,11 @@ * To help us fund GROMACS development, we humbly ask that you cite * the research papers on the package. Check out http://www.gromacs.org. */ -/*! \brief Declares the integrator type for mdrun +/*! \internal + * \brief Declares the integrator interface for mdrun * * \author David van der Spoel + * \author Mark Abraham * \ingroup module_mdrun */ #ifndef GMX_MDRUN_INTEGRATOR_H @@ -42,27 +44,27 @@ #include -#include "gromacs/gmxlib/nrnb.h" -#include "gromacs/mdlib/constr.h" -#include "gromacs/mdlib/vsite.h" -#include "gromacs/mdtypes/fcdata.h" -#include "gromacs/mdtypes/forcerec.h" -#include "gromacs/timing/wallcycle.h" -#include "gromacs/timing/walltime_accounting.h" #include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/real.h" class energyhistory_t; +struct gmx_constr; struct gmx_mtop_t; struct gmx_membed_t; struct gmx_multisim_t; struct gmx_output_env_t; +struct gmx_vsite_t; +struct gmx_wallcycle; +struct gmx_walltime_accounting; struct MdrunOptions; struct ObservablesHistory; struct ReplicaExchangeParameters; struct t_commrec; +struct t_fcdata; +struct t_forcerec; struct t_filenm; struct t_inputrec; +struct t_nrnb; class t_state; namespace gmx @@ -72,50 +74,105 @@ class IMDOutputProvider; class MDLogger; class MDAtoms; -/*! \brief Integrator algorithm implementation. +//! Function type for integrator code. +using IntegratorFunctionType = void(); + +/*! \internal + * \brief Struct to handle setting up and running the different "integrators". * - * \param[in] fplog Log file for output - * \param[in] cr Communication record - * \param[in] ms Handle to multi-simulation handler. - * \param[in] mdlog Log writer for important output - * \param[in] nfile Number of files - * \param[in] fnm Filename structure array - * \param[in] oenv Output information - * \param[in] mdrunOptions Options for mdrun - * \param[in] vsite Virtual site information - * \param[in] constr Constraint information - * \param[in] outputProvider Additional output provider - * \param[in] inputrec Input record with mdp options - * \param[in] top_global Molecular topology for the whole system - * \param[in] fcd Force and constraint data - * \param[in] state_global The state (x, v, f, box etc.) of the whole system - * \param[in] observablesHistory The observables statistics history - * \param[in] mdAtoms Atom information - * \param[in] nrnb Accounting for floating point operations - * \param[in] wcycle Wall cycle timing information - * \param[in] fr Force record with cut-off information and more - * \param[in] replExParams Parameters for the replica exchange algorithm - * \param[in] membed Membrane embedding data structure - * \param[in] walltime_accounting More timing information - */ -typedef double integrator_t (FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - ObservablesHistory *observablesHistory, - MDAtoms *mdatoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_membed_t gmx_unused * membed, - gmx_walltime_accounting_t walltime_accounting); + * This struct is a mere aggregate of parameters to pass to evaluate an + * energy, so that future changes to names and types of them consume + * less time when refactoring other code. + * + * Aggregate initialization is used, for which the chief risk is that + * if a member is added at the end and not all initializer lists are + * updated, then the member will be value initialized, which will + * typically mean initialization to zero. + * + * Having multiple integrators as member functions isn't a good + * design, and we definitely only intend one to be called, but the + * goal is to make it easy to change the names and types of members + * without having to make identical changes in several places in the + * code. Once many of them have become modules, we should change this + * approach. + * + * Note that the presence of const reference members means that the + * default constructor would be implicitly deleted. But we only want + * to make one of these when we know how to initialize these members, + * so that is perfect. To ensure this remains true even if we would + * remove those members, we explicitly delete this constructor. + * Other constructors, copies and moves are OK. */ +struct Integrator +{ + //! Handles logging. + FILE *fplog; + //! Handles communication. + t_commrec *cr; + //! Coordinates multi-simulations. + const gmx_multisim_t *ms; + //! Handles logging. + const MDLogger &mdlog; + //! Count of input file options. + int nfile; + //! Content of input file options. + const t_filenm *fnm; + //! Handles writing text output. + const gmx_output_env_t *oenv; + //! Contains command-line options to mdrun. + const MdrunOptions &mdrunOptions; + //! Handles virtual sites. + gmx_vsite_t *vsite; + //! Handles constraints. + gmx_constr *constr; + //! Handles writing output files. + IMDOutputProvider *outputProvider; + //! Contains user input mdp options. + t_inputrec *inputrec; + //! Full system topology. + gmx_mtop_t *top_global; + //! Helper struct for force calculations. + t_fcdata *fcd; + //! Full simulation state (only non-nullptr on master rank). + t_state *state_global; + //! History of simulation observables. + ObservablesHistory *observablesHistory; + //! Atom parameters for this domain. + MDAtoms *mdAtoms; + //! Manages flop accounting. + t_nrnb *nrnb; + //! Manages wall cycle accounting. + gmx_wallcycle *wcycle; + //! Parameters for force calculations. + t_forcerec *fr; + //! Parameters for replica exchange algorihtms. + const ReplicaExchangeParameters &replExParams; + //! Parameters for membrane embedding. + gmx_membed_t *membed; + //! Manages wall time accounting. + gmx_walltime_accounting *walltime_accounting; + //! Implements the normal MD integrators. + IntegratorFunctionType do_md; + //! Implements steepest descent EM. + IntegratorFunctionType do_steep; + //! Implements conjugate gradient energy minimization + IntegratorFunctionType do_cg; + //! Implements onjugate gradient energy minimization using the L-BFGS algorithm + IntegratorFunctionType do_lbfgs; + //! Implements normal mode analysis + IntegratorFunctionType do_nm; + //! Implements test particle insertion + IntegratorFunctionType do_tpi; + /*! \brief Function to run the correct IntegratorFunctionType, + * based on the .mdp integrator field. */ + void run(unsigned int ei); + //! We only intend to construct such objects with an initializer list. +#if __GNUC__ > 4 || (__GNUC__ == 4 && __GNUC_MINOR__ >= 9) + // Aspects of the C++11 spec changed after GCC 4.8.5, and + // compilation of the initializer list construction in runner.cpp + // fails in GCC 4.8.5. + Integrator() = delete; +#endif +}; } // namespace gmx diff --git a/src/gromacs/mdrun/md.cpp b/src/gromacs/mdrun/md.cpp index 4751be45d3..b48f2e229d 100644 --- a/src/gromacs/mdrun/md.cpp +++ b/src/gromacs/mdrun/md.cpp @@ -43,8 +43,6 @@ */ #include "gmxpre.h" -#include "md.h" - #include "config.h" #include @@ -136,6 +134,8 @@ #include "gromacs/utility/real.h" #include "gromacs/utility/smalloc.h" +#include "integrator.h" + #ifdef GMX_FAHCORE #include "corewrap.h" #endif @@ -260,46 +260,16 @@ static void prepareRerunState(const t_trxframe &rerunFrame, } } -/*! \libinternal - \copydoc integrator_t (FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - t_mdatoms *mdatoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_membed_t *membed, - gmx_walltime_accounting_t walltime_accounting) - */ -double gmx::do_md(FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *ir, - gmx_mtop_t *top_global, - t_fcdata *fcd, - t_state *state_global, - ObservablesHistory *observablesHistory, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_membed_t *membed, - gmx_walltime_accounting_t walltime_accounting) +void gmx::Integrator::do_md() { - gmx_mdoutf_t outf = nullptr; + // TODO Historically, the EM and MD "integrators" used different + // names for the t_inputrec *parameter, but these must have the + // same name, now that it's a member of a struct. We use this ir + // alias to avoid a large ripple of nearly useless changes. + // t_inputrec is being replaced by IMdpOptionsProvider, so this + // will go away eventually. + t_inputrec *ir = inputrec; + gmx_mdoutf *outf = nullptr; gmx_int64_t step, step_rel; double elapsed_time; double t, t0, lam0[efptNR]; @@ -2025,5 +1995,4 @@ double gmx::do_md(FILE *fplog, t_commrec *cr, mdAlgorithmsTearDownAtomData(fr->bonded_threading, top); fr->bonded_threading = nullptr; sfree(top); - return 0; } diff --git a/src/gromacs/mdrun/minimize.cpp b/src/gromacs/mdrun/minimize.cpp index e18f673949..0aef6dca81 100644 --- a/src/gromacs/mdrun/minimize.cpp +++ b/src/gromacs/mdrun/minimize.cpp @@ -44,8 +44,6 @@ */ #include "gmxpre.h" -#include "minimize.h" - #include "config.h" #include @@ -100,6 +98,8 @@ #include "gromacs/utility/logger.h" #include "gromacs/utility/smalloc.h" +#include "integrator.h" + //! Utility structure for manipulating states during EM typedef struct { //! Copy of the global state @@ -979,44 +979,8 @@ static real pr_beta(const t_commrec *cr, t_grpopts *opts, t_mdatoms *mdatoms, namespace gmx { -/*! \brief Do conjugate gradients minimization - \copydoc integrator_t(FILE *fplog, t_commrec *cr, - const gmx_multi_sim_t *, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - gmx_edsam_t ed, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) - */ -double do_cg(FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger gmx_unused &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t gmx_unused *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - ObservablesHistory *observablesHistory, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters gmx_unused &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) +void +Integrator::do_cg() { const char *CG = "Polak-Ribiere Conjugate Gradients"; @@ -1623,49 +1587,11 @@ double do_cg(FILE *fplog, t_commrec *cr, /* To print the actual number of steps we needed somewhere */ walltime_accounting_set_nsteps_done(walltime_accounting, step); +} - return 0; -} /* That's all folks */ - - -/*! \brief Do L-BFGS conjugate gradients minimization - \copydoc integrator_t(FILE *fplog, t_commrec *cr, - const gmx_multi_sim_t *, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - gmx_edsam_t ed, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) - */ -double do_lbfgs(FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger gmx_unused &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t gmx_unused *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - ObservablesHistory *observablesHistory, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters gmx_unused &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) + +void +Integrator::do_lbfgs() { static const char *LBFGS = "Low-Memory BFGS Minimizer"; em_state_t ems; @@ -2394,47 +2320,10 @@ double do_lbfgs(FILE *fplog, t_commrec *cr, /* To print the actual number of steps we needed somewhere */ walltime_accounting_set_nsteps_done(walltime_accounting, step); +} - return 0; -} /* That's all folks */ - -/*! \brief Do steepest descents minimization - \copydoc integrator_t(FILE *fplog, t_commrec *cr, - const gmx_multi_sim_t *, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - gmx_edsam_t ed, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_walltime_accounting_t walltime_accounting) - */ -double do_steep(FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger gmx_unused &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t gmx_unused *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - ObservablesHistory *observablesHistory, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters gmx_unused &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) +void +Integrator::do_steep() { const char *SD = "Steepest Descents"; gmx_localtop_t *top; @@ -2664,47 +2553,10 @@ double do_steep(FILE *fplog, t_commrec *cr, inputrec->nsteps = count; walltime_accounting_set_nsteps_done(walltime_accounting, count); +} - return 0; -} /* That's all folks */ - -/*! \brief Do normal modes analysis - \copydoc integrator_t(FILE *fplog, t_commrec *cr, - const gmx_multi_sim_t *, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - gmx_edsam_t ed, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_walltime_accounting_t walltime_accounting) - */ -double do_nm(FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t gmx_unused *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - ObservablesHistory gmx_unused *observablesHistory, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters gmx_unused &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) +void +Integrator::do_nm() { const char *NM = "Normal Mode Analysis"; gmx_mdoutf_t outf; @@ -2991,8 +2843,6 @@ double do_nm(FILE *fplog, t_commrec *cr, finish_em(cr, outf, walltime_accounting, wcycle); walltime_accounting_set_nsteps_done(walltime_accounting, atom_index.size()*2); - - return 0; } } // namespace gmx diff --git a/src/gromacs/mdrun/minimize.h b/src/gromacs/mdrun/minimize.h index 00ed45ffd7..61785396fc 100644 --- a/src/gromacs/mdrun/minimize.h +++ b/src/gromacs/mdrun/minimize.h @@ -45,18 +45,6 @@ namespace gmx { -//! Steepest descents energy minimization -integrator_t do_steep; - -//! Conjugate gradient energy minimization -integrator_t do_cg; - -//! Conjugate gradient energy minimization using the L-BFGS algorithm -integrator_t do_lbfgs; - -//! Normal mode analysis -integrator_t do_nm; - } // namespace gmx #endif // GMX_MDRUN_MINIMIZE_H diff --git a/src/gromacs/mdrun/runner.cpp b/src/gromacs/mdrun/runner.cpp index a35a94d96c..c385359706 100644 --- a/src/gromacs/mdrun/runner.cpp +++ b/src/gromacs/mdrun/runner.cpp @@ -63,6 +63,7 @@ #include "gromacs/fileio/oenv.h" #include "gromacs/fileio/tpxio.h" #include "gromacs/gmxlib/network.h" +#include "gromacs/gmxlib/nrnb.h" #include "gromacs/gpu_utils/gpu_utils.h" #include "gromacs/hardware/cpuinfo.h" #include "gromacs/hardware/detecthardware.h" @@ -124,9 +125,6 @@ #include "gromacs/utility/stringutil.h" #include "integrator.h" -#include "md.h" -#include "minimize.h" -#include "tpi.h" #ifdef GMX_FAHCORE #include "corewrap.h" @@ -341,43 +339,6 @@ static bool gpuAccelerationOfNonbondedIsUseful(const MDLogger &mdlog, return true; } -//! \brief Return the correct integrator function. -static integrator_t *my_integrator(unsigned int ei) -{ - switch (ei) - { - case eiMD: - case eiBD: - case eiSD1: - case eiVV: - case eiVVAK: - if (!EI_DYNAMICS(ei)) - { - GMX_THROW(APIError("do_md integrator would be called for a non-dynamical integrator")); - } - return do_md; - case eiSteep: - return do_steep; - case eiCG: - return do_cg; - case eiNM: - return do_nm; - case eiLBFGS: - return do_lbfgs; - case eiTPI: - case eiTPIC: - if (!EI_TPI(ei)) - { - GMX_THROW(APIError("do_tpi integrator would be called for a non-TPI integrator")); - } - return do_tpi; - case eiSD2_REMOVED: - GMX_THROW(NotImplementedError("SD2 integrator has been removed")); - default: - GMX_THROW(APIError("Non existing integrator selected")); - } -} - //! Initializes the logger for mdrun. static gmx::LoggerOwner buildLogger(FILE *fplog, const t_commrec *cr) { @@ -1339,20 +1300,22 @@ int Mdrunner::mdrunner() } /* Now do whatever the user wants us to do (how flexible...) */ - my_integrator(inputrec->eI) (fplog, cr, ms, mdlog, nfile, fnm, - oenv, - mdrunOptions, - vsite, constr, - mdModules->outputProvider(), - inputrec, &mtop, - fcd, - globalState.get(), - &observablesHistory, - mdAtoms.get(), nrnb, wcycle, fr, - replExParams, - membed, - walltime_accounting); - + Integrator integrator { + fplog, cr, ms, mdlog, nfile, fnm, + oenv, + mdrunOptions, + vsite, constr, + mdModules->outputProvider(), + inputrec, &mtop, + fcd, + globalState.get(), + &observablesHistory, + mdAtoms.get(), nrnb, wcycle, fr, + replExParams, + membed, + walltime_accounting + }; + integrator.run(inputrec->eI); if (inputrec->bRot) { finish_rot(inputrec->rot); diff --git a/src/gromacs/mdrun/tpi.cpp b/src/gromacs/mdrun/tpi.cpp index 54a202f500..dd4fd7d329 100644 --- a/src/gromacs/mdrun/tpi.cpp +++ b/src/gromacs/mdrun/tpi.cpp @@ -43,8 +43,6 @@ */ #include "gmxpre.h" -#include "tpi.h" - #include #include #include @@ -91,6 +89,8 @@ #include "gromacs/utility/gmxassert.h" #include "gromacs/utility/smalloc.h" +#include "integrator.h" + //! Global max algorithm static void global_max(t_commrec *cr, int *n) { @@ -126,44 +126,8 @@ static void realloc_bins(double **bin, int *nbin, int nbin_new) namespace gmx { -/*! \brief Do test particle insertion. - \copydoc integrator_t (FILE *fplog, t_commrec *cr, - const gmx_multi_sim_t *, - const gmx::MDLogger &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t *vsite, gmx_constr_t constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - gmx_edsam_t ed, - t_forcerec *fr, - const ReplicaExchangeParameters &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) - */ -double do_tpi(FILE *fplog, t_commrec *cr, - const gmx_multisim_t *ms, - const gmx::MDLogger gmx_unused &mdlog, - int nfile, const t_filenm fnm[], - const gmx_output_env_t *oenv, - const MdrunOptions &mdrunOptions, - gmx_vsite_t gmx_unused *vsite, gmx_constr_t gmx_unused constr, - gmx::IMDOutputProvider *outputProvider, - t_inputrec *inputrec, - gmx_mtop_t *top_global, t_fcdata *fcd, - t_state *state_global, - ObservablesHistory gmx_unused *observablesHistory, - gmx::MDAtoms *mdAtoms, - t_nrnb *nrnb, gmx_wallcycle_t wcycle, - t_forcerec *fr, - const ReplicaExchangeParameters gmx_unused &replExParams, - gmx_membed_t gmx_unused *membed, - gmx_walltime_accounting_t walltime_accounting) +void +Integrator::do_tpi() { gmx_localtop_t *top; gmx_groups_t *groups; @@ -892,8 +856,6 @@ double do_tpi(FILE *fplog, t_commrec *cr, sfree(sum_UgembU); walltime_accounting_set_nsteps_done(walltime_accounting, frame*inputrec->nsteps); - - return 0; } } // namespace gmx -- 2.11.4.GIT