]> www.ginac.de Git - ginac.git/blobdiff - ginac/inifcns_trans.cpp
[BUGFIX] Fix crash in parser.
[ginac.git] / ginac / inifcns_trans.cpp
index bf44d69833606e26fe2fb160342317f4ebe5e414..9a2a0b6c402bc6f666637ccadbd9c1331ec66b95 100644 (file)
@@ -4,7 +4,7 @@
  *  functions. */
 
 /*
- *  GiNaC Copyright (C) 1999-2014 Johannes Gutenberg University Mainz, Germany
+ *  GiNaC Copyright (C) 1999-2024 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
@@ -98,7 +98,7 @@ static ex exp_expand(const ex & arg, unsigned options)
                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 dynallocate<mul>(prodseq).setflag(status_flags::expanded);
        }
 
        return exp(exp_arg).hold();
@@ -128,13 +128,47 @@ static ex exp_conjugate(const ex & x)
        return exp(x.conjugate());
 }
 
+static ex exp_power(const ex & x, const ex & a)
+{
+       /*
+        * The power law (e^x)^a=e^(x*a) is used in two cases:
+        * a) a is an integer and x may be complex;
+        * b) both x and a are reals.
+        * Negative a is excluded to keep automatic simplifications like exp(x)/exp(x)=1.
+        */
+       if (a.info(info_flags::nonnegative)
+           && (a.info(info_flags::integer) || (x.info(info_flags::real) && a.info(info_flags::real))))
+               return exp(x*a);
+       else if (a.info(info_flags::negative)
+                && (a.info(info_flags::integer) || (x.info(info_flags::real) && a.info(info_flags::real))))
+               return power(exp(-x*a), _ex_1).hold();
+
+       return power(exp(x), a).hold();
+}
+
+static bool exp_info(const ex & x, unsigned inf)
+{
+       switch (inf) {
+       case info_flags::expanded:
+       case info_flags::real:
+               return x.info(inf);
+       case info_flags::positive:
+       case info_flags::nonnegative:
+               return x.info(info_flags::real);
+       default:
+               return false;
+       }
+}
+
 REGISTER_FUNCTION(exp, eval_func(exp_eval).
                        evalf_func(exp_evalf).
+                       info_func(exp_info).
                        expand_func(exp_expand).
                        derivative_func(exp_deriv).
                        real_part_func(exp_real_part).
                        imag_part_func(exp_imag_part).
                        conjugate_func(exp_conjugate).
+                       power_func(exp_power).
                        latex_name("\\exp"));
 
 //////////
@@ -197,7 +231,7 @@ static ex log_series(const ex &arg,
        // maybe substitution of rel into arg fails because of a pole
        try {
                arg_pt = arg.subs(rel, subs_options::no_pattern);
-       } catch (pole_error) {
+       } catch (pole_error &) {
                must_expand_arg = true;
        }
        // or we are at the branch point anyways
@@ -242,25 +276,19 @@ static ex log_series(const ex &arg,
                        // in this case n more (or less) terms are needed
                        // (sadly, to generate them, we have to start from the beginning)
                        if (n == 0 && coeff == 1) {
-                               epvector epv;
-                               ex acc = (new pseries(rel, epv))->setflag(status_flags::dynallocated);
-                               epv.reserve(2);
-                               epv.push_back(expair(-1, _ex0));
-                               epv.push_back(expair(Order(_ex1), order));
-                               ex rest = pseries(rel, epv).add_series(argser);
+                               ex rest = pseries(rel, epvector{expair(-1, _ex0), expair(Order(_ex1), order)}).add_series(argser);
+                               ex acc = dynallocate<pseries>(rel, epvector());
                                for (int i = order-1; i>0; --i) {
-                                       epvector cterm;
-                                       cterm.reserve(1);
-                                       cterm.push_back(expair(i%2 ? _ex1/i : _ex_1/i, _ex0));
-                                       acc = pseries(rel, cterm).add_series(ex_to<pseries>(acc));
+                                       epvector cterm { expair(i%2 ? _ex1/i : _ex_1/i, _ex0) };
+                                       acc = pseries(rel, std::move(cterm)).add_series(ex_to<pseries>(acc));
                                        acc = (ex_to<pseries>(rest)).mul_series(ex_to<pseries>(acc));
                                }
                                return acc;
                        }
                        const ex newarg = ex_to<pseries>((arg/coeff).series(rel, order+n, options)).shift_exponents(-n).convert_to_poly(true);
-                       return pseries(rel, seq).add_series(ex_to<pseries>(log(newarg).series(rel, order, options)));
+                       return pseries(rel, std::move(seq)).add_series(ex_to<pseries>(log(newarg).series(rel, order, options)));
                } else  // it was a monomial
-                       return pseries(rel, seq);
+                       return pseries(rel, std::move(seq));
        }
        if (!(options & series_options::suppress_branchcut) &&
             arg_pt.info(info_flags::negative)) {
@@ -272,9 +300,12 @@ static ex log_series(const ex &arg,
                const symbol foo;
                const ex replarg = series(log(arg), s==foo, order).subs(foo==point, subs_options::no_pattern);
                epvector seq;
-               seq.push_back(expair(-I*csgn(arg*I)*Pi, _ex0));
+               if (order > 0) {
+                       seq.reserve(2);
+                       seq.push_back(expair(-I*csgn(arg*I)*Pi, _ex0));
+               }
                seq.push_back(expair(Order(_ex1), order));
-               return series(replarg - I*Pi + pseries(rel, seq), rel, order);
+               return series(replarg - I*Pi + pseries(rel, std::move(seq)), rel, order);
        }
        throw do_taylor();  // caught by function::series()
 }
@@ -354,8 +385,21 @@ static ex log_conjugate(const ex & x)
        return conjugate_function(log(x)).hold();
 }
 
+static bool log_info(const ex & x, unsigned inf)
+{
+       switch (inf) {
+       case info_flags::expanded:
+               return x.info(inf);
+       case info_flags::real:
+               return x.info(info_flags::positive);
+       default:
+               return false;
+       }
+}
+
 REGISTER_FUNCTION(log, eval_func(log_eval).
                        evalf_func(log_evalf).
+                       info_func(log_info).
                        expand_func(log_expand).
                        derivative_func(log_deriv).
                        series_func(log_series).
@@ -463,8 +507,20 @@ static ex sin_conjugate(const ex & x)
        return sin(x.conjugate());
 }
 
+static bool trig_info(const ex & x, unsigned inf)
+{
+       switch (inf) {
+       case info_flags::expanded:
+       case info_flags::real:
+               return x.info(inf);
+       default:
+               return false;
+       }
+}
+
 REGISTER_FUNCTION(sin, eval_func(sin_eval).
                        evalf_func(sin_evalf).
+                       info_func(trig_info).
                        derivative_func(sin_deriv).
                        real_part_func(sin_real_part).
                        imag_part_func(sin_imag_part).
@@ -571,6 +627,7 @@ static ex cos_conjugate(const ex & x)
 }
 
 REGISTER_FUNCTION(cos, eval_func(cos_eval).
+                       info_func(trig_info).
                        evalf_func(cos_evalf).
                        derivative_func(cos_deriv).
                        real_part_func(cos_real_part).
@@ -696,6 +753,7 @@ static ex tan_conjugate(const ex & x)
 
 REGISTER_FUNCTION(tan, eval_func(tan_eval).
                        evalf_func(tan_evalf).
+                       info_func(trig_info).
                        derivative_func(tan_deriv).
                        series_func(tan_series).
                        real_part_func(tan_real_part).
@@ -770,8 +828,19 @@ static ex asin_conjugate(const ex & x)
        return conjugate_function(asin(x)).hold();
 }
 
+static bool asin_info(const ex & x, unsigned inf)
+{
+       switch (inf) {
+       case info_flags::expanded:
+               return x.info(inf);
+       default:
+               return false;
+       }
+}
+
 REGISTER_FUNCTION(asin, eval_func(asin_eval).
                         evalf_func(asin_evalf).
+                        info_func(asin_info).
                         derivative_func(asin_deriv).
                         conjugate_func(asin_conjugate).
                         latex_name("\\arcsin"));
@@ -845,6 +914,7 @@ static ex acos_conjugate(const ex & x)
 
 REGISTER_FUNCTION(acos, eval_func(acos_eval).
                         evalf_func(acos_evalf).
+                        info_func(asin_info). // Flags of acos are shared with asin functions
                         derivative_func(acos_deriv).
                         conjugate_func(acos_conjugate).
                         latex_name("\\arccos"));
@@ -936,9 +1006,12 @@ static ex atan_series(const ex &arg,
                else
                        Order0correction += log((I*arg_pt+_ex1)/(I*arg_pt+_ex_1))*I*_ex1_2;
                epvector seq;
-               seq.push_back(expair(Order0correction, _ex0));
+               if (order > 0) {
+                       seq.reserve(2);
+                       seq.push_back(expair(Order0correction, _ex0));
+               }
                seq.push_back(expair(Order(_ex1), order));
-               return series(replarg - pseries(rel, seq), rel, order);
+               return series(replarg - pseries(rel, std::move(seq)), rel, order);
        }
        throw do_taylor();
 }
@@ -959,8 +1032,24 @@ static ex atan_conjugate(const ex & x)
        return conjugate_function(atan(x)).hold();
 }
 
+static bool atan_info(const ex & x, unsigned inf)
+{
+       switch (inf) {
+       case info_flags::expanded:
+       case info_flags::real:
+               return x.info(inf);
+       case info_flags::positive:
+       case info_flags::negative:
+       case info_flags::nonnegative:
+               return x.info(info_flags::real) && x.info(inf);
+       default:
+               return false;
+       }
+}
+
 REGISTER_FUNCTION(atan, eval_func(atan_eval).
                         evalf_func(atan_evalf).
+                        info_func(atan_info).
                         derivative_func(atan_deriv).
                         series_func(atan_series).
                         conjugate_func(atan_conjugate).
@@ -1061,7 +1150,25 @@ static ex atan2_deriv(const ex & y, const ex & x, unsigned deriv_param)
        return -y*power(power(x,_ex2)+power(y,_ex2),_ex_1);
 }
 
+static bool atan2_info(const ex & y, const ex & x, unsigned inf)
+{
+       switch (inf) {
+       case info_flags::expanded:
+       case info_flags::real:
+               return y.info(inf) && x.info(inf);
+       case info_flags::positive:
+       case info_flags::negative:
+       case info_flags::nonnegative:
+               return y.info(info_flags::real) && x.info(info_flags::real)
+                       && y.info(inf);
+       default:
+               return false;
+       }
+}
+
 REGISTER_FUNCTION(atan2, eval_func(atan2_eval).
+                         evalf_func(atan2_evalf).
+                         info_func(atan2_info).
                          evalf_func(atan2_evalf).
                          derivative_func(atan2_deriv));
 
@@ -1143,6 +1250,7 @@ static ex sinh_conjugate(const ex & x)
 
 REGISTER_FUNCTION(sinh, eval_func(sinh_eval).
                         evalf_func(sinh_evalf).
+                        info_func(atan_info). // Flags of sinh are shared with atan functions
                         derivative_func(sinh_deriv).
                         real_part_func(sinh_real_part).
                         imag_part_func(sinh_imag_part).
@@ -1227,6 +1335,7 @@ static ex cosh_conjugate(const ex & x)
 
 REGISTER_FUNCTION(cosh, eval_func(cosh_eval).
                         evalf_func(cosh_evalf).
+                        info_func(exp_info). // Flags of cosh are shared with exp functions
                         derivative_func(cosh_deriv).
                         real_part_func(cosh_real_part).
                         imag_part_func(cosh_imag_part).
@@ -1331,6 +1440,7 @@ static ex tanh_conjugate(const ex & x)
 
 REGISTER_FUNCTION(tanh, eval_func(tanh_eval).
                         evalf_func(tanh_evalf).
+                        info_func(atan_info). // Flags of tanh are shared with atan functions
                         derivative_func(tanh_deriv).
                         series_func(tanh_series).
                         real_part_func(tanh_real_part).
@@ -1396,6 +1506,7 @@ static ex asinh_conjugate(const ex & x)
 
 REGISTER_FUNCTION(asinh, eval_func(asinh_eval).
                          evalf_func(asinh_evalf).
+                         info_func(atan_info). // Flags of asinh are shared with atan functions
                          derivative_func(asinh_deriv).
                          conjugate_func(asinh_conjugate));
 
@@ -1460,6 +1571,7 @@ static ex acosh_conjugate(const ex & x)
 
 REGISTER_FUNCTION(acosh, eval_func(acosh_eval).
                          evalf_func(acosh_evalf).
+                         info_func(asin_info). // Flags of acosh are shared with asin functions
                          derivative_func(acosh_deriv).
                          conjugate_func(acosh_conjugate));
 
@@ -1530,22 +1642,25 @@ static ex atanh_series(const ex &arg,
                return ((log(_ex1+arg)-log(_ex1-arg))*_ex1_2).series(rel, order, options);
        // ...and the branch cuts (the discontinuity at the cut being just I*Pi)
        if (!(options & series_options::suppress_branchcut)) {
-               // method:
-               // This is the branch cut: assemble the primitive series manually and
-               // then add the corresponding complex step function.
-               const symbol &s = ex_to<symbol>(rel.lhs());
-               const ex &point = rel.rhs();
-               const symbol foo;
-               const ex replarg = series(atanh(arg), s==foo, order).subs(foo==point, subs_options::no_pattern);
+               // method:
+               // This is the branch cut: assemble the primitive series manually and
+               // then add the corresponding complex step function.
+               const symbol &s = ex_to<symbol>(rel.lhs());
+               const ex &point = rel.rhs();
+               const symbol foo;
+               const ex replarg = series(atanh(arg), s==foo, order).subs(foo==point, subs_options::no_pattern);
                ex Order0correction = replarg.op(0)+csgn(I*arg)*Pi*I*_ex1_2;
                if (arg_pt<_ex0)
                        Order0correction += log((arg_pt+_ex_1)/(arg_pt+_ex1))*_ex1_2;
                else
                        Order0correction += log((arg_pt+_ex1)/(arg_pt+_ex_1))*_ex_1_2;
-               epvector seq;
-               seq.push_back(expair(Order0correction, _ex0));
-               seq.push_back(expair(Order(_ex1), order));
-               return series(replarg - pseries(rel, seq), rel, order);
+               epvector seq;
+               if (order > 0) {
+                       seq.reserve(2);
+                       seq.push_back(expair(Order0correction, _ex0));
+               }
+               seq.push_back(expair(Order(_ex1), order));
+               return series(replarg - pseries(rel, std::move(seq)), rel, order);
        }
        throw do_taylor();
 }
@@ -1563,6 +1678,7 @@ static ex atanh_conjugate(const ex & x)
 
 REGISTER_FUNCTION(atanh, eval_func(atanh_eval).
                          evalf_func(atanh_evalf).
+                         info_func(asin_info). // Flags of atanh are shared with asin functions
                          derivative_func(atanh_deriv).
                          series_func(atanh_series).
                          conjugate_func(atanh_conjugate));