Strassen's algorithm: Difference between revisions

Output is the same as the dynamically padded version.

<lang Julia>function Strassen(A::Matrix, B::Matrix)

@views A11 = A[1:n÷2, 1:n÷2]

@views A12 = A[1:n÷2, n÷2+1:n]

@views A21 = A[n÷2+1:n, 1:n÷2]

@views A22 = A[n÷2+1:n, n÷2+1:n]

@views B11 = B[1:n÷2, 1:n÷2]

@views B12 = B[1:n÷2, n÷2+1:n]

@views B21 = B[n÷2+1:n, 1:n÷2]

@views B22 = B[n÷2+1:n, n÷2+1:n]

P4 = Strassen(A11 + A12, Matrix(B22))

P5 = Strassen(Matrix(A11), B12 - B22)

P6 = Strassen(Matrix(A22), B21 - B11)

P7 = Strassen(A21 + A22, Matrix(B11))

end

const A = [[1, 2] [3, 4]]

const B = [[5, 6] [7, 8]]

const C = [[1, 1, 1, 1] [2, 4, 8, 16] [3, 9, 27, 81] [4, 16, 64, 256]]

const D = [[4, -3, 4/3, -1/4] [-13/3, 19/4, -7/3, 11/24] [3/2, -2, 7/6, -1/4] [-1/6, 1/4, -1/6, 1/24]]

const E = [[1, 2, 3, 4] [5, 6, 7, 8] [9, 10, 11, 12] [13, 14, 15, 16]]

const F = [[1, 0, 0, 0] [0, 1, 0, 0] [0, 0, 1, 0] [0, 0, 0, 1]]

intprint(s, mat) = println(s, map(x -> Int(round(x, digits=8)), mat)')

intprint("Regular multiply: ", A' * B')

intprint("Strassen multiply: ", Strassen(Matrix(A'), Matrix(B')))

intprint("Regular multiply: ", C * D)

intprint("Strassen multiply: ", Strassen(C, D))

intprint("Regular multiply: ", E * F)

intprint("Strassen multiply: ", Strassen(E, F))

const r = sqrt(2)/2

const R = [[r, r] [-r, r]]

intprint("Regular multiply: ", R * R)

intprint("Strassen multiply: ", Strassen(R,R))