* Implementation of GiNaC's symbolic exponentiation (basis^exponent). */
/*
- * GiNaC Copyright (C) 1999-2008 Johannes Gutenberg University Mainz, Germany
+ * GiNaC Copyright (C) 1999-2011 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
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
*/
-#include <vector>
-#include <iostream>
-#include <stdexcept>
-#include <limits>
-
#include "power.h"
#include "expairseq.h"
#include "add.h"
#include "relational.h"
#include "compiler.h"
+#include <iostream>
+#include <limits>
+#include <stdexcept>
+#include <vector>
+
namespace GiNaC {
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(power, basic,
// default constructor
//////////
-power::power() : inherited(&power::tinfo_static) { }
+power::power() { }
//////////
// other constructors
// archiving
//////////
-power::power(const archive_node &n, lst &sym_lst) : inherited(n, sym_lst)
+void power::read_archive(const archive_node &n, lst &sym_lst)
{
+ inherited::read_archive(n, sym_lst);
n.find_ex("basis", basis, sym_lst);
n.find_ex("exponent", exponent, sym_lst);
}
n.add_ex("exponent", exponent);
}
-DEFAULT_UNARCHIVE(power)
-
//////////
// functions overriding virtual functions from base classes
//////////
basis.info(inf);
case info_flags::expanded:
return (flags & status_flags::expanded);
+ case info_flags::positive:
+ return basis.info(info_flags::positive) && exponent.info(info_flags::real);
case info_flags::has_indices: {
if (flags & status_flags::has_indices)
return true;
const ex & ebasis = level==1 ? basis : basis.eval(level-1);
const ex & eexponent = level==1 ? exponent : exponent.eval(level-1);
- bool basis_is_numerical = false;
- bool exponent_is_numerical = false;
- const numeric *num_basis;
- const numeric *num_exponent;
+ const numeric *num_basis = NULL;
+ const numeric *num_exponent = NULL;
if (is_exactly_a<numeric>(ebasis)) {
- basis_is_numerical = true;
num_basis = &ex_to<numeric>(ebasis);
}
if (is_exactly_a<numeric>(eexponent)) {
- exponent_is_numerical = true;
num_exponent = &ex_to<numeric>(eexponent);
}
return ebasis;
// ^(0,c1) -> 0 or exception (depending on real value of c1)
- if (ebasis.is_zero() && exponent_is_numerical) {
+ if ( ebasis.is_zero() && num_exponent ) {
if ((num_exponent->real()).is_zero())
throw (std::domain_error("power::eval(): pow(0,I) is undefined"));
else if ((num_exponent->real()).is_negative())
if (is_exactly_a<power>(ebasis) && ebasis.op(0).info(info_flags::positive) && ebasis.op(1).info(info_flags::real))
return power(ebasis.op(0), ebasis.op(1) * eexponent);
- if (exponent_is_numerical) {
+ if ( num_exponent ) {
// ^(c1,c2) -> c1^c2 (c1, c2 numeric(),
// except if c1,c2 are rational, but c1^c2 is not)
- if (basis_is_numerical) {
+ if ( num_basis ) {
const bool basis_is_crational = num_basis->is_crational();
const bool exponent_is_crational = num_exponent->is_crational();
if (!basis_is_crational || !exponent_is_crational) {
if (num_coeff.is_positive()) {
mul *mulp = new mul(mulref);
mulp->overall_coeff = _ex1;
+ mulp->setflag(status_flags::dynallocated);
mulp->clearflag(status_flags::evaluated);
mulp->clearflag(status_flags::hash_calculated);
return (new mul(power(*mulp,exponent),
if (!num_coeff.is_equal(*_num_1_p)) {
mul *mulp = new mul(mulref);
mulp->overall_coeff = _ex_1;
+ mulp->setflag(status_flags::dynallocated);
mulp->clearflag(status_flags::evaluated);
mulp->clearflag(status_flags::hash_calculated);
return (new mul(power(*mulp,exponent),
ex power::conjugate() const
{
- ex newbasis = basis.conjugate();
- ex newexponent = exponent.conjugate();
- if (are_ex_trivially_equal(basis, newbasis) && are_ex_trivially_equal(exponent, newexponent)) {
- return *this;
+ // conjugate(pow(x,y))==pow(conjugate(x),conjugate(y)) unless on the
+ // branch cut which runs along the negative real axis.
+ if (basis.info(info_flags::positive)) {
+ ex newexponent = exponent.conjugate();
+ if (are_ex_trivially_equal(exponent, newexponent)) {
+ return *this;
+ }
+ return (new power(basis, newexponent))->setflag(status_flags::dynallocated);
}
- return (new power(newbasis, newexponent))->setflag(status_flags::dynallocated);
+ if (exponent.info(info_flags::integer)) {
+ ex newbasis = basis.conjugate();
+ if (are_ex_trivially_equal(basis, newbasis)) {
+ return *this;
+ }
+ return (new power(newbasis, exponent))->setflag(status_flags::dynallocated);
+ }
+ return conjugate_function(*this).hold();
}
ex power::real_part() const
ex b=basis.imag_part();
ex c=exponent.real_part();
ex d=exponent.imag_part();
- return
- power(abs(basis),c)*exp(-d*atan2(b,a))*sin(c*atan2(b,a)+d*log(abs(basis)));
+ return power(abs(basis),c)*exp(-d*atan2(b,a))*sin(c*atan2(b,a)+d*log(abs(basis)));
}
// protected
return result;
}
+GINAC_BIND_UNARCHIVER(power);
+
} // namespace GiNaC