X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?a=blobdiff_plain;ds=sidebyside;f=ginac%2Fpower.cpp;h=75fd2ff51a9d12bd09d867d30fdf50b860ec6240;hb=6c4ddbcabae766fc993d40de53cbbcdbd758d61d;hp=69e3dd7fc968c4545381f8890c95b4cec6badf23;hpb=5a8b8e3c4d882249db35b679ce3144a59a7012e8;p=ginac.git

diff --git a/ginac/power.cpp b/ginac/power.cpp
index 69e3dd7f..75fd2ff5 100644
--- a/ginac/power.cpp
+++ b/ginac/power.cpp
@@ -40,6 +40,7 @@
 #include "lst.h"
 #include "archive.h"
 #include "utils.h"
+#include "relational.h"
 
 namespace GiNaC {
 
@@ -258,6 +259,15 @@ ex power::map(map_function & f) const
 		return *this;
 }
 
+bool power::is_polynomial(const ex & var) const
+{
+	if (exponent.has(var))
+		return false;
+	if (!exponent.info(info_flags::nonnegint))
+		return false;
+	return basis.is_polynomial(var);
+}
+
 int power::degree(const ex & s) const
 {
 	if (is_equal(ex_to<basic>(s)))
@@ -605,13 +615,66 @@ ex power::conjugate() const
 	return (new power(newbasis, newexponent))->setflag(status_flags::dynallocated);
 }
 
+ex power::real_part() const
+{
+	if (exponent.info(info_flags::integer)) {
+		ex basis_real = basis.real_part();
+		if (basis_real == basis)
+			return *this;
+		realsymbol a("a"),b("b");
+		ex result;
+		if (exponent.info(info_flags::posint))
+			result = power(a+I*b,exponent);
+		else
+			result = power(a/(a*a+b*b)-I*b/(a*a+b*b),-exponent);
+		result = result.expand();
+		result = result.real_part();
+		result = result.subs(lst( a==basis_real, b==basis.imag_part() ));
+		return result;
+	}
+	
+	ex a = basis.real_part();
+	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))*cos(c*atan2(b,a)+d*log(abs(basis)));
+}
+
+ex power::imag_part() const
+{
+	if (exponent.info(info_flags::integer)) {
+		ex basis_real = basis.real_part();
+		if (basis_real == basis)
+			return 0;
+		realsymbol a("a"),b("b");
+		ex result;
+		if (exponent.info(info_flags::posint))
+			result = power(a+I*b,exponent);
+		else
+			result = power(a/(a*a+b*b)-I*b/(a*a+b*b),-exponent);
+		result = result.expand();
+		result = result.imag_part();
+		result = result.subs(lst( a==basis_real, b==basis.imag_part() ));
+		return result;
+	}
+	
+	ex a=basis.real_part();
+	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)));
+}
+
+// protected
+
 // protected
 
 /** Implementation of ex::diff() for a power.
  *  @see ex::diff */
 ex power::derivative(const symbol & s) const
 {
-	if (exponent.info(info_flags::real)) {
+	if (is_a<numeric>(exponent)) {
 		// D(b^r) = r * b^(r-1) * D(b) (faster than the formula below)
 		epvector newseq;
 		newseq.reserve(2);
@@ -643,7 +706,7 @@ unsigned power::return_type() const
 	return basis.return_type();
 }
 
-const basic* power::return_type_tinfo() const
+tinfo_t power::return_type_tinfo() const
 {
 	return basis.return_type_tinfo();
 }