263 SUBROUTINE zlalsa( ICOMPQ, SMLSIZ, N, NRHS, B, LDB, BX, LDBX, U,
264 $ LDU, VT, K, DIFL, DIFR, Z, POLES, GIVPTR,
265 $ GIVCOL, LDGCOL, PERM, GIVNUM, C, S, RWORK,
273 INTEGER ICOMPQ, INFO, LDB, LDBX, LDGCOL, LDU, N, NRHS,
277 INTEGER GIVCOL( LDGCOL, * ), GIVPTR( * ), IWORK( * ),
278 $ K( * ), PERM( LDGCOL, * )
279 DOUBLE PRECISION C( * ), DIFL( LDU, * ), DIFR( LDU, * ),
280 $ givnum( ldu, * ), poles( ldu, * ), rwork( * ),
281 $ s( * ), u( ldu, * ), vt( ldu, * ), z( ldu, * )
282 COMPLEX*16 B( LDB, * ), BX( LDBX, * )
288 DOUBLE PRECISION ZERO, ONE
289 PARAMETER ( ZERO = 0.0d0, one = 1.0d0 )
292 INTEGER I, I1, IC, IM1, INODE, J, JCOL, JIMAG, JREAL,
293 $ JROW, LF, LL, LVL, LVL2, ND, NDB1, NDIML,
294 $ NDIMR, NL, NLF, NLP1, NLVL, NR, NRF, NRP1, SQRE
300 INTRINSIC dble, dcmplx, dimag
308 IF( ( icompq.LT.0 ) .OR. ( icompq.GT.1 ) )
THEN
310 ELSE IF( smlsiz.LT.3 )
THEN
312 ELSE IF( n.LT.smlsiz )
THEN
314 ELSE IF( nrhs.LT.1 )
THEN
316 ELSE IF( ldb.LT.n )
THEN
318 ELSE IF( ldbx.LT.n )
THEN
320 ELSE IF( ldu.LT.n )
THEN
322 ELSE IF( ldgcol.LT.n )
THEN
326 CALL xerbla(
'ZLALSA', -info )
336 CALL dlasdt( n, nlvl, nd, iwork( inode ), iwork( ndiml ),
337 $ iwork( ndimr ), smlsiz )
342 IF( icompq.EQ.1 )
THEN
361 ic = iwork( inode+i1 )
362 nl = iwork( ndiml+i1 )
363 nr = iwork( ndimr+i1 )
375 DO 10 jrow = nlf, nlf + nl - 1
377 rwork( j ) = dble( b( jrow, jcol ) )
380 CALL dgemm(
'T',
'N', nl, nrhs, nl, one, u( nlf, 1 ), ldu,
381 $ rwork( 1+nl*nrhs*2 ), nl, zero, rwork( 1 ), nl )
384 DO 30 jrow = nlf, nlf + nl - 1
386 rwork( j ) = dimag( b( jrow, jcol ) )
389 CALL dgemm(
'T',
'N', nl, nrhs, nl, one, u( nlf, 1 ), ldu,
390 $ rwork( 1+nl*nrhs*2 ), nl, zero, rwork( 1+nl*nrhs ),
395 DO 50 jrow = nlf, nlf + nl - 1
398 bx( jrow, jcol ) = dcmplx( rwork( jreal ),
411 DO 70 jrow = nrf, nrf + nr - 1
413 rwork( j ) = dble( b( jrow, jcol ) )
416 CALL dgemm(
'T',
'N', nr, nrhs, nr, one, u( nrf, 1 ), ldu,
417 $ rwork( 1+nr*nrhs*2 ), nr, zero, rwork( 1 ), nr )
419 DO 100 jcol = 1, nrhs
420 DO 90 jrow = nrf, nrf + nr - 1
422 rwork( j ) = dimag( b( jrow, jcol ) )
425 CALL dgemm(
'T',
'N', nr, nrhs, nr, one, u( nrf, 1 ), ldu,
426 $ rwork( 1+nr*nrhs*2 ), nr, zero, rwork( 1+nr*nrhs ),
430 DO 120 jcol = 1, nrhs
431 DO 110 jrow = nrf, nrf + nr - 1
434 bx( jrow, jcol ) = dcmplx( rwork( jreal ),
445 ic = iwork( inode+i-1 )
446 CALL zcopy( nrhs, b( ic, 1 ), ldb, bx( ic, 1 ), ldbx )
455 DO 160 lvl = nlvl, 1, -1
470 ic = iwork( inode+im1 )
471 nl = iwork( ndiml+im1 )
472 nr = iwork( ndimr+im1 )
476 CALL zlals0( icompq, nl, nr, sqre, nrhs, bx( nlf, 1 ), ldbx,
477 $ b( nlf, 1 ), ldb, perm( nlf, lvl ),
478 $ givptr( j ), givcol( nlf, lvl2 ), ldgcol,
479 $ givnum( nlf, lvl2 ), ldu, poles( nlf, lvl2 ),
480 $ difl( nlf, lvl ), difr( nlf, lvl2 ),
481 $ z( nlf, lvl ), k( j ), c( j ), s( j ), rwork,
508 DO 180 i = ll, lf, -1
510 ic = iwork( inode+im1 )
511 nl = iwork( ndiml+im1 )
512 nr = iwork( ndimr+im1 )
521 CALL zlals0( icompq, nl, nr, sqre, nrhs, b( nlf, 1 ), ldb,
522 $ bx( nlf, 1 ), ldbx, perm( nlf, lvl ),
523 $ givptr( j ), givcol( nlf, lvl2 ), ldgcol,
524 $ givnum( nlf, lvl2 ), ldu, poles( nlf, lvl2 ),
525 $ difl( nlf, lvl ), difr( nlf, lvl2 ),
526 $ z( nlf, lvl ), k( j ), c( j ), s( j ), rwork,
538 ic = iwork( inode+i1 )
539 nl = iwork( ndiml+i1 )
540 nr = iwork( ndimr+i1 )
557 DO 210 jcol = 1, nrhs
558 DO 200 jrow = nlf, nlf + nlp1 - 1
560 rwork( j ) = dble( b( jrow, jcol ) )
563 CALL dgemm(
'T',
'N', nlp1, nrhs, nlp1, one, vt( nlf, 1 ), ldu,
564 $ rwork( 1+nlp1*nrhs*2 ), nlp1, zero, rwork( 1 ),
567 DO 230 jcol = 1, nrhs
568 DO 220 jrow = nlf, nlf + nlp1 - 1
570 rwork( j ) = dimag( b( jrow, jcol ) )
573 CALL dgemm(
'T',
'N', nlp1, nrhs, nlp1, one, vt( nlf, 1 ), ldu,
574 $ rwork( 1+nlp1*nrhs*2 ), nlp1, zero,
575 $ rwork( 1+nlp1*nrhs ), nlp1 )
578 DO 250 jcol = 1, nrhs
579 DO 240 jrow = nlf, nlf + nlp1 - 1
582 bx( jrow, jcol ) = dcmplx( rwork( jreal ),
594 DO 270 jcol = 1, nrhs
595 DO 260 jrow = nrf, nrf + nrp1 - 1
597 rwork( j ) = dble( b( jrow, jcol ) )
600 CALL dgemm(
'T',
'N', nrp1, nrhs, nrp1, one, vt( nrf, 1 ), ldu,
601 $ rwork( 1+nrp1*nrhs*2 ), nrp1, zero, rwork( 1 ),
604 DO 290 jcol = 1, nrhs
605 DO 280 jrow = nrf, nrf + nrp1 - 1
607 rwork( j ) = dimag( b( jrow, jcol ) )
610 CALL dgemm(
'T',
'N', nrp1, nrhs, nrp1, one, vt( nrf, 1 ), ldu,
611 $ rwork( 1+nrp1*nrhs*2 ), nrp1, zero,
612 $ rwork( 1+nrp1*nrhs ), nrp1 )
615 DO 310 jcol = 1, nrhs
616 DO 300 jrow = nrf, nrf + nrp1 - 1
619 bx( jrow, jcol ) = dcmplx( rwork( jreal ),
subroutine xerbla(srname, info)
subroutine zcopy(n, zx, incx, zy, incy)
ZCOPY
subroutine dgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
DGEMM
subroutine zlals0(icompq, nl, nr, sqre, nrhs, b, ldb, bx, ldbx, perm, givptr, givcol, ldgcol, givnum, ldgnum, poles, difl, difr, z, k, c, s, rwork, info)
ZLALS0 applies back multiplying factors in solving the least squares problem using divide and conquer...
subroutine zlalsa(icompq, smlsiz, n, nrhs, b, ldb, bx, ldbx, u, ldu, vt, k, difl, difr, z, poles, givptr, givcol, ldgcol, perm, givnum, c, s, rwork, iwork, info)
ZLALSA computes the SVD of the coefficient matrix in compact form. Used by sgelsd.
subroutine dlasdt(n, lvl, nd, inode, ndiml, ndimr, msub)
DLASDT creates a tree of subproblems for bidiagonal divide and conquer. Used by sbdsdc.