def FivePointsDiff(fx, x, h=0.001):
    """
    FivePointsDiff(@fx, x, h);
    Use five points difference to approximatee the derivative of function fx
    in points x, and with step length h
    The function fx must be defined as a function handle with input
    parameter x and the derivative as output parameter
    
    Parameters
    ----------
    fx : function
        A function defined as fx(x)
    x : float, list, numpy.ndarray
        The point(s) of function fx to compute the derivatives
    h : float, optional
        The step size
    
    Returns
    -------
    float, list, numpy.ndarray: The numerical derivatives of fx at x with 
        the same size as x and the type if from fx()
    """
    return (-fx(x+2*h) + 8*fx(x+h) - 8*fx(x-h) + fx(x-2*h)) / (12.0*h)