Expansion rules for exp(), log(), and abs().

exp(a + b) -> exp(a)*exp(b)
log(a*b) -> log(a) + log(b)
abs(z*w) -> abs(z)*abs(w)

log and exp are transformed only if expand_options::expand_transcendental
is given.

Signed-off-by: Vladimir V. Kisil <kisilv@maths.leeds.ac.uk>

diff --git a/check/exam_inifcns.cpp b/check/exam_inifcns.cpp
index 30d5b2128bac37df222fa1c8716561818d4dea0c..32368e6030a9d92d747ffccf2c2e84c0a1634ecd 100644 (file)
--- a/check/exam_inifcns.cpp
+++ b/check/exam_inifcns.cpp
@@ -243,6 +243,10 @@ static unsigned inifcns_consist_abs()
        if (!abs(pow(x+I*y,a+I*b)).eval().is_equal(abs(pow(x+I*y,a+I*b))))
                ++result;
 
+       // check expansion of abs
+       if (!abs(-7*z*a*p).expand(expand_options::expand_transcendental).is_equal(7*abs(z)*abs(a)*p))
+               ++result;
+
        if (!abs(z.conjugate()).eval().is_equal(abs(z)))
                ++result;
 
@@ -261,6 +265,65 @@ static unsigned inifcns_consist_abs()
        return result;
 }
 
+static unsigned inifcns_consist_exp()
+{
+       unsigned result = 0;
+       symbol a("a"), b("b");
+
+       if (!exp(a+b).expand(expand_options::expand_transcendental).is_equal(exp(a)*exp(b)))
+               ++result;
+
+       // shall not be expanded since the arg is not add
+       if (!exp(pow(a+b,2)).expand(expand_options::expand_transcendental).is_equal(exp(pow(a+b,2))))
+               ++result;
+
+       // expand now
+       if (!exp(pow(a+b,2)).expand(expand_options::expand_function_args | expand_options::expand_transcendental)
+               .is_equal(exp(a*a)*exp(b*b)*exp(2*a*b)))
+               ++result;
+
+       return result;
+}
+
+static unsigned inifcns_consist_log()
+{
+       unsigned result = 0;
+       symbol z("a"), w("b");
+       realsymbol a("a"), b("b");
+       possymbol p("p"), q("q");
+
+       // do not expand
+       if (!log(z*w).expand(expand_options::expand_transcendental).is_equal(log(z*w)))
+               ++result;
+
+       // do not expand
+       if (!log(a*b).expand(expand_options::expand_transcendental).is_equal(log(a*b)))
+               ++result;
+
+       // shall expand
+       if (!log(p*q).expand(expand_options::expand_transcendental).is_equal(log(p) + log(q)))
+               ++result;
+
+       // a bit more complicated
+       ex e1 = log(-7*p*pow(q,3)*a*pow(b,2)*z*w).expand(expand_options::expand_transcendental);
+       ex e2 = log(7)+log(p)+log(pow(q,3))+log(-z*a*w*pow(b,2));
+       if (!e1.is_equal(e2))
+               ++result;
+
+       if (!ex(log(pow(p,a))).is_equal(a*log(p)))
+               ++result;
+
+       // shall not do for non-real powers
+       if (ex(log(pow(p,z))).is_equal(z*log(p)))
+               ++result;
+
+       // shall not do for non-positive basis
+       if (ex(log(pow(a,b))).is_equal(b*log(a)))
+               ++result;
+
+       return result;
+}
+
 static unsigned inifcns_consist_various()
 {
        unsigned result = 0;
@@ -285,6 +348,8 @@ unsigned exam_inifcns()
        result += inifcns_consist_psi();  cout << '.' << flush;
        result += inifcns_consist_zeta();  cout << '.' << flush;
        result += inifcns_consist_abs();  cout << '.' << flush;
+       result += inifcns_consist_exp();  cout << '.' << flush;
+       result += inifcns_consist_log();  cout << '.' << flush;
        result += inifcns_consist_various();  cout << '.' << flush;
        
        return result;

diff --git a/doc/tutorial/ginac.texi b/doc/tutorial/ginac.texi
index 732487061da055bb2522796a5c95f91d96c34b43..d5bea846342325744ea73283293dc31c2815760d 100644 (file)
--- a/doc/tutorial/ginac.texi
+++ b/doc/tutorial/ginac.texi
@@ -5926,6 +5926,59 @@ towards -infinity.  These functions are continuous as the branch cut
 is approached coming around the finite endpoint of the cut in a
 counter clockwise direction.
 
+@c
+@subsection Expanding functions
+@cindex expand trancedent functions
+@cindex @code{expand_options::expand_transcendental}
+@cindex @code{expand_options::expand_function_args}
+GiNaC knows several expansion laws for trancedent functions, e.g.
+@tex
+$e^{a+b}=e^a e^b$,
+$|zw|=|z|\cdot |w|$
+@end tex
+@ifnottex
+@command{exp(a+b)=exp(a) exp(b), |zw|=|z| |w|}
+@end ifnottex
+or
+@tex
+$\log(c*d)=\log(c)+\log(d)$,
+@end tex
+@ifnottex
+@command{log(cd)=log(c)+log(d)}
+@end ifnottex
+(for positive
+@tex
+$c,\ d$
+@end tex
+@ifnottex
+@command{c, d}
+@end ifnottex
+). In order to use these rules you need to call @code{expand()} method
+with the option @code{expand_options::expand_transcendental}. Another
+relevant option is @code{expand_options::expand_function_args}. Their
+usage and interaction can be seen from the following example:
+@example
+@{
+       symbol x("x"),  y("y");
+       ex e=exp(pow(x+y,2));
+       cout << e.expand() << endl;
+       // -> exp((x+y)^2)
+       cout << e.expand(expand_options::expand_transcendental) << endl;
+       // -> exp((x+y)^2)
+       cout << e.expand(expand_options::expand_function_args) << endl;
+       // -> exp(2*x*y+x^2+y^2)
+       cout << e.expand(expand_options::expand_function_args
+                       | expand_options::expand_transcendental) << endl;
+       // -> exp(y^2)*exp(2*x*y)*exp(x^2)
+@}
+@end example
+If both flags are set (as in the last call), then GiNaC tries to get
+the maximal expansion. For example, for the exponent GiNaC firstly expands
+the argument and then the function. For the logarithm and absolute value,
+GiNaC uses the opposite order: firstly expands the function and then its
+argument. Of course, a user can fine-tune this behaviour by sequential
+calls of several @code{expand()} methods with desired flags.
+
 @node Multiple polylogarithms, Complex expressions, Built-in functions, Methods and functions
 @c    node-name, next, previous, up
 @subsection Multiple polylogarithms

diff --git a/ginac/inifcns.cpp b/ginac/inifcns.cpp
index 1144ce1997e14eed7911a6b6004a18fcc008b1a4..153a3c1e92e6417a6664fbed399d6492ec2af506 100644 (file)
--- a/ginac/inifcns.cpp
+++ b/ginac/inifcns.cpp
@@ -254,6 +254,27 @@ static ex abs_eval(const ex & arg)
        return abs(arg).hold();
 }
 
