def euler_bernoulli_beam(x, y, EI, f):
    """Euler-Bernoulli Beam Theory defining y'''' for Lu=f.

    This form is used because it is expected for Scipy's solve_ivp method.

    Keyword arguments:
    x -- independent variable
    y -- dependent variable
    EI -- EI(x) parameter
    f -- forcing function f(x)
    """
    return [y[1], y[2], y[3], -1*f/EI]