1: SUBROUTINE ZHER2(UPLO,N,ALPHA,X,INCX,Y,INCY,A,LDA)
2: * .. Scalar Arguments ..
3: DOUBLE COMPLEX ALPHA
4: INTEGER INCX,INCY,LDA,N
5: CHARACTER UPLO
6: * ..
7: * .. Array Arguments ..
8: DOUBLE COMPLEX A(LDA,*),X(*),Y(*)
9: * ..
10: *
11: * Purpose
12: * =======
13: *
14: * ZHER2 performs the hermitian rank 2 operation
15: *
16: * A := alpha*x*conjg( y' ) + conjg( alpha )*y*conjg( x' ) + A,
17: *
18: * where alpha is a scalar, x and y are n element vectors and A is an n
19: * by n hermitian matrix.
20: *
21: * Arguments
22: * ==========
23: *
24: * UPLO - CHARACTER*1.
25: * On entry, UPLO specifies whether the upper or lower
26: * triangular part of the array A is to be referenced as
27: * follows:
28: *
29: * UPLO = 'U' or 'u' Only the upper triangular part of A
30: * is to be referenced.
31: *
32: * UPLO = 'L' or 'l' Only the lower triangular part of A
33: * is to be referenced.
34: *
35: * Unchanged on exit.
36: *
37: * N - INTEGER.
38: * On entry, N specifies the order of the matrix A.
39: * N must be at least zero.
40: * Unchanged on exit.
41: *
42: * ALPHA - COMPLEX*16 .
43: * On entry, ALPHA specifies the scalar alpha.
44: * Unchanged on exit.
45: *
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
49: * element vector x.
50: * Unchanged on exit.
51: *
52: * INCX - INTEGER.
53: * On entry, INCX specifies the increment for the elements of
54: * X. INCX must not be zero.
55: * Unchanged on exit.
56: *
57: * Y - COMPLEX*16 array of dimension at least
58: * ( 1 + ( n - 1 )*abs( INCY ) ).
59: * Before entry, the incremented array Y must contain the n
60: * element vector y.
61: * Unchanged on exit.
62: *
63: * INCY - INTEGER.
64: * On entry, INCY specifies the increment for the elements of
65: * Y. INCY must not be zero.
66: * Unchanged on exit.
67: *
68: * A - COMPLEX*16 array of DIMENSION ( LDA, n ).
69: * Before entry with UPLO = 'U' or 'u', the leading n by n
70: * upper triangular part of the array A must contain the upper
71: * triangular part of the hermitian matrix and the strictly
72: * lower triangular part of A is not referenced. On exit, the
73: * upper triangular part of the array A is overwritten by the
74: * upper triangular part of the updated matrix.
75: * Before entry with UPLO = 'L' or 'l', the leading n by n
76: * lower triangular part of the array A must contain the lower
77: * triangular part of the hermitian matrix and the strictly
78: * upper triangular part of A is not referenced. On exit, the
79: * lower triangular part of the array A is overwritten by the
80: * lower triangular part of the updated matrix.
81: * Note that the imaginary parts of the diagonal elements need
82: * not be set, they are assumed to be zero, and on exit they
83: * are set to zero.
84: *
85: * LDA - INTEGER.
86: * On entry, LDA specifies the first dimension of A as declared
87: * in the calling (sub) program. LDA must be at least
88: * max( 1, n ).
89: * Unchanged on exit.
90: *
91: * Further Details
92: * ===============
93: *
94: * Level 2 Blas routine.
95: *
96: * -- Written on 22-October-1986.
97: * Jack Dongarra, Argonne National Lab.
98: * Jeremy Du Croz, Nag Central Office.
99: * Sven Hammarling, Nag Central Office.
100: * Richard Hanson, Sandia National Labs.
101: *
102: * =====================================================================
103: *
104: * .. Parameters ..
105: DOUBLE COMPLEX ZERO
106: PARAMETER (ZERO= (0.0D+0,0.0D+0))
107: * ..
108: * .. Local Scalars ..
109: DOUBLE COMPLEX TEMP1,TEMP2
110: INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY
111: * ..
112: * .. External Functions ..
113: LOGICAL LSAME
114: EXTERNAL LSAME
115: * ..
116: * .. External Subroutines ..
117: EXTERNAL XERBLA
118: * ..
119: * .. Intrinsic Functions ..
120: INTRINSIC DBLE,DCONJG,MAX
121: * ..
122: *
123: * Test the input parameters.
124: *
125: INFO = 0
126: IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
127: INFO = 1
128: ELSE IF (N.LT.0) THEN
129: INFO = 2
130: ELSE IF (INCX.EQ.0) THEN
131: INFO = 5
132: ELSE IF (INCY.EQ.0) THEN
133: INFO = 7
134: ELSE IF (LDA.LT.MAX(1,N)) THEN
135: INFO = 9
136: END IF
137: IF (INFO.NE.0) THEN
138: CALL XERBLA('ZHER2 ',INFO)
139: RETURN
140: END IF
141: *
142: * Quick return if possible.
143: *
144: IF ((N.EQ.0) .OR. (ALPHA.EQ.ZERO)) RETURN
145: *
146: * Set up the start points in X and Y if the increments are not both
147: * unity.
148: *
149: IF ((INCX.NE.1) .OR. (INCY.NE.1)) THEN
150: IF (INCX.GT.0) THEN
151: KX = 1
152: ELSE
153: KX = 1 - (N-1)*INCX
154: END IF
155: IF (INCY.GT.0) THEN
156: KY = 1
157: ELSE
158: KY = 1 - (N-1)*INCY
159: END IF
160: JX = KX
161: JY = KY
162: END IF
163: *
164: * Start the operations. In this version the elements of A are
165: * accessed sequentially with one pass through the triangular part
166: * of A.
167: *
168: IF (LSAME(UPLO,'U')) THEN
169: *
170: * Form A when A is stored in the upper triangle.
171: *
172: IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
173: DO 20 J = 1,N
174: IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
175: TEMP1 = ALPHA*DCONJG(Y(J))
176: TEMP2 = DCONJG(ALPHA*X(J))
177: DO 10 I = 1,J - 1
178: A(I,J) = A(I,J) + X(I)*TEMP1 + Y(I)*TEMP2
179: 10 CONTINUE
180: A(J,J) = DBLE(A(J,J)) +
181: + DBLE(X(J)*TEMP1+Y(J)*TEMP2)
182: ELSE
183: A(J,J) = DBLE(A(J,J))
184: END IF
185: 20 CONTINUE
186: ELSE
187: DO 40 J = 1,N
188: IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
189: TEMP1 = ALPHA*DCONJG(Y(JY))
190: TEMP2 = DCONJG(ALPHA*X(JX))
191: IX = KX
192: IY = KY
193: DO 30 I = 1,J - 1
194: A(I,J) = A(I,J) + X(IX)*TEMP1 + Y(IY)*TEMP2
195: IX = IX + INCX
196: IY = IY + INCY
197: 30 CONTINUE
198: A(J,J) = DBLE(A(J,J)) +
199: + DBLE(X(JX)*TEMP1+Y(JY)*TEMP2)
200: ELSE
201: A(J,J) = DBLE(A(J,J))
202: END IF
203: JX = JX + INCX
204: JY = JY + INCY
205: 40 CONTINUE
206: END IF
207: ELSE
208: *
209: * Form A when A is stored in the lower triangle.
210: *
211: IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
212: DO 60 J = 1,N
213: IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
214: TEMP1 = ALPHA*DCONJG(Y(J))
215: TEMP2 = DCONJG(ALPHA*X(J))
216: A(J,J) = DBLE(A(J,J)) +
217: + DBLE(X(J)*TEMP1+Y(J)*TEMP2)
218: DO 50 I = J + 1,N
219: A(I,J) = A(I,J) + X(I)*TEMP1 + Y(I)*TEMP2
220: 50 CONTINUE
221: ELSE
222: A(J,J) = DBLE(A(J,J))
223: END IF
224: 60 CONTINUE
225: ELSE
226: DO 80 J = 1,N
227: IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
228: TEMP1 = ALPHA*DCONJG(Y(JY))
229: TEMP2 = DCONJG(ALPHA*X(JX))
230: A(J,J) = DBLE(A(J,J)) +
231: + DBLE(X(JX)*TEMP1+Y(JY)*TEMP2)
232: IX = JX
233: IY = JY
234: DO 70 I = J + 1,N
235: IX = IX + INCX
236: IY = IY + INCY
237: A(I,J) = A(I,J) + X(IX)*TEMP1 + Y(IY)*TEMP2
238: 70 CONTINUE
239: ELSE
240: A(J,J) = DBLE(A(J,J))
241: END IF
242: JX = JX + INCX
243: JY = JY + INCY
244: 80 CONTINUE
245: END IF
246: END IF
247: *
248: RETURN
249: *
250: * End of ZHER2 .
251: *
252: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zher2.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / blas