LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine stfsm ( character  TRANSR,
character  SIDE,
character  UPLO,
character  TRANS,
character  DIAG,
integer  M,
integer  N,
real  ALPHA,
real, dimension( 0: * )  A,
real, dimension( 0: ldb-1, 0: * )  B,
integer  LDB 
)

STFSM solves a matrix equation (one operand is a triangular matrix in RFP format).

Download STFSM + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 Level 3 BLAS like routine for A in RFP Format.

 STFSM  solves the matrix equation

    op( A )*X = alpha*B  or  X*op( A ) = alpha*B

 where alpha is a scalar, X and B are m by n matrices, A is a unit, or
 non-unit,  upper or lower triangular matrix  and  op( A )  is one  of

    op( A ) = A   or   op( A ) = A**T.

 A is in Rectangular Full Packed (RFP) Format.

 The matrix X is overwritten on B.
Parameters
[in]TRANSR
          TRANSR is CHARACTER*1
          = 'N':  The Normal Form of RFP A is stored;
          = 'T':  The Transpose Form of RFP A is stored.
[in]SIDE
          SIDE is CHARACTER*1
           On entry, SIDE specifies whether op( A ) appears on the left
           or right of X as follows:

              SIDE = 'L' or 'l'   op( A )*X = alpha*B.

              SIDE = 'R' or 'r'   X*op( A ) = alpha*B.

           Unchanged on exit.
[in]UPLO
          UPLO is CHARACTER*1
           On entry, UPLO specifies whether the RFP matrix A came from
           an upper or lower triangular matrix as follows:
           UPLO = 'U' or 'u' RFP A came from an upper triangular matrix
           UPLO = 'L' or 'l' RFP A came from a  lower triangular matrix

           Unchanged on exit.
[in]TRANS
          TRANS is CHARACTER*1
           On entry, TRANS  specifies the form of op( A ) to be used
           in the matrix multiplication as follows:

              TRANS  = 'N' or 'n'   op( A ) = A.

              TRANS  = 'T' or 't'   op( A ) = A'.

           Unchanged on exit.
[in]DIAG
          DIAG is CHARACTER*1
           On entry, DIAG specifies whether or not RFP A is unit
           triangular as follows:

              DIAG = 'U' or 'u'   A is assumed to be unit triangular.

              DIAG = 'N' or 'n'   A is not assumed to be unit
                                  triangular.

           Unchanged on exit.
[in]M
          M is INTEGER
           On entry, M specifies the number of rows of B. M must be at
           least zero.
           Unchanged on exit.
[in]N
          N is INTEGER
           On entry, N specifies the number of columns of B.  N must be
           at least zero.
           Unchanged on exit.
[in]ALPHA
          ALPHA is REAL
           On entry,  ALPHA specifies the scalar  alpha. When  alpha is
           zero then  A is not referenced and  B need not be set before
           entry.
           Unchanged on exit.
[in]A
          A is REAL array, dimension (NT)
           NT = N*(N+1)/2. On entry, the matrix A in RFP Format.
           RFP Format is described by TRANSR, UPLO and N as follows:
           If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
           K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
           TRANSR = 'T' then RFP is the transpose of RFP A as
           defined when TRANSR = 'N'. The contents of RFP A are defined
           by UPLO as follows: If UPLO = 'U' the RFP A contains the NT
           elements of upper packed A either in normal or
           transpose Format. If UPLO = 'L' the RFP A contains
           the NT elements of lower packed A either in normal or
           transpose Format. The LDA of RFP A is (N+1)/2 when
           TRANSR = 'T'. When TRANSR is 'N' the LDA is N+1 when N is
           even and is N when is odd.
           See the Note below for more details. Unchanged on exit.
[in,out]B
          B is REAL array, DIMENSION (LDB,N)
           Before entry,  the leading  m by n part of the array  B must
           contain  the  right-hand  side  matrix  B,  and  on exit  is
           overwritten by the solution matrix  X.
[in]LDB
          LDB is INTEGER
           On entry, LDB specifies the first dimension of B as declared
           in  the  calling  (sub)  program.   LDB  must  be  at  least
           max( 1, m ).
           Unchanged on exit.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012
