1 *> \brief \b ZSPR performs the symmetrical rank-1 update of a complex symmetric packed matrix.
3 * =========== DOCUMENTATION ===========
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
9 *> Download ZSPR + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zspr.f">
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zspr.f">
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zspr.f">
21 * SUBROUTINE ZSPR( UPLO, N, ALPHA, X, INCX, AP )
23 * .. Scalar Arguments ..
28 * .. Array Arguments ..
29 * COMPLEX*16 AP( * ), X( * )
38 *> ZSPR performs the symmetric rank 1 operation
40 *> A := alpha*x*x**H + A,
42 *> where alpha is a complex scalar, x is an n element vector and A is an
43 *> n by n symmetric matrix, supplied in packed form.
51 *> UPLO is CHARACTER*1
52 *> On entry, UPLO specifies whether the upper or lower
53 *> triangular part of the matrix A is supplied in the packed
54 *> array AP as follows:
56 *> UPLO = 'U' or 'u' The upper triangular part of A is
59 *> UPLO = 'L' or 'l' The lower triangular part of A is
68 *> On entry, N specifies the order of the matrix A.
69 *> N must be at least zero.
75 *> ALPHA is COMPLEX*16
76 *> On entry, ALPHA specifies the scalar alpha.
82 *> X is COMPLEX*16 array, dimension at least
83 *> ( 1 + ( N - 1 )*abs( INCX ) ).
84 *> Before entry, the incremented array X must contain the N-
92 *> On entry, INCX specifies the increment for the elements of
93 *> X. INCX must not be zero.
99 *> AP is COMPLEX*16 array, dimension at least
100 *> ( ( N*( N + 1 ) )/2 ).
101 *> Before entry, with UPLO = 'U' or 'u', the array AP must
102 *> contain the upper triangular part of the symmetric matrix
103 *> packed sequentially, column by column, so that AP( 1 )
104 *> contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
105 *> and a( 2, 2 ) respectively, and so on. On exit, the array
106 *> AP is overwritten by the upper triangular part of the
108 *> Before entry, with UPLO = 'L' or 'l', the array AP must
109 *> contain the lower triangular part of the symmetric matrix
110 *> packed sequentially, column by column, so that AP( 1 )
111 *> contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
112 *> and a( 3, 1 ) respectively, and so on. On exit, the array
113 *> AP is overwritten by the lower triangular part of the
115 *> Note that the imaginary parts of the diagonal elements need
116 *> not be set, they are assumed to be zero, and on exit they
123 *> \author Univ. of Tennessee
124 *> \author Univ. of California Berkeley
125 *> \author Univ. of Colorado Denver
128 *> \date September 2012
130 *> \ingroup complex16OTHERauxiliary
132 * =====================================================================
133 SUBROUTINE ZSPR( UPLO, N, ALPHA, X, INCX, AP )
135 * -- LAPACK auxiliary routine (version 3.4.2) --
136 * -- LAPACK is a software package provided by Univ. of Tennessee, --
137 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
140 * .. Scalar Arguments ..
145 * .. Array Arguments ..
146 COMPLEX*16 AP( * ), X( * )
149 * =====================================================================
153 PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ) )
155 * .. Local Scalars ..
156 INTEGER I, INFO, IX, J, JX, K, KK, KX
159 * .. External Functions ..
163 * .. External Subroutines ..
166 * .. Executable Statements ..
168 * Test the input parameters.
171 IF( .NOT.LSAME( UPLO, 'U' ) .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
173 ELSE IF( N.LT.0 ) THEN
175 ELSE IF( INCX.EQ.0 ) THEN
179 CALL XERBLA( 'ZSPR ', INFO )
183 * Quick return if possible.
185 IF( ( N.EQ.0 ) .OR. ( ALPHA.EQ.ZERO ) )
188 * Set the start point in X if the increment is not unity.
191 KX = 1 - ( N-1 )*INCX
192 ELSE IF( INCX.NE.1 ) THEN
196 * Start the operations. In this version the elements of the array AP
197 * are accessed sequentially with one pass through AP.
200 IF( LSAME( UPLO, 'U' ) ) THEN
202 * Form A when upper triangle is stored in AP.
206 IF( X( J ).NE.ZERO ) THEN
210 AP( K ) = AP( K ) + X( I )*TEMP
213 AP( KK+J-1 ) = AP( KK+J-1 ) + X( J )*TEMP
215 AP( KK+J-1 ) = AP( KK+J-1 )
222 IF( X( JX ).NE.ZERO ) THEN
225 DO 30 K = KK, KK + J - 2
226 AP( K ) = AP( K ) + X( IX )*TEMP
229 AP( KK+J-1 ) = AP( KK+J-1 ) + X( JX )*TEMP
231 AP( KK+J-1 ) = AP( KK+J-1 )
239 * Form A when lower triangle is stored in AP.
243 IF( X( J ).NE.ZERO ) THEN
245 AP( KK ) = AP( KK ) + TEMP*X( J )
248 AP( K ) = AP( K ) + X( I )*TEMP
259 IF( X( JX ).NE.ZERO ) THEN
261 AP( KK ) = AP( KK ) + TEMP*X( JX )
263 DO 70 K = KK + 1, KK + N - J
265 AP( K ) = AP( K ) + X( IX )*TEMP