Minor accuracy improvement for statistics.correlation()

[rhettinger]

Perform sqrt(x * y) with only a single rounding. For example with x=0.1836988816815257 and y=0.9200609291953, the new computation gives 0.4111133222993899 instead of 0.41111332229938985.

Here's a script to demonstrate the improvement:

from math import sumprod, sqrt
from random import random
from decimal import Decimal, getcontext

getcontext().prec = 100

def sqrtprod_baseline(x: float, y: float) -> float:
    return sqrt(x * y)

def sqrtprod_decimal(x: float, y: float) -> float:
    x, y = map(Decimal, (x, y))
    return float((x * y).sqrt())

def sqrtprod_new(x: float, y: float) -> float:
    h = sqrt(x * y)
    x = sumprod((x, h), (y, -h))
    return h + x / (2.0 * h)

for i in range(100):
    x, y = random(), random()
    assert sqrtprod_new(x, y) == sqrtprod_decimal(x, y)
    if sqrtprod_baseline(x, y) != sqrtprod_decimal(x, y):
        print(i, x, y)
        print(sqrtprod_baseline(x, y))
        print(sqrtprod_decimal(x, y))
        print(sqrtprod_new(x, y))
        print()
    

Sample output: