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
315 IF( ( icompq.LT.0 ) .OR. ( icompq.GT.1 ) )
THEN
317 ELSE IF( nl.LT.1 )
THEN
319 ELSE IF( nr.LT.1 )
THEN
321 ELSE IF( ( sqre.LT.0 ) .OR. ( sqre.GT.1 ) )
THEN
329 ELSE IF( ldb.LT.n )
THEN
331 ELSE IF( ldbx.LT.n )
THEN
333 ELSE IF( givptr.LT.0 )
THEN
335 ELSE IF( ldgcol.LT.n )
THEN
337 ELSE IF( ldgnum.LT.n )
THEN
339 ELSE IF( k.LT.1 )
THEN
343 CALL xerbla(
'DLALS0', -info )
350 IF( icompq.EQ.0 )
THEN
357 CALL drot( nrhs, b( givcol( i, 2 ), 1 ), ldb,
358 $ b( givcol( i, 1 ), 1 ), ldb, givnum( i, 2 ),
364 CALL dcopy( nrhs, b( nlp1, 1 ), ldb, bx( 1, 1 ), ldbx )
366 CALL dcopy( nrhs, b( perm( i ), 1 ), ldb, bx( i, 1 ), ldbx )
373 CALL dcopy( nrhs, bx, ldbx, b, ldb )
374 IF( z( 1 ).LT.zero )
THEN
375 CALL dscal( nrhs, negone, b, ldb )
381 dsigj = -poles( j, 2 )
383 difrj = -difr( j, 1 )
384 dsigjp = -poles( j+1, 2 )
386 IF( ( z( j ).EQ.zero ) .OR. ( poles( j, 2 ).EQ.zero ) )
390 work( j ) = -poles( j, 2 )*z( j ) / diflj /
391 $ ( poles( j, 2 )+dj )
394 IF( ( z( i ).EQ.zero ) .OR.
395 $ ( poles( i, 2 ).EQ.zero ) )
THEN
398 work( i ) = poles( i, 2 )*z( i ) /
399 $ ( dlamc3( poles( i, 2 ), dsigj )-
400 $ diflj ) / ( poles( i, 2 )+dj )
404 IF( ( z( i ).EQ.zero ) .OR.
405 $ ( poles( i, 2 ).EQ.zero ) )
THEN
408 work( i ) = poles( i, 2 )*z( i ) /
409 $ ( dlamc3( poles( i, 2 ), dsigjp )+
410 $ difrj ) / ( poles( i, 2 )+dj )
414 temp = dnrm2( k, work, 1 )
415 CALL dgemv(
'T', k, nrhs, one, bx, ldbx, work, 1, zero,
417 CALL dlascl(
'G', 0, 0, temp, one, 1, nrhs, b( j, 1 ),
424 IF( k.LT.max( m, n ) )
425 $
CALL dlacpy(
'A', n-k, nrhs, bx( k+1, 1 ), ldbx,
435 CALL dcopy( nrhs, b, ldb, bx, ldbx )
438 dsigj = poles( j, 2 )
439 IF( z( j ).EQ.zero )
THEN
442 work( j ) = -z( j ) / difl( j ) /
443 $ ( dsigj+poles( j, 1 ) ) / difr( j, 2 )
446 IF( z( j ).EQ.zero )
THEN
449 work( i ) = z( j ) / ( dlamc3( dsigj, -poles( i+1,
450 $ 2 ) )-difr( i, 1 ) ) /
451 $ ( dsigj+poles( i, 1 ) ) / difr( i, 2 )
455 IF( z( j ).EQ.zero )
THEN
458 work( i ) = z( j ) / ( dlamc3( dsigj, -poles( i,
459 $ 2 ) )-difl( i ) ) /
460 $ ( dsigj+poles( i, 1 ) ) / difr( i, 2 )
463 CALL dgemv(
'T', k, nrhs, one, b, ldb, work, 1, zero,
472 CALL dcopy( nrhs, b( m, 1 ), ldb, bx( m, 1 ), ldbx )
473 CALL drot( nrhs, bx( 1, 1 ), ldbx, bx( m, 1 ), ldbx, c, s )
475 IF( k.LT.max( m, n ) )
476 $
CALL dlacpy(
'A', n-k, nrhs, b( k+1, 1 ), ldb, bx( k+1, 1 ),
481 CALL dcopy( nrhs, bx( 1, 1 ), ldbx, b( nlp1, 1 ), ldb )
483 CALL dcopy( nrhs, bx( m, 1 ), ldbx, b( m, 1 ), ldb )
486 CALL dcopy( nrhs, bx( i, 1 ), ldbx, b( perm( i ), 1 ), ldb )
491 DO 100 i = givptr, 1, -1
492 CALL drot( nrhs, b( givcol( i, 2 ), 1 ), ldb,
493 $ b( givcol( i, 1 ), 1 ), ldb, givnum( i, 2 ),
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
subroutine xerbla(SRNAME, INFO)
XERBLA
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 dscal(N, DA, DX, INCX)
DSCAL
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 drot(N, DX, INCX, DY, INCY, C, S)
DROT
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV