X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=blobdiff_plain;f=ginac%2Fnormal.cpp;h=bcb2a398bd8d57b8b8a983c2350963270e3ef3f5;hp=2a251a780704fdfb66ed96dc0befcf93567dc84b;hb=019faaa693e82ccccadacb185ddd999437409ba2;hpb=1602530f716ba1d425a0667b897182b99c374823 diff --git a/ginac/normal.cpp b/ginac/normal.cpp index 2a251a78..bcb2a398 100644 --- a/ginac/normal.cpp +++ b/ginac/normal.cpp @@ -6,7 +6,7 @@ * computation, square-free factorization and rational function normalization. */ /* - * GiNaC Copyright (C) 1999-2009 Johannes Gutenberg University Mainz, Germany + * GiNaC Copyright (C) 1999-2015 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 @@ -120,6 +120,11 @@ static bool get_first_symbol(const ex &e, ex &x) * * @see get_symbol_stats */ struct sym_desc { + /** Initialize symbol, leave other variables uninitialized */ + sym_desc(const ex& s) + : sym(s), deg_a(0), deg_b(0), ldeg_a(0), ldeg_b(0), max_deg(0), max_lcnops(0) + { } + /** Reference to symbol */ ex sym; @@ -157,15 +162,11 @@ typedef std::vector sym_desc_vec; // Add symbol the sym_desc_vec (used internally by get_symbol_stats()) static void add_symbol(const ex &s, sym_desc_vec &v) { - sym_desc_vec::const_iterator it = v.begin(), itend = v.end(); - while (it != itend) { - if (it->sym.is_equal(s)) // If it's already in there, don't add it a second time + for (auto & it : v) + if (it.sym.is_equal(s)) // If it's already in there, don't add it a second time return; - ++it; - } - sym_desc d; - d.sym = s; - v.push_back(d); + + v.push_back(sym_desc(s)); } // Collect all symbols of an expression (used internally by get_symbol_stats()) @@ -197,17 +198,15 @@ static void get_symbol_stats(const ex &a, const ex &b, sym_desc_vec &v) { collect_symbols(a.eval(), v); // eval() to expand assigned symbols collect_symbols(b.eval(), v); - sym_desc_vec::iterator it = v.begin(), itend = v.end(); - while (it != itend) { - int deg_a = a.degree(it->sym); - int deg_b = b.degree(it->sym); - it->deg_a = deg_a; - it->deg_b = deg_b; - it->max_deg = std::max(deg_a, deg_b); - it->max_lcnops = std::max(a.lcoeff(it->sym).nops(), b.lcoeff(it->sym).nops()); - it->ldeg_a = a.ldegree(it->sym); - it->ldeg_b = b.ldegree(it->sym); - ++it; + for (auto & it : v) { + int deg_a = a.degree(it.sym); + int deg_b = b.degree(it.sym); + it.deg_a = deg_a; + it.deg_b = deg_b; + it.max_deg = std::max(deg_a, deg_b); + it.max_lcnops = std::max(a.lcoeff(it.sym).nops(), b.lcoeff(it.sym).nops()); + it.ldeg_a = a.ldegree(it.sym); + it.ldeg_b = b.ldegree(it.sym); } std::sort(v.begin(), v.end()); @@ -321,15 +320,12 @@ numeric numeric::integer_content() const numeric add::integer_content() const { - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); numeric c = *_num0_p, l = *_num1_p; - while (it != itend) { - GINAC_ASSERT(!is_exactly_a(it->rest)); - GINAC_ASSERT(is_exactly_a(it->coeff)); - c = gcd(ex_to(it->coeff).numer(), c); - l = lcm(ex_to(it->coeff).denom(), l); - it++; + for (auto & it : seq) { + GINAC_ASSERT(!is_exactly_a(it.rest)); + GINAC_ASSERT(is_exactly_a(it.coeff)); + c = gcd(ex_to(it.coeff).numer(), c); + l = lcm(ex_to(it.coeff).denom(), l); } GINAC_ASSERT(is_exactly_a(overall_coeff)); c = gcd(ex_to(overall_coeff).numer(), c); @@ -340,11 +336,8 @@ numeric add::integer_content() const numeric mul::integer_content() const { #ifdef DO_GINAC_ASSERT - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); - while (it != itend) { - GINAC_ASSERT(!is_exactly_a(recombine_pair_to_ex(*it))); - ++it; + for (auto & it : seq) { + GINAC_ASSERT(!is_exactly_a(recombine_pair_to_ex(it))); } #endif // def DO_GINAC_ASSERT GINAC_ASSERT(is_exactly_a(overall_coeff)); @@ -796,10 +789,10 @@ static bool divide_in_z(const ex &a, const ex &b, ex &q, sym_desc_vec::const_ite if (is_exactly_a(b)) { ex qbar = a; - for (const_iterator itrb = b.begin(); itrb != b.end(); ++itrb) { + for (const auto & it : b) { sym_desc_vec sym_stats; - get_symbol_stats(a, *itrb, sym_stats); - if (!divide_in_z(qbar, *itrb, q, sym_stats.begin())) + get_symbol_stats(a, it, sym_stats); + if (!divide_in_z(qbar, it, q, sym_stats.begin())) return false; qbar = q; @@ -961,7 +954,7 @@ ex ex::content(const ex &x) const return lcoeff * c / lcoeff.unit(x); ex cont = _ex0; for (int i=ldeg; i<=deg; i++) - cont = gcd(r.coeff(x, i), cont, NULL, NULL, false); + cont = gcd(r.coeff(x, i), cont, nullptr, nullptr, false); return cont * c; } @@ -1101,7 +1094,7 @@ static ex sr_gcd(const ex &a, const ex &b, sym_desc_vec::const_iterator var) // Remove content from c and d, to be attached to GCD later ex cont_c = c.content(x); ex cont_d = d.content(x); - ex gamma = gcd(cont_c, cont_d, NULL, NULL, false); + ex gamma = gcd(cont_c, cont_d, nullptr, nullptr, false); if (ddeg == 0) return gamma; c = c.primpart(x, cont_c); @@ -1165,17 +1158,14 @@ numeric numeric::max_coefficient() const numeric add::max_coefficient() const { - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); GINAC_ASSERT(is_exactly_a(overall_coeff)); numeric cur_max = abs(ex_to(overall_coeff)); - while (it != itend) { + for (auto & it : seq) { numeric a; - GINAC_ASSERT(!is_exactly_a(it->rest)); - a = abs(ex_to(it->coeff)); + GINAC_ASSERT(!is_exactly_a(it.rest)); + a = abs(ex_to(it.coeff)); if (a > cur_max) cur_max = a; - it++; } return cur_max; } @@ -1183,11 +1173,8 @@ numeric add::max_coefficient() const numeric mul::max_coefficient() const { #ifdef DO_GINAC_ASSERT - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); - while (it != itend) { - GINAC_ASSERT(!is_exactly_a(recombine_pair_to_ex(*it))); - it++; + for (auto & it : seq) { + GINAC_ASSERT(!is_exactly_a(recombine_pair_to_ex(it))); } #endif // def DO_GINAC_ASSERT GINAC_ASSERT(is_exactly_a(overall_coeff)); @@ -1215,28 +1202,22 @@ ex add::smod(const numeric &xi) const { epvector newseq; newseq.reserve(seq.size()+1); - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); - while (it != itend) { - GINAC_ASSERT(!is_exactly_a(it->rest)); - numeric coeff = GiNaC::smod(ex_to(it->coeff), xi); + for (auto & it : seq) { + GINAC_ASSERT(!is_exactly_a(it.rest)); + numeric coeff = GiNaC::smod(ex_to(it.coeff), xi); if (!coeff.is_zero()) - newseq.push_back(expair(it->rest, coeff)); - it++; + newseq.push_back(expair(it.rest, coeff)); } GINAC_ASSERT(is_exactly_a(overall_coeff)); numeric coeff = GiNaC::smod(ex_to(overall_coeff), xi); - return (new add(newseq,coeff))->setflag(status_flags::dynallocated); + return (new add(std::move(newseq), coeff))->setflag(status_flags::dynallocated); } ex mul::smod(const numeric &xi) const { #ifdef DO_GINAC_ASSERT - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); - while (it != itend) { - GINAC_ASSERT(!is_exactly_a(recombine_pair_to_ex(*it))); - it++; + for (auto & it : seq) { + GINAC_ASSERT(!is_exactly_a(recombine_pair_to_ex(it))); } #endif // def DO_GINAC_ASSERT mul * mulcopyp = new mul(*this); @@ -1272,9 +1253,9 @@ class gcdheu_failed {}; * * @param a first integer multivariate polynomial (expanded) * @param b second integer multivariate polynomial (expanded) - * @param ca cofactor of polynomial a (returned), NULL to suppress + * @param ca cofactor of polynomial a (returned), nullptr to suppress * calculation of cofactor - * @param cb cofactor of polynomial b (returned), NULL to suppress + * @param cb cofactor of polynomial b (returned), nullptr to suppress * calculation of cofactor * @param var iterator to first element of vector of sym_desc structs * @param res the GCD (returned) @@ -1365,9 +1346,9 @@ static bool heur_gcd_z(ex& res, const ex &a, const ex &b, ex *ca, ex *cb, * * @param a first rational multivariate polynomial (expanded) * @param b second rational multivariate polynomial (expanded) - * @param ca cofactor of polynomial a (returned), NULL to suppress + * @param ca cofactor of polynomial a (returned), nullptr to suppress * calculation of cofactor - * @param cb cofactor of polynomial b (returned), NULL to suppress + * @param cb cofactor of polynomial b (returned), nullptr to suppress * calculation of cofactor * @param var iterator to first element of vector of sym_desc structs * @param res the GCD (returned) @@ -1431,8 +1412,8 @@ static ex gcd_pf_mul(const ex& a, const ex& b, ex* ca, ex* cb); * * @param a first multivariate polynomial * @param b second multivariate polynomial - * @param ca pointer to expression that will receive the cofactor of a, or NULL - * @param cb pointer to expression that will receive the cofactor of b, or NULL + * @param ca pointer to expression that will receive the cofactor of a, or nullptr + * @param cb pointer to expression that will receive the cofactor of b, or nullptr * @param check_args check whether a and b are polynomials with rational * coefficients (defaults to "true") * @return the GCD as a new expression */ @@ -1535,7 +1516,7 @@ ex gcd(const ex &a, const ex &b, ex *ca, ex *cb, bool check_args, unsigned optio if (ca) *ca = ex_to(aex)/g; if (cb) - *cb = bex/g; + *cb = bex/g; return g; } @@ -1863,11 +1844,8 @@ ex sqrfree(const ex &a, const lst &l) if (l.nops()==0) { sym_desc_vec sdv; get_symbol_stats(a, _ex0, sdv); - sym_desc_vec::const_iterator it = sdv.begin(), itend = sdv.end(); - while (it != itend) { - args.append(it->sym); - ++it; - } + for (auto & it : sdv) + args.append(it.sym); } else { args = l; } @@ -1890,18 +1868,15 @@ ex sqrfree(const ex &a, const lst &l) // recurse down the factors in remaining variables if (newargs.nops()>0) { - exvector::iterator i = factors.begin(); - while (i != factors.end()) { - *i = sqrfree(*i, newargs); - ++i; - } + for (auto & it : factors) + it = sqrfree(it, newargs); } // Done with recursion, now construct the final result ex result = _ex1; - exvector::const_iterator it = factors.begin(), itend = factors.end(); - for (int p = 1; it!=itend; ++it, ++p) - result *= power(*it, p); + int p = 1; + for (auto & it : factors) + result *= power(it, p++); // Yun's algorithm does not account for constant factors. (For univariate // polynomials it works only in the monic case.) We can correct this by @@ -1912,7 +1887,7 @@ ex sqrfree(const ex &a, const lst &l) else result *= quo(tmp, result, x); - // Put in the reational overall factor again and return + // Put in the rational overall factor again and return return result * lcm.inverse(); } @@ -2006,16 +1981,18 @@ ex sqrfree_parfrac(const ex & a, const symbol & x) * @see ex::normal */ static ex replace_with_symbol(const ex & e, exmap & repl, exmap & rev_lookup) { + // Since the repl contains replaced expressions we should search for them + ex e_replaced = e.subs(repl, subs_options::no_pattern); + // Expression already replaced? Then return the assigned symbol - exmap::const_iterator it = rev_lookup.find(e); + auto it = rev_lookup.find(e_replaced); if (it != rev_lookup.end()) return it->second; - + // Otherwise create new symbol and add to list, taking care that the // replacement expression doesn't itself contain symbols from repl, // because subs() is not recursive ex es = (new symbol)->setflag(status_flags::dynallocated); - ex e_replaced = e.subs(repl, subs_options::no_pattern); repl.insert(std::make_pair(es, e_replaced)); rev_lookup.insert(std::make_pair(e_replaced, es)); return es; @@ -2028,16 +2005,18 @@ static ex replace_with_symbol(const ex & e, exmap & repl, exmap & rev_lookup) * @see basic::to_polynomial */ static ex replace_with_symbol(const ex & e, exmap & repl) { + // Since the repl contains replaced expressions we should search for them + ex e_replaced = e.subs(repl, subs_options::no_pattern); + // Expression already replaced? Then return the assigned symbol - for (exmap::const_iterator it = repl.begin(); it != repl.end(); ++it) - if (it->second.is_equal(e)) - return it->first; - + for (auto & it : repl) + if (it.second.is_equal(e_replaced)) + return it.first; + // Otherwise create new symbol and add to list, taking care that the // replacement expression doesn't itself contain symbols from repl, // because subs() is not recursive ex es = (new symbol)->setflag(status_flags::dynallocated); - ex e_replaced = e.subs(repl, subs_options::no_pattern); repl.insert(std::make_pair(es, e_replaced)); return es; } @@ -2047,7 +2026,7 @@ static ex replace_with_symbol(const ex & e, exmap & repl) struct normal_map_function : public map_function { int level; normal_map_function(int l) : level(l) {} - ex operator()(const ex & e) { return normal(e, level); } + ex operator()(const ex & e) override { return normal(e, level); } }; /** Default implementation of ex::normal(). It normalizes the children and @@ -2177,12 +2156,10 @@ ex add::normal(exmap & repl, exmap & rev_lookup, int level) const exvector nums, dens; nums.reserve(seq.size()+1); dens.reserve(seq.size()+1); - epvector::const_iterator it = seq.begin(), itend = seq.end(); - while (it != itend) { - ex n = ex_to(recombine_pair_to_ex(*it)).normal(repl, rev_lookup, level-1); + for (auto & it : seq) { + ex n = ex_to(recombine_pair_to_ex(it)).normal(repl, rev_lookup, level-1); nums.push_back(n.op(0)); dens.push_back(n.op(1)); - it++; } ex n = ex_to(overall_coeff).normal(repl, rev_lookup, level-1); nums.push_back(n.op(0)); @@ -2194,8 +2171,8 @@ ex add::normal(exmap & repl, exmap & rev_lookup, int level) const //std::clog << "add::normal uses " << nums.size() << " summands:\n"; // Add fractions sequentially - exvector::const_iterator num_it = nums.begin(), num_itend = nums.end(); - exvector::const_iterator den_it = dens.begin(), den_itend = dens.end(); + auto num_it = nums.begin(), num_itend = nums.end(); + auto den_it = dens.begin(), den_itend = dens.end(); //std::clog << " num = " << *num_it << ", den = " << *den_it << std::endl; ex num = *num_it++, den = *den_it++; while (num_it != num_itend) { @@ -2208,7 +2185,7 @@ ex add::normal(exmap & repl, exmap & rev_lookup, int level) const num_it++; den_it++; } - // Additiion of two fractions, taking advantage of the fact that + // Addition of two fractions, taking advantage of the fact that // the heuristic GCD algorithm computes the cofactors at no extra cost ex co_den1, co_den2; ex g = gcd(den, next_den, &co_den1, &co_den2, false); @@ -2236,12 +2213,10 @@ ex mul::normal(exmap & repl, exmap & rev_lookup, int level) const exvector num; num.reserve(seq.size()); exvector den; den.reserve(seq.size()); ex n; - epvector::const_iterator it = seq.begin(), itend = seq.end(); - while (it != itend) { - n = ex_to(recombine_pair_to_ex(*it)).normal(repl, rev_lookup, level-1); + for (auto & it : seq) { + n = ex_to(recombine_pair_to_ex(it)).normal(repl, rev_lookup, level-1); num.push_back(n.op(0)); den.push_back(n.op(1)); - it++; } n = ex_to(overall_coeff).normal(repl, rev_lookup, level-1); num.push_back(n.op(0)); @@ -2315,14 +2290,12 @@ ex power::normal(exmap & repl, exmap & rev_lookup, int level) const ex pseries::normal(exmap & repl, exmap & rev_lookup, int level) const { epvector newseq; - epvector::const_iterator i = seq.begin(), end = seq.end(); - while (i != end) { - ex restexp = i->rest.normal(); + for (auto & it : seq) { + ex restexp = it.rest.normal(); if (!restexp.is_zero()) - newseq.push_back(expair(restexp, i->coeff)); - ++i; + newseq.push_back(expair(restexp, it.coeff)); } - ex n = pseries(relational(var,point), newseq); + ex n = pseries(relational(var,point), std::move(newseq)); return (new lst(replace_with_symbol(n, repl, rev_lookup), _ex1))->setflag(status_flags::dynallocated); } @@ -2394,7 +2367,7 @@ ex ex::denom() const return e.op(1).subs(repl, subs_options::no_pattern); } -/** Get numerator and denominator of an expression. If the expresison is not +/** Get numerator and denominator of an expression. If the expression is not * of the normal form "numerator/denominator", it is first converted to this * form and then a list [numerator, denominator] is returned. * @@ -2438,15 +2411,15 @@ ex ex::to_rational(lst & repl_lst) const { // Convert lst to exmap exmap m; - for (lst::const_iterator it = repl_lst.begin(); it != repl_lst.end(); ++it) - m.insert(std::make_pair(it->op(0), it->op(1))); + for (auto & it : repl_lst) + m.insert(std::make_pair(it.op(0), it.op(1))); ex ret = bp->to_rational(m); // Convert exmap back to lst repl_lst.remove_all(); - for (exmap::const_iterator it = m.begin(); it != m.end(); ++it) - repl_lst.append(it->first == it->second); + for (auto & it : m) + repl_lst.append(it.first == it.second); return ret; } @@ -2461,15 +2434,15 @@ ex ex::to_polynomial(lst & repl_lst) const { // Convert lst to exmap exmap m; - for (lst::const_iterator it = repl_lst.begin(); it != repl_lst.end(); ++it) - m.insert(std::make_pair(it->op(0), it->op(1))); + for (auto & it : repl_lst) + m.insert(std::make_pair(it.op(0), it.op(1))); ex ret = bp->to_polynomial(m); // Convert exmap back to lst repl_lst.remove_all(); - for (exmap::const_iterator it = m.begin(); it != m.end(); ++it) - repl_lst.append(it->first == it->second); + for (auto & it : m) + repl_lst.append(it.first == it.second); return ret; } @@ -2576,17 +2549,15 @@ ex expairseq::to_rational(exmap & repl) const { epvector s; s.reserve(seq.size()); - epvector::const_iterator i = seq.begin(), end = seq.end(); - while (i != end) { - s.push_back(split_ex_to_pair(recombine_pair_to_ex(*i).to_rational(repl))); - ++i; - } + for (auto & it : seq) + s.push_back(split_ex_to_pair(recombine_pair_to_ex(it).to_rational(repl))); + ex oc = overall_coeff.to_rational(repl); if (oc.info(info_flags::numeric)) - return thisexpairseq(s, overall_coeff); + return thisexpairseq(std::move(s), overall_coeff); else s.push_back(combine_ex_with_coeff_to_pair(oc, _ex1)); - return thisexpairseq(s, default_overall_coeff()); + return thisexpairseq(std::move(s), default_overall_coeff()); } /** Implementation of ex::to_polynomial() for expairseqs. */ @@ -2594,17 +2565,15 @@ ex expairseq::to_polynomial(exmap & repl) const { epvector s; s.reserve(seq.size()); - epvector::const_iterator i = seq.begin(), end = seq.end(); - while (i != end) { - s.push_back(split_ex_to_pair(recombine_pair_to_ex(*i).to_polynomial(repl))); - ++i; - } + for (auto & it : seq) + s.push_back(split_ex_to_pair(recombine_pair_to_ex(it).to_polynomial(repl))); + ex oc = overall_coeff.to_polynomial(repl); if (oc.info(info_flags::numeric)) - return thisexpairseq(s, overall_coeff); + return thisexpairseq(std::move(s), overall_coeff); else s.push_back(combine_ex_with_coeff_to_pair(oc, _ex1)); - return thisexpairseq(s, default_overall_coeff()); + return thisexpairseq(std::move(s), default_overall_coeff()); }