An efficient algorithm for computing the generalized null space decomposition

The generalized null space decomposition (GNSD) is a unitary reductionof a general matrix $A$ of order $n$ to a block upper triangular formthat reveals the structure of the Jordan blocks of $A$ correspondingto a zero eigenvalue. The reduction was introduced by Kublanovskaya.It was extended first by Ruhe and then by Golub and Wilkinson, whobased the reduction on the singular value decomposition.Unfortunately, if $A$ has large Jordan blocks, the complexity of thesealgorithms can approach the order of $n^4$. This paper presents analternative algorithm, based on repeated updates of a QR decompositionof $A$, that is guaranteed to be of order $n^{3}$. Numericalexperiments confirm the stability of this algorithm, which turns outto produce essentially the same form as that of Golub and Wilkinson.The effect of errors in $A$ on the ability to recover the originalstructure is investigated empirically. Several applications are discussed,including the computation of the Drazin inverse.