SUBROUTINE DTRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C, $ LDC, SCALE, INFO ) * * -- LAPACK routine (version 3.2) -- * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- * November 2006 * * .. Scalar Arguments .. CHARACTER TRANA, TRANB INTEGER INFO, ISGN, LDA, LDB, LDC, M, N DOUBLE PRECISION SCALE * .. * .. Array Arguments .. DOUBLE PRECISION A( LDA, * ), B( LDB, * ), C( LDC, * ) * .. * * Purpose * ======= * * DTRSYL solves the real Sylvester matrix equation: * * op(A)*X + X*op(B) = scale*C or * op(A)*X - X*op(B) = scale*C, * * where op(A) = A or A**T, and A and B are both upper quasi- * triangular. A is M-by-M and B is N-by-N; the right hand side C and * the solution X are M-by-N; and scale is an output scale factor, set * <= 1 to avoid overflow in X. * * A and B must be in Schur canonical form (as returned by DHSEQR), that * is, block upper triangular with 1-by-1 and 2-by-2 diagonal blocks; * each 2-by-2 diagonal block has its diagonal elements equal and its * off-diagonal elements of opposite sign. * * Arguments * ========= * * TRANA (input) CHARACTER*1 * Specifies the option op(A): * = 'N': op(A) = A (No transpose) * = 'T': op(A) = A**T (Transpose) * = 'C': op(A) = A**H (Conjugate transpose = Transpose) * * TRANB (input) CHARACTER*1 * Specifies the option op(B): * = 'N': op(B) = B (No transpose) * = 'T': op(B) = B**T (Transpose) * = 'C': op(B) = B**H (Conjugate transpose = Transpose) * * ISGN (input) INTEGER * Specifies the sign in the equation: * = +1: solve op(A)*X + X*op(B) = scale*C * = -1: solve op(A)*X - X*op(B) = scale*C * * M (input) INTEGER * The order of the matrix A, and the number of rows in the * matrices X and C. M >= 0. * * N (input) INTEGER * The order of the matrix B, and the number of columns in the * matrices X and C. N >= 0. * * A (input) DOUBLE PRECISION array, dimension (LDA,M) * The upper quasi-triangular matrix A, in Schur canonical form. * * LDA (input) INTEGER * The leading dimension of the array A. LDA >= max(1,M). * * B (input) DOUBLE PRECISION array, dimension (LDB,N) * The upper quasi-triangular matrix B, in Schur canonical form. * * LDB (input) INTEGER * The leading dimension of the array B. LDB >= max(1,N). * * C (input/output) DOUBLE PRECISION array, dimension (LDC,N) * On entry, the M-by-N right hand side matrix C. * On exit, C is overwritten by the solution matrix X. * * LDC (input) INTEGER * The leading dimension of the array C. LDC >= max(1,M) * * SCALE (output) DOUBLE PRECISION * The scale factor, scale, set <= 1 to avoid overflow in X. * * INFO (output) INTEGER * = 0: successful exit * < 0: if INFO = -i, the i-th argument had an illegal value * = 1: A and B have common or very close eigenvalues; perturbed * values were used to solve the equation (but the matrices * A and B are unchanged). * * ===================================================================== * * .. Parameters .. DOUBLE PRECISION ZERO, ONE PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) * .. * .. Local Scalars .. LOGICAL NOTRNA, NOTRNB INTEGER IERR, J, K, K1, K2, KNEXT, L, L1, L2, LNEXT DOUBLE PRECISION A11, BIGNUM, DA11, DB, EPS, SCALOC, SGN, SMIN, $ SMLNUM, SUML, SUMR, XNORM * .. * .. Local Arrays .. DOUBLE PRECISION DUM( 1 ), VEC( 2, 2 ), X( 2, 2 ) * .. * .. External Functions .. LOGICAL LSAME DOUBLE PRECISION DDOT, DLAMCH, DLANGE EXTERNAL LSAME, DDOT, DLAMCH, DLANGE * .. * .. External Subroutines .. EXTERNAL DLABAD, DLALN2, DLASY2, DSCAL, XERBLA * .. * .. Intrinsic Functions .. INTRINSIC ABS, DBLE, MAX, MIN * .. * .. Executable Statements .. * * Decode and Test input parameters * NOTRNA = LSAME( TRANA, 'N' ) NOTRNB = LSAME( TRANB, 'N' ) * INFO = 0 IF( .NOT.NOTRNA .AND. .NOT.LSAME( TRANA, 'T' ) .AND. .NOT. $ LSAME( TRANA, 'C' ) ) THEN INFO = -1 ELSE IF( .NOT.NOTRNB .AND. .NOT.LSAME( TRANB, 'T' ) .AND. .NOT. $ LSAME( TRANB, 'C' ) ) THEN INFO = -2 ELSE IF( ISGN.NE.1 .AND. ISGN.NE.-1 ) THEN INFO = -3 ELSE IF( M.LT.0 ) THEN INFO = -4 ELSE IF( N.LT.0 ) THEN INFO = -5 ELSE IF( LDA.LT.MAX( 1, M ) ) THEN INFO = -7 ELSE IF( LDB.LT.MAX( 1, N ) ) THEN INFO = -9 ELSE IF( LDC.LT.MAX( 1, M ) ) THEN INFO = -11 END IF IF( INFO.NE.0 ) THEN CALL XERBLA( 'DTRSYL', -INFO ) RETURN END IF * * Quick return if possible * SCALE = ONE IF( M.EQ.0 .OR. N.EQ.0 ) $ RETURN * * Set constants to control overflow * EPS = DLAMCH( 'P' ) SMLNUM = DLAMCH( 'S' ) BIGNUM = ONE / SMLNUM CALL DLABAD( SMLNUM, BIGNUM ) SMLNUM = SMLNUM*DBLE( M*N ) / EPS BIGNUM = ONE / SMLNUM * SMIN = MAX( SMLNUM, EPS*DLANGE( 'M', M, M, A, LDA, DUM ), $ EPS*DLANGE( 'M', N, N, B, LDB, DUM ) ) * SGN = ISGN * IF( NOTRNA .AND. NOTRNB ) THEN * * Solve A*X + ISGN*X*B = scale*C. * * The (K,L)th block of X is determined starting from * bottom-left corner column by column by * * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L) * * Where * M L-1 * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)]. * I=K+1 J=1 * * Start column loop (index = L) * L1 (L2) : column index of the first (first) row of X(K,L). * LNEXT = 1 DO 60 L = 1, N IF( L.LT.LNEXT ) $ GO TO 60 IF( L.EQ.N ) THEN L1 = L L2 = L ELSE IF( B( L+1, L ).NE.ZERO ) THEN L1 = L L2 = L + 1 LNEXT = L + 2 ELSE L1 = L L2 = L LNEXT = L + 1 END IF END IF * * Start row loop (index = K) * K1 (K2): row index of the first (last) row of X(K,L). * KNEXT = M DO 50 K = M, 1, -1 IF( K.GT.KNEXT ) $ GO TO 50 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 * 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( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) SCALOC = ONE * 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 * IF( SCALOC.NE.ONE ) THEN DO 10 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 10 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) * ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN * SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, $ C( MIN( K2+1, M ), L1 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, $ C( MIN( K2+1, M ), L1 ), 1 ) SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) * CALL DLALN2( .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 * IF( SCALOC.NE.ONE ) THEN DO 20 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 20 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) C( K2, L1 ) = X( 2, 1 ) * ELSE IF( L1.NE.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( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) ) * SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA, $ C( MIN( K1+1, M ), L2 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) ) * CALL DLALN2( .TRUE., 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 * IF( SCALOC.NE.ONE ) THEN DO 30 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 30 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) C( K1, L2 ) = X( 2, 1 ) * ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN * SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, $ C( MIN( K2+1, M ), L1 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, $ C( MIN( K2+1, M ), L2 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, $ C( MIN( K2+1, M ), L1 ), 1 ) SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, $ C( MIN( K2+1, M ), L2 ), 1 ) SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 ) VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR ) * CALL DLASY2( .FALSE., .FALSE., ISGN, 2, 2, $ A( K1, K1 ), LDA, B( L1, L1 ), LDB, VEC, $ 2, SCALOC, X, 2, XNORM, IERR ) IF( IERR.NE.0 ) $ INFO = 1 * IF( SCALOC.NE.ONE ) THEN DO 40 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 40 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 * 50 CONTINUE * 60 CONTINUE * ELSE IF( .NOT.NOTRNA .AND. NOTRNB ) THEN * * Solve A' *X + ISGN*X*B = scale*C. * * The (K,L)th block of X is determined starting from * upper-left corner column by column by * * A(K,K)'*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L) * * Where * K-1 L-1 * R(K,L) = SUM [A(I,K)'*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)] * I=1 J=1 * * Start column loop (index = L) * L1 (L2): column index of the first (last) row of X(K,L) * LNEXT = 1 DO 120 L = 1, N IF( L.LT.LNEXT ) $ GO TO 120 IF( L.EQ.N ) THEN L1 = L L2 = L ELSE IF( B( L+1, L ).NE.ZERO ) THEN L1 = L L2 = L + 1 LNEXT = L + 2 ELSE L1 = L L2 = L LNEXT = L + 1 END IF END IF * * Start row loop (index = K) * K1 (K2): row index of the first (last) row of X(K,L) * KNEXT = 1 DO 110 K = 1, M IF( K.LT.KNEXT ) $ GO TO 110 IF( K.EQ.M ) THEN K1 = K K2 = K ELSE IF( A( K+1, K ).NE.ZERO ) THEN K1 = K K2 = K + 1 KNEXT = K + 2 ELSE K1 = K K2 = K KNEXT = K + 1 END IF END IF * IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) SCALOC = ONE * 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 * IF( SCALOC.NE.ONE ) THEN DO 70 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 70 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) * ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 ) SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) * CALL DLALN2( .TRUE., 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 * IF( SCALOC.NE.ONE ) THEN DO 80 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 80 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) C( K2, L1 ) = X( 2, 1 ) * ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) ) * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) ) * CALL DLALN2( .TRUE., 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 * IF( SCALOC.NE.ONE ) THEN DO 90 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 90 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) C( K1, L2 ) = X( 2, 1 ) * ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 ) SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 ) SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) * SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 ) SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 ) VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR ) * CALL DLASY2( .TRUE., .FALSE., ISGN, 2, 2, A( K1, K1 ), $ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X, $ 2, XNORM, IERR ) IF( IERR.NE.0 ) $ INFO = 1 * IF( SCALOC.NE.ONE ) THEN DO 100 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 100 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 * 110 CONTINUE 120 CONTINUE * ELSE IF( .NOT.NOTRNA .AND. .NOT.NOTRNB ) THEN * * Solve A'*X + ISGN*X*B' = scale*C. * * The (K,L)th block of X is determined starting from * top-right corner column by column by * * A(K,K)'*X(K,L) + ISGN*X(K,L)*B(L,L)' = C(K,L) - R(K,L) * * Where * K-1 N * R(K,L) = SUM [A(I,K)'*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)']. * I=1 J=L+1 * * Start column loop (index = L) * L1 (L2): column index of the first (last) row of X(K,L) * LNEXT = N DO 180 L = N, 1, -1 IF( L.GT.LNEXT ) $ GO TO 180 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 * * Start row loop (index = K) * K1 (K2): row index of the first (last) row of X(K,L) * KNEXT = 1 DO 170 K = 1, M IF( K.LT.KNEXT ) $ GO TO 170 IF( K.EQ.M ) THEN K1 = K K2 = K ELSE IF( A( K+1, K ).NE.ZERO ) THEN K1 = K K2 = K + 1 KNEXT = K + 2 ELSE K1 = K K2 = K KNEXT = K + 1 END IF END IF * IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, 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 * 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 * IF( SCALOC.NE.ONE ) THEN DO 130 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 130 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) * ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, 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 ) * SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, 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 ) * CALL DLALN2( .TRUE., 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 * IF( SCALOC.NE.ONE ) THEN DO 140 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 140 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) C( K2, L1 ) = X( 2, 1 ) * ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, 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 ) ) * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, 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 ) ) * CALL DLALN2( .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 * IF( SCALOC.NE.ONE ) THEN DO 150 J = 1, N CALL DSCAL( M, SCALOC, C( 1, J ), 1 ) 150 CONTINUE SCALE = SCALE*SCALOC END IF C( K1, L1 ) = X( 1, 1 ) C( K1, L2 ) = X( 2, 1 ) * ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, 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 ) * SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, 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 ) * SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, 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 ) * 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 ) * CALL DLASY2( .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 * 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 * 170 CONTINUE 180 CONTINUE * ELSE IF( NOTRNA .AND. .NOT.NOTRNB ) THEN * * Solve A*X + ISGN*X*B' = scale*C. * * The (K,L)th block of X is determined starting from * bottom-right corner column by column by * * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L)' = C(K,L) - R(K,L) * * Where * M N * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)']. * I=K+1 J=L+1 * * Start column loop (index = L) * L1 (L2): column index of the first (last) row of X(K,L) * 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 * * Start row loop (index = K) * K1 (K2): row index of the first (last) row of X(K,L) * 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 * 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 * 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 * 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 ) * ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN * 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 ) * 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 ) * CALL DLALN2( .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 * 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 ) * ELSE IF( L1.NE.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-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 ) ) * 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 ) ) * CALL DLALN2( .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 * 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 ) * ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN * 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 ) * 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 ) * 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 ) * 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 ) * CALL DLASY2( .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 * 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 * 230 CONTINUE 240 CONTINUE * END IF * RETURN * * End of DTRSYL * END