?gghd3
?gghd3
Reduces a pair of matrices to generalized upper Hessenberg form.
Syntax
lapack_int
LAPACKE_sgghd3
(
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_dgghd3
(
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_cgghd3
(
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_zgghd3
(
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
?gghd3
reduces a pair of real or complex matrices (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 isA
*x
= λ
*B
*x
,and
B
is typically made upper triangular by computing its QR factorization and moving the orthogonal/unitary matrix Q
to the left side of the equation.This subroutine simultaneously reduces
A
to a Hessenberg matrix H
:Q
T
*A
*Z
= H
for real flavorsor
Q
T
*A
*Z
= H
for complex flavorsand transforms
B
to another upper triangular matrix T
:Q
T
*B
*Z
= T
for real flavorsor
Q
T
*B
*Z
= T
for complex flavorsin order to reduce the problem to its standard form
H
*y
= λ
*T
*y
where
y
= Z
T
*x
for real flavorsor
y
= Z
T
*x
for complex flavors.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 thatfor real flavors:
Q
1
* A
* Z
1
T
= (Q
1
*Q
) * H
* (Z
1
*Z
)T
Q
1
* B
* Z
1
T
= (Q
1
*Q
) * T
* (Z
1
*Z
)T
for complex flavors:
Q
1
* A
* Z
1
H
= (Q
1
*Q
) * H
* (Z
1
*Z
)T
Q
1
* B
* Z
1
T
= (Q
1
*Q
) * T
* (Z
1
*Z
)T
If
Q
1
is the orthogonal/unitary matrix from the QR factorization of B
in the original equation A
*x
= λ
*B
*x
, then ?gghd3
reduces the original problem to generalized Hessenberg form.This is a blocked variant of
?gghrd
, using matrix-matrix multiplications for parts of the computation to enhance performance.Input Parameters
- matrix_layout
- Specifies whether matrix storage layout is row major (LAPACK_ROW_MAJOR) or column major (LAPACK_COL_MAJOR).
- compq
- = 'N': do not computeq;= 'I':qis initialized to the unit matrix, and the orthogonal/unitary matrixQis returned;= 'V':qmust contain an orthogonal/unitary matrixQ1on entry, and the productQ1*qis returned.
- compz
- = 'N': do not computez;= 'I':zis initialized to the unit matrix, and the orthogonal/unitary matrixZis returned;= '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.iloandihiare normally set by a previous call to?ggbal; otherwise they should be set to 1 andn, respectively.1≤ilo≤ihi≤n, ifn> 0;ilo=1 andihi=0, ifn=0.
- a
- Array, size(.lda*n)On entry, then-by-ngeneral matrix to be reduced.
- lda
- The leading dimension of the arraya.lda≥max(1,n).
- b
- Array,(.ldb*n)On entry, then-by-nupper triangular matrixB.
- ldb
- The leading dimension of the arrayb.ldb≥max(1,n).
- q
- Array, size(.ldq*n)On entry, ifcompq= 'V', the orthogonal/unitary matrixQ1, typically from the QR factorization ofb.
- ldq
- The leading dimension of the arrayq.ldq≥nifcompq='V' or 'I';ldq≥1 otherwise.
- z
- Array, size(.ldz*n)On entry, ifcompz= 'V', the orthogonal/unitary matrixZ1.Not referenced ifcompz='N'.
- ldz
- The leading dimension of the arrayz.ldz≥nifcompz='V' or 'I';ldz≥1 otherwise.
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, the upper triangular matrixT=QTBZfor real flavors orT=QHBZfor complex flavors. The elements below the diagonal are set to zero.
- q
- On exit, ifcompq='I', the orthogonal/unitary matrixQ, and ifcompq= 'V', the productQ1*Q.Not referenced ifcompq='N'.
- z
- On exit, ifcompz='I', the orthogonal/unitary matrixZ, and ifcompz= 'V', the productZ1*Z.Not referenced ifcompz='N'.
Return Values
This function returns a value
info
.= 0: successful exit.
< 0: if
info
= -i
, the i
-th argument had an illegal value.Application Notes
This routine reduces
A
to Hessenberg form and maintains B
in using a blocked variant of Moler and Stewart's original algorithm, as described by Kagstrom, Kressner, Quintana-Orti, and Quintana-Orti (BIT 2008).