Error Bounds for the Generalized Singular Value Decomposition

The generalized (or quotient) singular value decomposition
of an ** m**-by-

where

is*U*-by-*m*,*m*is*V*-by-*p*,*p*is*Q*-by-*n*, and all three matrices are orthogonal. If*n*and*A*are complex, these matrices are unitary instead of orthogonal, and*B*should be replaced by*Q*^{T}in the pair of factorizations.*Q*^{H}is*R*-by-*r*, upper triangular and nonsingular.*r***[0,**is*R*]-by-*r*. The integer*n*is the rank of , and satisfies .*r*-
is
-by-*m*, is*r*-by-*p*, both are real, nonnegative and diagonal, and . Write and , where and lie in the interval from 0 to 1. The ratios are called the*r***generalized singular values**of the pair,*A*. If , then the generalized singular value is*B***infinite**. For details on the structure of , and, see section 2.3.5.3.*R*

The generalized singular value decomposition is
computed by driver routine xGGSVD (see section 2.3.5.3).
We will give error bounds for the generalized
singular values in the
common case where
has full
rank ** r=n**.
Let
and
be the values of
and ,
respectively,
computed by xGGSVD.
The approximate error
bound

Note that if is close to zero, then a true generalized singular value can differ greatly in magnitude from the computed generalized singular value , even if

Here is another way to interpret `SERRBD`:
if we think of
and
as representing the *subspace*
consisting of the straight line through the origin with slope
,
and similarly
and
representing the subspace ,
then
bounds the acute angle between
and .
Note that any two
lines through the origin with nearly vertical slopes
(very large
)
are close together in angle.
(This is related to the *chordal distance* in
section 4.10.1.)

`SERRBD` can be computed by the following code fragment,
which for simplicity assumes .
(The assumption ** r=n** implies only that .
Error bounds can also be computed when
,
with slightly more complicated code.)

EPSMCH = SLAMCH( 'E' ) * Compute generalized singular values of A and B CALL SGGSVD( 'N', 'N', 'N', M, N, P, K, L, A, LDA, B, $ LDB, ALPHA, BETA, U, LDU, V, LDV, Q, LDQ, $ WORK, IWORK, INFO ) * Compute rank of [A',B']' RANK = K+L IF( INFO.GT.0 ) THEN PRINT *,'SGGSVD did not converge' ELSE IF( RANK.LT.N ) THEN PRINT *,'[A**T,B**T]**T not full rank' ELSE IF ( M .GE. N .AND. N .GT. 0 ) THEN * Compute reciprocal condition number RCOND of R CALL STRCON( 'I', 'U', 'N', N, A, LDA, RCOND, WORK, IWORK, $ INFO ) RCOND = MAX( RCOND, EPSMCH ) SERRBD = EPSMCH / RCOND END IF

For example^{4.11}, if
,

then, to 4 decimal places,

, and the true errors are