konieczny.hpp

//
// libsemigroups - C++ library for semigroups and monoids
// Copyright (C) 2019-2025 Finn Smith
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program.  If not, see <http://www.gnu.org/licenses/>.
//
// This file contains a generic implementation of Konieczny's algorithm,
// originally for computing subsemigroups of the boolean matrix monoid.
 
// TODO(later):
// * exception safety!
// * expose iterators to relevant things in D-classes, in particular elements
 
#ifndef LIBSEMIGROUPS_KONIECZNY_HPP_
#define LIBSEMIGROUPS_KONIECZNY_HPP_
 
#include <algorithm>  // for binary_search
#include <cstddef>    // for size_t
#include <initializer_list>
#include <set>            // for set
#include <type_traits>    // for is_pointer
#include <unordered_map>  // for unordered_map
#include <unordered_set>  // for unordered_set
#include <utility>        // for pair, make_pair
#include <vector>         // for vector
 
#include "action.hpp"     // for LeftAction, RightAction
#include "adapters.hpp"   // for Lambda, etc
#include "constants.hpp"  // for UNDEFINED
#include "debug.hpp"      // for LIBSEMIGROUPS_ASSERT
#include "exception.hpp"  // for LIBSEMIGROUPS_EXCEPTION
#include "runner.hpp"     // for Runner
 
#include "detail/pool.hpp"    // for detail::Pool
#include "detail/report.hpp"  // for report_default
#include "detail/timer.hpp"   // for Timer
 
#include "detail/bruidhinn-traits.hpp"  // for BruidhinnTraits
 
namespace libsemigroups {


  template <typename Element>
  struct KoniecznyTraits {
    using element_type = typename detail::BruidhinnTraits<Element>::value_type;

    using const_element_type =
        typename detail::BruidhinnTraits<Element>::const_value_type;

    using lambda_value_type =
        typename ::libsemigroups::LambdaValue<element_type>::type;

    using rho_value_type =
        typename ::libsemigroups::RhoValue<element_type>::type;

    using rank_state_type = typename ::libsemigroups::RankState<element_type>;

    using lambda_orb_type
        = RightAction<element_type,
                      lambda_value_type,
                      ImageRightAction<element_type, lambda_value_type>>;

    using rho_orb_type
        = LeftAction<element_type,
                     rho_value_type,
                     ImageLeftAction<element_type, rho_value_type>>;

    using Lambda = ::libsemigroups::Lambda<element_type, lambda_value_type>;

    using Rho = ::libsemigroups::Rho<element_type, rho_value_type>;

    using Product = ::libsemigroups::Product<element_type>;

    using Rank = ::libsemigroups::Rank<element_type, rank_state_type>;

    using One = ::libsemigroups::One<element_type>;

    using ElementHash = ::libsemigroups::Hash<element_type>;

    using EqualTo = ::libsemigroups::EqualTo<element_type>;

    using Swap = ::libsemigroups::Swap<element_type>;

    using Less = ::libsemigroups::Less<element_type>;

    using Degree = ::libsemigroups::Degree<element_type>;
  };

  template <typename Element, typename Traits = KoniecznyTraits<Element>>
  class Konieczny : public Runner, private detail::BruidhinnTraits<Element> {
    // pointers are not currently supported
    static_assert(!std::is_pointer_v<Element>,
                  "Pointer types are not currently supported by Konieczny");

    // Konieczny - aliases - private
 
    using internal_element_type =
        typename detail::BruidhinnTraits<Element>::internal_value_type;
    using internal_const_element_type =
        typename detail::BruidhinnTraits<Element>::internal_const_value_type;
    using internal_const_reference =
        typename detail::BruidhinnTraits<Element>::internal_const_reference;
    using internal_reference =
        typename detail::BruidhinnTraits<Element>::internal_reference;
 
    using PairHash = Hash<std::pair<size_t, size_t>>;
 
   public:
    // Konieczny - aliases - public

    using element_type = typename Traits::element_type;

    using const_element_type = typename Traits::const_element_type;

    using const_reference =
        typename detail::BruidhinnTraits<Element>::const_reference;

    using D_class_index_type = size_t;

    using lambda_value_type = typename Traits::lambda_value_type;

    using lambda_orb_type = typename Traits::lambda_orb_type;

    using rho_value_type = typename Traits::rho_value_type;

    using rho_orb_type = typename Traits::rho_orb_type;

    using Degree = typename Traits::Degree;

    using EqualTo = typename Traits::EqualTo;

    using Lambda = typename Traits::Lambda;

    using Less = typename Traits::Less;

    using One = typename Traits::One;

    using Product = typename Traits::Product;

    using Rank = typename Traits::Rank;

    using Rho = typename Traits::Rho;

    using Swap = typename Traits::Swap;
 
   private:
    // Konieczny - aliases - private
 
    using lambda_orb_index_type    = typename lambda_orb_type::index_type;
    using rho_orb_index_type       = typename rho_orb_type::index_type;
    using rank_type                = size_t;
    using rank_state_type          = typename Traits::rank_state_type;
    using left_indices_index_type  = size_t;
    using right_indices_index_type = size_t;

    // Konieczny - internal structs - private
 
    struct InternalHash : private detail::BruidhinnTraits<Element> {
      size_t operator()(internal_const_element_type x) const {
        return Hash<Element>()(this->to_external_const(x));
      }
    };
 
    struct InternalEqualTo : private detail::BruidhinnTraits<Element> {
      bool operator()(internal_const_reference x,
                      internal_const_reference y) const {
        return EqualTo()(this->to_external_const(x),
                         this->to_external_const(y));
      }
    };
 
    struct InternalLess : private detail::BruidhinnTraits<Element> {
      bool operator()(internal_const_reference x,
                      internal_const_reference y) const {
        return Less()(this->to_external_const(x), this->to_external_const(y));
      }
    };
 
    struct InternalVecEqualTo : private detail::BruidhinnTraits<Element> {
      size_t operator()(std::vector<internal_element_type> const& x,
                        std::vector<internal_element_type> const& y) const {
        LIBSEMIGROUPS_ASSERT(x.size() == y.size());
        return std::equal(x.cbegin(), x.cend(), y.cbegin(), InternalEqualTo());
      }
    };
 
    struct InternalVecFree : private detail::BruidhinnTraits<Element> {
      void operator()(std::vector<internal_element_type> const& x) {
        for (auto it = x.cbegin(); it != x.cend(); ++it) {
          this->internal_free(*it);
        }
      }
    };
 
    struct InternalVecCopy : private detail::BruidhinnTraits<Element> {
      void operator()(std::vector<internal_element_type> const& source,
                      std::vector<internal_element_type>&       target) {
        InternalVecFree()(target);
        for (auto it = source.cbegin(); it != source.cend(); ++it) {
          target.push_back(this->internal_copy(*it));
        }
      }
    };
 
    struct OneParamLambda {
      lambda_value_type operator()(const_reference x) const {
        lambda_value_type lval;
        Lambda()(lval, x);
        return lval;
      }
    };
 
    struct OneParamRho {
      rho_value_type operator()(const_reference x) const {
        rho_value_type rval;
        Rho()(rval, x);
        return rval;
      }
    };
 
    struct InternalRank {
      template <typename SFINAE = size_t>
      auto operator()(void*, const_reference x)
          -> std::enable_if_t<std::is_void_v<typename rank_state_type::type>,
                              SFINAE> {
        return Rank()(x);
      }
 
      template <typename SFINAE = size_t>
      auto operator()(rank_state_type* state, const_reference x)
          -> std::enable_if_t<!std::is_void_v<typename rank_state_type::type>,
                              SFINAE> {
        return Rank()(*state, x);
      }
    };
 
    struct RepInfo {
      RepInfo(D_class_index_type    D_idx,
              internal_element_type elt,
              lambda_orb_index_type lambda_idx,
              rho_orb_index_type    rho_idx)
          : _D_idx(D_idx),
            _elt(elt),
            _lambda_idx(lambda_idx),
            _rho_idx(rho_idx) {}
 
      D_class_index_type    _D_idx;
      internal_element_type _elt;
      lambda_orb_index_type _lambda_idx;
      rho_orb_index_type    _rho_idx;
    };
 
   public:
    // Konieczny - constructor and destructor - public

    Konieczny();

    Konieczny& init();

    Konieczny(Konieczny const& that);

    Konieczny(Konieczny&& that);

    Konieczny& operator=(Konieczny const&);

    Konieczny& operator=(Konieczny&&);
 
    ~Konieczny();

    // Konieczny - forward class declarations - public/private

    class DClass : protected detail::BruidhinnTraits<Element> {
      // This friend is only here so that the virtual contains(x, lpos, rpos)
      // method and the cbegin_left_indices etc. methods can be private.
      friend class Konieczny<Element, Traits>;

