np.einsum matmul is 6-9x slower than @ unless optimize=True

np.einsum defaults to optimize=False, and then even a plain two-operand matrix product runs NumPy's own sum-of-products loop instead of BLAS. NumPy 2.5.4, OpenBLAS pinned to 1 thread, best of 5 on a loaded 8-core box: einsum('ij,jk->ik') took 0.22 / 5.46 / 60.2 ms for n = 100 / 300 / 600 vs 0.03 / 0.75 / 6.9 ms for A@B (6-9x); optimize=True brought it to 0.05 / 1.03 / 6.4 ms. Batched 200x64x64 'bij,bjk->bik': 15.1 ms vs 2.1 ms for matmul, 2.4 ms with optimize=True. Results agree to 2e-13. The opposite holds for contractions BLAS cannot express: row-wise dot einsum('ij,ij->i') was 0.09 ms vs 0.33 ms for (x*x).sum(1), because it avoids the temporary. With multithreaded BLAS the matmul gap would be larger; not measured here.