code: def run_reaction_production_check(self, model, solver, threshold,
                                        implicit_sinks=True):
          
          prob = solver.create_problem()
  
          # Create flux variables
          v = prob.namespace()
          for reaction_id in model.reactions:
              lower, upper = model.limits[reaction_id]
              v.define([reaction_id], lower=lower, upper=upper)
  
          # Build mass balance constraints
          massbalance_lhs = {compound: 0 for compound in model.compounds}
          for spec, value in iteritems(model.matrix):
              compound, reaction_id = spec
              massbalance_lhs[compound] += v(reaction_id) * value
  
          # Create production variables and apply constraints
          for compound, lhs in iteritems(massbalance_lhs):
              if implicit_sinks:
                  # The constraint is merely >0 meaning that we have implicit
                  # sinks for all compounds.
                  prob.add_linear_constraints(lhs >= 0)
              else:
                  prob.add_linear_constraints(lhs == 0)
  
          confirmed_production = set()
          for reaction in model.reactions:
              if all(c in confirmed_production for c, _ in
                     model.get_reaction_values(reaction)):
                  continue
  
              prob.set_objective(v(reaction))
              for sense in (lp.ObjectiveSense.Maximize,
                            lp.ObjectiveSense.Minimize):
                  try:
                      result = prob.solve(sense)
                  except lp.SolverError as e:
                      self.fail(
                          'Failed to solve for compound, reaction: {}, {}:'
                          ' {}'.format(compound, reaction, e))
  
                  flux = result.get_value(v(reaction))
                  for compound, value in model.get_reaction_values(reaction):
                      if compound in confirmed_production:
                          continue
  
                      production = 0
                      if sense == lp.ObjectiveSense.Maximize and flux > 0:
                          production = float(value) * flux
                      elif sense == lp.ObjectiveSense.Minimize and flux < 0:
                          production = float(value) * flux
  
                      if production >= threshold:
                          confirmed_production.add(compound)
  
          for compound in sorted(model.compounds):
              if compound not in confirmed_production:
                  yield compound
docs: Run reaction production check method.