[bugfix] chinrem_gcd: handle polynomials over rationals properly.

* extract_integer_content():
  integer_content() can also return a rational number, e.g if the expression
  is a polynomial over rationals (in this case expr/content is polynomial
  over integers with integer content being 1). Therefore check if
  integer_content() is really an integer (and if it's not just return 1,
  GCD for polynomials over a field is defined up to arbitrary element of
  the field).

  This fixes possible segfault when computing GCD of polynomials over
  rationals (this is not theoretical, see the added test case).

Thanks to Ernst Moritz Hahn for reporting this bug.

diff --git a/check/bugme_chinrem_gcd.cpp b/check/bugme_chinrem_gcd.cpp
index d3c0330ebbb7238a1ade185fbf7a4fbf8fa1a557..87c718c7114cf4ceff3ad617b092ab959c0454b3 100644 (file)
--- a/check/bugme_chinrem_gcd.cpp
+++ b/check/bugme_chinrem_gcd.cpp
@@ -38,13 +38,27 @@ static const std::string p1_srep("\
 static const std::string p2_srep("\
 1882371920+29943249139*x^9-21558061051*x^35+24497174109*x^69+3363043648*x^25+5186524611*x^98-17753230977*x^56+16461882236*x^72+11039962159*x^11-85814533599*x^43-12986831645*x^85+4813620791*x^29-2133682609*x^2+9141433582*x^88-14841292646*x^13+19494168654*x^51-426278523*x^15-18246235652*x^59-12424469513*x^16-14414308862*x^38-16262001424*x^75+1584505491*x^3+6616907060*x^46+9411879011*x^91+7056872613*x^20+29675566382*x^78-48441925739*x^54+12038293769*x^33-22463329638*x^65+21957440404*x^41+26252749471*x^81-28689993421*x^24+1190623161*x^94-3323333429*x^62+778770195*x^68-23673866858*x^49+10263027507*x^97+29115114125*x^28-34230002076*x^36-217623403*x^71-6344703601*x^84+2789684836*x^100-44973929094*x^57-6061422988*x^44+18964048763*x^87+3160532878*x^8-8039690791*x^19-5750277357*x^74+5544486596*x^23+1800283356*x^90-8174921854*x^52+2577247520*x^64-1088265300*x^10+18566882873*x^39+12678193001*x^77-7204610489*x^27+9980611956*x^60+12672890141*x^80+4463462090*x^12-30937311949*x^47-3883570155*x^93+7561875663*x^14-3850452489*x^55+20714527103*x^31-9973372012*x^34+440594870*x^67+10385086102*x^96-20693764726*x^63+11049483001*x^18-11578701274*x^70-5548876327*x^42+20393397488*x^83+20531692560*x^50+1309311388*x^99-7521320830*x-2750892631*x^37-6001481047*x^73+7149046134*x^22+10287410396*x^86+7332053562*x^58-1135211878*x^4-7420079075*x^45+9742202475*x^89-214559874*x^26+29530802947*x^30-2280622009*x^32-4237029618*x^5+13760397067*x^76+18073788685*x^53+2485463829*x^40+1889202286*x^79+8841198971*x^6+13562767020*x^92+110919258*x^7+6961020716*x^61+11700952246*x^17-13528331424*x^66-19801882818*x^21+25061328813*x^82+1553111326*x^48+4638169279*x^95");
 
-int main(int argc, char** argv)
+static void check_poly_cra()
 {
-       cout << "checking for bugs in poly_cra() " << flush;
        parser the_parser;
        ex p1 = the_parser(p1_srep);
        ex p2 = the_parser(p2_srep);
        ex g = chinrem_gcd(p1, p2);
+}
+
+static void check_extract_integer_content()
+{
+       parser readme;
+       ex A = readme("1/282901891126422365*(x + y)");
+       ex B = readme("165888/282901891126422365*(x - y)");
+       ex g = chinrem_gcd(A, B);
+}
+
+int main(int argc, char** argv)
+{
+       cout << "checking for bugs in poly_cra() and friends " << flush;
+       check_poly_cra();
+       check_extract_integer_content();
        cout << "not found.";
        return 0;
 }

diff --git a/ginac/polynomial/mgcd.cpp b/ginac/polynomial/mgcd.cpp
index 24da934690237e42e26a97caa35b09080eabc672..5b9da128223d8a6d9024b24e775398b2e3ed8b53 100644 (file)
--- a/ginac/polynomial/mgcd.cpp
+++ b/ginac/polynomial/mgcd.cpp
@@ -30,6 +30,9 @@
 #include <numeric> // std::accumulate
 
 #include <cln/integer.h>
+#include <cln/integer_ring.h>
+#include <cln/rational.h>
+#include <cln/rational_ring.h>
 
 namespace GiNaC {
 
@@ -37,6 +40,13 @@ static cln::cl_I extract_integer_content(ex& Apr, const ex& A)
 {
        static const cln::cl_I n1(1);
        const numeric icont_ = A.integer_content();
+       if (cln::instanceof(icont_.to_cl_N(), cln::cl_RA_ring)) {
+               Apr = (A/icont_).expand();
+               // A is a polynomail over rationals, so GCD is defined
+               // up to arbitrary rational number.
+               return n1;
+       }
+       GINAC_ASSERT(cln::instanceof(icont_.to_cl_N(), cln::cl_I_ring));
        const cln::cl_I icont = cln::the<cln::cl_I>(icont_.to_cl_N());
        if (icont != 1) {
                Apr = (A/icont_).expand();