1: *> \brief \b DLA_GEAMV computes a matrix-vector product using a general 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 DLA_GEAMV + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dla_geamv.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dla_geamv.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dla_geamv.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE DLA_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, TRANS
27: * ..
28: * .. Array Arguments ..
29: * DOUBLE PRECISION A( LDA, * ), X( * ), Y( * )
30: * ..
31: *
32: *
33: *> \par Purpose:
34: * =============
35: *>
36: *> \verbatim
37: *>
38: *> DLA_GEAMV performs one of the matrix-vector operations
39: *>
40: *> y := alpha*abs(A)*abs(x) + beta*abs(y),
41: *> or y := alpha*abs(A)**T*abs(x) + beta*abs(y),
42: *>
43: *> where alpha and beta are scalars, x and y are vectors and A is an
44: *> m by n 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] TRANS
60: *> \verbatim
61: *> TRANS is INTEGER
62: *> On entry, TRANS specifies the operation to be performed as
63: *> follows:
64: *>
65: *> BLAS_NO_TRANS y := alpha*abs(A)*abs(x) + beta*abs(y)
66: *> BLAS_TRANS y := alpha*abs(A**T)*abs(x) + beta*abs(y)
67: *> BLAS_CONJ_TRANS y := alpha*abs(A**T)*abs(x) + beta*abs(y)
68: *>
69: *> Unchanged on exit.
70: *> \endverbatim
71: *>
72: *> \param[in] M
73: *> \verbatim
74: *> M is INTEGER
75: *> On entry, M specifies the number of rows of the matrix A.
76: *> M must be at least zero.
77: *> Unchanged on exit.
78: *> \endverbatim
79: *>
80: *> \param[in] N
81: *> \verbatim
82: *> N is INTEGER
83: *> On entry, N specifies the number of columns of the matrix A.
84: *> N must be at least zero.
85: *> Unchanged on exit.
86: *> \endverbatim
87: *>
88: *> \param[in] ALPHA
89: *> \verbatim
90: *> ALPHA is DOUBLE PRECISION
91: *> On entry, ALPHA specifies the scalar alpha.
92: *> Unchanged on exit.
93: *> \endverbatim
94: *>
95: *> \param[in] A
96: *> \verbatim
97: *> A is DOUBLE PRECISION array, dimension ( LDA, n )
98: *> Before entry, the leading m by n part of the array A must
99: *> contain the matrix of coefficients.
100: *> Unchanged on exit.
101: *> \endverbatim
102: *>
103: *> \param[in] LDA
104: *> \verbatim
105: *> LDA is INTEGER
106: *> On entry, LDA specifies the first dimension of A as declared
107: *> in the calling (sub) program. LDA must be at least
108: *> max( 1, m ).
109: *> Unchanged on exit.
110: *> \endverbatim
111: *>
112: *> \param[in] X
113: *> \verbatim
114: *> X is DOUBLE PRECISION array, dimension
115: *> ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n'
116: *> and at least
117: *> ( 1 + ( m - 1 )*abs( INCX ) ) otherwise.
118: *> Before entry, the incremented array X must contain the
119: *> vector x.
120: *> Unchanged on exit.
121: *> \endverbatim
122: *>
123: *> \param[in] INCX
124: *> \verbatim
125: *> INCX is INTEGER
126: *> On entry, INCX specifies the increment for the elements of
127: *> X. INCX must not be zero.
128: *> Unchanged on exit.
129: *> \endverbatim
130: *>
131: *> \param[in] BETA
132: *> \verbatim
133: *> BETA is DOUBLE PRECISION
134: *> On entry, BETA specifies the scalar beta. When BETA is
135: *> supplied as zero then Y need not be set on input.
136: *> Unchanged on exit.
137: *> \endverbatim
138: *>
139: *> \param[in,out] Y
140: *> \verbatim
141: *> Y is DOUBLE PRECISION array,
142: *> dimension at least
143: *> ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n'
144: *> and at least
145: *> ( 1 + ( n - 1 )*abs( INCY ) ) otherwise.
146: *> Before entry with BETA non-zero, the incremented array Y
147: *> must contain the vector y. On exit, Y is overwritten by the
148: *> updated vector y.
149: *> \endverbatim
150: *>
151: *> \param[in] INCY
152: *> \verbatim
153: *> INCY is INTEGER
154: *> On entry, INCY specifies the increment for the elements of
155: *> Y. INCY must not be zero.
156: *> Unchanged on exit.
157: *>
158: *> Level 2 Blas routine.
159: *> \endverbatim
160: *
161: * Authors:
162: * ========
163: *
164: *> \author Univ. of Tennessee
165: *> \author Univ. of California Berkeley
166: *> \author Univ. of Colorado Denver
167: *> \author NAG Ltd.
168: *
169: *> \ingroup doubleGEcomputational
170: *
171: * =====================================================================
172: SUBROUTINE DLA_GEAMV ( TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA,
173: $ Y, INCY )
174: *
175: * -- LAPACK computational routine --
176: * -- LAPACK is a software package provided by Univ. of Tennessee, --
177: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
178: *
179: * .. Scalar Arguments ..
180: DOUBLE PRECISION ALPHA, BETA
181: INTEGER INCX, INCY, LDA, M, N, TRANS
182: * ..
183: * .. Array Arguments ..
184: DOUBLE PRECISION A( LDA, * ), X( * ), Y( * )
185: * ..
186: *
187: * =====================================================================
188: *
189: * .. Parameters ..
190: DOUBLE PRECISION ONE, ZERO
191: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
192: * ..
193: * .. Local Scalars ..
194: LOGICAL SYMB_ZERO
195: DOUBLE PRECISION TEMP, SAFE1
196: INTEGER I, INFO, IY, J, JX, KX, KY, LENX, LENY
197: * ..
198: * .. External Subroutines ..
199: EXTERNAL XERBLA, DLAMCH
200: DOUBLE PRECISION DLAMCH
201: * ..
202: * .. External Functions ..
203: EXTERNAL ILATRANS
204: INTEGER ILATRANS
205: * ..
206: * .. Intrinsic Functions ..
207: INTRINSIC MAX, ABS, SIGN
208: * ..
209: * .. Executable Statements ..
210: *
211: * Test the input parameters.
212: *
213: INFO = 0
214: IF ( .NOT.( ( TRANS.EQ.ILATRANS( 'N' ) )
215: $ .OR. ( TRANS.EQ.ILATRANS( 'T' ) )
216: $ .OR. ( TRANS.EQ.ILATRANS( 'C' )) ) ) THEN
217: INFO = 1
218: ELSE IF( M.LT.0 )THEN
219: INFO = 2
220: ELSE IF( N.LT.0 )THEN
221: INFO = 3
222: ELSE IF( LDA.LT.MAX( 1, M ) )THEN
223: INFO = 6
224: ELSE IF( INCX.EQ.0 )THEN
225: INFO = 8
226: ELSE IF( INCY.EQ.0 )THEN
227: INFO = 11
228: END IF
229: IF( INFO.NE.0 )THEN
230: CALL XERBLA( 'DLA_GEAMV ', INFO )
231: RETURN
232: END IF
233: *
234: * Quick return if possible.
235: *
236: IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR.
237: $ ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
238: $ RETURN
239: *
240: * Set LENX and LENY, the lengths of the vectors x and y, and set
241: * up the start points in X and Y.
242: *
243: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
244: LENX = N
245: LENY = M
246: ELSE
247: LENX = M
248: LENY = N
249: END IF
250: IF( INCX.GT.0 )THEN
251: KX = 1
252: ELSE
253: KX = 1 - ( LENX - 1 )*INCX
254: END IF
255: IF( INCY.GT.0 )THEN
256: KY = 1
257: ELSE
258: KY = 1 - ( LENY - 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(M*N) 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( TRANS.EQ.ILATRANS( 'N' ) )THEN
276: DO I = 1, LENY
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, LENX
288: TEMP = ABS( A( I, J ) )
289: SYMB_ZERO = SYMB_ZERO .AND.
290: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
291:
292: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
293: END DO
294: END IF
295:
296: IF ( .NOT.SYMB_ZERO )
297: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
298:
299: IY = IY + INCY
300: END DO
301: ELSE
302: DO I = 1, LENY
303: IF ( BETA .EQ. ZERO ) THEN
304: SYMB_ZERO = .TRUE.
305: Y( IY ) = 0.0D+0
306: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
307: SYMB_ZERO = .TRUE.
308: ELSE
309: SYMB_ZERO = .FALSE.
310: Y( IY ) = BETA * ABS( Y( IY ) )
311: END IF
312: IF ( ALPHA .NE. ZERO ) THEN
313: DO J = 1, LENX
314: TEMP = ABS( A( J, I ) )
315: SYMB_ZERO = SYMB_ZERO .AND.
316: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
317:
318: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
319: END DO
320: END IF
321:
322: IF ( .NOT.SYMB_ZERO )
323: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
324:
325: IY = IY + INCY
326: END DO
327: END IF
328: ELSE
329: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
330: DO I = 1, LENY
331: IF ( BETA .EQ. ZERO ) THEN
332: SYMB_ZERO = .TRUE.
333: Y( IY ) = 0.0D+0
334: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
335: SYMB_ZERO = .TRUE.
336: ELSE
337: SYMB_ZERO = .FALSE.
338: Y( IY ) = BETA * ABS( Y( IY ) )
339: END IF
340: IF ( ALPHA .NE. ZERO ) THEN
341: JX = KX
342: DO J = 1, LENX
343: TEMP = ABS( A( I, J ) )
344: SYMB_ZERO = SYMB_ZERO .AND.
345: $ ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
346:
347: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
348: JX = JX + INCX
349: END DO
350: END IF
351:
352: IF (.NOT.SYMB_ZERO)
353: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
354:
355: IY = IY + INCY
356: END DO
357: ELSE
358: DO I = 1, LENY
359: IF ( BETA .EQ. ZERO ) THEN
360: SYMB_ZERO = .TRUE.
361: Y( IY ) = 0.0D+0
362: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
363: SYMB_ZERO = .TRUE.
364: ELSE
365: SYMB_ZERO = .FALSE.
366: Y( IY ) = BETA * ABS( Y( IY ) )
367: END IF
368: IF ( ALPHA .NE. ZERO ) THEN
369: JX = KX
370: DO J = 1, LENX
371: TEMP = ABS( A( J, I ) )
372: SYMB_ZERO = SYMB_ZERO .AND.
373: $ ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
374:
375: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
376: JX = JX + INCX
377: END DO
378: END IF
379:
380: IF (.NOT.SYMB_ZERO)
381: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
382:
383: IY = IY + INCY
384: END DO
385: END IF
386:
387: END IF
388: *
389: RETURN
390: *
391: * End of DLA_GEAMV
392: *
393: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dla_geamv.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack