version 1.1.1.1, 2010/01/26 15:22:45
|
version 1.14, 2016/08/27 15:34:53
|
Line 1
|
Line 1
|
|
*> \brief \b ZHSEIN |
|
* |
|
* =========== DOCUMENTATION =========== |
|
* |
|
* Online html documentation available at |
|
* http://www.netlib.org/lapack/explore-html/ |
|
* |
|
*> \htmlonly |
|
*> Download ZHSEIN + dependencies |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhsein.f"> |
|
*> [TGZ]</a> |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhsein.f"> |
|
*> [ZIP]</a> |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhsein.f"> |
|
*> [TXT]</a> |
|
*> \endhtmlonly |
|
* |
|
* Definition: |
|
* =========== |
|
* |
|
* SUBROUTINE ZHSEIN( SIDE, EIGSRC, INITV, SELECT, N, H, LDH, W, VL, |
|
* LDVL, VR, LDVR, MM, M, WORK, RWORK, IFAILL, |
|
* IFAILR, INFO ) |
|
* |
|
* .. Scalar Arguments .. |
|
* CHARACTER EIGSRC, INITV, SIDE |
|
* INTEGER INFO, LDH, LDVL, LDVR, M, MM, N |
|
* .. |
|
* .. Array Arguments .. |
|
* LOGICAL SELECT( * ) |
|
* INTEGER IFAILL( * ), IFAILR( * ) |
|
* DOUBLE PRECISION RWORK( * ) |
|
* COMPLEX*16 H( LDH, * ), VL( LDVL, * ), VR( LDVR, * ), |
|
* $ W( * ), WORK( * ) |
|
* .. |
|
* |
|
* |
|
*> \par Purpose: |
|
* ============= |
|
*> |
|
*> \verbatim |
|
*> |
|
*> ZHSEIN uses inverse iteration to find specified right and/or left |
|
*> eigenvectors of a complex upper Hessenberg matrix H. |
|
*> |
|
*> The right eigenvector x and the left eigenvector y of the matrix H |
|
*> corresponding to an eigenvalue w are defined by: |
|
*> |
|
*> H * x = w * x, y**h * H = w * y**h |
|
*> |
|
*> where y**h denotes the conjugate transpose of the vector y. |
|
*> \endverbatim |
|
* |
|
* Arguments: |
|
* ========== |
|
* |
|
*> \param[in] SIDE |
|
*> \verbatim |
|
*> SIDE is CHARACTER*1 |
|
*> = 'R': compute right eigenvectors only; |
|
*> = 'L': compute left eigenvectors only; |
|
*> = 'B': compute both right and left eigenvectors. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] EIGSRC |
|
*> \verbatim |
|
*> EIGSRC is CHARACTER*1 |
|
*> Specifies the source of eigenvalues supplied in W: |
|
*> = 'Q': the eigenvalues were found using ZHSEQR; thus, if |
|
*> H has zero subdiagonal elements, and so is |
|
*> block-triangular, then the j-th eigenvalue can be |
|
*> assumed to be an eigenvalue of the block containing |
|
*> the j-th row/column. This property allows ZHSEIN to |
|
*> perform inverse iteration on just one diagonal block. |
|
*> = 'N': no assumptions are made on the correspondence |
|
*> between eigenvalues and diagonal blocks. In this |
|
*> case, ZHSEIN must always perform inverse iteration |
|
*> using the whole matrix H. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] INITV |
|
*> \verbatim |
|
*> INITV is CHARACTER*1 |
|
*> = 'N': no initial vectors are supplied; |
|
*> = 'U': user-supplied initial vectors are stored in the arrays |
|
*> VL and/or VR. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] SELECT |
|
*> \verbatim |
|
*> SELECT is LOGICAL array, dimension (N) |
|
*> Specifies the eigenvectors to be computed. To select the |
|
*> eigenvector corresponding to the eigenvalue W(j), |
|
*> SELECT(j) must be set to .TRUE.. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] N |
|
*> \verbatim |
|
*> N is INTEGER |
|
*> The order of the matrix H. N >= 0. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] H |
|
*> \verbatim |
|
*> H is COMPLEX*16 array, dimension (LDH,N) |
|
*> The upper Hessenberg matrix H. |
|
*> If a NaN is detected in H, the routine will return with INFO=-6. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDH |
|
*> \verbatim |
|
*> LDH is INTEGER |
|
*> The leading dimension of the array H. LDH >= max(1,N). |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] W |
|
*> \verbatim |
|
*> W is COMPLEX*16 array, dimension (N) |
|
*> On entry, the eigenvalues of H. |
|
*> On exit, the real parts of W may have been altered since |
|
*> close eigenvalues are perturbed slightly in searching for |
|
*> independent eigenvectors. |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] VL |
|
*> \verbatim |
|
*> VL is COMPLEX*16 array, dimension (LDVL,MM) |
|
*> On entry, if INITV = 'U' and SIDE = 'L' or 'B', VL must |
|
*> contain starting vectors for the inverse iteration for the |
|
*> left eigenvectors; the starting vector for each eigenvector |
|
*> must be in the same column in which the eigenvector will be |
|
*> stored. |
|
*> On exit, if SIDE = 'L' or 'B', the left eigenvectors |
|
*> specified by SELECT will be stored consecutively in the |
|
*> columns of VL, in the same order as their eigenvalues. |
|
*> If SIDE = 'R', VL is not referenced. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDVL |
|
*> \verbatim |
|
*> LDVL is INTEGER |
|
*> The leading dimension of the array VL. |
|
*> LDVL >= max(1,N) if SIDE = 'L' or 'B'; LDVL >= 1 otherwise. |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] VR |
|
*> \verbatim |
|
*> VR is COMPLEX*16 array, dimension (LDVR,MM) |
|
*> On entry, if INITV = 'U' and SIDE = 'R' or 'B', VR must |
|
*> contain starting vectors for the inverse iteration for the |
|
*> right eigenvectors; the starting vector for each eigenvector |
|
*> must be in the same column in which the eigenvector will be |
|
*> stored. |
|
*> On exit, if SIDE = 'R' or 'B', the right eigenvectors |
|
*> specified by SELECT will be stored consecutively in the |
|
*> columns of VR, in the same order as their eigenvalues. |
|
*> If SIDE = 'L', VR is not referenced. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDVR |
|
*> \verbatim |
|
*> LDVR is INTEGER |
|
*> The leading dimension of the array VR. |
|
*> LDVR >= max(1,N) if SIDE = 'R' or 'B'; LDVR >= 1 otherwise. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] MM |
|
*> \verbatim |
|
*> MM is INTEGER |
|
*> The number of columns in the arrays VL and/or VR. MM >= M. |
|
*> \endverbatim |
|
*> |
|
*> \param[out] M |
|
*> \verbatim |
|
*> M is INTEGER |
|
*> The number of columns in the arrays VL and/or VR required to |
|
*> store the eigenvectors (= the number of .TRUE. elements in |
|
*> SELECT). |
|
*> \endverbatim |
|
*> |
|
*> \param[out] WORK |
|
*> \verbatim |
|
*> WORK is COMPLEX*16 array, dimension (N*N) |
|
*> \endverbatim |
|
*> |
|
*> \param[out] RWORK |
|
*> \verbatim |
|
*> RWORK is DOUBLE PRECISION array, dimension (N) |
|
*> \endverbatim |
|
*> |
|
*> \param[out] IFAILL |
|
*> \verbatim |
|
*> IFAILL is INTEGER array, dimension (MM) |
|
*> If SIDE = 'L' or 'B', IFAILL(i) = j > 0 if the left |
|
*> eigenvector in the i-th column of VL (corresponding to the |
|
*> eigenvalue w(j)) failed to converge; IFAILL(i) = 0 if the |
|
*> eigenvector converged satisfactorily. |
|
*> If SIDE = 'R', IFAILL is not referenced. |
|
*> \endverbatim |
|
*> |
|
*> \param[out] IFAILR |
|
*> \verbatim |
|
*> IFAILR is INTEGER array, dimension (MM) |
|
*> If SIDE = 'R' or 'B', IFAILR(i) = j > 0 if the right |
|
*> eigenvector in the i-th column of VR (corresponding to the |
|
*> eigenvalue w(j)) failed to converge; IFAILR(i) = 0 if the |
|
*> eigenvector converged satisfactorily. |
|
*> If SIDE = 'L', IFAILR is not referenced. |
|
*> \endverbatim |
|
*> |
|
*> \param[out] INFO |
|
*> \verbatim |
|
*> INFO is INTEGER |
|
*> = 0: successful exit |
|
*> < 0: if INFO = -i, the i-th argument had an illegal value |
|
*> > 0: if INFO = i, i is the number of eigenvectors which |
|
*> failed to converge; see IFAILL and IFAILR for further |
|
*> details. |
|
*> \endverbatim |
|
* |
|
* Authors: |
|
* ======== |
|
* |
|
*> \author Univ. of Tennessee |
|
*> \author Univ. of California Berkeley |
|
*> \author Univ. of Colorado Denver |
|
*> \author NAG Ltd. |
|
* |
|
*> \date November 2013 |
|
* |
|
*> \ingroup complex16OTHERcomputational |
|
* |
|
*> \par Further Details: |
|
* ===================== |
|
*> |
|
*> \verbatim |
|
*> |
|
*> Each eigenvector is normalized so that the element of largest |
|
*> magnitude has magnitude 1; here the magnitude of a complex number |
|
*> (x,y) is taken to be |x|+|y|. |
|
*> \endverbatim |
|
*> |
|
* ===================================================================== |
SUBROUTINE ZHSEIN( SIDE, EIGSRC, INITV, SELECT, N, H, LDH, W, VL, |
SUBROUTINE ZHSEIN( SIDE, EIGSRC, INITV, SELECT, N, H, LDH, W, VL, |
$ LDVL, VR, LDVR, MM, M, WORK, RWORK, IFAILL, |
$ LDVL, VR, LDVR, MM, M, WORK, RWORK, IFAILL, |
$ IFAILR, INFO ) |
$ IFAILR, INFO ) |
* |
* |
* -- LAPACK routine (version 3.2) -- |
* -- LAPACK computational routine (version 3.5.0) -- |
* -- LAPACK is a software package provided by Univ. of Tennessee, -- |
* -- LAPACK is a software package provided by Univ. of Tennessee, -- |
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- |
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- |
* November 2006 |
* November 2013 |
* |
* |
* .. Scalar Arguments .. |
* .. Scalar Arguments .. |
CHARACTER EIGSRC, INITV, SIDE |
CHARACTER EIGSRC, INITV, SIDE |
Line 19
|
Line 262
|
$ W( * ), WORK( * ) |
$ W( * ), WORK( * ) |
* .. |
* .. |
* |
* |
* Purpose |
|
* ======= |
|
* |
|
* ZHSEIN uses inverse iteration to find specified right and/or left |
|
* eigenvectors of a complex upper Hessenberg matrix H. |
|
* |
|
* The right eigenvector x and the left eigenvector y of the matrix H |
|
* corresponding to an eigenvalue w are defined by: |
|
* |
|
* H * x = w * x, y**h * H = w * y**h |
|
* |
|
* where y**h denotes the conjugate transpose of the vector y. |
|
* |
|
* Arguments |
|
* ========= |
|
* |
|
* SIDE (input) CHARACTER*1 |
|
* = 'R': compute right eigenvectors only; |
|
* = 'L': compute left eigenvectors only; |
|
* = 'B': compute both right and left eigenvectors. |
|
* |
|
* EIGSRC (input) CHARACTER*1 |
|
* Specifies the source of eigenvalues supplied in W: |
|
* = 'Q': the eigenvalues were found using ZHSEQR; thus, if |
|
* H has zero subdiagonal elements, and so is |
|
* block-triangular, then the j-th eigenvalue can be |
|
* assumed to be an eigenvalue of the block containing |
|
* the j-th row/column. This property allows ZHSEIN to |
|
* perform inverse iteration on just one diagonal block. |
|
* = 'N': no assumptions are made on the correspondence |
|
* between eigenvalues and diagonal blocks. In this |
|
* case, ZHSEIN must always perform inverse iteration |
|
* using the whole matrix H. |
|
* |
|
* INITV (input) CHARACTER*1 |
|
* = 'N': no initial vectors are supplied; |
|
* = 'U': user-supplied initial vectors are stored in the arrays |
|
* VL and/or VR. |
|
* |
|
* SELECT (input) LOGICAL array, dimension (N) |
|
* Specifies the eigenvectors to be computed. To select the |
|
* eigenvector corresponding to the eigenvalue W(j), |
|
* SELECT(j) must be set to .TRUE.. |
|
* |
|
* N (input) INTEGER |
|
* The order of the matrix H. N >= 0. |
|
* |
|
* H (input) COMPLEX*16 array, dimension (LDH,N) |
|
* The upper Hessenberg matrix H. |
|
* |
|
* LDH (input) INTEGER |
|
* The leading dimension of the array H. LDH >= max(1,N). |
|
* |
|
* W (input/output) COMPLEX*16 array, dimension (N) |
|
* On entry, the eigenvalues of H. |
|
* On exit, the real parts of W may have been altered since |
|
* close eigenvalues are perturbed slightly in searching for |
|
* independent eigenvectors. |
|
* |
|
* VL (input/output) COMPLEX*16 array, dimension (LDVL,MM) |
|
* On entry, if INITV = 'U' and SIDE = 'L' or 'B', VL must |
|
* contain starting vectors for the inverse iteration for the |
|
* left eigenvectors; the starting vector for each eigenvector |
|
* must be in the same column in which the eigenvector will be |
|
* stored. |
|
* On exit, if SIDE = 'L' or 'B', the left eigenvectors |
|
* specified by SELECT will be stored consecutively in the |
|
* columns of VL, in the same order as their eigenvalues. |
|
* If SIDE = 'R', VL is not referenced. |
|
* |
|
* LDVL (input) INTEGER |
|
* The leading dimension of the array VL. |
|
* LDVL >= max(1,N) if SIDE = 'L' or 'B'; LDVL >= 1 otherwise. |
|
* |
|
* VR (input/output) COMPLEX*16 array, dimension (LDVR,MM) |
|
* On entry, if INITV = 'U' and SIDE = 'R' or 'B', VR must |
|
* contain starting vectors for the inverse iteration for the |
|
* right eigenvectors; the starting vector for each eigenvector |
|
* must be in the same column in which the eigenvector will be |
|
* stored. |
|
* On exit, if SIDE = 'R' or 'B', the right eigenvectors |
|
* specified by SELECT will be stored consecutively in the |
|
* columns of VR, in the same order as their eigenvalues. |
|
* If SIDE = 'L', VR is not referenced. |
|
* |
|
* LDVR (input) INTEGER |
|
* The leading dimension of the array VR. |
|
* LDVR >= max(1,N) if SIDE = 'R' or 'B'; LDVR >= 1 otherwise. |
|
* |
|
* MM (input) INTEGER |
|
* The number of columns in the arrays VL and/or VR. MM >= M. |
|
* |
|
* M (output) INTEGER |
|
* The number of columns in the arrays VL and/or VR required to |
|
* store the eigenvectors (= the number of .TRUE. elements in |
|
* SELECT). |
|
* |
|
* WORK (workspace) COMPLEX*16 array, dimension (N*N) |
|
* |
|
* RWORK (workspace) DOUBLE PRECISION array, dimension (N) |
|
* |
|
* IFAILL (output) INTEGER array, dimension (MM) |
|
* If SIDE = 'L' or 'B', IFAILL(i) = j > 0 if the left |
|
* eigenvector in the i-th column of VL (corresponding to the |
|
* eigenvalue w(j)) failed to converge; IFAILL(i) = 0 if the |
|
* eigenvector converged satisfactorily. |
|
* If SIDE = 'R', IFAILL is not referenced. |
|
* |
|
* IFAILR (output) INTEGER array, dimension (MM) |
|
* If SIDE = 'R' or 'B', IFAILR(i) = j > 0 if the right |
|
* eigenvector in the i-th column of VR (corresponding to the |
|
* eigenvalue w(j)) failed to converge; IFAILR(i) = 0 if the |
|
* eigenvector converged satisfactorily. |
|
* If SIDE = 'L', IFAILR is not referenced. |
|
* |
|
* INFO (output) INTEGER |
|
* = 0: successful exit |
|
* < 0: if INFO = -i, the i-th argument had an illegal value |
|
* > 0: if INFO = i, i is the number of eigenvectors which |
|
* failed to converge; see IFAILL and IFAILR for further |
|
* details. |
|
* |
|
* Further Details |
|
* =============== |
|
* |
|
* Each eigenvector is normalized so that the element of largest |
|
* magnitude has magnitude 1; here the magnitude of a complex number |
|
* (x,y) is taken to be |x|+|y|. |
|
* |
|
* ===================================================================== |
* ===================================================================== |
* |
* |
* .. Parameters .. |
* .. Parameters .. |
Line 163
|
Line 277
|
COMPLEX*16 CDUM, WK |
COMPLEX*16 CDUM, WK |
* .. |
* .. |
* .. External Functions .. |
* .. External Functions .. |
LOGICAL LSAME |
LOGICAL LSAME, DISNAN |
DOUBLE PRECISION DLAMCH, ZLANHS |
DOUBLE PRECISION DLAMCH, ZLANHS |
EXTERNAL LSAME, DLAMCH, ZLANHS |
EXTERNAL LSAME, DLAMCH, ZLANHS, DISNAN |
* .. |
* .. |
* .. External Subroutines .. |
* .. External Subroutines .. |
EXTERNAL XERBLA, ZLAEIN |
EXTERNAL XERBLA, ZLAEIN |
Line 286
|
Line 400
|
* has not ben computed before. |
* has not ben computed before. |
* |
* |
HNORM = ZLANHS( 'I', KR-KL+1, H( KL, KL ), LDH, RWORK ) |
HNORM = ZLANHS( 'I', KR-KL+1, H( KL, KL ), LDH, RWORK ) |
IF( HNORM.GT.RZERO ) THEN |
IF( DISNAN( HNORM ) ) THEN |
|
INFO = -6 |
|
RETURN |
|
ELSE IF( HNORM.GT.RZERO ) THEN |
EPS3 = HNORM*ULP |
EPS3 = HNORM*ULP |
ELSE |
ELSE |
EPS3 = SMLNUM |
EPS3 = SMLNUM |