Utilities/detail/CombinatoricIndices.ipp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 #include "Acts/Utilities/detail/CombinatoricIndices.hpp"
0012
0013 #include <format>
0014
0015 namespace Acts::detail {
0016
0017 template <std::size_t K>
0018 CombinatoricIndices<K>::CombinatoricIndices(const std::size_t N)
0019   requires(K > 0)
0020     : m_N{N} {
0021   if (N < K) {
0022     throw std::invalid_argument(
0023         std::format("CombinatoricIndices() - The set size {:} needs at least "
0024                     "to exceed the number of elements to draw {:}",
0025                     N, K));
0026   }
0027
0028   for (std::size_t slot = 0; slot < m_elementOccurance.size(); ++slot) {
0029     std::size_t maxCombNumb = 0;
0030     // Calculate in how many combinations the integer can occur
0031     const std::size_t nPrime = N - 1ul - slot;
0032     const std::size_t kPrime = K - 1ul - slot;
0033
0034     for (std::size_t leftCmb = nPrime; leftCmb >= kPrime; --leftCmb) {
0035       maxCombNumb += binomial(leftCmb, kPrime);
0036       m_elementOccurance[slot].push_back(maxCombNumb);
0037       if (leftCmb == 0ul) {
0038         break;
0039       }
0040     }
0041   }
0042 }
0043
0044 template <std::size_t K>
0045 std::size_t CombinatoricIndices<K>::size() const {
0046   return m_elementOccurance[0].back();
0047 }
0048
0049 template <std::size_t K>
0050 std::size_t CombinatoricIndices<K>::intervalSize() const {
0051   return m_N;
0052 }
0053
0054 template <std::size_t K>
0055 std::array<std::size_t, K> CombinatoricIndices<K>::draw(
0056     const std::size_t combination) const {
0057   if (combination >= size()) {
0058     throw std::invalid_argument(
0059         std::format("CombinatoricIndices::draw({:}) - You cannot draw more "
0060                     "than {:} combinations",
0061                     combination, size()));
0062   }
0063   std::array<std::size_t, K> indices{};
0064   for (std::size_t s = 0; s < indices.size(); ++s) {
0065     indices[s] = drawIndex(combination, s);
0066   }
0067   return indices;
0068 }
0069
0070 template <std::size_t K>
0071 std::size_t CombinatoricIndices<K>::drawIndex(const std::size_t combination,
0072                                               const std::size_t slot) const {
0073   std::size_t cmbPrime = combination;
0074   if (slot >= K) {
0075     throw std::runtime_error(
0076         std::format("CombinatoricIndices: The slot must not exceed {:} the "
0077                     "number of drawn elements {:}",
0078                     slot, K));
0079   }
0080   for (std::size_t s = 0; s <= slot; ++s) {
0081     const auto slotItr =
0082         std::ranges::upper_bound(m_elementOccurance[s], cmbPrime);
0083     if (s == slot) {
0084       return std::distance(m_elementOccurance[s].begin(), slotItr) + slot;
0085     }
0086     cmbPrime -= ((*slotItr) - m_elementOccurance[s].front());
0087   }
0088   return 0ul;
0089 }
0090
0091 template <std::size_t K>
0092 CombinatoricIndices<K>::iterator::iterator(const CombinatoricIndices* parent,
0093                                            const std::size_t combination)
0094     : m_parent{parent}, m_itr{combination} {
0095   updateArray();
0096 }
0097
0098 template <std::size_t K>
0099 CombinatoricIndices<K>::iterator&
0100 CombinatoricIndices<K>::iterator::operator++() {
0101   ++m_itr;
0102   updateArray();
0103   return *this;
0104 }
0105
0106 template <std::size_t K>
0107 bool CombinatoricIndices<K>::iterator::operator==(const iterator& other) const {
0108   return m_parent == other.m_parent && m_itr == other.m_itr;
0109 }
0110
0111 template <std::size_t K>
0112 bool CombinatoricIndices<K>::iterator::operator!=(const iterator& other) const {
0113   return m_parent != other.m_parent || m_itr != other.m_itr;
0114 }
0115
0116 template <std::size_t K>
0117 const CombinatoricIndices<K>::IndexArray&
0118 CombinatoricIndices<K>::iterator::operator*() const {
0119   return m_array;
0120 }
0121
0122 template <std::size_t K>
0123 void CombinatoricIndices<K>::iterator::updateArray() {
0124   if (m_parent == nullptr || m_itr >= m_parent->size()) {
0125     m_array =
0126         filledArray<std::size_t, K>(std::numeric_limits<std::size_t>::max());
0127   } else {
0128     for (std::size_t s = 0; s < m_array.size(); ++s) {
0129       m_array[s] = m_parent->drawIndex(m_itr, s);
0130     }
0131   }
0132 }
0133
0134 template <std::size_t K>
0135 CombinatoricIndices<K>::iterator CombinatoricIndices<K>::iterator::operator+(
0136     const std::size_t idx) const {
0137   return iterator{m_parent, m_itr + idx};
0138 }
0139
0140 }  // namespace Acts::detail