X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?a=blobdiff_plain;f=ginac%2Fnormal.cpp;h=5e6aa470eae2119d119f5daa1fe3b74f5f4b7c20;hb=0f617ea4ea1dc6a33cae5e2ce518b26d8055c177;hp=ee27868165e02e148960407369b8d18da282f180;hpb=9f410a18b08213185a6ae473295c9945ebdd9985;p=ginac.git diff --git a/ginac/normal.cpp b/ginac/normal.cpp index ee278681..5e6aa470 100644 --- a/ginac/normal.cpp +++ b/ginac/normal.cpp @@ -6,7 +6,7 @@ * computation, square-free factorization and rational function normalization. */ /* - * GiNaC Copyright (C) 1999-2001 Johannes Gutenberg University Mainz, Germany + * GiNaC Copyright (C) 1999-2002 Johannes Gutenberg University Mainz, Germany * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by @@ -381,7 +381,7 @@ ex quo(const ex &a, const ex &b, const symbol &x, bool check_args) int rdeg = r.degree(x); ex blcoeff = b.expand().coeff(x, bdeg); bool blcoeff_is_numeric = is_ex_exactly_of_type(blcoeff, numeric); - exvector v; v.reserve(rdeg - bdeg + 1); + exvector v; v.reserve(std::max(rdeg - bdeg + 1, 0)); while (rdeg >= bdeg) { ex term, rcoeff = r.coeff(x, rdeg); if (blcoeff_is_numeric) @@ -471,14 +471,14 @@ ex decomp_rational(const ex &a, const symbol &x) } -/** Pseudo-remainder of polynomials a(x) and b(x) in Z[x]. +/** Pseudo-remainder of polynomials a(x) and b(x) in Q[x]. * * @param a first polynomial in x (dividend) * @param b second polynomial in x (divisor) * @param x a and b are polynomials in x * @param check_args check whether a and b are polynomials with rational * coefficients (defaults to "true") - * @return pseudo-remainder of a(x) and b(x) in Z[x] */ + * @return pseudo-remainder of a(x) and b(x) in Q[x] */ ex prem(const ex &a, const ex &b, const symbol &x, bool check_args) { if (b.is_zero()) @@ -523,14 +523,14 @@ ex prem(const ex &a, const ex &b, const symbol &x, bool check_args) } -/** Sparse pseudo-remainder of polynomials a(x) and b(x) in Z[x]. +/** Sparse pseudo-remainder of polynomials a(x) and b(x) in Q[x]. * * @param a first polynomial in x (dividend) * @param b second polynomial in x (divisor) * @param x a and b are polynomials in x * @param check_args check whether a and b are polynomials with rational * coefficients (defaults to "true") - * @return sparse pseudo-remainder of a(x) and b(x) in Z[x] */ + * @return sparse pseudo-remainder of a(x) and b(x) in Q[x] */ ex sprem(const ex &a, const ex &b, const symbol &x, bool check_args) { if (b.is_zero()) @@ -581,14 +581,15 @@ ex sprem(const ex &a, const ex &b, const symbol &x, bool check_args) * @param check_args check whether a and b are polynomials with rational * coefficients (defaults to "true") * @return "true" when exact division succeeds (quotient returned in q), - * "false" otherwise */ + * "false" otherwise (q left untouched) */ bool divide(const ex &a, const ex &b, ex &q, bool check_args) { - q = _ex0; if (b.is_zero()) throw(std::overflow_error("divide: division by zero")); - if (a.is_zero()) + if (a.is_zero()) { + q = _ex0; return true; + } if (is_ex_exactly_of_type(b, numeric)) { q = a / b; return true; @@ -611,13 +612,15 @@ bool divide(const ex &a, const ex &b, ex &q, bool check_args) // Polynomial long division (recursive) ex r = a.expand(); - if (r.is_zero()) + if (r.is_zero()) { + q = _ex0; return true; + } int bdeg = b.degree(*x); int rdeg = r.degree(*x); ex blcoeff = b.expand().coeff(*x, bdeg); bool blcoeff_is_numeric = is_ex_exactly_of_type(blcoeff, numeric); - exvector v; v.reserve(rdeg - bdeg + 1); + exvector v; v.reserve(std::max(rdeg - bdeg + 1, 0)); while (rdeg >= bdeg) { ex term, rcoeff = r.coeff(*x, rdeg); if (blcoeff_is_numeric) @@ -778,7 +781,7 @@ static bool divide_in_z(const ex &a, const ex &b, ex &q, sym_desc_vec::const_ite int rdeg = adeg; ex eb = b.expand(); ex blcoeff = eb.coeff(*x, bdeg); - exvector v; v.reserve(rdeg - bdeg + 1); + exvector v; v.reserve(std::max(rdeg - bdeg + 1, 0)); while (rdeg >= bdeg) { ex term, rcoeff = r.coeff(*x, rdeg); if (!divide_in_z(rcoeff, blcoeff, term, var+1)) @@ -1738,11 +1741,41 @@ static exvector sqrfree_yun(const ex &a, const symbol &x) return res; } -/** Compute square-free factorization of 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 x lst of variables to factor in, may be left empty for autodetection - * @return polynomial a in square-free factored form. */ + * @return a square-free factorization of \p a. + * + * \note + * A polynomial \f$p(X) \in C[X]\f$ is said square-free + * if, whenever any two polynomials \f$q(X)\f$ and \f$r(X)\f$ + * are such that + * \f[ + * p(X) = q(X)^2 r(X), + * \f] + * we have \f$q(X) \in C\f$. + * This means that \f$p(X)\f$ has no repeated factors, apart + * eventually from constants. + * Given a polynomial \f$p(X) \in C[X]\f$, we say that the + * decomposition + * \f[ + * p(X) = b \cdot p_1(X)^{a_1} \cdot p_2(X)^{a_2} \cdots p_r(X)^{a_r} + * \f] + * is a square-free factorization of \f$p(X)\f$ if the + * following conditions hold: + * -# \f$b \in C\f$ and \f$b \neq 0\f$; + * -# \f$a_i\f$ is a positive integer for \f$i = 1, \ldots, r\f$; + * -# the degree of the polynomial \f$p_i\f$ is strictly positive + * for \f$i = 1, \ldots, r\f$; + * -# the polynomial \f$\Pi_{i=1}^r p_i(X)\f$ is square-free. + * + * Square-free factorizations need not be unique. For example, if + * \f$a_i\f$ is even, we could change the polynomial \f$p_i(X)\f$ + * into \f$-p_i(X)\f$. + * Observe also that the factors \f$p_i(X)\f$ need not be irreducible + * polynomials. + */ ex sqrfree(const ex &a, const lst &l) { if (is_a(a) || // algorithm does not trap a==0