1: DOUBLE PRECISION FUNCTION DLANTP( NORM, UPLO, DIAG, N, AP, WORK )
2: *
3: * -- LAPACK auxiliary 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 DIAG, NORM, UPLO
10: INTEGER N
11: * ..
12: * .. Array Arguments ..
13: DOUBLE PRECISION AP( * ), WORK( * )
14: * ..
15: *
16: * Purpose
17: * =======
18: *
19: * DLANTP returns the value of the one norm, or the Frobenius norm, or
20: * the infinity norm, or the element of largest absolute value of a
21: * triangular matrix A, supplied in packed form.
22: *
23: * Description
24: * ===========
25: *
26: * DLANTP returns the value
27: *
28: * DLANTP = ( max(abs(A(i,j))), NORM = 'M' or 'm'
29: * (
30: * ( norm1(A), NORM = '1', 'O' or 'o'
31: * (
32: * ( normI(A), NORM = 'I' or 'i'
33: * (
34: * ( normF(A), NORM = 'F', 'f', 'E' or 'e'
35: *
36: * where norm1 denotes the one norm of a matrix (maximum column sum),
37: * normI denotes the infinity norm of a matrix (maximum row sum) and
38: * normF denotes the Frobenius norm of a matrix (square root of sum of
39: * squares). Note that max(abs(A(i,j))) is not a consistent matrix norm.
40: *
41: * Arguments
42: * =========
43: *
44: * NORM (input) CHARACTER*1
45: * Specifies the value to be returned in DLANTP as described
46: * above.
47: *
48: * UPLO (input) CHARACTER*1
49: * Specifies whether the matrix A is upper or lower triangular.
50: * = 'U': Upper triangular
51: * = 'L': Lower triangular
52: *
53: * DIAG (input) CHARACTER*1
54: * Specifies whether or not the matrix A is unit triangular.
55: * = 'N': Non-unit triangular
56: * = 'U': Unit triangular
57: *
58: * N (input) INTEGER
59: * The order of the matrix A. N >= 0. When N = 0, DLANTP is
60: * set to zero.
61: *
62: * AP (input) DOUBLE PRECISION array, dimension (N*(N+1)/2)
63: * The upper or lower triangular matrix A, packed columnwise in
64: * a linear array. The j-th column of A is stored in the array
65: * AP as follows:
66: * if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
67: * if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n.
68: * Note that when DIAG = 'U', the elements of the array AP
69: * corresponding to the diagonal elements of the matrix A are
70: * not referenced, but are assumed to be one.
71: *
72: * WORK (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK)),
73: * where LWORK >= N when NORM = 'I'; otherwise, WORK is not
74: * referenced.
75: *
76: * =====================================================================
77: *
78: * .. Parameters ..
79: DOUBLE PRECISION ONE, ZERO
80: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
81: * ..
82: * .. Local Scalars ..
83: LOGICAL UDIAG
84: INTEGER I, J, K
85: DOUBLE PRECISION SCALE, SUM, VALUE
86: * ..
87: * .. External Subroutines ..
88: EXTERNAL DLASSQ
89: * ..
90: * .. External Functions ..
91: LOGICAL LSAME
92: EXTERNAL LSAME
93: * ..
94: * .. Intrinsic Functions ..
95: INTRINSIC ABS, MAX, SQRT
96: * ..
97: * .. Executable Statements ..
98: *
99: IF( N.EQ.0 ) THEN
100: VALUE = ZERO
101: ELSE IF( LSAME( NORM, 'M' ) ) THEN
102: *
103: * Find max(abs(A(i,j))).
104: *
105: K = 1
106: IF( LSAME( DIAG, 'U' ) ) THEN
107: VALUE = ONE
108: IF( LSAME( UPLO, 'U' ) ) THEN
109: DO 20 J = 1, N
110: DO 10 I = K, K + J - 2
111: VALUE = MAX( VALUE, ABS( AP( I ) ) )
112: 10 CONTINUE
113: K = K + J
114: 20 CONTINUE
115: ELSE
116: DO 40 J = 1, N
117: DO 30 I = K + 1, K + N - J
118: VALUE = MAX( VALUE, ABS( AP( I ) ) )
119: 30 CONTINUE
120: K = K + N - J + 1
121: 40 CONTINUE
122: END IF
123: ELSE
124: VALUE = ZERO
125: IF( LSAME( UPLO, 'U' ) ) THEN
126: DO 60 J = 1, N
127: DO 50 I = K, K + J - 1
128: VALUE = MAX( VALUE, ABS( AP( I ) ) )
129: 50 CONTINUE
130: K = K + J
131: 60 CONTINUE
132: ELSE
133: DO 80 J = 1, N
134: DO 70 I = K, K + N - J
135: VALUE = MAX( VALUE, ABS( AP( I ) ) )
136: 70 CONTINUE
137: K = K + N - J + 1
138: 80 CONTINUE
139: END IF
140: END IF
141: ELSE IF( ( LSAME( NORM, 'O' ) ) .OR. ( NORM.EQ.'1' ) ) THEN
142: *
143: * Find norm1(A).
144: *
145: VALUE = ZERO
146: K = 1
147: UDIAG = LSAME( DIAG, 'U' )
148: IF( LSAME( UPLO, 'U' ) ) THEN
149: DO 110 J = 1, N
150: IF( UDIAG ) THEN
151: SUM = ONE
152: DO 90 I = K, K + J - 2
153: SUM = SUM + ABS( AP( I ) )
154: 90 CONTINUE
155: ELSE
156: SUM = ZERO
157: DO 100 I = K, K + J - 1
158: SUM = SUM + ABS( AP( I ) )
159: 100 CONTINUE
160: END IF
161: K = K + J
162: VALUE = MAX( VALUE, SUM )
163: 110 CONTINUE
164: ELSE
165: DO 140 J = 1, N
166: IF( UDIAG ) THEN
167: SUM = ONE
168: DO 120 I = K + 1, K + N - J
169: SUM = SUM + ABS( AP( I ) )
170: 120 CONTINUE
171: ELSE
172: SUM = ZERO
173: DO 130 I = K, K + N - J
174: SUM = SUM + ABS( AP( I ) )
175: 130 CONTINUE
176: END IF
177: K = K + N - J + 1
178: VALUE = MAX( VALUE, SUM )
179: 140 CONTINUE
180: END IF
181: ELSE IF( LSAME( NORM, 'I' ) ) THEN
182: *
183: * Find normI(A).
184: *
185: K = 1
186: IF( LSAME( UPLO, 'U' ) ) THEN
187: IF( LSAME( DIAG, 'U' ) ) THEN
188: DO 150 I = 1, N
189: WORK( I ) = ONE
190: 150 CONTINUE
191: DO 170 J = 1, N
192: DO 160 I = 1, J - 1
193: WORK( I ) = WORK( I ) + ABS( AP( K ) )
194: K = K + 1
195: 160 CONTINUE
196: K = K + 1
197: 170 CONTINUE
198: ELSE
199: DO 180 I = 1, N
200: WORK( I ) = ZERO
201: 180 CONTINUE
202: DO 200 J = 1, N
203: DO 190 I = 1, J
204: WORK( I ) = WORK( I ) + ABS( AP( K ) )
205: K = K + 1
206: 190 CONTINUE
207: 200 CONTINUE
208: END IF
209: ELSE
210: IF( LSAME( DIAG, 'U' ) ) THEN
211: DO 210 I = 1, N
212: WORK( I ) = ONE
213: 210 CONTINUE
214: DO 230 J = 1, N
215: K = K + 1
216: DO 220 I = J + 1, N
217: WORK( I ) = WORK( I ) + ABS( AP( K ) )
218: K = K + 1
219: 220 CONTINUE
220: 230 CONTINUE
221: ELSE
222: DO 240 I = 1, N
223: WORK( I ) = ZERO
224: 240 CONTINUE
225: DO 260 J = 1, N
226: DO 250 I = J, N
227: WORK( I ) = WORK( I ) + ABS( AP( K ) )
228: K = K + 1
229: 250 CONTINUE
230: 260 CONTINUE
231: END IF
232: END IF
233: VALUE = ZERO
234: DO 270 I = 1, N
235: VALUE = MAX( VALUE, WORK( I ) )
236: 270 CONTINUE
237: ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
238: *
239: * Find normF(A).
240: *
241: IF( LSAME( UPLO, 'U' ) ) THEN
242: IF( LSAME( DIAG, 'U' ) ) THEN
243: SCALE = ONE
244: SUM = N
245: K = 2
246: DO 280 J = 2, N
247: CALL DLASSQ( J-1, AP( K ), 1, SCALE, SUM )
248: K = K + J
249: 280 CONTINUE
250: ELSE
251: SCALE = ZERO
252: SUM = ONE
253: K = 1
254: DO 290 J = 1, N
255: CALL DLASSQ( J, AP( K ), 1, SCALE, SUM )
256: K = K + J
257: 290 CONTINUE
258: END IF
259: ELSE
260: IF( LSAME( DIAG, 'U' ) ) THEN
261: SCALE = ONE
262: SUM = N
263: K = 2
264: DO 300 J = 1, N - 1
265: CALL DLASSQ( N-J, AP( K ), 1, SCALE, SUM )
266: K = K + N - J + 1
267: 300 CONTINUE
268: ELSE
269: SCALE = ZERO
270: SUM = ONE
271: K = 1
272: DO 310 J = 1, N
273: CALL DLASSQ( N-J+1, AP( K ), 1, SCALE, SUM )
274: K = K + N - J + 1
275: 310 CONTINUE
276: END IF
277: END IF
278: VALUE = SCALE*SQRT( SUM )
279: END IF
280: *
281: DLANTP = VALUE
282: RETURN
283: *
284: * End of DLANTP
285: *
286: END
CVSweb interface <joel.bertrand@systella.fr>