143 SUBROUTINE cpstrf( UPLO, N, A, LDA, PIV, RANK, TOL, WORK, INFO )
152 INTEGER INFO, LDA, N, RANK
165 parameter ( one = 1.0e+0, zero = 0.0e+0 )
167 parameter ( cone = ( 1.0e+0, 0.0e+0 ) )
171 REAL AJJ, SSTOP, STEMP
172 INTEGER I, ITEMP, J, JB, K, NB, PVT
178 LOGICAL LSAME, SISNAN
179 EXTERNAL slamch, ilaenv, lsame, sisnan
186 INTRINSIC conjg, max, min,
REAL, SQRT, MAXLOC
193 upper = lsame( uplo,
'U' )
194 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
196 ELSE IF( n.LT.0 )
THEN
198 ELSE IF( lda.LT.max( 1, n ) )
THEN
202 CALL xerbla(
'CPSTRF', -info )
213 nb = ilaenv( 1,
'CPOTRF', uplo, n, -1, -1, -1 )
214 IF( nb.LE.1 .OR. nb.GE.n )
THEN
218 CALL cpstf2( uplo, n, a( 1, 1 ), lda, piv, rank, tol, work,
233 work( i ) =
REAL( A( I, I ) )
235 pvt = maxloc( work( 1:n ), 1 )
236 ajj =
REAL( A( PVT, PVT ) )
237 IF( ajj.LE.zero.OR.sisnan( ajj ) )
THEN
245 IF( tol.LT.zero )
THEN
246 sstop = n * slamch(
'Epsilon' ) * ajj
260 jb = min( nb, n-k+1 )
269 DO 150 j = k, k + jb - 1
278 work( i ) = work( i ) +
279 $
REAL( CONJG( A( J-1, I ) )*
282 work( n+i ) =
REAL( A( I, I ) ) - WORK( i )
287 itemp = maxloc( work( (n+j):(2*n) ), 1 )
290 IF( ajj.LE.sstop.OR.sisnan( ajj ) )
THEN
300 a( pvt, pvt ) = a( j, j )
301 CALL cswap( j-1, a( 1, j ), 1, a( 1, pvt ), 1 )
303 $
CALL cswap( n-pvt, a( j, pvt+1 ), lda,
304 $ a( pvt, pvt+1 ), lda )
305 DO 140 i = j + 1, pvt - 1
306 ctemp = conjg( a( j, i ) )
307 a( j, i ) = conjg( a( i, pvt ) )
310 a( j, pvt ) = conjg( a( j, pvt ) )
315 work( j ) = work( pvt )
318 piv( pvt ) = piv( j )
328 CALL clacgv( j-1, a( 1, j ), 1 )
329 CALL cgemv(
'Trans', j-k, n-j, -cone, a( k, j+1 ),
330 $ lda, a( k, j ), 1, cone, a( j, j+1 ),
332 CALL clacgv( j-1, a( 1, j ), 1 )
333 CALL csscal( n-j, one / ajj, a( j, j+1 ), lda )
341 CALL cherk(
'Upper',
'Conj Trans', n-j+1, jb, -one,
342 $ a( k, j ), lda, one, a( j, j ), lda )
355 jb = min( nb, n-k+1 )
364 DO 200 j = k, k + jb - 1
373 work( i ) = work( i ) +
374 $
REAL( CONJG( A( I, J-1 ) )*
377 work( n+i ) =
REAL( A( I, I ) ) - WORK( i )
382 itemp = maxloc( work( (n+j):(2*n) ), 1 )
385 IF( ajj.LE.sstop.OR.sisnan( ajj ) )
THEN
395 a( pvt, pvt ) = a( j, j )
396 CALL cswap( j-1, a( j, 1 ), lda, a( pvt, 1 ), lda )
398 $
CALL cswap( n-pvt, a( pvt+1, j ), 1,
399 $ a( pvt+1, pvt ), 1 )
400 DO 190 i = j + 1, pvt - 1
401 ctemp = conjg( a( i, j ) )
402 a( i, j ) = conjg( a( pvt, i ) )
405 a( pvt, j ) = conjg( a( pvt, j ) )
410 work( j ) = work( pvt )
413 piv( pvt ) = piv( j )
423 CALL clacgv( j-1, a( j, 1 ), lda )
424 CALL cgemv(
'No Trans', n-j, j-k, -cone,
425 $ a( j+1, k ), lda, a( j, k ), lda, cone,
427 CALL clacgv( j-1, a( j, 1 ), lda )
428 CALL csscal( n-j, one / ajj, a( j+1, j ), 1 )
436 CALL cherk(
'Lower',
'No Trans', n-j+1, jb, -one,
437 $ a( j, k ), lda, one, a( j, j ), lda )
subroutine cpstf2(UPLO, N, A, LDA, PIV, RANK, TOL, WORK, INFO)
CPSTF2 computes the Cholesky factorization with complete pivoting of complex Hermitian positive semid...
subroutine cpstrf(UPLO, N, A, LDA, PIV, RANK, TOL, WORK, INFO)
CPSTRF computes the Cholesky factorization with complete pivoting of complex Hermitian positive semid...
subroutine cherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
CHERK
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine cgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CGEMV
subroutine clacgv(N, X, INCX)
CLACGV conjugates a complex vector.
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
subroutine csscal(N, SA, CX, INCX)
CSSCAL