def dynamics(state, time, config, intervention=None):
    """State dynamics for the Chen et al. opioid epidemic model.

    Parameters
    ----------
        state: np.ndarray or tuple or list
            Instantenous system state
        time: float
            Time to evaluate the derivative
        config: whynot.simulators.opioid.Config
            Configuration parameters for the dynamics
        intervention: whynot.simulators.opioid.Intervention
            (Optional) Parameters specifying when/how to update the dynamics

    Returns
    -------
        ds_dt: np.ndarray
            Derivative of the state with respect to time, evaluated at state and time.

    """
    if intervention:
        config = config.update(intervention)
    nonmedical_users, oud_users, illicit_users = state
    nonmedical_increase = config.nonmedical_incidence[time]
    nonmedical_decrease = nonmedical_users * (config.nonmedical_exit +
        config.nonmedical_overdose + config.nonmedical_to_oud + config.
        nonmedical_to_illicit)
    nonmedical_delta = nonmedical_increase - nonmedical_decrease
    oud_increase = config.nonmedical_to_oud * nonmedical_users
    oud_decrease = oud_users * (config.oud_exit[time] + config.oud_overdose
        [time] + config.oud_to_illicit[time])
    oud_delta = oud_increase - oud_decrease
    illicit_increase = config.illicit_incidence[time] + (config.
        oud_to_illicit[time] * oud_users + config.nonmedical_to_illicit *
        nonmedical_users)
    illicit_decrease = illicit_users * (config.illicit_exit + config.
        illicit_overdose[time])
    illicit_delta = illicit_increase - illicit_decrease
    return [nonmedical_delta, oud_delta, illicit_delta]
