Using QR algorithm to compute the SVD of a matrix How to use the QR algorithm to compute the SVD of a matrix $X\in R^{m\times n}$? Is there any algorithm for doing that?
 A: Quoting the wikipedia article on SVD:


*

*The left-singular vectors of $M$ (i.e. the columns of $U$) are eigenvectors of $MM^T$.

*The right-singular vectors of $M$ (i.e. the columns of $V$) are eigenvectors of $M^TM$.

*The non-zero singular values of $M$ (i.e. the non-zero diagonal elements of $\Sigma$) are the square roots of the non-zero eigenvalues of both $M^TM$ and $MM^T$.


The QR algorithm finds eigenvalues and eigenvectors of square matrices. $M^TM$ and $MM^T$ are square matrices. There might be a better way, I don't know, but this is the naive, obvious way.
A: The SVD can be obtained by computing the eigenvalue decomposition of the symmetric matrix
$$\begin{align}
\begin{bmatrix}
0&X\\X^T&0
\end{bmatrix}
&=
\begin{bmatrix}
U&0\\0&V
\end{bmatrix}⋅
\begin{bmatrix}
0&Σ\\Σ^T&0
\end{bmatrix}\cdot
\begin{bmatrix}
U&0\\0&V
\end{bmatrix}^T
\\&=
\frac1{\sqrt2}
\begin{bmatrix}
U&-U\\V&V
\end{bmatrix}⋅
\begin{bmatrix}
Σ&0\\0&-Σ
\end{bmatrix}⋅
\frac1{\sqrt2}
\begin{bmatrix}
U&-U\\V&V
\end{bmatrix}^T
\end{align}$$
The eigenvectors have the form $\begin{bmatrix}\pm u_k \\ v_k\end{bmatrix}$ with $u_k$ and $v_k$ being the left and right singular vectors, with eigenvalues $\pm\sigma_k$ and some zeros to fill the dimensions.
Since the Hessenberg form of symmetric matrices is tridiagonal, heavy simplifications are possible. These simplifications lead directly to the Golub-Kahan algorithm.
A: My answer to this is the following (see the MATLAB snippet below). You make QR-decomposition for A and then repetitively take the R matrix, transpose it and apply QR-decomposition to R'. Evidently, QR applied to the upper-triangular R gives the matrix R again producing nothing new but if you apply QR to the lower-triangular R' and keep doing it again and again you'll see that the resulting R will converge to a diagonal matrix with singular values on its diagonal. Whereas if you multiply all the orthogonal Q-matrices obtained in the process you'll get U and V matrices. See the code example:
   N = 15; % number of iterations
   R1 = A;
   U = eye(size(A));
   V = eye(size(A));
   for i = 1 : N
      [Q1, R1] = qr(R1);
      [Q2, R2] = qr(R1');
      R1 = R2';
      U = U * Q1;
      V = V * Q2;
   end
   S = R1;

This is a manifestation of the known QR algorithm for finding eigenvalues of a symmetric matrix but applied to a more general situation. Notice that here the matrix A need nor be symmetric and still we obtain both S and U and V together.
