3425
3426
3427
3428
3429
3430
3431
3432 INTEGER INCX, N
3433 DOUBLE PRECISION ERRBND, PREC, USCLR
3434
3435
3436 COMPLEX*16 X( * )
3437
3438
3439
3440
3441
3442
3443
3444
3445
3446
3447
3448
3449
3450
3451
3452
3453
3454
3455
3456
3457
3458
3459
3460
3461
3462
3463
3464
3465
3466
3467
3468
3469
3470
3471
3472
3473
3474
3475
3476
3477
3478
3479
3480
3481
3482
3483
3484
3485
3486
3487
3488 DOUBLE PRECISION ONE, TWO, ZERO
3489 parameter( one = 1.0d+0, two = 2.0d+0,
3490 $ zero = 0.0d+0 )
3491
3492
3493 INTEGER IX
3494 DOUBLE PRECISION ABSXI, ADDBND, FACT, SCALE, SSQ, SUMSCA, SUMSSQ
3495
3496
3497 INTRINSIC abs, dble, dimag
3498
3499
3500
3501 usclr = zero
3502 sumssq = one
3503 sumsca = zero
3504 addbnd = two * two * two * prec
3505 fact = one + two * ( ( one + prec )**3 - one )
3506
3507 scale = zero
3508 ssq = one
3509 DO 10 ix = 1, 1 + ( n - 1 )*incx, incx
3510 IF( dble( x( ix ) ).NE.zero ) THEN
3511 absxi = abs( dble( x( ix ) ) )
3512 IF( scale.LT.absxi )THEN
3513 sumssq = one + ( ssq*( scale/absxi )**2 ) * fact
3514 errbnd = addbnd * sumssq
3515 sumssq = sumssq + errbnd
3516 ssq = one + ssq*( scale/absxi )**2
3517 sumsca = absxi
3518 scale = absxi
3519 ELSE
3520 sumssq = ssq + ( ( absxi/scale )**2 ) * fact
3521 errbnd = addbnd * sumssq
3522 sumssq = sumssq + errbnd
3523 ssq = ssq + ( absxi/scale )**2
3524 END IF
3525 END IF
3526 IF( dimag( x( ix ) ).NE.zero ) THEN
3527 absxi = abs( dimag( x( ix ) ) )
3528 IF( scale.LT.absxi )THEN
3529 sumssq = one + ( ssq*( scale/absxi )**2 ) * fact
3530 errbnd = addbnd * sumssq
3531 sumssq = sumssq + errbnd
3532 ssq = one + ssq*( scale/absxi )**2
3533 sumsca = absxi
3534 scale = absxi
3535 ELSE
3536 sumssq = ssq + ( ( absxi/scale )**2 ) * fact
3537 errbnd = addbnd * sumssq
3538 sumssq = sumssq + errbnd
3539 ssq = ssq + ( absxi/scale )**2
3540 END IF
3541 END IF
3542 10 CONTINUE
3543
3544 usclr = scale * sqrt( ssq )
3545
3546
3547
3548 errbnd = sqrt( sumssq ) * ( one + two * ( 1.00001d+0 * prec ) )
3549
3550 errbnd = ( sumsca * errbnd ) - usclr
3551
3552 RETURN
3553
3554
3555