2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2009,2010,2011,2012,2013,2014, 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 * Implements neighborhood searching for analysis (from nbsearch.h).
40 * The grid implementation could still be optimized in several different ways:
41 * - Triclinic grid cells are not the most efficient shape, but make PBC
43 * - Precalculating the required PBC shift for a pair of cells outside the
44 * inner loop. After this is done, it should be quite straightforward to
45 * move to rectangular cells.
46 * - Pruning grid cells from the search list if they are completely outside
47 * the sphere that is being considered.
48 * - A better heuristic could be added for falling back to simple loops for a
49 * small number of reference particles.
50 * - A better heuristic for selecting the grid size.
51 * - A multi-level grid implementation could be used to be able to use small
52 * grids for short cutoffs with very inhomogeneous particle distributions
53 * without a memory cost.
55 * \author Teemu Murtola <teemu.murtola@gmail.com>
56 * \ingroup module_selection
58 #include "gromacs/selection/nbsearch.h"
65 #include "thread_mpi/mutex.h"
67 #include "gromacs/legacyheaders/typedefs.h"
69 #include "gromacs/math/vec.h"
70 #include "gromacs/pbcutil/pbc.h"
71 #include "gromacs/selection/position.h"
72 #include "gromacs/utility/arrayref.h"
73 #include "gromacs/utility/gmxassert.h"
74 #include "gromacs/utility/smalloc.h"
82 /********************************************************************
83 * Implementation class declarations
86 class AnalysisNeighborhoodSearchImpl
89 typedef AnalysisNeighborhoodPairSearch::ImplPointer
90 PairSearchImplPointer;
91 typedef std::vector<PairSearchImplPointer> PairSearchList;
92 typedef std::vector<std::vector<int> > CellList;
94 explicit AnalysisNeighborhoodSearchImpl(real cutoff);
95 ~AnalysisNeighborhoodSearchImpl();
98 * Initializes the search with a given box and reference positions.
100 * \param[in] mode Search mode to use.
101 * \param[in] pbc PBC information.
102 * \param[in] positions Set of reference positions.
104 void init(AnalysisNeighborhood::SearchMode mode,
106 const AnalysisNeighborhoodPositions &positions);
107 PairSearchImplPointer getPairSearch();
109 real cutoffSquared() const { return cutoff2_; }
110 bool usesGridSearch() const { return bGrid_; }
113 //! Calculates offsets to neighboring grid cells that should be considered.
114 void initGridCellNeighborList();
116 * Determines a suitable grid size and sets up the cells.
118 * \param[in] pbc Information about the box.
119 * \returns false if grid search is not suitable.
121 bool initGridCells(const t_pbc *pbc);
123 * Sets ua a search grid for a given box.
125 * \param[in] pbc Information about the box.
126 * \returns false if grid search is not suitable.
128 bool initGrid(const t_pbc *pbc);
130 * Maps a point into a grid cell.
132 * \param[in] x Point to map.
133 * \param[out] cell Indices of the grid cell in which \p x lies.
135 * \p x should be within the triclinic unit cell.
137 void mapPointToGridCell(const rvec x, ivec cell) const;
139 * Calculates linear index of a grid cell.
141 * \param[in] cell Cell indices.
142 * \returns Linear index of \p cell.
144 int getGridCellIndex(const ivec cell) const;
146 * Adds an index into a grid cell.
148 * \param[in] cell Cell into which \p i should be added.
149 * \param[in] i Index to add.
151 void addToGridCell(const ivec cell, int i);
153 //! Whether to try grid searching.
157 //! The cutoff squared.
160 //! Number of reference points for the current frame.
162 //! Reference point positions.
164 //! Reference position ids (NULL if not available).
169 //! Number of excluded reference positions for current test particle.
171 //! Exclusions for current test particle.
174 //! Whether grid searching is actually used for the current positions.
176 //! Array allocated for storing in-unit-cell reference positions.
178 //! Allocation count for xref_alloc.
180 //! false if the box is rectangular.
182 //! Box vectors of a single grid cell.
184 //! The reciprocal cell vectors as columns; the inverse of \p cellbox.
186 //! Number of cells along each dimension.
188 //! Data structure to hold the grid cell contents.
190 //! Number of neighboring cells to consider.
192 //! Offsets of the neighboring cells to consider.
194 //! Allocation count for \p gnboffs.
197 tMPI::mutex createPairSearchMutex_;
198 PairSearchList pairSearchList_;
200 friend class AnalysisNeighborhoodPairSearchImpl;
202 GMX_DISALLOW_COPY_AND_ASSIGN(AnalysisNeighborhoodSearchImpl);
205 class AnalysisNeighborhoodPairSearchImpl
208 explicit AnalysisNeighborhoodPairSearchImpl(const AnalysisNeighborhoodSearchImpl &search)
212 clear_ivec(testcell_);
216 //! Initializes a search to find reference positions neighboring \p x.
217 void startSearch(const AnalysisNeighborhoodPositions &positions);
218 //! Searches for the next neighbor.
219 template <class Action>
220 bool searchNext(Action action);
221 //! Initializes a pair representing the pair found by searchNext().
222 void initFoundPair(AnalysisNeighborhoodPair *pair) const;
223 //! Advances to the next test position, skipping any remaining pairs.
224 void nextTestPosition();
227 //! Clears the loop indices.
228 void reset(int testIndex);
229 //! Checks whether a reference positiong should be excluded.
230 bool isExcluded(int j);
232 //! Parent search object.
233 const AnalysisNeighborhoodSearchImpl &search_;
234 //! Reference to the test positions.
235 ConstArrayRef<rvec> testPositions_;
236 //! Index of the currently active test position in \p testPositions_.
238 //! Stores test position during a pair loop.
240 //! Stores the previous returned position during a pair loop.
242 //! Stores the current exclusion index during loops.
244 //! Stores the test particle cell index during loops.
246 //! Stores the current cell neighbor index during pair loops.
248 //! Stores the index within the current cell during pair loops.
251 GMX_DISALLOW_COPY_AND_ASSIGN(AnalysisNeighborhoodPairSearchImpl);
254 /********************************************************************
255 * AnalysisNeighborhoodSearchImpl
258 AnalysisNeighborhoodSearchImpl::AnalysisNeighborhoodSearchImpl(real cutoff)
264 cutoff_ = GMX_REAL_MAX;
267 cutoff2_ = sqr(cutoff_);
283 clear_mat(recipcell_);
284 clear_ivec(ncelldim_);
291 AnalysisNeighborhoodSearchImpl::~AnalysisNeighborhoodSearchImpl()
293 PairSearchList::const_iterator i;
294 for (i = pairSearchList_.begin(); i != pairSearchList_.end(); ++i)
296 GMX_RELEASE_ASSERT(i->unique(),
297 "Dangling AnalysisNeighborhoodPairSearch reference");
303 AnalysisNeighborhoodSearchImpl::PairSearchImplPointer
304 AnalysisNeighborhoodSearchImpl::getPairSearch()
306 tMPI::lock_guard<tMPI::mutex> lock(createPairSearchMutex_);
307 // TODO: Consider whether this needs to/can be faster, e.g., by keeping a
308 // separate pool of unused search objects.
309 PairSearchList::const_iterator i;
310 for (i = pairSearchList_.begin(); i != pairSearchList_.end(); ++i)
317 PairSearchImplPointer pairSearch(new AnalysisNeighborhoodPairSearchImpl(*this));
318 pairSearchList_.push_back(pairSearch);
322 void AnalysisNeighborhoodSearchImpl::initGridCellNeighborList()
324 int maxx, maxy, maxz;
327 /* Find the extent of the sphere in triclinic coordinates */
328 maxz = static_cast<int>(cutoff_ * recipcell_[ZZ][ZZ]) + 1;
329 rvnorm = sqrt(sqr(recipcell_[YY][YY]) + sqr(recipcell_[ZZ][YY]));
330 maxy = static_cast<int>(cutoff_ * rvnorm) + 1;
331 rvnorm = sqrt(sqr(recipcell_[XX][XX]) + sqr(recipcell_[YY][XX])
332 + sqr(recipcell_[ZZ][XX]));
333 maxx = static_cast<int>(cutoff_ * rvnorm) + 1;
335 /* Calculate the number of cells and reallocate if necessary */
336 ngridnb_ = (2 * maxx + 1) * (2 * maxy + 1) * (2 * maxz + 1);
337 if (gnboffs_nalloc_ < ngridnb_)
339 gnboffs_nalloc_ = ngridnb_;
340 srenew(gnboffs_, gnboffs_nalloc_);
343 /* Store the whole cube */
344 /* TODO: Prune off corners that are not needed */
346 for (int x = -maxx; x <= maxx; ++x)
348 for (int y = -maxy; y <= maxy; ++y)
350 for (int z = -maxz; z <= maxz; ++z)
361 bool AnalysisNeighborhoodSearchImpl::initGridCells(const t_pbc *pbc)
363 const real targetsize =
364 pow(pbc->box[XX][XX] * pbc->box[YY][YY] * pbc->box[ZZ][ZZ]
365 * 10 / nref_, static_cast<real>(1./3.));
368 for (int dd = 0; dd < DIM; ++dd)
370 ncelldim_[dd] = static_cast<int>(pbc->box[dd][dd] / targetsize);
371 cellCount *= ncelldim_[dd];
372 if (ncelldim_[dd] < 3)
377 // Never decrease the size of the cell vector to avoid reallocating
378 // memory for the nested vectors. The actual size of the vector is not
379 // used outside this function.
380 if (cells_.size() < static_cast<size_t>(cellCount))
382 cells_.resize(cellCount);
384 for (int ci = 0; ci < cellCount; ++ci)
391 bool AnalysisNeighborhoodSearchImpl::initGrid(const t_pbc *pbc)
393 /* TODO: This check could be improved. */
394 if (0.5*pbc->max_cutoff2 < cutoff2_)
399 if (!initGridCells(pbc))
404 bTric_ = TRICLINIC(pbc->box);
407 for (int dd = 0; dd < DIM; ++dd)
409 svmul(1.0 / ncelldim_[dd], pbc->box[dd], cellbox_[dd]);
411 m_inv_ur0(cellbox_, recipcell_);
415 for (int dd = 0; dd < DIM; ++dd)
417 cellbox_[dd][dd] = pbc->box[dd][dd] / ncelldim_[dd];
418 recipcell_[dd][dd] = 1.0 / cellbox_[dd][dd];
421 initGridCellNeighborList();
425 void AnalysisNeighborhoodSearchImpl::mapPointToGridCell(const rvec x,
432 tmvmul_ur0(recipcell_, x, tx);
433 for (int dd = 0; dd < DIM; ++dd)
435 cell[dd] = static_cast<int>(tx[dd]);
440 for (int dd = 0; dd < DIM; ++dd)
442 cell[dd] = static_cast<int>(x[dd] * recipcell_[dd][dd]);
447 int AnalysisNeighborhoodSearchImpl::getGridCellIndex(const ivec cell) const
449 GMX_ASSERT(cell[XX] >= 0 && cell[XX] < ncelldim_[XX],
450 "Grid cell X index out of range");
451 GMX_ASSERT(cell[YY] >= 0 && cell[YY] < ncelldim_[YY],
452 "Grid cell Y index out of range");
453 GMX_ASSERT(cell[ZZ] >= 0 && cell[ZZ] < ncelldim_[ZZ],
454 "Grid cell Z index out of range");
455 return cell[XX] + cell[YY] * ncelldim_[XX]
456 + cell[ZZ] * ncelldim_[XX] * ncelldim_[YY];
459 void AnalysisNeighborhoodSearchImpl::addToGridCell(const ivec cell, int i)
461 const int ci = getGridCellIndex(cell);
462 cells_[ci].push_back(i);
465 void AnalysisNeighborhoodSearchImpl::init(
466 AnalysisNeighborhood::SearchMode mode,
468 const AnalysisNeighborhoodPositions &positions)
470 GMX_RELEASE_ASSERT(positions.index_ == -1,
471 "Individual indexed positions not supported as reference");
472 pbc_ = const_cast<t_pbc *>(pbc);
473 nref_ = positions.count_;
474 // TODO: Consider whether it would be possible to support grid searching in
476 if (mode == AnalysisNeighborhood::eSearchMode_Simple
477 || pbc_ == NULL || pbc_->ePBC != epbcXYZ)
483 // TODO: Actually implement forcing eSearchMode_Grid
484 bGrid_ = initGrid(pbc_);
488 if (xref_nalloc_ < nref_)
490 srenew(xref_alloc_, nref_);
491 xref_nalloc_ = nref_;
495 for (int i = 0; i < nref_; ++i)
497 copy_rvec(positions.x_[i], xref_alloc_[i]);
499 put_atoms_in_triclinic_unitcell(ecenterTRIC, pbc_->box,
501 for (int i = 0; i < nref_; ++i)
505 mapPointToGridCell(xref_[i], refcell);
506 addToGridCell(refcell, i);
511 xref_ = positions.x_;
513 // TODO: Once exclusions are supported, this may need to be initialized.
519 * Sets the exclusions for the next neighborhood search.
521 * \param[in,out] d Neighborhood search data structure.
522 * \param[in] nexcl Number of reference positions to exclude from next
524 * \param[in] excl Indices of reference positions to exclude.
526 * The set exclusions remain in effect until the next call of this function.
529 gmx_ana_nbsearch_set_excl(gmx_ana_nbsearch_t *d, int nexcl, int excl[])
537 /********************************************************************
538 * AnalysisNeighborhoodPairSearchImpl
541 void AnalysisNeighborhoodPairSearchImpl::reset(int testIndex)
543 testIndex_ = testIndex;
544 if (testIndex_ >= 0 && testIndex_ < static_cast<int>(testPositions_.size()))
546 copy_rvec(testPositions_[testIndex_], xtest_);
549 put_atoms_in_triclinic_unitcell(ecenterTRIC, search_.pbc_->box,
551 search_.mapPointToGridCell(xtest_, testcell_);
560 void AnalysisNeighborhoodPairSearchImpl::nextTestPosition()
562 if (testIndex_ < static_cast<int>(testPositions_.size()))
569 bool AnalysisNeighborhoodPairSearchImpl::isExcluded(int j)
571 if (exclind_ < search_.nexcl_)
575 while (exclind_ < search_.nexcl_
576 && search_.excl_[exclind_] < search_.refid_[j])
580 if (exclind_ < search_.nexcl_
581 && search_.refid_[j] == search_.excl_[exclind_])
589 while (search_.bGrid_ && exclind_ < search_.nexcl_
590 && search_.excl_[exclind_] < j)
594 if (search_.excl_[exclind_] == j)
604 void AnalysisNeighborhoodPairSearchImpl::startSearch(
605 const AnalysisNeighborhoodPositions &positions)
607 if (positions.index_ < 0)
609 testPositions_ = ConstArrayRef<rvec>(positions.x_, positions.count_);
614 // Somewhat of a hack: setup the array such that only the last position
616 testPositions_ = ConstArrayRef<rvec>(positions.x_, positions.index_ + 1);
617 reset(positions.index_);
621 template <class Action>
622 bool AnalysisNeighborhoodPairSearchImpl::searchNext(Action action)
624 while (testIndex_ < static_cast<int>(testPositions_.size()))
629 int cai = prevcai_ + 1;
631 for (; nbi < search_.ngridnb_; ++nbi)
635 ivec_add(testcell_, search_.gnboffs_[nbi], cell);
636 cell[XX] = (cell[XX] + search_.ncelldim_[XX]) % search_.ncelldim_[XX];
637 cell[YY] = (cell[YY] + search_.ncelldim_[YY]) % search_.ncelldim_[YY];
638 cell[ZZ] = (cell[ZZ] + search_.ncelldim_[ZZ]) % search_.ncelldim_[ZZ];
640 const int ci = search_.getGridCellIndex(cell);
641 const int cellSize = static_cast<int>(search_.cells_[ci].size());
642 /* TODO: Calculate the required PBC shift outside the inner loop */
643 for (; cai < cellSize; ++cai)
645 const int i = search_.cells_[ci][cai];
651 pbc_dx_aiuc(search_.pbc_, xtest_, search_.xref_[i], dx);
652 const real r2 = norm2(dx);
653 if (r2 <= search_.cutoff2_)
670 for (int i = previ_ + 1; i < search_.nref_; ++i)
679 pbc_dx(search_.pbc_, xtest_, search_.xref_[i], dx);
683 rvec_sub(xtest_, search_.xref_[i], dx);
685 const real r2 = norm2(dx);
686 if (r2 <= search_.cutoff2_)
701 void AnalysisNeighborhoodPairSearchImpl::initFoundPair(
702 AnalysisNeighborhoodPair *pair) const
706 *pair = AnalysisNeighborhoodPair();
710 *pair = AnalysisNeighborhoodPair(previ_, testIndex_);
714 } // namespace internal
720 * Search action to find the next neighbor.
722 * Used as the action for AnalysisNeighborhoodPairSearchImpl::searchNext() to
723 * find the next neighbor.
725 * Simply breaks the loop on the first found neighbor.
727 bool withinAction(int /*i*/, real /*r2*/)
733 * Search action find the minimum distance.
735 * Used as the action for AnalysisNeighborhoodPairSearchImpl::searchNext() to
736 * find the nearest neighbor.
738 * With this action, AnalysisNeighborhoodPairSearchImpl::searchNext() always
739 * returns false, and the output is put into the variables passed by pointer
740 * into the constructor. If no neighbors are found, the output variables are
741 * not modified, i.e., the caller must initialize them.
747 * Initializes the action with given output locations.
749 * \param[out] closestPoint Index of the closest reference location.
750 * \param[out] minDist2 Minimum distance squared.
752 * The constructor call does not modify the pointed values, but only
753 * stores the pointers for later use.
754 * See the class description for additional semantics.
756 MindistAction(int *closestPoint, real *minDist2)
757 : closestPoint_(*closestPoint), minDist2_(*minDist2)
761 //! Processes a neighbor to find the nearest point.
762 bool operator()(int i, real r2)
776 GMX_DISALLOW_ASSIGN(MindistAction);
781 /********************************************************************
782 * AnalysisNeighborhood::Impl
785 class AnalysisNeighborhood::Impl
788 typedef AnalysisNeighborhoodSearch::ImplPointer SearchImplPointer;
789 typedef std::vector<SearchImplPointer> SearchList;
791 Impl() : cutoff_(0), mode_(eSearchMode_Automatic)
796 SearchList::const_iterator i;
797 for (i = searchList_.begin(); i != searchList_.end(); ++i)
799 GMX_RELEASE_ASSERT(i->unique(),
800 "Dangling AnalysisNeighborhoodSearch reference");
804 SearchImplPointer getSearch();
806 tMPI::mutex createSearchMutex_;
807 SearchList searchList_;
812 AnalysisNeighborhood::Impl::SearchImplPointer
813 AnalysisNeighborhood::Impl::getSearch()
815 tMPI::lock_guard<tMPI::mutex> lock(createSearchMutex_);
816 // TODO: Consider whether this needs to/can be faster, e.g., by keeping a
817 // separate pool of unused search objects.
818 SearchList::const_iterator i;
819 for (i = searchList_.begin(); i != searchList_.end(); ++i)
826 SearchImplPointer search(new internal::AnalysisNeighborhoodSearchImpl(cutoff_));
827 searchList_.push_back(search);
831 /********************************************************************
832 * AnalysisNeighborhood
835 AnalysisNeighborhood::AnalysisNeighborhood()
840 AnalysisNeighborhood::~AnalysisNeighborhood()
844 void AnalysisNeighborhood::setCutoff(real cutoff)
846 GMX_RELEASE_ASSERT(impl_->searchList_.empty(),
847 "Changing the cutoff after initSearch() not currently supported");
848 impl_->cutoff_ = cutoff;
851 void AnalysisNeighborhood::setMode(SearchMode mode)
856 AnalysisNeighborhood::SearchMode AnalysisNeighborhood::mode() const
861 AnalysisNeighborhoodSearch
862 AnalysisNeighborhood::initSearch(const t_pbc *pbc,
863 const AnalysisNeighborhoodPositions &positions)
865 Impl::SearchImplPointer search(impl_->getSearch());
866 search->init(mode(), pbc, positions);
867 return AnalysisNeighborhoodSearch(search);
870 /********************************************************************
871 * AnalysisNeighborhoodSearch
874 AnalysisNeighborhoodSearch::AnalysisNeighborhoodSearch()
878 AnalysisNeighborhoodSearch::AnalysisNeighborhoodSearch(const ImplPointer &impl)
883 void AnalysisNeighborhoodSearch::reset()
888 AnalysisNeighborhood::SearchMode AnalysisNeighborhoodSearch::mode() const
890 GMX_RELEASE_ASSERT(impl_, "Accessing an invalid search object");
891 return (impl_->usesGridSearch()
892 ? AnalysisNeighborhood::eSearchMode_Grid
893 : AnalysisNeighborhood::eSearchMode_Simple);
896 bool AnalysisNeighborhoodSearch::isWithin(
897 const AnalysisNeighborhoodPositions &positions) const
899 GMX_RELEASE_ASSERT(impl_, "Accessing an invalid search object");
900 internal::AnalysisNeighborhoodPairSearchImpl pairSearch(*impl_);
901 pairSearch.startSearch(positions);
902 return pairSearch.searchNext(&withinAction);
905 real AnalysisNeighborhoodSearch::minimumDistance(
906 const AnalysisNeighborhoodPositions &positions) const
908 GMX_RELEASE_ASSERT(impl_, "Accessing an invalid search object");
909 internal::AnalysisNeighborhoodPairSearchImpl pairSearch(*impl_);
910 pairSearch.startSearch(positions);
911 real minDist2 = impl_->cutoffSquared();
912 int closestPoint = -1;
913 MindistAction action(&closestPoint, &minDist2);
914 (void)pairSearch.searchNext(action);
915 return sqrt(minDist2);
918 AnalysisNeighborhoodPair
919 AnalysisNeighborhoodSearch::nearestPoint(
920 const AnalysisNeighborhoodPositions &positions) const
922 GMX_RELEASE_ASSERT(impl_, "Accessing an invalid search object");
923 internal::AnalysisNeighborhoodPairSearchImpl pairSearch(*impl_);
924 pairSearch.startSearch(positions);
925 real minDist2 = impl_->cutoffSquared();
926 int closestPoint = -1;
927 MindistAction action(&closestPoint, &minDist2);
928 (void)pairSearch.searchNext(action);
929 return AnalysisNeighborhoodPair(closestPoint, 0);
932 AnalysisNeighborhoodPairSearch
933 AnalysisNeighborhoodSearch::startPairSearch(
934 const AnalysisNeighborhoodPositions &positions) const
936 GMX_RELEASE_ASSERT(impl_, "Accessing an invalid search object");
937 Impl::PairSearchImplPointer pairSearch(impl_->getPairSearch());
938 pairSearch->startSearch(positions);
939 return AnalysisNeighborhoodPairSearch(pairSearch);
942 /********************************************************************
943 * AnalysisNeighborhoodPairSearch
946 AnalysisNeighborhoodPairSearch::AnalysisNeighborhoodPairSearch(
947 const ImplPointer &impl)
952 bool AnalysisNeighborhoodPairSearch::findNextPair(AnalysisNeighborhoodPair *pair)
954 bool bFound = impl_->searchNext(&withinAction);
955 impl_->initFoundPair(pair);
959 void AnalysisNeighborhoodPairSearch::skipRemainingPairsForTestPosition()
961 impl_->nextTestPosition();