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