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 of 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: *> \date December 2016
182: *
183: *> \ingroup complex16GBcomputational
184: *
185: * =====================================================================
186: SUBROUTINE ZLA_GBAMV( TRANS, M, N, KL, KU, ALPHA, AB, LDAB, X,
187: $ INCX, BETA, Y, INCY )
188: *
189: * -- LAPACK computational routine (version 3.7.0) --
190: * -- LAPACK is a software package provided by Univ. of Tennessee, --
191: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
192: * December 2016
193: *
194: * .. Scalar Arguments ..
195: DOUBLE PRECISION ALPHA, BETA
196: INTEGER INCX, INCY, LDAB, M, N, KL, KU, TRANS
197: * ..
198: * .. Array Arguments ..
199: COMPLEX*16 AB( LDAB, * ), X( * )
200: DOUBLE PRECISION Y( * )
201: * ..
202: *
203: * =====================================================================
204: *
205: * .. Parameters ..
206: COMPLEX*16 ONE, ZERO
207: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
208: * ..
209: * .. Local Scalars ..
210: LOGICAL SYMB_ZERO
211: DOUBLE PRECISION TEMP, SAFE1
212: INTEGER I, INFO, IY, J, JX, KX, KY, LENX, LENY, KD, KE
213: COMPLEX*16 CDUM
214: * ..
215: * .. External Subroutines ..
216: EXTERNAL XERBLA, DLAMCH
217: DOUBLE PRECISION DLAMCH
218: * ..
219: * .. External Functions ..
220: EXTERNAL ILATRANS
221: INTEGER ILATRANS
222: * ..
223: * .. Intrinsic Functions ..
224: INTRINSIC MAX, ABS, REAL, DIMAG, SIGN
225: * ..
226: * .. Statement Functions
227: DOUBLE PRECISION CABS1
228: * ..
229: * .. Statement Function Definitions ..
230: CABS1( CDUM ) = ABS( DBLE( CDUM ) ) + ABS( DIMAG( CDUM ) )
231: * ..
232: * .. Executable Statements ..
233: *
234: * Test the input parameters.
235: *
236: INFO = 0
237: IF ( .NOT.( ( TRANS.EQ.ILATRANS( 'N' ) )
238: $ .OR. ( TRANS.EQ.ILATRANS( 'T' ) )
239: $ .OR. ( TRANS.EQ.ILATRANS( 'C' ) ) ) ) THEN
240: INFO = 1
241: ELSE IF( M.LT.0 )THEN
242: INFO = 2
243: ELSE IF( N.LT.0 )THEN
244: INFO = 3
245: ELSE IF( KL.LT.0 .OR. KL.GT.M-1 ) THEN
246: INFO = 4
247: ELSE IF( KU.LT.0 .OR. KU.GT.N-1 ) THEN
248: INFO = 5
249: ELSE IF( LDAB.LT.KL+KU+1 )THEN
250: INFO = 6
251: ELSE IF( INCX.EQ.0 )THEN
252: INFO = 8
253: ELSE IF( INCY.EQ.0 )THEN
254: INFO = 11
255: END IF
256: IF( INFO.NE.0 )THEN
257: CALL XERBLA( 'ZLA_GBAMV ', INFO )
258: RETURN
259: END IF
260: *
261: * Quick return if possible.
262: *
263: IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR.
264: $ ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
265: $ RETURN
266: *
267: * Set LENX and LENY, the lengths of the vectors x and y, and set
268: * up the start points in X and Y.
269: *
270: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
271: LENX = N
272: LENY = M
273: ELSE
274: LENX = M
275: LENY = N
276: END IF
277: IF( INCX.GT.0 )THEN
278: KX = 1
279: ELSE
280: KX = 1 - ( LENX - 1 )*INCX
281: END IF
282: IF( INCY.GT.0 )THEN
283: KY = 1
284: ELSE
285: KY = 1 - ( LENY - 1 )*INCY
286: END IF
287: *
288: * Set SAFE1 essentially to be the underflow threshold times the
289: * number of additions in each row.
290: *
291: SAFE1 = DLAMCH( 'Safe minimum' )
292: SAFE1 = (N+1)*SAFE1
293: *
294: * Form y := alpha*abs(A)*abs(x) + beta*abs(y).
295: *
296: * The O(M*N) SYMB_ZERO tests could be replaced by O(N) queries to
297: * the inexact flag. Still doesn't help change the iteration order
298: * to per-column.
299: *
300: KD = KU + 1
301: KE = KL + 1
302: IY = KY
303: IF ( INCX.EQ.1 ) THEN
304: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
305: DO I = 1, LENY
306: IF ( BETA .EQ. 0.0D+0 ) THEN
307: SYMB_ZERO = .TRUE.
308: Y( IY ) = 0.0D+0
309: ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
310: SYMB_ZERO = .TRUE.
311: ELSE
312: SYMB_ZERO = .FALSE.
313: Y( IY ) = BETA * ABS( Y( IY ) )
314: END IF
315: IF ( ALPHA .NE. 0.0D+0 ) THEN
316: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
317: TEMP = CABS1( AB( KD+I-J, J ) )
318: SYMB_ZERO = SYMB_ZERO .AND.
319: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
320:
321: Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP
322: END DO
323: END IF
324:
325: IF ( .NOT.SYMB_ZERO)
326: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
327:
328: IY = IY + INCY
329: END DO
330: ELSE
331: DO I = 1, LENY
332: IF ( BETA .EQ. 0.0D+0 ) THEN
333: SYMB_ZERO = .TRUE.
334: Y( IY ) = 0.0D+0
335: ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
336: SYMB_ZERO = .TRUE.
337: ELSE
338: SYMB_ZERO = .FALSE.
339: Y( IY ) = BETA * ABS( Y( IY ) )
340: END IF
341: IF ( ALPHA .NE. 0.0D+0 ) THEN
342: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
343: TEMP = CABS1( AB( KE-I+J, I ) )
344: SYMB_ZERO = SYMB_ZERO .AND.
345: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
346:
347: Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP
348: END DO
349: END IF
350:
351: IF ( .NOT.SYMB_ZERO)
352: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
353:
354: IY = IY + INCY
355: END DO
356: END IF
357: ELSE
358: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
359: DO I = 1, LENY
360: IF ( BETA .EQ. 0.0D+0 ) THEN
361: SYMB_ZERO = .TRUE.
362: Y( IY ) = 0.0D+0
363: ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
364: SYMB_ZERO = .TRUE.
365: ELSE
366: SYMB_ZERO = .FALSE.
367: Y( IY ) = BETA * ABS( Y( IY ) )
368: END IF
369: IF ( ALPHA .NE. 0.0D+0 ) THEN
370: JX = KX
371: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
372: TEMP = CABS1( AB( KD+I-J, J ) )
373: SYMB_ZERO = SYMB_ZERO .AND.
374: $ ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
375:
376: Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP
377: JX = JX + INCX
378: END DO
379: END IF
380:
381: IF ( .NOT.SYMB_ZERO )
382: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
383:
384: IY = IY + INCY
385: END DO
386: ELSE
387: DO I = 1, LENY
388: IF ( BETA .EQ. 0.0D+0 ) THEN
389: SYMB_ZERO = .TRUE.
390: Y( IY ) = 0.0D+0
391: ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
392: SYMB_ZERO = .TRUE.
393: ELSE
394: SYMB_ZERO = .FALSE.
395: Y( IY ) = BETA * ABS( Y( IY ) )
396: END IF
397: IF ( ALPHA .NE. 0.0D+0 ) THEN
398: JX = KX
399: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
400: TEMP = CABS1( AB( KE-I+J, I ) )
401: SYMB_ZERO = SYMB_ZERO .AND.
402: $ ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
403:
404: Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP
405: JX = JX + INCX
406: END DO
407: END IF
408:
409: IF ( .NOT.SYMB_ZERO )
410: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
411:
412: IY = IY + INCY
413: END DO
414: END IF
415:
416: END IF
417: *
418: RETURN
419: *
420: * End of ZLA_GBAMV
421: *
422: 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