From b35147fdfe4c6a6a62f635691d871886d7783cea Mon Sep 17 00:00:00 2001 From: "M. Eric Irrgang" Date: Sat, 26 Sep 2020 12:54:32 +0300 Subject: [PATCH] Interface updates allowing client to provide MPI communicator. gmxapi::MpiContextManager is acquired earlier in API initialization and has a larger role in setting up the SimulationContext as a member of gmxapi::ContextImpl. The MpiContextManager can be acquired with appropriate default values for the various MPI-related GROMACS build variants, or it can be initialized with user-provided resources. gmxapi::Context creation method has been updated to accept optional non-default resources from the client. A template header helps MPI-enabled client code to provide resources to createContext(). The client requires minimal awareness of whether GROMACS was built for thread-MPI or an installed MPI library. Thread-MPI details are isolated from any translation unit that might include the client MPI tool chain. Refs #3644 --- api/gmxapi/CMakeLists.txt | 3 + api/gmxapi/cpp/CMakeLists.txt | 26 +-- api/gmxapi/include/gmxapi/context.h | 57 ++++++ api/gmxapi/include/gmxapi/mpi/gmxapi_mpi.h | 170 ++++++++++++++++++ .../multiprocessingresources.cmakein.h | 135 ++++++++++++++ .../2021/major/miscellaneous.rst | 9 + src/gromacs/mdrun/simulationcontext.h | 21 ++- 7 files changed, 402 insertions(+), 19 deletions(-) create mode 100644 api/gmxapi/include/gmxapi/mpi/gmxapi_mpi.h create mode 100644 api/gmxapi/include/multiprocessingresources.cmakein.h diff --git a/api/gmxapi/CMakeLists.txt b/api/gmxapi/CMakeLists.txt index 4c73464704..be182daa33 100644 --- a/api/gmxapi/CMakeLists.txt +++ b/api/gmxapi/CMakeLists.txt @@ -67,6 +67,7 @@ set_target_properties(gmxapi PROPERTIES # The include directory should be mostly empty so that we can use it internally as # the public interface include directory during build and testing. configure_file(include/gmxapiversion.h.in include/gmxapi/version.h) +configure_file(include/multiprocessingresources.cmakein.h include/gmxapi/mpi/multiprocessingresources.h) target_include_directories(gmxapi PUBLIC $ $ @@ -89,4 +90,6 @@ install(DIRECTORY include/gmxapi # Install "configured" files from the build tree. install(FILES ${CMAKE_CURRENT_BINARY_DIR}/include/gmxapi/version.h DESTINATION include/gmxapi) +install(FILES ${CMAKE_CURRENT_BINARY_DIR}/include/gmxapi/mpi/multiprocessingresources.h + DESTINATION include/gmxapi/mpi) diff --git a/api/gmxapi/cpp/CMakeLists.txt b/api/gmxapi/cpp/CMakeLists.txt index 813af99271..129ffb69fc 100644 --- a/api/gmxapi/cpp/CMakeLists.txt +++ b/api/gmxapi/cpp/CMakeLists.txt @@ -35,19 +35,19 @@ # This list file provides the Gromacs::gmxapi cmake target. target_sources(gmxapi PRIVATE - context.cpp - exceptions.cpp - gmxapi.cpp - md.cpp - mdmodule.cpp - mdsignals.cpp - session.cpp - status.cpp - system.cpp - version.cpp - workflow.cpp - tpr.cpp - ) + context.cpp + exceptions.cpp + gmxapi.cpp + md.cpp + mdmodule.cpp + mdsignals.cpp + session.cpp + status.cpp + system.cpp + version.cpp + workflow.cpp + tpr.cpp + ) set_target_properties(gmxapi PROPERTIES POSITION_INDEPENDENT_CODE ON) gmx_target_compile_options(gmxapi) target_compile_definitions(gmxapi PRIVATE HAVE_CONFIG_H) diff --git a/api/gmxapi/include/gmxapi/context.h b/api/gmxapi/include/gmxapi/context.h index c705d08a43..ec434059af 100644 --- a/api/gmxapi/include/gmxapi/context.h +++ b/api/gmxapi/include/gmxapi/context.h @@ -39,6 +39,24 @@ * * \author M. Eric Irrgang * \ingroup gmxapi + * + * # Notes on scope + * + * Client owns Context and MultiProcessingResources, which exist entirely in client scope. + * + * Note that MpiContextManager and Session live in library scope. + * + * # MPI-enabled clients + * + * When an MPI toolchain is compatible with the GROMACS build + * (see [#3672](https://gitlab.com/gromacs/gromacs/-/issues/3672)), + * multi-process clients may determine how processes are allocated to GROMACS + * simulations by providing Contexts with specific MPI Communicators to + * createContext(const MultiProcessingResources& resources) + * + * MultiProcessingResources is an opaque type of library-internal resources. + * Clients acquire such resources with the assignResource() template helper. + * See gmxapi_mpi.h */ #include @@ -168,6 +186,45 @@ private: std::shared_ptr impl_; }; +// Forward declaration for interoperation with the library. +// Client code implements a concrete ResourceAssignment indirectly by instantiating +// assignResource() through the gmxapi_mpi.h template header. +class ResourceAssignment; + +// See template header gmxapi_mpi.h +template +std::unique_ptr assignResource(CommT communicator); + +/*! + * \brief Initialize a new API Context to manage resources and software environment. + * + * The client is responsible for keeping the Context instance alive for at least + * as long as any API objects initialized from it. We allow this responsibility + * to be placed on the client (rather than using a global Singleton) because the + * library is theoretically reentrant, and multiple Context objects may exist. + * + * Client may explicitly assign run time computing resources for use within the Context. + * Otherwise, the library will try to automatically detect and allocate computing resources. + * + * Currently, only MPI communication scope can be assigned explicitly. + * Future development is needed before this interface can influence the subdivision + * of available GPUs or thread-based processor allocation. See #3650 and #3688. + * + * \return Initialized Context instance. + * + * \internal + * Use cases: + * + * 1. tMPI and client-provided comm: provide a place for safety checks, then construct a suitable + * dummy MpiContextManager. + * 2. tMPI and no client-provided comm: construct suitable dummy MpiContextManager + * 3. MPI and client-provided comm: use compatible comm. error if COMM_NULL + * 4. MPI and no client-provided comm: generate MpiContextManager with COMM_WORLD + * + */ +Context createContext(); +Context createContext(const ResourceAssignment& resources); + } // end namespace gmxapi #endif // header guard diff --git a/api/gmxapi/include/gmxapi/mpi/gmxapi_mpi.h b/api/gmxapi/include/gmxapi/mpi/gmxapi_mpi.h new file mode 100644 index 0000000000..b8f8db24d3 --- /dev/null +++ b/api/gmxapi/include/gmxapi/mpi/gmxapi_mpi.h @@ -0,0 +1,170 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2020, 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. + */ + +#ifndef GMXAPI_MPI_H +#define GMXAPI_MPI_H + +#include + +#include + +#include "gmxapi/context.h" +#include "gmxapi/exceptions.h" + +#include "gmxapi/mpi/multiprocessingresources.h" + +/*! \file + * \brief Provide details of any MPI implementation used when building the library. + * + * If the gmxapi library was built with MPI-enabled GROMACS, client software + * using a compatible MPI implementation may use this header to access additional + * functionality. + * + * Client software should use the CMake infrastructure to ensure a compatible + * MPI implementation. + * + * Use of an incompatible MPI implementation will produce linking errors. + * + * \todo With resolution of #3672, we can include additional checks. + * + */ + +namespace gmxapi +{ + +/*! + * \brief Convey resources assigned by the client to the new library Context. + * + * Only enabled if the required MPI library calls are available. + */ +template +class ResourceAssignmentImpl : public ResourceAssignment +{ +public: + ~ResourceAssignmentImpl() override = default; + + /*! + * \brief Initialize from borrowed communicator. + * + * \param communicator Client-provided communicator relayed through assignResource() + * + * Create an abstract wrapper for client-provided values with which to initialize + * simulation resources. When this wrapper is used, the client is responsible for providing + * a valid communicator that will remain valid for the life of the consumer. + */ + explicit ResourceAssignmentImpl(const CommT& communicator) : communicator_{ communicator } + { + if (communicator_ == MPI_COMM_NULL) + { + throw UsageError("Null communicator cannot be lent."); + } + int flag = 0; + MPI_Initialized(&flag); + if (!flag) + { + throw UsageError("Client has not initialized MPI context."); + } + } + + [[nodiscard]] int size() const override + { + assert(communicator_ != MPI_COMM_NULL && "Resource invariant implies a valid communicator."); + int size = 0; + MPI_Comm_size(communicator_, &size); + return size; + } + + [[nodiscard]] int rank() const override + { + assert(communicator_ != MPI_COMM_NULL && "Resource invariant implies a valid communicator."); + // This default value will never be read, but the compiler can't tell + // that it is initialized by the MPI call. + int rank = -1; + MPI_Comm_rank(communicator_, &rank); + return rank; + } + + void applyCommunicator(CommHandle* dst) const override + { + assert(communicator_ != MPI_COMM_NULL && "Resource invariant implies a valid communicator."); + offerComm(communicator_, dst); + } + + CommT communicator_; +}; + +/*! + * \brief Template header utility for connecting to MPI implementations. + * + * To use this helper function, the client software build environment must be + * configured for an MPI implementation compatible with the target GROMACS library. + * + * If the library and client software are built with compatible MPI implementations, + * the template will capture the communicator in an opaque container that can be + * passed to the ContextImpl creation function. Otherwise, the opaque container + * will carry a suitable NULL object. The template definition is necessarily + * provided by a generated template header, since the definition depends on the + * GROMACS installation. + * + * \tparam CommT Deduced type of client-provided communicator. + * \param communicator MPI (sub)communicator linking collaborating processes in the new Context. + * \return Opaque resource handle usable when setting up execution context. + * + * \see createContext(std::unique_ptr resources) + * + * The client provides a communicator for work executed within the scope of the Context. + * Client remains responsible for freeing the communicator and finalizing the MPI environment. + * + * The communicator resource type is a template parameter because MPI_Comm is commonly a C + * typedef that varies between implementations, so we do not want to couple our + * API to it, but we cannot forward-declare it. + * + * See also https://gitlab.com/gromacs/gromacs/-/issues/3650 + */ +template +std::unique_ptr assignResource(CommT communicator) +{ + if (communicator == MPI_COMM_NULL) + { + throw UsageError("Cannot assign a Null communicator."); + } + return std::make_unique>(communicator); +} + +// Also ref https://stackoverflow.com/questions/49259704/pybind11-possible-to-use-mpi4py + +} // end namespace gmxapi + +#endif // GMXAPI_MPI_H diff --git a/api/gmxapi/include/multiprocessingresources.cmakein.h b/api/gmxapi/include/multiprocessingresources.cmakein.h new file mode 100644 index 0000000000..8aa69d0906 --- /dev/null +++ b/api/gmxapi/include/multiprocessingresources.cmakein.h @@ -0,0 +1,135 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2020, 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. + */ + +#ifndef GMXAPI_MULTIPROCESSINGRESOURCES_H +#define GMXAPI_MULTIPROCESSINGRESOURCES_H + +/*! \file + * \brief Provide build-specific overloads for client-MPI-dependent stuff. + * + * Define the interfaces that a client must implement to generate MpiContextManager for the library. + * Client code should use the gmxapi_mpi.h template header to generate code supporting + * the required interface. (Client code should not need to include this header directly.) + * + * \note This is a generated header. Some definitions are determined when the GROMACS library + * build is configured. + * + * If the library is built with tMPI, CommHandle is empty and offerComm is not defined. + * + * If the library is built with an MPI library, CommHandle holds a native MPI_Comm object and + * the library-provided offerComm is used by the client MultiProcessingResources implementation + * to support passing a copy of the communicator. + * + * Note: The communicator is not `dup`ed in this call. If the library later duplicates + * the offered communicator, the library will be responsible for freeing the duplicate. + * However, the caller is responsible for keeping any MPI environment valid while the library is + * in use. For details, \see MpiContextManager. + * + * \author "M. Eric Irrgang " + */ + +#include +#include + +// The interface in this header is determined when the GROMACS library build is configured. +#cmakedefine01 GMX_LIB_MPI +#if GMX_LIB_MPI +# include +#endif + +namespace gmxapi +{ + +// Forward declaration for library resources. +// CommHandle is opaque to the public API, but allows extensions to implement +// required library interfaces by composing the behavior of helpers like offerComm(). +// The abstraction is important because the library implementations depend on +// the options with which GROMACS was built. +class CommHandle; + +/*! + * \brief Convey the resources that the client has directed the library to use within a Context. + * + * The object itself does not convey ownership of resources. However, the interfaces + * for declaring assigned resources must have well-specified (and documented) + * semantics for resource ownership. See assignResource(). + * (Note: the initial implementation only allows for assignment of an MPI Communicator.) + * + * The library and client share this interface definition. The implementation is + * provided by client code with the aid of a template header. + * Client developers should refer instead to gmxapi_mpi.h. + */ +class ResourceAssignment +{ +public: + virtual ~ResourceAssignment(); + [[nodiscard]] virtual int size() const = 0; + [[nodiscard]] virtual int rank() const = 0; + virtual void applyCommunicator(class CommHandle* dst) const; +}; + +/*! \brief Offer the client communicator to the library. + * + * Helper function allowing clients to provide the MPI communicator for the library. + * \param src client-provided communicator. + * \param dst library recipient, abstracted to hide library MPI type details. + * + */ +template +void offerComm([[maybe_unused]] CommT src, [[maybe_unused]] CommHandle* dst) +{ + // Default overload does nothing. +} + +#if GMX_LIB_MPI +/*! \brief Offer the client communicator to the library. + * + * Helper function allowing clients to provide communicator to library. + * \param src communicator offered by client + * \param dst opaque pointer to library resource destination + * + * This function is only available in MPI-enabled GROMACS. + * Used indirectly in template helpers for client code to implement library interfaces. + * Clients should use the assignResource() higher level function in explicit code, + * which will evaluable appropriately for the target GROMACS library. + * + * \todo Provide a stub for docs even when not available. + */ +void offerComm(MPI_Comm src, CommHandle* dst); +#endif + +} // end namespace gmxapi + +#endif // GMXAPI_MULTIPROCESSINGRESOURCES_H diff --git a/docs/release-notes/2021/major/miscellaneous.rst b/docs/release-notes/2021/major/miscellaneous.rst index 68c87e9f81..77cd830423 100644 --- a/docs/release-notes/2021/major/miscellaneous.rst +++ b/docs/release-notes/2021/major/miscellaneous.rst @@ -31,6 +31,15 @@ change outside of the users direct control we have removed the support for automatically setting booleans. GMX_BUILD_HELP and GMX_HWLOC are now disabled by default, while GMX_LOAD_PLUGINS is enabled by default. +gmxapi C++ interface +"""""""""""""""""""" + +``gmxapi::Context`` is now created with ``gmxapi::createContext()``, which allows +the client to provide an MPI communicator for the library to use instead of its default +(e.g MPI_COMM_WORLD). MPI-enabled clients may use the :file:`gmxapi/mpi/gmxapi_mpi.h` +template header and the ``assignResource()`` helper to generate the argument to +``createContext``. + Unification of several CUDA and OpenCL environment variables """""""""""""""""""""""""""""""""""""""""""""""""""""""""""" diff --git a/src/gromacs/mdrun/simulationcontext.h b/src/gromacs/mdrun/simulationcontext.h index 0fb81b5865..c57027311d 100644 --- a/src/gromacs/mdrun/simulationcontext.h +++ b/src/gromacs/mdrun/simulationcontext.h @@ -79,6 +79,15 @@ class ArrayRef; * The public interface of SimulationContext is not yet well-specified. * Client code can create an instance with gmx::createSimulationContext() * + * \warning The SimulationContext does not yet encapsulate the resources allocated + * to the simulator. See https://gitlab.com/gromacs/gromacs/-/issues/3650 + * + * \warning On thread-MPI code paths, the SimulationContext *does not* represent + * constraints on the number of simulation ranks, nor does it represent initialized + * communication resources. + * + * \todo Clarify and strengthen the invariant represented by SimulationContext. + * * \todo This class should also handle aspects of simulation * environment such as working directory and environment variables. * @@ -106,7 +115,7 @@ public: */ SimulationContext() = delete; /*! - * \brief Construct + * \brief Initialize from borrowed communicator. * * \param communicator Communicator for all collaborating simulation processes. * \param multiSimDirectoryNames Names of any directories used with -multidir @@ -114,8 +123,8 @@ public: * Caller is responsible for keeping *communicator* valid for the life of SimulationContext, * and then freeing the communicator, if appropriate. * - * With an external MPI library (GMX_LIB_MPI==1), the client promises that - * gmx::init() has called MPI_Init and the provided communicator is valid to use. + * With an external MPI library (non-thread-MPI chosen when configuring with CMake), + * the client promises that MPI has been initialized (such as by calling gmx::init()). * This communicator is "borrowed" (not duplicated) from the caller. * Additional communicators may be split from the provided communicator * during the life of the SimulationContext or its consumers. @@ -134,20 +143,20 @@ public: /*! * \brief MPI resources for the entire simulation context. * - * With an external MPI library (GMX_LIB_MPI==1), + * With an external MPI library (non-thread-MPI chosen when configuring with CMake), * gmx::init() has called MPI_Init and the provided communicator is valid to use. * The communicator is "borrowed" (not duplicated) from the caller. * * With thread-MPI, the communicator is set up later * during the process of spawning the threads that will be the MPI - * ranks. (Multi-simulation is not supported with thread-MPI.) + * ranks. */ MPI_Comm libraryWorldCommunicator_ = MPI_COMM_NULL; /*! * \brief MPI communicator object for this simulation. * - * With an external MPI library (GMX_LIB_MPI==1), + * With an external MPI library (non-thread-MPI chosen when configuring with CMake), * gmx::init() has called MPI_Init and the provided communicator is valid to use. * The communicator is "borrowed" (not duplicated) from the caller. * -- 2.22.0