![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE ZPOTRF( UPLO, N, A, LDA, INFO ) 2: * 3: * -- LAPACK routine (version 3.2) -- 4: * -- LAPACK is a software package provided by Univ. of Tennessee, -- 5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 6: * November 2006 7: * 8: * .. Scalar Arguments .. 9: CHARACTER UPLO 10: INTEGER INFO, LDA, N 11: * .. 12: * .. Array Arguments .. 13: COMPLEX*16 A( LDA, * ) 14: * .. 15: * 16: * Purpose 17: * ======= 18: * 19: * ZPOTRF computes the Cholesky factorization of a complex Hermitian 20: * positive definite matrix A. 21: * 22: * The factorization has the form 23: * A = U**H * U, if UPLO = 'U', or 24: * A = L * L**H, if UPLO = 'L', 25: * where U is an upper triangular matrix and L is lower triangular. 26: * 27: * This is the block version of the algorithm, calling Level 3 BLAS. 28: * 29: * Arguments 30: * ========= 31: * 32: * UPLO (input) CHARACTER*1 33: * = 'U': Upper triangle of A is stored; 34: * = 'L': Lower triangle of A is stored. 35: * 36: * N (input) INTEGER 37: * The order of the matrix A. N >= 0. 38: * 39: * A (input/output) COMPLEX*16 array, dimension (LDA,N) 40: * On entry, the Hermitian matrix A. If UPLO = 'U', the leading 41: * N-by-N upper triangular part of A contains the upper 42: * triangular part of the matrix A, and the strictly lower 43: * triangular part of A is not referenced. If UPLO = 'L', the 44: * leading N-by-N lower triangular part of A contains the lower 45: * triangular part of the matrix A, and the strictly upper 46: * triangular part of A is not referenced. 47: * 48: * On exit, if INFO = 0, the factor U or L from the Cholesky 49: * factorization A = U**H*U or A = L*L**H. 50: * 51: * LDA (input) INTEGER 52: * The leading dimension of the array A. LDA >= max(1,N). 53: * 54: * INFO (output) INTEGER 55: * = 0: successful exit 56: * < 0: if INFO = -i, the i-th argument had an illegal value 57: * > 0: if INFO = i, the leading minor of order i is not 58: * positive definite, and the factorization could not be 59: * completed. 60: * 61: * ===================================================================== 62: * 63: * .. Parameters .. 64: DOUBLE PRECISION ONE 65: COMPLEX*16 CONE 66: PARAMETER ( ONE = 1.0D+0, CONE = ( 1.0D+0, 0.0D+0 ) ) 67: * .. 68: * .. Local Scalars .. 69: LOGICAL UPPER 70: INTEGER J, JB, NB 71: * .. 72: * .. External Functions .. 73: LOGICAL LSAME 74: INTEGER ILAENV 75: EXTERNAL LSAME, ILAENV 76: * .. 77: * .. External Subroutines .. 78: EXTERNAL XERBLA, ZGEMM, ZHERK, ZPOTF2, ZTRSM 79: * .. 80: * .. Intrinsic Functions .. 81: INTRINSIC MAX, MIN 82: * .. 83: * .. Executable Statements .. 84: * 85: * Test the input parameters. 86: * 87: INFO = 0 88: UPPER = LSAME( UPLO, 'U' ) 89: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN 90: INFO = -1 91: ELSE IF( N.LT.0 ) THEN 92: INFO = -2 93: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 94: INFO = -4 95: END IF 96: IF( INFO.NE.0 ) THEN 97: CALL XERBLA( 'ZPOTRF', -INFO ) 98: RETURN 99: END IF 100: * 101: * Quick return if possible 102: * 103: IF( N.EQ.0 ) 104: $ RETURN 105: * 106: * Determine the block size for this environment. 107: * 108: NB = ILAENV( 1, 'ZPOTRF', UPLO, N, -1, -1, -1 ) 109: IF( NB.LE.1 .OR. NB.GE.N ) THEN 110: * 111: * Use unblocked code. 112: * 113: CALL ZPOTF2( UPLO, N, A, LDA, INFO ) 114: ELSE 115: * 116: * Use blocked code. 117: * 118: IF( UPPER ) THEN 119: * 120: * Compute the Cholesky factorization A = U'*U. 121: * 122: DO 10 J = 1, N, NB 123: * 124: * Update and factorize the current diagonal block and test 125: * for non-positive-definiteness. 126: * 127: JB = MIN( NB, N-J+1 ) 128: CALL ZHERK( 'Upper', 'Conjugate transpose', JB, J-1, 129: $ -ONE, A( 1, J ), LDA, ONE, A( J, J ), LDA ) 130: CALL ZPOTF2( 'Upper', JB, A( J, J ), LDA, INFO ) 131: IF( INFO.NE.0 ) 132: $ GO TO 30 133: IF( J+JB.LE.N ) THEN 134: * 135: * Compute the current block row. 136: * 137: CALL ZGEMM( 'Conjugate transpose', 'No transpose', JB, 138: $ N-J-JB+1, J-1, -CONE, A( 1, J ), LDA, 139: $ A( 1, J+JB ), LDA, CONE, A( J, J+JB ), 140: $ LDA ) 141: CALL ZTRSM( 'Left', 'Upper', 'Conjugate transpose', 142: $ 'Non-unit', JB, N-J-JB+1, CONE, A( J, J ), 143: $ LDA, A( J, J+JB ), LDA ) 144: END IF 145: 10 CONTINUE 146: * 147: ELSE 148: * 149: * Compute the Cholesky factorization A = L*L'. 150: * 151: DO 20 J = 1, N, NB 152: * 153: * Update and factorize the current diagonal block and test 154: * for non-positive-definiteness. 155: * 156: JB = MIN( NB, N-J+1 ) 157: CALL ZHERK( 'Lower', 'No transpose', JB, J-1, -ONE, 158: $ A( J, 1 ), LDA, ONE, A( J, J ), LDA ) 159: CALL ZPOTF2( 'Lower', JB, A( J, J ), LDA, INFO ) 160: IF( INFO.NE.0 ) 161: $ GO TO 30 162: IF( J+JB.LE.N ) THEN 163: * 164: * Compute the current block column. 165: * 166: CALL ZGEMM( 'No transpose', 'Conjugate transpose', 167: $ N-J-JB+1, JB, J-1, -CONE, A( J+JB, 1 ), 168: $ LDA, A( J, 1 ), LDA, CONE, A( J+JB, J ), 169: $ LDA ) 170: CALL ZTRSM( 'Right', 'Lower', 'Conjugate transpose', 171: $ 'Non-unit', N-J-JB+1, JB, CONE, A( J, J ), 172: $ LDA, A( J+JB, J ), LDA ) 173: END IF 174: 20 CONTINUE 175: END IF 176: END IF 177: GO TO 40 178: * 179: 30 CONTINUE 180: INFO = INFO + J - 1 181: * 182: 40 CONTINUE 183: RETURN 184: * 185: * End of ZPOTRF 186: * 187: END