 LAPACK 3.11.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ zunbdb()

 subroutine zunbdb ( character TRANS, character SIGNS, integer M, integer P, integer Q, complex*16, dimension( ldx11, * ) X11, integer LDX11, complex*16, dimension( ldx12, * ) X12, integer LDX12, complex*16, dimension( ldx21, * ) X21, integer LDX21, complex*16, dimension( ldx22, * ) X22, integer LDX22, double precision, dimension( * ) THETA, double precision, dimension( * ) PHI, complex*16, dimension( * ) TAUP1, complex*16, dimension( * ) TAUP2, complex*16, dimension( * ) TAUQ1, complex*16, dimension( * ) TAUQ2, complex*16, dimension( * ) WORK, integer LWORK, integer INFO )

ZUNBDB

Purpose:
``` ZUNBDB simultaneously bidiagonalizes the blocks of an M-by-M
partitioned unitary matrix X:

[ B11 | B12 0  0 ]
[ X11 | X12 ]   [ P1 |    ] [  0  |  0 -I  0 ] [ Q1 |    ]**H
X = [-----------] = [---------] [----------------] [---------]   .
[ X21 | X22 ]   [    | P2 ] [ B21 | B22 0  0 ] [    | Q2 ]
[  0  |  0  0  I ]

X11 is P-by-Q. Q must be no larger than P, M-P, or M-Q. (If this is
not the case, then X must be transposed and/or permuted. This can be
done in constant time using the TRANS and SIGNS options. See ZUNCSD
for details.)

The unitary matrices P1, P2, Q1, and Q2 are P-by-P, (M-P)-by-
(M-P), Q-by-Q, and (M-Q)-by-(M-Q), respectively. They are
represented implicitly by Householder vectors.

B11, B12, B21, and B22 are Q-by-Q bidiagonal matrices represented
implicitly by angles THETA, PHI.```
Parameters
 [in] TRANS ``` TRANS is CHARACTER = 'T': X, U1, U2, V1T, and V2T are stored in row-major order; otherwise: X, U1, U2, V1T, and V2T are stored in column- major order.``` [in] SIGNS ``` SIGNS is CHARACTER = 'O': The lower-left block is made nonpositive (the "other" convention); otherwise: The upper-right block is made nonpositive (the "default" convention).``` [in] M ``` M is INTEGER The number of rows and columns in X.``` [in] P ``` P is INTEGER The number of rows in X11 and X12. 0 <= P <= M.``` [in] Q ``` Q is INTEGER The number of columns in X11 and X21. 0 <= Q <= MIN(P,M-P,M-Q).``` [in,out] X11 ``` X11 is COMPLEX*16 array, dimension (LDX11,Q) On entry, the top-left block of the unitary matrix to be reduced. On exit, the form depends on TRANS: If TRANS = 'N', then the columns of tril(X11) specify reflectors for P1, the rows of triu(X11,1) specify reflectors for Q1; else TRANS = 'T', and the rows of triu(X11) specify reflectors for P1, the columns of tril(X11,-1) specify reflectors for Q1.``` [in] LDX11 ``` LDX11 is INTEGER The leading dimension of X11. If TRANS = 'N', then LDX11 >= P; else LDX11 >= Q.``` [in,out] X12 ``` X12 is COMPLEX*16 array, dimension (LDX12,M-Q) On entry, the top-right block of the unitary matrix to be reduced. On exit, the form depends on TRANS: If TRANS = 'N', then the rows of triu(X12) specify the first P reflectors for Q2; else TRANS = 'T', and the columns of tril(X12) specify the first P reflectors for Q2.``` [in] LDX12 ``` LDX12 is INTEGER The leading dimension of X12. If TRANS = 'N', then LDX12 >= P; else LDX11 >= M-Q.``` [in,out] X21 ``` X21 is COMPLEX*16 array, dimension (LDX21,Q) On entry, the bottom-left block of the unitary matrix to be reduced. On exit, the form depends on TRANS: If TRANS = 'N', then the columns of tril(X21) specify reflectors for P2; else TRANS = 'T', and the rows of triu(X21) specify reflectors for P2.``` [in] LDX21 ``` LDX21 is INTEGER The leading dimension of X21. If TRANS = 'N', then LDX21 >= M-P; else LDX21 >= Q.``` [in,out] X22 ``` X22 is COMPLEX*16 array, dimension (LDX22,M-Q) On entry, the bottom-right block of the unitary matrix to be reduced. On exit, the form depends on TRANS: If TRANS = 'N', then the rows of triu(X22(Q+1:M-P,P+1:M-Q)) specify the last M-P-Q reflectors for Q2, else TRANS = 'T', and the columns of tril(X22(P+1:M-Q,Q+1:M-P)) specify the last M-P-Q reflectors for P2.``` [in] LDX22 ``` LDX22 is INTEGER The leading dimension of X22. If TRANS = 'N', then LDX22 >= M-P; else LDX22 >= M-Q.``` [out] THETA ``` THETA is DOUBLE PRECISION array, dimension (Q) The entries of the bidiagonal blocks B11, B12, B21, B22 can be computed from the angles THETA and PHI. See Further Details.``` [out] PHI ``` PHI is DOUBLE PRECISION array, dimension (Q-1) The entries of the bidiagonal blocks B11, B12, B21, B22 can be computed from the angles THETA and PHI. See Further Details.``` [out] TAUP1 ``` TAUP1 is COMPLEX*16 array, dimension (P) The scalar factors of the elementary reflectors that define P1.``` [out] TAUP2 ``` TAUP2 is COMPLEX*16 array, dimension (M-P) The scalar factors of the elementary reflectors that define P2.``` [out] TAUQ1 ``` TAUQ1 is COMPLEX*16 array, dimension (Q) The scalar factors of the elementary reflectors that define Q1.``` [out] TAUQ2 ``` TAUQ2 is COMPLEX*16 array, dimension (M-Q) The scalar factors of the elementary reflectors that define Q2.``` [out] WORK ` WORK is COMPLEX*16 array, dimension (LWORK)` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= M-Q. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] INFO ``` INFO is INTEGER = 0: successful exit. < 0: if INFO = -i, the i-th argument had an illegal value.```
Further Details:
```  The bidiagonal blocks B11, B12, B21, and B22 are represented
implicitly by angles THETA(1), ..., THETA(Q) and PHI(1), ...,
PHI(Q-1). B11 and B21 are upper bidiagonal, while B21 and B22 are
lower bidiagonal. Every entry in each bidiagonal band is a product
of a sine or cosine of a THETA with a sine or cosine of a PHI. See
 or ZUNCSD for details.

P1, P2, Q1, and Q2 are represented as products of elementary
reflectors. See ZUNCSD for details on generating P1, P2, Q1, and Q2
using ZUNGQR and ZUNGLQ.```
References:
 Brian D. Sutton. Computing the complete CS decomposition. Numer. Algorithms, 50(1):33-65, 2009.

Definition at line 284 of file zunbdb.f.

287*
288* -- LAPACK computational routine --
289* -- LAPACK is a software package provided by Univ. of Tennessee, --
290* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
291*
292* .. Scalar Arguments ..
293 CHARACTER SIGNS, TRANS
294 INTEGER INFO, LDX11, LDX12, LDX21, LDX22, LWORK, M, P,
295 \$ Q
296* ..
297* .. Array Arguments ..
298 DOUBLE PRECISION PHI( * ), THETA( * )
299 COMPLEX*16 TAUP1( * ), TAUP2( * ), TAUQ1( * ), TAUQ2( * ),
300 \$ WORK( * ), X11( LDX11, * ), X12( LDX12, * ),
301 \$ X21( LDX21, * ), X22( LDX22, * )
302* ..
303*
304* ====================================================================
305*
306* .. Parameters ..
307 DOUBLE PRECISION REALONE
308 parameter( realone = 1.0d0 )
309 COMPLEX*16 ONE
310 parameter( one = (1.0d0,0.0d0) )
311* ..
312* .. Local Scalars ..
313 LOGICAL COLMAJOR, LQUERY
314 INTEGER I, LWORKMIN, LWORKOPT
315 DOUBLE PRECISION Z1, Z2, Z3, Z4
316* ..
317* .. External Subroutines ..
318 EXTERNAL zaxpy, zlarf, zlarfgp, zscal, xerbla
319 EXTERNAL zlacgv
320*
321* ..
322* .. External Functions ..
323 DOUBLE PRECISION DZNRM2
324 LOGICAL LSAME
325 EXTERNAL dznrm2, lsame
326* ..
327* .. Intrinsic Functions
328 INTRINSIC atan2, cos, max, min, sin
329 INTRINSIC dcmplx, dconjg
330* ..
331* .. Executable Statements ..
332*
333* Test input arguments
334*
335 info = 0
336 colmajor = .NOT. lsame( trans, 'T' )
337 IF( .NOT. lsame( signs, 'O' ) ) THEN
338 z1 = realone
339 z2 = realone
340 z3 = realone
341 z4 = realone
342 ELSE
343 z1 = realone
344 z2 = -realone
345 z3 = realone
346 z4 = -realone
347 END IF
348 lquery = lwork .EQ. -1
349*
350 IF( m .LT. 0 ) THEN
351 info = -3
352 ELSE IF( p .LT. 0 .OR. p .GT. m ) THEN
353 info = -4
354 ELSE IF( q .LT. 0 .OR. q .GT. p .OR. q .GT. m-p .OR.
355 \$ q .GT. m-q ) THEN
356 info = -5
357 ELSE IF( colmajor .AND. ldx11 .LT. max( 1, p ) ) THEN
358 info = -7
359 ELSE IF( .NOT.colmajor .AND. ldx11 .LT. max( 1, q ) ) THEN
360 info = -7
361 ELSE IF( colmajor .AND. ldx12 .LT. max( 1, p ) ) THEN
362 info = -9
363 ELSE IF( .NOT.colmajor .AND. ldx12 .LT. max( 1, m-q ) ) THEN
364 info = -9
365 ELSE IF( colmajor .AND. ldx21 .LT. max( 1, m-p ) ) THEN
366 info = -11
367 ELSE IF( .NOT.colmajor .AND. ldx21 .LT. max( 1, q ) ) THEN
368 info = -11
369 ELSE IF( colmajor .AND. ldx22 .LT. max( 1, m-p ) ) THEN
370 info = -13
371 ELSE IF( .NOT.colmajor .AND. ldx22 .LT. max( 1, m-q ) ) THEN
372 info = -13
373 END IF
374*
375* Compute workspace
376*
377 IF( info .EQ. 0 ) THEN
378 lworkopt = m - q
379 lworkmin = m - q
380 work(1) = lworkopt
381 IF( lwork .LT. lworkmin .AND. .NOT. lquery ) THEN
382 info = -21
383 END IF
384 END IF
385 IF( info .NE. 0 ) THEN
386 CALL xerbla( 'xORBDB', -info )
387 RETURN
388 ELSE IF( lquery ) THEN
389 RETURN
390 END IF
391*
392* Handle column-major and row-major separately
393*
394 IF( colmajor ) THEN
395*
396* Reduce columns 1, ..., Q of X11, X12, X21, and X22
397*
398 DO i = 1, q
399*
400 IF( i .EQ. 1 ) THEN
401 CALL zscal( p-i+1, dcmplx( z1, 0.0d0 ), x11(i,i), 1 )
402 ELSE
403 CALL zscal( p-i+1, dcmplx( z1*cos(phi(i-1)), 0.0d0 ),
404 \$ x11(i,i), 1 )
405 CALL zaxpy( p-i+1, dcmplx( -z1*z3*z4*sin(phi(i-1)),
406 \$ 0.0d0 ), x12(i,i-1), 1, x11(i,i), 1 )
407 END IF
408 IF( i .EQ. 1 ) THEN
409 CALL zscal( m-p-i+1, dcmplx( z2, 0.0d0 ), x21(i,i), 1 )
410 ELSE
411 CALL zscal( m-p-i+1, dcmplx( z2*cos(phi(i-1)), 0.0d0 ),
412 \$ x21(i,i), 1 )
413 CALL zaxpy( m-p-i+1, dcmplx( -z2*z3*z4*sin(phi(i-1)),
414 \$ 0.0d0 ), x22(i,i-1), 1, x21(i,i), 1 )
415 END IF
416*
417 theta(i) = atan2( dznrm2( m-p-i+1, x21(i,i), 1 ),
418 \$ dznrm2( p-i+1, x11(i,i), 1 ) )
419*
420 IF( p .GT. i ) THEN
421 CALL zlarfgp( p-i+1, x11(i,i), x11(i+1,i), 1, taup1(i) )
422 ELSE IF ( p .EQ. i ) THEN
423 CALL zlarfgp( p-i+1, x11(i,i), x11(i,i), 1, taup1(i) )
424 END IF
425 x11(i,i) = one
426 IF ( m-p .GT. i ) THEN
427 CALL zlarfgp( m-p-i+1, x21(i,i), x21(i+1,i), 1,
428 \$ taup2(i) )
429 ELSE IF ( m-p .EQ. i ) THEN
430 CALL zlarfgp( m-p-i+1, x21(i,i), x21(i,i), 1,
431 \$ taup2(i) )
432 END IF
433 x21(i,i) = one
434*
435 IF ( q .GT. i ) THEN
436 CALL zlarf( 'L', p-i+1, q-i, x11(i,i), 1,
437 \$ dconjg(taup1(i)), x11(i,i+1), ldx11, work )
438 CALL zlarf( 'L', m-p-i+1, q-i, x21(i,i), 1,
439 \$ dconjg(taup2(i)), x21(i,i+1), ldx21, work )
440 END IF
441 IF ( m-q+1 .GT. i ) THEN
442 CALL zlarf( 'L', p-i+1, m-q-i+1, x11(i,i), 1,
443 \$ dconjg(taup1(i)), x12(i,i), ldx12, work )
444 CALL zlarf( 'L', m-p-i+1, m-q-i+1, x21(i,i), 1,
445 \$ dconjg(taup2(i)), x22(i,i), ldx22, work )
446 END IF
447*
448 IF( i .LT. q ) THEN
449 CALL zscal( q-i, dcmplx( -z1*z3*sin(theta(i)), 0.0d0 ),
450 \$ x11(i,i+1), ldx11 )
451 CALL zaxpy( q-i, dcmplx( z2*z3*cos(theta(i)), 0.0d0 ),
452 \$ x21(i,i+1), ldx21, x11(i,i+1), ldx11 )
453 END IF
454 CALL zscal( m-q-i+1, dcmplx( -z1*z4*sin(theta(i)), 0.0d0 ),
455 \$ x12(i,i), ldx12 )
456 CALL zaxpy( m-q-i+1, dcmplx( z2*z4*cos(theta(i)), 0.0d0 ),
457 \$ x22(i,i), ldx22, x12(i,i), ldx12 )
458*
459 IF( i .LT. q )
460 \$ phi(i) = atan2( dznrm2( q-i, x11(i,i+1), ldx11 ),
461 \$ dznrm2( m-q-i+1, x12(i,i), ldx12 ) )
462*
463 IF( i .LT. q ) THEN
464 CALL zlacgv( q-i, x11(i,i+1), ldx11 )
465 IF ( i .EQ. q-1 ) THEN
466 CALL zlarfgp( q-i, x11(i,i+1), x11(i,i+1), ldx11,
467 \$ tauq1(i) )
468 ELSE
469 CALL zlarfgp( q-i, x11(i,i+1), x11(i,i+2), ldx11,
470 \$ tauq1(i) )
471 END IF
472 x11(i,i+1) = one
473 END IF
474 IF ( m-q+1 .GT. i ) THEN
475 CALL zlacgv( m-q-i+1, x12(i,i), ldx12 )
476 IF ( m-q .EQ. i ) THEN
477 CALL zlarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
478 \$ tauq2(i) )
479 ELSE
480 CALL zlarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
481 \$ tauq2(i) )
482 END IF
483 END IF
484 x12(i,i) = one
485*
486 IF( i .LT. q ) THEN
487 CALL zlarf( 'R', p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
488 \$ x11(i+1,i+1), ldx11, work )
489 CALL zlarf( 'R', m-p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
490 \$ x21(i+1,i+1), ldx21, work )
491 END IF
492 IF ( p .GT. i ) THEN
493 CALL zlarf( 'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
494 \$ x12(i+1,i), ldx12, work )
495 END IF
496 IF ( m-p .GT. i ) THEN
497 CALL zlarf( 'R', m-p-i, m-q-i+1, x12(i,i), ldx12,
498 \$ tauq2(i), x22(i+1,i), ldx22, work )
499 END IF
500*
501 IF( i .LT. q )
502 \$ CALL zlacgv( q-i, x11(i,i+1), ldx11 )
503 CALL zlacgv( m-q-i+1, x12(i,i), ldx12 )
504*
505 END DO
506*
507* Reduce columns Q + 1, ..., P of X12, X22
508*
509 DO i = q + 1, p
510*
511 CALL zscal( m-q-i+1, dcmplx( -z1*z4, 0.0d0 ), x12(i,i),
512 \$ ldx12 )
513 CALL zlacgv( m-q-i+1, x12(i,i), ldx12 )
514 IF ( i .GE. m-q ) THEN
515 CALL zlarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
516 \$ tauq2(i) )
517 ELSE
518 CALL zlarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
519 \$ tauq2(i) )
520 END IF
521 x12(i,i) = one
522*
523 IF ( p .GT. i ) THEN
524 CALL zlarf( 'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
525 \$ x12(i+1,i), ldx12, work )
526 END IF
527 IF( m-p-q .GE. 1 )
528 \$ CALL zlarf( 'R', m-p-q, m-q-i+1, x12(i,i), ldx12,
529 \$ tauq2(i), x22(q+1,i), ldx22, work )
530*
531 CALL zlacgv( m-q-i+1, x12(i,i), ldx12 )
532*
533 END DO
534*
535* Reduce columns P + 1, ..., M - Q of X12, X22
536*
537 DO i = 1, m - p - q
538*
539 CALL zscal( m-p-q-i+1, dcmplx( z2*z4, 0.0d0 ),
540 \$ x22(q+i,p+i), ldx22 )
541 CALL zlacgv( m-p-q-i+1, x22(q+i,p+i), ldx22 )
542 CALL zlarfgp( m-p-q-i+1, x22(q+i,p+i), x22(q+i,p+i+1),
543 \$ ldx22, tauq2(p+i) )
544 x22(q+i,p+i) = one
545 CALL zlarf( 'R', m-p-q-i, m-p-q-i+1, x22(q+i,p+i), ldx22,
546 \$ tauq2(p+i), x22(q+i+1,p+i), ldx22, work )
547*
548 CALL zlacgv( m-p-q-i+1, x22(q+i,p+i), ldx22 )
549*
550 END DO
551*
552 ELSE
553*
554* Reduce columns 1, ..., Q of X11, X12, X21, X22
555*
556 DO i = 1, q
557*
558 IF( i .EQ. 1 ) THEN
559 CALL zscal( p-i+1, dcmplx( z1, 0.0d0 ), x11(i,i),
560 \$ ldx11 )
561 ELSE
562 CALL zscal( p-i+1, dcmplx( z1*cos(phi(i-1)), 0.0d0 ),
563 \$ x11(i,i), ldx11 )
564 CALL zaxpy( p-i+1, dcmplx( -z1*z3*z4*sin(phi(i-1)),
565 \$ 0.0d0 ), x12(i-1,i), ldx12, x11(i,i), ldx11 )
566 END IF
567 IF( i .EQ. 1 ) THEN
568 CALL zscal( m-p-i+1, dcmplx( z2, 0.0d0 ), x21(i,i),
569 \$ ldx21 )
570 ELSE
571 CALL zscal( m-p-i+1, dcmplx( z2*cos(phi(i-1)), 0.0d0 ),
572 \$ x21(i,i), ldx21 )
573 CALL zaxpy( m-p-i+1, dcmplx( -z2*z3*z4*sin(phi(i-1)),
574 \$ 0.0d0 ), x22(i-1,i), ldx22, x21(i,i), ldx21 )
575 END IF
576*
577 theta(i) = atan2( dznrm2( m-p-i+1, x21(i,i), ldx21 ),
578 \$ dznrm2( p-i+1, x11(i,i), ldx11 ) )
579*
580 CALL zlacgv( p-i+1, x11(i,i), ldx11 )
581 CALL zlacgv( m-p-i+1, x21(i,i), ldx21 )
582*
583 CALL zlarfgp( p-i+1, x11(i,i), x11(i,i+1), ldx11, taup1(i) )
584 x11(i,i) = one
585 IF ( i .EQ. m-p ) THEN
586 CALL zlarfgp( m-p-i+1, x21(i,i), x21(i,i), ldx21,
587 \$ taup2(i) )
588 ELSE
589 CALL zlarfgp( m-p-i+1, x21(i,i), x21(i,i+1), ldx21,
590 \$ taup2(i) )
591 END IF
592 x21(i,i) = one
593*
594 CALL zlarf( 'R', q-i, p-i+1, x11(i,i), ldx11, taup1(i),
595 \$ x11(i+1,i), ldx11, work )
596 CALL zlarf( 'R', m-q-i+1, p-i+1, x11(i,i), ldx11, taup1(i),
597 \$ x12(i,i), ldx12, work )
598 CALL zlarf( 'R', q-i, m-p-i+1, x21(i,i), ldx21, taup2(i),
599 \$ x21(i+1,i), ldx21, work )
600 CALL zlarf( 'R', m-q-i+1, m-p-i+1, x21(i,i), ldx21,
601 \$ taup2(i), x22(i,i), ldx22, work )
602*
603 CALL zlacgv( p-i+1, x11(i,i), ldx11 )
604 CALL zlacgv( m-p-i+1, x21(i,i), ldx21 )
605*
606 IF( i .LT. q ) THEN
607 CALL zscal( q-i, dcmplx( -z1*z3*sin(theta(i)), 0.0d0 ),
608 \$ x11(i+1,i), 1 )
609 CALL zaxpy( q-i, dcmplx( z2*z3*cos(theta(i)), 0.0d0 ),
610 \$ x21(i+1,i), 1, x11(i+1,i), 1 )
611 END IF
612 CALL zscal( m-q-i+1, dcmplx( -z1*z4*sin(theta(i)), 0.0d0 ),
613 \$ x12(i,i), 1 )
614 CALL zaxpy( m-q-i+1, dcmplx( z2*z4*cos(theta(i)), 0.0d0 ),
615 \$ x22(i,i), 1, x12(i,i), 1 )
616*
617 IF( i .LT. q )
618 \$ phi(i) = atan2( dznrm2( q-i, x11(i+1,i), 1 ),
619 \$ dznrm2( m-q-i+1, x12(i,i), 1 ) )
620*
621 IF( i .LT. q ) THEN
622 CALL zlarfgp( q-i, x11(i+1,i), x11(i+2,i), 1, tauq1(i) )
623 x11(i+1,i) = one
624 END IF
625 CALL zlarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1, tauq2(i) )
626 x12(i,i) = one
627*
628 IF( i .LT. q ) THEN
629 CALL zlarf( 'L', q-i, p-i, x11(i+1,i), 1,
630 \$ dconjg(tauq1(i)), x11(i+1,i+1), ldx11, work )
631 CALL zlarf( 'L', q-i, m-p-i, x11(i+1,i), 1,
632 \$ dconjg(tauq1(i)), x21(i+1,i+1), ldx21, work )
633 END IF
634 CALL zlarf( 'L', m-q-i+1, p-i, x12(i,i), 1,
635 \$ dconjg(tauq2(i)), x12(i,i+1), ldx12, work )
636 IF ( m-p .GT. i ) THEN
637 CALL zlarf( 'L', m-q-i+1, m-p-i, x12(i,i), 1,
638 \$ dconjg(tauq2(i)), x22(i,i+1), ldx22, work )
639 END IF
640*
641 END DO
642*
643* Reduce columns Q + 1, ..., P of X12, X22
644*
645 DO i = q + 1, p
646*
647 CALL zscal( m-q-i+1, dcmplx( -z1*z4, 0.0d0 ), x12(i,i), 1 )
648 CALL zlarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1, tauq2(i) )
649 x12(i,i) = one
650*
651 IF ( p .GT. i ) THEN
652 CALL zlarf( 'L', m-q-i+1, p-i, x12(i,i), 1,
653 \$ dconjg(tauq2(i)), x12(i,i+1), ldx12, work )
654 END IF
655 IF( m-p-q .GE. 1 )
656 \$ CALL zlarf( 'L', m-q-i+1, m-p-q, x12(i,i), 1,
657 \$ dconjg(tauq2(i)), x22(i,q+1), ldx22, work )
658*
659 END DO
660*
661* Reduce columns P + 1, ..., M - Q of X12, X22
662*
663 DO i = 1, m - p - q
664*
665 CALL zscal( m-p-q-i+1, dcmplx( z2*z4, 0.0d0 ),
666 \$ x22(p+i,q+i), 1 )
667 CALL zlarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i+1,q+i), 1,
668 \$ tauq2(p+i) )
669 x22(p+i,q+i) = one
670*
671 IF ( m-p-q .NE. i ) THEN
672 CALL zlarf( 'L', m-p-q-i+1, m-p-q-i, x22(p+i,q+i), 1,
673 \$ dconjg(tauq2(p+i)), x22(p+i,q+i+1), ldx22,
674 \$ work )
675 END IF
676*
677 END DO
678*
679 END IF
680*
681 RETURN
682*
683* End of ZUNBDB
684*
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine zaxpy(N, ZA, ZX, INCX, ZY, INCY)
ZAXPY
Definition: zaxpy.f:88
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL
Definition: zscal.f:78
subroutine zlarfgp(N, ALPHA, X, INCX, TAU)
ZLARFGP generates an elementary reflector (Householder matrix) with non-negative beta.
Definition: zlarfgp.f:104
subroutine zlacgv(N, X, INCX)
ZLACGV conjugates a complex vector.
Definition: zlacgv.f:74
subroutine zlarf(SIDE, M, N, V, INCV, TAU, C, LDC, WORK)
ZLARF applies an elementary reflector to a general rectangular matrix.
Definition: zlarf.f:128
real(wp) function dznrm2(n, x, incx)
DZNRM2
Definition: dznrm2.f90:90
Here is the call graph for this function:
Here is the caller graph for this function: