docs: def get_residual(self, z, compute_norm=False):
          r'''Compute residual.
  
          For a given :math:`z\in\mathbb{C}^N`, the residual
  
          .. math::
  
            r = M M_l ( b - A z )
  
          is computed. If ``compute_norm == True``, then also the absolute
          residual norm
  
          .. math::
  
            \| M M_l (b-Az)\|_{M^{-1}}
  
          is computed.
  
          :param z: approximate solution with ``z.shape == (N, 1)``.
          :param compute_norm: (bool, optional) pass ``True`` if also the norm
            of the residual should be computed.
          '''
          if z is None:
              if compute_norm:
                  return self.MMlb, self.Mlb, self.MMlb_norm
              return self.MMlb, self.Mlb
          r = self.b - self.A*z
          Mlr = self.Ml*r
          MMlr = self.M*Mlr
          if compute_norm:
              return MMlr, Mlr, utils.norm(Mlr, MMlr, ip_B=self.ip_B)
          return MMlr, Mlr