def _run(self, cache_state=True):
    if (self._model is None):
        raise RuntimeError(f'{self.ID} was initialized without a suspended growth model.')
    else:
        isa = isinstance
        inf = self.ins[0]
        Q_in = inf.get_total_flow('m3/d')
        (eff, sludge) = self.outs
        eff.copy_like(inf)
        sludge.copy_like(inf)
        C_in = ((inf.mass / inf.F_vol) * 1000.0)
        cmps = self.components
        C = list(symbols(cmps.IDs))
        if (self._init_Vas is not None):
            V_0 = self._init_Vas
            C_0 = self._init_Cas
        else:
            V_0 = 0
            C_0 = C_in
        n = self._N_layer
        if (self._aeration.count(None) == len(self._aeration)):
            Vmax = self._V
            hj = (self._h / n)
        else:
            Vmax = (self._V * 0.75)
            hj = ((self._h * 0.75) / n)
        T_fill = ((Vmax - V_0) / Q_in)
        T = [(t / 24) for t in self._operation_cycle]
        if (T_fill <= T[0]):
            schedule = ([T_fill, (T[0] - T_fill)] + T[1:4])
            aer = ([self._aeration[0], self._aeration[0]] + list(self._aeration[1:4]))
            fill = ([True] + ([False] * 4))
            V_total = Vmax
        elif (T_fill <= (T[0] + T[1])):
            schedule = ([T[0], (T_fill - T[0]), ((T[0] + T[1]) - T_fill)] + T[2:4])
            aer = ((list(self._aeration[:2]) + [self._aeration[1]]) + list(self._aeration[2:4]))
            fill = (([True] * 2) + ([False] * 3))
            V_total = Vmax
        else:
            schedule = T[:4]
            aer = list(self._aeration[:4])
            fill = (([True] * 2) + ([False] * 2))
            V_total = ((Q_in * (T[0] + T[1])) + V_0)
            hj = (((V_total / self._V) * self._h) / n)
        for i in range(1, len(schedule)):
            if ((fill[(- i)] == fill[((- i) - 1)]) and (aer[(- i)] == aer[((- i) - 1)])):
                schedule[((- i) - 1)] += schedule[(- i)]
                schedule[(- i)] = 0
        t_arr = np.array([])
        y_mat = np.ndarray([])
        for i in range(len(schedule)):
            if (schedule[i] > 0):
                (dC_dt, J_func) = self._compile_dC_dt(V_0, Q_in, C_in, C, fill[i], aer[i])
                if isa(aer[i], (float, int)):
                    C_0[cmps.index(self._DO_ID)] = aer[i]
                sol = solve_ivp(dC_dt, (0, schedule[i]), C_0, method='BDF', jac=J_func)
                C_0 = sol.y.transpose()[(- 1)]
                V_0 += ((Q_in * schedule[i]) * fill[i])
                t_arr = np.concatenate((t_arr, (sol.t + t_arr[(- 1)])))
                y_mat = np.hstack((y_mat, sol.y))
        self._dynamic_composition = np.vstack((t_arr, y_mat)).transpose()
        eff.set_flow((C_0 * eff.F_vol), 'g/hr', self.components.IDs)
        X_0 = eff.get_TSS()
        X_min = (X_0 * self._fns)
        T_settle = T[4]

        def dX_dt(t, X):
            VX = [_settling_flux(x, self._v_max, self._v_max_p, X_min, self._rh, self._rp) for x in X]
            J = [(VX[j] if (X[(j + 1)] <= self._X_t) else min(VX[j], VX[(j + 1)])) for j in range((n - 1))]
            settle_out = np.array((J + [0]))
            settle_in = np.array(([0] + J))
            dXdt = ((settle_in - settle_out) / hj)
            return dXdt
        sol = solve_ivp(dX_dt, (0, T_settle), (np.ones(n) * X_0))
        X = sol.y.transpose()[(- 1)]
        V_eff = min((T[5] * self._Q_e), ((V_total * (n - 1)) / n))
        n_eff = (V_eff / V_total)
        w_eff = np.array((([1] * floor(n_eff)) + [(n_eff - floor(n_eff))]))
        X_eff = np.average(X[:ceil(n_eff)], weights=w_eff)
        eff_mass_flow = ((((((X_eff / X_0) * cmps.x) + (1 - cmps.x)) * C_0) * V_eff) / T[5])
        eff.set_flow(eff_mass_flow, 'g/d', cmps.IDs)
        V_was = min((T[6] * self._Q_WAS), (V_total - V_eff))
        X_as = (((V_total * X_0) - (V_eff * X_eff)) / (V_total - V_eff))
        C_as = ((((X_as / X_0) * cmps.x) + (1 - cmps.x)) * C_0)
        was_mass_flow = ((C_as * V_was) / T[6])
        sludge.set_flow(was_mass_flow, 'g/d', cmps.IDs)
        if self._cache_state:
            self._init_Vas = ((V_total - V_eff) - V_was)
            self._init_Cas = C_as
def _run(self, cache_state=True):
    if (self._model is None):
        raise RuntimeError(f'{self.ID} was initialized without a suspended growth model.')
    else:
        isa = isinstance
        inf = self.ins[0]
        Q_in = inf.get_total_flow('m3/d')
        (eff, sludge) = self.outs
        eff.copy_like(inf)
        sludge.copy_like(inf)
        C_in = ((inf.mass / inf.F_vol) * 1000.0)
        cmps = self.components
        C = list(symbols(cmps.IDs))
        if (self._init_Vas is not None):
            V_0 = self._init_Vas
            C_0 = self._init_Cas
        else:
            V_0 = 0
            C_0 = C_in
        n = self._N_layer
        if (self._aeration.count(None) == len(self._aeration)):
            Vmax = self._V
            hj = (self._h / n)
        else:
            Vmax = (self._V * 0.75)
            hj = ((self._h * 0.75) / n)
        T_fill = ((Vmax - V_0) / Q_in)
        T = [(t / 24) for t in self._operation_cycle]
        if (T_fill <= T[0]):
            schedule = ([T_fill, (T[0] - T_fill)] + T[1:4])
            aer = ([self._aeration[0], self._aeration[0]] + list(self._aeration[1:4]))
            fill = ([True] + ([False] * 4))
            V_total = Vmax
        elif (T_fill <= (T[0] + T[1])):
            schedule = ([T[0], (T_fill - T[0]), ((T[0] + T[1]) - T_fill)] + T[2:4])
            aer = ((list(self._aeration[:2]) + [self._aeration[1]]) + list(self._aeration[2:4]))
            fill = (([True] * 2) + ([False] * 3))
            V_total = Vmax
        else:
            schedule = T[:4]
            aer = list(self._aeration[:4])
            fill = (([True] * 2) + ([False] * 2))
            V_total = ((Q_in * (T[0] + T[1])) + V_0)
            hj = (((V_total / self._V) * self._h) / n)
        for i in range(1, len(schedule)):
            if ((fill[(- i)] == fill[((- i) - 1)]) and (aer[(- i)] == aer[((- i) - 1)])):
                schedule[((- i) - 1)] += schedule[(- i)]
                schedule[(- i)] = 0
        t_arr = np.array([])
        y_mat = np.ndarray([])
        for i in range(len(schedule)):
            if (schedule[i] > 0):
                (dC_dt, J_func) = self._compile_dC_dt(V_0, Q_in, C_in, C, fill[i], aer[i])
                if isa(aer[i], (float, int)):
                    C_0[cmps.index(self._DO_ID)] = aer[i]
                sol = solve_ivp(dC_dt, (0, schedule[i]), C_0, method='BDF', jac=J_func)
                C_0 = sol.y.transpose()[(- 1)]
                V_0 += ((Q_in * schedule[i]) * fill[i])
                t_arr = np.concatenate((t_arr, (sol.t + t_arr[(- 1)])))
                y_mat = np.hstack((y_mat, sol.y))
        self._dynamic_composition = np.vstack((t_arr, y_mat)).transpose()
        eff.set_flow((C_0 * eff.F_vol), 'g/hr', self.components.IDs)
        X_0 = eff.get_TSS()
        X_min = (X_0 * self._fns)
        T_settle = T[4]

        def dX_dt(t, X):
            VX = [_settling_flux(x, self._v_max, self._v_max_p, X_min, self._rh, self._rp) for x in X]
            J = [(VX[j] if (X[(j + 1)] <= self._X_t) else min(VX[j], VX[(j + 1)])) for j in range((n - 1))]
            settle_out = np.array((J + [0]))
            settle_in = np.array(([0] + J))
            dXdt = ((settle_in - settle_out) / hj)
            return dXdt
        sol = solve_ivp(dX_dt, (0, T_settle), (np.ones(n) * X_0))
        X = sol.y.transpose()[(- 1)]
        V_eff = min((T[5] * self._Q_e), ((V_total * (n - 1)) / n))
        n_eff = (V_eff / V_total)
        w_eff = np.array((([1] * floor(n_eff)) + [(n_eff - floor(n_eff))]))
        X_eff = np.average(X[:ceil(n_eff)], weights=w_eff)
        eff_mass_flow = ((((((X_eff / X_0) * cmps.x) + (1 - cmps.x)) * C_0) * V_eff) / T[5])
        eff.set_flow(eff_mass_flow, 'g/d', cmps.IDs)
        V_was = min((T[6] * self._Q_WAS), (V_total - V_eff))
        X_as = (((V_total * X_0) - (V_eff * X_eff)) / (V_total - V_eff))
        C_as = ((((X_as / X_0) * cmps.x) + (1 - cmps.x)) * C_0)
        was_mass_flow = ((C_as * V_was) / T[6])
        sludge.set_flow(was_mass_flow, 'g/d', cmps.IDs)
        if self._cache_state:
            self._init_Vas = ((V_total - V_eff) - V_was)
            self._init_Cas = C_as