code: def flux_variability(model, reactions, fixed, tfba, solver):
      
  
      fba = _get_fba_problem(model, tfba, solver)
  
      for reaction_id, value in iteritems(fixed):
          flux = fba.get_flux_var(reaction_id)
          fba.prob.add_linear_constraints(flux >= value)
  
      def min_max_solve(reaction_id):
          for direction in (-1, 1):
              yield fba.flux_bound(reaction_id, direction)
  
      # Solve for each reaction
      for reaction_id in reactions:
          yield reaction_id, tuple(min_max_solve(reaction_id))
docs: Find the variability of each reaction while fixing certain fluxes.
  
      Yields the reaction id, and a tuple of minimum and maximum value for each
      of the given reactions. The fixed reactions are given in a dictionary as
      a reaction id to value mapping.
  
      This is an implementation of flux variability analysis (FVA) as described
      in [Mahadevan03]_.
  
      Args:
          model: MetabolicModel to solve.
          reactions: Reactions on which to report variablity.
          fixed: dict of additional lower bounds on reaction fluxes.
          tfba: If True enable thermodynamic constraints.
          solver: LP solver instance to use.
  
      Returns:
          Iterator over pairs of reaction ID and bounds. Bounds are returned as
          pairs of lower and upper values.