      // DClass - aliases - protected
#ifndef LIBSEMIGROUPS_PARSED_BY_DOXYGEN
      using konieczny_type    = Konieczny<Element, Traits>;
      using internal_set_type = std::
          unordered_set<internal_element_type, InternalHash, InternalEqualTo>;

      // DClass - default constructor - private
 
     private:
      DClass();
 
      void clear();

      // DClass - constructors - protected
     protected:
      DClass(Konieczny* parent, internal_reference rep);
 
      DClass(DClass const& that);
 
      DClass& operator=(DClass const& that);
 
      DClass(DClass&&)            = delete;
      DClass& operator=(DClass&&) = delete;
 
#endif  // ndef LIBSEMIGROUPS_PARSED_BY_DOXYGEN
 
     public:
      // DClass - destructor - public
      virtual ~DClass();
      // DClass - member functions - public

      const_reference rep() const {
        return this->to_external_const(_rep);
      }

      size_t size() const {
        LIBSEMIGROUPS_ASSERT(this->class_computed());
        return number_of_L_classes() * number_of_R_classes() * size_H_class();
      }

      size_t number_of_L_classes() const {
        LIBSEMIGROUPS_ASSERT(_left_mults.size() > 0);
        LIBSEMIGROUPS_ASSERT(this->class_computed());
        return _left_mults.size();
      }

      size_t number_of_R_classes() const {
        // compute_right_mults();
        LIBSEMIGROUPS_ASSERT(_right_mults.size() > 0);
        LIBSEMIGROUPS_ASSERT(this->class_computed());
        return _right_mults.size();
      }

      size_t size_H_class() const {
        // compute_H_class();
        LIBSEMIGROUPS_ASSERT(_H_class.size() > 0);
        LIBSEMIGROUPS_ASSERT(this->class_computed());
        return _H_class.size();
      }

      bool is_regular_D_class() const noexcept {
        return _is_regular_D_class;
      }

      bool contains(const_reference x);

      virtual size_t number_of_idempotents() const {
        return 0;
      }
 
     protected:
#ifndef LIBSEMIGROUPS_PARSED_BY_DOXYGEN
      // DClass - iterators - protected
      using const_iterator =
          typename std::vector<internal_element_type>::const_iterator;
 
      const_iterator cbegin_left_reps() {
        compute_left_reps();
        return _left_reps.cbegin();
      }
 
      const_iterator cend_left_reps() {
        compute_left_reps();
        return _left_reps.cend();
      }
 
      const_iterator cbegin_right_reps() {
        compute_right_reps();
        return _right_reps.cbegin();
      }
 
      const_iterator cend_right_reps() {
        compute_right_reps();
        return _right_reps.cend();
      }
 
      const_iterator cbegin_left_mults() {
        compute_left_mults();
        return _left_mults.cbegin();
      }
 
      const_iterator cend_left_mults() {
        compute_left_mults();
        return _left_mults.cend();
      }
 
      const_iterator cbegin_right_mults() {
        compute_right_mults();
        return _right_mults.cbegin();
      }
 
      const_iterator cend_right_mults() {
        compute_right_mults();
        return _right_mults.cend();
      }
 
      const_iterator cbegin_H_class() {
        compute_H_class();
        return _H_class.cbegin();
      }
 
      const_iterator cend_H_class() {
        compute_H_class();
        return _H_class.cend();
      }
 
      internal_element_type left_mults_inv(size_t i) {
        compute_left_mults_inv();
        return _left_mults_inv[i];
      }
 
      internal_element_type right_mults_inv(size_t i) {
        compute_right_mults_inv();
        return _right_mults_inv[i];
      }
 
      internal_element_type H_class_no_checks(size_t i) const {
        return _H_class[i];
      }

      // DClass - initialisation member functions - protected
      virtual void compute_frame()           = 0;
      virtual void compute_left_indices()    = 0;
      virtual void compute_left_mults()      = 0;
      virtual void compute_left_mults_inv()  = 0;
      virtual void compute_left_reps()       = 0;
      virtual void compute_right_indices()   = 0;
      virtual void compute_right_mults()     = 0;
      virtual void compute_right_mults_inv() = 0;
      virtual void compute_right_reps()      = 0;
      virtual void compute_H_class()         = 0;

      // DClass - containment - protected
 
