LAPACK  3.9.0
LAPACK: Linear Algebra PACKage
clarfb.f
Go to the documentation of this file.
1 *> \brief \b CLARFB applies a block reflector or its conjugate-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 CLARFB + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/clarfb.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/clarfb.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/clarfb.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE CLARFB( 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 * COMPLEX C( LDC, * ), T( LDT, * ), V( LDV, * ),
30 * $ WORK( LDWORK, * )
31 * ..
32 *
33 *
34 *> \par Purpose:
35 * =============
36 *>
37 *> \verbatim
38 *>
39 *> CLARFB applies a complex block reflector H or its transpose H**H to a
40 *> complex 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**H from the Left
50 *> = 'R': apply H or H**H from the Right
51 *> \endverbatim
52 *>
53 *> \param[in] TRANS
54 *> \verbatim
55 *> TRANS is CHARACTER*1
56 *> = 'N': apply H (No transpose)
57 *> = 'C': apply H**H (Conjugate 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 COMPLEX 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 COMPLEX 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 COMPLEX array, dimension (LDC,N)
133 *> On entry, the M-by-N matrix C.
134 *> On exit, C is overwritten by H*C or H**H*C or C*H or C*H**H.
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 COMPLEX 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 complexOTHERauxiliary
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 clarfb( 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  COMPLEX C( LDC, * ), T( LDT, * ), V( LDV, * ),
211  $ work( ldwork, * )
212 * ..
213 *
214 * =====================================================================
215 *
216 * .. Parameters ..
217  COMPLEX ONE
218  parameter( one = ( 1.0e+0, 0.0e+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 ccopy, cgemm, clacgv, ctrmm
230 * ..
231 * .. Intrinsic Functions ..
232  INTRINSIC conjg
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 ccopy( n, c( j, 1 ), ldc, work( 1, j ), 1 )
266  CALL clacgv( n, work( 1, j ), 1 )
267  10 CONTINUE
268 *
269 * W := W * V1
270 *
271  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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 ) - conjg( 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 ccopy( m, c( 1, j ), 1, work( 1, j ), 1 )
321  40 CONTINUE
322 *
323 * W := W * V1
324 *
325  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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 ccopy( n, c( m-k+j, 1 ), ldc, work( 1, j ), 1 )
383  CALL clacgv( n, work( 1, j ), 1 )
384  70 CONTINUE
385 *
386 * W := W * V2
387 *
388  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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  $ conjg( 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 ccopy( m, c( 1, n-k+j ), 1, work( 1, j ), 1 )
440  100 CONTINUE
441 *
442 * W := W * V2
443 *
444  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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 ccopy( n, c( j, 1 ), ldc, work( 1, j ), 1 )
504  CALL clacgv( n, work( 1, j ), 1 )
505  130 CONTINUE
506 *
507 * W := W * V1**H
508 *
509  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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 ) - conjg( 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 ccopy( m, c( 1, j ), 1, work( 1, j ), 1 )
561  160 CONTINUE
562 *
563 * W := W * V1**H
564 *
565  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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 ccopy( n, c( m-k+j, 1 ), ldc, work( 1, j ), 1 )
623  CALL clacgv( n, work( 1, j ), 1 )
624  190 CONTINUE
625 *
626 * W := W * V2**H
627 *
628  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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  $ conjg( 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 ccopy( m, c( 1, n-k+j ), 1, work( 1, j ), 1 )
680  220 CONTINUE
681 *
682 * W := W * V2**H
683 *
684  CALL ctrmm( '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 cgemm( '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 ctrmm( '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 cgemm( '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 ctrmm( '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 CLARFB
732 *
733  END
cgemm
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
Definition: cgemm.f:189
clacgv
subroutine clacgv(N, X, INCX)
CLACGV conjugates a complex vector.
Definition: clacgv.f:76
clarfb
subroutine clarfb(SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV, T, LDT, C, LDC, WORK, LDWORK)
CLARFB applies a block reflector or its conjugate-transpose to a general rectangular matrix.
Definition: clarfb.f:199
ctrmm
subroutine ctrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
CTRMM
Definition: ctrmm.f:179
ccopy
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:83