cgebd2.f.html

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

HTML
275
字号
</span><span class="comment">*</span><span class="comment">  where d and e denote diagonal and off-diagonal elements of B, vi
</span><span class="comment">*</span><span class="comment">  denotes an element of the vector defining H(i), and ui an element of
</span><span class="comment">*</span><span class="comment">  the vector defining G(i).
</span><span class="comment">*</span><span class="comment">
</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>      COMPLEX            ZERO, ONE
      PARAMETER          ( ZERO = ( 0.0E+0, 0.0E+0 ),
     $                   ONE = ( 1.0E+0, 0.0E+0 ) )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Scalars ..
</span>      INTEGER            I
      COMPLEX            ALPHA
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           <a name="CLACGV.136"></a><a href="clacgv.f.html#CLACGV.1">CLACGV</a>, <a name="CLARF.136"></a><a href="clarf.f.html#CLARF.1">CLARF</a>, <a name="CLARFG.136"></a><a href="clarfg.f.html#CLARFG.1">CLARFG</a>, <a name="XERBLA.136"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          CONJG, 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">     Test the input parameters
</span><span class="comment">*</span><span class="comment">
</span>      INFO = 0
      IF( M.LT.0 ) THEN
         INFO = -1
      ELSE IF( N.LT.0 ) THEN
         INFO = -2
      ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
         INFO = -4
      END IF
      IF( INFO.LT.0 ) THEN
         CALL <a name="XERBLA.154"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="CGEBD2.154"></a><a href="cgebd2.f.html#CGEBD2.1">CGEBD2</a>'</span>, -INFO )
         RETURN
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( M.GE.N ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Reduce to upper bidiagonal form
</span><span class="comment">*</span><span class="comment">
</span>         DO 10 I = 1, N
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Generate elementary reflector H(i) to annihilate A(i+1:m,i)
</span><span class="comment">*</span><span class="comment">
</span>            ALPHA = A( I, I )
            CALL <a name="CLARFG.167"></a><a href="clarfg.f.html#CLARFG.1">CLARFG</a>( M-I+1, ALPHA, A( MIN( I+1, M ), I ), 1,
     $                   TAUQ( I ) )
            D( I ) = ALPHA
            A( I, I ) = ONE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Apply H(i)' to A(i:m,i+1:n) from the left
</span><span class="comment">*</span><span class="comment">
</span>            IF( I.LT.N )
     $         CALL <a name="CLARF.175"></a><a href="clarf.f.html#CLARF.1">CLARF</a>( <span class="string">'Left'</span>, M-I+1, N-I, A( I, I ), 1,
     $                     CONJG( TAUQ( I ) ), A( I, I+1 ), LDA, WORK )
            A( I, I ) = D( I )
<span class="comment">*</span><span class="comment">
</span>            IF( I.LT.N ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Generate elementary reflector G(i) to annihilate
</span><span class="comment">*</span><span class="comment">              A(i,i+2:n)
</span><span class="comment">*</span><span class="comment">
</span>               CALL <a name="CLACGV.184"></a><a href="clacgv.f.html#CLACGV.1">CLACGV</a>( N-I, A( I, I+1 ), LDA )
               ALPHA = A( I, I+1 )
               CALL <a name="CLARFG.186"></a><a href="clarfg.f.html#CLARFG.1">CLARFG</a>( N-I, ALPHA, A( I, MIN( I+2, N ) ),
     $                      LDA, TAUP( I ) )
               E( I ) = ALPHA
               A( I, I+1 ) = ONE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Apply G(i) to A(i+1:m,i+1:n) from the right
</span><span class="comment">*</span><span class="comment">
</span>               CALL <a name="CLARF.193"></a><a href="clarf.f.html#CLARF.1">CLARF</a>( <span class="string">'Right'</span>, M-I, N-I, A( I, I+1 ), LDA,
     $                     TAUP( I ), A( I+1, I+1 ), LDA, WORK )
               CALL <a name="CLACGV.195"></a><a href="clacgv.f.html#CLACGV.1">CLACGV</a>( N-I, A( I, I+1 ), LDA )
               A( I, I+1 ) = E( I )
            ELSE
               TAUP( I ) = ZERO
            END IF
   10    CONTINUE
      ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Reduce to lower bidiagonal form
</span><span class="comment">*</span><span class="comment">
</span>         DO 20 I = 1, M
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Generate elementary reflector G(i) to annihilate A(i,i+1:n)
</span><span class="comment">*</span><span class="comment">
</span>            CALL <a name="CLACGV.209"></a><a href="clacgv.f.html#CLACGV.1">CLACGV</a>( N-I+1, A( I, I ), LDA )
            ALPHA = A( I, I )
            CALL <a name="CLARFG.211"></a><a href="clarfg.f.html#CLARFG.1">CLARFG</a>( N-I+1, ALPHA, A( I, MIN( I+1, N ) ), LDA,
     $                   TAUP( I ) )
            D( I ) = ALPHA
            A( I, I ) = ONE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Apply G(i) to A(i+1:m,i:n) from the right
</span><span class="comment">*</span><span class="comment">
</span>            IF( I.LT.M )
     $         CALL <a name="CLARF.219"></a><a href="clarf.f.html#CLARF.1">CLARF</a>( <span class="string">'Right'</span>, M-I, N-I+1, A( I, I ), LDA,
     $                     TAUP( I ), A( I+1, I ), LDA, WORK )
            CALL <a name="CLACGV.221"></a><a href="clacgv.f.html#CLACGV.1">CLACGV</a>( N-I+1, A( I, I ), LDA )
            A( I, I ) = D( I )
<span class="comment">*</span><span class="comment">
</span>            IF( I.LT.M ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Generate elementary reflector H(i) to annihilate
</span><span class="comment">*</span><span class="comment">              A(i+2:m,i)
</span><span class="comment">*</span><span class="comment">
</span>               ALPHA = A( I+1, I )
               CALL <a name="CLARFG.230"></a><a href="clarfg.f.html#CLARFG.1">CLARFG</a>( M-I, ALPHA, A( MIN( I+2, M ), I ), 1,
     $                      TAUQ( I ) )
               E( I ) = ALPHA
               A( I+1, I ) = ONE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Apply H(i)' to A(i+1:m,i+1:n) from the left
</span><span class="comment">*</span><span class="comment">
</span>               CALL <a name="CLARF.237"></a><a href="clarf.f.html#CLARF.1">CLARF</a>( <span class="string">'Left'</span>, M-I, N-I, A( I+1, I ), 1,
     $                     CONJG( TAUQ( I ) ), A( I+1, I+1 ), LDA,
     $                     WORK )
               A( I+1, I ) = E( I )
            ELSE
               TAUQ( I ) = ZERO
            END IF
   20    CONTINUE
      END IF
      RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     End of <a name="CGEBD2.248"></a><a href="cgebd2.f.html#CGEBD2.1">CGEBD2</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>

⌨️ 快捷键说明

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