[GiNaC-devel] GiNaC 1.3.3 release

Chris Dams Chris.Dams at mi.infn.it
Mon Oct 24 17:57:57 CEST 2005


Dear developers,

On Mon, 24 Oct 2005, Richard B. Kreckel wrote:

> GiNaC 1.3.3 is out and available. The changes are:

Is the GiNaC development team making a habit out of ignoring patches? I
realize that because of adding a new overloaded version of
rename_dummy_indices_uniquely this probably can't go into 1.3.3 but the
separate symmetrization of different index types could have gone into it,
I think. Also, I would appreciate it if changes were added to the
1.4.0-branch. This spares me the trouble of managing a locally changed
version of GiNaC.

Here is the latest diff of my local copy of GiNaC with the CVS tree. New 
changes are:
(1) Made an automatic symplification that turns an indexed object without 
indices into its op(0). This is because it is sometimes useful to make an 
indexed object out of an exvector and an indexed object without indices is 
quite silly, really.

(2) Made the return_type and return_type_tinfo of an indexed object return
these values as they are of op(0). Reason: if the user wants to define his
own noncommutative functions/symbols and wants them to take indices, that
user probably wants these indexed objects to be noncummutative.

Best wishes,
Chris
-------------- next part --------------
Index: ginac/expairseq.cpp
===================================================================
RCS file: /home/cvs/GiNaC/ginac/expairseq.cpp,v
retrieving revision 1.77
diff -c -r1.77 expairseq.cpp
*** ginac/expairseq.cpp	1 May 2005 18:23:26 -0000	1.77
--- ginac/expairseq.cpp	24 Oct 2005 15:36:19 -0000
***************
*** 34,39 ****
--- 34,40 ----
  #include "archive.h"
  #include "operators.h"
  #include "utils.h"
+ #include "indexed.h"
  
  #if EXPAIRSEQ_USE_HASHTAB
  #include <cmath>
***************
*** 757,764 ****
  				construct_from_2_ex_via_exvector(lh,rh);
  			} else {
  #endif // EXPAIRSEQ_USE_HASHTAB
! 				construct_from_2_expairseq(ex_to<expairseq>(lh),
! 				                           ex_to<expairseq>(rh));
  #if EXPAIRSEQ_USE_HASHTAB
  			}
  #endif // EXPAIRSEQ_USE_HASHTAB
--- 758,772 ----
  				construct_from_2_ex_via_exvector(lh,rh);
  			} else {
  #endif // EXPAIRSEQ_USE_HASHTAB
! 				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));
  #if EXPAIRSEQ_USE_HASHTAB
  			}
  #endif // EXPAIRSEQ_USE_HASHTAB