Further Details:
  We first consider Rectangular Full Packed (RFP) Format when N is
  even. We give an example where N = 6.

      AP is Upper             AP is Lower

   00 01 02 03 04 05       00
      11 12 13 14 15       10 11
         22 23 24 25       20 21 22
            33 34 35       30 31 32 33
               44 45       40 41 42 43 44
                  55       50 51 52 53 54 55


  Let TRANSR = 'N'. RFP holds AP as follows:
  For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
  three columns of AP upper. The lower triangle A(4:6,0:2) consists of
  the transpose of the first three columns of AP upper.
  For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
  three columns of AP lower. The upper triangle A(0:2,0:2) consists of
  the transpose of the last three columns of AP lower.
  This covers the case N even and TRANSR = 'N'.

         RFP A                   RFP A

        03 04 05                33 43 53
        13 14 15                00 44 54
        23 24 25                10 11 55
        33 34 35                20 21 22
        00 44 45                30 31 32
        01 11 55                40 41 42
        02 12 22                50 51 52

  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
  transpose of RFP A above. One therefore gets:


           RFP A                   RFP A

     03 13 23 33 00 01 02    33 00 10 20 30 40 50
     04 14 24 34 44 11 12    43 44 11 21 31 41 51
     05 15 25 35 45 55 22    53 54 55 22 32 42 52


  We then consider Rectangular Full Packed (RFP) Format when N is
  odd. We give an example where N = 5.

     AP is Upper                 AP is Lower

   00 01 02 03 04              00
      11 12 13 14              10 11
         22 23 24              20 21 22
            33 34              30 31 32 33
               44              40 41 42 43 44


  Let TRANSR = 'N'. RFP holds AP as follows:
  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
  three columns of AP upper. The lower triangle A(3:4,0:1) consists of
  the transpose of the first two columns of AP upper.
  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
  three columns of AP lower. The upper triangle A(0:1,1:2) consists of
  the transpose of the last two columns of AP lower.
  This covers the case N odd and TRANSR = 'N'.

         RFP A                   RFP A

        02 03 04                00 33 43
        12 13 14                10 11 44
        22 23 24                20 21 22
        00 33 34                30 31 32
        01 11 44                40 41 42

  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
  transpose of RFP A above. One therefore gets:

           RFP A                   RFP A

     02 12 22 00 01             00 10 20 30 40 50
     03 13 23 33 11             33 11 21 31 41 51
     04 14 24 34 44             43 44 22 32 42 52

Definition at line 279 of file stfsm.f.

279 *
280 * -- LAPACK computational routine (version 3.4.2) --
281 * -- LAPACK is a software package provided by Univ. of Tennessee, --
282 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
283 * September 2012
284 *
285 * .. Scalar Arguments ..
286  CHARACTER transr, diag, side, trans, uplo
287  INTEGER ldb, m, n
288  REAL alpha
289 * ..
290 * .. Array Arguments ..
291  REAL a( 0: * ), b( 0: ldb-1, 0: * )
292 * ..
293 *
294 * =====================================================================
295 *
296 * ..
297 * .. Parameters ..
298  REAL one, zero
299  parameter ( one = 1.0e+0, zero = 0.0e+0 )
300 * ..
301 * .. Local Scalars ..
302  LOGICAL lower, lside, misodd, nisodd, normaltransr,
303  $ notrans
304  INTEGER m1, m2, n1, n2, k, info, i, j
305 * ..
306 * .. External Functions ..
307  LOGICAL lsame
308  EXTERNAL lsame
309 * ..
310 * .. External Subroutines ..
311  EXTERNAL sgemm, strsm, xerbla
312 * ..
313 * .. Intrinsic Functions ..
314  INTRINSIC max, mod
315 * ..
316 * .. Executable Statements ..
317 *
318 * Test the input parameters.
319 *
320  info = 0
321  normaltransr = lsame( transr, 'N' )
322  lside = lsame( side, 'L' )
323  lower = lsame( uplo, 'L' )
324  notrans = lsame( trans, 'N' )
325  IF( .NOT.normaltransr .AND. .NOT.lsame( transr, 'T' ) ) THEN
326  info = -1
327  ELSE IF( .NOT.lside .AND. .NOT.lsame( side, 'R' ) ) THEN
328  info = -2
329  ELSE IF( .NOT.lower .AND. .NOT.lsame( uplo, 'U' ) ) THEN
330  info = -3
331  ELSE IF( .NOT.notrans .AND. .NOT.lsame( trans, 'T' ) ) THEN
332  info = -4
333  ELSE IF( .NOT.lsame( diag, 'N' ) .AND. .NOT.lsame( diag, 'U' ) )
334  $ THEN
335  info = -5
336  ELSE IF( m.LT.0 ) THEN
337  info = -6
338  ELSE IF( n.LT.0 ) THEN
339  info = -7
340  ELSE IF( ldb.LT.max( 1, m ) ) THEN
341  info = -11
342  END IF
343  IF( info.NE.0 ) THEN
344  CALL xerbla( 'STFSM ', -info )
345  RETURN
346  END IF
347 *
348 * Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
349 *
350  IF( ( m.EQ.0 ) .OR. ( n.EQ.0 ) )
351  $ RETURN
352 *
353 * Quick return when ALPHA.EQ.(0D+0)
354 *
355  IF( alpha.EQ.zero ) THEN
356  DO 20 j = 0, n - 1
357  DO 10 i = 0, m - 1
358  b( i, j ) = zero
359  10 CONTINUE
360  20 CONTINUE
361  RETURN
362  END IF
363 *
364  IF( lside ) THEN
365 *
366 * SIDE = 'L'
367 *
368 * A is M-by-M.
369 * If M is odd, set NISODD = .TRUE., and M1 and M2.
370 * If M is even, NISODD = .FALSE., and M.
371 *
372  IF( mod( m, 2 ).EQ.0 ) THEN
373  misodd = .false.
374  k = m / 2
375  ELSE
376  misodd = .true.
377  IF( lower ) THEN
378  m2 = m / 2
379  m1 = m - m2
380  ELSE
381  m1 = m / 2
382  m2 = m - m1
383  END IF
384  END IF
385 *
386  IF( misodd ) THEN
387 *
388 * SIDE = 'L' and N is odd
389 *
390  IF( normaltransr ) THEN
391 *
392 * SIDE = 'L', N is odd, and TRANSR = 'N'
393 *
394  IF( lower ) THEN
395 *
396 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
397 *
398  IF( notrans ) THEN
399 *
400 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
401 * TRANS = 'N'
402 *
403  IF( m.EQ.1 ) THEN
404  CALL strsm( 'L', 'L', 'N', diag, m1, n, alpha,
405  $ a, m, b, ldb )
406  ELSE
407  CALL strsm( 'L', 'L', 'N', diag, m1, n, alpha,
408  $ a( 0 ), m, b, ldb )
409  CALL sgemm( 'N', 'N', m2, n, m1, -one, a( m1 ),
410  $ m, b, ldb, alpha, b( m1, 0 ), ldb )
411  CALL strsm( 'L', 'U', 'T', diag, m2, n, one,
412  $ a( m ), m, b( m1, 0 ), ldb )
413  END IF
414 *
415  ELSE
416 *
417 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
418 * TRANS = 'T'
419 *
420  IF( m.EQ.1 ) THEN
421  CALL strsm( 'L', 'L', 'T', diag, m1, n, alpha,
422  $ a( 0 ), m, b, ldb )
423  ELSE
424  CALL strsm( 'L', 'U', 'N', diag, m2, n, alpha,
425  $ a( m ), m, b( m1, 0 ), ldb )
426  CALL sgemm( 'T', 'N', m1, n, m2, -one, a( m1 ),
427  $ m, b( m1, 0 ), ldb, alpha, b, ldb )
428  CALL strsm( 'L', 'L', 'T', diag, m1, n, one,
429  $ a( 0 ), m, b, ldb )
430  END IF
431 *
432  END IF
433 *
434  ELSE
435 *
436 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
437 *
438  IF( .NOT.notrans ) THEN
439 *
440 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
441 * TRANS = 'N'
442 *
443  CALL strsm( 'L', 'L', 'N', diag, m1, n, alpha,
444  $ a( m2 ), m, b, ldb )
445  CALL sgemm( 'T', 'N', m2, n, m1, -one, a( 0 ), m,
446  $ b, ldb, alpha, b( m1, 0 ), ldb )
447  CALL strsm( 'L', 'U', 'T', diag, m2, n, one,
448  $ a( m1 ), m, b( m1, 0 ), ldb )
449 *
450  ELSE
451 *
452 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
453 * TRANS = 'T'
454 *
455  CALL strsm( 'L', 'U', 'N', diag, m2, n, alpha,
456  $ a( m1 ), m, b( m1, 0 ), ldb )
457  CALL sgemm( 'N', 'N', m1, n, m2, -one, a( 0 ), m,
458  $ b( m1, 0 ), ldb, alpha, b, ldb )
459  CALL strsm( 'L', 'L', 'T', diag, m1, n, one,
460  $ a( m2 ), m, b, ldb )
461 *
462  END IF
463 *
464  END IF
465 *
466  ELSE
467 *
468 * SIDE = 'L', N is odd, and TRANSR = 'T'
469 *
470  IF( lower ) THEN
471 *
472 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'L'
473 *
474  IF( notrans ) THEN
475 *
476 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
477 * TRANS = 'N'
478 *
479  IF( m.EQ.1 ) THEN
480  CALL strsm( 'L', 'U', 'T', diag, m1, n, alpha,
481  $ a( 0 ), m1, b, ldb )
482  ELSE
483  CALL strsm( 'L', 'U', 'T', diag, m1, n, alpha,
484  $ a( 0 ), m1, b, ldb )
485  CALL sgemm( 'T', 'N', m2, n, m1, -one,
486  $ a( m1*m1 ), m1, b, ldb, alpha,
487  $ b( m1, 0 ), ldb )
488  CALL strsm( 'L', 'L', 'N', diag, m2, n, one,
489  $ a( 1 ), m1, b( m1, 0 ), ldb )
490  END IF
491 *
492  ELSE
493 *
494 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
495 * TRANS = 'T'
496 *
497  IF( m.EQ.1 ) THEN
498  CALL strsm( 'L', 'U', 'N', diag, m1, n, alpha,
499  $ a( 0 ), m1, b, ldb )
500  ELSE
501  CALL strsm( 'L', 'L', 'T', diag, m2, n, alpha,
502  $ a( 1 ), m1, b( m1, 0 ), ldb )
503  CALL sgemm( 'N', 'N', m1, n, m2, -one,
504  $ a( m1*m1 ), m1, b( m1, 0 ), ldb,
505  $ alpha, b, ldb )
506  CALL strsm( 'L', 'U', 'N', diag, m1, n, one,
507  $ a( 0 ), m1, b, ldb )
508  END IF
509 *
510  END IF
511 *
512  ELSE
513 *
514 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'U'
515 *
516  IF( .NOT.notrans ) THEN
517 *
518 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
519 * TRANS = 'N'
520 *
521  CALL strsm( 'L', 'U', 'T', diag, m1, n, alpha,
522  $ a( m2*m2 ), m2, b, ldb )
523  CALL sgemm( 'N', 'N', m2, n, m1, -one, a( 0 ), m2,
524  $ b, ldb, alpha, b( m1, 0 ), ldb )
525  CALL strsm( 'L', 'L', 'N', diag, m2, n, one,
526  $ a( m1*m2 ), m2, b( m1, 0 ), ldb )
527 *
528  ELSE
529 *
530 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
531 * TRANS = 'T'
532 *
533  CALL strsm( 'L', 'L', 'T', diag, m2, n, alpha,
534  $ a( m1*m2 ), m2, b( m1, 0 ), ldb )
535  CALL sgemm( 'T', 'N', m1, n, m2, -one, a( 0 ), m2,
536  $ b( m1, 0 ), ldb, alpha, b, ldb )
537  CALL strsm( 'L', 'U', 'N', diag, m1, n, one,
538  $ a( m2*m2 ), m2, b, ldb )
539 *
540  END IF
541 *
542  END IF
543 *
544  END IF
545 *
546  ELSE
547 *
548 * SIDE = 'L' and N is even
549 *
550  IF( normaltransr ) THEN
551 *
552 * SIDE = 'L', N is even, and TRANSR = 'N'
553 *
554  IF( lower ) THEN
555 *
556 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L'
557 *
558  IF( notrans ) THEN
559 *
560 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
561 * and TRANS = 'N'
562 *
563  CALL strsm( 'L', 'L', 'N', diag, k, n, alpha,
564  $ a( 1 ), m+1, b, ldb )
565  CALL sgemm( 'N', 'N', k, n, k, -one, a( k+1 ),
566  $ m+1, b, ldb, alpha, b( k, 0 ), ldb )
567  CALL strsm( 'L', 'U', 'T', diag, k, n, one,
568  $ a( 0 ), m+1, b( k, 0 ), ldb )
569 *
570  ELSE
571 *
572 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
573 * and TRANS = 'T'
574 *
575  CALL strsm( 'L', 'U', 'N', diag, k, n, alpha,
576  $ a( 0 ), m+1, b( k, 0 ), ldb )
577  CALL sgemm( 'T', 'N', k, n, k, -one, a( k+1 ),
578  $ m+1, b( k, 0 ), ldb, alpha, b, ldb )
579  CALL strsm( 'L', 'L', 'T', diag, k, n, one,
580  $ a( 1 ), m+1, b, ldb )
581 *
582  END IF
583 *
584  ELSE
585 *
586 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U'
587 *
588  IF( .NOT.notrans ) THEN
589 *
590 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
591 * and TRANS = 'N'
592 *
593  CALL strsm( 'L', 'L', 'N', diag, k, n, alpha,
594  $ a( k+1 ), m+1, b, ldb )
595  CALL sgemm( 'T', 'N', k, n, k, -one, a( 0 ), m+1,
596  $ b, ldb, alpha, b( k, 0 ), ldb )
597  CALL strsm( 'L', 'U', 'T', diag, k, n, one,
598  $ a( k ), m+1, b( k, 0 ), ldb )
599 *
600  ELSE
601 *
602 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
603 * and TRANS = 'T'
604  CALL strsm( 'L', 'U', 'N', diag, k, n, alpha,
605  $ a( k ), m+1, b( k, 0 ), ldb )
606  CALL sgemm( 'N', 'N', k, n, k, -one, a( 0 ), m+1,
607  $ b( k, 0 ), ldb, alpha, b, ldb )
608  CALL strsm( 'L', 'L', 'T', diag, k, n, one,
609  $ a( k+1 ), m+1, b, ldb )
610 *
611  END IF
612 *
613  END IF
614 *
615  ELSE
616 *
617 * SIDE = 'L', N is even, and TRANSR = 'T'
618 *
619  IF( lower ) THEN
620 *
621 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'L'
622 *
623  IF( notrans ) THEN
624 *
625 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
626 * and TRANS = 'N'
627 *
628  CALL strsm( 'L', 'U', 'T', diag, k, n, alpha,
629  $ a( k ), k, b, ldb )
630  CALL sgemm( 'T', 'N', k, n, k, -one,
631  $ a( k*( k+1 ) ), k, b, ldb, alpha,
632  $ b( k, 0 ), ldb )
633  CALL strsm( 'L', 'L', 'N', diag, k, n, one,
634  $ a( 0 ), k, b( k, 0 ), ldb )
635 *
636  ELSE
637 *
638 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
639 * and TRANS = 'T'
640 *
641  CALL strsm( 'L', 'L', 'T', diag, k, n, alpha,
642  $ a( 0 ), k, b( k, 0 ), ldb )
643  CALL sgemm( 'N', 'N', k, n, k, -one,
644  $ a( k*( k+1 ) ), k, b( k, 0 ), ldb,
645  $ alpha, b, ldb )
646  CALL strsm( 'L', 'U', 'N', diag, k, n, one,
647  $ a( k ), k, b, ldb )
648 *
649  END IF
650 *
651  ELSE
652 *
653 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'U'
654 *
655  IF( .NOT.notrans ) THEN
656 *
657 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
658 * and TRANS = 'N'
659 *
660  CALL strsm( 'L', 'U', 'T', diag, k, n, alpha,
661  $ a( k*( k+1 ) ), k, b, ldb )
662  CALL sgemm( 'N', 'N', k, n, k, -one, a( 0 ), k, b,
663  $ ldb, alpha, b( k, 0 ), ldb )
664  CALL strsm( 'L', 'L', 'N', diag, k, n, one,
665  $ a( k*k ), k, b( k, 0 ), ldb )
666 *
667  ELSE
668 *
669 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
670 * and TRANS = 'T'
671 *
672  CALL strsm( 'L', 'L', 'T', diag, k, n, alpha,
673  $ a( k*k ), k, b( k, 0 ), ldb )
674  CALL sgemm( 'T', 'N', k, n, k, -one, a( 0 ), k,
675  $ b( k, 0 ), ldb, alpha, b, ldb )
676  CALL strsm( 'L', 'U', 'N', diag, k, n, one,
677  $ a( k*( k+1 ) ), k, b, ldb )
678 *
679  END IF
680 *
681  END IF
682 *
683  END IF
684 *
685  END IF
686 *
687  ELSE
688 *
689 * SIDE = 'R'
690 *
691 * A is N-by-N.
692 * If N is odd, set NISODD = .TRUE., and N1 and N2.
693 * If N is even, NISODD = .FALSE., and K.
694 *
695  IF( mod( n, 2 ).EQ.0 ) THEN
696  nisodd = .false.
697  k = n / 2
698  ELSE
699  nisodd = .true.
700  IF( lower ) THEN
701  n2 = n / 2
702  n1 = n - n2
703  ELSE
704  n1 = n / 2
705  n2 = n - n1
706  END IF
707  END IF
708 *
709  IF( nisodd ) THEN
710 *
711 * SIDE = 'R' and N is odd
712 *
713  IF( normaltransr ) THEN
714 *
715 * SIDE = 'R', N is odd, and TRANSR = 'N'
716 *
717  IF( lower ) THEN
718 *
719 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
720 *
721  IF( notrans ) THEN
722 *
723 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
724 * TRANS = 'N'
725 *
726  CALL strsm( 'R', 'U', 'T', diag, m, n2, alpha,
727  $ a( n ), n, b( 0, n1 ), ldb )
728  CALL sgemm( 'N', 'N', m, n1, n2, -one, b( 0, n1 ),
729  $ ldb, a( n1 ), n, alpha, b( 0, 0 ),
730  $ ldb )
731  CALL strsm( 'R', 'L', 'N', diag, m, n1, one,
732  $ a( 0 ), n, b( 0, 0 ), ldb )
733 *
734  ELSE
735 *
736 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
737 * TRANS = 'T'
738 *
739  CALL strsm( 'R', 'L', 'T', diag, m, n1, alpha,
740  $ a( 0 ), n, b( 0, 0 ), ldb )
741  CALL sgemm( 'N', 'T', m, n2, n1, -one, b( 0, 0 ),
742  $ ldb, a( n1 ), n, alpha, b( 0, n1 ),
743  $ ldb )
744  CALL strsm( 'R', 'U', 'N', diag, m, n2, one,
745  $ a( n ), n, b( 0, n1 ), ldb )
746 *
747  END IF
748 *
749  ELSE
750 *
751 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
752 *
753  IF( notrans ) THEN
754 *
755 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
756 * TRANS = 'N'
757 *
758  CALL strsm( 'R', 'L', 'T', diag, m, n1, alpha,
759  $ a( n2 ), n, b( 0, 0 ), ldb )
760  CALL sgemm( 'N', 'N', m, n2, n1, -one, b( 0, 0 ),
761  $ ldb, a( 0 ), n, alpha, b( 0, n1 ),
762  $ ldb )
763  CALL strsm( 'R', 'U', 'N', diag, m, n2, one,
764  $ a( n1 ), n, b( 0, n1 ), ldb )
765 *
766  ELSE
767 *
768 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
769 * TRANS = 'T'
770 *
771  CALL strsm( 'R', 'U', 'T', diag, m, n2, alpha,
772  $ a( n1 ), n, b( 0, n1 ), ldb )
773  CALL sgemm( 'N', 'T', m, n1, n2, -one, b( 0, n1 ),
774  $ ldb, a( 0 ), n, alpha, b( 0, 0 ), ldb )
775  CALL strsm( 'R', 'L', 'N', diag, m, n1, one,
776  $ a( n2 ), n, b( 0, 0 ), ldb )
777 *
778  END IF
779 *
780  END IF
781 *
782  ELSE
783 *
784 * SIDE = 'R', N is odd, and TRANSR = 'T'
785 *
786  IF( lower ) THEN
787 *
788 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'L'
789 *
790  IF( notrans ) THEN
791 *
792 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
793 * TRANS = 'N'
794 *
795  CALL strsm( 'R', 'L', 'N', diag, m, n2, alpha,
796  $ a( 1 ), n1, b( 0, n1 ), ldb )
797  CALL sgemm( 'N', 'T', m, n1, n2, -one, b( 0, n1 ),
798  $ ldb, a( n1*n1 ), n1, alpha, b( 0, 0 ),
799  $ ldb )
800  CALL strsm( 'R', 'U', 'T', diag, m, n1, one,
801  $ a( 0 ), n1, b( 0, 0 ), ldb )
802 *
803  ELSE
804 *
805 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
806 * TRANS = 'T'
807 *
808  CALL strsm( 'R', 'U', 'N', diag, m, n1, alpha,
809  $ a( 0 ), n1, b( 0, 0 ), ldb )
810  CALL sgemm( 'N', 'N', m, n2, n1, -one, b( 0, 0 ),
811  $ ldb, a( n1*n1 ), n1, alpha, b( 0, n1 ),
812  $ ldb )
813  CALL strsm( 'R', 'L', 'T', diag, m, n2, one,
814  $ a( 1 ), n1, b( 0, n1 ), ldb )
815 *
816  END IF
817 *
818  ELSE
819 *
820 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'U'
821 *
822  IF( notrans ) THEN
823 *
824 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
825 * TRANS = 'N'
826 *
827  CALL strsm( 'R', 'U', 'N', diag, m, n1, alpha,
828  $ a( n2*n2 ), n2, b( 0, 0 ), ldb )
829  CALL sgemm( 'N', 'T', m, n2, n1, -one, b( 0, 0 ),
830  $ ldb, a( 0 ), n2, alpha, b( 0, n1 ),
831  $ ldb )
832  CALL strsm( 'R', 'L', 'T', diag, m, n2, one,
833  $ a( n1*n2 ), n2, b( 0, n1 ), ldb )
834 *
835  ELSE
836 *
837 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
838 * TRANS = 'T'
839 *
840  CALL strsm( 'R', 'L', 'N', diag, m, n2, alpha,
841  $ a( n1*n2 ), n2, b( 0, n1 ), ldb )
842  CALL sgemm( 'N', 'N', m, n1, n2, -one, b( 0, n1 ),
843  $ ldb, a( 0 ), n2, alpha, b( 0, 0 ),
844  $ ldb )
845  CALL strsm( 'R', 'U', 'T', diag, m, n1, one,
846  $ a( n2*n2 ), n2, b( 0, 0 ), ldb )
847 *
848  END IF
849 *
850  END IF
851 *
852  END IF
853 *
854  ELSE
855 *
856 * SIDE = 'R' and N is even
857 *
858  IF( normaltransr ) THEN
859 *
860 * SIDE = 'R', N is even, and TRANSR = 'N'
861 *
862  IF( lower ) THEN
863 *
864 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L'
865 *
866  IF( notrans ) THEN
867 *
868 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
869 * and TRANS = 'N'
870 *
871  CALL strsm( 'R', 'U', 'T', diag, m, k, alpha,
872  $ a( 0 ), n+1, b( 0, k ), ldb )
873  CALL sgemm( 'N', 'N', m, k, k, -one, b( 0, k ),
874  $ ldb, a( k+1 ), n+1, alpha, b( 0, 0 ),
875  $ ldb )
876  CALL strsm( 'R', 'L', 'N', diag, m, k, one,
877  $ a( 1 ), n+1, b( 0, 0 ), ldb )
878 *
879  ELSE
880 *
881 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
882 * and TRANS = 'T'
883 *
884  CALL strsm( 'R', 'L', 'T', diag, m, k, alpha,
885  $ a( 1 ), n+1, b( 0, 0 ), ldb )
886  CALL sgemm( 'N', 'T', m, k, k, -one, b( 0, 0 ),
887  $ ldb, a( k+1 ), n+1, alpha, b( 0, k ),
888  $ ldb )
889  CALL strsm( 'R', 'U', 'N', diag, m, k, one,
890  $ a( 0 ), n+1, b( 0, k ), ldb )
891 *
892  END IF
893 *
894  ELSE
895 *
896 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U'
897 *
898  IF( notrans ) THEN
899 *
900 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
901 * and TRANS = 'N'
902 *
903  CALL strsm( 'R', 'L', 'T', diag, m, k, alpha,
904  $ a( k+1 ), n+1, b( 0, 0 ), ldb )
905  CALL sgemm( 'N', 'N', m, k, k, -one, b( 0, 0 ),
906  $ ldb, a( 0 ), n+1, alpha, b( 0, k ),
907  $ ldb )
908  CALL strsm( 'R', 'U', 'N', diag, m, k, one,
909  $ a( k ), n+1, b( 0, k ), ldb )
910 *
911  ELSE
912 *
913 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
914 * and TRANS = 'T'
915 *
916  CALL strsm( 'R', 'U', 'T', diag, m, k, alpha,
917  $ a( k ), n+1, b( 0, k ), ldb )
918  CALL sgemm( 'N', 'T', m, k, k, -one, b( 0, k ),
919  $ ldb, a( 0 ), n+1, alpha, b( 0, 0 ),
920  $ ldb )
921  CALL strsm( 'R', 'L', 'N', diag, m, k, one,
922  $ a( k+1 ), n+1, b( 0, 0 ), ldb )
923 *
924  END IF
925 *
926  END IF
927 *
928  ELSE
929 *
930 * SIDE = 'R', N is even, and TRANSR = 'T'
931 *
932  IF( lower ) THEN
933 *
934 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'L'
935 *
936  IF( notrans ) THEN
937 *
938 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
939 * and TRANS = 'N'
940 *
941  CALL strsm( 'R', 'L', 'N', diag, m, k, alpha,
942  $ a( 0 ), k, b( 0, k ), ldb )
943  CALL sgemm( 'N', 'T', m, k, k, -one, b( 0, k ),
944  $ ldb, a( ( k+1 )*k ), k, alpha,
945  $ b( 0, 0 ), ldb )
946  CALL strsm( 'R', 'U', 'T', diag, m, k, one,
947  $ a( k ), k, b( 0, 0 ), ldb )
948 *
949  ELSE
950 *
951 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
952 * and TRANS = 'T'
953 *
954  CALL strsm( 'R', 'U', 'N', diag, m, k, alpha,
955  $ a( k ), k, b( 0, 0 ), ldb )
956  CALL sgemm( 'N', 'N', m, k, k, -one, b( 0, 0 ),
957  $ ldb, a( ( k+1 )*k ), k, alpha,
958  $ b( 0, k ), ldb )
959  CALL strsm( 'R', 'L', 'T', diag, m, k, one,
960  $ a( 0 ), k, b( 0, k ), ldb )
961 *
962  END IF
963 *
964  ELSE
965 *
966 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'U'
967 *
968  IF( notrans ) THEN
969 *
970 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
971 * and TRANS = 'N'
972 *
973  CALL strsm( 'R', 'U', 'N', diag, m, k, alpha,
974  $ a( ( k+1 )*k ), k, b( 0, 0 ), ldb )
975  CALL sgemm( 'N', 'T', m, k, k, -one, b( 0, 0 ),
976  $ ldb, a( 0 ), k, alpha, b( 0, k ), ldb )
977  CALL strsm( 'R', 'L', 'T', diag, m, k, one,
978  $ a( k*k ), k, b( 0, k ), ldb )
979 *
980  ELSE
981 *
982 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
983 * and TRANS = 'T'
984 *
985  CALL strsm( 'R', 'L', 'N', diag, m, k, alpha,
986  $ a( k*k ), k, b( 0, k ), ldb )
987  CALL sgemm( 'N', 'N', m, k, k, -one, b( 0, k ),
988  $ ldb, a( 0 ), k, alpha, b( 0, 0 ), ldb )
989  CALL strsm( 'R', 'U', 'T', diag, m, k, one,
990  $ a( ( k+1 )*k ), k, b( 0, 0 ), ldb )
991 *
992  END IF
993 *
994  END IF
995 *
996  END IF
997 *
998  END IF
999  END IF
1000 *
1001  RETURN
1002 *
1003 * End of STFSM
1004 *
subroutine strsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRSM
Definition: strsm.f:183
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
Definition: sgemm.f:189
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function: