![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dla_syamv.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack|
Return to dla_syamv.f CVS log | Up to [local] / rpl / lapack / lapack |
1.8 bertrand 1: *> \brief \b DLA_SYAMV computes a matrix-vector product using a symmetric indefinite matrix to calculate error bounds.
1.5 bertrand 2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download DLA_SYAMV + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dla_syamv.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dla_syamv.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dla_syamv.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE DLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y,
22: * INCY )
23: *
24: * .. Scalar Arguments ..
25: * DOUBLE PRECISION ALPHA, BETA
26: * INTEGER INCX, INCY, LDA, N, UPLO
27: * ..
28: * .. Array Arguments ..
29: * DOUBLE PRECISION A( LDA, * ), X( * ), Y( * )
30: * ..
31: *
32: *
33: *> \par Purpose:
34: * =============
35: *>
36: *> \verbatim
37: *>
38: *> DLA_SYAMV performs the matrix-vector operation
39: *>
40: *> y := alpha*abs(A)*abs(x) + beta*abs(y),
41: *>
42: *> where alpha and beta are scalars, x and y are vectors and A is an
43: *> n by n symmetric matrix.
44: *>
45: *> This function is primarily used in calculating error bounds.
46: *> To protect against underflow during evaluation, components in
47: *> the resulting vector are perturbed away from zero by (N+1)
48: *> times the underflow threshold. To prevent unnecessarily large
49: *> errors for block-structure embedded in general matrices,
50: *> "symbolically" zero components are not perturbed. A zero
51: *> entry is considered "symbolic" if all multiplications involved
52: *> in computing that entry have at least one zero multiplicand.
53: *> \endverbatim
54: *
55: * Arguments:
56: * ==========
57: *
58: *> \param[in] UPLO
59: *> \verbatim
60: *> UPLO is INTEGER
61: *> On entry, UPLO specifies whether the upper or lower
62: *> triangular part of the array A is to be referenced as
63: *> follows:
64: *>
65: *> UPLO = BLAS_UPPER Only the upper triangular part of A
66: *> is to be referenced.
67: *>
68: *> UPLO = BLAS_LOWER Only the lower triangular part of A
69: *> is to be referenced.
70: *>
71: *> Unchanged on exit.
72: *> \endverbatim
73: *>
74: *> \param[in] N
75: *> \verbatim
76: *> N is INTEGER
77: *> On entry, N specifies the number of columns of the matrix A.
78: *> N must be at least zero.
79: *> Unchanged on exit.
80: *> \endverbatim
81: *>
82: *> \param[in] ALPHA
83: *> \verbatim
84: *> ALPHA is DOUBLE PRECISION .
85: *> On entry, ALPHA specifies the scalar alpha.
86: *> Unchanged on exit.
87: *> \endverbatim
88: *>
89: *> \param[in] A
90: *> \verbatim
91: *> A is DOUBLE PRECISION array of DIMENSION ( LDA, n ).
92: *> Before entry, the leading m by n part of the array A must
93: *> contain the matrix of coefficients.
94: *> Unchanged on exit.
95: *> \endverbatim
96: *>
97: *> \param[in] LDA
98: *> \verbatim
99: *> LDA is INTEGER
100: *> On entry, LDA specifies the first dimension of A as declared
101: *> in the calling (sub) program. LDA must be at least
102: *> max( 1, n ).
103: *> Unchanged on exit.
104: *> \endverbatim
105: *>
106: *> \param[in] X
107: *> \verbatim
108: *> X is DOUBLE PRECISION array, dimension
109: *> ( 1 + ( n - 1 )*abs( INCX ) )
110: *> Before entry, the incremented array X must contain the
111: *> vector x.
112: *> Unchanged on exit.
113: *> \endverbatim
114: *>
115: *> \param[in] INCX
116: *> \verbatim
117: *> INCX is INTEGER
118: *> On entry, INCX specifies the increment for the elements of
119: *> X. INCX must not be zero.
120: *> Unchanged on exit.
121: *> \endverbatim
122: *>
123: *> \param[in] BETA
124: *> \verbatim
125: *> BETA is DOUBLE PRECISION .
126: *> On entry, BETA specifies the scalar beta. When BETA is
127: *> supplied as zero then Y need not be set on input.
128: *> Unchanged on exit.
129: *> \endverbatim
130: *>
131: *> \param[in,out] Y
132: *> \verbatim
133: *> Y is DOUBLE PRECISION array, dimension
134: *> ( 1 + ( n - 1 )*abs( INCY ) )
135: *> Before entry with BETA non-zero, the incremented array Y
136: *> must contain the vector y. On exit, Y is overwritten by the
137: *> updated vector y.
138: *> \endverbatim
139: *>
140: *> \param[in] INCY
141: *> \verbatim
142: *> INCY is INTEGER
143: *> On entry, INCY specifies the increment for the elements of
144: *> Y. INCY must not be zero.
145: *> Unchanged on exit.
146: *> \endverbatim
147: *
148: * Authors:
149: * ========
150: *
151: *> \author Univ. of Tennessee
152: *> \author Univ. of California Berkeley
153: *> \author Univ. of Colorado Denver
154: *> \author NAG Ltd.
155: *
1.8 bertrand 156: *> \date September 2012
1.5 bertrand 157: *
158: *> \ingroup doubleSYcomputational
159: *
160: *> \par Further Details:
161: * =====================
162: *>
163: *> \verbatim
164: *>
165: *> Level 2 Blas routine.
166: *>
167: *> -- Written on 22-October-1986.
168: *> Jack Dongarra, Argonne National Lab.
169: *> Jeremy Du Croz, Nag Central Office.
170: *> Sven Hammarling, Nag Central Office.
171: *> Richard Hanson, Sandia National Labs.
172: *> -- Modified for the absolute-value product, April 2006
173: *> Jason Riedy, UC Berkeley
174: *> \endverbatim
175: *>
176: * =====================================================================
1.1 bertrand 177: SUBROUTINE DLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y,
178: $ INCY )
179: *
1.8 bertrand 180: * -- LAPACK computational routine (version 3.4.2) --
1.5 bertrand 181: * -- LAPACK is a software package provided by Univ. of Tennessee, --
182: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8 bertrand 183: * September 2012
1.1 bertrand 184: *
185: * .. Scalar Arguments ..
186: DOUBLE PRECISION ALPHA, BETA
187: INTEGER INCX, INCY, LDA, N, UPLO
188: * ..
189: * .. Array Arguments ..
190: DOUBLE PRECISION A( LDA, * ), X( * ), Y( * )
191: * ..
192: *
193: * =====================================================================
194: *
195: * .. Parameters ..
196: DOUBLE PRECISION ONE, ZERO
197: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
198: * ..
199: * .. Local Scalars ..
200: LOGICAL SYMB_ZERO
201: DOUBLE PRECISION TEMP, SAFE1
202: INTEGER I, INFO, IY, J, JX, KX, KY
203: * ..
204: * .. External Subroutines ..
205: EXTERNAL XERBLA, DLAMCH
206: DOUBLE PRECISION DLAMCH
207: * ..
208: * .. External Functions ..
209: EXTERNAL ILAUPLO
210: INTEGER ILAUPLO
211: * ..
212: * .. Intrinsic Functions ..
213: INTRINSIC MAX, ABS, SIGN
214: * ..
215: * .. Executable Statements ..
216: *
217: * Test the input parameters.
218: *
219: INFO = 0
220: IF ( UPLO.NE.ILAUPLO( 'U' ) .AND.
221: $ UPLO.NE.ILAUPLO( 'L' ) ) THEN
222: INFO = 1
223: ELSE IF( N.LT.0 )THEN
224: INFO = 2
225: ELSE IF( LDA.LT.MAX( 1, N ) )THEN
226: INFO = 5
227: ELSE IF( INCX.EQ.0 )THEN
228: INFO = 7
229: ELSE IF( INCY.EQ.0 )THEN
230: INFO = 10
231: END IF
232: IF( INFO.NE.0 )THEN
233: CALL XERBLA( 'DSYMV ', INFO )
234: RETURN
235: END IF
236: *
237: * Quick return if possible.
238: *
239: IF( ( N.EQ.0 ).OR.( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
240: $ RETURN
241: *
242: * Set up the start points in X and Y.
243: *
244: IF( INCX.GT.0 )THEN
245: KX = 1
246: ELSE
247: KX = 1 - ( N - 1 )*INCX
248: END IF
249: IF( INCY.GT.0 )THEN
250: KY = 1
251: ELSE
252: KY = 1 - ( N - 1 )*INCY
253: END IF
254: *
255: * Set SAFE1 essentially to be the underflow threshold times the
256: * number of additions in each row.
257: *
258: SAFE1 = DLAMCH( 'Safe minimum' )
259: SAFE1 = (N+1)*SAFE1
260: *
261: * Form y := alpha*abs(A)*abs(x) + beta*abs(y).
262: *
263: * The O(N^2) SYMB_ZERO tests could be replaced by O(N) queries to
264: * the inexact flag. Still doesn't help change the iteration order
265: * to per-column.
266: *
267: IY = KY
268: IF ( INCX.EQ.1 ) THEN
269: IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN
270: DO I = 1, N
271: IF ( BETA .EQ. ZERO ) THEN
272: SYMB_ZERO = .TRUE.
273: Y( IY ) = 0.0D+0
274: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
275: SYMB_ZERO = .TRUE.
276: ELSE
277: SYMB_ZERO = .FALSE.
278: Y( IY ) = BETA * ABS( Y( IY ) )
279: END IF
280: IF ( ALPHA .NE. ZERO ) THEN
281: DO J = 1, I
282: TEMP = ABS( A( J, I ) )
283: SYMB_ZERO = SYMB_ZERO .AND.
284: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
285:
286: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
287: END DO
288: DO J = I+1, N
289: TEMP = ABS( A( I, J ) )
290: SYMB_ZERO = SYMB_ZERO .AND.
291: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
292:
293: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
294: END DO
295: END IF
296:
297: IF ( .NOT.SYMB_ZERO )
298: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
299:
300: IY = IY + INCY
301: END DO
302: ELSE
303: DO I = 1, N
304: IF ( BETA .EQ. ZERO ) THEN
305: SYMB_ZERO = .TRUE.
306: Y( IY ) = 0.0D+0
307: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
308: SYMB_ZERO = .TRUE.
309: ELSE
310: SYMB_ZERO = .FALSE.
311: Y( IY ) = BETA * ABS( Y( IY ) )
312: END IF
313: IF ( ALPHA .NE. ZERO ) THEN
314: DO J = 1, I
315: TEMP = ABS( A( I, J ) )
316: SYMB_ZERO = SYMB_ZERO .AND.
317: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
318:
319: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
320: END DO
321: DO J = I+1, N
322: TEMP = ABS( A( J, I ) )
323: SYMB_ZERO = SYMB_ZERO .AND.
324: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
325:
326: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
327: END DO
328: END IF
329:
330: IF ( .NOT.SYMB_ZERO )
331: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
332:
333: IY = IY + INCY
334: END DO
335: END IF
336: ELSE
337: IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN
338: DO I = 1, N
339: IF ( BETA .EQ. ZERO ) THEN
340: SYMB_ZERO = .TRUE.
341: Y( IY ) = 0.0D+0
342: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
343: SYMB_ZERO = .TRUE.
344: ELSE
345: SYMB_ZERO = .FALSE.
346: Y( IY ) = BETA * ABS( Y( IY ) )
347: END IF
348: JX = KX
349: IF ( ALPHA .NE. ZERO ) THEN
350: DO J = 1, I
351: TEMP = ABS( A( J, I ) )
352: SYMB_ZERO = SYMB_ZERO .AND.
353: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
354:
355: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
356: JX = JX + INCX
357: END DO
358: DO J = I+1, N
359: TEMP = ABS( A( I, J ) )
360: SYMB_ZERO = SYMB_ZERO .AND.
361: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
362:
363: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
364: JX = JX + INCX
365: END DO
366: END IF
367:
368: IF ( .NOT.SYMB_ZERO )
369: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
370:
371: IY = IY + INCY
372: END DO
373: ELSE
374: DO I = 1, N
375: IF ( BETA .EQ. ZERO ) THEN
376: SYMB_ZERO = .TRUE.
377: Y( IY ) = 0.0D+0
378: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
379: SYMB_ZERO = .TRUE.
380: ELSE
381: SYMB_ZERO = .FALSE.
382: Y( IY ) = BETA * ABS( Y( IY ) )
383: END IF
384: JX = KX
385: IF ( ALPHA .NE. ZERO ) THEN
386: DO J = 1, I
387: TEMP = ABS( A( I, J ) )
388: SYMB_ZERO = SYMB_ZERO .AND.
389: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
390:
391: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
392: JX = JX + INCX
393: END DO
394: DO J = I+1, N
395: TEMP = ABS( A( J, I ) )
396: SYMB_ZERO = SYMB_ZERO .AND.
397: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
398:
399: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
400: JX = JX + INCX
401: END DO
402: END IF
403:
404: IF ( .NOT.SYMB_ZERO )
405: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
406:
407: IY = IY + INCY
408: END DO
409: END IF
410:
411: END IF
412: *
413: RETURN
414: *
415: * End of DLA_SYAMV
416: *
417: END