2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2018,2019, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
37 * Tests for TopologyInformation
39 * \author Mark Abraham <mark.j.abraham@gmail.com>
40 * \ingroup module_trajectoryanalysis
44 #include "gromacs/trajectoryanalysis/topologyinformation.h"
46 #include <gtest/gtest.h>
48 #include "gromacs/gmxpreprocess/grompp.h"
49 #include "gromacs/math/vectypes.h"
50 #include "gromacs/topology/topology.h"
51 #include "gromacs/utility/arrayref.h"
52 #include "gromacs/utility/exceptions.h"
53 #include "gromacs/utility/textwriter.h"
55 #include "testutils/cmdlinetest.h"
56 #include "testutils/testfilemanager.h"
58 #include "moduletest.h"
67 TEST(TopologyInformation, CantWorkWithoutReadingAFile)
69 TopologyInformation topInfo;
70 EXPECT_FALSE(topInfo.hasTopology());
71 EXPECT_FALSE(topInfo.hasFullTopology());
72 EXPECT_EQ(nullptr, topInfo.mtop());
73 EXPECT_EQ(nullptr, topInfo.expandedTopology());
74 auto atoms1 = topInfo.copyAtoms();
76 auto atoms2 = topInfo.copyAtoms();
78 EXPECT_NE(atoms1.get(), atoms2.get());
79 EXPECT_EQ(0, atoms1->nr);
80 EXPECT_EQ(-1, topInfo.ePBC());
81 EXPECT_THROW(topInfo.x().size(), gmx::APIError);
82 EXPECT_THROW(topInfo.v().size(), gmx::APIError);
85 EXPECT_EQ(0, box[XX][XX]);
86 EXPECT_EQ(0, box[XX][YY]);
87 EXPECT_EQ(0, box[XX][ZZ]);
88 EXPECT_EQ(0, box[YY][XX]);
89 EXPECT_EQ(0, box[YY][YY]);
90 EXPECT_EQ(0, box[YY][ZZ]);
91 EXPECT_EQ(0, box[ZZ][XX]);
92 EXPECT_EQ(0, box[ZZ][YY]);
93 EXPECT_EQ(0, box[ZZ][ZZ]);
94 EXPECT_FALSE(topInfo.name());
97 //! Common test code to reduce duplication
98 void runCommonTests(const TopologyInformation& topInfo, const int numAtoms)
100 EXPECT_TRUE(topInfo.hasTopology());
101 ASSERT_TRUE(topInfo.mtop());
102 EXPECT_EQ(numAtoms, topInfo.mtop()->natoms);
103 // TODO Dump mtop to refdata when that is possible
104 ASSERT_TRUE(topInfo.expandedTopology());
105 auto atoms1 = topInfo.copyAtoms();
107 auto atoms2 = topInfo.copyAtoms();
109 // Must be different pointer to a deep copy.
110 EXPECT_NE(atoms1.get(), atoms2.get());
111 auto atoms = topInfo.atoms();
112 // Must be a pointer to a different instance.
113 EXPECT_NE(atoms1.get(), atoms);
114 EXPECT_NE(atoms2.get(), atoms);
115 EXPECT_EQ(numAtoms, topInfo.x().size());
116 EXPECT_EQ(numAtoms, topInfo.v().size());
117 matrix box{ { -2 } };
119 EXPECT_FLOAT_EQ(5.9062, box[XX][XX]);
120 EXPECT_FLOAT_EQ(0, box[XX][YY]);
121 EXPECT_FLOAT_EQ(0, box[XX][ZZ]);
122 EXPECT_FLOAT_EQ(0, box[YY][XX]);
123 EXPECT_FLOAT_EQ(6.8451, box[YY][YY]);
124 EXPECT_FLOAT_EQ(0, box[YY][ZZ]);
125 EXPECT_FLOAT_EQ(0, box[ZZ][XX]);
126 EXPECT_FLOAT_EQ(0, box[ZZ][YY]);
127 EXPECT_FLOAT_EQ(3.0517, box[ZZ][ZZ]);
128 EXPECT_STREQ("First 10 residues from 1AKI", topInfo.name());
131 TEST(TopologyInformation, WorksWithGroFile)
133 const int numAtoms = 156;
134 TopologyInformation topInfo;
135 topInfo.fillFromInputFile(TestFileManager::getInputFilePath("lysozyme.gro"));
136 EXPECT_FALSE(topInfo.hasFullTopology());
137 runCommonTests(topInfo, numAtoms);
138 EXPECT_EQ(-1, topInfo.ePBC());
140 // Check the per-atom data
141 auto atoms = topInfo.copyAtoms();
142 ASSERT_EQ(numAtoms, atoms->nr);
143 EXPECT_TRUE(atoms->haveMass);
144 // TODO atommass.dat assumes united atom CA, which is probably not expected behaviour
145 EXPECT_FLOAT_EQ(13.019, atoms->atom[26].m);
146 EXPECT_FALSE(atoms->haveCharge);
147 EXPECT_FALSE(atoms->haveType);
148 EXPECT_EQ(0, atoms->atom[26].type);
149 EXPECT_EQ(0, atoms->atom[26].atomnumber);
150 EXPECT_EQ(1, atoms->atom[26].resind);
151 // gro files don't have the information that pdb files might
152 EXPECT_FALSE(atoms->havePdbInfo);
153 EXPECT_FALSE(atoms->pdbinfo);
154 EXPECT_EQ(10, atoms->nres);
155 ASSERT_TRUE(atoms->resinfo);
156 ASSERT_TRUE(atoms->resinfo[4].name);
157 EXPECT_STREQ("ARG", *atoms->resinfo[4].name);
158 EXPECT_EQ(5, atoms->resinfo[4].nr);
159 EXPECT_EQ(0, atoms->resinfo[4].chainnum);
160 EXPECT_EQ(' ', atoms->resinfo[4].chainid);
163 TEST(TopologyInformation, WorksWithPdbFile)
165 const int numAtoms = 156;
166 TopologyInformation topInfo;
167 topInfo.fillFromInputFile(TestFileManager::getInputFilePath("lysozyme.pdb"));
168 EXPECT_FALSE(topInfo.hasFullTopology());
169 runCommonTests(topInfo, numAtoms);
170 // TODO why does this differ from .gro?
171 EXPECT_EQ(0, topInfo.ePBC());
173 // Check the per-atom data
174 auto atoms = topInfo.copyAtoms();
175 ASSERT_EQ(numAtoms, atoms->nr);
176 EXPECT_TRUE(atoms->haveMass);
177 // TODO atommass.dat assumes united atom CA, which is probably not expected behaviour
178 EXPECT_FLOAT_EQ(13.019, atoms->atom[26].m);
179 EXPECT_FALSE(atoms->haveCharge);
180 EXPECT_FALSE(atoms->haveType);
181 EXPECT_EQ(0, atoms->atom[26].type);
182 EXPECT_EQ(0, atoms->atom[26].atomnumber);
183 EXPECT_EQ(1, atoms->atom[26].resind);
184 // pdb files can carry more information than gro
185 EXPECT_TRUE(atoms->havePdbInfo);
186 ASSERT_TRUE(atoms->pdbinfo);
187 EXPECT_EQ(10, atoms->nres);
188 ASSERT_TRUE(atoms->resinfo);
189 ASSERT_TRUE(atoms->resinfo[4].name);
190 EXPECT_STREQ("ARG", *atoms->resinfo[4].name);
191 EXPECT_EQ(5, atoms->resinfo[4].nr);
192 EXPECT_EQ(0, atoms->resinfo[4].chainnum);
193 EXPECT_EQ('B', atoms->resinfo[4].chainid);
196 TEST(TopologyInformation, WorksWithTprFromPdbFile)
198 TestFileManager fileManager;
200 // Make the tpr file to use
201 std::string name = "lysozyme";
202 const std::string mdpInputFileName = fileManager.getTemporaryFilePath(name + ".mdp");
203 // Ensure the seeds have a value so that the resulting .tpr dump
205 TextWriter::writeFileFromString(mdpInputFileName, "");
206 std::string tprName = fileManager.getTemporaryFilePath(name + ".tpr");
209 caller.append("grompp");
210 caller.addOption("-f", mdpInputFileName);
211 caller.addOption("-p", TestFileManager::getInputFilePath(name));
212 caller.addOption("-c", TestFileManager::getInputFilePath(name + ".pdb"));
213 caller.addOption("-o", tprName);
214 ASSERT_EQ(0, gmx_grompp(caller.argc(), caller.argv()));
217 const int numAtoms = 156;
218 TopologyInformation topInfo;
219 topInfo.fillFromInputFile(tprName);
220 EXPECT_TRUE(topInfo.hasFullTopology());
221 runCommonTests(topInfo, numAtoms);
222 // TODO why does this differ from .gro?
223 EXPECT_EQ(0, topInfo.ePBC());
225 // Check the per-atom data
226 auto atoms = topInfo.copyAtoms();
227 ASSERT_EQ(numAtoms, atoms->nr);
228 EXPECT_TRUE(atoms->haveMass);
229 EXPECT_FLOAT_EQ(12.011, atoms->atom[26].m);
230 EXPECT_TRUE(atoms->haveCharge);
231 EXPECT_TRUE(atoms->haveType);
232 EXPECT_EQ(2, atoms->atom[26].type);
233 EXPECT_EQ(6, atoms->atom[26].atomnumber);
234 EXPECT_EQ(1, atoms->atom[26].resind);
235 // tpr files also don't carry pdb information
236 EXPECT_FALSE(atoms->havePdbInfo);
237 EXPECT_FALSE(atoms->pdbinfo);
238 EXPECT_EQ(10, atoms->nres);
239 ASSERT_TRUE(atoms->resinfo);
240 ASSERT_TRUE(atoms->resinfo[4].name);
241 EXPECT_STREQ("ARG", *atoms->resinfo[4].name);
242 EXPECT_EQ(5, atoms->resinfo[4].nr);
243 EXPECT_EQ(0, atoms->resinfo[4].chainnum);
244 // In particular, chain ID does not get recorded in the .tpr file
245 EXPECT_EQ(0, atoms->resinfo[4].chainid);