Associativity and matrix multiplication

The matrix multiplication m1 @ m2 @ m3 can be done in two different ways: (m1 @ m2) @ m3 or m1 @ (m2 @ m3). Are these two orders equivalent or is there a better order?

import pprint
import numpy
import matplotlib.pyplot as plt
from pandas import DataFrame
from tqdm import tqdm
from teachcompute.ext_test_case import measure_time

First try

m1 = numpy.random.rand(100, 100)
m2 = numpy.random.rand(100, 10)
m3 = numpy.random.rand(10, 100)

m = m1 @ m2 @ m3

print(m.shape)

mm1 = (m1 @ m2) @ m3
mm2 = m1 @ (m2 @ m3)

print(mm1.shape, mm2.shape)

t1 = measure_time(lambda: (m1 @ m2) @ m3, context={}, number=50, repeat=50)
pprint.pprint(t1)

t2 = measure_time(lambda: m1 @ (m2 @ m3), context={}, number=50, repeat=50)
pprint.pprint(t2)

(100, 100)
(100, 100) (100, 100)
{'average': np.float64(3.261241159925703e-05),
 'context_size': 64,
 'deviation': np.float64(1.0855784817023806e-05),
 'max_exec': np.float64(9.926761995302514e-05),
 'min_exec': np.float64(2.255465995403938e-05),
 'number': 50,
 'repeat': 50,
 'ttime': np.float64(0.0016306205799628516),
 'warmup_time': 0.00016044600488385186}
{'average': np.float64(7.313718521909322e-05),
 'context_size': 64,
 'deviation': np.float64(1.0503636143782488e-05),
 'max_exec': np.float64(0.00013949059997685252),
 'min_exec': np.float64(6.590636010514572e-05),
 'number': 50,
 'repeat': 50,
 'ttime': np.float64(0.0036568592609546612),
 'warmup_time': 9.652199514675885e-05}

With different sizes

obs = []
for i in tqdm([50, 100, 125, 150, 175, 200]):
    m1 = numpy.random.rand(i, i)
    m2 = numpy.random.rand(i, 10)
    m3 = numpy.random.rand(10, i)

    t1 = measure_time(
        lambda m1=m1, m2=m2, m3=m3: (m1 @ m2) @ m3, context={}, number=50, repeat=50
    )
    t1["formula"] = "(m1 @ m2) @ m3"
    t1["size"] = i
    obs.append(t1)
    t2 = measure_time(
        lambda m1=m1, m2=m2, m3=m3: m1 @ (m2 @ m3), context={}, number=50, repeat=50
    )
    t2["formula"] = "m1 @ (m2 @ m3)"
    t2["size"] = i
    obs.append(t2)

df = DataFrame(obs)
piv = df.pivot(index="size", columns="formula", values="average")
piv

  0%|          | 0/6 [00:00<?, ?it/s]
 33%|███▎      | 2/6 [00:00<00:00,  6.20it/s]
 50%|█████     | 3/6 [00:00<00:00,  3.43it/s]
 67%|██████▋   | 4/6 [00:01<00:00,  2.21it/s]
 83%|████████▎ | 5/6 [00:03<00:00,  1.07it/s]
100%|██████████| 6/6 [00:07<00:00,  1.84s/it]
100%|██████████| 6/6 [00:07<00:00,  1.17s/it]

formula	(m1 @ m2) @ m3	m1 @ (m2 @ m3)
size
50	0.000009	0.000022
100	0.000024	0.000073
125	0.000045	0.000125
150	0.000057	0.000208
175	0.000084	0.000612
200	0.000102	0.001332

Graph

fig, ax = plt.subplots(1, 2, figsize=(12, 4))
piv.plot(
    logx=True,
    logy=True,
    ax=ax[0],
    title=f"{m1.shape!r} @ {m2.shape!r} @ {m3.shape!r}".replace("200", "size"),
)
piv["ratio"] = piv["m1 @ (m2 @ m3)"] / piv["(m1 @ m2) @ m3"]
piv[["ratio"]].plot(ax=ax[1])

<Axes: xlabel='size'>