154 SUBROUTINE clahef_aa( UPLO, J1, M, NB, A, LDA, IPIV,
155 $ h, ldh, work, info )
166 INTEGER M, NB, J1, LDA, LDH, INFO
170 COMPLEX A( lda, * ), H( ldh, * ), WORK( * )
176 parameter ( zero = (0.0e+0, 0.0e+0), one = (1.0e+0, 0.0e+0) )
179 INTEGER J, K, K1, I1, I2
184 INTEGER ICAMAX, ILAENV
185 EXTERNAL lsame, ilaenv, icamax
191 INTRINSIC REAL, CONJG, MAX
203 IF( lsame( uplo,
'U' ) )
THEN
210 IF ( j.GT.min(m, nb) )
231 CALL clacgv( j-k1, a( 1, j ), 1 )
232 CALL cgemv(
'No transpose', m-j+1, j-k1,
233 $ -one, h( j, k1 ), ldh,
235 $ one, h( j, j ), 1 )
236 CALL clacgv( j-k1, a( 1, j ), 1 )
241 CALL ccopy( m-j+1, h( j, j ), 1, work( 1 ), 1 )
248 alpha = -conjg( a( k-1, j ) )
249 CALL caxpy( m-j+1, alpha, a( k-2, j ), lda, work( 1 ), 1 )
254 a( k, j ) =
REAL( WORK( 1 ) )
263 CALL caxpy( m-j, alpha, a( k-1, j+1 ), lda,
269 i2 = icamax( m-j, work( 2 ), 1 ) + 1
274 IF( (i2.NE.2) .AND. (piv.NE.0) )
THEN
279 work( i2 ) = work( i1 )
286 CALL cswap( i2-i1-1, a( j1+i1-1, i1+1 ), lda,
287 $ a( j1+i1, i2 ), 1 )
288 CALL clacgv( i2-i1, a( j1+i1-1, i1+1 ), lda )
289 CALL clacgv( i2-i1-1, a( j1+i1, i2 ), 1 )
293 CALL cswap( m-i2, a( j1+i1-1, i2+1 ), lda,
294 $ a( j1+i2-1, i2+1 ), lda )
298 piv = a( i1+j1-1, i1 )
299 a( j1+i1-1, i1 ) = a( j1+i2-1, i2 )
300 a( j1+i2-1, i2 ) = piv
304 CALL cswap( i1-1, h( i1, 1 ), ldh, h( i2, 1 ), ldh )
307 IF( i1.GT.(k1-1) )
THEN
312 CALL cswap( i1-k1+1, a( 1, i1 ), 1,
321 a( k, j+1 ) = work( 2 )
322 IF( (a( k, j ).EQ.zero ) .AND.
323 $ ( (j.EQ.m) .OR. (a( k, j+1 ).EQ.zero)))
THEN
333 CALL ccopy( m-j, a( k+1, j+1 ), lda,
340 IF( a( k, j+1 ).NE.zero )
THEN
341 alpha = one / a( k, j+1 )
342 CALL ccopy( m-j-1, work( 3 ), 1, a( k, j+2 ), lda )
343 CALL cscal( m-j-1, alpha, a( k, j+2 ), lda )
345 CALL claset(
'Full', 1, m-j-1, zero, zero,
349 IF( (a( k, j ).EQ.zero) .AND. (info.EQ.0) )
THEN
364 IF( j.GT.min( m, nb ) )
385 CALL clacgv( j-k1, a( j, 1 ), lda )
386 CALL cgemv(
'No transpose', m-j+1, j-k1,
387 $ -one, h( j, k1 ), ldh,
389 $ one, h( j, j ), 1 )
390 CALL clacgv( j-k1, a( j, 1 ), lda )
395 CALL ccopy( m-j+1, h( j, j ), 1, work( 1 ), 1 )
402 alpha = -conjg( a( j, k-1 ) )
403 CALL caxpy( m-j+1, alpha, a( j, k-2 ), 1, work( 1 ), 1 )
408 a( j, k ) =
REAL( WORK( 1 ) )
417 CALL caxpy( m-j, alpha, a( j+1, k-1 ), 1,
423 i2 = icamax( m-j, work( 2 ), 1 ) + 1
428 IF( (i2.NE.2) .AND. (piv.NE.0) )
THEN
433 work( i2 ) = work( i1 )
440 CALL cswap( i2-i1-1, a( i1+1, j1+i1-1 ), 1,
441 $ a( i2, j1+i1 ), lda )
442 CALL clacgv( i2-i1, a( i1+1, j1+i1-1 ), 1 )
443 CALL clacgv( i2-i1-1, a( i2, j1+i1 ), lda )
447 CALL cswap( m-i2, a( i2+1, j1+i1-1 ), 1,
448 $ a( i2+1, j1+i2-1 ), 1 )
452 piv = a( i1, j1+i1-1 )
453 a( i1, j1+i1-1 ) = a( i2, j1+i2-1 )
454 a( i2, j1+i2-1 ) = piv
458 CALL cswap( i1-1, h( i1, 1 ), ldh, h( i2, 1 ), ldh )
461 IF( i1.GT.(k1-1) )
THEN
466 CALL cswap( i1-k1+1, a( i1, 1 ), lda,
475 a( j+1, k ) = work( 2 )
476 IF( (a( j, k ).EQ.zero) .AND.
477 $ ( (j.EQ.m) .OR. (a( j+1, k ).EQ.zero)) )
THEN
486 CALL ccopy( m-j, a( j+1, k+1 ), 1,
493 IF( a( j+1, k ).NE.zero )
THEN
494 alpha = one / a( j+1, k )
495 CALL ccopy( m-j-1, work( 3 ), 1, a( j+2, k ), 1 )
496 CALL cscal( m-j-1, alpha, a( j+2, k ), 1 )
498 CALL claset(
'Full', m-j-1, 1, zero, zero,
502 IF( (a( j, k ).EQ.zero) .AND. (j.EQ.m)
503 $ .AND. (info.EQ.0) ) info = j
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
subroutine cscal(N, CA, CX, INCX)
CSCAL
subroutine cgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CGEMV
subroutine clahef_aa(UPLO, J1, M, NB, A, LDA, IPIV, H, LDH, WORK, INFO)
CLAHEF_AA
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine clacgv(N, X, INCX)
CLACGV conjugates a complex vector.
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
subroutine caxpy(N, CA, CX, INCX, CY, INCY)
CAXPY