1

I am writing some software in Python that would allow the user to run simulations that involve simple additions and multiplications and interpolation on an arbitrary number of dimensions. In the process, the function has to repeatedly (and sequentially) add, multiply, and interpolate on an N-dimensional grid. I have to do this in a loop since in each iteration the results are dependent on the previous iteration (they're not parallelizable). In the 1 dimensional case, it looks like

T = 1000000000
b = 0
for t in range(T):
    b = b * f(b) + np.random.randn()
    b = b + f(b)

where f in the code above is an interpolant created by the data and grid supplied by the user. In each iteration, the same interpolant is used.

The first thing that comes to mind is to use numba to speed up the loops. However, to be able to interpolate on an N dimensional grid, I have to use scipy.interpolate.RegularGridInterpolator but using this would prevent me from running numba in nopython mode. I know that there is this package that allows faster interpolations, but I'm not sure if I can use this in the nopython mode. Also, I'm wondering if there are other options to get over this bottleneck? Any help is greatly appreciated.

hpaulj
  • 221,503
  • 14
  • 230
  • 353
tryingtosolve
  • 763
  • 5
  • 20
  • You need to name or tag that specialized numba package. Otherwise people who know about it won't see your question. It's too specialized for a generic numpy user. – hpaulj Apr 14 '18 at 02:34
  • I think you can read the source code of RegularGridInterpolator, and implement it in numba. – HYRY Apr 14 '18 at 03:00
  • If it is an option for you, you can wrap this specific piece of code in C++, you'd get a huge gain in time. Provided that you are willing to write a bit of C++, it is really easy with the use of PyBind, see [this answer](https://stackoverflow.com/questions/49582252/pybind-numpy-access-2d-nd-arrays/49693704#49693704) – Christian Apr 14 '18 at 09:35
  • @Christian Did you mean that I should write the function in C++, call ```scipy```'s interpolation function in C++ through Pybind, and write a Python wrapper for the entire function in C++? – tryingtosolve Apr 16 '18 at 14:15
  • @tryingtosolve No, sorry I wasn't very clear. When I said "provided that you are willing to write a bit of C++", that meant writing (or using a library) interpolation functions in C++. Just doing a loop in C++ but calling Python functions wouldn't increase your performances. Boost has a [cubic b-spline module](https://www.boost.org/doc/libs/1_65_0/libs/math/doc/html/math_toolkit/interpolate/cubic_b.html) or you can use [this C-code](http://bigwww.epfl.ch/thevenaz/interpolation/) written by a research lab to do spline interpolation. The latter may be easier to include. – Christian Apr 16 '18 at 14:27
  • @Christian do you have examples of boost in which multi-dimensional interpolations are implemented? – tryingtosolve Apr 16 '18 at 14:44
  • 1
    I never did it myself, no, maybe have a look at [this question](https://stackoverflow.com/questions/38306335/c-boost-store-objects-in-multidimensional-array-without-default-constructor) – Christian Apr 16 '18 at 14:49

0 Answers0