2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2010,2014,2015,2018, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
37 #ifndef GMX_TOPOLOGY_BLOCK_H
38 #define GMX_TOPOLOGY_BLOCK_H
46 #include "gromacs/utility/basedefinitions.h"
48 #include "gromacs/utility/gmxassert.h"
59 /*! \brief Division of a range of indices into consecutive blocks
61 * A range of consecutive indices 0 to full.range.end() is divided
62 * into numBlocks() consecutive blocks of consecutive indices.
63 * Block b contains indices i for which block(b).begin() <= i < block(b).end().
65 class RangePartitioning
68 /*! \brief Struct for returning the range of a block.
70 * Can be used in a range loop.
75 /*! \brief An iterator that loops over integers */
79 iterator(int value) : value_(value) {}
81 operator int () const { return value_; }
83 operator int &() { return value_; }
85 int operator* () const { return value_; }
86 //! Inequality comparison
87 bool operator!= (const iterator other) { return value_ != other; }
88 //! Increment operator
89 iterator &operator++() { ++value_; return *this; }
90 //! Increment operator
91 iterator operator++(int) { iterator tmp(*this); ++value_; return tmp; }
96 /*! \brief Constructor, constructs a range starting at 0 with 0 blocks */
104 /*! \brief Begin iterator/value */
105 const iterator begin() const { return begin_; };
106 /*! \brief End iterator/value */
107 const iterator end() const { return end_; };
109 /*! \brief The number of items in the block */
112 return end_ - begin_;
115 /*! \brief Returns whether \p index is within range of the block */
116 bool inRange(int index) const
118 return (begin_ <= index && index < end_);
122 const int begin_; /**< The start index of the block */
123 const int end_; /**< The end index of the block */
126 /*! \brief Returns the number of blocks */
127 int numBlocks() const
129 return index_.size() - 1;
132 /*! \brief Returns the size of the block with index \p blockIndex */
133 Block block(int blockIndex) const
135 return Block(index_[blockIndex], index_[blockIndex + 1]);
138 /*! \brief Returns the full range */
139 Block fullRange() const
141 return Block(index_.front(), index_.back());
144 /*! \brief Returns a range starting at \p blockIndexBegin and ending at \p blockIndexEnd */
145 Block subRange(int blockIndexBegin,
146 int blockIndexEnd) const
148 return Block(index_[blockIndexBegin], index_[blockIndexEnd]);
151 /*! \brief Returns true when all blocks have size 0 or numBlocks()=0 */
152 bool allBlocksHaveSizeOne() const
154 return (index_.back() == numBlocks());
157 /*! \brief Appends a block of size \p blockSize at the end of the range
159 * \note blocksize has to be >= 1
161 void appendBlock(int blockSize)
163 GMX_ASSERT(blockSize > 0, "block sizes should be >= 1");
164 index_.push_back(index_.back() + blockSize);
167 /*! \brief Removes all blocks */
173 /*! \brief Reduces the number of blocks to \p newNumBlocks
175 * \note \p newNumBlocks should be <= numBlocks().
177 void reduceNumBlocks(int newNumBlocks)
179 GMX_ASSERT(newNumBlocks <= numBlocks(), "Can only shrink to fewer blocks");
180 index_.resize(newNumBlocks + 1);
183 /*! \brief Sets the partitioning to \p numBlocks blocks each of size 1 */
184 void setAllBlocksSizeOne(int numBlocks);
186 /*! \brief Returns the raw block index array, avoid using this */
187 std::vector<int> &rawIndex()
193 std::vector<int> index_ = { 0 }; /**< The list of block begin/end indices */
197 #endif // __cplusplus
199 /* Deprecated, C-style version of BlockRanges */
200 typedef struct t_block
203 int blockSize(int blockIndex) const
205 GMX_ASSERT(blockIndex < nr, "blockIndex should be in range");
206 return index[blockIndex + 1] - index[blockIndex];
208 #endif // __cplusplus
210 int nr; /* The number of blocks */
211 int *index; /* Array of indices (dim: nr+1) */
212 int nalloc_index; /* The allocation size for index */
215 typedef struct t_blocka
217 int nr; /* The number of blocks */
218 int *index; /* Array of indices in a (dim: nr+1) */
219 int nra; /* The number of atoms */
220 int *a; /* Array of atom numbers in each group */
222 /* Block i (0<=i<nr) runs from */
223 /* index[i] to index[i+1]-1. There will */
224 /* allways be an extra entry in index */
225 /* to terminate the table */
226 int nalloc_index; /* The allocation size for index */
227 int nalloc_a; /* The allocation size for a */
230 void init_block(t_block *block);
231 void init_blocka(t_blocka *block);
232 t_blocka *new_blocka(void);
233 /* allocate new block */
235 void done_block(t_block *block);
236 void done_blocka(t_blocka *block);
238 void copy_blocka(const t_blocka *src, t_blocka *dest);
240 void stupid_fill_block(t_block *grp, int natom, gmx_bool bOneIndexGroup);
241 /* Fill a block structure with numbers identical to the index
242 * (0, 1, 2, .. natom-1)
243 * If bOneIndexGroup, then all atoms are lumped in one index group,
244 * otherwise there is one atom per index entry
247 void stupid_fill_blocka(t_blocka *grp, int natom);
248 /* Fill a block structure with numbers identical to the index
249 * (0, 1, 2, .. natom-1)
250 * There is one atom per index entry
253 void pr_block(FILE *fp, int indent, const char *title, const t_block *block, gmx_bool bShowNumbers);
254 void pr_blocka(FILE *fp, int indent, const char *title, const t_blocka *block, gmx_bool bShowNumbers);