[GiNaC-list] Yet another version of the patch
Chris Dams
Chris.Dams at mi.infn.it
Wed Oct 5 18:51:16 CEST 2005
Hi again!
I found yet another error. This time because
combine_ex_with_coeff_to_pair should be used to insert something into an
expairseq. The latest version of my patch is attached.
Best wishes,
Chris
-------------- next part --------------
Index: ginac/expairseq.cpp
===================================================================
RCS file: /home/cvs/GiNaC/ginac/expairseq.cpp,v
retrieving revision 1.77
diff -r1.77 expairseq.cpp
36a37
> #include "indexed.h"
760,761c761,769
< construct_from_2_expairseq(ex_to<expairseq>(lh),
< ex_to<expairseq>(rh));
---
> if(is_a<mul>(lh))
> {
> ex newrh=rename_dummy_indices_uniquely(lh, rh);
> construct_from_2_expairseq(ex_to<expairseq>(lh),
> ex_to<expairseq>(newrh));
> }
> else
> construct_from_2_expairseq(ex_to<expairseq>(lh),
> ex_to<expairseq>(rh));
1010a1019
> exvector dummy_indices;
1014,1017c1023,1038
< const expairseq &subseqref = ex_to<expairseq>(*cit);
< combine_overall_coeff(subseqref.overall_coeff);
< epvector::const_iterator cit_s = subseqref.seq.begin();
< while (cit_s!=subseqref.seq.end()) {
---
> const expairseq *subseqref;
> ex newfactor;
> if(is_a<mul>(*cit))
> {
> exvector dummies_of_factor = get_all_dummy_indices(*cit);
> newfactor = rename_dummy_indices_uniquely(dummy_indices, dummies_of_factor, *cit);
> subseqref = &(ex_to<expairseq>(newfactor));
> exvector new_dummy_indices;
> set_union(dummy_indices.begin(), dummy_indices.end(), dummies_of_factor.begin(), dummies_of_factor.end(), std::back_insert_iterator<exvector>(new_dummy_indices), ex_is_less());
> dummy_indices.swap(new_dummy_indices);
> }
> else
> subseqref = &ex_to<expairseq>(*cit);
> combine_overall_coeff(subseqref->overall_coeff);
> epvector::const_iterator cit_s = subseqref->seq.begin();
> while (cit_s!=subseqref->seq.end()) {
1581a1603,1651
> class safe_inserter
> {
> public:
> safe_inserter(const ex&);
> std::auto_ptr<epvector> getseq(){return epv;}
> void insert_old_pair(const expair &p)
> {
> epv->push_back(p);
> }
> void insert_new_pair(const expair &p);
> private:
> std::auto_ptr<epvector> epv;
> bool dodummies;
> exvector dummy_indices;
> void update_dummy_indices(const exvector&);
> };
>
> safe_inserter::safe_inserter(const ex&e):epv(new epvector)
> {
> epv->reserve(e.nops());
> bool dodummies=is_a<mul>(e);
> if(dodummies)
> dummy_indices = get_all_dummy_indices(e);
> }
>
> void safe_inserter::update_dummy_indices(const exvector &v)
> {
> exvector new_dummy_indices;
> set_union(dummy_indices.begin(), dummy_indices.end(), v.begin(), v.end(),
> std::back_insert_iterator<exvector>(new_dummy_indices), ex_is_less());
> dummy_indices.swap(new_dummy_indices);
> }
>
> void safe_inserter::insert_new_pair(const expair &p)
> {
> if(!dodummies)
> {
> epv->push_back(p);
> return;
> }
> exvector dummies_of_factor = get_all_dummy_indices(p.rest);
> if(dummies_of_factor.size() == 0)
> { epv->push_back(p);
> return;
> }
> ex newfactor = rename_dummy_indices_uniquely(dummy_indices, dummies_of_factor, p.rest);
> update_dummy_indices(dummies_of_factor);
> epv -> push_back(expair(newfactor, p.coeff));
> }
1617,1618c1687
< std::auto_ptr<epvector> s(new epvector);
< s->reserve(seq.size());
---
> safe_inserter s(*this);
1621c1690,1691
< s->insert(s->begin(), seq.begin(), cit);
---
> for(epvector::const_iterator i=seq.begin(); i!=cit; ++i)
> s.insert_old_pair(*i);
1624c1694
< s->push_back(split_ex_to_pair(subsed_ex));
---
> s.insert_new_pair(split_ex_to_pair(subsed_ex));
1629c1699,1704
< s->push_back(split_ex_to_pair(recombine_pair_to_ex(*cit).subs(m, options)));
---
> ex orig_ex = recombine_pair_to_ex(*cit);
> ex subsed_ex = orig_ex.subs(m, options);
> if(are_ex_trivially_equal(orig_ex, subsed_ex))
> s.insert_old_pair(*cit);
> else
> s.insert_new_pair(split_ex_to_pair(subsed_ex));
1632c1707
< return s;
---
> return s.getseq();
1648,1649c1723
< std::auto_ptr<epvector> s(new epvector);
< s->reserve(seq.size());
---
> safe_inserter s(*this);
1652c1726,1727
< s->insert(s->begin(), seq.begin(), cit);
---
> for(epvector::const_iterator i=seq.begin(); i!=cit; ++i)
> s.insert_old_pair(*i);
1655c1730
< s->push_back(combine_ex_with_coeff_to_pair(subsed_ex, cit->coeff));
---
> s.insert_new_pair(combine_ex_with_coeff_to_pair(subsed_ex, cit->coeff));
1660,1661c1735,1740
< s->push_back(combine_ex_with_coeff_to_pair(cit->rest.subs(m, options),
< cit->coeff));
---
> ex orig_ex = cit->rest;
> ex subsed_ex = orig_ex.subs(m, options);
> if(are_ex_trivially_equal(orig_ex, subsed_ex))
> s.insert_old_pair(*cit);
> else
> s.insert_new_pair(combine_ex_with_coeff_to_pair(subsed_ex, cit->coeff));
1664c1743
< return s;
---
> return s.getseq();
Index: ginac/indexed.cpp
===================================================================
RCS file: /home/cvs/GiNaC/ginac/indexed.cpp,v
retrieving revision 1.96
diff -r1.96 indexed.cpp
534a535,543
> template<class T> size_t number_of_type(const exvector&v)
> {
> size_t number = 0;
> for(exvector::const_iterator i=v.begin(); i!=v.end(); ++i)
> if(is_exactly_a<T>(*i))
> ++number;
> return number;
> }
>
543c552
< static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, exvector & local_dummy_indices)
---
> template<class T> static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, exvector & local_dummy_indices)
545,546c554,555
< size_t global_size = global_dummy_indices.size(),
< local_size = local_dummy_indices.size();
---
> size_t global_size = number_of_type<T>(global_dummy_indices),
> local_size = number_of_type<T>(local_dummy_indices);
560c569
< if (find_if(global_dummy_indices.begin(), global_dummy_indices.end(), bind2nd(op0_is_equal(), *it)) == global_dummy_indices.end()) {
---
> if (is_exactly_a<T>(*it) && find_if(global_dummy_indices.begin(), global_dummy_indices.end(), bind2nd(idx_is_equal_ignore_dim(), *it)) == global_dummy_indices.end()) {
578,579c587,589
< for (size_t i=0; i<local_size; i++)
< local_syms.push_back(local_dummy_indices[i].op(0));
---
> for (size_t i=0; local_syms.size()!=local_size; i++)
> if(is_exactly_a<T>(local_dummy_indices[i]))
> local_syms.push_back(local_dummy_indices[i].op(0));
581,582c591,593
< for (size_t i=0; i<local_size; i++) // don't use more global symbols than necessary
< global_syms.push_back(global_dummy_indices[i].op(0));
---
> for (size_t i=0; global_syms.size()!=local_size; i++) // don't use more global symbols than necessary
> if(is_exactly_a<T>(global_dummy_indices[i]))
> global_syms.push_back(global_dummy_indices[i].op(0));
706a718,729
> template<class T> ex idx_symmetrization(const ex& r,const exvector& local_dummy_indices)
> { exvector dummy_syms;
> dummy_syms.reserve(r.nops());
> for (exvector::const_iterator it = local_dummy_indices.begin(); it != local_dummy_indices.end(); ++it)
> if(is_exactly_a<T>(*it))
> dummy_syms.push_back(it->op(0));
> if(dummy_syms.size() < 2)
> return r;
> ex q=symmetrize(r, dummy_syms);
> return q;
> }
>
867,875c890,903
< if (local_dummy_indices.size() >= 2) {
< exvector dummy_syms;
< dummy_syms.reserve(local_dummy_indices.size());
< for (exvector::const_iterator it = local_dummy_indices.begin(); it != local_dummy_indices.end(); ++it)
< dummy_syms.push_back(it->op(0));
< if (symmetrize(r, dummy_syms).is_zero()) {
< free_indices.clear();
< return _ex0;
< }
---
> ex q = idx_symmetrization<idx>(r, local_dummy_indices);
> if (q.is_zero()) {
> free_indices.clear();
> return _ex0;
> }
> q = idx_symmetrization<varidx>(q, local_dummy_indices);
> if (q.is_zero()) {
> free_indices.clear();
> return _ex0;
> }
> q = idx_symmetrization<spinidx>(q, local_dummy_indices);
> if (q.is_zero()) {
> free_indices.clear();
> return _ex0;
879c907,909
< r = rename_dummy_indices(r, dummy_indices, local_dummy_indices);
---
> r = rename_dummy_indices<idx>(r, dummy_indices, local_dummy_indices);
> r = rename_dummy_indices<varidx>(r, dummy_indices, local_dummy_indices);
> r = rename_dummy_indices<spinidx>(r, dummy_indices, local_dummy_indices);
945a976,986
> bool hasindex(const ex &x, const ex &sym)
> {
> if(is_a<idx>(x) && x.op(0)==sym)
> return true;
> else
> for(size_t i=0; i<x.nops(); ++i)
> if(hasindex(x.op(i), sym))
> return true;
> return false;
> }
>
974c1015,1018
< return rename_dummy_indices(e_expanded, dummy_indices, local_dummy_indices);
---
> e_expanded = rename_dummy_indices<idx>(e_expanded, dummy_indices, local_dummy_indices);
> e_expanded = rename_dummy_indices<varidx>(e_expanded, dummy_indices, local_dummy_indices);
> e_expanded = rename_dummy_indices<spinidx>(e_expanded, dummy_indices, local_dummy_indices);
> return e_expanded;
1018,1023d1061
< // Yes, construct vector of all dummy index symbols
< exvector dummy_syms;
< dummy_syms.reserve(dummy_indices.size());
< for (exvector::const_iterator it = dummy_indices.begin(); it != dummy_indices.end(); ++it)
< dummy_syms.push_back(it->op(0));
<
1029c1067,1074
< ex term_symm = symmetrize(term, dummy_syms);
---
> exvector dummy_indices_of_term;
> dummy_indices_of_term.reserve(dummy_indices.size());
> for(exvector::iterator i=dummy_indices.begin(); i!=dummy_indices.end(); ++i)
> if(hasindex(term,i->op(0)))
> dummy_indices_of_term.push_back(*i);
> ex term_symm = idx_symmetrization<idx>(term, dummy_indices_of_term);
> term_symm = idx_symmetrization<varidx>(term_symm, dummy_indices_of_term);
> term_symm = idx_symmetrization<spinidx>(term_symm, dummy_indices_of_term);
1321c1366
< ex rename_dummy_indices_uniquely(const ex & a, const ex & b)
---
> ex rename_dummy_indices_uniquely(const exvector & va, const exvector & vb, const ex & b)
1323c1368
< exvector va = get_all_dummy_indices(a), vb = get_all_dummy_indices(b), common_indices;
---
> exvector common_indices;
1333,1336c1378,1394
< if (is_a<varidx>(*ip)) {
< varidx mu((new symbol)->setflag(status_flags::dynallocated), ex_to<varidx>(*ip).get_dim(), ex_to<varidx>(*ip).is_covariant());
< old_indices.push_back(*ip);
< new_indices.push_back(mu);
---
> ex newsym=(new symbol)->setflag(status_flags::dynallocated);
> ex newidx;
> if(is_exactly_a<spinidx>(*ip))
> newidx = (new spinidx(newsym, ex_to<spinidx>(*ip).get_dim(),
> ex_to<spinidx>(*ip).is_covariant(),
> ex_to<spinidx>(*ip).is_dotted()))
> -> setflag(status_flags::dynallocated);
> else if (is_exactly_a<varidx>(*ip))
> newidx = (new varidx(newsym, ex_to<varidx>(*ip).get_dim(),
> ex_to<varidx>(*ip).is_covariant()))
> -> setflag(status_flags::dynallocated);
> else
> newidx = (new idx(newsym, ex_to<idx>(*ip).get_dim()))
> -> setflag(status_flags::dynallocated);
> old_indices.push_back(*ip);
> new_indices.push_back(newidx);
> if(is_a<varidx>(*ip)) {
1338,1341c1396
< new_indices.push_back(mu.toggle_variance());
< } else {
< old_indices.push_back(*ip);
< new_indices.push_back(idx((new symbol)->setflag(status_flags::dynallocated), ex_to<varidx>(*ip).get_dim()));
---
> new_indices.push_back(ex_to<varidx>(newidx).toggle_variance());
1346a1402,1408
> }
>
> ex rename_dummy_indices_uniquely(const ex & a, const ex & b)
> {
> exvector va = get_all_dummy_indices(a);
> exvector vb = get_all_dummy_indices(b);
> return rename_dummy_indices_uniquely(va, vb, b);
Index: ginac/indexed.h
===================================================================
RCS file: /home/cvs/GiNaC/ginac/indexed.h,v
retrieving revision 1.52
diff -r1.52 indexed.h
258a259,261
> /** Same as above, where va and vb contain the indices of a and b and are sorted */
> ex rename_dummy_indices_uniquely(const exvector & va, const exvector & vb, const ex & b);
>
More information about the GiNaC-list
mailing list