Set the differentiation matrix as A to be N+2 x N+2
A = np.zeros((N+2,N+2))
Here we set them using a loop, which must be N+1 long (the size of the off diagonals)
for i in range(N+1):
A[i,i+1] = 1./(2*h)
A[i+1,i] = -1./(2*h)
Does @
mean multiply in this line below?
df_cent = A@f(x) # do the multiplication