--- rpl/lapack/lapack/dtrevc.f 2010/01/26 15:22:46 1.1 +++ rpl/lapack/lapack/dtrevc.f 2018/05/29 06:55:21 1.17 @@ -1,10 +1,231 @@ +*> \brief \b DTREVC +* +* =========== DOCUMENTATION =========== +* +* Online html documentation available at +* http://www.netlib.org/lapack/explore-html/ +* +*> \htmlonly +*> Download DTREVC + dependencies +*> +*> [TGZ] +*> +*> [ZIP] +*> +*> [TXT] +*> \endhtmlonly +* +* Definition: +* =========== +* +* SUBROUTINE DTREVC( SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, +* LDVR, MM, M, WORK, INFO ) +* +* .. Scalar Arguments .. +* CHARACTER HOWMNY, SIDE +* INTEGER INFO, LDT, LDVL, LDVR, M, MM, N +* .. +* .. Array Arguments .. +* LOGICAL SELECT( * ) +* DOUBLE PRECISION T( LDT, * ), VL( LDVL, * ), VR( LDVR, * ), +* $ WORK( * ) +* .. +* +* +*> \par Purpose: +* ============= +*> +*> \verbatim +*> +*> DTREVC computes some or all of the right and/or left eigenvectors of +*> a real upper quasi-triangular matrix T. +*> Matrices of this type are produced by the Schur factorization of +*> a real general matrix: A = Q*T*Q**T, as computed by DHSEQR. +*> +*> The right eigenvector x and the left eigenvector y of T corresponding +*> to an eigenvalue w are defined by: +*> +*> T*x = w*x, (y**H)*T = w*(y**H) +*> +*> where y**H denotes the conjugate transpose of y. +*> The eigenvalues are not input to this routine, but are read directly +*> from the diagonal blocks of T. +*> +*> This routine returns the matrices X and/or Y of right and left +*> eigenvectors of T, or the products Q*X and/or Q*Y, where Q is an +*> input matrix. If Q is the orthogonal factor that reduces a matrix +*> A to Schur form T, then Q*X and Q*Y are the matrices of right and +*> left eigenvectors of A. +*> \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] HOWMNY +*> \verbatim +*> HOWMNY is CHARACTER*1 +*> = 'A': compute all right and/or left eigenvectors; +*> = 'B': compute all right and/or left eigenvectors, +*> backtransformed by the matrices in VR and/or VL; +*> = 'S': compute selected right and/or left eigenvectors, +*> as indicated by the logical array SELECT. +*> \endverbatim +*> +*> \param[in,out] SELECT +*> \verbatim +*> SELECT is LOGICAL array, dimension (N) +*> If HOWMNY = 'S', SELECT specifies the eigenvectors to be +*> computed. +*> If w(j) is a real eigenvalue, the corresponding real +*> eigenvector is computed if SELECT(j) is .TRUE.. +*> If w(j) and w(j+1) are the real and imaginary parts of a +*> complex eigenvalue, the corresponding complex eigenvector is +*> computed if either SELECT(j) or SELECT(j+1) is .TRUE., and +*> on exit SELECT(j) is set to .TRUE. and SELECT(j+1) is set to +*> .FALSE.. +*> Not referenced if HOWMNY = 'A' or 'B'. +*> \endverbatim +*> +*> \param[in] N +*> \verbatim +*> N is INTEGER +*> The order of the matrix T. N >= 0. +*> \endverbatim +*> +*> \param[in] T +*> \verbatim +*> T is DOUBLE PRECISION array, dimension (LDT,N) +*> The upper quasi-triangular matrix T in Schur canonical form. +*> \endverbatim +*> +*> \param[in] LDT +*> \verbatim +*> LDT is INTEGER +*> The leading dimension of the array T. LDT >= max(1,N). +*> \endverbatim +*> +*> \param[in,out] VL +*> \verbatim +*> VL is DOUBLE PRECISION array, dimension (LDVL,MM) +*> On entry, if SIDE = 'L' or 'B' and HOWMNY = 'B', VL must +*> contain an N-by-N matrix Q (usually the orthogonal matrix Q +*> of Schur vectors returned by DHSEQR). +*> On exit, if SIDE = 'L' or 'B', VL contains: +*> if HOWMNY = 'A', the matrix Y of left eigenvectors of T; +*> if HOWMNY = 'B', the matrix Q*Y; +*> if HOWMNY = 'S', the left eigenvectors of T specified by +*> SELECT, stored consecutively in the columns +*> of VL, in the same order as their +*> eigenvalues. +*> A complex eigenvector corresponding to a complex eigenvalue +*> is stored in two consecutive columns, the first holding the +*> real part, and the second the imaginary part. +*> Not referenced if SIDE = 'R'. +*> \endverbatim +*> +*> \param[in] LDVL +*> \verbatim +*> LDVL is INTEGER +*> The leading dimension of the array VL. LDVL >= 1, and if +*> SIDE = 'L' or 'B', LDVL >= N. +*> \endverbatim +*> +*> \param[in,out] VR +*> \verbatim +*> VR is DOUBLE PRECISION array, dimension (LDVR,MM) +*> On entry, if SIDE = 'R' or 'B' and HOWMNY = 'B', VR must +*> contain an N-by-N matrix Q (usually the orthogonal matrix Q +*> of Schur vectors returned by DHSEQR). +*> On exit, if SIDE = 'R' or 'B', VR contains: +*> if HOWMNY = 'A', the matrix X of right eigenvectors of T; +*> if HOWMNY = 'B', the matrix Q*X; +*> if HOWMNY = 'S', the right eigenvectors of T specified by +*> SELECT, stored consecutively in the columns +*> of VR, in the same order as their +*> eigenvalues. +*> A complex eigenvector corresponding to a complex eigenvalue +*> is stored in two consecutive columns, the first holding the +*> real part and the second the imaginary part. +*> Not referenced if SIDE = 'L'. +*> \endverbatim +*> +*> \param[in] LDVR +*> \verbatim +*> LDVR is INTEGER +*> The leading dimension of the array VR. LDVR >= 1, and if +*> SIDE = 'R' or 'B', LDVR >= N. +*> \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 actually +*> used to store the eigenvectors. +*> If HOWMNY = 'A' or 'B', M is set to N. +*> Each selected real eigenvector occupies one column and each +*> selected complex eigenvector occupies two columns. +*> \endverbatim +*> +*> \param[out] WORK +*> \verbatim +*> WORK is DOUBLE PRECISION array, dimension (3*N) +*> \endverbatim +*> +*> \param[out] INFO +*> \verbatim +*> INFO is INTEGER +*> = 0: successful exit +*> < 0: if INFO = -i, the i-th argument had an illegal value +*> \endverbatim +* +* Authors: +* ======== +* +*> \author Univ. of Tennessee +*> \author Univ. of California Berkeley +*> \author Univ. of Colorado Denver +*> \author NAG Ltd. +* +*> \date November 2017 +* +*> \ingroup doubleOTHERcomputational +* +*> \par Further Details: +* ===================== +*> +*> \verbatim +*> +*> The algorithm used in this program is basically backward (forward) +*> substitution, with scaling to make the the code robust against +*> possible overflow. +*> +*> 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 DTREVC( SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, $ LDVR, MM, M, WORK, INFO ) * -* -- LAPACK routine (version 3.2) -- +* -- LAPACK computational routine (version 3.8.0) -- * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- -* November 2006 +* November 2017 * * .. Scalar Arguments .. CHARACTER HOWMNY, SIDE @@ -16,132 +237,6 @@ $ WORK( * ) * .. * -* Purpose -* ======= -* -* DTREVC computes some or all of the right and/or left eigenvectors of -* a real upper quasi-triangular matrix T. -* Matrices of this type are produced by the Schur factorization of -* a real general matrix: A = Q*T*Q**T, as computed by DHSEQR. -* -* The right eigenvector x and the left eigenvector y of T corresponding -* to an eigenvalue w are defined by: -* -* T*x = w*x, (y**H)*T = w*(y**H) -* -* where y**H denotes the conjugate transpose of y. -* The eigenvalues are not input to this routine, but are read directly -* from the diagonal blocks of T. -* -* This routine returns the matrices X and/or Y of right and left -* eigenvectors of T, or the products Q*X and/or Q*Y, where Q is an -* input matrix. If Q is the orthogonal factor that reduces a matrix -* A to Schur form T, then Q*X and Q*Y are the matrices of right and -* left eigenvectors of A. -* -* Arguments -* ========= -* -* SIDE (input) CHARACTER*1 -* = 'R': compute right eigenvectors only; -* = 'L': compute left eigenvectors only; -* = 'B': compute both right and left eigenvectors. -* -* HOWMNY (input) CHARACTER*1 -* = 'A': compute all right and/or left eigenvectors; -* = 'B': compute all right and/or left eigenvectors, -* backtransformed by the matrices in VR and/or VL; -* = 'S': compute selected right and/or left eigenvectors, -* as indicated by the logical array SELECT. -* -* SELECT (input/output) LOGICAL array, dimension (N) -* If HOWMNY = 'S', SELECT specifies the eigenvectors to be -* computed. -* If w(j) is a real eigenvalue, the corresponding real -* eigenvector is computed if SELECT(j) is .TRUE.. -* If w(j) and w(j+1) are the real and imaginary parts of a -* complex eigenvalue, the corresponding complex eigenvector is -* computed if either SELECT(j) or SELECT(j+1) is .TRUE., and -* on exit SELECT(j) is set to .TRUE. and SELECT(j+1) is set to -* .FALSE.. -* Not referenced if HOWMNY = 'A' or 'B'. -* -* N (input) INTEGER -* The order of the matrix T. N >= 0. -* -* T (input) DOUBLE PRECISION array, dimension (LDT,N) -* The upper quasi-triangular matrix T in Schur canonical form. -* -* LDT (input) INTEGER -* The leading dimension of the array T. LDT >= max(1,N). -* -* VL (input/output) DOUBLE PRECISION array, dimension (LDVL,MM) -* On entry, if SIDE = 'L' or 'B' and HOWMNY = 'B', VL must -* contain an N-by-N matrix Q (usually the orthogonal matrix Q -* of Schur vectors returned by DHSEQR). -* On exit, if SIDE = 'L' or 'B', VL contains: -* if HOWMNY = 'A', the matrix Y of left eigenvectors of T; -* if HOWMNY = 'B', the matrix Q*Y; -* if HOWMNY = 'S', the left eigenvectors of T specified by -* SELECT, stored consecutively in the columns -* of VL, in the same order as their -* eigenvalues. -* A complex eigenvector corresponding to a complex eigenvalue -* is stored in two consecutive columns, the first holding the -* real part, and the second the imaginary part. -* Not referenced if SIDE = 'R'. -* -* LDVL (input) INTEGER -* The leading dimension of the array VL. LDVL >= 1, and if -* SIDE = 'L' or 'B', LDVL >= N. -* -* VR (input/output) DOUBLE PRECISION array, dimension (LDVR,MM) -* On entry, if SIDE = 'R' or 'B' and HOWMNY = 'B', VR must -* contain an N-by-N matrix Q (usually the orthogonal matrix Q -* of Schur vectors returned by DHSEQR). -* On exit, if SIDE = 'R' or 'B', VR contains: -* if HOWMNY = 'A', the matrix X of right eigenvectors of T; -* if HOWMNY = 'B', the matrix Q*X; -* if HOWMNY = 'S', the right eigenvectors of T specified by -* SELECT, stored consecutively in the columns -* of VR, in the same order as their -* eigenvalues. -* A complex eigenvector corresponding to a complex eigenvalue -* is stored in two consecutive columns, the first holding the -* real part and the second the imaginary part. -* Not referenced if SIDE = 'L'. -* -* LDVR (input) INTEGER -* The leading dimension of the array VR. LDVR >= 1, and if -* SIDE = 'R' or 'B', LDVR >= N. -* -* 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 actually -* used to store the eigenvectors. -* If HOWMNY = 'A' or 'B', M is set to N. -* Each selected real eigenvector occupies one column and each -* selected complex eigenvector occupies two columns. -* -* WORK (workspace) DOUBLE PRECISION array, dimension (3*N) -* -* INFO (output) INTEGER -* = 0: successful exit -* < 0: if INFO = -i, the i-th argument had an illegal value -* -* Further Details -* =============== -* -* The algorithm used in this program is basically backward (forward) -* substitution, with scaling to make the the code robust against -* possible overflow. -* -* 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 .. @@ -162,7 +257,8 @@ EXTERNAL LSAME, IDAMAX, DDOT, DLAMCH * .. * .. External Subroutines .. - EXTERNAL DAXPY, DCOPY, DGEMV, DLALN2, DSCAL, XERBLA + EXTERNAL DLABAD, DAXPY, DCOPY, DGEMV, DLALN2, DSCAL, + $ XERBLA * .. * .. Intrinsic Functions .. INTRINSIC ABS, MAX, SQRT @@ -651,7 +747,7 @@ 160 CONTINUE * * Solve the quasi-triangular system: -* (T(KI+1:N,KI+1:N) - WR)'*X = SCALE*WORK +* (T(KI+1:N,KI+1:N) - WR)**T*X = SCALE*WORK * VMAX = ONE VCRIT = BIGNUM @@ -688,7 +784,7 @@ $ DDOT( J-KI-1, T( KI+1, J ), 1, $ WORK( KI+1+N ), 1 ) * -* Solve (T(J,J)-WR)'*X = WORK +* Solve (T(J,J)-WR)**T*X = WORK * CALL DLALN2( .FALSE., 1, 1, SMIN, ONE, T( J, J ), $ LDT, ONE, ONE, WORK( J+N ), N, WR, @@ -726,8 +822,8 @@ $ WORK( KI+1+N ), 1 ) * * Solve -* [T(J,J)-WR T(J,J+1) ]'* X = SCALE*( WORK1 ) -* [T(J+1,J) T(J+1,J+1)-WR] ( WORK2 ) +* [T(J,J)-WR T(J,J+1) ]**T * X = SCALE*( WORK1 ) +* [T(J+1,J) T(J+1,J+1)-WR] ( WORK2 ) * CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, T( J, J ), $ LDT, ONE, ONE, WORK( J+N ), N, WR, @@ -778,7 +874,7 @@ * Complex left eigenvector. * * Initial solve: -* ((T(KI,KI) T(KI,KI+1) )' - (WR - I* WI))*X = 0. +* ((T(KI,KI) T(KI,KI+1) )**T - (WR - I* WI))*X = 0. * ((T(KI+1,KI) T(KI+1,KI+1)) ) * IF( ABS( T( KI, KI+1 ) ).GE.ABS( T( KI+1, KI ) ) ) THEN @@ -891,8 +987,8 @@ $ WORK( KI+2+N2 ), 1 ) * * Solve 2-by-2 complex linear equation -* ([T(j,j) T(j,j+1) ]'-(wr-i*wi)*I)*X = SCALE*B -* ([T(j+1,j) T(j+1,j+1)] ) +* ([T(j,j) T(j,j+1) ]**T-(wr-i*wi)*I)*X = SCALE*B +* ([T(j+1,j) T(j+1,j+1)] ) * CALL DLALN2( .TRUE., 2, 2, SMIN, ONE, T( J, J ), $ LDT, ONE, ONE, WORK( J+N ), N, WR,