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

◆ sorbdb()

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

SORBDB

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

Purpose:
!>
!> SORBDB 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 SORCSD
!> 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
!>                       convention);
!>          otherwise:  The upper-right block is made nonpositive (the
!>                       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 REAL 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 REAL 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 REAL 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 REAL 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 REAL 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 REAL 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 REAL array, dimension (P)
!>          The scalar factors of the elementary reflectors that define
!>          P1.
!>
[out]TAUP2
!>          TAUP2 is REAL array, dimension (M-P)
!>          The scalar factors of the elementary reflectors that define
!>          P2.
!>
[out]TAUQ1
!>          TAUQ1 is REAL array, dimension (Q)
!>          The scalar factors of the elementary reflectors that define
!>          Q1.
!>
[out]TAUQ2
!>          TAUQ2 is REAL array, dimension (M-Q)
!>          The scalar factors of the elementary reflectors that define
!>          Q2.
!>
[out]WORK
!>          WORK is REAL 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 SORCSD for details.
!>
!>  P1, P2, Q1, and Q2 are represented as products of elementary
!>  reflectors. See SORCSD for details on generating P1, P2, Q1, and Q2
!>  using SORGQR and SORGLQ.
!>
References:
[1] Brian D. Sutton. Computing the complete CS decomposition. Numer. Algorithms, 50(1):33-65, 2009.

Definition at line 282 of file sorbdb.f.

286*
287* -- LAPACK computational routine --
288* -- LAPACK is a software package provided by Univ. of Tennessee, --
289* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
290*
291* .. Scalar Arguments ..
292 CHARACTER SIGNS, TRANS
293 INTEGER INFO, LDX11, LDX12, LDX21, LDX22, LWORK, M, P,
294 $ Q
295* ..
296* .. Array Arguments ..
297 REAL PHI( * ), THETA( * )
298 REAL TAUP1( * ), TAUP2( * ), TAUQ1( * ), TAUQ2( * ),
299 $ WORK( * ), X11( LDX11, * ), X12( LDX12, * ),
300 $ X21( LDX21, * ), X22( LDX22, * )
301* ..
302*
303* ====================================================================
304*
305* .. Parameters ..
306 REAL REALONE
307 parameter( realone = 1.0e0 )
308* ..
309* .. Local Scalars ..
310 LOGICAL COLMAJOR, LQUERY
311 INTEGER I, LWORKMIN, LWORKOPT
312 REAL Z1, Z2, Z3, Z4
313* ..
314* .. External Subroutines ..
315 EXTERNAL saxpy, slarf1f, slarfgp, sscal,
316 $ xerbla
317* ..
318* .. External Functions ..
319 REAL SNRM2
320 LOGICAL LSAME
321 EXTERNAL snrm2, lsame
322* ..
323* .. Intrinsic Functions
324 INTRINSIC atan2, cos, max, sin
325* ..
326* .. Executable Statements ..
327*
328* Test input arguments
329*
330 info = 0
331 colmajor = .NOT. lsame( trans, 'T' )
332 IF( .NOT. lsame( signs, 'O' ) ) THEN
333 z1 = realone
334 z2 = realone
335 z3 = realone
336 z4 = realone
337 ELSE
338 z1 = realone
339 z2 = -realone
340 z3 = realone
341 z4 = -realone
342 END IF
343 lquery = lwork .EQ. -1
344*
345 IF( m .LT. 0 ) THEN
346 info = -3
347 ELSE IF( p .LT. 0 .OR. p .GT. m ) THEN
348 info = -4
349 ELSE IF( q .LT. 0 .OR. q .GT. p .OR. q .GT. m-p .OR.
350 $ q .GT. m-q ) THEN
351 info = -5
352 ELSE IF( colmajor .AND. ldx11 .LT. max( 1, p ) ) THEN
353 info = -7
354 ELSE IF( .NOT.colmajor .AND. ldx11 .LT. max( 1, q ) ) THEN
355 info = -7
356 ELSE IF( colmajor .AND. ldx12 .LT. max( 1, p ) ) THEN
357 info = -9
358 ELSE IF( .NOT.colmajor .AND. ldx12 .LT. max( 1, m-q ) ) THEN
359 info = -9
360 ELSE IF( colmajor .AND. ldx21 .LT. max( 1, m-p ) ) THEN
361 info = -11
362 ELSE IF( .NOT.colmajor .AND. ldx21 .LT. max( 1, q ) ) THEN
363 info = -11
364 ELSE IF( colmajor .AND. ldx22 .LT. max( 1, m-p ) ) THEN
365 info = -13
366 ELSE IF( .NOT.colmajor .AND. ldx22 .LT. max( 1, m-q ) ) THEN
367 info = -13
368 END IF
369*
370* Compute workspace
371*
372 IF( info .EQ. 0 ) THEN
373 lworkopt = m - q
374 lworkmin = m - q
375 work(1) = real( lworkopt )
376 IF( lwork .LT. lworkmin .AND. .NOT. lquery ) THEN
377 info = -21
378 END IF
379 END IF
380 IF( info .NE. 0 ) THEN
381 CALL xerbla( 'xORBDB', -info )
382 RETURN
383 ELSE IF( lquery ) THEN
384 RETURN
385 END IF
386*
387* Handle column-major and row-major separately
388*
389 IF( colmajor ) THEN
390*
391* Reduce columns 1, ..., Q of X11, X12, X21, and X22
392*
393 DO i = 1, q
394*
395 IF( i .EQ. 1 ) THEN
396 CALL sscal( p-i+1, z1, x11(i,i), 1 )
397 ELSE
398 CALL sscal( p-i+1, z1*cos(phi(i-1)), x11(i,i), 1 )
399 CALL saxpy( p-i+1, -z1*z3*z4*sin(phi(i-1)), x12(i,
400 $ i-1),
401 $ 1, x11(i,i), 1 )
402 END IF
403 IF( i .EQ. 1 ) THEN
404 CALL sscal( m-p-i+1, z2, x21(i,i), 1 )
405 ELSE
406 CALL sscal( m-p-i+1, z2*cos(phi(i-1)), x21(i,i), 1 )
407 CALL saxpy( m-p-i+1, -z2*z3*z4*sin(phi(i-1)), x22(i,
408 $ i-1),
409 $ 1, x21(i,i), 1 )
410 END IF
411*
412 theta(i) = atan2( snrm2( m-p-i+1, x21(i,i), 1 ),
413 $ snrm2( p-i+1, x11(i,i), 1 ) )
414*
415 IF( p .GT. i ) THEN
416 CALL slarfgp( p-i+1, x11(i,i), x11(i+1,i), 1,
417 $ taup1(i) )
418 ELSE IF( p .EQ. i ) THEN
419 CALL slarfgp( p-i+1, x11(i,i), x11(i,i), 1, taup1(i) )
420 END IF
421 IF ( m-p .GT. i ) THEN
422 CALL slarfgp( 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 slarfgp( m-p-i+1, x21(i,i), x21(i,i), 1,
426 $ taup2(i) )
427 END IF
428*
429 IF ( q .GT. i ) THEN
430 CALL slarf1f( '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 slarf1f( 'L', p-i+1, m-q-i+1, x11(i,i), 1,
435 $ taup1(i), x12(i,i), ldx12, work )
436 END IF
437 IF ( q .GT. i ) THEN
438 CALL slarf1f( 'L', m-p-i+1, q-i, x21(i,i), 1,
439 $ taup2(i), x21(i,i+1), ldx21, work )
440 END IF
441 IF ( m-q+1 .GT. i ) THEN
442 CALL slarf1f( 'L', m-p-i+1, m-q-i+1, x21(i,i), 1,
443 $ taup2(i), x22(i,i), ldx22, work )
444 END IF
445*
446 IF( i .LT. q ) THEN
447 CALL sscal( q-i, -z1*z3*sin(theta(i)), x11(i,i+1),
448 $ ldx11 )
449 CALL saxpy( q-i, z2*z3*cos(theta(i)), x21(i,i+1),
450 $ ldx21,
451 $ x11(i,i+1), ldx11 )
452 END IF
453 CALL sscal( m-q-i+1, -z1*z4*sin(theta(i)), x12(i,i),
454 $ ldx12 )
455 CALL saxpy( m-q-i+1, z2*z4*cos(theta(i)), x22(i,i),
456 $ ldx22,
457 $ x12(i,i), ldx12 )
458*
459 IF( i .LT. q )
460 $ phi(i) = atan2( snrm2( q-i, x11(i,i+1), ldx11 ),
461 $ snrm2( m-q-i+1, x12(i,i), ldx12 ) )
462*
463 IF( i .LT. q ) THEN
464 IF ( q-i .EQ. 1 ) THEN
465 CALL slarfgp( q-i, x11(i,i+1), x11(i,i+1), ldx11,
466 $ tauq1(i) )
467 ELSE
468 CALL slarfgp( q-i, x11(i,i+1), x11(i,i+2), ldx11,
469 $ tauq1(i) )
470 END IF
471 END IF
472 IF ( q+i-1 .LT. m ) THEN
473 IF ( m-q .EQ. i ) THEN
474 CALL slarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
475 $ tauq2(i) )
476 ELSE
477 CALL slarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
478 $ tauq2(i) )
479 END IF
480 END IF
481*
482 IF( i .LT. q ) THEN
483 CALL slarf1f( 'R', p-i, q-i, x11(i,i+1), ldx11,
484 $ tauq1(i), x11(i+1,i+1), ldx11, work )
485 CALL slarf1f( 'R', m-p-i, q-i, x11(i,i+1), ldx11,
486 $ tauq1(i), x21(i+1,i+1), ldx21, work )
487 END IF
488 IF ( p .GT. i ) THEN
489 CALL slarf1f( 'R', p-i, m-q-i+1, x12(i,i), ldx12,
490 $ tauq2(i), x12(i+1,i), ldx12, work )
491 END IF
492 IF ( m-p .GT. i ) THEN
493 CALL slarf1f( 'R', m-p-i, m-q-i+1, x12(i,i), ldx12,
494 $ tauq2(i), x22(i+1,i), ldx22, work )
495 END IF
496*
497 END DO
498*
499* Reduce columns Q + 1, ..., P of X12, X22
500*
501 DO i = q + 1, p
502*
503 CALL sscal( m-q-i+1, -z1*z4, x12(i,i), ldx12 )
504 IF ( i .GE. m-q ) THEN
505 CALL slarfgp( m-q-i+1, x12(i,i), x12(i,i), ldx12,
506 $ tauq2(i) )
507 ELSE
508 CALL slarfgp( m-q-i+1, x12(i,i), x12(i,i+1), ldx12,
509 $ tauq2(i) )
510 END IF
511*
512 IF ( p .GT. i ) THEN
513 CALL slarf1f( 'R', p-i, m-q-i+1, x12(i,i), ldx12,
514 $ tauq2(i), x12(i+1,i), ldx12, work )
515 END IF
516 IF( m-p-q .GE. 1 )
517 $ CALL slarf1f( '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 sscal( m-p-q-i+1, z2*z4, x22(q+i,p+i), ldx22 )
527 IF ( i .EQ. m-p-q ) THEN
528 CALL slarfgp( m-p-q-i+1, x22(q+i,p+i), x22(q+i,p+i),
529 $ ldx22, tauq2(p+i) )
530 ELSE
531 CALL slarfgp( 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 IF ( i .LT. m-p-q ) THEN
535 CALL slarf1f( 'R', m-p-q-i, m-p-q-i+1, x22(q+i,p+i),
536 $ ldx22, tauq2(p+i), x22(q+i+1,p+i),
537 $ 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 sscal( p-i+1, z1, x11(i,i), ldx11 )
550 ELSE
551 CALL sscal( p-i+1, z1*cos(phi(i-1)), x11(i,i), ldx11 )
552 CALL saxpy( p-i+1, -z1*z3*z4*sin(phi(i-1)), x12(i-1,
553 $ i),
554 $ ldx12, x11(i,i), ldx11 )
555 END IF
556 IF( i .EQ. 1 ) THEN
557 CALL sscal( m-p-i+1, z2, x21(i,i), ldx21 )
558 ELSE
559 CALL sscal( m-p-i+1, z2*cos(phi(i-1)), x21(i,i),
560 $ ldx21 )
561 CALL saxpy( m-p-i+1, -z2*z3*z4*sin(phi(i-1)), x22(i-1,
562 $ i),
563 $ ldx22, x21(i,i), ldx21 )
564 END IF
565*
566 theta(i) = atan2( snrm2( m-p-i+1, x21(i,i), ldx21 ),
567 $ snrm2( p-i+1, x11(i,i), ldx11 ) )
568*
569 CALL slarfgp( p-i+1, x11(i,i), x11(i,i+1), ldx11,
570 $ taup1(i) )
571 IF ( i .EQ. m-p ) THEN
572 CALL slarfgp( m-p-i+1, x21(i,i), x21(i,i), ldx21,
573 $ taup2(i) )
574 ELSE
575 CALL slarfgp( m-p-i+1, x21(i,i), x21(i,i+1), ldx21,
576 $ taup2(i) )
577 END IF
578*
579 IF ( q .GT. i ) THEN
580 CALL slarf1f( 'R', q-i, p-i+1, x11(i,i), ldx11,
581 $ taup1(i), x11(i+1,i), ldx11, work )
582 END IF
583 IF ( m-q+1 .GT. i ) THEN
584 CALL slarf1f( 'R', m-q-i+1, p-i+1, x11(i,i), ldx11,
585 $ taup1(i), x12(i,i), ldx12, work )
586 END IF
587 IF ( q .GT. i ) THEN
588 CALL slarf1f( 'R', q-i, m-p-i+1, x21(i,i), ldx21,
589 $ taup2(i), x21(i+1,i), ldx21, work )
590 END IF
591 IF ( m-q+1 .GT. i ) THEN
592 CALL slarf1f( 'R', m-q-i+1, m-p-i+1, x21(i,i), ldx21,
593 $ taup2(i), x22(i,i), ldx22, work )
594 END IF
595*
596 IF( i .LT. q ) THEN
597 CALL sscal( q-i, -z1*z3*sin(theta(i)), x11(i+1,i), 1 )
598 CALL saxpy( q-i, z2*z3*cos(theta(i)), x21(i+1,i), 1,
599 $ x11(i+1,i), 1 )
600 END IF
601 CALL sscal( m-q-i+1, -z1*z4*sin(theta(i)), x12(i,i), 1 )
602 CALL saxpy( m-q-i+1, z2*z4*cos(theta(i)), x22(i,i), 1,
603 $ x12(i,i), 1 )
604*
605 IF( i .LT. q )
606 $ phi(i) = atan2( snrm2( q-i, x11(i+1,i), 1 ),
607 $ snrm2( m-q-i+1, x12(i,i), 1 ) )
608*
609 IF( i .LT. q ) THEN
610 IF ( q-i .EQ. 1) THEN
611 CALL slarfgp( q-i, x11(i+1,i), x11(i+1,i), 1,
612 $ tauq1(i) )
613 ELSE
614 CALL slarfgp( q-i, x11(i+1,i), x11(i+2,i), 1,
615 $ tauq1(i) )
616 END IF
617 END IF
618 IF ( m-q .GT. i ) THEN
619 CALL slarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1,
620 $ tauq2(i) )
621 ELSE
622 CALL slarfgp( m-q-i+1, x12(i,i), x12(i,i), 1,
623 $ tauq2(i) )
624 END IF
625*
626 IF( i .LT. q ) THEN
627 CALL slarf1f( 'L', q-i, p-i, x11(i+1,i), 1, tauq1(i),
628 $ x11(i+1,i+1), ldx11, work )
629 CALL slarf1f( 'L', q-i, m-p-i, x11(i+1,i), 1,
630 $ tauq1(i), x21(i+1,i+1), ldx21, work )
631 END IF
632 CALL slarf1f( '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 slarf1f( 'L', m-q-i+1, m-p-i, x12(i,i), 1,
636 $ tauq2(i), 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 sscal( m-q-i+1, -z1*z4, x12(i,i), 1 )
646 CALL slarfgp( m-q-i+1, x12(i,i), x12(i+1,i), 1,
647 $ tauq2(i) )
648*
649 IF ( p .GT. i ) THEN
650 CALL slarf1f( 'L', m-q-i+1, p-i, x12(i,i), 1,
651 $ tauq2(i), x12(i,i+1), ldx12, work )
652 END IF
653 IF( m-p-q .GE. 1 )
654 $ CALL slarf1f( 'L', m-q-i+1, m-p-q, x12(i,i), 1,
655 $ tauq2(i), 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 sscal( m-p-q-i+1, z2*z4, x22(p+i,q+i), 1 )
664 IF ( m-p-q .EQ. i ) THEN
665 CALL slarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i,q+i),
666 $ 1,
667 $ tauq2(p+i) )
668 ELSE
669 CALL slarfgp( m-p-q-i+1, x22(p+i,q+i), x22(p+i+1,q+i),
670 $ 1,
671 $ tauq2(p+i) )
672 CALL slarf1f( 'L', m-p-q-i+1, m-p-q-i, x22(p+i,q+i),
673 $ 1, tauq2(p+i), x22(p+i,q+i+1), ldx22,
674 $ work )
675 END IF
676*
677*
678 END DO
679*
680 END IF
681*
682 RETURN
683*
684* End of SORBDB
685*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
saxpy
subroutine saxpy(n, sa, sx, incx, sy, incy)
SAXPY
Definition saxpy.f:89
slarfgp
subroutine slarfgp(n, alpha, x, incx, tau)
SLARFGP generates an elementary reflector (Householder matrix) with non-negative beta.
Definition slarfgp.f:102
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
snrm2
real(wp) function snrm2(n, x, incx)
SNRM2
Definition snrm2.f90:89
sscal
subroutine sscal(n, sa, sx, incx)
SSCAL
Definition sscal.f:79
slarf1f
subroutine slarf1f(side, m, n, v, incv, tau, c, ldc, work)
SLARF1F applies an elementary reflector to a general rectangular
Definition slarf1f.f:123
Here is the call graph for this function:
Here is the caller graph for this function: