251 SUBROUTINE stprfb( SIDE, TRANS, DIRECT, STOREV, M, N, K, L,
252 $ v, ldv, t, ldt, a, lda, b, ldb, work, ldwork )
260 CHARACTER DIRECT, SIDE, STOREV, TRANS
261 INTEGER K, L, LDA, LDB, LDT, LDV, LDWORK, M, N
264 REAL A( lda, * ), B( ldb, * ), T( ldt, * ),
265 $ v( ldv, * ), work( ldwork, * )
272 parameter( one = 1.0, zero = 0.0 )
275 INTEGER I, J, MP, NP, KP
276 LOGICAL LEFT, FORWARD, COLUMN, RIGHT, BACKWARD, ROW
289 IF( m.LE.0 .OR. n.LE.0 .OR. k.LE.0 .OR. l.LT.0 )
RETURN
291 IF( lsame( storev,
'C' ) )
THEN
294 ELSE IF ( lsame( storev,
'R' ) )
THEN
302 IF( lsame( side,
'L' ) )
THEN
305 ELSE IF( lsame( side,
'R' ) )
THEN
313 IF( lsame( direct,
'F' ) )
THEN
316 ELSE IF( lsame( direct,
'B' ) )
THEN
326 IF( column .AND. forward .AND. left )
THEN
348 work( i, j ) = b( m-l+i, j )
351 CALL strmm(
'L',
'U',
'T',
'N', l, n, one, v( mp, 1 ), ldv,
353 CALL sgemm(
'T',
'N', l, n, m-l, one, v, ldv, b, ldb,
354 $ one, work, ldwork )
355 CALL sgemm(
'T',
'N', k-l, n, m, one, v( 1, kp ), ldv,
356 $ b, ldb, zero, work( kp, 1 ), ldwork )
360 work( i, j ) = work( i, j ) + a( i, j )
364 CALL strmm(
'L',
'U', trans,
'N', k, n, one, t, ldt,
369 a( i, j ) = a( i, j ) - work( i, j )
373 CALL sgemm(
'N',
'N', m-l, n, k, -one, v, ldv, work, ldwork,
375 CALL sgemm(
'N',
'N', l, n, k-l, -one, v( mp, kp ), ldv,
376 $ work( kp, 1 ), ldwork, one, b( mp, 1 ), ldb )
377 CALL strmm(
'L',
'U',
'N',
'N', l, n, one, v( mp, 1 ), ldv,
381 b( m-l+i, j ) = b( m-l+i, j ) - work( i, j )
387 ELSE IF( column .AND. forward .AND. right )
THEN
408 work( i, j ) = b( i, n-l+j )
411 CALL strmm(
'R',
'U',
'N',
'N', m, l, one, v( np, 1 ), ldv,
413 CALL sgemm(
'N',
'N', m, l, n-l, one, b, ldb,
414 $ v, ldv, one, work, ldwork )
415 CALL sgemm(
'N',
'N', m, k-l, n, one, b, ldb,
416 $ v( 1, kp ), ldv, zero, work( 1, kp ), ldwork )
420 work( i, j ) = work( i, j ) + a( i, j )
424 CALL strmm(
'R',
'U', trans,
'N', m, k, one, t, ldt,
429 a( i, j ) = a( i, j ) - work( i, j )
433 CALL sgemm(
'N',
'T', m, n-l, k, -one, work, ldwork,
434 $ v, ldv, one, b, ldb )
435 CALL sgemm(
'N',
'T', m, l, k-l, -one, work( 1, kp ), ldwork,
436 $ v( np, kp ), ldv, one, b( 1, np ), ldb )
437 CALL strmm(
'R',
'U',
'T',
'N', m, l, one, v( np, 1 ), ldv,
441 b( i, n-l+j ) = b( i, n-l+j ) - work( i, j )
447 ELSE IF( column .AND. backward .AND. left )
THEN
469 work( k-l+i, j ) = b( i, j )
473 CALL strmm(
'L',
'L',
'T',
'N', l, n, one, v( 1, kp ), ldv,
474 $ work( kp, 1 ), ldwork )
475 CALL sgemm(
'T',
'N', l, n, m-l, one, v( mp, kp ), ldv,
476 $ b( mp, 1 ), ldb, one, work( kp, 1 ), ldwork )
477 CALL sgemm(
'T',
'N', k-l, n, m, one, v, ldv,
478 $ b, ldb, zero, work, ldwork )
482 work( i, j ) = work( i, j ) + a( i, j )
486 CALL strmm(
'L',
'L', trans,
'N', k, n, one, t, ldt,
491 a( i, j ) = a( i, j ) - work( i, j )
495 CALL sgemm(
'N',
'N', m-l, n, k, -one, v( mp, 1 ), ldv,
496 $ work, ldwork, one, b( mp, 1 ), ldb )
497 CALL sgemm(
'N',
'N', l, n, k-l, -one, v, ldv,
498 $ work, ldwork, one, b, ldb )
499 CALL strmm(
'L',
'L',
'N',
'N', l, n, one, v( 1, kp ), ldv,
500 $ work( kp, 1 ), ldwork )
503 b( i, j ) = b( i, j ) - work( k-l+i, j )
509 ELSE IF( column .AND. backward .AND. right )
THEN
530 work( i, k-l+j ) = b( i, j )
533 CALL strmm(
'R',
'L',
'N',
'N', m, l, one, v( 1, kp ), ldv,
534 $ work( 1, kp ), ldwork )
535 CALL sgemm(
'N',
'N', m, l, n-l, one, b( 1, np ), ldb,
536 $ v( np, kp ), ldv, one, work( 1, kp ), ldwork )
537 CALL sgemm(
'N',
'N', m, k-l, n, one, b, ldb,
538 $ v, ldv, zero, work, ldwork )
542 work( i, j ) = work( i, j ) + a( i, j )
546 CALL strmm(
'R',
'L', trans,
'N', m, k, one, t, ldt,
551 a( i, j ) = a( i, j ) - work( i, j )
555 CALL sgemm(
'N',
'T', m, n-l, k, -one, work, ldwork,
556 $ v( np, 1 ), ldv, one, b( 1, np ), ldb )
557 CALL sgemm(
'N',
'T', m, l, k-l, -one, work, ldwork,
558 $ v, ldv, one, b, ldb )
559 CALL strmm(
'R',
'L',
'T',
'N', m, l, one, v( 1, kp ), ldv,
560 $ work( 1, kp ), ldwork )
563 b( i, j ) = b( i, j ) - work( i, k-l+j )
569 ELSE IF( row .AND. forward .AND. left )
THEN
590 work( i, j ) = b( m-l+i, j )
593 CALL strmm(
'L',
'L',
'N',
'N', l, n, one, v( 1, mp ), ldv,
595 CALL sgemm(
'N',
'N', l, n, m-l, one, v, ldv,b, ldb,
596 $ one, work, ldwork )
597 CALL sgemm(
'N',
'N', k-l, n, m, one, v( kp, 1 ), ldv,
598 $ b, ldb, zero, work( kp, 1 ), ldwork )
602 work( i, j ) = work( i, j ) + a( i, j )
606 CALL strmm(
'L',
'U', trans,
'N', k, n, one, t, ldt,
611 a( i, j ) = a( i, j ) - work( i, j )
615 CALL sgemm(
'T',
'N', m-l, n, k, -one, v, ldv, work, ldwork,
617 CALL sgemm(
'T',
'N', l, n, k-l, -one, v( kp, mp ), ldv,
618 $ work( kp, 1 ), ldwork, one, b( mp, 1 ), ldb )
619 CALL strmm(
'L',
'L',
'T',
'N', l, n, one, v( 1, mp ), ldv,
623 b( m-l+i, j ) = b( m-l+i, j ) - work( i, j )
629 ELSE IF( row .AND. forward .AND. right )
THEN
649 work( i, j ) = b( i, n-l+j )
652 CALL strmm(
'R',
'L',
'T',
'N', m, l, one, v( 1, np ), ldv,
654 CALL sgemm(
'N',
'T', m, l, n-l, one, b, ldb, v, ldv,
655 $ one, work, ldwork )
656 CALL sgemm(
'N',
'T', m, k-l, n, one, b, ldb,
657 $ v( kp, 1 ), ldv, zero, work( 1, kp ), ldwork )
661 work( i, j ) = work( i, j ) + a( i, j )
665 CALL strmm(
'R',
'U', trans,
'N', m, k, one, t, ldt,
670 a( i, j ) = a( i, j ) - work( i, j )
674 CALL sgemm(
'N',
'N', m, n-l, k, -one, work, ldwork,
675 $ v, ldv, one, b, ldb )
676 CALL sgemm(
'N',
'N', m, l, k-l, -one, work( 1, kp ), ldwork,
677 $ v( kp, np ), ldv, one, b( 1, np ), ldb )
678 CALL strmm(
'R',
'L',
'N',
'N', m, l, one, v( 1, np ), ldv,
682 b( i, n-l+j ) = b( i, n-l+j ) - work( i, j )
688 ELSE IF( row .AND. backward .AND. left )
THEN
709 work( k-l+i, j ) = b( i, j )
712 CALL strmm(
'L',
'U',
'N',
'N', l, n, one, v( kp, 1 ), ldv,
713 $ work( kp, 1 ), ldwork )
714 CALL sgemm(
'N',
'N', l, n, m-l, one, v( kp, mp ), ldv,
715 $ b( mp, 1 ), ldb, one, work( kp, 1 ), ldwork )
716 CALL sgemm(
'N',
'N', k-l, n, m, one, v, ldv, b, ldb,
717 $ zero, work, ldwork )
721 work( i, j ) = work( i, j ) + a( i, j )
725 CALL strmm(
'L',
'L ', trans,
'N', k, n, one, t, ldt,
730 a( i, j ) = a( i, j ) - work( i, j )
734 CALL sgemm(
'T',
'N', m-l, n, k, -one, v( 1, mp ), ldv,
735 $ work, ldwork, one, b( mp, 1 ), ldb )
736 CALL sgemm(
'T',
'N', l, n, k-l, -one, v, ldv,
737 $ work, ldwork, one, b, ldb )
738 CALL strmm(
'L',
'U',
'T',
'N', l, n, one, v( kp, 1 ), ldv,
739 $ work( kp, 1 ), ldwork )
742 b( i, j ) = b( i, j ) - work( k-l+i, j )
748 ELSE IF( row .AND. backward .AND. right )
THEN
768 work( i, k-l+j ) = b( i, j )
771 CALL strmm(
'R',
'U',
'T',
'N', m, l, one, v( kp, 1 ), ldv,
772 $ work( 1, kp ), ldwork )
773 CALL sgemm(
'N',
'T', m, l, n-l, one, b( 1, np ), ldb,
774 $ v( kp, np ), ldv, one, work( 1, kp ), ldwork )
775 CALL sgemm(
'N',
'T', m, k-l, n, one, b, ldb, v, ldv,
776 $ zero, work, ldwork )
780 work( i, j ) = work( i, j ) + a( i, j )
784 CALL strmm(
'R',
'L', trans,
'N', m, k, one, t, ldt,
789 a( i, j ) = a( i, j ) - work( i, j )
793 CALL sgemm(
'N',
'N', m, n-l, k, -one, work, ldwork,
794 $ v( 1, np ), ldv, one, b( 1, np ), ldb )
795 CALL sgemm(
'N',
'N', m, l, k-l , -one, work, ldwork,
796 $ v, ldv, one, b, ldb )
797 CALL strmm(
'R',
'U',
'N',
'N', m, l, one, v( kp, 1 ), ldv,
798 $ work( 1, kp ), ldwork )
801 b( i, j ) = b( i, j ) - work( i, k-l+j )
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
subroutine stprfb(SIDE, TRANS, DIRECT, STOREV, M, N, K, L, V, LDV, T, LDT, A, LDA, B, LDB, WORK, LDWORK)
STPRFB applies a real or complex "triangular-pentagonal" blocked reflector to a real or complex matri...
subroutine strmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRMM