LAPACK  3.9.0
LAPACK: Linear Algebra PACKage
dlarfb.f
Go to the documentation of this file.
1 *> \brief \b DLARFB applies a block reflector or its transpose to a general rectangular matrix.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download DLARFB + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlarfb.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlarfb.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlarfb.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE DLARFB( SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV,
22 * T, LDT, C, LDC, WORK, LDWORK )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER DIRECT, SIDE, STOREV, TRANS
26 * INTEGER K, LDC, LDT, LDV, LDWORK, M, N
27 * ..
28 * .. Array Arguments ..
29 * DOUBLE PRECISION C( LDC, * ), T( LDT, * ), V( LDV, * ),
30 * $ WORK( LDWORK, * )
31 * ..
32 *
33 *
34 *> \par Purpose:
35 * =============
36 *>
37 *> \verbatim
38 *>
39 *> DLARFB applies a real block reflector H or its transpose H**T to a
40 *> real m by n matrix C, from either the left or the right.
41 *> \endverbatim
42 *
43 * Arguments:
44 * ==========
45 *
46 *> \param[in] SIDE
47 *> \verbatim
48 *> SIDE is CHARACTER*1
49 *> = 'L': apply H or H**T from the Left
50 *> = 'R': apply H or H**T from the Right
51 *> \endverbatim
52 *>
53 *> \param[in] TRANS
54 *> \verbatim
55 *> TRANS is CHARACTER*1
56 *> = 'N': apply H (No transpose)
57 *> = 'T': apply H**T (Transpose)
58 *> \endverbatim
59 *>
60 *> \param[in] DIRECT
61 *> \verbatim
62 *> DIRECT is CHARACTER*1
63 *> Indicates how H is formed from a product of elementary
64 *> reflectors
65 *> = 'F': H = H(1) H(2) . . . H(k) (Forward)
66 *> = 'B': H = H(k) . . . H(2) H(1) (Backward)
67 *> \endverbatim
68 *>
69 *> \param[in] STOREV
70 *> \verbatim
71 *> STOREV is CHARACTER*1
72 *> Indicates how the vectors which define the elementary
73 *> reflectors are stored:
74 *> = 'C': Columnwise
75 *> = 'R': Rowwise
76 *> \endverbatim
77 *>
78 *> \param[in] M
79 *> \verbatim
80 *> M is INTEGER
81 *> The number of rows of the matrix C.
82 *> \endverbatim
83 *>
84 *> \param[in] N
85 *> \verbatim
86 *> N is INTEGER
87 *> The number of columns of the matrix C.
88 *> \endverbatim
89 *>
90 *> \param[in] K
91 *> \verbatim
92 *> K is INTEGER
93 *> The order of the matrix T (= the number of elementary
94 *> reflectors whose product defines the block reflector).
95 *> If SIDE = 'L', M >= K >= 0;
96 *> if SIDE = 'R', N >= K >= 0.
97 *> \endverbatim
98 *>
99 *> \param[in] V
100 *> \verbatim
101 *> V is DOUBLE PRECISION array, dimension
102 *> (LDV,K) if STOREV = 'C'
103 *> (LDV,M) if STOREV = 'R' and SIDE = 'L'
104 *> (LDV,N) if STOREV = 'R' and SIDE = 'R'
105 *> The matrix V. See Further Details.
106 *> \endverbatim
107 *>
108 *> \param[in] LDV
109 *> \verbatim
110 *> LDV is INTEGER
111 *> The leading dimension of the array V.
112 *> If STOREV = 'C' and SIDE = 'L', LDV >= max(1,M);
113 *> if STOREV = 'C' and SIDE = 'R', LDV >= max(1,N);
114 *> if STOREV = 'R', LDV >= K.
115 *> \endverbatim
116 *>
117 *> \param[in] T
118 *> \verbatim
119 *> T is DOUBLE PRECISION array, dimension (LDT,K)
120 *> The triangular k by k matrix T in the representation of the
121 *> block reflector.
122 *> \endverbatim
123 *>
124 *> \param[in] LDT
125 *> \verbatim
126 *> LDT is INTEGER
127 *> The leading dimension of the array T. LDT >= K.
128 *> \endverbatim
129 *>
130 *> \param[in,out] C
131 *> \verbatim
132 *> C is DOUBLE PRECISION array, dimension (LDC,N)
133 *> On entry, the m by n matrix C.
134 *> On exit, C is overwritten by H*C or H**T*C or C*H or C*H**T.
135 *> \endverbatim
136 *>
137 *> \param[in] LDC
138 *> \verbatim
139 *> LDC is INTEGER
140 *> The leading dimension of the array C. LDC >= max(1,M).
141 *> \endverbatim
142 *>
143 *> \param[out] WORK
144 *> \verbatim
145 *> WORK is DOUBLE PRECISION array, dimension (LDWORK,K)
146 *> \endverbatim
147 *>
148 *> \param[in] LDWORK
149 *> \verbatim
150 *> LDWORK is INTEGER
151 *> The leading dimension of the array WORK.
152 *> If SIDE = 'L', LDWORK >= max(1,N);
153 *> if SIDE = 'R', LDWORK >= max(1,M).
154 *> \endverbatim
155 *
156 * Authors:
157 * ========
158 *
159 *> \author Univ. of Tennessee
160 *> \author Univ. of California Berkeley
161 *> \author Univ. of Colorado Denver
162 *> \author NAG Ltd.
163 *
164 *> \date June 2013
165 *
166 *> \ingroup doubleOTHERauxiliary
167 *
168 *> \par Further Details:
169 * =====================
170 *>
171 *> \verbatim
172 *>
173 *> The shape of the matrix V and the storage of the vectors which define
174 *> the H(i) is best illustrated by the following example with n = 5 and
175 *> k = 3. The elements equal to 1 are not stored; the corresponding
176 *> array elements are modified but restored on exit. The rest of the
177 *> array is not used.
178 *>
179 *> DIRECT = 'F' and STOREV = 'C': DIRECT = 'F' and STOREV = 'R':
180 *>
181 *> V = ( 1 ) V = ( 1 v1 v1 v1 v1 )
182 *> ( v1 1 ) ( 1 v2 v2 v2 )
183 *> ( v1 v2 1 ) ( 1 v3 v3 )
184 *> ( v1 v2 v3 )
185 *> ( v1 v2 v3 )
186 *>
187 *> DIRECT = 'B' and STOREV = 'C': DIRECT = 'B' and STOREV = 'R':
188 *>
189 *> V = ( v1 v2 v3 ) V = ( v1 v1 1 )
190 *> ( v1 v2 v3 ) ( v2 v2 v2 1 )
191 *> ( 1 v2 v3 ) ( v3 v3 v3 v3 1 )
192 *> ( 1 v3 )
193 *> ( 1 )
194 *> \endverbatim
195 *>
196 * =====================================================================
197  SUBROUTINE dlarfb( SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV,
198  $ T, LDT, C, LDC, WORK, LDWORK )
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  DOUBLE PRECISION C( LDC, * ), T( LDT, * ), V( LDV, * ),
211  $ work( ldwork, * )
212 * ..
213 *
214 * =====================================================================
215 *
216 * .. Parameters ..
217  DOUBLE PRECISION ONE
218  parameter( one = 1.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 dcopy, dgemm, dtrmm
230 * ..
231 * .. Executable Statements ..
232 *
233 * Quick return if possible
234 *
235  IF( m.LE.0 .OR. n.LE.0 )
236  $ RETURN
237 *
238  IF( lsame( trans, 'N' ) ) THEN
239  transt = 'T'
240  ELSE
241  transt = 'N'
242  END IF
243 *
244  IF( lsame( storev, 'C' ) ) THEN
245 *
246  IF( lsame( direct, 'F' ) ) THEN
247 *
248 * Let V = ( V1 ) (first K rows)
249 * ( V2 )
250 * where V1 is unit lower triangular.
251 *
252  IF( lsame( side, 'L' ) ) THEN
253 *
254 * Form H * C or H**T * C where C = ( C1 )
255 * ( C2 )
256 *
257 * W := C**T * V = (C1**T * V1 + C2**T * V2) (stored in WORK)
258 *
259 * W := C1**T
260 *
261  DO 10 j = 1, k
262  CALL dcopy( n, c( j, 1 ), ldc, work( 1, j ), 1 )
263  10 CONTINUE
264 *
265 * W := W * V1
266 *
267  CALL dtrmm( 'Right', 'Lower', 'No transpose', 'Unit', n,
268  $ k, one, v, ldv, work, ldwork )
269  IF( m.GT.k ) THEN
270 *
271 * W := W + C2**T * V2
272 *
273  CALL dgemm( 'Transpose', 'No transpose', n, k, m-k,
274  $ one, c( k+1, 1 ), ldc, v( k+1, 1 ), ldv,
275  $ one, work, ldwork )
276  END IF
277 *
278 * W := W * T**T or W * T
279 *
280  CALL dtrmm( 'Right', 'Upper', transt, 'Non-unit', n, k,
281  $ one, t, ldt, work, ldwork )
282 *
283 * C := C - V * W**T
284 *
285  IF( m.GT.k ) THEN
286 *
287 * C2 := C2 - V2 * W**T
288 *
289  CALL dgemm( 'No transpose', 'Transpose', m-k, n, k,
290  $ -one, v( k+1, 1 ), ldv, work, ldwork, one,
291  $ c( k+1, 1 ), ldc )
292  END IF
293 *
294 * W := W * V1**T
295 *
296  CALL dtrmm( 'Right', 'Lower', 'Transpose', 'Unit', n, k,
297  $ one, v, ldv, work, ldwork )
298 *
299 * C1 := C1 - W**T
300 *
301  DO 30 j = 1, k
302  DO 20 i = 1, n
303  c( j, i ) = c( j, i ) - work( i, j )
304  20 CONTINUE
305  30 CONTINUE
306 *
307  ELSE IF( lsame( side, 'R' ) ) THEN
308 *
309 * Form C * H or C * H**T where C = ( C1 C2 )
310 *
311 * W := C * V = (C1*V1 + C2*V2) (stored in WORK)
312 *
313 * W := C1
314 *
315  DO 40 j = 1, k
316  CALL dcopy( m, c( 1, j ), 1, work( 1, j ), 1 )
317  40 CONTINUE
318 *
319 * W := W * V1
320 *
321  CALL dtrmm( 'Right', 'Lower', 'No transpose', 'Unit', m,
322  $ k, one, v, ldv, work, ldwork )
323  IF( n.GT.k ) THEN
324 *
325 * W := W + C2 * V2
326 *
327  CALL dgemm( 'No transpose', 'No transpose', m, k, n-k,
328  $ one, c( 1, k+1 ), ldc, v( k+1, 1 ), ldv,
329  $ one, work, ldwork )
330  END IF
331 *
332 * W := W * T or W * T**T
333 *
334  CALL dtrmm( 'Right', 'Upper', trans, 'Non-unit', m, k,
335  $ one, t, ldt, work, ldwork )
336 *
337 * C := C - W * V**T
338 *
339  IF( n.GT.k ) THEN
340 *
341 * C2 := C2 - W * V2**T
342 *
343  CALL dgemm( 'No transpose', 'Transpose', m, n-k, k,
344  $ -one, work, ldwork, v( k+1, 1 ), ldv, one,
345  $ c( 1, k+1 ), ldc )
346  END IF
347 *
348 * W := W * V1**T
349 *
350  CALL dtrmm( 'Right', 'Lower', 'Transpose', 'Unit', m, k,
351  $ one, v, ldv, work, ldwork )
352 *
353 * C1 := C1 - W
354 *
355  DO 60 j = 1, k
356  DO 50 i = 1, m
357  c( i, j ) = c( i, j ) - work( i, j )
358  50 CONTINUE
359  60 CONTINUE
360  END IF
361 *
362  ELSE
363 *
364 * Let V = ( V1 )
365 * ( V2 ) (last K rows)
366 * where V2 is unit upper triangular.
367 *
368  IF( lsame( side, 'L' ) ) THEN
369 *
370 * Form H * C or H**T * C where C = ( C1 )
371 * ( C2 )
372 *
373 * W := C**T * V = (C1**T * V1 + C2**T * V2) (stored in WORK)
374 *
375 * W := C2**T
376 *
377  DO 70 j = 1, k
378  CALL dcopy( n, c( m-k+j, 1 ), ldc, work( 1, j ), 1 )
379  70 CONTINUE
380 *
381 * W := W * V2
382 *
383  CALL dtrmm( 'Right', 'Upper', 'No transpose', 'Unit', n,
384  $ k, one, v( m-k+1, 1 ), ldv, work, ldwork )
385  IF( m.GT.k ) THEN
386 *
387 * W := W + C1**T * V1
388 *
389  CALL dgemm( 'Transpose', 'No transpose', n, k, m-k,
390  $ one, c, ldc, v, ldv, one, work, ldwork )
391  END IF
392 *
393 * W := W * T**T or W * T
394 *
395  CALL dtrmm( 'Right', 'Lower', transt, 'Non-unit', n, k,
396  $ one, t, ldt, work, ldwork )
397 *
398 * C := C - V * W**T
399 *
400  IF( m.GT.k ) THEN
401 *
402 * C1 := C1 - V1 * W**T
403 *
404  CALL dgemm( 'No transpose', 'Transpose', m-k, n, k,
405  $ -one, v, ldv, work, ldwork, one, c, ldc )
406  END IF
407 *
408 * W := W * V2**T
409 *
410  CALL dtrmm( 'Right', 'Upper', 'Transpose', 'Unit', n, k,
411  $ one, v( m-k+1, 1 ), ldv, work, ldwork )
412 *
413 * C2 := C2 - W**T
414 *
415  DO 90 j = 1, k
416  DO 80 i = 1, n
417  c( m-k+j, i ) = c( m-k+j, i ) - work( i, j )
418  80 CONTINUE
419  90 CONTINUE
420 *
421  ELSE IF( lsame( side, 'R' ) ) THEN
422 *
423 * Form C * H or C * H**T where C = ( C1 C2 )
424 *
425 * W := C * V = (C1*V1 + C2*V2) (stored in WORK)
426 *
427 * W := C2
428 *
429  DO 100 j = 1, k
430  CALL dcopy( m, c( 1, n-k+j ), 1, work( 1, j ), 1 )
431  100 CONTINUE
432 *
433 * W := W * V2
434 *
435  CALL dtrmm( 'Right', 'Upper', 'No transpose', 'Unit', m,
436  $ k, one, v( n-k+1, 1 ), ldv, work, ldwork )
437  IF( n.GT.k ) THEN
438 *
439 * W := W + C1 * V1
440 *
441  CALL dgemm( 'No transpose', 'No transpose', m, k, n-k,
442  $ one, c, ldc, v, ldv, one, work, ldwork )
443  END IF
444 *
445 * W := W * T or W * T**T
446 *
447  CALL dtrmm( 'Right', 'Lower', trans, 'Non-unit', m, k,
448  $ one, t, ldt, work, ldwork )
449 *
450 * C := C - W * V**T
451 *
452  IF( n.GT.k ) THEN
453 *
454 * C1 := C1 - W * V1**T
455 *
456  CALL dgemm( 'No transpose', 'Transpose', m, n-k, k,
457  $ -one, work, ldwork, v, ldv, one, c, ldc )
458  END IF
459 *
460 * W := W * V2**T
461 *
462  CALL dtrmm( 'Right', 'Upper', 'Transpose', 'Unit', m, k,
463  $ one, v( n-k+1, 1 ), ldv, work, ldwork )
464 *
465 * C2 := C2 - W
466 *
467  DO 120 j = 1, k
468  DO 110 i = 1, m
469  c( i, n-k+j ) = c( i, n-k+j ) - work( i, j )
470  110 CONTINUE
471  120 CONTINUE
472  END IF
473  END IF
474 *
475  ELSE IF( lsame( storev, 'R' ) ) THEN
476 *
477  IF( lsame( direct, 'F' ) ) THEN
478 *
479 * Let V = ( V1 V2 ) (V1: first K columns)
480 * where V1 is unit upper triangular.
481 *
482  IF( lsame( side, 'L' ) ) THEN
483 *
484 * Form H * C or H**T * C where C = ( C1 )
485 * ( C2 )
486 *
487 * W := C**T * V**T = (C1**T * V1**T + C2**T * V2**T) (stored in WORK)
488 *
489 * W := C1**T
490 *
491  DO 130 j = 1, k
492  CALL dcopy( n, c( j, 1 ), ldc, work( 1, j ), 1 )
493  130 CONTINUE
494 *
495 * W := W * V1**T
496 *
497  CALL dtrmm( 'Right', 'Upper', 'Transpose', 'Unit', n, k,
498  $ one, v, ldv, work, ldwork )
499  IF( m.GT.k ) THEN
500 *
501 * W := W + C2**T * V2**T
502 *
503  CALL dgemm( 'Transpose', 'Transpose', n, k, m-k, one,
504  $ c( k+1, 1 ), ldc, v( 1, k+1 ), ldv, one,
505  $ work, ldwork )
506  END IF
507 *
508 * W := W * T**T or W * T
509 *
510  CALL dtrmm( 'Right', 'Upper', transt, 'Non-unit', n, k,
511  $ one, t, ldt, work, ldwork )
512 *
513 * C := C - V**T * W**T
514 *
515  IF( m.GT.k ) THEN
516 *
517 * C2 := C2 - V2**T * W**T
518 *
519  CALL dgemm( 'Transpose', 'Transpose', m-k, n, k, -one,
520  $ v( 1, k+1 ), ldv, work, ldwork, one,
521  $ c( k+1, 1 ), ldc )
522  END IF
523 *
524 * W := W * V1
525 *
526  CALL dtrmm( 'Right', 'Upper', 'No transpose', 'Unit', n,
527  $ k, one, v, ldv, work, ldwork )
528 *
529 * C1 := C1 - W**T
530 *
531  DO 150 j = 1, k
532  DO 140 i = 1, n
533  c( j, i ) = c( j, i ) - work( i, j )
534  140 CONTINUE
535  150 CONTINUE
536 *
537  ELSE IF( lsame( side, 'R' ) ) THEN
538 *
539 * Form C * H or C * H**T where C = ( C1 C2 )
540 *
541 * W := C * V**T = (C1*V1**T + C2*V2**T) (stored in WORK)
542 *
543 * W := C1
544 *
545  DO 160 j = 1, k
546  CALL dcopy( m, c( 1, j ), 1, work( 1, j ), 1 )
547  160 CONTINUE
548 *
549 * W := W * V1**T
550 *
551  CALL dtrmm( 'Right', 'Upper', 'Transpose', 'Unit', m, k,
552  $ one, v, ldv, work, ldwork )
553  IF( n.GT.k ) THEN
554 *
555 * W := W + C2 * V2**T
556 *
557  CALL dgemm( 'No transpose', 'Transpose', m, k, n-k,
558  $ one, c( 1, k+1 ), ldc, v( 1, k+1 ), ldv,
559  $ one, work, ldwork )
560  END IF
561 *
562 * W := W * T or W * T**T
563 *
564  CALL dtrmm( 'Right', 'Upper', trans, 'Non-unit', m, k,
565  $ one, t, ldt, work, ldwork )
566 *
567 * C := C - W * V
568 *
569  IF( n.GT.k ) THEN
570 *
571 * C2 := C2 - W * V2
572 *
573  CALL dgemm( 'No transpose', 'No transpose', m, n-k, k,
574  $ -one, work, ldwork, v( 1, k+1 ), ldv, one,
575  $ c( 1, k+1 ), ldc )
576  END IF
577 *
578 * W := W * V1
579 *
580  CALL dtrmm( 'Right', 'Upper', 'No transpose', 'Unit', m,
581  $ k, one, v, ldv, work, ldwork )
582 *
583 * C1 := C1 - W
584 *
585  DO 180 j = 1, k
586  DO 170 i = 1, m
587  c( i, j ) = c( i, j ) - work( i, j )
588  170 CONTINUE
589  180 CONTINUE
590 *
591  END IF
592 *
593  ELSE
594 *
595 * Let V = ( V1 V2 ) (V2: last K columns)
596 * where V2 is unit lower triangular.
597 *
598  IF( lsame( side, 'L' ) ) THEN
599 *
600 * Form H * C or H**T * C where C = ( C1 )
601 * ( C2 )
602 *
603 * W := C**T * V**T = (C1**T * V1**T + C2**T * V2**T) (stored in WORK)
604 *
605 * W := C2**T
606 *
607  DO 190 j = 1, k
608  CALL dcopy( n, c( m-k+j, 1 ), ldc, work( 1, j ), 1 )
609  190 CONTINUE
610 *
611 * W := W * V2**T
612 *
613  CALL dtrmm( 'Right', 'Lower', 'Transpose', 'Unit', n, k,
614  $ one, v( 1, m-k+1 ), ldv, work, ldwork )
615  IF( m.GT.k ) THEN
616 *
617 * W := W + C1**T * V1**T
618 *
619  CALL dgemm( 'Transpose', 'Transpose', n, k, m-k, one,
620  $ c, ldc, v, ldv, one, work, ldwork )
621  END IF
622 *
623 * W := W * T**T or W * T
624 *
625  CALL dtrmm( 'Right', 'Lower', transt, 'Non-unit', n, k,
626  $ one, t, ldt, work, ldwork )
627 *
628 * C := C - V**T * W**T
629 *
630  IF( m.GT.k ) THEN
631 *
632 * C1 := C1 - V1**T * W**T
633 *
634  CALL dgemm( 'Transpose', 'Transpose', m-k, n, k, -one,
635  $ v, ldv, work, ldwork, one, c, ldc )
636  END IF
637 *
638 * W := W * V2
639 *
640  CALL dtrmm( 'Right', 'Lower', 'No transpose', 'Unit', n,
641  $ k, one, v( 1, m-k+1 ), ldv, work, ldwork )
642 *
643 * C2 := C2 - W**T
644 *
645  DO 210 j = 1, k
646  DO 200 i = 1, n
647  c( m-k+j, i ) = c( m-k+j, i ) - work( i, j )
648  200 CONTINUE
649  210 CONTINUE
650 *
651  ELSE IF( lsame( side, 'R' ) ) THEN
652 *
653 * Form C * H or C * H' where C = ( C1 C2 )
654 *
655 * W := C * V**T = (C1*V1**T + C2*V2**T) (stored in WORK)
656 *
657 * W := C2
658 *
659  DO 220 j = 1, k
660  CALL dcopy( m, c( 1, n-k+j ), 1, work( 1, j ), 1 )
661  220 CONTINUE
662 *
663 * W := W * V2**T
664 *
665  CALL dtrmm( 'Right', 'Lower', 'Transpose', 'Unit', m, k,
666  $ one, v( 1, n-k+1 ), ldv, work, ldwork )
667  IF( n.GT.k ) THEN
668 *
669 * W := W + C1 * V1**T
670 *
671  CALL dgemm( 'No transpose', 'Transpose', m, k, n-k,
672  $ one, c, ldc, v, ldv, one, work, ldwork )
673  END IF
674 *
675 * W := W * T or W * T**T
676 *
677  CALL dtrmm( 'Right', 'Lower', trans, 'Non-unit', m, k,
678  $ one, t, ldt, work, ldwork )
679 *
680 * C := C - W * V
681 *
682  IF( n.GT.k ) THEN
683 *
684 * C1 := C1 - W * V1
685 *
686  CALL dgemm( 'No transpose', 'No transpose', m, n-k, k,
687  $ -one, work, ldwork, v, ldv, one, c, ldc )
688  END IF
689 *
690 * W := W * V2
691 *
692  CALL dtrmm( 'Right', 'Lower', 'No transpose', 'Unit', m,
693  $ k, one, v( 1, n-k+1 ), ldv, work, ldwork )
694 *
695 * C1 := C1 - W
696 *
697  DO 240 j = 1, k
698  DO 230 i = 1, m
699  c( i, n-k+j ) = c( i, n-k+j ) - work( i, j )
700  230 CONTINUE
701  240 CONTINUE
702 *
703  END IF
704 *
705  END IF
706  END IF
707 *
708  RETURN
709 *
710 * End of DLARFB
711 *
712  END
dtrmm
subroutine dtrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRMM
Definition: dtrmm.f:179
dlarfb
subroutine dlarfb(SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV, T, LDT, C, LDC, WORK, LDWORK)
DLARFB applies a block reflector or its transpose to a general rectangular matrix.
Definition: dlarfb.f:199
dcopy
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:84
dgemm
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
Definition: dgemm.f:189