static unsigned check_normal(const ex &e, const ex &d)
{
ex en = e.normal();
- if (en.compare(d) != 0) {
+ if (!en.is_equal(d)) {
clog << "normal form of " << e << " erroneously returned "
<< en << " (should be " << d << ")" << endl;
return 1;
return result;
}
+/* Test content(), integer_content(), primpart(). */
+static unsigned check_content(const ex & e, const ex & x, const ex & ic, const ex & c, const ex & pp)
+{
+ unsigned result = 0;
+
+ ex r_ic = e.integer_content();
+ if (!r_ic.is_equal(ic)) {
+ clog << "integer_content(" << e << ") erroneously returned "
+ << r_ic << " instead of " << ic << endl;
+ ++result;
+ }
+
+ ex r_c = e.content(x);
+ if (!r_c.is_equal(c)) {
+ clog << "content(" << e << ", " << x << ") erroneously returned "
+ << r_c << " instead of " << c << endl;
+ ++result;
+ }
+
+ ex r_pp = e.primpart(x);
+ if (!r_pp.is_equal(pp)) {
+ clog << "primpart(" << e << ", " << x << ") erroneously returned "
+ << r_pp << " instead of " << pp << endl;
+ ++result;
+ }
+
+ ex r = r_c*r_pp*e.unit(x);
+ if (!(r - e).expand().is_zero()) {
+ clog << "product of unit, content, and primitive part of " << e << " yielded "
+ << r << " instead of " << e << endl;
+ ++result;
+ }
+
+ return result;
+}
+
+static unsigned exam_content()
+{
+ unsigned result = 0;
+ symbol x("x"), y("y");
+
+ result += check_content(ex(-3)/4, x, ex(3)/4, ex(3)/4, 1);
+ result += check_content(-x/4, x, ex(1)/4, ex(1)/4, x);
+ result += check_content(5*x-15, x, 5, 5, x-3);
+ result += check_content(5*x*y-15*y*y, x, 5, 5*y, x-3*y);
+ result += check_content(-15*x/2+ex(25)/3, x, ex(5)/6, ex(5)/6, 9*x-10);
+
+ return result;
+}
+
unsigned exam_normalization()
{
unsigned result = 0;
cout << "examining rational function normalization" << flush;
clog << "----------rational function normalization:" << endl;
- result += exam_normal1(); cout << '.' << flush;
- result += exam_normal2(); cout << '.' << flush;
- result += exam_normal3(); cout << '.' << flush;
- result += exam_normal4(); cout << '.' << flush;
+ result += exam_normal1(); cout << '.' << flush;
+ result += exam_normal2(); cout << '.' << flush;
+ result += exam_normal3(); cout << '.' << flush;
+ result += exam_normal4(); cout << '.' << flush;
+ result += exam_content(); cout << '.' << flush;
if (!result) {
cout << " passed " << endl;
/** Compute the integer content (= GCD of all numeric coefficients) of an
- * expanded polynomial.
+ * expanded polynomial. For a polynomial with rational coefficients, this
+ * returns g/l where g is the GCD of the coefficients' numerators and l
+ * is the LCM of the coefficients' denominators.
*
* @return integer content */
numeric ex::integer_content() const
{
epvector::const_iterator it = seq.begin();
epvector::const_iterator itend = seq.end();
- numeric c = _num0;
+ numeric c = _num0, l = _num1;
while (it != itend) {
GINAC_ASSERT(!is_exactly_a<numeric>(it->rest));
GINAC_ASSERT(is_exactly_a<numeric>(it->coeff));
- c = gcd(ex_to<numeric>(it->coeff), c);
+ c = gcd(ex_to<numeric>(it->coeff).numer(), c);
+ l = lcm(ex_to<numeric>(it->coeff).denom(), l);
it++;
}
GINAC_ASSERT(is_exactly_a<numeric>(overall_coeff));
- c = gcd(ex_to<numeric>(overall_coeff),c);
- return c;
+ c = gcd(ex_to<numeric>(overall_coeff).numer(), c);
+ l = lcm(ex_to<numeric>(overall_coeff).denom(), l);
+ return c/l;
}
numeric mul::integer_content() const
*/
/** Compute unit part (= sign of leading coefficient) of a multivariate
- * polynomial in Z[x]. The product of unit part, content part, and primitive
+ * polynomial in Q[x]. The product of unit part, content part, and primitive
* part is the polynomial itself.
*
* @param x variable in which to compute the unit part
/** Compute content part (= unit normal GCD of all coefficients) of a
- * multivariate polynomial in Z[x]. The product of unit part, content part,
+ * multivariate polynomial in Q[x]. The product of unit part, content part,
* and primitive part is the polynomial itself.
*
* @param x variable in which to compute the content part
if (e.is_zero())
return _ex0;
- // First, try the integer content
+ // First, divide out the integer content (which we can calculate very efficiently).
+ // If the leading coefficient of the quotient is an integer, we are done.
ex c = e.integer_content();
ex r = e / c;
ex lcoeff = r.lcoeff(x);
return c;
// GCD of all coefficients
- int deg = e.degree(x);
- int ldeg = e.ldegree(x);
+ int deg = r.degree(x);
+ int ldeg = r.ldegree(x);
if (deg == ldeg)
- return e.lcoeff(x) / e.unit(x);
- c = _ex0;
+ return lcoeff * c;
+ ex cont = _ex0;
for (int i=ldeg; i<=deg; i++)
- c = gcd(e.coeff(x, i), c, NULL, NULL, false);
- return c;
+ cont = gcd(r.coeff(x, i), cont, NULL, NULL, false);
+ return cont * c;
}
-/** Compute primitive part of a multivariate polynomial in Z[x].
+/** Compute primitive part of a multivariate polynomial in Q[x].
* The product of unit part, content part, and primitive part is the
* polynomial itself.
*
}
-/** Compute primitive part of a multivariate polynomial in Z[x] when the
+/** Compute primitive part of a multivariate polynomial in Q[x] when the
* content part is already known. This function is faster in computing the
* primitive part than the previous function.
*
.BI has( expression ", " pattern )
\- returns "1" if the first expression contains the pattern as a subexpression, "0" otherwise
.br
+.BI integer_content( expression )
+\- integer content of a polynomial
+.br
.BI inverse( matrix )
\- inverse of a matrix
.br
return found;
}
+static ex f_integer_content(const exprseq &e)
+{
+ return e[0].expand().integer_content();
+}
+
static ex f_inverse(const exprseq &e)
{
CHECK_ARG(0, matrix, inverse);
{"find", f_find, 2},
{"gcd", f_gcd, 2},
{"has", f_has, 2},
+ {"integer_content", f_integer_content, 1},
{"inverse", f_inverse, 1},
{"iprint", f_dummy, 0}, // for Tab-completion
{"is", f_is, 1},
projects / ginac.git / commitdiff