E = 0
for i = 1,...,N:
    E_new = | BS(x_i) - u_i |
    if E_new >= E
        E = E_new
    else:
        E = E
    end

# Compute the Error
E = np.zeros((N+1, M+1)) # so that E has NxM entries--one for each time-step and node in space
for m in range(1, M+1):
    for j in range(0, N+1):
        E_new = np.absolute(BS[j,m] - u[j,m])
        if (np.any(E_new >= E)):
            E = E_new
        else:
            E = E

print 'Error = ', E