]> www.ginac.de Git - ginac.git/blobdiff - ginac/mul.cpp
Finalize 1.7.6 release.
[ginac.git] / ginac / mul.cpp
index c899ee11f25248cc0edd547f8149bd1a556c512f..16bbc5b229ac9f21c26bcadbfacb1a26c5b6790f 100644 (file)
@@ -3,7 +3,7 @@
  *  Implementation of GiNaC's products of expressions. */
 
 /*
- *  GiNaC Copyright (C) 1999-2015 Johannes Gutenberg University Mainz, Germany
+ *  GiNaC Copyright (C) 1999-2019 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
@@ -89,6 +89,13 @@ mul::mul(const epvector & v, const ex & oc, bool do_index_renaming)
        GINAC_ASSERT(is_canonical());
 }
 
+mul::mul(epvector && vp)
+{
+       overall_coeff = _ex1;
+       construct_from_epvector(std::move(vp));
+       GINAC_ASSERT(is_canonical());
+}
+
 mul::mul(epvector && vp, const ex & oc, bool do_index_renaming)
 {
        overall_coeff = oc;
@@ -285,13 +292,6 @@ bool mul::info(unsigned inf) const
                                return true;
                        return overall_coeff.info(inf);
                }
-               case info_flags::algebraic: {
-                       for (auto & it : seq) {
-                               if (recombine_pair_to_ex(it).info(inf))
-                                       return true;
-                       }
-                       return false;
-               }
                case info_flags::positive:
                case info_flags::negative: {
                        if ((inf==info_flags::positive) && (flags & status_flags::is_positive))
@@ -435,7 +435,7 @@ ex mul::coeff(const ex & s, int n) const
                for (auto & it : seq)
                        coeffseq.push_back(recombine_pair_to_ex(it).coeff(s,n));
                coeffseq.push_back(overall_coeff);
-               return (new mul(coeffseq))->setflag(status_flags::dynallocated);
+               return dynallocate<mul>(coeffseq);
        }
        
        bool coeff_found = false;
@@ -451,7 +451,7 @@ ex mul::coeff(const ex & s, int n) const
        }
        if (coeff_found) {
                coeffseq.push_back(overall_coeff);
-               return (new mul(coeffseq))->setflag(status_flags::dynallocated);
+               return dynallocate<mul>(coeffseq);
        }
        
        return _ex0;
@@ -464,23 +464,21 @@ ex mul::coeff(const ex & s, int n) const
  *  - *(+(x1,x2,...);c) -> *(+(*(x1,c),*(x2,c),...))
  *  - *(x;1) -> x
  *  - *(;c) -> c
- *
- *  @param level cut-off in recursive evaluation */
-ex mul::eval(int level) const
+ */
+ex mul::eval() const
 {
-       epvector evaled = evalchildren(level);
-       if (unlikely(!evaled.empty())) {
-               // do more evaluation later
-               return (new mul(std::move(evaled), overall_coeff))->
-                       setflag(status_flags::dynallocated);
-       }
-
        if (flags & status_flags::evaluated) {
                GINAC_ASSERT(seq.size()>0);
                GINAC_ASSERT(seq.size()>1 || !overall_coeff.is_equal(_ex1));
                return *this;
        }
-       
+
+       const epvector evaled = evalchildren();
+       if (unlikely(!evaled.empty())) {
+               // start over evaluating a new object
+               return dynallocate<mul>(std::move(evaled), overall_coeff);
+       }
+
        size_t seq_size = seq.size();
        if (overall_coeff.is_zero()) {
                // *(...,x;0) -> 0
@@ -501,10 +499,9 @@ ex mul::eval(int level) const
                for (auto & it : addref.seq) {
                        distrseq.push_back(addref.combine_pair_with_coeff_to_pair(it, overall_coeff));
                }
-               return (new add(std::move(distrseq),
-                               ex_to<numeric>(addref.overall_coeff).
-                               mul_dyn(ex_to<numeric>(overall_coeff)))
-                      )->setflag(status_flags::dynallocated | status_flags::evaluated);
+               return dynallocate<add>(std::move(distrseq),
+                                       ex_to<numeric>(addref.overall_coeff).mul_dyn(ex_to<numeric>(overall_coeff)))
+                       .setflag(status_flags::evaluated);
        } else if ((seq_size >= 2) && (! (flags & status_flags::expanded))) {
                // Strip the content and the unit part from each term. Thus
                // things like (-x+a)*(3*x-3*a) automagically turn into - 3*(x-a)^2
@@ -554,14 +551,13 @@ ex mul::eval(int level) const
 
                        // divide add by the number in place to save at least 2 .eval() calls
                        const add& addref = ex_to<add>(i->rest);
-                       add* primitive = new add(addref);
-                       primitive->setflag(status_flags::dynallocated);
-                       primitive->clearflag(status_flags::hash_calculated);
-                       primitive->overall_coeff = ex_to<numeric>(primitive->overall_coeff).div_dyn(c);
-                       for (epvector::iterator ai = primitive->seq.begin(); ai != primitive->seq.end(); ++ai)
-                               ai->coeff = ex_to<numeric>(ai->coeff).div_dyn(c);
-                       
-                       s.push_back(expair(*primitive, _ex1));
+                       add & primitive = dynallocate<add>(addref);
+                       primitive.clearflag(status_flags::hash_calculated);
+                       primitive.overall_coeff = ex_to<numeric>(primitive.overall_coeff).div_dyn(c);
+                       for (auto & ai : primitive.seq)
+                               ai.coeff = ex_to<numeric>(ai.coeff).div_dyn(c);
+
+                       s.push_back(expair(primitive, _ex1));
 
                        ++i;
                        ++j;
@@ -571,31 +567,21 @@ ex mul::eval(int level) const
                                s.push_back(*j);
                                ++j;
                        }
-                       return (new mul(std::move(s), ex_to<numeric>(overall_coeff).mul_dyn(oc))
-                              )->setflag(status_flags::dynallocated);
+                       return dynallocate<mul>(std::move(s), ex_to<numeric>(overall_coeff).mul_dyn(oc));
                }
        }
 
        return this->hold();
 }
 
