docs: def implicit_step(self):
          """
          Integrate one step using trapezoidal method. Sets convergence and niter flags.
  
          Returns
          -------
          None
          """
          config = self.config
          system = self.system
          dae = self.system.dae
  
          # constant short names
          In = spdiag([1] * dae.n)
          h = self.h
  
          while self.err > config.tol and self.niter < config.maxit:
              if self.t - self.t_jac >= 5:
                  dae.rebuild = True
                  self.t_jac = self.t
              elif self.niter > 4:
                  dae.rebuild = True
              elif dae.factorize:
                  dae.rebuild = True
  
              # rebuild Jacobian
              if dae.rebuild:
                  exec(system.call.int)
                  dae.rebuild = False
              else:
                  exec(system.call.int_fg)
  
              # complete Jacobian matrix dae.Ac
              if config.method == 'euler':
                  dae.Ac = sparse(
                      [[In - h * dae.Fx, dae.Gx], [-h * dae.Fy, dae.Gy]],
                      'd')
                  dae.q = dae.x - self.x0 - h * dae.f
  
              elif config.method == 'trapezoidal':
                  dae.Ac = sparse([[In - h * 0.5 * dae.Fx, dae.Gx],
                                   [-h * 0.5 * dae.Fy, dae.Gy]], 'd')
                  dae.q = dae.x - self.x0 - h * 0.5 * (dae.f + self.f0)
  
              # windup limiters
              dae.reset_Ac()
  
              if dae.factorize:
                  self.F = self.solver.symbolic(dae.Ac)
                  dae.factorize = False
              self.inc = -matrix([dae.q, dae.g])
  
              try:
                  N = self.solver.numeric(dae.Ac, self.F)
                  self.solver.solve(dae.Ac, self.F, N, self.inc)
              except ArithmeticError:
                  logger.error('Singular matrix')
                  dae.check_diag(dae.Gy, 'unamey')
                  dae.check_diag(dae.Fx, 'unamex')
                  # force quit
                  self.niter = config.maxit + 1
                  break
              except ValueError:
                  logger.warning('Unexpected symbolic factorization')
                  dae.factorize = True
                  continue
              else:
                  inc_x = self.inc[:dae.n]
                  inc_y = self.inc[dae.n:dae.m + dae.n]
                  dae.x += inc_x
                  dae.y += inc_y
  
              self.err = max(abs(self.inc))
              if np.isnan(self.inc).any():
                  logger.error('Iteration error: NaN detected.')
                  self.niter = config.maxit + 1
                  break
  
              self.niter += 1
  
          if self.niter <= config.maxit:
              self.convergence = True