1
SUBROUTINE ZHPR(UPLO,N,ALPHA,X,INCX,AP)
2
* .. Scalar Arguments ..
7
* .. Array Arguments ..
8
DOUBLE COMPLEX AP(*),X(*)
14
* ZHPR performs the hermitian rank 1 operation
16
* A := alpha*x*conjg( x' ) + A,
18
* where alpha is a real scalar, x is an n element vector and A is an
19
* n by n hermitian matrix, supplied in packed form.
25
* On entry, UPLO specifies whether the upper or lower
26
* triangular part of the matrix A is supplied in the packed
27
* array AP as follows:
29
* UPLO = 'U' or 'u' The upper triangular part of A is
32
* UPLO = 'L' or 'l' The lower triangular part of A is
38
* On entry, N specifies the order of the matrix A.
39
* N must be at least zero.
42
* ALPHA - DOUBLE PRECISION.
43
* On entry, ALPHA specifies the scalar alpha.
46
* X - COMPLEX*16 array of dimension at least
47
* ( 1 + ( n - 1 )*abs( INCX ) ).
48
* Before entry, the incremented array X must contain the n
53
* On entry, INCX specifies the increment for the elements of
54
* X. INCX must not be zero.
57
* AP - COMPLEX*16 array of DIMENSION at least
58
* ( ( n*( n + 1 ) )/2 ).
59
* Before entry with UPLO = 'U' or 'u', the array AP must
60
* contain the upper triangular part of the hermitian matrix
61
* packed sequentially, column by column, so that AP( 1 )
62
* contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
63
* and a( 2, 2 ) respectively, and so on. On exit, the array
64
* AP is overwritten by the upper triangular part of the
66
* Before entry with UPLO = 'L' or 'l', the array AP must
67
* contain the lower triangular part of the hermitian matrix
68
* packed sequentially, column by column, so that AP( 1 )
69
* contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
70
* and a( 3, 1 ) respectively, and so on. On exit, the array
71
* AP is overwritten by the lower triangular part of the
73
* Note that the imaginary parts of the diagonal elements need
74
* not be set, they are assumed to be zero, and on exit they
78
* Level 2 Blas routine.
80
* -- Written on 22-October-1986.
81
* Jack Dongarra, Argonne National Lab.
82
* Jeremy Du Croz, Nag Central Office.
83
* Sven Hammarling, Nag Central Office.
84
* Richard Hanson, Sandia National Labs.
89
PARAMETER (ZERO= (0.0D+0,0.0D+0))
93
INTEGER I,INFO,IX,J,JX,K,KK,KX
95
* .. External Functions ..
99
* .. External Subroutines ..
102
* .. Intrinsic Functions ..
103
INTRINSIC DBLE,DCONJG
106
* Test the input parameters.
109
IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
111
ELSE IF (N.LT.0) THEN
113
ELSE IF (INCX.EQ.0) THEN
117
CALL XERBLA('ZHPR ',INFO)
121
* Quick return if possible.
123
IF ((N.EQ.0) .OR. (ALPHA.EQ.DBLE(ZERO))) RETURN
125
* Set the start point in X if the increment is not unity.
129
ELSE IF (INCX.NE.1) THEN
133
* Start the operations. In this version the elements of the array AP
134
* are accessed sequentially with one pass through AP.
137
IF (LSAME(UPLO,'U')) THEN
139
* Form A when upper triangle is stored in AP.
143
IF (X(J).NE.ZERO) THEN
144
TEMP = ALPHA*DCONJG(X(J))
147
AP(K) = AP(K) + X(I)*TEMP
150
AP(KK+J-1) = DBLE(AP(KK+J-1)) + DBLE(X(J)*TEMP)
152
AP(KK+J-1) = DBLE(AP(KK+J-1))
159
IF (X(JX).NE.ZERO) THEN
160
TEMP = ALPHA*DCONJG(X(JX))
162
DO 30 K = KK,KK + J - 2
163
AP(K) = AP(K) + X(IX)*TEMP
166
AP(KK+J-1) = DBLE(AP(KK+J-1)) + DBLE(X(JX)*TEMP)
168
AP(KK+J-1) = DBLE(AP(KK+J-1))
176
* Form A when lower triangle is stored in AP.
180
IF (X(J).NE.ZERO) THEN
181
TEMP = ALPHA*DCONJG(X(J))
182
AP(KK) = DBLE(AP(KK)) + DBLE(TEMP*X(J))
185
AP(K) = AP(K) + X(I)*TEMP
189
AP(KK) = DBLE(AP(KK))
196
IF (X(JX).NE.ZERO) THEN
197
TEMP = ALPHA*DCONJG(X(JX))
198
AP(KK) = DBLE(AP(KK)) + DBLE(TEMP*X(JX))
200
DO 70 K = KK + 1,KK + N - J
202
AP(K) = AP(K) + X(IX)*TEMP
205
AP(KK) = DBLE(AP(KK))