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*> \date June 2017
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.1) --
190*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
191*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
192*     June 2017
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
423