Sympy autowrap (cython):定义“助手"对于 sympy.Max, sympy.Heaviside [英] Sympy autowrap (cython): defining "helpers" for sympy.Max, sympy.Heaviside

问题描述

我有一个 sympy.Matrix，称为 J_sym，我想自动换行(最好使用 cython 后端)；相应的符号存储在列表 list_args 中.

I have a sympy.Matrix, called J_sym, that I would like to autowrap (preferably using the cython backend); the according symbols are stored in the list list_args.

然而，我遇到的问题是显然不支持某些 sympy 函数，在我的例子中特别是 sympy.Max 和 sympy.Heaviside.

However, the problem I run into is that apparently some sympy functions are not supported, in my case specifically sympy.Max and sympy.Heaviside.

具体来说，如果我尝试

J_num_autowrap = autowrap(J_sym, backend="cython", args=list_args)

我得到以下信息:

[...]
wrapped_code_10.c(4): warning C4013: 'Heaviside' undefined; assuming extern returning int
wrapped_code_10.c(4): warning C4013: 'Max' undefined; assuming extern returning int
[...]
: fatal error LNK1120: 2 unresolved externals
[...]
failed with exit status 1120

(当仅使用 numpy 进行lambdify 时会出现类似的问题；但可以轻松地为这些函数传递字典.)

(a similar problem arises when one lambdifies using only numpy; but there one can easily pass a dictionary for these functions.)

现在，似乎 autowrap 的helpers"参数应该提供一个解决方案；但是，我不知道如何正确使用它.文档字符串本身绝对是神秘的，我只找到了这个带有示例的链接:
https://github.com/sympy/sympy/issues/10572

Now, it seems that the "helpers" argument to autowrap should offer a solution; however, I can't work out how to properly use it. The docstring itself is absolutely cryptic, and I found only this link with an example:
https://github.com/sympy/sympy/issues/10572

我不确定如何正确使用助手 - 谁能帮忙?

I'm unsure how to properly use helpers here - can anyone help?

如果我按照以下方式尝试:

If I try something along these lines:

J_num_autowrap = autowrap(J_sym, backend="cython", args=list_args, helpers=("Max", max(x,y), [x,y]))

如果 x, y 未声明，我会收到 (a) 错误 - 我如何使用自由变量执行此操作?(b) 如果我在将 x 和 y 声明为 sympy 符号后做同样的事情，我只会得到错误关系的真值无法确定".

I get (a) an error if x, y are not declared - how can I do this with free variables? And (b) if I do the same after declaring x and y as sympy symbols, I just get the error "Truth value of relational can not be determined".

