3137
3138
3139
3140
3141
3142
3143
3144 INTEGER INCX, INCY, N
3145 DOUBLE PRECISION ERRBND, PREC
3146 COMPLEX*16 SCLR
3147
3148
3149 COMPLEX*16 X( * ), Y( * )
3150
3151
3152
3153
3154
3155
3156
3157
3158
3159
3160
3161
3162
3163
3164
3165
3166
3167
3168
3169
3170
3171
3172
3173
3174
3175
3176
3177
3178
3179
3180
3181
3182
3183
3184
3185
3186
3187
3188
3189
3190
3191
3192
3193
3194
3195
3196
3197
3198
3199
3200
3201
3202
3203
3204
3205
3206
3207
3208
3209
3210
3211 DOUBLE PRECISION ONE, TWO, ZERO
3212 parameter( one = 1.0d+0, two = 2.0d+0,
3213 $ zero = 0.0d+0 )
3214
3215
3216 INTEGER I, IX, IY
3217 DOUBLE PRECISION ADDBND, FACT, SUMINEG, SUMIPOS, SUMRNEG,
3218 $ SUMRPOS, TMP
3219
3220
3221 INTRINSIC abs, dble, dimag,
max
3222
3223
3224
3225 ix = 1
3226 iy = 1
3227 sclr = zero
3228 sumipos = zero
3229 sumineg = zero
3230 sumrpos = zero
3231 sumrneg = zero
3232 fact = two * ( one + prec )
3233 addbnd = two * two * two * prec
3234
3235 DO 10 i = 1, n
3236
3237 sclr = sclr + x( ix ) * y( iy )
3238
3239 tmp = dble( x( ix ) ) * dble( y( iy ) )
3240 IF( tmp.GE.zero ) THEN
3241 sumrpos = sumrpos + tmp * fact
3242 ELSE
3243 sumrneg = sumrneg - tmp * fact
3244 END IF
3245
3246 tmp = - dimag( x( ix ) ) * dimag( y( iy ) )
3247 IF( tmp.GE.zero ) THEN
3248 sumrpos = sumrpos + tmp * fact
3249 ELSE
3250 sumrneg = sumrneg - tmp * fact
3251 END IF
3252
3253 tmp = dimag( x( ix ) ) * dble( y( iy ) )
3254 IF( tmp.GE.zero ) THEN
3255 sumipos = sumipos + tmp * fact
3256 ELSE
3257 sumineg = sumineg - tmp * fact
3258 END IF
3259
3260 tmp = dble( x( ix ) ) * dimag( y( iy ) )
3261 IF( tmp.GE.zero ) THEN
3262 sumipos = sumipos + tmp * fact
3263 ELSE
3264 sumineg = sumineg - tmp * fact
3265 END IF
3266
3267 ix = ix + incx
3268 iy = iy + incy
3269
3270 10 CONTINUE
3271
3272 errbnd = addbnd *
max(
max( sumrpos, sumrneg ),
3273 $
max( sumipos, sumineg ) )
3274
3275 RETURN
3276
3277
3278