Almost the same patch for GiNaC 1.3, just in a case... The following transformation is done automagically now: (-x+y)*(x-y) -> -(x-y)^2 (3*x+y)/(6*x+2*y)^2 -> 1/4*(y+3*x)^(-1) --- ginac/mul.cpp | 55 ++++++++++++++++++++++++++++++++++++++++++++++++++++++- ginac/power.cpp | 22 ++++++++++++++++++++++ 2 files changed, 76 insertions(+), 1 deletions(-) diff --git a/ginac/mul.cpp b/ginac/mul.cpp index cbab7e7..8863949 100644 --- a/ginac/mul.cpp +++ b/ginac/mul.cpp @@ -440,7 +440,60 @@ ex mul::eval(int level) const ex_to<numeric>(addref.overall_coeff). mul_dyn(ex_to<numeric>(overall_coeff)))) ->setflag(status_flags::dynallocated | status_flags::evaluated); - } + } else if (seq_size >= 2) { + // Strip the content and the unit part from the each term. Thus + // things like (-x+a)*(3*x-3*a) automagically turn into - 3*(x-a)^2 + + epvector::const_iterator last = seq.end(); + epvector::const_iterator i = seq.begin(); + while (i!=last) { + if (! (is_a<add>(i->rest) && i->coeff.is_equal(_ex1))) { + // power::eval has such a rule, no need to handle powers here + ++i; + continue; + } + + // XXX: What is the best way to check if the polynomial is a primitive? + numeric c = i->rest.integer_content(); + const numeric& lead_coeff = + ex_to<numeric>(ex_to<add>(i->rest).seq.begin()->coeff).div_dyn(c); + const bool canonicalizable = lead_coeff.is_integer(); + + // XXX: The main variable is chosen in a random way, so this code + // does NOT transform the term into the canonical form (thus, in some + // very unlucky event it can even loop forever). Hopefully the main + // variable will be the same for all terms in *this + const bool unit_normal = lead_coeff.is_pos_integer(); + if ((c == *_num1_p) && ((! canonicalizable) || unit_normal)) { + ++i; + continue; + } + + std::auto_ptr<epvector> s(new epvector); + s->reserve(seq.size()); + + epvector::const_iterator j=seq.begin(); + while (j!=i) { + s->push_back(*j); + ++j; + } + + if (! unit_normal) + c = c.mul(*_num_1_p); + const ex primitive = (i->rest)/c; + s->push_back(expair(primitive, _ex1)); + ++j; + + while (j!=last) { + s->push_back(*j); + ++j; + } + return (new mul(s, + ex_to<numeric>(overall_coeff).mul_dyn(c)) + )->setflag(status_flags::dynallocated); + } + } + return this->hold(); } diff --git a/ginac/power.cpp b/ginac/power.cpp index a63556b..0573dc6 100644 --- a/ginac/power.cpp +++ b/ginac/power.cpp @@ -463,6 +463,28 @@ ex power::eval(int level) const return expand_mul(ex_to<mul>(ebasis), *num_exponent, 0); } + // (2*x + 6*y)^(-4) -> 1/16*(x + 3*y)^(-4) + if (num_exponent->is_integer() && is_exactly_a<add>(ebasis)) { + const numeric icont = ebasis.integer_content(); + const numeric& lead_coeff = + ex_to<numeric>(ex_to<add>(ebasis).seq.begin()->coeff).div_dyn(icont); + + const bool canonicalizable = lead_coeff.is_integer(); + const bool unit_normal = lead_coeff.is_pos_integer(); + + if (icont != *_num1_p) + return (new mul(power(ebasis/icont, *num_exponent), + power(icont, *num_exponent)))->setflag(status_flags::dynallocated); + + if (canonicalizable && (! unit_normal)) { + if (num_exponent->is_even()) + return power(-ebasis, *num_exponent); + else + return (new mul(power(-ebasis, *num_exponent), + *_num_1_p))->setflag(status_flags::dynallocated); + } + } + // ^(*(...,x;c1),c2) -> *(^(*(...,x;1),c2),c1^c2) (c1, c2 numeric(), c1>0) // ^(*(...,x;c1),c2) -> *(^(*(...,x;-1),c2),(-c1)^c2) (c1, c2 numeric(), c1<0) if (is_exactly_a<mul>(ebasis)) { -- 1.4.4.4 Best regards, Alexei -- All science is either physics or stamp collecting.