import sympy as sy
import numpy as np
from sympy.abc import x, y
z = ((x + y) ** 2).expand()
print(z)
X = np.arange(6).reshape(2, 3)
Y = np.arange(1, 7).reshape(2, 3)
f = sy.lambdify((x, y), z, 'numpy')
print(f(X, Y))
assert np.allclose(f(X, Y), X ** 2 + 2 * X * Y + Y ** 2)
How to substitute symbol for matrix using symPy and numPy