chseqr.f.html

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

      LQUERY = LWORK.EQ.-1
<span class="comment">*</span><span class="comment">
</span>      INFO = 0
      IF( .NOT.<a name="LSAME.265"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( JOB, <span class="string">'E'</span> ) .AND. .NOT.WANTT ) THEN
         INFO = -1
      ELSE IF( .NOT.<a name="LSAME.267"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( COMPZ, <span class="string">'N'</span> ) .AND. .NOT.WANTZ ) THEN
         INFO = -2
      ELSE IF( N.LT.0 ) THEN
         INFO = -3
      ELSE IF( ILO.LT.1 .OR. ILO.GT.MAX( 1, N ) ) THEN
         INFO = -4
      ELSE IF( IHI.LT.MIN( ILO, N ) .OR. IHI.GT.N ) THEN
         INFO = -5
      ELSE IF( LDH.LT.MAX( 1, N ) ) THEN
         INFO = -7
      ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.MAX( 1, N ) ) ) THEN
         INFO = -10
      ELSE IF( LWORK.LT.MAX( 1, N ) .AND. .NOT.LQUERY ) THEN
         INFO = -12
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( INFO.NE.0 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== Quick return in case of invalid argument. ====
</span><span class="comment">*</span><span class="comment">
</span>         CALL <a name="XERBLA.287"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="CHSEQR.287"></a><a href="chseqr.f.html#CHSEQR.1">CHSEQR</a>'</span>, -INFO )
         RETURN
<span class="comment">*</span><span class="comment">
</span>      ELSE IF( N.EQ.0 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== Quick return in case N = 0; nothing to do. ====
</span><span class="comment">*</span><span class="comment">
</span>         RETURN
<span class="comment">*</span><span class="comment">
</span>      ELSE IF( LQUERY ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== Quick return in case of a workspace query ====
</span><span class="comment">*</span><span class="comment">
</span>         CALL <a name="CLAQR0.300"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>( WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILO, IHI, Z,
     $                LDZ, WORK, LWORK, INFO )
<span class="comment">*</span><span class="comment">        ==== Ensure reported workspace size is backward-compatible with
</span><span class="comment">*</span><span class="comment">        .    previous LAPACK versions. ====
</span>         WORK( 1 ) = CMPLX( MAX( REAL( WORK( 1 ) ), REAL( MAX( 1,
     $               N ) ) ), RZERO )
         RETURN
<span class="comment">*</span><span class="comment">
</span>      ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== copy eigenvalues isolated by <a name="CGEBAL.310"></a><a href="cgebal.f.html#CGEBAL.1">CGEBAL</a> ====
</span><span class="comment">*</span><span class="comment">
</span>         IF( ILO.GT.1 )
     $      CALL CCOPY( ILO-1, H, LDH+1, W, 1 )
         IF( IHI.LT.N )
     $      CALL CCOPY( N-IHI, H( IHI+1, IHI+1 ), LDH+1, W( IHI+1 ), 1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== Initialize Z, if requested ====
</span><span class="comment">*</span><span class="comment">
</span>         IF( INITZ )
     $      CALL <a name="CLASET.320"></a><a href="claset.f.html#CLASET.1">CLASET</a>( <span class="string">'A'</span>, N, N, ZERO, ONE, Z, LDZ )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== Quick return if possible ====
</span><span class="comment">*</span><span class="comment">
</span>         IF( ILO.EQ.IHI ) THEN
            W( ILO ) = H( ILO, ILO )
            RETURN
         END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== <a name="CLAHQR.329"></a><a href="clahqr.f.html#CLAHQR.1">CLAHQR</a>/<a name="CLAQR0.329"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a> crossover point ====
</span><span class="comment">*</span><span class="comment">
</span>         NMIN = <a name="ILAENV.331"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="CHSEQR.331"></a><a href="chseqr.f.html#CHSEQR.1">CHSEQR</a>'</span>, JOB( : 1 ) // COMPZ( : 1 ), N, ILO,
     $          IHI, LWORK )
         NMIN = MAX( NTINY, NMIN )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== <a name="CLAQR0.335"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a> for big matrices; <a name="CLAHQR.335"></a><a href="clahqr.f.html#CLAHQR.1">CLAHQR</a> for small ones ====
</span><span class="comment">*</span><span class="comment">
</span>         IF( N.GT.NMIN ) THEN
            CALL <a name="CLAQR0.338"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>( WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILO, IHI,
     $                   Z, LDZ, WORK, LWORK, INFO )
         ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           ==== Small matrix ====
</span><span class="comment">*</span><span class="comment">
</span>            CALL <a name="CLAHQR.344"></a><a href="clahqr.f.html#CLAHQR.1">CLAHQR</a>( WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILO, IHI,
     $                   Z, LDZ, INFO )
<span class="comment">*</span><span class="comment">
</span>            IF( INFO.GT.0 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              ==== A rare <a name="CLAHQR.349"></a><a href="clahqr.f.html#CLAHQR.1">CLAHQR</a> failure!  <a name="CLAQR0.349"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a> sometimes succeeds
</span><span class="comment">*</span><span class="comment">              .    when <a name="CLAHQR.350"></a><a href="clahqr.f.html#CLAHQR.1">CLAHQR</a> fails. ====
</span><span class="comment">*</span><span class="comment">
</span>               KBOT = INFO
<span class="comment">*</span><span class="comment">
</span>               IF( N.GE.NL ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 ==== Larger matrices have enough subdiagonal scratch
</span><span class="comment">*</span><span class="comment">                 .    space to call <a name="CLAQR0.357"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a> directly. ====
</span><span class="comment">*</span><span class="comment">
</span>                  CALL <a name="CLAQR0.359"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>( WANTT, WANTZ, N, ILO, KBOT, H, LDH, W,
     $                         ILO, IHI, Z, LDZ, WORK, LWORK, INFO )
<span class="comment">*</span><span class="comment">
</span>               ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 ==== Tiny matrices don't have enough subdiagonal
</span><span class="comment">*</span><span class="comment">                 .    scratch space to benefit from <a name="CLAQR0.365"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>.  Hence,
</span><span class="comment">*</span><span class="comment">                 .    tiny matrices must be copied into a larger
</span><span class="comment">*</span><span class="comment">                 .    array before calling <a name="CLAQR0.367"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>. ====
</span><span class="comment">*</span><span class="comment">
</span>                  CALL <a name="CLACPY.369"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'A'</span>, N, N, H, LDH, HL, NL )
                  HL( N+1, N ) = ZERO
                  CALL <a name="CLASET.371"></a><a href="claset.f.html#CLASET.1">CLASET</a>( <span class="string">'A'</span>, NL, NL-N, ZERO, ZERO, HL( 1, N+1 ),
     $                         NL )
                  CALL <a name="CLAQR0.373"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>( WANTT, WANTZ, NL, ILO, KBOT, HL, NL, W,
     $                         ILO, IHI, Z, LDZ, WORKL, NL, INFO )
                  IF( WANTT .OR. INFO.NE.0 )
     $               CALL <a name="CLACPY.376"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'A'</span>, N, N, HL, NL, H, LDH )
               END IF
            END IF
         END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== Clear out the trash, if necessary. ====
</span><span class="comment">*</span><span class="comment">
</span>         IF( ( WANTT .OR. INFO.NE.0 ) .AND. N.GT.2 )
     $      CALL <a name="CLASET.384"></a><a href="claset.f.html#CLASET.1">CLASET</a>( <span class="string">'L'</span>, N-2, N-2, ZERO, ZERO, H( 3, 1 ), LDH )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== Ensure reported workspace size is backward-compatible with
</span><span class="comment">*</span><span class="comment">        .    previous LAPACK versions. ====
</span><span class="comment">*</span><span class="comment">
</span>         WORK( 1 ) = CMPLX( MAX( REAL( MAX( 1, N ) ),
     $               REAL( WORK( 1 ) ) ), RZERO )
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ==== End of <a name="CHSEQR.393"></a><a href="chseqr.f.html#CHSEQR.1">CHSEQR</a> ====
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>