Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How can you implement a C callable from Numba for efficient integration with nquad?

I need to do a numerical integration in 6D in python. Because the scipy.integrate.nquad function is slow I am currently trying to speed things up by defining the integrand as a scipy.LowLevelCallable with Numba.

I was able to do this in 1D with the scipy.integrate.quad by replicating the example given here:

import numpy as np
from numba import cfunc
from scipy import integrate

def integrand(t):
    return np.exp(-t) / t**2

nb_integrand = cfunc("float64(float64)")(integrand)

# regular integration
%timeit integrate.quad(integrand, 1, np.inf)

10000 loops, best of 3: 128 µs per loop

# integration with compiled function
%timeit integrate.quad(nb_integrand.ctypes, 1, np.inf)

100000 loops, best of 3: 7.08 µs per loop

When I want to do this now with nquad, the nquad documentation says:

If the user desires improved integration performance, then f may be a scipy.LowLevelCallable with one of the signatures:

double func(int n, double *xx)
double func(int n, double *xx, void *user_data)

where n is the number of extra parameters and args is an array of doubles of the additional parameters, the xx array contains the coordinates. The user_data is the data contained in the scipy.LowLevelCallable.

But the following code gives me an error:

from numba import cfunc
import ctypes

def func(n_arg,x):
    xe = x[0]
    xh = x[1]
    return np.sin(2*np.pi*xe)*np.sin(2*np.pi*xh)

nb_func = cfunc("float64(int64,CPointer(float64))")(func)

integrate.nquad(nb_func.ctypes, [[0,1],[0,1]], full_output=True)

error: quad: first argument is a ctypes function pointer with incorrect signature

Is it possible to compile a function with numba that can be used with nquad directly in the code and without defining the function in an external file?

Thank you very much in advance!

like image 696
dropther Avatar asked Aug 22 '17 16:08

dropther


People also ask

How does Numba work in Python?

Numba reads the Python bytecode for a decorated function and combines this with information about the types of the input arguments to the function. It analyzes and optimizes your code, and finally uses the LLVM compiler library to generate a machine code version of your function, tailored to your CPU capabilities.

What is Numba good for?

Built for Scientific Computing Numba is designed to be used with NumPy arrays and functions. Numba generates specialized code for different array data types and layouts to optimize performance. Special decorators can create universal functions that broadcast over NumPy arrays just like NumPy functions do.

Can I use Numba with NumPy?

NumPy arrays are directly supported in Numba. Access to Numpy arrays is very efficient, as indexing is lowered to direct memory accesses when possible. Numba is able to generate ufuncs and gufuncs.

Can I use scipy in Numba?

However numba and scipy are still not compatible. Yes, Scipy calls compiled C and Fortran, but it does so in a way that numba can't deal with.


1 Answers

Wrapping the function in a scipy.LowLevelCallable makes nquad happy:

si.nquad(sp.LowLevelCallable(nb_func.ctypes), [[0,1],[0,1]], full_output=True)
# (-2.3958561404687756e-19, 7.002641250699693e-15, {'neval': 1323})
like image 74
MB-F Avatar answered Sep 25 '22 21:09

MB-F