***************
*** 1008,1020 ****
  	seq.reserve(v.size()+noperands-nexpairseqs);
  	
  	// copy elements and split off numerical part
  	cit = v.begin();
  	while (cit!=v.end()) {
  		if (ex_to<basic>(*cit).tinfo()==this->tinfo()) {
! 			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()) {
  				seq.push_back(*cit_s);
  				++cit_s;
  			}
--- 1016,1041 ----
  	seq.reserve(v.size()+noperands-nexpairseqs);
  	
  	// copy elements and split off numerical part
+ 	exvector dummy_indices;
  	cit = v.begin();
  	while (cit!=v.end()) {
  		if (ex_to<basic>(*cit).tinfo()==this->tinfo()) {
! 			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()) {
  				seq.push_back(*cit_s);
  				++cit_s;
  			}
***************
*** 1579,1584 ****
--- 1600,1654 ----
  	return std::auto_ptr<epvector>(0); // signalling nothing has changed
  }
  
+ 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());
+ 	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));
+ }
  
  /** Member-wise substitute in this sequence.
   *
***************
*** 1614,1635 ****
  			if (!are_ex_trivially_equal(orig_ex, subsed_ex)) {
  
  				// Something changed, copy seq, subs and return it
! 				std::auto_ptr<epvector> s(new epvector);
! 				s->reserve(seq.size());
  
  				// Copy parts of seq which are known not to have changed
! 				s->insert(s->begin(), seq.begin(), cit);
  
  				// Copy first changed element
! 				s->push_back(split_ex_to_pair(subsed_ex));
  				++cit;
  
  				// Copy rest
  				while (cit != last) {
! 					s->push_back(split_ex_to_pair(recombine_pair_to_ex(*cit).subs(m, options)));
  					++cit;
  				}
! 				return s;
  			}
  
  			++cit;
--- 1684,1710 ----
  			if (!are_ex_trivially_equal(orig_ex, subsed_ex)) {
  
  				// Something changed, copy seq, subs and return it
! 				safe_inserter s(*this);
  
  				// Copy parts of seq which are known not to have changed
! 				for(epvector::const_iterator i=seq.begin(); i!=cit; ++i)
! 					s.insert_old_pair(*i);
  
  				// Copy first changed element
! 				s.insert_new_pair(split_ex_to_pair(subsed_ex));
  				++cit;
  
  				// Copy rest
  				while (cit != last) {
! 					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));
  					++cit;
  				}
! 				return s.getseq();
  			}
  
  			++cit;
***************
*** 1645,1667 ****
  			if (!are_ex_trivially_equal(cit->rest, subsed_ex)) {
  			
  				// Something changed, copy seq, subs and return it
! 				std::auto_ptr<epvector> s(new epvector);
! 				s->reserve(seq.size());
  
  				// Copy parts of seq which are known not to have changed
! 				s->insert(s->begin(), seq.begin(), cit);
  			
  				// Copy first changed element
! 				s->push_back(combine_ex_with_coeff_to_pair(subsed_ex, cit->coeff));
  				++cit;
  
  				// Copy rest
  				while (cit != last) {
! 					s->push_back(combine_ex_with_coeff_to_pair(cit->rest.subs(m, options),
! 					                                           cit->coeff));
  					++cit;
  				}
! 				return s;
  			}
  
  			++cit;
--- 1720,1746 ----
  			if (!are_ex_trivially_equal(cit->rest, subsed_ex)) {
  			
  				// Something changed, copy seq, subs and return it
! 				safe_inserter s(*this);
  
  				// Copy parts of seq which are known not to have changed
! 				for(epvector::const_iterator i=seq.begin(); i!=cit; ++i)
! 					s.insert_old_pair(*i);
  			
  				// Copy first changed element
! 				s.insert_new_pair(combine_ex_with_coeff_to_pair(subsed_ex, cit->coeff));
  				++cit;
  
  				// Copy rest
  				while (cit != last) {
! 					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));
  					++cit;
  				}
! 				return s.getseq();
  			}
  
  			++cit;
Index: ginac/indexed.cpp
===================================================================
RCS file: /home/cvs/GiNaC/ginac/indexed.cpp,v
retrieving revision 1.96
diff -c -r1.96 indexed.cpp
*** ginac/indexed.cpp	19 May 2005 14:10:40 -0000	1.96
--- ginac/indexed.cpp	24 Oct 2005 15:36:19 -0000
***************
*** 297,302 ****
--- 297,305 ----
  		return f * thiscontainer(v);
  	}
  
+ 	if(seq.size()==1)
+ 		return base;
+ 
  	// Canonicalize indices according to the symmetry properties
  	if (seq.size() > 2) {
  		exvector v = seq;
***************
*** 532,537 ****
--- 535,549 ----
  	return f.get_free_indices();
  }
  
+ 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;
+ }
+ 
  /** Rename dummy indices in an expression.
   *
   *  @param e Expression to work on
***************
*** 540,549 ****
   *  @param global_dummy_indices The set of dummy indices that have appeared
   *    before and which we would like to use in "e", too. This gets updated
   *    by the function */
