?hegv
?hegv
Computes all eigenvalues and, optionally, eigenvectors of a complex generalized Hermitian positive-definite eigenproblem.
Syntax
lapack_int LAPACKE_chegv
(
int
matrix_layout
,
lapack_int
itype
,
char
jobz
,
char
uplo
,
lapack_int
n
,
lapack_complex_float*
a
,
lapack_int
lda
,
lapack_complex_float*
b
,
lapack_int
ldb
,
float*
w
);
lapack_int LAPACKE_zhegv
(
int
matrix_layout
,
lapack_int
itype
,
char
jobz
,
char
uplo
,
lapack_int
n
,
lapack_complex_double*
a
,
lapack_int
lda
,
lapack_complex_double*
b
,
lapack_int
ldb
,
double*
w
);
Include Files
- mkl.h
Description
The routine computes all the eigenvalues, and optionally, the eigenvectors of a complex generalized Hermitian positive-definite eigenproblem, of the form
A
*x
= λ
*B
*x
, A
*B
*x
= λ
*x
B
*A
*x
= λ
*x
Here
A
and B
are assumed to be Hermitian and B
is also positive definite.Input Parameters
- matrix_layout
- Specifies whether matrix storage layout is row major (LAPACK_ROW_MAJOR) or column major (LAPACK_COL_MAJOR).
- itype
- Must be 1 or 2 or 3. Specifies the problem type to be solved:if, the problem type isitype= 1;A*x=lambda*B*xif, the problem type isitype= 2;A*B*x=lambda*xif, the problem type isitype= 3.B*A*x=lambda*x
- jobz
- Must be'N'or'V'.If, then compute eigenvalues only.jobz='N'If, then compute eigenvalues and eigenvectors.jobz='V'
- uplo
- Must be'U'or'L'.If, arraysuplo='U'aandbstore the upper triangles ofAandB;If, arraysuplo='L'aandbstore the lower triangles ofAandB.
- n
- The order of the matricesAandB().n≥0
- a,b
- Arrays:a(size at least max(1,contains the upper or lower triangle of the Hermitian matrixlda*n))A, as specified byuplo.b(size at least max(1,contains the upper or lower triangle of the Hermitian positive definite matrixldb*n))B, as specified byuplo.
- lda
- The leading dimension ofa; at least max(1,n).
- ldb
- The leading dimension ofb; at least max(1,n).
Output Parameters
- a
- On exit, if, then ifjobz='V',info= 0acontains the matrixZof eigenvectors. The eigenvectors are normalized as follows:ifor 2,itype= 1;Z*HB*Z= Iif,itype= 3;Z*inv(HB)*Z= IIf, then on exit the upper triangle (ifjobz='N') or the lower triangle (ifuplo='U') ofuplo='L'A, including the diagonal, is destroyed.
- b
- On exit, if, the part ofinfo≤nbcontaining the matrix is overwritten by the triangular factorUorLfrom the Cholesky factorizationorB=U*HU.B=L*LH
- w
- Array, size at least max(1,n).If, contains the eigenvalues in ascending order.info= 0
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 ,
info
> 0cpotrf
/zpotrf
or cheev
/zheev
return an error code:If ,
info
= i
≤
n
cheev
/zheev
fails to converge, and i
off-diagonal elements of an intermediate tridiagonal do not converge to zero;If , for
info
= n
+ i
1
, then the leading minor of order ≤
i
≤
n
i
of B
is not positive-definite. The factorization of B
can not be completed and no eigenvalues or eigenvectors are computed.