docs: def rdl_decomposition_nrev(T, k, norm='standard', ncv=None):
      r"""Compute the decomposition into left and right eigenvectors.
  
      Parameters
      ----------
      T : sparse matrix
          Transition matrix
      k : int
          Number of eigenvector/eigenvalue pairs
      norm: {'standard', 'reversible'}
          standard: (L'R) = Id, L[:,0] is a probability distribution,
              the stationary distribution mu of T. Right eigenvectors
              R have a 2-norm of 1.
          reversible: R and L are related via L=L[:,0]*R.
      ncv : int (optional)
          The number of Lanczos vectors generated, `ncv` must be greater than k;
          it is recommended that ncv > 2*k
  
      Returns
      -------
      R : (M, M) ndarray
          The normalized ("unit length") right eigenvectors, such that the
          column ``R[:,i]`` is the right eigenvector corresponding to the eigenvalue
          ``w[i]``, ``dot(T,R[:,i])``=``w[i]*R[:,i]``
      D : (M, M) ndarray
          A diagonal matrix containing the eigenvalues, each repeated
          according to its multiplicity
      L : (M, M) ndarray
          The normalized (with respect to `R`) left eigenvectors, such that the
          row ``L[i, :]`` is the left eigenvector corresponding to the eigenvalue
          ``w[i]``, ``dot(L[i, :], T)``=``w[i]*L[i, :]``
  
      """
      # Standard norm: Euclidean norm is 1 for r and LR = I.
      if norm == 'standard':
          v, R = scipy.sparse.linalg.eigs(T, k=k, which='LM', ncv=ncv)
          r, L = scipy.sparse.linalg.eigs(T.transpose(), k=k, which='LM', ncv=ncv)
  
          """Sort right eigenvectors"""
          ind = np.argsort(np.abs(v))[::-1]
          v = v[ind]
          R = R[:, ind]
  
          """Sort left eigenvectors"""
          ind = np.argsort(np.abs(r))[::-1]
          r = r[ind]
          L = L[:, ind]
  
          """l1-normalization of L[:, 0]"""
          L[:, 0] = L[:, 0] / np.sum(L[:, 0])
  
          """Standard normalization L'R=Id"""
          ov = np.diag(np.dot(np.transpose(L), R))
          R = R / ov[np.newaxis, :]
  
          """Diagonal matrix with eigenvalues"""
          D = np.diag(v)
  
          return R, D, np.transpose(L)
  
      # Reversible norm:
      elif norm == 'reversible':
          v, R = scipy.sparse.linalg.eigs(T, k=k, which='LM', ncv=ncv)
          mu = stationary_distribution(T)
  
          """Sort right eigenvectors"""
          ind = np.argsort(np.abs(v))[::-1]
          v = v[ind]
          R = R[:, ind]
  
          """Ensure that R[:,0] is positive"""
          R[:, 0] = R[:, 0] / np.sign(R[0, 0])
  
          """Diagonal matrix with eigenvalues"""
          D = np.diag(v)
  
          """Compute left eigenvectors from right ones"""
          L = mu[:, np.newaxis] * R
  
          """Compute overlap"""
          s = np.diag(np.dot(np.transpose(L), R))
  
          """Renormalize left-and right eigenvectors to ensure L'R=Id"""
          R = R / np.sqrt(s[np.newaxis, :])
          L = L / np.sqrt(s[np.newaxis, :])
  
          return R, D, np.transpose(L)
      else:
          raise ValueError("Keyword 'norm' has to be either 'standard' or 'reversible'")