dtrsyl.f.html
来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 938 行 · 第 1/4 页
HTML
938 行
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 )
SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
$ B( L2, MIN( L2+1, N ) ), LDB )
VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span> CALL <a name="DLASY2.693"></a><a href="dlasy2.f.html#DLASY2.1">DLASY2</a>( .TRUE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
$ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
$ 2, XNORM, IERR )
IF( IERR.NE.0 )
$ INFO = 1
<span class="comment">*</span><span class="comment">
</span> IF( SCALOC.NE.ONE ) THEN
DO 160 J = 1, N
CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
160 CONTINUE
SCALE = SCALE*SCALOC
END IF
C( K1, L1 ) = X( 1, 1 )
C( K1, L2 ) = X( 1, 2 )
C( K2, L1 ) = X( 2, 1 )
C( K2, L2 ) = X( 2, 2 )
END IF
<span class="comment">*</span><span class="comment">
</span> 170 CONTINUE
180 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' = 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"> bottom-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"> 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><span class="comment">*</span><span class="comment"> Start column loop (index = L)
</span><span class="comment">*</span><span class="comment"> L1 (L2): column index of the first (last) row of X(K,L)
</span><span class="comment">*</span><span class="comment">
</span> LNEXT = N
DO 240 L = N, 1, -1
IF( L.GT.LNEXT )
$ GO TO 240
IF( L.EQ.1 ) THEN
L1 = L
L2 = L
ELSE
IF( B( L, L-1 ).NE.ZERO ) THEN
L1 = L - 1
L2 = L
LNEXT = L - 2
ELSE
L1 = L
L2 = L
LNEXT = L - 1
END IF
END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Start row loop (index = K)
</span><span class="comment">*</span><span class="comment"> K1 (K2): row index of the first (last) row of X(K,L)
</span><span class="comment">*</span><span class="comment">
</span> KNEXT = M
DO 230 K = M, 1, -1
IF( K.GT.KNEXT )
$ GO TO 230
IF( K.EQ.1 ) THEN
K1 = K
K2 = K
ELSE
IF( A( K, K-1 ).NE.ZERO ) THEN
K1 = K - 1
K2 = K
KNEXT = K - 2
ELSE
K1 = K
K2 = K
KNEXT = K - 1
END IF
END IF
<span class="comment">*</span><span class="comment">
</span> IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
$ C( MIN( K1+1, M ), L1 ), 1 )
SUMR = DDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC,
$ B( L1, MIN( L1+1, N ) ), LDB )
VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
SCALOC = ONE
<span class="comment">*</span><span class="comment">
</span> A11 = A( K1, K1 ) + SGN*B( L1, L1 )
DA11 = ABS( A11 )
IF( DA11.LE.SMIN ) THEN
A11 = SMIN
DA11 = SMIN
INFO = 1
END IF
DB = ABS( VEC( 1, 1 ) )
IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
IF( DB.GT.BIGNUM*DA11 )
$ SCALOC = ONE / DB
END IF
X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
<span class="comment">*</span><span class="comment">
</span> IF( SCALOC.NE.ONE ) THEN
DO 190 J = 1, N
CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
190 CONTINUE
SCALE = SCALE*SCALOC
END IF
C( K1, L1 ) = X( 1, 1 )
<span class="comment">*</span><span class="comment">
</span> ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
$ C( MIN( K2+1, M ), L1 ), 1 )
SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
$ B( L1, MIN( L2+1, N ) ), LDB )
VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
$ C( MIN( K2+1, M ), L1 ), 1 )
SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
$ B( L1, MIN( L2+1, N ) ), LDB )
VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span> CALL <a name="DLALN2.816"></a><a href="dlaln2.f.html#DLALN2.1">DLALN2</a>( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ),
$ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
$ ZERO, X, 2, SCALOC, XNORM, IERR )
IF( IERR.NE.0 )
$ INFO = 1
<span class="comment">*</span><span class="comment">
</span> IF( SCALOC.NE.ONE ) THEN
DO 200 J = 1, N
CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
200 CONTINUE
SCALE = SCALE*SCALOC
END IF
C( K1, L1 ) = X( 1, 1 )
C( K2, L1 ) = X( 2, 1 )
<span class="comment">*</span><span class="comment">
</span> ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
$ C( MIN( K1+1, M ), L1 ), 1 )
SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
$ B( L1, MIN( L2+1, N ) ), LDB )
VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
$ C( MIN( K1+1, M ), L2 ), 1 )
SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
$ B( L2, MIN( L2+1, N ) ), LDB )
VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
<span class="comment">*</span><span class="comment">
</span> CALL <a name="DLALN2.845"></a><a href="dlaln2.f.html#DLALN2.1">DLALN2</a>( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ),
$ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
$ ZERO, X, 2, SCALOC, XNORM, IERR )
IF( IERR.NE.0 )
$ INFO = 1
<span class="comment">*</span><span class="comment">
</span> IF( SCALOC.NE.ONE ) THEN
DO 210 J = 1, N
CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
210 CONTINUE
SCALE = SCALE*SCALOC
END IF
C( K1, L1 ) = X( 1, 1 )
C( K1, L2 ) = X( 2, 1 )
<span class="comment">*</span><span class="comment">
</span> ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
$ C( MIN( K2+1, M ), L1 ), 1 )
SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
$ B( L1, MIN( L2+1, N ) ), LDB )
VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
$ C( MIN( K2+1, M ), L2 ), 1 )
SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
$ B( L2, MIN( L2+1, N ) ), LDB )
VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
$ C( MIN( K2+1, M ), L1 ), 1 )
SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
$ B( L1, MIN( L2+1, N ) ), LDB )
VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span> SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
$ C( MIN( K2+1, M ), L2 ), 1 )
SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
$ B( L2, MIN( L2+1, N ) ), LDB )
VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
<span class="comment">*</span><span class="comment">
</span> CALL <a name="DLASY2.886"></a><a href="dlasy2.f.html#DLASY2.1">DLASY2</a>( .FALSE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
$ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
$ 2, XNORM, IERR )
IF( IERR.NE.0 )
$ INFO = 1
<span class="comment">*</span><span class="comment">
</span> IF( SCALOC.NE.ONE ) THEN
DO 220 J = 1, N
CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
220 CONTINUE
SCALE = SCALE*SCALOC
END IF
C( K1, L1 ) = X( 1, 1 )
C( K1, L2 ) = X( 1, 2 )
C( K2, L1 ) = X( 2, 1 )
C( K2, L2 ) = X( 2, 2 )
END IF
<span class="comment">*</span><span class="comment">
</span> 230 CONTINUE
240 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="DTRSYL.911"></a><a href="dtrsyl.f.html#DTRSYL.1">DTRSYL</a>
</span><span class="comment">*</span><span class="comment">
</span> END
</pre>
</body>
</html>
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?