LAPACK 3.11.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ stprfb()

 subroutine stprfb ( character SIDE, character TRANS, character DIRECT, character STOREV, integer M, integer N, integer K, integer L, real, dimension( ldv, * ) V, integer LDV, real, dimension( ldt, * ) T, integer LDT, real, dimension( lda, * ) A, integer LDA, real, dimension( ldb, * ) B, integer LDB, real, dimension( ldwork, * ) WORK, integer LDWORK )

STPRFB applies a real "triangular-pentagonal" block reflector to a real matrix, which is composed of two blocks.

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

Purpose:
``` STPRFB applies a real "triangular-pentagonal" block reflector H or its
transpose H**T to a real matrix C, which is composed of two
blocks A and B, either from the left or right.```
Parameters
 [in] SIDE ``` SIDE is CHARACTER*1 = 'L': apply H or H**T from the Left = 'R': apply H or H**T from the Right``` [in] TRANS ``` TRANS is CHARACTER*1 = 'N': apply H (No transpose) = 'T': apply H**T (Transpose)``` [in] DIRECT ``` DIRECT is CHARACTER*1 Indicates how H is formed from a product of elementary reflectors = 'F': H = H(1) H(2) . . . H(k) (Forward) = 'B': H = H(k) . . . H(2) H(1) (Backward)``` [in] STOREV ``` STOREV is CHARACTER*1 Indicates how the vectors which define the elementary reflectors are stored: = 'C': Columns = 'R': Rows``` [in] M ``` M is INTEGER The number of rows of the matrix B. M >= 0.``` [in] N ``` N is INTEGER The number of columns of the matrix B. N >= 0.``` [in] K ``` K is INTEGER The order of the matrix T, i.e. the number of elementary reflectors whose product defines the block reflector. K >= 0.``` [in] L ``` L is INTEGER The order of the trapezoidal part of V. K >= L >= 0. See Further Details.``` [in] V ``` V is REAL array, dimension (LDV,K) if STOREV = 'C' (LDV,M) if STOREV = 'R' and SIDE = 'L' (LDV,N) if STOREV = 'R' and SIDE = 'R' The pentagonal matrix V, which contains the elementary reflectors H(1), H(2), ..., H(K). See Further Details.``` [in] LDV ``` LDV is INTEGER The leading dimension of the array V. If STOREV = 'C' and SIDE = 'L', LDV >= max(1,M); if STOREV = 'C' and SIDE = 'R', LDV >= max(1,N); if STOREV = 'R', LDV >= K.``` [in] T ``` T is REAL array, dimension (LDT,K) The triangular K-by-K matrix T in the representation of the block reflector.``` [in] LDT ``` LDT is INTEGER The leading dimension of the array T. LDT >= K.``` [in,out] A ``` A is REAL array, dimension (LDA,N) if SIDE = 'L' or (LDA,K) if SIDE = 'R' On entry, the K-by-N or M-by-K matrix A. On exit, A is overwritten by the corresponding block of H*C or H**T*C or C*H or C*H**T. See Further Details.``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. If SIDE = 'L', LDA >= max(1,K); If SIDE = 'R', LDA >= max(1,M).``` [in,out] B ``` B is REAL array, dimension (LDB,N) On entry, the M-by-N matrix B. On exit, B is overwritten by the corresponding block of H*C or H**T*C or C*H or C*H**T. See Further Details.``` [in] LDB ``` LDB is INTEGER The leading dimension of the array B. LDB >= max(1,M).``` [out] WORK ``` WORK is REAL array, dimension (LDWORK,N) if SIDE = 'L', (LDWORK,K) if SIDE = 'R'.``` [in] LDWORK ``` LDWORK is INTEGER The leading dimension of the array WORK. If SIDE = 'L', LDWORK >= K; if SIDE = 'R', LDWORK >= M.```
Further Details:
```  The matrix C is a composite matrix formed from blocks A and B.
The block B is of size M-by-N; if SIDE = 'R', A is of size M-by-K,
and if SIDE = 'L', A is of size K-by-N.

If SIDE = 'R' and DIRECT = 'F', C = [A B].

If SIDE = 'L' and DIRECT = 'F', C = [A]
[B].

If SIDE = 'R' and DIRECT = 'B', C = [B A].

If SIDE = 'L' and DIRECT = 'B', C = [B]
[A].

The pentagonal matrix V is composed of a rectangular block V1 and a
trapezoidal block V2.  The size of the trapezoidal block is determined by
the parameter L, where 0<=L<=K.  If L=K, the V2 block of V is triangular;
if L=0, there is no trapezoidal block, thus V = V1 is rectangular.

If DIRECT = 'F' and STOREV = 'C':  V = [V1]
[V2]
- V2 is upper trapezoidal (first L rows of K-by-K upper triangular)

If DIRECT = 'F' and STOREV = 'R':  V = [V1 V2]

- V2 is lower trapezoidal (first L columns of K-by-K lower triangular)

If DIRECT = 'B' and STOREV = 'C':  V = [V2]
[V1]
- V2 is lower trapezoidal (last L rows of K-by-K lower triangular)

If DIRECT = 'B' and STOREV = 'R':  V = [V2 V1]

- V2 is upper trapezoidal (last L columns of K-by-K upper triangular)

If STOREV = 'C' and SIDE = 'L', V is M-by-K with V2 L-by-K.

If STOREV = 'C' and SIDE = 'R', V is N-by-K with V2 L-by-K.

If STOREV = 'R' and SIDE = 'L', V is K-by-M with V2 K-by-L.

If STOREV = 'R' and SIDE = 'R', V is K-by-N with V2 K-by-L.```

Definition at line 249 of file stprfb.f.

251*
252* -- LAPACK auxiliary routine --
253* -- LAPACK is a software package provided by Univ. of Tennessee, --
254* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
255*
256* .. Scalar Arguments ..
257 CHARACTER DIRECT, SIDE, STOREV, TRANS
258 INTEGER K, L, LDA, LDB, LDT, LDV, LDWORK, M, N
259* ..
260* .. Array Arguments ..
261 REAL A( LDA, * ), B( LDB, * ), T( LDT, * ),
262 \$ V( LDV, * ), WORK( LDWORK, * )
263* ..
264*
265* ==========================================================================
266*
267* .. Parameters ..
268 REAL ONE, ZERO
269 parameter( one = 1.0, zero = 0.0 )
270* ..
271* .. Local Scalars ..
272 INTEGER I, J, MP, NP, KP
273 LOGICAL LEFT, FORWARD, COLUMN, RIGHT, BACKWARD, ROW
274* ..
275* .. External Functions ..
276 LOGICAL LSAME
277 EXTERNAL lsame
278* ..
279* .. External Subroutines ..
280 EXTERNAL sgemm, strmm
281* ..
282* .. Executable Statements ..
283*
284* Quick return if possible
285*
286 IF( m.LE.0 .OR. n.LE.0 .OR. k.LE.0 .OR. l.LT.0 ) RETURN
287*
288 IF( lsame( storev, 'C' ) ) THEN
289 column = .true.
290 row = .false.
291 ELSE IF ( lsame( storev, 'R' ) ) THEN
292 column = .false.
293 row = .true.
294 ELSE
295 column = .false.
296 row = .false.
297 END IF
298*
299 IF( lsame( side, 'L' ) ) THEN
300 left = .true.
301 right = .false.
302 ELSE IF( lsame( side, 'R' ) ) THEN
303 left = .false.
304 right = .true.
305 ELSE
306 left = .false.
307 right = .false.
308 END IF
309*
310 IF( lsame( direct, 'F' ) ) THEN
311 forward = .true.
312 backward = .false.
313 ELSE IF( lsame( direct, 'B' ) ) THEN
314 forward = .false.
315 backward = .true.
316 ELSE
317 forward = .false.
318 backward = .false.
319 END IF
320*
321* ---------------------------------------------------------------------------
322*
323 IF( column .AND. forward .AND. left ) THEN
324*
325* ---------------------------------------------------------------------------
326*
327* Let W = [ I ] (K-by-K)
328* [ V ] (M-by-K)
329*
330* Form H C or H**T C where C = [ A ] (K-by-N)
331* [ B ] (M-by-N)
332*
333* H = I - W T W**T or H**T = I - W T**T W**T
334*
335* A = A - T (A + V**T B) or A = A - T**T (A + V**T B)
336* B = B - V T (A + V**T B) or B = B - V T**T (A + V**T B)
337*
338* ---------------------------------------------------------------------------
339*
340 mp = min( m-l+1, m )
341 kp = min( l+1, k )
342*
343 DO j = 1, n
344 DO i = 1, l
345 work( i, j ) = b( m-l+i, j )
346 END DO
347 END DO
348 CALL strmm( 'L', 'U', 'T', 'N', l, n, one, v( mp, 1 ), ldv,
349 \$ work, ldwork )
350 CALL sgemm( 'T', 'N', l, n, m-l, one, v, ldv, b, ldb,
351 \$ one, work, ldwork )
352 CALL sgemm( 'T', 'N', k-l, n, m, one, v( 1, kp ), ldv,
353 \$ b, ldb, zero, work( kp, 1 ), ldwork )
354*
355 DO j = 1, n
356 DO i = 1, k
357 work( i, j ) = work( i, j ) + a( i, j )
358 END DO
359 END DO
360*
361 CALL strmm( 'L', 'U', trans, 'N', k, n, one, t, ldt,
362 \$ work, ldwork )
363*
364 DO j = 1, n
365 DO i = 1, k
366 a( i, j ) = a( i, j ) - work( i, j )
367 END DO
368 END DO
369*
370 CALL sgemm( 'N', 'N', m-l, n, k, -one, v, ldv, work, ldwork,
371 \$ one, b, ldb )
372 CALL sgemm( 'N', 'N', l, n, k-l, -one, v( mp, kp ), ldv,
373 \$ work( kp, 1 ), ldwork, one, b( mp, 1 ), ldb )
374 CALL strmm( 'L', 'U', 'N', 'N', l, n, one, v( mp, 1 ), ldv,
375 \$ work, ldwork )
376 DO j = 1, n
377 DO i = 1, l
378 b( m-l+i, j ) = b( m-l+i, j ) - work( i, j )
379 END DO
380 END DO
381*
382* ---------------------------------------------------------------------------
383*
384 ELSE IF( column .AND. forward .AND. right ) THEN
385*
386* ---------------------------------------------------------------------------
387*
388* Let W = [ I ] (K-by-K)
389* [ V ] (N-by-K)
390*
391* Form C H or C H**T where C = [ A B ] (A is M-by-K, B is M-by-N)
392*
393* H = I - W T W**T or H**T = I - W T**T W**T
394*
395* A = A - (A + B V) T or A = A - (A + B V) T**T
396* B = B - (A + B V) T V**T or B = B - (A + B V) T**T V**T
397*
398* ---------------------------------------------------------------------------
399*
400 np = min( n-l+1, n )
401 kp = min( l+1, k )
402*
403 DO j = 1, l
404 DO i = 1, m
405 work( i, j ) = b( i, n-l+j )
406 END DO
407 END DO
408 CALL strmm( 'R', 'U', 'N', 'N', m, l, one, v( np, 1 ), ldv,
409 \$ work, ldwork )
410 CALL sgemm( 'N', 'N', m, l, n-l, one, b, ldb,
411 \$ v, ldv, one, work, ldwork )
412 CALL sgemm( 'N', 'N', m, k-l, n, one, b, ldb,
413 \$ v( 1, kp ), ldv, zero, work( 1, kp ), ldwork )
414*
415 DO j = 1, k
416 DO i = 1, m
417 work( i, j ) = work( i, j ) + a( i, j )
418 END DO
419 END DO
420*
421 CALL strmm( 'R', 'U', trans, 'N', m, k, one, t, ldt,
422 \$ work, ldwork )
423*
424 DO j = 1, k
425 DO i = 1, m
426 a( i, j ) = a( i, j ) - work( i, j )
427 END DO
428 END DO
429*
430 CALL sgemm( 'N', 'T', m, n-l, k, -one, work, ldwork,
431 \$ v, ldv, one, b, ldb )
432 CALL sgemm( 'N', 'T', m, l, k-l, -one, work( 1, kp ), ldwork,
433 \$ v( np, kp ), ldv, one, b( 1, np ), ldb )
434 CALL strmm( 'R', 'U', 'T', 'N', m, l, one, v( np, 1 ), ldv,
435 \$ work, ldwork )
436 DO j = 1, l
437 DO i = 1, m
438 b( i, n-l+j ) = b( i, n-l+j ) - work( i, j )
439 END DO
440 END DO
441*
442* ---------------------------------------------------------------------------
443*
444 ELSE IF( column .AND. backward .AND. left ) THEN
445*
446* ---------------------------------------------------------------------------
447*
448* Let W = [ V ] (M-by-K)
449* [ I ] (K-by-K)
450*
451* Form H C or H**T C where C = [ B ] (M-by-N)
452* [ A ] (K-by-N)
453*
454* H = I - W T W**T or H**T = I - W T**T W**T
455*
456* A = A - T (A + V**T B) or A = A - T**T (A + V**T B)
457* B = B - V T (A + V**T B) or B = B - V T**T (A + V**T B)
458*
459* ---------------------------------------------------------------------------
460*
461 mp = min( l+1, m )
462 kp = min( k-l+1, k )
463*
464 DO j = 1, n
465 DO i = 1, l
466 work( k-l+i, j ) = b( i, j )
467 END DO
468 END DO
469*
470 CALL strmm( 'L', 'L', 'T', 'N', l, n, one, v( 1, kp ), ldv,
471 \$ work( kp, 1 ), ldwork )
472 CALL sgemm( 'T', 'N', l, n, m-l, one, v( mp, kp ), ldv,
473 \$ b( mp, 1 ), ldb, one, work( kp, 1 ), ldwork )
474 CALL sgemm( 'T', 'N', k-l, n, m, one, v, ldv,
475 \$ b, ldb, zero, work, ldwork )
476*
477 DO j = 1, n
478 DO i = 1, k
479 work( i, j ) = work( i, j ) + a( i, j )
480 END DO
481 END DO
482*
483 CALL strmm( 'L', 'L', trans, 'N', k, n, one, t, ldt,
484 \$ work, ldwork )
485*
486 DO j = 1, n
487 DO i = 1, k
488 a( i, j ) = a( i, j ) - work( i, j )
489 END DO
490 END DO
491*
492 CALL sgemm( 'N', 'N', m-l, n, k, -one, v( mp, 1 ), ldv,
493 \$ work, ldwork, one, b( mp, 1 ), ldb )
494 CALL sgemm( 'N', 'N', l, n, k-l, -one, v, ldv,
495 \$ work, ldwork, one, b, ldb )
496 CALL strmm( 'L', 'L', 'N', 'N', l, n, one, v( 1, kp ), ldv,
497 \$ work( kp, 1 ), ldwork )
498 DO j = 1, n
499 DO i = 1, l
500 b( i, j ) = b( i, j ) - work( k-l+i, j )
501 END DO
502 END DO
503*
504* ---------------------------------------------------------------------------
505*
506 ELSE IF( column .AND. backward .AND. right ) THEN
507*
508* ---------------------------------------------------------------------------
509*
510* Let W = [ V ] (N-by-K)
511* [ I ] (K-by-K)
512*
513* Form C H or C H**T where C = [ B A ] (B is M-by-N, A is M-by-K)
514*
515* H = I - W T W**T or H**T = I - W T**T W**T
516*
517* A = A - (A + B V) T or A = A - (A + B V) T**T
518* B = B - (A + B V) T V**T or B = B - (A + B V) T**T V**T
519*
520* ---------------------------------------------------------------------------
521*
522 np = min( l+1, n )
523 kp = min( k-l+1, k )
524*
525 DO j = 1, l
526 DO i = 1, m
527 work( i, k-l+j ) = b( i, j )
528 END DO
529 END DO
530 CALL strmm( 'R', 'L', 'N', 'N', m, l, one, v( 1, kp ), ldv,
531 \$ work( 1, kp ), ldwork )
532 CALL sgemm( 'N', 'N', m, l, n-l, one, b( 1, np ), ldb,
533 \$ v( np, kp ), ldv, one, work( 1, kp ), ldwork )
534 CALL sgemm( 'N', 'N', m, k-l, n, one, b, ldb,
535 \$ v, ldv, zero, work, ldwork )
536*
537 DO j = 1, k
538 DO i = 1, m
539 work( i, j ) = work( i, j ) + a( i, j )
540 END DO
541 END DO
542*
543 CALL strmm( 'R', 'L', trans, 'N', m, k, one, t, ldt,
544 \$ work, ldwork )
545*
546 DO j = 1, k
547 DO i = 1, m
548 a( i, j ) = a( i, j ) - work( i, j )
549 END DO
550 END DO
551*
552 CALL sgemm( 'N', 'T', m, n-l, k, -one, work, ldwork,
553 \$ v( np, 1 ), ldv, one, b( 1, np ), ldb )
554 CALL sgemm( 'N', 'T', m, l, k-l, -one, work, ldwork,
555 \$ v, ldv, one, b, ldb )
556 CALL strmm( 'R', 'L', 'T', 'N', m, l, one, v( 1, kp ), ldv,
557 \$ work( 1, kp ), ldwork )
558 DO j = 1, l
559 DO i = 1, m
560 b( i, j ) = b( i, j ) - work( i, k-l+j )
561 END DO
562 END DO
563*
564* ---------------------------------------------------------------------------
565*
566 ELSE IF( row .AND. forward .AND. left ) THEN
567*
568* ---------------------------------------------------------------------------
569*
570* Let W = [ I V ] ( I is K-by-K, V is K-by-M )
571*
572* Form H C or H**T C where C = [ A ] (K-by-N)
573* [ B ] (M-by-N)
574*
575* H = I - W**T T W or H**T = I - W**T T**T W
576*
577* A = A - T (A + V B) or A = A - T**T (A + V B)
578* B = B - V**T T (A + V B) or B = B - V**T T**T (A + V B)
579*
580* ---------------------------------------------------------------------------
581*
582 mp = min( m-l+1, m )
583 kp = min( l+1, k )
584*
585 DO j = 1, n
586 DO i = 1, l
587 work( i, j ) = b( m-l+i, j )
588 END DO
589 END DO
590 CALL strmm( 'L', 'L', 'N', 'N', l, n, one, v( 1, mp ), ldv,
591 \$ work, ldb )
592 CALL sgemm( 'N', 'N', l, n, m-l, one, v, ldv,b, ldb,
593 \$ one, work, ldwork )
594 CALL sgemm( 'N', 'N', k-l, n, m, one, v( kp, 1 ), ldv,
595 \$ b, ldb, zero, work( kp, 1 ), ldwork )
596*
597 DO j = 1, n
598 DO i = 1, k
599 work( i, j ) = work( i, j ) + a( i, j )
600 END DO
601 END DO
602*
603 CALL strmm( 'L', 'U', trans, 'N', k, n, one, t, ldt,
604 \$ work, ldwork )
605*
606 DO j = 1, n
607 DO i = 1, k
608 a( i, j ) = a( i, j ) - work( i, j )
609 END DO
610 END DO
611*
612 CALL sgemm( 'T', 'N', m-l, n, k, -one, v, ldv, work, ldwork,
613 \$ one, b, ldb )
614 CALL sgemm( 'T', 'N', l, n, k-l, -one, v( kp, mp ), ldv,
615 \$ work( kp, 1 ), ldwork, one, b( mp, 1 ), ldb )
616 CALL strmm( 'L', 'L', 'T', 'N', l, n, one, v( 1, mp ), ldv,
617 \$ work, ldwork )
618 DO j = 1, n
619 DO i = 1, l
620 b( m-l+i, j ) = b( m-l+i, j ) - work( i, j )
621 END DO
622 END DO
623*
624* ---------------------------------------------------------------------------
625*
626 ELSE IF( row .AND. forward .AND. right ) THEN
627*
628* ---------------------------------------------------------------------------
629*
630* Let W = [ I V ] ( I is K-by-K, V is K-by-N )
631*
632* Form C H or C H**T where C = [ A B ] (A is M-by-K, B is M-by-N)
633*
634* H = I - W**T T W or H**T = I - W**T T**T W
635*
636* A = A - (A + B V**T) T or A = A - (A + B V**T) T**T
637* B = B - (A + B V**T) T V or B = B - (A + B V**T) T**T V
638*
639* ---------------------------------------------------------------------------
640*
641 np = min( n-l+1, n )
642 kp = min( l+1, k )
643*
644 DO j = 1, l
645 DO i = 1, m
646 work( i, j ) = b( i, n-l+j )
647 END DO
648 END DO
649 CALL strmm( 'R', 'L', 'T', 'N', m, l, one, v( 1, np ), ldv,
650 \$ work, ldwork )
651 CALL sgemm( 'N', 'T', m, l, n-l, one, b, ldb, v, ldv,
652 \$ one, work, ldwork )
653 CALL sgemm( 'N', 'T', m, k-l, n, one, b, ldb,
654 \$ v( kp, 1 ), ldv, zero, work( 1, kp ), ldwork )
655*
656 DO j = 1, k
657 DO i = 1, m
658 work( i, j ) = work( i, j ) + a( i, j )
659 END DO
660 END DO
661*
662 CALL strmm( 'R', 'U', trans, 'N', m, k, one, t, ldt,
663 \$ work, ldwork )
664*
665 DO j = 1, k
666 DO i = 1, m
667 a( i, j ) = a( i, j ) - work( i, j )
668 END DO
669 END DO
670*
671 CALL sgemm( 'N', 'N', m, n-l, k, -one, work, ldwork,
672 \$ v, ldv, one, b, ldb )
673 CALL sgemm( 'N', 'N', m, l, k-l, -one, work( 1, kp ), ldwork,
674 \$ v( kp, np ), ldv, one, b( 1, np ), ldb )
675 CALL strmm( 'R', 'L', 'N', 'N', m, l, one, v( 1, np ), ldv,
676 \$ work, ldwork )
677 DO j = 1, l
678 DO i = 1, m
679 b( i, n-l+j ) = b( i, n-l+j ) - work( i, j )
680 END DO
681 END DO
682*
683* ---------------------------------------------------------------------------
684*
685 ELSE IF( row .AND. backward .AND. left ) THEN
686*
687* ---------------------------------------------------------------------------
688*
689* Let W = [ V I ] ( I is K-by-K, V is K-by-M )
690*
691* Form H C or H**T C where C = [ B ] (M-by-N)
692* [ A ] (K-by-N)
693*
694* H = I - W**T T W or H**T = I - W**T T**T W
695*
696* A = A - T (A + V B) or A = A - T**T (A + V B)
697* B = B - V**T T (A + V B) or B = B - V**T T**T (A + V B)
698*
699* ---------------------------------------------------------------------------
700*
701 mp = min( l+1, m )
702 kp = min( k-l+1, k )
703*
704 DO j = 1, n
705 DO i = 1, l
706 work( k-l+i, j ) = b( i, j )
707 END DO
708 END DO
709 CALL strmm( 'L', 'U', 'N', 'N', l, n, one, v( kp, 1 ), ldv,
710 \$ work( kp, 1 ), ldwork )
711 CALL sgemm( 'N', 'N', l, n, m-l, one, v( kp, mp ), ldv,
712 \$ b( mp, 1 ), ldb, one, work( kp, 1 ), ldwork )
713 CALL sgemm( 'N', 'N', k-l, n, m, one, v, ldv, b, ldb,
714 \$ zero, work, ldwork )
715*
716 DO j = 1, n
717 DO i = 1, k
718 work( i, j ) = work( i, j ) + a( i, j )
719 END DO
720 END DO
721*
722 CALL strmm( 'L', 'L ', trans, 'N', k, n, one, t, ldt,
723 \$ work, ldwork )
724*
725 DO j = 1, n
726 DO i = 1, k
727 a( i, j ) = a( i, j ) - work( i, j )
728 END DO
729 END DO
730*
731 CALL sgemm( 'T', 'N', m-l, n, k, -one, v( 1, mp ), ldv,
732 \$ work, ldwork, one, b( mp, 1 ), ldb )
733 CALL sgemm( 'T', 'N', l, n, k-l, -one, v, ldv,
734 \$ work, ldwork, one, b, ldb )
735 CALL strmm( 'L', 'U', 'T', 'N', l, n, one, v( kp, 1 ), ldv,
736 \$ work( kp, 1 ), ldwork )
737 DO j = 1, n
738 DO i = 1, l
739 b( i, j ) = b( i, j ) - work( k-l+i, j )
740 END DO
741 END DO
742*
743* ---------------------------------------------------------------------------
744*
745 ELSE IF( row .AND. backward .AND. right ) THEN
746*
747* ---------------------------------------------------------------------------
748*
749* Let W = [ V I ] ( I is K-by-K, V is K-by-N )
750*
751* Form C H or C H**T where C = [ B A ] (A is M-by-K, B is M-by-N)
752*
753* H = I - W**T T W or H**T = I - W**T T**T W
754*
755* A = A - (A + B V**T) T or A = A - (A + B V**T) T**T
756* B = B - (A + B V**T) T V or B = B - (A + B V**T) T**T V
757*
758* ---------------------------------------------------------------------------
759*
760 np = min( l+1, n )
761 kp = min( k-l+1, k )
762*
763 DO j = 1, l
764 DO i = 1, m
765 work( i, k-l+j ) = b( i, j )
766 END DO
767 END DO
768 CALL strmm( 'R', 'U', 'T', 'N', m, l, one, v( kp, 1 ), ldv,
769 \$ work( 1, kp ), ldwork )
770 CALL sgemm( 'N', 'T', m, l, n-l, one, b( 1, np ), ldb,
771 \$ v( kp, np ), ldv, one, work( 1, kp ), ldwork )
772 CALL sgemm( 'N', 'T', m, k-l, n, one, b, ldb, v, ldv,
773 \$ zero, work, ldwork )
774*
775 DO j = 1, k
776 DO i = 1, m
777 work( i, j ) = work( i, j ) + a( i, j )
778 END DO
779 END DO
780*
781 CALL strmm( 'R', 'L', trans, 'N', m, k, one, t, ldt,
782 \$ work, ldwork )
783*
784 DO j = 1, k
785 DO i = 1, m
786 a( i, j ) = a( i, j ) - work( i, j )
787 END DO
788 END DO
789*
790 CALL sgemm( 'N', 'N', m, n-l, k, -one, work, ldwork,
791 \$ v( 1, np ), ldv, one, b( 1, np ), ldb )
792 CALL sgemm( 'N', 'N', m, l, k-l , -one, work, ldwork,
793 \$ v, ldv, one, b, ldb )
794 CALL strmm( 'R', 'U', 'N', 'N', m, l, one, v( kp, 1 ), ldv,
795 \$ work( 1, kp ), ldwork )
796 DO j = 1, l
797 DO i = 1, m
798 b( i, j ) = b( i, j ) - work( i, k-l+j )
799 END DO
800 END DO
801*
802 END IF
803*
804 RETURN
805*
806* End of STPRFB
807*
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine strmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRMM
Definition: strmm.f:177
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
Definition: sgemm.f:187
Here is the call graph for this function:
Here is the caller graph for this function: