version 1.3, 2010/08/06 15:28:58
|
version 1.14, 2016/08/27 15:34:59
|
Line 1
|
Line 1
|
|
*> \brief \b ZLAQR4 computes the eigenvalues of a Hessenberg matrix, and optionally the matrices from the Schur decomposition. |
|
* |
|
* =========== DOCUMENTATION =========== |
|
* |
|
* Online html documentation available at |
|
* http://www.netlib.org/lapack/explore-html/ |
|
* |
|
*> \htmlonly |
|
*> Download ZLAQR4 + dependencies |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlaqr4.f"> |
|
*> [TGZ]</a> |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlaqr4.f"> |
|
*> [ZIP]</a> |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlaqr4.f"> |
|
*> [TXT]</a> |
|
*> \endhtmlonly |
|
* |
|
* Definition: |
|
* =========== |
|
* |
|
* SUBROUTINE ZLAQR4( WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, |
|
* IHIZ, Z, LDZ, WORK, LWORK, INFO ) |
|
* |
|
* .. Scalar Arguments .. |
|
* INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, LWORK, N |
|
* LOGICAL WANTT, WANTZ |
|
* .. |
|
* .. Array Arguments .. |
|
* COMPLEX*16 H( LDH, * ), W( * ), WORK( * ), Z( LDZ, * ) |
|
* .. |
|
* |
|
* |
|
*> \par Purpose: |
|
* ============= |
|
*> |
|
*> \verbatim |
|
*> |
|
*> ZLAQR4 implements one level of recursion for ZLAQR0. |
|
*> It is a complete implementation of the small bulge multi-shift |
|
*> QR algorithm. It may be called by ZLAQR0 and, for large enough |
|
*> deflation window size, it may be called by ZLAQR3. This |
|
*> subroutine is identical to ZLAQR0 except that it calls ZLAQR2 |
|
*> instead of ZLAQR3. |
|
*> |
|
*> ZLAQR4 computes the eigenvalues of a Hessenberg matrix H |
|
*> and, optionally, the matrices T and Z from the Schur decomposition |
|
*> H = Z T Z**H, where T is an upper triangular matrix (the |
|
*> Schur form), and Z is the unitary matrix of Schur vectors. |
|
*> |
|
*> Optionally Z may be postmultiplied into an input unitary |
|
*> matrix Q so that this routine can give the Schur factorization |
|
*> of a matrix A which has been reduced to the Hessenberg form H |
|
*> by the unitary matrix Q: A = Q*H*Q**H = (QZ)*H*(QZ)**H. |
|
*> \endverbatim |
|
* |
|
* Arguments: |
|
* ========== |
|
* |
|
*> \param[in] WANTT |
|
*> \verbatim |
|
*> WANTT is LOGICAL |
|
*> = .TRUE. : the full Schur form T is required; |
|
*> = .FALSE.: only eigenvalues are required. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] WANTZ |
|
*> \verbatim |
|
*> WANTZ is LOGICAL |
|
*> = .TRUE. : the matrix of Schur vectors Z is required; |
|
*> = .FALSE.: Schur vectors are not required. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] N |
|
*> \verbatim |
|
*> N is INTEGER |
|
*> The order of the matrix H. N .GE. 0. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] ILO |
|
*> \verbatim |
|
*> ILO is INTEGER |
|
*> \endverbatim |
|
*> |
|
*> \param[in] IHI |
|
*> \verbatim |
|
*> IHI is INTEGER |
|
*> It is assumed that H is already upper triangular in rows |
|
*> and columns 1:ILO-1 and IHI+1:N and, if ILO.GT.1, |
|
*> H(ILO,ILO-1) is zero. ILO and IHI are normally set by a |
|
*> previous call to ZGEBAL, and then passed to ZGEHRD when the |
|
*> matrix output by ZGEBAL is reduced to Hessenberg form. |
|
*> Otherwise, ILO and IHI should be set to 1 and N, |
|
*> respectively. If N.GT.0, then 1.LE.ILO.LE.IHI.LE.N. |
|
*> If N = 0, then ILO = 1 and IHI = 0. |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] H |
|
*> \verbatim |
|
*> H is COMPLEX*16 array, dimension (LDH,N) |
|
*> On entry, the upper Hessenberg matrix H. |
|
*> On exit, if INFO = 0 and WANTT is .TRUE., then H |
|
*> contains the upper triangular matrix T from the Schur |
|
*> decomposition (the Schur form). If INFO = 0 and WANT is |
|
*> .FALSE., then the contents of H are unspecified on exit. |
|
*> (The output value of H when INFO.GT.0 is given under the |
|
*> description of INFO below.) |
|
*> |
|
*> This subroutine may explicitly set H(i,j) = 0 for i.GT.j and |
|
*> j = 1, 2, ... ILO-1 or j = IHI+1, IHI+2, ... N. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDH |
|
*> \verbatim |
|
*> LDH is INTEGER |
|
*> The leading dimension of the array H. LDH .GE. max(1,N). |
|
*> \endverbatim |
|
*> |
|
*> \param[out] W |
|
*> \verbatim |
|
*> W is COMPLEX*16 array, dimension (N) |
|
*> The computed eigenvalues of H(ILO:IHI,ILO:IHI) are stored |
|
*> in W(ILO:IHI). If WANTT is .TRUE., then the eigenvalues are |
|
*> stored in the same order as on the diagonal of the Schur |
|
*> form returned in H, with W(i) = H(i,i). |
|
*> \endverbatim |
|
*> |
|
*> \param[in] ILOZ |
|
*> \verbatim |
|
*> ILOZ is INTEGER |
|
*> \endverbatim |
|
*> |
|
*> \param[in] IHIZ |
|
*> \verbatim |
|
*> IHIZ is INTEGER |
|
*> Specify the rows of Z to which transformations must be |
|
*> applied if WANTZ is .TRUE.. |
|
*> 1 .LE. ILOZ .LE. ILO; IHI .LE. IHIZ .LE. N. |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] Z |
|
*> \verbatim |
|
*> Z is COMPLEX*16 array, dimension (LDZ,IHI) |
|
*> If WANTZ is .FALSE., then Z is not referenced. |
|
*> If WANTZ is .TRUE., then Z(ILO:IHI,ILOZ:IHIZ) is |
|
*> replaced by Z(ILO:IHI,ILOZ:IHIZ)*U where U is the |
|
*> orthogonal Schur factor of H(ILO:IHI,ILO:IHI). |
|
*> (The output value of Z when INFO.GT.0 is given under |
|
*> the description of INFO below.) |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDZ |
|
*> \verbatim |
|
*> LDZ is INTEGER |
|
*> The leading dimension of the array Z. if WANTZ is .TRUE. |
|
*> then LDZ.GE.MAX(1,IHIZ). Otherwize, LDZ.GE.1. |
|
*> \endverbatim |
|
*> |
|
*> \param[out] WORK |
|
*> \verbatim |
|
*> WORK is COMPLEX*16 array, dimension LWORK |
|
*> On exit, if LWORK = -1, WORK(1) returns an estimate of |
|
*> the optimal value for LWORK. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LWORK |
|
*> \verbatim |
|
*> LWORK is INTEGER |
|
*> The dimension of the array WORK. LWORK .GE. max(1,N) |
|
*> is sufficient, but LWORK typically as large as 6*N may |
|
*> be required for optimal performance. A workspace query |
|
*> to determine the optimal workspace size is recommended. |
|
*> |
|
*> If LWORK = -1, then ZLAQR4 does a workspace query. |
|
*> In this case, ZLAQR4 checks the input parameters and |
|
*> estimates the optimal workspace size for the given |
|
*> values of N, ILO and IHI. The estimate is returned |
|
*> in WORK(1). No error message related to LWORK is |
|
*> issued by XERBLA. Neither H nor Z are accessed. |
|
*> \endverbatim |
|
*> |
|
*> \param[out] INFO |
|
*> \verbatim |
|
*> INFO is INTEGER |
|
*> = 0: successful exit |
|
*> .GT. 0: if INFO = i, ZLAQR4 failed to compute all of |
|
*> the eigenvalues. Elements 1:ilo-1 and i+1:n of WR |
|
*> and WI contain those eigenvalues which have been |
|
*> successfully computed. (Failures are rare.) |
|
*> |
|
*> If INFO .GT. 0 and WANT is .FALSE., then on exit, |
|
*> the remaining unconverged eigenvalues are the eigen- |
|
*> values of the upper Hessenberg matrix rows and |
|
*> columns ILO through INFO of the final, output |
|
*> value of H. |
|
*> |
|
*> If INFO .GT. 0 and WANTT is .TRUE., then on exit |
|
*> |
|
*> (*) (initial value of H)*U = U*(final value of H) |
|
*> |
|
*> where U is a unitary matrix. The final |
|
*> value of H is upper Hessenberg and triangular in |
|
*> rows and columns INFO+1 through IHI. |
|
*> |
|
*> If INFO .GT. 0 and WANTZ is .TRUE., then on exit |
|
*> |
|
*> (final value of Z(ILO:IHI,ILOZ:IHIZ) |
|
*> = (initial value of Z(ILO:IHI,ILOZ:IHIZ)*U |
|
*> |
|
*> where U is the unitary matrix in (*) (regard- |
|
*> less of the value of WANTT.) |
|
*> |
|
*> If INFO .GT. 0 and WANTZ is .FALSE., then Z is not |
|
*> accessed. |
|
*> \endverbatim |
|
* |
|
* Authors: |
|
* ======== |
|
* |
|
*> \author Univ. of Tennessee |
|
*> \author Univ. of California Berkeley |
|
*> \author Univ. of Colorado Denver |
|
*> \author NAG Ltd. |
|
* |
|
*> \date September 2012 |
|
* |
|
*> \ingroup complex16OTHERauxiliary |
|
* |
|
*> \par Contributors: |
|
* ================== |
|
*> |
|
*> Karen Braman and Ralph Byers, Department of Mathematics, |
|
*> University of Kansas, USA |
|
* |
|
*> \par References: |
|
* ================ |
|
*> |
|
*> K. Braman, R. Byers and R. Mathias, The Multi-Shift QR |
|
*> Algorithm Part I: Maintaining Well Focused Shifts, and Level 3 |
|
*> Performance, SIAM Journal of Matrix Analysis, volume 23, pages |
|
*> 929--947, 2002. |
|
*> \n |
|
*> K. Braman, R. Byers and R. Mathias, The Multi-Shift QR |
|
*> Algorithm Part II: Aggressive Early Deflation, SIAM Journal |
|
*> of Matrix Analysis, volume 23, pages 948--973, 2002. |
|
*> |
|
* ===================================================================== |
SUBROUTINE ZLAQR4( WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, |
SUBROUTINE ZLAQR4( WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, |
$ IHIZ, Z, LDZ, WORK, LWORK, INFO ) |
$ IHIZ, Z, LDZ, WORK, LWORK, INFO ) |
* |
* |
* -- LAPACK auxiliary routine (version 3.2) -- |
* -- LAPACK auxiliary routine (version 3.4.2) -- |
* Univ. of Tennessee, Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd.. |
* -- LAPACK is a software package provided by Univ. of Tennessee, -- |
* November 2006 |
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- |
|
* September 2012 |
* |
* |
* .. Scalar Arguments .. |
* .. Scalar Arguments .. |
INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, LWORK, N |
INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, LWORK, N |
Line 13
|
Line 260
|
COMPLEX*16 H( LDH, * ), W( * ), WORK( * ), Z( LDZ, * ) |
COMPLEX*16 H( LDH, * ), W( * ), WORK( * ), Z( LDZ, * ) |
* .. |
* .. |
* |
* |
* This subroutine implements one level of recursion for ZLAQR0. |
* ================================================================ |
* It is a complete implementation of the small bulge multi-shift |
|
* QR algorithm. It may be called by ZLAQR0 and, for large enough |
|
* deflation window size, it may be called by ZLAQR3. This |
|
* subroutine is identical to ZLAQR0 except that it calls ZLAQR2 |
|
* instead of ZLAQR3. |
|
* |
|
* Purpose |
|
* ======= |
|
* |
|
* ZLAQR4 computes the eigenvalues of a Hessenberg matrix H |
|
* and, optionally, the matrices T and Z from the Schur decomposition |
|
* H = Z T Z**H, where T is an upper triangular matrix (the |
|
* Schur form), and Z is the unitary matrix of Schur vectors. |
|
* |
|
* Optionally Z may be postmultiplied into an input unitary |
|
* matrix Q so that this routine can give the Schur factorization |
|
* of a matrix A which has been reduced to the Hessenberg form H |
|
* by the unitary matrix Q: A = Q*H*Q**H = (QZ)*H*(QZ)**H. |
|
* |
|
* Arguments |
|
* ========= |
|
* |
|
* WANTT (input) LOGICAL |
|
* = .TRUE. : the full Schur form T is required; |
|
* = .FALSE.: only eigenvalues are required. |
|
* |
|
* WANTZ (input) LOGICAL |
|
* = .TRUE. : the matrix of Schur vectors Z is required; |
|
* = .FALSE.: Schur vectors are not required. |
|
* |
|
* N (input) INTEGER |
|
* The order of the matrix H. N .GE. 0. |
|
* |
|
* ILO (input) INTEGER |
|
* IHI (input) INTEGER |
|
* It is assumed that H is already upper triangular in rows |
|
* and columns 1:ILO-1 and IHI+1:N and, if ILO.GT.1, |
|
* H(ILO,ILO-1) is zero. ILO and IHI are normally set by a |
|
* previous call to ZGEBAL, and then passed to ZGEHRD when the |
|
* matrix output by ZGEBAL is reduced to Hessenberg form. |
|
* Otherwise, ILO and IHI should be set to 1 and N, |
|
* respectively. If N.GT.0, then 1.LE.ILO.LE.IHI.LE.N. |
|
* If N = 0, then ILO = 1 and IHI = 0. |
|
* |
|
* H (input/output) COMPLEX*16 array, dimension (LDH,N) |
|
* On entry, the upper Hessenberg matrix H. |
|
* On exit, if INFO = 0 and WANTT is .TRUE., then H |
|
* contains the upper triangular matrix T from the Schur |
|
* decomposition (the Schur form). If INFO = 0 and WANT is |
|
* .FALSE., then the contents of H are unspecified on exit. |
|
* (The output value of H when INFO.GT.0 is given under the |
|
* description of INFO below.) |
|
* |
|
* This subroutine may explicitly set H(i,j) = 0 for i.GT.j and |
|
* j = 1, 2, ... ILO-1 or j = IHI+1, IHI+2, ... N. |
|
* |
|
* LDH (input) INTEGER |
|
* The leading dimension of the array H. LDH .GE. max(1,N). |
|
* |
|
* W (output) COMPLEX*16 array, dimension (N) |
|
* The computed eigenvalues of H(ILO:IHI,ILO:IHI) are stored |
|
* in W(ILO:IHI). If WANTT is .TRUE., then the eigenvalues are |
|
* stored in the same order as on the diagonal of the Schur |
|
* form returned in H, with W(i) = H(i,i). |
|
* |
|
* Z (input/output) COMPLEX*16 array, dimension (LDZ,IHI) |
|
* If WANTZ is .FALSE., then Z is not referenced. |
|
* If WANTZ is .TRUE., then Z(ILO:IHI,ILOZ:IHIZ) is |
|
* replaced by Z(ILO:IHI,ILOZ:IHIZ)*U where U is the |
|
* orthogonal Schur factor of H(ILO:IHI,ILO:IHI). |
|
* (The output value of Z when INFO.GT.0 is given under |
|
* the description of INFO below.) |
|
* |
|
* LDZ (input) INTEGER |
|
* The leading dimension of the array Z. if WANTZ is .TRUE. |
|
* then LDZ.GE.MAX(1,IHIZ). Otherwize, LDZ.GE.1. |
|
* |
|
* WORK (workspace/output) COMPLEX*16 array, dimension LWORK |
|
* On exit, if LWORK = -1, WORK(1) returns an estimate of |
|
* the optimal value for LWORK. |
|
* |
|
* LWORK (input) INTEGER |
|
* The dimension of the array WORK. LWORK .GE. max(1,N) |
|
* is sufficient, but LWORK typically as large as 6*N may |
|
* be required for optimal performance. A workspace query |
|
* to determine the optimal workspace size is recommended. |
|
* |
|
* If LWORK = -1, then ZLAQR4 does a workspace query. |
|
* In this case, ZLAQR4 checks the input parameters and |
|
* estimates the optimal workspace size for the given |
|
* values of N, ILO and IHI. The estimate is returned |
|
* in WORK(1). No error message related to LWORK is |
|
* issued by XERBLA. Neither H nor Z are accessed. |
|
* |
|
* |
|
* INFO (output) INTEGER |
|
* = 0: successful exit |
|
* .GT. 0: if INFO = i, ZLAQR4 failed to compute all of |
|
* the eigenvalues. Elements 1:ilo-1 and i+1:n of WR |
|
* and WI contain those eigenvalues which have been |
|
* successfully computed. (Failures are rare.) |
|
* |
|
* If INFO .GT. 0 and WANT is .FALSE., then on exit, |
|
* the remaining unconverged eigenvalues are the eigen- |
|
* values of the upper Hessenberg matrix rows and |
|
* columns ILO through INFO of the final, output |
|
* value of H. |
|
* |
|
* If INFO .GT. 0 and WANTT is .TRUE., then on exit |
|
* |
|
* (*) (initial value of H)*U = U*(final value of H) |
|
* |
|
* where U is a unitary matrix. The final |
|
* value of H is upper Hessenberg and triangular in |
|
* rows and columns INFO+1 through IHI. |
|
* |
|
* If INFO .GT. 0 and WANTZ is .TRUE., then on exit |
|
* |
|
* (final value of Z(ILO:IHI,ILOZ:IHIZ) |
|
* = (initial value of Z(ILO:IHI,ILOZ:IHIZ)*U |
|
* |
|
* where U is the unitary matrix in (*) (regard- |
|
* less of the value of WANTT.) |
|
* |
|
* If INFO .GT. 0 and WANTZ is .FALSE., then Z is not |
|
* accessed. |
|
* |
|
* ================================================================ |
|
* Based on contributions by |
|
* Karen Braman and Ralph Byers, Department of Mathematics, |
|
* University of Kansas, USA |
|
* |
|
* ================================================================ |
|
* References: |
|
* K. Braman, R. Byers and R. Mathias, The Multi-Shift QR |
|
* Algorithm Part I: Maintaining Well Focused Shifts, and Level 3 |
|
* Performance, SIAM Journal of Matrix Analysis, volume 23, pages |
|
* 929--947, 2002. |
|
* |
|
* K. Braman, R. Byers and R. Mathias, The Multi-Shift QR |
|
* Algorithm Part II: Aggressive Early Deflation, SIAM Journal |
|
* of Matrix Analysis, volume 23, pages 948--973, 2002. |
|
* |
* |
* ================================================================ |
|
* .. Parameters .. |
* .. Parameters .. |
* |
* |
* ==== Matrices of order NTINY or smaller must be processed by |
* ==== Matrices of order NTINY or smaller must be processed by |