      // Returns whether the element \p x belongs to this
      // \f$\mathscr{D}\f$-class.
      //
      // Given an element \p x of the semigroup represented by \c parent, this
      // function returns whether \p x is an element of the
      // \f$\mathscr{D}\f$-class represented by \c this. If \p x is not an
      // element of the semigroup, then the behaviour is undefined.
      // This member function involved computing most of the frame for
      // \c this, if it is not already known.
      bool contains_no_checks(internal_const_reference x);
 
      // Returns whether the element \p x belongs to this
      // \f$\mathscr{D}\f$-class.
      //
      // This overload of DClass::contains_no_checks is provided in order to
      // avoid recalculating the rank of \p x when it is already known.
      bool contains_no_checks(internal_const_reference x, size_t rank);
 
      // Returns whether the element \p x belongs to this
      // \f$\mathscr{D}\f$-class.
      //
      // This overload of DClass::contains_no_checks is provided in order to
      // avoid recalculating the rank, lambda value, and rho value of \p x
      // when they are already known.
      bool contains_no_checks(internal_const_reference x,
                              size_t                   rank,
                              lambda_orb_index_type    lpos,
                              rho_orb_index_type       rpos);
 
      // Returns whether the element \p x belongs to this
      // \f$\mathscr{D}\f$-class.
      //
      // This overload of DClass::contains_no_checks is provided in order to
      // avoid recalculating the lambda value and rho value of \p x  when they
      // are already known.
      virtual bool contains_no_checks(internal_const_reference x,
                                      lambda_orb_index_type    lpos,
                                      rho_orb_index_type       rpos)
          = 0;
 
      virtual bool contains(const_reference       x,
                            lambda_orb_index_type lpos,
                            rho_orb_index_type    rpos)
          = 0;

      // DClass - accessor member functions - protected
 
      size_t number_of_left_reps_no_checks() const noexcept {
        return _left_reps.size();
      }
 
      size_t number_of_right_reps_no_checks() const noexcept {
        return _right_reps.size();
      }
 
      size_t size_H_class_no_checks() const noexcept {
        return _H_class.size();
      }
 
      void push_left_mult(internal_const_reference x);
 
      void push_left_mult_inv(internal_const_reference x);
 
      void push_right_mult(internal_const_reference x);
 
      void push_right_mult_inv(internal_const_reference x);
 
      void push_left_rep(internal_const_reference x);
 
      void push_right_rep(internal_const_reference x);
 
      bool class_computed() const noexcept {
        return _class_computed;
      }
 
      bool mults_computed() const noexcept {
        return _mults_computed;
      }
 
      bool reps_computed() const noexcept {
        return _reps_computed;
      }
 
      bool H_class_computed() const noexcept {
        return _H_class_computed;
      }
 
      void set_class_computed(bool x) noexcept {
        _class_computed = x;
      }
 
      void set_mults_computed(bool x) noexcept {
        _mults_computed = x;
      }
 
      void set_reps_computed(bool x) noexcept {
        _reps_computed = x;
      }
 
      void set_H_class_computed(bool x) noexcept {
        _H_class_computed = x;
      }
 
      Konieczny* parent() const noexcept {
        return _parent;
      }
 
      void set_parent(Konieczny* x) noexcept {
        _parent = x;
      }
 
      // Watch out! Doesn't copy its argument
      void push_back_H_class(internal_element_type x) {
        _H_class.push_back(x);
      }
 
      std::vector<internal_element_type>& H_class() {
        return _H_class;
      }
 
      lambda_value_type& tmp_lambda_value() const noexcept {
        return _tmp_lambda_value;
      }
 
      rho_value_type& tmp_rho_value() const noexcept {
        return _tmp_rho_value;
      }
 
      rank_type rank() const noexcept {
        return _rank;
      }
 
      internal_set_type& internal_set() const noexcept {
        return _tmp_internal_set;
      }
 
      // Elements must be freed before next used
      std::vector<internal_element_type>& internal_vec() const noexcept {
        return _tmp_internal_vec;
      }
 
      internal_reference unsafe_rep() noexcept {
        return _rep;
      }
 
      std::vector<lambda_orb_index_type>& left_indices() {
        return _left_indices;
      }
 
      std::vector<rho_orb_index_type>& right_indices() {
        return _right_indices;
      }
 
     protected:
      // DClass - index iterators - protected
 
      typename std::vector<left_indices_index_type>::const_iterator
      cbegin_left_indices() {
        compute_left_indices();
        return _left_indices.cbegin();
      }
 
      typename std::vector<left_indices_index_type>::const_iterator
      cend_left_indices() {
        compute_left_indices();
        return _left_indices.cend();
      }
 
      typename std::vector<right_indices_index_type>::const_iterator
      cbegin_right_indices() {
        compute_right_indices();
        return _right_indices.cbegin();
      }
 
      typename std::vector<right_indices_index_type>::const_iterator
      cend_right_indices() {
        compute_right_indices();
        return _right_indices.cend();
      }
#endif
 
     private:
      // DClass - member functions - private
 
      // Returns a set of representatives of L- or R-classes covered by \c
      // this.
      //
      // The \f$\mathscr{D}\f$-classes of the parent semigroup are enumerated
      // either by finding representatives of all L-classes or all R-classes.
      // This member function returns the representatives obtainable by
      // multipliying the representatives  by generators on either the left or
      // right.
      std::vector<RepInfo>& covering_reps();

      // DClass - data - private
      bool                               _class_computed;
      std::vector<internal_element_type> _H_class;
      bool                               _H_class_computed;
      bool                               _is_regular_D_class;
      std::vector<lambda_orb_index_type> _left_indices;
      std::vector<internal_element_type> _left_mults;
      std::vector<internal_element_type> _left_mults_inv;
      std::vector<internal_element_type> _left_reps;
      bool                               _mults_computed;
      Konieczny*                         _parent;
      rank_type                          _rank;
      internal_element_type              _rep;
      bool                               _reps_computed;
      std::vector<rho_orb_index_type>    _right_indices;
      std::vector<internal_element_type> _right_mults;
      std::vector<internal_element_type> _right_mults_inv;
      std::vector<internal_element_type> _right_reps;
      mutable internal_set_type _tmp_internal_set;  // Does not own its elements
      mutable std::vector<RepInfo>
          _tmp_rep_info_vec;  // RepInfos do not own their rep
      mutable std::vector<internal_element_type> _tmp_internal_vec;  // Does not
      mutable lambda_value_type                  _tmp_lambda_value;
      mutable rho_value_type                     _tmp_rho_value;
    };
 
   private:
    class RegularDClass;
    class NonRegularDClass;
 
   public:
    // Konieczny - member functions - public

    size_t number_of_generators() const noexcept {
      return _gens.empty() ? 0 : _gens.size() - 1;
    }

    const_reference generator(size_t pos) const {
      if (pos >= number_of_generators()) {
        LIBSEMIGROUPS_EXCEPTION(
            "index out of bounds, expected value in [{}, {}) found {}",
            0,
            number_of_generators(),
            pos);
      }
      return this->to_external_const(_gens[pos]);
    }

    size_t number_of_D_classes() {
      run();
      return std::distance(cbegin_current_D_classes(),
                           cend_current_D_classes());
    }

    size_t current_number_of_D_classes() const {
      return std::distance(cbegin_current_D_classes(),
                           cend_current_D_classes());
    }

    size_t number_of_regular_D_classes() {
      run();
      return current_number_of_regular_D_classes();
    }

    size_t current_number_of_regular_D_classes() const {
      return std::distance(cbegin_current_regular_D_classes(),
                           cend_current_regular_D_classes());
    }

    size_t number_of_L_classes() {
      run();
      return current_number_of_L_classes();
    }

    size_t current_number_of_L_classes() const {
      size_t val = 0;
      std::for_each(
          cbegin_current_D_classes(),
          cend_current_D_classes(),
          [&val](DClass const& D) { val += D.number_of_L_classes(); });
      return val;
    }

    size_t number_of_regular_L_classes() {
      run();
      return current_number_of_regular_L_classes();
    }

    size_t current_number_of_regular_L_classes() const {
      size_t val = 0;
      std::for_each(
          cbegin_current_regular_D_classes(),
          cend_current_regular_D_classes(),
          [&val](DClass const& D) { val += D.number_of_L_classes(); });
      return val;
    }

    size_t number_of_R_classes() {
      run();
      return current_number_of_R_classes();
    }

    size_t current_number_of_R_classes() const {
      size_t val = 0;
      std::for_each(
          cbegin_current_D_classes(),
          cend_current_D_classes(),
          [&val](DClass const& D) { val += D.number_of_R_classes(); });
      return val;
    }

    size_t number_of_regular_R_classes() {
      run();
      return current_number_of_regular_R_classes();
    }

    size_t current_number_of_regular_R_classes() const {
      size_t val = 0;
      std::for_each(
          cbegin_current_regular_D_classes(),
          cend_current_regular_D_classes(),
          [&val](DClass const& D) { val += D.number_of_R_classes(); });
      return val;
    }

    size_t number_of_H_classes() {
      run();
      return current_number_of_H_classes();
    }

    size_t current_number_of_H_classes() const {
      size_t val = 0;
      std::for_each(cbegin_current_D_classes(),
                    cend_current_D_classes(),
                    [&val](DClass const& D) {
                      val += (D.number_of_R_classes()
                              * D.number_of_L_classes());
                    });
      return val;
    }

    size_t number_of_idempotents() {
      run();
      return current_number_of_idempotents();
    }

    size_t current_number_of_idempotents() const {
      size_t val = 0;
      std::for_each(
          cbegin_current_regular_D_classes(),
          cend_current_regular_D_classes(),
          [&val](RegularDClass const& D) { val += D.number_of_idempotents(); });
      return val;
    }

    size_t number_of_regular_elements() {
      run();
      return current_number_of_regular_elements();
    }

    size_t current_number_of_regular_elements() const {
      size_t val = 0;
      std::for_each(cbegin_current_regular_D_classes(),
                    cend_current_regular_D_classes(),
                    [&val](DClass const& D) { val += D.size(); });
      return val;
    }

    size_t size() {
      run();
      return current_size();
    }

    size_t current_size() const {
      size_t val = 0;
      std::for_each(cbegin_current_D_classes(),
                    cend_current_D_classes(),
                    [&val](DClass const& D) { val += D.size(); });
      return val;
    }

    size_t degree() const noexcept {
      return _degree;
    }

    bool contains(const_reference x) {
      return Degree()(x) == degree()
             && get_containing_D_class(this->to_internal_const(x), true)
                    != UNDEFINED;
    }

    DClass& D_class_of_element(const_reference x) {
      D_class_index_type i
          = get_containing_D_class(this->to_internal_const(x), true);
      if (i == UNDEFINED) {
        LIBSEMIGROUPS_EXCEPTION(
            "the argument does not belong to this semigroup!");
      }
      return *_D_classes[i];
    }

    bool is_regular_element(const_reference x) {
      return contains(x)
             && is_regular_element_no_checks(this->to_internal_const(x));
    }

    template <typename T>
    void add_generators(T const& first, T const& last);

    void add_generator(const_reference gen) {
      add_generators(&gen, &gen + 1);
    }

    // Konieczny - iterators - public

    using const_iterator
        = detail::BruidhinnConstIterator<element_type,
                                         std::vector<internal_element_type>>;

    //! \sa cend_generators
    const_iterator cbegin_generators() const noexcept {
      return const_iterator(_gens.cbegin());
    }

    //! \sa cbegin_generators
    const_iterator cend_generators() const noexcept {
      return const_iterator(_gens.cend() - 1);
    }
 
    // forward decl, defined in konieczny.tpp
    template <typename T>
    struct DClassIteratorTraits;

    using const_d_class_iterator
        = detail::ConstIteratorStateless<DClassIteratorTraits<DClass>>;

    //! \no_libsemigroups_except
    // not noexcept because operator++ isn't necessarily
    const_d_class_iterator cbegin_current_D_classes() const {
      auto it = _D_classes.cbegin();
      if (_run_initialised) {
        return const_d_class_iterator(it)
               + (_adjoined_identity_contained ? 0 : 1);
      } else {
        return const_d_class_iterator(it);
      }
    }

    //! \exceptions
    //! \no_libsemigroups_except
    const_d_class_iterator cend_current_D_classes() const noexcept {
      return const_d_class_iterator(_D_classes.cend());
    }

    //! \returns
    //! A value of type \c const_d_class_iterator.
    const_d_class_iterator cbegin_D_classes() {
      run();
      return cbegin_current_D_classes();
    }

    //! \returns
    //! A value of type \c const_d_class_iterator.
    const_d_class_iterator cend_D_classes() {
      run();
      return cend_current_D_classes();
    }

    using const_regular_d_class_iterator
        = detail::ConstIteratorStateless<DClassIteratorTraits<RegularDClass>>;

    //! \no_libsemigroups_except
    //!
    // not noexcept because operator++ isn't necessarily
    const_regular_d_class_iterator cbegin_current_regular_D_classes() const {
      auto it = _regular_D_classes.cbegin();
      if (_run_initialised) {
        return const_regular_d_class_iterator(it)
               + (_adjoined_identity_contained ? 0 : 1);
      } else {
        return const_regular_d_class_iterator(it);
      }
    }

    //! \exceptions
    //! \no_libsemigroups_except
    const_regular_d_class_iterator
    cend_current_regular_D_classes() const noexcept {
      return const_regular_d_class_iterator(_regular_D_classes.cend());
    }

    //! \exceptions
    //! \no_libsemigroups_except
    //!
    const_regular_d_class_iterator cbegin_regular_D_classes() {
      run();
      return cbegin_current_regular_D_classes();
    }

    //!
    //! \exceptions
    //! \no_libsemigroups_except
    const_regular_d_class_iterator cend_regular_D_classes() {
      run();
      return cend_current_regular_D_classes();
    }
 
   private:
    using PoolGuard = detail::PoolGuard<internal_element_type>;

    // Konieczny - utility methods - private
 
    // assumes its argument has valid lambda/rho values
    bool is_regular_element_no_checks(internal_const_reference x,
                                      lambda_orb_index_type    lpos = UNDEFINED,
                                      rho_orb_index_type rpos = UNDEFINED) {
      LIBSEMIGROUPS_ASSERT(_lambda_orb.finished() && _rho_orb.finished());
      lpos = lpos != UNDEFINED ? lpos : get_lpos(x);
      rpos = rpos != UNDEFINED ? rpos : get_rpos(x);
      return get_lambda_group_index(x, lpos, rpos) != UNDEFINED;
    }
 
    lambda_orb_index_type get_lpos(internal_const_reference x) const {
      Lambda()(_tmp_lambda_value1, this->to_external_const(x));
      return _lambda_orb.position(_tmp_lambda_value1);
    }
 
    rho_orb_index_type get_rpos(internal_const_reference x) const {
      Rho()(_tmp_rho_value1, this->to_external_const(x));
      return _rho_orb.position(_tmp_rho_value1);
    }
 
    // Returns a lambda orb index corresponding to a group H-class in the R-
    // class of \p x.
    // asserts its argument has lambda/rho values in the orbits.
    // modifies _tmp_lambda_value1
    // modifies _tmp_rho_value1
    lambda_orb_index_type get_lambda_group_index(internal_const_reference x,
                                                 lambda_orb_index_type    lpos
                                                 = UNDEFINED,
                                                 rho_orb_index_type rpos
                                                 = UNDEFINED);
 
    // Finds a group index of a H-class in the L-class of \p x.
    // modifies _tmp_lambda_value1
    // modifies _tmp_rho_value1
    rho_orb_index_type get_rho_group_index(internal_const_reference x,
                                           lambda_orb_index_type    lpos
                                           = UNDEFINED,
                                           rho_orb_index_type rpos = UNDEFINED);

    // TODO(later): it must be possible to do better than this
    void idem_in_H_class(internal_reference       res,
                         internal_const_reference x) const;

    // modifies _tmp_lambda_value1
    void make_idem(internal_reference x);

    void group_inverse(internal_element_type&   res,
                       internal_const_reference id,
                       internal_const_reference x) const;

    // modifies _tmp_lambda_value and _tmp_rho_value
    bool is_group_index(internal_const_reference x,
                        internal_const_reference y,
                        lambda_orb_index_type    lpos = UNDEFINED,
                        rho_orb_index_type       rpos = UNDEFINED) const;
 
    // pass full_check = true to use the contains method of the D-classes
    // instead of the contains_no_checks
    D_class_index_type get_containing_D_class(internal_const_reference x,
                                              bool const full_check = false);
 
    void add_to_D_maps(D_class_index_type d);

    // Konieczny - accessor member functions - private
 
    void add_D_class(RegularDClass* D);
    void add_D_class(NonRegularDClass* D);
 
    auto cbegin_internal_generators() const noexcept {
      return _gens.cbegin();
    }
 
    auto cend_internal_generators() const noexcept {
      return _gens.cend();
    }
 
    detail::Pool<internal_element_type>& element_pool() const {
      return _element_pool;
    }
 
    size_t max_rank() const noexcept {
      if (_ranks.empty()) {
        return UNDEFINED;
      }
      return *_ranks.rbegin();
    }

    // Konieczny - initialisation member functions - private
 
    void free_internals();
 
    void init_run();
 
    void init_data();
 
    void init_rank_state_and_rep_vecs();
 
    void compute_orbs();

    // Konieczny - validation member functions - private
 
    void throw_if_bad_degree(const_reference x) const {
      size_t const n = Degree()(x);
      if (degree() != UNDEFINED && n != degree()) {
        LIBSEMIGROUPS_EXCEPTION(
            "element has degree {} but should have degree {}", n, degree());
      }
    }
 
    template <typename Iterator>
    void throw_if_bad_degree(Iterator first, Iterator last) const;

    // Konieczny - Runner methods - private
 
    bool finished_impl() const override;
    void run_impl() override;
    void report_before_run();
    void report_progress();
    void report_after_run();

    // Konieczny - data - private
 
    bool                                         _adjoined_identity_contained;
    bool                                         _can_accept_generators;
    std::vector<DClass*>                         _D_classes;
    std::vector<std::vector<D_class_index_type>> _D_rels;
    bool                                         _data_initialised;
    size_t                                       _degree;
    mutable detail::Pool<internal_element_type>  _element_pool;
    std::vector<internal_element_type>           _gens;
    std::unordered_map<std::pair<rho_orb_index_type, lambda_orb_index_type>,
                       lambda_orb_index_type,
                       PairHash>
        _group_indices;
    std::unordered_map<std::pair<rho_orb_index_type, lambda_orb_index_type>,
                       rho_orb_index_type,
                       PairHash>
                    _group_indices_rev;
    lambda_orb_type _lambda_orb;
    std::unordered_map<lambda_orb_index_type, std::vector<D_class_index_type>>
                                      _lambda_to_D_map;
    std::vector<std::vector<RepInfo>> _nonregular_reps;
    internal_element_type             _one;
    rank_state_type*                  _rank_state;
    std::set<rank_type>               _ranks;
    std::vector<RegularDClass*>       _regular_D_classes;
    std::vector<std::vector<RepInfo>> _regular_reps;
    size_t                            _reps_processed;
    rho_orb_type                      _rho_orb;
    std::unordered_map<rho_orb_index_type, std::vector<D_class_index_type>>
                              _rho_to_D_map;
    bool                      _run_initialised;
    mutable lambda_value_type _tmp_lambda_value1;
    mutable lambda_value_type _tmp_lambda_value2;
    mutable rho_value_type    _tmp_rho_value1;
    mutable rho_value_type    _tmp_rho_value2;
  };  // class Konieczny
 
  template <typename Element, typename Traits>
  bool Konieczny<Element, Traits>::finished_impl() const {
    return _ranks.empty() && _run_initialised;
  }

  template <typename Element, typename Traits>
  [[nodiscard]] std::string
  to_human_readable_repr(Konieczny<Element, Traits> const& x);

  template <typename Element, typename Traits>
  [[nodiscard]] std::string
  to_human_readable_repr(typename Konieczny<Element, Traits>::DClass const& x);

  namespace konieczny {

    template <typename Element, typename Traits, typename T>
    void add_generators(Konieczny<Element, Traits> K, T const& coll) {
      static_assert(!std::is_pointer_v<T>,
                    "the template parameter T must not be a pointer");
      K.add_generators(coll.begin(), coll.end());
    }

    template <typename Element, typename Traits>
    void add_generators(
        Konieczny<Element, Traits>& K,
        std::initializer_list<
            typename Konieczny<Element, Traits>::const_element_type> coll) {
      K.add_generators(coll.begin(), coll.end());
    }
 
  }  // namespace konieczny


  // TODO(1) to tpp
  template <template <typename...> typename Thing,
            typename Container,
            typename Traits = KoniecznyTraits<typename Container::value_type>>
  [[nodiscard]] std::enable_if_t<
      std::is_same_v<Konieczny<typename Container::value_type, Traits>,
                     Thing<typename Container::value_type, Traits>>,
      Konieczny<typename Container::value_type, Traits>>
  make(Container const& gens) {
    if (gens.size() == 0) {
      LIBSEMIGROUPS_EXCEPTION(
          "expected a positive number of generators, but got 0");
    }
    Konieczny<typename Container::value_type, Traits> result;
    result.add_generators(std::begin(gens), std::end(gens));
    return result;
  }

  template <template <typename...> typename Thing,
            typename Element,
            typename Traits = KoniecznyTraits<Element>>
  [[nodiscard]] std::enable_if_t<
      std::is_same_v<Konieczny<Element, Traits>, Thing<Element, Traits>>,
      Konieczny<Element, Traits>>
  make(std::initializer_list<Element> const& gens) {
    return make<Konieczny, std::initializer_list<Element>>(gens);
  }
 
}  // namespace libsemigroups
#include "konieczny.tpp"
#endif  // LIBSEMIGROUPS_KONIECZNY_HPP_