dhseqr.f.html

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

HTML
432
字号
      WORK( 1 ) = DBLE( MAX( 1, N ) )
      LQUERY = LWORK.EQ.-1
<span class="comment">*</span><span class="comment">
</span>      INFO = 0
      IF( .NOT.<a name="LSAME.274"></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.276"></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 = -11
      ELSE IF( LWORK.LT.MAX( 1, N ) .AND. .NOT.LQUERY ) THEN
         INFO = -13
      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.296"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="DHSEQR.296"></a><a href="dhseqr.f.html#DHSEQR.1">DHSEQR</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="DLAQR0.309"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a>( WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI, 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 ) = MAX( DBLE( MAX( 1, N ) ), WORK( 1 ) )
         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="DGEBAL.318"></a><a href="dgebal.f.html#DGEBAL.1">DGEBAL</a> ====
</span><span class="comment">*</span><span class="comment">
</span>         DO 10 I = 1, ILO - 1
            WR( I ) = H( I, I )
            WI( I ) = ZERO
   10    CONTINUE
         DO 20 I = IHI + 1, N
            WR( I ) = H( I, I )
            WI( I ) = ZERO
   20    CONTINUE
<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="DLASET.332"></a><a href="dlaset.f.html#DLASET.1">DLASET</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
            WR( ILO ) = H( ILO, ILO )
            WI( ILO ) = ZERO
            RETURN
         END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        ==== <a name="DLAHQR.342"></a><a href="dlahqr.f.html#DLAHQR.1">DLAHQR</a>/<a name="DLAQR0.342"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a> crossover point ====
</span><span class="comment">*</span><span class="comment">
</span>         NMIN = <a name="ILAENV.344"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 12, <span class="string">'<a name="DHSEQR.344"></a><a href="dhseqr.f.html#DHSEQR.1">DHSEQR</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="DLAQR0.348"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a> for big matrices; <a name="DLAHQR.348"></a><a href="dlahqr.f.html#DLAHQR.1">DLAHQR</a> for small ones ====
</span><span class="comment">*</span><span class="comment">
</span>         IF( N.GT.NMIN ) THEN
            CALL <a name="DLAQR0.351"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a>( WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI, 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="DLAHQR.357"></a><a href="dlahqr.f.html#DLAHQR.1">DLAHQR</a>( WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI, 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="DLAHQR.362"></a><a href="dlahqr.f.html#DLAHQR.1">DLAHQR</a> failure!  <a name="DLAQR0.362"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a> sometimes succeeds
</span><span class="comment">*</span><span class="comment">              .    when <a name="DLAHQR.363"></a><a href="dlahqr.f.html#DLAHQR.1">DLAHQR</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="DLAQR0.370"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a> directly. ====
</span><span class="comment">*</span><span class="comment">
</span>                  CALL <a name="DLAQR0.372"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a>( WANTT, WANTZ, N, ILO, KBOT, H, LDH, WR,
     $                         WI, 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="DLAQR0.378"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</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="DLAQR0.380"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a>. ====
</span><span class="comment">*</span><span class="comment">
</span>                  CALL <a name="DLACPY.382"></a><a href="dlacpy.f.html#DLACPY.1">DLACPY</a>( <span class="string">'A'</span>, N, N, H, LDH, HL, NL )
                  HL( N+1, N ) = ZERO
                  CALL <a name="DLASET.384"></a><a href="dlaset.f.html#DLASET.1">DLASET</a>( <span class="string">'A'</span>, NL, NL-N, ZERO, ZERO, HL( 1, N+1 ),
     $                         NL )
                  CALL <a name="DLAQR0.386"></a><a href="dlaqr0.f.html#DLAQR0.1">DLAQR0</a>( WANTT, WANTZ, NL, ILO, KBOT, HL, NL, WR,
     $                         WI, ILO, IHI, Z, LDZ, WORKL, NL, INFO )
                  IF( WANTT .OR. INFO.NE.0 )
     $               CALL <a name="DLACPY.389"></a><a href="dlacpy.f.html#DLACPY.1">DLACPY</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="DLASET.397"></a><a href="dlaset.f.html#DLASET.1">DLASET</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 ) = MAX( DBLE( MAX( 1, N ) ), WORK( 1 ) )
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ==== End of <a name="DHSEQR.405"></a><a href="dhseqr.f.html#DHSEQR.1">DHSEQR</a> ====
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>

⌨️ 快捷键说明

复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?