Add evaluation of Order(x)^e -> Order(x^e) for positive integer e.

diff --git a/ginac/inifcns.cpp b/ginac/inifcns.cpp
index f8ec39a3035aabdf64e4c6516946de407d7c56f8..2202794dd33d784a0ab24b73a5c47e45089be60d 100644 (file)
--- a/ginac/inifcns.cpp
+++ b/ginac/inifcns.cpp
@@ -1021,6 +1021,20 @@ static ex Order_imag_part(const ex & x)
        return Order(x).hold();
 }
 
+static ex Order_power(const ex & x, const ex & e)
+{
+       // Order(x)^e -> Order(x^e) for positive integer e
+       if (is_exactly_a<numeric>(e) && e.info(info_flags::posint))
+               return Order(pow(x, e));
+       // NB: For negative exponents, the above could be wrong.
+       // This is because series() produces Order(x^n) to denote the order where
+       // it gave up. So, Order(x^n) can also be an x^(n+1) term if the x^n term
+       // vanishes. In this situation, 1/Order(x^n) can also be a x^(-n-1) term.
+       // Transforming it to Order(x^-n) would miss that.
+
+       return power(Order(x), e).hold();
+}
+
 static ex Order_expl_derivative(const ex & arg, const symbol & s)
 {
        return Order(arg.diff(s));
@@ -1032,7 +1046,8 @@ REGISTER_FUNCTION(Order, eval_func(Order_eval).
                          expl_derivative_func(Order_expl_derivative).
                          conjugate_func(Order_conjugate).
                          real_part_func(Order_real_part).
-                         imag_part_func(Order_imag_part));
+                         imag_part_func(Order_imag_part).
+                         power_func(Order_power));
 
 //////////
 // Solve linear system