def fd_8(x, v, w, t, a, b):
    """Implements the shooting method to solve linear second order BVPs

    Compute finite difference solution to the BVP but with the 8-order formula!


    t should be passed in as an n element array.   x, v, and w should be
    either n element arrays corresponding to x(t), v(t) and w(t) or
    scalars, in which case an n element array with the given value is
    generated for each of them.

    USAGE:
        u'' = x(t) + v(t) u + w(t) u'
        u = fd(x, v, w, t, a, b)

    INPUT:
        x,v,w - arrays containing x(t), v(t), and w(t) values.  May be
                specified as Python lists, NumPy arrays, or scalars.  In
                each case they are converted to NumPy arrays.
        t     - array of n time values to determine u at
        a     - solution value at the left boundary: a = u(t[0])
        b     - solution value at the right boundary: b = u(t[n-1])

    OUTPUT:
        u     - array of solution function values corresponding to the
                values in the supplied array t.
    """

    # Get the dimension of t and make sure that t is an n-element vector

    if type(t) != np.ndarray:
        if type(t) == list:
            t = np.array(t)
        else:
            t = np.array([ float( t ) ])

    n = len(t)

    # Make sure that x, v, and w are either scalars or n-element vectors.
    # If they are scalars then we create vectors with the scalar value in
    # each position.

    if type(x) == int or type(x) == float:
        x = np.array([float(x)] * n)

    if type(v) == int or type(v) == float:
        v = np.array([float(v)] * n)

    if type(w) == int or type(w) == float:
        w = np.array([float(w)] * n)

    # Compute the stepsize.  It is assumed that all elements in t are
    # equally spaced.

    h = t[1] - t[0];

    if (n > 7):
        # Construct tridiagonal system; boundary conditions appear as first and
        # last equations in system.

        #A = -( 1.0 + w[1:n] * h / 2.0 )

        A = - (490/(180*(h**2))) - v

        A[1:3] = -(2/h**2) - v[1:3]
        A[-3:-1] = -(2/h**2) - v[-3:-1]
        A[0] = A[-1] = 1.0

        B = np.zeros(n-1)
        for i in range(n-1):
            B[i] = 270/(180*(h**2))
        B[0:2] = B[-3:-1] = B[-1] = 1/h**2
        B[-1] = 0.0

        C = B[::-1].copy()

        D = np.zeros(n-2)
        for i in range(n-2):
            D[i] = - (27/(180*(h**2)))
        D[0] = 0.0
        D[-3:-1] = D[-1] = 0.0

        E = D[::-1].copy()

        F = np.zeros(n-3)
        for i in range(n-3):
            F[i] = 2/(180*(h**2))
        F[-3:-1] = F[-1] = 0.0

        G = F[::-1].copy()

        H = x
        H[0] = a
        H[-1] = b
        '''
        # Solve tridiagonal system?? How?