def _build_rhs(self):

        def ode_func(t, y):
            # TODO: add a cache mechanism
            self[self.settings["quantity"]] = y
            # TransientReactiveTransport._update_A_and_b(self)
            self.apply_callbacks()
            A = self.A.tocsc()
            b = self.b
            V = self.network[self.settings["pore_volume"]]
            return (-A.dot(y) + b) / V  # much faster than A*y

        return ode_func