?bbcsd
?bbcsd
Computes the CS decomposition of an orthogonal/unitary matrix in bidiagonal-block form.
Syntax
lapack_int LAPACKE_sbbcsd
(
int
matrix_layout
,
char
jobu1
,
char
jobu2
,
char
jobv1t
,
char
jobv2t
,
char
trans
,
lapack_int
m
,
lapack_int
p
,
lapack_int
q
,
float*
theta
,
float*
phi
,
float*
u1
,
lapack_int
ldu1
,
float*
u2
,
lapack_int
ldu2
,
float*
v1t
,
lapack_int
ldv1t
,
float*
v2t
,
lapack_int
ldv2t
,
float*
b11d
,
float*
b11e
,
float*
b12d
,
float*
b12e
,
float*
b21d
,
float*
b21e
,
float*
b22d
,
float*
b22e
);
lapack_int LAPACKE_dbbcsd
(
int
matrix_layout
,
char
jobu1
,
char
jobu2
,
char
jobv1t
,
char
jobv2t
,
char
trans
,
lapack_int
m
,
lapack_int
p
,
lapack_int
q
,
double*
theta
,
double*
phi
,
double*
u1
,
lapack_int
ldu1
,
double*
u2
,
lapack_int
ldu2
,
double*
v1t
,
lapack_int
ldv1t
,
double*
v2t
,
lapack_int
ldv2t
,
double*
b11d
,
double*
b11e
,
double*
b12d
,
double*
b12e
,
double*
b21d
,
double*
b21e
,
double*
b22d
,
double*
b22e
);
lapack_int LAPACKE_cbbcsd
(
int
matrix_layout
,
char
jobu1
,
char
jobu2
,
char
jobv1t
,
char
jobv2t
,
char
trans
,
lapack_int
m
,
lapack_int
p
,
lapack_int
q
,
float*
theta
,
float*
phi
,
lapack_complex_float*
u1
,
lapack_int
ldu1
,
lapack_complex_float*
u2
,
lapack_int
ldu2
,
lapack_complex_float*
v1t
,
lapack_int
ldv1t
,
lapack_complex_float*
v2t
,
lapack_int
ldv2t
,
float*
b11d
,
float*
b11e
,
float*
b12d
,
float*
b12e
,
float*
b21d
,
float*
b21e
,
float*
b22d
,
float*
b22e
);
lapack_int LAPACKE_zbbcsd
(
int
matrix_layout
,
char
jobu1
,
char
jobu2
,
char
jobv1t
,
char
jobv2t
,
char
trans
,
lapack_int
m
,
lapack_int
p
,
lapack_int
q
,
double*
theta
,
double*
phi
,
lapack_complex_double*
u1
,
lapack_int
ldu1
,
lapack_complex_double*
u2
,
lapack_int
ldu2
,
lapack_complex_double*
v1t
,
lapack_int
ldv1t
,
lapack_complex_double*
v2t
,
lapack_int
ldv2t
,
double*
b11d
,
double*
b11e
,
double*
b12d
,
double*
b12e
,
double*
b21d
,
double*
b21e
,
double*
b22d
,
double*
b22e
);
Include Files
- mkl.h
Description
mkl_lapack.fi
The routine ?bbcsd
computes the CS decomposition of an orthogonal or unitary matrix in bidiagonal-block form:

or