+static ex abs_expand(const ex & arg, unsigned options)
+{
+       if ((options & expand_options::expand_transcendental)
+               && is_exactly_a<mul>(arg)) {
+               exvector prodseq;
+               prodseq.reserve(arg.nops());
+               for (const_iterator i = arg.begin(); i != arg.end(); ++i) {
+                       if (options & expand_options::expand_function_args)
+                               prodseq.push_back(abs(i->expand(options)));
+                       else
+                               prodseq.push_back(abs(*i));
+               }
+               return (new mul(prodseq))->setflag(status_flags::dynallocated | status_flags::expanded);
+       }
+
+       if (options & expand_options::expand_function_args)
+               return abs(arg.expand(options)).hold();
+       else
+               return abs(arg).hold();
+}
+
 static void abs_print_latex(const ex & arg, const print_context & c)
 {
        c.s << "{|"; arg.print(c); c.s << "|}";
@@ -317,6 +338,7 @@ bool abs_info(const ex & arg, unsigned inf)
 
 REGISTER_FUNCTION(abs, eval_func(abs_eval).
                        evalf_func(abs_evalf).
+                       expand_func(abs_expand).
                        info_func(abs_info).
                        print_func<print_latex>(abs_print_latex).
                        print_func<print_csrc_float>(abs_print_csrc_float).

diff --git a/ginac/inifcns_trans.cpp b/ginac/inifcns_trans.cpp
index 10a36758f43001f6c02774531522c7fd2667c7c2..6ee0a23b0c388d1146123cb2c4654536f0442b3e 100644 (file)
--- a/ginac/inifcns_trans.cpp
+++ b/ginac/inifcns_trans.cpp
@@ -24,6 +24,8 @@
 #include "inifcns.h"
 #include "ex.h"
 #include "constant.h"
+#include "add.h"
+#include "mul.h"
 #include "numeric.h"
 #include "power.h"
 #include "operators.h"
@@ -81,6 +83,27 @@ static ex exp_eval(const ex & x)
        return exp(x).hold();
 }
 
+static ex exp_expand(const ex & arg, unsigned options)
+{
+       ex exp_arg;
+       if (options & expand_options::expand_function_args)
+               exp_arg = arg.expand(options);
+       else
+               exp_arg=arg;
+
+       if ((options & expand_options::expand_transcendental)
+               && is_exactly_a<add>(exp_arg)) {
+               exvector prodseq;
+               prodseq.reserve(exp_arg.nops());
+               for (const_iterator i = exp_arg.begin(); i != exp_arg.end(); ++i)
+                       prodseq.push_back(exp(*i));
+
+               return (new mul(prodseq))->setflag(status_flags::dynallocated | status_flags::expanded);
+       }
+
+       return exp(exp_arg).hold();
+}
+
 static ex exp_deriv(const ex & x, unsigned deriv_param)
 {
        GINAC_ASSERT(deriv_param==0);
@@ -107,6 +130,7 @@ static ex exp_conjugate(const ex & x)
 
 REGISTER_FUNCTION(exp, eval_func(exp_eval).
                        evalf_func(exp_evalf).
+                       expand_func(exp_expand).
                        derivative_func(exp_deriv).
                        real_part_func(exp_real_part).
                        imag_part_func(exp_imag_part).
@@ -265,6 +289,53 @@ static ex log_imag_part(const ex & x)
        return atan2(GiNaC::imag_part(x), GiNaC::real_part(x));
 }
 
+static ex log_expand(const ex & arg, unsigned options)
+{
+       if ((options & expand_options::expand_transcendental)
+               && is_exactly_a<mul>(arg) && !arg.info(info_flags::indefinite)) {
+               exvector sumseq;
+               exvector prodseq;
+               sumseq.reserve(arg.nops());
+               prodseq.reserve(arg.nops());
+               bool possign=true;
+
+               // searching for positive/negative factors
+               for (const_iterator i = arg.begin(); i != arg.end(); ++i) {
+                       ex e;
+                       if (options & expand_options::expand_function_args)
+                               e=i->expand(options);
+                       else
+                               e=*i;
+                       if (e.info(info_flags::positive))
+                               sumseq.push_back(log(e));
+                       else if (e.info(info_flags::negative)) {
+                               sumseq.push_back(log(-e));
+                               possign = !possign;
+                       } else
+                               prodseq.push_back(e);
+               }
+
+               if (sumseq.size() > 0) {
+                       ex newarg;
+                       if (options & expand_options::expand_function_args)
+                               newarg=((possign?_ex1:_ex_1)*mul(prodseq)).expand(options);
+                       else {
+                               newarg=(possign?_ex1:_ex_1)*mul(prodseq);
+                               ex_to<basic>(newarg).setflag(status_flags::purely_indefinite);
+                       }
+                       return add(sumseq)+log(newarg);
+               } else {
+                       if (!(options & expand_options::expand_function_args))
+                               ex_to<basic>(arg).setflag(status_flags::purely_indefinite);
+               }
+       }
+
+       if (options & expand_options::expand_function_args)
+               return log(arg.expand(options)).hold();
+       else
+               return log(arg).hold();
+}
+
 static ex log_conjugate(const ex & x)
 {
        // conjugate(log(x))==log(conjugate(x)) unless on the branch cut which
@@ -281,6 +352,7 @@ static ex log_conjugate(const ex & x)
 
 REGISTER_FUNCTION(log, eval_func(log_eval).
                        evalf_func(log_evalf).
+                       expand_func(log_expand).
                        derivative_func(log_deriv).
                        series_func(log_series).
                        real_part_func(log_real_part).