我正在尝试使用 sympy 中的 ccode() 将来自 sage 的大表达式转换为有效的 C 代码。然而,我的表达式有很多平方和立方项。由于 pow(x,2) 比 x*x 慢得多,因此我尝试在转换之前在表达式中扩展这些术语。
基于this https://groups.google.com/forum/#!topic/sympy/qaJGesRbX_0对话中,我写了以下代码:
from sympy import Symbol, Mul, Pow, pprint, Matrix, symbols
from sympy.core import numbers
def pow_to_mul(expr):
"""
Convert integer powers in an expression to Muls, like a**2 => a*a.
"""
pows = list(expr.atoms(Pow))
pows = [p for p in pows if p.as_base_exp()[1]>=0]
if any(not e.is_Integer for b, e in (i.as_base_exp() for i in pows)):
raise ValueError("A power contains a non-integer exponent")
repl = zip(pows, (Mul(*[b]*e,evaluate=False) for b,e in (i.as_base_exp() for i in pows)))
return expr.subs(repl)
它部分有效,但只要幂是乘法的参数就会失败:
>>>_=var('x')
>>>print pow_to_mul((x^3+2*x^2)._sympy_())
2*x**2 + x*x*x
>>>print pow_to_mul((x^2/(1+x^2)+(1-x^2)/(1+x^2))._sympy_())
x**2/(x*x + 1) - (x*x - 1)/(x*x + 1)
为什么?我怎样才能改变这一点?
非常感谢,
如果你编译-ffast-math
编译器会为你做这个优化。如果您使用的是古老的编译器或者不能影响构建过程中使用的优化级别,您可以将用户定义的函数传递给 ccode(使用 SymPy master 分支):
>>> ccode(x**97 + 4*x**7 + 5*x**3 + 3**pi, user_functions={'Pow': [
... (lambda b, e: e.is_Integer and e < 42, lambda b, e: '*'.join([b]*int(e))),
... (lambda b, e: not e.is_Integer, 'pow')]})
'pow(x, 97) + 4*x*x*x*x*x*x*x + 5*x*x*x + pow(3, M_PI)'
本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)