277 SUBROUTINE dtfsm( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
286 CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO
288 DOUBLE PRECISION ALPHA
291 DOUBLE PRECISION A( 0: * ), B( 0: ldb-1, 0: * )
298 DOUBLE PRECISION ONE, ZERO
299 parameter( one = 1.0d+0, zero = 0.0d+0 )
302 LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
304 INTEGER M1, M2, N1, N2, K, INFO, I, J
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
327 ELSE IF( .NOT.lside .AND. .NOT.lsame( side,
'R' ) )
THEN
329 ELSE IF( .NOT.lower .AND. .NOT.lsame( uplo,
'U' ) )
THEN
331 ELSE IF( .NOT.notrans .AND. .NOT.lsame( trans,
'T' ) )
THEN
333 ELSE IF( .NOT.lsame( diag,
'N' ) .AND. .NOT.lsame( diag,
'U' ) )
336 ELSE IF( m.LT.0 )
THEN
338 ELSE IF( n.LT.0 )
THEN
340 ELSE IF( ldb.LT.max( 1, m ) )
THEN
344 CALL xerbla(
'DTFSM ', -info )
350 IF( ( m.EQ.0 ) .OR. ( n.EQ.0 ) )
355 IF( alpha.EQ.zero )
THEN
372 IF( mod( m, 2 ).EQ.0 )
THEN
391 IF( normaltransr )
THEN
405 CALL dtrsm(
'L',
'L',
'N', diag, m1, n, alpha,
408 CALL dtrsm(
'L',
'L',
'N', diag, m1, n, alpha,
409 $ a( 0 ), m, b, ldb )
410 CALL dgemm(
'N',
'N', m2, n, m1, -one, a( m1 ),
411 $ m, b, ldb, alpha, b( m1, 0 ), ldb )
412 CALL dtrsm(
'L',
'U',
'T', diag, m2, n, one,
413 $ a( m ), m, b( m1, 0 ), ldb )
422 CALL dtrsm(
'L',
'L',
'T', diag, m1, n, alpha,
423 $ a( 0 ), m, b, ldb )
425 CALL dtrsm(
'L',
'U',
'N', diag, m2, n, alpha,
426 $ a( m ), m, b( m1, 0 ), ldb )
427 CALL dgemm(
'T',
'N', m1, n, m2, -one, a( m1 ),
428 $ m, b( m1, 0 ), ldb, alpha, b, ldb )
429 CALL dtrsm(
'L',
'L',
'T', diag, m1, n, one,
430 $ a( 0 ), m, b, ldb )
439 IF( .NOT.notrans )
THEN
444 CALL dtrsm(
'L',
'L',
'N', diag, m1, n, alpha,
445 $ a( m2 ), m, b, ldb )
446 CALL dgemm(
'T',
'N', m2, n, m1, -one, a( 0 ), m,
447 $ b, ldb, alpha, b( m1, 0 ), ldb )
448 CALL dtrsm(
'L',
'U',
'T', diag, m2, n, one,
449 $ a( m1 ), m, b( m1, 0 ), ldb )
456 CALL dtrsm(
'L',
'U',
'N', diag, m2, n, alpha,
457 $ a( m1 ), m, b( m1, 0 ), ldb )
458 CALL dgemm(
'N',
'N', m1, n, m2, -one, a( 0 ), m,
459 $ b( m1, 0 ), ldb, alpha, b, ldb )
460 CALL dtrsm(
'L',
'L',
'T', diag, m1, n, one,
461 $ a( m2 ), m, b, ldb )
481 CALL dtrsm(
'L',
'U',
'T', diag, m1, n, alpha,
482 $ a( 0 ), m1, b, ldb )
484 CALL dtrsm(
'L',
'U',
'T', diag, m1, n, alpha,
485 $ a( 0 ), m1, b, ldb )
486 CALL dgemm(
'T',
'N', m2, n, m1, -one,
487 $ a( m1*m1 ), m1, b, ldb, alpha,
489 CALL dtrsm(
'L',
'L',
'N', diag, m2, n, one,
490 $ a( 1 ), m1, b( m1, 0 ), ldb )
499 CALL dtrsm(
'L',
'U',
'N', diag, m1, n, alpha,
500 $ a( 0 ), m1, b, ldb )
502 CALL dtrsm(
'L',
'L',
'T', diag, m2, n, alpha,
503 $ a( 1 ), m1, b( m1, 0 ), ldb )
504 CALL dgemm(
'N',
'N', m1, n, m2, -one,
505 $ a( m1*m1 ), m1, b( m1, 0 ), ldb,
507 CALL dtrsm(
'L',
'U',
'N', diag, m1, n, one,
508 $ a( 0 ), m1, b, ldb )
517 IF( .NOT.notrans )
THEN
522 CALL dtrsm(
'L',
'U',
'T', diag, m1, n, alpha,
523 $ a( m2*m2 ), m2, b, ldb )
524 CALL dgemm(
'N',
'N', m2, n, m1, -one, a( 0 ), m2,
525 $ b, ldb, alpha, b( m1, 0 ), ldb )
526 CALL dtrsm(
'L',
'L',
'N', diag, m2, n, one,
527 $ a( m1*m2 ), m2, b( m1, 0 ), ldb )
534 CALL dtrsm(
'L',
'L',
'T', diag, m2, n, alpha,
535 $ a( m1*m2 ), m2, b( m1, 0 ), ldb )
536 CALL dgemm(
'T',
'N', m1, n, m2, -one, a( 0 ), m2,
537 $ b( m1, 0 ), ldb, alpha, b, ldb )
538 CALL dtrsm(
'L',
'U',
'N', diag, m1, n, one,
539 $ a( m2*m2 ), m2, b, ldb )
551 IF( normaltransr )
THEN
564 CALL dtrsm(
'L',
'L',
'N', diag, k, n, alpha,
565 $ a( 1 ), m+1, b, ldb )
566 CALL dgemm(
'N',
'N', k, n, k, -one, a( k+1 ),
567 $ m+1, b, ldb, alpha, b( k, 0 ), ldb )
568 CALL dtrsm(
'L',
'U',
'T', diag, k, n, one,
569 $ a( 0 ), m+1, b( k, 0 ), ldb )
576 CALL dtrsm(
'L',
'U',
'N', diag, k, n, alpha,
577 $ a( 0 ), m+1, b( k, 0 ), ldb )
578 CALL dgemm(
'T',
'N', k, n, k, -one, a( k+1 ),
579 $ m+1, b( k, 0 ), ldb, alpha, b, ldb )
580 CALL dtrsm(
'L',
'L',
'T', diag, k, n, one,
581 $ a( 1 ), m+1, b, ldb )
589 IF( .NOT.notrans )
THEN
594 CALL dtrsm(
'L',
'L',
'N', diag, k, n, alpha,
595 $ a( k+1 ), m+1, b, ldb )
596 CALL dgemm(
'T',
'N', k, n, k, -one, a( 0 ), m+1,
597 $ b, ldb, alpha, b( k, 0 ), ldb )
598 CALL dtrsm(
'L',
'U',
'T', diag, k, n, one,
599 $ a( k ), m+1, b( k, 0 ), ldb )
605 CALL dtrsm(
'L',
'U',
'N', diag, k, n, alpha,
606 $ a( k ), m+1, b( k, 0 ), ldb )
607 CALL dgemm(
'N',
'N', k, n, k, -one, a( 0 ), m+1,
608 $ b( k, 0 ), ldb, alpha, b, ldb )
609 CALL dtrsm(
'L',
'L',
'T', diag, k, n, one,
610 $ a( k+1 ), m+1, b, ldb )
629 CALL dtrsm(
'L',
'U',
'T', diag, k, n, alpha,
630 $ a( k ), k, b, ldb )
631 CALL dgemm(
'T',
'N', k, n, k, -one,
632 $ a( k*( k+1 ) ), k, b, ldb, alpha,
634 CALL dtrsm(
'L',
'L',
'N', diag, k, n, one,
635 $ a( 0 ), k, b( k, 0 ), ldb )
642 CALL dtrsm(
'L',
'L',
'T', diag, k, n, alpha,
643 $ a( 0 ), k, b( k, 0 ), ldb )
644 CALL dgemm(
'N',
'N', k, n, k, -one,
645 $ a( k*( k+1 ) ), k, b( k, 0 ), ldb,
647 CALL dtrsm(
'L',
'U',
'N', diag, k, n, one,
648 $ a( k ), k, b, ldb )
656 IF( .NOT.notrans )
THEN
661 CALL dtrsm(
'L',
'U',
'T', diag, k, n, alpha,
662 $ a( k*( k+1 ) ), k, b, ldb )
663 CALL dgemm(
'N',
'N', k, n, k, -one, a( 0 ), k, b,
664 $ ldb, alpha, b( k, 0 ), ldb )
665 CALL dtrsm(
'L',
'L',
'N', diag, k, n, one,
666 $ a( k*k ), k, b( k, 0 ), ldb )
673 CALL dtrsm(
'L',
'L',
'T', diag, k, n, alpha,
674 $ a( k*k ), k, b( k, 0 ), ldb )
675 CALL dgemm(
'T',
'N', k, n, k, -one, a( 0 ), k,
676 $ b( k, 0 ), ldb, alpha, b, ldb )
677 CALL dtrsm(
'L',
'U',
'N', diag, k, n, one,
678 $ a( k*( k+1 ) ), k, b, ldb )
696 IF( mod( n, 2 ).EQ.0 )
THEN
714 IF( normaltransr )
THEN
727 CALL dtrsm(
'R',
'U',
'T', diag, m, n2, alpha,
728 $ a( n ), n, b( 0, n1 ), ldb )
729 CALL dgemm(
'N',
'N', m, n1, n2, -one, b( 0, n1 ),
730 $ ldb, a( n1 ), n, alpha, b( 0, 0 ),
732 CALL dtrsm(
'R',
'L',
'N', diag, m, n1, one,
733 $ a( 0 ), n, b( 0, 0 ), ldb )
740 CALL dtrsm(
'R',
'L',
'T', diag, m, n1, alpha,
741 $ a( 0 ), n, b( 0, 0 ), ldb )
742 CALL dgemm(
'N',
'T', m, n2, n1, -one, b( 0, 0 ),
743 $ ldb, a( n1 ), n, alpha, b( 0, n1 ),
745 CALL dtrsm(
'R',
'U',
'N', diag, m, n2, one,
746 $ a( n ), n, b( 0, n1 ), ldb )
759 CALL dtrsm(
'R',
'L',
'T', diag, m, n1, alpha,
760 $ a( n2 ), n, b( 0, 0 ), ldb )
761 CALL dgemm(
'N',
'N', m, n2, n1, -one, b( 0, 0 ),
762 $ ldb, a( 0 ), n, alpha, b( 0, n1 ),
764 CALL dtrsm(
'R',
'U',
'N', diag, m, n2, one,
765 $ a( n1 ), n, b( 0, n1 ), ldb )
772 CALL dtrsm(
'R',
'U',
'T', diag, m, n2, alpha,
773 $ a( n1 ), n, b( 0, n1 ), ldb )
774 CALL dgemm(
'N',
'T', m, n1, n2, -one, b( 0, n1 ),
775 $ ldb, a( 0 ), n, alpha, b( 0, 0 ), ldb )
776 CALL dtrsm(
'R',
'L',
'N', diag, m, n1, one,
777 $ a( n2 ), n, b( 0, 0 ), ldb )
796 CALL dtrsm(
'R',
'L',
'N', diag, m, n2, alpha,
797 $ a( 1 ), n1, b( 0, n1 ), ldb )
798 CALL dgemm(
'N',
'T', m, n1, n2, -one, b( 0, n1 ),
799 $ ldb, a( n1*n1 ), n1, alpha, b( 0, 0 ),
801 CALL dtrsm(
'R',
'U',
'T', diag, m, n1, one,
802 $ a( 0 ), n1, b( 0, 0 ), ldb )
809 CALL dtrsm(
'R',
'U',
'N', diag, m, n1, alpha,
810 $ a( 0 ), n1, b( 0, 0 ), ldb )
811 CALL dgemm(
'N',
'N', m, n2, n1, -one, b( 0, 0 ),
812 $ ldb, a( n1*n1 ), n1, alpha, b( 0, n1 ),
814 CALL dtrsm(
'R',
'L',
'T', diag, m, n2, one,
815 $ a( 1 ), n1, b( 0, n1 ), ldb )
828 CALL dtrsm(
'R',
'U',
'N', diag, m, n1, alpha,
829 $ a( n2*n2 ), n2, b( 0, 0 ), ldb )
830 CALL dgemm(
'N',
'T', m, n2, n1, -one, b( 0, 0 ),
831 $ ldb, a( 0 ), n2, alpha, b( 0, n1 ),
833 CALL dtrsm(
'R',
'L',
'T', diag, m, n2, one,
834 $ a( n1*n2 ), n2, b( 0, n1 ), ldb )
841 CALL dtrsm(
'R',
'L',
'N', diag, m, n2, alpha,
842 $ a( n1*n2 ), n2, b( 0, n1 ), ldb )
843 CALL dgemm(
'N',
'N', m, n1, n2, -one, b( 0, n1 ),
844 $ ldb, a( 0 ), n2, alpha, b( 0, 0 ),
846 CALL dtrsm(
'R',
'U',
'T', diag, m, n1, one,
847 $ a( n2*n2 ), n2, b( 0, 0 ), ldb )
859 IF( normaltransr )
THEN
872 CALL dtrsm(
'R',
'U',
'T', diag, m, k, alpha,
873 $ a( 0 ), n+1, b( 0, k ), ldb )
874 CALL dgemm(
'N',
'N', m, k, k, -one, b( 0, k ),
875 $ ldb, a( k+1 ), n+1, alpha, b( 0, 0 ),
877 CALL dtrsm(
'R',
'L',
'N', diag, m, k, one,
878 $ a( 1 ), n+1, b( 0, 0 ), ldb )
885 CALL dtrsm(
'R',
'L',
'T', diag, m, k, alpha,
886 $ a( 1 ), n+1, b( 0, 0 ), ldb )
887 CALL dgemm(
'N',
'T', m, k, k, -one, b( 0, 0 ),
888 $ ldb, a( k+1 ), n+1, alpha, b( 0, k ),
890 CALL dtrsm(
'R',
'U',
'N', diag, m, k, one,
891 $ a( 0 ), n+1, b( 0, k ), ldb )
904 CALL dtrsm(
'R',
'L',
'T', diag, m, k, alpha,
905 $ a( k+1 ), n+1, b( 0, 0 ), ldb )
906 CALL dgemm(
'N',
'N', m, k, k, -one, b( 0, 0 ),
907 $ ldb, a( 0 ), n+1, alpha, b( 0, k ),
909 CALL dtrsm(
'R',
'U',
'N', diag, m, k, one,
910 $ a( k ), n+1, b( 0, k ), ldb )
917 CALL dtrsm(
'R',
'U',
'T', diag, m, k, alpha,
918 $ a( k ), n+1, b( 0, k ), ldb )
919 CALL dgemm(
'N',
'T', m, k, k, -one, b( 0, k ),
920 $ ldb, a( 0 ), n+1, alpha, b( 0, 0 ),
922 CALL dtrsm(
'R',
'L',
'N', diag, m, k, one,
923 $ a( k+1 ), n+1, b( 0, 0 ), ldb )
942 CALL dtrsm(
'R',
'L',
'N', diag, m, k, alpha,
943 $ a( 0 ), k, b( 0, k ), ldb )
944 CALL dgemm(
'N',
'T', m, k, k, -one, b( 0, k ),
945 $ ldb, a( ( k+1 )*k ), k, alpha,
947 CALL dtrsm(
'R',
'U',
'T', diag, m, k, one,
948 $ a( k ), k, b( 0, 0 ), ldb )
955 CALL dtrsm(
'R',
'U',
'N', diag, m, k, alpha,
956 $ a( k ), k, b( 0, 0 ), ldb )
957 CALL dgemm(
'N',
'N', m, k, k, -one, b( 0, 0 ),
958 $ ldb, a( ( k+1 )*k ), k, alpha,
960 CALL dtrsm(
'R',
'L',
'T', diag, m, k, one,
961 $ a( 0 ), k, b( 0, k ), ldb )
974 CALL dtrsm(
'R',
'U',
'N', diag, m, k, alpha,
975 $ a( ( k+1 )*k ), k, b( 0, 0 ), ldb )
976 CALL dgemm(
'N',
'T', m, k, k, -one, b( 0, 0 ),
977 $ ldb, a( 0 ), k, alpha, b( 0, k ), ldb )
978 CALL dtrsm(
'R',
'L',
'T', diag, m, k, one,
979 $ a( k*k ), k, b( 0, k ), ldb )
986 CALL dtrsm(
'R',
'L',
'N', diag, m, k, alpha,
987 $ a( k*k ), k, b( 0, k ), ldb )
988 CALL dgemm(
'N',
'N', m, k, k, -one, b( 0, k ),
989 $ ldb, a( 0 ), k, alpha, b( 0, 0 ), ldb )
990 CALL dtrsm(
'R',
'U',
'T', diag, m, k, one,
991 $ a( ( k+1 )*k ), k, b( 0, 0 ), ldb )
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine dtrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRSM
subroutine dtfsm(TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A, B, LDB)
DTFSM solves a matrix equation (one operand is a triangular matrix in RFP format).