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