LAPACK 3.3.0
|
00001 DOUBLE PRECISION FUNCTION ZLANTP( NORM, UPLO, DIAG, 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 DIAG, NORM, UPLO 00010 INTEGER N 00011 * .. 00012 * .. Array Arguments .. 00013 DOUBLE PRECISION WORK( * ) 00014 COMPLEX*16 AP( * ) 00015 * .. 00016 * 00017 * Purpose 00018 * ======= 00019 * 00020 * ZLANTP 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 * triangular matrix A, supplied in packed form. 00023 * 00024 * Description 00025 * =========== 00026 * 00027 * ZLANTP returns the value 00028 * 00029 * ZLANTP = ( 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 ZLANTP as described 00047 * above. 00048 * 00049 * UPLO (input) CHARACTER*1 00050 * Specifies whether the matrix A is upper or lower triangular. 00051 * = 'U': Upper triangular 00052 * = 'L': Lower triangular 00053 * 00054 * DIAG (input) CHARACTER*1 00055 * Specifies whether or not the matrix A is unit triangular. 00056 * = 'N': Non-unit triangular 00057 * = 'U': Unit triangular 00058 * 00059 * N (input) INTEGER 00060 * The order of the matrix A. N >= 0. When N = 0, ZLANTP is 00061 * set to zero. 00062 * 00063 * AP (input) COMPLEX*16 array, dimension (N*(N+1)/2) 00064 * The upper or lower triangular matrix A, packed columnwise in 00065 * a linear array. The j-th column of A is stored in the array 00066 * AP as follows: 00067 * if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; 00068 * if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. 00069 * Note that when DIAG = 'U', the elements of the array AP 00070 * corresponding to the diagonal elements of the matrix A are 00071 * not referenced, but are assumed to be one. 00072 * 00073 * WORK (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK)), 00074 * where LWORK >= N when NORM = 'I'; otherwise, WORK is not 00075 * referenced. 00076 * 00077 * ===================================================================== 00078 * 00079 * .. Parameters .. 00080 DOUBLE PRECISION ONE, ZERO 00081 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00082 * .. 00083 * .. Local Scalars .. 00084 LOGICAL UDIAG 00085 INTEGER I, J, K 00086 DOUBLE PRECISION SCALE, SUM, VALUE 00087 * .. 00088 * .. External Functions .. 00089 LOGICAL LSAME 00090 EXTERNAL LSAME 00091 * .. 00092 * .. External Subroutines .. 00093 EXTERNAL ZLASSQ 00094 * .. 00095 * .. Intrinsic Functions .. 00096 INTRINSIC ABS, MAX, SQRT 00097 * .. 00098 * .. Executable Statements .. 00099 * 00100 IF( N.EQ.0 ) THEN 00101 VALUE = ZERO 00102 ELSE IF( LSAME( NORM, 'M' ) ) THEN 00103 * 00104 * Find max(abs(A(i,j))). 00105 * 00106 K = 1 00107 IF( LSAME( DIAG, 'U' ) ) THEN 00108 VALUE = ONE 00109 IF( LSAME( UPLO, 'U' ) ) THEN 00110 DO 20 J = 1, N 00111 DO 10 I = K, K + J - 2 00112 VALUE = MAX( VALUE, ABS( AP( I ) ) ) 00113 10 CONTINUE 00114 K = K + J 00115 20 CONTINUE 00116 ELSE 00117 DO 40 J = 1, N 00118 DO 30 I = K + 1, K + N - J 00119 VALUE = MAX( VALUE, ABS( AP( I ) ) ) 00120 30 CONTINUE 00121 K = K + N - J + 1 00122 40 CONTINUE 00123 END IF 00124 ELSE 00125 VALUE = ZERO 00126 IF( LSAME( UPLO, 'U' ) ) THEN 00127 DO 60 J = 1, N 00128 DO 50 I = K, K + J - 1 00129 VALUE = MAX( VALUE, ABS( AP( I ) ) ) 00130 50 CONTINUE 00131 K = K + J 00132 60 CONTINUE 00133 ELSE 00134 DO 80 J = 1, N 00135 DO 70 I = K, K + N - J 00136 VALUE = MAX( VALUE, ABS( AP( I ) ) ) 00137 70 CONTINUE 00138 K = K + N - J + 1 00139 80 CONTINUE 00140 END IF 00141 END IF 00142 ELSE IF( ( LSAME( NORM, 'O' ) ) .OR. ( NORM.EQ.'1' ) ) THEN 00143 * 00144 * Find norm1(A). 00145 * 00146 VALUE = ZERO 00147 K = 1 00148 UDIAG = LSAME( DIAG, 'U' ) 00149 IF( LSAME( UPLO, 'U' ) ) THEN 00150 DO 110 J = 1, N 00151 IF( UDIAG ) THEN 00152 SUM = ONE 00153 DO 90 I = K, K + J - 2 00154 SUM = SUM + ABS( AP( I ) ) 00155 90 CONTINUE 00156 ELSE 00157 SUM = ZERO 00158 DO 100 I = K, K + J - 1 00159 SUM = SUM + ABS( AP( I ) ) 00160 100 CONTINUE 00161 END IF 00162 K = K + J 00163 VALUE = MAX( VALUE, SUM ) 00164 110 CONTINUE 00165 ELSE 00166 DO 140 J = 1, N 00167 IF( UDIAG ) THEN 00168 SUM = ONE 00169 DO 120 I = K + 1, K + N - J 00170 SUM = SUM + ABS( AP( I ) ) 00171 120 CONTINUE 00172 ELSE 00173 SUM = ZERO 00174 DO 130 I = K, K + N - J 00175 SUM = SUM + ABS( AP( I ) ) 00176 130 CONTINUE 00177 END IF 00178 K = K + N - J + 1 00179 VALUE = MAX( VALUE, SUM ) 00180 140 CONTINUE 00181 END IF 00182 ELSE IF( LSAME( NORM, 'I' ) ) THEN 00183 * 00184 * Find normI(A). 00185 * 00186 K = 1 00187 IF( LSAME( UPLO, 'U' ) ) THEN 00188 IF( LSAME( DIAG, 'U' ) ) THEN 00189 DO 150 I = 1, N 00190 WORK( I ) = ONE 00191 150 CONTINUE 00192 DO 170 J = 1, N 00193 DO 160 I = 1, J - 1 00194 WORK( I ) = WORK( I ) + ABS( AP( K ) ) 00195 K = K + 1 00196 160 CONTINUE 00197 K = K + 1 00198 170 CONTINUE 00199 ELSE 00200 DO 180 I = 1, N 00201 WORK( I ) = ZERO 00202 180 CONTINUE 00203 DO 200 J = 1, N 00204 DO 190 I = 1, J 00205 WORK( I ) = WORK( I ) + ABS( AP( K ) ) 00206 K = K + 1 00207 190 CONTINUE 00208 200 CONTINUE 00209 END IF 00210 ELSE 00211 IF( LSAME( DIAG, 'U' ) ) THEN 00212 DO 210 I = 1, N 00213 WORK( I ) = ONE 00214 210 CONTINUE 00215 DO 230 J = 1, N 00216 K = K + 1 00217 DO 220 I = J + 1, N 00218 WORK( I ) = WORK( I ) + ABS( AP( K ) ) 00219 K = K + 1 00220 220 CONTINUE 00221 230 CONTINUE 00222 ELSE 00223 DO 240 I = 1, N 00224 WORK( I ) = ZERO 00225 240 CONTINUE 00226 DO 260 J = 1, N 00227 DO 250 I = J, N 00228 WORK( I ) = WORK( I ) + ABS( AP( K ) ) 00229 K = K + 1 00230 250 CONTINUE 00231 260 CONTINUE 00232 END IF 00233 END IF 00234 VALUE = ZERO 00235 DO 270 I = 1, N 00236 VALUE = MAX( VALUE, WORK( I ) ) 00237 270 CONTINUE 00238 ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN 00239 * 00240 * Find normF(A). 00241 * 00242 IF( LSAME( UPLO, 'U' ) ) THEN 00243 IF( LSAME( DIAG, 'U' ) ) THEN 00244 SCALE = ONE 00245 SUM = N 00246 K = 2 00247 DO 280 J = 2, N 00248 CALL ZLASSQ( J-1, AP( K ), 1, SCALE, SUM ) 00249 K = K + J 00250 280 CONTINUE 00251 ELSE 00252 SCALE = ZERO 00253 SUM = ONE 00254 K = 1 00255 DO 290 J = 1, N 00256 CALL ZLASSQ( J, AP( K ), 1, SCALE, SUM ) 00257 K = K + J 00258 290 CONTINUE 00259 END IF 00260 ELSE 00261 IF( LSAME( DIAG, 'U' ) ) THEN 00262 SCALE = ONE 00263 SUM = N 00264 K = 2 00265 DO 300 J = 1, N - 1 00266 CALL ZLASSQ( N-J, AP( K ), 1, SCALE, SUM ) 00267 K = K + N - J + 1 00268 300 CONTINUE 00269 ELSE 00270 SCALE = ZERO 00271 SUM = ONE 00272 K = 1 00273 DO 310 J = 1, N 00274 CALL ZLASSQ( N-J+1, AP( K ), 1, SCALE, SUM ) 00275 K = K + N - J + 1 00276 310 CONTINUE 00277 END IF 00278 END IF 00279 VALUE = SCALE*SQRT( SUM ) 00280 END IF 00281 * 00282 ZLANTP = VALUE 00283 RETURN 00284 * 00285 * End of ZLANTP 00286 * 00287 END