 LAPACK  3.10.0 LAPACK: Linear Algebra PACKage

◆ slarfb()

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

SLARFB applies a block reflector or its transpose to a general rectangular matrix.

Purpose:
SLARFB applies a real block reflector H or its transpose H**T to a
real m by n matrix C, from either the left or the 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': 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 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 matrix V. 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] C C is REAL array, dimension (LDC,N) On entry, the m by n matrix C. On exit, C is overwritten by H*C or H**T*C or C*H or C*H**T. [in] LDC LDC is INTEGER The leading dimension of the array C. LDC >= max(1,M). [out] WORK WORK is REAL 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).
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 195 of file slarfb.f.

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