286 SUBROUTINE sorbdb( 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 REAL 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 )
316 LOGICAL colmajor, lquery
317 INTEGER i, lworkmin, lworkopt
329 INTRINSIC atan2, cos, max, sin
336 colmajor = .NOT.
lsame( trans,
'T' )
337 IF( .NOT.
lsame( signs,
'O' ) )
THEN
348 lquery = lwork .EQ. -1
352 ELSE IF( p .LT. 0 .OR. p .GT. m )
THEN
354 ELSE IF( q .LT. 0 .OR. q .GT. p .OR. q .GT. m-p .OR.
357 ELSE IF( colmajor .AND. ldx11 .LT. max( 1, p ) )
THEN
359 ELSE IF( .NOT.colmajor .AND. ldx11 .LT. max( 1, q ) )
THEN
361 ELSE IF( colmajor .AND. ldx12 .LT. max( 1, p ) )
THEN
363 ELSE IF( .NOT.colmajor .AND. ldx12 .LT. max( 1, m-q ) )
THEN
365 ELSE IF( colmajor .AND. ldx21 .LT. max( 1, m-p ) )
THEN
367 ELSE IF( .NOT.colmajor .AND. ldx21 .LT. max( 1, q ) )
THEN
369 ELSE IF( colmajor .AND. ldx22 .LT. max( 1, m-p ) )
THEN
371 ELSE IF( .NOT.colmajor .AND. ldx22 .LT. max( 1, m-q ) )
THEN
377 IF( info .EQ. 0 )
THEN
381 IF( lwork .LT. lworkmin .AND. .NOT. lquery )
THEN
385 IF( info .NE. 0 )
THEN
386 CALL
xerbla(
'xORBDB', -info )
388 ELSE IF( lquery )
THEN
401 CALL
sscal( p-i+1, z1, x11(i,i), 1 )
403 CALL
sscal( p-i+1, z1*cos(phi(i-1)), x11(i,i), 1 )
404 CALL
saxpy( p-i+1, -z1*z3*z4*sin(phi(i-1)), x12(i,i-1),
408 CALL
sscal( m-p-i+1, z2, x21(i,i), 1 )
410 CALL
sscal( m-p-i+1, z2*cos(phi(i-1)), x21(i,i), 1 )
411 CALL
saxpy( m-p-i+1, -z2*z3*z4*sin(phi(i-1)), x22(i,i-1),
415 theta(i) = atan2(
snrm2( m-p-i+1, x21(i,i), 1 ),
416 $
snrm2( p-i+1, x11(i,i), 1 ) )
419 CALL
slarfgp( p-i+1, x11(i,i), x11(i+1,i), 1, taup1(i) )
420 ELSE IF( p .EQ. i )
THEN
421 CALL
slarfgp( p-i+1, x11(i,i), x11(i,i), 1, taup1(i) )
424 IF ( m-p .GT. i )
THEN
425 CALL
slarfgp( m-p-i+1, x21(i,i), x21(i+1,i), 1,
427 ELSE IF ( m-p .EQ. i )
THEN
428 CALL
slarfgp( m-p-i+1, x21(i,i), x21(i,i), 1, taup2(i) )
433 CALL
slarf(
'L', p-i+1, q-i, x11(i,i), 1, taup1(i),
434 $ x11(i,i+1), ldx11, work )
436 IF ( m-q+1 .GT. i )
THEN
437 CALL
slarf(
'L', p-i+1, m-q-i+1, x11(i,i), 1, taup1(i),
438 $ x12(i,i), ldx12, work )
441 CALL
slarf(
'L', m-p-i+1, q-i, x21(i,i), 1, taup2(i),
442 $ x21(i,i+1), ldx21, work )
444 IF ( m-q+1 .GT. i )
THEN
445 CALL
slarf(
'L', m-p-i+1, m-q-i+1, x21(i,i), 1, taup2(i),
446 $ x22(i,i), ldx22, work )
450 CALL
sscal( q-i, -z1*z3*sin(theta(i)), x11(i,i+1),
452 CALL
saxpy( q-i, z2*z3*cos(theta(i)), x21(i,i+1), ldx21,
453 $ x11(i,i+1), ldx11 )
455 CALL
sscal( m-q-i+1, -z1*z4*sin(theta(i)), x12(i,i), ldx12 )
456 CALL
saxpy( m-q-i+1, z2*z4*cos(theta(i)), x22(i,i), ldx22,
460 $ phi(i) = atan2(
snrm2( q-i, x11(i,i+1), ldx11 ),
461 $
snrm2( m-q-i+1, x12(i,i), ldx12 ) )
464 IF ( q-i .EQ. 1 )
THEN
465 CALL
slarfgp( q-i, x11(i,i+1), x11(i,i+1), ldx11,
468 CALL
slarfgp( q-i, x11(i,i+1), x11(i,i+2), ldx11,
473 IF ( q+i-1 .LT. m )
THEN
474 IF ( m-q .EQ. i )
THEN
475 CALL
slarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
478 CALL
slarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
485 CALL
slarf(
'R', p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
486 $ x11(i+1,i+1), ldx11, work )
487 CALL
slarf(
'R', m-p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
488 $ x21(i+1,i+1), ldx21, work )
491 CALL
slarf(
'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
492 $ x12(i+1,i), ldx12, work )
494 IF ( m-p .GT. i )
THEN
495 CALL
slarf(
'R', m-p-i, m-q-i+1, x12(i,i), ldx12,
496 $ tauq2(i), x22(i+1,i), ldx22, work )
505 CALL
sscal( m-q-i+1, -z1*z4, x12(i,i), ldx12 )
506 IF ( i .GE. m-q )
THEN
507 CALL
slarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
510 CALL
slarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
516 CALL
slarf(
'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
517 $ x12(i+1,i), ldx12, work )
520 $ CALL
slarf(
'R', m-p-q, m-q-i+1, x12(i,i), ldx12,
521 $ tauq2(i), x22(q+1,i), ldx22, work )
529 CALL
sscal( m-p-q-i+1, z2*z4, x22(q+i,p+i), ldx22 )
530 IF ( i .EQ. m-p-q )
THEN
531 CALL
slarfgp( m-p-q-i+1, x22(q+i,p+i), x22(q+i,p+i),
532 $ ldx22, tauq2(p+i) )
534 CALL
slarfgp( m-p-q-i+1, x22(q+i,p+i), x22(q+i,p+i+1),
535 $ ldx22, tauq2(p+i) )
538 IF ( i .LT. m-p-q )
THEN
539 CALL
slarf(
'R', m-p-q-i, m-p-q-i+1, x22(q+i,p+i), ldx22,
540 $ tauq2(p+i), x22(q+i+1,p+i), ldx22, work )
552 CALL
sscal( p-i+1, z1, x11(i,i), ldx11 )
554 CALL
sscal( p-i+1, z1*cos(phi(i-1)), x11(i,i), ldx11 )
555 CALL
saxpy( p-i+1, -z1*z3*z4*sin(phi(i-1)), x12(i-1,i),
556 $ ldx12, x11(i,i), ldx11 )
559 CALL
sscal( m-p-i+1, z2, x21(i,i), ldx21 )
561 CALL
sscal( m-p-i+1, z2*cos(phi(i-1)), x21(i,i), ldx21 )
562 CALL
saxpy( m-p-i+1, -z2*z3*z4*sin(phi(i-1)), x22(i-1,i),
563 $ ldx22, x21(i,i), ldx21 )
566 theta(i) = atan2(
snrm2( m-p-i+1, x21(i,i), ldx21 ),
567 $
snrm2( p-i+1, x11(i,i), ldx11 ) )
569 CALL
slarfgp( p-i+1, x11(i,i), x11(i,i+1), ldx11, taup1(i) )
571 IF ( i .EQ. m-p )
THEN
572 CALL
slarfgp( m-p-i+1, x21(i,i), x21(i,i), ldx21,
575 CALL
slarfgp( m-p-i+1, x21(i,i), x21(i,i+1), ldx21,
581 CALL
slarf(
'R', q-i, p-i+1, x11(i,i), ldx11, taup1(i),
582 $ x11(i+1,i), ldx11, work )
584 IF ( m-q+1 .GT. i )
THEN
585 CALL
slarf(
'R', m-q-i+1, p-i+1, x11(i,i), ldx11,
586 $ taup1(i), x12(i,i), ldx12, work )
589 CALL
slarf(
'R', q-i, m-p-i+1, x21(i,i), ldx21, taup2(i),
590 $ x21(i+1,i), ldx21, work )
592 IF ( m-q+1 .GT. i )
THEN
593 CALL
slarf(
'R', m-q-i+1, m-p-i+1, x21(i,i), ldx21,
594 $ taup2(i), x22(i,i), ldx22, work )
598 CALL
sscal( q-i, -z1*z3*sin(theta(i)), x11(i+1,i), 1 )
599 CALL
saxpy( q-i, z2*z3*cos(theta(i)), x21(i+1,i), 1,
602 CALL
sscal( m-q-i+1, -z1*z4*sin(theta(i)), x12(i,i), 1 )
603 CALL
saxpy( m-q-i+1, z2*z4*cos(theta(i)), x22(i,i), 1,
607 $ phi(i) = atan2(
snrm2( q-i, x11(i+1,i), 1 ),
608 $
snrm2( m-q-i+1, x12(i,i), 1 ) )
611 IF ( q-i .EQ. 1)
THEN
612 CALL
slarfgp( q-i, x11(i+1,i), x11(i+1,i), 1,
615 CALL
slarfgp( q-i, x11(i+1,i), x11(i+2,i), 1,
620 IF ( m-q .GT. i )
THEN
621 CALL
slarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1,
624 CALL
slarfgp( m-q-i+1, x12(i,i), x12(i,i), 1,
630 CALL
slarf(
'L', q-i, p-i, x11(i+1,i), 1, tauq1(i),
631 $ x11(i+1,i+1), ldx11, work )
632 CALL
slarf(
'L', q-i, m-p-i, x11(i+1,i), 1, tauq1(i),
633 $ x21(i+1,i+1), ldx21, work )
635 CALL
slarf(
'L', m-q-i+1, p-i, x12(i,i), 1, tauq2(i),
636 $ x12(i,i+1), ldx12, work )
637 IF ( m-p-i .GT. 0 )
THEN
638 CALL
slarf(
'L', m-q-i+1, m-p-i, x12(i,i), 1, tauq2(i),
639 $ x22(i,i+1), ldx22, work )
648 CALL
sscal( m-q-i+1, -z1*z4, x12(i,i), 1 )
649 CALL
slarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1, tauq2(i) )
653 CALL
slarf(
'L', m-q-i+1, p-i, x12(i,i), 1, tauq2(i),
654 $ x12(i,i+1), ldx12, work )
657 $ CALL
slarf(
'L', m-q-i+1, m-p-q, x12(i,i), 1, tauq2(i),
658 $ x22(i,q+1), ldx22, work )
666 CALL
sscal( m-p-q-i+1, z2*z4, x22(p+i,q+i), 1 )
667 IF ( m-p-q .EQ. i )
THEN
668 CALL
slarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i,q+i), 1,
672 CALL
slarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i+1,q+i), 1,
675 CALL
slarf(
'L', m-p-q-i+1, m-p-q-i, x22(p+i,q+i), 1,
676 $ tauq2(p+i), x22(p+i,q+i+1), ldx22, work )
subroutine slarf(SIDE, M, N, V, INCV, TAU, C, LDC, WORK)
SLARF applies an elementary reflector to a general rectangular matrix.
subroutine slarfgp(N, ALPHA, X, INCX, TAU)
SLARFGP generates an elementary reflector (Householder matrix) with non-negatibe beta.
subroutine xerbla(SRNAME, INFO)
XERBLA
logical function lsame(CA, CB)
LSAME
subroutine saxpy(N, SA, SX, INCX, SY, INCY)
SAXPY
real function snrm2(N, X, INCX)
SNRM2
subroutine sorbdb(TRANS, SIGNS, M, P, Q, X11, LDX11, X12, LDX12, X21, LDX21, X22, LDX22, THETA, PHI, TAUP1, TAUP2, TAUQ1, TAUQ2, WORK, LWORK, INFO)
SORBDB
subroutine sscal(N, SA, SX, INCX)
SSCAL