Numpy Matrix Difference row by row into 3D tensor

Question

I am looking for a fast way to compute a 3D np.ndarray such that given X.shape == (m,p) and X2.shape == (n,p),

tau[i,j] = X[i] - X2[j] with tau.shape == (m,n,p)

I have the following iterative (read, slow) method thus far:

for i in range(X.shape[0]):
    for j in range(X2.shape[0]):
        tau[i,j] = X[i] - X2[j]

This works but I'd like to know if there is a faster way to do this using some broadcasting trick or something.

Thanks!