! static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, exvector & local_dummy_indices)
  {
! 	size_t global_size = global_dummy_indices.size(),
! 	       local_size = local_dummy_indices.size();
  
  	// Any local dummy indices at all?
  	if (local_size == 0)
--- 552,561 ----
   *  @param global_dummy_indices The set of dummy indices that have appeared
   *    before and which we would like to use in "e", too. This gets updated
   *    by the function */
! template<class T> static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, exvector & local_dummy_indices)
  {
! 	size_t global_size = number_of_type<T>(global_dummy_indices),
! 	       local_size = number_of_type<T>(local_dummy_indices);
  
  	// Any local dummy indices at all?
  	if (local_size == 0)
***************
*** 557,563 ****
  		int remaining = local_size - global_size;
  		exvector::const_iterator it = local_dummy_indices.begin(), itend = local_dummy_indices.end();
  		while (it != itend && remaining > 0) {
! 			if (find_if(global_dummy_indices.begin(), global_dummy_indices.end(), bind2nd(op0_is_equal(), *it)) == global_dummy_indices.end()) {
  				global_dummy_indices.push_back(*it);
  				global_size++;
  				remaining--;
--- 569,575 ----
  		int remaining = local_size - global_size;
  		exvector::const_iterator it = local_dummy_indices.begin(), itend = local_dummy_indices.end();
  		while (it != itend && remaining > 0) {
! 			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()) {
  				global_dummy_indices.push_back(*it);
  				global_size++;
  				remaining--;
***************
*** 575,585 ****
  	exvector local_syms, global_syms;
  	local_syms.reserve(local_size);
  	global_syms.reserve(local_size);
! 	for (size_t i=0; i<local_size; i++)
! 		local_syms.push_back(local_dummy_indices[i].op(0));
  	shaker_sort(local_syms.begin(), local_syms.end(), ex_is_less(), ex_swap());
! 	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));
  	shaker_sort(global_syms.begin(), global_syms.end(), ex_is_less(), ex_swap());
  
  	// Remove common indices
--- 587,599 ----
  	exvector local_syms, global_syms;
  	local_syms.reserve(local_size);
  	global_syms.reserve(local_size);
! 	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));
  	shaker_sort(local_syms.begin(), local_syms.end(), ex_is_less(), ex_swap());
! 	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));
  	shaker_sort(global_syms.begin(), global_syms.end(), ex_is_less(), ex_swap());
  
  	// Remove common indices
***************
*** 704,709 ****
--- 718,735 ----
  	}
  }
  
+ 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;
+ }
+ 
  /** Simplify product of indexed expressions (commutative, noncommutative and
   *  simple squares), return list of free indices. */
  ex simplify_indexed_product(const ex & e, exvector & free_indices, exvector & dummy_indices, const scalar_products & sp)
***************
*** 864,882 ****
  	// The result should be symmetric with respect to exchange of dummy
  	// indices, so if the symmetrization vanishes, the whole expression is
  	// zero. This detects things like eps.i.j.k * p.j * p.k = 0.
! 	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;
! 		}
  	}
  
  	// Dummy index renaming
! 	r = rename_dummy_indices(r, dummy_indices, local_dummy_indices);
  
  	// Product of indexed object with a scalar?
  	if (is_exactly_a<mul>(r) && r.nops() == 2
--- 890,915 ----
  	// The result should be symmetric with respect to exchange of dummy
  	// indices, so if the symmetrization vanishes, the whole expression is
  	// zero. This detects things like eps.i.j.k * p.j * p.k = 0.
! 	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;
  	}
  
  	// Dummy index renaming
! 	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);
  
  	// Product of indexed object with a scalar?
  	if (is_exactly_a<mul>(r) && r.nops() == 2
***************
*** 943,948 ****
--- 976,992 ----
  	}
  };
  
+ 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;
+ }
+ 
  /** Simplify indexed expression, return list of free indices. */
  ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indices, const scalar_products & sp)
  {
***************
*** 971,977 ****
  		}
  
  		// Rename the dummy indices
