 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 DOUBLE PRECISION FUNCTION ZLANHP( NORM, UPLO, N, AP, WORK ) 00002 * 00003 * -- LAPACK auxiliary routine (version 3.2) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 CHARACTER NORM, UPLO 00010 INTEGER N 00011 * .. 00012 * .. Array Arguments .. 00013 DOUBLE PRECISION WORK( * ) 00014 COMPLEX*16 AP( * ) 00015 * .. 00016 * 00017 * Purpose 00018 * ======= 00019 * 00020 * ZLANHP returns the value of the one norm, or the Frobenius norm, or 00021 * the infinity norm, or the element of largest absolute value of a 00022 * complex hermitian matrix A, supplied in packed form. 00023 * 00024 * Description 00025 * =========== 00026 * 00027 * ZLANHP returns the value 00028 * 00029 * ZLANHP = ( max(abs(A(i,j))), NORM = 'M' or 'm' 00030 * ( 00031 * ( norm1(A), NORM = '1', 'O' or 'o' 00032 * ( 00033 * ( normI(A), NORM = 'I' or 'i' 00034 * ( 00035 * ( normF(A), NORM = 'F', 'f', 'E' or 'e' 00036 * 00037 * where norm1 denotes the one norm of a matrix (maximum column sum), 00038 * normI denotes the infinity norm of a matrix (maximum row sum) and 00039 * normF denotes the Frobenius norm of a matrix (square root of sum of 00040 * squares). Note that max(abs(A(i,j))) is not a consistent matrix norm. 00041 * 00042 * Arguments 00043 * ========= 00044 * 00045 * NORM (input) CHARACTER*1 00046 * Specifies the value to be returned in ZLANHP as described 00047 * above. 00048 * 00049 * UPLO (input) CHARACTER*1 00050 * Specifies whether the upper or lower triangular part of the 00051 * hermitian matrix A is supplied. 00052 * = 'U': Upper triangular part of A is supplied 00053 * = 'L': Lower triangular part of A is supplied 00054 * 00055 * N (input) INTEGER 00056 * The order of the matrix A. N >= 0. When N = 0, ZLANHP is 00057 * set to zero. 00058 * 00059 * AP (input) COMPLEX*16 array, dimension (N*(N+1)/2) 00060 * The upper or lower triangle of the hermitian matrix A, packed 00061 * columnwise in a linear array. The j-th column of A is stored 00062 * in the array AP as follows: 00063 * if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; 00064 * if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. 00065 * Note that the imaginary parts of the diagonal elements need 00066 * not be set and are assumed to be zero. 00067 * 00068 * WORK (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK)), 00069 * where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise, 00070 * WORK is not referenced. 00071 * 00072 * ===================================================================== 00073 * 00074 * .. Parameters .. 00075 DOUBLE PRECISION ONE, ZERO 00076 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00077 * .. 00078 * .. Local Scalars .. 00079 INTEGER I, J, K 00080 DOUBLE PRECISION ABSA, SCALE, SUM, VALUE 00081 * .. 00082 * .. External Functions .. 00083 LOGICAL LSAME 00084 EXTERNAL LSAME 00085 * .. 00086 * .. External Subroutines .. 00087 EXTERNAL ZLASSQ 00088 * .. 00089 * .. Intrinsic Functions .. 00090 INTRINSIC ABS, DBLE, MAX, SQRT 00091 * .. 00092 * .. Executable Statements .. 00093 * 00094 IF( N.EQ.0 ) THEN 00095 VALUE = ZERO 00096 ELSE IF( LSAME( NORM, 'M' ) ) THEN 00097 * 00098 * Find max(abs(A(i,j))). 00099 * 00100 VALUE = ZERO 00101 IF( LSAME( UPLO, 'U' ) ) THEN 00102 K = 0 00103 DO 20 J = 1, N 00104 DO 10 I = K + 1, K + J - 1 00105 VALUE = MAX( VALUE, ABS( AP( I ) ) ) 00106 10 CONTINUE 00107 K = K + J 00108 VALUE = MAX( VALUE, ABS( DBLE( AP( K ) ) ) ) 00109 20 CONTINUE 00110 ELSE 00111 K = 1 00112 DO 40 J = 1, N 00113 VALUE = MAX( VALUE, ABS( DBLE( AP( K ) ) ) ) 00114 DO 30 I = K + 1, K + N - J 00115 VALUE = MAX( VALUE, ABS( AP( I ) ) ) 00116 30 CONTINUE 00117 K = K + N - J + 1 00118 40 CONTINUE 00119 END IF 00120 ELSE IF( ( LSAME( NORM, 'I' ) ) .OR. ( LSAME( NORM, 'O' ) ) .OR. 00121 $ ( NORM.EQ.'1' ) ) THEN 00122 * 00123 * Find normI(A) ( = norm1(A), since A is hermitian). 00124 * 00125 VALUE = ZERO 00126 K = 1 00127 IF( LSAME( UPLO, 'U' ) ) THEN 00128 DO 60 J = 1, N 00129 SUM = ZERO 00130 DO 50 I = 1, J - 1 00131 ABSA = ABS( AP( K ) ) 00132 SUM = SUM + ABSA 00133 WORK( I ) = WORK( I ) + ABSA 00134 K = K + 1 00135 50 CONTINUE 00136 WORK( J ) = SUM + ABS( DBLE( AP( K ) ) ) 00137 K = K + 1 00138 60 CONTINUE 00139 DO 70 I = 1, N 00140 VALUE = MAX( VALUE, WORK( I ) ) 00141 70 CONTINUE 00142 ELSE 00143 DO 80 I = 1, N 00144 WORK( I ) = ZERO 00145 80 CONTINUE 00146 DO 100 J = 1, N 00147 SUM = WORK( J ) + ABS( DBLE( AP( K ) ) ) 00148 K = K + 1 00149 DO 90 I = J + 1, N 00150 ABSA = ABS( AP( K ) ) 00151 SUM = SUM + ABSA 00152 WORK( I ) = WORK( I ) + ABSA 00153 K = K + 1 00154 90 CONTINUE 00155 VALUE = MAX( VALUE, SUM ) 00156 100 CONTINUE 00157 END IF 00158 ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN 00159 * 00160 * Find normF(A). 00161 * 00162 SCALE = ZERO 00163 SUM = ONE 00164 K = 2 00165 IF( LSAME( UPLO, 'U' ) ) THEN 00166 DO 110 J = 2, N 00167 CALL ZLASSQ( J-1, AP( K ), 1, SCALE, SUM ) 00168 K = K + J 00169 110 CONTINUE 00170 ELSE 00171 DO 120 J = 1, N - 1 00172 CALL ZLASSQ( N-J, AP( K ), 1, SCALE, SUM ) 00173 K = K + N - J + 1 00174 120 CONTINUE 00175 END IF 00176 SUM = 2*SUM 00177 K = 1 00178 DO 130 I = 1, N 00179 IF( DBLE( AP( K ) ).NE.ZERO ) THEN 00180 ABSA = ABS( DBLE( AP( K ) ) ) 00181 IF( SCALE.LT.ABSA ) THEN 00182 SUM = ONE + SUM*( SCALE / ABSA )**2 00183 SCALE = ABSA 00184 ELSE 00185 SUM = SUM + ( ABSA / SCALE )**2 00186 END IF 00187 END IF 00188 IF( LSAME( UPLO, 'U' ) ) THEN 00189 K = K + I + 1 00190 ELSE 00191 K = K + N - I + 1 00192 END IF 00193 130 CONTINUE 00194 VALUE = SCALE*SQRT( SUM ) 00195 END IF 00196 * 00197 ZLANHP = VALUE 00198 RETURN 00199 * 00200 * End of ZLANHP 00201 * 00202 END
1.7.3 
