Matrix Multiplication in NumPy: matmul, @ and np.dot

In NumPy, np.matmul(a, b) and the @ operator compute the matrix product of two arrays, while * and np.multiply multiply element by element. Mixing the two up is the most common cause of wrong shapes and silently wrong numbers. This guide shows which one to use, how shapes are handled, and how to invert a matrix with np.linalg.inv and check the result. All behaviour below is taken from the NumPy v2.5 reference manual.

What matrix multiplication is used for in research

Matrix products appear wherever data are stored as arrays: applying a weight or transformation matrix to many observations, combining a design matrix with coefficients in a linear model, or working with covariance matrices (see Mahalanobis distance and cluster analysis in Python). The three operations covered here are the matrix product, the element-wise product, and the matrix inverse.

Install or access

NumPy’s install page lists pip install numpy (ideally inside a virtual environment created with python -m venv my-env), conda install numpy for conda users, and uv pip install numpy for uv users. Source: NumPy install page, checked 2026-10-29. Then:

import numpy as np

For creating test arrays, see NumPy arrays and np.zeros.

Step-by-step: multiply matrices in NumPy

Which operation do you need?

You want	Write	What it does (NumPy docs)
Matrix product	a @ b or np.matmul(a, b)	@ implements the PEP 465 semantics and, on ndarrays, is shorthand for np.matmul. It uses optimized BLAS when possible.
Element-wise product	a * b or np.multiply(a, b)	Multiplies arguments element-wise; shapes must be broadcastable. * is shorthand for np.multiply on ndarrays.
Product of stacks of matrices	np.matmul(a, b)	Treats N-D arguments as a stack of matrices in the last two indices and broadcasts. np.dot does not behave this way.
Scalar times matrix	3 * a	np.matmul does not allow scalar arguments; use *.

1. Build two 2-D arrays and take the matrix product.

   a = np.array([[1, 0], [0, 1]])
   b = np.array([[4, 1], [2, 2]])
   np.matmul(a, b)   # array([[4, 1], [2, 2]])
   a @ b             # same result

2. Use a 1-D argument. If the first argument is 1-D, a 1 is prepended to its shape and removed after the product; if the second is 1-D, a 1 is appended and removed afterwards. With b = [1, 2] and a above, both np.matmul(a, b) and np.matmul(b, a) return [1, 2]. The result is a scalar only when both arguments are 1-D.

3. Multiply stacks of matrices.

   a = np.arange(16).reshape((2, 2, 4))
   b = np.arange(16).reshape((2, 4, 2))
   np.matmul(a, b).shape   # (2, 2, 2)
   np.matmul(a, b)[0, 1, 1]   # 98

4. Multiply element by element when that is what you mean.

   x1 = np.arange(9.0).reshape((3, 3))
   x2 = np.arange(3.0)
   np.multiply(x1, x2)   # [[0, 1, 4], [0, 4, 10], [0, 7, 16]]; same as x1 * x2

matmul versus dot for N-D arrays

For two stacks a = np.ones([9, 5, 7, 4]) and c = np.ones([9, 5, 4, 3]), the NumPy manual shows np.dot(a, c).shape is (9, 5, 7, 9, 5, 3) while np.matmul(a, c).shape is (9, 5, 7, 3). For batches of matrices, matmul is almost always what you want.

Do not use np.matrix to get * as a matrix product

For np.matrix, * is matrix multiplication, but for ordinary ndarrays it is element-wise. The NumPy manual states: “It is no longer recommended to use this class, even for linear algebra. Instead use regular arrays. The class may be removed in the future.” Use regular arrays and @.

Matrix inverse in NumPy

np.linalg.inv(a) takes an array of shape (..., M, M) and returns ainv such that a @ ainv = ainv @ a = eye(a.shape[0]). Stacks of matrices can be inverted in one call.

a = np.array([[1., 2.], [3., 4.]])
ainv = np.linalg.inv(a)    # [[-2., 1.], [1.5, -0.5]]
np.allclose(a @ ainv, np.eye(2))   # True

Checking a @ ainv against the identity with np.allclose is a quick way to confirm the inverse is usable.

Output and how to read it

• Shape rule: for 2-D inputs, an (m, k) array times a (k, n) array gives (m, n). The last dimension of the first argument must match the second-to-last dimension of the second.
• 1-D inputs: the promoted dimension is removed, so a matrix times a vector returns a 1-D result.
• Complex vectors: np.matmul([2j, 3j], [2j, 3j]) returns (-13+0j); there is no conjugation. x1 @ x2 gives the same value.
• Inverse accuracy: for a near-singular matrix, inv may return no error but very large entries. The manual’s example [[2, 4, 6], [2, 0, 2], [6, 8, 14]] returns entries around 1e15, and its condition number is about 8.66e+17 (the value may vary). The manual gives this rule of thumb: “if the condition number cond(a) = 10**k, then you may lose up to k digits of accuracy on top of what would be lost to the numerical method due to loss of precision from arithmetic methods.” Check np.linalg.cond(a) before trusting an inverse. See also np.linalg.norm for other linear algebra helpers.

Troubleshooting

ValueError: matmul: Input operand 1 does not have enough dimensions

You passed a scalar. np.matmul([1, 2], 3) raises this error because scalars are not allowed. Use * for scalar multiplication.

ValueError about mismatched dimensions

The last dimension of the first argument does not equal the second-to-last dimension of the second. Print both .shape values and transpose or reshape one argument.

Numbers look wrong but there is no error

Check whether you used * where you meant @. On ndarrays, * multiplies element-wise with broadcasting, so compatible shapes will not raise an error.

numpy.linalg.LinAlgError from inv

inv raises LinAlgError if the matrix is not square or inversion fails, including when it is singular. For an ill-conditioned matrix the error may or may not be raised and the result may be inaccurate, so inspect the condition number.

Alternatives

Function	Documented purpose
scipy.linalg.inv	“Compute the inverse of a matrix”; raises LinAlgError if the matrix is singular (SciPy reference). The NumPy inv page lists it under “See also”.
numpy.linalg.cond	Condition number, listed under “See also” on the NumPy inv page, for judging inversion accuracy.
numpy.linalg.solve	Computes the “exact” solution x of a well-determined linear system (NumPy reference).

FAQ

How do I multiply matrices in NumPy?

Use a @ b or np.matmul(a, b) for the matrix product. The inner dimensions must match.

What does numpy multiply do?

np.multiply(x1, x2) multiplies the arguments element-wise and is equivalent to x1 * x2. Shapes must be broadcastable.

How do I take the matrix inverse in NumPy?

Call np.linalg.inv(a) on a square array, then verify with np.allclose(a @ ainv, np.eye(n)). A singular matrix raises LinAlgError.

Sources

• NumPy v2.5 manual, numpy.matmul: https://numpy.org/doc/stable/reference/generated/numpy.matmul.html (checked 2026-10-29)
• NumPy v2.5 manual, numpy.multiply: https://numpy.org/doc/stable/reference/generated/numpy.multiply.html (checked 2026-10-29)
• NumPy v2.5 manual, numpy.linalg.inv: https://numpy.org/doc/stable/reference/generated/numpy.linalg.inv.html (checked 2026-10-29)
• NumPy v2.5 manual, numpy.matrix: https://numpy.org/doc/stable/reference/generated/numpy.matrix.html (checked 2026-10-29)
• NumPy install guide: https://numpy.org/install (checked 2026-10-29)
• SciPy reference, scipy.linalg.inv: https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.inv.html (checked 2026-10-29)