+/*
+ * 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 tpr reading
+ *
+ * \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/tpr.h"
+
+#include "gromacs/fileio/tpxio.h"
+#include "gromacs/gmxlib/network.h"
+#include "gromacs/listed_forces/listed_forces.h"
+#include "gromacs/mdtypes/forcerec.h"
+#include "gromacs/mdlib/mdatoms.h"
+#include "gromacs/mdlib/forcerec.h"
+#include "gromacs/mdtypes/commrec.h"
+#include "gromacs/mdtypes/iforceprovider.h"
+#include "gromacs/mdtypes/inputrec.h"
+#include "gromacs/mdtypes/mdatom.h"
+#include "gromacs/mdtypes/simulation_workload.h"
+#include "gromacs/mdtypes/state.h"
+#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/topology.h"
+#include "gromacs/utility/logger.h"
+#include "nblib/box.h"
+
+namespace nblib
+{
+
+TprReader::TprReader(std::string filename)
+{
+ t_inputrec inputRecord;
+ t_state globalState;
+ gmx_mtop_t molecularTopology;
+ gmx::SimulationWorkload simulationWorkload;
+
+ // If the file does not exist, this function will throw
+ PartialDeserializedTprFile partialDeserializedTpr =
+ read_tpx_state(filename.c_str(), &inputRecord, &globalState, &molecularTopology);
+
+ // init commrec
+ MPI_Comm simulationCommunicator = MPI_COMM_WORLD;
+ CommrecHandle crHandle = init_commrec(simulationCommunicator);
+ t_commrec* commrec = crHandle.get();
+ assert((commrec != nullptr) && "Must have valid commrec");
+
+ // init forcerec
+ t_forcerec forceRecord;
+ gmx::ForceProviders forceProviders;
+ forceRecord.forceProviders = &forceProviders;
+ init_forcerec(nullptr,
+ gmx::MDLogger(),
+ simulationWorkload,
+ &forceRecord,
+ inputRecord,
+ molecularTopology,
+ commrec,
+ globalState.box,
+ nullptr,
+ nullptr,
+ gmx::ArrayRef<const std::string>{},
+ -1);
+
+ nonbondedParameters_ = makeNonBondedParameterLists(molecularTopology.ffparams.atnr,
+ molecularTopology.ffparams.iparams,
+ forceRecord.haveBuckingham);
+
+ gmx_localtop_t localtop(molecularTopology.ffparams);
+ gmx_mtop_generate_local_top(molecularTopology, &localtop, false);
+ exclusionListElements_ = std::vector<int>(localtop.excls.elementsView().begin(),
+ localtop.excls.elementsView().end());
+ exclusionListRanges_ = std::vector<int>(localtop.excls.listRangesView().begin(),
+ localtop.excls.listRangesView().end());
+
+ int ntopatoms = molecularTopology.natoms;
+ std::unique_ptr<gmx::MDAtoms> mdAtoms =
+ gmx::makeMDAtoms(nullptr, molecularTopology, inputRecord, false);
+ atoms2md(molecularTopology, inputRecord, -1, {}, ntopatoms, mdAtoms.get());
+ update_mdatoms(mdAtoms->mdatoms(), inputRecord.fepvals->init_lambda);
+ auto numParticles = mdAtoms->mdatoms()->nr;
+ charges_.resize(numParticles);
+ particleTypeIdOfAllParticles_.resize(numParticles);
+ inverseMasses_.resize(numParticles);
+ for (int i = 0; i < numParticles; i++)
+ {
+ charges_[i] = mdAtoms->mdatoms()->chargeA[i];
+ particleTypeIdOfAllParticles_[i] = mdAtoms->mdatoms()->typeA[i];
+ inverseMasses_[i] = mdAtoms->mdatoms()->invmass[i];
+ }
+ particleInteractionFlags_ = forceRecord.atomInfo;
+
+ if (TRICLINIC(globalState.box))
+ {
+ throw InputException("Only rectangular unit-cells are supported here");
+ }
+ boxX_ = globalState.box[XX][XX];
+ boxY_ = globalState.box[YY][YY];
+ boxZ_ = globalState.box[ZZ][ZZ];
+ coordinates_.assign(globalState.x.begin(), globalState.x.end());
+ velocities_.assign(globalState.v.begin(), globalState.v.end());
+}
+
+Box TprReader::getBox() const
+{
+ return Box(boxX_, boxY_, boxZ_);
+}
+
+} // namespace nblib