Annotation of rpl/lapack/lapack/ztfsm.f, revision 1.17
1.10 bertrand 1: *> \brief \b ZTFSM solves a matrix equation (one operand is a triangular matrix in RFP format).
1.7 bertrand 2: *
3: * =========== DOCUMENTATION ===========
4: *
1.14 bertrand 5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
1.7 bertrand 7: *
8: *> \htmlonly
1.14 bertrand 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">
1.7 bertrand 15: *> [TXT]</a>
1.14 bertrand 16: *> \endhtmlonly
1.7 bertrand 17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
22: * B, LDB )
1.14 bertrand 23: *
1.7 bertrand 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: * ..
1.14 bertrand 32: *
1.7 bertrand 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: *
1.14 bertrand 180: *> \author Univ. of Tennessee
181: *> \author Univ. of California Berkeley
182: *> \author Univ. of Colorado Denver
183: *> \author NAG Ltd.
1.7 bertrand 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: * =====================================================================
1.1 bertrand 296: SUBROUTINE ZTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
1.6 bertrand 297: $ B, LDB )
1.1 bertrand 298: *
1.17 ! bertrand 299: * -- LAPACK computational routine --
1.1 bertrand 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: *
1.6 bertrand 312: * =====================================================================
1.1 bertrand 313: * ..
314: * .. Parameters ..
315: COMPLEX*16 CONE, CZERO
316: PARAMETER ( CONE = ( 1.0D+0, 0.0D+0 ),
1.6 bertrand 317: $ CZERO = ( 0.0D+0, 0.0D+0 ) )
1.1 bertrand 318: * ..
319: * .. Local Scalars ..
320: LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
1.6 bertrand 321: $ NOTRANS
1.1 bertrand 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' ) )
1.6 bertrand 352: $ THEN
1.1 bertrand 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 ) )
1.6 bertrand 369: $ RETURN
1.1 bertrand 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,
1.6 bertrand 423: $ A, M, B, LDB )
1.1 bertrand 424: ELSE
425: CALL ZTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
1.6 bertrand 426: $ A( 0 ), M, B, LDB )
1.1 bertrand 427: CALL ZGEMM( 'N', 'N', M2, N, M1, -CONE, A( M1 ),
1.6 bertrand 428: $ M, B, LDB, ALPHA, B( M1, 0 ), LDB )
1.1 bertrand 429: CALL ZTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
1.6 bertrand 430: $ A( M ), M, B( M1, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 440: $ A( 0 ), M, B, LDB )
1.1 bertrand 441: ELSE
442: CALL ZTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
1.6 bertrand 443: $ A( M ), M, B( M1, 0 ), LDB )
1.1 bertrand 444: CALL ZGEMM( 'C', 'N', M1, N, M2, -CONE, A( M1 ),
1.6 bertrand 445: $ M, B( M1, 0 ), LDB, ALPHA, B, LDB )
1.1 bertrand 446: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
1.6 bertrand 447: $ A( 0 ), M, B, LDB )
1.1 bertrand 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,
1.6 bertrand 462: $ A( M2 ), M, B, LDB )
1.1 bertrand 463: CALL ZGEMM( 'C', 'N', M2, N, M1, -CONE, A( 0 ), M,
1.6 bertrand 464: $ B, LDB, ALPHA, B( M1, 0 ), LDB )
1.1 bertrand 465: CALL ZTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
1.6 bertrand 466: $ A( M1 ), M, B( M1, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 474: $ A( M1 ), M, B( M1, 0 ), LDB )
1.1 bertrand 475: CALL ZGEMM( 'N', 'N', M1, N, M2, -CONE, A( 0 ), M,
1.6 bertrand 476: $ B( M1, 0 ), LDB, ALPHA, B, LDB )
1.1 bertrand 477: CALL ZTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
1.6 bertrand 478: $ A( M2 ), M, B, LDB )
1.1 bertrand 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,
1.6 bertrand 499: $ A( 0 ), M1, B, LDB )
1.1 bertrand 500: ELSE
501: CALL ZTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
1.6 bertrand 502: $ A( 0 ), M1, B, LDB )
1.1 bertrand 503: CALL ZGEMM( 'C', 'N', M2, N, M1, -CONE,
1.6 bertrand 504: $ A( M1*M1 ), M1, B, LDB, ALPHA,
505: $ B( M1, 0 ), LDB )
1.1 bertrand 506: CALL ZTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
1.6 bertrand 507: $ A( 1 ), M1, B( M1, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 517: $ A( 0 ), M1, B, LDB )
1.1 bertrand 518: ELSE
519: CALL ZTRSM( 'L', 'L', 'C', DIAG, M2, N, ALPHA,
1.6 bertrand 520: $ A( 1 ), M1, B( M1, 0 ), LDB )
1.1 bertrand 521: CALL ZGEMM( 'N', 'N', M1, N, M2, -CONE,
1.6 bertrand 522: $ A( M1*M1 ), M1, B( M1, 0 ), LDB,
523: $ ALPHA, B, LDB )
1.1 bertrand 524: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
1.6 bertrand 525: $ A( 0 ), M1, B, LDB )
1.1 bertrand 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,
1.6 bertrand 540: $ A( M2*M2 ), M2, B, LDB )
1.1 bertrand 541: CALL ZGEMM( 'N', 'N', M2, N, M1, -CONE, A( 0 ), M2,
1.6 bertrand 542: $ B, LDB, ALPHA, B( M1, 0 ), LDB )
1.1 bertrand 543: CALL ZTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
1.6 bertrand 544: $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 552: $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
1.1 bertrand 553: CALL ZGEMM( 'C', 'N', M1, N, M2, -CONE, A( 0 ), M2,
1.6 bertrand 554: $ B( M1, 0 ), LDB, ALPHA, B, LDB )
1.1 bertrand 555: CALL ZTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
1.6 bertrand 556: $ A( M2*M2 ), M2, B, LDB )
1.1 bertrand 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,
1.6 bertrand 582: $ A( 1 ), M+1, B, LDB )
1.1 bertrand 583: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( K+1 ),
1.6 bertrand 584: $ M+1, B, LDB, ALPHA, B( K, 0 ), LDB )
1.1 bertrand 585: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
1.6 bertrand 586: $ A( 0 ), M+1, B( K, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 594: $ A( 0 ), M+1, B( K, 0 ), LDB )
1.1 bertrand 595: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( K+1 ),
1.6 bertrand 596: $ M+1, B( K, 0 ), LDB, ALPHA, B, LDB )
1.1 bertrand 597: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
1.6 bertrand 598: $ A( 1 ), M+1, B, LDB )
1.1 bertrand 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,
1.6 bertrand 612: $ A( K+1 ), M+1, B, LDB )
1.1 bertrand 613: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), M+1,
1.6 bertrand 614: $ B, LDB, ALPHA, B( K, 0 ), LDB )
1.1 bertrand 615: CALL ZTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
1.6 bertrand 616: $ A( K ), M+1, B( K, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 623: $ A( K ), M+1, B( K, 0 ), LDB )
1.1 bertrand 624: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), M+1,
1.6 bertrand 625: $ B( K, 0 ), LDB, ALPHA, B, LDB )
1.1 bertrand 626: CALL ZTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
1.6 bertrand 627: $ A( K+1 ), M+1, B, LDB )
1.1 bertrand 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,
1.6 bertrand 647: $ A( K ), K, B, LDB )
1.1 bertrand 648: CALL ZGEMM( 'C', 'N', K, N, K, -CONE,
1.6 bertrand 649: $ A( K*( K+1 ) ), K, B, LDB, ALPHA,
650: $ B( K, 0 ), LDB )
1.1 bertrand 651: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
1.6 bertrand 652: $ A( 0 ), K, B( K, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 660: $ A( 0 ), K, B( K, 0 ), LDB )
1.1 bertrand 661: CALL ZGEMM( 'N', 'N', K, N, K, -CONE,
1.6 bertrand 662: $ A( K*( K+1 ) ), K, B( K, 0 ), LDB,
663: $ ALPHA, B, LDB )
1.1 bertrand 664: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
1.6 bertrand 665: $ A( K ), K, B, LDB )
1.1 bertrand 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,
1.6 bertrand 679: $ A( K*( K+1 ) ), K, B, LDB )
1.1 bertrand 680: CALL ZGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), K, B,
1.6 bertrand 681: $ LDB, ALPHA, B( K, 0 ), LDB )
1.1 bertrand 682: CALL ZTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
1.6 bertrand 683: $ A( K*K ), K, B( K, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 691: $ A( K*K ), K, B( K, 0 ), LDB )
1.1 bertrand 692: CALL ZGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), K,
1.6 bertrand 693: $ B( K, 0 ), LDB, ALPHA, B, LDB )
1.1 bertrand 694: CALL ZTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
1.6 bertrand 695: $ A( K*( K+1 ) ), K, B, LDB )
1.1 bertrand 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,
1.6 bertrand 745: $ A( N ), N, B( 0, N1 ), LDB )
1.1 bertrand 746: CALL ZGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
1.6 bertrand 747: $ LDB, A( N1 ), N, ALPHA, B( 0, 0 ),
748: $ LDB )
1.1 bertrand 749: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
1.6 bertrand 750: $ A( 0 ), N, B( 0, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 758: $ A( 0 ), N, B( 0, 0 ), LDB )
1.1 bertrand 759: CALL ZGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
1.6 bertrand 760: $ LDB, A( N1 ), N, ALPHA, B( 0, N1 ),
761: $ LDB )
1.1 bertrand 762: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
1.6 bertrand 763: $ A( N ), N, B( 0, N1 ), LDB )
1.1 bertrand 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,
1.6 bertrand 777: $ A( N2 ), N, B( 0, 0 ), LDB )
1.1 bertrand 778: CALL ZGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
1.6 bertrand 779: $ LDB, A( 0 ), N, ALPHA, B( 0, N1 ),
780: $ LDB )
1.1 bertrand 781: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
1.6 bertrand 782: $ A( N1 ), N, B( 0, N1 ), LDB )
1.1 bertrand 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,
1.6 bertrand 790: $ A( N1 ), N, B( 0, N1 ), LDB )
1.1 bertrand 791: CALL ZGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
1.6 bertrand 792: $ LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB )
1.1 bertrand 793: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
1.6 bertrand 794: $ A( N2 ), N, B( 0, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 814: $ A( 1 ), N1, B( 0, N1 ), LDB )
1.1 bertrand 815: CALL ZGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
1.6 bertrand 816: $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ),
817: $ LDB )
1.1 bertrand 818: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
1.6 bertrand 819: $ A( 0 ), N1, B( 0, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 827: $ A( 0 ), N1, B( 0, 0 ), LDB )
1.1 bertrand 828: CALL ZGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
1.6 bertrand 829: $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ),
830: $ LDB )
1.1 bertrand 831: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
1.6 bertrand 832: $ A( 1 ), N1, B( 0, N1 ), LDB )
1.1 bertrand 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,
1.6 bertrand 846: $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
1.1 bertrand 847: CALL ZGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
1.6 bertrand 848: $ LDB, A( 0 ), N2, ALPHA, B( 0, N1 ),
849: $ LDB )
1.1 bertrand 850: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
1.6 bertrand 851: $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
1.1 bertrand 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,
1.6 bertrand 859: $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
1.1 bertrand 860: CALL ZGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
1.6 bertrand 861: $ LDB, A( 0 ), N2, ALPHA, B( 0, 0 ),
862: $ LDB )
1.1 bertrand 863: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
1.6 bertrand 864: $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 890: $ A( 0 ), N+1, B( 0, K ), LDB )
1.1 bertrand 891: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
1.6 bertrand 892: $ LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ),
893: $ LDB )
1.1 bertrand 894: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
1.6 bertrand 895: $ A( 1 ), N+1, B( 0, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 903: $ A( 1 ), N+1, B( 0, 0 ), LDB )
1.1 bertrand 904: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
1.6 bertrand 905: $ LDB, A( K+1 ), N+1, ALPHA, B( 0, K ),
906: $ LDB )
1.1 bertrand 907: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
1.6 bertrand 908: $ A( 0 ), N+1, B( 0, K ), LDB )
1.1 bertrand 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,
1.6 bertrand 922: $ A( K+1 ), N+1, B( 0, 0 ), LDB )
1.1 bertrand 923: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
1.6 bertrand 924: $ LDB, A( 0 ), N+1, ALPHA, B( 0, K ),
925: $ LDB )
1.1 bertrand 926: CALL ZTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
1.6 bertrand 927: $ A( K ), N+1, B( 0, K ), LDB )
1.1 bertrand 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,
1.6 bertrand 935: $ A( K ), N+1, B( 0, K ), LDB )
1.1 bertrand 936: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
1.6 bertrand 937: $ LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ),
938: $ LDB )
1.1 bertrand 939: CALL ZTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
1.6 bertrand 940: $ A( K+1 ), N+1, B( 0, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 960: $ A( 0 ), K, B( 0, K ), LDB )
1.1 bertrand 961: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
1.6 bertrand 962: $ LDB, A( ( K+1 )*K ), K, ALPHA,
963: $ B( 0, 0 ), LDB )
1.1 bertrand 964: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
1.6 bertrand 965: $ A( K ), K, B( 0, 0 ), LDB )
1.1 bertrand 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,
1.6 bertrand 973: $ A( K ), K, B( 0, 0 ), LDB )
1.1 bertrand 974: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
1.6 bertrand 975: $ LDB, A( ( K+1 )*K ), K, ALPHA,
976: $ B( 0, K ), LDB )
1.1 bertrand 977: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
1.6 bertrand 978: $ A( 0 ), K, B( 0, K ), LDB )
1.1 bertrand 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,
1.6 bertrand 992: $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
1.1 bertrand 993: CALL ZGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
1.6 bertrand 994: $ LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB )
1.1 bertrand 995: CALL ZTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
1.6 bertrand 996: $ A( K*K ), K, B( 0, K ), LDB )
1.1 bertrand 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,
1.6 bertrand 1004: $ A( K*K ), K, B( 0, K ), LDB )
1.1 bertrand 1005: CALL ZGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
1.6 bertrand 1006: $ LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB )
1.1 bertrand 1007: CALL ZTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
1.6 bertrand 1008: $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
1.1 bertrand 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>