1: SUBROUTINE ZTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
2: + B, LDB )
3: *
4: * -- LAPACK routine (version 3.3.0) --
5: *
6: * -- Contributed by Fred Gustavson of the IBM Watson Research Center --
7: * November 2010
8: *
9: * -- LAPACK is a software package provided by Univ. of Tennessee, --
10: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
11: *
12: * ..
13: * .. Scalar Arguments ..
14: CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO
15: INTEGER LDB, M, N
16: COMPLEX*16 ALPHA
17: * ..
18: * .. Array Arguments ..
19: COMPLEX*16 A( 0: * ), B( 0: LDB-1, 0: * )
20: * ..
21: *
22: * Purpose
23: * =======
24: *
25: * Level 3 BLAS like routine for A in RFP Format.
26: *
27: * ZTFSM solves the matrix equation
28: *
29: * op( A )*X = alpha*B or X*op( A ) = alpha*B
30: *
31: * where alpha is a scalar, X and B are m by n matrices, A is a unit, or
32: * non-unit, upper or lower triangular matrix and op( A ) is one of
33: *
34: * op( A ) = A or op( A ) = conjg( A' ).
35: *
36: * A is in Rectangular Full Packed (RFP) Format.
37: *
38: * The matrix X is overwritten on B.
39: *
40: * Arguments
41: * ==========
42: *
43: * TRANSR (input) CHARACTER*1
44: * = 'N': The Normal Form of RFP A is stored;
45: * = 'C': The Conjugate-transpose Form of RFP A is stored.
46: *
47: * SIDE (input) CHARACTER*1
48: * On entry, SIDE specifies whether op( A ) appears on the left
49: * or right of X as follows:
50: *
51: * SIDE = 'L' or 'l' op( A )*X = alpha*B.
52: *
53: * SIDE = 'R' or 'r' X*op( A ) = alpha*B.
54: *
55: * Unchanged on exit.
56: *
57: * UPLO (input) CHARACTER*1
58: * On entry, UPLO specifies whether the RFP matrix A came from
59: * an upper or lower triangular matrix as follows:
60: * UPLO = 'U' or 'u' RFP A came from an upper triangular matrix
61: * UPLO = 'L' or 'l' RFP A came from a lower triangular matrix
62: *
63: * Unchanged on exit.
64: *
65: * TRANS (input) CHARACTER*1
66: * On entry, TRANS specifies the form of op( A ) to be used
67: * in the matrix multiplication as follows:
68: *
69: * TRANS = 'N' or 'n' op( A ) = A.
70: *
71: * TRANS = 'C' or 'c' op( A ) = conjg( A' ).
72: *
73: * Unchanged on exit.
74: *
75: * DIAG (input) CHARACTER*1
76: * On entry, DIAG specifies whether or not RFP A is unit
77: * triangular as follows:
78: *
79: * DIAG = 'U' or 'u' A is assumed to be unit triangular.
80: *
81: * DIAG = 'N' or 'n' A is not assumed to be unit
82: * triangular.
83: *
84: * Unchanged on exit.
85: *
86: * M (input) INTEGER
87: * On entry, M specifies the number of rows of B. M must be at
88: * least zero.
89: * Unchanged on exit.
90: *
91: * N (input) INTEGER
92: * On entry, N specifies the number of columns of B. N must be
93: * at least zero.
94: * Unchanged on exit.
95: *
96: * ALPHA (input) COMPLEX*16
97: * On entry, ALPHA specifies the scalar alpha. When alpha is
98: * zero then A is not referenced and B need not be set before
99: * entry.
100: * Unchanged on exit.
101: *
102: * A (input) COMPLEX*16 array, dimension (N*(N+1)/2)
103: * NT = N*(N+1)/2. On entry, the matrix A in RFP Format.
104: * RFP Format is described by TRANSR, UPLO and N as follows:
105: * If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
106: * K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
107: * TRANSR = 'C' then RFP is the Conjugate-transpose of RFP A as
108: * defined when TRANSR = 'N'. The contents of RFP A are defined
109: * by UPLO as follows: If UPLO = 'U' the RFP A contains the NT
110: * elements of upper packed A either in normal or
111: * conjugate-transpose Format. If UPLO = 'L' the RFP A contains
112: * the NT elements of lower packed A either in normal or
113: * conjugate-transpose Format. The LDA of RFP A is (N+1)/2 when
114: * TRANSR = 'C'. When TRANSR is 'N' the LDA is N+1 when N is
115: * even and is N when is odd.
116: * See the Note below for more details. Unchanged on exit.
117: *
118: * B (input/output) COMPLEX*16 array, dimension (LDB,N)
119: * Before entry, the leading m by n part of the array B must
120: * contain the right-hand side matrix B, and on exit is
121: * overwritten by the solution matrix X.
122: *
123: * LDB (input) INTEGER
124: * On entry, LDB specifies the first dimension of B as declared
125: * in the calling (sub) program. LDB must be at least
126: * max( 1, m ).
127: * Unchanged on exit.
128: *
129: * Further Details
130: * ===============
131: *
132: * We first consider Standard Packed Format when N is even.
133: * We give an example where N = 6.
134: *
135: * AP is Upper AP is Lower
136: *
137: * 00 01 02 03 04 05 00
138: * 11 12 13 14 15 10 11
139: * 22 23 24 25 20 21 22
140: * 33 34 35 30 31 32 33
141: * 44 45 40 41 42 43 44
142: * 55 50 51 52 53 54 55
143: *
144: *
145: * Let TRANSR = 'N'. RFP holds AP as follows:
146: * For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
147: * three columns of AP upper. The lower triangle A(4:6,0:2) consists of
148: * conjugate-transpose of the first three columns of AP upper.
149: * For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
150: * three columns of AP lower. The upper triangle A(0:2,0:2) consists of
151: * conjugate-transpose of the last three columns of AP lower.
152: * To denote conjugate we place -- above the element. This covers the
153: * case N even and TRANSR = 'N'.
154: *
155: * RFP A RFP A
156: *
157: * -- -- --
158: * 03 04 05 33 43 53
159: * -- --
160: * 13 14 15 00 44 54
161: * --
162: * 23 24 25 10 11 55
163: *
164: * 33 34 35 20 21 22
165: * --
166: * 00 44 45 30 31 32
167: * -- --
168: * 01 11 55 40 41 42
169: * -- -- --
170: * 02 12 22 50 51 52
171: *
172: * Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
173: * transpose of RFP A above. One therefore gets:
174: *
175: *
176: * RFP A RFP A
177: *
178: * -- -- -- -- -- -- -- -- -- --
179: * 03 13 23 33 00 01 02 33 00 10 20 30 40 50
180: * -- -- -- -- -- -- -- -- -- --
181: * 04 14 24 34 44 11 12 43 44 11 21 31 41 51
182: * -- -- -- -- -- -- -- -- -- --
183: * 05 15 25 35 45 55 22 53 54 55 22 32 42 52
184: *
185: *
186: * We next consider Standard Packed Format when N is odd.
187: * We give an example where N = 5.
188: *
189: * AP is Upper AP is Lower
190: *
191: * 00 01 02 03 04 00
192: * 11 12 13 14 10 11
193: * 22 23 24 20 21 22
194: * 33 34 30 31 32 33
195: * 44 40 41 42 43 44
196: *
197: *
198: * Let TRANSR = 'N'. RFP holds AP as follows:
199: * For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
200: * three columns of AP upper. The lower triangle A(3:4,0:1) consists of
201: * conjugate-transpose of the first two columns of AP upper.
202: * For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
203: * three columns of AP lower. The upper triangle A(0:1,1:2) consists of
204: * conjugate-transpose of the last two columns of AP lower.
205: * To denote conjugate we place -- above the element. This covers the
206: * case N odd and TRANSR = 'N'.
207: *
208: * RFP A RFP A
209: *
210: * -- --
211: * 02 03 04 00 33 43
212: * --
213: * 12 13 14 10 11 44
214: *
215: * 22 23 24 20 21 22
216: * --
217: * 00 33 34 30 31 32
218: * -- --
219: * 01 11 44 40 41 42
220: *
221: * Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
222: * transpose of RFP A above. One therefore gets:
223: *
224: *
225: * RFP A RFP A
226: *
227: * -- -- -- -- -- -- -- -- --
228: * 02 12 22 00 01 00 10 20 30 40 50
229: * -- -- -- -- -- -- -- -- --
230: * 03 13 23 33 11 33 11 21 31 41 51
231: * -- -- -- -- -- -- -- -- --
232: * 04 14 24 34 44 43 44 22 32 42 52
233: *
234: * ..
235: * .. Parameters ..
236: COMPLEX*16 CONE, CZERO
237: PARAMETER ( CONE = ( 1.0D+0, 0.0D+0 ),
238: + CZERO = ( 0.0D+0, 0.0D+0 ) )
239: * ..
240: * .. Local Scalars ..
241: LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
242: + NOTRANS
243: INTEGER M1, M2, N1, N2, K, INFO, I, J
244: * ..
245: * .. External Functions ..
246: LOGICAL LSAME
247: EXTERNAL LSAME
248: * ..
249: * .. External Subroutines ..
250: EXTERNAL XERBLA, ZGEMM, ZTRSM
251: * ..
252: * .. Intrinsic Functions ..
253: INTRINSIC MAX, MOD
254: * ..
255: * .. Executable Statements ..
256: *
257: * Test the input parameters.
258: *
259: INFO = 0
260: NORMALTRANSR = LSAME( TRANSR, 'N' )
261: LSIDE = LSAME( SIDE, 'L' )
262: LOWER = LSAME( UPLO, 'L' )
263: NOTRANS = LSAME( TRANS, 'N' )
264: IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'C' ) ) THEN
265: INFO = -1
266: ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
267: INFO = -2
268: ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
269: INFO = -3
270: ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN
271: INFO = -4
272: ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) )
273: + THEN
274: INFO = -5
275: ELSE IF( M.LT.0 ) THEN
276: INFO = -6
277: ELSE IF( N.LT.0 ) THEN
278: INFO = -7
279: ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
280: INFO = -11
281: END IF
282: IF( INFO.NE.0 ) THEN
283: CALL XERBLA( 'ZTFSM ', -INFO )
284: RETURN
285: END IF
286: *
287: * Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
288: *
289: IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
290: + RETURN
291: *
292: * Quick return when ALPHA.EQ.(0D+0,0D+0)
293: *
294: IF( ALPHA.EQ.CZERO ) THEN
295: DO 20 J = 0, N - 1
296: DO 10 I = 0, M - 1
297: B( I, J ) = CZERO
298: 10 CONTINUE
299: 20 CONTINUE
300: RETURN
301: END IF
302: *
303: IF( LSIDE ) THEN
304: *
305: * SIDE = 'L'
306: *
307: * A is M-by-M.
308: * If M is odd, set NISODD = .TRUE., and M1 and M2.
309: * If M is even, NISODD = .FALSE., and M.
310: *
311: IF( MOD( M, 2 ).EQ.0 ) THEN
312: MISODD = .FALSE.
313: K = M / 2
314: ELSE
315: MISODD = .TRUE.
316: IF( LOWER ) THEN
317: M2 = M / 2
318: M1 = M - M2
319: ELSE
320: M1 = M / 2
321: M2 = M - M1
322: END IF
323: END IF
324: *
325: IF( MISODD ) THEN
326: *
327: * SIDE = 'L' and N is odd
328: *
329: IF( NORMALTRANSR ) THEN
330: *
331: * SIDE = 'L', N is odd, and TRANSR = 'N'
332: *
333: IF( LOWER ) THEN
334: *
335: * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
336: *
337: IF( NOTRANS ) THEN
338: *
339: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
340: * TRANS = 'N'
341: *
342: IF( M.EQ.1 ) THEN
343: CALL ZTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
344: + A, M, B, LDB )
345: ELSE
346: CALL ZTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
347: + A( 0 ), M, B, LDB )
348: CALL ZGEMM( 'N', 'N', M2, N, M1, -CONE, A( M1 ),
349: + M, B, LDB, ALPHA, B( M1, 0 ), LDB )
350: CALL ZTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
351: + A( M ), M, B( M1, 0 ), LDB )
352: END IF
353: *
354: ELSE
355: *
356: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
357: * TRANS = 'C'
358: *
359: IF( M.EQ.1 ) THEN
360: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, ALPHA,
361: + A( 0 ), M, B, LDB )
362: ELSE
363: CALL ZTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
364: + A( M ), M, B( M1, 0 ), LDB )
365: CALL ZGEMM( 'C', 'N', M1, N, M2, -CONE, A( M1 ),
366: + M, B( M1, 0 ), LDB, ALPHA, B, LDB )
367: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
368: + A( 0 ), M, B, LDB )
369: END IF
370: *
371: END IF
372: *
373: ELSE
374: *
375: * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
376: *
377: IF( .NOT.NOTRANS ) THEN
378: *
379: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
380: * TRANS = 'N'
381: *
382: CALL ZTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
383: + A( M2 ), M, B, LDB )
384: CALL ZGEMM( 'C', 'N', M2, N, M1, -CONE, A( 0 ), M,
385: + B, LDB, ALPHA, B( M1, 0 ), LDB )
386: CALL ZTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
387: + A( M1 ), M, B( M1, 0 ), LDB )
388: *
389: ELSE
390: *
391: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
392: * TRANS = 'C'
393: *
394: CALL ZTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
395: + A( M1 ), M, B( M1, 0 ), LDB )
396: CALL ZGEMM( 'N', 'N', M1, N, M2, -CONE, A( 0 ), M,
397: + B( M1, 0 ), LDB, ALPHA, B, LDB )
398: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
399: + A( M2 ), M, B, LDB )
400: *
401: END IF
402: *
403: END IF
404: *
405: ELSE
406: *
407: * SIDE = 'L', N is odd, and TRANSR = 'C'
408: *
409: IF( LOWER ) THEN
410: *
411: * SIDE ='L', N is odd, TRANSR = 'C', and UPLO = 'L'
412: *
413: IF( NOTRANS ) THEN
414: *
415: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'L', and
416: * TRANS = 'N'
417: *
418: IF( M.EQ.1 ) THEN
419: CALL ZTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
420: + A( 0 ), M1, B, LDB )
421: ELSE
422: CALL ZTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
423: + A( 0 ), M1, B, LDB )
424: CALL ZGEMM( 'C', 'N', M2, N, M1, -CONE,
425: + A( M1*M1 ), M1, B, LDB, ALPHA,
426: + B( M1, 0 ), LDB )
427: CALL ZTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
428: + A( 1 ), M1, B( M1, 0 ), LDB )
429: END IF
430: *
431: ELSE
432: *
433: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'L', and
434: * TRANS = 'C'
435: *
436: IF( M.EQ.1 ) THEN
437: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA,
438: + A( 0 ), M1, B, LDB )
439: ELSE
440: CALL ZTRSM( 'L', 'L', 'C', DIAG, M2, N, ALPHA,
441: + A( 1 ), M1, B( M1, 0 ), LDB )
442: CALL ZGEMM( 'N', 'N', M1, N, M2, -CONE,
443: + A( M1*M1 ), M1, B( M1, 0 ), LDB,
444: + ALPHA, B, LDB )
445: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
446: + A( 0 ), M1, B, LDB )
447: END IF
448: *
449: END IF
450: *
451: ELSE
452: *
453: * SIDE ='L', N is odd, TRANSR = 'C', and UPLO = 'U'
454: *
455: IF( .NOT.NOTRANS ) THEN
456: *
457: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'U', and
458: * TRANS = 'N'
459: *
460: CALL ZTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
461: + A( M2*M2 ), M2, B, LDB )
462: CALL ZGEMM( 'N', 'N', M2, N, M1, -CONE, A( 0 ), M2,
463: + B, LDB, ALPHA, B( M1, 0 ), LDB )
464: CALL ZTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
465: + A( M1*M2 ), M2, B( M1, 0 ), LDB )
466: *
467: ELSE
468: *
469: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'U', and
470: * TRANS = 'C'
471: *
472: CALL ZTRSM( 'L', 'L', 'C', DIAG, M2, N, ALPHA,
473: + A( M1*M2 ), M2, B( M1, 0 ), LDB )
474: CALL ZGEMM( 'C', 'N', M1, N, M2, -CONE, A( 0 ), M2,
475: + B( M1, 0 ), LDB, ALPHA, B, LDB )
476: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
477: + A( M2*M2 ), M2, B, LDB )
478: *
479: END IF
480: *
481: END IF
482: *
483: END IF
484: *
485: ELSE
486: *
487: * SIDE = 'L' and N is even
488: *
489: IF( NORMALTRANSR ) THEN
490: *
491: * SIDE = 'L', N is even, and TRANSR = 'N'
492: *
493: IF( LOWER ) THEN
494: *
495: * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L'
496: *
497: IF( NOTRANS ) THEN
498: *
499: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
500: * and TRANS = 'N'
501: *
502: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
503: + A( 1 ), M+1, B, LDB )
504: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( K+1 ),
505: + M+1, B, LDB, ALPHA, B( K, 0 ), LDB )
506: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
507: + A( 0 ), M+1, B( K, 0 ), LDB )
508: *
509: ELSE
510: *
511: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
512: * and TRANS = 'C'
513: *
514: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
515: + A( 0 ), M+1, B( K, 0 ), LDB )
516: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( K+1 ),
517: + M+1, B( K, 0 ), LDB, ALPHA, B, LDB )
518: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
519: + A( 1 ), M+1, B, LDB )
520: *
521: END IF
522: *
523: ELSE
524: *
525: * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U'
526: *
527: IF( .NOT.NOTRANS ) THEN
528: *
529: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
530: * and TRANS = 'N'
531: *
532: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
533: + A( K+1 ), M+1, B, LDB )
534: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), M+1,
535: + B, LDB, ALPHA, B( K, 0 ), LDB )
536: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
537: + A( K ), M+1, B( K, 0 ), LDB )
538: *
539: ELSE
540: *
541: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
542: * and TRANS = 'C'
543: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
544: + A( K ), M+1, B( K, 0 ), LDB )
545: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), M+1,
546: + B( K, 0 ), LDB, ALPHA, B, LDB )
547: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
548: + A( K+1 ), M+1, B, LDB )
549: *
550: END IF
551: *
552: END IF
553: *
554: ELSE
555: *
556: * SIDE = 'L', N is even, and TRANSR = 'C'
557: *
558: IF( LOWER ) THEN
559: *
560: * SIDE ='L', N is even, TRANSR = 'C', and UPLO = 'L'
561: *
562: IF( NOTRANS ) THEN
563: *
564: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'L',
565: * and TRANS = 'N'
566: *
567: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, ALPHA,
568: + A( K ), K, B, LDB )
569: CALL ZGEMM( 'C', 'N', K, N, K, -CONE,
570: + A( K*( K+1 ) ), K, B, LDB, ALPHA,
571: + B( K, 0 ), LDB )
572: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
573: + A( 0 ), K, B( K, 0 ), LDB )
574: *
575: ELSE
576: *
577: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'L',
578: * and TRANS = 'C'
579: *
580: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, ALPHA,
581: + A( 0 ), K, B( K, 0 ), LDB )
582: CALL ZGEMM( 'N', 'N', K, N, K, -CONE,
583: + A( K*( K+1 ) ), K, B( K, 0 ), LDB,
584: + ALPHA, B, LDB )
585: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
586: + A( K ), K, B, LDB )
587: *
588: END IF
589: *
590: ELSE
591: *
592: * SIDE ='L', N is even, TRANSR = 'C', and UPLO = 'U'
593: *
594: IF( .NOT.NOTRANS ) THEN
595: *
596: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'U',
597: * and TRANS = 'N'
598: *
599: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, ALPHA,
600: + A( K*( K+1 ) ), K, B, LDB )
601: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), K, B,
602: + LDB, ALPHA, B( K, 0 ), LDB )
603: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
604: + A( K*K ), K, B( K, 0 ), LDB )
605: *
606: ELSE
607: *
608: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'U',
609: * and TRANS = 'C'
610: *
611: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, ALPHA,
612: + A( K*K ), K, B( K, 0 ), LDB )
613: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), K,
614: + B( K, 0 ), LDB, ALPHA, B, LDB )
615: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
616: + A( K*( K+1 ) ), K, B, LDB )
617: *
618: END IF
619: *
620: END IF
621: *
622: END IF
623: *
624: END IF
625: *
626: ELSE
627: *
628: * SIDE = 'R'
629: *
630: * A is N-by-N.
631: * If N is odd, set NISODD = .TRUE., and N1 and N2.
632: * If N is even, NISODD = .FALSE., and K.
633: *
634: IF( MOD( N, 2 ).EQ.0 ) THEN
635: NISODD = .FALSE.
636: K = N / 2
637: ELSE
638: NISODD = .TRUE.
639: IF( LOWER ) THEN
640: N2 = N / 2
641: N1 = N - N2
642: ELSE
643: N1 = N / 2
644: N2 = N - N1
645: END IF
646: END IF
647: *
648: IF( NISODD ) THEN
649: *
650: * SIDE = 'R' and N is odd
651: *
652: IF( NORMALTRANSR ) THEN
653: *
654: * SIDE = 'R', N is odd, and TRANSR = 'N'
655: *
656: IF( LOWER ) THEN
657: *
658: * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
659: *
660: IF( NOTRANS ) THEN
661: *
662: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
663: * TRANS = 'N'
664: *
665: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N2, ALPHA,
666: + A( N ), N, B( 0, N1 ), LDB )
667: CALL ZGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
668: + LDB, A( N1 ), N, ALPHA, B( 0, 0 ),
669: + LDB )
670: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
671: + A( 0 ), N, B( 0, 0 ), LDB )
672: *
673: ELSE
674: *
675: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
676: * TRANS = 'C'
677: *
678: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N1, ALPHA,
679: + A( 0 ), N, B( 0, 0 ), LDB )
680: CALL ZGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
681: + LDB, A( N1 ), N, ALPHA, B( 0, N1 ),
682: + LDB )
683: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
684: + A( N ), N, B( 0, N1 ), LDB )
685: *
686: END IF
687: *
688: ELSE
689: *
690: * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
691: *
692: IF( NOTRANS ) THEN
693: *
694: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
695: * TRANS = 'N'
696: *
697: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N1, ALPHA,
698: + A( N2 ), N, B( 0, 0 ), LDB )
699: CALL ZGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
700: + LDB, A( 0 ), N, ALPHA, B( 0, N1 ),
701: + LDB )
702: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
703: + A( N1 ), N, B( 0, N1 ), LDB )
704: *
705: ELSE
706: *
707: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
708: * TRANS = 'C'
709: *
710: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N2, ALPHA,
711: + A( N1 ), N, B( 0, N1 ), LDB )
712: CALL ZGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
713: + LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB )
714: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
715: + A( N2 ), N, B( 0, 0 ), LDB )
716: *
717: END IF
718: *
719: END IF
720: *
721: ELSE
722: *
723: * SIDE = 'R', N is odd, and TRANSR = 'C'
724: *
725: IF( LOWER ) THEN
726: *
727: * SIDE ='R', N is odd, TRANSR = 'C', and UPLO = 'L'
728: *
729: IF( NOTRANS ) THEN
730: *
731: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'L', and
732: * TRANS = 'N'
733: *
734: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
735: + A( 1 ), N1, B( 0, N1 ), LDB )
736: CALL ZGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
737: + LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ),
738: + LDB )
739: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
740: + A( 0 ), N1, B( 0, 0 ), LDB )
741: *
742: ELSE
743: *
744: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'L', and
745: * TRANS = 'C'
746: *
747: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
748: + A( 0 ), N1, B( 0, 0 ), LDB )
749: CALL ZGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
750: + LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ),
751: + LDB )
752: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
753: + A( 1 ), N1, B( 0, N1 ), LDB )
754: *
755: END IF
756: *
757: ELSE
758: *
759: * SIDE ='R', N is odd, TRANSR = 'C', and UPLO = 'U'
760: *
761: IF( NOTRANS ) THEN
762: *
763: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'U', and
764: * TRANS = 'N'
765: *
766: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
767: + A( N2*N2 ), N2, B( 0, 0 ), LDB )
768: CALL ZGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
769: + LDB, A( 0 ), N2, ALPHA, B( 0, N1 ),
770: + LDB )
771: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
772: + A( N1*N2 ), N2, B( 0, N1 ), LDB )
773: *
774: ELSE
775: *
776: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'U', and
777: * TRANS = 'C'
778: *
779: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
780: + A( N1*N2 ), N2, B( 0, N1 ), LDB )
781: CALL ZGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
782: + LDB, A( 0 ), N2, ALPHA, B( 0, 0 ),
783: + LDB )
784: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
785: + A( N2*N2 ), N2, B( 0, 0 ), LDB )
786: *
787: END IF
788: *
789: END IF
790: *
791: END IF
792: *
793: ELSE
794: *
795: * SIDE = 'R' and N is even
796: *
797: IF( NORMALTRANSR ) THEN
798: *
799: * SIDE = 'R', N is even, and TRANSR = 'N'
800: *
801: IF( LOWER ) THEN
802: *
803: * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L'
804: *
805: IF( NOTRANS ) THEN
806: *
807: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
808: * and TRANS = 'N'
809: *
810: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, ALPHA,
811: + A( 0 ), N+1, B( 0, K ), LDB )
812: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
813: + LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ),
814: + LDB )
815: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
816: + A( 1 ), N+1, B( 0, 0 ), LDB )
817: *
818: ELSE
819: *
820: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
821: * and TRANS = 'C'
822: *
823: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, ALPHA,
824: + A( 1 ), N+1, B( 0, 0 ), LDB )
825: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
826: + LDB, A( K+1 ), N+1, ALPHA, B( 0, K ),
827: + LDB )
828: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
829: + A( 0 ), N+1, B( 0, K ), LDB )
830: *
831: END IF
832: *
833: ELSE
834: *
835: * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U'
836: *
837: IF( NOTRANS ) THEN
838: *
839: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
840: * and TRANS = 'N'
841: *
842: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, ALPHA,
843: + A( K+1 ), N+1, B( 0, 0 ), LDB )
844: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
845: + LDB, A( 0 ), N+1, ALPHA, B( 0, K ),
846: + LDB )
847: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
848: + A( K ), N+1, B( 0, K ), LDB )
849: *
850: ELSE
851: *
852: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
853: * and TRANS = 'C'
854: *
855: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, ALPHA,
856: + A( K ), N+1, B( 0, K ), LDB )
857: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
858: + LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ),
859: + LDB )
860: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
861: + A( K+1 ), N+1, B( 0, 0 ), LDB )
862: *
863: END IF
864: *
865: END IF
866: *
867: ELSE
868: *
869: * SIDE = 'R', N is even, and TRANSR = 'C'
870: *
871: IF( LOWER ) THEN
872: *
873: * SIDE ='R', N is even, TRANSR = 'C', and UPLO = 'L'
874: *
875: IF( NOTRANS ) THEN
876: *
877: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'L',
878: * and TRANS = 'N'
879: *
880: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
881: + A( 0 ), K, B( 0, K ), LDB )
882: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
883: + LDB, A( ( K+1 )*K ), K, ALPHA,
884: + B( 0, 0 ), LDB )
885: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
886: + A( K ), K, B( 0, 0 ), LDB )
887: *
888: ELSE
889: *
890: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'L',
891: * and TRANS = 'C'
892: *
893: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
894: + A( K ), K, B( 0, 0 ), LDB )
895: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
896: + LDB, A( ( K+1 )*K ), K, ALPHA,
897: + B( 0, K ), LDB )
898: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
899: + A( 0 ), K, B( 0, K ), LDB )
900: *
901: END IF
902: *
903: ELSE
904: *
905: * SIDE ='R', N is even, TRANSR = 'C', and UPLO = 'U'
906: *
907: IF( NOTRANS ) THEN
908: *
909: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'U',
910: * and TRANS = 'N'
911: *
912: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
913: + A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
914: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
915: + LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB )
916: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
917: + A( K*K ), K, B( 0, K ), LDB )
918: *
919: ELSE
920: *
921: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'U',
922: * and TRANS = 'C'
923: *
924: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
925: + A( K*K ), K, B( 0, K ), LDB )
926: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
927: + LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB )
928: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
929: + A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
930: *
931: END IF
932: *
933: END IF
934: *
935: END IF
936: *
937: END IF
938: END IF
939: *
940: RETURN
941: *
942: * End of ZTFSM
943: *
944: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to ztfsm.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack