import numpy as np

def f(x):
    return x**2 + 1

dx = 0.1
x_array = np.arange(-2, 2+dx, dx)
y_array = f(x_array)  # Vectorized calculation of y