sympy.MatrixSymbol : 微分问题答案

【问题标题】：sympy.MatrixSymbol : Problems with differentiationsympy.MatrixSymbol : 微分问题
【发布时间】：2017-06-05 21:22:41
【问题描述】：

我正在为 sympy.MatrixSymbol 的元素似乎与 sympy 的微分例程无法很好地交互这一事实而苦苦挣扎。

我正在尝试使用 sympy.MatrixSymbol 元素而不是“普通” sympy 符号的事实是因为我想自动包装一个大函数，这似乎是克服参数限制并启用的唯一方法单个数组的输入。

为了让读者了解可能的解决方案的限制，我将从概述我的意图开始；但是，仓促的读者不妨跳到下面的代码块，这说明了我的问题。

声明某种变量的向量或数组。
利用前者的元素构建一些表达式；这些表达式构成所述向量的向量值函数的分量。除了这个功能，我想获得雅可比 w.r.t.向量。
使用自动换行（带有 cython 后端）来获得向量函数及其雅可比行列式的数值实现。这对前面的步骤施加了一些限制：(a) 希望函数的输入以向量的形式给出，而不是符号列表。（这既是因为 autwrapped 函数的输入数量似乎是有限的，也是为了便于以后与 scipy 交互，即避免经常将 numpy 向量解包到列表中）。

在旅途中，我遇到了 2 个问题：

Cython 似乎不喜欢某些 sympy 函数，其中我非常依赖 sympy.Max。 autowrap 的“助手”kwarg 似乎无法同时处理多个助手。
这本身没什么大不了的，因为我学会了使用 abs() 或 sign() 来规避它，cython 很容易理解。

（另见上文this question）

如前所述，autowrap/cython 不接受超过 509 个符号形式的参数，至少在我的编译器设置中不接受。（另见here）
我的最新问题从这里开始：我意识到 sympy.MatrixSymbol 的元素在许多方面不像“其他” sympy 符号。必须单独分配属性 real 和 commutative ，但这似乎可以正常工作。然而，当我试图得到雅可比行列时，我真正的麻烦开始了。 sympy 似乎无法立即获得元素的衍生物：

import sympy

X= sympy.MatrixSymbol("X",10,1)

for element in X:
    element._assumptions.update({"real":True, "commutative":True})

X[0].diff(X[0])
Out[2]: Derivative(X[0, 0], X[0, 0])


X[1].diff(X[0])
Out[15]: Derivative(X[1, 0], X[0, 0])

以下块是我想做的一个最小示例，但这里使用普通符号：（我认为它包含了我需要的所有内容，如果我忘记了什么我会在稍后添加。）

import sympy
from sympy.utilities.autowrap import autowrap

X = sympy.symbols("X:2", real = True)
expr0 = X[1]*( (X[0] - abs(X[0]) ) /2)**2
expr1 = X[0]*( (X[1] - abs(X[1]) ) /2)**2
F = sympy.Matrix([expr0, expr1])
J = F.jacobian([X[0],X[1]])
J_num = autowrap(J, args = [X[0],X[1]], backend="cython")

这是我（目前）使用 sympy.MatrixSymbol 的最佳猜测，然后当然会失败，因为 J 中的 Derivative-表达式：

X= sympy.MatrixSymbol("X",2,1)
for element in X:
    element._assumptions.update({"real":True, "commutative":True, "complex":False})
expr0 = X[1]*( (X[0] - abs(X[0]) ) /2)**2
expr1 = X[0]*( (X[1] - abs(X[1]) ) /2)**2
F = sympy.Matrix([expr0, expr1])
J = F.jacobian([X[0],X[1]])
J_num = autowrap(J, args = [X], backend="cython")

这是J运行上述代码后的样子：

J
Out[50]: 
Matrix([
[(1 - Derivative(X[0, 0], X[0, 0])*X[0, 0]/Abs(X[0, 0]))*(-Abs(X[0, 0])/2 + X[0, 0]/2)*X[1, 0], (-Abs(X[0, 0])/2 + X[0, 0]/2)**2],
[(-Abs(X[1, 0])/2 + X[1, 0]/2)**2, (1 - Derivative(X[1, 0], X[1, 0])*X[1, 0]/Abs(X[1, 0]))*(-Abs(X[1, 0])/2 + X[1, 0]/2)*X[0, 0]]])

不出所料，自动换行不喜欢：

[...]
wrapped_code_2.c(4): warning C4013: 'Derivative' undefined; assuming extern returning int
[...]
wrapped_code_2.obj : error LNK2001: unresolved external symbol Derivative

我如何告诉 sympy X[0].diff(X[0])=1 和 X[0].diff(X[1])=0？甚至可能是abs(X[0]).diff(X[0]) = sign(X[0])。

或者有什么办法可以使用 sympy.MatrixSymbol 和 still get a cythonized function，其中输入是单个向量而不是符号列表？

对于任何输入都非常有用，很可能是上述过程的任何步骤的解决方法。感谢阅读！

编辑：
一个简短的评论：我自己想出的一个解决方案是：使用普通符号构造F 和J；然后用一些 sympy.MatrixSymbol 的元素替换两个表达式中的符号。这似乎完成了工作，但更换需要相当长的时间，因为J 可以达到~1000x1000 及以上的尺寸。因此，我宁愿避免这种方法。

【问题讨论】：

标签： python cython sympy

【解决方案1】：

经过更广泛的研究，我上面描述的问题似乎已经在 development/github 版本中得到修复。相应更新后，所有涉及MatrixElement的Derivative术语均已正确解析！

请参阅here 以供参考。

【讨论】：