Determining the binary operation from a Cayley Table

Question

I have tried a lot of things with this Cayley Table (several Julia/Python scripts which iterate over various functions, symbolic regression, semi-manually trying various permutation groups, octonions, etc) but they haven't really turned out to give me an answer yet.

I have this Cayley Table which appears to be non-associative (so not a group?). I am not 'really a mathematician', so perhaps this is simpler than I it seems to me. This problem is from a recent paper on 'grokking' if that matters to anyone.

Here is the Cayley Table:

★	a	b	c	d	e
a	a	d	?	c	d
b	c	d	d	a	c
c	?	e	d	b	d
d	a	?	?	b	c
e	b	b	c	?	a

Am I missing a typical way to approach this problem other than trying symbolic regression on it? Unfortunately - at least the way I have tried - this requires assuming the elements are some set of arbitrarily chosen integers, or other data structure - which is obviously not ideal. I tried iterating through the 'properties' of the system (i.e. non-associative) but it left me empty handed since I don't know what the general ontology I'm traversing is. Apologies for completely misusing the terminology here - I still have a lot to learn since mathematics is not directly my field of study.

Answer 1

There are $5^{25}$=298 quadrillion possible binary operations $[0..4]^2\to[0..4]$, so you probably need to narrow things down somehow. All possible functions are attainable as polynomials, so that doesn't limit anything. (For proof, make "indicator" functions like $I_3(k) := 4(k-0)(k-1)(k-2)(k-4)$ that vanish except at 3, then put $f(x, y) = \sum_{a,b} f(a,b)I_a(x)I_b(y)$.)

I tried this Python script, which tries most of the operations listed Appendix A in that paper, and shuffles among all possible assignments of a, b, c, d, and e:

ops = [
    "(x+y)%5",
    "(x-y)%5",
    "x*pow(y,-1,5)%5 if y else 0",
    "x*pow(y,-1,5)%5 if y&1 else (x-y)%5",
    "(x*x+y*y)%5",
    "(x*x+x*y+y*y)%5",
    "(x**2+x*y+y**2+x)%5",
    "(x**3+x*y)%5",
    "(x**3+x*y**2+y)%5"]

opfuncs = [(s, lambda x, y, s=s:eval(s, {'x':x, 'y':y}))
           for s in ops]

import itertools


data = "aaa abd adc aed bac bbd bcd bda bec cbe ccd cdb ced daa ddb dec eab ebb ecc eea".split()

for s, opfunc in opfuncs:
    for tup in itertools.permutations(range(5)):
        D = dict(zip("abcde", tup))
        inverse = {v:k for k, v in D.items()}
        for x,y,z in data:
            if opfunc(D[x], D[y]) != D[z]:
                break
        else:
            # consistent
            print(s, D)
            for x, y in "ac ca db dc ed".split():
                print(x, y, "-->", inverse[opfunc(D[x], D[y])])

The result:

(x**2+x*y+y**2+x)%5 {'a': 0, 'b': 2, 'c': 1, 'd': 4, 'e': 3}
a c --> c
c a --> b
d b --> b
d c --> a
e d --> a

So in other words, $a=0$, $b=2$, $c=1$, $d=4$, and $e=3$, with $a\star b = x^2 + xy + y^2 + x \pmod 5$.