# Range of a matrix from its complete orthogonal decomposition

In this StackOverflow answer, @Gokul has shown how to get a basis of the kernel of a matrix with the help of the 'Eigen' function CompleteOrthogonalDecomposition. The complete orthogonal decomposition of a matrix is related to the so-called URV decomposition of this matrix. I don't know this topic. I read somewhere that the URV decomposition also allows to get the image (range, or span) of the matrix. I would like to know how to use Eigen::CompleteOrthogonalDecomposition to get a basis of the image of a matrix.

• Basically, the image of the decomposed matrix is spanned by the rows of U that correspond to non-vanishing singular values. But, as you write, if you don't know the SVD, you should seriously read it up. It's the single most important decomposition in numerical linear algebra. – davidhigh Nov 17 '20 at 19:12
• @davidhigh OP is asking about URV, not SVD. It is a faster decomposition. – Federico Poloni Nov 18 '20 at 14:32
• Ok, thanks for pointing. Still, I'd be interested in the use case where the higher performance of the URV decomposition really pays for the much more complex implementation in the answer below. My first choice here would always be the SVD. I just mention it because I'm not sure if this is transparent to the OP, as he was asking explicitly about the range and not about performance. – davidhigh Nov 18 '20 at 23:19
• @davidhigh My underlying motivation was to get both the kernel and the range from a single shot. – Stéphane Laurent Nov 18 '20 at 23:30
• Then habe a look at the SVD: it is perfectly appropriate for your task and offers much less complex coding as compared to your answer. Given the decomposition UDV^+, for the range just take the rows of U corresponding to non-vanishing singular values, and for the Kernel the rows of V corresponding to vanishing singular values. Or better read it up anywhere over the net. – davidhigh Nov 19 '20 at 21:59

I get it.

template <typename Number> // e.g. double or complex
Eigen::Matrix<Number, Eigen::Dynamic, Eigen::Dynamic> image_COD(
const Eigen::Matrix<Number, Eigen::Dynamic, Eigen::Dynamic>& M) {
Eigen::CompleteOrthogonalDecomposition<
Eigen::Matrix<Number, Eigen::Dynamic, Eigen::Dynamic>> cod(M);
const Eigen::Matrix<Number, Eigen::Dynamic, Eigen::Dynamic> Q =
cod.householderQ();
return Q.leftCols(cod.rank());
}


It is slightly slower than the QR method, but has the advantage that one can get both the kernel and the range from this decomposition. The LU method is faster but does not return orthonormal bases.