-ex mul::evalf(int level) const
+ex mul::evalf() const
 {
-       if (level==1)
-               return mul(seq,overall_coeff);
-       
-       if (level==-max_recursion_level)
-               throw(std::runtime_error("max recursion level reached"));
-       
        epvector s;
        s.reserve(seq.size());
 
-       --level;
-       for (auto & it : seq) {
-               s.push_back(combine_ex_with_coeff_to_pair(it.rest.evalf(level),
-                                                         it.coeff));
-       }
-       return mul(std::move(s), overall_coeff.evalf(level));
+       for (auto & it : seq)
+               s.push_back(expair(it.rest.evalf(), it.coeff));
+       return dynallocate<mul>(std::move(s), overall_coeff.evalf());
 }
 
 void mul::find_real_imag(ex & rp, ex & ip) const
@@ -664,11 +650,11 @@ ex mul::evalm() const
                // into that matrix.
                matrix m = ex_to<matrix>(the_matrix->rest);
                s.erase(the_matrix);
-               ex scalar = (new mul(std::move(s), overall_coeff))->setflag(status_flags::dynallocated);
+               ex scalar = dynallocate<mul>(std::move(s), overall_coeff);
                return m.mul_scalar(scalar);
 
        } else
-               return (new mul(std::move(s), overall_coeff))->setflag(status_flags::dynallocated);
+               return dynallocate<mul>(std::move(s), overall_coeff);
 }
 
 ex mul::eval_ncmul(const exvector & v) const
@@ -803,10 +789,10 @@ retry1:
                                        subsed[j] = true;
                        ex subsed_pattern
                                = it.first.subs(repls, subs_options::no_pattern);
-                       divide_by *= power(subsed_pattern, nummatches);
+                       divide_by *= pow(subsed_pattern, nummatches);
                        ex subsed_result
                                = it.second.subs(repls, subs_options::no_pattern);
-                       multiply_by *= power(subsed_result, nummatches);
+                       multiply_by *= pow(subsed_result, nummatches);
                        goto retry1;
 
                } else {
@@ -818,10 +804,10 @@ retry1:
                                        subsed[j] = true;
                                        ex subsed_pattern
                                                = it.first.subs(repls, subs_options::no_pattern);
-                                       divide_by *= power(subsed_pattern, nummatches);
+                                       divide_by *= pow(subsed_pattern, nummatches);
                                        ex subsed_result
                                                = it.second.subs(repls, subs_options::no_pattern);
-                                       multiply_by *= power(subsed_result, nummatches);
+                                       multiply_by *= pow(subsed_result, nummatches);
                                }
                        }
                }
@@ -885,14 +871,14 @@ ex mul::derivative(const symbol & s) const
        auto i = seq.begin(), end = seq.end();
        auto i2 = mulseq.begin();
        while (i != end) {
-               expair ep = split_ex_to_pair(power(i->rest, i->coeff - _ex1) *
+               expair ep = split_ex_to_pair(pow(i->rest, i->coeff - _ex1) *
                                             i->rest.diff(s));
                ep.swap(*i2);
-               addseq.push_back((new mul(mulseq, overall_coeff * i->coeff))->setflag(status_flags::dynallocated));
+               addseq.push_back(dynallocate<mul>(mulseq, overall_coeff * i->coeff));
                ep.swap(*i2);
                ++i; ++i2;
        }
-       return (new add(addseq))->setflag(status_flags::dynallocated);
+       return dynallocate<add>(addseq);
 }
 
 int mul::compare_same_type(const basic & other) const
@@ -949,12 +935,12 @@ return_type_t mul::return_type_tinfo() const
 
 ex mul::thisexpairseq(const epvector & v, const ex & oc, bool do_index_renaming) const
 {
-       return (new mul(v, oc, do_index_renaming))->setflag(status_flags::dynallocated);
+       return dynallocate<mul>(v, oc, do_index_renaming);
 }
 
 ex mul::thisexpairseq(epvector && vp, const ex & oc, bool do_index_renaming) const
 {
-       return (new mul(std::move(vp), oc, do_index_renaming))->setflag(status_flags::dynallocated);
+       return dynallocate<mul>(std::move(vp), oc, do_index_renaming);
 }
 
 expair mul::split_ex_to_pair(const ex & e) const
@@ -970,35 +956,52 @@ expair mul::split_ex_to_pair(const ex & e) const
 expair mul::combine_ex_with_coeff_to_pair(const ex & e,
                                           const ex & c) const
 {
+       GINAC_ASSERT(is_exactly_a<numeric>(c));
+
+       // First, try a common shortcut:
+       if (is_exactly_a<symbol>(e))
+               return expair(e, c);
+
+       // trivial case: exponent 1
+       if (c.is_equal(_ex1))
+               return split_ex_to_pair(e);
+
        // to avoid duplication of power simplification rules,
        // we create a temporary power object
        // otherwise it would be hard to correctly evaluate
        // expression like (4^(1/3))^(3/2)
-       if (c.is_equal(_ex1))
-               return split_ex_to_pair(e);
-
-       return split_ex_to_pair(power(e,c));
+       return split_ex_to_pair(pow(e,c));
 }
 
 expair mul::combine_pair_with_coeff_to_pair(const expair & p,
                                             const ex & c) const
 {
+       GINAC_ASSERT(is_exactly_a<numeric>(p.coeff));
+       GINAC_ASSERT(is_exactly_a<numeric>(c));
+
+       // First, try a common shortcut:
+       if (is_exactly_a<symbol>(p.rest))
+               return expair(p.rest, p.coeff * c);
+
+       // trivial case: exponent 1
+       if (c.is_equal(_ex1))
+               return p;
+       if (p.coeff.is_equal(_ex1))
+               return expair(p.rest, c);
+
        // to avoid duplication of power simplification rules,
        // we create a temporary power object
        // otherwise it would be hard to correctly evaluate
        // expression like (4^(1/3))^(3/2)
-       if (c.is_equal(_ex1))
-               return p;
-
-       return split_ex_to_pair(power(recombine_pair_to_ex(p),c));
+       return split_ex_to_pair(pow(recombine_pair_to_ex(p),c));
 }
 
 ex mul::recombine_pair_to_ex(const expair & p) const
 {
-       if (ex_to<numeric>(p.coeff).is_equal(*_num1_p)) 
+       if (p.coeff.is_equal(_ex1))
                return p.rest;
        else
-               return (new power(p.rest,p.coeff))->setflag(status_flags::dynallocated);
+               return dynallocate<power>(p.rest, p.coeff);
 }
 
 bool mul::expair_needs_further_processing(epp it)
