From 7b27cc8b0e4da74e995b7fd59ce7e778b63fb8cb Mon Sep 17 00:00:00 2001 From: Christian Bauer Date: Mon, 19 May 2003 18:21:43 +0000 Subject: [PATCH] lst now provides (read-only) iterator types, the use of which makes iterating over the elements of a lst an operation of order O(N), instead of O(N^2) when using nops()/op(); this should speed up subs() --- NEWS | 4 +++ configure.ac | 2 +- doc/tutorial/ginac.texi | 45 ++++++++++++++++++++++++++++++ ginac/basic.cpp | 27 ++++++++++-------- ginac/clifford.cpp | 5 ++-- ginac/container.pl | 11 ++++++++ ginac/expairseq.cpp | 10 +++---- ginac/idx.cpp | 11 ++++---- ginac/indexed.cpp | 11 ++------ ginac/matrix.cpp | 62 ++++++++++++++++++++++++++--------------- ginac/mul.cpp | 15 +++++----- ginac/ncmul.cpp | 2 +- ginac/normal.cpp | 13 +++++---- ginac/power.cpp | 8 ++++-- ginac/symbol.cpp | 6 ++-- ginac/symmetry.cpp | 15 ++-------- 16 files changed, 160 insertions(+), 87 deletions(-) diff --git a/NEWS b/NEWS index 34b333b6..59ba650a 100644 --- a/NEWS +++ b/NEWS @@ -1,5 +1,9 @@ This file records noteworthy changes. +1.1.1 () +* lst (and exprseq) provide iterators for read-only element access. For + sequential access this is one order faster than using op(). + 1.1.0 (3 April 2003) * Removed deprecated macros is_ex_a, is_ex_exactly_a and friends for good. * The scalar_products mechanism allows the specification of an index dimension. diff --git a/configure.ac b/configure.ac index 46d636ef..763348d0 100644 --- a/configure.ac +++ b/configure.ac @@ -20,7 +20,7 @@ dnl (don't we all *love* M4?)... GINACLIB_MAJOR_VERSION=1 GINACLIB_MINOR_VERSION=1 -GINACLIB_MICRO_VERSION=0 +GINACLIB_MICRO_VERSION=1 GINACLIB_INTERFACE_AGE=0 GINACLIB_BINARY_AGE=0 GINACLIB_VERSION=$GINACLIB_MAJOR_VERSION.$GINACLIB_MINOR_VERSION.$GINACLIB_MICRO_VERSION diff --git a/doc/tutorial/ginac.texi b/doc/tutorial/ginac.texi index 04cedf43..358ea911 100644 --- a/doc/tutorial/ginac.texi +++ b/doc/tutorial/ginac.texi @@ -1304,6 +1304,51 @@ individual elements: ... @end example +As with the standard @code{list} container, accessing random elements of a +@code{lst} is generally an operation of order @math{O(N)}. Faster read-only +sequential access to the elements of a list is possible with the +iterator types provided by the @code{lst} class: + +@example +typedef ... lst::const_iterator; +typedef ... lst::const_reverse_iterator; +lst::const_iterator lst::begin() const; +lst::const_iterator lst::end() const; +lst::const_reverse_iterator lst::rbegin() const; +lst::const_reverse_iterator lst::rend() const; +@end example + +For example, to print the elements of a list individually you can use: + +@example + ... + // O(N) + for (lst::const_iterator i = l.begin(); i != l.end(); ++i) + cout << *i << endl; + ... +@end example + +which is one order faster than + +@example + ... + // O(N^2) + for (size_t i = 0; i < l.nops(); ++i) + cout << l.op(i) << endl; + ... +@end example + +These iterators also allow you to use some of the algorithms provided by +the C++ standard library: + +@example + ... + // sum up the elements of the list (requires #include ) + ex sum = accumulate(l.begin(), l.end(), ex(0)); + cout << sum << endl; // prints '2+2*x+2*y' + ... +@end example + @code{lst} is one of the few GiNaC classes that allow in-place modifications (the only other one is @code{matrix}). You can modify single elements: diff --git a/ginac/basic.cpp b/ginac/basic.cpp index 85294cac..6fb3952f 100644 --- a/ginac/basic.cpp +++ b/ginac/basic.cpp @@ -484,9 +484,9 @@ bool basic::match(const ex & pattern, lst & repl_lst) const // Wildcard matches anything, but check whether we already have found // a match for that wildcard first (if so, it the earlier match must // be the same expression) - for (size_t i=0; i(repl_lst.op(i).op(1))); + for (lst::const_iterator it = repl_lst.begin(); it != repl_lst.end(); ++it) { + if (it->op(0).is_equal(pattern)) + return is_equal(ex_to(it->op(1))); } repl_lst.append(pattern == *this); return true; @@ -526,16 +526,18 @@ ex basic::subs(const lst & ls, const lst & lr, unsigned options) const { GINAC_ASSERT(ls.nops() == lr.nops()); + lst::const_iterator its, itr; + if (options & subs_options::subs_no_pattern) { - for (size_t i=0; i(ls.op(i)))) - return lr.op(i); + for (its = ls.begin(), itr = lr.begin(); its != ls.end(); ++its, ++itr) { + if (is_equal(ex_to(*its))) + return *itr; } } else { - for (size_t i=0; i(ls.op(i)), repl_lst)) - return lr.op(i).subs(repl_lst, options | subs_options::subs_no_pattern); // avoid infinite recursion when re-substituting the wildcards + if (match(ex_to(*its), repl_lst)) + return itr->subs(repl_lst, options | subs_options::subs_no_pattern); // avoid infinite recursion when re-substituting the wildcards } } @@ -711,10 +713,13 @@ ex basic::subs(const ex & e, unsigned options) const if (!e.info(info_flags::list)) { throw(std::invalid_argument("basic::subs(ex): argument must be a list")); } + + // Split list into two lst ls; lst lr; - for (size_t i=0; i(e)); + for (lst::const_iterator it = ex_to(e).begin(); it != ex_to(e).end(); ++it) { + ex r = *it; if (!r.info(info_flags::relation_equal)) { throw(std::invalid_argument("basic::subs(ex): argument must be a list of equations")); } diff --git a/ginac/clifford.cpp b/ginac/clifford.cpp index a8b87343..8ea924b1 100644 --- a/ginac/clifford.cpp +++ b/ginac/clifford.cpp @@ -695,8 +695,9 @@ ex canonicalize_clifford(const ex & e) ex aux = e.to_rational(srl); for (size_t i=0; i(rhs) && rhs.return_type() == return_types::noncommutative diff --git a/ginac/container.pl b/ginac/container.pl index 7dcd0be6..66bef262 100755 --- a/ginac/container.pl +++ b/ginac/container.pl @@ -244,6 +244,10 @@ class ${CONTAINER} : public basic { GINAC_DECLARE_REGISTERED_CLASS(${CONTAINER}, basic) +public: + typedef ${STLT}::const_iterator const_iterator; + typedef ${STLT}::const_reverse_iterator const_reverse_iterator; + public: ${CONTAINER}(${STLT} const & s, bool discardable = false); ${CONTAINER}(${STLT} * vp); // vp will be deleted @@ -276,6 +280,13 @@ protected: virtual ex this${CONTAINER}(${STLT} const & v) const; virtual ex this${CONTAINER}(${STLT} * vp) const; + // non-virtual functions in this class +public: + const_iterator begin() const {return seq.begin();} + const_iterator end() const {return seq.end();} + const_reverse_iterator rbegin() const {return seq.rbegin();} + const_reverse_iterator rend() const {return seq.rend();} + protected: bool is_canonical() const; ${STLT} evalchildren(int level) const; diff --git a/ginac/expairseq.cpp b/ginac/expairseq.cpp index c41c694b..041994af 100644 --- a/ginac/expairseq.cpp +++ b/ginac/expairseq.cpp @@ -379,9 +379,9 @@ found: ; for (size_t i=0; ipush_back(split_ex_to_pair(ops[i])); ex rest = thisexpairseq(vp, default_overall_coeff()); - for (size_t i=0; iop(0).is_equal(global_wildcard)) + return rest.is_equal(it->op(1)); } repl_lst.append(global_wildcard == rest); return true; @@ -1564,8 +1564,8 @@ epvector * expairseq::subschildren(const lst &ls, const lst &lr, unsigned option // is a product or power. In this case we have to recombine the pairs // because the numeric coefficients may be part of the search pattern. bool complex_subs = false; - for (size_t i=0; i(ls.op(i)) || is_exactly_a(ls.op(i))) { + for (lst::const_iterator it = ls.begin(); it != ls.end(); ++it) { + if (is_exactly_a(*it) || is_exactly_a(*it)) { complex_subs = true; break; } diff --git a/ginac/idx.cpp b/ginac/idx.cpp index 55475568..35952871 100644 --- a/ginac/idx.cpp +++ b/ginac/idx.cpp @@ -376,16 +376,17 @@ ex idx::subs(const lst & ls, const lst & lr, unsigned options) const GINAC_ASSERT(ls.nops() == lr.nops()); // First look for index substitutions - for (size_t i=0; i(ls.op(i)))) { + lst::const_iterator its, itr; + for (its = ls.begin(), itr = lr.begin(); its != ls.end(); ++its, ++itr) { + if (is_equal(ex_to(*its))) { // Substitution index->index - if (is_a(lr.op(i))) - return lr.op(i); + if (is_a(*itr)) + return *itr; // Otherwise substitute value idx *i_copy = static_cast(duplicate()); - i_copy->value = lr.op(i); + i_copy->value = *itr; i_copy->clearflag(status_flags::hash_calculated); return i_copy->setflag(status_flags::dynallocated); } diff --git a/ginac/indexed.cpp b/ginac/indexed.cpp index 695768f6..de2e8c96 100644 --- a/ginac/indexed.cpp +++ b/ginac/indexed.cpp @@ -1224,14 +1224,9 @@ void scalar_products::add(const ex & v1, const ex & v2, const ex & dim, const ex void scalar_products::add_vectors(const lst & l, const ex & dim) { // Add all possible pairs of products - size_t num = l.nops(); - for (size_t i=0; i= r) break; // matrix smaller than list: throw away excessive elements - m[y*c+x] = l.op(i); + m[y*c+x] = *it; } } @@ -880,6 +882,7 @@ ex matrix::charpoly(const symbol & lambda) const // trapped and we use Leverrier's algorithm which goes as row^3 for // every coefficient. The expensive part is the matrix multiplication. if (numeric_flag) { + matrix B(*this); ex c = B.trace(); ex poly = power(lambda,row)-c*power(lambda,row-1); @@ -887,20 +890,22 @@ ex matrix::charpoly(const symbol & lambda) const for (unsigned j=0; jmul(B); - c = B.trace()/ex(i+1); + c = B.trace() / ex(i+1); poly -= c*power(lambda,row-i-1); } if (row%2) return -poly; else return poly; - } + + } else { - matrix M(*this); - for (unsigned r=0; r cols) - cols = l.op(i).nops(); + size_t rows = l.nops(), cols = 0; + for (itr = l.begin(); itr != l.end(); ++itr) { + if (!is_a(*itr)) + throw (std::invalid_argument("lst_to_matrix: argument must be a list of lists")); + if (itr->nops() > cols) + cols = itr->nops(); + } // Allocate and fill matrix matrix &M = *new matrix(rows, cols); M.setflag(status_flags::dynallocated); - for (i=0; i j) - M(i, j) = l.op(i).op(j); - else - M(i, j) = _ex0; + + unsigned i; + for (itr = l.begin(), i = 0; itr != l.end(); ++itr, ++i) { + unsigned j; + for (itc = ex_to(*itr).begin(), j = 0; itc != ex_to(*itr).end(); ++itc, ++j) + M(i, j) = *itc; + } + return M; } ex diag_matrix(const lst & l) { + lst::const_iterator it; size_t dim = l.nops(); - matrix &m = *new matrix(dim, dim); - m.setflag(status_flags::dynallocated); - for (size_t i=0; i subsed(seq.size(), false); exvector subsresult(seq.size()); - for (size_t i=0; i(ls.op(i))) { + if (is_exactly_a(*its)) { int nummatches = std::numeric_limits::max(); std::vector currsubsed(seq.size(), false); bool succeed = true; lst repls; - for (size_t j=0; jnops(); j++) { bool found=false; for (size_t k=0; kop(j), nummatches, repls)) { currsubsed[k] = true; found = true; break; @@ -601,7 +602,7 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con subsresult[j] = op(j); else { foundfirstsubsedfactor = true; - subsresult[j] = op(j) * power(lr.op(i).subs(ex(repls), subs_options::subs_no_pattern) / ls.op(i).subs(ex(repls), subs_options::subs_no_pattern), nummatches); + subsresult[j] = op(j) * power(itr->subs(ex(repls), subs_options::subs_no_pattern) / its->subs(ex(repls), subs_options::subs_no_pattern), nummatches); } subsed[j] = true; } @@ -613,9 +614,9 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con lst repls; for (size_t j=0; jnops(); j++) { - if (!subsed[j] && tryfactsubs(op(j), ls.op(i), nummatches, repls)) { + if (!subsed[j] && tryfactsubs(op(j), *its, nummatches, repls)) { subsed[j] = true; - subsresult[j] = op(j) * power(lr.op(i).subs(ex(repls), subs_options::subs_no_pattern) / ls.op(i).subs(ex(repls), subs_options::subs_no_pattern), nummatches); + subsresult[j] = op(j) * power(itr->subs(ex(repls), subs_options::subs_no_pattern) / its->subs(ex(repls), subs_options::subs_no_pattern), nummatches); } } } diff --git a/ginac/ncmul.cpp b/ginac/ncmul.cpp index b0fbb93d..110afa0d 100644 --- a/ginac/ncmul.cpp +++ b/ginac/ncmul.cpp @@ -269,7 +269,7 @@ void ncmul::append_factors(exvector & v, const ex & e) const if ((is_exactly_a(e)&&(e.return_type()!=return_types::commutative))|| (is_exactly_a(e))) { for (size_t i=0; iis_equal(e)) + return *its; // Otherwise create new symbol and add to list, taking care that the // replacement expression doesn't contain symbols from the sym_lst @@ -1712,9 +1713,9 @@ static ex replace_with_symbol(const ex &e, lst &sym_lst, lst &repl_lst) static ex replace_with_symbol(const ex &e, lst &repl_lst) { // Expression already in repl_lst? Then return the assigned symbol - for (size_t i=0; iop(1).is_equal(e)) + return it->op(0); // Otherwise create new symbol and add to list, taking care that the // replacement expression doesn't contain symbols from the sym_lst diff --git a/ginac/power.cpp b/ginac/power.cpp index 07a52ebe..2cae8a88 100644 --- a/ginac/power.cpp +++ b/ginac/power.cpp @@ -524,6 +524,7 @@ ex power::evalm(void) const return (new power(ebasis, eexponent))->setflag(status_flags::dynallocated); } +// from mul.cpp extern bool tryfactsubs(const ex &, const ex &, int &, lst &); ex power::subs(const lst & ls, const lst & lr, unsigned options) const @@ -538,11 +539,12 @@ ex power::subs(const lst & ls, const lst & lr, unsigned options) const if(!(options & subs_options::subs_algebraic)) return basic::subs(ls, lr, options); - for (size_t i=0; i::max(); lst repls; - if (tryfactsubs(*this, ls.op(i), nummatches, repls)) - return (ex_to((*this) * power(lr.op(i).subs(ex(repls), subs_options::subs_no_pattern) / ls.op(i).subs(ex(repls), subs_options::subs_no_pattern), nummatches))).basic::subs(ls, lr, options); + if (tryfactsubs(*this, *its, nummatches, repls)) + return (ex_to((*this) * power(itr->subs(ex(repls), subs_options::subs_no_pattern) / its->subs(ex(repls), subs_options::subs_no_pattern), nummatches))).basic::subs(ls, lr, options); } ex result=basic::subs(ls, lr, options); diff --git a/ginac/symbol.cpp b/ginac/symbol.cpp index 43e3684c..062a2890 100644 --- a/ginac/symbol.cpp +++ b/ginac/symbol.cpp @@ -110,9 +110,9 @@ ex symbol::unarchive(const archive_node &n, lst &sym_lst) ex s = (new symbol(n, sym_lst))->setflag(status_flags::dynallocated); // If symbol is in sym_lst, return the existing symbol - for (size_t i=0; i(sym_lst.op(i)) && (ex_to(sym_lst.op(i)).name == ex_to(s).name)) - return sym_lst.op(i); + for (lst::const_iterator it = sym_lst.begin(); it != sym_lst.end(); ++it) { + if (is_a(*it) && (ex_to(*it).name == ex_to(s).name)) + return *it; } // Otherwise add new symbol to list and return it diff --git a/ginac/symmetry.cpp b/ginac/symmetry.cpp index 2058a279..7934d673 100644 --- a/ginac/symmetry.cpp +++ b/ginac/symmetry.cpp @@ -426,20 +426,14 @@ ex symmetrize_cyclic(const ex & e, exvector::const_iterator first, exvector::con /** Symmetrize expression over a list of objects (symbols, indices). */ ex ex::symmetrize(const lst & l) const { - exvector v; - v.reserve(l.nops()); - for (size_t i=0; i