@@ -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
- *            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)
@@ -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);
+       if (g.is_zero()) {
+               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);
+               if (w.is_zero()) {
+                       return res;
+               }
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].
*
- *  @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.
*
@@ -1833,8 +1839,8 @@ static exvector sqrfree_yun(const ex &a, const symbol &x)
*/
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