2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2013, 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 the index group handling in the selection engine.
40 * Tests for other functions, at least the set operations.
42 * \author Teemu Murtola <teemu.murtola@gmail.com>
43 * \ingroup module_selection
45 #include <gtest/gtest.h>
47 #include "gromacs/legacyheaders/typedefs.h"
49 #include "gromacs/selection/indexutil.h"
51 #include "testutils/refdata.h"
58 /********************************************************************
62 class IndexBlockTest : public ::testing::Test
68 void setGroup(int count, const int atoms[]);
70 void setGroup(const int (&atoms)[count])
72 setGroup(count, atoms);
77 gmx::test::TestReferenceData data_;
78 gmx::test::TopologyManager topManager_;
83 IndexBlockTest::IndexBlockTest()
87 blocka_.nalloc_index = 0;
91 gmx_ana_index_clear(&g_);
94 IndexBlockTest::~IndexBlockTest()
96 done_blocka(&blocka_);
99 void IndexBlockTest::setGroup(int count, const int atoms[])
102 g_.index = const_cast<int *>(atoms);
105 void IndexBlockTest::checkBlocka()
107 gmx::test::TestReferenceChecker compound(
108 data_.rootChecker().checkCompound("BlockAtoms", "Block"));
109 compound.checkSequenceArray(g_.isize, g_.index, "Input");
110 compound.checkInteger(blocka_.nr, "Count");
111 for (int i = 0; i < blocka_.nr; ++i)
113 gmx::test::TestReferenceChecker blockCompound(
114 compound.checkCompound("Block", NULL));
115 blockCompound.checkSequence(&blocka_.a[blocka_.index[i]],
116 &blocka_.a[blocka_.index[i+1]],
121 /********************************************************************
122 * gmx_ana_index_make_block() tests
125 TEST_F(IndexBlockTest, CreatesUnknownBlock)
127 gmx_ana_index_make_block(&blocka_, NULL, NULL, INDEX_UNKNOWN, false);
129 done_blocka(&blocka_);
130 gmx_ana_index_make_block(&blocka_, NULL, NULL, INDEX_UNKNOWN, false);
134 TEST_F(IndexBlockTest, CreatesAtomBlock)
136 const int group[] = { 0, 1, 3, 4, 6 };
138 gmx_ana_index_make_block(&blocka_, NULL, &g_, INDEX_ATOM, false);
140 done_blocka(&blocka_);
141 gmx_ana_index_make_block(&blocka_, NULL, &g_, INDEX_ATOM, true);
145 TEST_F(IndexBlockTest, CreatesResidueBlock)
147 const int group[] = { 0, 1, 3, 6, 7 };
148 topManager_.initAtoms(9);
149 topManager_.initUniformResidues(3);
151 gmx_ana_index_make_block(&blocka_, topManager_.topology(), &g_,
156 TEST_F(IndexBlockTest, CreatesMoleculeBlock)
158 const int group[] = { 3, 4, 7, 8, 13 };
159 topManager_.initAtoms(18);
160 topManager_.initUniformMolecules(3);
162 gmx_ana_index_make_block(&blocka_, topManager_.topology(), &g_,
167 TEST_F(IndexBlockTest, CreatesResidueBlockWithCompletion)
169 const int group[] = { 3, 4, 7, 8, 13 };
170 topManager_.initAtoms(18);
171 topManager_.initUniformResidues(3);
173 gmx_ana_index_make_block(&blocka_, topManager_.topology(), &g_,
178 TEST_F(IndexBlockTest, CreatesMoleculeBlockWithCompletion)
180 const int group[] = { 3, 4, 7, 8, 13 };
181 topManager_.initAtoms(18);
182 topManager_.initUniformMolecules(3);
184 gmx_ana_index_make_block(&blocka_, topManager_.topology(), &g_,
189 TEST_F(IndexBlockTest, CreatesSingleBlock)
191 const int group[] = { 0, 1, 3, 4, 6 };
193 gmx_ana_index_make_block(&blocka_, NULL, &g_, INDEX_ALL, false);
195 done_blocka(&blocka_);
196 gmx_ana_index_make_block(&blocka_, NULL, &g_, INDEX_ALL, true);
200 /********************************************************************
201 * gmx_ana_index_has_full_ablocks() tests
204 TEST_F(IndexBlockTest, ChecksGroupForFullBlocksPositive)
206 const int maxGroup[] = {
207 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17
209 const int testGroup[] = { 3, 4, 5, 6, 7, 8, 12, 13, 14 };
210 topManager_.initAtoms(18);
211 topManager_.initUniformResidues(3);
213 gmx_ana_index_make_block(&blocka_, topManager_.topology(), &g_,
216 EXPECT_TRUE(gmx_ana_index_has_full_ablocks(&g_, &blocka_));
219 TEST_F(IndexBlockTest, ChecksOutOfOrderGroupForFullBlocksPositive)
221 const int maxGroup[] = {
222 15, 16, 17, 2, 1, 0, 12, 13, 14, 5, 4, 3, 9, 10, 11, 8, 7, 6
224 const int testGroup[] = { 2, 1, 0, 5, 4, 3, 8, 7, 6, };
225 topManager_.initAtoms(18);
226 topManager_.initUniformResidues(3);
228 gmx_ana_index_make_block(&blocka_, topManager_.topology(), &g_,
231 EXPECT_TRUE(gmx_ana_index_has_full_ablocks(&g_, &blocka_));
234 TEST_F(IndexBlockTest, ChecksGroupForFullBlocksNegative)
236 const int maxGroup[] = {
237 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17
239 const int testGroup1[] = { 3, 4, 5, 6, 7, 8, 12, 13 };
240 const int testGroup2[] = { 3, 4, 5, 6, 7, 12, 13, 14 };
241 const int testGroup3[] = { 4, 5, 6, 7, 8, 12, 13, 14 };
243 topManager_.initAtoms(18);
244 topManager_.initUniformResidues(3);
246 gmx_ana_index_make_block(&blocka_, topManager_.topology(), &g_,
249 setGroup(testGroup1);
250 EXPECT_FALSE(gmx_ana_index_has_full_ablocks(&g_, &blocka_));
252 setGroup(testGroup2);
253 EXPECT_FALSE(gmx_ana_index_has_full_ablocks(&g_, &blocka_));
255 setGroup(testGroup3);
256 EXPECT_FALSE(gmx_ana_index_has_full_ablocks(&g_, &blocka_));
259 /********************************************************************
260 * gmx_ana_index_has_complete_elems() tests
263 TEST_F(IndexBlockTest, ChecksGroupForCompleteElementsTrivial)
265 const int group[] = { 0, 1, 2 };
267 EXPECT_TRUE(gmx_ana_index_has_complete_elems(&g_, INDEX_ATOM, NULL));
268 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_ALL, NULL));
269 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_UNKNOWN, NULL));
272 TEST_F(IndexBlockTest, ChecksGroupForCompleteResiduesPositive)
274 const int group1[] = { 0, 1, 2, 6, 7, 8, 12, 13, 14 };
275 const int group2[] = { 3, 4, 5, 6, 7, 8 };
277 topManager_.initAtoms(15);
278 topManager_.initUniformResidues(3);
279 t_topology *top = topManager_.topology();
282 EXPECT_TRUE(gmx_ana_index_has_complete_elems(&g_, INDEX_RES, top));
285 EXPECT_TRUE(gmx_ana_index_has_complete_elems(&g_, INDEX_RES, top));
288 TEST_F(IndexBlockTest, ChecksGroupForCompleteResiduesNegative)
290 const int group1[] = { 3, 4, 5, 6, 7, 8, 12, 13 };
291 const int group2[] = { 3, 4, 5, 6, 7, 12, 13, 14 };
292 const int group3[] = { 4, 5, 6, 7, 8, 12, 13, 14 };
294 topManager_.initAtoms(18);
295 topManager_.initUniformResidues(3);
296 t_topology *top = topManager_.topology();
299 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_RES, top));
302 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_RES, top));
305 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_RES, top));
308 TEST_F(IndexBlockTest, ChecksGroupForCompleteMoleculesPositive)
310 const int group[] = { 0, 1, 2, 6, 7, 8, 12, 13, 14 };
312 topManager_.initAtoms(15);
313 topManager_.initUniformMolecules(3);
314 t_topology *top = topManager_.topology();
317 EXPECT_TRUE(gmx_ana_index_has_complete_elems(&g_, INDEX_MOL, top));
320 TEST_F(IndexBlockTest, ChecksGroupForCompleteMoleculesNegative)
322 const int group1[] = { 3, 4, 5, 6, 7, 8, 12, 13 };
323 const int group2[] = { 3, 4, 5, 6, 7, 12, 13, 14 };
324 const int group3[] = { 4, 5, 6, 7, 8, 12, 13, 14 };
326 topManager_.initAtoms(18);
327 topManager_.initUniformMolecules(3);
328 t_topology *top = topManager_.topology();
331 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_MOL, top));
334 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_MOL, top));
337 EXPECT_FALSE(gmx_ana_index_has_complete_elems(&g_, INDEX_MOL, top));
340 /********************************************************************
344 class IndexMapTest : public ::testing::Test
350 void testInit(int atomCount, const int atoms[], e_index_t type);
351 void testUpdate(int atomCount, const int atoms[], bool bMaskOnly,
354 void testInit(const int (&atoms)[count], e_index_t type)
356 testInit(count, atoms, type);
359 void testUpdate(const int (&atoms)[count], bool bMaskOnly,
362 testUpdate(count, atoms, bMaskOnly, name);
365 void checkMapping(int atomCount, const int atoms[], const char *name);
367 gmx::test::TestReferenceData data_;
368 gmx::test::TestReferenceChecker checker_;
369 gmx::test::TopologyManager topManager_;
370 gmx_ana_indexmap_t map_;
373 gmx_ana_index_t initGroup_;
376 IndexMapTest::IndexMapTest()
377 : checker_(data_.rootChecker())
379 gmx_ana_indexmap_clear(&map_);
380 gmx_ana_index_clear(&initGroup_);
383 IndexMapTest::~IndexMapTest()
385 gmx_ana_indexmap_deinit(&map_);
388 void IndexMapTest::testInit(int atomCount, const int atoms[], e_index_t type)
390 initGroup_.isize = atomCount;
391 initGroup_.index = const_cast<int *>(atoms);
392 gmx_ana_indexmap_init(&map_, &initGroup_, topManager_.topology(), type);
393 EXPECT_EQ(type, map_.type);
394 checkMapping(atomCount, atoms, "Initialized");
397 void IndexMapTest::testUpdate(int atomCount, const int atoms[], bool bMaskOnly,
402 g.index = const_cast<int *>(atoms);
403 gmx_ana_indexmap_update(&map_, &g, bMaskOnly);
410 checkMapping(initGroup_.isize, initGroup_.index, name);
414 checkMapping(atomCount, atoms, name);
418 void IndexMapTest::checkMapping(int atomCount, const int atoms[],
421 gmx::test::TestReferenceChecker compound(
422 checker_.checkCompound("IndexMapping", name));
423 compound.checkSequenceArray(atomCount, atoms, "Input");
424 compound.checkInteger(map_.mapb.nr, "Count");
425 for (int i = 0; i < map_.mapb.nr; ++i)
427 gmx::test::TestReferenceChecker blockCompound(
428 compound.checkCompound("Block", NULL));
429 blockCompound.checkSequence(&atoms[map_.mapb.index[i]],
430 &atoms[map_.mapb.index[i+1]],
432 blockCompound.checkInteger(map_.refid[i], "RefId");
433 blockCompound.checkInteger(map_.mapid[i], "MapId");
434 int originalIdIndex = (map_.refid[i] != -1 ? map_.refid[i] : i);
435 EXPECT_EQ(map_.orgid[originalIdIndex], map_.mapid[i]);
439 /********************************************************************
440 * gmx_ana_indexmap_t tests
443 TEST_F(IndexMapTest, InitializesAtomBlock)
445 const int maxGroup[] = { 1, 2, 4, 5 };
446 testInit(maxGroup, INDEX_ATOM);
449 TEST_F(IndexMapTest, InitializesMoleculeBlock)
451 const int maxGroup[] = { 3, 4, 7, 8, 13 };
452 topManager_.initAtoms(18);
453 topManager_.initUniformMolecules(3);
454 testInit(maxGroup, INDEX_MOL);
457 TEST_F(IndexMapTest, MapsSingleBlock)
459 const int maxGroup[] = { 0, 1, 2, 3 };
460 const int evalGroup[] = { 0, 2 };
461 testInit(maxGroup, INDEX_ALL);
462 testUpdate(evalGroup, false, NULL);
465 TEST_F(IndexMapTest, MapsResidueBlocks)
467 const int maxGroup[] = {
468 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17
470 const int evalGroup[] = { 3, 4, 7, 8, 13 };
471 topManager_.initAtoms(18);
472 topManager_.initUniformResidues(3);
473 testInit(maxGroup, INDEX_RES);
474 testUpdate(evalGroup, false, NULL);
477 TEST_F(IndexMapTest, MapsResidueBlocksWithMask)
479 const int maxGroup[] = {
480 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17
482 const int evalGroup[] = { 3, 4, 7, 8, 13 };
483 topManager_.initAtoms(18);
484 topManager_.initUniformResidues(3);
485 testInit(maxGroup, INDEX_RES);
486 testUpdate(evalGroup, true, NULL);
489 TEST_F(IndexMapTest, HandlesMultipleRequests)
491 const int maxGroup[] = {
492 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17
494 const int evalGroup[] = { 3, 4, 7, 8, 13 };
495 topManager_.initAtoms(18);
496 topManager_.initUniformResidues(3);
497 testInit(maxGroup, INDEX_RES);
498 testUpdate(evalGroup, false, "EvaluatedNoMask");
499 testUpdate(evalGroup, true, "EvaluatedMask");
500 testUpdate(maxGroup, true, "Initialized");
501 testUpdate(evalGroup, true, "EvaluatedMask");
502 testUpdate(evalGroup, false, "EvaluatedNoMask");
503 testUpdate(maxGroup, false, "Initialized");