The GSVD algorithm used in LAPACK ([83,10,8]) is backward stable:

Let the computed GSVD ofandAbe and . This is nearly the exact GSVD ofBandA+Ein the following sense.B+FandEare small:F

there exist small , , and such that , , and are exactly orthogonal (or unitary):

and

is the exact GSVD ofandA+E. HereB+Fis a modestly growing function ofp(n), and we takenin the above code fragment.p(n)=1

Let and be the square roots of the diagonal entries of the exact and , and let and the square roots of the diagonal entries of the computed and . Let

Then providedand have full rankG, one can show [96,82] thatn

In the code fragment we approximate the numerator of the last expression by and approximate the denominator by in order to computeSERRBD; STRCON returns an approximationRCONDto .

We assume that the rank ** r** of

then

The reason the code fragment assumes
is that in this case
is
stored overwritten on ** A**, and can be passed to STRCON in order to compute