1: *> \brief \b ZTRSYL3
2: *
3: * Definition:
4: * ===========
5: *
6: *
7: *> \par Purpose
8: * =============
9: *>
10: *> \verbatim
11: *>
12: *> ZTRSYL3 solves the complex Sylvester matrix equation:
13: *>
14: *> op(A)*X + X*op(B) = scale*C or
15: *> op(A)*X - X*op(B) = scale*C,
16: *>
17: *> where op(A) = A or A**H, and A and B are both upper triangular. A is
18: *> M-by-M and B is N-by-N; the right hand side C and the solution X are
19: *> M-by-N; and scale is an output scale factor, set <= 1 to avoid
20: *> overflow in X.
21: *>
22: *> This is the block version of the algorithm.
23: *> \endverbatim
24: *
25: * Arguments
26: * =========
27: *
28: *> \param[in] TRANA
29: *> \verbatim
30: *> TRANA is CHARACTER*1
31: *> Specifies the option op(A):
32: *> = 'N': op(A) = A (No transpose)
33: *> = 'C': op(A) = A**H (Conjugate transpose)
34: *> \endverbatim
35: *>
36: *> \param[in] TRANB
37: *> \verbatim
38: *> TRANB is CHARACTER*1
39: *> Specifies the option op(B):
40: *> = 'N': op(B) = B (No transpose)
41: *> = 'C': op(B) = B**H (Conjugate transpose)
42: *> \endverbatim
43: *>
44: *> \param[in] ISGN
45: *> \verbatim
46: *> ISGN is INTEGER
47: *> Specifies the sign in the equation:
48: *> = +1: solve op(A)*X + X*op(B) = scale*C
49: *> = -1: solve op(A)*X - X*op(B) = scale*C
50: *> \endverbatim
51: *>
52: *> \param[in] M
53: *> \verbatim
54: *> M is INTEGER
55: *> The order of the matrix A, and the number of rows in the
56: *> matrices X and C. M >= 0.
57: *> \endverbatim
58: *>
59: *> \param[in] N
60: *> \verbatim
61: *> N is INTEGER
62: *> The order of the matrix B, and the number of columns in the
63: *> matrices X and C. N >= 0.
64: *> \endverbatim
65: *>
66: *> \param[in] A
67: *> \verbatim
68: *> A is COMPLEX*16 array, dimension (LDA,M)
69: *> The upper triangular matrix A.
70: *> \endverbatim
71: *>
72: *> \param[in] LDA
73: *> \verbatim
74: *> LDA is INTEGER
75: *> The leading dimension of the array A. LDA >= max(1,M).
76: *> \endverbatim
77: *>
78: *> \param[in] B
79: *> \verbatim
80: *> B is COMPLEX*16 array, dimension (LDB,N)
81: *> The upper triangular matrix B.
82: *> \endverbatim
83: *>
84: *> \param[in] LDB
85: *> \verbatim
86: *> LDB is INTEGER
87: *> The leading dimension of the array B. LDB >= max(1,N).
88: *> \endverbatim
89: *>
90: *> \param[in,out] C
91: *> \verbatim
92: *> C is COMPLEX*16 array, dimension (LDC,N)
93: *> On entry, the M-by-N right hand side matrix C.
94: *> On exit, C is overwritten by the solution matrix X.
95: *> \endverbatim
96: *>
97: *> \param[in] LDC
98: *> \verbatim
99: *> LDC is INTEGER
100: *> The leading dimension of the array C. LDC >= max(1,M)
101: *> \endverbatim
102: *>
103: *> \param[out] SCALE
104: *> \verbatim
105: *> SCALE is DOUBLE PRECISION
106: *> The scale factor, scale, set <= 1 to avoid overflow in X.
107: *> \endverbatim
108: *>
109: *> \param[out] SWORK
110: *> \verbatim
111: *> SWORK is DOUBLE PRECISION array, dimension (MAX(2, ROWS),
112: *> MAX(1,COLS)).
113: *> On exit, if INFO = 0, SWORK(1) returns the optimal value ROWS
114: *> and SWORK(2) returns the optimal COLS.
115: *> \endverbatim
116: *>
117: *> \param[in] LDSWORK
118: *> \verbatim
119: *> LDSWORK is INTEGER
120: *> LDSWORK >= MAX(2,ROWS), where ROWS = ((M + NB - 1) / NB + 1)
121: *> and NB is the optimal block size.
122: *>
123: *> If LDSWORK = -1, then a workspace query is assumed; the routine
124: *> only calculates the optimal dimensions of the SWORK matrix,
125: *> returns these values as the first and second entry of the SWORK
126: *> matrix, and no error message related LWORK is issued by XERBLA.
127: *> \endverbatim
128: *>
129: *> \param[out] INFO
130: *> \verbatim
131: *> INFO is INTEGER
132: *> = 0: successful exit
133: *> < 0: if INFO = -i, the i-th argument had an illegal value
134: *> = 1: A and B have common or very close eigenvalues; perturbed
135: *> values were used to solve the equation (but the matrices
136: *> A and B are unchanged).
137: *> \endverbatim
138: *
139: *> \ingroup complex16SYcomputational
140: *
141: * =====================================================================
142: * References:
143: * E. S. Quintana-Orti and R. A. Van De Geijn (2003). Formal derivation of
144: * algorithms: The triangular Sylvester equation, ACM Transactions
145: * on Mathematical Software (TOMS), volume 29, pages 218--243.
146: *
147: * A. Schwarz and C. C. Kjelgaard Mikkelsen (2020). Robust Task-Parallel
148: * Solution of the Triangular Sylvester Equation. Lecture Notes in
149: * Computer Science, vol 12043, pages 82--92, Springer.
150: *
151: * Contributor:
152: * Angelika Schwarz, Umea University, Sweden.
153: *
154: * =====================================================================
155: SUBROUTINE ZTRSYL3( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C,
156: $ LDC, SCALE, SWORK, LDSWORK, INFO )
157: IMPLICIT NONE
158: *
159: * .. Scalar Arguments ..
160: CHARACTER TRANA, TRANB
161: INTEGER INFO, ISGN, LDA, LDB, LDC, LDSWORK, M, N
162: DOUBLE PRECISION SCALE
163: * ..
164: * .. Array Arguments ..
165: COMPLEX*16 A( LDA, * ), B( LDB, * ), C( LDC, * )
166: DOUBLE PRECISION SWORK( LDSWORK, * )
167: * ..
168: * .. Parameters ..
169: DOUBLE PRECISION ZERO, ONE
170: PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
171: COMPLEX*16 CONE
172: PARAMETER ( CONE = ( 1.0D0, 0.0D0 ) )
173: * ..
174: * .. Local Scalars ..
175: LOGICAL NOTRNA, NOTRNB, LQUERY
176: INTEGER AWRK, BWRK, I, I1, I2, IINFO, J, J1, J2, JJ,
177: $ K, K1, K2, L, L1, L2, LL, NBA, NB, NBB
178: DOUBLE PRECISION ANRM, BIGNUM, BNRM, CNRM, SCAL, SCALOC,
179: $ SCAMIN, SGN, XNRM, BUF, SMLNUM
180: COMPLEX*16 CSGN
181: * ..
182: * .. Local Arrays ..
183: DOUBLE PRECISION WNRM( MAX( M, N ) )
184: * ..
185: * .. External Functions ..
186: LOGICAL LSAME
187: INTEGER ILAENV
188: DOUBLE PRECISION DLAMCH, DLARMM, ZLANGE
189: EXTERNAL DLAMCH, DLARMM, ILAENV, LSAME, ZLANGE
190: * ..
191: * .. External Subroutines ..
192: EXTERNAL XERBLA, ZDSCAL, ZGEMM, ZLASCL, ZTRSYL
193: * ..
194: * .. Intrinsic Functions ..
195: INTRINSIC ABS, DBLE, DIMAG, EXPONENT, MAX, MIN
196: * ..
197: * .. Executable Statements ..
198: *
199: * Decode and Test input parameters
200: *
201: NOTRNA = LSAME( TRANA, 'N' )
202: NOTRNB = LSAME( TRANB, 'N' )
203: *
204: * Use the same block size for all matrices.
205: *
206: NB = MAX( 8, ILAENV( 1, 'ZTRSYL', '', M, N, -1, -1) )
207: *
208: * Compute number of blocks in A and B
209: *
210: NBA = MAX( 1, (M + NB - 1) / NB )
211: NBB = MAX( 1, (N + NB - 1) / NB )
212: *
213: * Compute workspace
214: *
215: INFO = 0
216: LQUERY = ( LDSWORK.EQ.-1 )
217: IF( LQUERY ) THEN
218: LDSWORK = 2
219: SWORK(1,1) = MAX( NBA, NBB )
220: SWORK(2,1) = 2 * NBB + NBA
221: END IF
222: *
223: * Test the input arguments
224: *
225: IF( .NOT.NOTRNA .AND. .NOT. LSAME( TRANA, 'C' ) ) THEN
226: INFO = -1
227: ELSE IF( .NOT.NOTRNB .AND. .NOT. LSAME( TRANB, 'C' ) ) THEN
228: INFO = -2
229: ELSE IF( ISGN.NE.1 .AND. ISGN.NE.-1 ) THEN
230: INFO = -3
231: ELSE IF( M.LT.0 ) THEN
232: INFO = -4
233: ELSE IF( N.LT.0 ) THEN
234: INFO = -5
235: ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
236: INFO = -7
237: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
238: INFO = -9
239: ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
240: INFO = -11
241: END IF
242: IF( INFO.NE.0 ) THEN
243: CALL XERBLA( 'ZTRSYL3', -INFO )
244: RETURN
245: ELSE IF( LQUERY ) THEN
246: RETURN
247: END IF
248: *
249: * Quick return if possible
250: *
251: SCALE = ONE
252: IF( M.EQ.0 .OR. N.EQ.0 )
253: $ RETURN
254: *
255: * Use unblocked code for small problems or if insufficient
256: * workspace is provided
257: *
258: IF( MIN( NBA, NBB ).EQ.1 .OR. LDSWORK.LT.MAX( NBA, NBB ) ) THEN
259: CALL ZTRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB,
260: $ C, LDC, SCALE, INFO )
261: RETURN
262: END IF
263: *
264: * Set constants to control overflow
265: *
266: SMLNUM = DLAMCH( 'S' )
267: BIGNUM = ONE / SMLNUM
268: *
269: * Set local scaling factors.
270: *
271: DO L = 1, NBB
272: DO K = 1, NBA
273: SWORK( K, L ) = ONE
274: END DO
275: END DO
276: *
277: * Fallback scaling factor to prevent flushing of SWORK( K, L ) to zero.
278: * This scaling is to ensure compatibility with TRSYL and may get flushed.
279: *
280: BUF = ONE
281: *
282: * Compute upper bounds of blocks of A and B
283: *
284: AWRK = NBB
285: DO K = 1, NBA
286: K1 = (K - 1) * NB + 1
287: K2 = MIN( K * NB, M ) + 1
288: DO L = K, NBA
289: L1 = (L - 1) * NB + 1
290: L2 = MIN( L * NB, M ) + 1
291: IF( NOTRNA ) THEN
292: SWORK( K, AWRK + L ) = ZLANGE( 'I', K2-K1, L2-L1,
293: $ A( K1, L1 ), LDA, WNRM )
294: ELSE
295: SWORK( L, AWRK + K ) = ZLANGE( '1', K2-K1, L2-L1,
296: $ A( K1, L1 ), LDA, WNRM )
297: END IF
298: END DO
299: END DO
300: BWRK = NBB + NBA
301: DO K = 1, NBB
302: K1 = (K - 1) * NB + 1
303: K2 = MIN( K * NB, N ) + 1
304: DO L = K, NBB
305: L1 = (L - 1) * NB + 1
306: L2 = MIN( L * NB, N ) + 1
307: IF( NOTRNB ) THEN
308: SWORK( K, BWRK + L ) = ZLANGE( 'I', K2-K1, L2-L1,
309: $ B( K1, L1 ), LDB, WNRM )
310: ELSE
311: SWORK( L, BWRK + K ) = ZLANGE( '1', K2-K1, L2-L1,
312: $ B( K1, L1 ), LDB, WNRM )
313: END IF
314: END DO
315: END DO
316: *
317: SGN = DBLE( ISGN )
318: CSGN = DCMPLX( SGN, ZERO )
319: *
320: IF( NOTRNA .AND. NOTRNB ) THEN
321: *
322: * Solve A*X + ISGN*X*B = scale*C.
323: *
324: * The (K,L)th block of X is determined starting from
325: * bottom-left corner column by column by
326: *
327: * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
328: *
329: * Where
330: * M L-1
331: * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)].
332: * I=K+1 J=1
333: *
334: * Start loop over block rows (index = K) and block columns (index = L)
335: *
336: DO K = NBA, 1, -1
337: *
338: * K1: row index of the first row in X( K, L )
339: * K2: row index of the first row in X( K+1, L )
340: * so the K2 - K1 is the column count of the block X( K, L )
341: *
342: K1 = (K - 1) * NB + 1
343: K2 = MIN( K * NB, M ) + 1
344: DO L = 1, NBB
345: *
346: * L1: column index of the first column in X( K, L )
347: * L2: column index of the first column in X( K, L + 1)
348: * so that L2 - L1 is the row count of the block X( K, L )
349: *
350: L1 = (L - 1) * NB + 1
351: L2 = MIN( L * NB, N ) + 1
352: *
353: CALL ZTRSYL( TRANA, TRANB, ISGN, K2-K1, L2-L1,
354: $ A( K1, K1 ), LDA,
355: $ B( L1, L1 ), LDB,
356: $ C( K1, L1 ), LDC, SCALOC, IINFO )
357: INFO = MAX( INFO, IINFO )
358: *
359: IF( SCALOC * SWORK( K, L ) .EQ. ZERO ) THEN
360: IF( SCALOC .EQ. ZERO ) THEN
361: * The magnitude of the largest entry of X(K1:K2-1, L1:L2-1)
362: * is larger than the product of BIGNUM**2 and cannot be
363: * represented in the form (1/SCALE)*X(K1:K2-1, L1:L2-1).
364: * Mark the computation as pointless.
365: BUF = ZERO
366: ELSE
367: BUF = BUF*2.D0**EXPONENT( SCALOC )
368: END IF
369: DO JJ = 1, NBB
370: DO LL = 1, NBA
371: * Bound by BIGNUM to not introduce Inf. The value
372: * is irrelevant; corresponding entries of the
373: * solution will be flushed in consistency scaling.
374: SWORK( LL, JJ ) = MIN( BIGNUM,
375: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
376: END DO
377: END DO
378: END IF
379: SWORK( K, L ) = SCALOC * SWORK( K, L )
380: XNRM = ZLANGE( 'I', K2-K1, L2-L1, C( K1, L1 ), LDC,
381: $ WNRM )
382: *
383: DO I = K - 1, 1, -1
384: *
385: * C( I, L ) := C( I, L ) - A( I, K ) * C( K, L )
386: *
387: I1 = (I - 1) * NB + 1
388: I2 = MIN( I * NB, M ) + 1
389: *
390: * Compute scaling factor to survive the linear update
391: * simulating consistent scaling.
392: *
393: CNRM = ZLANGE( 'I', I2-I1, L2-L1, C( I1, L1 ),
394: $ LDC, WNRM )
395: SCAMIN = MIN( SWORK( I, L ), SWORK( K, L ) )
396: CNRM = CNRM * ( SCAMIN / SWORK( I, L ) )
397: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
398: ANRM = SWORK( I, AWRK + K )
399: SCALOC = DLARMM( ANRM, XNRM, CNRM )
400: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
401: * Use second scaling factor to prevent flushing to zero.
402: BUF = BUF*2.D0**EXPONENT( SCALOC )
403: DO JJ = 1, NBB
404: DO LL = 1, NBA
405: SWORK( LL, JJ ) = MIN( BIGNUM,
406: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
407: END DO
408: END DO
409: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
410: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
411: END IF
412: CNRM = CNRM * SCALOC
413: XNRM = XNRM * SCALOC
414: *
415: * Simultaneously apply the robust update factor and the
416: * consistency scaling factor to C( I, L ) and C( K, L ).
417: *
418: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
419: IF( SCAL .NE. ONE ) THEN
420: DO JJ = L1, L2-1
421: CALL ZDSCAL( K2-K1, SCAL, C( K1, JJ ), 1)
422: END DO
423: ENDIF
424: *
425: SCAL = ( SCAMIN / SWORK( I, L ) ) * SCALOC
426: IF( SCAL .NE. ONE ) THEN
427: DO LL = L1, L2-1
428: CALL ZDSCAL( I2-I1, SCAL, C( I1, LL ), 1)
429: END DO
430: ENDIF
431: *
432: * Record current scaling factor
433: *
434: SWORK( K, L ) = SCAMIN * SCALOC
435: SWORK( I, L ) = SCAMIN * SCALOC
436: *
437: CALL ZGEMM( 'N', 'N', I2-I1, L2-L1, K2-K1, -CONE,
438: $ A( I1, K1 ), LDA, C( K1, L1 ), LDC,
439: $ CONE, C( I1, L1 ), LDC )
440: *
441: END DO
442: *
443: DO J = L + 1, NBB
444: *
445: * C( K, J ) := C( K, J ) - SGN * C( K, L ) * B( L, J )
446: *
447: J1 = (J - 1) * NB + 1
448: J2 = MIN( J * NB, N ) + 1
449: *
450: * Compute scaling factor to survive the linear update
451: * simulating consistent scaling.
452: *
453: CNRM = ZLANGE( 'I', K2-K1, J2-J1, C( K1, J1 ),
454: $ LDC, WNRM )
455: SCAMIN = MIN( SWORK( K, J ), SWORK( K, L ) )
456: CNRM = CNRM * ( SCAMIN / SWORK( K, J ) )
457: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
458: BNRM = SWORK(L, BWRK + J)
459: SCALOC = DLARMM( BNRM, XNRM, CNRM )
460: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
461: * Use second scaling factor to prevent flushing to zero.
462: BUF = BUF*2.D0**EXPONENT( SCALOC )
463: DO JJ = 1, NBB
464: DO LL = 1, NBA
465: SWORK( LL, JJ ) = MIN( BIGNUM,
466: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
467: END DO
468: END DO
469: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
470: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
471: END IF
472: CNRM = CNRM * SCALOC
473: XNRM = XNRM * SCALOC
474: *
475: * Simultaneously apply the robust update factor and the
476: * consistency scaling factor to C( K, J ) and C( K, L).
477: *
478: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
479: IF( SCAL .NE. ONE ) THEN
480: DO LL = L1, L2-1
481: CALL ZDSCAL( K2-K1, SCAL, C( K1, LL ), 1 )
482: END DO
483: ENDIF
484: *
485: SCAL = ( SCAMIN / SWORK( K, J ) ) * SCALOC
486: IF( SCAL .NE. ONE ) THEN
487: DO JJ = J1, J2-1
488: CALL ZDSCAL( K2-K1, SCAL, C( K1, JJ ), 1 )
489: END DO
490: ENDIF
491: *
492: * Record current scaling factor
493: *
494: SWORK( K, L ) = SCAMIN * SCALOC
495: SWORK( K, J ) = SCAMIN * SCALOC
496: *
497: CALL ZGEMM( 'N', 'N', K2-K1, J2-J1, L2-L1, -CSGN,
498: $ C( K1, L1 ), LDC, B( L1, J1 ), LDB,
499: $ CONE, C( K1, J1 ), LDC )
500: END DO
501: END DO
502: END DO
503: ELSE IF( .NOT.NOTRNA .AND. NOTRNB ) THEN
504: *
505: * Solve A**H *X + ISGN*X*B = scale*C.
506: *
507: * The (K,L)th block of X is determined starting from
508: * upper-left corner column by column by
509: *
510: * A(K,K)**H*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
511: *
512: * Where
513: * K-1 L-1
514: * R(K,L) = SUM [A(I,K)**H*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)]
515: * I=1 J=1
516: *
517: * Start loop over block rows (index = K) and block columns (index = L)
518: *
519: DO K = 1, NBA
520: *
521: * K1: row index of the first row in X( K, L )
522: * K2: row index of the first row in X( K+1, L )
523: * so the K2 - K1 is the column count of the block X( K, L )
524: *
525: K1 = (K - 1) * NB + 1
526: K2 = MIN( K * NB, M ) + 1
527: DO L = 1, NBB
528: *
529: * L1: column index of the first column in X( K, L )
530: * L2: column index of the first column in X( K, L + 1)
531: * so that L2 - L1 is the row count of the block X( K, L )
532: *
533: L1 = (L - 1) * NB + 1
534: L2 = MIN( L * NB, N ) + 1
535: *
536: CALL ZTRSYL( TRANA, TRANB, ISGN, K2-K1, L2-L1,
537: $ A( K1, K1 ), LDA,
538: $ B( L1, L1 ), LDB,
539: $ C( K1, L1 ), LDC, SCALOC, IINFO )
540: INFO = MAX( INFO, IINFO )
541: *
542: IF( SCALOC * SWORK( K, L ) .EQ. ZERO ) THEN
543: IF( SCALOC .EQ. ZERO ) THEN
544: * The magnitude of the largest entry of X(K1:K2-1, L1:L2-1)
545: * is larger than the product of BIGNUM**2 and cannot be
546: * represented in the form (1/SCALE)*X(K1:K2-1, L1:L2-1).
547: * Mark the computation as pointless.
548: BUF = ZERO
549: ELSE
550: * Use second scaling factor to prevent flushing to zero.
551: BUF = BUF*2.D0**EXPONENT( SCALOC )
552: END IF
553: DO JJ = 1, NBB
554: DO LL = 1, NBA
555: * Bound by BIGNUM to not introduce Inf. The value
556: * is irrelevant; corresponding entries of the
557: * solution will be flushed in consistency scaling.
558: SWORK( LL, JJ ) = MIN( BIGNUM,
559: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
560: END DO
561: END DO
562: END IF
563: SWORK( K, L ) = SCALOC * SWORK( K, L )
564: XNRM = ZLANGE( 'I', K2-K1, L2-L1, C( K1, L1 ), LDC,
565: $ WNRM )
566: *
567: DO I = K + 1, NBA
568: *
569: * C( I, L ) := C( I, L ) - A( K, I )**H * C( K, L )
570: *
571: I1 = (I - 1) * NB + 1
572: I2 = MIN( I * NB, M ) + 1
573: *
574: * Compute scaling factor to survive the linear update
575: * simulating consistent scaling.
576: *
577: CNRM = ZLANGE( 'I', I2-I1, L2-L1, C( I1, L1 ),
578: $ LDC, WNRM )
579: SCAMIN = MIN( SWORK( I, L ), SWORK( K, L ) )
580: CNRM = CNRM * ( SCAMIN / SWORK( I, L ) )
581: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
582: ANRM = SWORK( I, AWRK + K )
583: SCALOC = DLARMM( ANRM, XNRM, CNRM )
584: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
585: * Use second scaling factor to prevent flushing to zero.
586: BUF = BUF*2.D0**EXPONENT( SCALOC )
587: DO JJ = 1, NBB
588: DO LL = 1, NBA
589: SWORK( LL, JJ ) = MIN( BIGNUM,
590: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
591: END DO
592: END DO
593: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
594: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
595: END IF
596: CNRM = CNRM * SCALOC
597: XNRM = XNRM * SCALOC
598: *
599: * Simultaneously apply the robust update factor and the
600: * consistency scaling factor to to C( I, L ) and C( K, L).
601: *
602: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
603: IF( SCAL .NE. ONE ) THEN
604: DO LL = L1, L2-1
605: CALL ZDSCAL( K2-K1, SCAL, C( K1, LL ), 1 )
606: END DO
607: ENDIF
608: *
609: SCAL = ( SCAMIN / SWORK( I, L ) ) * SCALOC
610: IF( SCAL .NE. ONE ) THEN
611: DO LL = L1, L2-1
612: CALL ZDSCAL( I2-I1, SCAL, C( I1, LL ), 1 )
613: END DO
614: ENDIF
615: *
616: * Record current scaling factor
617: *
618: SWORK( K, L ) = SCAMIN * SCALOC
619: SWORK( I, L ) = SCAMIN * SCALOC
620: *
621: CALL ZGEMM( 'C', 'N', I2-I1, L2-L1, K2-K1, -CONE,
622: $ A( K1, I1 ), LDA, C( K1, L1 ), LDC,
623: $ CONE, C( I1, L1 ), LDC )
624: END DO
625: *
626: DO J = L + 1, NBB
627: *
628: * C( K, J ) := C( K, J ) - SGN * C( K, L ) * B( L, J )
629: *
630: J1 = (J - 1) * NB + 1
631: J2 = MIN( J * NB, N ) + 1
632: *
633: * Compute scaling factor to survive the linear update
634: * simulating consistent scaling.
635: *
636: CNRM = ZLANGE( 'I', K2-K1, J2-J1, C( K1, J1 ),
637: $ LDC, WNRM )
638: SCAMIN = MIN( SWORK( K, J ), SWORK( K, L ) )
639: CNRM = CNRM * ( SCAMIN / SWORK( K, J ) )
640: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
641: BNRM = SWORK( L, BWRK + J )
642: SCALOC = DLARMM( BNRM, XNRM, CNRM )
643: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
644: * Use second scaling factor to prevent flushing to zero.
645: BUF = BUF*2.D0**EXPONENT( SCALOC )
646: DO JJ = 1, NBB
647: DO LL = 1, NBA
648: SWORK( LL, JJ ) = MIN( BIGNUM,
649: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
650: END DO
651: END DO
652: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
653: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
654: END IF
655: CNRM = CNRM * SCALOC
656: XNRM = XNRM * SCALOC
657: *
658: * Simultaneously apply the robust update factor and the
659: * consistency scaling factor to to C( K, J ) and C( K, L).
660: *
661: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
662: IF( SCAL .NE. ONE ) THEN
663: DO LL = L1, L2-1
664: CALL ZDSCAL( K2-K1, SCAL, C( K1, LL ), 1 )
665: END DO
666: ENDIF
667: *
668: SCAL = ( SCAMIN / SWORK( K, J ) ) * SCALOC
669: IF( SCAL .NE. ONE ) THEN
670: DO JJ = J1, J2-1
671: CALL ZDSCAL( K2-K1, SCAL, C( K1, JJ ), 1 )
672: END DO
673: ENDIF
674: *
675: * Record current scaling factor
676: *
677: SWORK( K, L ) = SCAMIN * SCALOC
678: SWORK( K, J ) = SCAMIN * SCALOC
679: *
680: CALL ZGEMM( 'N', 'N', K2-K1, J2-J1, L2-L1, -CSGN,
681: $ C( K1, L1 ), LDC, B( L1, J1 ), LDB,
682: $ CONE, C( K1, J1 ), LDC )
683: END DO
684: END DO
685: END DO
686: ELSE IF( .NOT.NOTRNA .AND. .NOT.NOTRNB ) THEN
687: *
688: * Solve A**H *X + ISGN*X*B**H = scale*C.
689: *
690: * The (K,L)th block of X is determined starting from
691: * top-right corner column by column by
692: *
693: * A(K,K)**H*X(K,L) + ISGN*X(K,L)*B(L,L)**H = C(K,L) - R(K,L)
694: *
695: * Where
696: * K-1 N
697: * R(K,L) = SUM [A(I,K)**H*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**H].
698: * I=1 J=L+1
699: *
700: * Start loop over block rows (index = K) and block columns (index = L)
701: *
702: DO K = 1, NBA
703: *
704: * K1: row index of the first row in X( K, L )
705: * K2: row index of the first row in X( K+1, L )
706: * so the K2 - K1 is the column count of the block X( K, L )
707: *
708: K1 = (K - 1) * NB + 1
709: K2 = MIN( K * NB, M ) + 1
710: DO L = NBB, 1, -1
711: *
712: * L1: column index of the first column in X( K, L )
713: * L2: column index of the first column in X( K, L + 1)
714: * so that L2 - L1 is the row count of the block X( K, L )
715: *
716: L1 = (L - 1) * NB + 1
717: L2 = MIN( L * NB, N ) + 1
718: *
719: CALL ZTRSYL( TRANA, TRANB, ISGN, K2-K1, L2-L1,
720: $ A( K1, K1 ), LDA,
721: $ B( L1, L1 ), LDB,
722: $ C( K1, L1 ), LDC, SCALOC, IINFO )
723: INFO = MAX( INFO, IINFO )
724: *
725: IF( SCALOC * SWORK( K, L ) .EQ. ZERO ) THEN
726: IF( SCALOC .EQ. ZERO ) THEN
727: * The magnitude of the largest entry of X(K1:K2-1, L1:L2-1)
728: * is larger than the product of BIGNUM**2 and cannot be
729: * represented in the form (1/SCALE)*X(K1:K2-1, L1:L2-1).
730: * Mark the computation as pointless.
731: BUF = ZERO
732: ELSE
733: * Use second scaling factor to prevent flushing to zero.
734: BUF = BUF*2.D0**EXPONENT( SCALOC )
735: END IF
736: DO JJ = 1, NBB
737: DO LL = 1, NBA
738: * Bound by BIGNUM to not introduce Inf. The value
739: * is irrelevant; corresponding entries of the
740: * solution will be flushed in consistency scaling.
741: SWORK( LL, JJ ) = MIN( BIGNUM,
742: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
743: END DO
744: END DO
745: END IF
746: SWORK( K, L ) = SCALOC * SWORK( K, L )
747: XNRM = ZLANGE( 'I', K2-K1, L2-L1, C( K1, L1 ), LDC,
748: $ WNRM )
749: *
750: DO I = K + 1, NBA
751: *
752: * C( I, L ) := C( I, L ) - A( K, I )**H * C( K, L )
753: *
754: I1 = (I - 1) * NB + 1
755: I2 = MIN( I * NB, M ) + 1
756: *
757: * Compute scaling factor to survive the linear update
758: * simulating consistent scaling.
759: *
760: CNRM = ZLANGE( 'I', I2-I1, L2-L1, C( I1, L1 ),
761: $ LDC, WNRM )
762: SCAMIN = MIN( SWORK( I, L ), SWORK( K, L ) )
763: CNRM = CNRM * ( SCAMIN / SWORK( I, L ) )
764: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
765: ANRM = SWORK( I, AWRK + K )
766: SCALOC = DLARMM( ANRM, XNRM, CNRM )
767: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
768: * Use second scaling factor to prevent flushing to zero.
769: BUF = BUF*2.D0**EXPONENT( SCALOC )
770: DO JJ = 1, NBB
771: DO LL = 1, NBA
772: SWORK( LL, JJ ) = MIN( BIGNUM,
773: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
774: END DO
775: END DO
776: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
777: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
778: END IF
779: CNRM = CNRM * SCALOC
780: XNRM = XNRM * SCALOC
781: *
782: * Simultaneously apply the robust update factor and the
783: * consistency scaling factor to C( I, L ) and C( K, L).
784: *
785: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
786: IF( SCAL .NE. ONE ) THEN
787: DO LL = L1, L2-1
788: CALL ZDSCAL( K2-K1, SCAL, C( K1, LL ), 1 )
789: END DO
790: ENDIF
791: *
792: SCAL = ( SCAMIN / SWORK( I, L ) ) * SCALOC
793: IF( SCAL .NE. ONE ) THEN
794: DO LL = L1, L2-1
795: CALL ZDSCAL( I2-I1, SCAL, C( I1, LL ), 1 )
796: END DO
797: ENDIF
798: *
799: * Record current scaling factor
800: *
801: SWORK( K, L ) = SCAMIN * SCALOC
802: SWORK( I, L ) = SCAMIN * SCALOC
803: *
804: CALL ZGEMM( 'C', 'N', I2-I1, L2-L1, K2-K1, -CONE,
805: $ A( K1, I1 ), LDA, C( K1, L1 ), LDC,
806: $ CONE, C( I1, L1 ), LDC )
807: END DO
808: *
809: DO J = 1, L - 1
810: *
811: * C( K, J ) := C( K, J ) - SGN * C( K, L ) * B( J, L )**H
812: *
813: J1 = (J - 1) * NB + 1
814: J2 = MIN( J * NB, N ) + 1
815: *
816: * Compute scaling factor to survive the linear update
817: * simulating consistent scaling.
818: *
819: CNRM = ZLANGE( 'I', K2-K1, J2-J1, C( K1, J1 ),
820: $ LDC, WNRM )
821: SCAMIN = MIN( SWORK( K, J ), SWORK( K, L ) )
822: CNRM = CNRM * ( SCAMIN / SWORK( K, J ) )
823: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
824: BNRM = SWORK( L, BWRK + J )
825: SCALOC = DLARMM( BNRM, XNRM, CNRM )
826: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
827: * Use second scaling factor to prevent flushing to zero.
828: BUF = BUF*2.D0**EXPONENT( SCALOC )
829: DO JJ = 1, NBB
830: DO LL = 1, NBA
831: SWORK( LL, JJ ) = MIN( BIGNUM,
832: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
833: END DO
834: END DO
835: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
836: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
837: END IF
838: CNRM = CNRM * SCALOC
839: XNRM = XNRM * SCALOC
840: *
841: * Simultaneously apply the robust update factor and the
842: * consistency scaling factor to C( K, J ) and C( K, L).
843: *
844: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
845: IF( SCAL .NE. ONE ) THEN
846: DO LL = L1, L2-1
847: CALL ZDSCAL( K2-K1, SCAL, C( K1, LL ), 1)
848: END DO
849: ENDIF
850: *
851: SCAL = ( SCAMIN / SWORK( K, J ) ) * SCALOC
852: IF( SCAL .NE. ONE ) THEN
853: DO JJ = J1, J2-1
854: CALL ZDSCAL( K2-K1, SCAL, C( K1, JJ ), 1 )
855: END DO
856: ENDIF
857: *
858: * Record current scaling factor
859: *
860: SWORK( K, L ) = SCAMIN * SCALOC
861: SWORK( K, J ) = SCAMIN * SCALOC
862: *
863: CALL ZGEMM( 'N', 'C', K2-K1, J2-J1, L2-L1, -CSGN,
864: $ C( K1, L1 ), LDC, B( J1, L1 ), LDB,
865: $ CONE, C( K1, J1 ), LDC )
866: END DO
867: END DO
868: END DO
869: ELSE IF( NOTRNA .AND. .NOT.NOTRNB ) THEN
870: *
871: * Solve A*X + ISGN*X*B**H = scale*C.
872: *
873: * The (K,L)th block of X is determined starting from
874: * bottom-right corner column by column by
875: *
876: * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L)**H = C(K,L) - R(K,L)
877: *
878: * Where
879: * M N
880: * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**H].
881: * I=K+1 J=L+1
882: *
883: * Start loop over block rows (index = K) and block columns (index = L)
884: *
885: DO K = NBA, 1, -1
886: *
887: * K1: row index of the first row in X( K, L )
888: * K2: row index of the first row in X( K+1, L )
889: * so the K2 - K1 is the column count of the block X( K, L )
890: *
891: K1 = (K - 1) * NB + 1
892: K2 = MIN( K * NB, M ) + 1
893: DO L = NBB, 1, -1
894: *
895: * L1: column index of the first column in X( K, L )
896: * L2: column index of the first column in X( K, L + 1)
897: * so that L2 - L1 is the row count of the block X( K, L )
898: *
899: L1 = (L - 1) * NB + 1
900: L2 = MIN( L * NB, N ) + 1
901: *
902: CALL ZTRSYL( TRANA, TRANB, ISGN, K2-K1, L2-L1,
903: $ A( K1, K1 ), LDA,
904: $ B( L1, L1 ), LDB,
905: $ C( K1, L1 ), LDC, SCALOC, IINFO )
906: INFO = MAX( INFO, IINFO )
907: *
908: IF( SCALOC * SWORK( K, L ) .EQ. ZERO ) THEN
909: IF( SCALOC .EQ. ZERO ) THEN
910: * The magnitude of the largest entry of X(K1:K2-1, L1:L2-1)
911: * is larger than the product of BIGNUM**2 and cannot be
912: * represented in the form (1/SCALE)*X(K1:K2-1, L1:L2-1).
913: * Mark the computation as pointless.
914: BUF = ZERO
915: ELSE
916: * Use second scaling factor to prevent flushing to zero.
917: BUF = BUF*2.D0**EXPONENT( SCALOC )
918: END IF
919: DO JJ = 1, NBB
920: DO LL = 1, NBA
921: * Bound by BIGNUM to not introduce Inf. The value
922: * is irrelevant; corresponding entries of the
923: * solution will be flushed in consistency scaling.
924: SWORK( LL, JJ ) = MIN( BIGNUM,
925: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
926: END DO
927: END DO
928: END IF
929: SWORK( K, L ) = SCALOC * SWORK( K, L )
930: XNRM = ZLANGE( 'I', K2-K1, L2-L1, C( K1, L1 ), LDC,
931: $ WNRM )
932: *
933: DO I = 1, K - 1
934: *
935: * C( I, L ) := C( I, L ) - A( I, K ) * C( K, L )
936: *
937: I1 = (I - 1) * NB + 1
938: I2 = MIN( I * NB, M ) + 1
939: *
940: * Compute scaling factor to survive the linear update
941: * simulating consistent scaling.
942: *
943: CNRM = ZLANGE( 'I', I2-I1, L2-L1, C( I1, L1 ),
944: $ LDC, WNRM )
945: SCAMIN = MIN( SWORK( I, L ), SWORK( K, L ) )
946: CNRM = CNRM * ( SCAMIN / SWORK( I, L ) )
947: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
948: ANRM = SWORK( I, AWRK + K )
949: SCALOC = DLARMM( ANRM, XNRM, CNRM )
950: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
951: * Use second scaling factor to prevent flushing to zero.
952: BUF = BUF*2.D0**EXPONENT( SCALOC )
953: DO JJ = 1, NBB
954: DO LL = 1, NBA
955: SWORK( LL, JJ ) = MIN( BIGNUM,
956: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
957: END DO
958: END DO
959: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
960: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
961: END IF
962: CNRM = CNRM * SCALOC
963: XNRM = XNRM * SCALOC
964: *
965: * Simultaneously apply the robust update factor and the
966: * consistency scaling factor to C( I, L ) and C( K, L).
967: *
968: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
969: IF( SCAL .NE. ONE ) THEN
970: DO LL = L1, L2-1
971: CALL ZDSCAL( K2-K1, SCAL, C( K1, LL ), 1 )
972: END DO
973: ENDIF
974: *
975: SCAL = ( SCAMIN / SWORK( I, L ) ) * SCALOC
976: IF( SCAL .NE. ONE ) THEN
977: DO LL = L1, L2-1
978: CALL ZDSCAL( I2-I1, SCAL, C( I1, LL ), 1 )
979: END DO
980: ENDIF
981: *
982: * Record current scaling factor
983: *
984: SWORK( K, L ) = SCAMIN * SCALOC
985: SWORK( I, L ) = SCAMIN * SCALOC
986: *
987: CALL ZGEMM( 'N', 'N', I2-I1, L2-L1, K2-K1, -CONE,
988: $ A( I1, K1 ), LDA, C( K1, L1 ), LDC,
989: $ CONE, C( I1, L1 ), LDC )
990: *
991: END DO
992: *
993: DO J = 1, L - 1
994: *
995: * C( K, J ) := C( K, J ) - SGN * C( K, L ) * B( J, L )**H
996: *
997: J1 = (J - 1) * NB + 1
998: J2 = MIN( J * NB, N ) + 1
999: *
1000: * Compute scaling factor to survive the linear update
1001: * simulating consistent scaling.
1002: *
1003: CNRM = ZLANGE( 'I', K2-K1, J2-J1, C( K1, J1 ),
1004: $ LDC, WNRM )
1005: SCAMIN = MIN( SWORK( K, J ), SWORK( K, L ) )
1006: CNRM = CNRM * ( SCAMIN / SWORK( K, J ) )
1007: XNRM = XNRM * ( SCAMIN / SWORK( K, L ) )
1008: BNRM = SWORK( L, BWRK + J )
1009: SCALOC = DLARMM( BNRM, XNRM, CNRM )
1010: IF( SCALOC * SCAMIN .EQ. ZERO ) THEN
1011: * Use second scaling factor to prevent flushing to zero.
1012: BUF = BUF*2.D0**EXPONENT( SCALOC )
1013: DO JJ = 1, NBB
1014: DO LL = 1, NBA
1015: SWORK( LL, JJ ) = MIN( BIGNUM,
1016: $ SWORK( LL, JJ ) / 2.D0**EXPONENT( SCALOC ) )
1017: END DO
1018: END DO
1019: SCAMIN = SCAMIN / 2.D0**EXPONENT( SCALOC )
1020: SCALOC = SCALOC / 2.D0**EXPONENT( SCALOC )
1021: END IF
1022: CNRM = CNRM * SCALOC
1023: XNRM = XNRM * SCALOC
1024: *
1025: * Simultaneously apply the robust update factor and the
1026: * consistency scaling factor to C( K, J ) and C( K, L).
1027: *
1028: SCAL = ( SCAMIN / SWORK( K, L ) ) * SCALOC
1029: IF( SCAL .NE. ONE ) THEN
1030: DO JJ = L1, L2-1
1031: CALL ZDSCAL( K2-K1, SCAL, C( K1, JJ ), 1 )
1032: END DO
1033: ENDIF
1034: *
1035: SCAL = ( SCAMIN / SWORK( K, J ) ) * SCALOC
1036: IF( SCAL .NE. ONE ) THEN
1037: DO JJ = J1, J2-1
1038: CALL ZDSCAL( K2-K1, SCAL, C( K1, JJ ), 1 )
1039: END DO
1040: ENDIF
1041: *
1042: * Record current scaling factor
1043: *
1044: SWORK( K, L ) = SCAMIN * SCALOC
1045: SWORK( K, J ) = SCAMIN * SCALOC
1046: *
1047: CALL ZGEMM( 'N', 'C', K2-K1, J2-J1, L2-L1, -CSGN,
1048: $ C( K1, L1 ), LDC, B( J1, L1 ), LDB,
1049: $ CONE, C( K1, J1 ), LDC )
1050: END DO
1051: END DO
1052: END DO
1053: *
1054: END IF
1055: *
1056: * Reduce local scaling factors
1057: *
1058: SCALE = SWORK( 1, 1 )
1059: DO K = 1, NBA
1060: DO L = 1, NBB
1061: SCALE = MIN( SCALE, SWORK( K, L ) )
1062: END DO
1063: END DO
1064: IF( SCALE .EQ. ZERO ) THEN
1065: *
1066: * The magnitude of the largest entry of the solution is larger
1067: * than the product of BIGNUM**2 and cannot be represented in the
1068: * form (1/SCALE)*X if SCALE is DOUBLE PRECISION. Set SCALE to
1069: * zero and give up.
1070: *
1071: SWORK(1,1) = MAX( NBA, NBB )
1072: SWORK(2,1) = 2 * NBB + NBA
1073: RETURN
1074: END IF
1075: *
1076: * Realize consistent scaling
1077: *
1078: DO K = 1, NBA
1079: K1 = (K - 1) * NB + 1
1080: K2 = MIN( K * NB, M ) + 1
1081: DO L = 1, NBB
1082: L1 = (L - 1) * NB + 1
1083: L2 = MIN( L * NB, N ) + 1
1084: SCAL = SCALE / SWORK( K, L )
1085: IF( SCAL .NE. ONE ) THEN
1086: DO LL = L1, L2-1
1087: CALL ZDSCAL( K2-K1, SCAL, C( K1, LL ), 1 )
1088: END DO
1089: ENDIF
1090: END DO
1091: END DO
1092: *
1093: IF( BUF .NE. ONE .AND. BUF.GT.ZERO ) THEN
1094: *
1095: * Decrease SCALE as much as possible.
1096: *
1097: SCALOC = MIN( SCALE / SMLNUM, ONE / BUF )
1098: BUF = BUF * SCALOC
1099: SCALE = SCALE / SCALOC
1100: END IF
1101: *
1102: IF( BUF.NE.ONE .AND. BUF.GT.ZERO ) THEN
1103: *
1104: * In case of overly aggressive scaling during the computation,
1105: * flushing of the global scale factor may be prevented by
1106: * undoing some of the scaling. This step is to ensure that
1107: * this routine flushes only scale factors that TRSYL also
1108: * flushes and be usable as a drop-in replacement.
1109: *
1110: * How much can the normwise largest entry be upscaled?
1111: *
1112: SCAL = MAX( ABS( DBLE( C( 1, 1 ) ) ),
1113: $ ABS( DIMAG( C ( 1, 1 ) ) ) )
1114: DO K = 1, M
1115: DO L = 1, N
1116: SCAL = MAX( SCAL, ABS( DBLE ( C( K, L ) ) ),
1117: $ ABS( DIMAG ( C( K, L ) ) ) )
1118: END DO
1119: END DO
1120: *
1121: * Increase BUF as close to 1 as possible and apply scaling.
1122: *
1123: SCALOC = MIN( BIGNUM / SCAL, ONE / BUF )
1124: BUF = BUF * SCALOC
1125: CALL ZLASCL( 'G', -1, -1, ONE, SCALOC, M, N, C, LDC, IINFO )
1126: END IF
1127: *
1128: * Combine with buffer scaling factor. SCALE will be flushed if
1129: * BUF is less than one here.
1130: *
1131: SCALE = SCALE * BUF
1132: *
1133: * Restore workspace dimensions
1134: *
1135: SWORK(1,1) = MAX( NBA, NBB )
1136: SWORK(2,1) = 2 * NBB + NBA
1137: *
1138: RETURN
1139: *
1140: * End of ZTRSYL3
1141: *
1142: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to ztrsyl3.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack