gh-35915: Is fully commutative for Coxeter group elements
    
<!-- Please provide a concise, informative and self-explanatory title.
-->
<!-- Don't put issue numbers in the title. Put it in the Description
below. -->
<!-- For example, instead of "Fixes #12345", use "Add a new method to
multiply two integers" -->

### 📚 Description

This is adding a method "is_fully_commutative" to all Coxeter group
elements.

On the way, introduce iterators over the sets of reduced words of a
given element.

<!-- Describe your changes here in detail. -->
<!-- Why is this change required? What problem does it solve? -->
<!-- If this PR resolves an open issue, please link to it here. For
example "Fixes #12345". -->
<!-- If your change requires a documentation PR, please link it
appropriately. -->

### 📝 Checklist

<!-- Put an `x` in all the boxes that apply. It should be `[x]` not `[x
]`. -->

- [x] The title is concise, informative, and self-explanatory.
- [x] The description explains in detail what this PR is about.
- [ ] I have linked a relevant issue or discussion.
- [x] I have created tests covering the changes.
- [x] I have updated the documentation accordingly.

### ⌛ Dependencies

<!-- List all open PRs that this PR logically depends on
- #12345: short description why this is a dependency
- #34567: ...
-->

<!-- If you're unsure about any of these, don't hesitate to ask. We're
here to help! -->
    
URL: #35915
Reported by: Frédéric Chapoton
Reviewer(s): Travis Scrimshaw

Author: Release Manager

src/sage/categories/coxeter_groups.py

@@ -251,6 +251,53 @@ def braid_group_as_finitely_presented_group(self):
             rels = self.braid_relations()
             return F / [prod(S[I.index(i)] for i in l) * prod(S[I.index(i)]**-1 for i in reversed(r)) for l, r in rels]
 
+        def braid_orbit_iter(self, word):
+            r"""
+            Iterate over the braid orbit of a word ``word`` of indices.
+
+            The input word does not need to be a reduced expression of
+            an element.
+
+            INPUT:
+
+            - ``word`` -- a list (or iterable) of indices in
+              ``self.index_set()``
+
+            OUTPUT:
+
+            all lists that can be obtained from
+            ``word`` by replacements of braid relations
+
+            EXAMPLES::
+
+                sage: W = CoxeterGroups().example()
+                sage: sorted(W.braid_orbit_iter([0, 1, 2, 1]))                # optional - sage.combinat sage.groups
+                [[0, 1, 2, 1], [0, 2, 1, 2], [2, 0, 1, 2]]
+            """
+            word = list(word)
+            from sage.combinat.root_system.braid_orbit import BraidOrbit
+
+            braid_rels = self.braid_relations()
+            I = self.index_set()
+
+            from sage.rings.integer_ring import ZZ
+            be_careful = any(i not in ZZ for i in I)
+
+            if be_careful:
+                Iinv = {i: j for j, i in enumerate(I)}
+                word = [Iinv[i] for i in word]
+                braid_rels = [[[Iinv[i] for i in l],
+                               [Iinv[i] for i in r]] for l, r in braid_rels]
+
+            orb = BraidOrbit(word, braid_rels)
+
+            if be_careful:
+                for word in orb:
+                    yield [I[i] for i in word]
+            else:
+                for I in orb:
+                    yield list(I)
+
         def braid_orbit(self, word):
             r"""
             Return the braid orbit of a word ``word`` of indices.
@@ -314,26 +361,7 @@ def braid_orbit(self, word):
 
                 :meth:`.reduced_words`
             """
-            word = list(word)
-            from sage.combinat.root_system.braid_orbit import BraidOrbit
-
-            braid_rels = self.braid_relations()
-            I = self.index_set()
-
-            from sage.rings.integer_ring import ZZ
-            be_careful = any(i not in ZZ for i in I)
-
-            if be_careful:
-                Iinv = {i: j for j, i in enumerate(I)}
-                word = [Iinv[i] for i in word]
-                braid_rels = [[[Iinv[i] for i in l],
-                               [Iinv[i] for i in r]] for l, r in braid_rels]
-
-            orb = BraidOrbit(word, braid_rels)
-
-            if be_careful:
-                return [[I[i] for i in word] for word in orb]
-            return [list(I) for I in orb]
+            return list(self.braid_orbit_iter(word))
 
         def __iter__(self):
             r"""
@@ -1493,7 +1521,7 @@ def descents(self, side='right', index_set=None, positive=False):
             return [i for i in index_set if self.has_descent(i, side=side,
                                                              positive=positive)]
 
-        def is_grassmannian(self, side="right"):
+        def is_grassmannian(self, side="right") -> bool:
             """
             Return whether ``self`` is Grassmannian.
 
