src/gromacs/math/multidimarray.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 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.
   8  *
   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.
  13  *
  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.
  18  *
  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.
  23  *
  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.
  31  *
  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.
  34  */
  35 /*! \libinternal
  36  * \file
  37  * \brief Declares MultiDimArray.
  38  *
  39  * \author Christian Blau <cblau@gwdg.de>
  40  * \ingroup module_math
  41  * \ingroup module_mdspan
  42  * \inlibraryapi
  43  */
  44
  45 #ifndef GMX_MATH_MULTIDIMARRAY_H_
  46 #define GMX_MATH_MULTIDIMARRAY_H_
  47
  48 #include "gromacs/mdspan/mdspan.h"
  49 #include "gromacs/utility/arrayref.h"
  50
  51 namespace gmx
  52 {
  53
  54 namespace detail
  55 {
  56 //! Same as std::void_t from C++17
  57 template<class...>
  58 using void_t = void;
  59
  60 template<typename T, typename = void>
  61 struct is_resizable : std::false_type
  62 {
  63 };
  64
  65 template<typename T>
  66 struct is_resizable<T, void_t<decltype(std::declval<T>().resize(size_t()))>> : std::true_type
  67 {
  68 };
  69
  70 //! Type has a resize member function callable with size_t argument
  71 template<typename T>
  72 constexpr bool is_resizable_v = is_resizable<T>::value;
  73 } // namespace detail
  74
  75 /*! \libinternal \brief
  76  * Multidimensional array that manages its own memory.
  77  *
  78  * \note No bounds checking when accessing memory
  79  *
  80  * \note That the view holds a valid pointer to the data is a class invariant.
  81  *
  82  * The Container type that stores the data may be resizable (std::vector or similar)
  83  * or static (std::array or similar). Copy and move assignment routines as well as
  84  * swapping are designed to yield good performances in both cases, notably
  85  * foregoing the copy-swap idiom due to the bad performance in swapping std::array.
  86  *
  87  * This class avoids throwing exeptions, apart from the ones that might be thrown
  88  * from the containers during resizing an allocation. (bad_malloc from std::vector
  89  * is a likely candidate)
  90  *
  91  *
  92  * Holds as many elements as required by a multidimensional view.
  93  * \tparam TContainer   Data type container for the data to be stored
  94  *                      as MultiDimArray with random element access and
  95  *                      value_type, refence and const_reference exposed
  96  * \tparam Extents      An extents class describing the array dimensions
  97  *                      as used in module_mdspan
  98  * \tparam LayoutPolicy The data layout as in module_mdspan describes
  99  *                      translation of indices to memory offset.
 100  *                      Defaults to right aligned, so that the right-most index
 101  *                      is contiguous in memory.
 102  */
 103 template<class TContainer, class Extents, class LayoutPolicy = layout_right>
 104 class MultiDimArray
 105 {
 106 public:
 107     //! the type of values that are stored
 108     using value_type = typename TContainer::value_type;
 109     //! reference type to the stored values
 110     using reference = typename TContainer::reference;
 111     //! const reference type to the stored values
 112     using const_reference = typename TContainer::const_reference;
 113     //! the view used to access the data
 114     using view_type = basic_mdspan<value_type, Extents, LayoutPolicy>;
 115     //! const view on the data
 116     using const_view_type = basic_mdspan<const value_type, Extents, LayoutPolicy>;
 117     /*! \brief Iterator type for contiguous iteration over the stored data.
 118      * Used, e.g., in free begin and end functions
 119      */
 120     using iterator = typename ArrayRef<value_type>::iterator;
 121     /*! \brief Const iterator type for contiguous iteration over the stored data.
 122      *  used, e.g., in free begin and end functions
 123      */
 124     using const_iterator = const typename ArrayRef<const value_type>::const_iterator;
 125
 126     static_assert(detail::is_resizable_v<TContainer> == (Extents::rank_dynamic() > 0),
 127                   "Resizable container (e.g. std::vector) requires at least one dynamic rank. "
 128                   "Non-resizable container (e.g. std::array) requires zero dynamic ranks.");
 129
 130     /*! \brief
 131      * Allocate dynamic array data and set view with the dynamic extents.
 132      *
 133      * \param[in] dynamicExtent A parameter pack that describes the dynamic
 134      *                          size of the array. Empty if purely static.
 135      *
 136      * \tparam IndexType        Parameter pack type holding the dynamic
 137      *                          extents of the multidimensional array
 138      */
 139     template<class... IndexType, typename T = TContainer, typename = typename std::enable_if_t<detail::is_resizable_v<T>>>
 140     MultiDimArray(IndexType... dynamicExtent)
 141     {
 142         resize(dynamicExtent...);
 143     }
 144     /*! \brief
 145      * Construction from fixed sized arrays if the array size is static and
 146      * layout policy allows compile time determination of the container size.
 147      *
 148      * Enables the expected initialization
 149      * MultiDimArray<std::array<float, 9>, extents<3,3>> arr = {{1,2...}}
 150      * \tparam T Template parameter for activation via SFINAE.
 151      */
 152     // SFINAE required because std::vector::size isn't constexpr and is_constexpr doesn't exist.
 153     template<typename T = TContainer, typename = typename std::enable_if_t<!detail::is_resizable_v<T>>>
 154     constexpr MultiDimArray(const TContainer& data = {}) noexcept : data_(data), view_(data_.data())
 155     {
 156         static_assert(TContainer().size() == typename view_type::mapping_type().required_span_size(),
 157                       "Non-resizable container type size must match static MultiDimArray size.");
 158     }
 159     //! Copy constructor
 160     constexpr MultiDimArray(const MultiDimArray& o) :
 161         data_(o.data_),
 162         view_(data_.data(), o.view_.extents())
 163     {
 164     }
 165     //! Move constructor
 166     MultiDimArray(MultiDimArray&& o) noexcept :
 167         data_(std::move(o.data_)),
 168         view_(data_.data(), o.view_.extents())
 169     {
 170     }
 171     //! Copy assignment
 172     MultiDimArray& operator=(const MultiDimArray& o)
 173     {
 174         data_ = o.data_;
 175         view_ = view_type(data_.data(), o.view_.extents());
 176         return *this;
 177     }
 178     //! Move assignment
 179     MultiDimArray& operator=(MultiDimArray&& o) noexcept
 180     {
 181         data_ = std::move(o.data_);
 182         view_ = view_type(data_.data(), o.view_.extents());
 183         return *this;
 184     }
 185     //! Swaps content with other
 186     void swap(MultiDimArray& o)
 187     {
 188         using std::swap;
 189         swap(data_, o.data_);
 190         // swap(view_, o.view_) also swaps the pointer to the data and thus does not work
 191         // instead, restore the view as class invariant after the data swapping operation
 192         o.view_ = view_type(o.data_.data(), view_.extents());
 193         view_   = view_type(data_.data(), o.view_.extents());
 194     }
 195     /*! \brief
 196      * Resize the dynamic extents of the array if any and set container size
 197      * accordingly.
 198      *
 199      * Invalidates data and views of this array.
 200      *
 201      * \param[in] dynamicExtent A parameter pack that describes the dynamic
 202      *                          size of the array. Empty if purely static.
 203      * \tparam IndexType        Parameter pack type holding the dynamic
 204      *                          extents of the multidimensional array
 205      */
 206     template<class... IndexType>
 207     void resize(IndexType... dynamicExtent)
 208     {
 209         // use a mapping object to determine the required span size;
 210         layout_right::mapping<Extents> map{ Extents{ dynamicExtent... } };
 211         data_.resize(map.required_span_size());
 212         // to construct a valid view on the data, the container has to be resized before
 213         // the assignment, so that data_.data() is valid
 214         view_ = view_type(data_.data(), dynamicExtent...);
 215     }
 216     /*! \brief Data access via multidimensional indices.
 217      * This allows referencing rank R array elements as array(x_0,x_1,x_2, .., x_R)
 218      *
 219      * \param[in] index multidimensional indices as parameter pack
 220      *                  the number of parameters must match the rank of the array.
 221      *
 222      * \returns reference to array element
 223      */
 224     template<class... IndexType>
 225     reference operator()(IndexType... index) noexcept
 226     {
 227         return view_(index...);
 228     }
 229     /*! \brief Const data access via multidimensional indices.
 230      * This allows referencing rank R array elements as array(x_0,x_1,x_2, .., x_R)
 231      *
 232      * \param[in] index multidimensional indices as parameter pack
 233      *                  the number of parameters must match the rank of the array.
 234      *
 235      * \returns const reference to array element
 236      */
 237     template<class... IndexType>
 238     constexpr const_reference operator()(IndexType... index) const noexcept
 239     {
 240         return view_(index...);
 241     }
 242     /*! \brief Contiguous access to the data.
 243      * \returns ArrayRef to stored data.
 244      */
 245     ArrayRef<value_type> toArrayRef() { return { data_.data(), data_.data() + data_.size() }; }
 246     /*! \brief Contiguous const access to the data.
 247      * \returns ArrayRef to stored data.
 248      */
 249     constexpr ArrayRef<const value_type> toArrayRef() const
 250     {
 251         return { data_.data(), data_.data() + data_.size() };
 252     }
 253     /*! \brief Return the extent.
 254      * \param[in] k dimension to query for extent
 255      * \returns extent along specified dimension
 256      */
 257     constexpr typename view_type::index_type extent(int k) const noexcept
 258     {
 259         return view_.extent(k);
 260     }
 261     //! Conversion to multidimensional view on the data
 262     constexpr view_type asView() noexcept { return view_; }
 263     //! Conversion to const multidimensional view on the data
 264     constexpr const_view_type asConstView() const noexcept
 265     {
 266         return { data_.data(), view_.mapping() };
 267     }
 268
 269 private:
 270     //! The contiguous data that is equipped with multidimensional indexing in this class
 271     TContainer data_;
 272     //! Multidimensional view into data_.
 273     view_type view_;
 274 };
 275
 276 //! Free MultiDimArray begin function addressing its contiguous memory.
 277 template<class TContainer, class Extents>
 278 constexpr typename MultiDimArray<TContainer, Extents>::const_iterator
 279 begin(const MultiDimArray<TContainer, Extents>& multiDimArray)
 280 {
 281     return multiDimArray.toArrayRef().begin();
 282 }
 283
 284 //! Free MultiDimArray begin function addressing its contiguous memory.
 285 template<class TContainer, class Extents>
 286 constexpr typename MultiDimArray<TContainer, Extents>::iterator begin(MultiDimArray<TContainer, Extents>& multiDimArray)
 287 {
 288     return multiDimArray.toArrayRef().begin();
 289 }
 290
 291 //! Free MultiDimArray end function addressing its contiguous memory.
 292 template<class TContainer, class Extents>
 293 constexpr typename MultiDimArray<TContainer, Extents>::const_iterator
 294 end(const MultiDimArray<TContainer, Extents>& multiDimArray)
 295 {
 296     return multiDimArray.toArrayRef().end();
 297 }
 298
 299 //! Free MultiDimArray end function addressing its contiguous memory.
 300 template<class TContainer, class Extents>
 301 constexpr typename MultiDimArray<TContainer, Extents>::iterator end(MultiDimArray<TContainer, Extents>& multiDimArray)
 302 {
 303     return multiDimArray.toArrayRef().end();
 304 }
 305
 306 //! Swap function
 307 template<class TContainer, class Extents>
 308 void swap(MultiDimArray<TContainer, Extents>& a, MultiDimArray<TContainer, Extents>& b) noexcept
 309 {
 310     a.swap(b);
 311 }
 312
 313 } // namespace gmx
 314
 315 #endif // GMX_MATH_MULTIDIMARRAY_H_