How to compile C extension for Python where C function uses LAPACK library?

Question

I wrote a C extesion for Python and the module is compiled in a .so file successfully. However when I am trying to use the wrapped C function in Python side (a test code in python that calls the wrapped C function) I get the following ImportError

ImportError: /home/username/newModule.cpython-36m-x86_64-linux-gnu.so: undefined symbol: dgetri_

I am pretty sure that undefined symbol: dgetri_ in the import error is because the generated .so file did not find link to LAPACK library. So my question is as following,

How do I compile c extension code for python when the wrapped C function depends on LAPACK library to generate module in .so format?

Currently I am compiling the C code using python's utils.core module. I think I need to compile the C code from command line to link LAPACK but do not exactly know which are the appropriate commands to use?

Any help is appreciated.

Overcompensation

Answer 1

You may be interested in using scipy.linalg.cython_lapack. It provides access to LAPACK's function dgetri among others. And the good news is:

This makes it possible to use SciPy's BLAS and LAPACK from any 3rd party Cython module without explicitely linking with the libraries. This means that projects like scikit-learn and statsmodels do not need to maintain a separate build dependency on BLAS and LAPACK.

An exemple using dger is available at Calling BLAS / LAPACK directly using the SciPy interface and Cython . See also Improving Cython Lapack performance with internal array definitions? I detailed how to use cython_blas in my answer to MPI python-Open-MPI , so here is how it can be adapted to dgetri:

1. The critical part of the code are written in Cython, in a dedicated file myinverse.pyx.

2. This file is turned into a myinverse.c file by Cython

3. This c file is compiled by your favorite c compiler gcc to build a shared library myinverse.so

4. The optimized function can be used in your program after import myinverse.

Here is a cython module, to be placed in the .pyx file:

import numpy

cimport numpy
cimport scipy.linalg.cython_lapack
ctypedef numpy.float64_t DTYPE_t
cimport cython
from libc.stdlib cimport malloc, free

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
def invert(numpy.ndarray[DTYPE_t, ndim=2] array):

    cdef  int rows = array.shape[0]
    cdef   int cols = array.shape[1]
    cdef  int info = 0
    if cols !=rows:
        return array,1,"not a square matrix"

    cdef int* ipiv = <int *> malloc(rows * sizeof(int))
    if not ipiv:
        raise MemoryError()

    scipy.linalg.cython_lapack.dgetrf(&cols,&rows,&array[0,0],&rows,ipiv,&info)
    if info !=0:
        free(ipiv)
        return array,info,"dgetrf failed, INFO="+str(info)
    #workspace query
    cdef double workl
    cdef int lwork=-1
    scipy.linalg.cython_lapack.dgetri(&cols,&array[0,0],&rows,ipiv,&workl,&lwork,&info)
    if info !=0:
        free(ipiv)
        return array,info,"dgetri failed, workspace query, INFO="+str(info)
    #allocation workspace
    lwork= int(workl)
    cdef double* work = <double *> malloc(lwork * sizeof(double))
    if not work:
        raise MemoryError()

    scipy.linalg.cython_lapack.dgetri(&cols,&array[0,0],&rows,ipiv,work,&lwork,&info)
    if info !=0:
        free(ipiv)
        free(work)
        return array,info,"dgetri failed, INFO="+str(info)

    free(ipiv)
    free(work)

    return array,info,""

To cythonize and compile the .pyx file, the following makefile can be used (I hope you are using Linux...)

all: myinverse myinverseb


myinverse: myinverse.pyx
    cython -a myinverse.pyx

myinverseb: myinverse.c
    gcc -shared -pthread -fPIC -fwrapv -O2 -Wall -fno-strict-aliasing -I/usr/include/python2.7 -o myinverse.so myinverse.c

The new python myinverse function, chainng LAPACK's dgetrf() and dgetri(), is called in the main python file:

import numpy as np

import myinverse
n=42

#A=np.zeros((n,n))
#for i in range(n):
#    A[i,i]=10
A=np.random.rand(n,n)
#A=np.zeros((n,n))
Am,info,string=myinverse.invert(A.copy())
if info==0:
    print np.linalg.norm(A.dot(Am)-np.identity(n), np.inf)
else :
    print "inversion failed, info=",info, string

Lafontaine

as additional advantage: scipy has some compatibility code for different versions of BLAS/LAPACK libraries that are installed, e.g. whether scipy was build with MKL, OpenBLAS or others. This means when using the cython wrappers from scipy, the calling code doesn't need to worry about which linalg library is used. – Mccloud

Thank you for all the details. They are interesting reads. However I got my code worked by using the following two lines while compiling in linux – Overcompensation

Answer 2

Alternatively I got my code worked by using the following two lines while compiling in linux

$gcc -DNDEBUG -Wall -Wstrict-prototypes -fPIC -I/home/username/anaconda3/include/python3.6m -c stackDoc.cpp -o mydemo.o
$gcc -shared mydemo.o -o mydemo.so

The following link I found useful, https://docs.python.org/2/extending/building.html

Overcompensation

Answer 3

Alternatively you can make the module to be linked against lapack.

e.g.

from distutils.core import setup
from distutils.extension import Extension

setup(
    name='MyExtension',
    version='0.1',
    ext_modules=[
        Extension('lib_name', ['lib_name.cpp'], extra_link_args=['-lopenblas']),
    ],
    scripts=['lib_name.cpp', '__init__.py'],
)

Infield