[ivansavz]

Since y'all are code-literate, here the SymPy function for finding the CR-decomposition of any matrix A:

  def crd(A):
      """
      Computes the CR decomposition of the matrix A.
      """
      rrefA, licols = A.rref()  # compute RREF(A)
      C = A[:, licols]          # linearly indep. cols of A
      r = len(licols)           # = rank(A)
      R = rrefA[0:r, :]         # non-zero rows in RREF(A)
      return C, R

Test to check it works: https://live.sympy.org/?evaluate=A%20%3D%20Matrix(%5B%0A%20%...