[Python-ideas] Yet another sum function (fractions.sum)

Oscar Benjamin oscar.j.benjamin at gmail.com
Tue Aug 20 18:39:59 CEST 2013


On 19 August 2013 10:35, Stefan Behnel <stefan_ml at behnel.de> wrote:
> David Mertz, 16.08.2013 22:06:
>> My intuition was that one might do better than mergesum() by binning
>> numbers.  I.e. this implementation:
>>
>>   def binsum(seq):
>>       bins = dict()
>>       for f in seq:
>>           bins[f.denominator] = bins.get(f.denominator, 0) + f
>>       return mergesum(list(bins.values()))

Good thinking David. Actually for some distributions of inputs this
massively outperforms the other implementations. One common use case I
have for Fractions is to check the accuracy of floating point
computation by repeating the computation with Fractions. In this case
I would be working with Fractions whose denominators are always powers
of 2 (and usually only a handful of different values). And as Stefan
says if we're binning on the denominator then we can make it really
fast by adding the numerators with int.__add__.

> Simply sorting by denominator gives me a visible advantage over the above:
>
> def sortsum(seq):
>     def key(f):
>         return f.denominator if isinstance(f, F) else 1
>     seq = sorted(seq, key=key)
>     if len(seq) < 3:
>         return sum(seq)
>     return mergesum(seq)

Interesting.

I've developed an iterator version of mergesum:

def imergesum(iterable):
    stack = []
    it = iter(iterable)
    nextpow2 = 1
    while True:
        for n2, val in enumerate(it, 1):
            if n2 == nextpow2:
                stack.append(val)
                nextpow2 *= 2
                break
            elif n2 % 2:
                stack[-1] += val
            else:
                ishift = -1
                while not n2 % 2:
                    val, stack[ishift] = stack[ishift], val
                    ishift -= 1
                    n2 >>= 1
                stack[ishift] += val
        else:
            return mergesum(stack)  # Could just use sum here

It uses log(N) storage so you'd run out of space on an iterator of say
2 ** sys.maxint Fractions but that's unlikely to bother anyone.

I combined the different approaches to make a rationalsum() function
that is scalable and tries to take advantage of the binsum() and
sortsum() improvements where possible:

def rationalsum(iterable, tablesize=1000):
    def binsumgen(iterable):
        iterator = iter(iterable)
        bins = defaultdict(int)
        finished = False
        denom = itemgetter(0)
        while not finished:
            finished = True
            for n, num in zip(range(tablesize), iterator):
                bins[num.denominator] += num.numerator
                finished = False
            if len(bins) >= tablesize or finished:
                dn = sorted(bins.items(), key=denom)
                yield mergesum([F(n, d) for d, n in dn])
                bins.clear()
    return imergesum(binsumgen(iterable))

I tested this on a suite of different rational sequences (full script
at the bottom):
1) float - Fractions made from floats with power of 2 denominators.
2) d1000 - the original numbers with denominators uniform in [1, 1000].
3) d10 - like above but [1, 10].
4) bigint - big ints rather than Fractions (just for comparison).
5) e - the sequence 1 + 1 + 1/2 + 1/3! + 1/4! + ... (this is slow).

In different situations the binsum(), sortedsum() and mergesum()
functions parform differently (in some cases wildly so). The
rationalsum() function is never the best but always has the same order
of magnitude as the best. Apart from the bigint case builtins.sum() is
always among the worst performing.

Here's the benchmark output (Python 3.3 on Windows XP 32-bit - the
script takes a while):

Running benchmarks (times in microseconds)...


--------------
float
---------------
     n | sum    | merges | imerge | binsum | sortsu | ration
------------------------------------------------------------
    10 | 282    | 269    | 298    | 238    | 284    | 249
   100 | 2804   | 2717   | 3045   | 766    | 2689   | 780
  1000 | 29597  | 27468  | 31275  | 2093   | 27151  | 2224
 10000 | 2.9e+5 | 2.7e+5 | 3.1e+5 | 12908  | 2.7e+5 | 13122
100000 | 3.1e+6 | 2.8e+6 | 3.1e+6 | 1.2e+5 | 2.8e+6 | 1.2e+5


--------------
d1000
---------------
     n | sum    | merges | imerge | binsum | sortsu | ration
------------------------------------------------------------
    10 | 214    | 190    | 203    | 318    | 206    | 343
   100 | 5879   | 2233   | 2578   | 3178   | 2338   | 3197
  1000 | 2.6e+5 | 24866  | 29015  | 21929  | 23344  | 21901
 10000 | 6.1e+6 | 2.7e+5 | 3e+5   | 48111  | 1.8e+5 | 48814
100000 | 6.6e+7 | 2.7e+6 | 3e+6   | 1.5e+5 | 1.8e+6 | 3.5e+5


--------------
bigint
---------------
     n | sum    | merges | imerge | binsum | sortsu | ration
------------------------------------------------------------
    10 | 1      | 14     | 19     | 19     | 19     | 36
   100 | 15     | 54     | 140    | 59     | 72     | 81
  1000 | 158    | 329    | 1294   | 461    | 499    | 550
 10000 | 1655   | 2957   | 13102  | 4477   | 4497   | 5242
100000 | 16539  | 35929  | 1.3e+5 | 44905  | 55101  | 52732


--------------
e
---------------
     n | sum    | merges | imerge | binsum | sortsu | ration
------------------------------------------------------------
    10 | 155    | 156    | 164    | 262    | 166    | 274
   100 | 9464   | 2777   | 3190   | 9282   | 2893   | 4079
  1000 | 1.1e+7 | 8.7e+5 | 9.4e+5 | 1.1e+7 | 8.8e+5 | 9e+5


--------------
d10
---------------
     n | sum    | merges | imerge | binsum | sortsu | ration
------------------------------------------------------------
    10 | 142    | 144    | 156    | 139    | 157    | 153
   100 | 1456   | 1470   | 1822   | 343    | 1541   | 347
  1000 | 15080  | 14978  | 16955  | 1278   | 15557  | 1241
 10000 | 1.5e+5 | 1.5e+5 | 1.6e+5 | 9845   | 1.5e+5 | 10343
100000 | 1.5e+6 | 1.5e+6 | 1.7e+6 | 96015  | 1.5e+6 | 99856


And here's the script:

# tmpsum.py

from __future__ import print_function

def mergesum(seq):
    while len(seq) > 1:
        new = [a + b for a, b in zip(seq[:-1:2], seq[1::2])]
        if len(seq) % 2:
            new.append(seq[-1])
        seq = new
    return seq[0]

def imergesum(iterable):
    stack = []
    it = iter(iterable)
    nextpow2 = 1
    while True:
        for n2, val in enumerate(it, 1):
            if n2 == nextpow2:
                stack.append(val)
                nextpow2 *= 2
                break
            elif n2 % 2:
                stack[-1] += val
            else:
                ishift = -1
                while not n2 % 2:
                    val, stack[ishift] = stack[ishift], val
                    ishift -= 1
                    n2 >>= 1
                stack[ishift] += val
        else:
            return mergesum(stack)

from collections import defaultdict

def binsum(iterable):
    bins = defaultdict(int)
    for num in iterable:
        bins[num.denominator] += num.numerator
    return mergesum([F(n, d) for d, n in bins.items()])

from operator import attrgetter

def sortsum(seq):
    seq = sorted(seq, key=attrgetter('denominator'))
    if len(seq) < 3:
        return sum(seq)
    return mergesum(seq)

from operator import itemgetter

def rationalsum(iterable, tablesize=1000):
    def binsumgen(iterable):
        iterator = iter(iterable)
        bins = defaultdict(int)
        finished = False
        denom = itemgetter(0)
        while not finished:
            finished = True
            for n, num in zip(range(tablesize), iterator):
                bins[num.denominator] += num.numerator
                finished = False
            if len(bins) >= tablesize or finished:
                dn = sorted(bins.items(), key=denom)
                yield mergesum([F(n, d) for d, n in dn])
                bins.clear()
    return imergesum(binsumgen(iterable))

sumfuncs = sum, mergesum, imergesum, binsum, sortsum, rationalsum

# Just a quick test
if True:
    print('testing', end=' ')
    from fractions import Fraction as F
    nums = [F(n, 2*n +1) for n in range(3000)]
    for n in 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 20, 30, 40, 50, 100, 2005:
        ns = nums[:n]
        result = sum(ns)
        for func in sumfuncs:
            assert func(ns) == result, func.__name__
        print('.', end='')
    print(' passed!')

if True:
    print('generating data...')
    from random import randint, gauss, seed
    from fractions import Fraction as F
    seed(123456789)
    nmax = 10 ** 5
    nums_d1000 = [F(randint(-1000, 1000), randint(1, 1000)) for _ in
range(nmax)]
    nums_d10 = [F(randint(-10, 10), randint(1, 10)) for _ in range(nmax)]
    nums_float = [F(*gauss(0, 1).as_integer_ratio()) for _ in range(nmax)]
    nums_e = [F(1)]
    for n in range(1, 1000):
        nums_e.append(nums_e[-1] / n)
    nums_bigint = [10**65 + randint(1, 100) for n in range(nmax)]
    nums_all = {
        'd1000': nums_d1000,
        'd10': nums_d10,
        'float': nums_float,
        'e': nums_e,
        'bigint': nums_bigint,
    }

if True:
    print('\nRunning benchmarks (times in microseconds)...')

    import timeit

    def mytimeit(stmt, setup):
        n = 10
        time = lambda : timeit.timeit(stmt, setup, number=n) / n
        t = time()
        if t * n < 1e-1:
            while t * n < 1e-2:
                n *= 10
            if n == 10:
                ts = [t, time(), time()]
            else:
                ts = [time(), time(), time()]
            t = min(ts)
        return 1e6 * t

    def fmtnum(n):
        s = str(int(n))
        if len(s) > 5:
            s = '%1.2g' % n
            if s[-4:-1] == 'e+0':
                s = s[:-2] + s[-1]
        return '%-6s' % s

    gapline = '-' * 60
    fnames = [f.__name__ for f in sumfuncs]
    header = '     n | ' + ' | '.join('%-6s' % name[:6] for name in fnames)
    sum = sum
    setup = 'from __main__ import ' + ', '.join(['nums'] + fnames)
    setup += '; nums=nums[:%s]'

    for name, nums in nums_all.items():
        print('\n\n--------------\n%s\n---------------' % name)
        print(header)
        print(gapline)
        for n in 10, 100, 1000, 10000, 100000:
            if n > len(nums): break
            times = []
            for func in sumfuncs:
                stmt = '%s(nums)' % (func.__name__,)
                times.append(mytimeit(stmt, setup % n))
            print(('%6s | ' % n), end='')
            print(' | '.join(fmtnum(t) for t in times))


Oscar


More information about the Python-ideas mailing list