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

◆ ctprfb()

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

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

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

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