Add ToddCoxeterBand method

doc/freeband.xml

@@ -190,3 +190,23 @@ gap> ContentOfFreeBandElementCollection([x, y]);
     </Description>
   </ManSection>
 <#/GAPDoc>
+
+<#GAPDoc Label="ToddCoxeterBand">
+  <ManSection>
+    <Oper Name="ToddCoxeterBand"  Arg="n, R"/>
+    <Description>
+      This operation takes band presentation, where <A>n</A> is the size
+      of alphabet <C>A = [1 .. n]</C> and <A>R</A> is a list of lists of
+      words over <C>A</C>, representing the relations. It computes the
+      band defined by this band presentation via a band-specific version
+      of the Todd-Coxeter algorithm. If <A>R</A> is the empty list, then
+      the free band over <C>A</C> is computed.
+
+      <!--
+      TODO: I think we should write a little more about the algorithm here
+      once we've finalised it. We should also add some examples, once we've
+      decided on exactly what the output should be.
+      -->
+    </Description>
+  </ManSection>
+<#/GAPDoc>

doc/z-chap10.xml

@@ -154,6 +154,7 @@ x1x2x2^-1x1^-1x1x2]]></Example>
     <#Include Label = "IsFreeBandElementCollection">
     <#Include Label = "IsFreeBandSubsemigroup">
     <#Include Label = "ContentOfFreeBandElement">
+    <#Include Label = "ToddCoxeterBand">
 
   </Section>
 

gap/fp/freeband.gd

@@ -18,3 +18,5 @@ DeclareGlobalFunction("FreeBand");
 DeclareAttribute("ContentOfFreeBandElement", IsFreeBandElement);
 DeclareAttribute("ContentOfFreeBandElementCollection",
                  IsFreeBandElementCollection);
+
+DeclareOperation("ToddCoxeterBand", [IsPosInt, IsList]);

gap/fp/freeband.gi

@@ -550,3 +550,336 @@ function(x, hashlen)
   return rec(func := SEMIGROUPS.HashFunctionForFreeBandElements,
              data := hashlen);
 end);
+
+SEMIGROUPS.PrefixTupleOfFreeBandElement := function(word, n)
+  local lookup, distinct, L, i;
+  # word is a list of pos ints. n is the content size.
+  # Returns a list: first element is largest prefix of word which has content
+  # size n - 1.
+  # Second element is nth character to make its first appearance.
+
+  # Start with some argument checks and easy outputs (shouldn't be used).
+  if not (IsList(word) and not IsEmpty(word)) then
+    ErrorNoReturn("expected a non-empty list of pos ints as first argument");
+  fi;
+  if not IsPosInt(n) then
+    ErrorNoReturn("expected a positive integer as second argument");
+  fi;
+  if n > Length(word) then
+    ErrorNoReturn("second argument (content size) must be smaller than length ",
+                  "of first argument (word)");
+  fi;
+  if n = 1 then
+    return [[], word[1]];
+  fi;
+
+  # If the arguments are nontrivial, start scanning the word.
+  lookup          := [];
+  lookup[word[1]] := 1;
+  distinct        := 1;
+  L               := Length(word);
+  for i in [2 .. L] do
+    if not IsBound(lookup[word[i]]) then
+      distinct := distinct + 1;
+      if distinct = n then
+        # in this case the number of distinct letters has hit the max
+        return [word{[1 .. i - 1]}, word[i]];
+      fi;
+      lookup[word[i]] := distinct;
+    fi;
+  od;
+
+  # if distinct never reached n, then n was larger than content of word.
+  ErrorNoReturn("second argument (content size) must be at most equal to the ",
+                "number of distinct entries in the first argument");
+end;
+
+SEMIGROUPS.ShortCanonicalFormOfFreeBandElement := function(word)
+  local L, lookup, distinct, l1, l2, preftup, sufftup, pref, suff,
+  n1, n2, n, char, i;
+  # Returns shortest representation of the free band equivalence class of
+  # the input.
+  if not IsList(word) then
+    ErrorNoReturn("Expected a list of pos ints as input");
+  fi;
+  for char in word do
+    if not IsPosInt(char) then
+      ErrorNoReturn("Expected a list of pos ints as input");
+    fi;
+  od;
+
+  # Check if word is short, in which case output is easy.
+  L := Length(word);
+  if L = 0 then
+    return [];
+  elif L = 1 then
+    return word;
+  fi;
+
+  # If word is longer: first we need size of the alph.
+  lookup          := [];
+  lookup[word[1]] := 1;
+  distinct        := 1;
+  for i in [2 .. L] do
+    if not IsBound(lookup[word[i]]) then
+      distinct        := distinct + 1;
+      lookup[word[i]] := distinct;
+    fi;
+  od;
+
+  # If the  size of the alph is 1 or 2, output is easy.
+  if distinct = 1 then
+    # only one letter. Delete copies.
+    return [word[1]];
+  fi;
+  if distinct = 2 then
+    # then only the first and last letters matter.
+    if word[1] <> word[L] then
+      return [word[1], word[L]];
+    fi;
+    # otherwise the first and last are same. Need to know middle.
+    l1 := word[1];
+    for char in word do
+      if char <> l1 then
+        l2 := char;
+        return [l1, l2, l1];
+      fi;
+    od;
+  fi;
+
+  # If we've made it this far then the content is at least of size 3.
+  # Need to find prefix and suffix. Note suffix will be reversed in the output.
+  preftup := SEMIGROUPS.PrefixTupleOfFreeBandElement(word, distinct);
+  sufftup := SEMIGROUPS.PrefixTupleOfFreeBandElement(Reversed(word), distinct);
+
+  # Turn the two tuples into the two halves of the output.
+  pref := SEMIGROUPS.ShortCanonicalFormOfFreeBandElement(preftup[1]);
+  Add(pref, preftup[2]);
+
+  suff := SEMIGROUPS.ShortCanonicalFormOfFreeBandElement(sufftup[1]);
+  Add(suff, sufftup[2]);
+  suff := Reversed(suff);
+
+  # Now see if any cancellations are possible between suffixes of pref and
+  # prefixes of suff. Give priority to longest possible.
+  n1 := Length(pref);
+  n2 := Length(suff);
+  n  := Minimum(n1, n2);
+  for i in [0 .. n - 1] do
+    if pref{[n1 - n + i + 1 .. n1]} = suff{[1 .. n - i]} then
+      # a cancellation is possible. Cancel and return.
+      return Concatenation(pref, suff{[n - i + 1 .. n2]});
+    fi;
+  od;
+
+  # If we've made it this far then no cancellations are possible.
+  return Concatenation(pref, suff);
+end;
+
+InstallMethod(ToddCoxeterBand, "for a pos int and list of lists of words",
+[IsPosInt, IsList],
+function(N, R)
+  local new_coset, tauf, canon, push_relation, process_coincidences,
+  A, k, active_cosets, table, coincidences, words, n, word,
+  pair, char, coset, i, tau;
+
+  new_coset := function(coset, char)
+    local new_word, target, cosetm, charm, pword;
+    # intelligently creates a new coset ONLY if (non-forced) tau(coset, char)
+    # is not defined.
+    if table[coset][char] = 0 then
+      new_word := canon(Concatenation(words[coset], [char]));
+      target   := tau(1, new_word);
+
+      if target = 0 then
+        # in this case following new_word from empty word does not lead us the
+        # full way and we need to define more cosets.
+        # Need to follow word again to see how far we got before undefined.
+        cosetm := 1;
+        pword  := [];  # partial word, add letters to it each time
+        for charm in new_word do
+          # extend partial word
+          Add(pword, charm);
+          if table[cosetm][charm] = 0 then
+            # edge is undefined, define a new one.
+            table[cosetm][charm] := k;
+            active_cosets[k]     := true;
+            Add(table, ListWithIdenticalEntries(Length(A), 0));
+            Add(words, ShallowCopy(pword));
+            # we need to re-canonicalise pword at the moment because canon does
+            # not always output the shortlex-least word.
+            # TODO remove this comment soon
+            k := k + 1;
+          fi;
+          cosetm := table[cosetm][charm];
+        od;
+
+        # now we've defined all the intermediate words, get the original
+        # request to point in the right place.
+        table[coset][char] := k - 1;  # k had been incremented 1 too many
+
+      else
+        # in this case following new_word led us somewhere and we should
+        # point there
+        table[coset][char] := target;
+      fi;
+    fi;
+  end;
+
+  tauf := function(coset, word)
+    local char;
+    # forced tau. This creates new cosets as necessary.
+    if Length(word) = 0 then
+      return coset;
+    fi;
+    for char in word do
+      if table[coset][char] = 0 then
+        # if the product is undefined, define it, and start coset back up
+        # at the newly defined value (k-1).
+        new_coset(coset, char);
+      fi;
+      coset := table[coset][char];
+    od;
+    return coset;
+  end;
+
+  tau := function(coset, word)
+    local char;
+    # non-forced tau, checks whether you can get the whole way.
+    # for use in new-coset.
+    if Length(word) = 0 then
+      return coset;
+    fi;
+    for char in word do
+      if table[coset][char] = 0 then
+        return 0;
+      fi;
+      coset := table[coset][char];
+    od;
+    return coset;
+  end;
+
+  canon := function(word)
+    # expresses a word in free band-canonical form.
+    # NOTE: it is essential to the validity of this algorithm that the canonical
+    # form returned is shortlex minimal. Otherwise new_coset doesn't work
+    # properly.
+    return SEMIGROUPS.ShortCanonicalFormOfFreeBandElement(word);
+  end;
+
+  push_relation := function(coset, u, v)
+    local ut, vt;
+    ut := tauf(coset, u);
+    vt := tauf(coset, v);
+    if ut <> vt then
+      Add(coincidences, [ut, vt]);
+    fi;
+  end;
+
+  process_coincidences := function()
+    # changed to depth-first.
+    local i, j, char, coset, pair, current, counter;
+    if Length(coincidences) = 0 then
+      return;
+    fi;
+    while Length(coincidences) <> 0 do
+      # current := Length(coincidences);
+      current := 1;
+      i       := Minimum(coincidences[current]);
+      j       := Maximum(coincidences[current]);
+      if i = j then
+      fi;
+      counter := 0;
+      if i <> j then
+        for char in A do
+          if table[j][char] <> 0 then
+            if table[i][char] = 0 then
+              table[i][char] := table[j][char];
+            elif table[i][char] <> 0 then
+              counter := counter + 1;
+              Add(coincidences, [table[i][char], table[j][char]]);
+            fi;
+          fi;
+        od;
+        # for coset in ListBlist([1 .. k - 1], active_cosets) do
+        for coset in [1 .. k - 1] do
+          for char in A do
+            if table[coset][char] = j then
+              table[coset][char] := i;
+            fi;
+          od;
+        od;
+        for pair in coincidences do
+          if pair[1] = j then
+            pair[1] := i;
+          fi;
+          if pair[2] = j then
+            pair[2] := i;
+          fi;
+        od;
+        active_cosets[j] := false;
+      fi;
+      Remove(coincidences, current);
+      # Unbind(parents[j]);
+      # Unbind(edges[j]);
+    od;
+  end;
+
+  A             := [1 .. N];
+  # F             := FreeBand(N);               # obsolete since new canon func
+  # G             := GeneratorsOfSemigroup(F);  # obsolete since new canon func
+  k             := 2;
+  active_cosets := [true];
+  table         := [[]];
+  coincidences  := [];
+  words         := [[]];
+  for char in A do
+    table[1][char] := 0;
+  od;
+  n := 0;
+  repeat
+
+    n  := n + 1;
+
+    # only do anything if the current coset is active
+    if active_cosets[n] then
+
+      # populate the current line of the table with new cosets if need be
+      for char in A do
+        new_coset(n, char);
+      od;
+
+      # push the coset n through every explicit relation
+      for pair in R do
+        push_relation(n, pair[1], pair[2]);
+      od;
+
+      # push the current coset through every known implicit relation
+      for word in words do
+        pair := [Concatenation(word, word), word];
+        push_relation(n, pair[1], pair[2]);  # word is already canonical
+      od;
+
+      # push every previous coset through the current implicit relation
+      word := words[n];
+      pair := [Concatenation(word, word), word];
+      for i in ListBlist([1 .. n - 1], active_cosets{[1 .. n - 1]}) do
+        push_relation(i, pair[1], pair[2]);
+      od;
+
+    fi;
+
+    process_coincidences();
+
+  until n = k - 1;
+
+  for pair in R do
+    if Length(pair[1]) = 0 or Length(pair[2]) = 0 then
+      # then one of these must be a monoid presentation.
+      return Length(ListBlist([1 .. k - 1], active_cosets));
+    fi;
+  od;
+
+  # if no relations have the empty word then this is not a monoid presentation.
+  return Length(ListBlist([1 .. k - 1], active_cosets)) - 1;
+end);