[bugfix] fix sqrfree(poly) for zero polynomials in disguise.

Yun's algorithm now handles polynomials which would become zero after
expanding.

diff --git a/check/exam_paranoia.cpp b/check/exam_paranoia.cpp
index 579ddc57087b25c06078fd1f566dd207cce9a10b..5c8ef03ccffbd555256362d19a32f185fee1bef8 100644 (file)
--- a/check/exam_paranoia.cpp
+++ b/check/exam_paranoia.cpp
@@ -577,6 +577,52 @@ static unsigned exam_paranoia23()
        return result;
 }
 
        return result;
 }
 

+// Bug in sqrfree_yun (fixed 2016-02-02).
+static unsigned exam_paranoia24()
+{
+       unsigned result = 0;
+       symbol x("x");
+       ex e;
+
+       e = (x-1)*(x+1) - x*x + 1;  // an unexpanded 0...
+       try {
+               ex f = sqrfree(e);
+               if (!f.is_zero()) {
+                       clog << "sqrfree(" << e << ") returns " << f << " instead of 0\n";
+                       ++result;
+               }
+       } catch (const exception &err) {
+               clog << "sqrfree(" << e << ") throws " << err.what() << endl;
+               ++result;
+       }
+
+       e = pow(x-1,3) - expand(pow(x-1,3));  // ...still after differentiating...
+       try {
+               ex f = sqrfree(e);
+               if (!f.is_zero()) {
+                       clog << "sqrfree(" << e << ") returns " << f << " instead of 0\n";
+                       ++result;
+               }
+       } catch (const exception &err) {
+               clog << "sqrfree(" << e << ") throws " << err.what() << endl;
+               ++result;
+       }
+
+       e = pow(x-1,4) - expand(pow(x-1,4));  // ...and after differentiating twice.
+       try {
+               ex f = sqrfree(e);
+               if (!f.is_zero()) {
+                       clog << "sqrfree(" << e << ") returns " << f << " instead of 0\n";
+                       ++result;
+               }
+       } catch (const exception &err) {
+               clog << "sqrfree(" << e << ") throws " << err.what() << endl;
+               ++result;
+       }
+
+       return result;
+}
+

 unsigned exam_paranoia()
 {
        unsigned result = 0;
 unsigned exam_paranoia()
 {
        unsigned result = 0;


@@ -607,6 +653,7 @@ unsigned exam_paranoia()
        result += exam_paranoia21();  cout << '.' << flush;
        result += exam_paranoia22();  cout << '.' << flush;
        result += exam_paranoia23();  cout << '.' << flush;
        result += exam_paranoia21();  cout << '.' << flush;
        result += exam_paranoia22();  cout << '.' << flush;
        result += exam_paranoia23();  cout << '.' << flush;

+       result += exam_paranoia24();  cout << '.' << flush;

        
        return result;
 }
        
        return result;
 }

diff --git a/ginac/normal.cpp b/ginac/normal.cpp
index 439a644606388ad1f70fc08d1d537f3a7189a7be..a8a1e1eb1bd3f7282418eff3160f2ba1a25114e6 100644 (file)
--- a/ginac/normal.cpp
+++ b/ginac/normal.cpp
@@ -1771,7 +1771,7 @@ ex lcm(const ex &a, const ex &b, bool check_args)
  *  Yun's algorithm.  Used internally by sqrfree().
  *
  *  @param a  multivariate polynomial over Z[X], treated here as univariate
  *  Yun's algorithm.  Used internally by sqrfree().
  *
  *  @param a  multivariate polynomial over Z[X], treated here as univariate

- *            polynomial in x.
+ *            polynomial in x (needs not be expanded).

  *  @param x  variable to factor in
  *  @return   vector of factors sorted in ascending degree */
 static exvector sqrfree_yun(const ex &a, const symbol &x)
  *  @param x  variable to factor in
  *  @return   vector of factors sorted in ascending degree */
 static exvector sqrfree_yun(const ex &a, const symbol &x)


@@ -1780,6 +1780,9 @@ static exvector sqrfree_yun(const ex &a, const symbol &x)
        ex w = a;
        ex z = w.diff(x);
        ex g = gcd(w, z);
        ex w = a;
        ex z = w.diff(x);
        ex g = gcd(w, z);

+       if (g.is_zero()) {
+               return res;
+       }

        if (g.is_equal(_ex1)) {
                res.push_back(a);
                return res;
        if (g.is_equal(_ex1)) {
                res.push_back(a);
                return res;


@@ -1787,6 +1790,9 @@ static exvector sqrfree_yun(const ex &a, const symbol &x)
        ex y;
        do {
                w = quo(w, g, x);
        ex y;
        do {
                w = quo(w, g, x);

+               if (w.is_zero()) {
+                       return res;
+               }

                y = quo(z, g, x);
                z = y - w.diff(x);
                g = gcd(w, z);
                y = quo(z, g, x);
                z = y - w.diff(x);
                g = gcd(w, z);


@@ -1798,7 +1804,7 @@ static exvector sqrfree_yun(const ex &a, const symbol &x)
 
 /** Compute a square-free factorization of a multivariate polynomial in Q[X].
  *
 
 /** Compute a square-free factorization of a multivariate polynomial in Q[X].
  *

- *  @param a  multivariate polynomial over Q[X]
+ *  @param a  multivariate polynomial over Q[X] (needs not be expanded)

  *  @param l  lst of variables to factor in, may be left empty for autodetection
  *  @return   a square-free factorization of \p a.
  *
  *  @param l  lst of variables to factor in, may be left empty for autodetection
  *  @return   a square-free factorization of \p a.
  *


@@ -1833,8 +1839,8 @@ static exvector sqrfree_yun(const ex &a, const symbol &x)
  */
 ex sqrfree(const ex &a, const lst &l)
 {
  */
 ex sqrfree(const ex &a, const lst &l)
 {

-       if (is_exactly_a<numeric>(a) ||     // algorithm does not trap a==0
-           is_a<symbol>(a))        // shortcut
+       if (is_exactly_a<numeric>(a) ||
+           is_a<symbol>(a))        // shortcuts

                return a;
 
        // If no lst of variables to factorize in was specified we have to
                return a;
 
        // If no lst of variables to factorize in was specified we have to