From df3b598f6d15d1e2e70e1c6c3bf5debc49069c47 Mon Sep 17 00:00:00 2001 From: "M. Eric Irrgang" Date: Sun, 14 Oct 2018 15:33:48 +0300 Subject: [PATCH] Preparation for de-duplicating mdrun setup Recent work introducing the API extracted MdFilenames and then essentially duplicated everything else. This commit expands MdFilenames into LegacyMdrunOptions, and clarified the plan that all such code needs to become either module-based option handling, or aspects of SimulationContext suitable for both CLI and API. This commit introduces the Legacy MdrunOptions struct, so that a future commit can make it available to the API setup for an mdrun session. The layout of the changes is intended to minimize the textual diff for the CLI. Code motion for clarity and also use by the API will follow. Removed MdFilenames, made its data a std::vector in LegacyMdrunOptions, and used it via an ArrayRef. This clarifies that it is read-only input data supplied by the client. Change-Id: I3438cfe1ab52eeb2fffbd16f4bd5254b2a74f03f --- src/api/cpp/context-impl.h | 7 + src/api/cpp/context.cpp | 6 +- src/gromacs/mdrun/legacyoptions.h | 278 ++++++++++++++++++++++++++++++ src/gromacs/mdrun/runner.cpp | 47 +++--- src/gromacs/mdrun/runner.h | 12 +- src/gromacs/mdrun/simulationcontext.h | 8 +- src/programs/legacymodules.cpp | 2 +- src/programs/mdrun/mdrun.cpp | 306 +++++++++------------------------- src/programs/mdrun/mdrun_main.h | 7 +- 9 files changed, 403 insertions(+), 270 deletions(-) create mode 100644 src/gromacs/mdrun/legacyoptions.h diff --git a/src/api/cpp/context-impl.h b/src/api/cpp/context-impl.h index a3f4a82d3f..c34ebfe33c 100644 --- a/src/api/cpp/context-impl.h +++ b/src/api/cpp/context-impl.h @@ -44,6 +44,8 @@ #include #include +#include "gromacs/mdrun/mdfilenames.h" + #include "gmxapi/context.h" #include "gmxapi/session.h" @@ -138,6 +140,11 @@ class ContextImpl final : public std::enable_shared_from_this * evolves. */ MDArgs mdArgs_; + + /*! + * \brief Ensure lifetime of MD filenames options object exceeds that of session. + */ + gmx::MdFilenames mdFilenames_; }; } // end namespace gmxapi diff --git a/src/api/cpp/context.cpp b/src/api/cpp/context.cpp index 3ae3e637ff..09a2f03da3 100644 --- a/src/api/cpp/context.cpp +++ b/src/api/cpp/context.cpp @@ -58,6 +58,7 @@ #include "gromacs/gmxlib/network.h" #include "gromacs/mdlib/stophandler.h" #include "gromacs/mdrun/logging.h" +#include "gromacs/mdrun/mdfilenames.h" #include "gromacs/mdrun/multisim.h" #include "gromacs/mdrun/runner.h" #include "gromacs/mdrunutility/handlerestart.h" @@ -178,7 +179,8 @@ std::shared_ptr ContextImpl::launch(const Workflow &work) ReplicaExchangeParameters replExParams; // Filenames and properties from command-line argument values. - auto filenames = gmx::MdFilenames(); + mdFilenames_ = gmx::MdFilenames(); + auto &filenames = mdFilenames_; // Print a warning if any force is larger than this (in kJ/mol nm). real pforce = -1; @@ -491,7 +493,7 @@ std::shared_ptr ContextImpl::launch(const Workflow &work) // Need to establish run-time values from various inputs to provide a resource handle to Mdrunner builder.addHardwareOptions(hw_opt); // \todo File names are parameters that should be managed modularly through further factoring. - builder.addFilenames(filenames); + builder.addFilenames(filenames()); // Note: The gmx_output_env_t life time is not managed after the call to parse_common_args. // \todo Implement lifetime management for gmx_output_env_t. // \todo Output environment should be configured outside of Mdrunner and provided as a resource. diff --git a/src/gromacs/mdrun/legacyoptions.h b/src/gromacs/mdrun/legacyoptions.h new file mode 100644 index 0000000000..742e04a873 --- /dev/null +++ b/src/gromacs/mdrun/legacyoptions.h @@ -0,0 +1,278 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 1991-2000, University of Groningen, The Netherlands. + * Copyright (c) 2001-2004, The GROMACS development team. + * Copyright (c) 2011,2012,2013,2014,2015,2016,2017,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. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \libinternal \file + * + * \brief This file declares helper functionality for legacy option handling for mdrun + * + * \author Berk Hess + * \author David van der Spoel + * \author Erik Lindahl + * \author Mark Abraham + * + * \ingroup module_mdrun + * \inlibraryapi + */ +#include "gromacs/commandline/filenm.h" +#include "gromacs/commandline/pargs.h" +#include "gromacs/domdec/domdec.h" +#include "gromacs/hardware/hw_info.h" +#include "gromacs/mdlib/mdrun.h" +#include "gromacs/mdrun/logging.h" + +#include "replicaexchange.h" + +struct gmx_multisim_t; + +namespace gmx +{ + +/*! \libinternal + * \brief This class provides the same command-line option + * functionality to both CLI and API sessions. + * + * This class should not exist, but is necessary now to introduce + * support for the CLI and API without duplicating code. It should be + * eliminated following the TODOs below. + * + * \todo Modules in mdrun should acquire proper option handling so + * that all of these declarations and defaults are local to the + * modules. + * + * \todo Contextual aspects, such as working directory, MPI + * environment, and environment variable handling are more properly + * the role of SimulationContext, and should be moved there */ +class LegacyMdrunOptions +{ + public: + //! Ongoing collection of mdrun options + MdrunOptions mdrunOptions; + //! Options for the domain decomposition. + DomdecOptions domdecOptions; + //! Parallelism-related user options. + gmx_hw_opt_t hw_opt; + //! Command-line override for the duration of a neighbor list with the Verlet scheme. + int nstlist_cmdline = 0; + //! Parameters for replica-exchange simulations. + ReplicaExchangeParameters replExParams; + + //! Filename options to fill from command-line argument values. + std::vector filenames = + {{{ efTPR, nullptr, nullptr, ffREAD }, + { efTRN, "-o", nullptr, ffWRITE }, + { efCOMPRESSED, "-x", nullptr, ffOPTWR }, + { efCPT, "-cpi", nullptr, ffOPTRD | ffALLOW_MISSING }, + { efCPT, "-cpo", nullptr, ffOPTWR }, + { efSTO, "-c", "confout", ffWRITE }, + { efEDR, "-e", "ener", ffWRITE }, + { efLOG, "-g", "md", ffWRITE }, + { efXVG, "-dhdl", "dhdl", ffOPTWR }, + { efXVG, "-field", "field", ffOPTWR }, + { efXVG, "-table", "table", ffOPTRD }, + { efXVG, "-tablep", "tablep", ffOPTRD }, + { efXVG, "-tableb", "table", ffOPTRDMULT }, + { efTRX, "-rerun", "rerun", ffOPTRD }, + { efXVG, "-tpi", "tpi", ffOPTWR }, + { efXVG, "-tpid", "tpidist", ffOPTWR }, + { efEDI, "-ei", "sam", ffOPTRD }, + { efXVG, "-eo", "edsam", ffOPTWR }, + { efXVG, "-devout", "deviatie", ffOPTWR }, + { efXVG, "-runav", "runaver", ffOPTWR }, + { efXVG, "-px", "pullx", ffOPTWR }, + { efXVG, "-pf", "pullf", ffOPTWR }, + { efXVG, "-ro", "rotation", ffOPTWR }, + { efLOG, "-ra", "rotangles", ffOPTWR }, + { efLOG, "-rs", "rotslabs", ffOPTWR }, + { efLOG, "-rt", "rottorque", ffOPTWR }, + { efMTX, "-mtx", "nm", ffOPTWR }, + { efRND, "-multidir", nullptr, ffOPTRDMULT}, + { efXVG, "-awh", "awhinit", ffOPTRD }, + { efDAT, "-membed", "membed", ffOPTRD }, + { efTOP, "-mp", "membed", ffOPTRD }, + { efNDX, "-mn", "membed", ffOPTRD }, + { efXVG, "-if", "imdforces", ffOPTWR }, + { efXVG, "-swap", "swapions", ffOPTWR }}}; + + //! Print a warning if any force is larger than this (in kJ/mol nm). + real pforce = -1; + + /*! \brief Output context for writing text files + * + * \todo Clarify initialization, ownership, and lifetime. */ + gmx_output_env_t *oenv = nullptr; + + //! Handle to file used for logging. + LogFilePtr logFileGuard = nullptr; + + /*! \brief Command line options, defaults, docs and storage for them to fill. */ + /*! \{ */ + rvec realddxyz = {0, 0, 0}; + const char *ddrank_opt_choices[static_cast(DdRankOrder::nr)+1] = + { nullptr, "interleave", "pp_pme", "cartesian", nullptr }; + const char *dddlb_opt_choices[static_cast(DlbOption::nr)+1] = + { nullptr, "auto", "no", "yes", nullptr }; + const char *thread_aff_opt_choices[threadaffNR+1] = + { nullptr, "auto", "on", "off", nullptr }; + const char *nbpu_opt_choices[5] = + { nullptr, "auto", "cpu", "gpu", nullptr }; + const char *pme_opt_choices[5] = + { nullptr, "auto", "cpu", "gpu", nullptr }; + const char *pme_fft_opt_choices[5] = + { nullptr, "auto", "cpu", "gpu", nullptr }; + gmx_bool bTryToAppendFiles = TRUE; + const char *gpuIdsAvailable = ""; + const char *userGpuTaskAssignment = ""; + + ImdOptions &imdOptions = mdrunOptions.imdOptions; + + t_pargs pa[47] = { + + { "-dd", FALSE, etRVEC, {&realddxyz}, + "Domain decomposition grid, 0 is optimize" }, + { "-ddorder", FALSE, etENUM, {ddrank_opt_choices}, + "DD rank order" }, + { "-npme", FALSE, etINT, {&domdecOptions.numPmeRanks}, + "Number of separate ranks to be used for PME, -1 is guess" }, + { "-nt", FALSE, etINT, {&hw_opt.nthreads_tot}, + "Total number of threads to start (0 is guess)" }, + { "-ntmpi", FALSE, etINT, {&hw_opt.nthreads_tmpi}, + "Number of thread-MPI ranks to start (0 is guess)" }, + { "-ntomp", FALSE, etINT, {&hw_opt.nthreads_omp}, + "Number of OpenMP threads per MPI rank to start (0 is guess)" }, + { "-ntomp_pme", FALSE, etINT, {&hw_opt.nthreads_omp_pme}, + "Number of OpenMP threads per MPI rank to start (0 is -ntomp)" }, + { "-pin", FALSE, etENUM, {thread_aff_opt_choices}, + "Whether mdrun should try to set thread affinities" }, + { "-pinoffset", FALSE, etINT, {&hw_opt.core_pinning_offset}, + "The lowest logical core number to which mdrun should pin the first thread" }, + { "-pinstride", FALSE, etINT, {&hw_opt.core_pinning_stride}, + "Pinning distance in logical cores for threads, use 0 to minimize the number of threads per physical core" }, + { "-gpu_id", FALSE, etSTR, {&gpuIdsAvailable}, + "List of unique GPU device IDs available to use" }, + { "-gputasks", FALSE, etSTR, {&userGpuTaskAssignment}, + "List of GPU device IDs, mapping each PP task on each node to a device" }, + { "-ddcheck", FALSE, etBOOL, {&domdecOptions.checkBondedInteractions}, + "Check for all bonded interactions with DD" }, + { "-ddbondcomm", FALSE, etBOOL, {&domdecOptions.useBondedCommunication}, + "HIDDENUse special bonded atom communication when [TT]-rdd[tt] > cut-off" }, + { "-rdd", FALSE, etREAL, {&domdecOptions.minimumCommunicationRange}, + "The maximum distance for bonded interactions with DD (nm), 0 is determine from initial coordinates" }, + { "-rcon", FALSE, etREAL, {&domdecOptions.constraintCommunicationRange}, + "Maximum distance for P-LINCS (nm), 0 is estimate" }, + { "-dlb", FALSE, etENUM, {dddlb_opt_choices}, + "Dynamic load balancing (with DD)" }, + { "-dds", FALSE, etREAL, {&domdecOptions.dlbScaling}, + "Fraction in (0,1) by whose reciprocal the initial DD cell size will be increased in order to " + "provide a margin in which dynamic load balancing can act while preserving the minimum cell size." }, + { "-ddcsx", FALSE, etSTR, {&domdecOptions.cellSizeX}, + "HIDDENA string containing a vector of the relative sizes in the x " + "direction of the corresponding DD cells. Only effective with static " + "load balancing." }, + { "-ddcsy", FALSE, etSTR, {&domdecOptions.cellSizeY}, + "HIDDENA string containing a vector of the relative sizes in the y " + "direction of the corresponding DD cells. Only effective with static " + "load balancing." }, + { "-ddcsz", FALSE, etSTR, {&domdecOptions.cellSizeZ}, + "HIDDENA string containing a vector of the relative sizes in the z " + "direction of the corresponding DD cells. Only effective with static " + "load balancing." }, + { "-gcom", FALSE, etINT, {&mdrunOptions.globalCommunicationInterval}, + "Global communication frequency" }, + { "-nb", FALSE, etENUM, {nbpu_opt_choices}, + "Calculate non-bonded interactions on" }, + { "-nstlist", FALSE, etINT, {&nstlist_cmdline}, + "Set nstlist when using a Verlet buffer tolerance (0 is guess)" }, + { "-tunepme", FALSE, etBOOL, {&mdrunOptions.tunePme}, + "Optimize PME load between PP/PME ranks or GPU/CPU (only with the Verlet cut-off scheme)" }, + { "-pme", FALSE, etENUM, {pme_opt_choices}, + "Perform PME calculations on" }, + { "-pmefft", FALSE, etENUM, {pme_fft_opt_choices}, + "Perform PME FFT calculations on" }, + { "-v", FALSE, etBOOL, {&mdrunOptions.verbose}, + "Be loud and noisy" }, + { "-pforce", FALSE, etREAL, {&pforce}, + "Print all forces larger than this (kJ/mol nm)" }, + { "-reprod", FALSE, etBOOL, {&mdrunOptions.reproducible}, + "Try to avoid optimizations that affect binary reproducibility" }, + { "-cpt", FALSE, etREAL, {&mdrunOptions.checkpointOptions.period}, + "Checkpoint interval (minutes)" }, + { "-cpnum", FALSE, etBOOL, {&mdrunOptions.checkpointOptions.keepAndNumberCheckpointFiles}, + "Keep and number checkpoint files" }, + { "-append", FALSE, etBOOL, {&bTryToAppendFiles}, + "Append to previous output files when continuing from checkpoint instead of adding the simulation part number to all file names" }, + { "-nsteps", FALSE, etINT64, {&mdrunOptions.numStepsCommandline}, + "Run this number of steps, overrides .mdp file option (-1 means infinite, -2 means use mdp option, smaller is invalid)" }, + { "-maxh", FALSE, etREAL, {&mdrunOptions.maximumHoursToRun}, + "Terminate after 0.99 times this time (hours)" }, + { "-replex", FALSE, etINT, {&replExParams.exchangeInterval}, + "Attempt replica exchange periodically with this period (steps)" }, + { "-nex", FALSE, etINT, {&replExParams.numExchanges}, + "Number of random exchanges to carry out each exchange interval (N^3 is one suggestion). -nex zero or not specified gives neighbor replica exchange." }, + { "-reseed", FALSE, etINT, {&replExParams.randomSeed}, + "Seed for replica exchange, -1 is generate a seed" }, + { "-imdport", FALSE, etINT, {&imdOptions.port}, + "HIDDENIMD listening port" }, + { "-imdwait", FALSE, etBOOL, {&imdOptions.wait}, + "HIDDENPause the simulation while no IMD client is connected" }, + { "-imdterm", FALSE, etBOOL, {&imdOptions.terminatable}, + "HIDDENAllow termination of the simulation from IMD client" }, + { "-imdpull", FALSE, etBOOL, {&imdOptions.pull}, + "HIDDENAllow pulling in the simulation from IMD client" }, + { "-rerunvsite", FALSE, etBOOL, {&mdrunOptions.rerunConstructVsites}, + "HIDDENRecalculate virtual site coordinates with [TT]-rerun[tt]" }, + { "-confout", FALSE, etBOOL, {&mdrunOptions.writeConfout}, + "HIDDENWrite the last configuration with [TT]-c[tt] and force checkpointing at the last step" }, + { "-stepout", FALSE, etINT, {&mdrunOptions.verboseStepPrintInterval}, + "HIDDENFrequency of writing the remaining wall clock time for the run" }, + { "-resetstep", FALSE, etINT, {&mdrunOptions.timingOptions.resetStep}, + "HIDDENReset cycle counters after these many time steps" }, + { "-resethway", FALSE, etBOOL, {&mdrunOptions.timingOptions.resetHalfway}, + "HIDDENReset the cycle counters after half the number of steps or halfway [TT]-maxh[tt]" } + }; + /*! \} */ + + //! Handle to communication object. + t_commrec *cr = nullptr; + //! Multi-simulation object. + gmx_multisim_t *ms = nullptr; + + //! Parses the command-line input and prepares to start mdrun. + int updateFromCommandLine(int argc, char **argv, ArrayRef desc); + + ~LegacyMdrunOptions(); +}; + +} // end namespace gmx diff --git a/src/gromacs/mdrun/runner.cpp b/src/gromacs/mdrun/runner.cpp index f7080cca45..8634958ced 100644 --- a/src/gromacs/mdrun/runner.cpp +++ b/src/gromacs/mdrun/runner.cpp @@ -96,6 +96,7 @@ #include "gromacs/mdlib/sighandler.h" #include "gromacs/mdlib/sim_util.h" #include "gromacs/mdlib/stophandler.h" +#include "gromacs/mdrun/legacyoptions.h" #include "gromacs/mdrun/logging.h" #include "gromacs/mdrun/multisim.h" #include "gromacs/mdrun/simulationcontext.h" @@ -437,7 +438,7 @@ int Mdrunner::mdrunner() t_inputrec *inputrec = &inputrecInstance; gmx_mtop_t mtop; - bool doMembed = opt2bSet("-membed", filenames().size(), filenames().data()); + bool doMembed = opt2bSet("-membed", filenames.size(), filenames.data()); bool doRerun = mdrunOptions.rerun; // Handle task-assignment related user options. @@ -550,7 +551,7 @@ int Mdrunner::mdrunner() globalState = compat::make_unique(); /* Read (nearly) all data required for the simulation */ - read_tpx_state(ftp2fn(efTPR, filenames().size(), filenames().data()), inputrec, globalState.get(), &mtop); + read_tpx_state(ftp2fn(efTPR, filenames.size(), filenames.data()), inputrec, globalState.get(), &mtop); /* In rerun, set velocities to zero if present */ if (doRerun && ((globalState->flags & (1 << estV)) != 0)) @@ -834,7 +835,7 @@ int Mdrunner::mdrunner() */ gmx_bool bReadEkin; - load_checkpoint(opt2fn_master("-cpi", filenames().size(), filenames().data(), cr), + load_checkpoint(opt2fn_master("-cpi", filenames.size(), filenames.data(), cr), logFileHandle, cr, domdecOptions.numCells, inputrec, globalState.get(), @@ -1177,7 +1178,7 @@ int Mdrunner::mdrunner() /* Note that membed cannot work in parallel because mtop is * changed here. Fix this if we ever want to make it run with * multiple ranks. */ - membed = init_membed(fplog, filenames().size(), filenames().data(), &mtop, inputrec, globalState.get(), cr, + membed = init_membed(fplog, filenames.size(), filenames.data(), &mtop, inputrec, globalState.get(), cr, &mdrunOptions .checkpointOptions.period); } @@ -1193,9 +1194,9 @@ int Mdrunner::mdrunner() fr->forceProviders = mdModules->initForceProviders(); init_forcerec(fplog, mdlog, fr, fcd, inputrec, &mtop, cr, box, - opt2fn("-table", filenames().size(), filenames().data()), - opt2fn("-tablep", filenames().size(), filenames().data()), - opt2fns("-tableb", filenames().size(), filenames().data()), + opt2fn("-table", filenames.size(), filenames.data()), + opt2fn("-tablep", filenames.size(), filenames.data()), + opt2fns("-tableb", filenames.size(), filenames.data()), *hwinfo, nonbondedDeviceInfo, FALSE, pforce); @@ -1334,7 +1335,7 @@ int Mdrunner::mdrunner() if (EI_DYNAMICS(inputrec->eI) && MASTER(cr)) { init_pull_output_files(inputrec->pull_work, - filenames().size(), filenames().data(), oenv, + filenames.size(), filenames.data(), oenv, continuationOptions); } } @@ -1345,8 +1346,8 @@ int Mdrunner::mdrunner() /* Initialize enforced rotation code */ enforcedRotation = init_rot(fplog, inputrec, - filenames().size(), - filenames().data(), + filenames.size(), + filenames.data(), cr, &atomSets, globalState.get(), @@ -1358,7 +1359,7 @@ int Mdrunner::mdrunner() if (inputrec->eSwapCoords != eswapNO) { /* Initialize ion swapping code */ - init_swapcoords(fplog, inputrec, opt2fn_master("-swap", filenames().size(), filenames().data(), cr), + init_swapcoords(fplog, inputrec, opt2fn_master("-swap", filenames.size(), filenames.data(), cr), &mtop, globalState.get(), &observablesHistory, cr, &atomSets, oenv, mdrunOptions); } @@ -1366,7 +1367,7 @@ int Mdrunner::mdrunner() /* Let makeConstraints know whether we have essential dynamics constraints. * TODO: inputrec should tell us whether we use an algorithm, not a file option or the checkpoint */ - bool doEssentialDynamics = (opt2fn_null("-ei", filenames().size(), filenames().data()) != nullptr + bool doEssentialDynamics = (opt2fn_null("-ei", filenames.size(), filenames.data()) != nullptr || observablesHistory.edsamHistory); auto constr = makeConstraints(mtop, *inputrec, doEssentialDynamics, fplog, *mdAtoms->mdatoms(), @@ -1386,7 +1387,7 @@ int Mdrunner::mdrunner() GMX_ASSERT(stopHandlerBuilder_, "Runner must provide StopHandlerBuilder to integrator."); /* Now do whatever the user wants us to do (how flexible...) */ Integrator integrator { - fplog, cr, ms, mdlog, static_cast(filenames().size()), filenames().data(), + fplog, cr, ms, mdlog, static_cast(filenames.size()), filenames.data(), oenv, mdrunOptions, vsite.get(), constr.get(), @@ -1550,7 +1551,7 @@ class Mdrunner::BuilderImplementation void addHardwareOptions(const gmx_hw_opt_t &hardwareOptions); - void addFilenames(const MdFilenames &filenames); + void addFilenames(ArrayRef filenames); void addOutputEnvironment(gmx_output_env_t* outputEnvironment); @@ -1598,7 +1599,7 @@ class Mdrunner::BuilderImplementation gmx_hw_opt_t hardwareOptions_; //! filename options for simulation. - std::unique_ptr filenames_ = nullptr; + ArrayRef filenames_; /*! \brief Handle to output environment. * @@ -1676,15 +1677,7 @@ Mdrunner Mdrunner::BuilderImplementation::build() newRunner.replExParams = replicaExchangeParameters_; - // \todo determine an invariant to checkout or establish that all MdFilenames objects are valid - if (filenames_) - { - newRunner.filenames = *filenames_; - } - else - { - GMX_THROW(gmx::APIError("MdrunnerBuilder::addFilenames() is required before build()")); - } + newRunner.filenames = filenames_; GMX_ASSERT(context_->communicationRecord_, "SimulationContext communications not initialized."); newRunner.cr = context_->communicationRecord_; @@ -1763,9 +1756,9 @@ void Mdrunner::BuilderImplementation::addHardwareOptions(const gmx_hw_opt_t &har hardwareOptions_ = hardwareOptions; } -void Mdrunner::BuilderImplementation::addFilenames(const MdFilenames &filenames) +void Mdrunner::BuilderImplementation::addFilenames(ArrayRef filenames) { - filenames_ = compat::make_unique(filenames); + filenames_ = filenames; } void Mdrunner::BuilderImplementation::addOutputEnvironment(gmx_output_env_t* outputEnvironment) @@ -1855,7 +1848,7 @@ MdrunnerBuilder &MdrunnerBuilder::addHardwareOptions(const gmx_hw_opt_t &hardwar return *this; } -MdrunnerBuilder &MdrunnerBuilder::addFilenames(const MdFilenames &filenames) +MdrunnerBuilder &MdrunnerBuilder::addFilenames(ArrayRef filenames) { impl_->addFilenames(filenames); return *this; diff --git a/src/gromacs/mdrun/runner.h b/src/gromacs/mdrun/runner.h index 4c981e6d24..75ee32f51b 100644 --- a/src/gromacs/mdrun/runner.h +++ b/src/gromacs/mdrun/runner.h @@ -47,12 +47,12 @@ #include #include +#include "gromacs/commandline/filenm.h" #include "gromacs/compat/pointers.h" #include "gromacs/domdec/domdec.h" #include "gromacs/hardware/hw_info.h" #include "gromacs/math/vec.h" #include "gromacs/mdlib/mdrun.h" -#include "gromacs/mdrun/mdfilenames.h" #include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/real.h" @@ -200,7 +200,7 @@ class Mdrunner gmx_hw_opt_t hw_opt; //! Filenames and properties from command-line argument values. - MdFilenames filenames; + ArrayRef filenames; /*! \brief Output context for writing text files * @@ -484,17 +484,17 @@ class MdrunnerBuilder final MdrunnerBuilder &addHardwareOptions(const gmx_hw_opt_t &hardwareOptions); /*! - * \brief Provide the filenames options structure. + * \brief Provide the filenames options structure with option values chosen * - * Required. The validity of a default-constructed MdFilenames object is unspecified. - * The object is assumed to have been updated by parse_common_args. + * Required. The object is assumed to have been updated by + * parse_common_args or equivalent. * * \param filenames Filenames and properties from command-line argument values or defaults. * * \internal * \todo Modules should manage their own filename options and defaults. */ - MdrunnerBuilder &addFilenames(const MdFilenames &filenames); + MdrunnerBuilder &addFilenames(ArrayRef filenames); /*! * \brief Provide parameters for setting up output environment. diff --git a/src/gromacs/mdrun/simulationcontext.h b/src/gromacs/mdrun/simulationcontext.h index 0f635be71c..f80033113c 100644 --- a/src/gromacs/mdrun/simulationcontext.h +++ b/src/gromacs/mdrun/simulationcontext.h @@ -50,7 +50,6 @@ #include #include "gromacs/hardware/hw_info.h" -#include "gromacs/mdrun/mdfilenames.h" struct t_filenm; struct t_commrec; @@ -70,6 +69,9 @@ namespace gmx * The public interface of SimulationContext is not yet well-specified. * Client code can create an instance with gmx::createSimulationContext() * + * \todo This class should also handle aspects of simulation + * environment such as working directory and environment variables. + * * \ingroup module_mdrun * \inlibraryapi * @@ -78,7 +80,7 @@ namespace gmx * Interfaces for different API levels are not yet final, but also depend on * additional development of t_commrec and other resources. * \todo Impose sensible access restrictions. - * Either the Context should be passed to the Runner as logically constant or + * Either the SimulationContext should be passed to the Mdrunner as logically constant or * a separate handle class can provide access to resources that have been * allocated by (negotiated with) the client for the current simulation * (or simulation segment). @@ -143,8 +145,6 @@ class SimulationContext final * * \param simulationCommunicator Handle to communication data structure. * - * \todo move filenames to separate builder method. - * * \ingroup module_mdrun */ SimulationContext createSimulationContext(t_commrec* simulationCommunicator); diff --git a/src/programs/legacymodules.cpp b/src/programs/legacymodules.cpp index b8ec535f0a..ead9b6d4e5 100644 --- a/src/programs/legacymodules.cpp +++ b/src/programs/legacymodules.cpp @@ -182,7 +182,7 @@ void registerLegacyModules(gmx::CommandLineModuleManager *manager) registerModule(manager, &gmx_x2top, "x2top", "Generate a primitive topology from coordinates"); - registerModuleNoNice(manager, &gmx_mdrun, "mdrun", + registerModuleNoNice(manager, &gmx::gmx_mdrun, "mdrun", "Perform a simulation, do a normal mode analysis or an energy minimization"); gmx::ICommandLineOptionsModule::registerModuleFactory( diff --git a/src/programs/mdrun/mdrun.cpp b/src/programs/mdrun/mdrun.cpp index 40f9ce0405..7c1def5201 100644 --- a/src/programs/mdrun/mdrun.cpp +++ b/src/programs/mdrun/mdrun.cpp @@ -64,6 +64,7 @@ #include "gromacs/fileio/gmxfio.h" #include "gromacs/gmxlib/network.h" #include "gromacs/mdlib/mdrun.h" +#include "gromacs/mdrun/legacyoptions.h" #include "gromacs/mdrun/logging.h" #include "gromacs/mdrun/multisim.h" #include "gromacs/mdrun/replicaexchange.h" @@ -71,12 +72,16 @@ #include "gromacs/mdrun/simulationcontext.h" #include "gromacs/mdrunutility/handlerestart.h" #include "gromacs/mdtypes/commrec.h" +#include "gromacs/utility/arrayref.h" #include "gromacs/utility/arraysize.h" #include "gromacs/utility/fatalerror.h" #include "gromacs/utility/smalloc.h" #include "mdrun_main.h" +namespace gmx +{ + /*! \brief Return whether the command-line parameter that * will trigger a multi-simulation is set */ static bool is_multisim_option_set(int argc, const char *const argv[]) @@ -91,11 +96,14 @@ static bool is_multisim_option_set(int argc, const char *const argv[]) return false; } +//! Helper function. +static int launchMdrun(LegacyMdrunOptions &&options, + compat::not_null context); + //! Implements C-style main function for mdrun int gmx_mdrun(int argc, char *argv[]) { - using namespace gmx; - const char *desc[] = { + std::vectordesc = { "[THISMODULE] is the main computational chemistry engine", "within GROMACS. Obviously, it performs Molecular Dynamics simulations,", "but it can also perform Stochastic Dynamics, Energy Minimization,", @@ -241,159 +249,43 @@ int gmx_mdrun(int argc, char *argv[]) "When [TT]mdrun[tt] is started with MPI, it does not run niced by default." }; - // Ongoing collection of mdrun options - MdrunOptions mdrunOptions; - // Options for the domain decomposition. - DomdecOptions domdecOptions; - // Parallelism-related user options. - gmx_hw_opt_t hw_opt; - // Command-line override for the duration of a neighbor list with the Verlet scheme. - int nstlist_cmdline = 0; - // Parameters for replica-exchange simulations. - ReplicaExchangeParameters replExParams; - - // Filenames and properties from command-line argument values. - auto filenames = gmx::MdFilenames(); - - // Print a warning if any force is larger than this (in kJ/mol nm). - real pforce = -1; - - // Output context for writing text files - // \todo Clarify initialization, ownership, and lifetime. - gmx_output_env_t *oenv = nullptr; - - /* Command line options */ - rvec realddxyz = {0, 0, 0}; - const char *ddrank_opt_choices[static_cast(DdRankOrder::nr)+1] = - { nullptr, "interleave", "pp_pme", "cartesian", nullptr }; - const char *dddlb_opt_choices[static_cast(DlbOption::nr)+1] = - { nullptr, "auto", "no", "yes", nullptr }; - const char *thread_aff_opt_choices[threadaffNR+1] = - { nullptr, "auto", "on", "off", nullptr }; - const char *nbpu_opt_choices[] = - { nullptr, "auto", "cpu", "gpu", nullptr }; - const char *pme_opt_choices[] = - { nullptr, "auto", "cpu", "gpu", nullptr }; - const char *pme_fft_opt_choices[] = - { nullptr, "auto", "cpu", "gpu", nullptr }; - gmx_bool bTryToAppendFiles = TRUE; - const char *gpuIdsAvailable = ""; - const char *userGpuTaskAssignment = ""; - - ImdOptions &imdOptions = mdrunOptions.imdOptions; - - t_pargs pa[] = { - - { "-dd", FALSE, etRVEC, {&realddxyz}, - "Domain decomposition grid, 0 is optimize" }, - { "-ddorder", FALSE, etENUM, {ddrank_opt_choices}, - "DD rank order" }, - { "-npme", FALSE, etINT, {&domdecOptions.numPmeRanks}, - "Number of separate ranks to be used for PME, -1 is guess" }, - { "-nt", FALSE, etINT, {&hw_opt.nthreads_tot}, - "Total number of threads to start (0 is guess)" }, - { "-ntmpi", FALSE, etINT, {&hw_opt.nthreads_tmpi}, - "Number of thread-MPI ranks to start (0 is guess)" }, - { "-ntomp", FALSE, etINT, {&hw_opt.nthreads_omp}, - "Number of OpenMP threads per MPI rank to start (0 is guess)" }, - { "-ntomp_pme", FALSE, etINT, {&hw_opt.nthreads_omp_pme}, - "Number of OpenMP threads per MPI rank to start (0 is -ntomp)" }, - { "-pin", FALSE, etENUM, {thread_aff_opt_choices}, - "Whether mdrun should try to set thread affinities" }, - { "-pinoffset", FALSE, etINT, {&hw_opt.core_pinning_offset}, - "The lowest logical core number to which mdrun should pin the first thread" }, - { "-pinstride", FALSE, etINT, {&hw_opt.core_pinning_stride}, - "Pinning distance in logical cores for threads, use 0 to minimize the number of threads per physical core" }, - { "-gpu_id", FALSE, etSTR, {&gpuIdsAvailable}, - "List of unique GPU device IDs available to use" }, - { "-gputasks", FALSE, etSTR, {&userGpuTaskAssignment}, - "List of GPU device IDs, mapping each PP task on each node to a device" }, - { "-ddcheck", FALSE, etBOOL, {&domdecOptions.checkBondedInteractions}, - "Check for all bonded interactions with DD" }, - { "-ddbondcomm", FALSE, etBOOL, {&domdecOptions.useBondedCommunication}, - "HIDDENUse special bonded atom communication when [TT]-rdd[tt] > cut-off" }, - { "-rdd", FALSE, etREAL, {&domdecOptions.minimumCommunicationRange}, - "The maximum distance for bonded interactions with DD (nm), 0 is determine from initial coordinates" }, - { "-rcon", FALSE, etREAL, {&domdecOptions.constraintCommunicationRange}, - "Maximum distance for P-LINCS (nm), 0 is estimate" }, - { "-dlb", FALSE, etENUM, {dddlb_opt_choices}, - "Dynamic load balancing (with DD)" }, - { "-dds", FALSE, etREAL, {&domdecOptions.dlbScaling}, - "Fraction in (0,1) by whose reciprocal the initial DD cell size will be increased in order to " - "provide a margin in which dynamic load balancing can act while preserving the minimum cell size." }, - { "-ddcsx", FALSE, etSTR, {&domdecOptions.cellSizeX}, - "HIDDENA string containing a vector of the relative sizes in the x " - "direction of the corresponding DD cells. Only effective with static " - "load balancing." }, - { "-ddcsy", FALSE, etSTR, {&domdecOptions.cellSizeY}, - "HIDDENA string containing a vector of the relative sizes in the y " - "direction of the corresponding DD cells. Only effective with static " - "load balancing." }, - { "-ddcsz", FALSE, etSTR, {&domdecOptions.cellSizeZ}, - "HIDDENA string containing a vector of the relative sizes in the z " - "direction of the corresponding DD cells. Only effective with static " - "load balancing." }, - { "-gcom", FALSE, etINT, {&mdrunOptions.globalCommunicationInterval}, - "Global communication frequency (deprecated)" }, - { "-nb", FALSE, etENUM, {nbpu_opt_choices}, - "Calculate non-bonded interactions on" }, - { "-nstlist", FALSE, etINT, {&nstlist_cmdline}, - "Set nstlist when using a Verlet buffer tolerance (0 is guess)" }, - { "-tunepme", FALSE, etBOOL, {&mdrunOptions.tunePme}, - "Optimize PME load between PP/PME ranks or GPU/CPU (only with the Verlet cut-off scheme)" }, - { "-pme", FALSE, etENUM, {pme_opt_choices}, - "Perform PME calculations on" }, - { "-pmefft", FALSE, etENUM, {pme_fft_opt_choices}, - "Perform PME FFT calculations on" }, - { "-v", FALSE, etBOOL, {&mdrunOptions.verbose}, - "Be loud and noisy" }, - { "-pforce", FALSE, etREAL, {&pforce}, - "Print all forces larger than this (kJ/mol nm)" }, - { "-reprod", FALSE, etBOOL, {&mdrunOptions.reproducible}, - "Try to avoid optimizations that affect binary reproducibility" }, - { "-cpt", FALSE, etREAL, {&mdrunOptions.checkpointOptions.period}, - "Checkpoint interval (minutes)" }, - { "-cpnum", FALSE, etBOOL, {&mdrunOptions.checkpointOptions.keepAndNumberCheckpointFiles}, - "Keep and number checkpoint files" }, - { "-append", FALSE, etBOOL, {&bTryToAppendFiles}, - "Append to previous output files when continuing from checkpoint instead of adding the simulation part number to all file names" }, - { "-nsteps", FALSE, etINT64, {&mdrunOptions.numStepsCommandline}, - "Run this number of steps, overrides .mdp file option (-1 means infinite, -2 means use mdp option, smaller is invalid) (deprecated)" }, - { "-maxh", FALSE, etREAL, {&mdrunOptions.maximumHoursToRun}, - "Terminate after 0.99 times this time (hours)" }, - { "-replex", FALSE, etINT, {&replExParams.exchangeInterval}, - "Attempt replica exchange periodically with this period (steps)" }, - { "-nex", FALSE, etINT, {&replExParams.numExchanges}, - "Number of random exchanges to carry out each exchange interval (N^3 is one suggestion). -nex zero or not specified gives neighbor replica exchange." }, - { "-reseed", FALSE, etINT, {&replExParams.randomSeed}, - "Seed for replica exchange, -1 is generate a seed" }, - { "-imdport", FALSE, etINT, {&imdOptions.port}, - "HIDDENIMD listening port" }, - { "-imdwait", FALSE, etBOOL, {&imdOptions.wait}, - "HIDDENPause the simulation while no IMD client is connected" }, - { "-imdterm", FALSE, etBOOL, {&imdOptions.terminatable}, - "HIDDENAllow termination of the simulation from IMD client" }, - { "-imdpull", FALSE, etBOOL, {&imdOptions.pull}, - "HIDDENAllow pulling in the simulation from IMD client" }, - { "-rerunvsite", FALSE, etBOOL, {&mdrunOptions.rerunConstructVsites}, - "HIDDENRecalculate virtual site coordinates with [TT]-rerun[tt]" }, - { "-confout", FALSE, etBOOL, {&mdrunOptions.writeConfout}, - "HIDDENWrite the last configuration with [TT]-c[tt] and force checkpointing at the last step (deprecated)" }, - { "-stepout", FALSE, etINT, {&mdrunOptions.verboseStepPrintInterval}, - "HIDDENFrequency of writing the remaining wall clock time for the run (deprecated)" }, - { "-resetstep", FALSE, etINT, {&mdrunOptions.timingOptions.resetStep}, - "HIDDENReset cycle counters after these many time steps (deprecated)" }, - { "-resethway", FALSE, etBOOL, {&mdrunOptions.timingOptions.resetHalfway}, - "HIDDENReset the cycle counters after half the number of steps or halfway [TT]-maxh[tt] (deprecated)" } - }; - int rc; + LegacyMdrunOptions options; // pointer-to-t_commrec is the de facto handle type for communications record. - // Complete shared / borrowed ownership requires a reference to this stack variable - // (or pointer-to-pointer-to-t_commrec) since borrowing code may update the pointer. // \todo Define the ownership and lifetime management semantics for a communication record, handle or value type. - t_commrec * cr = init_commrec(); + options.cr = init_commrec(); + + if (options.updateFromCommandLine(argc, argv, desc) == 0) + { + return 0; + } + + if (MASTER(options.cr)) + { + options.logFileGuard = openLogFile(ftp2fn(efLOG, + options.filenames.size(), + options.filenames.data()), + options.mdrunOptions.continuationOptions.appendFiles, + options.cr->nodeid, + options.cr->nnodes); + } + + /* The SimulationContext is a resource owned by the client code. + * A more complete design should address handles to resources with appropriate + * lifetimes and invariants for the resources allocated to the client, + * to the current simulation and to scheduled tasks within the simulation. + * + * \todo Clarify Context lifetime-management requirements and reconcile with scoped ownership. + * + * \todo Take ownership of and responsibility for communications record (cr). + */ + auto context = createSimulationContext(options.cr); + return launchMdrun(std::move(options), compat::not_null(&context)); +} + +int LegacyMdrunOptions::updateFromCommandLine(int argc, char **argv, ArrayRef desc) +{ unsigned long PCA_Flags = PCA_CAN_SET_DEFFNM; // With -multidir, the working directory still needs to be // changed, so we can't check for the existence of files during @@ -405,10 +297,9 @@ int gmx_mdrun(int argc, char *argv[]) } if (!parse_common_args(&argc, argv, PCA_Flags, - static_cast(filenames().size()), filenames().data(), asize(pa), pa, - asize(desc), desc, 0, nullptr, &oenv)) + static_cast(filenames.size()), filenames.data(), asize(pa), pa, + static_cast(desc.size()), desc.data(), 0, nullptr, &oenv)) { - sfree(cr); return 0; } @@ -454,9 +345,9 @@ int gmx_mdrun(int argc, char *argv[]) hw_opt.thread_affinity = nenum(thread_aff_opt_choices); // now check for a multi-simulation - gmx::ArrayRef multidir = opt2fnsIfOptionSet("-multidir", - static_cast(filenames().size()), - filenames().data()); + ArrayRef multidir = opt2fnsIfOptionSet("-multidir", + static_cast(filenames.size()), + filenames.data()); if (replExParams.exchangeInterval != 0 && multidir.size() < 2) { @@ -468,7 +359,7 @@ int gmx_mdrun(int argc, char *argv[]) gmx_fatal(FARGS, "Replica exchange number of exchanges needs to be positive"); } - gmx_multisim_t* ms = init_multisystem(MPI_COMM_WORLD, multidir); + ms = init_multisystem(MPI_COMM_WORLD, multidir); /* Prepare the intra-simulation communication */ // TODO consolidate this with init_commrec, after changing the @@ -485,7 +376,7 @@ int gmx_mdrun(int argc, char *argv[]) #endif if (!opt2bSet("-cpi", - static_cast(filenames().size()), filenames().data())) + static_cast(filenames.size()), filenames.data())) { // If we are not starting from a checkpoint we never allow files to be appended // to, since that has caused a ton of strange behaviour and bugs in the past. @@ -506,75 +397,28 @@ int gmx_mdrun(int argc, char *argv[]) continuationOptions.appendFilesOptionSet = opt2parg_bSet("-append", asize(pa), pa); handleRestart(cr, ms, bTryToAppendFiles, - static_cast(filenames().size()), - filenames().data(), + static_cast(filenames.size()), + filenames.data(), &continuationOptions.appendFiles, &continuationOptions.startedFromCheckpoint); mdrunOptions.rerun = opt2bSet("-rerun", - static_cast(filenames().size()), - filenames().data()); + static_cast(filenames.size()), + filenames.data()); mdrunOptions.ntompOptionIsSet = opt2parg_bSet("-ntomp", asize(pa), pa); - LogFilePtr logFileGuard(nullptr); - if (MASTER(cr)) - { - logFileGuard = openLogFile(ftp2fn(efLOG, - filenames().size(), - filenames().data()), - continuationOptions.appendFiles, - cr->nodeid, - cr->nnodes); - } - domdecOptions.rankOrder = static_cast(nenum(ddrank_opt_choices)); domdecOptions.dlbOption = static_cast(nenum(dddlb_opt_choices)); domdecOptions.numCells[XX] = roundToInt(realddxyz[XX]); domdecOptions.numCells[YY] = roundToInt(realddxyz[YY]); domdecOptions.numCells[ZZ] = roundToInt(realddxyz[ZZ]); - if (logFileGuard != nullptr && !mdrunOptions.continuationOptions.appendFiles) - { - FILE *fplog = gmx_fio_getfp(logFileGuard.get()); - if (opt2parg_bSet("-gcom", asize(pa), pa)) - { - fprintf(fplog, "Note that gmx mdrun -gcom is deprecated, and will be removed\n" - "in a future version of GROMACS.\n\n"); - } - if (opt2bSet("-table", filenames().size(), filenames().data()) || - opt2bSet("-tablep", filenames().size(), filenames().data()) || - opt2bSet("-tableb", filenames().size(), filenames().data())) - { - fprintf(fplog, "Note that gmx mdrun -table options are deprecated. Table reading\n" - "will likely be supported by grompp in a future version of GROMACS.\n\n"); - } - if (opt2parg_bSet("-confout", asize(pa), pa) || - opt2parg_bSet("-resetstep", asize(pa), pa) || - opt2parg_bSet("-resethway", asize(pa), pa)) - { - fprintf(fplog, "Note that gmx mdrun options intended for benchmarking, including\n" - "-confout, -resetstep, and -resethway are deprecated. They will\n" - "likely be supported by gmx benchmark in a future version of GROMACS.\n\n"); - } - if (opt2parg_bSet("-nsteps", asize(pa), pa)) - { - fprintf(fplog, "Note that gmx mdrun -nsteps is deprecated. Use gmx convert-tpr\n" - "or gmx grompp to change the intended number of steps for your\n" - "simulation.\n\n"); - } - } - - /* The Context is a shared resource owned by the client code. - * A more complete design should address handles to resources with appropriate - * lifetimes and invariants for the resources allocated to the client, - * to the current simulation and to scheduled tasks within the simulation. - * - * \todo Clarify Context lifetime-management requirements and reconcile with scoped ownership. - * - * \todo Take ownership of and responsibility for communications record (cr). - */ - auto context = gmx::createSimulationContext(cr); + return 1; +} +int launchMdrun(LegacyMdrunOptions &&options, + compat::not_null context) +{ /* The named components for the builder exposed here are descriptive of the * state of mdrun at implementation and are not intended to be prescriptive * of future design. (Note the ICommandLineOptions... framework used elsewhere.) @@ -586,36 +430,40 @@ int gmx_mdrun(int argc, char *argv[]) * We would prefer to rebuild resources only as necessary, but we defer such * details to future optimizations. */ - auto builder = gmx::MdrunnerBuilder(compat::not_null(&context)); - builder.addSimulationMethod(mdrunOptions, pforce); - builder.addDomainDecomposition(domdecOptions); + auto builder = MdrunnerBuilder(context); + builder.addSimulationMethod(options.mdrunOptions, options.pforce); + builder.addDomainDecomposition(options.domdecOptions); // \todo pass by value - builder.addNonBonded(nbpu_opt_choices[0]); + builder.addNonBonded(options.nbpu_opt_choices[0]); // \todo pass by value - builder.addElectrostatics(pme_opt_choices[0], pme_fft_opt_choices[0]); - builder.addNeighborList(nstlist_cmdline); - builder.addReplicaExchange(replExParams); + builder.addElectrostatics(options.pme_opt_choices[0], options.pme_fft_opt_choices[0]); + builder.addNeighborList(options.nstlist_cmdline); + builder.addReplicaExchange(options.replExParams); // \todo take ownership of multisim resources (ms) - builder.addMultiSim(ms); + builder.addMultiSim(options.ms); // \todo Provide parallelism resources through SimulationContext. // Need to establish run-time values from various inputs to provide a resource handle to Mdrunner - builder.addHardwareOptions(hw_opt); + builder.addHardwareOptions(options.hw_opt); // \todo File names are parameters that should be managed modularly through further factoring. - builder.addFilenames(filenames); + builder.addFilenames(options.filenames); // Note: The gmx_output_env_t life time is not managed after the call to parse_common_args. // \todo Implement lifetime management for gmx_output_env_t. // \todo Output environment should be configured outside of Mdrunner and provided as a resource. - builder.addOutputEnvironment(oenv); - builder.addLogFile(logFileGuard.get()); + builder.addOutputEnvironment(options.oenv); + builder.addLogFile(options.logFileGuard.get()); auto runner = builder.build(); - rc = runner.mdrunner(); + return runner.mdrunner(); +} +LegacyMdrunOptions::~LegacyMdrunOptions() +{ if (GMX_LIB_MPI) { done_commrec(cr); } done_multisim(ms); - return rc; +} + } diff --git a/src/programs/mdrun/mdrun_main.h b/src/programs/mdrun/mdrun_main.h index 8b0799f138..38f899b66f 100644 --- a/src/programs/mdrun/mdrun_main.h +++ b/src/programs/mdrun/mdrun_main.h @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2013, by the GROMACS development team, led by + * Copyright (c) 2013,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. @@ -35,6 +35,11 @@ #ifndef GMX_PROGRAMS_MDRUN_MDRUN_H #define GMX_PROGRAMS_MDRUN_MDRUN_H +namespace gmx +{ + int gmx_mdrun(int argc, char *argv[]); +} // namespace gmx + #endif -- 2.11.4.GIT