Thin QR Matrices in Eigen

Alec Jacobson

October 18, 2024

weblog/

The QR Module of Eigen is full-featured and usually you don't really need/want to realize the $Q$, $R$ matrices (e.g., for solving linear systems). Nevertheless, it's pretty frustraing how confusing it is to extract these as matrices. Here's a small example of how to implement matlab's [Q,R,P] = qr(A,'econ') "thin"/"economic"/"reduced" QR decomposition in Eigen.

I specifically use fixed sized matrices because it's extra-annoying to convert the internal QR types to the right sizes.

#include <Eigen/Core>
#include <Eigen/QR>
#include <igl/matlab_format.h>
#include <iostream>
int main()
{
  Eigen::Matrix<double,4,3> A(4,3);
  A<<
    0,3,8,
    9,2,3,
    5,7,3,
    9,1,5;
  auto qr = A.colPivHouseholderQr();
  Eigen::Matrix<double,4,3> Q = Eigen::Matrix<double,4,4>(qr.householderQ()).leftCols(A.cols());
  Eigen::Matrix<double,3,3> R = Eigen::Matrix<double,4,3>(qr.matrixQR().triangularView<Eigen::Upper>()).topLeftCorner(A.cols(),A.cols()).eval();
  Eigen::Matrix<double,3,3> P = Eigen::Matrix<double,3,3>(qr.colsPermutation());

  std::cout<<igl::matlab_format(A,"A")<<std::endl;
  std::cout<<igl::matlab_format( Q ,"Q")<<std::endl;
  std::cout<<igl::matlab_format( R ,"R")<<std::endl;
  std::cout<<igl::matlab_format( P ,"P")<<std::endl;
  std::cout<<igl::matlab_format( (Q*R*P.transpose()) ,"QRPᵀ")<<std::endl;
}

This should output:

A = [
  0 3 8
  9 2 3
  5 7 3
  9 1 5
];
Q = [
                   0  0.980845491305081 0.0474282849225397
  -0.658145181714418  -0.14555327477923 0.0891720659621116
  -0.365636212063565 0.0826113181179414 -0.904869663552425
  -0.658145181714418 0.0996580980470403  0.413533302678124
];
R = [
  -13.6747943311773 -6.36207008990604 -4.53388902958821
                  0  8.15622855069198  3.32936724922941
                  0                 0 -5.59992535549701
];
P = [
  1 0 0
  0 0 1
  0 1 0
];
QRPᵀ = [
  0 3 8
  9 2 3
  5 7 3
  9 1 5
];