sbdsqr.f.html

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

HTML
767
字号
            E( LL ) = H*OLDSN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Update singular vectors
</span><span class="comment">*</span><span class="comment">
</span>            IF( NCVT.GT.0 )
     $         CALL <a name="SLASR.563"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'L'</span>, <span class="string">'V'</span>, <span class="string">'B'</span>, M-LL+1, NCVT, WORK( NM12+1 ),
     $                     WORK( NM13+1 ), VT( LL, 1 ), LDVT )
            IF( NRU.GT.0 )
     $         CALL <a name="SLASR.566"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'R'</span>, <span class="string">'V'</span>, <span class="string">'B'</span>, NRU, M-LL+1, WORK( 1 ),
     $                     WORK( N ), U( 1, LL ), LDU )
            IF( NCC.GT.0 )
     $         CALL <a name="SLASR.569"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'L'</span>, <span class="string">'V'</span>, <span class="string">'B'</span>, M-LL+1, NCC, WORK( 1 ),
     $                     WORK( N ), C( LL, 1 ), LDC )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Test convergence
</span><span class="comment">*</span><span class="comment">
</span>            IF( ABS( E( LL ) ).LE.THRESH )
     $         E( LL ) = ZERO
         END IF
      ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Use nonzero shift
</span><span class="comment">*</span><span class="comment">
</span>         IF( IDIR.EQ.1 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Chase bulge from top to bottom
</span><span class="comment">*</span><span class="comment">           Save cosines and sines for later singular vector updates
</span><span class="comment">*</span><span class="comment">
</span>            F = ( ABS( D( LL ) )-SHIFT )*
     $          ( SIGN( ONE, D( LL ) )+SHIFT / D( LL ) )
            G = E( LL )
            DO 140 I = LL, M - 1
               CALL <a name="SLARTG.590"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( F, G, COSR, SINR, R )
               IF( I.GT.LL )
     $            E( I-1 ) = R
               F = COSR*D( I ) + SINR*E( I )
               E( I ) = COSR*E( I ) - SINR*D( I )
               G = SINR*D( I+1 )
               D( I+1 ) = COSR*D( I+1 )
               CALL <a name="SLARTG.597"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( F, G, COSL, SINL, R )
               D( I ) = R
               F = COSL*E( I ) + SINL*D( I+1 )
               D( I+1 ) = COSL*D( I+1 ) - SINL*E( I )
               IF( I.LT.M-1 ) THEN
                  G = SINL*E( I+1 )
                  E( I+1 ) = COSL*E( I+1 )
               END IF
               WORK( I-LL+1 ) = COSR
               WORK( I-LL+1+NM1 ) = SINR
               WORK( I-LL+1+NM12 ) = COSL
               WORK( I-LL+1+NM13 ) = SINL
  140       CONTINUE
            E( M-1 ) = F
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Update singular vectors
</span><span class="comment">*</span><span class="comment">
</span>            IF( NCVT.GT.0 )
     $         CALL <a name="SLASR.615"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'L'</span>, <span class="string">'V'</span>, <span class="string">'F'</span>, M-LL+1, NCVT, WORK( 1 ),
     $                     WORK( N ), VT( LL, 1 ), LDVT )
            IF( NRU.GT.0 )
     $         CALL <a name="SLASR.618"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'R'</span>, <span class="string">'V'</span>, <span class="string">'F'</span>, NRU, M-LL+1, WORK( NM12+1 ),
     $                     WORK( NM13+1 ), U( 1, LL ), LDU )
            IF( NCC.GT.0 )
     $         CALL <a name="SLASR.621"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'L'</span>, <span class="string">'V'</span>, <span class="string">'F'</span>, M-LL+1, NCC, WORK( NM12+1 ),
     $                     WORK( NM13+1 ), C( LL, 1 ), LDC )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Test convergence
</span><span class="comment">*</span><span class="comment">
</span>            IF( ABS( E( M-1 ) ).LE.THRESH )
     $         E( M-1 ) = ZERO
<span class="comment">*</span><span class="comment">
</span>         ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Chase bulge from bottom to top
</span><span class="comment">*</span><span class="comment">           Save cosines and sines for later singular vector updates
</span><span class="comment">*</span><span class="comment">
</span>            F = ( ABS( D( M ) )-SHIFT )*( SIGN( ONE, D( M ) )+SHIFT /
     $          D( M ) )
            G = E( M-1 )
            DO 150 I = M, LL + 1, -1
               CALL <a name="SLARTG.638"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( F, G, COSR, SINR, R )
               IF( I.LT.M )
     $            E( I ) = R
               F = COSR*D( I ) + SINR*E( I-1 )
               E( I-1 ) = COSR*E( I-1 ) - SINR*D( I )
               G = SINR*D( I-1 )
               D( I-1 ) = COSR*D( I-1 )
               CALL <a name="SLARTG.645"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( F, G, COSL, SINL, R )
               D( I ) = R
               F = COSL*E( I-1 ) + SINL*D( I-1 )
               D( I-1 ) = COSL*D( I-1 ) - SINL*E( I-1 )
               IF( I.GT.LL+1 ) THEN
                  G = SINL*E( I-2 )
                  E( I-2 ) = COSL*E( I-2 )
               END IF
               WORK( I-LL ) = COSR
               WORK( I-LL+NM1 ) = -SINR
               WORK( I-LL+NM12 ) = COSL
               WORK( I-LL+NM13 ) = -SINL
  150       CONTINUE
            E( LL ) = F
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Test convergence
</span><span class="comment">*</span><span class="comment">
</span>            IF( ABS( E( LL ) ).LE.THRESH )
     $         E( LL ) = ZERO
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Update singular vectors if desired
</span><span class="comment">*</span><span class="comment">
</span>            IF( NCVT.GT.0 )
     $         CALL <a name="SLASR.668"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'L'</span>, <span class="string">'V'</span>, <span class="string">'B'</span>, M-LL+1, NCVT, WORK( NM12+1 ),
     $                     WORK( NM13+1 ), VT( LL, 1 ), LDVT )
            IF( NRU.GT.0 )
     $         CALL <a name="SLASR.671"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'R'</span>, <span class="string">'V'</span>, <span class="string">'B'</span>, NRU, M-LL+1, WORK( 1 ),
     $                     WORK( N ), U( 1, LL ), LDU )
            IF( NCC.GT.0 )
     $         CALL <a name="SLASR.674"></a><a href="slasr.f.html#SLASR.1">SLASR</a>( <span class="string">'L'</span>, <span class="string">'V'</span>, <span class="string">'B'</span>, M-LL+1, NCC, WORK( 1 ),
     $                     WORK( N ), C( LL, 1 ), LDC )
         END IF
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     QR iteration finished, go back and check convergence
</span><span class="comment">*</span><span class="comment">
</span>      GO TO 60
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     All singular values converged, so make them positive
</span><span class="comment">*</span><span class="comment">
</span>  160 CONTINUE
      DO 170 I = 1, N
         IF( D( I ).LT.ZERO ) THEN
            D( I ) = -D( I )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Change sign of singular vectors, if desired
</span><span class="comment">*</span><span class="comment">
</span>            IF( NCVT.GT.0 )
     $         CALL SSCAL( NCVT, NEGONE, VT( I, 1 ), LDVT )
         END IF
  170 CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Sort the singular values into decreasing order (insertion sort on
</span><span class="comment">*</span><span class="comment">     singular values, but only one transposition per singular vector)
</span><span class="comment">*</span><span class="comment">
</span>      DO 190 I = 1, N - 1
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Scan for smallest D(I)
</span><span class="comment">*</span><span class="comment">
</span>         ISUB = 1
         SMIN = D( 1 )
         DO 180 J = 2, N + 1 - I
            IF( D( J ).LE.SMIN ) THEN
               ISUB = J
               SMIN = D( J )
            END IF
  180    CONTINUE
         IF( ISUB.NE.N+1-I ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Swap singular values and vectors
</span><span class="comment">*</span><span class="comment">
</span>            D( ISUB ) = D( N+1-I )
            D( N+1-I ) = SMIN
            IF( NCVT.GT.0 )
     $         CALL SSWAP( NCVT, VT( ISUB, 1 ), LDVT, VT( N+1-I, 1 ),
     $                     LDVT )
            IF( NRU.GT.0 )
     $         CALL SSWAP( NRU, U( 1, ISUB ), 1, U( 1, N+1-I ), 1 )
            IF( NCC.GT.0 )
     $         CALL SSWAP( NCC, C( ISUB, 1 ), LDC, C( N+1-I, 1 ), LDC )
         END IF
  190 CONTINUE
      GO TO 220
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Maximum number of iterations exceeded, failure to converge
</span><span class="comment">*</span><span class="comment">
</span>  200 CONTINUE
      INFO = 0
      DO 210 I = 1, N - 1
         IF( E( I ).NE.ZERO )
     $      INFO = INFO + 1
  210 CONTINUE
  220 CONTINUE
      RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     End of <a name="SBDSQR.740"></a><a href="sbdsqr.f.html#SBDSQR.1">SBDSQR</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>

⌨️ 快捷键说明

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