141 SUBROUTINE dpstf2( UPLO, N, A, LDA, PIV, RANK, TOL, WORK, INFO )
150 INTEGER info, lda, n, rank
154 DOUBLE PRECISION a( lda, * ), work( 2*n )
161 DOUBLE PRECISION one, zero
162 parameter( one = 1.0d+0, zero = 0.0d+0 )
165 DOUBLE PRECISION ajj, dstop, dtemp
166 INTEGER i, itemp,
j, pvt
178 INTRINSIC max, sqrt, maxloc
185 upper =
lsame( uplo,
'U' )
186 IF( .NOT.upper .AND. .NOT.
lsame( uplo,
'L' ) )
THEN
188 ELSE IF( n.LT.0 )
THEN
190 ELSE IF( lda.LT.max( 1, n ) )
THEN
194 CALL
xerbla(
'DPSTF2', -info )
214 IF( a( i, i ).GT.ajj )
THEN
219 IF( ajj.EQ.zero.OR.
disnan( ajj ) )
THEN
227 IF( tol.LT.zero )
THEN
228 dstop = n *
dlamch(
'Epsilon' ) * ajj
252 work( i ) = work( i ) + a(
j-1, i )**2
254 work( n+i ) = a( i, i ) - work( i )
259 itemp = maxloc( work( (n+
j):(2*n) ), 1 )
262 IF( ajj.LE.dstop.OR.
disnan( ajj ) )
THEN
272 a( pvt, pvt ) = a(
j,
j )
273 CALL
dswap(
j-1, a( 1,
j ), 1, a( 1, pvt ), 1 )
275 $ CALL
dswap( n-pvt, a(
j, pvt+1 ), lda,
276 $ a( pvt, pvt+1 ), lda )
277 CALL
dswap( pvt-
j-1, a(
j,
j+1 ), lda, a(
j+1, pvt ), 1 )
282 work(
j ) = work( pvt )
285 piv( pvt ) = piv(
j )
295 CALL
dgemv(
'Trans',
j-1, n-
j, -one, a( 1,
j+1 ), lda,
296 $ a( 1,
j ), 1, one, a(
j,
j+1 ), lda )
297 CALL
dscal( n-
j, one / ajj, a(
j,
j+1 ), lda )
315 work( i ) = work( i ) + a( i,
j-1 )**2
317 work( n+i ) = a( i, i ) - work( i )
322 itemp = maxloc( work( (n+
j):(2*n) ), 1 )
325 IF( ajj.LE.dstop.OR.
disnan( ajj ) )
THEN
335 a( pvt, pvt ) = a(
j,
j )
336 CALL
dswap(
j-1, a(
j, 1 ), lda, a( pvt, 1 ), lda )
338 $ CALL
dswap( n-pvt, a( pvt+1,
j ), 1, a( pvt+1, pvt ),
340 CALL
dswap( pvt-
j-1, a(
j+1,
j ), 1, a( pvt,
j+1 ), lda )
345 work(
j ) = work( pvt )
348 piv( pvt ) = piv(
j )
358 CALL
dgemv(
'No Trans', n-
j,
j-1, -one, a(
j+1, 1 ), lda,
359 $ a(
j, 1 ), lda, one, a(
j+1,
j ), 1 )
360 CALL
dscal( n-
j, one / ajj, a(
j+1,
j ), 1 )