|
2329
|
1 SUBROUTINE DTRMV ( UPLO, TRANS, DIAG, N, A, LDA, X, INCX ) |
|
|
2 * .. Scalar Arguments .. |
|
|
3 INTEGER INCX, LDA, N |
|
|
4 CHARACTER*1 DIAG, TRANS, UPLO |
|
|
5 * .. Array Arguments .. |
|
|
6 DOUBLE PRECISION A( LDA, * ), X( * ) |
|
|
7 * .. |
|
|
8 * |
|
|
9 * Purpose |
|
|
10 * ======= |
|
|
11 * |
|
|
12 * DTRMV performs one of the matrix-vector operations |
|
|
13 * |
|
|
14 * x := A*x, or x := A'*x, |
|
|
15 * |
|
|
16 * where x is an n element vector and A is an n by n unit, or non-unit, |
|
|
17 * upper or lower triangular matrix. |
|
|
18 * |
|
|
19 * Parameters |
|
|
20 * ========== |
|
|
21 * |
|
|
22 * UPLO - CHARACTER*1. |
|
|
23 * On entry, UPLO specifies whether the matrix is an upper or |
|
|
24 * lower triangular matrix as follows: |
|
|
25 * |
|
|
26 * UPLO = 'U' or 'u' A is an upper triangular matrix. |
|
|
27 * |
|
|
28 * UPLO = 'L' or 'l' A is a lower triangular matrix. |
|
|
29 * |
|
|
30 * Unchanged on exit. |
|
|
31 * |
|
|
32 * TRANS - CHARACTER*1. |
|
|
33 * On entry, TRANS specifies the operation to be performed as |
|
|
34 * follows: |
|
|
35 * |
|
|
36 * TRANS = 'N' or 'n' x := A*x. |
|
|
37 * |
|
|
38 * TRANS = 'T' or 't' x := A'*x. |
|
|
39 * |
|
|
40 * TRANS = 'C' or 'c' x := A'*x. |
|
|
41 * |
|
|
42 * Unchanged on exit. |
|
|
43 * |
|
|
44 * DIAG - CHARACTER*1. |
|
|
45 * On entry, DIAG specifies whether or not A is unit |
|
|
46 * triangular as follows: |
|
|
47 * |
|
|
48 * DIAG = 'U' or 'u' A is assumed to be unit triangular. |
|
|
49 * |
|
|
50 * DIAG = 'N' or 'n' A is not assumed to be unit |
|
|
51 * triangular. |
|
|
52 * |
|
|
53 * Unchanged on exit. |
|
|
54 * |
|
|
55 * N - INTEGER. |
|
|
56 * On entry, N specifies the order of the matrix A. |
|
|
57 * N must be at least zero. |
|
|
58 * Unchanged on exit. |
|
|
59 * |
|
|
60 * A - DOUBLE PRECISION array of DIMENSION ( LDA, n ). |
|
|
61 * Before entry with UPLO = 'U' or 'u', the leading n by n |
|
|
62 * upper triangular part of the array A must contain the upper |
|
|
63 * triangular matrix and the strictly lower triangular part of |
|
|
64 * A is not referenced. |
|
|
65 * Before entry with UPLO = 'L' or 'l', the leading n by n |
|
|
66 * lower triangular part of the array A must contain the lower |
|
|
67 * triangular matrix and the strictly upper triangular part of |
|
|
68 * A is not referenced. |
|
|
69 * Note that when DIAG = 'U' or 'u', the diagonal elements of |
|
|
70 * A are not referenced either, but are assumed to be unity. |
|
|
71 * Unchanged on exit. |
|
|
72 * |
|
|
73 * LDA - INTEGER. |
|
|
74 * On entry, LDA specifies the first dimension of A as declared |
|
|
75 * in the calling (sub) program. LDA must be at least |
|
|
76 * max( 1, n ). |
|
|
77 * Unchanged on exit. |
|
|
78 * |
|
|
79 * X - DOUBLE PRECISION array of dimension at least |
|
|
80 * ( 1 + ( n - 1 )*abs( INCX ) ). |
|
|
81 * Before entry, the incremented array X must contain the n |
|
|
82 * element vector x. On exit, X is overwritten with the |
|
|
83 * tranformed vector x. |
|
|
84 * |
|
|
85 * INCX - INTEGER. |
|
|
86 * On entry, INCX specifies the increment for the elements of |
|
|
87 * X. INCX must not be zero. |
|
|
88 * Unchanged on exit. |
|
|
89 * |
|
|
90 * |
|
|
91 * Level 2 Blas routine. |
|
|
92 * |
|
|
93 * -- Written on 22-October-1986. |
|
|
94 * Jack Dongarra, Argonne National Lab. |
|
|
95 * Jeremy Du Croz, Nag Central Office. |
|
|
96 * Sven Hammarling, Nag Central Office. |
|
|
97 * Richard Hanson, Sandia National Labs. |
|
|
98 * |
|
|
99 * |
|
|
100 * .. Parameters .. |
|
|
101 DOUBLE PRECISION ZERO |
|
|
102 PARAMETER ( ZERO = 0.0D+0 ) |
|
|
103 * .. Local Scalars .. |
|
|
104 DOUBLE PRECISION TEMP |
|
|
105 INTEGER I, INFO, IX, J, JX, KX |
|
|
106 LOGICAL NOUNIT |
|
|
107 * .. External Functions .. |
|
|
108 LOGICAL LSAME |
|
|
109 EXTERNAL LSAME |
|
|
110 * .. External Subroutines .. |
|
|
111 EXTERNAL XERBLA |
|
|
112 * .. Intrinsic Functions .. |
|
|
113 INTRINSIC MAX |
|
|
114 * .. |
|
|
115 * .. Executable Statements .. |
|
|
116 * |
|
|
117 * Test the input parameters. |
|
|
118 * |
|
|
119 INFO = 0 |
|
|
120 IF ( .NOT.LSAME( UPLO , 'U' ).AND. |
|
|
121 $ .NOT.LSAME( UPLO , 'L' ) )THEN |
|
|
122 INFO = 1 |
|
|
123 ELSE IF( .NOT.LSAME( TRANS, 'N' ).AND. |
|
|
124 $ .NOT.LSAME( TRANS, 'T' ).AND. |
|
|
125 $ .NOT.LSAME( TRANS, 'C' ) )THEN |
|
|
126 INFO = 2 |
|
|
127 ELSE IF( .NOT.LSAME( DIAG , 'U' ).AND. |
|
|
128 $ .NOT.LSAME( DIAG , 'N' ) )THEN |
|
|
129 INFO = 3 |
|
|
130 ELSE IF( N.LT.0 )THEN |
|
|
131 INFO = 4 |
|
|
132 ELSE IF( LDA.LT.MAX( 1, N ) )THEN |
|
|
133 INFO = 6 |
|
|
134 ELSE IF( INCX.EQ.0 )THEN |
|
|
135 INFO = 8 |
|
|
136 END IF |
|
|
137 IF( INFO.NE.0 )THEN |
|
|
138 CALL XERBLA( 'DTRMV ', INFO ) |
|
|
139 RETURN |
|
|
140 END IF |
|
|
141 * |
|
|
142 * Quick return if possible. |
|
|
143 * |
|
|
144 IF( N.EQ.0 ) |
|
|
145 $ RETURN |
|
|
146 * |
|
|
147 NOUNIT = LSAME( DIAG, 'N' ) |
|
|
148 * |
|
|
149 * Set up the start point in X if the increment is not unity. This |
|
|
150 * will be ( N - 1 )*INCX too small for descending loops. |
|
|
151 * |
|
|
152 IF( INCX.LE.0 )THEN |
|
|
153 KX = 1 - ( N - 1 )*INCX |
|
|
154 ELSE IF( INCX.NE.1 )THEN |
|
|
155 KX = 1 |
|
|
156 END IF |
|
|
157 * |
|
|
158 * Start the operations. In this version the elements of A are |
|
|
159 * accessed sequentially with one pass through A. |
|
|
160 * |
|
|
161 IF( LSAME( TRANS, 'N' ) )THEN |
|
|
162 * |
|
|
163 * Form x := A*x. |
|
|
164 * |
|
|
165 IF( LSAME( UPLO, 'U' ) )THEN |
|
|
166 IF( INCX.EQ.1 )THEN |
|
|
167 DO 20, J = 1, N |
|
|
168 IF( X( J ).NE.ZERO )THEN |
|
|
169 TEMP = X( J ) |
|
|
170 DO 10, I = 1, J - 1 |
|
|
171 X( I ) = X( I ) + TEMP*A( I, J ) |
|
|
172 10 CONTINUE |
|
|
173 IF( NOUNIT ) |
|
|
174 $ X( J ) = X( J )*A( J, J ) |
|
|
175 END IF |
|
|
176 20 CONTINUE |
|
|
177 ELSE |
|
|
178 JX = KX |
|
|
179 DO 40, J = 1, N |
|
|
180 IF( X( JX ).NE.ZERO )THEN |
|
|
181 TEMP = X( JX ) |
|
|
182 IX = KX |
|
|
183 DO 30, I = 1, J - 1 |
|
|
184 X( IX ) = X( IX ) + TEMP*A( I, J ) |
|
|
185 IX = IX + INCX |
|
|
186 30 CONTINUE |
|
|
187 IF( NOUNIT ) |
|
|
188 $ X( JX ) = X( JX )*A( J, J ) |
|
|
189 END IF |
|
|
190 JX = JX + INCX |
|
|
191 40 CONTINUE |
|
|
192 END IF |
|
|
193 ELSE |
|
|
194 IF( INCX.EQ.1 )THEN |
|
|
195 DO 60, J = N, 1, -1 |
|
|
196 IF( X( J ).NE.ZERO )THEN |
|
|
197 TEMP = X( J ) |
|
|
198 DO 50, I = N, J + 1, -1 |
|
|
199 X( I ) = X( I ) + TEMP*A( I, J ) |
|
|
200 50 CONTINUE |
|
|
201 IF( NOUNIT ) |
|
|
202 $ X( J ) = X( J )*A( J, J ) |
|
|
203 END IF |
|
|
204 60 CONTINUE |
|
|
205 ELSE |
|
|
206 KX = KX + ( N - 1 )*INCX |
|
|
207 JX = KX |
|
|
208 DO 80, J = N, 1, -1 |
|
|
209 IF( X( JX ).NE.ZERO )THEN |
|
|
210 TEMP = X( JX ) |
|
|
211 IX = KX |
|
|
212 DO 70, I = N, J + 1, -1 |
|
|
213 X( IX ) = X( IX ) + TEMP*A( I, J ) |
|
|
214 IX = IX - INCX |
|
|
215 70 CONTINUE |
|
|
216 IF( NOUNIT ) |
|
|
217 $ X( JX ) = X( JX )*A( J, J ) |
|
|
218 END IF |
|
|
219 JX = JX - INCX |
|
|
220 80 CONTINUE |
|
|
221 END IF |
|
|
222 END IF |
|
|
223 ELSE |
|
|
224 * |
|
|
225 * Form x := A'*x. |
|
|
226 * |
|
|
227 IF( LSAME( UPLO, 'U' ) )THEN |
|
|
228 IF( INCX.EQ.1 )THEN |
|
|
229 DO 100, J = N, 1, -1 |
|
|
230 TEMP = X( J ) |
|
|
231 IF( NOUNIT ) |
|
|
232 $ TEMP = TEMP*A( J, J ) |
|
|
233 DO 90, I = J - 1, 1, -1 |
|
|
234 TEMP = TEMP + A( I, J )*X( I ) |
|
|
235 90 CONTINUE |
|
|
236 X( J ) = TEMP |
|
|
237 100 CONTINUE |
|
|
238 ELSE |
|
|
239 JX = KX + ( N - 1 )*INCX |
|
|
240 DO 120, J = N, 1, -1 |
|
|
241 TEMP = X( JX ) |
|
|
242 IX = JX |
|
|
243 IF( NOUNIT ) |
|
|
244 $ TEMP = TEMP*A( J, J ) |
|
|
245 DO 110, I = J - 1, 1, -1 |
|
|
246 IX = IX - INCX |
|
|
247 TEMP = TEMP + A( I, J )*X( IX ) |
|
|
248 110 CONTINUE |
|
|
249 X( JX ) = TEMP |
|
|
250 JX = JX - INCX |
|
|
251 120 CONTINUE |
|
|
252 END IF |
|
|
253 ELSE |
|
|
254 IF( INCX.EQ.1 )THEN |
|
|
255 DO 140, J = 1, N |
|
|
256 TEMP = X( J ) |
|
|
257 IF( NOUNIT ) |
|
|
258 $ TEMP = TEMP*A( J, J ) |
|
|
259 DO 130, I = J + 1, N |
|
|
260 TEMP = TEMP + A( I, J )*X( I ) |
|
|
261 130 CONTINUE |
|
|
262 X( J ) = TEMP |
|
|
263 140 CONTINUE |
|
|
264 ELSE |
|
|
265 JX = KX |
|
|
266 DO 160, J = 1, N |
|
|
267 TEMP = X( JX ) |
|
|
268 IX = JX |
|
|
269 IF( NOUNIT ) |
|
|
270 $ TEMP = TEMP*A( J, J ) |
|
|
271 DO 150, I = J + 1, N |
|
|
272 IX = IX + INCX |
|
|
273 TEMP = TEMP + A( I, J )*X( IX ) |
|
|
274 150 CONTINUE |
|
|
275 X( JX ) = TEMP |
|
|
276 JX = JX + INCX |
|
|
277 160 CONTINUE |
|
|
278 END IF |
|
|
279 END IF |
|
|
280 END IF |
|
|
281 * |
|
|
282 RETURN |
|
|
283 * |
|
|
284 * End of DTRMV . |
|
|
285 * |
|
|
286 END |
