def carbonic_ode(S, t, pars):
    """
    Equations from Scott's Thesis page 53
    """
    k = pars[0]  # rate constant / s^-1
    alpha = pars[1]  # H2(16)O / H2(18)O ratio

    S44, S46, S48 = S[0], S[1], S[2]

    dS44 = k * (-2 / 3 * (1 - alpha) * S44 + 1 / 3 * alpha * S46)
    dS46 = k * (2 / 3 * (1 - alpha) * S44 - 1 / 3 * S46 + 2 / 3 * alpha * S48)
    dS48 = k * (1 / 3 * (1 - alpha) * S46 - 2 / 3 * alpha * S48)

    return [dS44, dS46, dS48]