ztrsyl.f.html

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

<span class="comment">*</span><span class="comment">
</span>               SUML = ZDOTU( M-K, A( K, MIN( K+1, M ) ), LDA,
     $                C( MIN( K+1, M ), L ), 1 )
               SUMR = ZDOTU( L-1, C( K, 1 ), LDC, B( 1, L ), 1 )
               VEC = C( K, L ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span>               SCALOC = ONE
               A11 = A( K, K ) + SGN*B( L, L )
               DA11 = ABS( DBLE( A11 ) ) + ABS( DIMAG( A11 ) )
               IF( DA11.LE.SMIN ) THEN
                  A11 = SMIN
                  DA11 = SMIN
                  INFO = 1
               END IF
               DB = ABS( DBLE( VEC ) ) + ABS( DIMAG( VEC ) )
               IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
                  IF( DB.GT.BIGNUM*DA11 )
     $               SCALOC = ONE / DB
               END IF
               X11 = <a name="ZLADIV.196"></a><a href="zladiv.f.html#ZLADIV.1">ZLADIV</a>( VEC*DCMPLX( SCALOC ), A11 )
<span class="comment">*</span><span class="comment">
</span>               IF( SCALOC.NE.ONE ) THEN
                  DO 10 J = 1, N
                     CALL ZDSCAL( M, SCALOC, C( 1, J ), 1 )
   10             CONTINUE
                  SCALE = SCALE*SCALOC
               END IF
               C( K, L ) = X11
<span class="comment">*</span><span class="comment">
</span>   20       CONTINUE
   30    CONTINUE
<span class="comment">*</span><span class="comment">
</span>      ELSE IF( .NOT.NOTRNA .AND. NOTRNB ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Solve    A' *X + ISGN*X*B = scale*C.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        The (K,L)th block of X is determined starting from
</span><span class="comment">*</span><span class="comment">        upper-left corner column by column by
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">            A'(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Where
</span><span class="comment">*</span><span class="comment">                   K-1                         L-1
</span><span class="comment">*</span><span class="comment">          R(K,L) = SUM [A'(I,K)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)]
</span><span class="comment">*</span><span class="comment">                   I=1                         J=1
</span><span class="comment">*</span><span class="comment">
</span>         DO 60 L = 1, N
            DO 50 K = 1, M
<span class="comment">*</span><span class="comment">
</span>               SUML = ZDOTC( K-1, A( 1, K ), 1, C( 1, L ), 1 )
               SUMR = ZDOTU( L-1, C( K, 1 ), LDC, B( 1, L ), 1 )
               VEC = C( K, L ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span>               SCALOC = ONE
               A11 = DCONJG( A( K, K ) ) + SGN*B( L, L )
               DA11 = ABS( DBLE( A11 ) ) + ABS( DIMAG( A11 ) )
               IF( DA11.LE.SMIN ) THEN
                  A11 = SMIN
                  DA11 = SMIN
                  INFO = 1
               END IF
               DB = ABS( DBLE( VEC ) ) + ABS( DIMAG( VEC ) )
               IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
                  IF( DB.GT.BIGNUM*DA11 )
     $               SCALOC = ONE / DB
               END IF
<span class="comment">*</span><span class="comment">
</span>               X11 = <a name="ZLADIV.244"></a><a href="zladiv.f.html#ZLADIV.1">ZLADIV</a>( VEC*DCMPLX( SCALOC ), A11 )
<span class="comment">*</span><span class="comment">
</span>               IF( SCALOC.NE.ONE ) THEN
                  DO 40 J = 1, N
                     CALL ZDSCAL( M, SCALOC, C( 1, J ), 1 )
   40             CONTINUE
                  SCALE = SCALE*SCALOC
               END IF
               C( K, L ) = X11
<span class="comment">*</span><span class="comment">
</span>   50       CONTINUE
   60    CONTINUE
<span class="comment">*</span><span class="comment">
</span>      ELSE IF( .NOT.NOTRNA .AND. .NOT.NOTRNB ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Solve    A'*X + ISGN*X*B' = C.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        The (K,L)th block of X is determined starting from
</span><span class="comment">*</span><span class="comment">        upper-right corner column by column by
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">            A'(K,K)*X(K,L) + ISGN*X(K,L)*B'(L,L) = C(K,L) - R(K,L)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Where
</span><span class="comment">*</span><span class="comment">                    K-1
</span><span class="comment">*</span><span class="comment">           R(K,L) = SUM [A'(I,K)*X(I,L)] +
</span><span class="comment">*</span><span class="comment">                    I=1
</span><span class="comment">*</span><span class="comment">                           N
</span><span class="comment">*</span><span class="comment">                     ISGN*SUM [X(K,J)*B'(L,J)].
</span><span class="comment">*</span><span class="comment">                          J=L+1
</span><span class="comment">*</span><span class="comment">
</span>         DO 90 L = N, 1, -1
            DO 80 K = 1, M
<span class="comment">*</span><span class="comment">
</span>               SUML = ZDOTC( K-1, A( 1, K ), 1, C( 1, L ), 1 )
               SUMR = ZDOTC( N-L, C( K, MIN( L+1, N ) ), LDC,
     $                B( L, MIN( L+1, N ) ), LDB )
               VEC = C( K, L ) - ( SUML+SGN*DCONJG( SUMR ) )
<span class="comment">*</span><span class="comment">
</span>               SCALOC = ONE
               A11 = DCONJG( A( K, K )+SGN*B( L, L ) )
               DA11 = ABS( DBLE( A11 ) ) + ABS( DIMAG( A11 ) )
               IF( DA11.LE.SMIN ) THEN
                  A11 = SMIN
                  DA11 = SMIN
                  INFO = 1
               END IF
               DB = ABS( DBLE( VEC ) ) + ABS( DIMAG( VEC ) )
               IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
                  IF( DB.GT.BIGNUM*DA11 )
     $               SCALOC = ONE / DB
               END IF
<span class="comment">*</span><span class="comment">
</span>               X11 = <a name="ZLADIV.296"></a><a href="zladiv.f.html#ZLADIV.1">ZLADIV</a>( VEC*DCMPLX( SCALOC ), A11 )
<span class="comment">*</span><span class="comment">
</span>               IF( SCALOC.NE.ONE ) THEN
                  DO 70 J = 1, N
                     CALL ZDSCAL( M, SCALOC, C( 1, J ), 1 )
   70             CONTINUE
                  SCALE = SCALE*SCALOC
               END IF
               C( K, L ) = X11
<span class="comment">*</span><span class="comment">
</span>   80       CONTINUE
   90    CONTINUE
<span class="comment">*</span><span class="comment">
</span>      ELSE IF( NOTRNA .AND. .NOT.NOTRNB ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Solve    A*X + ISGN*X*B' = C.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        The (K,L)th block of X is determined starting from
</span><span class="comment">*</span><span class="comment">        bottom-left corner column by column by
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           A(K,K)*X(K,L) + ISGN*X(K,L)*B'(L,L) = C(K,L) - R(K,L)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Where
</span><span class="comment">*</span><span class="comment">                    M                          N
</span><span class="comment">*</span><span class="comment">          R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B'(L,J)]
</span><span class="comment">*</span><span class="comment">                  I=K+1                      J=L+1
</span><span class="comment">*</span><span class="comment">
</span>         DO 120 L = N, 1, -1
            DO 110 K = M, 1, -1
<span class="comment">*</span><span class="comment">
</span>               SUML = ZDOTU( M-K, A( K, MIN( K+1, M ) ), LDA,
     $                C( MIN( K+1, M ), L ), 1 )
               SUMR = ZDOTC( N-L, C( K, MIN( L+1, N ) ), LDC,
     $                B( L, MIN( L+1, N ) ), LDB )
               VEC = C( K, L ) - ( SUML+SGN*DCONJG( SUMR ) )
<span class="comment">*</span><span class="comment">
</span>               SCALOC = ONE
               A11 = A( K, K ) + SGN*DCONJG( B( L, L ) )
               DA11 = ABS( DBLE( A11 ) ) + ABS( DIMAG( A11 ) )
               IF( DA11.LE.SMIN ) THEN
                  A11 = SMIN
                  DA11 = SMIN
                  INFO = 1
               END IF
               DB = ABS( DBLE( VEC ) ) + ABS( DIMAG( VEC ) )
               IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
                  IF( DB.GT.BIGNUM*DA11 )
     $               SCALOC = ONE / DB
               END IF
<span class="comment">*</span><span class="comment">
</span>               X11 = <a name="ZLADIV.346"></a><a href="zladiv.f.html#ZLADIV.1">ZLADIV</a>( VEC*DCMPLX( SCALOC ), A11 )
<span class="comment">*</span><span class="comment">
</span>               IF( SCALOC.NE.ONE ) THEN
                  DO 100 J = 1, N
                     CALL ZDSCAL( M, SCALOC, C( 1, J ), 1 )
  100             CONTINUE
                  SCALE = SCALE*SCALOC
               END IF
               C( K, L ) = X11
<span class="comment">*</span><span class="comment">
</span>  110       CONTINUE
  120    CONTINUE
<span class="comment">*</span><span class="comment">
</span>      END IF
<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="ZTRSYL.363"></a><a href="ztrsyl.f.html#ZTRSYL.1">ZTRSYL</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>