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