(上面的链接提到无论如何都不可能传递一个以上的帮助程序 - 尽管在我看来这将在 1.0 中修复.我仍然收到ValueError:没有足够的值来解包(预期为 3, got 2)" 尝试传递 2 个(胡言乱语)辅助元组时.不确定状态 - 也许我在这里运气不好).

(The above link mentions that it is impossible to pass more than one helper anyway - even though it seems to me that this was to be fixed in 1.0. I still get the "ValueError: not enough values to unpack (expected 3, got 2)" when trying to pass 2 (gibberish) helper tuples. Not sure about the status - maybe I'm just out of luck here).

如果有任何其他方法可以从我的矩阵中生成 C 代码，我也会很感激在这个方向上的任何提示.

If there is any alternative way to generate C code from my matrix, I would appreciate any hints in that direction as well.

推荐答案

作为 Max 问题的解决方法，您可以使用

As a workaround for Max issue, you can use

helpers=("Max", (abs(x+y) + abs(x-y))/2, [x, y])

表达式 (abs(x+y) + abs(xy))/2 在数学上等价于 max(x, y) 并且不会产生任何抱怨来自 autowrap.一个完整的例子:

The expression (abs(x+y) + abs(x-y))/2 is mathematically equivalent to max(x, y) and does not generate any complaints from autowrap. A complete example:

from sympy import *
from sympy.utilities.autowrap import autowrap
x, y = symbols('x y')
f = autowrap(Max(2*x, y+1), args=[x, y], backend="cython", helpers=("Max", (abs(x+y) + abs(x-y))/2, [x, y]))
print([f(5, 6), f(1, 2)])    # outputs 10 and 3 

类似地，Heaviside(x) 与 (x + abs(x))/(2*x) 相同——除了后者表达式为 NaN 时x=0.更丑但安全的版本是 (x + abs(x))/(2*abs(x) + 1e-300)，其中添加的 1e-300 几乎不会改变结果.

Similarly, Heaviside(x) is the same as (x + abs(x))/(2*x) -- except the latter expression is NaN when x=0. The uglier but safe version is (x + abs(x))/(2*abs(x) + 1e-300) where the added 1e-300 will almost never change the result.

您需要 Max 和 Heaviside 的助手.这就是您遇到未决问题的地方:autowrap 中存在一个错误，导致不可能使用多个助手.文档表明该格式将就像

You want helpers for both Max and Heaviside. And this is where you hit an open issue: there is a bug in autowrap which makes it impossible to use more than one helper. Documentation suggests that the format would be like

helpers=[("Max", ..., [x, y]), ("Heaviside", ..., [x])]

但是在这一行 autowrap做一些对一个助手很方便的事情(我们不必把它放在一个列表中)，但对一个以上的人来说却是致命的(额外的一层包装):

But on this line autowrap does something that is convenient for one helper (we don't have to put it in a list), but fatal for more than one (extra layer of wrapping):

helpers = [helpers] if helpers else ()

自然地，后续解包for name_h, expr_h, args_h in helpers 失败.

Naturally, subsequent unpacking for name_h, expr_h, args_h in helpers fails.

ufuncify 正确处理 helpers 参数.不幸的是，ufuncify 最终会为 NumPy 以外的所有后端调用 autowrap，因此错误仍然出现在 autowrap 中.但如果你愿意使用 NumPy 后端，这是一个解决方案:

The ufuncify handles helpers argument correctly. Unfortunately, ufuncify eventually calls autowrap for all backends other than NumPy, so the error still occurs in autowrap. But if you are willing to use NumPy backend, this is a solution:

from sympy.utilities.autowrap import ufuncify
x, y = symbols('x y', real=True)
my_helpers = [("Max", abs(x+y)/2 + abs(x-y)/2, [x, y]), ("Heaviside", (x + abs(x)) / (2*abs(x) + 1e-300), [x])]
f = ufuncify([x,y], Max(2*x, y+1) + Heaviside(x-4), backend="numpy", helpers=my_helpers)
print([f(5, 6), f(1, 2)])   $ outputs 11  and  3

减少到一个帮手

如果您可以更改表达式以消除 Max 和 Heaviside 之一(或什至两者)，则可以使用 Cython 后端.例如，也许您只需要 Max(x, 0)，因此您定义了一个 Python 函数正部分":

Reducing to one helper

If you can change the expression to eliminate one of Max and Heaviside (or even both), Cython backend can be used. For example, maybe you only need Max(x, 0) and so you define a Python function "positive part":

pos = lambda x: (x+abs(x))/2

那么autowrap就没有pos的问题了，你只需要帮助Heaviside:

Then the autowrap has no problem with pos, and you only need to help with Heaviside:

f = autowrap(pos(9-x-y) + Heaviside(x-4), args = [x, y], backend = "cython", helpers=("Heaviside", (x + abs(x)) / (2*abs(x) + 1e-300), [x]))
print([f(5, 6), f(1, 2)])   #  outputs 1 and 6

这种替换的实用性取决于您的符号表达式是如何获得的.

How practical such replacements are depends on how your symbolic expression is obtained.

这篇关于Sympy autowrap (cython):定义“助手"对于 sympy.Max, sympy.Heaviside的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！