Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <algorithm>
       7             : #include <array>
       8             : #include <cstddef>
       9             : 
      10             : #include "Utilities/Algorithm.hpp"
      11             : #include "Utilities/Array.hpp"
      12             : #include "Utilities/ConstantExpressions.hpp"
      13             : #include "Utilities/Numeric.hpp"
      14             : 
      15             : /*!
      16             :  * \ingroup DataStructuresGroup
      17             :  *
      18             :  * \brief Iterate over all nonzero index permutations for a Levi-Civita
      19             :  * symbol.
      20             :  *
      21             :  * \details This class provides an iterator that allows you to loop over
      22             :  * only the nonzero index permutations of a Levi-Civita symbol of dimension
      23             :  * `dimension`. Inside the loop, the operator `()`
      24             :  * returns an `std::array` containing an ordered list of the indices of this
      25             :  * permutation, the operator `[]` returns a specific index from the same
      26             :  * `std::array`, and the function `sign()` returns the sign of the
      27             :  * Levi-Civita symbol for this permutation.
      28             :  *
      29             :  * Note that the Levi-Civita symbol in odd dimensions (like 3) is invariant
      30             :  * under cyclic permutations of its indices, but in even dimensions (like 4)
      31             :  * it changes sign under cyclic permutations.
      32             :  *
      33             :  * \example
      34             :  * \snippet Test_LeviCivitaIterator.cpp levi_civita_iterator_example
      35             :  */
      36             : template <size_t Dim>
      37           1 : class LeviCivitaIterator {
      38             :  public:
      39           0 :   explicit LeviCivitaIterator() = default;
      40             : 
      41             :   /// \cond HIDDEN_SYMBOLS
      42             :   ~LeviCivitaIterator() = default;
      43             :   /// @{
      44             :   /// No copy or move semantics
      45             :   LeviCivitaIterator(const LeviCivitaIterator<Dim>&) = delete;
      46             :   LeviCivitaIterator(LeviCivitaIterator<Dim>&&) = delete;
      47             :   LeviCivitaIterator<Dim>& operator=(const LeviCivitaIterator<Dim>&) = delete;
      48             :   LeviCivitaIterator<Dim>& operator=(LeviCivitaIterator<Dim>&&) = delete;
      49             :   /// @}
      50             :   /// \endcond
      51             : 
      52             :   // Return false if the end of the loop is reached
      53           0 :   explicit operator bool() const { return valid_; }
      54             : 
      55             :   // Increment the current permutation
      56           0 :   LeviCivitaIterator& operator++() {
      57             :     ++permutation_;
      58             :     if (permutation_ >= signs_.size()) {
      59             :       valid_ = false;
      60             :     }
      61             :     return *this;
      62             :   }
      63             : 
      64             :   /// Return a `std::array` containing the current multi-index, an ordered list
      65             :   /// of indices for the current permutation.
      66           1 :   const std::array<size_t, Dim> operator()() const {
      67             :     return static_cast<std::array<size_t, Dim>>(indexes_[permutation_]);
      68             :   }
      69             : 
      70             :   /// Return a specific index from the multi-index of the current permutation.
      71           1 :   const size_t& operator[](const size_t i) const {
      72             :     return indexes_[permutation_][i];
      73             :   }
      74             : 
      75             :   /// Return the sign of the Levi-Civita symbol for the current permutation.
      76           1 :   int sign() const { return signs_[permutation_]; };
      77             : 
      78             :  private:
      79             :   static constexpr cpp20::array<cpp20::array<size_t, Dim>, factorial(Dim)>
      80           0 :   indexes() {
      81             :     cpp20::array<cpp20::array<size_t, Dim>, factorial(Dim)> indexes{};
      82             :     cpp20::array<size_t, Dim> index{};
      83             :     cpp2b::iota(index.begin(), index.end(), size_t(0));
      84             : 
      85             :     indexes[0] = index;
      86             : 
      87             :     // cpp20::next_permutation generates the different permutations of index
      88             :     // loop over them to fill in the rest of the permutations in indexes
      89             :     size_t permutation = 1;
      90             :     while (cpp20::next_permutation(index.begin(), index.end())) {
      91             :       indexes[permutation] = index;
      92             :       ++permutation;
      93             :     }
      94             :     return indexes;
      95             :   };
      96             : 
      97           0 :   static constexpr cpp20::array<int, factorial(Dim)> signs() {
      98             :     cpp20::array<int, factorial(Dim)> signs{};
      99             :     // By construction, the sign of the first permutation is +1
     100             :     signs[0] = 1;
     101             : 
     102             :     // How do you know whether the corresponding Levi Civita symbol is
     103             :     // +1 or -1? To find out, compute the number, in the factoradic number
     104             :     // system, corresponding to each permutation, and sum the digits. If the sum
     105             :     // is even (odd), the corresponding Levi-Civita symbol will be +1 (-1). This
     106             :     // works because there are Dim! unique permutations of the Levi Civita
     107             :     // symbol indices, so each one can be represented by a
     108             :     // (`Dim`)-digit factorial-number-system number. For more on the
     109             :     // factoradic number system, see
     110             :     // https://en.wikipedia.org/wiki/Factorial_number_system
     111             :     auto factoradic_counter = cpp20::array<size_t, Dim>();
     112             :     for (size_t permutation = 1; permutation < factorial(Dim); ++permutation) {
     113             :       for (size_t i = 0; i < Dim; ++i) {
     114             :         factoradic_counter[i]++;
     115             :         if (factoradic_counter[i] < i + 1) {
     116             :           break;
     117             :         } else {
     118             :           factoradic_counter[i] = 0;
     119             :         }
     120             :       }
     121             :       signs[permutation] =
     122             :           (cpp2b::accumulate(factoradic_counter.begin(),
     123             :                              factoradic_counter.end(), size_t(0)) %
     124             :            2) == 0
     125             :               ? 1
     126             :               : -1;
     127             :     }
     128             :     return signs;
     129             :   };
     130             : 
     131             :   // Note: here and throughout, use cpp20::array, which is constexpr
     132             :   // (unlike std::array), to enable constexpr signs_ and indexes_.
     133             :   static constexpr cpp20::array<cpp20::array<size_t, Dim>, factorial(Dim)>
     134           0 :       indexes_ = indexes();
     135           0 :   static constexpr cpp20::array<int, factorial(Dim)> signs_{signs()};
     136           0 :   size_t permutation_{0};
     137           0 :   bool valid_{true};
     138             : };