LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zla_syamv.f
Go to the documentation of this file.
1*> \brief \b ZLA_SYAMV computes a matrix-vector product using a symmetric indefinite 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 ZLA_SYAMV + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zla_syamv.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zla_syamv.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zla_syamv.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE ZLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y,
22* INCY )
23*
24* .. Scalar Arguments ..
25* DOUBLE PRECISION ALPHA, BETA
26* INTEGER INCX, INCY, LDA, N
27* INTEGER UPLO
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_SYAMV performs the matrix-vector operation
41*>
42*> y := alpha*abs(A)*abs(x) + beta*abs(y),
43*>
44*> where alpha and beta are scalars, x and y are vectors and A is an
45*> n by n symmetric 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] UPLO
61*> \verbatim
62*> UPLO is INTEGER
63*> On entry, UPLO specifies whether the upper or lower
64*> triangular part of the array A is to be referenced as
65*> follows:
66*>
67*> UPLO = BLAS_UPPER Only the upper triangular part of A
68*> is to be referenced.
69*>
70*> UPLO = BLAS_LOWER Only the lower triangular part of A
71*> is to be referenced.
72*>
73*> Unchanged on exit.
74*> \endverbatim
75*>
76*> \param[in] N
77*> \verbatim
78*> N is INTEGER
79*> On entry, N specifies the number of columns of the matrix A.
80*> N must be at least zero.
81*> Unchanged on exit.
82*> \endverbatim
83*>
84*> \param[in] ALPHA
85*> \verbatim
86*> ALPHA is DOUBLE PRECISION .
87*> On entry, ALPHA specifies the scalar alpha.
88*> Unchanged on exit.
89*> \endverbatim
90*>
91*> \param[in] A
92*> \verbatim
93*> A is COMPLEX*16 array, dimension ( LDA, n ).
94*> Before entry, the leading m by n part of the array A must
95*> contain the matrix of coefficients.
96*> Unchanged on exit.
97*> \endverbatim
98*>
99*> \param[in] LDA
100*> \verbatim
101*> LDA is INTEGER
102*> On entry, LDA specifies the first dimension of A as declared
103*> in the calling (sub) program. LDA must be at least
104*> max( 1, n ).
105*> Unchanged on exit.
106*> \endverbatim
107*>
108*> \param[in] X
109*> \verbatim
110*> X is COMPLEX*16 array, dimension at least
111*> ( 1 + ( n - 1 )*abs( INCX ) )
112*> Before entry, the incremented array X must contain the
113*> vector x.
114*> Unchanged on exit.
115*> \endverbatim
116*>
117*> \param[in] INCX
118*> \verbatim
119*> INCX is INTEGER
120*> On entry, INCX specifies the increment for the elements of
121*> X. INCX must not be zero.
122*> Unchanged on exit.
123*> \endverbatim
124*>
125*> \param[in] BETA
126*> \verbatim
127*> BETA is DOUBLE PRECISION .
128*> On entry, BETA specifies the scalar beta. When BETA is
129*> supplied as zero then Y need not be set on input.
130*> Unchanged on exit.
131*> \endverbatim
132*>
133*> \param[in,out] Y
134*> \verbatim
135*> Y is DOUBLE PRECISION array, dimension
136*> ( 1 + ( n - 1 )*abs( INCY ) )
137*> Before entry with BETA non-zero, the incremented array Y
138*> must contain the vector y. On exit, Y is overwritten by the
139*> updated vector y.
140*> \endverbatim
141*>
142*> \param[in] INCY
143*> \verbatim
144*> INCY is INTEGER
145*> On entry, INCY specifies the increment for the elements of
146*> Y. INCY must not be zero.
147*> Unchanged on exit.
148*> \endverbatim
149*
150* Authors:
151* ========
152*
153*> \author Univ. of Tennessee
154*> \author Univ. of California Berkeley
155*> \author Univ. of Colorado Denver
156*> \author NAG Ltd.
157*
158*> \ingroup la_heamv
159*
160*> \par Further Details:
161* =====================
162*>
163*> \verbatim
164*>
165*> Level 2 Blas routine.
166*>
167*> -- Written on 22-October-1986.
168*> Jack Dongarra, Argonne National Lab.
169*> Jeremy Du Croz, Nag Central Office.
170*> Sven Hammarling, Nag Central Office.
171*> Richard Hanson, Sandia National Labs.
172*> -- Modified for the absolute-value product, April 2006
173*> Jason Riedy, UC Berkeley
174*> \endverbatim
175*>
176* =====================================================================
177 SUBROUTINE zla_syamv( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y,
178 $ INCY )
179*
180* -- LAPACK computational routine --
181* -- LAPACK is a software package provided by Univ. of Tennessee, --
182* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
183*
184* .. Scalar Arguments ..
185 DOUBLE PRECISION ALPHA, BETA
186 INTEGER INCX, INCY, LDA, N
187 INTEGER UPLO
188* ..
189* .. Array Arguments ..
190 COMPLEX*16 A( LDA, * ), X( * )
191 DOUBLE PRECISION Y( * )
192* ..
193*
194* =====================================================================
195*
196* .. Parameters ..
197 DOUBLE PRECISION ONE, ZERO
198 parameter( one = 1.0d+0, zero = 0.0d+0 )
199* ..
200* .. Local Scalars ..
201 LOGICAL SYMB_ZERO
202 DOUBLE PRECISION TEMP, SAFE1
203 INTEGER I, INFO, IY, J, JX, KX, KY
204 COMPLEX*16 ZDUM
205* ..
206* .. External Subroutines ..
207 EXTERNAL xerbla, dlamch
208 DOUBLE PRECISION DLAMCH
209* ..
210* .. External Functions ..
211 EXTERNAL ilauplo
212 INTEGER ILAUPLO
213* ..
214* .. Intrinsic Functions ..
215 INTRINSIC max, abs, sign, real, dimag
216* ..
217* .. Statement Functions ..
218 DOUBLE PRECISION CABS1
219* ..
220* .. Statement Function Definitions ..
221 cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
222* ..
223* .. Executable Statements ..
224*
225* Test the input parameters.
226*
227 info = 0
228 IF ( uplo.NE.ilauplo( 'U' ) .AND.
229 $ uplo.NE.ilauplo( 'L' ) )THEN
230 info = 1
231 ELSE IF( n.LT.0 )THEN
232 info = 2
233 ELSE IF( lda.LT.max( 1, n ) )THEN
234 info = 5
235 ELSE IF( incx.EQ.0 )THEN
236 info = 7
237 ELSE IF( incy.EQ.0 )THEN
238 info = 10
239 END IF
240 IF( info.NE.0 )THEN
241 CALL xerbla( 'ZLA_SYAMV', info )
242 RETURN
243 END IF
244*
245* Quick return if possible.
246*
247 IF( ( n.EQ.0 ).OR.( ( alpha.EQ.zero ).AND.( beta.EQ.one ) ) )
248 $ RETURN
249*
250* Set up the start points in X and Y.
251*
252 IF( incx.GT.0 )THEN
253 kx = 1
254 ELSE
255 kx = 1 - ( n - 1 )*incx
256 END IF
257 IF( incy.GT.0 )THEN
258 ky = 1
259 ELSE
260 ky = 1 - ( n - 1 )*incy
261 END IF
262*
263* Set SAFE1 essentially to be the underflow threshold times the
264* number of additions in each row.
265*
266 safe1 = dlamch( 'Safe minimum' )
267 safe1 = (n+1)*safe1
268*
269* Form y := alpha*abs(A)*abs(x) + beta*abs(y).
270*
271* The O(N^2) SYMB_ZERO tests could be replaced by O(N) queries to
272* the inexact flag. Still doesn't help change the iteration order
273* to per-column.
274*
275 iy = ky
276 IF ( incx.EQ.1 ) THEN
277 IF ( uplo .EQ. ilauplo( 'U' ) ) THEN
278 DO i = 1, n
279 IF ( beta .EQ. zero ) THEN
280 symb_zero = .true.
281 y( iy ) = 0.0d+0
282 ELSE IF ( y( iy ) .EQ. zero ) THEN
283 symb_zero = .true.
284 ELSE
285 symb_zero = .false.
286 y( iy ) = beta * abs( y( iy ) )
287 END IF
288 IF ( alpha .NE. zero ) THEN
289 DO j = 1, i
290 temp = cabs1( a( j, i ) )
291 symb_zero = symb_zero .AND.
292 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
293
294 y( iy ) = y( iy ) + alpha*cabs1( x( j ) )*temp
295 END DO
296 DO j = i+1, n
297 temp = cabs1( a( i, j ) )
298 symb_zero = symb_zero .AND.
299 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
300
301 y( iy ) = y( iy ) + alpha*cabs1( x( j ) )*temp
302 END DO
303 END IF
304
305 IF ( .NOT.symb_zero )
306 $ y( iy ) = y( iy ) + sign( safe1, y( iy ) )
307
308 iy = iy + incy
309 END DO
310 ELSE
311 DO i = 1, n
312 IF ( beta .EQ. zero ) THEN
313 symb_zero = .true.
314 y( iy ) = 0.0d+0
315 ELSE IF ( y( iy ) .EQ. zero ) THEN
316 symb_zero = .true.
317 ELSE
318 symb_zero = .false.
319 y( iy ) = beta * abs( y( iy ) )
320 END IF
321 IF ( alpha .NE. zero ) THEN
322 DO j = 1, i
323 temp = cabs1( a( i, j ) )
324 symb_zero = symb_zero .AND.
325 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
326
327 y( iy ) = y( iy ) + alpha*cabs1( x( j ) )*temp
328 END DO
329 DO j = i+1, n
330 temp = cabs1( a( j, i ) )
331 symb_zero = symb_zero .AND.
332 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
333
334 y( iy ) = y( iy ) + alpha*cabs1( x( j ) )*temp
335 END DO
336 END IF
337
338 IF ( .NOT.symb_zero )
339 $ y( iy ) = y( iy ) + sign( safe1, y( iy ) )
340
341 iy = iy + incy
342 END DO
343 END IF
344 ELSE
345 IF ( uplo .EQ. ilauplo( 'U' ) ) THEN
346 DO i = 1, n
347 IF ( beta .EQ. zero ) THEN
348 symb_zero = .true.
349 y( iy ) = 0.0d+0
350 ELSE IF ( y( iy ) .EQ. zero ) THEN
351 symb_zero = .true.
352 ELSE
353 symb_zero = .false.
354 y( iy ) = beta * abs( y( iy ) )
355 END IF
356 jx = kx
357 IF ( alpha .NE. zero ) THEN
358 DO j = 1, i
359 temp = cabs1( a( j, i ) )
360 symb_zero = symb_zero .AND.
361 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
362
363 y( iy ) = y( iy ) + alpha*cabs1( x( jx ) )*temp
364 jx = jx + incx
365 END DO
366 DO j = i+1, n
367 temp = cabs1( a( i, j ) )
368 symb_zero = symb_zero .AND.
369 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
370
371 y( iy ) = y( iy ) + alpha*cabs1( x( jx ) )*temp
372 jx = jx + incx
373 END DO
374 END IF
375
376 IF ( .NOT.symb_zero )
377 $ y( iy ) = y( iy ) + sign( safe1, y( iy ) )
378
379 iy = iy + incy
380 END DO
381 ELSE
382 DO i = 1, n
383 IF ( beta .EQ. zero ) THEN
384 symb_zero = .true.
385 y( iy ) = 0.0d+0
386 ELSE IF ( y( iy ) .EQ. zero ) THEN
387 symb_zero = .true.
388 ELSE
389 symb_zero = .false.
390 y( iy ) = beta * abs( y( iy ) )
391 END IF
392 jx = kx
393 IF ( alpha .NE. zero ) THEN
394 DO j = 1, i
395 temp = cabs1( a( i, j ) )
396 symb_zero = symb_zero .AND.
397 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
398
399 y( iy ) = y( iy ) + alpha*cabs1( x( jx ) )*temp
400 jx = jx + incx
401 END DO
402 DO j = i+1, n
403 temp = cabs1( a( j, i ) )
404 symb_zero = symb_zero .AND.
405 $ ( x( j ) .EQ. zero .OR. temp .EQ. zero )
406
407 y( iy ) = y( iy ) + alpha*cabs1( x( jx ) )*temp
408 jx = jx + incx
409 END DO
410 END IF
411
412 IF ( .NOT.symb_zero )
413 $ y( iy ) = y( iy ) + sign( safe1, y( iy ) )
414
415 iy = iy + incy
416 END DO
417 END IF
418
419 END IF
420*
421 RETURN
422*
423* End of ZLA_SYAMV
424*
425 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
integer function ilauplo(uplo)
ILAUPLO
Definition ilauplo.f:58
subroutine zla_syamv(uplo, n, alpha, a, lda, x, incx, beta, y, incy)
ZLA_SYAMV computes a matrix-vector product using a symmetric indefinite matrix to calculate error bou...
Definition zla_syamv.f:179
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69