267 SUBROUTINE dlals0( ICOMPQ, NL, NR, SQRE, NRHS, B, LDB, BX, LDBX,
268 $ perm, givptr, givcol, ldgcol, givnum, ldgnum,
269 $ poles, difl, difr, z, k, c, s, work, info )
277 INTEGER GIVPTR, ICOMPQ, INFO, K, LDB, LDBX, LDGCOL,
278 $ ldgnum, nl, nr, nrhs, sqre
279 DOUBLE PRECISION C, S
282 INTEGER GIVCOL( ldgcol, * ), PERM( * )
283 DOUBLE PRECISION B( ldb, * ), BX( ldbx, * ), DIFL( * ),
284 $ difr( ldgnum, * ), givnum( ldgnum, * ),
285 $ poles( ldgnum, * ), work( * ), z( * )
291 DOUBLE PRECISION ONE, ZERO, NEGONE
292 parameter ( one = 1.0d0, zero = 0.0d0, negone = -1.0d0 )
295 INTEGER I, J, M, N, NLP1
296 DOUBLE PRECISION DIFLJ, DIFRJ, DJ, DSIGJ, DSIGJP, TEMP
303 DOUBLE PRECISION DLAMC3, DNRM2
304 EXTERNAL dlamc3, dnrm2
316 IF( ( icompq.LT.0 ) .OR. ( icompq.GT.1 ) )
THEN
318 ELSE IF( nl.LT.1 )
THEN
320 ELSE IF( nr.LT.1 )
THEN
322 ELSE IF( ( sqre.LT.0 ) .OR. ( sqre.GT.1 ) )
THEN
324 ELSE IF( nrhs.LT.1 )
THEN
326 ELSE IF( ldb.LT.n )
THEN
328 ELSE IF( ldbx.LT.n )
THEN
330 ELSE IF( givptr.LT.0 )
THEN
332 ELSE IF( ldgcol.LT.n )
THEN
334 ELSE IF( ldgnum.LT.n )
THEN
336 ELSE IF( k.LT.1 )
THEN
340 CALL xerbla(
'DLALS0', -info )
347 IF( icompq.EQ.0 )
THEN
354 CALL drot( nrhs, b( givcol( i, 2 ), 1 ), ldb,
355 $ b( givcol( i, 1 ), 1 ), ldb, givnum( i, 2 ),
361 CALL dcopy( nrhs, b( nlp1, 1 ), ldb, bx( 1, 1 ), ldbx )
363 CALL dcopy( nrhs, b( perm( i ), 1 ), ldb, bx( i, 1 ), ldbx )
370 CALL dcopy( nrhs, bx, ldbx, b, ldb )
371 IF( z( 1 ).LT.zero )
THEN
372 CALL dscal( nrhs, negone, b, ldb )
378 dsigj = -poles( j, 2 )
380 difrj = -difr( j, 1 )
381 dsigjp = -poles( j+1, 2 )
383 IF( ( z( j ).EQ.zero ) .OR. ( poles( j, 2 ).EQ.zero ) )
387 work( j ) = -poles( j, 2 )*z( j ) / diflj /
388 $ ( poles( j, 2 )+dj )
391 IF( ( z( i ).EQ.zero ) .OR.
392 $ ( poles( i, 2 ).EQ.zero ) )
THEN
395 work( i ) = poles( i, 2 )*z( i ) /
396 $ ( dlamc3( poles( i, 2 ), dsigj )-
397 $ diflj ) / ( poles( i, 2 )+dj )
401 IF( ( z( i ).EQ.zero ) .OR.
402 $ ( poles( i, 2 ).EQ.zero ) )
THEN
405 work( i ) = poles( i, 2 )*z( i ) /
406 $ ( dlamc3( poles( i, 2 ), dsigjp )+
407 $ difrj ) / ( poles( i, 2 )+dj )
411 temp = dnrm2( k, work, 1 )
412 CALL dgemv(
'T', k, nrhs, one, bx, ldbx, work, 1, zero,
414 CALL dlascl(
'G', 0, 0, temp, one, 1, nrhs, b( j, 1 ),
421 IF( k.LT.max( m, n ) )
422 $
CALL dlacpy(
'A', n-k, nrhs, bx( k+1, 1 ), ldbx,
432 CALL dcopy( nrhs, b, ldb, bx, ldbx )
435 dsigj = poles( j, 2 )
436 IF( z( j ).EQ.zero )
THEN
439 work( j ) = -z( j ) / difl( j ) /
440 $ ( dsigj+poles( j, 1 ) ) / difr( j, 2 )
443 IF( z( j ).EQ.zero )
THEN
446 work( i ) = z( j ) / ( dlamc3( dsigj, -poles( i+1,
447 $ 2 ) )-difr( i, 1 ) ) /
448 $ ( dsigj+poles( i, 1 ) ) / difr( i, 2 )
452 IF( z( j ).EQ.zero )
THEN
455 work( i ) = z( j ) / ( dlamc3( dsigj, -poles( i,
456 $ 2 ) )-difl( i ) ) /
457 $ ( dsigj+poles( i, 1 ) ) / difr( i, 2 )
460 CALL dgemv(
'T', k, nrhs, one, b, ldb, work, 1, zero,
469 CALL dcopy( nrhs, b( m, 1 ), ldb, bx( m, 1 ), ldbx )
470 CALL drot( nrhs, bx( 1, 1 ), ldbx, bx( m, 1 ), ldbx, c, s )
472 IF( k.LT.max( m, n ) )
473 $
CALL dlacpy(
'A', n-k, nrhs, b( k+1, 1 ), ldb, bx( k+1, 1 ),
478 CALL dcopy( nrhs, bx( 1, 1 ), ldbx, b( nlp1, 1 ), ldb )
480 CALL dcopy( nrhs, bx( m, 1 ), ldbx, b( m, 1 ), ldb )
483 CALL dcopy( nrhs, bx( i, 1 ), ldbx, b( perm( i ), 1 ), ldb )
488 DO 100 i = givptr, 1, -1
489 CALL drot( nrhs, b( givcol( i, 2 ), 1 ), ldb,
490 $ b( givcol( i, 1 ), 1 ), ldb, givnum( i, 2 ),
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
subroutine dlals0(ICOMPQ, NL, NR, SQRE, NRHS, B, LDB, BX, LDBX, PERM, GIVPTR, GIVCOL, LDGCOL, GIVNUM, LDGNUM, POLES, DIFL, DIFR, Z, K, C, S, WORK, INFO)
DLALS0 applies back multiplying factors in solving the least squares problem using divide and conquer...
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine drot(N, DX, INCX, DY, INCY, C, S)
DROT
subroutine dlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine dscal(N, DA, DX, INCX)
DSCAL