@@ -1068,20 +1071,22 @@ bool mul::can_be_further_expanded(const ex & e)
 
 ex mul::expand(unsigned options) const
 {
-       {
-       // trivial case: expanding the monomial (~ 30% of all calls)
-               epvector::const_iterator i = seq.begin(), seq_end = seq.end();
-               while ((i != seq.end()) &&  is_a<symbol>(i->rest) && i->coeff.info(info_flags::integer))
-                       ++i;
-               if (i == seq_end) {
-                       setflag(status_flags::expanded);
-                       return *this;
+       // Check for trivial case: expanding the monomial (~ 30% of all calls)
+       bool monomial_case = true;
+       for (const auto & i : seq) {
+               if (!is_a<symbol>(i.rest) || !i.coeff.info(info_flags::integer)) {
+                       monomial_case = false;
+                       break;
                }
        }
+       if (monomial_case) {
+               setflag(status_flags::expanded);
+               return *this;
+       }
 
        // do not rename indices if the object has no indices at all
        if ((!(options & expand_options::expand_rename_idx)) && 
-                       this->info(info_flags::has_indices))
+           this->info(info_flags::has_indices))
                options |= expand_options::expand_rename_idx;
 
        const bool skip_idx_rename = !(options & expand_options::expand_rename_idx);
@@ -1134,7 +1139,7 @@ ex mul::expand(unsigned options) const
                                }
 
                                // Compute the new overall coefficient and put it together:
-                               ex tmp_accu = (new add(distrseq, add1.overall_coeff*add2.overall_coeff))->setflag(status_flags::dynallocated);
+                               ex tmp_accu = dynallocate<add>(distrseq, add1.overall_coeff*add2.overall_coeff);
 
                                exvector add1_dummy_indices, add2_dummy_indices, add_indices;
                                lst dummy_subs;
@@ -1168,14 +1173,14 @@ ex mul::expand(unsigned options) const
                                        for (const auto & i1 : add1.seq) {
                                                // Don't push_back expairs which might have a rest that evaluates to a numeric,
                                                // since that would violate an invariant of expairseq:
-                                               const ex rest = (new mul(i1.rest, i2_new))->setflag(status_flags::dynallocated);
+                                               const ex rest = dynallocate<mul>(i1.rest, i2_new);
                                                if (is_exactly_a<numeric>(rest)) {
                                                        oc += ex_to<numeric>(rest).mul(ex_to<numeric>(i1.coeff).mul(ex_to<numeric>(i2.coeff)));
                                                } else {
                                                        distrseq2.push_back(expair(rest, ex_to<numeric>(i1.coeff).mul_dyn(ex_to<numeric>(i2.coeff))));
                                                }
                                        }
-                                       tmp_accu += (new add(distrseq2, oc))->setflag(status_flags::dynallocated);
+                                       tmp_accu += dynallocate<add>(std::move(distrseq2), oc);
                                }
                                last_expanded = tmp_accu;
                        } else {
@@ -1207,7 +1212,7 @@ ex mul::expand(unsigned options) const
                                factors.push_back(split_ex_to_pair(last_expanded.op(i)));
                        else
                                factors.push_back(split_ex_to_pair(rename_dummy_indices_uniquely(va, last_expanded.op(i))));
-                       ex term = (new mul(factors, overall_coeff))->setflag(status_flags::dynallocated);
+                       ex term = dynallocate<mul>(factors, overall_coeff);
                        if (can_be_further_expanded(term)) {
                                distrseq.push_back(term.expand());
                        } else {
@@ -1217,12 +1222,11 @@ ex mul::expand(unsigned options) const
                        }
                }
 
-               return ((new add(distrseq))->
-                       setflag(status_flags::dynallocated | (options == 0 ? status_flags::expanded : 0)));
+               return dynallocate<add>(distrseq).setflag(options == 0 ? status_flags::expanded : 0);
        }
 
        non_adds.push_back(split_ex_to_pair(last_expanded));
-       ex result = (new mul(non_adds, overall_coeff))->setflag(status_flags::dynallocated);
+       ex result = dynallocate<mul>(non_adds, overall_coeff);
        if (can_be_further_expanded(result)) {
                return result.expand();
        } else {