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