Implement Strassen's algorithm [closed]

Question

Strassen's algorithm was the first method of matrix multiplication in subcubic time complexity, namely O(n**log2(7)) for a pair of n*n matrices (assuming the numbers therein are large enough that their O(n*log2(n)) exact multiplication has become the constraint on your performance, and any nested structures and function calls are negligible in comparison).

For a 2*2 matrix, it is defined as

lambda a,b: (lambda a,b,c,d,e,f,g: (a+d-e+g,c+e,b+d,a-b+c+f))((a[0]+a[3])*(b[0]+b[3]),(a[2]+a[3])*b[0],a[0]*(b[1]-b[3]),a[3]*(b[2]-b[0]),(a[0]+a[1])*b[3],(a[2]-a[0])*(b[0]+b[1]),(a[1]-a[3])*(b[2]+b[3]))

And for larger square ones, as splitting them into quarters and then calling this but with numbers' addition, negation and multiplication methods substituted for matrices' (notably using itself for the latter, its seven self-calls per doubling of the width and height being the reason for the exponent).

For simplicity (so you won't need to deal with implementing standard matrix multiplication and subdividing into Strassen-able ones), your inputs will be two 2**n*2**n matrices of integers, represented as length-2**(2*n) tuples (or lists or your language's equivalent), encoded with elements in reading order, and you will return another such tuple. For instance, when inputted with these two

(5,2,0,0,
 4,0,5,2,
 3,4,5,0,
 3,1,4,2)

(7,4,5,3,
 4,0,2,7,
 1,4,1,1,
 0,5,3,5)

, it should return

(43,20,29,29,
 33,46,31,27,
 42,32,28,42,
 29,38,27,30)

Answer 1

One-line Python (1027 bytes)

(explained in this gist, would be greatly reducible if you were to be less pedantic about it all being one line)

from functools import reduce
from math import isqrt
from random import random
lap=(lambda f,*a: list(map(f,*a)))
strassen=(lambda a,b: tuple((lambda f: f(a,b) if len(a)==4 else (lambda m: (lambda m,d: ((m[y*2+x][d*j+i] for y in range(2) for j in range(d) for x in range(2) for i in range(d))))(m,isqrt(len(m[0]))))(f(*(map((lambda m: (lambda m,d: tuple((tuple((m[i*d+j] for i in range(d//2*y,d//(2-y)) for j in range(d//2*x,d//(2-x)))) for y in range(2) for x in range(2))))(m,isqrt(len(m)))),(a,b))))))((lambda a,b: (lambda s,n,o,m: (lambda a,b,c,d,e,f,g: lap(s,((a,d,n(e),g),(c,e),(b,d),(a,n(b),c,f))))(*map(m,*map(lambda x: map(o,x),(((a[0],a[3]),(a[2],a[3]),(a[0],),(a[3],),(a[0],a[1]),(a[2],n(a[0])),(a[1],n(a[3]))),((b[0],b[3]),(b[0],),(b[1],n(b[3])),(b[2],n(b[0])),(b[3],),(b[0],b[1]),(b[2],b[3])))))))(*((sum,int.__neg__,sum,int.__mul__) if type(a[0])==int else ((lambda x: reduce(lambda a,b: lap(int.__add__,a,b),x)),(lambda x: map(int.__neg__,x)),(lambda t: t[0] if len(t)==1 else lap(int.__add__,*t)),strassen)))))))