! 		return rename_dummy_indices(e_expanded, dummy_indices, local_dummy_indices);
  	}
  
  	// Simplification of sum = sum of simplifications, check consistency of
--- 1015,1024 ----
  		}
  
  		// Rename the 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;
  	}
  
  	// Simplification of sum = sum of simplifications, check consistency of
***************
*** 1015,1032 ****
  		if (num_terms_orig < 2 || dummy_indices.size() < 2)
  			return sum;
  
- 		// 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));
- 
  		// Chop the sum into terms and symmetrize each one over the dummy
  		// indices
  		std::vector<terminfo> terms;
  		for (size_t i=0; i<sum.nops(); i++) {
  			const ex & term = sum.op(i);
! 			ex term_symm = symmetrize(term, dummy_syms);
  			if (term_symm.is_zero())
  				continue;
  			terms.push_back(terminfo(term, term_symm));
--- 1062,1080 ----
  		if (num_terms_orig < 2 || dummy_indices.size() < 2)
  			return sum;
  
  		// Chop the sum into terms and symmetrize each one over the dummy
  		// indices
  		std::vector<terminfo> terms;
  		for (size_t i=0; i<sum.nops(); i++) {
  			const ex & term = sum.op(i);
! 			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);
  			if (term_symm.is_zero())
  				continue;
  			terms.push_back(terminfo(term, term_symm));
***************
*** 1318,1326 ****
  	return v;
  }
  
! ex rename_dummy_indices_uniquely(const ex & a, const ex & b)
  {
! 	exvector va = get_all_dummy_indices(a), vb = get_all_dummy_indices(b), common_indices;
  	set_intersection(va.begin(), va.end(), vb.begin(), vb.end(), std::back_insert_iterator<exvector>(common_indices), ex_is_less());
  	if (common_indices.empty()) {
  		return b;
--- 1366,1374 ----
  	return v;
  }
  
! ex rename_dummy_indices_uniquely(const exvector & va, const exvector & vb, const ex & b)
  {
! 	exvector common_indices;
  	set_intersection(va.begin(), va.end(), vb.begin(), vb.end(), std::back_insert_iterator<exvector>(common_indices), ex_is_less());
  	if (common_indices.empty()) {
  		return b;
***************
*** 1330,1349 ****
  		new_indices.reserve(2*common_indices.size());
  		exvector::const_iterator ip = common_indices.begin(), ipend = common_indices.end();
  		while (ip != ipend) {
! 			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);
  				old_indices.push_back(ex_to<varidx>(*ip).toggle_variance());
! 				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()));
  			}
  			++ip;
  		}
  		return b.subs(lst(old_indices.begin(), old_indices.end()), lst(new_indices.begin(), new_indices.end()), subs_options::no_pattern);
  	}
  }
  
  ex expand_dummy_sum(const ex & e, bool subs_idx)
--- 1378,1414 ----
  		new_indices.reserve(2*common_indices.size());
  		exvector::const_iterator ip = common_indices.begin(), ipend = common_indices.end();
  		while (ip != ipend) {
! 			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)) {
  				old_indices.push_back(ex_to<varidx>(*ip).toggle_variance());
! 				new_indices.push_back(ex_to<varidx>(newidx).toggle_variance());
  			}
  			++ip;
  		}
  		return b.subs(lst(old_indices.begin(), old_indices.end()), lst(new_indices.begin(), new_indices.end()), subs_options::no_pattern);
  	}
+ }
+ 
+ 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);
  }
  
  ex expand_dummy_sum(const ex & e, bool subs_idx)