@@ -1529,6 +1557,49 @@ def is_grassmannian(self, side="right"):
             """
             return len(self.descents(side=side)) <= 1
 
+        def is_fully_commutative(self) -> bool:
+            r"""
+            Check if ``self`` is a fully-commutative element.
+
+            We use the characterization that an element `w` in a Coxeter
+            system `(W,S)` is fully-commutative if and only if for every pair
+            of generators `s,t \in S` for which `m(s,t)>2`, no reduced
+            word of `w` contains the 'braid' word `sts...` of length
+            `m(s,t)` as a contiguous subword. See [Ste1996]_.
+
+            EXAMPLES::
+
+                sage: W = CoxeterGroup(['A', 3])
+                sage: len([1 for w in W if w.is_fully_commutative()])
+                14
+                sage: W = CoxeterGroup(['B', 3])
+                sage: len([1 for w in W if w.is_fully_commutative()])
+                24
+
+            TESTS::
+
+                sage: W = CoxeterGroup(matrix(2,2,[1,7,7,1]),index_set='ab')
+                sage: len([1 for w in W if w.is_fully_commutative()])
+                13
+            """
+            word = self.reduced_word()
+            from sage.combinat.root_system.braid_orbit import is_fully_commutative as is_fully_comm
+
+            group = self.parent()
+            braid_rels = group.braid_relations()
+            I = group.index_set()
+
+            from sage.rings.integer_ring import ZZ
+            be_careful = any(i not in ZZ for i in I)
+
+            if be_careful:
+                Iinv = {i: j for j, i in enumerate(I)}
+                word = [Iinv[i] for i in word]
+                braid_rels = [[[Iinv[i] for i in l],
+                               [Iinv[i] for i in r]] for l, r in braid_rels]
+
+            return is_fully_comm(word, braid_rels)
+
         def reduced_word_reverse_iterator(self):
             """
             Return a reverse iterator on a reduced word for ``self``.
@@ -1589,6 +1660,31 @@ def reduced_word(self):
             result = list(self.reduced_word_reverse_iterator())
             return list(reversed(result))
 
+        def reduced_words_iter(self):
+            r"""
+            Iterate over all reduced words for ``self``.
+
+            See :meth:`reduced_word` for the definition of a reduced
+            word.
+
+            The algorithm uses the Matsumoto property that any two
+            reduced expressions are related by braid relations, see
+            Theorem 3.3.1(ii) in [BB2005]_.
+
+            .. SEEALSO::
+
+                :meth:`braid_orbit_iter`
+
+            EXAMPLES::
+
+                sage: W = CoxeterGroups().example()
+                sage: s = W.simple_reflections()
+                sage: w = s[0] * s[2]
+                sage: sorted(w.reduced_words_iter())                                    # optional - sage.combinat
+                [[0, 2], [2, 0]]
+            """
+            return self.parent().braid_orbit_iter(self.reduced_word())
+
         def reduced_words(self):
             r"""
             Return all reduced words for ``self``.
@@ -1648,7 +1744,7 @@ def reduced_words(self):
                 :meth:`.reduced_word`, :meth:`.reduced_word_reverse_iterator`,
                 :meth:`length`, :meth:`reduced_word_graph`
             """
-            return self.parent().braid_orbit(self.reduced_word())
+            return list(self.reduced_words_iter())
 
         def support(self):
             r"""

src/sage/combinat/affine_permutation.py

@@ -217,7 +217,7 @@ def __call__(self, i):
         """
         return self.value(i)
 
-    def is_i_grassmannian(self, i=0, side="right"):
+    def is_i_grassmannian(self, i=0, side="right") -> bool:
         r"""
         Test whether ``self`` is `i`-grassmannian, i.e., either is the
         identity or has ``i`` as the sole descent.
@@ -277,7 +277,7 @@ def lower_covers(self,side="right"):
         S = self.descents(side)
         return [self.apply_simple_reflection(i, side) for i in S]
 
-    def is_one(self):
+    def is_one(self) -> bool:
         r"""
         Tests whether the affine permutation is the identity.
 
@@ -876,18 +876,21 @@ def to_lehmer_code(self, typ='decreasing', side='right'):
                         code[i] += (b-i) // (self.k+1) + 1
         return Composition(code)
 
-    def is_fully_commutative(self):
+    def is_fully_commutative(self) -> bool:
         r"""
-        Determine whether ``self`` is fully commutative, i.e., has no
-        reduced words with a braid.
+        Determine whether ``self`` is fully commutative.
+
+        This means that it has no reduced word with a braid.
+
+        This uses a specific algorithm.
 
         EXAMPLES::
 
             sage: A = AffinePermutationGroup(['A',7,1])
-            sage: p=A([3, -1, 0, 6, 5, 4, 10, 9])
+            sage: p = A([3, -1, 0, 6, 5, 4, 10, 9])
             sage: p.is_fully_commutative()
             False
-            sage: q=A([-3, -2, 0, 7, 9, 2, 11, 12])
+            sage: q = A([-3, -2, 0, 7, 9, 2, 11, 12])
             sage: q.is_fully_commutative()
             True
         """
@@ -900,14 +903,12 @@ def is_fully_commutative(self):
             if c[i] > 0:
                 if firstnonzero is None:
                     firstnonzero = i
-                if m != -1 and c[i] - (i-m) >= c[m]:
+                if m != -1 and c[i] - (i - m) >= c[m]:
                     return False
                 m = i
-        #now check m (the last non-zero) against firstnonzero.
-        d = self.n-(m-firstnonzero)
-        if c[firstnonzero]-d >= c[m]:
-            return False
-        return True
+        # now check m (the last non-zero) against first non-zero.
+        d = self.n - (m - firstnonzero)
+        return not c[firstnonzero] - d >= c[m]
 
     def to_bounded_partition(self, typ='decreasing', side='right'):
         r"""

src/sage/combinat/fully_commutative_elements.py

@@ -153,44 +153,23 @@ def is_fully_commutative(self):
             sage: x = FC.element_class(FC, [1, 2, 1], check=False); x.is_fully_commutative()
             False
         """
-        matrix = self.parent().coxeter_group().coxeter_matrix()
-        w = tuple(self)
-
-        # The following function detects 'braid' words.
-        def contains_long_braid(w):
-            for i in range(len(w) - 2):
-                a = w[i]
-                b = w[i + 1]
-                m = matrix[a, b]
-                if m > 2 and i + m <= len(w):
-                    ab_braid = (a, b) * (m // 2) + ((a,) if m % 2 else ())
-                    if w[i:i + m] == ab_braid:
-                        return True
-            return False
+        word = list(self)
+        from sage.combinat.root_system.braid_orbit import is_fully_commutative as is_fully_comm
 
-        # The following function applies a commutation relation on a word.
-        def commute_once(word, i):
-            return word[:i] + (word[i + 1], word[i]) + word[i + 2:]
+        group = self.parent().coxeter_group()
+        braid_rels = group.braid_relations()
+        I = group.index_set()
 
-        # A word is the reduced word of an FC element iff no sequence of
-        # commutation relations on it yields a word with a 'braid' word:
-        if contains_long_braid(w):
-            return False
-        else:
-            l, checked, queue = len(w), {w}, deque([w])
-            while queue:
-                word = queue.pop()
-                for i in range(l - 1):
-                    a, b = word[i], word[i + 1]
-                    if matrix[a, b] == 2:
-                        new_word = commute_once(word, i)
-                        if new_word not in checked:
-                            if contains_long_braid(new_word):
-                                return False
-                            else:
-                                checked.add(new_word)
-                                queue.appendleft(new_word)
-            return True
+        from sage.rings.integer_ring import ZZ
+        be_careful = any(i not in ZZ for i in I)
+
+        if be_careful:
+            Iinv = {i: j for j, i in enumerate(I)}
+            word = [Iinv[i] for i in word]
+            braid_rels = [[[Iinv[i] for i in l],
+                           [Iinv[i] for i in r]] for l, r in braid_rels]
+
+        return is_fully_comm(word, braid_rels)
 
     # Representing FC elements: Heaps
     def heap(self, **kargs):

src/sage/combinat/root_system/braid_orbit.pyx

@@ -74,6 +74,61 @@ cpdef set BraidOrbit(list word, list rels):
     return words
 
 
+cpdef bint is_fully_commutative(list word, list rels):
+    r"""
+    Check if the braid orbit of ``word`` is using a braid relation.
+
+    INPUT:
+
+    - ``word`` -- list of integers
+
+    - ``rels`` -- list of pairs ``(A, B)``, where ``A`` and ``B`` are
+      lists of integers the same length
+
+    EXAMPLES::
+
+        sage: from sage.combinat.root_system.braid_orbit import is_fully_commutative
+        sage: rels = [[[2, 1, 2], [1, 2, 1]], [[3, 1], [1, 3]], [[3, 2, 3], [2, 3, 2]]]
+        sage: word = [1,2,1,3,2,1]
+        sage: is_fully_commutative(word, rels)
+        False
+        sage: word = [1,2,3]
+        sage: is_fully_commutative(word, rels)
+        True
+    """
+    cdef int i, l, rel_l, loop_ind, list_len
+    cdef tuple left, right, test_word, new_word
+    cdef list rel
+
+    l = len(word)
+    cdef set words = set( [tuple(word)] )
+    cdef list test_words = [ tuple(word) ]
+
+    rels = rels + [[b, a] for a, b in rels]
+    rels = [[tuple(a), tuple(b), len(a)] for a, b in rels]
+
+    loop_ind = 0
+    list_len = 1
+    while loop_ind < list_len:
+        sig_check()
+        test_word = <tuple> test_words[loop_ind]
+        loop_ind += 1
+        for rel in rels:
+            left = <tuple> rel[0]
+            right = <tuple> rel[1]
+            rel_l = <int> rel[2]
+            for i in range(l-rel_l+1):
+                if pattern_match(test_word, i, left, rel_l):
+                    if rel_l > 2:
+                        return False
+                    new_word = test_word[:i] + right + test_word[i+rel_l:]
+                    if new_word not in words:
+                        words.add(new_word)
+                        test_words.append(new_word)
+                        list_len += 1
+    return True
+
+
 cdef inline bint pattern_match(tuple L, int i, tuple X, int l):
     r"""
     Return ``True`` if ``L[i:i+l] == X``.