docs: def _time_dependent_to_qutip(
          op, full_space=None, time_symbol=symbols("t", real=True),
          convert_as='pyfunc'):
      """Convert a possiblty time-dependent operator into the nested-list
      structure required by QuTiP"""
      if full_space is None:
          full_space = op.space
      if time_symbol in op.free_symbols:
          op = op.expand()
          if isinstance(op, OperatorPlus):
              result = []
              for o in op.operands:
                  if time_symbol not in o.free_symbols:
                      if len(result) == 0:
                          result.append(convert_to_qutip(o,
                                                         full_space=full_space))
                      else:
                          result[0] += convert_to_qutip(o, full_space=full_space)
              for o in op.operands:
                  if time_symbol in o.free_symbols:
                      result.append(_time_dependent_to_qutip(o, full_space,
                                    time_symbol, convert_as))
              return result
          elif (
                  isinstance(op, ScalarTimesOperator) and
                  isinstance(op.coeff, ScalarValue)):
              if convert_as == 'pyfunc':
                  func_no_args = lambdify(time_symbol, op.coeff.val)
                  if {time_symbol, } == op.coeff.free_symbols:
                      def func(t, args):
                          # args are ignored for increased efficiency, since we
                          # know there are no free symbols except t
                          return func_no_args(t)
                  else:
                      def func(t, args):
                          return func_no_args(t).subs(args)
                  coeff = func
              elif convert_as == 'str':
                  # a bit of a hack to replace imaginary unit
                  # TODO: we can probably use one of the sympy code generation
                  # routines, or lambdify with 'numexpr' to implement this in a
                  # more robust way
                  coeff = re.sub("I", "(1.0j)", str(op.coeff.val))
              else:
                  raise ValueError(("Invalid value '%s' for `convert_as`, must "
                                    "be one of 'str', 'pyfunc'") % convert_as)
              return [convert_to_qutip(op.term, full_space), coeff]
          else:
              raise ValueError("op cannot be expressed in qutip. It must have "
                               "the structure op = sum_i f_i(t) * op_i")
      else:
          return convert_to_qutip(op, full_space=full_space)