Index: ginac/indexed.h
===================================================================
RCS file: /home/cvs/GiNaC/ginac/indexed.h,v
retrieving revision 1.52
diff -c -r1.52 indexed.h
*** ginac/indexed.h	19 May 2005 14:10:40 -0000	1.52
--- ginac/indexed.h	24 Oct 2005 15:36:19 -0000
***************
*** 153,159 ****
  	ex derivative(const symbol & s) const;
  	ex thiscontainer(const exvector & v) const;
  	ex thiscontainer(std::auto_ptr<exvector> vp) const;
! 	unsigned return_type() const { return return_types::commutative; }
  	ex expand(unsigned options = 0) const;
  
  	// new virtual functions which can be overridden by derived classes
--- 153,160 ----
  	ex derivative(const symbol & s) const;
  	ex thiscontainer(const exvector & v) const;
  	ex thiscontainer(std::auto_ptr<exvector> vp) const;
! 	unsigned return_type() const { return op(0).return_type(); }
! 	unsigned return_type_tinfo() const { return op(0).return_type_tinfo(); }
  	ex expand(unsigned options = 0) const;
  
  	// new virtual functions which can be overridden by derived classes
***************
*** 255,260 ****
--- 256,264 ----
  
  /** Returns b with all dummy indices, which are common with a, renamed */
  ex rename_dummy_indices_uniquely(const ex & a, const ex & b);
+ 
+ /** 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);
  
  /** This function returns the given expression with expanded sums
   *  for all dummy index summations, where the dimensionality of 
Index: ginac/operators.cpp
===================================================================
RCS file: /home/cvs/GiNaC/ginac/operators.cpp,v
retrieving revision 1.33
diff -c -r1.33 operators.cpp
*** ginac/operators.cpp	29 Jun 2005 14:09:24 -0000	1.33
--- ginac/operators.cpp	24 Oct 2005 15:36:19 -0000
***************
*** 47,52 ****
--- 47,55 ----
  	// Check if we are constructing a mul object or a ncmul object.  Due to
  	// ncmul::eval()'s rule to pull out commutative elements we need to check
  	// only one of the elements.
+ 	std::cout << "exmul: multiplying " << lh << " and " << rh << std::endl;
+ 	std::cout << "lhs is commutative " << int(lh.return_type()==return_types::commutative) << std::endl;
+ 	std::cout << "rhs is commutative " << int(rh.return_type()==return_types::commutative) << std::endl;
  	if (rh.return_type()==return_types::commutative ||
  	    lh.return_type()==return_types::commutative)
  		return (new mul(lh,rh))->setflag(status_flags::dynallocated);
Index: ginac/tensor.h
===================================================================
RCS file: /home/cvs/GiNaC/ginac/tensor.h,v
retrieving revision 1.23
diff -c -r1.23 tensor.h
*** ginac/tensor.h	1 May 2005 18:23:27 -0000	1.23
--- ginac/tensor.h	24 Oct 2005 15:36:19 -0000
***************
*** 65,70 ****
--- 65,71 ----
  
  	// non-virtual functions in this class
  protected:
+ 	unsigned return_type() const { return return_types::commutative; }
  	void do_print(const print_context & c, unsigned level) const;
  	void do_print_latex(const print_latex & c, unsigned level) const;
  };
***************
*** 84,89 ****
--- 85,91 ----
  
  	// non-virtual functions in this class
  protected:
+ 	unsigned return_type() const { return return_types::commutative; }
  	void do_print(const print_context & c, unsigned level) const;
  };
  
***************
*** 106,111 ****
--- 108,114 ----
  
  	// non-virtual functions in this class
  protected:
+ 	unsigned return_type() const { return return_types::commutative; }
  	void do_print(const print_context & c, unsigned level) const;
  	void do_print_latex(const print_latex & c, unsigned level) const;
  
***************
*** 153,158 ****
--- 156,162 ----
  
  	// non-virtual functions in this class
  protected:
+ 	unsigned return_type() const { return return_types::commutative; }
  	void do_print(const print_context & c, unsigned level) const;
  	void do_print_latex(const print_latex & c, unsigned level) const;
  


More information about the GiNaC-devel mailing list