[GiNaC-devel] gcd

[Alexei Sheplyakov]

Hello!

On Tue, Apr 22, 2008 at 03:52:09PM +0200, Jens Vollinga wrote:

> there must be some bug in gcd_heur():
> 
> "gcd(5/4+n, 25/8+n^2+15/4*n)" causes an infinite loop.

I think heur_gcd() *itself* is OK.
 
> Any ideas?

heur_gcd() handles only polynomials over Z.  However, gcd() passes a polynomial
over Q, and ... everything can happen. In particular, if the first polynomial
divides the second one...


1326	// Apply evaluation homomorphism and calculate GCD
1327	ex cp, cq;
1328	ex gamma = heur_gcd(p.subs(x == xi, subs_options::no_pattern), q.subs(x == xi, subs_options::no_pattern), &cp, &cq, var+1).expand();
1329	if (!is_exactly_a<fail>(gamma)) {
1330
1331		// Reconstruct polynomial from GCD of mapped polynomials
1332		ex g = interpolate(gamma, xi, x, maxdeg);
1333 

... we have gamma = 1 here, and (for your polynomials) xi = 9/2. So, interpolate()
will loop forever:

1253	ex e = gamma;
1254	numeric rxi = xi.inverse();
1255	for (int i=0; !e.is_zero(); i++) {
1256		ex gi = e.smod(xi);
1257		g.push_back(gi * power(x, i));
1258		e = (e - gi) * rxi;
1259	}

because gi will be always zero.

So, we need to convert the polynomial from Q[X] into Z[X] somewhere. (Sorry,
no patch yet).

Best regards,
	Alexei

-- 
All science is either physics or stamp collecting.

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 827 bytes
Desc: Digital signature
URL: <http://www.ginac.de/pipermail/ginac-devel/attachments/20080422/4d119043/attachment.sig>