LAPACK 3.3.1
Linear Algebra PACKage

dtrsyl.f

Go to the documentation of this file.
00001       SUBROUTINE DTRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C,
00002      $                   LDC, SCALE, INFO )
00003 *
00004 *  -- LAPACK routine (version 3.3.1) --
00005 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00006 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00007 *  -- April 2011                                                      --
00008 *
00009 *     .. Scalar Arguments ..
00010       CHARACTER          TRANA, TRANB
00011       INTEGER            INFO, ISGN, LDA, LDB, LDC, M, N
00012       DOUBLE PRECISION   SCALE
00013 *     ..
00014 *     .. Array Arguments ..
00015       DOUBLE PRECISION   A( LDA, * ), B( LDB, * ), C( LDC, * )
00016 *     ..
00017 *
00018 *  Purpose
00019 *  =======
00020 *
00021 *  DTRSYL solves the real Sylvester matrix equation:
00022 *
00023 *     op(A)*X + X*op(B) = scale*C or
00024 *     op(A)*X - X*op(B) = scale*C,
00025 *
00026 *  where op(A) = A or A**T, and  A and B are both upper quasi-
00027 *  triangular. A is M-by-M and B is N-by-N; the right hand side C and
00028 *  the solution X are M-by-N; and scale is an output scale factor, set
00029 *  <= 1 to avoid overflow in X.
00030 *
00031 *  A and B must be in Schur canonical form (as returned by DHSEQR), that
00032 *  is, block upper triangular with 1-by-1 and 2-by-2 diagonal blocks;
00033 *  each 2-by-2 diagonal block has its diagonal elements equal and its
00034 *  off-diagonal elements of opposite sign.
00035 *
00036 *  Arguments
00037 *  =========
00038 *
00039 *  TRANA   (input) CHARACTER*1
00040 *          Specifies the option op(A):
00041 *          = 'N': op(A) = A    (No transpose)
00042 *          = 'T': op(A) = A**T (Transpose)
00043 *          = 'C': op(A) = A**H (Conjugate transpose = Transpose)
00044 *
00045 *  TRANB   (input) CHARACTER*1
00046 *          Specifies the option op(B):
00047 *          = 'N': op(B) = B    (No transpose)
00048 *          = 'T': op(B) = B**T (Transpose)
00049 *          = 'C': op(B) = B**H (Conjugate transpose = Transpose)
00050 *
00051 *  ISGN    (input) INTEGER
00052 *          Specifies the sign in the equation:
00053 *          = +1: solve op(A)*X + X*op(B) = scale*C
00054 *          = -1: solve op(A)*X - X*op(B) = scale*C
00055 *
00056 *  M       (input) INTEGER
00057 *          The order of the matrix A, and the number of rows in the
00058 *          matrices X and C. M >= 0.
00059 *
00060 *  N       (input) INTEGER
00061 *          The order of the matrix B, and the number of columns in the
00062 *          matrices X and C. N >= 0.
00063 *
00064 *  A       (input) DOUBLE PRECISION array, dimension (LDA,M)
00065 *          The upper quasi-triangular matrix A, in Schur canonical form.
00066 *
00067 *  LDA     (input) INTEGER
00068 *          The leading dimension of the array A. LDA >= max(1,M).
00069 *
00070 *  B       (input) DOUBLE PRECISION array, dimension (LDB,N)
00071 *          The upper quasi-triangular matrix B, in Schur canonical form.
00072 *
00073 *  LDB     (input) INTEGER
00074 *          The leading dimension of the array B. LDB >= max(1,N).
00075 *
00076 *  C       (input/output) DOUBLE PRECISION array, dimension (LDC,N)
00077 *          On entry, the M-by-N right hand side matrix C.
00078 *          On exit, C is overwritten by the solution matrix X.
00079 *
00080 *  LDC     (input) INTEGER
00081 *          The leading dimension of the array C. LDC >= max(1,M)
00082 *
00083 *  SCALE   (output) DOUBLE PRECISION
00084 *          The scale factor, scale, set <= 1 to avoid overflow in X.
00085 *
00086 *  INFO    (output) INTEGER
00087 *          = 0: successful exit
00088 *          < 0: if INFO = -i, the i-th argument had an illegal value
00089 *          = 1: A and B have common or very close eigenvalues; perturbed
00090 *               values were used to solve the equation (but the matrices
00091 *               A and B are unchanged).
00092 *
00093 *  =====================================================================
00094 *
00095 *     .. Parameters ..
00096       DOUBLE PRECISION   ZERO, ONE
00097       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
00098 *     ..
00099 *     .. Local Scalars ..
00100       LOGICAL            NOTRNA, NOTRNB
00101       INTEGER            IERR, J, K, K1, K2, KNEXT, L, L1, L2, LNEXT
00102       DOUBLE PRECISION   A11, BIGNUM, DA11, DB, EPS, SCALOC, SGN, SMIN,
00103      $                   SMLNUM, SUML, SUMR, XNORM
00104 *     ..
00105 *     .. Local Arrays ..
00106       DOUBLE PRECISION   DUM( 1 ), VEC( 2, 2 ), X( 2, 2 )
00107 *     ..
00108 *     .. External Functions ..
00109       LOGICAL            LSAME
00110       DOUBLE PRECISION   DDOT, DLAMCH, DLANGE
00111       EXTERNAL           LSAME, DDOT, DLAMCH, DLANGE
00112 *     ..
00113 *     .. External Subroutines ..
00114       EXTERNAL           DLABAD, DLALN2, DLASY2, DSCAL, XERBLA
00115 *     ..
00116 *     .. Intrinsic Functions ..
00117       INTRINSIC          ABS, DBLE, MAX, MIN
00118 *     ..
00119 *     .. Executable Statements ..
00120 *
00121 *     Decode and Test input parameters
00122 *
00123       NOTRNA = LSAME( TRANA, 'N' )
00124       NOTRNB = LSAME( TRANB, 'N' )
00125 *
00126       INFO = 0
00127       IF( .NOT.NOTRNA .AND. .NOT.LSAME( TRANA, 'T' ) .AND. .NOT.
00128      $    LSAME( TRANA, 'C' ) ) THEN
00129          INFO = -1
00130       ELSE IF( .NOT.NOTRNB .AND. .NOT.LSAME( TRANB, 'T' ) .AND. .NOT.
00131      $         LSAME( TRANB, 'C' ) ) THEN
00132          INFO = -2
00133       ELSE IF( ISGN.NE.1 .AND. ISGN.NE.-1 ) THEN
00134          INFO = -3
00135       ELSE IF( M.LT.0 ) THEN
00136          INFO = -4
00137       ELSE IF( N.LT.0 ) THEN
00138          INFO = -5
00139       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
00140          INFO = -7
00141       ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
00142          INFO = -9
00143       ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
00144          INFO = -11
00145       END IF
00146       IF( INFO.NE.0 ) THEN
00147          CALL XERBLA( 'DTRSYL', -INFO )
00148          RETURN
00149       END IF
00150 *
00151 *     Quick return if possible
00152 *
00153       SCALE = ONE
00154       IF( M.EQ.0 .OR. N.EQ.0 )
00155      $   RETURN
00156 *
00157 *     Set constants to control overflow
00158 *
00159       EPS = DLAMCH( 'P' )
00160       SMLNUM = DLAMCH( 'S' )
00161       BIGNUM = ONE / SMLNUM
00162       CALL DLABAD( SMLNUM, BIGNUM )
00163       SMLNUM = SMLNUM*DBLE( M*N ) / EPS
00164       BIGNUM = ONE / SMLNUM
00165 *
00166       SMIN = MAX( SMLNUM, EPS*DLANGE( 'M', M, M, A, LDA, DUM ),
00167      $       EPS*DLANGE( 'M', N, N, B, LDB, DUM ) )
00168 *
00169       SGN = ISGN
00170 *
00171       IF( NOTRNA .AND. NOTRNB ) THEN
00172 *
00173 *        Solve    A*X + ISGN*X*B = scale*C.
00174 *
00175 *        The (K,L)th block of X is determined starting from
00176 *        bottom-left corner column by column by
00177 *
00178 *         A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
00179 *
00180 *        Where
00181 *                  M                         L-1
00182 *        R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)].
00183 *                I=K+1                       J=1
00184 *
00185 *        Start column loop (index = L)
00186 *        L1 (L2) : column index of the first (first) row of X(K,L).
00187 *
00188          LNEXT = 1
00189          DO 60 L = 1, N
00190             IF( L.LT.LNEXT )
00191      $         GO TO 60
00192             IF( L.EQ.N ) THEN
00193                L1 = L
00194                L2 = L
00195             ELSE
00196                IF( B( L+1, L ).NE.ZERO ) THEN
00197                   L1 = L
00198                   L2 = L + 1
00199                   LNEXT = L + 2
00200                ELSE
00201                   L1 = L
00202                   L2 = L
00203                   LNEXT = L + 1
00204                END IF
00205             END IF
00206 *
00207 *           Start row loop (index = K)
00208 *           K1 (K2): row index of the first (last) row of X(K,L).
00209 *
00210             KNEXT = M
00211             DO 50 K = M, 1, -1
00212                IF( K.GT.KNEXT )
00213      $            GO TO 50
00214                IF( K.EQ.1 ) THEN
00215                   K1 = K
00216                   K2 = K
00217                ELSE
00218                   IF( A( K, K-1 ).NE.ZERO ) THEN
00219                      K1 = K - 1
00220                      K2 = K
00221                      KNEXT = K - 2
00222                   ELSE
00223                      K1 = K
00224                      K2 = K
00225                      KNEXT = K - 1
00226                   END IF
00227                END IF
00228 *
00229                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
00230                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
00231      $                   C( MIN( K1+1, M ), L1 ), 1 )
00232                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00233                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00234                   SCALOC = ONE
00235 *
00236                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
00237                   DA11 = ABS( A11 )
00238                   IF( DA11.LE.SMIN ) THEN
00239                      A11 = SMIN
00240                      DA11 = SMIN
00241                      INFO = 1
00242                   END IF
00243                   DB = ABS( VEC( 1, 1 ) )
00244                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
00245                      IF( DB.GT.BIGNUM*DA11 )
00246      $                  SCALOC = ONE / DB
00247                   END IF
00248                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
00249 *
00250                   IF( SCALOC.NE.ONE ) THEN
00251                      DO 10 J = 1, N
00252                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00253    10                CONTINUE
00254                      SCALE = SCALE*SCALOC
00255                   END IF
00256                   C( K1, L1 ) = X( 1, 1 )
00257 *
00258                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
00259 *
00260                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
00261      $                   C( MIN( K2+1, M ), L1 ), 1 )
00262                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00263                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00264 *
00265                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
00266      $                   C( MIN( K2+1, M ), L1 ), 1 )
00267                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
00268                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00269 *
00270                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ),
00271      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
00272      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00273                   IF( IERR.NE.0 )
00274      $               INFO = 1
00275 *
00276                   IF( SCALOC.NE.ONE ) THEN
00277                      DO 20 J = 1, N
00278                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00279    20                CONTINUE
00280                      SCALE = SCALE*SCALOC
00281                   END IF
00282                   C( K1, L1 ) = X( 1, 1 )
00283                   C( K2, L1 ) = X( 2, 1 )
00284 *
00285                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
00286 *
00287                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
00288      $                   C( MIN( K1+1, M ), L1 ), 1 )
00289                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00290                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
00291 *
00292                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
00293      $                   C( MIN( K1+1, M ), L2 ), 1 )
00294                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
00295                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
00296 *
00297                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ),
00298      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
00299      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00300                   IF( IERR.NE.0 )
00301      $               INFO = 1
00302 *
00303                   IF( SCALOC.NE.ONE ) THEN
00304                      DO 30 J = 1, N
00305                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00306    30                CONTINUE
00307                      SCALE = SCALE*SCALOC
00308                   END IF
00309                   C( K1, L1 ) = X( 1, 1 )
00310                   C( K1, L2 ) = X( 2, 1 )
00311 *
00312                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
00313 *
00314                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
00315      $                   C( MIN( K2+1, M ), L1 ), 1 )
00316                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00317                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00318 *
00319                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
00320      $                   C( MIN( K2+1, M ), L2 ), 1 )
00321                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
00322                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
00323 *
00324                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
00325      $                   C( MIN( K2+1, M ), L1 ), 1 )
00326                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
00327                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00328 *
00329                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
00330      $                   C( MIN( K2+1, M ), L2 ), 1 )
00331                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 )
00332                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
00333 *
00334                   CALL DLASY2( .FALSE., .FALSE., ISGN, 2, 2,
00335      $                         A( K1, K1 ), LDA, B( L1, L1 ), LDB, VEC,
00336      $                         2, SCALOC, X, 2, XNORM, IERR )
00337                   IF( IERR.NE.0 )
00338      $               INFO = 1
00339 *
00340                   IF( SCALOC.NE.ONE ) THEN
00341                      DO 40 J = 1, N
00342                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00343    40                CONTINUE
00344                      SCALE = SCALE*SCALOC
00345                   END IF
00346                   C( K1, L1 ) = X( 1, 1 )
00347                   C( K1, L2 ) = X( 1, 2 )
00348                   C( K2, L1 ) = X( 2, 1 )
00349                   C( K2, L2 ) = X( 2, 2 )
00350                END IF
00351 *
00352    50       CONTINUE
00353 *
00354    60    CONTINUE
00355 *
00356       ELSE IF( .NOT.NOTRNA .AND. NOTRNB ) THEN
00357 *
00358 *        Solve    A**T *X + ISGN*X*B = scale*C.
00359 *
00360 *        The (K,L)th block of X is determined starting from
00361 *        upper-left corner column by column by
00362 *
00363 *          A(K,K)**T*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
00364 *
00365 *        Where
00366 *                   K-1                          L-1
00367 *          R(K,L) = SUM [A(I,K)**T*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)]
00368 *                   I=1                          J=1
00369 *
00370 *        Start column loop (index = L)
00371 *        L1 (L2): column index of the first (last) row of X(K,L)
00372 *
00373          LNEXT = 1
00374          DO 120 L = 1, N
00375             IF( L.LT.LNEXT )
00376      $         GO TO 120
00377             IF( L.EQ.N ) THEN
00378                L1 = L
00379                L2 = L
00380             ELSE
00381                IF( B( L+1, L ).NE.ZERO ) THEN
00382                   L1 = L
00383                   L2 = L + 1
00384                   LNEXT = L + 2
00385                ELSE
00386                   L1 = L
00387                   L2 = L
00388                   LNEXT = L + 1
00389                END IF
00390             END IF
00391 *
00392 *           Start row loop (index = K)
00393 *           K1 (K2): row index of the first (last) row of X(K,L)
00394 *
00395             KNEXT = 1
00396             DO 110 K = 1, M
00397                IF( K.LT.KNEXT )
00398      $            GO TO 110
00399                IF( K.EQ.M ) THEN
00400                   K1 = K
00401                   K2 = K
00402                ELSE
00403                   IF( A( K+1, K ).NE.ZERO ) THEN
00404                      K1 = K
00405                      K2 = K + 1
00406                      KNEXT = K + 2
00407                   ELSE
00408                      K1 = K
00409                      K2 = K
00410                      KNEXT = K + 1
00411                   END IF
00412                END IF
00413 *
00414                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
00415                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00416                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00417                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00418                   SCALOC = ONE
00419 *
00420                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
00421                   DA11 = ABS( A11 )
00422                   IF( DA11.LE.SMIN ) THEN
00423                      A11 = SMIN
00424                      DA11 = SMIN
00425                      INFO = 1
00426                   END IF
00427                   DB = ABS( VEC( 1, 1 ) )
00428                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
00429                      IF( DB.GT.BIGNUM*DA11 )
00430      $                  SCALOC = ONE / DB
00431                   END IF
00432                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
00433 *
00434                   IF( SCALOC.NE.ONE ) THEN
00435                      DO 70 J = 1, N
00436                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00437    70                CONTINUE
00438                      SCALE = SCALE*SCALOC
00439                   END IF
00440                   C( K1, L1 ) = X( 1, 1 )
00441 *
00442                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
00443 *
00444                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00445                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00446                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00447 *
00448                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
00449                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
00450                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00451 *
00452                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ),
00453      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
00454      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00455                   IF( IERR.NE.0 )
00456      $               INFO = 1
00457 *
00458                   IF( SCALOC.NE.ONE ) THEN
00459                      DO 80 J = 1, N
00460                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00461    80                CONTINUE
00462                      SCALE = SCALE*SCALOC
00463                   END IF
00464                   C( K1, L1 ) = X( 1, 1 )
00465                   C( K2, L1 ) = X( 2, 1 )
00466 *
00467                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
00468 *
00469                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00470                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00471                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
00472 *
00473                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
00474                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
00475                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
00476 *
00477                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ),
00478      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
00479      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00480                   IF( IERR.NE.0 )
00481      $               INFO = 1
00482 *
00483                   IF( SCALOC.NE.ONE ) THEN
00484                      DO 90 J = 1, N
00485                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00486    90                CONTINUE
00487                      SCALE = SCALE*SCALOC
00488                   END IF
00489                   C( K1, L1 ) = X( 1, 1 )
00490                   C( K1, L2 ) = X( 2, 1 )
00491 *
00492                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
00493 *
00494                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00495                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
00496                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00497 *
00498                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
00499                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
00500                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
00501 *
00502                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
00503                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
00504                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00505 *
00506                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 )
00507                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 )
00508                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
00509 *
00510                   CALL DLASY2( .TRUE., .FALSE., ISGN, 2, 2, A( K1, K1 ),
00511      $                         LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
00512      $                         2, XNORM, IERR )
00513                   IF( IERR.NE.0 )
00514      $               INFO = 1
00515 *
00516                   IF( SCALOC.NE.ONE ) THEN
00517                      DO 100 J = 1, N
00518                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00519   100                CONTINUE
00520                      SCALE = SCALE*SCALOC
00521                   END IF
00522                   C( K1, L1 ) = X( 1, 1 )
00523                   C( K1, L2 ) = X( 1, 2 )
00524                   C( K2, L1 ) = X( 2, 1 )
00525                   C( K2, L2 ) = X( 2, 2 )
00526                END IF
00527 *
00528   110       CONTINUE
00529   120    CONTINUE
00530 *
00531       ELSE IF( .NOT.NOTRNA .AND. .NOT.NOTRNB ) THEN
00532 *
00533 *        Solve    A**T*X + ISGN*X*B**T = scale*C.
00534 *
00535 *        The (K,L)th block of X is determined starting from
00536 *        top-right corner column by column by
00537 *
00538 *           A(K,K)**T*X(K,L) + ISGN*X(K,L)*B(L,L)**T = C(K,L) - R(K,L)
00539 *
00540 *        Where
00541 *                     K-1                            N
00542 *            R(K,L) = SUM [A(I,K)**T*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**T].
00543 *                     I=1                          J=L+1
00544 *
00545 *        Start column loop (index = L)
00546 *        L1 (L2): column index of the first (last) row of X(K,L)
00547 *
00548          LNEXT = N
00549          DO 180 L = N, 1, -1
00550             IF( L.GT.LNEXT )
00551      $         GO TO 180
00552             IF( L.EQ.1 ) THEN
00553                L1 = L
00554                L2 = L
00555             ELSE
00556                IF( B( L, L-1 ).NE.ZERO ) THEN
00557                   L1 = L - 1
00558                   L2 = L
00559                   LNEXT = L - 2
00560                ELSE
00561                   L1 = L
00562                   L2 = L
00563                   LNEXT = L - 1
00564                END IF
00565             END IF
00566 *
00567 *           Start row loop (index = K)
00568 *           K1 (K2): row index of the first (last) row of X(K,L)
00569 *
00570             KNEXT = 1
00571             DO 170 K = 1, M
00572                IF( K.LT.KNEXT )
00573      $            GO TO 170
00574                IF( K.EQ.M ) THEN
00575                   K1 = K
00576                   K2 = K
00577                ELSE
00578                   IF( A( K+1, K ).NE.ZERO ) THEN
00579                      K1 = K
00580                      K2 = K + 1
00581                      KNEXT = K + 2
00582                   ELSE
00583                      K1 = K
00584                      K2 = K
00585                      KNEXT = K + 1
00586                   END IF
00587                END IF
00588 *
00589                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
00590                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00591                   SUMR = DDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC,
00592      $                   B( L1, MIN( L1+1, N ) ), LDB )
00593                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00594                   SCALOC = ONE
00595 *
00596                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
00597                   DA11 = ABS( A11 )
00598                   IF( DA11.LE.SMIN ) THEN
00599                      A11 = SMIN
00600                      DA11 = SMIN
00601                      INFO = 1
00602                   END IF
00603                   DB = ABS( VEC( 1, 1 ) )
00604                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
00605                      IF( DB.GT.BIGNUM*DA11 )
00606      $                  SCALOC = ONE / DB
00607                   END IF
00608                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
00609 *
00610                   IF( SCALOC.NE.ONE ) THEN
00611                      DO 130 J = 1, N
00612                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00613   130                CONTINUE
00614                      SCALE = SCALE*SCALOC
00615                   END IF
00616                   C( K1, L1 ) = X( 1, 1 )
00617 *
00618                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
00619 *
00620                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00621                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00622      $                   B( L1, MIN( L2+1, N ) ), LDB )
00623                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00624 *
00625                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
00626                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
00627      $                   B( L1, MIN( L2+1, N ) ), LDB )
00628                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00629 *
00630                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ),
00631      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
00632      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00633                   IF( IERR.NE.0 )
00634      $               INFO = 1
00635 *
00636                   IF( SCALOC.NE.ONE ) THEN
00637                      DO 140 J = 1, N
00638                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00639   140                CONTINUE
00640                      SCALE = SCALE*SCALOC
00641                   END IF
00642                   C( K1, L1 ) = X( 1, 1 )
00643                   C( K2, L1 ) = X( 2, 1 )
00644 *
00645                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
00646 *
00647                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00648                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00649      $                   B( L1, MIN( L2+1, N ) ), LDB )
00650                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
00651 *
00652                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
00653                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00654      $                   B( L2, MIN( L2+1, N ) ), LDB )
00655                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
00656 *
00657                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ),
00658      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
00659      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00660                   IF( IERR.NE.0 )
00661      $               INFO = 1
00662 *
00663                   IF( SCALOC.NE.ONE ) THEN
00664                      DO 150 J = 1, N
00665                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00666   150                CONTINUE
00667                      SCALE = SCALE*SCALOC
00668                   END IF
00669                   C( K1, L1 ) = X( 1, 1 )
00670                   C( K1, L2 ) = X( 2, 1 )
00671 *
00672                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
00673 *
00674                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
00675                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00676      $                   B( L1, MIN( L2+1, N ) ), LDB )
00677                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00678 *
00679                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
00680                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00681      $                   B( L2, MIN( L2+1, N ) ), LDB )
00682                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
00683 *
00684                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
00685                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
00686      $                   B( L1, MIN( L2+1, N ) ), LDB )
00687                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00688 *
00689                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 )
00690                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
00691      $                   B( L2, MIN( L2+1, N ) ), LDB )
00692                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
00693 *
00694                   CALL DLASY2( .TRUE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
00695      $                         LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
00696      $                         2, XNORM, IERR )
00697                   IF( IERR.NE.0 )
00698      $               INFO = 1
00699 *
00700                   IF( SCALOC.NE.ONE ) THEN
00701                      DO 160 J = 1, N
00702                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00703   160                CONTINUE
00704                      SCALE = SCALE*SCALOC
00705                   END IF
00706                   C( K1, L1 ) = X( 1, 1 )
00707                   C( K1, L2 ) = X( 1, 2 )
00708                   C( K2, L1 ) = X( 2, 1 )
00709                   C( K2, L2 ) = X( 2, 2 )
00710                END IF
00711 *
00712   170       CONTINUE
00713   180    CONTINUE
00714 *
00715       ELSE IF( NOTRNA .AND. .NOT.NOTRNB ) THEN
00716 *
00717 *        Solve    A*X + ISGN*X*B**T = scale*C.
00718 *
00719 *        The (K,L)th block of X is determined starting from
00720 *        bottom-right corner column by column by
00721 *
00722 *            A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L)**T = C(K,L) - R(K,L)
00723 *
00724 *        Where
00725 *                      M                          N
00726 *            R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**T].
00727 *                    I=K+1                      J=L+1
00728 *
00729 *        Start column loop (index = L)
00730 *        L1 (L2): column index of the first (last) row of X(K,L)
00731 *
00732          LNEXT = N
00733          DO 240 L = N, 1, -1
00734             IF( L.GT.LNEXT )
00735      $         GO TO 240
00736             IF( L.EQ.1 ) THEN
00737                L1 = L
00738                L2 = L
00739             ELSE
00740                IF( B( L, L-1 ).NE.ZERO ) THEN
00741                   L1 = L - 1
00742                   L2 = L
00743                   LNEXT = L - 2
00744                ELSE
00745                   L1 = L
00746                   L2 = L
00747                   LNEXT = L - 1
00748                END IF
00749             END IF
00750 *
00751 *           Start row loop (index = K)
00752 *           K1 (K2): row index of the first (last) row of X(K,L)
00753 *
00754             KNEXT = M
00755             DO 230 K = M, 1, -1
00756                IF( K.GT.KNEXT )
00757      $            GO TO 230
00758                IF( K.EQ.1 ) THEN
00759                   K1 = K
00760                   K2 = K
00761                ELSE
00762                   IF( A( K, K-1 ).NE.ZERO ) THEN
00763                      K1 = K - 1
00764                      K2 = K
00765                      KNEXT = K - 2
00766                   ELSE
00767                      K1 = K
00768                      K2 = K
00769                      KNEXT = K - 1
00770                   END IF
00771                END IF
00772 *
00773                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
00774                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
00775      $                   C( MIN( K1+1, M ), L1 ), 1 )
00776                   SUMR = DDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC,
00777      $                   B( L1, MIN( L1+1, N ) ), LDB )
00778                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00779                   SCALOC = ONE
00780 *
00781                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
00782                   DA11 = ABS( A11 )
00783                   IF( DA11.LE.SMIN ) THEN
00784                      A11 = SMIN
00785                      DA11 = SMIN
00786                      INFO = 1
00787                   END IF
00788                   DB = ABS( VEC( 1, 1 ) )
00789                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
00790                      IF( DB.GT.BIGNUM*DA11 )
00791      $                  SCALOC = ONE / DB
00792                   END IF
00793                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
00794 *
00795                   IF( SCALOC.NE.ONE ) THEN
00796                      DO 190 J = 1, N
00797                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00798   190                CONTINUE
00799                      SCALE = SCALE*SCALOC
00800                   END IF
00801                   C( K1, L1 ) = X( 1, 1 )
00802 *
00803                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
00804 *
00805                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
00806      $                   C( MIN( K2+1, M ), L1 ), 1 )
00807                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00808      $                   B( L1, MIN( L2+1, N ) ), LDB )
00809                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00810 *
00811                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
00812      $                   C( MIN( K2+1, M ), L1 ), 1 )
00813                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
00814      $                   B( L1, MIN( L2+1, N ) ), LDB )
00815                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00816 *
00817                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ),
00818      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
00819      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00820                   IF( IERR.NE.0 )
00821      $               INFO = 1
00822 *
00823                   IF( SCALOC.NE.ONE ) THEN
00824                      DO 200 J = 1, N
00825                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00826   200                CONTINUE
00827                      SCALE = SCALE*SCALOC
00828                   END IF
00829                   C( K1, L1 ) = X( 1, 1 )
00830                   C( K2, L1 ) = X( 2, 1 )
00831 *
00832                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
00833 *
00834                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
00835      $                   C( MIN( K1+1, M ), L1 ), 1 )
00836                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00837      $                   B( L1, MIN( L2+1, N ) ), LDB )
00838                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
00839 *
00840                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
00841      $                   C( MIN( K1+1, M ), L2 ), 1 )
00842                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00843      $                   B( L2, MIN( L2+1, N ) ), LDB )
00844                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
00845 *
00846                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ),
00847      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
00848      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
00849                   IF( IERR.NE.0 )
00850      $               INFO = 1
00851 *
00852                   IF( SCALOC.NE.ONE ) THEN
00853                      DO 210 J = 1, N
00854                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00855   210                CONTINUE
00856                      SCALE = SCALE*SCALOC
00857                   END IF
00858                   C( K1, L1 ) = X( 1, 1 )
00859                   C( K1, L2 ) = X( 2, 1 )
00860 *
00861                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
00862 *
00863                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
00864      $                   C( MIN( K2+1, M ), L1 ), 1 )
00865                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00866      $                   B( L1, MIN( L2+1, N ) ), LDB )
00867                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
00868 *
00869                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
00870      $                   C( MIN( K2+1, M ), L2 ), 1 )
00871                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
00872      $                   B( L2, MIN( L2+1, N ) ), LDB )
00873                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
00874 *
00875                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
00876      $                   C( MIN( K2+1, M ), L1 ), 1 )
00877                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
00878      $                   B( L1, MIN( L2+1, N ) ), LDB )
00879                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
00880 *
00881                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
00882      $                   C( MIN( K2+1, M ), L2 ), 1 )
00883                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
00884      $                   B( L2, MIN( L2+1, N ) ), LDB )
00885                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
00886 *
00887                   CALL DLASY2( .FALSE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
00888      $                         LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
00889      $                         2, XNORM, IERR )
00890                   IF( IERR.NE.0 )
00891      $               INFO = 1
00892 *
00893                   IF( SCALOC.NE.ONE ) THEN
00894                      DO 220 J = 1, N
00895                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
00896   220                CONTINUE
00897                      SCALE = SCALE*SCALOC
00898                   END IF
00899                   C( K1, L1 ) = X( 1, 1 )
00900                   C( K1, L2 ) = X( 1, 2 )
00901                   C( K2, L1 ) = X( 2, 1 )
00902                   C( K2, L2 ) = X( 2, 2 )
00903                END IF
00904 *
00905   230       CONTINUE
00906   240    CONTINUE
00907 *
00908       END IF
00909 *
00910       RETURN
00911 *
00912 *     End of DTRSYL
00913 *
00914       END
 All Files Functions