from sympy import Matrix, Poly, pprint
from sympy.abc import a, t

n = 4

M = Matrix.zeros(n + 1, n + 1)

for j in range(n + 1):
    p = t**j
    new_p = p.subs(t, t - a).expand()
    coefficients = list(reversed(Poly(new_p, t).all_coeffs()))
    for i in range(len(coefficients)):
        M[i, j] = coefficients[i]

print()
pprint(M)
print()