In [1]:

    

import numpy as np
from pint.core import roundfloat as rf


    

In [2]:

    

def fasttwosum(a, b):
    x = a + b
    tmp = x - a
    y = b - tmp
    return x, y


    

In [3]:

    

def twosum(a, b):
    x = np.zeros(a.shape)
    y = np.zeros(a.shape)
    # aBb: a is bigger than b
    aBb = abs(a) > abs(b)
    ta = a[aBb]
    tb = b[aBb]
    tx, ty = fasttwosum(ta, tb)
    x[aBb] = tx
    y[aBb] = ty
    # 実質 abs(a) <= abs(b)
    bBa = ~aBb
    ta = a[bBa]
    tb = b[bBa]
    tx, ty = fasttwosum(tb, ta)
    x[bBa] = tx
    y[bBa] = ty
    return x, y


    

In [4]:

    

pntwosum = rf.twosum
def twosum_old(a, b):
    x = np.zeros(a.shape)
    y = np.zeros(a.shape)
    lenx = len(x)
    for j in range(0, lenx):
        for i in range(0, lenx):
            x[i][j], y[i][j] = pntwosum(a[i][j], b[i][j])
    return x, y


    

In [5]:

    

size = 1000
a = np.random.rand(size, size)
b = np.random.rand(size, size)


    

In [6]:

    

x, y = twosum(a, b)
old_x, old_y = twosum_old(a, b)


    

In [7]:

    

# もし、違う項があれば検出される
print(x[x != old_x])
print(y[y != old_y])


    

    




[]
[]

In [8]:

    

%timeit twosum(a, b)


    

    




10 loops, best of 3: 80.7 ms per loop

In [9]:

    

%timeit twosum_old(a, b)


    

    




1 loop, best of 3: 2 s per loop