LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ stfsm()

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.
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 275 of file stfsm.f.

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