From 045eab7e7622d81549b5c9e65d02929352184c32 Mon Sep 17 00:00:00 2001
From: Erich Neun <erich-9@gmx.de>
Date: Mon, 1 Oct 2018 20:03:17 +0200
Subject: [PATCH] Add HighLevelGroebnerBasis as an alternative to
 GBNPGroebnerBasis.

  This is a high-level implementation and therefore certainly not as fast
  as GBNPGroebnerBasis. However, GBNPGroebnerBasis sometimes results in
  unpredictable bugs such that having an alternative seems reasonable.
---
 init.g              |   1 +
 lib/gbhighlevel.gd  |   4 ++
 lib/gbhighlevel.gi  | 113 ++++++++++++++++++++++++++++++++++++++++++++
 lib/pathalgideal.gi |   2 +-
 read.g              |   1 +
 5 files changed, 120 insertions(+), 1 deletion(-)
 create mode 100644 lib/gbhighlevel.gd
 create mode 100644 lib/gbhighlevel.gi

diff --git a/init.g b/init.g
index 0de4557..eef74d7 100644
--- a/init.g
+++ b/init.g
@@ -23,6 +23,7 @@ ReadPackage("QPA", "lib/projpathalgmodule.gd");
 ReadPackage("QPA", "lib/moduledecomp.gd");
 ReadPackage("QPA", "lib/idempotent.gd");
 ReadPackage("QPA", "lib/gbnp.gd");
+ReadPackage("QPA", "lib/gbhighlevel.gd");
 ReadPackage("QPA", "lib/moduleprojres.gd");
 ReadPackage("QPA", "lib/pathalgtensor.gd");
 ReadPackage("QPA", "lib/pathalgpredef.gd");
diff --git a/lib/gbhighlevel.gd b/lib/gbhighlevel.gd
new file mode 100644
index 0000000..ab2781b
--- /dev/null
+++ b/lib/gbhighlevel.gd
@@ -0,0 +1,4 @@
+DeclareOperation( "HighLevelGroebnerBasis", [ IsList, IsPathAlgebra ] );
+DeclareOperation( "RemainderOfDivision", [ IsElementOfMagmaRingModuloRelations, IsList, IsPathAlgebra ] );
+DeclareOperation( "LeftmostOccurrence", [ IsList, IsList ] );
+DeclareSynonym( "TipWalk", x -> WalkOfPath(TipMonomial(x)) );
diff --git a/lib/gbhighlevel.gi b/lib/gbhighlevel.gi
new file mode 100644
index 0000000..a77132b
--- /dev/null
+++ b/lib/gbhighlevel.gi
@@ -0,0 +1,113 @@
+InstallMethod( HighLevelGroebnerBasis,
+  "compute a complete Groebner Basis",
+  [ IsList, IsPathAlgebra ],
+  function(els, A)
+    local gb, el, el_tip,
+          n, i, j, x, y, k, l, r, b, c,
+          overlap, remainder;
+
+    if not QPA_InArrowIdeal(els, A) then
+      Error("elements do not belong to the arrow ideal of the path algebra");
+    fi;
+
+    els := MakeUniform(els);
+
+    gb := [];
+
+    while Length(els) > 0 do
+      for el in els do
+        el_tip := Tip(el);
+        Add(gb, el/TipCoefficient(el_tip));
+      od;
+
+      n := Length(gb);
+      els := [];
+
+      for i in [1..n] do
+        x := TipWalk(gb[i]);
+        k := Length(x);
+
+        for j in [1..n] do
+          y := TipWalk(gb[j]);
+          l := Length(y);
+
+          for r in [Maximum(0, k-l)..k-1] do
+            if x{[r+1..k]} = y{[1..k-r]} then
+              b := x{[1..r]};
+              c := y{[k-r+1..l]};
+
+              overlap := gb[i]*Product(c, One(A)) - Product(b, One(A))*gb[j];
+              remainder := RemainderOfDivision(overlap, gb, A);
+
+              if not IsZero(remainder) then
+                AddSet(els, remainder);
+              fi;
+            fi;
+          od;
+        od;
+      od;
+    od;
+
+    return gb;
+  end
+);
+
+
+InstallMethod( RemainderOfDivision,
+  "for a path-algebra element and a list of monic path-algebra elements",
+  [ IsElementOfMagmaRingModuloRelations, IsList, IsPathAlgebra ],
+  function(y, X, A)
+    local r, n, y_tip, y_wtip, divided, i, p, u, v;
+
+    r := Zero(A);
+    n := Length(X);
+
+    while not IsZero(y) do
+      y_tip := Tip(y);
+      y_wtip := TipWalk(y_tip);
+
+      divided := false;
+
+      for i in [1..n] do
+        p := LeftmostOccurrence(y_wtip, TipWalk(X[i]));
+
+        if p <> fail then
+          u := Product(y_wtip{[1..p[1]-1]}, One(A));
+          v := Product(y_wtip{[p[2]+1..Length(y_wtip)]}, One(A));
+
+          y := y - TipCoefficient(y_tip) * u*X[i]*v;
+
+          divided := true;
+          break;
+        fi;
+      od;
+
+      if not divided then
+        r := r + y_tip;
+        y := y - y_tip;
+      fi;
+    od;
+
+    return r;
+  end
+);
+
+
+InstallMethod( LeftmostOccurrence,
+  "find second list as sublist of first list",
+  [ IsList, IsList ],
+  function(b, c)
+    local lb, lc, i;
+
+    lb := Length(b);
+    lc := Length(c);
+
+    for i in [1..lb-lc+1] do
+      if b{[i..i+lc-1]} = c then
+        return [i, i+lc-1];
+      fi;
+    od;
+
+    return fail;
+  end
+);
diff --git a/lib/pathalgideal.gi b/lib/pathalgideal.gi
index 92fb0e6..96147fc 100644
--- a/lib/pathalgideal.gi
+++ b/lib/pathalgideal.gi
@@ -65,7 +65,7 @@ InstallTrueMethod( IsIdealInPathAlgebra,
 #O  \in(<elt>, <I>)
 ##
 ##  This function performs a membership test for an ideal in path algebra using 
-##  Groebner Bases machinery.
+##  (by default GBNP) Groebner Bases machinery.
 ##  It returns true if <elt> is a member of an ideal <I>, false otherwise.
 ##  For the efficiency reasons, it computes Groebner basis
 ##  for <I> only if it has not been computed. Similarly, it performs
diff --git a/read.g b/read.g
index 5d4507e..fe3eb62 100644
--- a/read.g
+++ b/read.g
@@ -23,6 +23,7 @@ ReadPackage("QPA", "lib/projpathalgmodule.gi");
 ReadPackage("QPA", "lib/moduledecomp.gi");
 ReadPackage("QPA", "lib/idempotent.gi");
 ReadPackage("QPA", "lib/gbnp.gi");
+ReadPackage("QPA", "lib/gbhighlevel.gi");
 ReadPackage("QPA", "lib/moduleprojres.gi");
 ReadPackage("QPA", "lib/pathalgtensor.gi");
 ReadPackage("QPA", "lib/pathalgpredef.gi");