docs: def run(self):
          """Run GapFill command"""
  
          # Load compound information
          def compound_name(id):
              if id not in self._model.compounds:
                  return id
              return self._model.compounds[id].properties.get('name', id)
  
          # Calculate penalty if penalty file exists
          penalties = {}
          if self._args.penalty is not None:
              for line in self._args.penalty:
                  line, _, comment = line.partition('#')
                  line = line.strip()
                  if line == '':
                      continue
                  rxnid, penalty = line.split(None, 1)
                  penalties[rxnid] = float(penalty)
  
          core = set(self._mm.reactions)
  
          solver = self._get_solver(integer=True)
          default_comp = self._model.default_compartment
          epsilon = self._args.epsilon
          v_max = float(self._model.default_flux_limit)
  
          blocked = set()
          for compound in self._args.compound:
              if compound.compartment is None:
                  compound = compound.in_compartment(default_comp)
              blocked.add(compound)
  
          if len(blocked) > 0:
              logger.info('Unblocking compounds: {}...'.format(
                  ', '.join(text_type(c) for c in sorted(blocked))))
          else:
              logger.info(
                  'Unblocking all compounds in model. Use --compound option to'
                  ' unblock specific compounds.')
              blocked = set(self._mm.compounds)
  
          exclude = set()
          if self._model.biomass_reaction is not None:
              exclude.add(self._model.biomass_reaction)
  
          # Add exchange and transport reactions to database
          model_complete, weights = create_extended_model(
              self._model,
              db_penalty=self._args.db_penalty,
              ex_penalty=self._args.ex_penalty,
              tp_penalty=self._args.tp_penalty,
              penalties=penalties)
  
          implicit_sinks = not self._args.no_implicit_sinks
  
          logger.info('Searching for reactions to fill gaps')
          try:
              added_reactions, no_bounds_reactions = gapfill(
                  model_complete, core, blocked, exclude, solver=solver,
                  epsilon=epsilon, v_max=v_max, weights=weights,
                  implicit_sinks=implicit_sinks,
                  allow_bounds_expansion=self._args.allow_bounds_expansion)
          except GapFillError as e:
              self._log_epsilon_and_fail(epsilon, e)
  
          for reaction_id in sorted(self._mm.reactions):
              rx = self._mm.get_reaction(reaction_id)
              rxt = rx.translated_compounds(compound_name)
              print('{}\t{}\t{}\t{}'.format(reaction_id, 'Model', 0, rxt))
  
          for rxnid in sorted(added_reactions):
              rx = model_complete.get_reaction(rxnid)
              rxt = rx.translated_compounds(compound_name)
              print('{}\t{}\t{}\t{}'.format(
                  rxnid, 'Add', weights.get(rxnid, 1), rxt))
  
          for rxnid in sorted(no_bounds_reactions):
              rx = model_complete.get_reaction(rxnid)
              rxt = rx.translated_compounds(compound_name)
              print('{}\t{}\t{}\t{}'.format(
                  rxnid, 'Remove bounds', weights.get(rxnid, 1), rxt))