Resons for potentially bad performance of numpy?

I have been running into performance issues with numpy. I tested different ways to compute a dot product, but similar performance issues have been seen with other numpy functions.

import timeit
import numpy as np
import tensorflow as tf

def np_dot(U, X):
    tic = timeit.default_timer()
    X_bar = np.dot(U.T, X)
    toc = timeit.default_timer()
    print("np.dot took {:.4f} s".format(toc-tic))
    return X_bar

def np_at(U, X):
    tic = timeit.default_timer()
    X_bar = U.T @ X
    toc = timeit.default_timer()
    print("np @ took {:.4f} s".format(toc-tic))
    return X_bar

def np_matmul(U, X):
    tic = timeit.default_timer()
    X_bar = np.matmul(U.T, X)
    toc = timeit.default_timer()
    print("np.matmul took {:.4f} s".format(toc-tic))
    return X_bar

def np_einsum(U, X):
    tic = timeit.default_timer()
    X_bar = np.einsum('ij,jk', U.T, X)
    toc = timeit.default_timer()
    print("np.einsum took {:.4f} s".format(toc-tic))
    return X_bar

def tf_matmul(U, X):
    tic = timeit.default_timer()
    X_bar = tf.matmul(U.T, X)
    toc = timeit.default_timer()
    print("tf.matmul took {:.4f} s".format(toc-tic))
    return X_bar


if __name__ == "__main__":
    print("is_gpu_available?", tf.test.is_gpu_available(), tf.test.is_built_with_cuda())
    print("numpy version:", np.__version__)
    print("tensorflow version:", tf.__version__)
    
    M = N = 1000
    X = np.random.rand(M, N)
    U = np.random.rand(M, N)
    X_bar = np_dot(U, X)
    X_bar0 = np_at(U, X)
    X_bar1 = np_matmul(U, X)
    X_bar2 = np_einsum(U, X)
    X_bar3 = tf_matmul(U, X)
    print(np.allclose(X_bar, X_bar0))
    print(np.allclose(X_bar, X_bar1))
    print(np.allclose(X_bar, X_bar2))
    print(np.allclose(X_bar, X_bar3))

The output on my laptop (with Windows 10 and python 3.8) is:

is_gpu_available? False False
numpy version: 1.19.2
tensorflow version: 2.2.0
np.dot took 12.5368 s
np @ took 14.3303 s
np.matmul took 14.0634 s
np.einsum took 0.5937 s
tf.matmul took 0.0263 s
True
True
True
True

meaning the standard @ is 545 times slower than than it could be! On a friends laptop the result look more like expected (with Ubuntu and python 3.6.9):

is_gpu_available? False False
numpy version: 1.19.5
tensorflow version: 1.14.0
np.dot took 0.0192 s
np @ took 0.0196 s
np.matmul took 0.0208 s
np.einsum took 0.3893 s
tf.matmul took 0.0280 s
True
True
True
True

Interestingly, numpys einsum does not seem to be affected by whatever is going wrong on my laptop and neither is tensorflow.

What is going on here? What are the potential reasons for such a huge performance difference?

edit:

Always make sure numpy is using some BLAS library. Read more about it here: Boosting numpy: Why BLAS Matters.

In my case numpy.show_config() returns the following:

blas_mkl_info:
  NOT AVAILABLE
blis_info:
  NOT AVAILABLE
openblas_info:
  NOT AVAILABLE
atlas_3_10_blas_threads_info:
  NOT AVAILABLE
atlas_3_10_blas_info:
  NOT AVAILABLE
atlas_blas_threads_info:
  NOT AVAILABLE
atlas_blas_info:
  NOT AVAILABLE
accelerate_info:
  NOT AVAILABLE
blas_info:
  NOT AVAILABLE
blas_src_info:
  NOT AVAILABLE
blas_opt_info:
  NOT AVAILABLE
lapack_mkl_info:
  NOT AVAILABLE
openblas_lapack_info:
  NOT AVAILABLE
openblas_clapack_info:
  NOT AVAILABLE
flame_info:
  NOT AVAILABLE
atlas_3_10_threads_info:
  NOT AVAILABLE
atlas_3_10_info:
  NOT AVAILABLE
atlas_threads_info:
  NOT AVAILABLE
atlas_info:
  NOT AVAILABLE
lapack_info:
  NOT AVAILABLE
lapack_src_info:
  NOT AVAILABLE
lapack_opt_info:
  NOT AVAILABLE
numpy_linalg_lapack_lite:
    language = c
    define_macros = [('HAVE_BLAS_ILP64', None), ('BLAS_SYMBOL_SUFFIX', '64_')]