001:       SUBROUTINE STFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
002:      +                  B, LDB )
003: *
004: *  -- LAPACK routine (version 3.2)                                    --
005: *
006: *  -- Contributed by Fred Gustavson of the IBM Watson Research Center --
007: *  -- November 2008                                                   --
008: *
009: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
010: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
011: *
012: *     ..
013: *     .. Scalar Arguments ..
014:       CHARACTER          TRANSR, DIAG, SIDE, TRANS, UPLO
015:       INTEGER            LDB, M, N
016:       REAL               ALPHA
017: *     ..
018: *     .. Array Arguments ..
019:       REAL               A( 0: * ), B( 0: LDB-1, 0: * )
020: *     ..
021: *
022: *  Purpose
023: *  =======
024: *
025: *  Level 3 BLAS like routine for A in RFP Format.
026: *
027: *  STFSM  solves the matrix equation
028: *
029: *     op( A )*X = alpha*B  or  X*op( A ) = alpha*B
030: *
031: *  where alpha is a scalar, X and B are m by n matrices, A is a unit, or
032: *  non-unit,  upper or lower triangular matrix  and  op( A )  is one  of
033: *
034: *     op( A ) = A   or   op( A ) = A'.
035: *
036: *  A is in Rectangular Full Packed (RFP) Format.
037: *
038: *  The matrix X is overwritten on B.
039: *
040: *  Arguments
041: *  ==========
042: *
043: *  TRANSR - (input) CHARACTER
044: *          = 'N':  The Normal Form of RFP A is stored;
045: *          = 'T':  The Transpose Form of RFP A is stored.
046: *
047: *  SIDE   - (input) CHARACTER
048: *           On entry, SIDE specifies whether op( A ) appears on the left
049: *           or right of X as follows:
050: *
051: *              SIDE = 'L' or 'l'   op( A )*X = alpha*B.
052: *
053: *              SIDE = 'R' or 'r'   X*op( A ) = alpha*B.
054: *
055: *           Unchanged on exit.
056: *
057: *  UPLO   - (input) CHARACTER
058: *           On entry, UPLO specifies whether the RFP matrix A came from
059: *           an upper or lower triangular matrix as follows:
060: *           UPLO = 'U' or 'u' RFP A came from an upper triangular matrix
061: *           UPLO = 'L' or 'l' RFP A came from a  lower triangular matrix
062: *
063: *           Unchanged on exit.
064: *
065: *  TRANS  - (input) CHARACTER
066: *           On entry, TRANS  specifies the form of op( A ) to be used
067: *           in the matrix multiplication as follows:
068: *
069: *              TRANS  = 'N' or 'n'   op( A ) = A.
070: *
071: *              TRANS  = 'T' or 't'   op( A ) = A'.
072: *
073: *           Unchanged on exit.
074: *
075: *  DIAG   - (input) CHARACTER
076: *           On entry, DIAG specifies whether or not RFP A is unit
077: *           triangular as follows:
078: *
079: *              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
080: *
081: *              DIAG = 'N' or 'n'   A is not assumed to be unit
082: *                                  triangular.
083: *
084: *           Unchanged on exit.
085: *
086: *  M      - (input) INTEGER.
087: *           On entry, M specifies the number of rows of B. M must be at
088: *           least zero.
089: *           Unchanged on exit.
090: *
091: *  N      - (input) INTEGER.
092: *           On entry, N specifies the number of columns of B.  N must be
093: *           at least zero.
094: *           Unchanged on exit.
095: *
096: *  ALPHA  - (input) REAL.
097: *           On entry,  ALPHA specifies the scalar  alpha. When  alpha is
098: *           zero then  A is not referenced and  B need not be set before
099: *           entry.
100: *           Unchanged on exit.
101: *
102: *  A      - (input) REAL array, dimension (NT);
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 = 'T' then RFP is the 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: *           transpose Format. If UPLO = 'L' the RFP A contains
112: *           the NT elements of lower packed A either in normal or
113: *           transpose Format. The LDA of RFP A is (N+1)/2 when
114: *           TRANSR = 'T'. 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/ouptut) REAL 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: *  Notes
130: *  =====
131: *
132: *  We first consider Rectangular Full Packed (RFP) Format when N is
133: *  even. 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: *  the 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: *  the transpose of the last three columns of AP lower.
152: *  This covers the case N even and TRANSR = 'N'.
153: *
154: *         RFP A                   RFP A
155: *
156: *        03 04 05                33 43 53
157: *        13 14 15                00 44 54
158: *        23 24 25                10 11 55
159: *        33 34 35                20 21 22
160: *        00 44 45                30 31 32
161: *        01 11 55                40 41 42
162: *        02 12 22                50 51 52
163: *
164: *  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
165: *  transpose of RFP A above. One therefore gets:
166: *
167: *
168: *           RFP A                   RFP A
169: *
170: *     03 13 23 33 00 01 02    33 00 10 20 30 40 50
171: *     04 14 24 34 44 11 12    43 44 11 21 31 41 51
172: *     05 15 25 35 45 55 22    53 54 55 22 32 42 52
173: *
174: *
175: *  We first consider Rectangular Full Packed (RFP) Format when N is
176: *  odd. We give an example where N = 5.
177: *
178: *     AP is Upper                 AP is Lower
179: *
180: *   00 01 02 03 04              00
181: *      11 12 13 14              10 11
182: *         22 23 24              20 21 22
183: *            33 34              30 31 32 33
184: *               44              40 41 42 43 44
185: *
186: *
187: *  Let TRANSR = 'N'. RFP holds AP as follows:
188: *  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
189: *  three columns of AP upper. The lower triangle A(3:4,0:1) consists of
190: *  the transpose of the first two columns of AP upper.
191: *  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
192: *  three columns of AP lower. The upper triangle A(0:1,1:2) consists of
193: *  the transpose of the last two columns of AP lower.
194: *  This covers the case N odd and TRANSR = 'N'.
195: *
196: *         RFP A                   RFP A
197: *
198: *        02 03 04                00 33 43
199: *        12 13 14                10 11 44
200: *        22 23 24                20 21 22
201: *        00 33 34                30 31 32
202: *        01 11 44                40 41 42
203: *
204: *  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
205: *  transpose of RFP A above. One therefore gets:
206: *
207: *           RFP A                   RFP A
208: *
209: *     02 12 22 00 01             00 10 20 30 40 50
210: *     03 13 23 33 11             33 11 21 31 41 51
211: *     04 14 24 34 44             43 44 22 32 42 52
212: *
213: *  Reference
214: *  =========
215: *
216: *  =====================================================================
217: *
218: *     ..
219: *     .. Parameters ..
220:       REAL               ONE, ZERO
221:       PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
222: *     ..
223: *     .. Local Scalars ..
224:       LOGICAL            LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
225:      +                   NOTRANS
226:       INTEGER            M1, M2, N1, N2, K, INFO, I, J
227: *     ..
228: *     .. External Functions ..
229:       LOGICAL            LSAME
230:       EXTERNAL           LSAME
231: *     ..
232: *     .. External Subroutines ..
233:       EXTERNAL           SGEMM, STRSM, XERBLA
234: *     ..
235: *     .. Intrinsic Functions ..
236:       INTRINSIC          MAX, MOD
237: *     ..
238: *     .. Executable Statements ..
239: *
240: *     Test the input parameters.
241: *
242:       INFO = 0
243:       NORMALTRANSR = LSAME( TRANSR, 'N' )
244:       LSIDE = LSAME( SIDE, 'L' )
245:       LOWER = LSAME( UPLO, 'L' )
246:       NOTRANS = LSAME( TRANS, 'N' )
247:       IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN
248:          INFO = -1
249:       ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
250:          INFO = -2
251:       ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
252:          INFO = -3
253:       ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN
254:          INFO = -4
255:       ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) )
256:      +         THEN
257:          INFO = -5
258:       ELSE IF( M.LT.0 ) THEN
259:          INFO = -6
260:       ELSE IF( N.LT.0 ) THEN
261:          INFO = -7
262:       ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
263:          INFO = -11
264:       END IF
265:       IF( INFO.NE.0 ) THEN
266:          CALL XERBLA( 'STFSM ', -INFO )
267:          RETURN
268:       END IF
269: *
270: *     Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
271: *
272:       IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
273:      +   RETURN
274: *
275: *     Quick return when ALPHA.EQ.(0D+0)
276: *
277:       IF( ALPHA.EQ.ZERO ) THEN
278:          DO 20 J = 0, N - 1
279:             DO 10 I = 0, M - 1
280:                B( I, J ) = ZERO
281:    10       CONTINUE
282:    20    CONTINUE
283:          RETURN
284:       END IF
285: *
286:       IF( LSIDE ) THEN
287: *
288: *        SIDE = 'L'
289: *
290: *        A is M-by-M.
291: *        If M is odd, set NISODD = .TRUE., and M1 and M2.
292: *        If M is even, NISODD = .FALSE., and M.
293: *
294:          IF( MOD( M, 2 ).EQ.0 ) THEN
295:             MISODD = .FALSE.
296:             K = M / 2
297:          ELSE
298:             MISODD = .TRUE.
299:             IF( LOWER ) THEN
300:                M2 = M / 2
301:                M1 = M - M2
302:             ELSE
303:                M1 = M / 2
304:                M2 = M - M1
305:             END IF
306:          END IF
307: *
308:          IF( MISODD ) THEN
309: *
310: *           SIDE = 'L' and N is odd
311: *
312:             IF( NORMALTRANSR ) THEN
313: *
314: *              SIDE = 'L', N is odd, and TRANSR = 'N'
315: *
316:                IF( LOWER ) THEN
317: *
318: *                 SIDE  ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
319: *
320:                   IF( NOTRANS ) THEN
321: *
322: *                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
323: *                    TRANS = 'N'
324: *
325:                      CALL STRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
326:      +                           A( 0 ), M, B, LDB )
327:                      CALL SGEMM( 'N', 'N', M2, N, M1, -ONE, A( M1 ), M,
328:      +                           B, LDB, ALPHA, B( M1, 0 ), LDB )
329:                      CALL STRSM( 'L', 'U', 'T', DIAG, M2, N, ONE,
330:      +                           A( M ), M, B( M1, 0 ), LDB )
331: *
332:                   ELSE
333: *
334: *                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
335: *                    TRANS = 'T'
336: *
337:                      CALL STRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
338:      +                           A( M ), M, B( M1, 0 ), LDB )
339:                      CALL SGEMM( 'T', 'N', M1, N, M2, -ONE, A( M1 ), M,
340:      +                           B( M1, 0 ), LDB, ALPHA, B, LDB )
341:                      CALL STRSM( 'L', 'L', 'T', DIAG, M1, N, ONE,
342:      +                           A( 0 ), M, B, LDB )
343: *
344:                   END IF
345: *
346:                ELSE
347: *
348: *                 SIDE  ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
349: *
350:                   IF( .NOT.NOTRANS ) THEN
351: *
352: *                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
353: *                    TRANS = 'N'
354: *
355:                      CALL STRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
356:      +                           A( M2 ), M, B, LDB )
357:                      CALL SGEMM( 'T', 'N', M2, N, M1, -ONE, A( 0 ), M,
358:      +                           B, LDB, ALPHA, B( M1, 0 ), LDB )
359:                      CALL STRSM( 'L', 'U', 'T', DIAG, M2, N, ONE,
360:      +                           A( M1 ), M, B( M1, 0 ), LDB )
361: *
362:                   ELSE
363: *
364: *                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
365: *                    TRANS = 'T'
366: *
367:                      CALL STRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
368:      +                           A( M1 ), M, B( M1, 0 ), LDB )
369:                      CALL SGEMM( 'N', 'N', M1, N, M2, -ONE, A( 0 ), M,
370:      +                           B( M1, 0 ), LDB, ALPHA, B, LDB )
371:                      CALL STRSM( 'L', 'L', 'T', DIAG, M1, N, ONE,
372:      +                           A( M2 ), M, B, LDB )
373: *
374:                   END IF
375: *
376:                END IF
377: *
378:             ELSE
379: *
380: *              SIDE = 'L', N is odd, and TRANSR = 'T'
381: *
382:                IF( LOWER ) THEN
383: *
384: *                 SIDE  ='L', N is odd, TRANSR = 'T', and UPLO = 'L'
385: *
386:                   IF( NOTRANS ) THEN
387: *
388: *                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
389: *                    TRANS = 'N'
390: *
391:                      CALL STRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
392:      +                           A( 0 ), M1, B, LDB )
393:                      CALL SGEMM( 'T', 'N', M2, N, M1, -ONE, A( M1*M1 ),
394:      +                           M1, B, LDB, ALPHA, B( M1, 0 ), LDB )
395:                      CALL STRSM( 'L', 'L', 'N', DIAG, M2, N, ONE,
396:      +                           A( 1 ), M1, B( M1, 0 ), LDB )
397: *
398:                   ELSE
399: *
400: *                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
401: *                    TRANS = 'T'
402: *
403:                      CALL STRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA,
404:      +                           A( 1 ), M1, B( M1, 0 ), LDB )
405:                      CALL SGEMM( 'N', 'N', M1, N, M2, -ONE, A( M1*M1 ),
406:      +                           M1, B( M1, 0 ), LDB, ALPHA, B, LDB )
407:                      CALL STRSM( 'L', 'U', 'N', DIAG, M1, N, ONE,
408:      +                           A( 0 ), M1, B, LDB )
409: *
410:                   END IF
411: *
412:                ELSE
413: *
414: *                 SIDE  ='L', N is odd, TRANSR = 'T', and UPLO = 'U'
415: *
416:                   IF( .NOT.NOTRANS ) THEN
417: *
418: *                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
419: *                    TRANS = 'N'
420: *
421:                      CALL STRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
422:      +                           A( M2*M2 ), M2, B, LDB )
423:                      CALL SGEMM( 'N', 'N', M2, N, M1, -ONE, A( 0 ), M2,
424:      +                           B, LDB, ALPHA, B( M1, 0 ), LDB )
425:                      CALL STRSM( 'L', 'L', 'N', DIAG, M2, N, ONE,
426:      +                           A( M1*M2 ), M2, B( M1, 0 ), LDB )
427: *
428:                   ELSE
429: *
430: *                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
431: *                    TRANS = 'T'
432: *
433:                      CALL STRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA,
434:      +                           A( M1*M2 ), M2, B( M1, 0 ), LDB )
435:                      CALL SGEMM( 'T', 'N', M1, N, M2, -ONE, A( 0 ), M2,
436:      +                           B( M1, 0 ), LDB, ALPHA, B, LDB )
437:                      CALL STRSM( 'L', 'U', 'N', DIAG, M1, N, ONE,
438:      +                           A( M2*M2 ), M2, B, LDB )
439: *
440:                   END IF
441: *
442:                END IF
443: *
444:             END IF
445: *
446:          ELSE
447: *
448: *           SIDE = 'L' and N is even
449: *
450:             IF( NORMALTRANSR ) THEN
451: *
452: *              SIDE = 'L', N is even, and TRANSR = 'N'
453: *
454:                IF( LOWER ) THEN
455: *
456: *                 SIDE  ='L', N is even, TRANSR = 'N', and UPLO = 'L'
457: *
458:                   IF( NOTRANS ) THEN
459: *
460: *                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'L',
461: *                    and TRANS = 'N'
462: *
463:                      CALL STRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
464:      +                           A( 1 ), M+1, B, LDB )
465:                      CALL SGEMM( 'N', 'N', K, N, K, -ONE, A( K+1 ),
466:      +                           M+1, B, LDB, ALPHA, B( K, 0 ), LDB )
467:                      CALL STRSM( 'L', 'U', 'T', DIAG, K, N, ONE,
468:      +                           A( 0 ), M+1, B( K, 0 ), LDB )
469: *
470:                   ELSE
471: *
472: *                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'L',
473: *                    and TRANS = 'T'
474: *
475:                      CALL STRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
476:      +                           A( 0 ), M+1, B( K, 0 ), LDB )
477:                      CALL SGEMM( 'T', 'N', K, N, K, -ONE, A( K+1 ),
478:      +                           M+1, B( K, 0 ), LDB, ALPHA, B, LDB )
479:                      CALL STRSM( 'L', 'L', 'T', DIAG, K, N, ONE,
480:      +                           A( 1 ), M+1, B, LDB )
481: *
482:                   END IF
483: *
484:                ELSE
485: *
486: *                 SIDE  ='L', N is even, TRANSR = 'N', and UPLO = 'U'
487: *
488:                   IF( .NOT.NOTRANS ) THEN
489: *
490: *                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'U',
491: *                    and TRANS = 'N'
492: *
493:                      CALL STRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
494:      +                           A( K+1 ), M+1, B, LDB )
495:                      CALL SGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), M+1,
496:      +                           B, LDB, ALPHA, B( K, 0 ), LDB )
497:                      CALL STRSM( 'L', 'U', 'T', DIAG, K, N, ONE,
498:      +                           A( K ), M+1, B( K, 0 ), LDB )
499: *
500:                   ELSE
501: *
502: *                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'U',
503: *                    and TRANS = 'T'
504:                      CALL STRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
505:      +                           A( K ), M+1, B( K, 0 ), LDB )
506:                      CALL SGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), M+1,
507:      +                           B( K, 0 ), LDB, ALPHA, B, LDB )
508:                      CALL STRSM( 'L', 'L', 'T', DIAG, K, N, ONE,
509:      +                           A( K+1 ), M+1, B, LDB )
510: *
511:                   END IF
512: *
513:                END IF
514: *
515:             ELSE
516: *
517: *              SIDE = 'L', N is even, and TRANSR = 'T'
518: *
519:                IF( LOWER ) THEN
520: *
521: *                 SIDE  ='L', N is even, TRANSR = 'T', and UPLO = 'L'
522: *
523:                   IF( NOTRANS ) THEN
524: *
525: *                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'L',
526: *                    and TRANS = 'N'
527: *
528:                      CALL STRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA,
529:      +                           A( K ), K, B, LDB )
530:                      CALL SGEMM( 'T', 'N', K, N, K, -ONE,
531:      +                           A( K*( K+1 ) ), K, B, LDB, ALPHA,
532:      +                           B( K, 0 ), LDB )
533:                      CALL STRSM( 'L', 'L', 'N', DIAG, K, N, ONE,
534:      +                           A( 0 ), K, B( K, 0 ), LDB )
535: *
536:                   ELSE
537: *
538: *                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'L',
539: *                    and TRANS = 'T'
540: *
541:                      CALL STRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA,
542:      +                           A( 0 ), K, B( K, 0 ), LDB )
543:                      CALL SGEMM( 'N', 'N', K, N, K, -ONE,
544:      +                           A( K*( K+1 ) ), K, B( K, 0 ), LDB,
545:      +                           ALPHA, B, LDB )
546:                      CALL STRSM( 'L', 'U', 'N', DIAG, K, N, ONE,
547:      +                           A( K ), K, B, LDB )
548: *
549:                   END IF
550: *
551:                ELSE
552: *
553: *                 SIDE  ='L', N is even, TRANSR = 'T', and UPLO = 'U'
554: *
555:                   IF( .NOT.NOTRANS ) THEN
556: *
557: *                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'U',
558: *                    and TRANS = 'N'
559: *
560:                      CALL STRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA,
561:      +                           A( K*( K+1 ) ), K, B, LDB )
562:                      CALL SGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), K, B,
563:      +                           LDB, ALPHA, B( K, 0 ), LDB )
564:                      CALL STRSM( 'L', 'L', 'N', DIAG, K, N, ONE,
565:      +                           A( K*K ), K, B( K, 0 ), LDB )
566: *
567:                   ELSE
568: *
569: *                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'U',
570: *                    and TRANS = 'T'
571: *
572:                      CALL STRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA,
573:      +                           A( K*K ), K, B( K, 0 ), LDB )
574:                      CALL SGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), K,
575:      +                           B( K, 0 ), LDB, ALPHA, B, LDB )
576:                      CALL STRSM( 'L', 'U', 'N', DIAG, K, N, ONE,
577:      +                           A( K*( K+1 ) ), K, B, LDB )
578: *
579:                   END IF
580: *
581:                END IF
582: *
583:             END IF
584: *
585:          END IF
586: *
587:       ELSE
588: *
589: *        SIDE = 'R'
590: *
591: *        A is N-by-N.
592: *        If N is odd, set NISODD = .TRUE., and N1 and N2.
593: *        If N is even, NISODD = .FALSE., and K.
594: *
595:          IF( MOD( N, 2 ).EQ.0 ) THEN
596:             NISODD = .FALSE.
597:             K = N / 2
598:          ELSE
599:             NISODD = .TRUE.
600:             IF( LOWER ) THEN
601:                N2 = N / 2
602:                N1 = N - N2
603:             ELSE
604:                N1 = N / 2
605:                N2 = N - N1
606:             END IF
607:          END IF
608: *
609:          IF( NISODD ) THEN
610: *
611: *           SIDE = 'R' and N is odd
612: *
613:             IF( NORMALTRANSR ) THEN
614: *
615: *              SIDE = 'R', N is odd, and TRANSR = 'N'
616: *
617:                IF( LOWER ) THEN
618: *
619: *                 SIDE  ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
620: *
621:                   IF( NOTRANS ) THEN
622: *
623: *                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
624: *                    TRANS = 'N'
625: *
626:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA,
627:      +                           A( N ), N, B( 0, N1 ), LDB )
628:                      CALL SGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ),
629:      +                           LDB, A( N1 ), N, ALPHA, B( 0, 0 ),
630:      +                           LDB )
631:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, N1, ONE,
632:      +                           A( 0 ), N, B( 0, 0 ), LDB )
633: *
634:                   ELSE
635: *
636: *                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
637: *                    TRANS = 'T'
638: *
639:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA,
640:      +                           A( 0 ), N, B( 0, 0 ), LDB )
641:                      CALL SGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ),
642:      +                           LDB, A( N1 ), N, ALPHA, B( 0, N1 ),
643:      +                           LDB )
644:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, N2, ONE,
645:      +                           A( N ), N, B( 0, N1 ), LDB )
646: *
647:                   END IF
648: *
649:                ELSE
650: *
651: *                 SIDE  ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
652: *
653:                   IF( NOTRANS ) THEN
654: *
655: *                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
656: *                    TRANS = 'N'
657: *
658:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA,
659:      +                           A( N2 ), N, B( 0, 0 ), LDB )
660:                      CALL SGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ),
661:      +                           LDB, A( 0 ), N, ALPHA, B( 0, N1 ),
662:      +                           LDB )
663:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, N2, ONE,
664:      +                           A( N1 ), N, B( 0, N1 ), LDB )
665: *
666:                   ELSE
667: *
668: *                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
669: *                    TRANS = 'T'
670: *
671:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA,
672:      +                           A( N1 ), N, B( 0, N1 ), LDB )
673:                      CALL SGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ),
674:      +                           LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB )
675:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, N1, ONE,
676:      +                           A( N2 ), N, B( 0, 0 ), LDB )
677: *
678:                   END IF
679: *
680:                END IF
681: *
682:             ELSE
683: *
684: *              SIDE = 'R', N is odd, and TRANSR = 'T'
685: *
686:                IF( LOWER ) THEN
687: *
688: *                 SIDE  ='R', N is odd, TRANSR = 'T', and UPLO = 'L'
689: *
690:                   IF( NOTRANS ) THEN
691: *
692: *                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
693: *                    TRANS = 'N'
694: *
695:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
696:      +                           A( 1 ), N1, B( 0, N1 ), LDB )
697:                      CALL SGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ),
698:      +                           LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ),
699:      +                           LDB )
700:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, N1, ONE,
701:      +                           A( 0 ), N1, B( 0, 0 ), LDB )
702: *
703:                   ELSE
704: *
705: *                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
706: *                    TRANS = 'T'
707: *
708:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
709:      +                           A( 0 ), N1, B( 0, 0 ), LDB )
710:                      CALL SGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ),
711:      +                           LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ),
712:      +                           LDB )
713:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, N2, ONE,
714:      +                           A( 1 ), N1, B( 0, N1 ), LDB )
715: *
716:                   END IF
717: *
718:                ELSE
719: *
720: *                 SIDE  ='R', N is odd, TRANSR = 'T', and UPLO = 'U'
721: *
722:                   IF( NOTRANS ) THEN
723: *
724: *                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
725: *                    TRANS = 'N'
726: *
727:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
728:      +                           A( N2*N2 ), N2, B( 0, 0 ), LDB )
729:                      CALL SGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ),
730:      +                           LDB, A( 0 ), N2, ALPHA, B( 0, N1 ),
731:      +                           LDB )
732:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, N2, ONE,
733:      +                           A( N1*N2 ), N2, B( 0, N1 ), LDB )
734: *
735:                   ELSE
736: *
737: *                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
738: *                    TRANS = 'T'
739: *
740:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
741:      +                           A( N1*N2 ), N2, B( 0, N1 ), LDB )
742:                      CALL SGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ),
743:      +                           LDB, A( 0 ), N2, ALPHA, B( 0, 0 ),
744:      +                           LDB )
745:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, N1, ONE,
746:      +                           A( N2*N2 ), N2, B( 0, 0 ), LDB )
747: *
748:                   END IF
749: *
750:                END IF
751: *
752:             END IF
753: *
754:          ELSE
755: *
756: *           SIDE = 'R' and N is even
757: *
758:             IF( NORMALTRANSR ) THEN
759: *
760: *              SIDE = 'R', N is even, and TRANSR = 'N'
761: *
762:                IF( LOWER ) THEN
763: *
764: *                 SIDE  ='R', N is even, TRANSR = 'N', and UPLO = 'L'
765: *
766:                   IF( NOTRANS ) THEN
767: *
768: *                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'L',
769: *                    and TRANS = 'N'
770: *
771:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA,
772:      +                           A( 0 ), N+1, B( 0, K ), LDB )
773:                      CALL SGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ),
774:      +                           LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ),
775:      +                           LDB )
776:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, K, ONE,
777:      +                           A( 1 ), N+1, B( 0, 0 ), LDB )
778: *
779:                   ELSE
780: *
781: *                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'L',
782: *                    and TRANS = 'T'
783: *
784:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA,
785:      +                           A( 1 ), N+1, B( 0, 0 ), LDB )
786:                      CALL SGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ),
787:      +                           LDB, A( K+1 ), N+1, ALPHA, B( 0, K ),
788:      +                           LDB )
789:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, K, ONE,
790:      +                           A( 0 ), N+1, B( 0, K ), LDB )
791: *
792:                   END IF
793: *
794:                ELSE
795: *
796: *                 SIDE  ='R', N is even, TRANSR = 'N', and UPLO = 'U'
797: *
798:                   IF( NOTRANS ) THEN
799: *
800: *                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'U',
801: *                    and TRANS = 'N'
802: *
803:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA,
804:      +                           A( K+1 ), N+1, B( 0, 0 ), LDB )
805:                      CALL SGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ),
806:      +                           LDB, A( 0 ), N+1, ALPHA, B( 0, K ),
807:      +                           LDB )
808:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, K, ONE,
809:      +                           A( K ), N+1, B( 0, K ), LDB )
810: *
811:                   ELSE
812: *
813: *                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'U',
814: *                    and TRANS = 'T'
815: *
816:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA,
817:      +                           A( K ), N+1, B( 0, K ), LDB )
818:                      CALL SGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ),
819:      +                           LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ),
820:      +                           LDB )
821:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, K, ONE,
822:      +                           A( K+1 ), N+1, B( 0, 0 ), LDB )
823: *
824:                   END IF
825: *
826:                END IF
827: *
828:             ELSE
829: *
830: *              SIDE = 'R', N is even, and TRANSR = 'T'
831: *
832:                IF( LOWER ) THEN
833: *
834: *                 SIDE  ='R', N is even, TRANSR = 'T', and UPLO = 'L'
835: *
836:                   IF( NOTRANS ) THEN
837: *
838: *                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'L',
839: *                    and TRANS = 'N'
840: *
841:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
842:      +                           A( 0 ), K, B( 0, K ), LDB )
843:                      CALL SGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ),
844:      +                           LDB, A( ( K+1 )*K ), K, ALPHA,
845:      +                           B( 0, 0 ), LDB )
846:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, K, ONE,
847:      +                           A( K ), K, B( 0, 0 ), LDB )
848: *
849:                   ELSE
850: *
851: *                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'L',
852: *                    and TRANS = 'T'
853: *
854:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
855:      +                           A( K ), K, B( 0, 0 ), LDB )
856:                      CALL SGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ),
857:      +                           LDB, A( ( K+1 )*K ), K, ALPHA,
858:      +                           B( 0, K ), LDB )
859:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, K, ONE,
860:      +                           A( 0 ), K, B( 0, K ), LDB )
861: *
862:                   END IF
863: *
864:                ELSE
865: *
866: *                 SIDE  ='R', N is even, TRANSR = 'T', and UPLO = 'U'
867: *
868:                   IF( NOTRANS ) THEN
869: *
870: *                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'U',
871: *                    and TRANS = 'N'
872: *
873:                      CALL STRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
874:      +                           A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
875:                      CALL SGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ),
876:      +                           LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB )
877:                      CALL STRSM( 'R', 'L', 'T', DIAG, M, K, ONE,
878:      +                           A( K*K ), K, B( 0, K ), LDB )
879: *
880:                   ELSE
881: *
882: *                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'U',
883: *                    and TRANS = 'T'
884: *
885:                      CALL STRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
886:      +                           A( K*K ), K, B( 0, K ), LDB )
887:                      CALL SGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ),
888:      +                           LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB )
889:                      CALL STRSM( 'R', 'U', 'T', DIAG, M, K, ONE,
890:      +                           A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
891: *
892:                   END IF
893: *
894:                END IF
895: *
896:             END IF
897: *
898:          END IF
899:       END IF
900: *
901:       RETURN
902: *
903: *     End of STFSM
904: *
905:       END
906: