zhseqr.f.html

来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 420 行 · 第 1/3 页

</span><span class="comment">*</span><span class="comment">                columns ILO through INFO of the final, output
</span><span class="comment">*</span><span class="comment">                value of H.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                If INFO .GT. 0 and JOB   = 'S', then on exit
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           (*)  (initial value of H)*U  = U*(final value of H)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                where U is a unitary matrix.  The final
</span><span class="comment">*</span><span class="comment">                value of  H is upper Hessenberg and triangular in
</span><span class="comment">*</span><span class="comment">                rows and columns INFO+1 through IHI.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                If INFO .GT. 0 and COMPZ = 'V', then on exit
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                  (final value of Z)  =  (initial value of Z)*U
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                where U is the unitary matrix in (*) (regard-
</span><span class="comment">*</span><span class="comment">                less of the value of JOB.)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                If INFO .GT. 0 and COMPZ = 'I', then on exit
</span><span class="comment">*</span><span class="comment">                      (final value of Z)  = U
</span><span class="comment">*</span><span class="comment">                where U is the unitary matrix in (*) (regard-
</span><span class="comment">*</span><span class="comment">                less of the value of JOB.)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                If INFO .GT. 0 and COMPZ = 'N', then Z is not
</span><span class="comment">*</span><span class="comment">                accessed.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ================================================================
</span><span class="comment">*</span><span class="comment">             Default values supplied by
</span><span class="comment">*</span><span class="comment">             <a name="ILAENV.150"></a><a href="hfy-index.html#ILAENV">ILAENV</a>(ISPEC,'<a name="ZHSEQR.150"></a><a href="zhseqr.f.html#ZHSEQR.1">ZHSEQR</a>',JOB(:1)//COMPZ(:1),N,ILO,IHI,LWORK).
</span><span class="comment">*</span><span class="comment">             It is suggested that these defaults be adjusted in order
</span><span class="comment">*</span><span class="comment">             to attain best performance in each particular
</span><span class="comment">*</span><span class="comment">             computational environment.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">            ISPEC=1:  The <a name="ZLAHQR.155"></a><a href="zlahqr.f.html#ZLAHQR.1">ZLAHQR</a> vs <a name="ZLAQR0.155"></a><a href="zlaqr0.f.html#ZLAQR0.1">ZLAQR0</a> crossover point.
</span><span class="comment">*</span><span class="comment">                      Default: 75. (Must be at least 11.)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">            ISPEC=2:  Recommended deflation window size.
</span><span class="comment">*</span><span class="comment">                      This depends on ILO, IHI and NS.  NS is the
</span><span class="comment">*</span><span class="comment">                      number of simultaneous shifts returned
</span><span class="comment">*</span><span class="comment">                      by <a name="ILAENV.161"></a><a href="hfy-index.html#ILAENV">ILAENV</a>(ISPEC=4).  (See ISPEC=4 below.)
</span><span class="comment">*</span><span class="comment">                      The default for (IHI-ILO+1).LE.500 is NS.
</span><span class="comment">*</span><span class="comment">                      The default for (IHI-ILO+1).GT.500 is 3*NS/2.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">            ISPEC=3:  Nibble crossover point. (See <a name="ILAENV.165"></a><a href="hfy-index.html#ILAENV">ILAENV</a> for
</span><span class="comment">*</span><span class="comment">                      details.)  Default: 14% of deflation window
</span><span class="comment">*</span><span class="comment">                      size.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">            ISPEC=4:  Number of simultaneous shifts, NS, in
</span><span class="comment">*</span><span class="comment">                      a multi-shift QR iteration.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                      If IHI-ILO+1 is ...
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                      greater than      ...but less    ... the
</span><span class="comment">*</span><span class="comment">                      or equal to ...      than        default is
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                           1               30          NS -   2(+)
</span><span class="comment">*</span><span class="comment">                          30               60          NS -   4(+)
</span><span class="comment">*</span><span class="comment">                          60              150          NS =  10(+)
</span><span class="comment">*</span><span class="comment">                         150              590          NS =  **
</span><span class="comment">*</span><span class="comment">                         590             3000          NS =  64
</span><span class="comment">*</span><span class="comment">                        3000             6000          NS = 128
</span><span class="comment">*</span><span class="comment">                        6000             infinity      NS = 256
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                  (+)  By default some or all matrices of this order 
</span><span class="comment">*</span><span class="comment">                       are passed to the implicit double shift routine
</span><span class="comment">*</span><span class="comment">                       <a name="ZLAHQR.187"></a><a href="zlahqr.f.html#ZLAHQR.1">ZLAHQR</a> and NS is ignored.  See ISPEC=1 above 
</span><span class="comment">*</span><span class="comment">                       and comments in IPARM for details.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                       The asterisks (**) indicate an ad-hoc
</span><span class="comment">*</span><span class="comment">                       function of N increasing from 10 to 64.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">            ISPEC=5:  Select structured matrix multiply.
</span><span class="comment">*</span><span class="comment">                      (See <a name="ILAENV.194"></a><a href="hfy-index.html#ILAENV">ILAENV</a> for details.) Default: 3.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ================================================================
</span><span class="comment">*</span><span class="comment">     Based on contributions by
</span><span class="comment">*</span><span class="comment">        Karen Braman and Ralph Byers, Department of Mathematics,
</span><span class="comment">*</span><span class="comment">        University of Kansas, USA
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ================================================================
</span><span class="comment">*</span><span class="comment">     References:
</span><span class="comment">*</span><span class="comment">       K. Braman, R. Byers and R. Mathias, The Multi-Shift QR
</span><span class="comment">*</span><span class="comment">       Algorithm Part I: Maintaining Well Focused Shifts, and Level 3
</span><span class="comment">*</span><span class="comment">       Performance, SIAM Journal of Matrix Analysis, volume 23, pages
</span><span class="comment">*</span><span class="comment">       929--947, 2002.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">       K. Braman, R. Byers and R. Mathias, The Multi-Shift QR
</span><span class="comment">*</span><span class="comment">       Algorithm Part II: Aggressive Early Deflation, SIAM Journal
</span><span class="comment">*</span><span class="comment">       of Matrix Analysis, volume 23, pages 948--973, 2002.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ================================================================
</span><span class="comment">*</span><span class="comment">     .. Parameters ..
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ==== Matrices of order NTINY or smaller must be processed by
</span><span class="comment">*</span><span class="comment">     .    <a name="ZLAHQR.216"></a><a href="zlahqr.f.html#ZLAHQR.1">ZLAHQR</a> because of insufficient subdiagonal scratch space.
</span><span class="comment">*</span><span class="comment">     .    (This is a hard limit.) ====
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ==== NL allocates some local workspace to help small matrices
</span><span class="comment">*</span><span class="comment">     .    through a rare <a name="ZLAHQR.220"></a><a href="zlahqr.f.html#ZLAHQR.1">ZLAHQR</a> failure.  NL .GT. NTINY = 11 is
</span><span class="comment">*</span><span class="comment">     .    required and NL .LE. NMIN = <a name="ILAENV.221"></a><a href="hfy-index.html#ILAENV">ILAENV</a>(ISPEC=1,...) is recom-
</span><span class="comment">*</span><span class="comment">     .    mended.  (The default value of NMIN is 75.)  Using NL = 49
</span><span class="comment">*</span><span class="comment">     .    allows up to six simultaneous shifts and a 16-by-16
</span><span class="comment">*</span><span class="comment">     .    deflation window.  ====
</span><span class="comment">*</span><span class="comment">
</span>      INTEGER            NTINY
      PARAMETER          ( NTINY = 11 )
      INTEGER            NL
      PARAMETER          ( NL = 49 )
      COMPLEX*16         ZERO, ONE
      PARAMETER          ( ZERO = ( 0.0d0, 0.0d0 ),
     $                   ONE = ( 1.0d0, 0.0d0 ) )
      DOUBLE PRECISION   RZERO
      PARAMETER          ( RZERO = 0.0d0 )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Arrays ..
