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