286 SUBROUTINE cunbdb( TRANS, SIGNS, M, P, Q, X11, LDX11, X12, LDX12,
287 $ x21, ldx21, x22, ldx22, theta, phi, taup1,
288 $ taup2, tauq1, tauq2, work, lwork, info )
296 CHARACTER signs, trans
297 INTEGER info, ldx11, ldx12, ldx21, ldx22, lwork, m, p,
301 REAL phi( * ), theta( * )
302 COMPLEX taup1( * ), taup2( * ), tauq1( * ), tauq2( * ),
303 $ work( * ), x11( ldx11, * ), x12( ldx12, * ),
304 $ x21( ldx21, * ), x22( ldx22, * )
311 parameter( realone = 1.0e0 )
313 parameter( one = (1.0e0,0.0e0) )
316 LOGICAL colmajor, lquery
317 INTEGER i, lworkmin, lworkopt
331 INTRINSIC atan2, cos, max, min, sin
332 INTRINSIC cmplx, conjg
339 colmajor = .NOT.
lsame( trans,
'T' )
340 IF( .NOT.
lsame( signs,
'O' ) )
THEN
351 lquery = lwork .EQ. -1
355 ELSE IF( p .LT. 0 .OR. p .GT. m )
THEN
357 ELSE IF( q .LT. 0 .OR. q .GT. p .OR. q .GT. m-p .OR.
360 ELSE IF( colmajor .AND. ldx11 .LT. max( 1, p ) )
THEN
362 ELSE IF( .NOT.colmajor .AND. ldx11 .LT. max( 1, q ) )
THEN
364 ELSE IF( colmajor .AND. ldx12 .LT. max( 1, p ) )
THEN
366 ELSE IF( .NOT.colmajor .AND. ldx12 .LT. max( 1, m-q ) )
THEN
368 ELSE IF( colmajor .AND. ldx21 .LT. max( 1, m-p ) )
THEN
370 ELSE IF( .NOT.colmajor .AND. ldx21 .LT. max( 1, q ) )
THEN
372 ELSE IF( colmajor .AND. ldx22 .LT. max( 1, m-p ) )
THEN
374 ELSE IF( .NOT.colmajor .AND. ldx22 .LT. max( 1, m-q ) )
THEN
380 IF( info .EQ. 0 )
THEN
384 IF( lwork .LT. lworkmin .AND. .NOT. lquery )
THEN
388 IF( info .NE. 0 )
THEN
389 CALL
xerbla(
'xORBDB', -info )
391 ELSE IF( lquery )
THEN
404 CALL
cscal( p-i+1, cmplx( z1, 0.0e0 ), x11(i,i), 1 )
406 CALL
cscal( p-i+1, cmplx( z1*cos(phi(i-1)), 0.0e0 ),
408 CALL
caxpy( p-i+1, cmplx( -z1*z3*z4*sin(phi(i-1)),
409 $ 0.0e0 ), x12(i,i-1), 1, x11(i,i), 1 )
412 CALL
cscal( m-p-i+1, cmplx( z2, 0.0e0 ), x21(i,i), 1 )
414 CALL
cscal( m-p-i+1, cmplx( z2*cos(phi(i-1)), 0.0e0 ),
416 CALL
caxpy( m-p-i+1, cmplx( -z2*z3*z4*sin(phi(i-1)),
417 $ 0.0e0 ), x22(i,i-1), 1, x21(i,i), 1 )
420 theta(i) = atan2(
scnrm2( m-p-i+1, x21(i,i), 1 ),
421 $
scnrm2( p-i+1, x11(i,i), 1 ) )
424 CALL
clarfgp( p-i+1, x11(i,i), x11(i+1,i), 1, taup1(i) )
425 ELSE IF ( p .EQ. i )
THEN
426 CALL
clarfgp( p-i+1, x11(i,i), x11(i,i), 1, taup1(i) )
429 IF ( m-p .GT. i )
THEN
430 CALL
clarfgp( m-p-i+1, x21(i,i), x21(i+1,i), 1,
432 ELSE IF ( m-p .EQ. i )
THEN
433 CALL
clarfgp( m-p-i+1, x21(i,i), x21(i,i), 1,
439 CALL
clarf(
'L', p-i+1, q-i, x11(i,i), 1,
440 $ conjg(taup1(i)), x11(i,i+1), ldx11, work )
441 CALL
clarf(
'L', m-p-i+1, q-i, x21(i,i), 1,
442 $ conjg(taup2(i)), x21(i,i+1), ldx21, work )
444 IF ( m-q+1 .GT. i )
THEN
445 CALL
clarf(
'L', p-i+1, m-q-i+1, x11(i,i), 1,
446 $ conjg(taup1(i)), x12(i,i), ldx12, work )
447 CALL
clarf(
'L', m-p-i+1, m-q-i+1, x21(i,i), 1,
448 $ conjg(taup2(i)), x22(i,i), ldx22, work )
452 CALL
cscal( q-i, cmplx( -z1*z3*sin(theta(i)), 0.0e0 ),
453 $ x11(i,i+1), ldx11 )
454 CALL
caxpy( q-i, cmplx( z2*z3*cos(theta(i)), 0.0e0 ),
455 $ x21(i,i+1), ldx21, x11(i,i+1), ldx11 )
457 CALL
cscal( m-q-i+1, cmplx( -z1*z4*sin(theta(i)), 0.0e0 ),
459 CALL
caxpy( m-q-i+1, cmplx( z2*z4*cos(theta(i)), 0.0e0 ),
460 $ x22(i,i), ldx22, x12(i,i), ldx12 )
463 $ phi(i) = atan2(
scnrm2( q-i, x11(i,i+1), ldx11 ),
464 $
scnrm2( m-q-i+1, x12(i,i), ldx12 ) )
467 CALL
clacgv( q-i, x11(i,i+1), ldx11 )
468 IF ( i .EQ. q-1 )
THEN
469 CALL
clarfgp( q-i, x11(i,i+1), x11(i,i+1), ldx11,
472 CALL
clarfgp( q-i, x11(i,i+1), x11(i,i+2), ldx11,
477 IF ( m-q+1 .GT. i )
THEN
478 CALL
clacgv( m-q-i+1, x12(i,i), ldx12 )
479 IF ( m-q .EQ. i )
THEN
480 CALL
clarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
483 CALL
clarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
490 CALL
clarf(
'R', p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
491 $ x11(i+1,i+1), ldx11, work )
492 CALL
clarf(
'R', m-p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
493 $ x21(i+1,i+1), ldx21, work )
496 CALL
clarf(
'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
497 $ x12(i+1,i), ldx12, work )
499 IF ( m-p .GT. i )
THEN
500 CALL
clarf(
'R', m-p-i, m-q-i+1, x12(i,i), ldx12,
501 $ tauq2(i), x22(i+1,i), ldx22, work )
505 $ CALL
clacgv( q-i, x11(i,i+1), ldx11 )
506 CALL
clacgv( m-q-i+1, x12(i,i), ldx12 )
514 CALL
cscal( m-q-i+1, cmplx( -z1*z4, 0.0e0 ), x12(i,i),
516 CALL
clacgv( m-q-i+1, x12(i,i), ldx12 )
517 IF ( i .GE. m-q )
THEN
518 CALL
clarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
521 CALL
clarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
527 CALL
clarf(
'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
528 $ x12(i+1,i), ldx12, work )
531 $ CALL
clarf(
'R', m-p-q, m-q-i+1, x12(i,i), ldx12,
532 $ tauq2(i), x22(q+1,i), ldx22, work )
534 CALL
clacgv( m-q-i+1, x12(i,i), ldx12 )
542 CALL
cscal( m-p-q-i+1, cmplx( z2*z4, 0.0e0 ),
543 $ x22(q+i,p+i), ldx22 )
544 CALL
clacgv( m-p-q-i+1, x22(q+i,p+i), ldx22 )
545 CALL
clarfgp( m-p-q-i+1, x22(q+i,p+i), x22(q+i,p+i+1),
546 $ ldx22, tauq2(p+i) )
548 CALL
clarf(
'R', m-p-q-i, m-p-q-i+1, x22(q+i,p+i), ldx22,
549 $ tauq2(p+i), x22(q+i+1,p+i), ldx22, work )
551 CALL
clacgv( m-p-q-i+1, x22(q+i,p+i), ldx22 )
562 CALL
cscal( p-i+1, cmplx( z1, 0.0e0 ), x11(i,i),
565 CALL
cscal( p-i+1, cmplx( z1*cos(phi(i-1)), 0.0e0 ),
567 CALL
caxpy( p-i+1, cmplx( -z1*z3*z4*sin(phi(i-1)),
568 $ 0.0e0 ), x12(i-1,i), ldx12, x11(i,i), ldx11 )
571 CALL
cscal( m-p-i+1, cmplx( z2, 0.0e0 ), x21(i,i),
574 CALL
cscal( m-p-i+1, cmplx( z2*cos(phi(i-1)), 0.0e0 ),
576 CALL
caxpy( m-p-i+1, cmplx( -z2*z3*z4*sin(phi(i-1)),
577 $ 0.0e0 ), x22(i-1,i), ldx22, x21(i,i), ldx21 )
580 theta(i) = atan2(
scnrm2( m-p-i+1, x21(i,i), ldx21 ),
581 $
scnrm2( p-i+1, x11(i,i), ldx11 ) )
583 CALL
clacgv( p-i+1, x11(i,i), ldx11 )
584 CALL
clacgv( m-p-i+1, x21(i,i), ldx21 )
586 CALL
clarfgp( p-i+1, x11(i,i), x11(i,i+1), ldx11, taup1(i) )
588 IF ( i .EQ. m-p )
THEN
589 CALL
clarfgp( m-p-i+1, x21(i,i), x21(i,i), ldx21,
592 CALL
clarfgp( m-p-i+1, x21(i,i), x21(i,i+1), ldx21,
597 CALL
clarf(
'R', q-i, p-i+1, x11(i,i), ldx11, taup1(i),
598 $ x11(i+1,i), ldx11, work )
599 CALL
clarf(
'R', m-q-i+1, p-i+1, x11(i,i), ldx11, taup1(i),
600 $ x12(i,i), ldx12, work )
601 CALL
clarf(
'R', q-i, m-p-i+1, x21(i,i), ldx21, taup2(i),
602 $ x21(i+1,i), ldx21, work )
603 CALL
clarf(
'R', m-q-i+1, m-p-i+1, x21(i,i), ldx21,
604 $ taup2(i), x22(i,i), ldx22, work )
606 CALL
clacgv( p-i+1, x11(i,i), ldx11 )
607 CALL
clacgv( m-p-i+1, x21(i,i), ldx21 )
610 CALL
cscal( q-i, cmplx( -z1*z3*sin(theta(i)), 0.0e0 ),
612 CALL
caxpy( q-i, cmplx( z2*z3*cos(theta(i)), 0.0e0 ),
613 $ x21(i+1,i), 1, x11(i+1,i), 1 )
615 CALL
cscal( m-q-i+1, cmplx( -z1*z4*sin(theta(i)), 0.0e0 ),
617 CALL
caxpy( m-q-i+1, cmplx( z2*z4*cos(theta(i)), 0.0e0 ),
618 $ x22(i,i), 1, x12(i,i), 1 )
621 $ phi(i) = atan2(
scnrm2( q-i, x11(i+1,i), 1 ),
622 $
scnrm2( m-q-i+1, x12(i,i), 1 ) )
625 CALL
clarfgp( q-i, x11(i+1,i), x11(i+2,i), 1, tauq1(i) )
628 CALL
clarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1, tauq2(i) )
632 CALL
clarf(
'L', q-i, p-i, x11(i+1,i), 1,
633 $ conjg(tauq1(i)), x11(i+1,i+1), ldx11, work )
634 CALL
clarf(
'L', q-i, m-p-i, x11(i+1,i), 1,
635 $ conjg(tauq1(i)), x21(i+1,i+1), ldx21, work )
637 CALL
clarf(
'L', m-q-i+1, p-i, x12(i,i), 1, conjg(tauq2(i)),
638 $ x12(i,i+1), ldx12, work )
640 IF ( m-p .GT. i )
THEN
641 CALL
clarf(
'L', m-q-i+1, m-p-i, x12(i,i), 1,
642 $ conjg(tauq2(i)), x22(i,i+1), ldx22, work )
650 CALL
cscal( m-q-i+1, cmplx( -z1*z4, 0.0e0 ), x12(i,i), 1 )
651 CALL
clarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1, tauq2(i) )
655 CALL
clarf(
'L', m-q-i+1, p-i, x12(i,i), 1,
656 $ conjg(tauq2(i)), x12(i,i+1), ldx12, work )
659 $ CALL
clarf(
'L', m-q-i+1, m-p-q, x12(i,i), 1,
660 $ conjg(tauq2(i)), x22(i,q+1), ldx22, work )
668 CALL
cscal( m-p-q-i+1, cmplx( z2*z4, 0.0e0 ),
670 CALL
clarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i+1,q+i), 1,
673 IF ( m-p-q .NE. i )
THEN
674 CALL
clarf(
'L', m-p-q-i+1, m-p-q-i, x22(p+i,q+i), 1,
675 $ conjg(tauq2(p+i)), x22(p+i,q+i+1), ldx22,
LOGICAL function lsame(CA, CB)
LSAME
subroutine cscal(N, CA, CX, INCX)
CSCAL
subroutine cunbdb(TRANS, SIGNS, M, P, Q, X11, LDX11, X12, LDX12, X21, LDX21, X22, LDX22, THETA, PHI, TAUP1, TAUP2, TAUQ1, TAUQ2, WORK, LWORK, INFO)
CUNBDB
subroutine caxpy(N, CA, CX, INCX, CY, INCY)
CAXPY
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine clacgv(N, X, INCX)
CLACGV conjugates a complex vector.
subroutine clarfgp(N, ALPHA, X, INCX, TAU)
CLARFGP generates an elementary reflector (Householder matrix) with non-negatibe beta.
REAL function scnrm2(N, X, INCX)
SCNRM2
subroutine clarf(SIDE, M, N, V, INCV, TAU, C, LDC, WORK)
CLARF applies an elementary reflector to a general rectangular matrix.