243 SUBROUTINE chetrd_he2hb( UPLO, N, KD, A, LDA, AB, LDAB, TAU,
244 $ work, lwork, info )
255 INTEGER INFO, LDA, LDAB, LWORK, N, KD
258 COMPLEX A( lda, * ), AB( ldab, * ),
259 $ tau( * ), work( * )
266 COMPLEX ZERO, ONE, HALF
267 parameter ( rone = 1.0e+0,
268 $ zero = ( 0.0e+0, 0.0e+0 ),
269 $ one = ( 1.0e+0, 0.0e+0 ),
270 $ half = ( 0.5e+0, 0.0e+0 ) )
273 LOGICAL LQUERY, UPPER
274 INTEGER I, J, IINFO, LWMIN, PN, PK, LK,
275 $ ldt, ldw, lds2, lds1,
277 $ tpos, wpos, s2pos, s1pos
289 EXTERNAL lsame, ilaenv
297 upper = lsame( uplo,
'U' )
298 lquery = ( lwork.EQ.-1 )
299 lwmin = ilaenv( 20,
'CHETRD_HE2HB',
'', n, kd, -1, -1 )
301 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
303 ELSE IF( n.LT.0 )
THEN
305 ELSE IF( kd.LT.0 )
THEN
307 ELSE IF( lda.LT.max( 1, n ) )
THEN
309 ELSE IF( ldab.LT.max( 1, kd+1 ) )
THEN
311 ELSE IF( lwork.LT.lwmin .AND. .NOT.lquery )
THEN
316 CALL xerbla(
'CHETRD_HE2HB', -info )
318 ELSE IF( lquery )
THEN
330 CALL ccopy( lk, a( i-lk+1, i ), 1,
331 $ ab( kd+1-lk+1, i ), 1 )
335 lk = min( kd+1, n-i+1 )
336 CALL ccopy( lk, a( i, i ), 1, ab( 1, i ), 1 )
350 ls2 = lwmin - lt - lw - ls1
368 CALL claset(
"A", ldt, kd, zero, zero, work( tpos ), ldt )
371 DO 10 i = 1, n - kd, kd
373 pk = min( n-i-kd+1, kd )
377 CALL cgelqf( kd, pn, a( i, i+kd ), lda,
378 $ tau( i ), work( s2pos ), ls2, iinfo )
383 lk = min( kd, n-j ) + 1
384 CALL ccopy( lk, a( j, j ), lda, ab( kd+1, j ), ldab-1 )
387 CALL claset(
'Lower', pk, pk, zero, one,
388 $ a( i, i+kd ), lda )
392 CALL clarft(
'Forward',
'Rowwise', pn, pk,
393 $ a( i, i+kd ), lda, tau( i ),
394 $ work( tpos ), ldt )
398 CALL cgemm(
'Conjugate',
'No transpose', pk, pn, pk,
399 $ one, work( tpos ), ldt,
401 $ zero, work( s2pos ), lds2 )
403 CALL chemm(
'Right', uplo, pk, pn,
404 $ one, a( i+kd, i+kd ), lda,
405 $ work( s2pos ), lds2,
406 $ zero, work( wpos ), ldw )
408 CALL cgemm(
'No transpose',
'Conjugate', pk, pk, pn,
409 $ one, work( wpos ), ldw,
410 $ work( s2pos ), lds2,
411 $ zero, work( s1pos ), lds1 )
413 CALL cgemm(
'No transpose',
'No transpose', pk, pn, pk,
414 $ -half, work( s1pos ), lds1,
416 $ one, work( wpos ), ldw )
422 CALL cher2k( uplo,
'Conjugate', pn, pk,
423 $ -one, a( i, i+kd ), lda,
425 $ rone, a( i+kd, i+kd ), lda )
431 lk = min(kd, n-j) + 1
432 CALL ccopy( lk, a( j, j ), lda, ab( kd+1, j ), ldab-1 )
439 DO 40 i = 1, n - kd, kd
441 pk = min( n-i-kd+1, kd )
445 CALL cgeqrf( pn, kd, a( i+kd, i ), lda,
446 $ tau( i ), work( s2pos ), ls2, iinfo )
451 lk = min( kd, n-j ) + 1
452 CALL ccopy( lk, a( j, j ), 1, ab( 1, j ), 1 )
455 CALL claset(
'Upper', pk, pk, zero, one,
456 $ a( i+kd, i ), lda )
460 CALL clarft(
'Forward',
'Columnwise', pn, pk,
461 $ a( i+kd, i ), lda, tau( i ),
462 $ work( tpos ), ldt )
466 CALL cgemm(
'No transpose',
'No transpose', pn, pk, pk,
467 $ one, a( i+kd, i ), lda,
469 $ zero, work( s2pos ), lds2 )
471 CALL chemm(
'Left', uplo, pn, pk,
472 $ one, a( i+kd, i+kd ), lda,
473 $ work( s2pos ), lds2,
474 $ zero, work( wpos ), ldw )
476 CALL cgemm(
'Conjugate',
'No transpose', pk, pk, pn,
477 $ one, work( s2pos ), lds2,
479 $ zero, work( s1pos ), lds1 )
481 CALL cgemm(
'No transpose',
'No transpose', pn, pk, pk,
482 $ -half, a( i+kd, i ), lda,
483 $ work( s1pos ), lds1,
484 $ one, work( wpos ), ldw )
490 CALL cher2k( uplo,
'No transpose', pn, pk,
491 $ -one, a( i+kd, i ), lda,
493 $ rone, a( i+kd, i+kd ), lda )
506 lk = min(kd, n-j) + 1
507 CALL ccopy( lk, a( j, j ), 1, ab( 1, j ), 1 )
subroutine clarft(DIRECT, STOREV, N, K, V, LDV, TAU, T, LDT)
CLARFT forms the triangular factor T of a block reflector H = I - vtvH
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 chemm(SIDE, UPLO, M, N, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CHEMM
subroutine chetrd_he2hb(UPLO, N, KD, A, LDA, AB, LDAB, TAU, WORK, LWORK, INFO)
CHETRD_HE2HB
subroutine cgelqf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGELQF
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
subroutine cher2k(UPLO, TRANS, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CHER2K
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM