LAPACK 3.3.0
|
00001 SUBROUTINE ZSPMV( UPLO, N, ALPHA, AP, X, INCX, BETA, Y, INCY ) 00002 * 00003 * -- LAPACK auxiliary routine (version 3.2) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 CHARACTER UPLO 00010 INTEGER INCX, INCY, N 00011 COMPLEX*16 ALPHA, BETA 00012 * .. 00013 * .. Array Arguments .. 00014 COMPLEX*16 AP( * ), X( * ), Y( * ) 00015 * .. 00016 * 00017 * Purpose 00018 * ======= 00019 * 00020 * ZSPMV performs the matrix-vector operation 00021 * 00022 * y := alpha*A*x + beta*y, 00023 * 00024 * where alpha and beta are scalars, x and y are n element vectors and 00025 * A is an n by n symmetric matrix, supplied in packed form. 00026 * 00027 * Arguments 00028 * ========== 00029 * 00030 * UPLO (input) CHARACTER*1 00031 * On entry, UPLO specifies whether the upper or lower 00032 * triangular part of the matrix A is supplied in the packed 00033 * array AP as follows: 00034 * 00035 * UPLO = 'U' or 'u' The upper triangular part of A is 00036 * supplied in AP. 00037 * 00038 * UPLO = 'L' or 'l' The lower triangular part of A is 00039 * supplied in AP. 00040 * 00041 * Unchanged on exit. 00042 * 00043 * N (input) INTEGER 00044 * On entry, N specifies the order of the matrix A. 00045 * N must be at least zero. 00046 * Unchanged on exit. 00047 * 00048 * ALPHA (input) COMPLEX*16 00049 * On entry, ALPHA specifies the scalar alpha. 00050 * Unchanged on exit. 00051 * 00052 * AP (input) COMPLEX*16 array, dimension at least 00053 * ( ( N*( N + 1 ) )/2 ). 00054 * Before entry, with UPLO = 'U' or 'u', the array AP must 00055 * contain the upper triangular part of the symmetric matrix 00056 * packed sequentially, column by column, so that AP( 1 ) 00057 * contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 ) 00058 * and a( 2, 2 ) respectively, and so on. 00059 * Before entry, with UPLO = 'L' or 'l', the array AP must 00060 * contain the lower triangular part of the symmetric matrix 00061 * packed sequentially, column by column, so that AP( 1 ) 00062 * contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 ) 00063 * and a( 3, 1 ) respectively, and so on. 00064 * Unchanged on exit. 00065 * 00066 * X (input) COMPLEX*16 array, dimension at least 00067 * ( 1 + ( N - 1 )*abs( INCX ) ). 00068 * Before entry, the incremented array X must contain the N- 00069 * element vector x. 00070 * Unchanged on exit. 00071 * 00072 * INCX (input) INTEGER 00073 * On entry, INCX specifies the increment for the elements of 00074 * X. INCX must not be zero. 00075 * Unchanged on exit. 00076 * 00077 * BETA (input) COMPLEX*16 00078 * On entry, BETA specifies the scalar beta. When BETA is 00079 * supplied as zero then Y need not be set on input. 00080 * Unchanged on exit. 00081 * 00082 * Y (input/output) COMPLEX*16 array, dimension at least 00083 * ( 1 + ( N - 1 )*abs( INCY ) ). 00084 * Before entry, the incremented array Y must contain the n 00085 * element vector y. On exit, Y is overwritten by the updated 00086 * vector y. 00087 * 00088 * INCY (input) INTEGER 00089 * On entry, INCY specifies the increment for the elements of 00090 * Y. INCY must not be zero. 00091 * Unchanged on exit. 00092 * 00093 * ===================================================================== 00094 * 00095 * .. Parameters .. 00096 COMPLEX*16 ONE 00097 PARAMETER ( ONE = ( 1.0D+0, 0.0D+0 ) ) 00098 COMPLEX*16 ZERO 00099 PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ) ) 00100 * .. 00101 * .. Local Scalars .. 00102 INTEGER I, INFO, IX, IY, J, JX, JY, K, KK, KX, KY 00103 COMPLEX*16 TEMP1, TEMP2 00104 * .. 00105 * .. External Functions .. 00106 LOGICAL LSAME 00107 EXTERNAL LSAME 00108 * .. 00109 * .. External Subroutines .. 00110 EXTERNAL XERBLA 00111 * .. 00112 * .. Executable Statements .. 00113 * 00114 * Test the input parameters. 00115 * 00116 INFO = 0 00117 IF( .NOT.LSAME( UPLO, 'U' ) .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN 00118 INFO = 1 00119 ELSE IF( N.LT.0 ) THEN 00120 INFO = 2 00121 ELSE IF( INCX.EQ.0 ) THEN 00122 INFO = 6 00123 ELSE IF( INCY.EQ.0 ) THEN 00124 INFO = 9 00125 END IF 00126 IF( INFO.NE.0 ) THEN 00127 CALL XERBLA( 'ZSPMV ', INFO ) 00128 RETURN 00129 END IF 00130 * 00131 * Quick return if possible. 00132 * 00133 IF( ( N.EQ.0 ) .OR. ( ( ALPHA.EQ.ZERO ) .AND. ( BETA.EQ.ONE ) ) ) 00134 $ RETURN 00135 * 00136 * Set up the start points in X and Y. 00137 * 00138 IF( INCX.GT.0 ) THEN 00139 KX = 1 00140 ELSE 00141 KX = 1 - ( N-1 )*INCX 00142 END IF 00143 IF( INCY.GT.0 ) THEN 00144 KY = 1 00145 ELSE 00146 KY = 1 - ( N-1 )*INCY 00147 END IF 00148 * 00149 * Start the operations. In this version the elements of the array AP 00150 * are accessed sequentially with one pass through AP. 00151 * 00152 * First form y := beta*y. 00153 * 00154 IF( BETA.NE.ONE ) THEN 00155 IF( INCY.EQ.1 ) THEN 00156 IF( BETA.EQ.ZERO ) THEN 00157 DO 10 I = 1, N 00158 Y( I ) = ZERO 00159 10 CONTINUE 00160 ELSE 00161 DO 20 I = 1, N 00162 Y( I ) = BETA*Y( I ) 00163 20 CONTINUE 00164 END IF 00165 ELSE 00166 IY = KY 00167 IF( BETA.EQ.ZERO ) THEN 00168 DO 30 I = 1, N 00169 Y( IY ) = ZERO 00170 IY = IY + INCY 00171 30 CONTINUE 00172 ELSE 00173 DO 40 I = 1, N 00174 Y( IY ) = BETA*Y( IY ) 00175 IY = IY + INCY 00176 40 CONTINUE 00177 END IF 00178 END IF 00179 END IF 00180 IF( ALPHA.EQ.ZERO ) 00181 $ RETURN 00182 KK = 1 00183 IF( LSAME( UPLO, 'U' ) ) THEN 00184 * 00185 * Form y when AP contains the upper triangle. 00186 * 00187 IF( ( INCX.EQ.1 ) .AND. ( INCY.EQ.1 ) ) THEN 00188 DO 60 J = 1, N 00189 TEMP1 = ALPHA*X( J ) 00190 TEMP2 = ZERO 00191 K = KK 00192 DO 50 I = 1, J - 1 00193 Y( I ) = Y( I ) + TEMP1*AP( K ) 00194 TEMP2 = TEMP2 + AP( K )*X( I ) 00195 K = K + 1 00196 50 CONTINUE 00197 Y( J ) = Y( J ) + TEMP1*AP( KK+J-1 ) + ALPHA*TEMP2 00198 KK = KK + J 00199 60 CONTINUE 00200 ELSE 00201 JX = KX 00202 JY = KY 00203 DO 80 J = 1, N 00204 TEMP1 = ALPHA*X( JX ) 00205 TEMP2 = ZERO 00206 IX = KX 00207 IY = KY 00208 DO 70 K = KK, KK + J - 2 00209 Y( IY ) = Y( IY ) + TEMP1*AP( K ) 00210 TEMP2 = TEMP2 + AP( K )*X( IX ) 00211 IX = IX + INCX 00212 IY = IY + INCY 00213 70 CONTINUE 00214 Y( JY ) = Y( JY ) + TEMP1*AP( KK+J-1 ) + ALPHA*TEMP2 00215 JX = JX + INCX 00216 JY = JY + INCY 00217 KK = KK + J 00218 80 CONTINUE 00219 END IF 00220 ELSE 00221 * 00222 * Form y when AP contains the lower triangle. 00223 * 00224 IF( ( INCX.EQ.1 ) .AND. ( INCY.EQ.1 ) ) THEN 00225 DO 100 J = 1, N 00226 TEMP1 = ALPHA*X( J ) 00227 TEMP2 = ZERO 00228 Y( J ) = Y( J ) + TEMP1*AP( KK ) 00229 K = KK + 1 00230 DO 90 I = J + 1, N 00231 Y( I ) = Y( I ) + TEMP1*AP( K ) 00232 TEMP2 = TEMP2 + AP( K )*X( I ) 00233 K = K + 1 00234 90 CONTINUE 00235 Y( J ) = Y( J ) + ALPHA*TEMP2 00236 KK = KK + ( N-J+1 ) 00237 100 CONTINUE 00238 ELSE 00239 JX = KX 00240 JY = KY 00241 DO 120 J = 1, N 00242 TEMP1 = ALPHA*X( JX ) 00243 TEMP2 = ZERO 00244 Y( JY ) = Y( JY ) + TEMP1*AP( KK ) 00245 IX = JX 00246 IY = JY 00247 DO 110 K = KK + 1, KK + N - J 00248 IX = IX + INCX 00249 IY = IY + INCY 00250 Y( IY ) = Y( IY ) + TEMP1*AP( K ) 00251 TEMP2 = TEMP2 + AP( K )*X( IX ) 00252 110 CONTINUE 00253 Y( JY ) = Y( JY ) + ALPHA*TEMP2 00254 JX = JX + INCX 00255 JY = JY + INCY 00256 KK = KK + ( N-J+1 ) 00257 120 CONTINUE 00258 END IF 00259 END IF 00260 * 00261 RETURN 00262 * 00263 * End of ZSPMV 00264 * 00265 END