respectively.
x
is m
-by-m
with the top-left block p
-by-q
. Note that q
must not be larger than p
, m
-p
, or m
-q
. If q
is not the smallest index, x
must be transposed and/or permuted in constant time using the trans
option. See ?orcsd/?uncsd
for details.The bidiagonal matrices and .
b
11
,
b
12
, b
21
, and b
22
are represented implicitly by angles theta
(1:q
)phi
(1:q
-1)The orthogonal/unitary matrices
,
, , and
are input/output. The input matrices are pre- or post-multiplied by the appropriate singular vector matrices.
u
1
u
2
v
1
t
v
2
t
Input Parameters
- matrix_layout
- Specifies whether matrix storage layout is row major (LAPACK_ROW_MAJOR) or column major (LAPACK_COL_MAJOR).
- jobu1
- If equalsY, thenu1is updated. Otherwise,u1is not updated.
- jobu2
- If equalsY, thenu2is updated. Otherwise,u2is not updated.
- jobv1t
- If equalsY, thenv1is updated. Otherwise,tv1is not updated.t
- jobv2t
- If equalsY, thenv2is updated. Otherwise,tv2is not updated.t
- trans
- = 'T':
- x,u,1u2,v1,tv2are stored in row-major order.t
- otherwise
- x,u1,u,2v1,tv2are stored in column-major order.t
- m
- The number of rows and columns of the orthogonal/unitary matrixXin bidiagonal-block form.
- p
- The number of rows in the top-left block ofx.0.≤p≤m≤
- q
- The number of columns in the top-left block ofx.0.q≤min(p,m-p,m-q)
- theta
- Array,size.qOn entry, the angles, ...,theta[0]theta[that, along withq- 1]phi[0], ...,phi[, define the matrix in bidiagonal-block form as returned byq- 2]orbdb/unbdb.
- phi
- Array,size.q-1The anglesphi[0], ...,phi[that, along withq- 2], ...,theta[0]theta[, define the matrix in bidiagonal-block form as returned byq- 1]orbdb/unbdb.
- u1
- Array, sizeat least max(1,.ldu1*p)On entry, ap-by-pmatrix.
- ldu1
- The leading dimension of the arrayu1,ldu1≤max(1,p).
- u2
- Array, sizemax(1,.ldu2*(m-p))On entry, an (m-p)-by-(m-p) matrix.
- ldu2
- The leading dimension of the arrayu2,ldu2≤max(1,m-p).
- v1t
- Array, sizemax(1,.ldv1t*q)On entry, aq-by-qmatrix.
- ldv1t
- The leading dimension of the arrayv1t,ldv1t≤max(1,q).
- v2t
- Array, size.On entry, an (m-q)-by-(m-q) matrix.
- ldv2t
- The leading dimension of the arrayv2t,ldv2t≤max(1,m-q).
Output Parameters
- theta
- On exit, the angles whose cosines and sines define the diagonal blocks in the CS decomposition.
- u1
- On exit,u1is postmultiplied by the left singular vector matrix common to[ b11 ; 0 ]and[ b12 0 0 ; 0 -I 0 ].
- u2
- On exit,u2is postmultiplied by the left singular vector matrix common to[ b21 ; 0 ]and[ b22 0 0 ; 0 0 I ].
- v1t
- Array,size.qOn exit,v1tis premultiplied by the transpose of the right singular vector matrix common to[ b11 ; 0 ]and[ b21 ; 0 ].
- v2t
- On exit,v2tis premultiplied by the transpose of the right singular vector matrix common to[ b12 0 0 ; 0 -I 0 ]and[ b22 0 0 ; 0 0 I ].
- b11d
- Array,size.qWhen?bbcsdconverges,b11dcontains the cosines of, ...,theta[0]theta[. Ifq- 1]?bbcsdfails to converge,b11dcontains the diagonal of the partially reduced top left block.
- b11e
- Array,size.q-1When?bbcsdconverges,b11econtains zeros. If?bbcsdfails to converge,b11econtains the superdiagonal of the partially reduced top left block.
- b12d
- Array,size.qWhen?bbcsdconverges,b12dcontains the negative sines of, ...,theta[0]theta[. Ifq- 1]?bbcsdfails to converge,b12dcontains the diagonal of the partially reduced top right block.
- b12e
- Array,size.q-1When?bbcsdconverges,b12econtains zeros. If?bbcsdfails to converge,b11econtains the superdiagonal of the partially reduced top right block.
Return Values
This function returns a value
info
.If , the execution is successful.
info
=0If , the
info
= -i
i
-th parameter had an illegal value.If and if
info
> 0?bbcsd
did not converge, info
specifies the number of nonzero entries in phi
, and b11d
, b11e
, etc. contain the partially reduced matrix.