<p><strong>您可能对使用<a href="https://docs.scipy.org/doc/scipy/reference/linalg.cython_lapack.html" rel="nofollow noreferrer">scipy.linalg.cython_lapack</a>感兴趣。</strong>
它提供对LAPACK函数<a href="http://www.netlib.org/lapack/explore-html/df/da4/dgetri_8f_source.html" rel="nofollow noreferrer">^{<cd1>}</a>的访问。而<a href="https://groups.google.com/forum/#!topic/cython-users/aQMy_4APZlY" rel="nofollow noreferrer">the good news</a>是:</p>
<blockquote>
<p>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. </p>
</blockquote>
<p>使用<code>dger</code>的示例在<a href="https://stackoverflow.com/questions/44710838/calling-blas-lapack-directly-using-the-scipy-interface-and-cython">Calling BLAS / LAPACK directly using the SciPy interface and Cython</a>上提供。另请参见<a href="https://stackoverflow.com/questions/31994879/improving-cython-lapack-performance-with-internal-array-definitions">Improving Cython Lapack performance with internal array definitions?</a>
在我对<a href="https://stackoverflow.com/questions/40363111/mpi-python-open-mpi/40410216#40410216">MPI python-Open-MPI</a>的回答中,我详细介绍了如何使用cythonãblas,因此下面是如何将其改编为dgetri的方法:</p>
<ol>
<li><p>代码的关键部分写在<a href="http://cython.org/" rel="nofollow noreferrer">^{<cd3>}</a>,在一个专用文件<code>myinverse.pyx</code>中。</p></li>
<li><p>这个文件被Cython转换成<code>myinverse.c</code>文件</p></li>
<li><p>这个c文件是由您最喜欢的c编译器<code>gcc</code>编译来构建一个共享库<code>myinverse.so</code></p></li>
<li><p>优化后的函数可以在程序中使用<code>import myinverse</code>。</p></li>
</ol>
<p>下面是一个cython模块,它将放在.pyx文件中:</p>
<pre><code>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,""
</code></pre>
<p>要对.pyx文件进行cythonize和编译,可以使用以下makefile(我希望您使用的是Linux…)</p>
^{pr2}$
<p>新的python myinverse函数chainng LAPACK的<code>dgetrf()</code>和{<cd10>},在python主文件中调用:</p>
^{3}$