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