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