From c0d698dfb0dc87af2f03e19860c832ad262d6b8d Mon Sep 17 00:00:00 2001 From: Joe Jordan Date: Tue, 29 Sep 2020 13:48:33 +0000 Subject: [PATCH] Add nblib backend: Part 2 of 2 --- api/nblib/CMakeLists.txt | 5 + api/nblib/forcecalculator.cpp | 69 ++++ api/nblib/gmxcalculator.cpp | 105 ++++++ api/nblib/gmxcalculator.h | 116 +++++++ api/nblib/gmxsetup.cpp | 327 ++++++++++++++++++ api/nblib/gmxsetup.h | 127 +++++++ api/nblib/integrator.cpp | 76 ++++ api/nblib/simulationstate.cpp | 146 ++++++++ api/nblib/simulationstateimpl.h | 87 +++++ api/nblib/tests/CMakeLists.txt | 32 ++ api/nblib/tests/gmxcalculator.cpp | 81 +++++ api/nblib/tests/integrator.cpp | 137 ++++++++ api/nblib/tests/nbkernelsystem.cpp | 207 +++++++++++ api/nblib/tests/nbnxnsetup.cpp | 65 ++++ .../NBlibTest_ArgonForcesAreCorrect.xml | 67 ++++ .../NBlibTest_SpcMethanolForcesAreCorrect.xml | 37 ++ api/nblib/tests/simstate.cpp | 154 +++++++++ api/nblib/tests/testsystems.cpp | 107 ++++++ api/nblib/tests/testsystems.h | 60 ++++ 19 files changed, 2005 insertions(+) create mode 100644 api/nblib/forcecalculator.cpp create mode 100644 api/nblib/gmxcalculator.cpp create mode 100644 api/nblib/gmxcalculator.h create mode 100644 api/nblib/gmxsetup.cpp create mode 100644 api/nblib/gmxsetup.h create mode 100644 api/nblib/integrator.cpp create mode 100644 api/nblib/simulationstate.cpp create mode 100644 api/nblib/simulationstateimpl.h create mode 100644 api/nblib/tests/gmxcalculator.cpp create mode 100644 api/nblib/tests/integrator.cpp create mode 100644 api/nblib/tests/nbkernelsystem.cpp create mode 100644 api/nblib/tests/nbnxnsetup.cpp create mode 100644 api/nblib/tests/refdata/NBlibTest_ArgonForcesAreCorrect.xml create mode 100644 api/nblib/tests/refdata/NBlibTest_SpcMethanolForcesAreCorrect.xml create mode 100644 api/nblib/tests/simstate.cpp diff --git a/api/nblib/CMakeLists.txt b/api/nblib/CMakeLists.txt index 3d307976af..6adf27a36b 100644 --- a/api/nblib/CMakeLists.txt +++ b/api/nblib/CMakeLists.txt @@ -88,9 +88,14 @@ add_library(nblib SHARED "") target_sources(nblib PRIVATE box.cpp + forcecalculator.cpp + gmxcalculator.cpp + gmxsetup.cpp + integrator.cpp interactions.cpp molecules.cpp particletype.cpp + simulationstate.cpp topologyhelpers.cpp topology.cpp ) diff --git a/api/nblib/forcecalculator.cpp b/api/nblib/forcecalculator.cpp new file mode 100644 index 0000000000..61dc95af31 --- /dev/null +++ b/api/nblib/forcecalculator.cpp @@ -0,0 +1,69 @@ +/* + * 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. + */ +/*! \internal \file + * \brief Implements nblib ForceCalculator + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ +#include "nblib/forcecalculator.h" +#include "nblib/gmxcalculator.h" +#include "nblib/gmxsetup.h" + +namespace nblib +{ + +ForceCalculator::~ForceCalculator() = default; + +ForceCalculator::ForceCalculator(const SimulationState& system, const NBKernelOptions& options) +{ + gmxForceCalculator_ = nblib::GmxSetupDirector::setupGmxForceCalculator(system, options); +} + +void ForceCalculator::compute(gmx::ArrayRef coordinates, gmx::ArrayRef forces) +{ + return gmxForceCalculator_->compute(coordinates, forces); +} + +void ForceCalculator::updatePairList(gmx::ArrayRef particleInfoAllVdW, + gmx::ArrayRef coordinates, + const Box& box) +{ + gmxForceCalculator_->setParticlesOnGrid(particleInfoAllVdW, coordinates, box); +} + +} // namespace nblib diff --git a/api/nblib/gmxcalculator.cpp b/api/nblib/gmxcalculator.cpp new file mode 100644 index 0000000000..30e2594e40 --- /dev/null +++ b/api/nblib/gmxcalculator.cpp @@ -0,0 +1,105 @@ +/* + * 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. + */ +/*! \internal \file + * \brief Implements a force calculator based on GROMACS data structures. + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ +#include "nblib/gmxcalculator.h" +#include "gromacs/ewald/ewald_utils.h" +#include "gromacs/gmxlib/nrnb.h" +#include "gromacs/mdlib/rf_util.h" +#include "gromacs/mdtypes/enerdata.h" +#include "gromacs/mdtypes/forcerec.h" +#include "gromacs/mdtypes/interaction_const.h" +#include "gromacs/mdtypes/simulation_workload.h" +#include "gromacs/nbnxm/nbnxm.h" +#include "gromacs/utility/range.h" +#include "nblib/exception.h" +#include "nblib/simulationstate.h" + +namespace nblib +{ + +GmxForceCalculator::GmxForceCalculator() +{ + enerd_ = std::make_unique(1, 0); + forcerec_ = std::make_unique(); + interactionConst_ = std::make_unique(); + stepWork_ = std::make_unique(); + nrnb_ = std::make_unique(); +} + +GmxForceCalculator::~GmxForceCalculator() = default; + +void GmxForceCalculator::compute(gmx::ArrayRef coordinateInput, + gmx::ArrayRef forceOutput) +{ + // update the coordinates in the backend + nbv_->convertCoordinates(gmx::AtomLocality::Local, false, coordinateInput); + + // set forces to zero + std::fill(forceOutput.begin(), forceOutput.end(), gmx::RVec{ 0, 0, 0 }); + + nbv_->dispatchNonbondedKernel(gmx::InteractionLocality::Local, *interactionConst_, *stepWork_, + enbvClearFYes, *forcerec_, enerd_.get(), nrnb_.get()); + + nbv_->atomdata_add_nbat_f_to_f(gmx::AtomLocality::All, forceOutput); +} + +void GmxForceCalculator::setParticlesOnGrid(gmx::ArrayRef particleInfoAllVdw, + gmx::ArrayRef coordinates, + const Box& box) +{ + auto legacyBox = box.legacyMatrix(); + + if (TRICLINIC(legacyBox)) + { + throw InputException("Only rectangular unit-cells are supported here"); + } + const rvec lowerCorner = { 0, 0, 0 }; + const rvec upperCorner = { legacyBox[dimX][dimX], legacyBox[dimY][dimY], legacyBox[dimZ][dimZ] }; + + const real particleDensity = coordinates.size() / det(legacyBox); + + nbnxn_put_on_grid(nbv_.get(), legacyBox, 0, lowerCorner, upperCorner, nullptr, + { 0, int(coordinates.size()) }, particleDensity, particleInfoAllVdw, + coordinates, 0, nullptr); +} + +} // namespace nblib diff --git a/api/nblib/gmxcalculator.h b/api/nblib/gmxcalculator.h new file mode 100644 index 0000000000..610983ccf7 --- /dev/null +++ b/api/nblib/gmxcalculator.h @@ -0,0 +1,116 @@ +/* + * 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. + */ +/*! \libinternal \file + * \brief + * Implements a force calculator based on GROMACS data structures. + * + * Intended for internal use inside the ForceCalculator. + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ + +#ifndef NBLIB_GMXCALCULATOR_H +#define NBLIB_GMXCALCULATOR_H + +#include + +#include "nblib/vector.h" + +struct nonbonded_verlet_t; +struct t_forcerec; +struct t_nrnb; +struct interaction_const_t; +struct gmx_enerdata_t; + +namespace gmx +{ +template +class ArrayRef; +class StepWorkload; +} // namespace gmx + +namespace nblib +{ +class Box; +class NbvSetupUtil; +class SimulationState; +struct NBKernelOptions; + +class GmxForceCalculator final +{ +public: + GmxForceCalculator(); + + ~GmxForceCalculator(); + + //! Compute forces and return + void compute(gmx::ArrayRef coordinateInput, gmx::ArrayRef forceOutput); + + //! Puts particles on a grid based on bounds specified by the box (for every NS step) + void setParticlesOnGrid(gmx::ArrayRef particleInfoAllVdw, + gmx::ArrayRef coordinates, + const Box& box); + +private: + friend class NbvSetupUtil; + + //! Non-Bonded Verlet object for force calculation + std::unique_ptr nbv_; + + //! Only nbfp and shift_vec are used + std::unique_ptr forcerec_; + + //! Parameters for various interactions in the system + std::unique_ptr interactionConst_; + + //! Tasks to perform in an MD Step + std::unique_ptr stepWork_; + + //! Energies of different interaction types; currently only needed as an argument for dispatchNonbondedKernel + std::unique_ptr enerd_; + + //! Non-bonded flop counter; currently only needed as an argument for dispatchNonbondedKernel + std::unique_ptr nrnb_; + + //! Legacy matrix for box + matrix box_{ { 0 } }; +}; + +} // namespace nblib + +#endif // NBLIB_GMXCALCULATOR_H diff --git a/api/nblib/gmxsetup.cpp b/api/nblib/gmxsetup.cpp new file mode 100644 index 0000000000..66afb08add --- /dev/null +++ b/api/nblib/gmxsetup.cpp @@ -0,0 +1,327 @@ +/* + * 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. + */ +/*! \internal \file + * \brief Translation layer to GROMACS data structures for force calculations. + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ +#include "nblib/gmxsetup.h" +#include "gromacs/ewald/ewald_utils.h" +#include "gromacs/gmxlib/nrnb.h" +#include "gromacs/mdlib/forcerec.h" +#include "gromacs/mdlib/gmx_omp_nthreads.h" +#include "gromacs/mdlib/rf_util.h" +#include "gromacs/mdtypes/forcerec.h" +#include "gromacs/mdtypes/interaction_const.h" +#include "gromacs/mdtypes/simulation_workload.h" +#include "gromacs/nbnxm/atomdata.h" +#include "gromacs/nbnxm/nbnxm.h" +#include "gromacs/nbnxm/nbnxm_simd.h" +#include "gromacs/nbnxm/pairlistset.h" +#include "gromacs/nbnxm/pairlistsets.h" +#include "gromacs/nbnxm/pairsearch.h" +#include "gromacs/pbcutil/pbc.h" +#include "gromacs/utility/logger.h" +#include "gromacs/utility/smalloc.h" +#include "nblib/exception.h" +#include "nblib/kerneloptions.h" +#include "nblib/particletype.h" +#include "nblib/simulationstate.h" + +namespace nblib +{ + +//! Helper to translate between the different enumeration values. +static Nbnxm::KernelType translateBenchmarkEnum(const SimdKernels& kernel) +{ + int kernelInt = static_cast(kernel); + return static_cast(kernelInt); +} + +/*! \brief Checks the kernel setup + * + * Returns an error string when the kernel is not available. + */ +static void checkKernelSetup(const NBKernelOptions& options) +{ + if (options.nbnxmSimd >= SimdKernels::Count && options.nbnxmSimd == SimdKernels::SimdAuto) + { + throw InputException("Need a valid kernel SIMD type"); + } + // Check SIMD support + if ((options.nbnxmSimd != SimdKernels::SimdNo && !GMX_SIMD) +#ifndef GMX_NBNXN_SIMD_4XN + || options.nbnxmSimd == SimdKernels::Simd4XM +#endif +#ifndef GMX_NBNXN_SIMD_2XNN + || options.nbnxmSimd == SimdKernels::Simd2XMM +#endif + ) + { + throw InputException("The requested SIMD kernel was not set up at configuration time"); + } +} + +NbvSetupUtil::NbvSetupUtil() : gmxForceCalculator_(std::make_unique()) {} + +void NbvSetupUtil::setExecutionContext(const NBKernelOptions& options) +{ + // Todo: find a more general way to initialize hardware + gmx_omp_nthreads_set(emntPairsearch, options.numOpenMPThreads); + gmx_omp_nthreads_set(emntNonbonded, options.numOpenMPThreads); +} + +Nbnxm::KernelSetup NbvSetupUtil::getKernelSetup(const NBKernelOptions& options) +{ + checkKernelSetup(options); + + Nbnxm::KernelSetup kernelSetup; + + // The int enum options.nbnxnSimd is set up to match Nbnxm::KernelType + 1 + kernelSetup.kernelType = translateBenchmarkEnum(options.nbnxmSimd); + // The plain-C kernel does not support analytical ewald correction + if (kernelSetup.kernelType == Nbnxm::KernelType::Cpu4x4_PlainC) + { + kernelSetup.ewaldExclusionType = Nbnxm::EwaldExclusionType::Table; + } + else + { + kernelSetup.ewaldExclusionType = options.useTabulatedEwaldCorr + ? Nbnxm::EwaldExclusionType::Table + : Nbnxm::EwaldExclusionType::Analytical; + } + + return kernelSetup; +} + +void NbvSetupUtil::setParticleInfoAllVdv(const size_t numParticles) + +{ + particleInfoAllVdw_.resize(numParticles); + for (size_t particleI = 0; particleI < numParticles; particleI++) + { + SET_CGINFO_HAS_VDW(particleInfoAllVdw_[particleI]); + SET_CGINFO_HAS_Q(particleInfoAllVdw_[particleI]); + } +} + +void NbvSetupUtil::setNonBondedParameters(const std::vector& particleTypes, + const NonBondedInteractionMap& nonBondedInteractionMap) +{ + /* Todo: Refactor nbnxm to take nonbondedParameters_ directly + * + * initial self-handling of combination rules + * size: 2*(numParticleTypes^2) + */ + nonbondedParameters_.reserve(2 * particleTypes.size() * particleTypes.size()); + + constexpr real c6factor = 6.0; + constexpr real c12factor = 12.0; + + for (const ParticleType& particleType1 : particleTypes) + { + for (const ParticleType& particleType2 : particleTypes) + { + nonbondedParameters_.push_back( + nonBondedInteractionMap.getC6(particleType1.name(), particleType2.name()) * c6factor); + nonbondedParameters_.push_back( + nonBondedInteractionMap.getC12(particleType1.name(), particleType2.name()) * c12factor); + } + } +} + +void NbvSetupUtil::setAtomProperties(const std::vector& particleTypeIdOfAllParticles, + const std::vector& charges) +{ + gmxForceCalculator_->nbv_->setAtomProperties(particleTypeIdOfAllParticles, charges, particleInfoAllVdw_); +} + +//! Sets up and returns a Nbnxm object for the given options and system +void NbvSetupUtil::setupNbnxmInstance(const size_t numParticleTypes, const NBKernelOptions& options) +{ + const auto pinPolicy = (options.useGpu ? gmx::PinningPolicy::PinnedIfSupported + : gmx::PinningPolicy::CannotBePinned); + const int numThreads = options.numOpenMPThreads; + // Note: the options and Nbnxm combination rule enums values should match + const int combinationRule = static_cast(options.ljCombinationRule); + + checkKernelSetup(options); // throws exception is setup is invalid + + Nbnxm::KernelSetup kernelSetup = getKernelSetup(options); + + PairlistParams pairlistParams(kernelSetup.kernelType, false, options.pairlistCutoff, false); + Nbnxm::GridSet gridSet(PbcType::Xyz, false, nullptr, nullptr, pairlistParams.pairlistType, + false, numThreads, pinPolicy); + auto pairlistSets = std::make_unique(pairlistParams, false, 0); + auto pairSearch = + std::make_unique(PbcType::Xyz, false, nullptr, nullptr, + pairlistParams.pairlistType, false, numThreads, pinPolicy); + + auto atomData = std::make_unique(pinPolicy); + + // Put everything together + auto nbv = std::make_unique(std::move(pairlistSets), std::move(pairSearch), + std::move(atomData), kernelSetup, nullptr, nullptr); + + // Needs to be called with the number of unique ParticleTypes + nbnxn_atomdata_init(gmx::MDLogger(), nbv->nbat.get(), kernelSetup.kernelType, combinationRule, + numParticleTypes, nonbondedParameters_, 1, numThreads); + + gmxForceCalculator_->nbv_ = std::move(nbv); +} + +//! Computes the Ewald splitting coefficient for Coulomb +static real ewaldCoeff(const real ewald_rtol, const real pairlistCutoff) +{ + return calc_ewaldcoeff_q(pairlistCutoff, ewald_rtol); +} + +void NbvSetupUtil::setupStepWorkload(const NBKernelOptions& options) +{ + gmxForceCalculator_->stepWork_->computeForces = true; + gmxForceCalculator_->stepWork_->computeNonbondedForces = true; + + if (options.computeVirialAndEnergy) + { + gmxForceCalculator_->stepWork_->computeVirial = true; + gmxForceCalculator_->stepWork_->computeEnergy = true; + } +} + +void NbvSetupUtil::setupInteractionConst(const NBKernelOptions& options) +{ + gmxForceCalculator_->interactionConst_->vdwtype = evdwCUT; + gmxForceCalculator_->interactionConst_->vdw_modifier = eintmodPOTSHIFT; + gmxForceCalculator_->interactionConst_->rvdw = options.pairlistCutoff; + + switch (options.coulombType) + { + case CoulombType::Pme: gmxForceCalculator_->interactionConst_->eeltype = eelPME; break; + case CoulombType::Cutoff: gmxForceCalculator_->interactionConst_->eeltype = eelCUT; break; + case CoulombType::ReactionField: + gmxForceCalculator_->interactionConst_->eeltype = eelRF; + break; + case CoulombType::Count: throw InputException("Unsupported electrostatic interaction"); + } + gmxForceCalculator_->interactionConst_->coulomb_modifier = eintmodPOTSHIFT; + gmxForceCalculator_->interactionConst_->rcoulomb = options.pairlistCutoff; + // Note: values correspond to ic->coulomb_modifier = eintmodPOTSHIFT + gmxForceCalculator_->interactionConst_->dispersion_shift.cpot = + -1.0 / gmx::power6(gmxForceCalculator_->interactionConst_->rvdw); + gmxForceCalculator_->interactionConst_->repulsion_shift.cpot = + -1.0 / gmx::power12(gmxForceCalculator_->interactionConst_->rvdw); + + // These are the initialized values but we leave them here so that later + // these can become options. + gmxForceCalculator_->interactionConst_->epsilon_r = 1.0; + gmxForceCalculator_->interactionConst_->epsilon_rf = 1.0; + + /* Set the Coulomb energy conversion factor */ + if (gmxForceCalculator_->interactionConst_->epsilon_r != 0) + { + gmxForceCalculator_->interactionConst_->epsfac = + ONE_4PI_EPS0 / gmxForceCalculator_->interactionConst_->epsilon_r; + } + else + { + /* eps = 0 is infinite dieletric: no Coulomb interactions */ + gmxForceCalculator_->interactionConst_->epsfac = 0; + } + + calc_rffac(nullptr, gmxForceCalculator_->interactionConst_->epsilon_r, + gmxForceCalculator_->interactionConst_->epsilon_rf, + gmxForceCalculator_->interactionConst_->rcoulomb, + &gmxForceCalculator_->interactionConst_->k_rf, + &gmxForceCalculator_->interactionConst_->c_rf); + + if (EEL_PME_EWALD(gmxForceCalculator_->interactionConst_->eeltype)) + { + // Ewald coefficients, we ignore the potential shift + gmxForceCalculator_->interactionConst_->ewaldcoeff_q = ewaldCoeff(1e-5, options.pairlistCutoff); + if (gmxForceCalculator_->interactionConst_->ewaldcoeff_q <= 0) + { + throw InputException("Ewald coefficient should be > 0"); + } + gmxForceCalculator_->interactionConst_->coulombEwaldTables = + std::make_unique(); + init_interaction_const_tables(nullptr, gmxForceCalculator_->interactionConst_.get(), 0); + } +} + +void NbvSetupUtil::setupForceRec(const matrix& box) +{ + assert((gmxForceCalculator_->forcerec_ && "Forcerec not initialized")); + gmxForceCalculator_->forcerec_->nbfp = nonbondedParameters_; + snew(gmxForceCalculator_->forcerec_->shift_vec, numShiftVectors); + calc_shifts(box, gmxForceCalculator_->forcerec_->shift_vec); +} + +void NbvSetupUtil::setParticlesOnGrid(const std::vector& coordinates, const Box& box) +{ + gmxForceCalculator_->setParticlesOnGrid(particleInfoAllVdw_, coordinates, box); +} + +void NbvSetupUtil::constructPairList(const gmx::ListOfLists& exclusions) +{ + gmxForceCalculator_->nbv_->constructPairlist(gmx::InteractionLocality::Local, exclusions, 0, + gmxForceCalculator_->nrnb_.get()); +} + + +std::unique_ptr GmxSetupDirector::setupGmxForceCalculator(const SimulationState& system, + const NBKernelOptions& options) +{ + NbvSetupUtil nbvSetupUtil; + nbvSetupUtil.setExecutionContext(options); + nbvSetupUtil.setNonBondedParameters(system.topology().getParticleTypes(), + system.topology().getNonBondedInteractionMap()); + nbvSetupUtil.setParticleInfoAllVdv(system.topology().numParticles()); + + nbvSetupUtil.setupInteractionConst(options); + nbvSetupUtil.setupStepWorkload(options); + nbvSetupUtil.setupNbnxmInstance(system.topology().getParticleTypes().size(), options); + nbvSetupUtil.setParticlesOnGrid(system.coordinates(), system.box()); + nbvSetupUtil.constructPairList(system.topology().getGmxExclusions()); + nbvSetupUtil.setAtomProperties(system.topology().getParticleTypeIdOfAllParticles(), + system.topology().getCharges()); + nbvSetupUtil.setupForceRec(system.box().legacyMatrix()); + + return nbvSetupUtil.getGmxForceCalculator(); +} + +} // namespace nblib diff --git a/api/nblib/gmxsetup.h b/api/nblib/gmxsetup.h new file mode 100644 index 0000000000..628172ca99 --- /dev/null +++ b/api/nblib/gmxsetup.h @@ -0,0 +1,127 @@ +/* + * 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. + */ +/*! \libinternal \file + * \brief Translation layer to GROMACS data structures for force calculations. + * + * Implements the translation layer between the user scope and + * GROMACS data structures for force calculations. Sets up the + * non-bonded verlet. + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ + +#ifndef NBLIB_GMXSETUP_H +#define NBLIB_GMXSETUP_H + +#include "nblib/gmxcalculator.h" +#include "nblib/simulationstate.h" + +namespace Nbnxm +{ +struct KernelSetup; +} + +namespace nblib +{ + +class NbvSetupUtil final +{ +public: + NbvSetupUtil(); + + //! Sets hardware params from the execution context + void setExecutionContext(const NBKernelOptions& options); + + //! Sets non-bonded parameters to be used to build GMX data structures + void setNonBondedParameters(const std::vector& particleTypes, + const NonBondedInteractionMap& nonBondedInteractionMap); + + //! Marks particles to have Van der Waals interactions + void setParticleInfoAllVdv(size_t numParticles); + + //! Returns the kernel setup + Nbnxm::KernelSetup getKernelSetup(const NBKernelOptions& options); + + //! Set up StepWorkload data + void setupStepWorkload(const NBKernelOptions& options); + + //! Return an interaction constants struct with members set appropriately + void setupInteractionConst(const NBKernelOptions& options); + + //! Sets Particle Types and Charges and VdW params + void setAtomProperties(const std::vector& particleTypeIdOfAllParticles, + const std::vector& charges); + + //! Sets up non-bonded verlet on the GmxForceCalculator + void setupNbnxmInstance(size_t numParticleTypes, const NBKernelOptions& options); + + //! Puts particles on a grid based on bounds specified by the box + void setParticlesOnGrid(const std::vector& coordinates, const Box& box); + + //! Constructs pair lists + void constructPairList(const gmx::ListOfLists& exclusions); + + //! Sets up t_forcerec object on the GmxForceCalculator + void setupForceRec(const matrix& box); + + std::unique_ptr getGmxForceCalculator() + { + return std::move(gmxForceCalculator_); + } + +private: + //! Storage for parameters for short range interactions. + std::vector nonbondedParameters_; + + //! Particle info where all particles are marked to have Van der Waals interactions + std::vector particleInfoAllVdw_; + + //! GROMACS force calculator to compute forces + std::unique_ptr gmxForceCalculator_; +}; + +class GmxSetupDirector +{ +public: + //! Sets up and returns a GmxForceCalculator + static std::unique_ptr setupGmxForceCalculator(const SimulationState& system, + const NBKernelOptions& options); +}; + +} // namespace nblib +#endif // NBLIB_GMXSETUP_H diff --git a/api/nblib/integrator.cpp b/api/nblib/integrator.cpp new file mode 100644 index 0000000000..970eea0c6a --- /dev/null +++ b/api/nblib/integrator.cpp @@ -0,0 +1,76 @@ +/* + * 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. + */ +/*! \internal \file + * \brief + * Implements nblib integrator + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + * \author Artem Zhmurov + */ +#include "nblib/integrator.h" +#include "gromacs/pbcutil/pbc.h" +#include "gromacs/utility/arrayref.h" +#include "nblib/topology.h" + +namespace nblib +{ +// NOLINTNEXTLINE(performance-unnecessary-value-param) +LeapFrog::LeapFrog(const Topology& topology, const Box& box) : box_(box) +{ + inverseMasses_.resize(topology.numParticles()); + for (int i = 0; i < topology.numParticles(); i++) + { + int typeIndex = topology.getParticleTypeIdOfAllParticles()[i]; + inverseMasses_[i] = 1.0 / topology.getParticleTypes()[typeIndex].mass(); + } +} + +void LeapFrog::integrate(const real dt, gmx::ArrayRef x, gmx::ArrayRef v, gmx::ArrayRef f) +{ + for (size_t i = 0; i < x.size(); i++) + { + for (int dim = 0; dim < dimSize; dim++) + { + v[i][dim] += f[i][dim] * dt * inverseMasses_[i]; + x[i][dim] += v[i][dim] * dt; + } + } + put_atoms_in_box(PbcType::Xyz, box_.legacyMatrix(), x); +} + +} // namespace nblib diff --git a/api/nblib/simulationstate.cpp b/api/nblib/simulationstate.cpp new file mode 100644 index 0000000000..9e38433644 --- /dev/null +++ b/api/nblib/simulationstate.cpp @@ -0,0 +1,146 @@ +/* + * 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. + */ +/*! \internal \file + * \brief + * Implements nblib SimulationState + * + * \author Berk Hess + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + * \author Artem Zhmurov + */ +#include + +#include "gromacs/pbcutil/pbc.h" +#include "nblib/exception.h" +#include "nblib/simulationstate.h" +#include "nblib/simulationstateimpl.h" + +namespace nblib +{ + +SimulationState::SimulationState(const std::vector& coordinates, + const std::vector& velocities, + const std::vector& forces, + Box box, + Topology topology) : + simulationStatePtr_(std::make_shared(coordinates, velocities, forces, box, topology)) +{ +} + +SimulationState::Impl::Impl(const std::vector& coordinates, + const std::vector& velocities, + const std::vector& forces, + const Box& box, + Topology topology) : + box_(box), + topology_(std::move(topology)) +{ + if (!checkNumericValues(coordinates)) + { + throw InputException("Input coordinates has at least one NaN"); + } + coordinates_ = coordinates; + if (!checkNumericValues(velocities)) + { + throw InputException("Input velocities has at least one NaN"); + } + + velocities_ = velocities; + + forces_ = forces; + + // Ensure that the coordinates are in a box following PBC + put_atoms_in_box(PbcType::Xyz, box.legacyMatrix(), coordinates_); +} + +const Topology& SimulationState::Impl::topology() const +{ + return topology_; +} + +Box SimulationState::Impl::box() const +{ + return box_; +} + +std::vector& SimulationState::Impl::coordinates() +{ + return coordinates_; +} + +std::vector& SimulationState::Impl::velocities() +{ + return velocities_; +} + +std::vector& SimulationState::Impl::forces() +{ + return forces_; +} + +const Topology& SimulationState::topology() const +{ + return simulationStatePtr_->topology(); +} + +Box SimulationState::box() const +{ + return simulationStatePtr_->box(); +} + +std::vector& SimulationState::coordinates() +{ + return simulationStatePtr_->coordinates(); +} + +const std::vector& SimulationState::coordinates() const +{ + return simulationStatePtr_->coordinates(); +} + +std::vector& SimulationState::velocities() +{ + return simulationStatePtr_->velocities(); +} + +std::vector& SimulationState::forces() +{ + return simulationStatePtr_->forces(); +} + +} // namespace nblib diff --git a/api/nblib/simulationstateimpl.h b/api/nblib/simulationstateimpl.h new file mode 100644 index 0000000000..8ed604bb8a --- /dev/null +++ b/api/nblib/simulationstateimpl.h @@ -0,0 +1,87 @@ +/* + * 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. + */ +/*! \inpublicapi \file + * \brief + * Implements nblib SimulationState + * + * \author Berk Hess + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + * \author Artem Zhmurov + */ +#ifndef NBLIB_SIMULATIONSTATE_IMPL_H +#define NBLIB_SIMULATIONSTATE_IMPL_H + +namespace nblib +{ + +class SimulationState::Impl final +{ +public: + //! Constructor + Impl(const std::vector& coordinates, + const std::vector& velocities, + const std::vector& forces, + const Box& box, + Topology topology); + + //! Returns topology of the current state + const Topology& topology() const; + + //! Returns the box + Box box() const; + + //! Returns a vector of particle coordinates + std::vector& coordinates(); + + //! Returns a vector of particle velocities + std::vector& velocities(); + + //! Returns a vector of forces + std::vector& forces(); + +private: + std::vector coordinates_ = {}; + std::vector velocities_ = {}; + std::vector forces_ = {}; + Box box_; + Topology topology_; +}; + +} // namespace nblib + +#endif // NBLIB_SIMULATIONSTATE_IMPL_H diff --git a/api/nblib/tests/CMakeLists.txt b/api/nblib/tests/CMakeLists.txt index f7b17dd26d..366027339c 100644 --- a/api/nblib/tests/CMakeLists.txt +++ b/api/nblib/tests/CMakeLists.txt @@ -64,3 +64,35 @@ target_link_libraries(${exename} PRIVATE nblib_test_infrastructure nblib) target_include_directories(${exename} PRIVATE ${PROJECT_SOURCE_DIR}/api) gmx_register_gtest_test(${testname} ${exename} INTEGRATION_TEST) add_dependencies(check-nblib ${exename}) + +set(testname "NbLibIntegrationTests") +set(exename "nblib-integration-test") + +gmx_add_gtest_executable( + ${exename} + CPP_SOURCE_FILES + # files with code for tests + gmxcalculator.cpp + nbkernelsystem.cpp + nbnxnsetup.cpp + simstate.cpp + ) +target_link_libraries(${exename} PRIVATE nblib_test_infrastructure nblib) +target_include_directories(${exename} PRIVATE ${PROJECT_SOURCE_DIR}/api) +gmx_register_gtest_test(${testname} ${exename} INTEGRATION_TEST) +add_dependencies(check-nblib ${exename}) + +# The integrator sometimes times out on TSAN so it gets its own test harness +set(testname "NbLibIntegratorTests") +set(exename "nblib-integrator-test") + +gmx_add_gtest_executable( + ${exename} + CPP_SOURCE_FILES + integrator.cpp +) +target_link_libraries(${exename} PRIVATE nblib_test_infrastructure nblib) +target_include_directories(${exename} PRIVATE ${PROJECT_SOURCE_DIR}/api) +gmx_register_gtest_test(${testname} ${exename} INTEGRATION_TEST) +add_dependencies(check-nblib ${exename}) + diff --git a/api/nblib/tests/gmxcalculator.cpp b/api/nblib/tests/gmxcalculator.cpp new file mode 100644 index 0000000000..c50b2c83d8 --- /dev/null +++ b/api/nblib/tests/gmxcalculator.cpp @@ -0,0 +1,81 @@ +/* + * 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. + */ +/*! \internal \file + * \brief + * This implements basic nblib utility tests + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ +#include + +#include "nblib/gmxsetup.h" +#include "nblib/kerneloptions.h" +#include "nblib/simulationstate.h" +#include "nblib/tests/testsystems.h" + +namespace nblib +{ +namespace test +{ +namespace +{ +TEST(NBlibTest, GmxForceCalculatorCanCompute) +{ + ArgonSimulationStateBuilder argonSystemBuilder; + SimulationState simState = argonSystemBuilder.setupSimulationState(); + NBKernelOptions options = NBKernelOptions(); + options.nbnxmSimd = SimdKernels::SimdNo; + std::unique_ptr gmxForceCalculator = + nblib::GmxSetupDirector::setupGmxForceCalculator(simState, options); + EXPECT_NO_THROW(gmxForceCalculator->compute(simState.coordinates(), simState.forces())); +} + +TEST(NBlibTest, CanSetupStepWorkload) +{ + NBKernelOptions options; + EXPECT_NO_THROW(NbvSetupUtil{}.setupStepWorkload(options)); +} + +TEST(NBlibTest, GmxForceCalculatorCanSetupInteractionConst) +{ + NBKernelOptions options; + EXPECT_NO_THROW(NbvSetupUtil{}.setupInteractionConst(options)); +} +} // namespace +} // namespace test +} // namespace nblib diff --git a/api/nblib/tests/integrator.cpp b/api/nblib/tests/integrator.cpp new file mode 100644 index 0000000000..a1ee57fe0f --- /dev/null +++ b/api/nblib/tests/integrator.cpp @@ -0,0 +1,137 @@ +/* + * 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. + */ +/*! \internal \file + * \brief + * This implements molecule setup tests + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + * \author Artem Zhmurov + */ +#include "nblib/integrator.h" +#include "gromacs/pbcutil/pbc.h" +#include "nblib/molecules.h" +#include "nblib/particletype.h" +#include "nblib/simulationstate.h" +#include "nblib/topology.h" +#include "nblib/util/internal.h" + +#include "testutils/testasserts.h" + +namespace nblib +{ +namespace test +{ +namespace +{ + +TEST(NBlibTest, IntegratorWorks) +{ + int numAtoms = 1; + int numSteps = 100; + real dt = 0.001; + + ParticleType particleType(ParticleTypeName("H"), Mass(1.0)); + Molecule molecule(MoleculeName("SomeMolecule")); + molecule.addParticle(ParticleName("SomeAtom"), particleType); + + ParticleTypesInteractions interactions; + interactions.add(particleType.name(), C6{ 0 }, C12{ 0 }); + + TopologyBuilder topologyBuilder; + topologyBuilder.addMolecule(molecule, numAtoms); + topologyBuilder.addParticleTypesInteractions(interactions); + Topology topology = topologyBuilder.buildTopology(); + + std::vector x(numAtoms); + std::vector v(numAtoms); + std::vector f(numAtoms); + + f[0][XX] = 1.0; + f[0][YY] = 2.0; + f[0][ZZ] = 0.0; + + Box box(100); + + std::vector x0(x); + std::vector v0(v); + + SimulationState simulationState(x, v, f, box, topology); + put_atoms_in_box(PbcType::Xyz, box.legacyMatrix(), x0); + + LeapFrog integrator(simulationState.topology(), simulationState.box()); + + gmx::test::FloatingPointTolerance tolerance = gmx::test::absoluteTolerance(numSteps * 0.000005); + for (int step = 0; step < numSteps; step++) + { + real totalTime = step * dt; + + Vec3 xAnalytical; + Vec3 vAnalytical; + + for (int i = 0; i < numAtoms; i++) + { + for (int d = 0; d < dimSize; d++) + { + // Analytical solution for constant-force particle movement + int typeIndex = simulationState.topology().getParticleTypeIdOfAllParticles()[i]; + real im = 1.0 / simulationState.topology().getParticleTypes()[typeIndex].mass(); + xAnalytical[d] = + x0[i][d] + v0[i][d] * totalTime + 0.5 * f[i][d] * totalTime * totalTime * im; + vAnalytical[d] = v0[i][d] + f[i][d] * totalTime * im; + + EXPECT_REAL_EQ_TOL(xAnalytical[d], simulationState.coordinates()[i][d], tolerance) + << formatString( + "Coordinate {} of atom {} is different from analytical solution " + "at step {}.", + d, i, step); + + EXPECT_REAL_EQ_TOL(vAnalytical[d], simulationState.velocities()[i][d], tolerance) + << formatString( + "Velocity component {} of atom {} is different from analytical " + "solution at step {}.", + d, i, step); + } + integrator.integrate(dt, simulationState.coordinates(), simulationState.velocities(), + simulationState.forces()); + } + } +} + +} // namespace +} // namespace test +} // namespace nblib diff --git a/api/nblib/tests/nbkernelsystem.cpp b/api/nblib/tests/nbkernelsystem.cpp new file mode 100644 index 0000000000..1d022f5694 --- /dev/null +++ b/api/nblib/tests/nbkernelsystem.cpp @@ -0,0 +1,207 @@ +/* + * 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. + */ +/*! \internal \file + * \brief + * This implements topology setup tests + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + * \author Artem Zhmurov + */ +#include + +#include "gromacs/topology/exclusionblocks.h" +#include "nblib/forcecalculator.h" +#include "nblib/gmxsetup.h" +#include "nblib/integrator.h" +#include "nblib/tests/testhelpers.h" +#include "nblib/tests/testsystems.h" +#include "nblib/topology.h" + +namespace nblib +{ +namespace test +{ +namespace +{ + +// This is defined in src/gromacs/mdtypes/forcerec.h but there is also a +// legacy C6 macro defined there that conflicts with the nblib C6 type. +// Todo: Once that C6 has been refactored into a regular function, this +// file can just include forcerec.h +//! Macro to set Van der Waals interactions to atoms +#define SET_CGINFO_HAS_VDW(cgi) (cgi) = ((cgi) | (1 << 23)) + +TEST(NBlibTest, SpcMethanolForcesAreCorrect) +{ + auto options = NBKernelOptions(); + options.nbnxmSimd = SimdKernels::SimdNo; + options.coulombType = CoulombType::Cutoff; + + SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder; + + auto simState = spcMethanolSystemBuilder.setupSimulationState(); + auto forceCalculator = ForceCalculator(simState, options); + + gmx::ArrayRef forces(simState.forces()); + ASSERT_NO_THROW(forceCalculator.compute(simState.coordinates(), forces)); + + Vector3DTest forcesOutputTest; + forcesOutputTest.testVectors(forces, "SPC-methanol forces"); +} + +TEST(NBlibTest, ExpectedNumberOfForces) +{ + auto options = NBKernelOptions(); + options.nbnxmSimd = SimdKernels::SimdNo; + + SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder; + + auto simState = spcMethanolSystemBuilder.setupSimulationState(); + auto forceCalculator = ForceCalculator(simState, options); + + gmx::ArrayRef forces(simState.forces()); + forceCalculator.compute(simState.coordinates(), forces); + EXPECT_EQ(simState.topology().numParticles(), forces.size()); +} + +TEST(NBlibTest, CanIntegrateSystem) +{ + auto options = NBKernelOptions(); + options.nbnxmSimd = SimdKernels::SimdNo; + options.numIterations = 1; + + SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder; + + auto simState = spcMethanolSystemBuilder.setupSimulationState(); + auto forceCalculator = ForceCalculator(simState, options); + + LeapFrog integrator(simState.topology(), simState.box()); + + for (int iter = 0; iter < options.numIterations; iter++) + { + gmx::ArrayRef forces(simState.forces()); + forceCalculator.compute(simState.coordinates(), simState.forces()); + EXPECT_NO_THROW(integrator.integrate(1.0, simState.coordinates(), simState.velocities(), + simState.forces())); + } +} + +/*! + * Check if the following aspects of the ForceCalculator and + * LeapFrog (integrator) work as expected: + * + * 1. Calling the ForceCalculator::compute() function makes no change + * to the internal representation of the system. Calling it repeatedly + * without update should return the same vector of forces. + * + * 2. Once the LeapFrog objects integrates for the given time using the + * forces, there the coordinates in SimulationState must change. + * Calling the compute() function must now generate a new set of forces. + * + */ +TEST(NBlibTest, UpdateChangesForces) +{ + auto options = NBKernelOptions(); + options.nbnxmSimd = SimdKernels::SimdNo; + options.numIterations = 1; + + SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder; + + auto simState = spcMethanolSystemBuilder.setupSimulationState(); + auto forceCalculator = ForceCalculator(simState, options); + + LeapFrog integrator(simState.topology(), simState.box()); + + // step 1 + gmx::ArrayRef forces(simState.forces()); + forceCalculator.compute(simState.coordinates(), simState.forces()); + + std::vector forces_1(forces.size()); + std::copy(forces.begin(), forces.end(), begin(forces_1)); + + // check if forces change without update step + forceCalculator.compute(simState.coordinates(), forces); + + // check if forces change without update + for (size_t i = 0; i < forces_1.size(); i++) + { + for (int j = 0; j < dimSize; j++) + { + EXPECT_EQ(forces[i][j], forces_1[i][j]); + } + } + + // update + integrator.integrate(1.0, simState.coordinates(), simState.velocities(), simState.forces()); + + // step 2 + forceCalculator.compute(simState.coordinates(), simState.forces()); + std::vector forces_2(forces.size()); + std::copy(forces.begin(), forces.end(), begin(forces_2)); + + // check if forces change after update + for (size_t i = 0; i < forces_1.size(); i++) + { + for (int j = 0; j < dimSize; j++) + { + EXPECT_NE(forces_1[i][j], forces_2[i][j]); + } + } +} + +TEST(NBlibTest, ArgonForcesAreCorrect) +{ + auto options = NBKernelOptions(); + options.nbnxmSimd = SimdKernels::SimdNo; + options.coulombType = CoulombType::Cutoff; + + ArgonSimulationStateBuilder argonSystemBuilder; + + auto simState = argonSystemBuilder.setupSimulationState(); + auto forceCalculator = ForceCalculator(simState, options); + + gmx::ArrayRef testForces(simState.forces()); + forceCalculator.compute(simState.coordinates(), simState.forces()); + + Vector3DTest forcesOutputTest; + forcesOutputTest.testVectors(testForces, "Argon forces"); +} + +} // namespace +} // namespace test +} // namespace nblib diff --git a/api/nblib/tests/nbnxnsetup.cpp b/api/nblib/tests/nbnxnsetup.cpp new file mode 100644 index 0000000000..ada24d5044 --- /dev/null +++ b/api/nblib/tests/nbnxnsetup.cpp @@ -0,0 +1,65 @@ +/* + * 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. + */ +/*! \internal \file + * \brief + * This implements nbnxn setup tests + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ +#include "nblib/gmxsetup.h" +#include "nblib/tests/testhelpers.h" +#include "nblib/tests/testsystems.h" + +#include "testutils/testasserts.h" + +namespace nblib +{ + +namespace test +{ + +TEST(NBlibTest, CanConstructNbvSetupUtil) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + SimulationState system = argonSimulationStateBuilder.setupSimulationState(); + EXPECT_NO_THROW(NbvSetupUtil()); +} + +} // namespace test + +} // namespace nblib diff --git a/api/nblib/tests/refdata/NBlibTest_ArgonForcesAreCorrect.xml b/api/nblib/tests/refdata/NBlibTest_ArgonForcesAreCorrect.xml new file mode 100644 index 0000000000..ddaa587c23 --- /dev/null +++ b/api/nblib/tests/refdata/NBlibTest_ArgonForcesAreCorrect.xml @@ -0,0 +1,67 @@ + + + + + 12 + + -0.41298821868424429 + -1.0982427445010643 + -0.11318936363938557 + + + 0 + 0 + 0 + + + 0.41298821868424429 + 1.0982427445010643 + 0.11318936363938557 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + 0 + 0 + 0 + + + diff --git a/api/nblib/tests/refdata/NBlibTest_SpcMethanolForcesAreCorrect.xml b/api/nblib/tests/refdata/NBlibTest_SpcMethanolForcesAreCorrect.xml new file mode 100644 index 0000000000..9f812d2ad1 --- /dev/null +++ b/api/nblib/tests/refdata/NBlibTest_SpcMethanolForcesAreCorrect.xml @@ -0,0 +1,37 @@ + + + + + 6 + + -0.38183025207307253 + 0.87922789140299651 + -6.1406494302583212 + + + 8.3033346966856811 + -7.3388265216240605 + 27.783758884487042 + + + -9.4221671563754512 + 6.2920799125341684 + -33.986135630631082 + + + 27.494322535207331 + 8.3914914141170893 + 39.693768750211909 + + + -19.009852871290562 + -5.3975006967248795 + -15.545639878395072 + + + -6.983806952153925 + -2.8264719997053138 + -11.805102695414478 + + + diff --git a/api/nblib/tests/simstate.cpp b/api/nblib/tests/simstate.cpp new file mode 100644 index 0000000000..01a0bf5432 --- /dev/null +++ b/api/nblib/tests/simstate.cpp @@ -0,0 +1,154 @@ +/* + * 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. + */ +/*! \internal \file + * \brief + * This implements SimulationState tests + * + * \author Victor Holanda + * \author Joe Jordan + * \author Prashanth Kanduri + * \author Sebastian Keller + */ +#include + +#include "nblib/box.h" +#include "nblib/exception.h" +#include "nblib/simulationstate.h" +#include "nblib/simulationstateimpl.h" +#include "nblib/tests/testhelpers.h" +#include "nblib/tests/testsystems.h" +#include "nblib/topology.h" + +#include "testutils/testasserts.h" + +namespace nblib +{ +namespace test +{ +namespace +{ + +//! Utility function to compare 2 std::vectors of gmx::RVec used to compare cartesians +void compareValues(const std::vector& ref, const std::vector& test) +{ + for (size_t i = 0; i < ref.size(); i++) + { + for (int j = 0; j < dimSize; j++) + { + EXPECT_EQ(ref[i][j], test.at(i)[j]); + } + } +} + +TEST(NBlibTest, CanConstructSimulationState) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + EXPECT_NO_THROW(argonSimulationStateBuilder.setupSimulationState()); +} + +TEST(NBlibTest, SimulationStateThrowsCoordinateNAN) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + argonSimulationStateBuilder.setCoordinate(2, 0, NAN); + EXPECT_THROW(argonSimulationStateBuilder.setupSimulationState(), InputException); +} + +TEST(NBlibTest, SimulationStateThrowsCoordinateINF) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + argonSimulationStateBuilder.setCoordinate(2, 0, INFINITY); + EXPECT_THROW(argonSimulationStateBuilder.setupSimulationState(), InputException); +} + +TEST(NBlibTest, SimulationStateThrowsVelocityNAN) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + argonSimulationStateBuilder.setVelocity(2, 0, NAN); + EXPECT_THROW(argonSimulationStateBuilder.setupSimulationState(), InputException); +} + +TEST(NBlibTest, SimulationStateThrowsVelocityINF) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + argonSimulationStateBuilder.setVelocity(2, 0, INFINITY); + EXPECT_THROW(argonSimulationStateBuilder.setupSimulationState(), InputException); +} + +TEST(NBlibTest, SimulationStateCanMove) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + SimulationState simState = argonSimulationStateBuilder.setupSimulationState(); + EXPECT_NO_THROW(SimulationState movedSimState = std::move(simState)); +} + +TEST(NBlibTest, SimulationStateCanAssign) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + SimulationState simState = argonSimulationStateBuilder.setupSimulationState(); + EXPECT_NO_THROW(const SimulationState& gmx_unused AssignedSimState = simState); +} + +TEST(NBlibTest, SimulationStateHasBox) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + SimulationState simState = argonSimulationStateBuilder.setupSimulationState(); + const Box& testBox = simState.box(); + const Box& refBox = argonSimulationStateBuilder.box(); + // GTEST does not like the comparison operator in a different namespace + const bool compare = (refBox == testBox); + EXPECT_TRUE(compare); +} + +TEST(NBlibTest, SimulationStateHasCorrectCoordinates) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + SimulationState simState = argonSimulationStateBuilder.setupSimulationState(); + std::vector test = simState.coordinates(); + std::vector ref = argonSimulationStateBuilder.coordinates(); + compareValues(ref, test); +} + +TEST(NBlibTest, SimulationStateHasCorrectVelocities) +{ + ArgonSimulationStateBuilder argonSimulationStateBuilder; + SimulationState simState = argonSimulationStateBuilder.setupSimulationState(); + std::vector test = simState.velocities(); + std::vector ref = argonSimulationStateBuilder.velocities(); + compareValues(ref, test); +} + +} // namespace +} // namespace test +} // namespace nblib diff --git a/api/nblib/tests/testsystems.cpp b/api/nblib/tests/testsystems.cpp index 02abe3e20f..70e86dd813 100644 --- a/api/nblib/tests/testsystems.cpp +++ b/api/nblib/tests/testsystems.cpp @@ -241,5 +241,112 @@ Topology ArgonTopologyBuilder::argonTopology() return topologyBuilder_.buildTopology(); } +ArgonSimulationStateBuilder::ArgonSimulationStateBuilder() : + box_(6.05449), + topology_(ArgonTopologyBuilder(12).argonTopology()) +{ + + coordinates_ = { + { 0.794, 1.439, 0.610 }, { 1.397, 0.673, 1.916 }, { 0.659, 1.080, 0.573 }, + { 1.105, 0.090, 3.431 }, { 1.741, 1.291, 3.432 }, { 1.936, 1.441, 5.873 }, + { 0.960, 2.246, 1.659 }, { 0.382, 3.023, 2.793 }, { 0.053, 4.857, 4.242 }, + { 2.655, 5.057, 2.211 }, { 4.114, 0.737, 0.614 }, { 5.977, 5.104, 5.217 }, + }; + + velocities_ = { + { 0.0055, -0.1400, 0.2127 }, { 0.0930, -0.0160, -0.0086 }, { 0.1678, 0.2476, -0.0660 }, + { 0.1591, -0.0934, -0.0835 }, { -0.0317, 0.0573, 0.1453 }, { 0.0597, 0.0013, -0.0462 }, + { 0.0484, -0.0357, 0.0168 }, { 0.0530, 0.0295, -0.2694 }, { -0.0550, -0.0896, 0.0494 }, + { -0.0799, -0.2534, -0.0079 }, { 0.0436, -0.1557, 0.1849 }, { -0.0214, 0.0446, 0.0758 }, + }; + forces_ = { + { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, + { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, + { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, + { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, { 0.0000, 0.0000, 0.0000 }, + }; +} + +void ArgonSimulationStateBuilder::setCoordinate(int particleNum, int dimension, real value) +{ + if (dimension < 0 or dimension > 2) + { + throw InputException("Must provide a valid dimension\n"); + } + coordinates_.at(particleNum)[dimension] = value; +} + +void ArgonSimulationStateBuilder::setVelocity(int particleNum, int dimension, real value) +{ + if (dimension < 0 or dimension > 2) + { + throw InputException("Must provide a valid dimension\n"); + } + velocities_.at(particleNum)[dimension] = value; +} + +SimulationState ArgonSimulationStateBuilder::setupSimulationState() +{ + return SimulationState(coordinates_, velocities_, forces_, box_, topology_); +} + +const Topology& ArgonSimulationStateBuilder::topology() const +{ + return topology_; +} + +Box& ArgonSimulationStateBuilder::box() +{ + return box_; +} + +std::vector& ArgonSimulationStateBuilder::coordinates() +{ + return coordinates_; +} + +std::vector& ArgonSimulationStateBuilder::velocities() +{ + return velocities_; +} + +SpcMethanolSimulationStateBuilder::SpcMethanolSimulationStateBuilder() : + box_(3.01000), + topology_(SpcMethanolTopologyBuilder().buildTopology(1, 1)) +{ + coordinates_ = { + { 1.970, 1.460, 1.209 }, // Me1 + { 1.978, 1.415, 1.082 }, // O2 + { 1.905, 1.460, 1.030 }, // H3 + { 1.555, 1.511, 0.703 }, // Ow + { 1.498, 1.495, 0.784 }, // Hw1 + { 1.496, 1.521, 0.623 }, // Hw2 + }; + + velocities_ = { + { -0.8587, -0.1344, -0.0643 }, { 0.0623, -0.1787, 0.0036 }, { -0.5020, -0.9564, 0.0997 }, + { 0.869, 1.245, 1.665 }, { 0.169, 0.275, 1.565 }, { 0.269, 2.275, 1.465 }, + }; + + forces_ = { + { 0.000, 0.000, 0.000 }, { 0.000, 0.000, 0.000 }, { 0.000, 0.000, 0.000 }, + { 0.000, 0.000, 0.000 }, { 0.000, 0.000, 0.000 }, { 0.000, 0.000, 0.000 }, + }; +} + +SimulationState SpcMethanolSimulationStateBuilder::setupSimulationState() +{ + return SimulationState(coordinates_, velocities_, forces_, box_, topology_); +} + +std::vector& SpcMethanolSimulationStateBuilder::coordinates() +{ + return coordinates_; +} + +std::vector& SpcMethanolSimulationStateBuilder::velocities() +{ + return velocities_; +} } // namespace nblib diff --git a/api/nblib/tests/testsystems.h b/api/nblib/tests/testsystems.h index e2f7cfe666..f4d4dace63 100644 --- a/api/nblib/tests/testsystems.h +++ b/api/nblib/tests/testsystems.h @@ -158,5 +158,65 @@ private: TopologyBuilder topologyBuilder_; }; +//! \internal \brief Build simulation state for the argon example +class ArgonSimulationStateBuilder +{ +public: + ArgonSimulationStateBuilder(); + + //! Set coordinates of particles in the defined system + void setCoordinate(int particleNum, int dimension, real value); + + //! Set particle velocities + void setVelocity(int particleNum, int dimension, real value); + + //! Setup simulation state + SimulationState setupSimulationState(); + + //! Get the topology + const Topology& topology() const; + + //! Get the box bounding the system + Box& box(); + + //! Get current coordinates + std::vector& coordinates(); + + //! Get current velocities + std::vector& velocities(); + +private: + std::vector coordinates_; + std::vector velocities_; + std::vector forces_; + + Box box_; + Topology topology_; +}; + +//! \internal \brief Build simulation state for the SPC-Methanol example +class SpcMethanolSimulationStateBuilder +{ +public: + SpcMethanolSimulationStateBuilder(); + + //! Setup simulation state + SimulationState setupSimulationState(); + + //! Get current coordinates + std::vector& coordinates(); + + //! Get current velocities + std::vector& velocities(); + +private: + std::vector coordinates_; + std::vector velocities_; + std::vector forces_; + + Box box_; + Topology topology_; +}; + } // namespace nblib #endif // NBLIB_TESTSYSTEMS_H -- 2.22.0