From 5a4d1c6f54f63e605232902b3f21196ad3f6bd3f Mon Sep 17 00:00:00 2001
From: Luke Kenneth Casson Leighton <lkcl@lkcl.net>
Date: Sun, 6 Mar 2022 09:36:06 +0000
Subject: [PATCH] add quick demo gf2.py

---
 openpower/sv/bitmanip.mdwn |  3 ++
 openpower/sv/gf2.py        | 95 ++++++++++++++++++++++++++++++++++++++
 2 files changed, 98 insertions(+)
 create mode 100644 openpower/sv/gf2.py

diff --git a/openpower/sv/bitmanip.mdwn b/openpower/sv/bitmanip.mdwn
index d5ebbd6ea..0e3de5985 100644
--- a/openpower/sv/bitmanip.mdwn
+++ b/openpower/sv/bitmanip.mdwn
@@ -663,6 +663,9 @@ GFMOD is a pseudo-op where RA=0
 
 ## gf invert
 
+note this function is incorrect / incomplete, as-is
+https://stackoverflow.com/questions/45442396/
+
 ```
 def gf_degree(a) :
   res = 0
diff --git a/openpower/sv/gf2.py b/openpower/sv/gf2.py
new file mode 100644
index 000000000..8a52a2b99
--- /dev/null
+++ b/openpower/sv/gf2.py
@@ -0,0 +1,95 @@
+from functools import reduce
+
+def gf_degree(a) :
+  res = 0
+  a >>= 1
+  while (a != 0) :
+    a >>= 1;
+    res += 1;
+  return res
+
+# constants used in the multGF2 function
+mask1 = mask2 = polyred = None
+
+def setGF2(irPoly):
+    """Define parameters of binary finite field GF(2^m)/g(x)
+       - irPoly: coefficients of irreducible polynomial g(x)
+    """
+    # degree: extension degree of binary field
+    degree = gf_degree(irPoly)
+
+    def i2P(sInt):
+        """Convert an integer into a polynomial"""
+        return [(sInt >> i) & 1
+                for i in reversed(range(sInt.bit_length()))]    
+    
+    global mask1, mask2, polyred
+    mask1 = mask2 = 1 << degree
+    mask2 -= 1
+    polyred = reduce(lambda x, y: (x << 1) + y, i2P(irPoly)[1:])
+        
+
+def multGF2(p1, p2):
+    """Multiply two polynomials in GF(2^m)/g(x)"""
+    p = 0
+    while p2:
+        if p2 & 1:
+            p ^= p1
+        p1 <<= 1
+        if p1 & mask1:
+            p1 ^= polyred
+        p2 >>= 1
+    return p & mask2
+
+
+def divmodGF2(f, v):
+    fDegree, vDegree = gf_degree(f), gf_degree(v)
+    res, rem = 0, f
+    i = fDegree
+    mask = 1 << i
+    while (i >= vDegree):
+        if (mask & rem):
+            res ^= (1 << (i - vDegree))
+            rem ^= ( v << (i - vDegree))
+        i -= 1
+        mask >>= 1
+    return (res, rem)
+
+
+def xgcd(a, b):
+    """return (g, x, y) such that a*x + b*y = g = gcd(a, b)"""
+    x0, x1, y0, y1 = 0, 1, 1, 0
+    while a != 0:
+        (q, a), b = divmod(b, a), a
+        y0, y1 = y1, y0 - q * y1
+        x0, x1 = x1, x0 - q * x1
+    return b, x0, y0
+
+
+if __name__ == "__main__":
+  
+    # Define binary field GF(2^3)/x^3 + x + 1
+    setGF2(0b1011) # degree 3
+
+    # Evaluate the product (x^2 + x + 1)(x^2 + 1)
+    x = multGF2(0b111, 0b101)
+    print("%02x" % x)
+    
+    # Define binary field GF(2^8)/x^8 + x^4 + x^3 + x + 1
+    # (used in the Advanced Encryption Standard-AES)
+    setGF2(0b100011011) # degree 8
+    
+    # Evaluate the product (x^7 + x^2)(x^7 + x + 1)
+    x = 0b10000100
+    y = 0b10000011
+    z = multGF2(x, y)
+    print("%02x * %02x = %02x" % (x, y, z))
+
+    # divide z by y into result/remainder
+    res, rem = divmodGF2(z, y)
+    print("%02x / %02x = (%02x, %02x)" % (z, y, res, rem))
+
+    # reconstruct x by multiplying divided result by y and adding the remainder
+    x1 = multGF2(res, y)
+    print("%02x == %02x" % (z, x1 ^ rem))
+
-- 
2.30.2