version 1.7, 2015/11/26 11:44:23
|
version 1.13, 2023/08/07 08:39:25
|
Line 2
|
Line 2
|
* |
* |
* =========== DOCUMENTATION =========== |
* =========== DOCUMENTATION =========== |
* |
* |
* Online html documentation available at |
* Online html documentation available at |
* http://www.netlib.org/lapack/explore-html/ |
* http://www.netlib.org/lapack/explore-html/ |
* |
* |
*> \htmlonly |
*> \htmlonly |
*> Download ZHETRI2X + dependencies |
*> Download ZHETRI2X + dependencies |
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhetri2x.f"> |
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhetri2x.f"> |
*> [TGZ]</a> |
*> [TGZ]</a> |
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhetri2x.f"> |
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhetri2x.f"> |
*> [ZIP]</a> |
*> [ZIP]</a> |
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhetri2x.f"> |
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhetri2x.f"> |
*> [TXT]</a> |
*> [TXT]</a> |
*> \endhtmlonly |
*> \endhtmlonly |
* |
* |
* Definition: |
* Definition: |
* =========== |
* =========== |
* |
* |
* SUBROUTINE ZHETRI2X( UPLO, N, A, LDA, IPIV, WORK, NB, INFO ) |
* SUBROUTINE ZHETRI2X( UPLO, N, A, LDA, IPIV, WORK, NB, INFO ) |
* |
* |
* .. Scalar Arguments .. |
* .. Scalar Arguments .. |
* CHARACTER UPLO |
* CHARACTER UPLO |
* INTEGER INFO, LDA, N, NB |
* INTEGER INFO, LDA, N, NB |
Line 28
|
Line 28
|
* INTEGER IPIV( * ) |
* INTEGER IPIV( * ) |
* COMPLEX*16 A( LDA, * ), WORK( N+NB+1,* ) |
* COMPLEX*16 A( LDA, * ), WORK( N+NB+1,* ) |
* .. |
* .. |
* |
* |
* |
* |
*> \par Purpose: |
*> \par Purpose: |
* ============= |
* ============= |
Line 108
|
Line 108
|
* Authors: |
* Authors: |
* ======== |
* ======== |
* |
* |
*> \author Univ. of Tennessee |
*> \author Univ. of Tennessee |
*> \author Univ. of California Berkeley |
*> \author Univ. of California Berkeley |
*> \author Univ. of Colorado Denver |
*> \author Univ. of Colorado Denver |
*> \author NAG Ltd. |
*> \author NAG Ltd. |
* |
|
*> \date November 2015 |
|
* |
* |
*> \ingroup complex16HEcomputational |
*> \ingroup complex16HEcomputational |
* |
* |
* ===================================================================== |
* ===================================================================== |
SUBROUTINE ZHETRI2X( UPLO, N, A, LDA, IPIV, WORK, NB, INFO ) |
SUBROUTINE ZHETRI2X( UPLO, N, A, LDA, IPIV, WORK, NB, INFO ) |
* |
* |
* -- LAPACK computational routine (version 3.6.0) -- |
* -- LAPACK computational routine -- |
* -- LAPACK is a software package provided by Univ. of Tennessee, -- |
* -- LAPACK is a software package provided by Univ. of Tennessee, -- |
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- |
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- |
* November 2015 |
|
* |
* |
* .. Scalar Arguments .. |
* .. Scalar Arguments .. |
CHARACTER UPLO |
CHARACTER UPLO |
Line 137
|
Line 134
|
* ===================================================================== |
* ===================================================================== |
* |
* |
* .. Parameters .. |
* .. Parameters .. |
REAL ONE |
DOUBLE PRECISION ONE |
COMPLEX*16 CONE, ZERO |
COMPLEX*16 CONE, ZERO |
PARAMETER ( ONE = 1.0D+0, |
PARAMETER ( ONE = 1.0D+0, |
$ CONE = ( 1.0D+0, 0.0D+0 ), |
$ CONE = ( 1.0D+0, 0.0D+0 ), |
Line 215
|
Line 212
|
INFO = 0 |
INFO = 0 |
* |
* |
* Splitting Workspace |
* Splitting Workspace |
* U01 is a block (N,NB+1) |
* U01 is a block (N,NB+1) |
* The first element of U01 is in WORK(1,1) |
* The first element of U01 is in WORK(1,1) |
* U11 is a block (NB+1,NB+1) |
* U11 is a block (NB+1,NB+1) |
* The first element of U11 is in WORK(N+1,1) |
* The first element of U11 is in WORK(N+1,1) |
Line 231
|
Line 228
|
CALL ZTRTRI( UPLO, 'U', N, A, LDA, INFO ) |
CALL ZTRTRI( UPLO, 'U', N, A, LDA, INFO ) |
* |
* |
* inv(D) and inv(D)*inv(U) |
* inv(D) and inv(D)*inv(U) |
* |
* |
K=1 |
K=1 |
DO WHILE ( K .LE. N ) |
DO WHILE ( K .LE. N ) |
IF( IPIV( K ).GT.0 ) THEN |
IF( IPIV( K ).GT.0 ) THEN |
Line 242
|
Line 239
|
ELSE |
ELSE |
* 2 x 2 diagonal NNB |
* 2 x 2 diagonal NNB |
T = ABS ( WORK(K+1,1) ) |
T = ABS ( WORK(K+1,1) ) |
AK = REAL ( A( K, K ) ) / T |
AK = DBLE ( A( K, K ) ) / T |
AKP1 = REAL ( A( K+1, K+1 ) ) / T |
AKP1 = DBLE ( A( K+1, K+1 ) ) / T |
AKKP1 = WORK(K+1,1) / T |
AKKP1 = WORK(K+1,1) / T |
D = T*( AK*AKP1-ONE ) |
D = T*( AK*AKP1-ONE ) |
WORK(K,INVD) = AKP1 / D |
WORK(K,INVD) = AKP1 / D |
WORK(K+1,INVD+1) = AK / D |
WORK(K+1,INVD+1) = AK / D |
WORK(K,INVD+1) = -AKKP1 / D |
WORK(K,INVD+1) = -AKKP1 / D |
WORK(K+1,INVD) = DCONJG (WORK(K,INVD+1) ) |
WORK(K+1,INVD) = DCONJG (WORK(K,INVD+1) ) |
K=K+2 |
K=K+2 |
END IF |
END IF |
Line 265
|
Line 262
|
NNB=CUT |
NNB=CUT |
ELSE |
ELSE |
COUNT = 0 |
COUNT = 0 |
* count negative elements, |
* count negative elements, |
DO I=CUT+1-NNB,CUT |
DO I=CUT+1-NNB,CUT |
IF (IPIV(I) .LT. 0) COUNT=COUNT+1 |
IF (IPIV(I) .LT. 0) COUNT=COUNT+1 |
END DO |
END DO |
Line 275
|
Line 272
|
|
|
CUT=CUT-NNB |
CUT=CUT-NNB |
* |
* |
* U01 Block |
* U01 Block |
* |
* |
DO I=1,CUT |
DO I=1,CUT |
DO J=1,NNB |
DO J=1,NNB |
Line 338
|
Line 335
|
I=I+2 |
I=I+2 |
END IF |
END IF |
END DO |
END DO |
* |
* |
* U11**H*invD1*U11->U11 |
* U11**H*invD1*U11->U11 |
* |
* |
CALL ZTRMM('L','U','C','U',NNB, NNB, |
CALL ZTRMM('L','U','C','U',NNB, NNB, |
Line 382
|
Line 379
|
END DO |
END DO |
* |
* |
* Apply PERMUTATIONS P and P**H: P * inv(U**H)*inv(D)*inv(U) *P**H |
* Apply PERMUTATIONS P and P**H: P * inv(U**H)*inv(D)*inv(U) *P**H |
* |
* |
I=1 |
I=1 |
DO WHILE ( I .LE. N ) |
DO WHILE ( I .LE. N ) |
IF( IPIV(I) .GT. 0 ) THEN |
IF( IPIV(I) .GT. 0 ) THEN |
Line 392
|
Line 389
|
ELSE |
ELSE |
IP=-IPIV(I) |
IP=-IPIV(I) |
I=I+1 |
I=I+1 |
IF ( (I-1) .LT. IP) |
IF ( (I-1) .LT. IP) |
$ CALL ZHESWAPR( UPLO, N, A, LDA, I-1 ,IP ) |
$ CALL ZHESWAPR( UPLO, N, A, LDA, I-1 ,IP ) |
IF ( (I-1) .GT. IP) |
IF ( (I-1) .GT. IP) |
$ CALL ZHESWAPR( UPLO, N, A, LDA, IP ,I-1 ) |
$ CALL ZHESWAPR( UPLO, N, A, LDA, IP ,I-1 ) |
ENDIF |
ENDIF |
I=I+1 |
I=I+1 |
Line 408
|
Line 405
|
CALL ZTRTRI( UPLO, 'U', N, A, LDA, INFO ) |
CALL ZTRTRI( UPLO, 'U', N, A, LDA, INFO ) |
* |
* |
* inv(D) and inv(D)*inv(U) |
* inv(D) and inv(D)*inv(U) |
* |
* |
K=N |
K=N |
DO WHILE ( K .GE. 1 ) |
DO WHILE ( K .GE. 1 ) |
IF( IPIV( K ).GT.0 ) THEN |
IF( IPIV( K ).GT.0 ) THEN |
Line 419
|
Line 416
|
ELSE |
ELSE |
* 2 x 2 diagonal NNB |
* 2 x 2 diagonal NNB |
T = ABS ( WORK(K-1,1) ) |
T = ABS ( WORK(K-1,1) ) |
AK = REAL ( A( K-1, K-1 ) ) / T |
AK = DBLE ( A( K-1, K-1 ) ) / T |
AKP1 = REAL ( A( K, K ) ) / T |
AKP1 = DBLE ( A( K, K ) ) / T |
AKKP1 = WORK(K-1,1) / T |
AKKP1 = WORK(K-1,1) / T |
D = T*( AK*AKP1-ONE ) |
D = T*( AK*AKP1-ONE ) |
WORK(K-1,INVD) = AKP1 / D |
WORK(K-1,INVD) = AKP1 / D |
WORK(K,INVD) = AK / D |
WORK(K,INVD) = AK / D |
WORK(K,INVD+1) = -AKKP1 / D |
WORK(K,INVD+1) = -AKKP1 / D |
WORK(K-1,INVD+1) = DCONJG (WORK(K,INVD+1) ) |
WORK(K-1,INVD+1) = DCONJG (WORK(K,INVD+1) ) |
K=K-2 |
K=K-2 |
END IF |
END IF |
Line 442
|
Line 439
|
NNB=N-CUT |
NNB=N-CUT |
ELSE |
ELSE |
COUNT = 0 |
COUNT = 0 |
* count negative elements, |
* count negative elements, |
DO I=CUT+1,CUT+NNB |
DO I=CUT+1,CUT+NNB |
IF (IPIV(I) .LT. 0) COUNT=COUNT+1 |
IF (IPIV(I) .LT. 0) COUNT=COUNT+1 |
END DO |
END DO |
Line 509
|
Line 506
|
I=I-2 |
I=I-2 |
END IF |
END IF |
END DO |
END DO |
* |
* |
* L11**H*invD1*L11->L11 |
* L11**H*invD1*L11->L11 |
* |
* |
CALL ZTRMM('L',UPLO,'C','U',NNB, NNB, |
CALL ZTRMM('L',UPLO,'C','U',NNB, NNB, |
Line 527
|
Line 524
|
* |
* |
CALL ZGEMM('C','N',NNB,NNB,N-NNB-CUT,CONE,A(CUT+NNB+1,CUT+1) |
CALL ZGEMM('C','N',NNB,NNB,N-NNB-CUT,CONE,A(CUT+NNB+1,CUT+1) |
$ ,LDA,WORK,N+NB+1, ZERO, WORK(U11+1,1), N+NB+1) |
$ ,LDA,WORK,N+NB+1, ZERO, WORK(U11+1,1), N+NB+1) |
|
|
* |
* |
* L11 = L11**H*invD1*L11 + U01**H*invD*U01 |
* L11 = L11**H*invD1*L11 + U01**H*invD*U01 |
* |
* |
Line 565
|
Line 562
|
END DO |
END DO |
* |
* |
* Apply PERMUTATIONS P and P**H: P * inv(U**H)*inv(D)*inv(U) *P**H |
* Apply PERMUTATIONS P and P**H: P * inv(U**H)*inv(D)*inv(U) *P**H |
* |
* |
I=N |
I=N |
DO WHILE ( I .GE. 1 ) |
DO WHILE ( I .GE. 1 ) |
IF( IPIV(I) .GT. 0 ) THEN |
IF( IPIV(I) .GT. 0 ) THEN |