Add 'src/nmutil/openpower_sv_bitmanip_in_wiki/' from commit '7c40cd7bd4f5a42684c223c6...

git-subtree-dir: src/nmutil/openpower_sv_bitmanip_in_wiki
git-subtree-mainline: 7a3ab7ecce390b2003ff61befca59cd4ad6b7428
git-subtree-split: 7c40cd7bd4f5a42684c223c6b6289f392a4dec01

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/__init__.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/__init__.py

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/cldivrem.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..33a528da422e69243edfaac62f21411c9273ef0b
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/cldivrem.py
@@@ -1,0 -1,0 +1,30 @@@
++from .log2 import floor_log2
++
++
++def cldivrem(n, d, width):
++    """ Carry-less Division and Remainder.
++        `n` and `d` are integers, `width` is the number of bits needed to hold
++        each input/output.
++        Returns a tuple `q, r` of the quotient and remainder.
++    """
++    assert d != 0, "TODO: decide what happens on division by zero"
++    assert 0 <= n < 1 << width, f"bad n (doesn't fit in {width}-bit uint)"
++    assert 0 <= d < 1 << width, f"bad d (doesn't fit in {width}-bit uint)"
++    r = n
++    q = 0
++    d <<= width
++    for _ in range(width):
++        d >>= 1
++        q <<= 1
++        if degree(d) == degree(r):
++            r ^= d
++            q |= 1
++    return q, r
++
++
++def degree(v):
++    """the degree of the GF(2) polynomial `v`. `v` is a non-negative integer.
++    """
++    if v == 0:
++        return -1
++    return floor_log2(v)

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/clmul.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..7451c3cc8c77dbb5792535ff8d4aa9c80c4dac32
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/clmul.py
@@@ -1,0 -1,0 +1,8 @@@
++def clmul(a, b):
++    x = 0
++    i = 0
++    while b >> i != 0:
++        if (b >> i) & 1:
++            x ^= a << i
++        i += 1
++    return x

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/clmulh.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..b170aca3b798ac8d6180fec1a5486b3664732239
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/clmulh.py
@@@ -1,0 -1,0 +1,5 @@@
++from .clmul import clmul
++
++
++def clmulh(a, b, XLEN):
++    return clmul(a, b) >> XLEN

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/clmulr.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..5b155b5c59bbdbe1b9531c6f49e4d087eef33cce
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/clmulr.py
@@@ -1,0 -1,0 +1,5 @@@
++from .clmul import clmul
++
++
++def clmulh(a, b, XLEN):
++    return clmul(a, b) >> (XLEN - 1)

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/decode_reducing_polynomial.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..ac6623b5f04c5789d3f84fd5dd1353f4c2fa0533
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/decode_reducing_polynomial.py
@@@ -1,0 -1,0 +1,19 @@@
++from .state import ST
++
++
++def decode_reducing_polynomial():
++    """returns the decoded reducing polynomial as an integer.
++        Note: the returned integer is `XLEN + 1` bits wide.
++    """
++    v = ST.GFBREDPOLY & ((1 << ST.XLEN) - 1)  # mask to XLEN bits
++    if v == 0 or v == 2:  # GF(2)
++        return 0b10  # degree = 1, poly = x
++    if (v & 1) == 0:
++        # all reducing polynomials of degree > 1 must have the LSB set,
++        # because they must be irreducible polynomials (meaning they
++        # can't be factored), if the LSB was clear, then they would
++        # have `x` as a factor. Therefore, we can reuse the LSB clear
++        # to instead mean the polynomial has degree XLEN.
++        v |= 1 << ST.XLEN
++        v |= 1  # LSB must be set
++    return v

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfbinv.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..57aa5a86e9175a69938db55300c4f1ae5e098224
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfbinv.py
@@@ -1,0 -1,0 +1,40 @@@
++from .decode_reducing_polynomial import decode_reducing_polynomial
++from .cldivrem import degree
++
++
++def gfbinv(a):
++    """compute the GF(2^m) inverse of `a`."""
++    # Derived from Algorithm 3, from [7] in:
++    # https://ftp.libre-soc.org/ARITH18_Kobayashi.pdf
++
++    s = decode_reducing_polynomial()
++    m = degree(s)
++    assert a >> m == 0, "`a` is out-of-range"
++    r = a
++    v = 0
++    u = 1
++    delta = 0
++
++    for _ in range(2 * m):
++        # could use count-leading-zeros here to skip ahead
++        if r >> m == 0:  # if the MSB of `r` is zero
++            r <<= 1
++            u <<= 1
++            delta += 1
++        else:
++            if s >> m != 0:  # if the MSB of `s` isn't zero
++                s ^= r
++                v ^= u
++            s <<= 1
++            if delta == 0:
++                r, s = s, r  # swap r and s
++                u, v = v << 1, u  # shift v and swap
++                delta = 1
++            else:
++                u >>= 1
++                delta -= 1
++    if a == 0:
++        # we specifically choose 0 as the result of inverting 0, rather than an
++        # error or undefined, since that's what Rijndael needs.
++        return 0
++    return u

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfbmadd.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..df826c4f0edf9a651e6bb3fe8936fe4162aff681
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfbmadd.py
@@@ -1,0 -1,0 +1,11 @@@
++from .state import ST
++from .decode_reducing_polynomial import decode_reducing_polynomial
++from .clmul import clmul
++from .cldivrem import cldivrem
++
++
++def gfbmadd(a, b, c):
++    v = clmul(a, b) ^ c
++    red_poly = decode_reducing_polynomial()
++    q, r = cldivrem(v, red_poly, width=ST.XLEN + 1)
++    return r

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfbmul.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..6295a29d8bd086972da9ca7deb324cfdf4d26379
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfbmul.py
@@@ -1,0 -1,0 +1,11 @@@
++from .state import ST
++from .decode_reducing_polynomial import decode_reducing_polynomial
++from .clmul import clmul
++from .cldivrem import cldivrem
++
++
++def gfbmul(a, b):
++    product = clmul(a, b)
++    red_poly = decode_reducing_polynomial()
++    q, r = cldivrem(product, red_poly, width=ST.XLEN + 1)
++    return r

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfbredpoly.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..18129679494d9cc7d663cc2c992ea509fa10ca87
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfbredpoly.py
@@@ -1,0 -1,0 +1,6 @@@
++from .state import ST
++
++
++def gfbredpoly(immed):
++    # TODO: figure out how `immed` should be encoded
++    ST.GFBREDPOLY = immed

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfpadd.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..f00bf8c694aa52de8991331793f01960858cd2e6
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfpadd.py
@@@ -1,0 -1,0 +1,5 @@@
++from .state import ST
++
++
++def gfpadd(a, b):
++    return (a + b) % ST.GFPRIME

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfpinv.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..45b6dbbdd4cb19bdbfa0538fb9dfcb79981655fe
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfpinv.py
@@@ -1,0 -1,0 +1,51 @@@
++from .state import ST
++
++
++def gfpinv(a):
++    # based on Algorithm ExtEucdInv from:
++    # https://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.90.5233&rep=rep1&type=pdf
++    p = ST.GFPRIME
++    assert p >= 2, "GFPRIME isn't a prime"
++    assert a != 0, "TODO: decide what happens for division by zero"
++    assert isinstance(a, int) and 0 < a < p, "a out of range"
++    if p == 2:
++        return 1  # the only value possible
++
++    u = p
++    v = a
++    r = 0
++    s = 1
++    while v > 0:
++        # implementations could use count-zeros on
++        # both u and r to save cycles
++        if u & 1 == 0:
++            u >>= 1
++            if r & 1 == 0:
++                r >>= 1
++            else:
++                r = (r + p) >> 1
++        # implementations could use count-zeros on
++        # both v and s to save cycles
++        elif v & 1 == 0:
++            v >>= 1
++            if s & 1 == 0:
++                s >>= 1
++            else:
++                s = (s + p) >> 1
++        else:
++            x = u - v
++            if x > 0:
++                u = x
++                r -= s
++                if r < 0:
++                    r += p
++            else:
++                v = -x
++                s -= r
++                if s < 0:
++                    s += p
++    if r > p:
++        r -= p
++    if r < 0:
++        r += p
++    return r

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmadd.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..e7159a2da4cfc1a1080e55c3e533628b5857b104
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmadd.py
@@@ -1,0 -1,0 +1,5 @@@
++from .state import ST
++
++
++def gfpmadd(a, b, c):
++    return (a * b + c) % ST.GFPRIME

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmsub.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..fd0912480015ea3a49f8a2ddac033240d8136399
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmsub.py
@@@ -1,0 -1,0 +1,5 @@@
++from .state import ST
++
++
++def gfpmsub(a, b, c):
++    return (a * b - c) % ST.GFPRIME

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmsubr.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..2d349f85803fced2dd462f5197d77367abff9e26
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmsubr.py
@@@ -1,0 -1,0 +1,5 @@@
++from .state import ST
++
++
++def gfpmsubr(a, b, c):
++    return (c - a * b) % ST.GFPRIME

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmul.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..42c43a23313e79fb1254a81e4559d958bb802e76
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfpmul.py
@@@ -1,0 -1,0 +1,5 @@@
++from .state import ST
++
++
++def gfpmul(a, b):
++    return (a * b) % ST.GFPRIME

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/gfpsub.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..44e5798d7a174b4381f47c3004949e1efc4136e7
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/gfpsub.py
@@@ -1,0 -1,0 +1,5 @@@
++from .state import ST
++
++
++def gfpsub(a, b):
++    return (a - b) % ST.GFPRIME

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/log2.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..0defcc89f899645b152aef02555505220966827a
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/log2.py
@@@ -1,0 -1,0 +1,12 @@@
++def floor_log2(v):
++    """return floor(log2(v))."""
++    assert isinstance(v, int)
++    assert v > 0
++    return v.bit_length() - 1
++
++
++def ceil_log2(v):
++    """return ceil(log2(v))."""
++    assert isinstance(v, int)
++    assert v > 0
++    return (v - 1).bit_length()

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/pack_poly.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..f5ab64467950ed22f29729e3fd8b13dbb3898943
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/pack_poly.py
@@@ -1,0 -1,0 +1,19 @@@
++"""Polynomials with GF(2) coefficients."""
++
++
++def pack_poly(poly):
++    """`poly` is a list where `poly[i]` is the coefficient for `x ** i`"""
++    retval = 0
++    for i, v in enumerate(poly):
++        retval |= v << i
++    return retval
++
++
++def unpack_poly(v):
++    """returns a list `poly`, where `poly[i]` is the coefficient for `x ** i`.
++    """
++    poly = []
++    while v != 0:
++        poly.append(v & 1)
++        v >>= 1
++    return poly

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/state.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..91cb5426f29c51314a671a46639649c9ee8ccedc
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/state.py
@@@ -1,0 -1,0 +1,19 @@@
++from .log2 import floor_log2
++from threading import local
++
++
++class State(local):
++    # thread local so unit tests can be run in parallel without breaking
++    def __init__(self, *, XLEN=64, GFBREDPOLY=0, GFPRIME=31):
++        assert isinstance(XLEN, int) and 2 ** floor_log2(XLEN) == XLEN
++        assert isinstance(GFBREDPOLY, int) and 0 <= GFBREDPOLY < 2 ** 64
++        assert isinstance(GFPRIME, int) and 0 <= GFPRIME < 2 ** 64
++        self.XLEN = XLEN
++        self.GFBREDPOLY = GFBREDPOLY
++        self.GFPRIME = GFPRIME
++
++    def reinit(self, *, XLEN=64, GFBREDPOLY=0, GFPRIME=31):
++        self.__init__(XLEN=XLEN, GFBREDPOLY=GFBREDPOLY, GFPRIME=GFPRIME)
++
++
++ST = State()

diff --cc src/nmutil/openpower_sv_bitmanip_in_wiki/test_cl_gfb_gfp.py
index 0000000000000000000000000000000000000000,0000000000000000000000000000000000000000..34133cbe3a4a65a7c281f02b1e39a0d458adb2e5
new file mode 100644 (file)
--- /dev/null
--- /dev/null
+++ b/src/nmutil/openpower_sv_bitmanip_in_wiki/test_cl_gfb_gfp.py
@@@ -1,0 -1,0 +1,701 @@@
++from .state import ST
++from .cldivrem import cldivrem
++from .clmul import clmul
++from .gfbmul import gfbmul
++from .gfbmadd import gfbmadd
++from .gfbinv import gfbinv
++from .gfpadd import gfpadd
++from .gfpsub import gfpsub
++from .gfpmul import gfpmul
++from .gfpinv import gfpinv
++from .gfpmadd import gfpmadd
++from .gfpmsub import gfpmsub
++from .gfpmsubr import gfpmsubr
++from .pack_poly import pack_poly, unpack_poly
++import unittest
++
++
++class GF2Poly:
++    """Polynomial with GF(2) coefficients.
++
++    `self.coefficients`: a list where `coefficients[-1] != 0`.
++    `coefficients[i]` is the coefficient for `x ** i`.
++    """
++
++    def __init__(self, coefficients=None):
++        self.coefficients = []
++        if coefficients is not None:
++            if not isinstance(coefficients, (tuple, list)):
++                coefficients = list(coefficients)
++            # reversed to resize self.coefficients once
++            for i in reversed(range(len(coefficients))):
++                self[i] = coefficients[i]
++
++    def __len__(self):
++        return len(self.coefficients)
++
++    @property
++    def degree(self):
++        return len(self) - 1
++
++    @property
++    def lc(self):
++        """leading coefficient."""
++        return 0 if len(self) == 0 else self.coefficients[-1]
++
++    def __getitem__(self, key):
++        assert key >= 0
++        if key < len(self):
++            return self.coefficients[key]
++        return 0
++
++    def __setitem__(self, key, value):
++        assert key >= 0
++        assert value == 0 or value == 1
++        if key < len(self):
++            self.coefficients[key] = value
++            while len(self) and self.coefficients[-1] == 0:
++                self.coefficients.pop()
++        elif value != 0:
++            self.coefficients += [0] * (key + 1 - len(self))
++            self.coefficients[key] = value
++
++    def __repr__(self):
++        return f"GF2Poly({self.coefficients})"
++
++    def __iadd__(self, rhs):
++        for i in range(max(len(self), len(rhs))):
++            self[i] ^= rhs[i]
++        return self
++
++    def __add__(self, rhs):
++        return GF2Poly(self).__iadd__(rhs)
++
++    def __isub__(self, rhs):
++        return self.__iadd__(rhs)
++
++    def __sub__(self, rhs):
++        return self.__add__(rhs)
++
++    def __iter__(self):
++        return iter(self.coefficients)
++
++    def __mul__(self, rhs):
++        retval = GF2Poly()
++        # reversed to resize retval.coefficients once
++        for i in reversed(range(len(self))):
++            if self[i]:
++                for j in reversed(range(len(rhs))):
++                    retval[i + j] ^= rhs[j]
++        return retval
++
++    def __ilshift__(self, amount):
++        """multiplies `self` by the polynomial `x**amount`"""
++        if len(self) != 0:
++            self.coefficients[:0] = [0] * amount
++        return self
++
++    def __lshift__(self, amount):
++        """returns the polynomial `self * x**amount`"""
++        return GF2Poly(self).__ilshift__(amount)
++
++    def __irshift__(self, amount):
++        """divides `self` by the polynomial `x**amount`, discarding the
++            remainder.
++        """
++        if amount < len(self):
++            del self.coefficients[:amount]
++        else:
++            del self.coefficients[:]
++        return self
++
++    def __rshift__(self, amount):
++        """divides `self` by the polynomial `x**amount`, discarding the
++            remainder.
++        """
++        return GF2Poly(self).__irshift__(amount)
++
++    def __divmod__(self, divisor):
++        # based on https://en.wikipedia.org/wiki/Polynomial_greatest_common_divisor#Euclidean_division
++        assert isinstance(divisor, GF2Poly)
++        if len(divisor) == 0:
++            raise ZeroDivisionError
++        q = GF2Poly()
++        r = GF2Poly(self)
++        while r.degree >= divisor.degree:
++            shift = r.degree - divisor.degree
++            q[shift] ^= 1
++            r -= divisor << shift
++        return q, r
++
++    def __floordiv__(self, divisor):
++        q, r = divmod(self, divisor)
++        return q
++
++    def __mod__(self, divisor):
++        q, r = divmod(self, divisor)
++        return r
++
++    def __pow__(self, exponent, modulus=None):
++        assert isinstance(exponent, int) and exponent >= 0
++        assert modulus is None or isinstance(modulus, GF2Poly)
++        retval = GF2Poly([1])
++        pow2 = GF2Poly(self)
++        while exponent != 0:
++            if exponent & 1:
++                retval *= pow2
++                if modulus is not None:
++                    retval %= modulus
++                exponent &= ~1
++            else:
++                pow2 *= pow2
++                if modulus is not None:
++                    pow2 %= modulus
++                exponent >>= 1
++        return retval
++
++    def __eq__(self, rhs):
++        if isinstance(rhs, GF2Poly):
++            return self.coefficients == rhs.coefficients
++        return NotImplemented
++
++
++class TestGF2Poly(unittest.TestCase):
++    def test_add(self):
++        a = GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1])
++        b = GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1, 0, 0])
++        c = a + b
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(b, GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1]))
++        self.assertEqual(c, GF2Poly([1, 0, 1, 1, 1, 0, 0, 1, 0, 0, 1]))
++        c = b + a
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(b, GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1]))
++        self.assertEqual(c, GF2Poly([1, 0, 1, 1, 1, 0, 0, 1, 0, 0, 1]))
++        a = GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1])
++        b = GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1])
++        c = a + b
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(b, GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly([1, 0, 1, 1, 1, 0, 0, 1]))
++        c = a - b
++        self.assertEqual(c, GF2Poly([1, 0, 1, 1, 1, 0, 0, 1]))
++        c = b - a
++        self.assertEqual(c, GF2Poly([1, 0, 1, 1, 1, 0, 0, 1]))
++
++    def test_shift(self):
++        a = GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1])
++        c = a << 0
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        c = a << 5
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly(
++            [0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        c = a << 10
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly(
++            [0] * 10 + [1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        c = a >> 0
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        c = a >> 5
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly([1, 1, 0, 1, 0, 1]))
++        c = a >> 10
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly([1]))
++        c = a >> 11
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly([]))
++        c = a >> 100
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(c, GF2Poly([]))
++
++    def test_mul(self):
++        a = GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1])
++        b = GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1])
++        c = a * b
++        expected = GF2Poly([0, 0, 1, 0, 1, 0, 1, 1, 1, 0,
++                            0, 1, 0, 1, 1, 0, 0, 1, 1])
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(b, GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1]))
++        self.assertEqual(c, expected)
++
++    def test_divmod(self):
++        a = GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1])
++        b = GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1])
++        q, r = divmod(a, b)
++        self.assertEqual(a, GF2Poly([1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1]))
++        self.assertEqual(b, GF2Poly([0, 0, 1, 0, 1, 1, 1, 1, 1]))
++        self.assertEqual(q, GF2Poly([1, 1, 1]))
++        self.assertEqual(r, GF2Poly([1, 0, 1, 0, 0, 1, 0, 1]))
++        q = a // b
++        self.assertEqual(q, GF2Poly([1, 1, 1]))
++        r = a % b
++        self.assertEqual(r, GF2Poly([1, 0, 1, 0, 0, 1, 0, 1]))
++
++    def test_pow(self):
++        b = GF2Poly([0, 1])
++        for e in range(8):
++            expected = GF2Poly([0] * e + [1])
++            with self.subTest(b=str(b), e=e, expected=str(expected)):
++                v = b ** e
++                self.assertEqual(b, GF2Poly([0, 1]))
++                self.assertEqual(v, expected)
++
++        # AES's finite field reducing polynomial
++        m = GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1])
++        period = 2 ** m.degree - 1
++        b = GF2Poly([1, 1, 0, 0, 1, 0, 1])
++        e = period - 1
++        expected = GF2Poly([0, 1, 0, 1, 0, 0, 1, 1])
++        v = pow(b, e, m)
++        self.assertEqual(m, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(b, GF2Poly([1, 1, 0, 0, 1, 0, 1]))
++        self.assertEqual(v, expected)
++
++        # test that pow doesn't take inordinately long when given a modulus.
++        # adding a multiple of `period` should leave results unchanged.
++        e += period * 10 ** 15
++        v = pow(b, e, m)
++        self.assertEqual(m, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(b, GF2Poly([1, 1, 0, 0, 1, 0, 1]))
++        self.assertEqual(v, expected)
++
++
++class GFB:
++    def __init__(self, value, red_poly=None):
++        if isinstance(value, GFB):
++            # copy value
++            assert red_poly is None
++            self.red_poly = GF2Poly(value.red_poly)
++            self.value = GF2Poly(value.value)
++            return
++        assert isinstance(value, GF2Poly)
++        assert isinstance(red_poly, GF2Poly)
++        assert red_poly.degree > 0
++        self.value = value % red_poly
++        self.red_poly = red_poly
++
++    def __repr__(self):
++        return f"GFB({self.value}, {self.red_poly})"
++
++    def __add__(self, rhs):
++        assert isinstance(rhs, GFB)
++        assert self.red_poly == rhs.red_poly
++        return GFB((self.value + rhs.value) % self.red_poly, self.red_poly)
++
++    def __sub__(self, rhs):
++        return self.__add__(rhs)
++
++    def __eq__(self, rhs):
++        if isinstance(rhs, GFB):
++            return self.value == rhs.value and self.red_poly == rhs.red_poly
++        return NotImplemented
++
++    def __mul__(self, rhs):
++        assert isinstance(rhs, GFB)
++        assert self.red_poly == rhs.red_poly
++        return GFB((self.value * rhs.value) % self.red_poly, self.red_poly)
++
++    def __div__(self, rhs):
++        assert isinstance(rhs, GFB)
++        assert self.red_poly == rhs.red_poly
++        return self * rhs ** -1
++
++    @property
++    def __pow_period(self):
++        period = (1 << self.red_poly.degree) - 1
++        assert period > 0, "internal logic error"
++        return period
++
++    def __pow__(self, exponent):
++        assert isinstance(exponent, int)
++        if len(self.value) == 0:
++            if exponent < 0:
++                raise ZeroDivisionError
++            else:
++                return GFB(self)
++        exponent %= self.__pow_period
++        return GFB(pow(self.value, exponent, self.red_poly), self.red_poly)
++
++
++class TestGFBClass(unittest.TestCase):
++    def test_add(self):
++        # AES's finite field reducing polynomial
++        red_poly = GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1])
++        a = GFB(GF2Poly([0, 1, 0, 1]), red_poly)
++        b = GFB(GF2Poly([0, 0, 0, 0, 0, 0, 1, 1]), red_poly)
++        expected = GFB(GF2Poly([0, 1, 0, 1, 0, 0, 1, 1]), red_poly)
++        c = a + b
++        self.assertEqual(red_poly, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(a, GFB(GF2Poly([0, 1, 0, 1]), red_poly))
++        self.assertEqual(b, GFB(GF2Poly([0, 0, 0, 0, 0, 0, 1, 1]), red_poly))
++        self.assertEqual(c, expected)
++        c = a - b
++        self.assertEqual(red_poly, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(a, GFB(GF2Poly([0, 1, 0, 1]), red_poly))
++        self.assertEqual(b, GFB(GF2Poly([0, 0, 0, 0, 0, 0, 1, 1]), red_poly))
++        self.assertEqual(c, expected)
++        c = b - a
++        self.assertEqual(red_poly, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(a, GFB(GF2Poly([0, 1, 0, 1]), red_poly))
++        self.assertEqual(b, GFB(GF2Poly([0, 0, 0, 0, 0, 0, 1, 1]), red_poly))
++        self.assertEqual(c, expected)
++
++    def test_mul(self):
++        # AES's finite field reducing polynomial
++        red_poly = GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1])
++        a = GFB(GF2Poly([0, 1, 0, 1, 0, 0, 1, 1]), red_poly)
++        b = GFB(GF2Poly([1, 1, 0, 0, 1, 0, 1]), red_poly)
++        expected = GFB(GF2Poly([1]), red_poly)
++        c = a * b
++        self.assertEqual(red_poly, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(a, GFB(GF2Poly([0, 1, 0, 1, 0, 0, 1, 1]), red_poly))
++        self.assertEqual(b, GFB(GF2Poly([1, 1, 0, 0, 1, 0, 1]), red_poly))
++        self.assertEqual(c, expected)
++
++    def test_pow(self):
++        # AES's finite field reducing polynomial
++        red_poly = GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1])
++        period = 2 ** red_poly.degree - 1
++        b = GFB(GF2Poly([1, 1, 0, 0, 1, 0, 1]), red_poly)
++        e = period - 1
++        expected = GFB(GF2Poly([0, 1, 0, 1, 0, 0, 1, 1]), red_poly)
++        v = b ** e
++        self.assertEqual(red_poly, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(b, GFB(GF2Poly([1, 1, 0, 0, 1, 0, 1]), red_poly))
++        self.assertEqual(v, expected)
++        e = -1
++        v = b ** e
++        self.assertEqual(red_poly, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(b, GFB(GF2Poly([1, 1, 0, 0, 1, 0, 1]), red_poly))
++        self.assertEqual(v, expected)
++
++        # test that pow doesn't take inordinately long when given a modulus.
++        # adding a multiple of `period` should leave results unchanged.
++        e += period * 10 ** 15
++        v = b ** e
++        self.assertEqual(red_poly, GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1]))
++        self.assertEqual(b, GFB(GF2Poly([1, 1, 0, 0, 1, 0, 1]), red_poly))
++        self.assertEqual(v, expected)
++
++
++class GFP:
++    def __init__(self, value, size):
++        assert isinstance(value, int)
++        assert isinstance(size, int) and size >= 2, "size is not a prime"
++        self.value = value % size
++        self.size = size
++
++    def __repr__(self):
++        return f"GFP({self.value}, {self.size})"
++
++    def __eq__(self, rhs):
++        if isinstance(rhs, GFP):
++            return self.value == rhs.value and self.size == rhs.size
++        return NotImplemented
++
++    def __add__(self, rhs):
++        assert isinstance(rhs, GFP)
++        assert self.size == rhs.size
++        return GFP((self.value + rhs.value) % self.size, self.size)
++
++    def __sub__(self, rhs):
++        assert isinstance(rhs, GFP)
++        assert self.size == rhs.size
++        return GFP((self.value - rhs.value) % self.size, self.size)
++
++    def __mul__(self, rhs):
++        assert isinstance(rhs, GFP)
++        assert self.size == rhs.size
++        return GFP((self.value * rhs.value) % self.size, self.size)
++
++    def __div__(self, rhs):
++        assert isinstance(rhs, GFP)
++        assert self.size == rhs.size
++        return self * rhs ** -1
++
++    @property
++    def __pow_period(self):
++        period = self.size - 1
++        assert period > 0, "internal logic error"
++        return period
++
++    def __pow__(self, exponent):
++        assert isinstance(exponent, int)
++        if self.value == 0:
++            if exponent < 0:
++                raise ZeroDivisionError
++            else:
++                return GFP(self.value, self.size)
++        exponent %= self.__pow_period
++        return GFP(pow(self.value, exponent, self.size), self.size)
++
++
++PRIMES = 2, 3, 5, 7, 11, 13, 17, 19
++"""handy list of small primes for testing"""
++
++
++class TestGFPClass(unittest.TestCase):
++    def test_add_sub(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    with self.subTest(av=av, bv=bv, prime=prime):
++                        a = GFP(av, prime)
++                        b = GFP(bv, prime)
++                        expected = GFP((av + bv) % prime, prime)
++                        c = a + b
++                        self.assertEqual(a, GFP(av, prime))
++                        self.assertEqual(b, GFP(bv, prime))
++                        self.assertEqual(c, expected)
++                        c = b + a
++                        self.assertEqual(a, GFP(av, prime))
++                        self.assertEqual(b, GFP(bv, prime))
++                        self.assertEqual(c, expected)
++                        expected = GFP((av - bv) % prime, prime)
++                        c = a - b
++                        self.assertEqual(a, GFP(av, prime))
++                        self.assertEqual(b, GFP(bv, prime))
++                        self.assertEqual(c, expected)
++                        expected = GFP((bv - av) % prime, prime)
++                        c = b - a
++                        self.assertEqual(a, GFP(av, prime))
++                        self.assertEqual(b, GFP(bv, prime))
++                        self.assertEqual(c, expected)
++
++    def test_mul(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    with self.subTest(av=av, bv=bv, prime=prime):
++                        a = GFP(av, prime)
++                        b = GFP(bv, prime)
++                        expected = GFP((av * bv) % prime, prime)
++                        c = a * b
++                        self.assertEqual(a, GFP(av, prime))
++                        self.assertEqual(b, GFP(bv, prime))
++                        self.assertEqual(c, expected)
++                        c = b * a
++                        self.assertEqual(a, GFP(av, prime))
++                        self.assertEqual(b, GFP(bv, prime))
++                        self.assertEqual(c, expected)
++
++    def test_pow(self):
++        for prime in PRIMES:
++            for bv in range(prime):
++                with self.subTest(bv=bv, prime=prime):
++                    b = GFP(bv, prime)
++                    period = prime - 1
++                    e = period - 1
++                    expected = GFP(pow(bv, e, prime) if bv != 0 else 0, prime)
++                    v = b ** e
++                    self.assertEqual(b, GFP(bv, prime))
++                    self.assertEqual(v, expected)
++                    e = -1
++                    if bv != 0:
++                        v = b ** e
++                        self.assertEqual(b, GFP(bv, prime))
++                        self.assertEqual(v, expected)
++
++                    # test that pow doesn't take inordinately long when given
++                    # a modulus. adding a multiple of `period` should leave
++                    # results unchanged.
++                    e += period * 10 ** 15
++                    v = b ** e
++                    self.assertEqual(b, GFP(bv, prime))
++                    self.assertEqual(v, expected)
++
++
++class TestCL(unittest.TestCase):
++    def test_cldivrem(self):
++        n_width = 8
++        d_width = 4
++        width = max(n_width, d_width)
++        for nv in range(2 ** n_width):
++            n = GF2Poly(unpack_poly(nv))
++            for dv in range(1, 2 ** d_width):
++                d = GF2Poly(unpack_poly(dv))
++                with self.subTest(n=str(n), nv=nv, d=str(d), dv=dv):
++                    q_expected, r_expected = divmod(n, d)
++                    self.assertEqual(q_expected * d + r_expected, n)
++                    q, r = cldivrem(nv, dv, width)
++                    q_expected = pack_poly(q_expected.coefficients)
++                    r_expected = pack_poly(r_expected.coefficients)
++                    self.assertEqual((q, r), (q_expected, r_expected))
++
++    def test_clmul(self):
++        a_width = 5
++        b_width = 5
++        for av in range(2 ** a_width):
++            a = GF2Poly(unpack_poly(av))
++            for bv in range(2 ** b_width):
++                b = GF2Poly(unpack_poly(bv))
++                with self.subTest(a=str(a), av=av, b=str(b), bv=bv):
++                    expected = a * b
++                    product = clmul(av, bv)
++                    expected = pack_poly(expected.coefficients)
++                    self.assertEqual(product, expected)
++
++
++class TestGFBInstructions(unittest.TestCase):
++    @staticmethod
++    def init_aes_red_poly():
++        # AES's finite field reducing polynomial
++        red_poly = GF2Poly([1, 1, 0, 1, 1, 0, 0, 0, 1])
++        ST.reinit(GFBREDPOLY=pack_poly(red_poly.coefficients))
++        return red_poly
++
++    def test_gfbmul(self):
++        # AES's finite field reducing polynomial
++        red_poly = self.init_aes_red_poly()
++        a_width = 8
++        b_width = 4
++        for av in range(2 ** a_width):
++            a = GFB(GF2Poly(unpack_poly(av)), red_poly)
++            for bv in range(2 ** b_width):
++                b = GFB(GF2Poly(unpack_poly(bv)), red_poly)
++                expected = a * b
++                with self.subTest(a=str(a), av=av, b=str(b), bv=bv, expected=str(expected)):
++                    product = gfbmul(av, bv)
++                    expectedv = pack_poly(expected.value.coefficients)
++                    self.assertEqual(product, expectedv)
++
++    def test_gfbmadd(self):
++        # AES's finite field reducing polynomial
++        red_poly = self.init_aes_red_poly()
++        a_width = 5
++        b_width = 4
++        c_width = 4
++        for av in range(2 ** a_width):
++            a = GFB(GF2Poly(unpack_poly(av)), red_poly)
++            for bv in range(2 ** b_width):
++                b = GFB(GF2Poly(unpack_poly(bv)), red_poly)
++                for cv in range(2 ** c_width):
++                    c = GFB(GF2Poly(unpack_poly(cv)), red_poly)
++                    expected = a * b + c
++                    with self.subTest(a=str(a), av=av,
++                                      b=str(b), bv=bv,
++                                      c=str(c), cv=cv,
++                                      expected=str(expected)):
++                        result = gfbmadd(av, bv, cv)
++                        expectedv = pack_poly(expected.value.coefficients)
++                        self.assertEqual(result, expectedv)
++
++    def test_gfbinv(self):
++        # AES's finite field reducing polynomial
++        red_poly = self.init_aes_red_poly()
++        width = 8
++        for av in range(2 ** width):
++            a = GFB(GF2Poly(unpack_poly(av)), red_poly)
++            expected = a ** -1 if av != 0 else GFB(GF2Poly(), red_poly)
++            with self.subTest(a=str(a), av=av, expected=str(expected)):
++                result = gfbinv(av)
++                expectedv = pack_poly(expected.value.coefficients)
++                self.assertEqual(result, expectedv)
++
++
++class TestGFPInstructions(unittest.TestCase):
++    def test_gfpadd(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    a = GFP(av, prime)
++                    b = GFP(bv, prime)
++                    expected = a + b
++                    with self.subTest(a=str(a), b=str(b),
++                                      expected=str(expected)):
++                        ST.reinit(GFPRIME=prime)
++                        v = gfpadd(av, bv)
++                        self.assertEqual(v, expected.value)
++
++    def test_gfpsub(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    a = GFP(av, prime)
++                    b = GFP(bv, prime)
++                    expected = a - b
++                    with self.subTest(a=str(a), b=str(b),
++                                      expected=str(expected)):
++                        ST.reinit(GFPRIME=prime)
++                        v = gfpsub(av, bv)
++                        self.assertEqual(v, expected.value)
++
++    def test_gfpmul(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    a = GFP(av, prime)
++                    b = GFP(bv, prime)
++                    expected = a * b
++                    with self.subTest(a=str(a), b=str(b),
++                                      expected=str(expected)):
++                        ST.reinit(GFPRIME=prime)
++                        v = gfpmul(av, bv)
++                        self.assertEqual(v, expected.value)
++
++    def test_gfpinv(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                a = GFP(av, prime)
++                if av == 0:
++                    # TODO: determine what's expected for division by zero
++                    continue
++                else:
++                    expected = a ** -1
++                with self.subTest(a=str(a), expected=str(expected)):
++                    ST.reinit(GFPRIME=prime)
++                    v = gfpinv(av)
++                    self.assertEqual(v, expected.value)
++
++    def test_gfpmadd(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    for cv in range(prime):
++                        a = GFP(av, prime)
++                        b = GFP(bv, prime)
++                        c = GFP(cv, prime)
++                        expected = a * b + c
++                        with self.subTest(a=str(a), b=str(b), c=str(c),
++                                          expected=str(expected)):
++                            ST.reinit(GFPRIME=prime)
++                            v = gfpmadd(av, bv, cv)
++                            self.assertEqual(v, expected.value)
++
++    def test_gfpmsub(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    for cv in range(prime):
++                        a = GFP(av, prime)
++                        b = GFP(bv, prime)
++                        c = GFP(cv, prime)
++                        expected = a * b - c
++                        with self.subTest(a=str(a), b=str(b), c=str(c),
++                                          expected=str(expected)):
++                            ST.reinit(GFPRIME=prime)
++                            v = gfpmsub(av, bv, cv)
++                            self.assertEqual(v, expected.value)
++
++    def test_gfpmsubr(self):
++        for prime in PRIMES:
++            for av in range(prime):
++                for bv in range(prime):
++                    for cv in range(prime):
++                        a = GFP(av, prime)
++                        b = GFP(bv, prime)
++                        c = GFP(cv, prime)
++                        expected = c - a * b
++                        with self.subTest(a=str(a), b=str(b), c=str(c),
++                                          expected=str(expected)):
++                            ST.reinit(GFPRIME=prime)
++                            v = gfpmsubr(av, bv, cv)
++                            self.assertEqual(v, expected.value)
++
++
++if __name__ == "__main__":
++    unittest.main()