# Ensure that "make tests" builds all nblib tests so the top-level
# "make check" will work.
if (BUILD_TESTING)
- add_dependencies(tests nblib-tests)
+ add_dependencies(tests nblib-tests)
- # this allows all nblib tests to be run with "make check-nblib"
- add_custom_target(check-nblib
- COMMAND ${CMAKE_CTEST_COMMAND} --output-on-failure -R NbLib
- COMMENT "Running nblib tests"
- USES_TERMINAL VERBATIM)
- add_dependencies(check-nblib nblib-tests)
+ # this allows all nblib tests to be run with "make check-nblib"
+ add_custom_target(check-nblib
+ COMMAND ${CMAKE_CTEST_COMMAND} --output-on-failure -R NbLib
+ COMMENT "Running nblib tests"
+ USES_TERMINAL VERBATIM)
+ add_dependencies(check-nblib nblib-tests)
endif()
set(NBLIB_MAJOR 0)
target_sources(nblib
PRIVATE
box.cpp
- gmxcalculator.cpp
- gmxsetup.cpp
+ gmxcalculatorcpu.cpp
integrator.cpp
interactions.cpp
molecules.cpp
basicdefinitions.h
box.h
exception.h
- gmxcalculator.h
+ gmxcalculatorcpu.h
integrator.h
interactions.h
molecules.h
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2020,2021, by the GROMACS development team, led by
+ * Copyright (c) 2021, 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.
* 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
+/*! \inpublicapi \file
* \brief
- * Implements a force calculator based on GROMACS data structures.
- *
- * Intended for internal use inside the ForceCalculator.
+ * Implements nblib simulation box
*
* \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
+#ifndef NBLIB_GMXBACKENDDATA_H
+#define NBLIB_GMXBACKENDDATA_H
+
+#include "gromacs/gmxlib/nrnb.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/timing/wallcycle.h"
+#include "gromacs/utility/listoflists.h"
+#include "gromacs/utility/range.h"
+#include "nblib/kerneloptions.h"
+#include "nblib/nbnxmsetuphelpers.h"
namespace nblib
{
-class Box;
-class NbvSetupUtil;
-class SimulationState;
-struct NBKernelOptions;
/*! \brief GROMACS non-bonded force calculation backend
*
* interactions (PP) on a single node CPU.
*
*/
-
-class GmxForceCalculator final
+class GmxBackendData
{
public:
- GmxForceCalculator();
+ GmxBackendData() = default;
+ GmxBackendData(const NBKernelOptions& options,
+ int numEnergyGroups,
+ gmx::ArrayRef<int> exclusionRanges,
+ gmx::ArrayRef<int> exclusionElements)
+ {
+ // Set hardware params from the execution context
+ setGmxNonBondedNThreads(options.numOpenMPThreads);
+
+ // Set interaction constants struct
+ interactionConst_ = createInteractionConst(options);
- ~GmxForceCalculator();
+ // Set up StepWorkload data
+ stepWork_ = createStepWorkload();
- //! Compute forces and return
- void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput, gmx::ArrayRef<gmx::RVec> forceOutput);
+ // Set up gmx_enerdata_t (holds energy information)
+ enerd_ = gmx_enerdata_t{ numEnergyGroups, 0 };
- //! Puts particles on a grid based on bounds specified by the box (for every NS step)
- void setParticlesOnGrid(gmx::ArrayRef<const int64_t> particleInfoAllVdw,
- gmx::ArrayRef<const gmx::RVec> coordinates,
- const Box& box);
+ // Construct pair lists
+ std::vector<int> exclusionRanges_(exclusionRanges.begin(), exclusionRanges.end());
+ std::vector<int> exclusionElements_(exclusionElements.begin(), exclusionElements.end());
+ exclusions_ = gmx::ListOfLists<int>(std::move(exclusionRanges_), std::move(exclusionElements_));
+ }
-private:
- //! Friend to allow setting up private members in this class
- friend class NbvSetupUtil;
+ //! exclusions in gmx format
+ gmx::ListOfLists<int> exclusions_;
//! 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_;
+ //! Only shift_vec is used
+ t_forcerec forcerec_;
//! Parameters for various interactions in the system
- std::unique_ptr<interaction_const_t> interactionConst_;
+ interaction_const_t interactionConst_;
//! Tasks to perform in an MD Step
- std::unique_ptr<gmx::StepWorkload> stepWork_;
+ gmx::StepWorkload stepWork_;
+
+ gmx::SimulationWorkload simulationWork_;
//! Energies of different interaction types; currently only needed as an argument for dispatchNonbondedKernel
- std::unique_ptr<gmx_enerdata_t> enerd_;
+ gmx_enerdata_t enerd_{ 1, 0 };
//! Non-bonded flop counter; currently only needed as an argument for dispatchNonbondedKernel
- std::unique_ptr<t_nrnb> nrnb_;
+ t_nrnb nrnb_;
- //! Legacy matrix for box
- matrix box_{ { 0 } };
+ //! Keep track of whether updatePairlist has been called at least once
+ bool updatePairlistCalled{ false };
};
} // namespace nblib
-
-#endif // NBLIB_GMXCALCULATOR_H
+#endif // NBLIB_GMXBACKENDDATA_H
+++ /dev/null
-/*
- * This file is part of the GROMACS molecular simulation package.
- *
- * Copyright (c) 2020,2021, 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, coordinateInput);
-
- nbv_->dispatchNonbondedKernel(
- gmx::InteractionLocality::Local,
- *interactionConst_,
- *stepWork_,
- enbvClearFYes,
- forcerec_->shift_vec,
- enerd_->grpp.energyGroupPairTerms[forcerec_->haveBuckingham ? NonBondedEnergyTerms::BuckinghamSR
- : NonBondedEnergyTerms::LJSR],
- enerd_->grpp.energyGroupPairTerms[NonBondedEnergyTerms::CoulombSR],
- nrnb_.get());
-
- nbv_->atomdata_add_nbat_f_to_f(gmx::AtomLocality::All, forceOutput);
-}
-
-void GmxForceCalculator::setParticlesOnGrid(gmx::ArrayRef<const int64_t> 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,2021, 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 "gromacs/ewald/ewald_utils.h"
+#include "gromacs/mdtypes/enerdata.h"
+#include "gromacs/mdtypes/interaction_const.h"
+#include "gromacs/nbnxm/atomdata.h"
+#include "gromacs/nbnxm/nbnxm.h"
+#include "gromacs/nbnxm/pairlistset.h"
+#include "gromacs/nbnxm/pairlistsets.h"
+#include "gromacs/nbnxm/pairsearch.h"
+#include "gromacs/utility/listoflists.h"
+#include "gromacs/utility/range.h"
+#include "nblib/exception.h"
+#include "nblib/gmxbackenddata.h"
+#include "nblib/gmxcalculatorcpu.h"
+#include "nblib/nbnxmsetuphelpers.h"
+#include "nblib/pbc.hpp"
+#include "nblib/systemdescription.h"
+#include "nblib/topology.h"
+#include "nblib/virials.h"
+
+namespace nblib
+{
+
+class GmxNBForceCalculatorCpu::CpuImpl final
+{
+public:
+ CpuImpl(gmx::ArrayRef<int> particleTypeIdOfAllParticles,
+ gmx::ArrayRef<real> nonBondedParams,
+ gmx::ArrayRef<real> charges,
+ gmx::ArrayRef<int64_t> particleInteractionFlags,
+ gmx::ArrayRef<int> exclusionRanges,
+ gmx::ArrayRef<int> exclusionElements,
+ const NBKernelOptions& options);
+
+ //! calculates a new pair list based on new coordinates (for every NS step)
+ void updatePairlist(gmx::ArrayRef<const gmx::RVec> coordinates, const Box& box);
+
+ //! Compute forces and return
+ void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput);
+
+ //! Compute forces and virial tensor
+ void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput);
+
+ //! Compute forces, virial tensor and potential energies
+ void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput,
+ gmx::ArrayRef<real> energyOutput);
+
+private:
+ //! \brief client-side provided system description data
+ SystemDescription system_;
+
+ //! \brief Gmx backend objects, employed for calculating the forces
+ GmxBackendData backend_;
+};
+
+GmxNBForceCalculatorCpu::CpuImpl::CpuImpl(gmx::ArrayRef<int> particleTypeIdOfAllParticles,
+ gmx::ArrayRef<real> nonBondedParams,
+ gmx::ArrayRef<real> charges,
+ gmx::ArrayRef<int64_t> particleInteractionFlags,
+ gmx::ArrayRef<int> exclusionRanges,
+ gmx::ArrayRef<int> exclusionElements,
+ const NBKernelOptions& options) :
+ system_(SystemDescription(particleTypeIdOfAllParticles, nonBondedParams, charges, particleInteractionFlags)),
+ backend_(GmxBackendData(options, findNumEnergyGroups(particleInteractionFlags), exclusionRanges, exclusionElements))
+{
+ // Set up non-bonded verlet in the backend
+ backend_.nbv_ = createNbnxmCPU(system_.numParticleTypes_,
+ options,
+ findNumEnergyGroups(particleInteractionFlags),
+ system_.nonBondedParams_);
+}
+
+void GmxNBForceCalculatorCpu::CpuImpl::updatePairlist(gmx::ArrayRef<const gmx::RVec> coordinates,
+ const Box& box)
+{
+ if (coordinates.size() != system_.numParticles_)
+ {
+ throw InputException(
+ "Coordinate array containing different number of entries than particles in the "
+ "system");
+ }
+
+ const auto* legacyBox = box.legacyMatrix();
+ system_.box_ = box;
+ updateForcerec(&backend_.forcerec_, 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 = static_cast<real>(coordinates.size()) / det(legacyBox);
+
+ // Put particles on a grid based on bounds specified by the box
+ nbnxn_put_on_grid(backend_.nbv_.get(),
+ legacyBox,
+ 0,
+ lowerCorner,
+ upperCorner,
+ nullptr,
+ { 0, int(coordinates.size()) },
+ particleDensity,
+ system_.particleInfo_,
+ coordinates,
+ 0,
+ nullptr);
+
+ backend_.nbv_->constructPairlist(
+ gmx::InteractionLocality::Local, backend_.exclusions_, 0, &backend_.nrnb_);
+
+ // Set Particle Types and Charges and VdW params
+ backend_.nbv_->setAtomProperties(
+ system_.particleTypeIdOfAllParticles_, system_.charges_, system_.particleInfo_);
+ backend_.updatePairlistCalled = true;
+}
+
+void GmxNBForceCalculatorCpu::CpuImpl::compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput,
+ gmx::ArrayRef<real> energyOutput)
+{
+ if (coordinateInput.size() != forceOutput.size())
+ {
+ throw InputException("coordinate array and force buffer size mismatch");
+ }
+
+ if (!backend_.updatePairlistCalled)
+ {
+ throw InputException("compute called without updating pairlist at least once");
+ }
+
+ // update the box if changed
+ if (!(system_.box_ == box))
+ {
+ system_.box_ = box;
+ updateForcerec(&backend_.forcerec_, box.legacyMatrix());
+ }
+
+ bool computeVirial = !virialOutput.empty();
+ bool computeEnergies = !energyOutput.empty();
+ backend_.stepWork_.computeVirial = computeVirial;
+ backend_.stepWork_.computeEnergy = computeEnergies;
+
+ // update the coordinates in the backend
+ backend_.nbv_->convertCoordinates(gmx::AtomLocality::Local, coordinateInput);
+
+ backend_.nbv_->dispatchNonbondedKernel(
+ gmx::InteractionLocality::Local,
+ backend_.interactionConst_,
+ backend_.stepWork_,
+ enbvClearFYes,
+ backend_.forcerec_.shift_vec,
+ backend_.enerd_.grpp.energyGroupPairTerms[backend_.forcerec_.haveBuckingham ? NonBondedEnergyTerms::BuckinghamSR
+ : NonBondedEnergyTerms::LJSR],
+ backend_.enerd_.grpp.energyGroupPairTerms[NonBondedEnergyTerms::CoulombSR],
+ &backend_.nrnb_);
+
+ backend_.nbv_->atomdata_add_nbat_f_to_f(gmx::AtomLocality::All, forceOutput);
+
+ if (computeVirial)
+ {
+ // calculate shift forces and turn into an array ref
+ std::vector<Vec3> shiftForcesVector(gmx::c_numShiftVectors, Vec3(0.0, 0.0, 0.0));
+ nbnxn_atomdata_add_nbat_fshift_to_fshift(*backend_.nbv_->nbat, shiftForcesVector);
+ auto shiftForcesRef = constArrayRefFromArray(shiftForcesVector.data(), shiftForcesVector.size());
+
+ std::vector<Vec3> shiftVectorsArray(gmx::c_numShiftVectors);
+
+ // copy shift vectors from ForceRec
+ std::copy(backend_.forcerec_.shift_vec.begin(),
+ backend_.forcerec_.shift_vec.end(),
+ shiftVectorsArray.begin());
+
+ computeVirialTensor(
+ coordinateInput, forceOutput, shiftVectorsArray, shiftForcesRef, box, virialOutput);
+ }
+
+ // extract term energies (per interaction type)
+ if (computeEnergies)
+ {
+ int nGroupPairs = backend_.enerd_.grpp.nener;
+ if (int(energyOutput.size()) != int(NonBondedEnergyTerms::Count) * nGroupPairs)
+ {
+ throw InputException("Array size for energy output is wrong\n");
+ }
+
+ for (int eg = 0; eg < int(NonBondedEnergyTerms::Count); ++eg)
+ {
+ std::copy(begin(backend_.enerd_.grpp.energyGroupPairTerms[eg]),
+ end(backend_.enerd_.grpp.energyGroupPairTerms[eg]),
+ energyOutput.begin() + eg * nGroupPairs);
+ }
+ }
+}
+
+void GmxNBForceCalculatorCpu::CpuImpl::compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput)
+{
+ // compute forces and fill in force buffer
+ compute(coordinateInput, box, forceOutput, gmx::ArrayRef<real>{}, gmx::ArrayRef<real>{});
+}
+
+void GmxNBForceCalculatorCpu::CpuImpl::compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput)
+{
+ // compute forces and fill in force buffer
+ compute(coordinateInput, box, forceOutput, virialOutput, gmx::ArrayRef<real>{});
+}
+
+GmxNBForceCalculatorCpu::GmxNBForceCalculatorCpu(gmx::ArrayRef<int> particleTypeIdOfAllParticles,
+ gmx::ArrayRef<real> nonBondedParams,
+ gmx::ArrayRef<real> charges,
+ gmx::ArrayRef<int64_t> particleInteractionFlags,
+ gmx::ArrayRef<int> exclusionRanges,
+ gmx::ArrayRef<int> exclusionElements,
+ const NBKernelOptions& options)
+{
+ if (options.useGpu)
+ {
+ throw InputException("Use GmxNBForceCalculatorGpu for GPU support");
+ }
+
+ impl_ = std::make_unique<CpuImpl>(particleTypeIdOfAllParticles,
+ nonBondedParams,
+ charges,
+ particleInteractionFlags,
+ exclusionRanges,
+ exclusionElements,
+ options);
+}
+
+GmxNBForceCalculatorCpu::~GmxNBForceCalculatorCpu() = default;
+
+//! calculates a new pair list based on new coordinates (for every NS step)
+void GmxNBForceCalculatorCpu::updatePairlist(gmx::ArrayRef<const gmx::RVec> coordinates, const Box& box)
+{
+ impl_->updatePairlist(coordinates, box);
+}
+
+//! Compute forces and return
+void GmxNBForceCalculatorCpu::compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput)
+{
+ impl_->compute(coordinateInput, box, forceOutput);
+}
+
+//! Compute forces and virial tensor
+void GmxNBForceCalculatorCpu::compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput)
+{
+ impl_->compute(coordinateInput, box, forceOutput, virialOutput);
+}
+
+//! Compute forces, virial tensor and potential energies
+void GmxNBForceCalculatorCpu::compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput,
+ gmx::ArrayRef<real> energyOutput)
+{
+ impl_->compute(coordinateInput, box, forceOutput, virialOutput, energyOutput);
+}
+
+std::unique_ptr<GmxNBForceCalculatorCpu> setupGmxForceCalculatorCpu(const Topology& topology,
+ const NBKernelOptions& options)
+{
+ std::vector<real> nonBondedParameters = createNonBondedParameters(
+ topology.getParticleTypes(), topology.getNonBondedInteractionMap());
+
+ std::vector<int64_t> particleInteractionFlags = createParticleInfoAllVdw(topology.numParticles());
+
+ return std::make_unique<GmxNBForceCalculatorCpu>(topology.getParticleTypeIdOfAllParticles(),
+ nonBondedParameters,
+ topology.getCharges(),
+ particleInteractionFlags,
+ topology.exclusionLists().ListRanges,
+ topology.exclusionLists().ListElements,
+ options);
+}
+
+} // namespace nblib
--- /dev/null
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2020,2021, 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_GMXCALCULATORCPU_H
+#define NBLIB_GMXCALCULATORCPU_H
+
+#include <memory>
+#include <vector>
+
+#include "nblib/box.h"
+#include "nblib/vector.h"
+
+namespace gmx
+{
+template<typename T>
+class ArrayRef;
+} // namespace gmx
+
+namespace nblib
+{
+struct NBKernelOptions;
+class Topology;
+
+class GmxNBForceCalculatorCpu final
+{
+public:
+ GmxNBForceCalculatorCpu(gmx::ArrayRef<int> particleTypeIdOfAllParticles,
+ gmx::ArrayRef<real> nonBondedParams,
+ gmx::ArrayRef<real> charges,
+ gmx::ArrayRef<int64_t> particleInteractionFlags,
+ gmx::ArrayRef<int> exclusionRanges,
+ gmx::ArrayRef<int> exclusionElements,
+ const NBKernelOptions& options);
+
+ ~GmxNBForceCalculatorCpu();
+
+ //! calculates a new pair list based on new coordinates (for every NS step)
+ void updatePairlist(gmx::ArrayRef<const gmx::RVec> coordinates, const Box& box);
+
+ //! Compute forces and return
+ void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput);
+
+ //! Compute forces and virial tensor
+ void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput);
+
+ //! Compute forces, virial tensor and potential energies
+ void compute(gmx::ArrayRef<const gmx::RVec> coordinateInput,
+ const Box& box,
+ gmx::ArrayRef<gmx::RVec> forceOutput,
+ gmx::ArrayRef<real> virialOutput,
+ gmx::ArrayRef<real> energyOutput);
+
+private:
+ //! Private implementation
+ class CpuImpl;
+ std::unique_ptr<CpuImpl> impl_;
+};
+
+std::unique_ptr<GmxNBForceCalculatorCpu> setupGmxForceCalculatorCpu(const Topology& topology,
+ const NBKernelOptions& options);
+
+} // namespace nblib
+
+#endif // NBLIB_GMXCALCULATORCPU_H
#include "nblib/basicdefinitions.h"
#include "nblib/box.h"
-#include "nblib/gmxcalculator.h"
-#include "nblib/gmxsetup.h"
+#include "nblib/gmxcalculatorcpu.h"
#include "nblib/integrator.h"
#include "nblib/interactions.h"
#include "nblib/kerneloptions.h"
#include "nblib/listed_forces/definitions.h"
#include "nblib/listed_forces/calculator.h"
#include "nblib/molecules.h"
+#include "nblib/nbnxmsetuphelpers.h"
#include "nblib/particlesequencer.h"
#include "nblib/particletype.h"
#include "nblib/simulationstate.h"
return nonbondedParameters;
}
-gmx::StepWorkload createStepWorkload([[maybe_unused]] const NBKernelOptions& options)
+gmx::StepWorkload createStepWorkload()
{
gmx::StepWorkload stepWorkload;
stepWorkload.computeForces = true;
const NonBondedInteractionMap& nonBondedInteractionMap);
//! Create a step work object
-gmx::StepWorkload createStepWorkload(const NBKernelOptions& options);
+gmx::StepWorkload createStepWorkload();
//! Computes the Ewald splitting coefficient for Coulomb
real ewaldCoeff(real ewald_rtol, real pairlistCutoff);
// Some performance flags can be set a run time
options.nbnxmSimd = nblib::SimdKernels::SimdNo;
// The force calculator contains all the data needed to compute forces.
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
+ auto forceCalculator = nblib::setupGmxForceCalculatorCpu(simState.topology(), options);
+ // build the pairlist
+ forceCalculator->updatePairlist(simState.coordinates(), simState.box());
// Integration requires masses, positions, and forces
nblib::LeapFrog integrator(simState.topology(), simState.box());
// Print some diagnostic info
printf("initial forces on particle 0: x %4f y %4f z %4f\n", forces[0][0], forces[0][1], forces[0][2]);
// The forces are computed for the user
gmx::ArrayRef<nblib::Vec3> userForces(simState.forces());
- forceCalculator->compute(simState.coordinates(), userForces);
+ forceCalculator->compute(simState.coordinates(), simState.box(), userForces);
// Print some diagnostic info
printf(" final forces on particle 0: x %4f y %4f z %4f\n",
userForces[0][0],
SimulationState simulationState(coordinates, velocities, forces, box, topology);
// The non-bonded force calculator contains all the data needed to compute forces
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simulationState.topology(), simulationState.coordinates(), simulationState.box(), options);
+ auto forceCalculator = setupGmxForceCalculatorCpu(simulationState.topology(), options);
+
+ // build the pair list
+ forceCalculator->updatePairlist(simulationState.coordinates(), simulationState.box());
// The listed force calculator is also initialized with the required arguments
ListedForceCalculator listedForceCalculator(
// Integrator is initialized with an array of inverse masses (constructed from topology) and
// the bounding box
- LeapFrog integrator(topology, box);
+ LeapFrog integrator(simulationState.topology(), simulationState.box());
// Print some diagnostic info
printf("initial position of particle 0: x %4f y %4f z %4f\n",
{
zeroCartesianArray(simulationState.forces());
- forceCalculator->compute(simulationState.coordinates(), simulationState.forces());
+ forceCalculator->compute(
+ simulationState.coordinates(), simulationState.box(), simulationState.forces());
listedForceCalculator.compute(simulationState.coordinates(), simulationState.forces());
--- /dev/null
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2021, 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 simulation box
+ *
+ * \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_SYSTEMDESCRIPTION_H
+#define NBLIB_SYSTEMDESCRIPTION_H
+
+#include "gromacs/ewald/ewald_utils.h"
+#include "gromacs/gmxlib/nrnb.h"
+#include "gromacs/gpu_utils/device_stream_manager.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/gpu_data_mgmt.h"
+#include "gromacs/nbnxm/nbnxm_gpu.h"
+#include "gromacs/nbnxm/atomdata.h"
+#include "gromacs/nbnxm/nbnxm.h"
+#include "gromacs/nbnxm/pairlistset.h"
+#include "gromacs/nbnxm/pairlistsets.h"
+#include "gromacs/nbnxm/pairsearch.h"
+#include "gromacs/timing/wallcycle.h"
+#include "gromacs/utility/listoflists.h"
+#include "gromacs/utility/range.h"
+
+namespace nblib
+{
+
+//! \brief client-side provided system description data
+struct SystemDescription
+{
+ SystemDescription() = default;
+
+ SystemDescription(gmx::ArrayRef<int> particleTypeIdOfAllParticles,
+ gmx::ArrayRef<real> nonBondedParams,
+ gmx::ArrayRef<real> charges,
+ gmx::ArrayRef<int64_t> particleInteractionFlags)
+ {
+ std::array inputSizes{ particleTypeIdOfAllParticles.size(),
+ charges.size(),
+ particleInteractionFlags.size() };
+ if (static_cast<unsigned long>(std::count(begin(inputSizes), end(inputSizes), inputSizes[0]))
+ != inputSizes.size())
+ {
+ throw InputException("input array size inconsistent");
+ }
+
+ int numParticleTypes = int(std::round(std::sqrt(nonBondedParams.size() / 2)));
+ if (2 * numParticleTypes * numParticleTypes != int(nonBondedParams.size()))
+ {
+ throw InputException("Wrong size of nonBondedParams");
+ }
+
+ numParticles_ = particleTypeIdOfAllParticles.size();
+ numParticleTypes_ = numParticleTypes;
+
+ particleTypeIdOfAllParticles_ = std::vector<int>(particleTypeIdOfAllParticles.begin(),
+ particleTypeIdOfAllParticles.end());
+
+ nonBondedParams_ = std::vector<real>(nonBondedParams.begin(), nonBondedParams.end());
+ charges_ = std::vector<real>(charges.begin(), charges.end());
+ particleInfo_ =
+ std::vector<int64_t>(particleInteractionFlags.begin(), particleInteractionFlags.end());
+ }
+
+ //! number of particles
+ size_t numParticles_{ 0 };
+
+ //! number of particle types
+ size_t numParticleTypes_{ 0 };
+
+ //! particle type id of all particles
+ std::vector<int> particleTypeIdOfAllParticles_;
+
+ //! Storage for parameters for short range interactions.
+ std::vector<real> nonBondedParams_;
+
+ //! electrostatic charges
+ std::vector<real> charges_;
+
+ //! flag for each particle to set LJ and Q interactions
+ std::vector<int64_t> particleInfo_;
+
+ //! Legacy matrix for box
+ Box box_{ 0 };
+};
+
+
+} // namespace nblib
+#endif // NBLIB_SYSTEMDESCRIPTION_H
*/
#include <gtest/gtest.h>
-#include "nblib/gmxsetup.h"
+#include "nblib/gmxcalculatorcpu.h"
#include "nblib/kerneloptions.h"
#include "nblib/simulationstate.h"
+#include "nblib/tests/testhelpers.h"
#include "nblib/tests/testsystems.h"
#include "gromacs/utility/arrayref.h"
SimulationState simState = argonSystemBuilder.setupSimulationState();
NBKernelOptions options = NBKernelOptions();
options.nbnxmSimd = SimdKernels::SimdNo;
- std::unique_ptr<GmxForceCalculator> gmxForceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
- ;
- EXPECT_NO_THROW(gmxForceCalculator->compute(simState.coordinates(), simState.forces()));
+
+ std::unique_ptr<GmxNBForceCalculatorCpu> gmxForceCalculator =
+ setupGmxForceCalculatorCpu(simState.topology(), options);
+ gmxForceCalculator->updatePairlist(simState.coordinates(), simState.box());
+
+ EXPECT_NO_THROW(gmxForceCalculator->compute(simState.coordinates(), simState.box(), simState.forces()));
}
-TEST(NBlibTest, CanSetupStepWorkload)
+TEST(NBlibTest, ArgonVirialsAreCorrect)
{
- NBKernelOptions options;
- EXPECT_NO_THROW(NbvSetupUtil{}.setupStepWorkload(options));
+ ArgonSimulationStateBuilder argonSystemBuilder(fftypes::OPLSA);
+ SimulationState simState = argonSystemBuilder.setupSimulationState();
+ NBKernelOptions options = NBKernelOptions();
+ options.nbnxmSimd = SimdKernels::SimdNo;
+ std::unique_ptr<GmxNBForceCalculatorCpu> gmxForceCalculator =
+ setupGmxForceCalculatorCpu(simState.topology(), options);
+ gmxForceCalculator->updatePairlist(simState.coordinates(), simState.box());
+
+ std::vector<real> virialArray(9, 0.0);
+
+ gmxForceCalculator->compute(simState.coordinates(), simState.box(), simState.forces(), virialArray);
+
+ RefDataChecker virialsOutputTest(1e-7);
+ virialsOutputTest.testArrays<real>(virialArray, "Virials");
}
-TEST(NBlibTest, GmxForceCalculatorCanSetupInteractionConst)
+TEST(NBlibTest, ArgonEnergiesAreCorrect)
{
- NBKernelOptions options;
- EXPECT_NO_THROW(NbvSetupUtil{}.setupInteractionConst(options));
+ ArgonSimulationStateBuilder argonSystemBuilder(fftypes::OPLSA);
+ SimulationState simState = argonSystemBuilder.setupSimulationState();
+ NBKernelOptions options = NBKernelOptions();
+ options.nbnxmSimd = SimdKernels::SimdNo;
+ std::unique_ptr<GmxNBForceCalculatorCpu> gmxForceCalculator =
+ setupGmxForceCalculatorCpu(simState.topology(), options);
+ gmxForceCalculator->updatePairlist(simState.coordinates(), simState.box());
+
+ // number of energy kinds is 5: COULSR, LJSR, BHAMSR, COUL14, LJ14,
+ std::vector<real> energies(5, 0.0);
+
+ gmxForceCalculator->compute(
+ simState.coordinates(), simState.box(), simState.forces(), gmx::ArrayRef<real>{}, energies);
+
+ RefDataChecker energiesOutputTest(5e-5);
+ energiesOutputTest.testArrays<real>(energies, "Argon energies");
+}
+
+TEST(NBlibTest, SpcMethanolEnergiesAreCorrect)
+{
+ SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder;
+ SimulationState simState = spcMethanolSystemBuilder.setupSimulationState();
+ NBKernelOptions options = NBKernelOptions();
+ options.nbnxmSimd = SimdKernels::SimdNo;
+ std::unique_ptr<GmxNBForceCalculatorCpu> gmxForceCalculator =
+ setupGmxForceCalculatorCpu(simState.topology(), options);
+ gmxForceCalculator->updatePairlist(simState.coordinates(), simState.box());
+
+ // number of energy kinds is 5: COULSR, LJSR, BHAMSR, COUL14, LJ14,
+ std::vector<real> energies(5, 0.0);
+
+ gmxForceCalculator->compute(
+ simState.coordinates(), simState.box(), simState.forces(), gmx::ArrayRef<real>{}, energies);
+
+ RefDataChecker energiesOutputTest(5e-5);
+ energiesOutputTest.testArrays<real>(energies, "SPC-methanol energies");
}
+
} // namespace
} // namespace test
} // namespace nblib
#include <gtest/gtest.h>
#include "gromacs/topology/exclusionblocks.h"
-#include "nblib/gmxcalculator.h"
-#include "nblib/gmxsetup.h"
+#include "nblib/gmxcalculatorcpu.h"
#include "nblib/integrator.h"
#include "nblib/tests/testhelpers.h"
#include "nblib/tests/testsystems.h"
SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder;
- auto simState = spcMethanolSystemBuilder.setupSimulationState();
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
+ auto simState = spcMethanolSystemBuilder.setupSimulationState();
+
+ auto forceCalculator = setupGmxForceCalculatorCpu(simState.topology(), options);
+
+ forceCalculator->updatePairlist(simState.coordinates(), simState.box());
gmx::ArrayRef<Vec3> forces(simState.forces());
- ASSERT_NO_THROW(forceCalculator->compute(simState.coordinates(), forces));
+ ASSERT_NO_THROW(forceCalculator->compute(simState.coordinates(), simState.box(), forces));
RefDataChecker forcesOutputTest(5e-5);
forcesOutputTest.testArrays<Vec3>(forces, "SPC-methanol forces");
SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder;
- auto simState = spcMethanolSystemBuilder.setupSimulationState();
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
+ auto simState = spcMethanolSystemBuilder.setupSimulationState();
+
+ auto forceCalculator = setupGmxForceCalculatorCpu(simState.topology(), options);
+
+ forceCalculator->updatePairlist(simState.coordinates(), simState.box());
gmx::ArrayRef<Vec3> forces(simState.forces());
- forceCalculator->compute(simState.coordinates(), forces);
+ forceCalculator->compute(simState.coordinates(), simState.box(), forces);
EXPECT_EQ(simState.topology().numParticles(), forces.size());
}
SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder;
- auto simState = spcMethanolSystemBuilder.setupSimulationState();
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
+ auto simState = spcMethanolSystemBuilder.setupSimulationState();
+
+ auto forceCalculator = setupGmxForceCalculatorCpu(simState.topology(), options);
+
+ forceCalculator->updatePairlist(simState.coordinates(), simState.box());
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()));
+ forceCalculator->compute(simState.coordinates(), simState.box(), forces);
+ EXPECT_NO_THROW(integrator.integrate(1.0, simState.coordinates(), simState.velocities(), forces));
}
}
SpcMethanolSimulationStateBuilder spcMethanolSystemBuilder;
- auto simState = spcMethanolSystemBuilder.setupSimulationState();
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
+ auto simState = spcMethanolSystemBuilder.setupSimulationState();
+
+ auto forceCalculator = setupGmxForceCalculatorCpu(simState.topology(), options);
+
+ forceCalculator->updatePairlist(simState.coordinates(), simState.box());
LeapFrog integrator(simState.topology(), simState.box());
// step 1
gmx::ArrayRef<Vec3> forces(simState.forces());
- forceCalculator->compute(simState.coordinates(), simState.forces());
+
+ forceCalculator->compute(simState.coordinates(), simState.box(), simState.forces());
// copy computed forces to another array
std::vector<Vec3> forces_1(forces.size());
zeroCartesianArray(forces);
// check if forces change without update step
- forceCalculator->compute(simState.coordinates(), forces);
+ forceCalculator->compute(simState.coordinates(), simState.box(), forces);
// check if forces change without update
for (size_t i = 0; i < forces_1.size(); i++)
zeroCartesianArray(forces);
// step 2
- forceCalculator->compute(simState.coordinates(), forces);
+ forceCalculator->compute(simState.coordinates(), simState.box(), forces);
+
std::vector<Vec3> forces_2(forces.size());
std::copy(forces.begin(), forces.end(), begin(forces_2));
ArgonSimulationStateBuilder argonSystemBuilder(fftypes::OPLSA);
- auto simState = argonSystemBuilder.setupSimulationState();
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
+ auto simState = argonSystemBuilder.setupSimulationState();
+
+ auto forceCalculator = setupGmxForceCalculatorCpu(simState.topology(), options);
+
+ forceCalculator->updatePairlist(simState.coordinates(), simState.box());
gmx::ArrayRef<Vec3> testForces(simState.forces());
- forceCalculator->compute(simState.coordinates(), simState.forces());
+ forceCalculator->compute(simState.coordinates(), simState.box(), simState.forces());
RefDataChecker forcesOutputTest(1e-7);
forcesOutputTest.testArrays<Vec3>(testForces, "Argon forces");
ArgonSimulationStateBuilder argonSystemBuilder(fftypes::GROMOS43A1);
- auto simState = argonSystemBuilder.setupSimulationState();
- auto forceCalculator = nblib::setupGmxForceCalculator(
- simState.topology(), simState.coordinates(), simState.box(), options);
+ auto simState = argonSystemBuilder.setupSimulationState();
+
+ auto forceCalculator = setupGmxForceCalculatorCpu(simState.topology(), options);
+
+ forceCalculator->updatePairlist(simState.coordinates(), simState.box());
gmx::ArrayRef<Vec3> testForces(simState.forces());
- forceCalculator->compute(simState.coordinates(), simState.forces());
+ forceCalculator->compute(simState.coordinates(), simState.box(), simState.forces());
- RefDataChecker forcesOutputTest(1e-7);
+ RefDataChecker forcesOutputTest;
forcesOutputTest.testArrays<Vec3>(testForces, "Argon forces");
}
--- /dev/null
+<?xml version="1.0"?>
+<?xml-stylesheet type="text/xsl" href="referencedata.xsl"?>
+<ReferenceData>
+ <Sequence Name="Argon energies">
+ <Int Name="Length">5</Int>
+ <Real>0</Real>
+ <Real>-1.0175924</Real>
+ <Real>0</Real>
+ <Real>0</Real>
+ <Real>0</Real>
+ </Sequence>
+</ReferenceData>
--- /dev/null
+<?xml version="1.0"?>
+<?xml-stylesheet type="text/xsl" href="referencedata.xsl"?>
+<ReferenceData>
+ <Sequence Name="Virials">
+ <Int Name="Length">9</Int>
+ <Real>0.094908442751139188</Real>
+ <Real>0.25238615516784413</Real>
+ <Real>0.02601194356883077</Real>
+ <Real>0.25238615516784413</Real>
+ <Real>0.67116022003893372</Real>
+ <Real>0.069172501786742557</Real>
+ <Real>0.026011943568830798</Real>
+ <Real>0.069172501786742391</Real>
+ <Real>0.0071291993484943611</Real>
+ </Sequence>
+</ReferenceData>
\ No newline at end of file