</span>      COMPLEX*16         HL( NL, NL ), WORKL( NL )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Scalars ..
</span>      INTEGER            KBOT, NMIN
      LOGICAL            INITZ, LQUERY, WANTT, WANTZ
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Functions ..
</span>      INTEGER            <a name="ILAENV.244"></a><a href="hfy-index.html#ILAENV">ILAENV</a>
      LOGICAL            <a name="LSAME.245"></a><a href="lsame.f.html#LSAME.1">LSAME</a>
      EXTERNAL           <a name="ILAENV.246"></a><a href="hfy-index.html#ILAENV">ILAENV</a>, <a name="LSAME.246"></a><a href="lsame.f.html#LSAME.1">LSAME</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           <a name="XERBLA.249"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>, ZCOPY, <a name="ZLACPY.249"></a><a href="zlacpy.f.html#ZLACPY.1">ZLACPY</a>, <a name="ZLAHQR.249"></a><a href="zlahqr.f.html#ZLAHQR.1">ZLAHQR</a>, <a name="ZLAQR0.249"></a><a href="zlaqr0.f.html#ZLAQR0.1">ZLAQR0</a>, <a name="ZLASET.249"></a><a href="zlaset.f.html#ZLASET.1">ZLASET</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          DBLE, DCMPLX, MAX, MIN
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Executable Statements ..
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ==== Decode and check the input parameters. ====
</span><span class="comment">*</span><span class="comment">
</span>      WANTT = <a name="LSAME.258"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( JOB, <span class="string">'S'</span> )
      INITZ = <a name="LSAME.259"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( COMPZ, <span class="string">'I'</span> )
      WANTZ = INITZ .OR. <a name="LSAME.260"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( COMPZ, <span class="string">'V'</span> )
      WORK( 1 ) = DCMPLX( DBLE( MAX( 1, N ) ), RZERO )