From: Richard Kreckel Date: Fri, 21 Aug 2020 15:47:19 +0000 (+0200) Subject: Avoid unnecessary expansion in sqrfree_yun(). X-Git-Tag: release_1-8-0~26 X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=commitdiff_plain;h=e65cbada7281938475aa0cf31a3d56405ac91d4e;ds=sidebyside Avoid unnecessary expansion in sqrfree_yun(). Don't make the polynomial primitive inside Yun's algorithm. Computing the primitive part will incur an extra expansion of the polynomial and this may distroy a pre-factored structure in the other variables. In the multi- variate case, working with primitive polynomials isn't sufficient anyways to avoid a final division to discover lost factors. --- diff --git a/ginac/normal.cpp b/ginac/normal.cpp index 6ad61137..f8615499 100644 --- a/ginac/normal.cpp +++ b/ginac/normal.cpp @@ -1797,8 +1797,7 @@ ex lcm(const ex &a, const ex &b, bool check_args) * @return vector of expairs (factor, exponent), sorted by exponents */ static epvector sqrfree_yun(const ex &a, const symbol &x) { - // make polynomial primitive w.r.t. x - ex w = a.primpart(x); + ex w = a; ex z = w.diff(x); ex g = gcd(w, z); if (g.is_zero()) { @@ -1814,6 +1813,10 @@ static epvector sqrfree_yun(const ex &a, const symbol &x) ex i = 0; // exponent do { w = quo(w, g, x); + if (w.is_zero()) { + // hidden zero + break; + } z = quo(z, g, x) - w.diff(x); i += 1; if (w.is_equal(x)) { @@ -1829,20 +1832,18 @@ static epvector sqrfree_yun(const ex &a, const symbol &x) } while (!z.is_zero()); // correct for lost factor - // (NB: This factor is not necessarily the unit and content lost when - // we made the polynomial primitive since Yun's algorithm only finds - // factors up to a unit.) + // (being based on GCDs, Yun's algorithm only finds factors up to a unit) const ex lost_factor = quo(a, mul{factors}, x); if (lost_factor.is_equal(_ex1)) { // trivial lost factor return factors; } - if (factors[0].coeff.is_equal(1)) { + if (!factors.empty() && factors[0].coeff.is_equal(1)) { // multiply factor^1 with lost_factor factors[0].rest *= lost_factor; return factors; } - // only factors with higher exponent: prepend lost_factor^1 to the results + // no factor^1: prepend lost_factor^1 to the results epvector results = {expair(lost_factor, 1)}; std::move(factors.begin(), factors.end(), std::back_inserter(results)); return results;