001:       SUBROUTINE DTRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C,
002:      $                   LDC, SCALE, INFO )
003: *
004: *  -- LAPACK routine (version 3.2) --
005: *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
006: *     November 2006
007: *
008: *     .. Scalar Arguments ..
009:       CHARACTER          TRANA, TRANB
010:       INTEGER            INFO, ISGN, LDA, LDB, LDC, M, N
011:       DOUBLE PRECISION   SCALE
012: *     ..
013: *     .. Array Arguments ..
014:       DOUBLE PRECISION   A( LDA, * ), B( LDB, * ), C( LDC, * )
015: *     ..
016: *
017: *  Purpose
018: *  =======
019: *
020: *  DTRSYL solves the real Sylvester matrix equation:
021: *
022: *     op(A)*X + X*op(B) = scale*C or
023: *     op(A)*X - X*op(B) = scale*C,
024: *
025: *  where op(A) = A or A**T, and  A and B are both upper quasi-
026: *  triangular. A is M-by-M and B is N-by-N; the right hand side C and
027: *  the solution X are M-by-N; and scale is an output scale factor, set
028: *  <= 1 to avoid overflow in X.
029: *
030: *  A and B must be in Schur canonical form (as returned by DHSEQR), that
031: *  is, block upper triangular with 1-by-1 and 2-by-2 diagonal blocks;
032: *  each 2-by-2 diagonal block has its diagonal elements equal and its
033: *  off-diagonal elements of opposite sign.
034: *
035: *  Arguments
036: *  =========
037: *
038: *  TRANA   (input) CHARACTER*1
039: *          Specifies the option op(A):
040: *          = 'N': op(A) = A    (No transpose)
041: *          = 'T': op(A) = A**T (Transpose)
042: *          = 'C': op(A) = A**H (Conjugate transpose = Transpose)
043: *
044: *  TRANB   (input) CHARACTER*1
045: *          Specifies the option op(B):
046: *          = 'N': op(B) = B    (No transpose)
047: *          = 'T': op(B) = B**T (Transpose)
048: *          = 'C': op(B) = B**H (Conjugate transpose = Transpose)
049: *
050: *  ISGN    (input) INTEGER
051: *          Specifies the sign in the equation:
052: *          = +1: solve op(A)*X + X*op(B) = scale*C
053: *          = -1: solve op(A)*X - X*op(B) = scale*C
054: *
055: *  M       (input) INTEGER
056: *          The order of the matrix A, and the number of rows in the
057: *          matrices X and C. M >= 0.
058: *
059: *  N       (input) INTEGER
060: *          The order of the matrix B, and the number of columns in the
061: *          matrices X and C. N >= 0.
062: *
063: *  A       (input) DOUBLE PRECISION array, dimension (LDA,M)
064: *          The upper quasi-triangular matrix A, in Schur canonical form.
065: *
066: *  LDA     (input) INTEGER
067: *          The leading dimension of the array A. LDA >= max(1,M).
068: *
069: *  B       (input) DOUBLE PRECISION array, dimension (LDB,N)
070: *          The upper quasi-triangular matrix B, in Schur canonical form.
071: *
072: *  LDB     (input) INTEGER
073: *          The leading dimension of the array B. LDB >= max(1,N).
074: *
075: *  C       (input/output) DOUBLE PRECISION array, dimension (LDC,N)
076: *          On entry, the M-by-N right hand side matrix C.
077: *          On exit, C is overwritten by the solution matrix X.
078: *
079: *  LDC     (input) INTEGER
080: *          The leading dimension of the array C. LDC >= max(1,M)
081: *
082: *  SCALE   (output) DOUBLE PRECISION
083: *          The scale factor, scale, set <= 1 to avoid overflow in X.
084: *
085: *  INFO    (output) INTEGER
086: *          = 0: successful exit
087: *          < 0: if INFO = -i, the i-th argument had an illegal value
088: *          = 1: A and B have common or very close eigenvalues; perturbed
089: *               values were used to solve the equation (but the matrices
090: *               A and B are unchanged).
091: *
092: *  =====================================================================
093: *
094: *     .. Parameters ..
095:       DOUBLE PRECISION   ZERO, ONE
096:       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
097: *     ..
098: *     .. Local Scalars ..
099:       LOGICAL            NOTRNA, NOTRNB
100:       INTEGER            IERR, J, K, K1, K2, KNEXT, L, L1, L2, LNEXT
101:       DOUBLE PRECISION   A11, BIGNUM, DA11, DB, EPS, SCALOC, SGN, SMIN,
102:      $                   SMLNUM, SUML, SUMR, XNORM
103: *     ..
104: *     .. Local Arrays ..
105:       DOUBLE PRECISION   DUM( 1 ), VEC( 2, 2 ), X( 2, 2 )
106: *     ..
107: *     .. External Functions ..
108:       LOGICAL            LSAME
109:       DOUBLE PRECISION   DDOT, DLAMCH, DLANGE
110:       EXTERNAL           LSAME, DDOT, DLAMCH, DLANGE
111: *     ..
112: *     .. External Subroutines ..
113:       EXTERNAL           DLABAD, DLALN2, DLASY2, DSCAL, XERBLA
114: *     ..
115: *     .. Intrinsic Functions ..
116:       INTRINSIC          ABS, DBLE, MAX, MIN
117: *     ..
118: *     .. Executable Statements ..
119: *
120: *     Decode and Test input parameters
121: *
122:       NOTRNA = LSAME( TRANA, 'N' )
123:       NOTRNB = LSAME( TRANB, 'N' )
124: *
125:       INFO = 0
126:       IF( .NOT.NOTRNA .AND. .NOT.LSAME( TRANA, 'T' ) .AND. .NOT.
127:      $    LSAME( TRANA, 'C' ) ) THEN
128:          INFO = -1
129:       ELSE IF( .NOT.NOTRNB .AND. .NOT.LSAME( TRANB, 'T' ) .AND. .NOT.
130:      $         LSAME( TRANB, 'C' ) ) THEN
131:          INFO = -2
132:       ELSE IF( ISGN.NE.1 .AND. ISGN.NE.-1 ) THEN
133:          INFO = -3
134:       ELSE IF( M.LT.0 ) THEN
135:          INFO = -4
136:       ELSE IF( N.LT.0 ) THEN
137:          INFO = -5
138:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
139:          INFO = -7
140:       ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
141:          INFO = -9
142:       ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
143:          INFO = -11
144:       END IF
145:       IF( INFO.NE.0 ) THEN
146:          CALL XERBLA( 'DTRSYL', -INFO )
147:          RETURN
148:       END IF
149: *
150: *     Quick return if possible
151: *
152:       IF( M.EQ.0 .OR. N.EQ.0 )
153:      $   RETURN
154: *
155: *     Set constants to control overflow
156: *
157:       EPS = DLAMCH( 'P' )
158:       SMLNUM = DLAMCH( 'S' )
159:       BIGNUM = ONE / SMLNUM
160:       CALL DLABAD( SMLNUM, BIGNUM )
161:       SMLNUM = SMLNUM*DBLE( M*N ) / EPS
162:       BIGNUM = ONE / SMLNUM
163: *
164:       SMIN = MAX( SMLNUM, EPS*DLANGE( 'M', M, M, A, LDA, DUM ),
165:      $       EPS*DLANGE( 'M', N, N, B, LDB, DUM ) )
166: *
167:       SCALE = ONE
168:       SGN = ISGN
169: *
170:       IF( NOTRNA .AND. NOTRNB ) THEN
171: *
172: *        Solve    A*X + ISGN*X*B = scale*C.
173: *
174: *        The (K,L)th block of X is determined starting from
175: *        bottom-left corner column by column by
176: *
177: *         A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
178: *
179: *        Where
180: *                  M                         L-1
181: *        R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)].
182: *                I=K+1                       J=1
183: *
184: *        Start column loop (index = L)
185: *        L1 (L2) : column index of the first (first) row of X(K,L).
186: *
187:          LNEXT = 1
188:          DO 60 L = 1, N
189:             IF( L.LT.LNEXT )
190:      $         GO TO 60
191:             IF( L.EQ.N ) THEN
192:                L1 = L
193:                L2 = L
194:             ELSE
195:                IF( B( L+1, L ).NE.ZERO ) THEN
196:                   L1 = L
197:                   L2 = L + 1
198:                   LNEXT = L + 2
199:                ELSE
200:                   L1 = L
201:                   L2 = L
202:                   LNEXT = L + 1
203:                END IF
204:             END IF
205: *
206: *           Start row loop (index = K)
207: *           K1 (K2): row index of the first (last) row of X(K,L).
208: *
209:             KNEXT = M
210:             DO 50 K = M, 1, -1
211:                IF( K.GT.KNEXT )
212:      $            GO TO 50
213:                IF( K.EQ.1 ) THEN
214:                   K1 = K
215:                   K2 = K
216:                ELSE
217:                   IF( A( K, K-1 ).NE.ZERO ) THEN
218:                      K1 = K - 1
219:                      K2 = K
220:                      KNEXT = K - 2
221:                   ELSE
222:                      K1 = K
223:                      K2 = K
224:                      KNEXT = K - 1
225:                   END IF
226:                END IF
227: *
228:                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
229:                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
230:      $                   C( MIN( K1+1, M ), L1 ), 1 )
231:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
232:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
233:                   SCALOC = ONE
234: *
235:                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
236:                   DA11 = ABS( A11 )
237:                   IF( DA11.LE.SMIN ) THEN
238:                      A11 = SMIN
239:                      DA11 = SMIN
240:                      INFO = 1
241:                   END IF
242:                   DB = ABS( VEC( 1, 1 ) )
243:                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
244:                      IF( DB.GT.BIGNUM*DA11 )
245:      $                  SCALOC = ONE / DB
246:                   END IF
247:                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
248: *
249:                   IF( SCALOC.NE.ONE ) THEN
250:                      DO 10 J = 1, N
251:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
252:    10                CONTINUE
253:                      SCALE = SCALE*SCALOC
254:                   END IF
255:                   C( K1, L1 ) = X( 1, 1 )
256: *
257:                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
258: *
259:                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
260:      $                   C( MIN( K2+1, M ), L1 ), 1 )
261:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
262:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
263: *
264:                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
265:      $                   C( MIN( K2+1, M ), L1 ), 1 )
266:                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
267:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
268: *
269:                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ),
270:      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
271:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
272:                   IF( IERR.NE.0 )
273:      $               INFO = 1
274: *
275:                   IF( SCALOC.NE.ONE ) THEN
276:                      DO 20 J = 1, N
277:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
278:    20                CONTINUE
279:                      SCALE = SCALE*SCALOC
280:                   END IF
281:                   C( K1, L1 ) = X( 1, 1 )
282:                   C( K2, L1 ) = X( 2, 1 )
283: *
284:                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
285: *
286:                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
287:      $                   C( MIN( K1+1, M ), L1 ), 1 )
288:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
289:                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
290: *
291:                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
292:      $                   C( MIN( K1+1, M ), L2 ), 1 )
293:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
294:                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
295: *
296:                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ),
297:      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
298:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
299:                   IF( IERR.NE.0 )
300:      $               INFO = 1
301: *
302:                   IF( SCALOC.NE.ONE ) THEN
303:                      DO 30 J = 1, N
304:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
305:    30                CONTINUE
306:                      SCALE = SCALE*SCALOC
307:                   END IF
308:                   C( K1, L1 ) = X( 1, 1 )
309:                   C( K1, L2 ) = X( 2, 1 )
310: *
311:                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
312: *
313:                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
314:      $                   C( MIN( K2+1, M ), L1 ), 1 )
315:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
316:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
317: *
318:                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
319:      $                   C( MIN( K2+1, M ), L2 ), 1 )
320:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
321:                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
322: *
323:                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
324:      $                   C( MIN( K2+1, M ), L1 ), 1 )
325:                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
326:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
327: *
328:                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
329:      $                   C( MIN( K2+1, M ), L2 ), 1 )
330:                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 )
331:                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
332: *
333:                   CALL DLASY2( .FALSE., .FALSE., ISGN, 2, 2,
334:      $                         A( K1, K1 ), LDA, B( L1, L1 ), LDB, VEC,
335:      $                         2, SCALOC, X, 2, XNORM, IERR )
336:                   IF( IERR.NE.0 )
337:      $               INFO = 1
338: *
339:                   IF( SCALOC.NE.ONE ) THEN
340:                      DO 40 J = 1, N
341:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
342:    40                CONTINUE
343:                      SCALE = SCALE*SCALOC
344:                   END IF
345:                   C( K1, L1 ) = X( 1, 1 )
346:                   C( K1, L2 ) = X( 1, 2 )
347:                   C( K2, L1 ) = X( 2, 1 )
348:                   C( K2, L2 ) = X( 2, 2 )
349:                END IF
350: *
351:    50       CONTINUE
352: *
353:    60    CONTINUE
354: *
355:       ELSE IF( .NOT.NOTRNA .AND. NOTRNB ) THEN
356: *
357: *        Solve    A' *X + ISGN*X*B = scale*C.
358: *
359: *        The (K,L)th block of X is determined starting from
360: *        upper-left corner column by column by
361: *
362: *          A(K,K)'*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
363: *
364: *        Where
365: *                   K-1                        L-1
366: *          R(K,L) = SUM [A(I,K)'*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)]
367: *                   I=1                        J=1
368: *
369: *        Start column loop (index = L)
370: *        L1 (L2): column index of the first (last) row of X(K,L)
371: *
372:          LNEXT = 1
373:          DO 120 L = 1, N
374:             IF( L.LT.LNEXT )
375:      $         GO TO 120
376:             IF( L.EQ.N ) THEN
377:                L1 = L
378:                L2 = L
379:             ELSE
380:                IF( B( L+1, L ).NE.ZERO ) THEN
381:                   L1 = L
382:                   L2 = L + 1
383:                   LNEXT = L + 2
384:                ELSE
385:                   L1 = L
386:                   L2 = L
387:                   LNEXT = L + 1
388:                END IF
389:             END IF
390: *
391: *           Start row loop (index = K)
392: *           K1 (K2): row index of the first (last) row of X(K,L)
393: *
394:             KNEXT = 1
395:             DO 110 K = 1, M
396:                IF( K.LT.KNEXT )
397:      $            GO TO 110
398:                IF( K.EQ.M ) THEN
399:                   K1 = K
400:                   K2 = K
401:                ELSE
402:                   IF( A( K+1, K ).NE.ZERO ) THEN
403:                      K1 = K
404:                      K2 = K + 1
405:                      KNEXT = K + 2
406:                   ELSE
407:                      K1 = K
408:                      K2 = K
409:                      KNEXT = K + 1
410:                   END IF
411:                END IF
412: *
413:                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
414:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
415:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
416:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
417:                   SCALOC = ONE
418: *
419:                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
420:                   DA11 = ABS( A11 )
421:                   IF( DA11.LE.SMIN ) THEN
422:                      A11 = SMIN
423:                      DA11 = SMIN
424:                      INFO = 1
425:                   END IF
426:                   DB = ABS( VEC( 1, 1 ) )
427:                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
428:                      IF( DB.GT.BIGNUM*DA11 )
429:      $                  SCALOC = ONE / DB
430:                   END IF
431:                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
432: *
433:                   IF( SCALOC.NE.ONE ) THEN
434:                      DO 70 J = 1, N
435:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
436:    70                CONTINUE
437:                      SCALE = SCALE*SCALOC
438:                   END IF
439:                   C( K1, L1 ) = X( 1, 1 )
440: *
441:                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
442: *
443:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
444:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
445:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
446: *
447:                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
448:                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
449:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
450: *
451:                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ),
452:      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
453:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
454:                   IF( IERR.NE.0 )
455:      $               INFO = 1
456: *
457:                   IF( SCALOC.NE.ONE ) THEN
458:                      DO 80 J = 1, N
459:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
460:    80                CONTINUE
461:                      SCALE = SCALE*SCALOC
462:                   END IF
463:                   C( K1, L1 ) = X( 1, 1 )
464:                   C( K2, L1 ) = X( 2, 1 )
465: *
466:                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
467: *
468:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
469:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
470:                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
471: *
472:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
473:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
474:                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
475: *
476:                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ),
477:      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
478:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
479:                   IF( IERR.NE.0 )
480:      $               INFO = 1
481: *
482:                   IF( SCALOC.NE.ONE ) THEN
483:                      DO 90 J = 1, N
484:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
485:    90                CONTINUE
486:                      SCALE = SCALE*SCALOC
487:                   END IF
488:                   C( K1, L1 ) = X( 1, 1 )
489:                   C( K1, L2 ) = X( 2, 1 )
490: *
491:                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
492: *
493:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
494:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 )
495:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
496: *
497:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
498:                   SUMR = DDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 )
499:                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
500: *
501:                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
502:                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 )
503:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
504: *
505:                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 )
506:                   SUMR = DDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 )
507:                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
508: *
509:                   CALL DLASY2( .TRUE., .FALSE., ISGN, 2, 2, A( K1, K1 ),
510:      $                         LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
511:      $                         2, XNORM, IERR )
512:                   IF( IERR.NE.0 )
513:      $               INFO = 1
514: *
515:                   IF( SCALOC.NE.ONE ) THEN
516:                      DO 100 J = 1, N
517:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
518:   100                CONTINUE
519:                      SCALE = SCALE*SCALOC
520:                   END IF
521:                   C( K1, L1 ) = X( 1, 1 )
522:                   C( K1, L2 ) = X( 1, 2 )
523:                   C( K2, L1 ) = X( 2, 1 )
524:                   C( K2, L2 ) = X( 2, 2 )
525:                END IF
526: *
527:   110       CONTINUE
528:   120    CONTINUE
529: *
530:       ELSE IF( .NOT.NOTRNA .AND. .NOT.NOTRNB ) THEN
531: *
532: *        Solve    A'*X + ISGN*X*B' = scale*C.
533: *
534: *        The (K,L)th block of X is determined starting from
535: *        top-right corner column by column by
536: *
537: *           A(K,K)'*X(K,L) + ISGN*X(K,L)*B(L,L)' = C(K,L) - R(K,L)
538: *
539: *        Where
540: *                     K-1                          N
541: *            R(K,L) = SUM [A(I,K)'*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)'].
542: *                     I=1                        J=L+1
543: *
544: *        Start column loop (index = L)
545: *        L1 (L2): column index of the first (last) row of X(K,L)
546: *
547:          LNEXT = N
548:          DO 180 L = N, 1, -1
549:             IF( L.GT.LNEXT )
550:      $         GO TO 180
551:             IF( L.EQ.1 ) THEN
552:                L1 = L
553:                L2 = L
554:             ELSE
555:                IF( B( L, L-1 ).NE.ZERO ) THEN
556:                   L1 = L - 1
557:                   L2 = L
558:                   LNEXT = L - 2
559:                ELSE
560:                   L1 = L
561:                   L2 = L
562:                   LNEXT = L - 1
563:                END IF
564:             END IF
565: *
566: *           Start row loop (index = K)
567: *           K1 (K2): row index of the first (last) row of X(K,L)
568: *
569:             KNEXT = 1
570:             DO 170 K = 1, M
571:                IF( K.LT.KNEXT )
572:      $            GO TO 170
573:                IF( K.EQ.M ) THEN
574:                   K1 = K
575:                   K2 = K
576:                ELSE
577:                   IF( A( K+1, K ).NE.ZERO ) THEN
578:                      K1 = K
579:                      K2 = K + 1
580:                      KNEXT = K + 2
581:                   ELSE
582:                      K1 = K
583:                      K2 = K
584:                      KNEXT = K + 1
585:                   END IF
586:                END IF
587: *
588:                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
589:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
590:                   SUMR = DDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC,
591:      $                   B( L1, MIN( L1+1, N ) ), LDB )
592:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
593:                   SCALOC = ONE
594: *
595:                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
596:                   DA11 = ABS( A11 )
597:                   IF( DA11.LE.SMIN ) THEN
598:                      A11 = SMIN
599:                      DA11 = SMIN
600:                      INFO = 1
601:                   END IF
602:                   DB = ABS( VEC( 1, 1 ) )
603:                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
604:                      IF( DB.GT.BIGNUM*DA11 )
605:      $                  SCALOC = ONE / DB
606:                   END IF
607:                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
608: *
609:                   IF( SCALOC.NE.ONE ) THEN
610:                      DO 130 J = 1, N
611:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
612:   130                CONTINUE
613:                      SCALE = SCALE*SCALOC
614:                   END IF
615:                   C( K1, L1 ) = X( 1, 1 )
616: *
617:                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
618: *
619:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
620:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
621:      $                   B( L1, MIN( L2+1, N ) ), LDB )
622:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
623: *
624:                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
625:                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
626:      $                   B( L1, MIN( L2+1, N ) ), LDB )
627:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
628: *
629:                   CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ),
630:      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
631:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
632:                   IF( IERR.NE.0 )
633:      $               INFO = 1
634: *
635:                   IF( SCALOC.NE.ONE ) THEN
636:                      DO 140 J = 1, N
637:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
638:   140                CONTINUE
639:                      SCALE = SCALE*SCALOC
640:                   END IF
641:                   C( K1, L1 ) = X( 1, 1 )
642:                   C( K2, L1 ) = X( 2, 1 )
643: *
644:                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
645: *
646:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
647:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
648:      $                   B( L1, MIN( L2+1, N ) ), LDB )
649:                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
650: *
651:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
652:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
653:      $                   B( L2, MIN( L2+1, N ) ), LDB )
654:                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
655: *
656:                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ),
657:      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
658:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
659:                   IF( IERR.NE.0 )
660:      $               INFO = 1
661: *
662:                   IF( SCALOC.NE.ONE ) THEN
663:                      DO 150 J = 1, N
664:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
665:   150                CONTINUE
666:                      SCALE = SCALE*SCALOC
667:                   END IF
668:                   C( K1, L1 ) = X( 1, 1 )
669:                   C( K1, L2 ) = X( 2, 1 )
670: *
671:                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
672: *
673:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
674:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
675:      $                   B( L1, MIN( L2+1, N ) ), LDB )
676:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
677: *
678:                   SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
679:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
680:      $                   B( L2, MIN( L2+1, N ) ), LDB )
681:                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
682: *
683:                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
684:                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
685:      $                   B( L1, MIN( L2+1, N ) ), LDB )
686:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
687: *
688:                   SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 )
689:                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
690:      $                   B( L2, MIN( L2+1, N ) ), LDB )
691:                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
692: *
693:                   CALL DLASY2( .TRUE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
694:      $                         LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
695:      $                         2, XNORM, IERR )
696:                   IF( IERR.NE.0 )
697:      $               INFO = 1
698: *
699:                   IF( SCALOC.NE.ONE ) THEN
700:                      DO 160 J = 1, N
701:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
702:   160                CONTINUE
703:                      SCALE = SCALE*SCALOC
704:                   END IF
705:                   C( K1, L1 ) = X( 1, 1 )
706:                   C( K1, L2 ) = X( 1, 2 )
707:                   C( K2, L1 ) = X( 2, 1 )
708:                   C( K2, L2 ) = X( 2, 2 )
709:                END IF
710: *
711:   170       CONTINUE
712:   180    CONTINUE
713: *
714:       ELSE IF( NOTRNA .AND. .NOT.NOTRNB ) THEN
715: *
716: *        Solve    A*X + ISGN*X*B' = scale*C.
717: *
718: *        The (K,L)th block of X is determined starting from
719: *        bottom-right corner column by column by
720: *
721: *            A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L)' = C(K,L) - R(K,L)
722: *
723: *        Where
724: *                      M                          N
725: *            R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)'].
726: *                    I=K+1                      J=L+1
727: *
728: *        Start column loop (index = L)
729: *        L1 (L2): column index of the first (last) row of X(K,L)
730: *
731:          LNEXT = N
732:          DO 240 L = N, 1, -1
733:             IF( L.GT.LNEXT )
734:      $         GO TO 240
735:             IF( L.EQ.1 ) THEN
736:                L1 = L
737:                L2 = L
738:             ELSE
739:                IF( B( L, L-1 ).NE.ZERO ) THEN
740:                   L1 = L - 1
741:                   L2 = L
742:                   LNEXT = L - 2
743:                ELSE
744:                   L1 = L
745:                   L2 = L
746:                   LNEXT = L - 1
747:                END IF
748:             END IF
749: *
750: *           Start row loop (index = K)
751: *           K1 (K2): row index of the first (last) row of X(K,L)
752: *
753:             KNEXT = M
754:             DO 230 K = M, 1, -1
755:                IF( K.GT.KNEXT )
756:      $            GO TO 230
757:                IF( K.EQ.1 ) THEN
758:                   K1 = K
759:                   K2 = K
760:                ELSE
761:                   IF( A( K, K-1 ).NE.ZERO ) THEN
762:                      K1 = K - 1
763:                      K2 = K
764:                      KNEXT = K - 2
765:                   ELSE
766:                      K1 = K
767:                      K2 = K
768:                      KNEXT = K - 1
769:                   END IF
770:                END IF
771: *
772:                IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
773:                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
774:      $                   C( MIN( K1+1, M ), L1 ), 1 )
775:                   SUMR = DDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC,
776:      $                   B( L1, MIN( L1+1, N ) ), LDB )
777:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
778:                   SCALOC = ONE
779: *
780:                   A11 = A( K1, K1 ) + SGN*B( L1, L1 )
781:                   DA11 = ABS( A11 )
782:                   IF( DA11.LE.SMIN ) THEN
783:                      A11 = SMIN
784:                      DA11 = SMIN
785:                      INFO = 1
786:                   END IF
787:                   DB = ABS( VEC( 1, 1 ) )
788:                   IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
789:                      IF( DB.GT.BIGNUM*DA11 )
790:      $                  SCALOC = ONE / DB
791:                   END IF
792:                   X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
793: *
794:                   IF( SCALOC.NE.ONE ) THEN
795:                      DO 190 J = 1, N
796:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
797:   190                CONTINUE
798:                      SCALE = SCALE*SCALOC
799:                   END IF
800:                   C( K1, L1 ) = X( 1, 1 )
801: *
802:                ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
803: *
804:                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
805:      $                   C( MIN( K2+1, M ), L1 ), 1 )
806:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
807:      $                   B( L1, MIN( L2+1, N ) ), LDB )
808:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
809: *
810:                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
811:      $                   C( MIN( K2+1, M ), L1 ), 1 )
812:                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
813:      $                   B( L1, MIN( L2+1, N ) ), LDB )
814:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
815: *
816:                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ),
817:      $                         LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
818:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
819:                   IF( IERR.NE.0 )
820:      $               INFO = 1
821: *
822:                   IF( SCALOC.NE.ONE ) THEN
823:                      DO 200 J = 1, N
824:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
825:   200                CONTINUE
826:                      SCALE = SCALE*SCALOC
827:                   END IF
828:                   C( K1, L1 ) = X( 1, 1 )
829:                   C( K2, L1 ) = X( 2, 1 )
830: *
831:                ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
832: *
833:                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
834:      $                   C( MIN( K1+1, M ), L1 ), 1 )
835:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
836:      $                   B( L1, MIN( L2+1, N ) ), LDB )
837:                   VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
838: *
839:                   SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
840:      $                   C( MIN( K1+1, M ), L2 ), 1 )
841:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
842:      $                   B( L2, MIN( L2+1, N ) ), LDB )
843:                   VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
844: *
845:                   CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ),
846:      $                         LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
847:      $                         ZERO, X, 2, SCALOC, XNORM, IERR )
848:                   IF( IERR.NE.0 )
849:      $               INFO = 1
850: *
851:                   IF( SCALOC.NE.ONE ) THEN
852:                      DO 210 J = 1, N
853:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
854:   210                CONTINUE
855:                      SCALE = SCALE*SCALOC
856:                   END IF
857:                   C( K1, L1 ) = X( 1, 1 )
858:                   C( K1, L2 ) = X( 2, 1 )
859: *
860:                ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
861: *
862:                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
863:      $                   C( MIN( K2+1, M ), L1 ), 1 )
864:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
865:      $                   B( L1, MIN( L2+1, N ) ), LDB )
866:                   VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
867: *
868:                   SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
869:      $                   C( MIN( K2+1, M ), L2 ), 1 )
870:                   SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
871:      $                   B( L2, MIN( L2+1, N ) ), LDB )
872:                   VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
873: *
874:                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
875:      $                   C( MIN( K2+1, M ), L1 ), 1 )
876:                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
877:      $                   B( L1, MIN( L2+1, N ) ), LDB )
878:                   VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
879: *
880:                   SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
881:      $                   C( MIN( K2+1, M ), L2 ), 1 )
882:                   SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
883:      $                   B( L2, MIN( L2+1, N ) ), LDB )
884:                   VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
885: *
886:                   CALL DLASY2( .FALSE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
887:      $                         LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
888:      $                         2, XNORM, IERR )
889:                   IF( IERR.NE.0 )
890:      $               INFO = 1
891: *
892:                   IF( SCALOC.NE.ONE ) THEN
893:                      DO 220 J = 1, N
894:                         CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
895:   220                CONTINUE
896:                      SCALE = SCALE*SCALOC
897:                   END IF
898:                   C( K1, L1 ) = X( 1, 1 )
899:                   C( K1, L2 ) = X( 1, 2 )
900:                   C( K2, L1 ) = X( 2, 1 )
901:                   C( K2, L2 ) = X( 2, 2 )
902:                END IF
903: *
904:   230       CONTINUE
905:   240    CONTINUE
906: *
907:       END IF
908: *
909:       RETURN
910: *
911: *     End of DTRSYL
912: *
913:       END
914: