LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ zlavhp()

subroutine zlavhp ( character  UPLO,
character  TRANS,
character  DIAG,
integer  N,
integer  NRHS,
complex*16, dimension( * )  A,
integer, dimension( * )  IPIV,
complex*16, dimension( ldb, * )  B,
integer  LDB,
integer  INFO 
)

ZLAVHP

Purpose:
    ZLAVHP  performs one of the matrix-vector operations
       x := A*x  or  x := A^H*x,
    where x is an N element vector and  A is one of the factors
    from the symmetric factorization computed by ZHPTRF.
    ZHPTRF produces a factorization of the form
         U * D * U^H     or     L * D * L^H,
    where U (or L) is a product of permutation and unit upper (lower)
    triangular matrices, U^H (or L^H) is the conjugate transpose of
    U (or L), and D is Hermitian and block diagonal with 1 x 1 and
    2 x 2 diagonal blocks.  The multipliers for the transformations
    and the upper or lower triangular parts of the diagonal blocks
    are stored columnwise in packed format in the linear array A.

    If TRANS = 'N' or 'n', ZLAVHP multiplies either by U or U * D
    (or L or L * D).
    If TRANS = 'C' or 'c', ZLAVHP multiplies either by U^H or D * U^H
    (or L^H or D * L^H ).
  UPLO   - CHARACTER*1
           On entry, UPLO specifies whether the triangular matrix
           stored in A is upper or lower triangular.
              UPLO = 'U' or 'u'   The matrix is upper triangular.
              UPLO = 'L' or 'l'   The matrix is lower triangular.
           Unchanged on exit.

  TRANS  - CHARACTER*1
           On entry, TRANS specifies the operation to be performed as
           follows:
              TRANS = 'N' or 'n'   x := A*x.
              TRANS = 'C' or 'c'   x := A^H*x.
           Unchanged on exit.

  DIAG   - CHARACTER*1
           On entry, DIAG specifies whether the diagonal blocks are
           assumed to be unit matrices, as follows:
              DIAG = 'U' or 'u'   Diagonal blocks are unit matrices.
              DIAG = 'N' or 'n'   Diagonal blocks are non-unit.
           Unchanged on exit.

  N      - INTEGER
           On entry, N specifies the order of the matrix A.
           N must be at least zero.
           Unchanged on exit.

  NRHS   - INTEGER
           On entry, NRHS specifies the number of right hand sides,
           i.e., the number of vectors x to be multiplied by A.
           NRHS must be at least zero.
           Unchanged on exit.

  A      - COMPLEX*16 array, dimension( N*(N+1)/2 )
           On entry, A contains a block diagonal matrix and the
           multipliers of the transformations used to obtain it,
           stored as a packed triangular matrix.
           Unchanged on exit.

  IPIV   - INTEGER array, dimension( N )
           On entry, IPIV contains the vector of pivot indices as
           determined by ZSPTRF or ZHPTRF.
           If IPIV( K ) = K, no interchange was done.
           If IPIV( K ) <> K but IPIV( K ) > 0, then row K was inter-
           changed with row IPIV( K ) and a 1 x 1 pivot block was used.
           If IPIV( K ) < 0 and UPLO = 'U', then row K-1 was exchanged
           with row | IPIV( K ) | and a 2 x 2 pivot block was used.
           If IPIV( K ) < 0 and UPLO = 'L', then row K+1 was exchanged
           with row | IPIV( K ) | and a 2 x 2 pivot block was used.

  B      - COMPLEX*16 array, dimension( LDB, NRHS )
           On entry, B contains NRHS vectors of length N.
           On exit, B is overwritten with the product A * B.

  LDB    - INTEGER
           On entry, LDB contains the leading dimension of B as
           declared in the calling program.  LDB must be at least
           max( 1, N ).
           Unchanged on exit.

  INFO   - INTEGER
           INFO is the error flag.
           On exit, a value of 0 indicates a successful exit.
           A negative value, say -K, indicates that the K-th argument
           has an illegal value.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
December 2016

Definition at line 133 of file zlavhp.f.

133 *
134 * -- LAPACK test routine (version 3.7.0) --
135 * -- LAPACK is a software package provided by Univ. of Tennessee, --
136 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
137 * December 2016
138 *
139 * .. Scalar Arguments ..
140  CHARACTER diag, trans, uplo
141  INTEGER info, ldb, n, nrhs
142 * ..
143 * .. Array Arguments ..
144  INTEGER ipiv( * )
145  COMPLEX*16 a( * ), b( ldb, * )
146 * ..
147 *
148 * =====================================================================
149 *
150 * .. Parameters ..
151  COMPLEX*16 one
152  parameter( one = ( 1.0d+0, 0.0d+0 ) )
153 * ..
154 * .. Local Scalars ..
155  LOGICAL nounit
156  INTEGER j, k, kc, kcnext, kp
157  COMPLEX*16 d11, d12, d21, d22, t1, t2
158 * ..
159 * .. External Functions ..
160  LOGICAL lsame
161  EXTERNAL lsame
162 * ..
163 * .. External Subroutines ..
164  EXTERNAL xerbla, zgemv, zgeru, zlacgv, zscal, zswap
165 * ..
166 * .. Intrinsic Functions ..
167  INTRINSIC abs, dconjg, max
168 * ..
169 * .. Executable Statements ..
170 *
171 * Test the input parameters.
172 *
173  info = 0
174  IF( .NOT.lsame( uplo, 'U' ) .AND. .NOT.lsame( uplo, 'L' ) ) THEN
175  info = -1
176  ELSE IF( .NOT.lsame( trans, 'N' ) .AND. .NOT.lsame( trans, 'C' ) )
177  $ THEN
178  info = -2
179  ELSE IF( .NOT.lsame( diag, 'U' ) .AND. .NOT.lsame( diag, 'N' ) )
180  $ THEN
181  info = -3
182  ELSE IF( n.LT.0 ) THEN
183  info = -4
184  ELSE IF( ldb.LT.max( 1, n ) ) THEN
185  info = -8
186  END IF
187  IF( info.NE.0 ) THEN
188  CALL xerbla( 'ZLAVHP ', -info )
189  RETURN
190  END IF
191 *
192 * Quick return if possible.
193 *
194  IF( n.EQ.0 )
195  $ RETURN
196 *
197  nounit = lsame( diag, 'N' )
198 *------------------------------------------
199 *
200 * Compute B := A * B (No transpose)
201 *
202 *------------------------------------------
203  IF( lsame( trans, 'N' ) ) THEN
204 *
205 * Compute B := U*B
206 * where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1))
207 *
208  IF( lsame( uplo, 'U' ) ) THEN
209 *
210 * Loop forward applying the transformations.
211 *
212  k = 1
213  kc = 1
214  10 CONTINUE
215  IF( k.GT.n )
216  $ GO TO 30
217 *
218 * 1 x 1 pivot block
219 *
220  IF( ipiv( k ).GT.0 ) THEN
221 *
222 * Multiply by the diagonal element if forming U * D.
223 *
224  IF( nounit )
225  $ CALL zscal( nrhs, a( kc+k-1 ), b( k, 1 ), ldb )
226 *
227 * Multiply by P(K) * inv(U(K)) if K > 1.
228 *
229  IF( k.GT.1 ) THEN
230 *
231 * Apply the transformation.
232 *
233  CALL zgeru( k-1, nrhs, one, a( kc ), 1, b( k, 1 ),
234  $ ldb, b( 1, 1 ), ldb )
235 *
236 * Interchange if P(K) != I.
237 *
238  kp = ipiv( k )
239  IF( kp.NE.k )
240  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
241  END IF
242  kc = kc + k
243  k = k + 1
244  ELSE
245 *
246 * 2 x 2 pivot block
247 *
248  kcnext = kc + k
249 *
250 * Multiply by the diagonal block if forming U * D.
251 *
252  IF( nounit ) THEN
253  d11 = a( kcnext-1 )
254  d22 = a( kcnext+k )
255  d12 = a( kcnext+k-1 )
256  d21 = dconjg( d12 )
257  DO 20 j = 1, nrhs
258  t1 = b( k, j )
259  t2 = b( k+1, j )
260  b( k, j ) = d11*t1 + d12*t2
261  b( k+1, j ) = d21*t1 + d22*t2
262  20 CONTINUE
263  END IF
264 *
265 * Multiply by P(K) * inv(U(K)) if K > 1.
266 *
267  IF( k.GT.1 ) THEN
268 *
269 * Apply the transformations.
270 *
271  CALL zgeru( k-1, nrhs, one, a( kc ), 1, b( k, 1 ),
272  $ ldb, b( 1, 1 ), ldb )
273  CALL zgeru( k-1, nrhs, one, a( kcnext ), 1,
274  $ b( k+1, 1 ), ldb, b( 1, 1 ), ldb )
275 *
276 * Interchange if P(K) != I.
277 *
278  kp = abs( ipiv( k ) )
279  IF( kp.NE.k )
280  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
281  END IF
282  kc = kcnext + k + 1
283  k = k + 2
284  END IF
285  GO TO 10
286  30 CONTINUE
287 *
288 * Compute B := L*B
289 * where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m)) .
290 *
291  ELSE
292 *
293 * Loop backward applying the transformations to B.
294 *
295  k = n
296  kc = n*( n+1 ) / 2 + 1
297  40 CONTINUE
298  IF( k.LT.1 )
299  $ GO TO 60
300  kc = kc - ( n-k+1 )
301 *
302 * Test the pivot index. If greater than zero, a 1 x 1
303 * pivot was used, otherwise a 2 x 2 pivot was used.
304 *
305  IF( ipiv( k ).GT.0 ) THEN
306 *
307 * 1 x 1 pivot block:
308 *
309 * Multiply by the diagonal element if forming L * D.
310 *
311  IF( nounit )
312  $ CALL zscal( nrhs, a( kc ), b( k, 1 ), ldb )
313 *
314 * Multiply by P(K) * inv(L(K)) if K < N.
315 *
316  IF( k.NE.n ) THEN
317  kp = ipiv( k )
318 *
319 * Apply the transformation.
320 *
321  CALL zgeru( n-k, nrhs, one, a( kc+1 ), 1, b( k, 1 ),
322  $ ldb, b( k+1, 1 ), ldb )
323 *
324 * Interchange if a permutation was applied at the
325 * K-th step of the factorization.
326 *
327  IF( kp.NE.k )
328  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
329  END IF
330  k = k - 1
331 *
332  ELSE
333 *
334 * 2 x 2 pivot block:
335 *
336  kcnext = kc - ( n-k+2 )
337 *
338 * Multiply by the diagonal block if forming L * D.
339 *
340  IF( nounit ) THEN
341  d11 = a( kcnext )
342  d22 = a( kc )
343  d21 = a( kcnext+1 )
344  d12 = dconjg( d21 )
345  DO 50 j = 1, nrhs
346  t1 = b( k-1, j )
347  t2 = b( k, j )
348  b( k-1, j ) = d11*t1 + d12*t2
349  b( k, j ) = d21*t1 + d22*t2
350  50 CONTINUE
351  END IF
352 *
353 * Multiply by P(K) * inv(L(K)) if K < N.
354 *
355  IF( k.NE.n ) THEN
356 *
357 * Apply the transformation.
358 *
359  CALL zgeru( n-k, nrhs, one, a( kc+1 ), 1, b( k, 1 ),
360  $ ldb, b( k+1, 1 ), ldb )
361  CALL zgeru( n-k, nrhs, one, a( kcnext+2 ), 1,
362  $ b( k-1, 1 ), ldb, b( k+1, 1 ), ldb )
363 *
364 * Interchange if a permutation was applied at the
365 * K-th step of the factorization.
366 *
367  kp = abs( ipiv( k ) )
368  IF( kp.NE.k )
369  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
370  END IF
371  kc = kcnext
372  k = k - 2
373  END IF
374  GO TO 40
375  60 CONTINUE
376  END IF
377 *-------------------------------------------------
378 *
379 * Compute B := A^H * B (conjugate transpose)
380 *
381 *-------------------------------------------------
382  ELSE
383 *
384 * Form B := U^H*B
385 * where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1))
386 * and U^H = inv(U^H(1))*P(1)* ... *inv(U^H(m))*P(m)
387 *
388  IF( lsame( uplo, 'U' ) ) THEN
389 *
390 * Loop backward applying the transformations.
391 *
392  k = n
393  kc = n*( n+1 ) / 2 + 1
394  70 CONTINUE
395  IF( k.LT.1 )
396  $ GO TO 90
397  kc = kc - k
398 *
399 * 1 x 1 pivot block.
400 *
401  IF( ipiv( k ).GT.0 ) THEN
402  IF( k.GT.1 ) THEN
403 *
404 * Interchange if P(K) != I.
405 *
406  kp = ipiv( k )
407  IF( kp.NE.k )
408  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
409 *
410 * Apply the transformation:
411 * y := y - B' * conjg(x)
412 * where x is a column of A and y is a row of B.
413 *
414  CALL zlacgv( nrhs, b( k, 1 ), ldb )
415  CALL zgemv( 'Conjugate', k-1, nrhs, one, b, ldb,
416  $ a( kc ), 1, one, b( k, 1 ), ldb )
417  CALL zlacgv( nrhs, b( k, 1 ), ldb )
418  END IF
419  IF( nounit )
420  $ CALL zscal( nrhs, a( kc+k-1 ), b( k, 1 ), ldb )
421  k = k - 1
422 *
423 * 2 x 2 pivot block.
424 *
425  ELSE
426  kcnext = kc - ( k-1 )
427  IF( k.GT.2 ) THEN
428 *
429 * Interchange if P(K) != I.
430 *
431  kp = abs( ipiv( k ) )
432  IF( kp.NE.k-1 )
433  $ CALL zswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ),
434  $ ldb )
435 *
436 * Apply the transformations.
437 *
438  CALL zlacgv( nrhs, b( k, 1 ), ldb )
439  CALL zgemv( 'Conjugate', k-2, nrhs, one, b, ldb,
440  $ a( kc ), 1, one, b( k, 1 ), ldb )
441  CALL zlacgv( nrhs, b( k, 1 ), ldb )
442 *
443  CALL zlacgv( nrhs, b( k-1, 1 ), ldb )
444  CALL zgemv( 'Conjugate', k-2, nrhs, one, b, ldb,
445  $ a( kcnext ), 1, one, b( k-1, 1 ), ldb )
446  CALL zlacgv( nrhs, b( k-1, 1 ), ldb )
447  END IF
448 *
449 * Multiply by the diagonal block if non-unit.
450 *
451  IF( nounit ) THEN
452  d11 = a( kc-1 )
453  d22 = a( kc+k-1 )
454  d12 = a( kc+k-2 )
455  d21 = dconjg( d12 )
456  DO 80 j = 1, nrhs
457  t1 = b( k-1, j )
458  t2 = b( k, j )
459  b( k-1, j ) = d11*t1 + d12*t2
460  b( k, j ) = d21*t1 + d22*t2
461  80 CONTINUE
462  END IF
463  kc = kcnext
464  k = k - 2
465  END IF
466  GO TO 70
467  90 CONTINUE
468 *
469 * Form B := L^H*B
470 * where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m))
471 * and L^H = inv(L(m))*P(m)* ... *inv(L(1))*P(1)
472 *
473  ELSE
474 *
475 * Loop forward applying the L-transformations.
476 *
477  k = 1
478  kc = 1
479  100 CONTINUE
480  IF( k.GT.n )
481  $ GO TO 120
482 *
483 * 1 x 1 pivot block
484 *
485  IF( ipiv( k ).GT.0 ) THEN
486  IF( k.LT.n ) THEN
487 *
488 * Interchange if P(K) != I.
489 *
490  kp = ipiv( k )
491  IF( kp.NE.k )
492  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
493 *
494 * Apply the transformation
495 *
496  CALL zlacgv( nrhs, b( k, 1 ), ldb )
497  CALL zgemv( 'Conjugate', n-k, nrhs, one, b( k+1, 1 ),
498  $ ldb, a( kc+1 ), 1, one, b( k, 1 ), ldb )
499  CALL zlacgv( nrhs, b( k, 1 ), ldb )
500  END IF
501  IF( nounit )
502  $ CALL zscal( nrhs, a( kc ), b( k, 1 ), ldb )
503  kc = kc + n - k + 1
504  k = k + 1
505 *
506 * 2 x 2 pivot block.
507 *
508  ELSE
509  kcnext = kc + n - k + 1
510  IF( k.LT.n-1 ) THEN
511 *
512 * Interchange if P(K) != I.
513 *
514  kp = abs( ipiv( k ) )
515  IF( kp.NE.k+1 )
516  $ CALL zswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ),
517  $ ldb )
518 *
519 * Apply the transformation
520 *
521  CALL zlacgv( nrhs, b( k+1, 1 ), ldb )
522  CALL zgemv( 'Conjugate', n-k-1, nrhs, one,
523  $ b( k+2, 1 ), ldb, a( kcnext+1 ), 1, one,
524  $ b( k+1, 1 ), ldb )
525  CALL zlacgv( nrhs, b( k+1, 1 ), ldb )
526 *
527  CALL zlacgv( nrhs, b( k, 1 ), ldb )
528  CALL zgemv( 'Conjugate', n-k-1, nrhs, one,
529  $ b( k+2, 1 ), ldb, a( kc+2 ), 1, one,
530  $ b( k, 1 ), ldb )
531  CALL zlacgv( nrhs, b( k, 1 ), ldb )
532  END IF
533 *
534 * Multiply by the diagonal block if non-unit.
535 *
536  IF( nounit ) THEN
537  d11 = a( kc )
538  d22 = a( kcnext )
539  d21 = a( kc+1 )
540  d12 = dconjg( d21 )
541  DO 110 j = 1, nrhs
542  t1 = b( k, j )
543  t2 = b( k+1, j )
544  b( k, j ) = d11*t1 + d12*t2
545  b( k+1, j ) = d21*t1 + d22*t2
546  110 CONTINUE
547  END IF
548  kc = kcnext + ( n-k )
549  k = k + 2
550  END IF
551  GO TO 100
552  120 CONTINUE
553  END IF
554 *
555  END IF
556  RETURN
557 *
558 * End of ZLAVHP
559 *
subroutine zlacgv(N, X, INCX)
ZLACGV conjugates a complex vector.
Definition: zlacgv.f:76
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
Definition: zswap.f:83
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL
Definition: zscal.f:80
subroutine zgeru(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
ZGERU
Definition: zgeru.f:132
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZGEMV
Definition: zgemv.f:160
Here is the call graph for this function:
Here is the caller graph for this function: