code: def flux_consumers(F, rtol=1e-05, atol=1e-12):
      r
      # can be used with sparse or dense
      n = np.shape(F)[0]
      influxes = np.array(np.sum(F, axis=0)).flatten()  # all that flows in
      outfluxes = np.array(np.sum(F, axis=1)).flatten()  # all that flows out
      # net in flux absolute
      con_abs = np.maximum(influxes - outfluxes, np.zeros(n))
      # net in flux relative
      con_rel = con_abs / (np.maximum(outfluxes, influxes))
      # return all indexes that are produces in terms of absolute and relative tolerance
      return list(np.where((con_abs > atol) * (con_rel > rtol))[0])
docs: r"""Return indexes of states that are net flux producers.
  
      Parameters
      ----------
      F : (n, n) ndarray
          Matrix of flux values between pairs of states.
      rtol : float
          relative tolerance. fulfilled if max(outflux-influx, 0) / max(outflux,influx) < rtol
      atol : float
          absolute tolerance. fulfilled if max(outflux-influx, 0) < atol
  
      Returns
      -------
      producers : (n) ndarray of int
          indexes of states that are net flux producers. May include "dirty" producers, i.e.
          states that have influx but still produce more outflux and thereby violate flux
          conservation.