Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

Chapter Contents
Chapter Introduction
NAG Toolbox

# NAG Toolbox: nag_lapack_dggsvp (f08ve)

## Purpose

nag_lapack_dggsvp (f08ve) uses orthogonal transformations to simultaneously reduce the $m$ by $n$ matrix $A$ and the $p$ by $n$ matrix $B$ to upper triangular form. This factorization is usually used as a preprocessing step for computing the generalized singular value decomposition (GSVD).

## Syntax

[a, b, k, l, u, v, q, info] = f08ve(jobu, jobv, jobq, a, b, tola, tolb, 'm', m, 'p', p, 'n', n)
[a, b, k, l, u, v, q, info] = nag_lapack_dggsvp(jobu, jobv, jobq, a, b, tola, tolb, 'm', m, 'p', p, 'n', n)

## Description

nag_lapack_dggsvp (f08ve) computes orthogonal matrices $U$, $V$ and $Q$ such that
where the $k$ by $k$ matrix ${A}_{12}$ and $l$ by $l$ matrix ${B}_{13}$ are nonsingular upper triangular; ${A}_{23}$ is $l$ by $l$ upper triangular if $m-k-l\ge 0$ and is $\left(m-k\right)$ by $l$ upper trapezoidal otherwise. $\left(k+l\right)$ is the effective numerical rank of the $\left(m+p\right)$ by $n$ matrix ${\left(\begin{array}{cc}{A}^{\mathrm{T}}& {B}^{\mathrm{T}}\end{array}\right)}^{\mathrm{T}}$.
This decomposition is usually used as the preprocessing step for computing the Generalized Singular Value Decomposition (GSVD), see function nag_lapack_dggsvd (f08va).

## References

Anderson E, Bai Z, Bischof C, Blackford S, Demmel J, Dongarra J J, Du Croz J J, Greenbaum A, Hammarling S, McKenney A and Sorensen D (1999) LAPACK Users' Guide (3rd Edition) SIAM, Philadelphia http://www.netlib.org/lapack/lug
Golub G H and Van Loan C F (1996) Matrix Computations (3rd Edition) Johns Hopkins University Press, Baltimore

## Parameters

### Compulsory Input Parameters

1:     $\mathrm{jobu}$ – string (length ≥ 1)
If ${\mathbf{jobu}}=\text{'U'}$, the orthogonal matrix $U$ is computed.
If ${\mathbf{jobu}}=\text{'N'}$, $U$ is not computed.
Constraint: ${\mathbf{jobu}}=\text{'U'}$ or $\text{'N'}$.
2:     $\mathrm{jobv}$ – string (length ≥ 1)
If ${\mathbf{jobv}}=\text{'V'}$, the orthogonal matrix $V$ is computed.
If ${\mathbf{jobv}}=\text{'N'}$, $V$ is not computed.
Constraint: ${\mathbf{jobv}}=\text{'V'}$ or $\text{'N'}$.
3:     $\mathrm{jobq}$ – string (length ≥ 1)
If ${\mathbf{jobq}}=\text{'Q'}$, the orthogonal matrix $Q$ is computed.
If ${\mathbf{jobq}}=\text{'N'}$, $Q$ is not computed.
Constraint: ${\mathbf{jobq}}=\text{'Q'}$ or $\text{'N'}$.
4:     $\mathrm{a}\left(\mathit{lda},:\right)$ – double array
The first dimension of the array a must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{m}}\right)$.
The second dimension of the array a must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$.
The $m$ by $n$ matrix $A$.
5:     $\mathrm{b}\left(\mathit{ldb},:\right)$ – double array
The first dimension of the array b must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{p}}\right)$.
The second dimension of the array b must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$.
The $p$ by $n$ matrix $B$.
6:     $\mathrm{tola}$ – double scalar
7:     $\mathrm{tolb}$ – double scalar
tola and tolb are the thresholds to determine the effective numerical rank of matrix $B$ and a subblock of $A$. Generally, they are set to
 $tola=maxm,nAε, tolb=maxp,nBε,$
where $\epsilon$ is the machine precision.
The size of tola and tolb may affect the size of backward errors of the decomposition.

### Optional Input Parameters

1:     $\mathrm{m}$int64int32nag_int scalar
Default: the first dimension of the array a.
$m$, the number of rows of the matrix $A$.
Constraint: ${\mathbf{m}}\ge 0$.
2:     $\mathrm{p}$int64int32nag_int scalar
Default: the first dimension of the array b.
$p$, the number of rows of the matrix $B$.
Constraint: ${\mathbf{p}}\ge 0$.
3:     $\mathrm{n}$int64int32nag_int scalar
Default: the second dimension of the array a.
$n$, the number of columns of the matrices $A$ and $B$.
Constraint: ${\mathbf{n}}\ge 0$.

### Output Parameters

1:     $\mathrm{a}\left(\mathit{lda},:\right)$ – double array
The first dimension of the array a will be $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{m}}\right)$.
The second dimension of the array a will be $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$.
Contains the triangular (or trapezoidal) matrix described in Description.
2:     $\mathrm{b}\left(\mathit{ldb},:\right)$ – double array
The first dimension of the array b will be $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{p}}\right)$.
The second dimension of the array b will be $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$.
Contains the triangular matrix described in Description.
3:     $\mathrm{k}$int64int32nag_int scalar
4:     $\mathrm{l}$int64int32nag_int scalar
k and l specify the dimension of the subblocks $k$ and $l$ as described in Description; $\left(k+l\right)$ is the effective numerical rank of ${\left(\begin{array}{cc}{{\mathbf{a}}}^{\mathrm{T}}& {{\mathbf{b}}}^{\mathrm{T}}\end{array}\right)}^{\mathrm{T}}$.
5:     $\mathrm{u}\left(\mathit{ldu},:\right)$ – double array
The first dimension, $\mathit{ldu}$, of the array u will be
• if ${\mathbf{jobu}}=\text{'U'}$, $\mathit{ldu}=\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{m}}\right)$;
• otherwise $\mathit{ldu}=1$.
The second dimension of the array u will be $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{m}}\right)$ if ${\mathbf{jobu}}=\text{'U'}$ and $1$ otherwise.
If ${\mathbf{jobu}}=\text{'U'}$, u contains the orthogonal matrix $U$.
If ${\mathbf{jobu}}=\text{'N'}$, u is not referenced.
6:     $\mathrm{v}\left(\mathit{ldv},:\right)$ – double array
The first dimension, $\mathit{ldv}$, of the array v will be
• if ${\mathbf{jobv}}=\text{'V'}$, $\mathit{ldv}=\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{p}}\right)$;
• otherwise $\mathit{ldv}=1$.
The second dimension of the array v will be $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{p}}\right)$ if ${\mathbf{jobv}}=\text{'V'}$ and $1$ otherwise.
If ${\mathbf{jobv}}=\text{'V'}$, v contains the orthogonal matrix $V$.
If ${\mathbf{jobv}}=\text{'N'}$, v is not referenced.
7:     $\mathrm{q}\left(\mathit{ldq},:\right)$ – double array
The first dimension, $\mathit{ldq}$, of the array q will be
• if ${\mathbf{jobq}}=\text{'Q'}$, $\mathit{ldq}=\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$;
• otherwise $\mathit{ldq}=1$.
The second dimension of the array q will be $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$ if ${\mathbf{jobq}}=\text{'Q'}$ and $1$ otherwise.
If ${\mathbf{jobq}}=\text{'Q'}$, q contains the orthogonal matrix $Q$.
If ${\mathbf{jobq}}=\text{'N'}$, q is not referenced.
8:     $\mathrm{info}$int64int32nag_int scalar
${\mathbf{info}}=0$ unless the function detects an error (see Error Indicators and Warnings).

## Error Indicators and Warnings

${\mathbf{info}}=-i$
If ${\mathbf{info}}=-i$, parameter $i$ had an illegal value on entry. The parameters are numbered as follows:
1: jobu, 2: jobv, 3: jobq, 4: m, 5: p, 6: n, 7: a, 8: lda, 9: b, 10: ldb, 11: tola, 12: tolb, 13: k, 14: l, 15: u, 16: ldu, 17: v, 18: ldv, 19: q, 20: ldq, 21: iwork, 22: tau, 23: work, 24: info.
It is possible that info refers to a parameter that is omitted from the MATLAB interface. This usually indicates that an error in one of the other input parameters has caused an incorrect value to be inferred.

## Accuracy

The computed factorization is nearly the exact factorization for nearby matrices $\left(A+E\right)$ and $\left(B+F\right)$, where
 $E2 = OεA2 and F2= OεB2,$
and $\epsilon$ is the machine precision.

The complex analogue of this function is nag_lapack_zggsvp (f08vs).

## Example

This example finds the generalized factorization
 $A = UΣ1 0 S QT , B= VΣ2 0 T QT ,$
of the matrix pair $\left(\begin{array}{cc}A& B\end{array}\right)$, where
 $A = 123 321 456 788 and B= -2-33 465 .$
```function f08ve_example

fprintf('f08ve example results\n\n');

% Generalized singular values of (A,B) where:
m = 4;
n = 3;
a = [ 1  2  3;
3  2  1;
4  5  6;
7  8  8];
p = 2;
b = [-2 -3  3;
4  6  5];

% Reduce A and B to upper triangular form S = U^T A Q, T = V^T B Q
tola = max(m,n)*norm(a,1)*x02aj;
tolb = max(p,n)*norm(a,1)*x02aj;
[S, T, k, l, U, V, Q, info] = ...
f08ve( ...
'U', 'V', 'Q', a, b, tola, tolb);

% Compute singular values
[S, T, alpha, beta, U, V, Q, ncycle, info] = ...
f08ye( ...
'U', 'V', 'Q', k, l, S, T, tola, tolb, U, V, Q);

fprintf('Number of infinite generalized singular values = %3d\n',k);
fprintf('Number of   finite generalized singular values = %3d\n',l);
fprintf('Effective rank of the matrix pair  (A^T B^T)^T = %3d\n',k+l);
fprintf('Number of cycles of the Kogbetliantz method    = %3d\n\n',ncycle);
disp('Finite generalized singular values');
disp(alpha(k+1:k+l)./beta(k+1:k+l));

```
```f08ve example results

Number of infinite generalized singular values =   1
Number of   finite generalized singular values =   2
Effective rank of the matrix pair  (A^T B^T)^T =   3
Number of cycles of the Kogbetliantz method    =   2

Finite generalized singular values
1.3151
0.0802

```