179 SUBROUTINE slalsd( UPLO, SMLSIZ, N, NRHS, D, E, B, LDB, RCOND,
180 $ rank, work, iwork, info )
189 INTEGER info, ldb, n, nrhs, rank, smlsiz
194 REAL b( ldb, * ), d( * ), e( * ), work( * )
201 parameter( zero = 0.0e0, one = 1.0e0, two = 2.0e0 )
204 INTEGER bx, bxst, c, difl, difr, givcol, givnum,
205 $ givptr, i, icmpq1, icmpq2, iwk,
j, k, nlvl,
206 $ nm1, nsize, nsub, nwork, perm, poles, s, sizei,
207 $ smlszp, sqre, st, st1, u, vt, z
208 REAL cs, eps, orgnrm, r, rcnd, sn, tol
220 INTRINSIC abs, int, log,
REAL, sign
230 ELSE IF( nrhs.LT.1 )
THEN
232 ELSE IF( ( ldb.LT.1 ) .OR. ( ldb.LT.n ) )
THEN
236 CALL
xerbla(
'SLALSD', -info )
244 IF( ( rcond.LE.zero ) .OR. ( rcond.GE.one ) )
THEN
256 ELSE IF( n.EQ.1 )
THEN
257 IF( d( 1 ).EQ.zero )
THEN
258 CALL
slaset(
'A', 1, nrhs, zero, zero,
b, ldb )
261 CALL
slascl(
'G', 0, 0, d( 1 ), one, 1, nrhs,
b, ldb, info )
262 d( 1 ) = abs( d( 1 ) )
269 IF( uplo.EQ.
'L' )
THEN
271 CALL
slartg( d( i ), e( i ), cs, sn, r )
274 d( i+1 ) = cs*d( i+1 )
276 CALL
srot( 1,
b( i, 1 ), 1,
b( i+1, 1 ), 1, cs, sn )
287 CALL
srot( 1,
b(
j, i ), 1,
b(
j+1, i ), 1, cs, sn )
296 orgnrm =
slanst(
'M', n, d, e )
297 IF( orgnrm.EQ.zero )
THEN
298 CALL
slaset(
'A', n, nrhs, zero, zero,
b, ldb )
302 CALL
slascl(
'G', 0, 0, orgnrm, one, n, 1, d, n, info )
303 CALL
slascl(
'G', 0, 0, orgnrm, one, nm1, 1, e, nm1, info )
308 IF( n.LE.smlsiz )
THEN
310 CALL
slaset(
'A', n, n, zero, one, work, n )
311 CALL
slasdq(
'U', 0, n, n, 0, nrhs, d, e, work, n, work, n,
b,
312 $ ldb, work( nwork ), info )
316 tol = rcnd*abs( d(
isamax( n, d, 1 ) ) )
318 IF( d( i ).LE.tol )
THEN
319 CALL
slaset(
'A', 1, nrhs, zero, zero,
b( i, 1 ), ldb )
321 CALL
slascl(
'G', 0, 0, d( i ), one, 1, nrhs,
b( i, 1 ),
326 CALL
sgemm(
'T',
'N', n, nrhs, n, one, work, n,
b, ldb, zero,
328 CALL
slacpy(
'A', n, nrhs, work( nwork ), n,
b, ldb )
332 CALL
slascl(
'G', 0, 0, one, orgnrm, n, 1, d, n, info )
333 CALL
slasrt(
'D', n, d, info )
334 CALL
slascl(
'G', 0, 0, orgnrm, one, n, nrhs,
b, ldb, info )
341 nlvl = int( log(
REAL( N ) /
REAL( SMLSIZ+1 ) ) / log( two ) ) + 1
353 givnum = poles + 2*nlvl*n
354 bx = givnum + 2*nlvl*n
361 givcol = perm + nlvl*n
362 iwk = givcol + nlvl*n*2
371 IF( abs( d( i ) ).LT.eps )
THEN
372 d( i ) = sign( eps, d( i ) )
377 IF( ( abs( e( i ) ).LT.eps ) .OR. ( i.EQ.nm1 ) )
THEN
389 iwork( sizei+nsub-1 ) = nsize
390 ELSE IF( abs( e( i ) ).GE.eps )
THEN
395 iwork( sizei+nsub-1 ) = nsize
403 iwork( sizei+nsub-1 ) = nsize
406 iwork( sizei+nsub-1 ) = 1
407 CALL
scopy( nrhs,
b( n, 1 ), ldb, work( bx+nm1 ), n )
410 IF( nsize.EQ.1 )
THEN
415 CALL
scopy( nrhs,
b( st, 1 ), ldb, work( bx+st1 ), n )
416 ELSE IF( nsize.LE.smlsiz )
THEN
420 CALL
slaset(
'A', nsize, nsize, zero, one,
421 $ work( vt+st1 ), n )
422 CALL
slasdq(
'U', 0, nsize, nsize, 0, nrhs, d( st ),
423 $ e( st ), work( vt+st1 ), n, work( nwork ),
424 $ n,
b( st, 1 ), ldb, work( nwork ), info )
428 CALL
slacpy(
'A', nsize, nrhs,
b( st, 1 ), ldb,
429 $ work( bx+st1 ), n )
434 CALL
slasda( icmpq1, smlsiz, nsize, sqre, d( st ),
435 $ e( st ), work( u+st1 ), n, work( vt+st1 ),
436 $ iwork( k+st1 ), work( difl+st1 ),
437 $ work( difr+st1 ), work( z+st1 ),
438 $ work( poles+st1 ), iwork( givptr+st1 ),
439 $ iwork( givcol+st1 ), n, iwork( perm+st1 ),
440 $ work( givnum+st1 ), work( c+st1 ),
441 $ work( s+st1 ), work( nwork ), iwork( iwk ),
447 CALL
slalsa( icmpq2, smlsiz, nsize, nrhs,
b( st, 1 ),
448 $ ldb, work( bxst ), n, work( u+st1 ), n,
449 $ work( vt+st1 ), iwork( k+st1 ),
450 $ work( difl+st1 ), work( difr+st1 ),
451 $ work( z+st1 ), work( poles+st1 ),
452 $ iwork( givptr+st1 ), iwork( givcol+st1 ), n,
453 $ iwork( perm+st1 ), work( givnum+st1 ),
454 $ work( c+st1 ), work( s+st1 ), work( nwork ),
455 $ iwork( iwk ), info )
466 tol = rcnd*abs( d(
isamax( n, d, 1 ) ) )
473 IF( abs( d( i ) ).LE.tol )
THEN
474 CALL
slaset(
'A', 1, nrhs, zero, zero, work( bx+i-1 ), n )
477 CALL
slascl(
'G', 0, 0, d( i ), one, 1, nrhs,
478 $ work( bx+i-1 ), n, info )
480 d( i ) = abs( d( i ) )
489 nsize = iwork( sizei+i-1 )
491 IF( nsize.EQ.1 )
THEN
492 CALL
scopy( nrhs, work( bxst ), n,
b( st, 1 ), ldb )
493 ELSE IF( nsize.LE.smlsiz )
THEN
494 CALL
sgemm(
'T',
'N', nsize, nrhs, nsize, one,
495 $ work( vt+st1 ), n, work( bxst ), n, zero,
498 CALL
slalsa( icmpq2, smlsiz, nsize, nrhs, work( bxst ), n,
499 $
b( st, 1 ), ldb, work( u+st1 ), n,
500 $ work( vt+st1 ), iwork( k+st1 ),
501 $ work( difl+st1 ), work( difr+st1 ),
502 $ work( z+st1 ), work( poles+st1 ),
503 $ iwork( givptr+st1 ), iwork( givcol+st1 ), n,
504 $ iwork( perm+st1 ), work( givnum+st1 ),
505 $ work( c+st1 ), work( s+st1 ), work( nwork ),
506 $ iwork( iwk ), info )
515 CALL
slascl(
'G', 0, 0, one, orgnrm, n, 1, d, n, info )
516 CALL
slasrt(
'D', n, d, info )
517 CALL
slascl(
'G', 0, 0, orgnrm, one, n, nrhs,
b, ldb, info )
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
INTEGER function isamax(N, SX, INCX)
ISAMAX
subroutine slalsd(UPLO, SMLSIZ, N, NRHS, D, E, B, LDB, RCOND, RANK, WORK, IWORK, INFO)
SLALSD uses the singular value decomposition of A to solve the least squares problem.
REAL function slamch(CMACH)
SLAMCH
subroutine slalsa(ICOMPQ, SMLSIZ, N, NRHS, B, LDB, BX, LDBX, U, LDU, VT, K, DIFL, DIFR, Z, POLES, GIVPTR, GIVCOL, LDGCOL, PERM, GIVNUM, C, S, WORK, IWORK, INFO)
SLALSA computes the SVD of the coefficient matrix in compact form. Used by sgelsd.
REAL function slanst(NORM, N, D, E)
SLANST returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real symmetric tridiagonal matrix.
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine slasdq(UPLO, SQRE, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, LDU, C, LDC, WORK, INFO)
SLASDQ computes the SVD of a real bidiagonal matrix with diagonal d and off-diagonal e...
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real b(3) integer i
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine slasda(ICOMPQ, SMLSIZ, N, SQRE, D, E, U, LDU, VT, K, DIFL, DIFR, Z, POLES, GIVPTR, GIVCOL, LDGCOL, PERM, GIVNUM, C, S, WORK, IWORK, INFO)
SLASDA computes the singular value decomposition (SVD) of a real upper bidiagonal matrix with diagona...
subroutine slasrt(ID, N, D, INFO)
SLASRT sorts numbers in increasing or decreasing order.
subroutine slartg(F, G, CS, SN, R)
SLARTG generates a plane rotation with real cosine and real sine.
subroutine srot(N, SX, INCX, SY, INCY, C, S)
SROT