schreier-sims.hpp

//
// libsemigroups - C++ library for semigroups and monoids
// Copyright (C) 2019-2026 James D. Mitchell
//
// 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/>.
//
 
// TODO(later)
//
// 0. check use of:
//    * element_type []
//    * const_value_type []
//    * reference        []
//    * const_element_reference  []
//    * internal_element_type []
//    * internal_const_element_type []
//    * internal_reference []
//    * internal_const_reference []
// 1. iterator to the elements (Finn requires this)
// 2. stabilizer member func
// 3. change base
// 4. random version
// 5. try it with Digraphs
 
// TODO:
// * code coverage (can't do it on MacBook Pro at present)
 
#ifndef LIBSEMIGROUPS_SCHREIER_SIMS_HPP_
#define LIBSEMIGROUPS_SCHREIER_SIMS_HPP_
 
#include <array>          // for array
#include <cstddef>        // for size_t
#include <cstdint>        // for uint64_t
#include <iterator>       // for distance
#include <limits>         // for numeric_limits
#include <memory>         // for make_unique
#include <string>         // for operator+, basic_string
#include <string_view>    // for string_view
#include <type_traits>    // for is_same
#include <unordered_set>  // for unordered_set
 
#include "adapters.hpp"   // for action, degree, inverse
#include "config.hpp"     // for LIBSEMIGROUPS_HPCOMBI_ENABLED
#include "debug.hpp"      // for LIBSEMIGROUPS_ASSERT
#include "exception.hpp"  // for LIBSEMIGROUPS_EXCEPTION
#include "transf.hpp"     // for Perm
#include "types.hpp"      // for SmallestInteger
 
#include "detail/bruidhinn-traits.hpp"  // for BruidhinnTraits
#include "detail/containers.hpp"        // for Array2, StaticTriVector2
#include "detail/fmt.hpp"               // for format
#include "detail/int-range.hpp"         // for IntRange
#include "detail/stl.hpp"               // for has_call_operator_size_t_v
 
namespace libsemigroups {


  template <size_t N, typename Point, typename Element>
  struct SchreierSimsTraits {
    using index_type = size_t;

    using domain_type = detail::IntRange<Point>;

    using point_type = Point;

    using element_type = Element;

    using Action = libsemigroups::ImageRightAction<element_type, Point>;

    using Degree = libsemigroups::Degree<element_type>;

    using EqualTo = libsemigroups::EqualTo<element_type>;

    using Inverse = libsemigroups::Inverse<element_type>;

    using One = libsemigroups::One<element_type>;

    using Product = libsemigroups::Product<element_type>;

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

  template <size_t N,
            typename Point   = typename SmallestInteger<N>::type,
            typename Element = LeastPerm<N>,
            typename Traits  = SchreierSimsTraits<N, Point, Element>>
  class SchreierSims : private detail::BruidhinnTraits<Element> {
    static_assert(std::is_same_v<Point, typename Traits::point_type>,
                  "incompatible point types, Traits::point_type and Point "
                  "must be the same");
    static_assert(
        std::is_same_v<Element, typename Traits::element_type>,
        "incompatible element types, Traits::element_type and Element "
        "must be the same");
    static_assert(N <= std::numeric_limits<Point>::max(),
                  "the first template parameter N cannot be expressed as value "
                  "of type equal to the second template parameter Point");
 
    using internal_element_type =
        typename detail::BruidhinnTraits<Element>::internal_value_type;
    using internal_reference =
        typename detail::BruidhinnTraits<Element>::internal_reference;
    using internal_const_element_type =
        typename detail::BruidhinnTraits<Element>::internal_const_value_type;
 
   public:
    using element_type = typename detail::BruidhinnTraits<Element>::value_type;

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

    using element_reference =
        typename detail::BruidhinnTraits<Element>::reference;

    using point_type = Point;

    using domain_type = typename Traits::domain_type;

    using index_type = typename Traits::index_type;

    using Action = typename Traits::Action;

    using Degree = typename Traits::Degree;

    using EqualTo = typename Traits::EqualTo;

    using Inverse = typename Traits::Inverse;

    using One = typename Traits::One;

    using Product = typename Traits::Product;

    using Swap = typename Traits::Swap;
 
    static_assert(detail::has_call_operator_size_t_v<One>,
                  "Traits::One must implement operator()(size_t)!");
 
   private:
    // TODO replace Array2 by TriArray2 everywhere below
    std::array<point_type, N>                          _base;
    index_type                                         _base_size;
    domain_type                                        _domain;
    bool                                               _finished;
    internal_element_type                              _one;
    detail::StaticTriVector2<point_type, N>            _orbits;
    detail::Array2<bool, N>                            _orbits_lookup;
    detail::StaticTriVector2<internal_element_type, N> _strong_gens;
    mutable internal_element_type                      _tmp_element1;
    mutable internal_element_type                      _tmp_element2;
    detail::Array2<internal_element_type, N>           _transversal;
    detail::Array2<internal_element_type, N>           _inversal;
 
   public:
    SchreierSims();

    SchreierSims& init();
 
    ~SchreierSims();

    SchreierSims(SchreierSims&&);

    // TODO this requires more tests
    SchreierSims(SchreierSims const& that);

    // TODO fix in the same way as the move constructor above
    SchreierSims& operator=(SchreierSims&&);

    // TODO this requires more tests
    SchreierSims& operator=(SchreierSims const& that);

    bool add_generator_no_checks(const_element_reference x);

    bool add_generator(const_element_reference x);

    // Not noexcept because strong_generator isn't
    [[nodiscard]] const_element_reference generator(index_type index) const {
      return strong_generator(0, index);
    }

    // Not noexcept because strong_generator_no_checks isn't
    [[nodiscard]] const_element_reference
    generator_no_checks(index_type index) const {
      return strong_generator_no_checks(0, index);
    }

    // Not noexcept because number_of_strong_generators_no_checks is not.
    [[nodiscard]] size_t number_of_generators() const;

    [[nodiscard]] size_t number_of_strong_generators(index_type depth) const {
      throw_if_bad_depth(depth);
      return number_of_strong_generators_no_checks(depth);
    }

    // Not noexcept because StaticTriVector2::size(size_t) is not.
    [[nodiscard]] size_t
    number_of_strong_generators_no_checks(index_type depth) const {
      return _strong_gens.size(depth);
    }

    [[nodiscard]] const_element_reference
    strong_generator(index_type depth, index_type index) const;

    // not noexcept because StaticTriVector2 isn't
    [[nodiscard]] const_element_reference
    strong_generator_no_checks(index_type depth, index_type index) const {
      return this->to_external_const(_strong_gens.at(depth, index));
    }

    // not noexcept because Array2::operator[] isn't
    [[nodiscard]] const_element_reference
    transversal_element_no_checks(index_type depth, point_type pt) const {
      return this->to_external_const(_transversal[depth][pt]);
    }

    // Not noexcept because throws
    [[nodiscard]] const_element_reference
    transversal_element(index_type depth, point_type pt) const;

    // Not noexcept because std::array::operator[] isn't
    [[nodiscard]] const_element_reference
    inverse_transversal_element_no_checks(index_type depth,
                                          point_type pt) const {
      return this->to_external_const(_inversal[depth][pt]);
    }

    [[nodiscard]] const_element_reference
    inverse_transversal_element(index_type depth, point_type pt) const;

    // Not noexcept because std::array::operator[] isn't
    [[nodiscard]] bool orbit_lookup_no_checks(index_type depth,
                                              point_type pt) const {
      return _orbits_lookup[depth][pt];
    }

    [[nodiscard]] bool orbit_lookup(index_type depth, point_type pt) const {
      throw_if_bad_depth(depth);
      throw_if_point_gt_degree(pt);
      return orbit_lookup_no_checks(depth, pt);
    }

    // Not noexcept because StaticTriVector2::size isn't
    [[nodiscard]] bool empty() const {
      return _strong_gens.size(0) == 0;
    }

    // not noexcept because run isn't
    [[nodiscard]] uint64_t size();

    [[nodiscard]] uint64_t current_size() const;

    // TODO tests
    // not noexcept because internal_sift isn't
    void sift_inplace_no_checks(element_reference x) const {
      // changes x in place, and uses _tmp_element1
      internal_sift(this->to_internal(x));
    }

    // TODO tests
    // not noexcept because can throw
    void sift_inplace(element_reference x) const {
      throw_if_bad_degree(x);
      sift_inplace_no_checks(x);
    }

    // not noexcept because internal_sift isn't
    const_element_reference sift_no_checks(const_element_reference x) const {
      this->to_external(_tmp_element2) = x;
      // changes _tmp_element2 in place, and uses _tmp_element1
      internal_sift(_tmp_element2);
      return this->to_external_const(_tmp_element2);
    }

    // not noexcept because can throw
    [[nodiscard]] const_element_reference
    sift(const_element_reference x) const {
      throw_if_bad_degree(x);
      return sift_no_checks(x);
    }

    // At present JDM thinks this will always return false if finished()
    // isn't true, because we either run to the end or haven't done anything at
    // all.
    // not noexcept because sift_no_checks isn't
    [[nodiscard]] bool currently_contains(const_element_reference x) const;

    [[nodiscard]] bool contains(const_element_reference x);

    [[nodiscard]] const_element_reference one() const
        noexcept(noexcept(this->to_external_const(_one))) {
      return this->to_external_const(_one);
    }

    [[nodiscard]] bool finished() const noexcept {
      return _finished;
    }

    // TODO rename to add_basepoint?
    // not noexcept because can throw
    SchreierSims& add_base_point(point_type pt);

    // Not noexcept because std::array::operator[] isn't
    [[nodiscard]] point_type base_no_checks(index_type index) const {
      return _base[index];
    }

    // not noexcept because can throw
    [[nodiscard]] point_type base(index_type index) const {
      throw_if_bad_depth(index);
      return base_no_checks(index);
    }

    [[nodiscard]] size_t base_size() const noexcept {
      return _base_size;
    }

    // not noexcept because it can call mem fns that aren't
    void run();
 
   private:
    // SchreierSims - validation - private
 
    bool is_valid_degree(point_type x) const noexcept {
      return
#ifdef LIBSEMIGROUPS_HPCOMBI_ENABLED
          std::is_same_v<HPCombi::Perm16, element_type> ||
#endif
          x == N;
    }
 
    void throw_if_bad_degree(const_element_reference x,
                             std::string_view        arg_pos = "1st") const;
 
    void throw_if_bad_depth(size_t           depth,
                            std::string_view arg_pos = "1st") const;
 
    void throw_if_point_gt_degree(point_type       pt,
                                  std::string_view arg_pos = "1st") const;
 
    void throw_if_point_not_in_orbit(index_type       depth,
                                     point_type       pt,
                                     std::string_view depth_arg_pos = "1st",
                                     std::string_view pt_arg_pos = "2nd") const;

    // SchreierSims - member functions - private
 
    bool internal_equal_to(internal_const_element_type x,
                           internal_const_element_type y) const
        noexcept(noexcept(
            EqualTo()(this->to_external_const(x),
                      this->to_external_const(
                          y)) && noexcept(this->to_external_const(x)))) {
      return EqualTo()(this->to_external_const(x), this->to_external_const(y));
    }
 
    // Used by copy constructor and assignment operator.
    void init_strong_gens_traversal_inversal(SchreierSims const& that);
    void free_strong_gens_traversal_inversal();
 
    void internal_add_base_point(point_type pt);
    void orbit_enumerate(index_type depth, index_type first = 0);
    void orbit_add_gen(index_type depth, internal_element_type gen);
    void orbit_add_point(index_type            depth,
                         internal_element_type x,
                         point_type            pt);
    // Changes _tmp_element2 in-place, and returns the depth reached in the
    // sifting.
    index_type internal_sift(internal_reference x) const;
 
    typename domain_type::const_iterator
    first_non_fixed_point(internal_const_element_type x) const {
      auto const& y = this->to_external_const(x);
      return std::find_if(_domain.cbegin(), _domain.cend(), [&y](auto pt) {
        return pt != Action()(pt, y);
      });
    }
  };

  namespace schreier_sims {

    // TODO(later) example
    // TODO (from RC):
    // 1. Implement orbit refinement heuristic for intersection.
    // 2. Make the Screier-Sims object during runtime, since we compute the
    //    stabilizers of the intersection already.
    // 3. Refactor for more generality (i.e. so the template parameters N don't
    //    all have to be the same.
    //
    // TODO (from JDM):
    // * use the no_checks mem fns of SchreierSims now that they exist
    template <size_t N, typename Point, typename Element, typename Traits>
    void intersection(SchreierSims<N, Point, Element, Traits>& T,
                      SchreierSims<N, Point, Element, Traits>& S1,
                      SchreierSims<N, Point, Element, Traits>& S2);
  }  // namespace schreier_sims

  template <size_t N, typename Point, typename Element, typename Traits>
  [[nodiscard]] std::string
  to_human_readable_repr(SchreierSims<N, Point, Element, Traits> const& S,
                         size_t max_width = 72);
}  // namespace libsemigroups
 
#include "schreier-sims.tpp"
 
#endif  // LIBSEMIGROUPS_SCHREIER_SIMS_HPP_