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
)
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+#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<const Vec3> coordinates, gmx::ArrayRef<Vec3> forces)
+{
+ return gmxForceCalculator_->compute(coordinates, forces);
+}
+
+void ForceCalculator::updatePairList(gmx::ArrayRef<const int> particleInfoAllVdW,
+ gmx::ArrayRef<Vec3> coordinates,
+ const Box& box)
+{
+ gmxForceCalculator_->setParticlesOnGrid(particleInfoAllVdW, coordinates, box);
+}
+
+} // namespace nblib
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+#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<gmx_enerdata_t>(1, 0);
+ forcerec_ = std::make_unique<t_forcerec>();
+ interactionConst_ = std::make_unique<interaction_const_t>();
+ stepWork_ = std::make_unique<gmx::StepWorkload>();
+ nrnb_ = std::make_unique<t_nrnb>();
+}
+
+GmxForceCalculator::~GmxForceCalculator() = default;
+
+void GmxForceCalculator::compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ gmx::ArrayRef<gmx::RVec> 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<const int> particleInfoAllVdw,
+ gmx::ArrayRef<const gmx::RVec> 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
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+
+#ifndef NBLIB_GMXCALCULATOR_H
+#define NBLIB_GMXCALCULATOR_H
+
+#include <memory>
+
+#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<typename T>
+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<const gmx::RVec> coordinateInput, gmx::ArrayRef<gmx::RVec> forceOutput);
+
+ //! Puts particles on a grid based on bounds specified by the box (for every NS step)
+ void setParticlesOnGrid(gmx::ArrayRef<const int> particleInfoAllVdw,
+ gmx::ArrayRef<const gmx::RVec> coordinates,
+ const Box& box);
+
+private:
+ friend class NbvSetupUtil;
+
+ //! Non-Bonded Verlet object for force calculation
+ std::unique_ptr<nonbonded_verlet_t> nbv_;
+
+ //! Only nbfp and shift_vec are used
+ std::unique_ptr<t_forcerec> forcerec_;
+
+ //! Parameters for various interactions in the system
+ std::unique_ptr<interaction_const_t> interactionConst_;
+
+ //! Tasks to perform in an MD Step
+ std::unique_ptr<gmx::StepWorkload> stepWork_;
+
+ //! Energies of different interaction types; currently only needed as an argument for dispatchNonbondedKernel
+ std::unique_ptr<gmx_enerdata_t> enerd_;
+
+ //! Non-bonded flop counter; currently only needed as an argument for dispatchNonbondedKernel
+ std::unique_ptr<t_nrnb> nrnb_;
+
+ //! Legacy matrix for box
+ matrix box_{ { 0 } };
+};
+
+} // namespace nblib
+
+#endif // NBLIB_GMXCALCULATOR_H
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+#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<int>(kernel);
+ return static_cast<Nbnxm::KernelType>(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<GmxForceCalculator>()) {}
+
+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<ParticleType>& 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<int>& particleTypeIdOfAllParticles,
+ const std::vector<real>& 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<int>(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<PairlistSets>(pairlistParams, false, 0);
+ auto pairSearch =
+ std::make_unique<PairSearch>(PbcType::Xyz, false, nullptr, nullptr,
+ pairlistParams.pairlistType, false, numThreads, pinPolicy);
+
+ auto atomData = std::make_unique<nbnxn_atomdata_t>(pinPolicy);
+
+ // Put everything together
+ auto nbv = std::make_unique<nonbonded_verlet_t>(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<EwaldCorrectionTables>();
+ 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<Vec3>& coordinates, const Box& box)
+{
+ gmxForceCalculator_->setParticlesOnGrid(particleInfoAllVdw_, coordinates, box);
+}
+
+void NbvSetupUtil::constructPairList(const gmx::ListOfLists<int>& exclusions)
+{
+ gmxForceCalculator_->nbv_->constructPairlist(gmx::InteractionLocality::Local, exclusions, 0,
+ gmxForceCalculator_->nrnb_.get());
+}
+
+
+std::unique_ptr<GmxForceCalculator> 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
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+
+#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<ParticleType>& 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<int>& particleTypeIdOfAllParticles,
+ const std::vector<real>& 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<Vec3>& coordinates, const Box& box);
+
+ //! Constructs pair lists
+ void constructPairList(const gmx::ListOfLists<int>& exclusions);
+
+ //! Sets up t_forcerec object on the GmxForceCalculator
+ void setupForceRec(const matrix& box);
+
+ std::unique_ptr<GmxForceCalculator> getGmxForceCalculator()
+ {
+ return std::move(gmxForceCalculator_);
+ }
+
+private:
+ //! Storage for parameters for short range interactions.
+ std::vector<real> nonbondedParameters_;
+
+ //! Particle info where all particles are marked to have Van der Waals interactions
+ std::vector<int> particleInfoAllVdw_;
+
+ //! GROMACS force calculator to compute forces
+ std::unique_ptr<GmxForceCalculator> gmxForceCalculator_;
+};
+
+class GmxSetupDirector
+{
+public:
+ //! Sets up and returns a GmxForceCalculator
+ static std::unique_ptr<GmxForceCalculator> setupGmxForceCalculator(const SimulationState& system,
+ const NBKernelOptions& options);
+};
+
+} // namespace nblib
+#endif // NBLIB_GMXSETUP_H
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ */
+#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<Vec3> x, gmx::ArrayRef<Vec3> v, gmx::ArrayRef<const Vec3> 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
--- /dev/null
+/*
+ * 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 <hess@kth.se>
+ * \author Victor Holanda <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ */
+#include <vector>
+
+#include "gromacs/pbcutil/pbc.h"
+#include "nblib/exception.h"
+#include "nblib/simulationstate.h"
+#include "nblib/simulationstateimpl.h"
+
+namespace nblib
+{
+
+SimulationState::SimulationState(const std::vector<Vec3>& coordinates,
+ const std::vector<Vec3>& velocities,
+ const std::vector<Vec3>& forces,
+ Box box,
+ Topology topology) :
+ simulationStatePtr_(std::make_shared<Impl>(coordinates, velocities, forces, box, topology))
+{
+}
+
+SimulationState::Impl::Impl(const std::vector<Vec3>& coordinates,
+ const std::vector<Vec3>& velocities,
+ const std::vector<Vec3>& 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<Vec3>& SimulationState::Impl::coordinates()
+{
+ return coordinates_;
+}
+
+std::vector<Vec3>& SimulationState::Impl::velocities()
+{
+ return velocities_;
+}
+
+std::vector<Vec3>& SimulationState::Impl::forces()
+{
+ return forces_;
+}
+
+const Topology& SimulationState::topology() const
+{
+ return simulationStatePtr_->topology();
+}
+
+Box SimulationState::box() const
+{
+ return simulationStatePtr_->box();
+}
+
+std::vector<Vec3>& SimulationState::coordinates()
+{
+ return simulationStatePtr_->coordinates();
+}
+
+const std::vector<Vec3>& SimulationState::coordinates() const
+{
+ return simulationStatePtr_->coordinates();
+}
+
+std::vector<Vec3>& SimulationState::velocities()
+{
+ return simulationStatePtr_->velocities();
+}
+
+std::vector<Vec3>& SimulationState::forces()
+{
+ return simulationStatePtr_->forces();
+}
+
+} // namespace nblib
--- /dev/null
+/*
+ * 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 <hess@kth.se>
+ * \author Victor Holanda <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ */
+#ifndef NBLIB_SIMULATIONSTATE_IMPL_H
+#define NBLIB_SIMULATIONSTATE_IMPL_H
+
+namespace nblib
+{
+
+class SimulationState::Impl final
+{
+public:
+ //! Constructor
+ Impl(const std::vector<Vec3>& coordinates,
+ const std::vector<Vec3>& velocities,
+ const std::vector<Vec3>& 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<Vec3>& coordinates();
+
+ //! Returns a vector of particle velocities
+ std::vector<Vec3>& velocities();
+
+ //! Returns a vector of forces
+ std::vector<Vec3>& forces();
+
+private:
+ std::vector<Vec3> coordinates_ = {};
+ std::vector<Vec3> velocities_ = {};
+ std::vector<Vec3> forces_ = {};
+ Box box_;
+ Topology topology_;
+};
+
+} // namespace nblib
+
+#endif // NBLIB_SIMULATIONSTATE_IMPL_H
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})
+
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+#include <gtest/gtest.h>
+
+#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> 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
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ */
+#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<Vec3> x(numAtoms);
+ std::vector<Vec3> v(numAtoms);
+ std::vector<Vec3> f(numAtoms);
+
+ f[0][XX] = 1.0;
+ f[0][YY] = 2.0;
+ f[0][ZZ] = 0.0;
+
+ Box box(100);
+
+ std::vector<Vec3> x0(x);
+ std::vector<Vec3> 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
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ */
+#include <gtest/gtest.h>
+
+#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<Vec3> 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<Vec3> 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<Vec3> 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<Vec3> forces(simState.forces());
+ forceCalculator.compute(simState.coordinates(), simState.forces());
+
+ std::vector<Vec3> 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<Vec3> 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<Vec3> testForces(simState.forces());
+ forceCalculator.compute(simState.coordinates(), simState.forces());
+
+ Vector3DTest forcesOutputTest;
+ forcesOutputTest.testVectors(testForces, "Argon forces");
+}
+
+} // namespace
+} // namespace test
+} // namespace nblib
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+#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
--- /dev/null
+<?xml version="1.0"?>
+<?xml-stylesheet type="text/xsl" href="referencedata.xsl"?>
+<ReferenceData>
+ <Sequence Name="Argon forces">
+ <Int Name="Length">12</Int>
+ <Vector>
+ <Real Name="X">-0.41298821868424429</Real>
+ <Real Name="Y">-1.0982427445010643</Real>
+ <Real Name="Z">-0.11318936363938557</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0.41298821868424429</Real>
+ <Real Name="Y">1.0982427445010643</Real>
+ <Real Name="Z">0.11318936363938557</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">0</Real>
+ <Real Name="Y">0</Real>
+ <Real Name="Z">0</Real>
+ </Vector>
+ </Sequence>
+</ReferenceData>
--- /dev/null
+<?xml version="1.0"?>
+<?xml-stylesheet type="text/xsl" href="referencedata.xsl"?>
+<ReferenceData>
+ <Sequence Name="SPC-methanol forces">
+ <Int Name="Length">6</Int>
+ <Vector>
+ <Real Name="X">-0.38183025207307253</Real>
+ <Real Name="Y">0.87922789140299651</Real>
+ <Real Name="Z">-6.1406494302583212</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">8.3033346966856811</Real>
+ <Real Name="Y">-7.3388265216240605</Real>
+ <Real Name="Z">27.783758884487042</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">-9.4221671563754512</Real>
+ <Real Name="Y">6.2920799125341684</Real>
+ <Real Name="Z">-33.986135630631082</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">27.494322535207331</Real>
+ <Real Name="Y">8.3914914141170893</Real>
+ <Real Name="Z">39.693768750211909</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">-19.009852871290562</Real>
+ <Real Name="Y">-5.3975006967248795</Real>
+ <Real Name="Z">-15.545639878395072</Real>
+ </Vector>
+ <Vector>
+ <Real Name="X">-6.983806952153925</Real>
+ <Real Name="Y">-2.8264719997053138</Real>
+ <Real Name="Z">-11.805102695414478</Real>
+ </Vector>
+ </Sequence>
+</ReferenceData>
--- /dev/null
+/*
+ * 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 <victor.holanda@cscs.ch>
+ * \author Joe Jordan <ejjordan@kth.se>
+ * \author Prashanth Kanduri <kanduri@cscs.ch>
+ * \author Sebastian Keller <keller@cscs.ch>
+ */
+#include <vector>
+
+#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<Vec3>& ref, const std::vector<Vec3>& 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<Vec3> test = simState.coordinates();
+ std::vector<Vec3> ref = argonSimulationStateBuilder.coordinates();
+ compareValues(ref, test);
+}
+
+TEST(NBlibTest, SimulationStateHasCorrectVelocities)
+{
+ ArgonSimulationStateBuilder argonSimulationStateBuilder;
+ SimulationState simState = argonSimulationStateBuilder.setupSimulationState();
+ std::vector<Vec3> test = simState.velocities();
+ std::vector<Vec3> ref = argonSimulationStateBuilder.velocities();
+ compareValues(ref, test);
+}
+
+} // namespace
+} // namespace test
+} // namespace nblib
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<Vec3>& ArgonSimulationStateBuilder::coordinates()
+{
+ return coordinates_;
+}
+
+std::vector<Vec3>& 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<Vec3>& SpcMethanolSimulationStateBuilder::coordinates()
+{
+ return coordinates_;
+}
+
+std::vector<Vec3>& SpcMethanolSimulationStateBuilder::velocities()
+{
+ return velocities_;
+}
} // namespace nblib
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<Vec3>& coordinates();
+
+ //! Get current velocities
+ std::vector<Vec3>& velocities();
+
+private:
+ std::vector<Vec3> coordinates_;
+ std::vector<Vec3> velocities_;
+ std::vector<Vec3> 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<Vec3>& coordinates();
+
+ //! Get current velocities
+ std::vector<Vec3>& velocities();
+
+private:
+ std::vector<Vec3> coordinates_;
+ std::vector<Vec3> velocities_;
+ std::vector<Vec3> forces_;
+
+ Box box_;
+ Topology topology_;
+};
+
} // namespace nblib
#endif // NBLIB_TESTSYSTEMS_H
biod.pnpi.spb.ru / alexxy / gromacs.git / commitdiff