slasdq.f.html

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

</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Test the input parameters.
</span><span class="comment">*</span><span class="comment">
</span>      INFO = 0
      IUPLO = 0
      IF( <a name="LSAME.158"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( UPLO, <span class="string">'U'</span> ) )
     $   IUPLO = 1
      IF( <a name="LSAME.160"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( UPLO, <span class="string">'L'</span> ) )
     $   IUPLO = 2
      IF( IUPLO.EQ.0 ) THEN
         INFO = -1
      ELSE IF( ( SQRE.LT.0 ) .OR. ( SQRE.GT.1 ) ) THEN
         INFO = -2
      ELSE IF( N.LT.0 ) THEN
         INFO = -3
      ELSE IF( NCVT.LT.0 ) THEN
         INFO = -4
      ELSE IF( NRU.LT.0 ) THEN
         INFO = -5
      ELSE IF( NCC.LT.0 ) THEN
         INFO = -6
      ELSE IF( ( NCVT.EQ.0 .AND. LDVT.LT.1 ) .OR.
     $         ( NCVT.GT.0 .AND. LDVT.LT.MAX( 1, N ) ) ) THEN
         INFO = -10
      ELSE IF( LDU.LT.MAX( 1, NRU ) ) THEN
         INFO = -12
      ELSE IF( ( NCC.EQ.0 .AND. LDC.LT.1 ) .OR.
     $         ( NCC.GT.0 .AND. LDC.LT.MAX( 1, N ) ) ) THEN
         INFO = -14
      END IF
      IF( INFO.NE.0 ) THEN
         CALL <a name="XERBLA.184"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="SLASDQ.184"></a><a href="slasdq.f.html#SLASDQ.1">SLASDQ</a>'</span>, -INFO )
         RETURN
      END IF
      IF( N.EQ.0 )
     $   RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ROTATE is true if any singular vectors desired, false otherwise
</span><span class="comment">*</span><span class="comment">
</span>      ROTATE = ( NCVT.GT.0 ) .OR. ( NRU.GT.0 ) .OR. ( NCC.GT.0 )
      NP1 = N + 1
      SQRE1 = SQRE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     If matrix non-square upper bidiagonal, rotate to be lower
</span><span class="comment">*</span><span class="comment">     bidiagonal.  The rotations are on the right.
</span><span class="comment">*</span><span class="comment">
</span>      IF( ( IUPLO.EQ.1 ) .AND. ( SQRE1.EQ.1 ) ) THEN
         DO 10 I = 1, N - 1
            CALL <a name="SLARTG.201"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( D( I ), E( I ), CS, SN, R )
            D( I ) = R
            E( I ) = SN*D( I+1 )
            D( I+1 ) = CS*D( I+1 )
            IF( ROTATE ) THEN
               WORK( I ) = CS
               WORK( N+I ) = SN
            END IF
   10    CONTINUE
         CALL <a name="SLARTG.210"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( D( N ), E( N ), CS, SN, R )
         D( N ) = R
         E( N ) = ZERO
         IF( ROTATE ) THEN
            WORK( N ) = CS
            WORK( N+N ) = SN
         END IF
         IUPLO = 2
         SQRE1 = 0
<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.223"></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>, NP1, NCVT, WORK( 1 ),
     $                  WORK( NP1 ), VT, LDVT )
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     If matrix lower bidiagonal, rotate to be upper bidiagonal
</span><span class="comment">*</span><span class="comment">     by applying Givens rotations on the left.
</span><span class="comment">*</span><span class="comment">
</span>      IF( IUPLO.EQ.2 ) THEN
         DO 20 I = 1, N - 1
            CALL <a name="SLARTG.232"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( D( I ), E( I ), CS, SN, R )
            D( I ) = R
            E( I ) = SN*D( I+1 )
            D( I+1 ) = CS*D( I+1 )
            IF( ROTATE ) THEN
               WORK( I ) = CS
               WORK( N+I ) = SN
            END IF
   20    CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        If matrix (N+1)-by-N lower bidiagonal, one additional
</span><span class="comment">*</span><span class="comment">        rotation is needed.
</span><span class="comment">*</span><span class="comment">
</span>         IF( SQRE1.EQ.1 ) THEN
            CALL <a name="SLARTG.246"></a><a href="slartg.f.html#SLARTG.1">SLARTG</a>( D( N ), E( N ), CS, SN, R )
            D( N ) = R
            IF( ROTATE ) THEN
               WORK( N ) = CS
               WORK( N+N ) = SN
            END IF
         END IF
<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( NRU.GT.0 ) THEN
            IF( SQRE1.EQ.0 ) THEN
               CALL <a name="SLASR.258"></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, N, WORK( 1 ),
     $                     WORK( NP1 ), U, LDU )
            ELSE
               CALL <a name="SLASR.261"></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, NP1, WORK( 1 ),
     $                     WORK( NP1 ), U, LDU )
            END IF
         END IF
         IF( NCC.GT.0 ) THEN
            IF( SQRE1.EQ.0 ) THEN
               CALL <a name="SLASR.267"></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>, N, NCC, WORK( 1 ),
     $                     WORK( NP1 ), C, LDC )
            ELSE
               CALL <a name="SLASR.270"></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>, NP1, NCC, WORK( 1 ),
     $                     WORK( NP1 ), C, LDC )
            END IF
         END IF
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Call <a name="SBDSQR.276"></a><a href="sbdsqr.f.html#SBDSQR.1">SBDSQR</a> to compute the SVD of the reduced real
</span><span class="comment">*</span><span class="comment">     N-by-N upper bidiagonal matrix.
</span><span class="comment">*</span><span class="comment">
</span>      CALL <a name="SBDSQR.279"></a><a href="sbdsqr.f.html#SBDSQR.1">SBDSQR</a>( <span class="string">'U'</span>, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, LDU, C,
     $             LDC, WORK, INFO )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Sort the singular values into ascending 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 40 I = 1, N
<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 = I
         SMIN = D( I )
         DO 30 J = I + 1, N
            IF( D( J ).LT.SMIN ) THEN
               ISUB = J
               SMIN = D( J )
            END IF
   30    CONTINUE
         IF( ISUB.NE.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( I )
            D( I ) = SMIN
            IF( NCVT.GT.0 )
     $         CALL SSWAP( NCVT, VT( ISUB, 1 ), LDVT, VT( I, 1 ), LDVT )
            IF( NRU.GT.0 )
     $         CALL SSWAP( NRU, U( 1, ISUB ), 1, U( 1, I ), 1 )
            IF( NCC.GT.0 )
     $         CALL SSWAP( NCC, C( ISUB, 1 ), LDC, C( I, 1 ), LDC )
         END IF
   40 CONTINUE
<span class="comment">*</span><span class="comment">
</span>      RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     End of <a name="SLASDQ.314"></a><a href="slasdq.f.html#SLASDQ.1">SLASDQ</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>