LAPACK  3.9.0
LAPACK: Linear Algebra PACKage

◆ zlarfb()

subroutine zlarfb ( character  SIDE,
character  TRANS,
character  DIRECT,
character  STOREV,
integer  M,
integer  N,
integer  K,
complex*16, dimension( ldv, * )  V,
integer  LDV,
complex*16, dimension( ldt, * )  T,
integer  LDT,
complex*16, dimension( ldc, * )  C,
integer  LDC,
complex*16, dimension( ldwork, * )  WORK,
integer  LDWORK 
)

ZLARFB applies a block reflector or its conjugate-transpose to a general rectangular matrix.

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

Purpose:
 ZLARFB applies a complex block reflector H or its transpose H**H to a
 complex M-by-N matrix C, from either the left or the 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': Columnwise
          = 'R': Rowwise
[in]M
          M is INTEGER
          The number of rows of the matrix C.
[in]N
          N is INTEGER
          The number of columns of the matrix C.
[in]K
          K is INTEGER
          The order of the matrix T (= the number of elementary
          reflectors whose product defines the block reflector).
          If SIDE = 'L', M >= K >= 0;
          if SIDE = 'R', N >= K >= 0.
[in]V
          V is COMPLEX*16 array, dimension
                                (LDV,K) if STOREV = 'C'
                                (LDV,M) if STOREV = 'R' and SIDE = 'L'
                                (LDV,N) if STOREV = 'R' and SIDE = 'R'
          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*16 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]C
          C is COMPLEX*16 array, dimension (LDC,N)
          On entry, the M-by-N matrix C.
          On exit, C is overwritten by H*C or H**H*C or C*H or C*H**H.
[in]LDC
          LDC is INTEGER
          The leading dimension of the array C. LDC >= max(1,M).
[out]WORK
          WORK is COMPLEX*16 array, dimension (LDWORK,K)
[in]LDWORK
          LDWORK is INTEGER
          The leading dimension of the array WORK.
          If SIDE = 'L', LDWORK >= max(1,N);
          if SIDE = 'R', LDWORK >= max(1,M).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
June 2013
Further Details:
  The shape of the matrix V and the storage of the vectors which define
  the H(i) is best illustrated by the following example with n = 5 and
  k = 3. The elements equal to 1 are not stored; the corresponding
  array elements are modified but restored on exit. The rest of the
  array is not used.

  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':

               V = (  1       )                 V = (  1 v1 v1 v1 v1 )
                   ( v1  1    )                     (     1 v2 v2 v2 )
                   ( v1 v2  1 )                     (        1 v3 v3 )
                   ( v1 v2 v3 )
                   ( v1 v2 v3 )

  DIRECT = 'B' and STOREV = 'C':         DIRECT = 'B' and STOREV = 'R':

               V = ( v1 v2 v3 )                 V = ( v1 v1  1       )
                   ( v1 v2 v3 )                     ( v2 v2 v2  1    )
                   (  1 v2 v3 )                     ( v3 v3 v3 v3  1 )
                   (     1 v3 )
                   (        1 )

Definition at line 199 of file zlarfb.f.

199 *
200 * -- LAPACK auxiliary routine (version 3.7.0) --
201 * -- LAPACK is a software package provided by Univ. of Tennessee, --
202 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
203 * June 2013
204 *
205 * .. Scalar Arguments ..
206  CHARACTER DIRECT, SIDE, STOREV, TRANS
207  INTEGER K, LDC, LDT, LDV, LDWORK, M, N
208 * ..
209 * .. Array Arguments ..
210  COMPLEX*16 C( LDC, * ), T( LDT, * ), V( LDV, * ),
211  $ WORK( LDWORK, * )
212 * ..
213 *
214 * =====================================================================
215 *
216 * .. Parameters ..
217  COMPLEX*16 ONE
218  parameter( one = ( 1.0d+0, 0.0d+0 ) )
219 * ..
220 * .. Local Scalars ..
221  CHARACTER TRANST
222  INTEGER I, J
223 * ..
224 * .. External Functions ..
225  LOGICAL LSAME
226  EXTERNAL lsame
227 * ..
228 * .. External Subroutines ..
229  EXTERNAL zcopy, zgemm, zlacgv, ztrmm
230 * ..
231 * .. Intrinsic Functions ..
232  INTRINSIC dconjg
233 * ..
234 * .. Executable Statements ..
235 *
236 * Quick return if possible
237 *
238  IF( m.LE.0 .OR. n.LE.0 )
239  $ RETURN
240 *
241  IF( lsame( trans, 'N' ) ) THEN
242  transt = 'C'
243  ELSE
244  transt = 'N'
245  END IF
246 *
247  IF( lsame( storev, 'C' ) ) THEN
248 *
249  IF( lsame( direct, 'F' ) ) THEN
250 *
251 * Let V = ( V1 ) (first K rows)
252 * ( V2 )
253 * where V1 is unit lower triangular.
254 *
255  IF( lsame( side, 'L' ) ) THEN
256 *
257 * Form H * C or H**H * C where C = ( C1 )
258 * ( C2 )
259 *
260 * W := C**H * V = (C1**H * V1 + C2**H * V2) (stored in WORK)
261 *
262 * W := C1**H
263 *
264  DO 10 j = 1, k
265  CALL zcopy( n, c( j, 1 ), ldc, work( 1, j ), 1 )
266  CALL zlacgv( n, work( 1, j ), 1 )
267  10 CONTINUE
268 *
269 * W := W * V1
270 *
271  CALL ztrmm( 'Right', 'Lower', 'No transpose', 'Unit', n,
272  $ k, one, v, ldv, work, ldwork )
273  IF( m.GT.k ) THEN
274 *
275 * W := W + C2**H * V2
276 *
277  CALL zgemm( 'Conjugate transpose', 'No transpose', n,
278  $ k, m-k, one, c( k+1, 1 ), ldc,
279  $ v( k+1, 1 ), ldv, one, work, ldwork )
280  END IF
281 *
282 * W := W * T**H or W * T
283 *
284  CALL ztrmm( 'Right', 'Upper', transt, 'Non-unit', n, k,
285  $ one, t, ldt, work, ldwork )
286 *
287 * C := C - V * W**H
288 *
289  IF( m.GT.k ) THEN
290 *
291 * C2 := C2 - V2 * W**H
292 *
293  CALL zgemm( 'No transpose', 'Conjugate transpose',
294  $ m-k, n, k, -one, v( k+1, 1 ), ldv, work,
295  $ ldwork, one, c( k+1, 1 ), ldc )
296  END IF
297 *
298 * W := W * V1**H
299 *
300  CALL ztrmm( 'Right', 'Lower', 'Conjugate transpose',
301  $ 'Unit', n, k, one, v, ldv, work, ldwork )
302 *
303 * C1 := C1 - W**H
304 *
305  DO 30 j = 1, k
306  DO 20 i = 1, n
307  c( j, i ) = c( j, i ) - dconjg( work( i, j ) )
308  20 CONTINUE
309  30 CONTINUE
310 *
311  ELSE IF( lsame( side, 'R' ) ) THEN
312 *
313 * Form C * H or C * H**H where C = ( C1 C2 )
314 *
315 * W := C * V = (C1*V1 + C2*V2) (stored in WORK)
316 *
317 * W := C1
318 *
319  DO 40 j = 1, k
320  CALL zcopy( m, c( 1, j ), 1, work( 1, j ), 1 )
321  40 CONTINUE
322 *
323 * W := W * V1
324 *
325  CALL ztrmm( 'Right', 'Lower', 'No transpose', 'Unit', m,
326  $ k, one, v, ldv, work, ldwork )
327  IF( n.GT.k ) THEN
328 *
329 * W := W + C2 * V2
330 *
331  CALL zgemm( 'No transpose', 'No transpose', m, k, n-k,
332  $ one, c( 1, k+1 ), ldc, v( k+1, 1 ), ldv,
333  $ one, work, ldwork )
334  END IF
335 *
336 * W := W * T or W * T**H
337 *
338  CALL ztrmm( 'Right', 'Upper', trans, 'Non-unit', m, k,
339  $ one, t, ldt, work, ldwork )
340 *
341 * C := C - W * V**H
342 *
343  IF( n.GT.k ) THEN
344 *
345 * C2 := C2 - W * V2**H
346 *
347  CALL zgemm( 'No transpose', 'Conjugate transpose', m,
348  $ n-k, k, -one, work, ldwork, v( k+1, 1 ),
349  $ ldv, one, c( 1, k+1 ), ldc )
350  END IF
351 *
352 * W := W * V1**H
353 *
354  CALL ztrmm( 'Right', 'Lower', 'Conjugate transpose',
355  $ 'Unit', m, k, one, v, ldv, work, ldwork )
356 *
357 * C1 := C1 - W
358 *
359  DO 60 j = 1, k
360  DO 50 i = 1, m
361  c( i, j ) = c( i, j ) - work( i, j )
362  50 CONTINUE
363  60 CONTINUE
364  END IF
365 *
366  ELSE
367 *
368 * Let V = ( V1 )
369 * ( V2 ) (last K rows)
370 * where V2 is unit upper triangular.
371 *
372  IF( lsame( side, 'L' ) ) THEN
373 *
374 * Form H * C or H**H * C where C = ( C1 )
375 * ( C2 )
376 *
377 * W := C**H * V = (C1**H * V1 + C2**H * V2) (stored in WORK)
378 *
379 * W := C2**H
380 *
381  DO 70 j = 1, k
382  CALL zcopy( n, c( m-k+j, 1 ), ldc, work( 1, j ), 1 )
383  CALL zlacgv( n, work( 1, j ), 1 )
384  70 CONTINUE
385 *
386 * W := W * V2
387 *
388  CALL ztrmm( 'Right', 'Upper', 'No transpose', 'Unit', n,
389  $ k, one, v( m-k+1, 1 ), ldv, work, ldwork )
390  IF( m.GT.k ) THEN
391 *
392 * W := W + C1**H * V1
393 *
394  CALL zgemm( 'Conjugate transpose', 'No transpose', n,
395  $ k, m-k, one, c, ldc, v, ldv, one, work,
396  $ ldwork )
397  END IF
398 *
399 * W := W * T**H or W * T
400 *
401  CALL ztrmm( 'Right', 'Lower', transt, 'Non-unit', n, k,
402  $ one, t, ldt, work, ldwork )
403 *
404 * C := C - V * W**H
405 *
406  IF( m.GT.k ) THEN
407 *
408 * C1 := C1 - V1 * W**H
409 *
410  CALL zgemm( 'No transpose', 'Conjugate transpose',
411  $ m-k, n, k, -one, v, ldv, work, ldwork,
412  $ one, c, ldc )
413  END IF
414 *
415 * W := W * V2**H
416 *
417  CALL ztrmm( 'Right', 'Upper', 'Conjugate transpose',
418  $ 'Unit', n, k, one, v( m-k+1, 1 ), ldv, work,
419  $ ldwork )
420 *
421 * C2 := C2 - W**H
422 *
423  DO 90 j = 1, k
424  DO 80 i = 1, n
425  c( m-k+j, i ) = c( m-k+j, i ) -
426  $ dconjg( work( i, j ) )
427  80 CONTINUE
428  90 CONTINUE
429 *
430  ELSE IF( lsame( side, 'R' ) ) THEN
431 *
432 * Form C * H or C * H**H where C = ( C1 C2 )
433 *
434 * W := C * V = (C1*V1 + C2*V2) (stored in WORK)
435 *
436 * W := C2
437 *
438  DO 100 j = 1, k
439  CALL zcopy( m, c( 1, n-k+j ), 1, work( 1, j ), 1 )
440  100 CONTINUE
441 *
442 * W := W * V2
443 *
444  CALL ztrmm( 'Right', 'Upper', 'No transpose', 'Unit', m,
445  $ k, one, v( n-k+1, 1 ), ldv, work, ldwork )
446  IF( n.GT.k ) THEN
447 *
448 * W := W + C1 * V1
449 *
450  CALL zgemm( 'No transpose', 'No transpose', m, k, n-k,
451  $ one, c, ldc, v, ldv, one, work, ldwork )
452  END IF
453 *
454 * W := W * T or W * T**H
455 *
456  CALL ztrmm( 'Right', 'Lower', trans, 'Non-unit', m, k,
457  $ one, t, ldt, work, ldwork )
458 *
459 * C := C - W * V**H
460 *
461  IF( n.GT.k ) THEN
462 *
463 * C1 := C1 - W * V1**H
464 *
465  CALL zgemm( 'No transpose', 'Conjugate transpose', m,
466  $ n-k, k, -one, work, ldwork, v, ldv, one,
467  $ c, ldc )
468  END IF
469 *
470 * W := W * V2**H
471 *
472  CALL ztrmm( 'Right', 'Upper', 'Conjugate transpose',
473  $ 'Unit', m, k, one, v( n-k+1, 1 ), ldv, work,
474  $ ldwork )
475 *
476 * C2 := C2 - W
477 *
478  DO 120 j = 1, k
479  DO 110 i = 1, m
480  c( i, n-k+j ) = c( i, n-k+j ) - work( i, j )
481  110 CONTINUE
482  120 CONTINUE
483  END IF
484  END IF
485 *
486  ELSE IF( lsame( storev, 'R' ) ) THEN
487 *
488  IF( lsame( direct, 'F' ) ) THEN
489 *
490 * Let V = ( V1 V2 ) (V1: first K columns)
491 * where V1 is unit upper triangular.
492 *
493  IF( lsame( side, 'L' ) ) THEN
494 *
495 * Form H * C or H**H * C where C = ( C1 )
496 * ( C2 )
497 *
498 * W := C**H * V**H = (C1**H * V1**H + C2**H * V2**H) (stored in WORK)
499 *
500 * W := C1**H
501 *
502  DO 130 j = 1, k
503  CALL zcopy( n, c( j, 1 ), ldc, work( 1, j ), 1 )
504  CALL zlacgv( n, work( 1, j ), 1 )
505  130 CONTINUE
506 *
507 * W := W * V1**H
508 *
509  CALL ztrmm( 'Right', 'Upper', 'Conjugate transpose',
510  $ 'Unit', n, k, one, v, ldv, work, ldwork )
511  IF( m.GT.k ) THEN
512 *
513 * W := W + C2**H * V2**H
514 *
515  CALL zgemm( 'Conjugate transpose',
516  $ 'Conjugate transpose', n, k, m-k, one,
517  $ c( k+1, 1 ), ldc, v( 1, k+1 ), ldv, one,
518  $ work, ldwork )
519  END IF
520 *
521 * W := W * T**H or W * T
522 *
523  CALL ztrmm( 'Right', 'Upper', transt, 'Non-unit', n, k,
524  $ one, t, ldt, work, ldwork )
525 *
526 * C := C - V**H * W**H
527 *
528  IF( m.GT.k ) THEN
529 *
530 * C2 := C2 - V2**H * W**H
531 *
532  CALL zgemm( 'Conjugate transpose',
533  $ 'Conjugate transpose', m-k, n, k, -one,
534  $ v( 1, k+1 ), ldv, work, ldwork, one,
535  $ c( k+1, 1 ), ldc )
536  END IF
537 *
538 * W := W * V1
539 *
540  CALL ztrmm( 'Right', 'Upper', 'No transpose', 'Unit', n,
541  $ k, one, v, ldv, work, ldwork )
542 *
543 * C1 := C1 - W**H
544 *
545  DO 150 j = 1, k
546  DO 140 i = 1, n
547  c( j, i ) = c( j, i ) - dconjg( work( i, j ) )
548  140 CONTINUE
549  150 CONTINUE
550 *
551  ELSE IF( lsame( side, 'R' ) ) THEN
552 *
553 * Form C * H or C * H**H where C = ( C1 C2 )
554 *
555 * W := C * V**H = (C1*V1**H + C2*V2**H) (stored in WORK)
556 *
557 * W := C1
558 *
559  DO 160 j = 1, k
560  CALL zcopy( m, c( 1, j ), 1, work( 1, j ), 1 )
561  160 CONTINUE
562 *
563 * W := W * V1**H
564 *
565  CALL ztrmm( 'Right', 'Upper', 'Conjugate transpose',
566  $ 'Unit', m, k, one, v, ldv, work, ldwork )
567  IF( n.GT.k ) THEN
568 *
569 * W := W + C2 * V2**H
570 *
571  CALL zgemm( 'No transpose', 'Conjugate transpose', m,
572  $ k, n-k, one, c( 1, k+1 ), ldc,
573  $ v( 1, k+1 ), ldv, one, work, ldwork )
574  END IF
575 *
576 * W := W * T or W * T**H
577 *
578  CALL ztrmm( 'Right', 'Upper', trans, 'Non-unit', m, k,
579  $ one, t, ldt, work, ldwork )
580 *
581 * C := C - W * V
582 *
583  IF( n.GT.k ) THEN
584 *
585 * C2 := C2 - W * V2
586 *
587  CALL zgemm( 'No transpose', 'No transpose', m, n-k, k,
588  $ -one, work, ldwork, v( 1, k+1 ), ldv, one,
589  $ c( 1, k+1 ), ldc )
590  END IF
591 *
592 * W := W * V1
593 *
594  CALL ztrmm( 'Right', 'Upper', 'No transpose', 'Unit', m,
595  $ k, one, v, ldv, work, ldwork )
596 *
597 * C1 := C1 - W
598 *
599  DO 180 j = 1, k
600  DO 170 i = 1, m
601  c( i, j ) = c( i, j ) - work( i, j )
602  170 CONTINUE
603  180 CONTINUE
604 *
605  END IF
606 *
607  ELSE
608 *
609 * Let V = ( V1 V2 ) (V2: last K columns)
610 * where V2 is unit lower triangular.
611 *
612  IF( lsame( side, 'L' ) ) THEN
613 *
614 * Form H * C or H**H * C where C = ( C1 )
615 * ( C2 )
616 *
617 * W := C**H * V**H = (C1**H * V1**H + C2**H * V2**H) (stored in WORK)
618 *
619 * W := C2**H
620 *
621  DO 190 j = 1, k
622  CALL zcopy( n, c( m-k+j, 1 ), ldc, work( 1, j ), 1 )
623  CALL zlacgv( n, work( 1, j ), 1 )
624  190 CONTINUE
625 *
626 * W := W * V2**H
627 *
628  CALL ztrmm( 'Right', 'Lower', 'Conjugate transpose',
629  $ 'Unit', n, k, one, v( 1, m-k+1 ), ldv, work,
630  $ ldwork )
631  IF( m.GT.k ) THEN
632 *
633 * W := W + C1**H * V1**H
634 *
635  CALL zgemm( 'Conjugate transpose',
636  $ 'Conjugate transpose', n, k, m-k, one, c,
637  $ ldc, v, ldv, one, work, ldwork )
638  END IF
639 *
640 * W := W * T**H or W * T
641 *
642  CALL ztrmm( 'Right', 'Lower', transt, 'Non-unit', n, k,
643  $ one, t, ldt, work, ldwork )
644 *
645 * C := C - V**H * W**H
646 *
647  IF( m.GT.k ) THEN
648 *
649 * C1 := C1 - V1**H * W**H
650 *
651  CALL zgemm( 'Conjugate transpose',
652  $ 'Conjugate transpose', m-k, n, k, -one, v,
653  $ ldv, work, ldwork, one, c, ldc )
654  END IF
655 *
656 * W := W * V2
657 *
658  CALL ztrmm( 'Right', 'Lower', 'No transpose', 'Unit', n,
659  $ k, one, v( 1, m-k+1 ), ldv, work, ldwork )
660 *
661 * C2 := C2 - W**H
662 *
663  DO 210 j = 1, k
664  DO 200 i = 1, n
665  c( m-k+j, i ) = c( m-k+j, i ) -
666  $ dconjg( work( i, j ) )
667  200 CONTINUE
668  210 CONTINUE
669 *
670  ELSE IF( lsame( side, 'R' ) ) THEN
671 *
672 * Form C * H or C * H**H where C = ( C1 C2 )
673 *
674 * W := C * V**H = (C1*V1**H + C2*V2**H) (stored in WORK)
675 *
676 * W := C2
677 *
678  DO 220 j = 1, k
679  CALL zcopy( m, c( 1, n-k+j ), 1, work( 1, j ), 1 )
680  220 CONTINUE
681 *
682 * W := W * V2**H
683 *
684  CALL ztrmm( 'Right', 'Lower', 'Conjugate transpose',
685  $ 'Unit', m, k, one, v( 1, n-k+1 ), ldv, work,
686  $ ldwork )
687  IF( n.GT.k ) THEN
688 *
689 * W := W + C1 * V1**H
690 *
691  CALL zgemm( 'No transpose', 'Conjugate transpose', m,
692  $ k, n-k, one, c, ldc, v, ldv, one, work,
693  $ ldwork )
694  END IF
695 *
696 * W := W * T or W * T**H
697 *
698  CALL ztrmm( 'Right', 'Lower', trans, 'Non-unit', m, k,
699  $ one, t, ldt, work, ldwork )
700 *
701 * C := C - W * V
702 *
703  IF( n.GT.k ) THEN
704 *
705 * C1 := C1 - W * V1
706 *
707  CALL zgemm( 'No transpose', 'No transpose', m, n-k, k,
708  $ -one, work, ldwork, v, ldv, one, c, ldc )
709  END IF
710 *
711 * W := W * V2
712 *
713  CALL ztrmm( 'Right', 'Lower', 'No transpose', 'Unit', m,
714  $ k, one, v( 1, n-k+1 ), ldv, work, ldwork )
715 *
716 * C1 := C1 - W
717 *
718  DO 240 j = 1, k
719  DO 230 i = 1, m
720  c( i, n-k+j ) = c( i, n-k+j ) - work( i, j )
721  230 CONTINUE
722  240 CONTINUE
723 *
724  END IF
725 *
726  END IF
727  END IF
728 *
729  RETURN
730 *
731 * End of ZLARFB
732 *
Here is the call graph for this function:
Here is the caller graph for this function:
zlacgv
subroutine zlacgv(N, X, INCX)
ZLACGV conjugates a complex vector.
Definition: zlacgv.f:76
zcopy
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:83
zgemm
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:189
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
ztrmm
subroutine ztrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
ZTRMM
Definition: ztrmm.f:179