LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ dorbdb()

subroutine dorbdb ( character  trans,
character  signs,
integer  m,
integer  p,
integer  q,
double precision, dimension( ldx11, * )  x11,
integer  ldx11,
double precision, dimension( ldx12, * )  x12,
integer  ldx12,
double precision, dimension( ldx21, * )  x21,
integer  ldx21,
double precision, dimension( ldx22, * )  x22,
integer  ldx22,
double precision, dimension( * )  theta,
double precision, dimension( * )  phi,
double precision, dimension( * )  taup1,
double precision, dimension( * )  taup2,
double precision, dimension( * )  tauq1,
double precision, dimension( * )  tauq2,
double precision, dimension( * )  work,
integer  lwork,
integer  info 
)

DORBDB

Download DORBDB + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DORBDB simultaneously bidiagonalizes the blocks of an M-by-M
 partitioned orthogonal matrix X:

                                 [ B11 | B12 0  0 ]
     [ X11 | X12 ]   [ P1 |    ] [  0  |  0 -I  0 ] [ Q1 |    ]**T
 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 DORCSD
 for details.)

 The orthogonal 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 DOUBLE PRECISION array, dimension (LDX11,Q)
          On entry, the top-left block of the orthogonal 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 DOUBLE PRECISION array, dimension (LDX12,M-Q)
          On entry, the top-right block of the orthogonal 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 DOUBLE PRECISION array, dimension (LDX21,Q)
          On entry, the bottom-left block of the orthogonal 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 DOUBLE PRECISION array, dimension (LDX22,M-Q)
          On entry, the bottom-right block of the orthogonal 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 DOUBLE PRECISION array, dimension (P)
          The scalar factors of the elementary reflectors that define
          P1.
[out]TAUP2
          TAUP2 is DOUBLE PRECISION array, dimension (M-P)
          The scalar factors of the elementary reflectors that define
          P2.
[out]TAUQ1
          TAUQ1 is DOUBLE PRECISION array, dimension (Q)
          The scalar factors of the elementary reflectors that define
          Q1.
[out]TAUQ2
          TAUQ2 is DOUBLE PRECISION array, dimension (M-Q)
          The scalar factors of the elementary reflectors that define
          Q2.
[out]WORK
          WORK is DOUBLE PRECISION 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.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
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
  [1] or DORCSD for details.

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

Definition at line 284 of file dorbdb.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 DOUBLE PRECISION 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 DOUBLE PRECISION ONE
310 parameter( one = 1.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 daxpy, dlarf, dlarfgp, dscal, xerbla
319* ..
320* .. External Functions ..
321 DOUBLE PRECISION DNRM2
322 LOGICAL LSAME
323 EXTERNAL dnrm2, lsame
324* ..
325* .. Intrinsic Functions
326 INTRINSIC atan2, cos, max, sin
327* ..
328* .. Executable Statements ..
329*
330* Test input arguments
331*
332 info = 0
333 colmajor = .NOT. lsame( trans, 'T' )
334 IF( .NOT. lsame( signs, 'O' ) ) THEN
335 z1 = realone
336 z2 = realone
337 z3 = realone
338 z4 = realone
339 ELSE
340 z1 = realone
341 z2 = -realone
342 z3 = realone
343 z4 = -realone
344 END IF
345 lquery = lwork .EQ. -1
346*
347 IF( m .LT. 0 ) THEN
348 info = -3
349 ELSE IF( p .LT. 0 .OR. p .GT. m ) THEN
350 info = -4
351 ELSE IF( q .LT. 0 .OR. q .GT. p .OR. q .GT. m-p .OR.
352 $ q .GT. m-q ) THEN
353 info = -5
354 ELSE IF( colmajor .AND. ldx11 .LT. max( 1, p ) ) THEN
355 info = -7
356 ELSE IF( .NOT.colmajor .AND. ldx11 .LT. max( 1, q ) ) THEN
357 info = -7
358 ELSE IF( colmajor .AND. ldx12 .LT. max( 1, p ) ) THEN
359 info = -9
360 ELSE IF( .NOT.colmajor .AND. ldx12 .LT. max( 1, m-q ) ) THEN
361 info = -9
362 ELSE IF( colmajor .AND. ldx21 .LT. max( 1, m-p ) ) THEN
363 info = -11
364 ELSE IF( .NOT.colmajor .AND. ldx21 .LT. max( 1, q ) ) THEN
365 info = -11
366 ELSE IF( colmajor .AND. ldx22 .LT. max( 1, m-p ) ) THEN
367 info = -13
368 ELSE IF( .NOT.colmajor .AND. ldx22 .LT. max( 1, m-q ) ) THEN
369 info = -13
370 END IF
371*
372* Compute workspace
373*
374 IF( info .EQ. 0 ) THEN
375 lworkopt = m - q
376 lworkmin = m - q
377 work(1) = lworkopt
378 IF( lwork .LT. lworkmin .AND. .NOT. lquery ) THEN
379 info = -21
380 END IF
381 END IF
382 IF( info .NE. 0 ) THEN
383 CALL xerbla( 'xORBDB', -info )
384 RETURN
385 ELSE IF( lquery ) THEN
386 RETURN
387 END IF
388*
389* Handle column-major and row-major separately
390*
391 IF( colmajor ) THEN
392*
393* Reduce columns 1, ..., Q of X11, X12, X21, and X22
394*
395 DO i = 1, q
396*
397 IF( i .EQ. 1 ) THEN
398 CALL dscal( p-i+1, z1, x11(i,i), 1 )
399 ELSE
400 CALL dscal( p-i+1, z1*cos(phi(i-1)), x11(i,i), 1 )
401 CALL daxpy( p-i+1, -z1*z3*z4*sin(phi(i-1)), x12(i,i-1),
402 $ 1, x11(i,i), 1 )
403 END IF
404 IF( i .EQ. 1 ) THEN
405 CALL dscal( m-p-i+1, z2, x21(i,i), 1 )
406 ELSE
407 CALL dscal( m-p-i+1, z2*cos(phi(i-1)), x21(i,i), 1 )
408 CALL daxpy( m-p-i+1, -z2*z3*z4*sin(phi(i-1)), x22(i,i-1),
409 $ 1, x21(i,i), 1 )
410 END IF
411*
412 theta(i) = atan2( dnrm2( m-p-i+1, x21(i,i), 1 ),
413 $ dnrm2( p-i+1, x11(i,i), 1 ) )
414*
415 IF( p .GT. i ) THEN
416 CALL dlarfgp( p-i+1, x11(i,i), x11(i+1,i), 1, taup1(i) )
417 ELSE IF( p .EQ. i ) THEN
418 CALL dlarfgp( p-i+1, x11(i,i), x11(i,i), 1, taup1(i) )
419 END IF
420 x11(i,i) = one
421 IF ( m-p .GT. i ) THEN
422 CALL dlarfgp( m-p-i+1, x21(i,i), x21(i+1,i), 1,
423 $ taup2(i) )
424 ELSE IF ( m-p .EQ. i ) THEN
425 CALL dlarfgp( m-p-i+1, x21(i,i), x21(i,i), 1, taup2(i) )
426 END IF
427 x21(i,i) = one
428*
429 IF ( q .GT. i ) THEN
430 CALL dlarf( 'L', p-i+1, q-i, x11(i,i), 1, taup1(i),
431 $ x11(i,i+1), ldx11, work )
432 END IF
433 IF ( m-q+1 .GT. i ) THEN
434 CALL dlarf( 'L', p-i+1, m-q-i+1, x11(i,i), 1, taup1(i),
435 $ x12(i,i), ldx12, work )
436 END IF
437 IF ( q .GT. i ) THEN
438 CALL dlarf( 'L', m-p-i+1, q-i, x21(i,i), 1, taup2(i),
439 $ x21(i,i+1), ldx21, work )
440 END IF
441 IF ( m-q+1 .GT. i ) THEN
442 CALL dlarf( 'L', m-p-i+1, m-q-i+1, x21(i,i), 1, taup2(i),
443 $ x22(i,i), ldx22, work )
444 END IF
445*
446 IF( i .LT. q ) THEN
447 CALL dscal( q-i, -z1*z3*sin(theta(i)), x11(i,i+1),
448 $ ldx11 )
449 CALL daxpy( q-i, z2*z3*cos(theta(i)), x21(i,i+1), ldx21,
450 $ x11(i,i+1), ldx11 )
451 END IF
452 CALL dscal( m-q-i+1, -z1*z4*sin(theta(i)), x12(i,i), ldx12 )
453 CALL daxpy( m-q-i+1, z2*z4*cos(theta(i)), x22(i,i), ldx22,
454 $ x12(i,i), ldx12 )
455*
456 IF( i .LT. q )
457 $ phi(i) = atan2( dnrm2( q-i, x11(i,i+1), ldx11 ),
458 $ dnrm2( m-q-i+1, x12(i,i), ldx12 ) )
459*
460 IF( i .LT. q ) THEN
461 IF ( q-i .EQ. 1 ) THEN
462 CALL dlarfgp( q-i, x11(i,i+1), x11(i,i+1), ldx11,
463 $ tauq1(i) )
464 ELSE
465 CALL dlarfgp( q-i, x11(i,i+1), x11(i,i+2), ldx11,
466 $ tauq1(i) )
467 END IF
468 x11(i,i+1) = one
469 END IF
470 IF ( q+i-1 .LT. m ) THEN
471 IF ( m-q .EQ. i ) THEN
472 CALL dlarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
473 $ tauq2(i) )
474 ELSE
475 CALL dlarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
476 $ tauq2(i) )
477 END IF
478 END IF
479 x12(i,i) = one
480*
481 IF( i .LT. q ) THEN
482 CALL dlarf( 'R', p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
483 $ x11(i+1,i+1), ldx11, work )
484 CALL dlarf( 'R', m-p-i, q-i, x11(i,i+1), ldx11, tauq1(i),
485 $ x21(i+1,i+1), ldx21, work )
486 END IF
487 IF ( p .GT. i ) THEN
488 CALL dlarf( 'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
489 $ x12(i+1,i), ldx12, work )
490 END IF
491 IF ( m-p .GT. i ) THEN
492 CALL dlarf( 'R', m-p-i, m-q-i+1, x12(i,i), ldx12,
493 $ tauq2(i), x22(i+1,i), ldx22, work )
494 END IF
495*
496 END DO
497*
498* Reduce columns Q + 1, ..., P of X12, X22
499*
500 DO i = q + 1, p
501*
502 CALL dscal( m-q-i+1, -z1*z4, x12(i,i), ldx12 )
503 IF ( i .GE. m-q ) THEN
504 CALL dlarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
505 $ tauq2(i) )
506 ELSE
507 CALL dlarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
508 $ tauq2(i) )
509 END IF
510 x12(i,i) = one
511*
512 IF ( p .GT. i ) THEN
513 CALL dlarf( 'R', p-i, m-q-i+1, x12(i,i), ldx12, tauq2(i),
514 $ x12(i+1,i), ldx12, work )
515 END IF
516 IF( m-p-q .GE. 1 )
517 $ CALL dlarf( 'R', m-p-q, m-q-i+1, x12(i,i), ldx12,
518 $ tauq2(i), x22(q+1,i), ldx22, work )
519*
520 END DO
521*
522* Reduce columns P + 1, ..., M - Q of X12, X22
523*
524 DO i = 1, m - p - q
525*
526 CALL dscal( m-p-q-i+1, z2*z4, x22(q+i,p+i), ldx22 )
527 IF ( i .EQ. m-p-q ) THEN
528 CALL dlarfgp( m-p-q-i+1, x22(q+i,p+i), x22(q+i,p+i),
529 $ ldx22, tauq2(p+i) )
530 ELSE
531 CALL dlarfgp( m-p-q-i+1, x22(q+i,p+i), x22(q+i,p+i+1),
532 $ ldx22, tauq2(p+i) )
533 END IF
534 x22(q+i,p+i) = one
535 IF ( i .LT. m-p-q ) THEN
536 CALL dlarf( 'R', m-p-q-i, m-p-q-i+1, x22(q+i,p+i), ldx22,
537 $ tauq2(p+i), x22(q+i+1,p+i), ldx22, work )
538 END IF
539*
540 END DO
541*
542 ELSE
543*
544* Reduce columns 1, ..., Q of X11, X12, X21, X22
545*
546 DO i = 1, q
547*
548 IF( i .EQ. 1 ) THEN
549 CALL dscal( p-i+1, z1, x11(i,i), ldx11 )
550 ELSE
551 CALL dscal( p-i+1, z1*cos(phi(i-1)), x11(i,i), ldx11 )
552 CALL daxpy( p-i+1, -z1*z3*z4*sin(phi(i-1)), x12(i-1,i),
553 $ ldx12, x11(i,i), ldx11 )
554 END IF
555 IF( i .EQ. 1 ) THEN
556 CALL dscal( m-p-i+1, z2, x21(i,i), ldx21 )
557 ELSE
558 CALL dscal( m-p-i+1, z2*cos(phi(i-1)), x21(i,i), ldx21 )
559 CALL daxpy( m-p-i+1, -z2*z3*z4*sin(phi(i-1)), x22(i-1,i),
560 $ ldx22, x21(i,i), ldx21 )
561 END IF
562*
563 theta(i) = atan2( dnrm2( m-p-i+1, x21(i,i), ldx21 ),
564 $ dnrm2( p-i+1, x11(i,i), ldx11 ) )
565*
566 CALL dlarfgp( p-i+1, x11(i,i), x11(i,i+1), ldx11, taup1(i) )
567 x11(i,i) = one
568 IF ( i .EQ. m-p ) THEN
569 CALL dlarfgp( m-p-i+1, x21(i,i), x21(i,i), ldx21,
570 $ taup2(i) )
571 ELSE
572 CALL dlarfgp( m-p-i+1, x21(i,i), x21(i,i+1), ldx21,
573 $ taup2(i) )
574 END IF
575 x21(i,i) = one
576*
577 IF ( q .GT. i ) THEN
578 CALL dlarf( 'R', q-i, p-i+1, x11(i,i), ldx11, taup1(i),
579 $ x11(i+1,i), ldx11, work )
580 END IF
581 IF ( m-q+1 .GT. i ) THEN
582 CALL dlarf( 'R', m-q-i+1, p-i+1, x11(i,i), ldx11,
583 $ taup1(i), x12(i,i), ldx12, work )
584 END IF
585 IF ( q .GT. i ) THEN
586 CALL dlarf( 'R', q-i, m-p-i+1, x21(i,i), ldx21, taup2(i),
587 $ x21(i+1,i), ldx21, work )
588 END IF
589 IF ( m-q+1 .GT. i ) THEN
590 CALL dlarf( 'R', m-q-i+1, m-p-i+1, x21(i,i), ldx21,
591 $ taup2(i), x22(i,i), ldx22, work )
592 END IF
593*
594 IF( i .LT. q ) THEN
595 CALL dscal( q-i, -z1*z3*sin(theta(i)), x11(i+1,i), 1 )
596 CALL daxpy( q-i, z2*z3*cos(theta(i)), x21(i+1,i), 1,
597 $ x11(i+1,i), 1 )
598 END IF
599 CALL dscal( m-q-i+1, -z1*z4*sin(theta(i)), x12(i,i), 1 )
600 CALL daxpy( m-q-i+1, z2*z4*cos(theta(i)), x22(i,i), 1,
601 $ x12(i,i), 1 )
602*
603 IF( i .LT. q )
604 $ phi(i) = atan2( dnrm2( q-i, x11(i+1,i), 1 ),
605 $ dnrm2( m-q-i+1, x12(i,i), 1 ) )
606*
607 IF( i .LT. q ) THEN
608 IF ( q-i .EQ. 1) THEN
609 CALL dlarfgp( q-i, x11(i+1,i), x11(i+1,i), 1,
610 $ tauq1(i) )
611 ELSE
612 CALL dlarfgp( q-i, x11(i+1,i), x11(i+2,i), 1,
613 $ tauq1(i) )
614 END IF
615 x11(i+1,i) = one
616 END IF
617 IF ( m-q .GT. i ) THEN
618 CALL dlarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1,
619 $ tauq2(i) )
620 ELSE
621 CALL dlarfgp( m-q-i+1, x12(i,i), x12(i,i), 1,
622 $ tauq2(i) )
623 END IF
624 x12(i,i) = one
625*
626 IF( i .LT. q ) THEN
627 CALL dlarf( 'L', q-i, p-i, x11(i+1,i), 1, tauq1(i),
628 $ x11(i+1,i+1), ldx11, work )
629 CALL dlarf( 'L', q-i, m-p-i, x11(i+1,i), 1, tauq1(i),
630 $ x21(i+1,i+1), ldx21, work )
631 END IF
632 CALL dlarf( 'L', m-q-i+1, p-i, x12(i,i), 1, tauq2(i),
633 $ x12(i,i+1), ldx12, work )
634 IF ( m-p-i .GT. 0 ) THEN
635 CALL dlarf( 'L', m-q-i+1, m-p-i, x12(i,i), 1, tauq2(i),
636 $ x22(i,i+1), ldx22, work )
637 END IF
638*
639 END DO
640*
641* Reduce columns Q + 1, ..., P of X12, X22
642*
643 DO i = q + 1, p
644*
645 CALL dscal( m-q-i+1, -z1*z4, x12(i,i), 1 )
646 CALL dlarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1, tauq2(i) )
647 x12(i,i) = one
648*
649 IF ( p .GT. i ) THEN
650 CALL dlarf( 'L', m-q-i+1, p-i, x12(i,i), 1, tauq2(i),
651 $ x12(i,i+1), ldx12, work )
652 END IF
653 IF( m-p-q .GE. 1 )
654 $ CALL dlarf( 'L', m-q-i+1, m-p-q, x12(i,i), 1, tauq2(i),
655 $ x22(i,q+1), ldx22, work )
656*
657 END DO
658*
659* Reduce columns P + 1, ..., M - Q of X12, X22
660*
661 DO i = 1, m - p - q
662*
663 CALL dscal( m-p-q-i+1, z2*z4, x22(p+i,q+i), 1 )
664 IF ( m-p-q .EQ. i ) THEN
665 CALL dlarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i,q+i), 1,
666 $ tauq2(p+i) )
667 ELSE
668 CALL dlarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i+1,q+i), 1,
669 $ tauq2(p+i) )
670 CALL dlarf( 'L', m-p-q-i+1, m-p-q-i, x22(p+i,q+i), 1,
671 $ tauq2(p+i), x22(p+i,q+i+1), ldx22, work )
672 END IF
673 x22(p+i,q+i) = one
674*
675 END DO
676*
677 END IF
678*
679 RETURN
680*
681* End of DORBDB
682*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
daxpy
subroutine daxpy(n, da, dx, incx, dy, incy)
DAXPY
Definition daxpy.f:89
dlarf
subroutine dlarf(side, m, n, v, incv, tau, c, ldc, work)
DLARF applies an elementary reflector to a general rectangular matrix.
Definition dlarf.f:124
dlarfgp
subroutine dlarfgp(n, alpha, x, incx, tau)
DLARFGP generates an elementary reflector (Householder matrix) with non-negative beta.
Definition dlarfgp.f:104
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
dnrm2
real(wp) function dnrm2(n, x, incx)
DNRM2
Definition dnrm2.f90:89
dscal
subroutine dscal(n, da, dx, incx)
DSCAL
Definition dscal.f:79
Here is the call graph for this function:
Here is the caller graph for this function: