Let *A* be a general real *m*-by-*n* matrix. The **singular value
decomposition (SVD)** of *A* is the factorization , where
*U* and *V* are orthogonal and
, ,
with . If *A* is complex,
its SVD is , where *U* and *V* are unitary
and is as before with real
diagonal elements.
The are called the **singular values** ,
the first *r* columns of *V*
the **right singular vectors** , and
the first *r* columns of *U*
the **left singular vectors** .

The routines described in this section, and listed in table 3.10, are used to compute this decomposition. The computation proceeds in the following stages:

- The matrix
*A*is reduced to bidiagonal form: if*A*is real ( if*A*is complex), where and are orthogonal (unitary if*A*is complex), and*B*is real and upper-bidiagonal when and lower bidiagonal when*m*<*n*, so that*B*is nonzero only on the main diagonal and either on the first superdiagonal (if ) or the first subdiagonal (if*m*<*n*). - The SVD of the bidiagonal matrix
*B*is computed: , where and are orthogonal and is diagonal as described above. The singular vectors of*A*are then and .

The reduction to bidiagonal form is performed by the subroutine PxGEBRD
and the SVD of *B* is computed using the LAPACK routine xBDSQR.

The routine PxGEBRD represents
and in factored form as products of elementary reflectors,
as described in section 3.4.
If *A* is real,
the matrices and may be multiplied by other matrices without
forming and using routine PxORMBR .
If *A* is complex, one instead uses PxUNMBR .

If , it may be more efficient to first perform a *QR* factorization
of *A*, using the routine
PxGEQRF ,
and then to compute the SVD of the *n*-by-*n* matrix *R*, since
if *A* = *QR* and , then the SVD of *A* is given by
.
Similarly, if , it may be more efficient to first perform an
*LQ* factorization of *A*, using xGELQF. These preliminary *QR* and *LQ*
factorizations are performed by the driver PxGESVD.

The SVD may be used to find a minimum norm solution to a (possibly)
rank-deficient linear least squares
problem (3.1). The effective rank, *k*, of *A* can be determined as the
number of singular values which exceed a suitable threshold.
Let be the leading *k*-by-*k* submatrix of , and
be the matrix consisting of the first *k* columns of *V*.
Then the solution is given by

where consists of the first *k* elements of . can be computed by using PxORMBR.

**Table 3.10:** Computational routines for the singular value decomposition

Tue May 13 09:21:01 EDT 1997