Well, the thing that bugs me the most is that you can do the same things with what would appear to be rational functions:
In [7]: Poly(1/x) Out[7]: Poly(1/x, 1/x, domain='ZZ') In [8]: Poly(1/x)*x Out[8]: Poly(x*1/x, 1/x, domain='ZZ[x]') In [9]: Poly(1/x)*x - 1 Out[9]: Poly(x*1/x - 1, 1/x, domain='ZZ[x]') In [10]: (Poly(1/x)*x - 1).is_zero Out[10]: False This was the original problem in issue 2032. On Jun 16, 2011, at 2:51 PM, Mateusz Paprocki wrote: > Hi, > > On 16 June 2011 13:28, Aaron Meurer <asmeu...@gmail.com> wrote: > This is related to > http://code.google.com/p/sympy/issues/detail?id=2032. The polys > pretend that they can work in K[x, 1/x], but they actually do not > implement things properly, which can lead to wrong results: > > In [3]: Poly(exp(-x)) > Out[3]: Poly(exp(-x), exp(-x), domain='ZZ') > > In [4]: Poly(exp(-x))*exp(x) > Out[4]: Poly(exp(x)*exp(-x), exp(-x), domain='ZZ[exp(x)]') > > In [5]: Poly(exp(-x))*exp(x) - 1 > Out[5]: Poly(exp(x)*exp(-x) - 1, exp(-x), domain='ZZ[exp(x)]') > > In [6]: (Poly(exp(-x))*exp(x) - 1).is_zero > Out[6]: False > > And I could go further, as an incorrect is_zero result can easily be > exacerbated to larger wrong results, but you get the idea. > > Result [6] is actually a valid result, because generators of polynomial > expressions are assumed to be algebraically independent, so in [6] you can > write (Poly(u)*v - 1).is_zero (which of course is False) as well. If you want > [6] to give True, then you would have to imply an algebraic relation between > exp(-x) and exp(x) (really u and v). A reasonable solution I see to this > problem is to implement LaurentPoly which would allow negative exponents, > which in turn would mean that exp(-x) and exp(x) would result in the same > generator, with exponent -1 and 1 respectively. > Another solution would be to add postprocessor to Poly, which would apply > known algebraic relations of generators, after performing the main > computation. This is a great idea. You see, there is actually a complete algorithm to find algebraic relations on elementary functions, called the Risch structure theorems. The transcendental Risch algorithm requires that each elementary extension be transcendental over the others, or else certain results (namely nonelementary results) will not be correct. The algorithm not only determines if a function is transcendental over a differential extension, but provides the algebraic relation when it is not. So it was necessary to implement this algorithm. So far, only the algorithm for exponentials and logarithms has been implemented, but there is also an algorithm for tangents, arctangents, and nonelementary extensions (functions defined in terms of a nonelementary integral of a transcendental function, like erf()). You can handle cosine and sine by first converting them to tangents. If you want to play around with this, the easiest way is to import DifferentialExtension from sympy.integrals.risch in my integration3 branch and create a differential extension out of the expression. For example, you can see it handles the above example with ease In [1]: from sympy.integrals.risch import DifferentialExtension In [2]: DifferentialExtension(exp(x) + exp(-x), x) Out[2]: (Poly(_t0**2 + 1, _t0, domain='ZZ'), Poly(_t0, _t0, domain='ZZ'), [Poly(1, x, domain='ZZ'), Poly(_t0, _t0, domain='ZZ')], [x, _t0], [Lambda(_i, exp(_i))], [], [1], [x], [], []) The output consists of a lot of stuff, but the main thing is that the first two Polys are the numerator and denominator of the resulting expression, and the list of Lambdas near the end gives the function definition of t0, t1, t2, etc. So in this case, we see that it has rewritten exp(x) + exp(-x) as (exp(x)**2 + 1)/(exp(x)). Note that I have implemented heuristics to take the exponentials in such a way to avoid algebraic relations, like In [3]: DifferentialExtension(exp(x) + exp(x/2) + exp(2*x), x) Out[3]: (Poly(_t0**4 + _t0**2 + _t0, _t0, domain='ZZ'), Poly(1, _t0, domain='ZZ'), [Poly(1, x, domain='ZZ'), Poly(1/2*_t0, _t0, domain='QQ')], [x, _t0], [Lambda(_i, exp(_i/2))], [], [1], [x/2], [], []) Notice that t0 is exp(x/2). If we had made t0 exp(x), then exp(x/2) would have been sqrt(t0), which is no good for the Risch algorithm (and also in the polys), because we can't handle algebraic relations yet. Note that the algorithms here are very powerful, and can handle much more nontrivial cases than any heuristic would In [4]: DifferentialExtension(log(exp(x) + 1) + log(exp(-x) + 1), x) Out[4]: (Poly(2*_t1 - x, _t1, domain='ZZ[x]'), Poly(1, _t1, domain='ZZ'), [Poly(1, x, domain='ZZ'), Poly(_t0, _t0, domain='ZZ'), Poly(_t0/(1 + _t0), _t1, domain='ZZ(_t0)')], [x, _t0, _t1], [Lambda(_i, exp(_i) ), Lambda(_i, log(1 + _t0))], [], [1], [x], [2], [1 + _t0]) (rewrote log(exp(x) + 1) + log(exp(-x) + 1) as 2*log(exp(x) + 1) - x. The nature of the algorithms is to always rewrite the expression in terms of some part of the rest. Sometimes, it is necessary to rewrite things in a more canonical form before integrating, for example In [5]: DifferentialExtension(2**x, x) Out[5]: (Poly(_t0, _t0, domain='ZZ'), Poly(1, _t0, domain='ZZ'), [Poly(1, x, domain='ZZ'), Poly(log(2)*_t0, _t0, domain='EX')], [x, _t0], [Lambda(_i, exp(_i*log(2)))], [(exp(x*log(2)), 2**x)], [1], [x*log(2)] , [], []) 2**x is rewritten as exp(x*log(2)). But the object contains an attribute backsubs which puts it back to the way it was after integration In [6]: DifferentialExtension(2**x, x).backsubs Out[6]: ⎡⎛ x⋅log(2) x⎞⎤ ⎣⎝ℯ , 2 ⎠⎦ Finally, I want to note that the algorithm for sin, cos, and tan is very nontrivial. It's described in Bronstein's paper, "Simplification of Real Elementary Functions," and involves rewriting all sin's and cos's as tangents (using the half angle formulas), then the dependence on the tangents is given by Flory's theorem (rewriting tangents in terms of symmetric polynomials). Right now, this code is all written with integration in mind, but it would not be too difficult to adapt it to other uses, like finding algebraic dependancies in Poly). By the way, this all applies only to the univariate case. I don't know if any of it can be extended to the multivariate case. Aaron Meurer > > > Aaron Meurer > > On Thu, Jun 16, 2011 at 12:15 PM, Mateusz Paprocki <matt...@gmail.com> wrote: > > Hi, > > > > On 16 June 2011 20:09, smichr <smi...@gmail.com> wrote: > >> > >> The exp(x)*exp(-x) term in the Poly should cancel, shouldn't it? > >> >>> Poly(exp(x) + exp(-x) - y)*exp(x) > >> Poly(-y*exp(x) + exp(-x)*exp(x) + exp(x)**2, y, exp(-x), exp(x), > >> domain='ZZ') > > > > Polynomials can't contain negative exponents, so exp(x) and exp(-x) = > > 1/exp(x) are treated as two different generators. > > > >> > >> -- > >> You received this message because you are subscribed to the Google Groups > >> "sympy" group. > >> To post to this group, send email to sympy@googlegroups.com. > >> To unsubscribe from this group, send email to > >> sympy+unsubscr...@googlegroups.com. > >> For more options, visit this group at > >> http://groups.google.com/group/sympy?hl=en. > >> > > > > Mateusz > > > > -- > > You received this message because you are subscribed to the Google Groups > > "sympy" group. > > To post to this group, send email to sympy@googlegroups.com. > > To unsubscribe from this group, send email to > > sympy+unsubscr...@googlegroups.com. > > For more options, visit this group at > > http://groups.google.com/group/sympy?hl=en. > > > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To post to this group, send email to sympy@googlegroups.com. > To unsubscribe from this group, send email to > sympy+unsubscr...@googlegroups.com. > For more options, visit this group at > http://groups.google.com/group/sympy?hl=en. > > > Mateusz > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To post to this group, send email to sympy@googlegroups.com. > To unsubscribe from this group, send email to > sympy+unsubscr...@googlegroups.com. > For more options, visit this group at > http://groups.google.com/group/sympy?hl=en. -- You received this message because you are subscribed to the Google Groups "sympy" group. To post to this group, send email to sympy@googlegroups.com. To unsubscribe from this group, send email to sympy+unsubscr...@googlegroups.com. For more options, visit this group at http://groups.google.com/group/sympy?hl=en.
