Commits

Volker Braun committed dd427ea

finished the infinities rewrite

  • Participants
  • Parent commits 2e11565

Comments (0)

Files changed (8)

File cleanup.patch

-# HG changeset patch
-# Parent d388e126049584546003d4835dfe636730a3b619
-
-diff --git a/ginac/expairseq.h b/ginac/expairseq.h
---- a/ginac/expairseq.h
-+++ b/ginac/expairseq.h
-@@ -75,6 +75,7 @@
- 	expairseq(const exvector & v);
- 	expairseq(const epvector & v, const ex & oc, bool do_index_renaming = false);
- 	expairseq(std::auto_ptr<epvector>, const ex & oc, bool do_index_renaming = false);
-+	~expairseq() { if (!seq_sorted) delete seq_sorted; }
- 	
- 	// functions overriding virtual functions from base classes
- public:
-diff --git a/ginac/order.cpp b/ginac/order.cpp
---- a/ginac/order.cpp
-+++ b/ginac/order.cpp
-@@ -30,24 +30,18 @@
- 
- namespace GiNaC {
- 
--/** What SAGE does for printing:
-- To print multivariate polynomials, SAGE uses "reversed" (bigger terms first) 'degrevlex' order
-- i.e. "reversed" graded reversed lexicographic order, the lexicographic order is defined below.
-- The variables are ordered according to their "creation" order
-- i.e PR.<x,a,t> = PolynomialRing(QQ) gives x > a > t
-- and x+t+a is printed x + a + t
-+/** What Sage does for printing:
-+ To print multivariate polynomials, SAGE uses "reversed" (bigger terms first)
-+ 'degrevlex' order i.e. "reversed" graded reversed lexicographic order, the
-+ lexicographic order is defined below.  The variables are ordered according to
-+ their "creation" order i.e PR.<x,a,t> = PolynomialRing(QQ) gives x > a > t and
-+ x+t+a is printed x + a + t
- */
- 
--/** Returns true if lhex > rhex for 'degrevlex' of SAGE
-+/** Returns true if lhex > rhex for 'degrevlex' of Sage
-  i.e. graded reversed lexicographic order
- */
- bool ex_is_greater_degrevlex::operator() (const ex &lhex, const ex &rhex) const {
--        //std::cout<<"in ex_is_greater_degrevlex::operator()"<<std::endl;
--	//std::cout<<"lh: ";
--	//lhex.dbgprint();
--	//std::cout<<std::endl<<"rh: ";
--	//rhex.dbgprint();
--	//std::cout<<std::endl;
- 	const basic *lh = get_pointer(lhex.bp);
- 	const basic *rh = get_pointer(rhex.bp);
- 	return compare(lh, rh) == 1;
-@@ -71,14 +65,8 @@
- 	}
- 
- 	// compare coeffs which are numerics
--	//std::cout << "comparing coeffs" << std::endl;
--	//lhex.coeff.bp->dbgprinttree();
--	//rhex.coeff.bp->dbgprinttree();
- 	numeric lh_coeff = ex_to<numeric>(lhex.coeff);
- 	numeric rh_coeff = ex_to<numeric>(rhex.coeff);
--	// strange...
--	//std::cout<<lh->compare_same_type(rh)<<std::endl;
--	//return lh->compare_same_type(rh) == 1;
- 	double lh_deg = 1;
- 	double rh_deg = 1;
- 		if (lh_coeff.is_real()) {
-@@ -98,10 +86,10 @@
- }
- 
- /** Comparison functions:
-- They should implement 'degrevlex' of SAGE
-+ They should implement 'degrevlex' of Sage
-  i.e. graded reversed lexicographic order
-- The lexicographic order should depend on the "creation" order of variables ?
-- Or should it be a "natural" one on strings : a > b > ... > x > y > ... ?
-+ The lexicographic order is the "natural" one on strings:
-+ a > b > ... > x > y > ... ?
- */
- 
- /** Return values for comparison functions :

File expair_is_greater.patch

-# HG changeset patch
-# Parent 9fe973db3ada1cab96f0233df1160819744322aa
-
-diff --git a/ginac/add.cpp b/ginac/add.cpp
---- a/ginac/add.cpp
-+++ b/ginac/add.cpp
-@@ -30,6 +30,7 @@
- #include "ncmul.h"
- #include "constant.h"
- #include "compiler.h"
-+#include "order.h"
- 
- #include <sstream>
- #include <iostream>
-@@ -707,4 +708,15 @@
- 	return (new add(vp, overall_coeff))->setflag(status_flags::dynallocated | (options == 0 ? status_flags::expanded : 0));
- }
- 
-+epvector* add::get_sorted_seq() const
-+{
-+	if (!this->seq_sorted) {
-+		seq_sorted = new epvector(seq.size());
-+		partial_sort_copy(seq.begin(), seq.end(),
-+				seq_sorted->begin(), seq_sorted->end(),
-+				expair_is_greater_degrevlex_mul());
-+	}
-+	return seq_sorted;
-+}
-+
- } // namespace GiNaC
-diff --git a/ginac/add.h b/ginac/add.h
---- a/ginac/add.h
-+++ b/ginac/add.h
-@@ -63,6 +63,7 @@
- 	ex imag_part() const;
- 	exvector get_free_indices() const;
- 	ex eval_ncmul(const exvector & v) const;
-+	epvector* get_sorted_seq() const;
- protected:
- 	ex derivative(const symbol & s) const;
- 	unsigned return_type() const;
-diff --git a/ginac/ex.h b/ginac/ex.h
---- a/ginac/ex.h
-+++ b/ginac/ex.h
-@@ -84,6 +84,7 @@
- 	template<class T> friend inline bool is_exactly_a(const ex &);
- 	
- 	friend struct ex_is_greater_degrevlex;
-+	friend struct ex_is_greater_degrevlex_mul;
- 	friend struct expair_is_greater_degrevlex;
- 	// default constructor, copy constructor and assignment operator
- public:
-diff --git a/ginac/expairseq.cpp b/ginac/expairseq.cpp
---- a/ginac/expairseq.cpp
-+++ b/ginac/expairseq.cpp
-@@ -32,7 +32,6 @@
- #include "utils.h"
- #include "indexed.h"
- #include "constant.h"
--#include "order.h"
- 
- #include <iostream>
- #include <algorithm>
-@@ -1171,17 +1170,6 @@
- 	std::sort(seq.begin(), seq.end(), expair_rest_is_less());
- }
- 
--epvector* expairseq::get_sorted_seq() const
--{
--	if (!this->seq_sorted) {
--		seq_sorted = new epvector(seq.size());
--		partial_sort_copy(seq.begin(), seq.end(),
--				seq_sorted->begin(), seq_sorted->end(),
--				  expair_is_greater_degrevlex::in_type(tinfo()));
--	}
--	return seq_sorted;
--}
--
- 
- /** Compact a presorted expairseq by combining all matching expairs to one
-  *  each.  On an add object, this is responsible for 2*x+3*x+y -> 5*x+y, for
-diff --git a/ginac/expairseq.h b/ginac/expairseq.h
---- a/ginac/expairseq.h
-+++ b/ginac/expairseq.h
-@@ -134,7 +134,6 @@
- 	void make_flat(const epvector & v, bool do_index_renaming = false);
- 	void canonicalize();
- 	void combine_same_terms_sorted_seq();
--	epvector* get_sorted_seq() const;
- #if EXPAIRSEQ_USE_HASHTAB
- 	void combine_same_terms();
- 	unsigned calc_hashtabsize(unsigned sz) const;
-diff --git a/ginac/mul.cpp b/ginac/mul.cpp
---- a/ginac/mul.cpp
-+++ b/ginac/mul.cpp
-@@ -1620,6 +1620,17 @@
- 	}
- }
- 
-+epvector* mul::get_sorted_seq() const
-+{
-+	if (!this->seq_sorted) {
-+		seq_sorted = new epvector(seq.size());
-+		partial_sort_copy(seq.begin(), seq.end(),
-+				seq_sorted->begin(), seq_sorted->end(),
-+				expair_is_greater_degrevlex_mul());
-+	}
-+	return seq_sorted;
-+}
-+
-   
- //////////
- // new virtual functions which can be overridden by derived classes
-diff --git a/ginac/mul.h b/ginac/mul.h
---- a/ginac/mul.h
-+++ b/ginac/mul.h
-@@ -95,6 +95,7 @@
- public:
- 	ex algebraic_subs_mul(const exmap & m, unsigned options) const;
- 	double total_degree() const;
-+	epvector* get_sorted_seq() const;
- 	//int compare_symbol(const symbol &other) const;
- 	//int compare_pow(const power &other) const;
- protected:
-diff --git a/ginac/order.cpp b/ginac/order.cpp
---- a/ginac/order.cpp
-+++ b/ginac/order.cpp
-@@ -30,6 +30,15 @@
- 
- namespace GiNaC {
- 
-+double numeric_to_double(const numeric& exp)
-+{
-+	if (exp.is_real())
-+		return exp.to_double();
-+	else 
-+		return std::sqrt(std::pow(exp.real().to_double(), 2) + 
-+				std::pow(exp.imag().to_double(), 2));
-+}
-+
- /** What Sage does for printing:
-  To print multivariate polynomials, SAGE uses "reversed" (bigger terms first)
-  'degrevlex' order i.e. "reversed" graded reversed lexicographic order, the
-@@ -53,34 +62,34 @@
- 	return compare(lh, rh);
- }
- 
-+bool expair_is_greater_degrevlex_mul::operator() (const expair &lhex, const expair &rhex) const
-+{
-+	int cmpval = ex_is_greater_degrevlex_mul().compare(lhex.rest, rhex.rest);
-+	if (cmpval != 0) {
-+		return cmpval == 1;
-+	}
-+	return compare_degrees(lhex, rhex);
-+}
-+
- bool expair_is_greater_degrevlex::operator() (const expair &lhex, const expair &rhex) const
- {
- 	const basic *lh = get_pointer(lhex.rest.bp);
- 	const basic *rh = get_pointer(rhex.rest.bp);
- 
- 	// compare rests
--        int cmpval = ex_is_greater_degrevlex::in_type(in_type_id).compare(lh, rh);
-+	int cmpval = ex_is_greater_degrevlex().compare(lh, rh);
- 	if (cmpval != 0) {
--	        return cmpval == 1;
-+		return cmpval == 1;
- 	}
-+	return compare_degrees(lhex, rhex);
-+}
- 
-+bool expair_is_greater_degrevlex::compare_degrees(const expair &lhex,
-+		const expair &rhex) const
-+{
- 	// compare coeffs which are numerics
--	numeric lh_coeff = ex_to<numeric>(lhex.coeff);
--	numeric rh_coeff = ex_to<numeric>(rhex.coeff);
--	double lh_deg = 1;
--	double rh_deg = 1;
--		if (lh_coeff.is_real()) {
--			lh_deg = lh_coeff.to_double();
--		} else {
--			lh_deg = std::sqrt(std::pow(lh_coeff.real().to_double(), 2) + 
--					std::pow(lh_coeff.imag().to_double(), 2));
--		}
--		if (rh_coeff.is_real()) {
--			rh_deg = rh_coeff.to_double();
--		} else {
--			rh_deg = std::sqrt(std::pow(rh_coeff.real().to_double(), 2) + 
--					std::pow(rh_coeff.imag().to_double(), 2));
--		}
-+	double lh_deg = numeric_to_double(ex_to<numeric>(lhex.coeff));
-+	double rh_deg = numeric_to_double(ex_to<numeric>(rhex.coeff));
- 
- 	return lh_deg > rh_deg;	
- }
-@@ -130,23 +139,23 @@
- 					static_cast<const fderivative*>(lh),
- 					static_cast<const fderivative*>(rh));
- 		} else {
--		        // using GiNaC functions by default
-+			// using GiNaC functions by default
- 			return lh->compare_same_type(*rh);
- 		}
- 	// at present numerics are combined into overall_coefficient
- 	// even when hold parameter is used
- 	} else if (typeid_lh == numeric_id) {
- 	 	//print numerics after anything else
--	        return -1;
-+		return -1;
- 	} else if (typeid_rh == numeric_id) {
- 	 	//print numerics after anything else
--	        return 1;
-+		return 1;
- 	} else if (typeid_lh == constant_id) {
- 	 	//print constants after anything else
--	        return -1;
-+		return -1;
- 	} else if (typeid_rh == constant_id) {
- 	 	//print constants after anything else
--	        return 1;
-+		return 1;
- 	} else if (typeid_lh == fderivative_id) {
- 		//print fderivatives after everything else
- 		return -1;
-@@ -155,10 +164,10 @@
- 		return 1;
- 	} else if (typeid_lh == function_id) {
- 		//print functions before fderivatives, after anything else
--	        return -1;
-+		return -1;
- 	} else if (typeid_rh == function_id) {
- 		//print functions before fderivatives, after anything else
--	        return 1;
-+		return 1;
- 	} else if (typeid_lh == mul_id) {
- 		if (typeid_rh == power_id) {
- 			return compare_mul_power(
-@@ -215,9 +224,8 @@
- 					static_cast<const power*>(rh),
- 					static_cast<const symbol*>(lh));
- 		}
--        }
--	std::cout<<"comparing typeid's"<<std::endl;
--	return (typeid_lh<typeid_rh ? -1 : 1);
-+	}
-+	throw (std::runtime_error("comparing typeid's"));
- }
- 
- // compare a mul and a symbol objects
-@@ -227,16 +235,13 @@
- int ex_is_greater_degrevlex::compare_mul_symbol(const mul *lh,
- 		const symbol *rh) const
- {
--  //std::cout<<"comparing mul and symbol"<<std::endl;
--  //lh->dbgprint();
--  //rh->dbgprint();
- 	int cmpval;
- 
- 	double tdeg;
- 	tdeg = lh->total_degree();
- 
- 	if (tdeg != 1)
--	        return tdeg > 1 ? 1 : -1;
-+		return tdeg > 1 ? 1 : -1;
- 	
- 	const expair smallest_item = lh->get_sorted_seq()->back();
- 
-@@ -259,9 +264,8 @@
- 	return 1;
- }
- 
--
- // compare a mul and a pow objects
--// same behavior wihtin mul and add objects:
-+// same behavior within mul and add objects:
- // first we compare total degrees
- // if equal we compare the smallest basis in the sequence to the basis in other
- // then their exponents
-@@ -274,13 +278,7 @@
- 	double rh_deg = 1;
- 	numeric rh_exp;
- 	if (is_a<numeric>(rh->exponent)) {
--		rh_exp = ex_to<numeric>(rh->exponent);
--		if (rh_exp.is_real()) {
--			rh_deg = rh_exp.to_double();
--		} else {
--			rh_deg = std::sqrt(std::pow(rh_exp.real().to_double(), 2) + 
--					std::pow(rh_exp.imag().to_double(), 2));
--		}
-+		rh_deg = numeric_to_double(ex_to<numeric>(rh->exponent));
- 	}
- 	if (rh_deg != lh_deg)
- 		return lh_deg < rh_deg ? -1 : 1;
-@@ -298,7 +296,7 @@
- 
- 	// compare exponents
- 	cmpval = -compare(get_pointer(smallest_item.coeff.bp),
--		        get_pointer(rh->exponent.bp));
-+			get_pointer(rh->exponent.bp));
- 	if (cmpval != 0) {
- 		return cmpval;
- 	}
-@@ -319,15 +317,11 @@
- int ex_is_greater_degrevlex::compare_same_type_mul(const mul *lh,
- 		const mul *rh) const
- {
--  //std::cout<<"comparing mul and mul"<<std::endl;
--  //lh->dbgprint();
--  //rh->dbgprint();
- 	int cmpval;
- 	
- 	// compare total degrees
- 	double lh_deg = lh->total_degree();
- 	double rh_deg = rh->total_degree();
--        //std::cout<<"degree "<<lh_deg<<" and "<<rh_deg<<std::endl;
- 	if (lh_deg != rh_deg)
- 		return lh_deg < rh_deg ? -1 : 1;
- 	
-@@ -341,13 +335,15 @@
- 
- 	for (; (cit1!=last1)&&(cit2!=last2); ++cit1, ++cit2) {
- 		// compare bases
--		cmpval = compare(get_pointer(cit1->rest.bp),get_pointer(cit2->rest.bp));
-+		cmpval = compare(get_pointer(cit1->rest.bp),
-+				get_pointer(cit2->rest.bp));
- 		if (cmpval != 0) {
- 			return cmpval;
- 		}
- 
- 		// compare exponents
--	        cmpval = -compare(get_pointer(cit1->coeff.bp),get_pointer(cit2->coeff.bp));
-+		cmpval = -compare(get_pointer(cit1->coeff.bp),
-+				get_pointer(cit2->coeff.bp));
- 		if (cmpval != 0) {
- 			return cmpval;
- 		}
-@@ -372,15 +368,12 @@
- int ex_is_greater_degrevlex::compare_add_symbol(const add *lh,
- 		const symbol *rh) const
- {
--  //std::cout<<"comparing add and symbol"<<std::endl;
--  //lh->dbgprint();
--  //rh->dbgprint();
- 	int cmpval;
- 
- 	const expair biggest_item = lh->get_sorted_seq()->front();
- 
- 	// compare bases
--	cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(biggest_item.rest.bp), rh);
-+	cmpval = ex_is_greater_degrevlex().compare(get_pointer(biggest_item.rest.bp), rh);
- 	if (cmpval != 0) {
- 		return cmpval;
- 	}
-@@ -404,11 +397,11 @@
- int ex_is_greater_degrevlex::compare_add_mul(const add *lh,
- 		const mul *rh) const
- {
--        int cmpval;
-+	int cmpval;
- 	const expair biggest_item = lh->get_sorted_seq()->front();
- 
- 	// compare bases
--	cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(biggest_item.rest.bp), rh);
-+	cmpval = ex_is_greater_degrevlex().compare(get_pointer(biggest_item.rest.bp), rh);
- 	if (cmpval != 0) {
- 		return cmpval;
- 	}
-@@ -432,11 +425,11 @@
- int ex_is_greater_degrevlex::compare_add_power(const add *lh,
- 		const power *rh) const
- {
--        int cmpval;
-+	int cmpval;
- 	const expair biggest_item = lh->get_sorted_seq()->front();
- 
- 	// compare bases
--	cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(biggest_item.rest.bp), rh);
-+	cmpval = ex_is_greater_degrevlex().compare(get_pointer(biggest_item.rest.bp), rh);
- 	if (cmpval != 0) {
- 		return cmpval;
- 	}
-@@ -473,13 +466,13 @@
- 
- 	for (; (cit1!=last1)&&(cit2!=last2); ++cit1, ++cit2) {
- 		// compare bases
--		cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(cit1->rest.bp),get_pointer(cit2->rest.bp));
-+		cmpval = ex_is_greater_degrevlex().compare(get_pointer(cit1->rest.bp),get_pointer(cit2->rest.bp));
- 		if (cmpval != 0) {
- 			return cmpval;
- 		}
- 
- 		// compare coefficients
--	        cmpval = compare(get_pointer(cit1->coeff.bp),get_pointer(cit2->coeff.bp));
-+		cmpval = compare(get_pointer(cit1->coeff.bp),get_pointer(cit2->coeff.bp));
- 		if (cmpval != 0) {
- 			return cmpval;
- 		}
-@@ -506,32 +499,27 @@
- // in add object:
- // first exponents
- // then bases
-+int ex_is_greater_degrevlex_mul::compare_power_symbol(const power *lh,
-+		const symbol *rh) const
-+{
-+	int cmpval;
-+	cmpval = compare(get_pointer(lh->basis.bp), rh);
-+	if (cmpval != 0)
-+		  return cmpval;
-+	
-+	cmpval = compare(get_pointer(lh->exponent.bp), _num1_p);
-+	
-+	return cmpval;
-+}
-+
- int ex_is_greater_degrevlex::compare_power_symbol(const power *lh,
- 		const symbol *rh) const
- {
--        int cmpval;
--
--	//std::cout<<"in compare_power_symbol"<<std::endl;
--        if (in_type_id == &mul::tinfo_static) {
--	        cmpval = compare(get_pointer(lh->basis.bp), rh);
--		if (cmpval != 0)
--        		  return cmpval;
--		
--		cmpval = compare(get_pointer(lh->exponent.bp), _num1_p);
--		
--		return cmpval;
--	}
-+	int cmpval;
- 
- 	// We are in an add object
- 	if (is_a<numeric>(lh->exponent)) {
--		numeric lh_exp = ex_to<numeric>(lh->exponent);
--		double lh_deg;
--		if (lh_exp.is_real()) {
--			lh_deg = lh_exp.to_double();
--		} else {
--			lh_deg = std::sqrt(std::pow(lh_exp.real().to_double(), 2) + 
--					std::pow(lh_exp.imag().to_double(), 2));
--		}
-+		double lh_deg = numeric_to_double(ex_to<numeric>(lh->exponent));
- 		if (lh_deg != 1)
- 			return lh_deg < 1 ? -1 : 1;
- 	}
-@@ -549,50 +537,41 @@
- // in add object:
- // first exponents
- // then bases
-+int ex_is_greater_degrevlex_mul::compare_same_type_power(const power *lh,
-+		const power *rh) const
-+{
-+	int cmpval;
-+	cmpval = compare(get_pointer(lh->basis.bp), get_pointer(rh->basis.bp));
-+	if (cmpval != 0)
-+		  return cmpval;
-+	cmpval = compare(get_pointer(lh->exponent.bp), get_pointer(rh->exponent.bp));
-+	
-+	return cmpval;
-+}
- int ex_is_greater_degrevlex::compare_same_type_power(const power *lh,
--	        const power *rh) const
-+		const power *rh) const
- {
--        int cmpval;
--        if (in_type_id == &mul::tinfo_static) {
--	        cmpval = compare(get_pointer(lh->basis.bp), get_pointer(rh->basis.bp));
--		if (cmpval != 0)
--        		  return cmpval;
--		cmpval = compare(get_pointer(lh->exponent.bp), get_pointer(rh->exponent.bp));
--		
--		return cmpval;
--	}
-+	int cmpval;
- 
- 	double lh_deg = 1;
- 	double rh_deg = 1;
- 	if (is_a<numeric>(lh->exponent)) {
--		numeric lh_exp = ex_to<numeric>(lh->exponent);
--		if (lh_exp.is_real()) {
--			lh_deg = lh_exp.to_double();
--		} else {
--			lh_deg = std::sqrt(std::pow(lh_exp.real().to_double(), 2) + 
--					std::pow(lh_exp.imag().to_double(), 2));
--		}
-+		lh_deg = numeric_to_double(ex_to<numeric>(lh->exponent));
- 	}
- 	if (is_a<numeric>(rh->exponent)) {
--		numeric rh_exp = ex_to<numeric>(rh->exponent);
--		if (rh_exp.is_real()) {
--			rh_deg = rh_exp.to_double();
--		} else {
--			rh_deg = std::sqrt(std::pow(rh_exp.real().to_double(), 2) + 
--					std::pow(rh_exp.imag().to_double(), 2));
--		}
-+		rh_deg = numeric_to_double(ex_to<numeric>(rh->exponent));
- 	}
- 	if (rh_deg != lh_deg)
- 		return lh_deg < rh_deg ? -1 : 1;
- 
--	cmpval = compare(get_pointer(lh->basis.bp),
--			 get_pointer(rh->basis.bp));
-+	cmpval = compare(get_pointer(lh->basis.bp), get_pointer(rh->basis.bp));
- 	if (cmpval != 0)
--        	return cmpval;
-+		return cmpval;
- 
- 	if (is_a<numeric>(lh->exponent) && is_a<numeric>(rh->exponent))
--	        return 0;
--	return compare(get_pointer(lh->exponent.bp), get_pointer(rh->exponent.bp));
-+		return 0;
-+	return compare(get_pointer(lh->exponent.bp),
-+			get_pointer(rh->exponent.bp));
- }
- 
- // compare two symbol objects
-@@ -601,26 +580,11 @@
- int ex_is_greater_degrevlex::compare_same_type_symbol(const symbol *lh,
- 		const symbol *rh) const
- {
--        //std::cout<<"in compare_same symbol"<<std::endl;
--	//std::cout<<"lh: ";
--	//lh->dbgprint();
--	//std::cout<<"rh: ";
--	//rh->dbgprint();
--	//std::cout<<std::endl;
--
--        // Weird because Sage sorts based on creation order.
--        // SAGE/Pynac: Sorting based on creation order doesn't work for Sage.
--        // instead we sort on variable name. -- William Stein
--
--	//std::cout<<"comparing names: "<<(lh->name < rh->name)<<std::endl;
--	/* Reversed ordering on char encoding (i.e. "a" < "b") then length (i.e. "abc" < "abcd")
-+	/* Reversed ordering on char encoding (i.e. "a" < "b") then length
-+	 * (i.e. "abc" < "abcd")
- 	   i.e. "x" > "y" and "xyz" > "xyzt" */
- 	if (lh->serial==rh->serial) return 0;
--	//std::cout<<"after if"<<std::endl;
- 	return lh->name < rh->name ? 1 : -1;
--
--	//return lh->serial < rh->serial ? 1 : -1;
--	//return -lh->compare_same_type(*rh);
- }
- 
- // compare two containers of the same type
-@@ -628,11 +592,11 @@
- int ex_is_greater_degrevlex::compare_same_type_container(const container<C> *lh,
- 							 const container<C> *rh) const
- {
--        typename C<ex>::const_iterator it1 = lh->seq.begin(), it1end = lh->seq.end(),
--	                      it2 = rh->seq.begin(), it2end = rh->seq.end();
-+	typename C<ex>::const_iterator it1 = lh->seq.begin(), it1end = lh->seq.end(),
-+			      it2 = rh->seq.begin(), it2end = rh->seq.end();
- 
- 	while (it1 != it1end && it2 != it2end) {
--	        int cmpval = compare(get_pointer(it1->bp), get_pointer(it2->bp));
-+		int cmpval = compare(get_pointer(it1->bp), get_pointer(it2->bp));
- 		if (cmpval)
- 			return cmpval;
- 		++it1; ++it2;
-@@ -648,16 +612,10 @@
- 		const function *rh) const
- {
- 
--	//std::cout<<"comparing names: "<<(lh->get_name() < rh->get_name())<<std::endl;
--	/* Reversed ordering on char encoding (i.e. "a" < "b") then length (i.e. "abc" < "abcd")
--	   i.e. "x" > "y" and "xyz" > "xyzt" */
--        if (lh->serial==rh->serial) //return 0;
--	        return compare_same_type_container(lh,rh);
--	//std::cout<<"after if"<<std::endl;
-+	if (lh->serial==rh->serial)
-+		return compare_same_type_container(lh,rh);
- 	return lh->get_name() < rh->get_name() ? 1 : -1;	
- 
--	//return lh->serial < rh->serial ? 1 : -1;
--	//return -lh->compare_same_type(*rh);
- }
- 
- // compare two fderivative objects
-@@ -666,11 +624,11 @@
- int ex_is_greater_degrevlex::compare_same_type_fderivative(const fderivative *lh,
- 		const fderivative *rh) const
- {
--        int cmpval = compare_same_type_function(lh, rh);
-+	int cmpval = compare_same_type_function(lh, rh);
- 	if (cmpval != 0)
--	        return cmpval;
-+		return cmpval;
- 	if (lh->parameter_set != rh->parameter_set)
--	        return lh->parameter_set < rh->parameter_set ? -1 : 1;
-+		return lh->parameter_set < rh->parameter_set ? -1 : 1;
- 	return 0;
- }
- 
-diff --git a/ginac/order.h b/ginac/order.h
---- a/ginac/order.h
-+++ b/ginac/order.h
-@@ -35,7 +35,7 @@
- 
- namespace GiNaC {
- 
--struct ex_is_greater_degrevlex : public std::binary_function<ex, ex, bool> {
-+class ex_is_greater_degrevlex : public std::binary_function<ex, ex, bool> {
- 	const tinfo_t function_id;// = find_tinfo_key("function");
- 	const tinfo_t fderivative_id;// = find_tinfo_key("fderivative");
- 	const tinfo_t power_id;// = find_tinfo_key("power");
-@@ -44,9 +44,9 @@
- 	const tinfo_t add_id;// = find_tinfo_key("add");
- 	const tinfo_t numeric_id;// = find_tinfo_key("numeric");
- 	const tinfo_t constant_id;// = find_tinfo_key("constant");
--	const tinfo_t in_type_id;
- 
--	ex_is_greater_degrevlex(const tinfo_t type_id):
-+	public:
-+	ex_is_greater_degrevlex():
- 		function_id(find_tinfo_key("function")),
- 		fderivative_id(find_tinfo_key("fderivative")),
- 		power_id(find_tinfo_key("power")),
-@@ -54,8 +54,7 @@
- 		mul_id(find_tinfo_key("mul")),
- 		add_id(find_tinfo_key("add")),
- 		numeric_id(find_tinfo_key("numeric")),
--		constant_id(find_tinfo_key("constant")),
--		in_type_id(type_id) {};
-+		constant_id(find_tinfo_key("constant")) {};
- 
- 	bool operator() (const ex &lh, const ex &rh) const;
- 	int compare(const ex &lh, const ex &rh) const;
-@@ -81,51 +80,26 @@
- 	int compare_same_type_function(const function *lh, const function *rh) const;
- 	// fderivative objects
- 	int compare_same_type_fderivative(const fderivative *lh, const fderivative *rh) const;
-+};
- 
--	static const ex_is_greater_degrevlex& in_type(tinfo_t in_type_id) {
--	  static ex_is_greater_degrevlex in_type[2] = {ex_is_greater_degrevlex(&add::tinfo_static),
--						       ex_is_greater_degrevlex(&mul::tinfo_static)};
--	        if (in_type_id == &mul::tinfo_static)
--		        return in_type[1];
--        	return in_type[0];
--	}
-+class ex_is_greater_degrevlex_mul : public ex_is_greater_degrevlex {
-+	int compare_power_symbol(const power *lh, const symbol *rh) const;
-+	int compare_same_type_power(const power *lh, const power *rh) const;
- };
- 
- // We have to define the following class to sort held expressions
- // E.g. 3*x+2*x which does not get simplified to 5*x.
--struct expair_is_greater_degrevlex : public std::binary_function<expair, expair, bool>
-+struct expair_is_greater_degrevlex : \
-+		public std::binary_function<expair, expair, bool>
- {
--        const tinfo_t in_type_id;
--        expair_is_greater_degrevlex(tinfo_t in_type):
--	        in_type_id(in_type) {};
- 	bool operator() (const expair &lh, const expair &rh) const;
--
--	static const expair_is_greater_degrevlex& in_type(tinfo_t in_type_id) {
--        	static expair_is_greater_degrevlex in_type[2] = {expair_is_greater_degrevlex(&add::tinfo_static),expair_is_greater_degrevlex(&add::tinfo_static)};
--		if (in_type_id == &mul::tinfo_static)
--		        return in_type[1];
--        	return in_type[0];
--	}
-+	bool compare_degrees(const expair &lhex, const expair &rhex) const;
- };
- 
--struct expair_rest_is_greater_degrevlex : public std::binary_function<expair, expair, bool>
-+struct expair_is_greater_degrevlex_mul : public expair_is_greater_degrevlex
- {
--        const tinfo_t in_type_id;
--        expair_rest_is_greater_degrevlex(tinfo_t in_type):
--	        in_type_id(in_type) {};
--	bool operator() (const expair &lh, const expair &rh) const {
--	        return ex_is_greater_degrevlex::in_type(in_type_id)(lh.rest, rh.rest);
--	}
--
--	static const expair_rest_is_greater_degrevlex& in_type(tinfo_t in_type_id) {
--        	static expair_rest_is_greater_degrevlex in_type[2] = {expair_rest_is_greater_degrevlex(&add::tinfo_static),
--					expair_rest_is_greater_degrevlex(&mul::tinfo_static)};
--		if (in_type_id == &mul::tinfo_static)
--		        return in_type[1];
--        	return in_type[0];
--	}
-+	bool operator() (const expair &lh, const expair &rh) const;
- };
- 
- } // namespace GiNaC
--
- #endif // ndef __GINAC_ORDER_H__
-diff --git a/ginac/power.h b/ginac/power.h
---- a/ginac/power.h
-+++ b/ginac/power.h
-@@ -41,6 +41,7 @@
- 	GINAC_DECLARE_REGISTERED_CLASS(power, basic)
- 	
- 	friend struct ex_is_greater_degrevlex;
-+	friend struct ex_is_greater_degrevlex_mul;
- 	friend class mul;
- 	
- // member functions

File num_symbol.patch

-# HG changeset patch
-# Parent 7c20443cf59b8ec0d0b58b50a5fd5a69324666d0
-Add method to ex for counting number of symbols in an expression.
-
-diff --git a/ginac/ex.cpp b/ginac/ex.cpp
---- a/ginac/ex.cpp
-+++ b/ginac/ex.cpp
-@@ -21,6 +21,7 @@
-  */
- 
- #include "ex.h"
-+#include "symbol.h"
- #include "add.h"
- #include "mul.h"
- #include "ncmul.h"
-@@ -265,6 +266,18 @@
- 	}
- }
- 
-+size_t ex::nsymbols() const
-+{
-+	int res = 0;
-+	if (is_a<symbol>(*this)) {
-+		res=1;
-+	} else {
-+		for (size_t i=0; i < nops(); i++)
-+			res += op(i).nsymbols();
-+	}
-+	return res;
-+}
-+
- // private
- 
- /** Make this ex writable (if more than one ex handle the same basic) by 
-diff --git a/ginac/ex.h b/ginac/ex.h
---- a/ginac/ex.h
-+++ b/ginac/ex.h
-@@ -143,6 +143,7 @@
- 
- 	// operand access
- 	size_t nops() const { return bp->nops(); }
-+	size_t nsymbols() const;
- 	ex op(size_t i) const { return bp->op(i); }
- 	ex operator[](const ex & index) const { return (*bp)[index]; }
- 	ex operator[](size_t i) const { return (*bp)[i]; }

File numeric_constant_order.patch

-# HG changeset patch
-# Parent 7f070564fc144f791fe3c5208212b523586d4093
-Make numerics and constants print first again, as in sqrt(2)*x.
-
-diff --git a/ginac/order.cpp b/ginac/order.cpp
---- a/ginac/order.cpp
-+++ b/ginac/order.cpp
-@@ -145,17 +145,17 @@
- 	// at present numerics are combined into overall_coefficient
- 	// even when hold parameter is used
- 	} else if (typeid_lh == numeric_id) {
--	 	//print numerics after anything else
-+	 	//print numerics before anything else
-+		return 1;
-+	} else if (typeid_rh == numeric_id) {
-+	 	//print numerics before anything else
- 		return -1;
--	} else if (typeid_rh == numeric_id) {
--	 	//print numerics after anything else
-+	} else if (typeid_lh == constant_id) {
-+	 	//print constants before anything else (but numerics)
- 		return 1;
--	} else if (typeid_lh == constant_id) {
--	 	//print constants after anything else
-+	} else if (typeid_rh == constant_id) {
-+	 	//print constants before anything else (but numerics)
- 		return -1;
--	} else if (typeid_rh == constant_id) {
--	 	//print constants after anything else
--		return 1;
- 	} else if (typeid_lh == fderivative_id) {
- 		//print fderivatives after everything else
- 		return -1;
 trac9880_pynac_order_burcin_original.patch
 trac_9880-pynac_order_jp_new-p2.take2.patch
-cleanup.patch
-expair_is_greater.patch
+trac_9880-cleanup.patch
 num_symbol.patch
-numeric_constant_order.patch
 stable_op.patch
 trac_9880_ginac_infinities_rewrite.patch

File stable_op.patch

-# HG changeset patch
-# Parent 6425f0c03d0211d3a58a7d259045bf9ad6fdeaf3
-Add ex::sorted_op() function which indexes sequence sorted with printing order.
-
-diff --git a/ginac/ex.cpp b/ginac/ex.cpp
---- a/ginac/ex.cpp
-+++ b/ginac/ex.cpp
-@@ -278,6 +278,13 @@
- 	return res;
- }
- 
-+ex ex::sorted_op(size_t i) const
-+{
-+	if (is_a<expairseq>(*this))
-+		return dynamic_cast<const expairseq&>(*bp).stable_op(i);
-+	else
-+		return bp->op(i);
-+}
- // private
- 
- /** Make this ex writable (if more than one ex handle the same basic) by 
-diff --git a/ginac/ex.h b/ginac/ex.h
---- a/ginac/ex.h
-+++ b/ginac/ex.h
-@@ -145,6 +145,7 @@
- 	size_t nops() const { return bp->nops(); }
- 	size_t nsymbols() const;
- 	ex op(size_t i) const { return bp->op(i); }
-+	ex sorted_op(size_t i) const;
- 	ex operator[](const ex & index) const { return (*bp)[index]; }
- 	ex operator[](size_t i) const { return (*bp)[i]; }
- 	ex & let_op(size_t i);
-diff --git a/ginac/expairseq.cpp b/ginac/expairseq.cpp
---- a/ginac/expairseq.cpp
-+++ b/ginac/expairseq.cpp
-@@ -308,6 +308,16 @@
- 	return overall_coeff;
- }
- 
-+ex expairseq::stable_op(size_t i) const
-+{
-+	if (i < seq.size()) {
-+		const epvector* sorted_seq = get_sorted_seq();
-+		return recombine_pair_to_ex((*sorted_seq)[i]);
-+	}
-+	GINAC_ASSERT(!overall_coeff.is_equal(default_overall_coeff()));
-+	return overall_coeff;
-+}
-+
- ex expairseq::map(map_function &f) const
- {
- 	std::auto_ptr<epvector> v(new epvector);
-diff --git a/ginac/expairseq.h b/ginac/expairseq.h
---- a/ginac/expairseq.h
-+++ b/ginac/expairseq.h
-@@ -83,6 +83,7 @@
- 	bool info(unsigned inf) const;
- 	size_t nops() const;
- 	ex op(size_t i) const;
-+	virtual ex stable_op(size_t i) const;
- 	ex map(map_function & f) const;
- 	ex eval(int level=0) const;
- 	ex to_rational(exmap & repl) const;
-@@ -91,6 +92,7 @@
- 	ex subs(const exmap & m, unsigned options = 0) const;
- 	ex conjugate() const;
- 	numeric calc_total_degree() const;
-+	virtual const epvector* get_sorted_seq() const { return &seq; }
- protected:
- 	bool is_equal_same_type(const basic & other) const;
- 	unsigned return_type() const;

File trac_9880-cleanup.patch

+# HG changeset patch
+# Parent 9f0de2356c31320bc122ffe6e95766f10bc54fff
+Clean up new order functions. Rename to print_order_*.
+Make numerics and constants print first again, as in sqrt(2)*x.
+
+diff --git a/ginac/add.cpp b/ginac/add.cpp
+--- a/ginac/add.cpp
++++ b/ginac/add.cpp
+@@ -30,6 +30,7 @@
+ #include "ncmul.h"
+ #include "constant.h"
+ #include "compiler.h"
++#include "order.h"
+ 
+ #include <sstream>
+ #include <iostream>
+@@ -707,4 +708,15 @@
+ 	return (new add(vp, overall_coeff))->setflag(status_flags::dynallocated | (options == 0 ? status_flags::expanded : 0));
+ }
+ 
++epvector* add::get_sorted_seq() const
++{
++	if (!this->seq_sorted) {
++		seq_sorted = new epvector(seq.size());
++		partial_sort_copy(seq.begin(), seq.end(),
++				seq_sorted->begin(), seq_sorted->end(),
++				print_order_pair());
++	}
++	return seq_sorted;
++}
++
+ } // namespace GiNaC
+diff --git a/ginac/add.h b/ginac/add.h
+--- a/ginac/add.h
++++ b/ginac/add.h
+@@ -63,6 +63,7 @@
+ 	ex imag_part() const;
+ 	exvector get_free_indices() const;
+ 	ex eval_ncmul(const exvector & v) const;
++	epvector* get_sorted_seq() const;
+ protected:
+ 	ex derivative(const symbol & s) const;
+ 	unsigned return_type() const;
+diff --git a/ginac/basic.h b/ginac/basic.h
+--- a/ginac/basic.h
++++ b/ginac/basic.h
+@@ -108,8 +108,8 @@
+ 	GINAC_DECLARE_REGISTERED_CLASS_NO_CTORS(basic, void)
+ 	
+ 	friend class ex;
+-	friend struct ex_is_greater_degrevlex;
+-	friend struct expair_is_greater_degrevlex;
++	friend struct print_order;
++	friend struct print_order_pair;
+ 	// default constructor, destructor, copy constructor and assignment operator
+ protected:
+ 	basic() : tinfo_key(&tinfo_static), flags(0) {}
+diff --git a/ginac/container.h b/ginac/container.h
+--- a/ginac/container.h
++++ b/ginac/container.h
+@@ -129,7 +129,7 @@
+ class container : public basic, public container_storage<C> {
+ 	GINAC_DECLARE_REGISTERED_CLASS(container, basic)
+ 
+-	friend struct ex_is_greater_degrevlex;
++	friend struct print_order;
+ 
+ protected:
+ 	typedef typename container_storage<C>::STLT STLT;
+diff --git a/ginac/ex.h b/ginac/ex.h
+--- a/ginac/ex.h
++++ b/ginac/ex.h
+@@ -83,8 +83,9 @@
+ 	template<class T> friend inline bool is_a(const ex &);
+ 	template<class T> friend inline bool is_exactly_a(const ex &);
+ 	
+-	friend struct ex_is_greater_degrevlex;
+-	friend struct expair_is_greater_degrevlex;
++	friend struct print_order;
++	friend struct print_order_mul;
++	friend struct print_order_pair;
+ 	// default constructor, copy constructor and assignment operator
+ public:
+ 	ex() throw();
+diff --git a/ginac/expairseq.cpp b/ginac/expairseq.cpp
+--- a/ginac/expairseq.cpp
++++ b/ginac/expairseq.cpp
+@@ -32,7 +32,6 @@
+ #include "utils.h"
+ #include "indexed.h"
+ #include "constant.h"
+-#include "order.h"
+ 
+ #include <iostream>
+ #include <algorithm>
+@@ -1171,17 +1170,6 @@
+ 	std::sort(seq.begin(), seq.end(), expair_rest_is_less());
+ }
+ 
+-epvector* expairseq::get_sorted_seq() const
+-{
+-	if (!this->seq_sorted) {
+-		seq_sorted = new epvector(seq.size());
+-		partial_sort_copy(seq.begin(), seq.end(),
+-				seq_sorted->begin(), seq_sorted->end(),
+-				  expair_is_greater_degrevlex::in_type(tinfo()));
+-	}
+-	return seq_sorted;
+-}
+-
+ 
+ /** Compact a presorted expairseq by combining all matching expairs to one
+  *  each.  On an add object, this is responsible for 2*x+3*x+y -> 5*x+y, for
+diff --git a/ginac/expairseq.h b/ginac/expairseq.h
+--- a/ginac/expairseq.h
++++ b/ginac/expairseq.h
+@@ -68,13 +68,14 @@
+ {
+ 	GINAC_DECLARE_REGISTERED_CLASS(expairseq, basic)
+ 
+-	friend struct ex_is_greater_degrevlex;
++	friend struct print_order;
+ 	// other constructors
+ public:
+ 	expairseq(const ex & lh, const ex & rh);
+ 	expairseq(const exvector & v);
+ 	expairseq(const epvector & v, const ex & oc, bool do_index_renaming = false);
+ 	expairseq(std::auto_ptr<epvector>, const ex & oc, bool do_index_renaming = false);
++	~expairseq() { if (!seq_sorted) delete seq_sorted; }
+ 	
+ 	// functions overriding virtual functions from base classes
+ public:
+@@ -133,7 +134,6 @@
+ 	void make_flat(const epvector & v, bool do_index_renaming = false);
+ 	void canonicalize();
+ 	void combine_same_terms_sorted_seq();
+-	epvector* get_sorted_seq() const;
+ #if EXPAIRSEQ_USE_HASHTAB
+ 	void combine_same_terms();
+ 	unsigned calc_hashtabsize(unsigned sz) const;
+diff --git a/ginac/fderivative.h b/ginac/fderivative.h
+--- a/ginac/fderivative.h
++++ b/ginac/fderivative.h
+@@ -40,7 +40,7 @@
+ {
+ 	GINAC_DECLARE_REGISTERED_CLASS(fderivative, function)
+ 
+-	friend struct ex_is_greater_degrevlex;
++	friend struct print_order;
+ 	// other constructors
+ public:
+ 	/** Construct derivative with respect to one parameter.
+diff --git a/ginac/function.h b/ginac/function.h
+--- a/ginac/function.h
++++ b/ginac/function.h
+@@ -310,7 +310,7 @@
+ {
+ 	GINAC_DECLARE_REGISTERED_CLASS(function, exprseq)
+ 
+-	friend struct ex_is_greater_degrevlex;
++	friend struct print_order;
+ 	// CINT has a linking problem
+ #ifndef __MAKECINT__
+ 	friend void ginsh_get_ginac_functions();
+diff --git a/ginac/mul.cpp b/ginac/mul.cpp
+--- a/ginac/mul.cpp
++++ b/ginac/mul.cpp
+@@ -1620,6 +1620,17 @@
+ 	}
+ }
+ 
++epvector* mul::get_sorted_seq() const
++{
++	if (!this->seq_sorted) {
++		seq_sorted = new epvector(seq.size());
++		partial_sort_copy(seq.begin(), seq.end(),
++				seq_sorted->begin(), seq_sorted->end(),
++				print_order_pair_mul());
++	}
++	return seq_sorted;
++}
++
+   
+ //////////
+ // new virtual functions which can be overridden by derived classes
+diff --git a/ginac/mul.h b/ginac/mul.h
+--- a/ginac/mul.h
++++ b/ginac/mul.h
+@@ -33,7 +33,7 @@
+ {
+ 	GINAC_DECLARE_REGISTERED_CLASS(mul, expairseq)
+ 	
+-	friend struct ex_is_greater_degrevlex;
++	friend struct print_order;
+ 	friend class add;
+ 	friend class ncmul;
+ 	friend class power;
+@@ -95,6 +95,7 @@
+ public:
+ 	ex algebraic_subs_mul(const exmap & m, unsigned options) const;
+ 	double total_degree() const;
++	epvector* get_sorted_seq() const;
+ 	//int compare_symbol(const symbol &other) const;
+ 	//int compare_pow(const power &other) const;
+ protected:
+diff --git a/ginac/numeric.h b/ginac/numeric.h
+--- a/ginac/numeric.h
++++ b/ginac/numeric.h
+@@ -207,7 +207,6 @@
+   class numeric : public basic
+   {
+     GINAC_DECLARE_REGISTERED_CLASS(numeric, basic)
+-    //friend struct expair_is_greater_degrevlex;
+     // member functions
+ 	
+     // other constructors
+diff --git a/ginac/order.cpp b/ginac/order.cpp
+--- a/ginac/order.cpp
++++ b/ginac/order.cpp
+@@ -30,78 +30,75 @@
+ 
+ namespace GiNaC {
+ 
+-/** What SAGE does for printing:
+- To print multivariate polynomials, SAGE uses "reversed" (bigger terms first) 'degrevlex' order
+- i.e. "reversed" graded reversed lexicographic order, the lexicographic order is defined below.
+- The variables are ordered according to their "creation" order
+- i.e PR.<x,a,t> = PolynomialRing(QQ) gives x > a > t
+- and x+t+a is printed x + a + t
++double numeric_to_double(const numeric& exp)
++{
++	if (exp.is_real())
++		return exp.to_double();
++	else 
++		return std::sqrt(std::pow(exp.real().to_double(), 2) + 
++				std::pow(exp.imag().to_double(), 2));
++}
++
++/** What Sage does for printing:
++ To print multivariate polynomials, SAGE uses "reversed" (bigger terms first)
++ 'degrevlex' order i.e. "reversed" graded reversed lexicographic order, the
++ lexicographic order is defined below.  The variables are ordered according to
++ their "creation" order i.e PR.<x,a,t> = PolynomialRing(QQ) gives x > a > t and
++ x+t+a is printed x + a + t
+ */
+ 
+-/** Returns true if lhex > rhex for 'degrevlex' of SAGE
++/** Returns true if lhex > rhex for 'degrevlex' of Sage
+  i.e. graded reversed lexicographic order
+ */
+-bool ex_is_greater_degrevlex::operator() (const ex &lhex, const ex &rhex) const {
+-        //std::cout<<"in ex_is_greater_degrevlex::operator()"<<std::endl;
+-	//std::cout<<"lh: ";
+-	//lhex.dbgprint();
+-	//std::cout<<std::endl<<"rh: ";
+-	//rhex.dbgprint();
+-	//std::cout<<std::endl;
++bool print_order::operator() (const ex &lhex, const ex &rhex) const {
+ 	const basic *lh = get_pointer(lhex.bp);
+ 	const basic *rh = get_pointer(rhex.bp);
+ 	return compare(lh, rh) == 1;
+ }
+ 
+-int ex_is_greater_degrevlex::compare (const ex &lhex, const ex &rhex) const {
++int print_order::compare (const ex &lhex, const ex &rhex) const {
+ 	const basic *lh = get_pointer(lhex.bp);
+ 	const basic *rh = get_pointer(rhex.bp);
+ 	return compare(lh, rh);
+ }
+ 
+-bool expair_is_greater_degrevlex::operator() (const expair &lhex, const expair &rhex) const
++bool print_order_pair_mul::operator() (const expair &lhex, const expair &rhex) const
++{
++	int cmpval = print_order_mul().compare(lhex.rest, rhex.rest);
++	if (cmpval != 0) {
++		return cmpval == 1;
++	}
++	return compare_degrees(lhex, rhex);
++}
++
++bool print_order_pair::operator() (const expair &lhex, const expair &rhex) const
+ {
+ 	const basic *lh = get_pointer(lhex.rest.bp);
+ 	const basic *rh = get_pointer(rhex.rest.bp);
+ 
+ 	// compare rests
+-        int cmpval = ex_is_greater_degrevlex::in_type(in_type_id).compare(lh, rh);
++	int cmpval = print_order().compare(lh, rh);
+ 	if (cmpval != 0) {
+-	        return cmpval == 1;
++		return cmpval == 1;
+ 	}
++	return compare_degrees(lhex, rhex);
++}
+ 
++bool print_order_pair::compare_degrees(const expair &lhex,
++		const expair &rhex) const
++{
+ 	// compare coeffs which are numerics
+-	//std::cout << "comparing coeffs" << std::endl;
+-	//lhex.coeff.bp->dbgprinttree();
+-	//rhex.coeff.bp->dbgprinttree();
+-	numeric lh_coeff = ex_to<numeric>(lhex.coeff);
+-	numeric rh_coeff = ex_to<numeric>(rhex.coeff);
+-	// strange...
+-	//std::cout<<lh->compare_same_type(rh)<<std::endl;
+-	//return lh->compare_same_type(rh) == 1;
+-	double lh_deg = 1;
+-	double rh_deg = 1;
+-		if (lh_coeff.is_real()) {
+-			lh_deg = lh_coeff.to_double();
+-		} else {
+-			lh_deg = std::sqrt(std::pow(lh_coeff.real().to_double(), 2) + 
+-					std::pow(lh_coeff.imag().to_double(), 2));
+-		}
+-		if (rh_coeff.is_real()) {
+-			rh_deg = rh_coeff.to_double();
+-		} else {
+-			rh_deg = std::sqrt(std::pow(rh_coeff.real().to_double(), 2) + 
+-					std::pow(rh_coeff.imag().to_double(), 2));
+-		}
++	double lh_deg = numeric_to_double(ex_to<numeric>(lhex.coeff));
++	double rh_deg = numeric_to_double(ex_to<numeric>(rhex.coeff));
+ 
+ 	return lh_deg > rh_deg;	
+ }
+ 
+ /** Comparison functions:
+- They should implement 'degrevlex' of SAGE
++ They should implement 'degrevlex' of Sage
+  i.e. graded reversed lexicographic order
+- The lexicographic order should depend on the "creation" order of variables ?
+- Or should it be a "natural" one on strings : a > b > ... > x > y > ... ?
++ The lexicographic order is the "natural" one on strings:
++ a > b > ... > x > y > ... ?
+ */
+ 
+ /** Return values for comparison functions :
+@@ -112,7 +109,7 @@
+  as <=> in Perl and GiNaC internal functions
+ */
+ 
+-int ex_is_greater_degrevlex::compare(const basic *lh, const basic *rh) const {
++int print_order::compare(const basic *lh, const basic *rh) const {
+ 	const tinfo_t typeid_lh = lh->tinfo();
+ 	const tinfo_t typeid_rh = rh->tinfo();
+ 
+@@ -142,23 +139,23 @@
+ 					static_cast<const fderivative*>(lh),
+ 					static_cast<const fderivative*>(rh));
+ 		} else {
+-		        // using GiNaC functions by default
++			// using GiNaC functions by default
+ 			return lh->compare_same_type(*rh);
+ 		}
+ 	// at present numerics are combined into overall_coefficient
+ 	// even when hold parameter is used
+ 	} else if (typeid_lh == numeric_id) {
+-	 	//print numerics after anything else
+-	        return -1;
++	 	//print numerics before anything else
++		return 1;
+ 	} else if (typeid_rh == numeric_id) {
+-	 	//print numerics after anything else
+-	        return 1;
++	 	//print numerics before anything else
++		return -1;
+ 	} else if (typeid_lh == constant_id) {
+-	 	//print constants after anything else
+-	        return -1;
++	 	//print constants before anything else (but numerics)
++		return 1;
+ 	} else if (typeid_rh == constant_id) {
+-	 	//print constants after anything else
+-	        return 1;
++	 	//print constants before anything else (but numerics)
++		return -1;
+ 	} else if (typeid_lh == fderivative_id) {
+ 		//print fderivatives after everything else
+ 		return -1;
+@@ -167,10 +164,10 @@
+ 		return 1;
+ 	} else if (typeid_lh == function_id) {
+ 		//print functions before fderivatives, after anything else
+-	        return -1;
++		return -1;
+ 	} else if (typeid_rh == function_id) {
+ 		//print functions before fderivatives, after anything else
+-	        return 1;
++		return 1;
+ 	} else if (typeid_lh == mul_id) {
+ 		if (typeid_rh == power_id) {
+ 			return compare_mul_power(
+@@ -227,28 +224,24 @@
+ 					static_cast<const power*>(rh),
+ 					static_cast<const symbol*>(lh));
+ 		}
+-        }
+-	std::cout<<"comparing typeid's"<<std::endl;
+-	return (typeid_lh<typeid_rh ? -1 : 1);
++	}
++	throw (std::runtime_error("comparing typeid's"));
+ }
+ 
+ // compare a mul and a symbol objects
+ // same behavior within mul and add objects:
+ // the total degree of the symbol is 1
+ // check total degree of mul then compare smallest item to symbol
+-int ex_is_greater_degrevlex::compare_mul_symbol(const mul *lh,
++int print_order::compare_mul_symbol(const mul *lh,
+ 		const symbol *rh) const
+ {
+-  //std::cout<<"comparing mul and symbol"<<std::endl;
+-  //lh->dbgprint();
+-  //rh->dbgprint();
+ 	int cmpval;
+ 
+ 	double tdeg;
+ 	tdeg = lh->total_degree();
+ 
+ 	if (tdeg != 1)
+-	        return tdeg > 1 ? 1 : -1;
++		return tdeg > 1 ? 1 : -1;
+ 	
+ 	const expair smallest_item = lh->get_sorted_seq()->back();
+ 
+@@ -271,13 +264,12 @@
+ 	return 1;
+ }
+ 
+-
+ // compare a mul and a pow objects
+-// same behavior wihtin mul and add objects:
++// same behavior within mul and add objects:
+ // first we compare total degrees
+ // if equal we compare the smallest basis in the sequence to the basis in other
+ // then their exponents
+-int ex_is_greater_degrevlex::compare_mul_power(const mul *lh,
++int print_order::compare_mul_power(const mul *lh,
+ 		const power *rh) const
+ {
+ 	int cmpval;
+@@ -286,13 +278,7 @@
+ 	double rh_deg = 1;
+ 	numeric rh_exp;
+ 	if (is_a<numeric>(rh->exponent)) {
+-		rh_exp = ex_to<numeric>(rh->exponent);
+-		if (rh_exp.is_real()) {
+-			rh_deg = rh_exp.to_double();
+-		} else {
+-			rh_deg = std::sqrt(std::pow(rh_exp.real().to_double(), 2) + 
+-					std::pow(rh_exp.imag().to_double(), 2));
+-		}
++		rh_deg = numeric_to_double(ex_to<numeric>(rh->exponent));
+ 	}
+ 	if (rh_deg != lh_deg)
+ 		return lh_deg < rh_deg ? -1 : 1;
+@@ -310,7 +296,7 @@
+ 
+ 	// compare exponents
+ 	cmpval = -compare(get_pointer(smallest_item.coeff.bp),
+-		        get_pointer(rh->exponent.bp));
++			get_pointer(rh->exponent.bp));
+ 	if (cmpval != 0) {
+ 		return cmpval;
+ 	}
+@@ -328,18 +314,14 @@
+ // if equal we compare the basis of the smallest items
+ // then their exponents
+ // and so on
+-int ex_is_greater_degrevlex::compare_same_type_mul(const mul *lh,
++int print_order::compare_same_type_mul(const mul *lh,
+ 		const mul *rh) const
+ {
+-  //std::cout<<"comparing mul and mul"<<std::endl;
+-  //lh->dbgprint();
+-  //rh->dbgprint();
+ 	int cmpval;
+ 	
+ 	// compare total degrees
+ 	double lh_deg = lh->total_degree();
+ 	double rh_deg = rh->total_degree();
+-        //std::cout<<"degree "<<lh_deg<<" and "<<rh_deg<<std::endl;
+ 	if (lh_deg != rh_deg)
+ 		return lh_deg < rh_deg ? -1 : 1;
+ 	
+@@ -353,13 +335,15 @@
+ 
+ 	for (; (cit1!=last1)&&(cit2!=last2); ++cit1, ++cit2) {
+ 		// compare bases
+-		cmpval = compare(get_pointer(cit1->rest.bp),get_pointer(cit2->rest.bp));
++		cmpval = compare(get_pointer(cit1->rest.bp),
++				get_pointer(cit2->rest.bp));
+ 		if (cmpval != 0) {
+ 			return cmpval;
+ 		}
+ 
+ 		// compare exponents
+-	        cmpval = -compare(get_pointer(cit1->coeff.bp),get_pointer(cit2->coeff.bp));
++		cmpval = -compare(get_pointer(cit1->coeff.bp),
++				get_pointer(cit2->coeff.bp));
+ 		if (cmpval != 0) {
+ 			return cmpval;
+ 		}
+@@ -381,18 +365,15 @@
+ // compare an add and a symbol objects
+ // same behavior wihtin mul and add objects:
+ // the coefficient of the symbol is 1
+-int ex_is_greater_degrevlex::compare_add_symbol(const add *lh,
++int print_order::compare_add_symbol(const add *lh,
+ 		const symbol *rh) const
+ {
+-  //std::cout<<"comparing add and symbol"<<std::endl;
+-  //lh->dbgprint();
+-  //rh->dbgprint();
+ 	int cmpval;
+ 
+ 	const expair biggest_item = lh->get_sorted_seq()->front();
+ 
+ 	// compare bases
+-	cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(biggest_item.rest.bp), rh);
++	cmpval = print_order().compare(get_pointer(biggest_item.rest.bp), rh);
+ 	if (cmpval != 0) {
+ 		return cmpval;
+ 	}
+@@ -413,14 +394,14 @@
+ // compare an add and a mul objects
+ // same behavior within mul and add objects:
+ // the coefficient of the mul object is 1
+-int ex_is_greater_degrevlex::compare_add_mul(const add *lh,
++int print_order::compare_add_mul(const add *lh,
+ 		const mul *rh) const
+ {
+-        int cmpval;
++	int cmpval;
+ 	const expair biggest_item = lh->get_sorted_seq()->front();
+ 
+ 	// compare bases
+-	cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(biggest_item.rest.bp), rh);
++	cmpval = print_order().compare(get_pointer(biggest_item.rest.bp), rh);
+ 	if (cmpval != 0) {
+ 		return cmpval;
+ 	}
+@@ -441,14 +422,14 @@
+ // compare an add and a pow objects
+ // same behavior wihtin mul and add objects:
+ // the coefficient of the power object is 1
+-int ex_is_greater_degrevlex::compare_add_power(const add *lh,
++int print_order::compare_add_power(const add *lh,
+ 		const power *rh) const
+ {
+-        int cmpval;
++	int cmpval;
+ 	const expair biggest_item = lh->get_sorted_seq()->front();
+ 
+ 	// compare bases
+-	cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(biggest_item.rest.bp), rh);
++	cmpval = print_order().compare(get_pointer(biggest_item.rest.bp), rh);
+ 	if (cmpval != 0) {
+ 		return cmpval;
+ 	}
+@@ -471,7 +452,7 @@
+ // first we compare the basis of the biggest items
+ // then their coefficients
+ // and so on
+-int ex_is_greater_degrevlex::compare_same_type_add(const add *lh,
++int print_order::compare_same_type_add(const add *lh,
+ 		const add *rh) const
+ {
+ 	int cmpval;
+@@ -485,13 +466,13 @@
+ 
+ 	for (; (cit1!=last1)&&(cit2!=last2); ++cit1, ++cit2) {
+ 		// compare bases
+-		cmpval = ex_is_greater_degrevlex::in_type(&add::tinfo_static).compare(get_pointer(cit1->rest.bp),get_pointer(cit2->rest.bp));
++		cmpval = print_order().compare(get_pointer(cit1->rest.bp),get_pointer(cit2->rest.bp));
+ 		if (cmpval != 0) {
+ 			return cmpval;
+ 		}
+ 
+ 		// compare coefficients
+-	        cmpval = compare(get_pointer(cit1->coeff.bp),get_pointer(cit2->coeff.bp));
++		cmpval = compare(get_pointer(cit1->coeff.bp),get_pointer(cit2->coeff.bp));
+ 		if (cmpval != 0) {
+ 			return cmpval;
+ 		}
+@@ -518,32 +499,27 @@
+ // in add object:
+ // first exponents
+ // then bases
+-int ex_is_greater_degrevlex::compare_power_symbol(const power *lh,
++int print_order_mul::compare_power_symbol(const power *lh,
+ 		const symbol *rh) const
+ {
+-        int cmpval;
++	int cmpval;
++	cmpval = compare(get_pointer(lh->basis.bp), rh);
++	if (cmpval != 0)
++		  return cmpval;
++	
++	cmpval = compare(get_pointer(lh->exponent.bp), _num1_p);
++	
++	return cmpval;
++}
+ 
+-	//std::cout<<"in compare_power_symbol"<<std::endl;
+-        if (in_type_id == &mul::tinfo_static) {
+-	        cmpval = compare(get_pointer(lh->basis.bp), rh);
+-		if (cmpval != 0)
+-        		  return cmpval;
+-		
+-		cmpval = compare(get_pointer(lh->exponent.bp), _num1_p);
+-		
+-		return cmpval;
+-	}
++int print_order::compare_power_symbol(const power *lh,
++		const symbol *rh) const
++{
++	int cmpval;
+ 
+ 	// We are in an add object
+ 	if (is_a<numeric>(lh->exponent)) {
+-		numeric lh_exp = ex_to<numeric>(lh->exponent);
+-		double lh_deg;
+-		if (lh_exp.is_real()) {
+-			lh_deg = lh_exp.to_double();
+-		} else {
+-			lh_deg = std::sqrt(std::pow(lh_exp.real().to_double(), 2) + 
+-					std::pow(lh_exp.imag().to_double(), 2));
+-		}
++		double lh_deg = numeric_to_double(ex_to<numeric>(lh->exponent));
+ 		if (lh_deg != 1)
+ 			return lh_deg < 1 ? -1 : 1;
+ 	}
+@@ -561,90 +537,66 @@
+ // in add object:
+ // first exponents
+ // then bases
+-int ex_is_greater_degrevlex::compare_same_type_power(const power *lh,
+-	        const power *rh) const
++int print_order_mul::compare_same_type_power(const power *lh,
++		const power *rh) const
+ {
+-        int cmpval;
+-        if (in_type_id == &mul::tinfo_static) {
+-	        cmpval = compare(get_pointer(lh->basis.bp), get_pointer(rh->basis.bp));
+-		if (cmpval != 0)
+-        		  return cmpval;
+-		cmpval = compare(get_pointer(lh->exponent.bp), get_pointer(rh->exponent.bp));
+-		
+-		return cmpval;
+-	}
++	int cmpval;
++	cmpval = compare(get_pointer(lh->basis.bp), get_pointer(rh->basis.bp));
++	if (cmpval != 0)
++		  return cmpval;
++	cmpval = compare(get_pointer(lh->exponent.bp), get_pointer(rh->exponent.bp));
++	
++	return cmpval;
++}
++int print_order::compare_same_type_power(const power *lh,
++		const power *rh) const
++{
++	int cmpval;
+ 
+ 	double lh_deg = 1;
+ 	double rh_deg = 1;
+ 	if (is_a<numeric>(lh->exponent)) {
+-		numeric lh_exp = ex_to<numeric>(lh->exponent);
+-		if (lh_exp.is_real()) {
+-			lh_deg = lh_exp.to_double();
+-		} else {
+-			lh_deg = std::sqrt(std::pow(lh_exp.real().to_double(), 2) + 
+-					std::pow(lh_exp.imag().to_double(), 2));
+-		}
++		lh_deg = numeric_to_double(ex_to<numeric>(lh->exponent));
+ 	}
+ 	if (is_a<numeric>(rh->exponent)) {
+-		numeric rh_exp = ex_to<numeric>(rh->exponent);
+-		if (rh_exp.is_real()) {
+-			rh_deg = rh_exp.to_double();
+-		} else {
+-			rh_deg = std::sqrt(std::pow(rh_exp.real().to_double(), 2) + 
+-					std::pow(rh_exp.imag().to_double(), 2));
+-		}
++		rh_deg = numeric_to_double(ex_to<numeric>(rh->exponent));
+ 	}
+ 	if (rh_deg != lh_deg)
+ 		return lh_deg < rh_deg ? -1 : 1;
+ 
+-	cmpval = compare(get_pointer(lh->basis.bp),
+-			 get_pointer(rh->basis.bp));
++	cmpval = compare(get_pointer(lh->basis.bp), get_pointer(rh->basis.bp));
+ 	if (cmpval != 0)
+-        	return cmpval;
++		return cmpval;
+ 
+ 	if (is_a<numeric>(lh->exponent) && is_a<numeric>(rh->exponent))
+-	        return 0;
+-	return compare(get_pointer(lh->exponent.bp), get_pointer(rh->exponent.bp));
++		return 0;
++	return compare(get_pointer(lh->exponent.bp),
++			get_pointer(rh->exponent.bp));
+ }
+ 
+ // compare two symbol objects
+ // same behavior wihtin mul and add objects:
+ // we compare names
+-int ex_is_greater_degrevlex::compare_same_type_symbol(const symbol *lh,
++int print_order::compare_same_type_symbol(const symbol *lh,
+ 		const symbol *rh) const
+ {
+-        //std::cout<<"in compare_same symbol"<<std::endl;
+-	//std::cout<<"lh: ";
+-	//lh->dbgprint();
+-	//std::cout<<"rh: ";
+-	//rh->dbgprint();
+-	//std::cout<<std::endl;
+-
+-        // Weird because Sage sorts based on creation order.
+-        // SAGE/Pynac: Sorting based on creation order doesn't work for Sage.
+-        // instead we sort on variable name. -- William Stein
+-
+-	//std::cout<<"comparing names: "<<(lh->name < rh->name)<<std::endl;
+-	/* Reversed ordering on char encoding (i.e. "a" < "b") then length (i.e. "abc" < "abcd")
++	/* Reversed ordering on char encoding (i.e. "a" < "b") then length
++	 * (i.e. "abc" < "abcd")
+ 	   i.e. "x" > "y" and "xyz" > "xyzt" */
+ 	if (lh->serial==rh->serial) return 0;
+-	//std::cout<<"after if"<<std::endl;
+ 	return lh->name < rh->name ? 1 : -1;
+-
+-	//return lh->serial < rh->serial ? 1 : -1;
+-	//return -lh->compare_same_type(*rh);
+ }
+ 
+ // compare two containers of the same type
+ template <template <class T, class = std::allocator<T> > class C>
+-int ex_is_greater_degrevlex::compare_same_type_container(const container<C> *lh,
++int print_order::compare_same_type_container(const container<C> *lh,
+ 							 const container<C> *rh) const
+ {
+-        typename C<ex>::const_iterator it1 = lh->seq.begin(), it1end = lh->seq.end(),
+-	                      it2 = rh->seq.begin(), it2end = rh->seq.end();
++	typename C<ex>::const_iterator it1 = lh->seq.begin(), it1end = lh->seq.end(),
++			      it2 = rh->seq.begin(), it2end = rh->seq.end();
+ 
+ 	while (it1 != it1end && it2 != it2end) {
+-	        int cmpval = compare(get_pointer(it1->bp), get_pointer(it2->bp));
++		int cmpval = compare(get_pointer(it1->bp), get_pointer(it2->bp));
+ 		if (cmpval)
+ 			return cmpval;
+ 		++it1; ++it2;
+@@ -656,33 +608,27 @@
+ // compare two function objects
+ // same behavior wihtin mul and add objects:
+ // we compare names
+-int ex_is_greater_degrevlex::compare_same_type_function(const function *lh,
++int print_order::compare_same_type_function(const function *lh,
+ 		const function *rh) const
+ {
+ 
+-	//std::cout<<"comparing names: "<<(lh->get_name() < rh->get_name())<<std::endl;
+-	/* Reversed ordering on char encoding (i.e. "a" < "b") then length (i.e. "abc" < "abcd")
+-	   i.e. "x" > "y" and "xyz" > "xyzt" */
+-        if (lh->serial==rh->serial) //return 0;
+-	        return compare_same_type_container(lh,rh);
+-	//std::cout<<"after if"<<std::endl;
++	if (lh->serial==rh->serial)
++		return compare_same_type_container(lh,rh);
+ 	return lh->get_name() < rh->get_name() ? 1 : -1;	
+ 
+-	//return lh->serial < rh->serial ? 1 : -1;
+-	//return -lh->compare_same_type(*rh);
+ }
+ 
+ // compare two fderivative objects
+ // same behavior wihtin mul and add objects:
+ // we compare names
+-int ex_is_greater_degrevlex::compare_same_type_fderivative(const fderivative *lh,
++int print_order::compare_same_type_fderivative(const fderivative *lh,
+ 		const fderivative *rh) const
+ {
+-        int cmpval = compare_same_type_function(lh, rh);
++	int cmpval = compare_same_type_function(lh, rh);
+ 	if (cmpval != 0)
+-	        return cmpval;
++		return cmpval;
+ 	if (lh->parameter_set != rh->parameter_set)
+-	        return lh->parameter_set < rh->parameter_set ? -1 : 1;
++		return lh->parameter_set < rh->parameter_set ? 1 : -1;
+ 	return 0;
+ }
+ 
+diff --git a/ginac/order.h b/ginac/order.h
+--- a/ginac/order.h
++++ b/ginac/order.h
+@@ -35,7 +35,7 @@
+ 
+ namespace GiNaC {
+ 
+-struct ex_is_greater_degrevlex : public std::binary_function<ex, ex, bool> {
++class print_order : public std::binary_function<ex, ex, bool> {
+ 	const tinfo_t function_id;// = find_tinfo_key("function");
+ 	const tinfo_t fderivative_id;// = find_tinfo_key("fderivative");
+ 	const tinfo_t power_id;// = find_tinfo_key("power");
+@@ -44,9 +44,9 @@
+ 	const tinfo_t add_id;// = find_tinfo_key("add");
+ 	const tinfo_t numeric_id;// = find_tinfo_key("numeric");
+ 	const tinfo_t constant_id;// = find_tinfo_key("constant");
+-	const tinfo_t in_type_id;
+ 
+-	ex_is_greater_degrevlex(const tinfo_t type_id):
++	public:
++	print_order():
+ 		function_id(find_tinfo_key("function")),
+ 		fderivative_id(find_tinfo_key("fderivative")),
+ 		power_id(find_tinfo_key("power")),
+@@ -54,8 +54,7 @@
+ 		mul_id(find_tinfo_key("mul")),
+ 		add_id(find_tinfo_key("add")),
+ 		numeric_id(find_tinfo_key("numeric")),
+-		constant_id(find_tinfo_key("constant")),
+-		in_type_id(type_id) {};
++		constant_id(find_tinfo_key("constant")) {};
+ 
+ 	bool operator() (const ex &lh, const ex &rh) const;
+ 	int compare(const ex &lh, const ex &rh) const;
+@@ -70,8 +69,8 @@
+ 	int compare_add_mul(const add *lh, const mul *rh) const;
+ 	int compare_same_type_add(const add *lh, const add *rh) const;
+ 	// power objects
+-	int compare_power_symbol(const power *lh, const symbol *rh) const;
+-	int compare_same_type_power(const power *lh, const power *rh) const;
++	virtual int compare_power_symbol(const power *lh, const symbol *rh) const;
++	virtual int compare_same_type_power(const power *lh, const power *rh) const;
+ 	// symbol objects
+ 	int compare_same_type_symbol(const symbol *lh, const symbol *rh) const;
+ 	// container objects
+@@ -81,51 +80,26 @@
+ 	int compare_same_type_function(const function *lh, const function *rh) const;
+ 	// fderivative objects
+ 	int compare_same_type_fderivative(const fderivative *lh, const fderivative *rh) const;
++};
+ 
+-	static const ex_is_greater_degrevlex& in_type(tinfo_t in_type_id) {
+-	  static ex_is_greater_degrevlex in_type[2] = {ex_is_greater_degrevlex(&add::tinfo_static),
+-						       ex_is_greater_degrevlex(&mul::tinfo_static)};
+-	        if (in_type_id == &mul::tinfo_static)
+-		        return in_type[1];
+-        	return in_type[0];
+-	}
++class print_order_mul : public print_order {
++	int compare_power_symbol(const power *lh, const symbol *rh) const;
++	int compare_same_type_power(const power *lh, const power *rh) const;
+ };
+ 
+ // We have to define the following class to sort held expressions
+ // E.g. 3*x+2*x which does not get simplified to 5*x.
+-struct expair_is_greater_degrevlex : public std::binary_function<expair, expair, bool>
++struct print_order_pair : \
++		public std::binary_function<expair, expair, bool>
+ {
+-        const tinfo_t in_type_id;
+-        expair_is_greater_degrevlex(tinfo_t in_type):
+-	        in_type_id(in_type) {};
+ 	bool operator() (const expair &lh, const expair &rh) const;
+-
+-	static const expair_is_greater_degrevlex& in_type(tinfo_t in_type_id) {
+-        	static expair_is_greater_degrevlex in_type[2] = {expair_is_greater_degrevlex(&add::tinfo_static),expair_is_greater_degrevlex(&add::tinfo_static)};
+-		if (in_type_id == &mul::tinfo_static)
+-		        return in_type[1];
+-        	return in_type[0];
+-	}
++	bool compare_degrees(const expair &lhex, const expair &rhex) const;
+ };
+ 
+-struct expair_rest_is_greater_degrevlex : public std::binary_function<expair, expair, bool>
++struct print_order_pair_mul : public print_order_pair
+ {
+-        const tinfo_t in_type_id;
+-        expair_rest_is_greater_degrevlex(tinfo_t in_type):
+-	        in_type_id(in_type) {};
+-	bool operator() (const expair &lh, const expair &rh) const {
+-	        return ex_is_greater_degrevlex::in_type(in_type_id)(lh.rest, rh.rest);
+-	}
+-
+-	static const expair_rest_is_greater_degrevlex& in_type(tinfo_t in_type_id) {
+-        	static expair_rest_is_greater_degrevlex in_type[2] = {expair_rest_is_greater_degrevlex(&add::tinfo_static),
+-					expair_rest_is_greater_degrevlex(&mul::tinfo_static)};
+-		if (in_type_id == &mul::tinfo_static)
+-		        return in_type[1];
+-        	return in_type[0];
+-	}
++	bool operator() (const expair &lh, const expair &rh) const;
+ };
+ 
+ } // namespace GiNaC
+-
+ #endif // ndef __GINAC_ORDER_H__
+diff --git a/ginac/power.h b/ginac/power.h
+--- a/ginac/power.h
++++ b/ginac/power.h
+@@ -40,7 +40,8 @@
+ {
+ 	GINAC_DECLARE_REGISTERED_CLASS(power, basic)
+ 	
+-	friend struct ex_is_greater_degrevlex;
++	friend struct print_order;
++	friend struct print_order_mul;
+ 	friend class mul;
+ 	
+ // member functions
+diff --git a/ginac/symbol.h b/ginac/symbol.h
+--- a/ginac/symbol.h
++++ b/ginac/symbol.h
+@@ -43,7 +43,7 @@
+ 
+ 	friend class realsymbol;
+ 	friend class possymbol;
+-	friend struct ex_is_greater_degrevlex;
++	friend struct print_order;
+ 
+ // types
+ 	

File trac_9880_ginac_infinities_rewrite.patch

+# HG changeset patch
+# Parent aed6fef803be8eb290dcbb6e0e34e51abdb1b863
 # HG changeset patch
 # Parent aed6fef803be8eb290dcbb6e0e34e51abdb1b863
 
+Trac #9880: The pynac C++ part of the infinities rewrite
 diff --git a/.hgignore b/.hgignore
 --- a/.hgignore
 +++ b/.hgignore
  	}
 -	
 +		
-+        // handle infinity
-+        for (epvector::const_iterator i = seq.begin(); i != seq.end(); i++)
++	// handle infinity
++	for (epvector::const_iterator i = seq.begin(); i != seq.end(); i++)
 +		if (unlikely(is_exactly_a<infinity>(i->rest)))
 +			return eval_infinity(i);
 +
 new file mode 100644
 --- /dev/null
 +++ b/ginac/infinity.cpp
-@@ -0,0 +1,340 @@
+@@ -0,0 +1,342 @@
 +/** @file infinity.cpp
 + *
 + *  Implementation of PyNaC's "infinity". */
 +{
 +	if (inf==info_flags::infinity)
 +		return true;
-+		
-+	// if (inf == info_flags::real)
-+	// 	return domain==domain::real || domain==domain::positive ;
-+	// if (inf==info_flags::positive || inf==info_flags::nonnegative)
-+	// 	return domain == domain::positive;
-+
++	if (inf == info_flags::real || 
++	    inf == info_flags::positive ||
++	    inf == info_flags::negative)
++		return direction.info(inf);
++	if (inf == info_flags::nonnegative)
++		return direction.info(info_flags::positive);
++	if (inf == info_flags::nonnegative)
++		return direction.info(info_flags::positive);
 +	return inherited::info(inf);
 +}
 +
 +	else if (rhs.is_zero())
 +		throw(std::runtime_error("indeterminate expression: "
 +					 "0 * infinity encountered."));
-+	else if (rhs.info(info_flags::positive))
++	else if (rhs.info(info_flags::positive)) {
 +		return *this;
-+	else if (rhs.info(info_flags::negative)) {
++	} else if (rhs.info(info_flags::negative)) {
 +		direction = mul(-1, direction);
 +		return *this;
-+	} else if (is_a<numeric>(rhs)) {
++	} else if (rhs.nsymbols()==0) {
 +		set_direction(mul(direction, rhs));
 +		return *this;
 +	}
 +	throw(std::runtime_error("indeterminate expression: "
-+				 "infinity * (expression of unknown sign) encountered."));
++				 "infinity * f(x) encountered."));
 +}
 +
 +
  		/* end paranoia */
  		++i;
  	}
-@@ -628,52 +631,29 @@
+@@ -627,53 +630,30 @@
+ 		GINAC_ASSERT(seq.size()>1 || !overall_coeff.is_equal(_ex1));
  		return *this;
  	}
- 	
+-	
++
 +	// handle infinity and handle exp(a)*exp(b) -> exp(a+b) and
 +	unsigned exp_count = 0;
 +	for (epvector::const_iterator i = seq.begin(); i != seq.end(); i++) {
  			if (s->empty()) {
  				return ex_to<numeric>(overall_coeff).mul_dyn(oc);
  			}
-@@ -904,9 +747,72 @@
+@@ -904,9 +747,53 @@
  			       )->setflag(status_flags::dynallocated);
  		}
  	}
 +
 +ex mul::eval_exponentials() const
 +{
-+	// if (seq_size==1) {
-+	// 	// even though exp defines a power simplification rule,
-+	// 	// it is possible to end up with i->coeff != 1, e.g.,
-+	// 	// as a result of construct_from_2_ex
-+	// 	const numeric & coeff = ex_to<numeric>(seg.begin()->coeff);
-+	// 	const ex & rest = seg.begin()->rest;
-+	// 	GINAC_ASSERT(is_ex_the_function(rest, exp));
-+	// 	GINAC_ASSERT(coeff.is_integer());
-+	// 	GINAC_ASSERT(not coeff.is_equal(_ex1));
-+	// 	ex new_exp = pow(rest, coeff);
-+	// 	if (is_a<numeric>(new_exp))
-+	// 		return ex_to<numeric>(overall_coeff).mul(
-+	// 				ex_to<numeric>(new_exp));
-+	// 	return (new mul(new_exp, ex_to<numeric>(overall_coeff)))
-+	// 		->setflag(status_flags::dynallocated);	
-+	// }
-+
 +	ex exp_arg = _ex0;
 +	numeric oc = *_num1_p;
 +	std::auto_ptr<epvector> s(new epvector);
 +{
 +	GINAC_ASSERT(is_exactly_a<infinity>(infinity_iter->rest));
 +	GINAC_ASSERT(infinity_iter->coeff.is_equal(_ex1));
-+	infinity result = ex_to<infinity>(infinity_iter->rest);
++	infinity result = ex_to<infinity>(recombine_pair_to_ex(*infinity_iter));
 +	result *= overall_coeff;
 +
 +        for (epvector::const_iterator i = seq.begin(); i != seq.end(); i++) {
-+                if (not is_exactly_a<infinity>(i->rest)) continue;
 +                if (i == infinity_iter) continue;
-+		GINAC_ASSERT(i->coeff.is_equal(_ex1));
-+		result *= ex_to<infinity>(i->rest);
++		result *= recombine_pair_to_ex(*i);
 +        }
 +	return result;
 +}
  #include "operators.h"
  #include "inifcns.h" // for log() in power::derivative() and exp for printing
  #include "matrix.h"
+@@ -466,53 +467,38 @@
+ 	}
+ 
+ 	// ^(\infty, x)
+-	// error if x is not numeric and real
+-	// -> 0 if x < 0
+-	// -> error if x == 0
+-	// -> Infinity if \infty is NegInfinity and x is even
+-	// -> \infty otherwise
+-	if (ebasis.info(info_flags::infinity)) {
+-		if (exponent_is_numerical) {
+-			if (!num_exponent->is_real()) {
+-				throw(std::domain_error("power::eval(): pow(Infinity, x) is not defined for complex x."));
+-			} else if (num_exponent->csgn() == -1)
+-				return _ex0;
+-			else if (num_exponent->is_zero())
+-				throw(std::domain_error("power::eval(): pow(Infinity, 0) is undefined."));
+-			else if (ebasis.is_equal(NegInfinity) && 
+-					num_exponent->is_even())
++	if (is_a<infinity>(ebasis)) {
++		const infinity & basis_inf = ex_to<infinity>(ebasis);
++		if (eexponent.nsymbols()>0)
++			throw(std::domain_error("power::eval(): pow(Infinity, f(x)) is not defined."));
++		if (eexponent.is_zero())
++			throw(std::domain_error("power::eval(): pow(Infinity, 0) is undefined."));
++		if (eexponent.info(info_flags::negative))
++			return _ex0;
++		if (eexponent.info(info_flags::positive))
++			if (basis_inf.is_plus_infinity())
+ 				return Infinity;
+ 			else
+-				return ebasis;
+-		} else
+-			throw(std::domain_error("power::eval(): pow(Infinity, x) for non numeric x is not defined."));
++				return UnsignedInfinity;
++		throw(std::domain_error("power::eval(): pow(Infinity, c)"
++					" for constant of undetermined sign is not defined."));
+ 	}
++
+ 	// ^(x, \infty)
+-	// error if x is not numeric
+-	// error if \infty is UnsignedInfinity
+-	// error if x \in {0, 1, -1}
+-	// 0 if \infty is NegInfinity
+-	// UnsignedInfinity if x < 0
+-	// Infinity otherwise
+-	if (eexponent.info(info_flags::infinity)) {
+-		if (!basis_is_numerical) {
+-			throw(std::domain_error("power::eval(): pow(x, Infinity) for non numeric x is not defined."));
+-		} else if (!num_basis->is_real()) {
+-			throw(std::domain_error("power::eval(): pow(x, Infinity) for non real x is not defined."));
+-		} else if (ebasis.is_equal(UnsignedInfinity)) {
+-			throw(std::domain_error("power::eval(): pow(x, UnsignedInfinity) is not defined."));
+-		} else if (num_basis->is_zero()) {
+-			throw(std::domain_error("power::eval(): pow(0, Infinity) is not defined."));
+-		} else if (num_basis->is_equal(*_num1_p) || 
+-				num_basis->is_equal(*_num_1_p)) {
++	if (is_a<infinity>(eexponent)) {
++		const infinity & exp_inf = ex_to<infinity>(eexponent);
++		if (exp_inf.is_unsigned_infinity())
++			throw(std::domain_error("power::eval(): pow(x, unsigned_infinity) is not defined."));
++		if (ebasis.nsymbols()>0) 
++			throw(std::domain_error("power::eval(): pow(f(x), infinity) is not defined."));
++		// x^(c*oo) --> (x^c)^(+oo)
++		const ex abs_base = abs(pow(ebasis, exp_inf.get_direction()));
++		if (abs_base > _ex1) return Infinity;
++		if (abs_base < _ex1) return _ex0;
++		if (abs_base == _ex1)
+ 			throw(std::domain_error("power::eval(): pow(1, Infinity) is not defined."));
+-		} else if (eexponent.is_equal(NegInfinity)) {
+-			return _ex0;
+-		} else if (num_basis->csgn() == -1) {
+-			return UnsignedInfinity;
+-		} else {
+-			return Infinity;
+-		}
++		throw(std::domain_error("power::eval(): pow(c, Infinity)"
++					" for unknown magnitude |c| is not defined."));
+ 	}
+ 	
+ 	// ^(x,0) -> 1  (0^0 also handled here)
+diff --git a/ginac/relational.cpp b/ginac/relational.cpp
+--- a/ginac/relational.cpp
++++ b/ginac/relational.cpp
+@@ -20,11 +20,13 @@
+  *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
+  */
+ 
++#include "compiler.h"
+ #include "relational.h"
+ #include "operators.h"
+ #include "numeric.h"
+ #include "archive.h"
+ #include "utils.h"
++#include "infinity.h"
+ 
+ #include <iostream>
+ #include <stdexcept>
+@@ -381,6 +383,27 @@
+  *  unequal or undecidable). */