Calling scipy's nonnegative least squares from nopython block

Hey guys,

I have been a huge fan of numba for some time now. I am mostly using it to speed up the computation of cost functions that are optimized by scipy then.

Now my cost function itself calls scipy’s nonnegative least squares (https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.nnls.html) itself in each iteration. I was wondering if it was possible to circumvent the python interpreter here and call the nnls function directly from numba. Unfortunately, it is wrapped by f2py, so I am not sure if it is as painfree as with C functions. Could you give me a hint here?

Best

Hi @dschmitz89,

Unfortunately I don’t think this is possible from the nnls bindings that SciPy provides. Whilst the routine is Fortran, it’s wrapped in a C extension and the symbol to the Fortran routine is local:

$ nm `python -c 'import scipy.optimize.__nnls as x; print(x.__file__)'`|grep 't nnls_\|T nnls_'
00000000000057b0 t nnls_

as a result, there’s no reasonable way of calling that function directly.

Were SciPy to elect to make this symbol global or make it e.g. a Cython export then Numba could bind to it and this would be just a case of wiring data.

Assuming SciPy doesn’t make the symbol global, I think your options are:

  1. Compile and ship the Fortran and bind to it, or find another package which ships this routine in a library with the symbol exported and bind to that.
  2. Translate the Fortran to Python and let Numba compile it. Given the routine makes use of things that are available in BLAS/LAPACK (Givens rotations and Householder transforms) and Numba can bind to BLAS/LAPACK, this might be the way to go.

Hope this helps?


References:
NNLS C extension import: https://github.com/scipy/scipy/blob/v1.5.4/scipy/optimize/_nnls.py#L1
NNLS Fortran routine: https://github.com/scipy/scipy/blob/v1.5.4/scipy/optimize/__nnls/nnls.f#L1

Hi @stuartarchibald,

thank you for the extensive answer. I got the same feeling that it would be hard to directly bind to the scipy version. Fortunately, i found a f2c translation of the code online: https://tinyurl.com/y5v9qrac .

After compiling this to a shared library, linking numba to a ctypes wrapper of the function should be doable with help of the documentation. Would complicate the build process of the Python application but is unavoidable. Speedups for small NNLs problems should be definitely possible as scipy does some sanity checks on the inputs which can be skipped.

Would this be of general interest for numba?

Glad you have found a work-around.

https://github.com/numba/numba-scipy is the place for SciPy functionality in Numba, however it is a pure python library so unfortunately a piece of C or Fortran code would unlikely be accepted. Were the NNLs routine translated into Numba JIT compiled Python, I think that could be added.

Just so to mark this question as solved: I ended up wrapping a slightly newer version of the Fortran code via f2py myself. In my wrapper there is no overhead caused by sanity checks anymore, so it is faster than scipy’s version.

Going the Fortran->C->Cython->Numba road turned out very hard for someone like me without proper experience in shipping compiled code. Reimplementing the whole thing also seems too much work at the moment.

Thanks to Stuart again for the answers though!

@dschmitz89 thanks for explaining what you did, glad you found something that worked :slight_smile: