1: DOUBLE PRECISION FUNCTION DLANSY( NORM, UPLO, N, A, LDA, 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 NORM, UPLO
10: INTEGER LDA, N
11: * ..
12: * .. Array Arguments ..
13: DOUBLE PRECISION A( LDA, * ), WORK( * )
14: * ..
15: *
16: * Purpose
17: * =======
18: *
19: * DLANSY 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: * real symmetric matrix A.
22: *
23: * Description
24: * ===========
25: *
26: * DLANSY returns the value
27: *
28: * DLANSY = ( 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 DLANSY as described
46: * above.
47: *
48: * UPLO (input) CHARACTER*1
49: * Specifies whether the upper or lower triangular part of the
50: * symmetric matrix A is to be referenced.
51: * = 'U': Upper triangular part of A is referenced
52: * = 'L': Lower triangular part of A is referenced
53: *
54: * N (input) INTEGER
55: * The order of the matrix A. N >= 0. When N = 0, DLANSY is
56: * set to zero.
57: *
58: * A (input) DOUBLE PRECISION array, dimension (LDA,N)
59: * The symmetric matrix A. If UPLO = 'U', the leading n by n
60: * upper triangular part of A contains the upper triangular part
61: * of the matrix A, and the strictly lower triangular part of A
62: * is not referenced. If UPLO = 'L', the leading n by n lower
63: * triangular part of A contains the lower triangular part of
64: * the matrix A, and the strictly upper triangular part of A is
65: * not referenced.
66: *
67: * LDA (input) INTEGER
68: * The leading dimension of the array A. LDA >= max(N,1).
69: *
70: * WORK (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK)),
71: * where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise,
72: * WORK is not referenced.
73: *
74: * =====================================================================
75: *
76: * .. Parameters ..
77: DOUBLE PRECISION ONE, ZERO
78: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
79: * ..
80: * .. Local Scalars ..
81: INTEGER I, J
82: DOUBLE PRECISION ABSA, SCALE, SUM, VALUE
83: * ..
84: * .. External Subroutines ..
85: EXTERNAL DLASSQ
86: * ..
87: * .. External Functions ..
88: LOGICAL LSAME
89: EXTERNAL LSAME
90: * ..
91: * .. Intrinsic Functions ..
92: INTRINSIC ABS, MAX, SQRT
93: * ..
94: * .. Executable Statements ..
95: *
96: IF( N.EQ.0 ) THEN
97: VALUE = ZERO
98: ELSE IF( LSAME( NORM, 'M' ) ) THEN
99: *
100: * Find max(abs(A(i,j))).
101: *
102: VALUE = ZERO
103: IF( LSAME( UPLO, 'U' ) ) THEN
104: DO 20 J = 1, N
105: DO 10 I = 1, J
106: VALUE = MAX( VALUE, ABS( A( I, J ) ) )
107: 10 CONTINUE
108: 20 CONTINUE
109: ELSE
110: DO 40 J = 1, N
111: DO 30 I = J, N
112: VALUE = MAX( VALUE, ABS( A( I, J ) ) )
113: 30 CONTINUE
114: 40 CONTINUE
115: END IF
116: ELSE IF( ( LSAME( NORM, 'I' ) ) .OR. ( LSAME( NORM, 'O' ) ) .OR.
117: $ ( NORM.EQ.'1' ) ) THEN
118: *
119: * Find normI(A) ( = norm1(A), since A is symmetric).
120: *
121: VALUE = ZERO
122: IF( LSAME( UPLO, 'U' ) ) THEN
123: DO 60 J = 1, N
124: SUM = ZERO
125: DO 50 I = 1, J - 1
126: ABSA = ABS( A( I, J ) )
127: SUM = SUM + ABSA
128: WORK( I ) = WORK( I ) + ABSA
129: 50 CONTINUE
130: WORK( J ) = SUM + ABS( A( J, J ) )
131: 60 CONTINUE
132: DO 70 I = 1, N
133: VALUE = MAX( VALUE, WORK( I ) )
134: 70 CONTINUE
135: ELSE
136: DO 80 I = 1, N
137: WORK( I ) = ZERO
138: 80 CONTINUE
139: DO 100 J = 1, N
140: SUM = WORK( J ) + ABS( A( J, J ) )
141: DO 90 I = J + 1, N
142: ABSA = ABS( A( I, J ) )
143: SUM = SUM + ABSA
144: WORK( I ) = WORK( I ) + ABSA
145: 90 CONTINUE
146: VALUE = MAX( VALUE, SUM )
147: 100 CONTINUE
148: END IF
149: ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
150: *
151: * Find normF(A).
152: *
153: SCALE = ZERO
154: SUM = ONE
155: IF( LSAME( UPLO, 'U' ) ) THEN
156: DO 110 J = 2, N
157: CALL DLASSQ( J-1, A( 1, J ), 1, SCALE, SUM )
158: 110 CONTINUE
159: ELSE
160: DO 120 J = 1, N - 1
161: CALL DLASSQ( N-J, A( J+1, J ), 1, SCALE, SUM )
162: 120 CONTINUE
163: END IF
164: SUM = 2*SUM
165: CALL DLASSQ( N, A, LDA+1, SCALE, SUM )
166: VALUE = SCALE*SQRT( SUM )
167: END IF
168: *
169: DLANSY = VALUE
170: RETURN
171: *
172: * End of DLANSY
173: *
174: END
CVSweb interface <joel.bertrand@systella.fr>