1: SUBROUTINE ZTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
2: $ B, LDB )
3: *
4: * -- LAPACK routine (version 3.3.1) --
5: *
6: * -- Contributed by Fred Gustavson of the IBM Watson Research Center --
7: * -- April 2011 --
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 ) = A**H.
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: * ..
236: * .. Parameters ..
237: COMPLEX*16 CONE, CZERO
238: PARAMETER ( CONE = ( 1.0D+0, 0.0D+0 ),
239: $ CZERO = ( 0.0D+0, 0.0D+0 ) )
240: * ..
241: * .. Local Scalars ..
242: LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
243: $ NOTRANS
244: INTEGER M1, M2, N1, N2, K, INFO, I, J
245: * ..
246: * .. External Functions ..
247: LOGICAL LSAME
248: EXTERNAL LSAME
249: * ..
250: * .. External Subroutines ..
251: EXTERNAL XERBLA, ZGEMM, ZTRSM
252: * ..
253: * .. Intrinsic Functions ..
254: INTRINSIC MAX, MOD
255: * ..
256: * .. Executable Statements ..
257: *
258: * Test the input parameters.
259: *
260: INFO = 0
261: NORMALTRANSR = LSAME( TRANSR, 'N' )
262: LSIDE = LSAME( SIDE, 'L' )
263: LOWER = LSAME( UPLO, 'L' )
264: NOTRANS = LSAME( TRANS, 'N' )
265: IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'C' ) ) THEN
266: INFO = -1
267: ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
268: INFO = -2
269: ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
270: INFO = -3
271: ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN
272: INFO = -4
273: ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) )
274: $ THEN
275: INFO = -5
276: ELSE IF( M.LT.0 ) THEN
277: INFO = -6
278: ELSE IF( N.LT.0 ) THEN
279: INFO = -7
280: ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
281: INFO = -11
282: END IF
283: IF( INFO.NE.0 ) THEN
284: CALL XERBLA( 'ZTFSM ', -INFO )
285: RETURN
286: END IF
287: *
288: * Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
289: *
290: IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
291: $ RETURN
292: *
293: * Quick return when ALPHA.EQ.(0D+0,0D+0)
294: *
295: IF( ALPHA.EQ.CZERO ) THEN
296: DO 20 J = 0, N - 1
297: DO 10 I = 0, M - 1
298: B( I, J ) = CZERO
299: 10 CONTINUE
300: 20 CONTINUE
301: RETURN
302: END IF
303: *
304: IF( LSIDE ) THEN
305: *
306: * SIDE = 'L'
307: *
308: * A is M-by-M.
309: * If M is odd, set NISODD = .TRUE., and M1 and M2.
310: * If M is even, NISODD = .FALSE., and M.
311: *
312: IF( MOD( M, 2 ).EQ.0 ) THEN
313: MISODD = .FALSE.
314: K = M / 2
315: ELSE
316: MISODD = .TRUE.
317: IF( LOWER ) THEN
318: M2 = M / 2
319: M1 = M - M2
320: ELSE
321: M1 = M / 2
322: M2 = M - M1
323: END IF
324: END IF
325: *
326: IF( MISODD ) THEN
327: *
328: * SIDE = 'L' and N is odd
329: *
330: IF( NORMALTRANSR ) THEN
331: *
332: * SIDE = 'L', N is odd, and TRANSR = 'N'
333: *
334: IF( LOWER ) THEN
335: *
336: * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
337: *
338: IF( NOTRANS ) THEN
339: *
340: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
341: * TRANS = 'N'
342: *
343: IF( M.EQ.1 ) THEN
344: CALL ZTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
345: $ A, M, B, LDB )
346: ELSE
347: CALL ZTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
348: $ A( 0 ), M, B, LDB )
349: CALL ZGEMM( 'N', 'N', M2, N, M1, -CONE, A( M1 ),
350: $ M, B, LDB, ALPHA, B( M1, 0 ), LDB )
351: CALL ZTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
352: $ A( M ), M, B( M1, 0 ), LDB )
353: END IF
354: *
355: ELSE
356: *
357: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
358: * TRANS = 'C'
359: *
360: IF( M.EQ.1 ) THEN
361: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, ALPHA,
362: $ A( 0 ), M, B, LDB )
363: ELSE
364: CALL ZTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
365: $ A( M ), M, B( M1, 0 ), LDB )
366: CALL ZGEMM( 'C', 'N', M1, N, M2, -CONE, A( M1 ),
367: $ M, B( M1, 0 ), LDB, ALPHA, B, LDB )
368: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
369: $ A( 0 ), M, B, LDB )
370: END IF
371: *
372: END IF
373: *
374: ELSE
375: *
376: * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
377: *
378: IF( .NOT.NOTRANS ) THEN
379: *
380: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
381: * TRANS = 'N'
382: *
383: CALL ZTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
384: $ A( M2 ), M, B, LDB )
385: CALL ZGEMM( 'C', 'N', M2, N, M1, -CONE, A( 0 ), M,
386: $ B, LDB, ALPHA, B( M1, 0 ), LDB )
387: CALL ZTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
388: $ A( M1 ), M, B( M1, 0 ), LDB )
389: *
390: ELSE
391: *
392: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
393: * TRANS = 'C'
394: *
395: CALL ZTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
396: $ A( M1 ), M, B( M1, 0 ), LDB )
397: CALL ZGEMM( 'N', 'N', M1, N, M2, -CONE, A( 0 ), M,
398: $ B( M1, 0 ), LDB, ALPHA, B, LDB )
399: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
400: $ A( M2 ), M, B, LDB )
401: *
402: END IF
403: *
404: END IF
405: *
406: ELSE
407: *
408: * SIDE = 'L', N is odd, and TRANSR = 'C'
409: *
410: IF( LOWER ) THEN
411: *
412: * SIDE ='L', N is odd, TRANSR = 'C', and UPLO = 'L'
413: *
414: IF( NOTRANS ) THEN
415: *
416: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'L', and
417: * TRANS = 'N'
418: *
419: IF( M.EQ.1 ) THEN
420: CALL ZTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
421: $ A( 0 ), M1, B, LDB )
422: ELSE
423: CALL ZTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
424: $ A( 0 ), M1, B, LDB )
425: CALL ZGEMM( 'C', 'N', M2, N, M1, -CONE,
426: $ A( M1*M1 ), M1, B, LDB, ALPHA,
427: $ B( M1, 0 ), LDB )
428: CALL ZTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
429: $ A( 1 ), M1, B( M1, 0 ), LDB )
430: END IF
431: *
432: ELSE
433: *
434: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'L', and
435: * TRANS = 'C'
436: *
437: IF( M.EQ.1 ) THEN
438: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA,
439: $ A( 0 ), M1, B, LDB )
440: ELSE
441: CALL ZTRSM( 'L', 'L', 'C', DIAG, M2, N, ALPHA,
442: $ A( 1 ), M1, B( M1, 0 ), LDB )
443: CALL ZGEMM( 'N', 'N', M1, N, M2, -CONE,
444: $ A( M1*M1 ), M1, B( M1, 0 ), LDB,
445: $ ALPHA, B, LDB )
446: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
447: $ A( 0 ), M1, B, LDB )
448: END IF
449: *
450: END IF
451: *
452: ELSE
453: *
454: * SIDE ='L', N is odd, TRANSR = 'C', and UPLO = 'U'
455: *
456: IF( .NOT.NOTRANS ) THEN
457: *
458: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'U', and
459: * TRANS = 'N'
460: *
461: CALL ZTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
462: $ A( M2*M2 ), M2, B, LDB )
463: CALL ZGEMM( 'N', 'N', M2, N, M1, -CONE, A( 0 ), M2,
464: $ B, LDB, ALPHA, B( M1, 0 ), LDB )
465: CALL ZTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
466: $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
467: *
468: ELSE
469: *
470: * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'U', and
471: * TRANS = 'C'
472: *
473: CALL ZTRSM( 'L', 'L', 'C', DIAG, M2, N, ALPHA,
474: $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
475: CALL ZGEMM( 'C', 'N', M1, N, M2, -CONE, A( 0 ), M2,
476: $ B( M1, 0 ), LDB, ALPHA, B, LDB )
477: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
478: $ A( M2*M2 ), M2, B, LDB )
479: *
480: END IF
481: *
482: END IF
483: *
484: END IF
485: *
486: ELSE
487: *
488: * SIDE = 'L' and N is even
489: *
490: IF( NORMALTRANSR ) THEN
491: *
492: * SIDE = 'L', N is even, and TRANSR = 'N'
493: *
494: IF( LOWER ) THEN
495: *
496: * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L'
497: *
498: IF( NOTRANS ) THEN
499: *
500: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
501: * and TRANS = 'N'
502: *
503: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
504: $ A( 1 ), M+1, B, LDB )
505: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( K+1 ),
506: $ M+1, B, LDB, ALPHA, B( K, 0 ), LDB )
507: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
508: $ A( 0 ), M+1, B( K, 0 ), LDB )
509: *
510: ELSE
511: *
512: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
513: * and TRANS = 'C'
514: *
515: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
516: $ A( 0 ), M+1, B( K, 0 ), LDB )
517: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( K+1 ),
518: $ M+1, B( K, 0 ), LDB, ALPHA, B, LDB )
519: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
520: $ A( 1 ), M+1, B, LDB )
521: *
522: END IF
523: *
524: ELSE
525: *
526: * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U'
527: *
528: IF( .NOT.NOTRANS ) THEN
529: *
530: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
531: * and TRANS = 'N'
532: *
533: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
534: $ A( K+1 ), M+1, B, LDB )
535: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), M+1,
536: $ B, LDB, ALPHA, B( K, 0 ), LDB )
537: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
538: $ A( K ), M+1, B( K, 0 ), LDB )
539: *
540: ELSE
541: *
542: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
543: * and TRANS = 'C'
544: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
545: $ A( K ), M+1, B( K, 0 ), LDB )
546: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), M+1,
547: $ B( K, 0 ), LDB, ALPHA, B, LDB )
548: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
549: $ A( K+1 ), M+1, B, LDB )
550: *
551: END IF
552: *
553: END IF
554: *
555: ELSE
556: *
557: * SIDE = 'L', N is even, and TRANSR = 'C'
558: *
559: IF( LOWER ) THEN
560: *
561: * SIDE ='L', N is even, TRANSR = 'C', and UPLO = 'L'
562: *
563: IF( NOTRANS ) THEN
564: *
565: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'L',
566: * and TRANS = 'N'
567: *
568: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, ALPHA,
569: $ A( K ), K, B, LDB )
570: CALL ZGEMM( 'C', 'N', K, N, K, -CONE,
571: $ A( K*( K+1 ) ), K, B, LDB, ALPHA,
572: $ B( K, 0 ), LDB )
573: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
574: $ A( 0 ), K, B( K, 0 ), LDB )
575: *
576: ELSE
577: *
578: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'L',
579: * and TRANS = 'C'
580: *
581: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, ALPHA,
582: $ A( 0 ), K, B( K, 0 ), LDB )
583: CALL ZGEMM( 'N', 'N', K, N, K, -CONE,
584: $ A( K*( K+1 ) ), K, B( K, 0 ), LDB,
585: $ ALPHA, B, LDB )
586: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
587: $ A( K ), K, B, LDB )
588: *
589: END IF
590: *
591: ELSE
592: *
593: * SIDE ='L', N is even, TRANSR = 'C', and UPLO = 'U'
594: *
595: IF( .NOT.NOTRANS ) THEN
596: *
597: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'U',
598: * and TRANS = 'N'
599: *
600: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, ALPHA,
601: $ A( K*( K+1 ) ), K, B, LDB )
602: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), K, B,
603: $ LDB, ALPHA, B( K, 0 ), LDB )
604: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
605: $ A( K*K ), K, B( K, 0 ), LDB )
606: *
607: ELSE
608: *
609: * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'U',
610: * and TRANS = 'C'
611: *
612: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, ALPHA,
613: $ A( K*K ), K, B( K, 0 ), LDB )
614: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), K,
615: $ B( K, 0 ), LDB, ALPHA, B, LDB )
616: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
617: $ A( K*( K+1 ) ), K, B, LDB )
618: *
619: END IF
620: *
621: END IF
622: *
623: END IF
624: *
625: END IF
626: *
627: ELSE
628: *
629: * SIDE = 'R'
630: *
631: * A is N-by-N.
632: * If N is odd, set NISODD = .TRUE., and N1 and N2.
633: * If N is even, NISODD = .FALSE., and K.
634: *
635: IF( MOD( N, 2 ).EQ.0 ) THEN
636: NISODD = .FALSE.
637: K = N / 2
638: ELSE
639: NISODD = .TRUE.
640: IF( LOWER ) THEN
641: N2 = N / 2
642: N1 = N - N2
643: ELSE
644: N1 = N / 2
645: N2 = N - N1
646: END IF
647: END IF
648: *
649: IF( NISODD ) THEN
650: *
651: * SIDE = 'R' and N is odd
652: *
653: IF( NORMALTRANSR ) THEN
654: *
655: * SIDE = 'R', N is odd, and TRANSR = 'N'
656: *
657: IF( LOWER ) THEN
658: *
659: * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
660: *
661: IF( NOTRANS ) THEN
662: *
663: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
664: * TRANS = 'N'
665: *
666: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N2, ALPHA,
667: $ A( N ), N, B( 0, N1 ), LDB )
668: CALL ZGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
669: $ LDB, A( N1 ), N, ALPHA, B( 0, 0 ),
670: $ LDB )
671: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
672: $ A( 0 ), N, B( 0, 0 ), LDB )
673: *
674: ELSE
675: *
676: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
677: * TRANS = 'C'
678: *
679: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N1, ALPHA,
680: $ A( 0 ), N, B( 0, 0 ), LDB )
681: CALL ZGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
682: $ LDB, A( N1 ), N, ALPHA, B( 0, N1 ),
683: $ LDB )
684: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
685: $ A( N ), N, B( 0, N1 ), LDB )
686: *
687: END IF
688: *
689: ELSE
690: *
691: * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
692: *
693: IF( NOTRANS ) THEN
694: *
695: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
696: * TRANS = 'N'
697: *
698: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N1, ALPHA,
699: $ A( N2 ), N, B( 0, 0 ), LDB )
700: CALL ZGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
701: $ LDB, A( 0 ), N, ALPHA, B( 0, N1 ),
702: $ LDB )
703: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
704: $ A( N1 ), N, B( 0, N1 ), LDB )
705: *
706: ELSE
707: *
708: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
709: * TRANS = 'C'
710: *
711: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N2, ALPHA,
712: $ A( N1 ), N, B( 0, N1 ), LDB )
713: CALL ZGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
714: $ LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB )
715: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
716: $ A( N2 ), N, B( 0, 0 ), LDB )
717: *
718: END IF
719: *
720: END IF
721: *
722: ELSE
723: *
724: * SIDE = 'R', N is odd, and TRANSR = 'C'
725: *
726: IF( LOWER ) THEN
727: *
728: * SIDE ='R', N is odd, TRANSR = 'C', and UPLO = 'L'
729: *
730: IF( NOTRANS ) THEN
731: *
732: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'L', and
733: * TRANS = 'N'
734: *
735: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
736: $ A( 1 ), N1, B( 0, N1 ), LDB )
737: CALL ZGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
738: $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ),
739: $ LDB )
740: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
741: $ A( 0 ), N1, B( 0, 0 ), LDB )
742: *
743: ELSE
744: *
745: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'L', and
746: * TRANS = 'C'
747: *
748: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
749: $ A( 0 ), N1, B( 0, 0 ), LDB )
750: CALL ZGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
751: $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ),
752: $ LDB )
753: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
754: $ A( 1 ), N1, B( 0, N1 ), LDB )
755: *
756: END IF
757: *
758: ELSE
759: *
760: * SIDE ='R', N is odd, TRANSR = 'C', and UPLO = 'U'
761: *
762: IF( NOTRANS ) THEN
763: *
764: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'U', and
765: * TRANS = 'N'
766: *
767: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
768: $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
769: CALL ZGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
770: $ LDB, A( 0 ), N2, ALPHA, B( 0, N1 ),
771: $ LDB )
772: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
773: $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
774: *
775: ELSE
776: *
777: * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'U', and
778: * TRANS = 'C'
779: *
780: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
781: $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
782: CALL ZGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
783: $ LDB, A( 0 ), N2, ALPHA, B( 0, 0 ),
784: $ LDB )
785: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
786: $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
787: *
788: END IF
789: *
790: END IF
791: *
792: END IF
793: *
794: ELSE
795: *
796: * SIDE = 'R' and N is even
797: *
798: IF( NORMALTRANSR ) THEN
799: *
800: * SIDE = 'R', N is even, and TRANSR = 'N'
801: *
802: IF( LOWER ) THEN
803: *
804: * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L'
805: *
806: IF( NOTRANS ) THEN
807: *
808: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
809: * and TRANS = 'N'
810: *
811: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, ALPHA,
812: $ A( 0 ), N+1, B( 0, K ), LDB )
813: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
814: $ LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ),
815: $ LDB )
816: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
817: $ A( 1 ), N+1, B( 0, 0 ), LDB )
818: *
819: ELSE
820: *
821: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
822: * and TRANS = 'C'
823: *
824: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, ALPHA,
825: $ A( 1 ), N+1, B( 0, 0 ), LDB )
826: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
827: $ LDB, A( K+1 ), N+1, ALPHA, B( 0, K ),
828: $ LDB )
829: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
830: $ A( 0 ), N+1, B( 0, K ), LDB )
831: *
832: END IF
833: *
834: ELSE
835: *
836: * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U'
837: *
838: IF( NOTRANS ) THEN
839: *
840: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
841: * and TRANS = 'N'
842: *
843: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, ALPHA,
844: $ A( K+1 ), N+1, B( 0, 0 ), LDB )
845: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
846: $ LDB, A( 0 ), N+1, ALPHA, B( 0, K ),
847: $ LDB )
848: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
849: $ A( K ), N+1, B( 0, K ), LDB )
850: *
851: ELSE
852: *
853: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
854: * and TRANS = 'C'
855: *
856: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, ALPHA,
857: $ A( K ), N+1, B( 0, K ), LDB )
858: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
859: $ LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ),
860: $ LDB )
861: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
862: $ A( K+1 ), N+1, B( 0, 0 ), LDB )
863: *
864: END IF
865: *
866: END IF
867: *
868: ELSE
869: *
870: * SIDE = 'R', N is even, and TRANSR = 'C'
871: *
872: IF( LOWER ) THEN
873: *
874: * SIDE ='R', N is even, TRANSR = 'C', and UPLO = 'L'
875: *
876: IF( NOTRANS ) THEN
877: *
878: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'L',
879: * and TRANS = 'N'
880: *
881: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
882: $ A( 0 ), K, B( 0, K ), LDB )
883: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
884: $ LDB, A( ( K+1 )*K ), K, ALPHA,
885: $ B( 0, 0 ), LDB )
886: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
887: $ A( K ), K, B( 0, 0 ), LDB )
888: *
889: ELSE
890: *
891: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'L',
892: * and TRANS = 'C'
893: *
894: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
895: $ A( K ), K, B( 0, 0 ), LDB )
896: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
897: $ LDB, A( ( K+1 )*K ), K, ALPHA,
898: $ B( 0, K ), LDB )
899: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
900: $ A( 0 ), K, B( 0, K ), LDB )
901: *
902: END IF
903: *
904: ELSE
905: *
906: * SIDE ='R', N is even, TRANSR = 'C', and UPLO = 'U'
907: *
908: IF( NOTRANS ) THEN
909: *
910: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'U',
911: * and TRANS = 'N'
912: *
913: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
914: $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
915: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
916: $ LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB )
917: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
918: $ A( K*K ), K, B( 0, K ), LDB )
919: *
920: ELSE
921: *
922: * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'U',
923: * and TRANS = 'C'
924: *
925: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
926: $ A( K*K ), K, B( 0, K ), LDB )
927: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
928: $ LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB )
929: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
930: $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
931: *
932: END IF
933: *
934: END IF
935: *
936: END IF
937: *
938: END IF
939: END IF
940: *
941: RETURN
942: *
943: * End of ZTFSM
944: *
945: 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