可以设置 "ufuncify" 函数来接受复杂的输入吗？ "ufunc 'wrapper_module_0' not supported for the input types" 错误
Can a "ufuncify" function be set up to take complex input? "ufunc 'wrapper_module_0' not supported for the input types" error
我正在学习如何使用 sympy 进行符号变量操作。我现在有一个表达式，我想在优化方案中对其进行多次评估，因此我希望它很快达到 运行。 sympy documentation on Numeric Computation 描述了一些方法：subs/evalf；羔羊化； lambdify-numpy；函数化； Theano.
到目前为止，我已经让 lambdify 开始工作，但它对我来说仍然不够快。进一步阅读，看起来 lambdify 以 python 速度运行，而 ufuncify 可能会在后台将代码重新编译为 C 代码，所以我现在正在研究这个选项。
到目前为止，我已经能够对我的表达式进行 ufuncify，但是每当我将复杂的输入作为参数传递时它都会抛出错误。
我修改了 this link 中的 ufuncify 示例以创建我的 MWE。当我 运行 这个：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# Create an example expression
a, b, c = sp.symbols('a, b, c')
expr = a + b*c
# Create a binary (compiled) function that broadcasts it's arguments
func = ufuncify((a, b, c), expr)
b = func(np.arange(5), 2.0+1j, 3.0)
print(b)
我收到此错误：TypeError：输入类型不支持 ufunc 'wrapper_module_0'，并且根据转换规则“'safe'”无法将输入安全地强制转换为任何受支持的类型
如果我更改代码以仅删除第二个参数的虚部，运行没问题：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# Create an example expression
a, b, c = sp.symbols('a, b, c')
expr = a + b*c
# Create a binary (compiled) function that broadcasts it's arguments
func = ufuncify((a, b, c), expr)
b = func(np.arange(5), 2.0, 3.0)
print(b)
returns: [6. 7. 8. 9. 10.]
理想情况下，我想使用 cython 或 f2py 作为后端，因为我希望它们是最快的...但是我遇到了类似的错误：
func = ufuncify((a, b, c), expr, backend='cython')
returns 类型错误：参数 'b_5226316' 的类型不正确（预期 numpy.ndarray，变得复杂）
可能还有另一个答案，但我通常会避免在 sympy 中封装复数值，因此您可以在一对 ar 和 ai 中设置每个值（对于实数和虚数部分）。
您还可以定义小表达式 (a = ar + sp.I*ai) 以简化代码编写。还要预先选择要允许复数的位置。
代码将是：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# Create an example expression
ar, ai, br, bi, cr, ci = sp.symbols('ar, ai, br, bi, cr, ci', real=True)
a = ar + sp.I*ai
b = br + sp.I*bi
c = cr + sp.I*ci
expr = a + b*c
# Create a binary (compiled) function that broadcasts it's arguments
func = ufuncify((ar, ai, br, bi, cr, ci), expr)
b = func(np.arange(5), 0, 2.0, 1, 3.0, 0)
print(b)
我没有编译它，因为我没有安装这些东西，但它至少开始编译了，所以我很确定它应该可以工作。
Update 没有编译是个错误。 ufuncify 函数目前不允许复杂的输入（它在 CodeGenerator 中也被称为里程碑）。 ufuncify 使用 Codegen 模块从给定表达式创建 C-code。在上面的例子中，它看起来像：
double autofunc0(double ar, double ai, double br, double bi, double cr, double ci) {
double autofunc0_result;
autofunc0_result = I*ai + ar + (I*bi + br)*(I*ci + cr);
return autofunc0_result;
}
它目前不使用任何复杂的扩展，尽管这些扩展在 C99 中可用。我想这将是对 Sympy 代码生成器的一个很好的扩展。
在文档中，他们还建议可以编写自己的本机代码。
1) If you are really concerned about speed or memory optimizations,
you will probably get better results by working directly with the
wrapper tools and the low level code. However, the files generated
by this utility may provide a useful starting point and reference
code. Temporary files will be left intact if you supply the keyword
tempdir="path/to/files/".
我想这将是一个替代选项，即使用生成的文件作为起点，然后使用 header 来实现复杂的值。不过，这不是一个令人满意的答案。
In [251]: a, b, c = symbols('a, b, c')
...: expr = a + b*c
...:
In [252]: f = lambdify((a,b,c), expr)
In [254]: print(f.__doc__)
Created with lambdify. Signature:
func(a, b, c)
Expression:
a + b*c
Source code:
def _lambdifygenerated(a, b, c):
return (a + b*c)
Imported modules:
对于这个简单的expr，lambdify numpy 代码看起来不错——利用整个数组编译运算符：
In [255]: f(np.arange(5), 2.0, 3.0)
Out[255]: array([ 6., 7., 8., 9., 10.])
In [256]: timeit f(np.arange(5), 2.0, 3.0)
5.55 µs ± 16 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
In [257]: timeit np.arange(5) + 2.0 * 3.0
5.04 µs ± 16.7 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
如果我使用 numba 会更快:
In [258]: import numba
In [259]: @numba.njit
...: def foo(a,b,c):
...: return a+b*c
...:
In [260]: foo(np.arange(5),2.0,3.0)
Out[260]: array([ 6., 7., 8., 9., 10.])
In [261]: timeit foo(np.arange(5),2.0,3.0)
2.68 µs ± 69.1 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
它有点 hacky，它可能只适用于 f2py，但我认为这就是你想要做的：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# This seems to be False by default and prevents the f2py codegen from
# creating complex variables in fortran
sp.utilities.codegen.COMPLEX_ALLOWED = True
# Create an example expression, mark them as Complex so that the codegen
# tool doesn't assume they're just floats. May not even be necessary.
a, b, c = sp.symbols('a, b, c', complex=True)
expr = a + b*c
# Create a binary (compiled) function that broadcasts its arguments
func = ufuncify((a, b, c), expr, backend="f2py")
# The f2py backend will want identically-sized arrays
b = func(np.arange(5), np.ones(5) * (2.0 + 1.0j), np.ones(5) * 3)
print(b)
您可以确认生成的 Fortran 代码需要复杂变量
如果您在 ufuncify 上设置 tempdir="." 并打开 wrapped_code_0.f90 文件。
!******************************************************************************
!* Code generated with sympy 1.6.2 *
!* *
!* See http://www.sympy.org/ for more information. *
!* *
!* This file is part of 'autowrap' *
!******************************************************************************
subroutine autofunc(y_1038050, a_1038054, b_1038055, c_1038056, &
m_1038051)
implicit none
INTEGER*4, intent(in) :: m_1038051
COMPLEX*16, intent(out), dimension(1:m_1038051) :: y_1038050
COMPLEX*16, intent(in), dimension(1:m_1038051) :: a_1038054
COMPLEX*16, intent(in), dimension(1:m_1038051) :: b_1038055
COMPLEX*16, intent(in), dimension(1:m_1038051) :: c_1038056
INTEGER*4 :: i_1038052
do i_1038052 = 1, m_1038051
y_1038050(i_1038052) = a_1038054(i_1038052) + b_1038055(i_1038052)* &
c_1038056(i_1038052)
end do
end subroutine
我正在学习如何使用 sympy 进行符号变量操作。我现在有一个表达式，我想在优化方案中对其进行多次评估，因此我希望它很快达到 运行。 sympy documentation on Numeric Computation 描述了一些方法：subs/evalf；羔羊化； lambdify-numpy；函数化； Theano.
到目前为止，我已经让 lambdify 开始工作，但它对我来说仍然不够快。进一步阅读，看起来 lambdify 以 python 速度运行，而 ufuncify 可能会在后台将代码重新编译为 C 代码，所以我现在正在研究这个选项。
到目前为止，我已经能够对我的表达式进行 ufuncify，但是每当我将复杂的输入作为参数传递时它都会抛出错误。
我修改了 this link 中的 ufuncify 示例以创建我的 MWE。当我 运行 这个：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# Create an example expression
a, b, c = sp.symbols('a, b, c')
expr = a + b*c
# Create a binary (compiled) function that broadcasts it's arguments
func = ufuncify((a, b, c), expr)
b = func(np.arange(5), 2.0+1j, 3.0)
print(b)
我收到此错误：TypeError：输入类型不支持 ufunc 'wrapper_module_0'，并且根据转换规则“'safe'”无法将输入安全地强制转换为任何受支持的类型
如果我更改代码以仅删除第二个参数的虚部，运行没问题：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# Create an example expression
a, b, c = sp.symbols('a, b, c')
expr = a + b*c
# Create a binary (compiled) function that broadcasts it's arguments
func = ufuncify((a, b, c), expr)
b = func(np.arange(5), 2.0, 3.0)
print(b)
returns: [6. 7. 8. 9. 10.]
理想情况下，我想使用 cython 或 f2py 作为后端，因为我希望它们是最快的...但是我遇到了类似的错误：
func = ufuncify((a, b, c), expr, backend='cython')
returns 类型错误：参数 'b_5226316' 的类型不正确（预期 numpy.ndarray，变得复杂）
可能还有另一个答案，但我通常会避免在 sympy 中封装复数值，因此您可以在一对 ar 和 ai 中设置每个值（对于实数和虚数部分）。
您还可以定义小表达式 (a = ar + sp.I*ai) 以简化代码编写。还要预先选择要允许复数的位置。
代码将是：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# Create an example expression
ar, ai, br, bi, cr, ci = sp.symbols('ar, ai, br, bi, cr, ci', real=True)
a = ar + sp.I*ai
b = br + sp.I*bi
c = cr + sp.I*ci
expr = a + b*c
# Create a binary (compiled) function that broadcasts it's arguments
func = ufuncify((ar, ai, br, bi, cr, ci), expr)
b = func(np.arange(5), 0, 2.0, 1, 3.0, 0)
print(b)
我没有编译它，因为我没有安装这些东西，但它至少开始编译了，所以我很确定它应该可以工作。
Update 没有编译是个错误。 ufuncify 函数目前不允许复杂的输入（它在 CodeGenerator 中也被称为里程碑）。 ufuncify 使用 Codegen 模块从给定表达式创建 C-code。在上面的例子中，它看起来像：
double autofunc0(double ar, double ai, double br, double bi, double cr, double ci) {
double autofunc0_result;
autofunc0_result = I*ai + ar + (I*bi + br)*(I*ci + cr);
return autofunc0_result;
}
它目前不使用任何复杂的扩展，尽管这些扩展在 C99 中可用。我想这将是对 Sympy 代码生成器的一个很好的扩展。 在文档中，他们还建议可以编写自己的本机代码。
1) If you are really concerned about speed or memory optimizations, you will probably get better results by working directly with the wrapper tools and the low level code. However, the files generated by this utility may provide a useful starting point and reference code. Temporary files will be left intact if you supply the keyword tempdir="path/to/files/".
我想这将是一个替代选项，即使用生成的文件作为起点，然后使用
In [251]: a, b, c = symbols('a, b, c')
...: expr = a + b*c
...:
In [252]: f = lambdify((a,b,c), expr)
In [254]: print(f.__doc__)
Created with lambdify. Signature:
func(a, b, c)
Expression:
a + b*c
Source code:
def _lambdifygenerated(a, b, c):
return (a + b*c)
Imported modules:
对于这个简单的expr，lambdify numpy 代码看起来不错——利用整个数组编译运算符：
In [255]: f(np.arange(5), 2.0, 3.0)
Out[255]: array([ 6., 7., 8., 9., 10.])
In [256]: timeit f(np.arange(5), 2.0, 3.0)
5.55 µs ± 16 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
In [257]: timeit np.arange(5) + 2.0 * 3.0
5.04 µs ± 16.7 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
如果我使用 numba 会更快:
In [258]: import numba
In [259]: @numba.njit
...: def foo(a,b,c):
...: return a+b*c
...:
In [260]: foo(np.arange(5),2.0,3.0)
Out[260]: array([ 6., 7., 8., 9., 10.])
In [261]: timeit foo(np.arange(5),2.0,3.0)
2.68 µs ± 69.1 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
它有点 hacky，它可能只适用于 f2py，但我认为这就是你想要做的：
import sympy as sp
from sympy.utilities.autowrap import ufuncify
import numpy as np
# This seems to be False by default and prevents the f2py codegen from
# creating complex variables in fortran
sp.utilities.codegen.COMPLEX_ALLOWED = True
# Create an example expression, mark them as Complex so that the codegen
# tool doesn't assume they're just floats. May not even be necessary.
a, b, c = sp.symbols('a, b, c', complex=True)
expr = a + b*c
# Create a binary (compiled) function that broadcasts its arguments
func = ufuncify((a, b, c), expr, backend="f2py")
# The f2py backend will want identically-sized arrays
b = func(np.arange(5), np.ones(5) * (2.0 + 1.0j), np.ones(5) * 3)
print(b)
您可以确认生成的 Fortran 代码需要复杂变量
如果您在 ufuncify 上设置 tempdir="." 并打开 wrapped_code_0.f90 文件。
!******************************************************************************
!* Code generated with sympy 1.6.2 *
!* *
!* See http://www.sympy.org/ for more information. *
!* *
!* This file is part of 'autowrap' *
!******************************************************************************
subroutine autofunc(y_1038050, a_1038054, b_1038055, c_1038056, &
m_1038051)
implicit none
INTEGER*4, intent(in) :: m_1038051
COMPLEX*16, intent(out), dimension(1:m_1038051) :: y_1038050
COMPLEX*16, intent(in), dimension(1:m_1038051) :: a_1038054
COMPLEX*16, intent(in), dimension(1:m_1038051) :: b_1038055
COMPLEX*16, intent(in), dimension(1:m_1038051) :: c_1038056
INTEGER*4 :: i_1038052
do i_1038052 = 1, m_1038051
y_1038050(i_1038052) = a_1038054(i_1038052) + b_1038055(i_1038052)* &
c_1038056(i_1038052)
end do
end subroutine