code: def peo(A, p):
      
  
      n = A.size[0]
      assert type(A) == spmatrix, "A must be a sparse matrix"
      assert A.size[1] == n, "A must be a square matrix"
      assert len(p) == n, "length of p must be equal to the order of A" 
      if isinstance(p, list): p = matrix(p)
      
      As = symmetrize(A)
      cp,ri,_ = As.CCS
  
      # compute inverse permutation array
      ip = matrix(0,(n,1))
      ip[p] = matrix(range(n),(n,1))
  
      # test set inclusion
      for k in range(n):
          v = p[k]  # next vertex to be eliminated
  
          # indices of neighbors that correspond to strictly lower triangular elements in reordered pattern
          r = set([rj for rj in ri[cp[v]:cp[v+1]] if ip[rj] > k])  
  
          for rj in r:
              if not r.issubset(set(ri[cp[rj]:cp[rj+1]])): return False
              
      return True
docs: Checks whether an ordering is a perfect elmimination order.
  
      Returns `True` if the permutation :math:`p` is a perfect elimination order
      for a Cholesky factorization :math:`PAP^T = LL^T`. Only the lower
      triangular part of :math:`A` is accessed.
      
      :param A:   :py:class:`spmatrix`
      
      :param p:   :py:class:`matrix` or :class:`list` of length `A.size[0]`