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