Search code examples
pythonnumpypytablesnumexpr

Capping a sub-expression in numexpr


How do I efficiently express the following using numexpr?

z = min(x-y, 1.0) / (x+y)

Here, x and y are some large NumPy arrays of the same shape.

In other words, I am trying to cap x-y to 1.0 before dividing it by x+y.

I would like to do this using a single numexpr expression (x and y are huge, and I don't want to have to iterate over them more than once).


Solution

  • Maybe something like this would work?

    In [11]: import numpy as np
    In [12]: import numexpr as ne    
    In [13]:     
    In [13]: x = np.linspace(0.02, 5.0, 1e7)
    In [14]: y = np.sin(x)
    In [15]:     
    In [15]: timeit z0 = ((x-y) - ((x-y) > 1) * (x-y - 1))/(x+y)
    1 loops, best of 3: 1.02 s per loop
    In [16]: timeit z1 = ne.evaluate("((x-y) - ((x-y) > 1.) * ((x-y) - 1.))/(x+y)")
    10 loops, best of 3: 120 ms per loop    
    In [17]: timeit z2 = ne.evaluate("((x-y)/(x+y))")
    10 loops, best of 3: 103 ms per loop
    

    There's a penalty for the capping above the division, but it's not too bad. Unfortunately when I tried it for some larger arrays it segfaulted. :-/

    Update: this is much prettier, and a little faster too:

    In [40]: timeit w0 = ne.evaluate("where(x-y>1,1,x-y)/(x+y)")
    10 loops, best of 3: 114 ms per loop