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