Developer Reference

Contents

?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 (
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
, and
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
H
*
y
=
λ
*
T
*
y
, where
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
Q
1
is the orthogonal/unitary matrix from the
QR
factorization of
B
in the original equation
A
*
x
=
λ
*
B
*
x
, then the routine
?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
compq
=
'N'
, matrix
Q
is not computed.
If
compq
=
'I'
,
Q
is initialized to the unit matrix, and the orthogonal/unitary matrix
Q
is returned;
If
compq
=
'V'
,
Q
must contain an orthogonal/unitary matrix
Q
1
on entry, and the product
Q
1
*
Q
is returned.
compz
Must be
'N'
,
'I'
, or
'V'
.
If
compz
=
'N'
, matrix
Z
is not computed.
If
compz
=
'I'
,
Z
is initialized to the unit matrix, and the orthogonal/unitary matrix
Z
is returned;
If
compz
=
'V'
,
Z
must contain an orthogonal/unitary matrix
Z
1
on entry, and the product
Z
1
*
Z
is returned.
n
The order of the matrices
A
and
B
(
n
0).
ilo
,
ihi
ilo
and
ihi
mark the rows and columns of
A
which are to be reduced. It is assumed that
A
is already upper triangular in rows and columns 1:
ilo
-1 and
ihi
+1:
n
. Values of
ilo
and
ihi
are normally set by a previous call to ggbal; otherwise they should be set to 1 and
n
respectively.
Constraint:
If
n
> 0
, then 1
ilo
ihi
n
;
if
n
= 0
, then
ilo
= 1
and
ihi
= 0
.
a
,
b
,
q
,
z
Arrays:
a
(size max(1,
lda
*
n
))
contains the
n
-by-
n
general matrix
A
.
b
(size max(1,
ldb
*
n
))
contains the
n
-by-
n
upper triangular matrix
B
.
q
(size max(1,
ldq
*
n
))
If
compq
=
'N'
, then
q
is not referenced.
If
compq
=
'V'
, then
q
must contain the orthogonal/unitary matrix
Q
1
, typically from the
QR
factorization of
B
.
z
(size max(1,
ldz
*
n
))
If
compz
=
'N'
, then
z
is not referenced.
If
compz
=
'V'
, then
z
must contain the orthogonal/unitary matrix
Z
1
.
lda
The leading dimension of
a
; at least max(1,
n
).
ldb
The leading dimension of
b
; at least max(1,
n
).
ldq
The leading dimension of
q
;
If
compq
=
'N'
, then
ldq
1.
If
compq
=
'I'
or
'V'
, then
ldq
max(1,
n
).
ldz
The leading dimension of
z
;
If
compz
=
'N'
, then
ldz
1.
If
compz
=
'I'
or
'V'
, then
ldz
max(1,
n
).
Output Parameters
a
On exit, the upper triangle and the first subdiagonal of
A
are overwritten with the upper Hessenberg matrix
H
, and the rest is set to zero.
b
On exit, overwritten by the upper triangular matrix
T
=
Q
H
*
B*
Z
. The elements below the diagonal are set to zero.
q
If
compq
=
'I'
, then
q
contains the orthogonal/unitary matrix
Q
, ;
If
compq
=
'V'
, then
q
is overwritten by the product
Q
1
*
Q
.
z
If
compz
=
'I'
, then
z
contains the orthogonal/unitary matrix
Z
;
If
compz
=
'V'
, then
z
is overwritten by the product
Z
1
*
Z
.
Return Values
This function returns a value
info
.
If
info
=0
, the execution is successful.
If
info
=
-i
, the
i
-th parameter had an illegal value.

Product and Performance Information

1

Intel's compilers may or may not optimize to the same degree for non-Intel microprocessors for optimizations that are not unique to Intel microprocessors. These optimizations include SSE2, SSE3, and SSSE3 instruction sets and other optimizations. Intel does not guarantee the availability, functionality, or effectiveness of any optimization on microprocessors not manufactured by Intel. Microprocessor-dependent optimizations in this product are intended for use with Intel microprocessors. Certain optimizations not specific to Intel microarchitecture are reserved for Intel microprocessors. Please refer to the applicable product User and Reference Guides for more information regarding the specific instruction sets covered by this notice.

Notice revision #20110804