* Implementation of GiNaC's symbolic exponentiation (basis^exponent). */
/*
- * GiNaC Copyright (C) 1999-2000 Johannes Gutenberg University Mainz, Germany
+ * GiNaC Copyright (C) 1999-2001 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
power::~power()
{
debugmsg("power destructor",LOGLEVEL_DESTRUCT);
- destroy(0);
+ destroy(false);
}
power::power(const power & other)
{
debugmsg("power operator=",LOGLEVEL_ASSIGNMENT);
if (this != &other) {
- destroy(1);
+ destroy(true);
copy(other);
}
return *this;
debugmsg("power print csrc", LOGLEVEL_PRINT);
// Integer powers of symbols are printed in a special, optimized way
- if (exponent.info(info_flags::integer) &&
- (is_ex_exactly_of_type(basis, symbol) ||
- is_ex_exactly_of_type(basis, constant))) {
+ if (exponent.info(info_flags::integer)
+ && (is_ex_exactly_of_type(basis, symbol) || is_ex_exactly_of_type(basis, constant))) {
int exp = ex_to_numeric(exponent).to_int();
if (exp > 0)
os << "(";
int power::degree(const symbol & s) const
{
if (is_exactly_of_type(*exponent.bp,numeric)) {
- if ((*basis.bp).compare(s)==0)
- return ex_to_numeric(exponent).to_int();
- else
+ if ((*basis.bp).compare(s)==0) {
+ if (ex_to_numeric(exponent).is_integer())
+ return ex_to_numeric(exponent).to_int();
+ else
+ return 0;
+ } else
return basis.degree(s) * ex_to_numeric(exponent).to_int();
}
return 0;
int power::ldegree(const symbol & s) const
{
if (is_exactly_of_type(*exponent.bp,numeric)) {
- if ((*basis.bp).compare(s)==0)
- return ex_to_numeric(exponent).to_int();
- else
+ if ((*basis.bp).compare(s)==0) {
+ if (ex_to_numeric(exponent).is_integer())
+ return ex_to_numeric(exponent).to_int();
+ else
+ return 0;
+ } else
return basis.ldegree(s) * ex_to_numeric(exponent).to_int();
}
return 0;
{
if ((*basis.bp).compare(s)!=0) {
// basis not equal to s
- if (n==0) {
+ if (n == 0)
return *this;
- } else {
+ else
return _ex0();
+ } else {
+ // basis equal to s
+ if (is_exactly_of_type(*exponent.bp, numeric) && ex_to_numeric(exponent).is_integer()) {
+ // integer exponent
+ int int_exp = ex_to_numeric(exponent).to_int();
+ if (n == int_exp)
+ return _ex1();
+ else
+ return _ex0();
+ } else {
+ // non-integer exponents are treated as zero
+ if (n == 0)
+ return *this;
+ else
+ return _ex0();
}
- } else if (is_exactly_of_type(*exponent.bp,numeric)&&
- (static_cast<const numeric &>(*exponent.bp).compare(numeric(n))==0)) {
- return _ex1();
}
-
- return _ex0();
}
ex power::eval(int level) const
mulp->clearflag(status_flags::evaluated);
mulp->clearflag(status_flags::hash_calculated);
return (new mul(power(*mulp,exponent),
- power(num_coeff,*num_exponent)))->
- setflag(status_flags::dynallocated);
+ power(num_coeff,*num_exponent)))->setflag(status_flags::dynallocated);
} else {
GINAC_ASSERT(num_coeff.compare(_num0())<0);
if (num_coeff.compare(_num_1())!=0) {
mulp->clearflag(status_flags::evaluated);
mulp->clearflag(status_flags::hash_calculated);
return (new mul(power(*mulp,exponent),
- power(abs(num_coeff),*num_exponent)))->
- setflag(status_flags::dynallocated);
+ power(abs(num_coeff),*num_exponent)))->setflag(status_flags::dynallocated);
}
}
}
return mul(newseq, exponent);
} else {
// D(b^e) = b^e * (D(e)*ln(b) + e*D(b)/b)
- return mul(power(basis, exponent),
- add(mul(exponent.diff(s), log(basis)),
- mul(mul(exponent, basis.diff(s)), power(basis, -1))));
+ return mul(*this,
+ add(mul(exponent.diff(s), log(basis)),
+ mul(mul(exponent, basis.diff(s)), power(basis, _ex_1()))));
}
}
return *this;
ex expanded_basis = basis.expand(options);
-
- if (!is_ex_exactly_of_type(exponent,numeric) ||
- !ex_to_numeric(exponent).is_integer()) {
- if (are_ex_trivially_equal(basis,expanded_basis)) {
+ ex expanded_exponent = exponent.expand(options);
+
+ // x^(a+b) -> x^a * x^b
+ if (is_ex_exactly_of_type(expanded_exponent, add)) {
+ const add &a = ex_to_add(expanded_exponent);
+ exvector distrseq;
+ distrseq.reserve(a.seq.size() + 1);
+ epvector::const_iterator last = a.seq.end();
+ epvector::const_iterator cit = a.seq.begin();
+ while (cit!=last) {
+ distrseq.push_back(power(expanded_basis, a.recombine_pair_to_ex(*cit)));
+ cit++;
+ }
+
+ // Make sure that e.g. (x+y)^(2+a) expands the (x+y)^2 factor
+ if (ex_to_numeric(a.overall_coeff).is_integer()) {
+ const numeric &num_exponent = ex_to_numeric(a.overall_coeff);
+ int int_exponent = num_exponent.to_int();
+ if (int_exponent > 0 && is_ex_exactly_of_type(expanded_basis, add))
+ distrseq.push_back(expand_add(ex_to_add(expanded_basis), int_exponent));
+ else
+ distrseq.push_back(power(expanded_basis, a.overall_coeff));
+ } else
+ distrseq.push_back(power(expanded_basis, a.overall_coeff));
+
+ // Make sure that e.g. (x+y)^(1+a) -> x*(x+y)^a + y*(x+y)^a
+ ex r = (new mul(distrseq))->setflag(status_flags::dynallocated);
+ return r.expand();
+ }
+
+ if (!is_ex_exactly_of_type(expanded_exponent, numeric) ||
+ !ex_to_numeric(expanded_exponent).is_integer()) {
+ if (are_ex_trivially_equal(basis,expanded_basis) && are_ex_trivially_equal(exponent,expanded_exponent)) {
return this->hold();
} else {
- return (new power(expanded_basis,exponent))->
- setflag(status_flags::dynallocated |
- status_flags::expanded);
+ return (new power(expanded_basis,expanded_exponent))->setflag(status_flags::dynallocated | status_flags::expanded);
}
}
// integer numeric exponent
- const numeric & num_exponent = ex_to_numeric(exponent);
+ const numeric & num_exponent = ex_to_numeric(expanded_exponent);
int int_exponent = num_exponent.to_int();
+ // (x+y)^n, n>0
if (int_exponent > 0 && is_ex_exactly_of_type(expanded_basis,add)) {
return expand_add(ex_to_add(expanded_basis), int_exponent);
}
+ // (x*y)^n -> x^n * y^n
if (is_ex_exactly_of_type(expanded_basis,mul)) {
return expand_mul(ex_to_mul(expanded_basis), num_exponent);
}
// cannot expand further
- if (are_ex_trivially_equal(basis,expanded_basis)) {
+ if (are_ex_trivially_equal(basis,expanded_basis) && are_ex_trivially_equal(exponent,expanded_exponent)) {
return this->hold();
} else {
- return (new power(expanded_basis,exponent))->
- setflag(status_flags::dynallocated |
- status_flags::expanded);
+ return (new power(expanded_basis,expanded_exponent))->setflag(status_flags::dynallocated | status_flags::expanded);
}
}
for (l=0; l<m-1; l++) {
const ex & b = a.op(l);
GINAC_ASSERT(!is_ex_exactly_of_type(b,add));
- GINAC_ASSERT(!is_ex_exactly_of_type(b,power)||
- !is_ex_exactly_of_type(ex_to_power(b).exponent,numeric)||
- !ex_to_numeric(ex_to_power(b).exponent).is_pos_integer()||
- !is_ex_exactly_of_type(ex_to_power(b).basis,add)||
- !is_ex_exactly_of_type(ex_to_power(b).basis,mul)||
- !is_ex_exactly_of_type(ex_to_power(b).basis,power));
+ GINAC_ASSERT(!is_ex_exactly_of_type(b,power) ||
+ !is_ex_exactly_of_type(ex_to_power(b).exponent,numeric) ||
+ !ex_to_numeric(ex_to_power(b).exponent).is_pos_integer() ||
+ !is_ex_exactly_of_type(ex_to_power(b).basis,add) ||
+ !is_ex_exactly_of_type(ex_to_power(b).basis,mul) ||
+ !is_ex_exactly_of_type(ex_to_power(b).basis,power));
if (is_ex_exactly_of_type(b,mul)) {
term.push_back(expand_mul(ex_to_mul(b),numeric(k[l])));
} else {
const ex & b = a.op(l);
GINAC_ASSERT(!is_ex_exactly_of_type(b,add));
- GINAC_ASSERT(!is_ex_exactly_of_type(b,power)||
- !is_ex_exactly_of_type(ex_to_power(b).exponent,numeric)||
- !ex_to_numeric(ex_to_power(b).exponent).is_pos_integer()||
- !is_ex_exactly_of_type(ex_to_power(b).basis,add)||
- !is_ex_exactly_of_type(ex_to_power(b).basis,mul)||
- !is_ex_exactly_of_type(ex_to_power(b).basis,power));
+ GINAC_ASSERT(!is_ex_exactly_of_type(b,power) ||
+ !is_ex_exactly_of_type(ex_to_power(b).exponent,numeric) ||
+ !ex_to_numeric(ex_to_power(b).exponent).is_pos_integer() ||
+ !is_ex_exactly_of_type(ex_to_power(b).basis,add) ||
+ !is_ex_exactly_of_type(ex_to_power(b).basis,mul) ||
+ !is_ex_exactly_of_type(ex_to_power(b).basis,power));
if (is_ex_exactly_of_type(b,mul)) {
term.push_back(expand_mul(ex_to_mul(b),numeric(n-k_cum[m-2])));
} else {
const ex & c=(*cit0).coeff;
GINAC_ASSERT(!is_ex_exactly_of_type(r,add));
- GINAC_ASSERT(!is_ex_exactly_of_type(r,power)||
- !is_ex_exactly_of_type(ex_to_power(r).exponent,numeric)||
- !ex_to_numeric(ex_to_power(r).exponent).is_pos_integer()||
- !is_ex_exactly_of_type(ex_to_power(r).basis,add)||
- !is_ex_exactly_of_type(ex_to_power(r).basis,mul)||
- !is_ex_exactly_of_type(ex_to_power(r).basis,power));
+ GINAC_ASSERT(!is_ex_exactly_of_type(r,power) ||
+ !is_ex_exactly_of_type(ex_to_power(r).exponent,numeric) ||
+ !ex_to_numeric(ex_to_power(r).exponent).is_pos_integer() ||
+ !is_ex_exactly_of_type(ex_to_power(r).basis,add) ||
+ !is_ex_exactly_of_type(ex_to_power(r).basis,mul) ||
+ !is_ex_exactly_of_type(ex_to_power(r).basis,power));
if (are_ex_trivially_equal(c,_ex1())) {
if (is_ex_exactly_of_type(r,mul)) {
- sum.push_back(expair(expand_mul(ex_to_mul(r),_num2()),_ex1()));
+ sum.push_back(expair(expand_mul(ex_to_mul(r),_num2()),
+ _ex1()));
} else {
sum.push_back(expair((new power(r,_ex2()))->setflag(status_flags::dynallocated),
- _ex1()));
+ _ex1()));
}
} else {
if (is_ex_exactly_of_type(r,mul)) {
sum.push_back(expair(expand_mul(ex_to_mul(r),_num2()),
- ex_to_numeric(c).power_dyn(_num2())));
+ ex_to_numeric(c).power_dyn(_num2())));
} else {
sum.push_back(expair((new power(r,_ex2()))->setflag(status_flags::dynallocated),
- ex_to_numeric(c).power_dyn(_num2())));
+ ex_to_numeric(c).power_dyn(_num2())));
}
}
const ex & r1=(*cit1).rest;
const ex & c1=(*cit1).coeff;
sum.push_back(a.combine_ex_with_coeff_to_pair((new mul(r,r1))->setflag(status_flags::dynallocated),
- _num2().mul(ex_to_numeric(c)).mul_dyn(ex_to_numeric(c1))));
+ _num2().mul(ex_to_numeric(c)).mul_dyn(ex_to_numeric(c1))));
}
}
GINAC_ASSERT(sum.size()==(a_nops*(a_nops+1))/2);
- return (new add(sum))->setflag(status_flags::dynallocated |
- status_flags::expanded );
+ return (new add(sum))->setflag(status_flags::dynallocated | status_flags::expanded);
}
/** Expand factors of m in m^n where m is a mul and n is and integer
} else {
// it is safe not to call mul::combine_pair_with_coeff_to_pair()
// since n is an integer
- distrseq.push_back(expair((*cit).rest,
- ex_to_numeric((*cit).coeff).mul(n)));
+ distrseq.push_back(expair((*cit).rest, ex_to_numeric((*cit).coeff).mul(n)));
}
++cit;
}
- return (new mul(distrseq,ex_to_numeric(m.overall_coeff).power_dyn(n)))
- ->setflag(status_flags::dynallocated);
+ return (new mul(distrseq,ex_to_numeric(m.overall_coeff).power_dyn(n)))->setflag(status_flags::dynallocated);
}
/*
ex power::expand_commutative_3(const ex & basis, const numeric & exponent,
- unsigned options) const
+ unsigned options) const
{
// obsolete
int n=exponent.to_int();
for (int k=0; k<=n; k++) {
- distrseq.push_back(binomial(n,k)*power(first_operands,numeric(k))*
- power(last_operand,numeric(n-k)));
+ distrseq.push_back(binomial(n,k) * power(first_operands,numeric(k))
+ * power(last_operand,numeric(n-k)));
}
- return ex((new add(distrseq))->setflag(status_flags::expanded |
- status_flags::dynallocated )).
- expand(options);
+ return ex((new add(distrseq))->setflag(status_flags::expanded | status_flags::dynallocated)).expand(options);
}
*/
ex power::expand_noncommutative(const ex & basis, const numeric & exponent,
unsigned options) const
{
- ex rest_power=ex(power(basis,exponent.add(_num_1()))).
- expand(options | expand_options::internal_do_not_expand_power_operands);
+ ex rest_power = ex(power(basis,exponent.add(_num_1()))).
+ expand(options | expand_options::internal_do_not_expand_power_operands);
return ex(mul(rest_power,basis),0).
- expand(options | expand_options::internal_do_not_expand_mul_operands);
+ expand(options | expand_options::internal_do_not_expand_mul_operands);
}
*/
//////////
const power some_power;
-const type_info & typeid_power=typeid(some_power);
+const std::type_info & typeid_power=typeid(some_power);
// helper function