?gghrd
?gghrd
Reduces a pair of matrices to generalized upper Hessenberg form using orthogonal/unitary transformations.
Syntax
lapack_int
LAPACKE_sgghrd
(
int
matrix_layout
,
char
compq
,
char
compz
,
lapack_int
n
,
lapack_int
ilo
,
lapack_int
ihi
,
float
*
a
,
lapack_int
lda
,
float
*
b
,
lapack_int
ldb
,
float
*
q
,
lapack_int
ldq
,
float
*
z
,
lapack_int
ldz
);
lapack_int
LAPACKE_dgghrd
(
int
matrix_layout
,
char
compq
,
char
compz
,
lapack_int
n
,
lapack_int
ilo
,
lapack_int
ihi
,
double
*
a
,
lapack_int
lda
,
double
*
b
,
lapack_int
ldb
,
double
*
q
,
lapack_int
ldq
,
double
*
z
,
lapack_int
ldz
);
lapack_int
LAPACKE_cgghrd
(
int
matrix_layout
,
char
compq
,
char
compz
,
lapack_int
n
,
lapack_int
ilo
,
lapack_int
ihi
,
lapack_complex_float
*
a
,
lapack_int
lda
,
lapack_complex_float
*
b
,
lapack_int
ldb
,
lapack_complex_float
*
q
,
lapack_int
ldq
,
lapack_complex_float
*
z
,
lapack_int
ldz
);
lapack_int
LAPACKE_zgghrd
(
int
matrix_layout
,
char
compq
,
char
compz
,
lapack_int
n
,
lapack_int
ilo
,
lapack_int
ihi
,
lapack_complex_double
*
a
,
lapack_int
lda
,
lapack_complex_double
*
b
,
lapack_int
ldb
,
lapack_complex_double
*
q
,
lapack_int
ldq
,
lapack_complex_double
*
z
,
lapack_int
ldz
);
Include Files
- mkl.h
Description
The routine reduces a pair of real/complex matrices (, and
A
,B
) to generalized upper Hessenberg form using orthogonal/unitary transformations, where A
is a general matrix and B
is upper triangular. The form of the generalized eigenvalue problem is A
*x
= λ
*B
*x
B
is typically made upper triangular by computing its QR
factorization and moving the orthogonal matrix Q
to the left side of the equation. This routine simultaneously reduces
A
to a Hessenberg matrix H
:Q
H
A
*Z
= H
and transforms
B
to another upper triangular matrix T
:Q
H
B
*Z
= T
in order to reduce the problem to its standard form , where *.
H
*y
= λ
*T
*y
y
= Z
H
x
The orthogonal/unitary matrices
Q
and Z
are determined as products of Givens rotations. They may either be formed explicitly, or they may be postmultiplied into input matrices Q
1
and Z
1
, so that Q
1
*A
*Z
1
H
Q
1
*Q
)*H
*(Z
1
*Z
)H
Q
1
*B*
Z
1
H
Q
1
*Q
)*T*
(Z
1
*Z
)H
If , then the routine
Q
1
is the orthogonal/unitary matrix from the QR
factorization of B
in the original equation A
*x
= λ
*B
*x
?gghrd
reduces the original problem to generalized Hessenberg form.Input Parameters
- matrix_layout
- Specifies whether matrix storage layout is row major (LAPACK_ROW_MAJOR) or column major (LAPACK_COL_MAJOR).
- compq
- Must be'N','I', or'V'.If, matrixcompq='N'Qis not computed.If,compq='I'Qis initialized to the unit matrix, and the orthogonal/unitary matrixQis returned;If,compq='V'Qmust contain an orthogonal/unitary matrixQ1on entry, and the productQ1*Qis returned.
- compz
- Must be'N','I', or'V'.If, matrixcompz='N'Zis not computed.If,compz='I'Zis initialized to the unit matrix, and the orthogonal/unitary matrixZis returned;If,compz='V'Zmust contain an orthogonal/unitary matrixZ1on entry, and the productZ1*Zis returned.
- n
- The order of the matricesAandB(n≥0).
- ilo,ihi
- iloandihimark the rows and columns ofAwhich are to be reduced. It is assumed thatAis already upper triangular in rows and columns 1:ilo-1 andihi+1:n. Values ofiloandihiare normally set by a previous call to ggbal; otherwise they should be set to 1 andnrespectively.Constraint:If, then 1n> 0≤ilo≤ihi≤n;if, thenn= 0andilo= 1.ihi= 0
- a,b,q,z
- Arrays:a(size max(1,contains thelda*n))n-by-ngeneral matrixA.b(size max(1,contains theldb*n))n-by-nupper triangular matrixB.q(size max(1,ldq*n))If, thencompq='N'qis not referenced.If, thencompq='V'qmust contain the orthogonal/unitary matrixQ1, typically from theQRfactorization ofB.z(size max(1,ldz*n))If, thencompz='N'zis not referenced.If, thencompz='V'zmust contain the orthogonal/unitary matrixZ1.
- lda
- The leading dimension ofa; at least max(1,n).
- ldb
- The leading dimension ofb; at least max(1,n).
- ldq
- The leading dimension ofq;If, thencompq='N'ldq≥1.Iforcompq='I''V', thenldq≥max(1,n).
- ldz
- The leading dimension ofz;If, thencompz='N'ldz≥1.Iforcompz='I''V', thenldz≥max(1,n).
Output Parameters
- a
- On exit, the upper triangle and the first subdiagonal ofAare overwritten with the upper Hessenberg matrixH, and the rest is set to zero.
- b
- On exit, overwritten by the upper triangular matrixT=Q*HB*Z. The elements below the diagonal are set to zero.
- q
- If, thencompq='I'qcontains the orthogonal/unitary matrixQ, ;If, thencompq='V'qis overwritten by the productQ1*Q.
- z
- If, thencompz='I'zcontains the orthogonal/unitary matrixZ;If, thencompz='V'zis overwritten by the productZ1*Z.
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.