if (divide_in_z(p, g, ca ? *ca : dummy, var) && divide_in_z(q, g, cb ? *cb : dummy, var)) {
g *= gc;
ex lc = g.lcoeff(*x);
- if (is_ex_exactly_of_type(lc, numeric) && lc.compare(_ex0()) < 0)
+ if (is_ex_exactly_of_type(lc, numeric) && ex_to_numeric(lc).is_negative())
return -g;
else
return g;
ex gcd(const ex &a, const ex &b, ex *ca, ex *cb, bool check_args)
{
+//clog << "gcd(" << a << "," << b << ")\n";
+
// Partially factored cases (to avoid expanding large expressions)
if (is_ex_exactly_of_type(a, mul)) {
if (is_ex_exactly_of_type(b, mul) && b.nops() > a.nops())
// as defined by get_first_symbol() is made positive)
const symbol *x;
if (get_first_symbol(den, x)) {
- if (den.unit(*x).compare(_ex0()) < 0) {
+ GINAC_ASSERT(is_ex_exactly_of_type(den.unit(*x),numeric));
+ if (ex_to_numeric(den.unit(*x)).is_negative()) {
num *= _ex_1();
den *= _ex_1();
}
}
// Return result as list
+//clog << " returns num = " << num << ", den = " << den << ", pre_factor = " << pre_factor << endl;
return (new lst(num * pre_factor.numer(), den * pre_factor.denom()))->setflag(status_flags::dynallocated);
}
// split it into numerator and denominator
GINAC_ASSERT(ex_to_numeric(ex_to_add(n.op(0)).overall_coeff).is_rational());
numeric overall = ex_to_numeric(ex_to_add(n.op(0)).overall_coeff);
- o.push_back((new lst(overall.numer(), overall.denom()))->setflag(status_flags::dynallocated));
+ o.push_back((new lst(overall.numer(), overall.denom() * n.op(1)))->setflag(status_flags::dynallocated));
} else
o.push_back(n);
it++;
// Determine common denominator
ex den = _ex1();
exvector::const_iterator ait = o.begin(), aitend = o.end();
+//clog << "add::normal uses the following summands:\n";
while (ait != aitend) {
+//clog << " num = " << ait->op(0) << ", den = " << ait->op(1) << endl;
den = lcm(ait->op(1), den, false);
ait++;
}
+//clog << " common denominator = " << den << endl;
// Add fractions
if (den.is_equal(_ex1())) {
// should not happen
throw(std::runtime_error("invalid expression in add::normal, division failed"));
}
- num_seq.push_back(ait->op(0) * q);
+ num_seq.push_back((ait->op(0) * q).expand());
}
ex num = (new add(num_seq))->setflag(status_flags::dynallocated);
* @see ex::normal */
ex power::normal(lst &sym_lst, lst &repl_lst, int level) const
{
- if (exponent.info(info_flags::posint)) {
- // Integer powers are distributed
- ex n = basis.bp->normal(sym_lst, repl_lst, level-1);
- return (new lst(power(n.op(0), exponent), power(n.op(1), exponent)))->setflag(status_flags::dynallocated);
- } else if (exponent.info(info_flags::negint)) {
- // Integer powers are distributed
- ex n = basis.bp->normal(sym_lst, repl_lst, level-1);
- return (new lst(power(n.op(1), -exponent), power(n.op(0), -exponent)))->setflag(status_flags::dynallocated);
- } else {
- // Non-integer powers are replaced by temporary symbol (after normalizing basis)
- ex n = basis.bp->normal(sym_lst, repl_lst, level-1);
- return (new lst(replace_with_symbol(power(n.op(0) / n.op(1), exponent), sym_lst, repl_lst), _ex1()))->setflag(status_flags::dynallocated);
+ // Normalize basis
+ ex n = basis.bp->normal(sym_lst, repl_lst, level-1);
+
+ if (exponent.info(info_flags::integer)) {
+
+ if (exponent.info(info_flags::positive)) {
+
+ // (a/b)^n -> {a^n, b^n}
+ return (new lst(power(n.op(0), exponent), power(n.op(1), exponent)))->setflag(status_flags::dynallocated);
+
+ } else if (exponent.info(info_flags::negative)) {
+
+ // (a/b)^-n -> {b^n, a^n}
+ return (new lst(power(n.op(1), -exponent), power(n.op(0), -exponent)))->setflag(status_flags::dynallocated);
+ }
+
+ } else {
+
+ if (exponent.info(info_flags::positive)) {
+
+ // (a/b)^x -> {sym((a/b)^x), 1}
+ return (new lst(replace_with_symbol(power(n.op(0) / n.op(1), exponent), sym_lst, repl_lst), _ex1()))->setflag(status_flags::dynallocated);
+
+ } else if (exponent.info(info_flags::negative)) {
+
+ if (n.op(1).is_equal(_ex1())) {
+
+ // a^-x -> {1, sym(a^x)}
+ return (new lst(_ex1(), replace_with_symbol(power(n.op(0), -exponent), sym_lst, repl_lst)))->setflag(status_flags::dynallocated);
+
+ } else {
+
+ // (a/b)^-x -> {sym((b/a)^x), 1}
+ return (new lst(replace_with_symbol(power(n.op(1) / n.op(0), -exponent), sym_lst, repl_lst), _ex1()))->setflag(status_flags::dynallocated);
+ }
+
+ } else { // exponent not numeric
+
+ // (a/b)^x -> {sym((a/b)^x, 1}
+ return (new lst(replace_with_symbol(power(n.op(0) / n.op(1), exponent), sym_lst, repl_lst), _ex1()))->setflag(status_flags::dynallocated);
+ }
}
}
return e.op(0) / e.op(1);
}
+/** Numerator of an expression. If the expression is not of the normal form
+ * "numerator/denominator", it is first converted to this form and then the
+ * numerator is returned.
+ *
+ * @see ex::normal
+ * @return numerator */
+ex ex::numer(void) const
+{
+ lst sym_lst, repl_lst;
+
+ ex e = bp->normal(sym_lst, repl_lst, 0);
+ GINAC_ASSERT(is_ex_of_type(e, lst));
+
+ // Re-insert replaced symbols
+ if (sym_lst.nops() > 0)
+ return e.op(0).subs(sym_lst, repl_lst);
+ else
+ return e.op(0);
+}
+
+/** Denominator of an expression. If the expression is not of the normal form
+ * "numerator/denominator", it is first converted to this form and then the
+ * denominator is returned.
+ *
+ * @see ex::normal
+ * @return denominator */
+ex ex::denom(void) const
+{
+ lst sym_lst, repl_lst;
+
+ ex e = bp->normal(sym_lst, repl_lst, 0);
+ GINAC_ASSERT(is_ex_of_type(e, lst));
+
+ // Re-insert replaced symbols
+ if (sym_lst.nops() > 0)
+ return e.op(1).subs(sym_lst, repl_lst);
+ else
+ return e.op(1);
+}
+
#ifndef NO_NAMESPACE_GINAC
} // namespace GiNaC
#endif // ndef NO_NAMESPACE_GINAC