From dcc848e23e4d793b2b1095012a8951f5245a68d5 Mon Sep 17 00:00:00 2001
From: Christian Bauer <Christian.Bauer@uni-mainz.de>
Date: Tue, 4 Mar 2003 23:59:06 +0000
Subject: [PATCH] return type and parameter type of nops() and op()/let_op()
 (respectively) is now a size_t

---
 doc/tutorial/ginac.texi | 20 +++++++-------
 ginac/add.cpp           |  2 +-
 ginac/basic.cpp         | 57 +++++++++++++++++++++------------------
 ginac/basic.h           | 13 ++++-----
 ginac/clifford.cpp      | 28 +++++++++----------
 ginac/color.cpp         |  6 ++---
 ginac/container.pl      | 18 ++++++-------
 ginac/ex.cpp            |  6 ++---
 ginac/ex.h              | 14 +++++-----
 ginac/expairseq.cpp     | 28 +++++++++----------
 ginac/expairseq.h       |  4 +--
 ginac/exprseq_suppl.cpp |  4 +--
 ginac/fderivative.cpp   |  4 +--
 ginac/function.pl       | 16 +++++------
 ginac/idx.cpp           |  6 ++---
 ginac/idx.h             |  8 +++---
 ginac/indexed.cpp       | 48 ++++++++++++++++-----------------
 ginac/inifcns.cpp       | 14 +++++-----
 ginac/matrix.cpp        | 22 +++++++--------
 ginac/matrix.h          |  8 +++---
 ginac/mul.cpp           | 45 +++++++++----------------------
 ginac/ncmul.cpp         | 44 +++++++++++++++---------------
 ginac/ncmul.h           |  2 +-
 ginac/normal.cpp        | 60 ++++++++++++++++++++---------------------
 ginac/power.cpp         | 17 ++++++------
 ginac/power.h           |  4 +--
 ginac/pseries.cpp       | 14 +++++-----
 ginac/pseries.h         |  4 +--
 ginac/relational.cpp    |  5 ++--
 ginac/relational.h      |  4 +--
 ginac/remember.cpp      |  4 +--
 ginac/symbol.cpp        |  2 +-
 ginac/symmetry.cpp      | 10 +++----
 ginac/tensor.cpp        | 14 +++++-----
 ginsh/ginsh_parser.yy   |  4 +--
 35 files changed, 269 insertions(+), 290 deletions(-)

diff --git a/doc/tutorial/ginac.texi b/doc/tutorial/ginac.texi
index 7ffc1d0a..f0369865 100644
--- a/doc/tutorial/ginac.texi
+++ b/doc/tutorial/ginac.texi
@@ -2984,8 +2984,8 @@ for an explanation of these.
 GiNaC provides the two methods
 
 @example
-unsigned ex::nops();
-ex ex::op(unsigned i);
+size_t ex::nops();
+ex ex::op(size_t i);
 @end example
 
 for accessing the subexpressions in the container-like GiNaC classes like
@@ -3450,7 +3450,7 @@ ex calc_trace(ex e)
         return ex_to<matrix>(e).trace();
     else if (is_a<add>(e)) @{
         ex sum = 0;
-        for (unsigned i=0; i<e.nops(); i++)
+        for (size_t i=0; i<e.nops(); i++)
             sum += calc_trace(e.op(i));
         return sum;
     @} else if (is_a<mul>)(e)) @{
@@ -4491,9 +4491,9 @@ static void my_print(const ex & e)
     else
         cout << e.bp->class_name();
     cout << "(";
-    unsigned n = e.nops();
+    size_t n = e.nops();
     if (n)
-        for (unsigned i=0; i<n; i++) @{
+        for (size_t i=0; i<n; i++) @{
             my_print(e.op(i));
             if (i != n-1)
                 cout << ",";
@@ -4678,8 +4678,8 @@ static void my_print2(const archive_node & n)
     archive_node::propinfovector p;
     n.get_properties(p);
 
-    unsigned num = p.size();
-    for (unsigned i=0; i<num; i++) @{
+    size_t num = p.size();
+    for (size_t i=0; i<num; i++) @{
         const string &name = p[i].name;
         if (name == "class")
             continue;
@@ -5548,9 +5548,9 @@ If your class stores sub-expressions you will probably want to override
 
 @cindex @code{let_op()}
 @example
-unsigned nops() cont;
-ex op(int i) const;
-ex & let_op(int i);
+size_t nops() cont;
+ex op(size_t i) const;
+ex & let_op(size_t i);
 ex map(map_function & f) const;
 ex subs(const lst & ls, const lst & lr, unsigned options = 0) const;
 @end example
diff --git a/ginac/add.cpp b/ginac/add.cpp
index b409513c..6a6bbf6e 100644
--- a/ginac/add.cpp
+++ b/ginac/add.cpp
@@ -164,7 +164,7 @@ void add::print(const print_context & c, unsigned level) const
 
 		c.s << class_name() << '(';
 		op(0).print(c);
-		for (unsigned i=1; i<nops(); ++i) {
+		for (size_t i=1; i<nops(); ++i) {
 			c.s << ',';
 			op(i).print(c);
 		}
diff --git a/ginac/basic.cpp b/ginac/basic.cpp
index 70237930..d2b624a6 100644
--- a/ginac/basic.cpp
+++ b/ginac/basic.cpp
@@ -116,7 +116,7 @@ void basic::print(const print_context & c, unsigned level) const
 		    << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 		    << ", nops=" << nops()
 		    << std::endl;
-		for (unsigned i=0; i<nops(); ++i)
+		for (size_t i=0; i<nops(); ++i)
 			op(i).print(c, level + static_cast<const print_tree &>(c).delta_indent);
 
 	} else
@@ -168,7 +168,7 @@ bool basic::info(unsigned inf) const
 }
 
 /** Number of operands/members. */
-unsigned basic::nops() const
+size_t basic::nops() const
 {
 	// iterating from 0 to nops() on atomic objects should be an empty loop,
 	// and accessing their elements is a range error.  Container objects should
@@ -177,27 +177,27 @@ unsigned basic::nops() const
 }
 
 /** Return operand/member at position i. */
-ex basic::op(int i) const
+ex basic::op(size_t i) const
 {
-	throw(std::runtime_error(std::string("op() not implemented by ") + class_name()));
+	throw(std::range_error(std::string("basic::op(): ") + class_name() + std::string(" has no operands")));
 }
 
 /** Return modifyable operand/member at position i. */
-ex & basic::let_op(int i)
+ex & basic::let_op(size_t i)
 {
 	ensure_if_modifiable();
-	throw(std::runtime_error(std::string("let_op() not implemented by ") + class_name()));
+	throw(std::range_error(std::string("basic::let_op(): ") + class_name() + std::string(" has no operands")));
 }
 
 ex basic::operator[](const ex & index) const
 {
 	if (is_exactly_a<numeric>(index))
-		return op(ex_to<numeric>(index).to_int());
+		return op(static_cast<size_t>(ex_to<numeric>(index).to_int()));
 
-	throw(std::invalid_argument("non-numeric indices not supported by this type"));
+	throw(std::invalid_argument(std::string("non-numeric indices not supported by ") + class_name()));
 }
 
-ex basic::operator[](int i) const
+ex basic::operator[](size_t i) const
 {
 	return op(i);
 }
@@ -207,10 +207,10 @@ ex & basic::operator[](const ex & index)
 	if (is_exactly_a<numeric>(index))
 		return let_op(ex_to<numeric>(index).to_int());
 
-	throw(std::invalid_argument("non-numeric indices not supported by this type"));
+	throw(std::invalid_argument(std::string("non-numeric indices not supported by ") + class_name()));
 }
 
-ex & basic::operator[](int i)
+ex & basic::operator[](size_t i)
 {
 	return let_op(i);
 }
@@ -224,7 +224,7 @@ bool basic::has(const ex & pattern) const
 	lst repl_lst;
 	if (match(pattern, repl_lst))
 		return true;
-	for (unsigned i=0; i<nops(); i++)
+	for (size_t i=0; i<nops(); i++)
 		if (op(i).has(pattern))
 			return true;
 	
@@ -235,7 +235,7 @@ bool basic::has(const ex & pattern) const
  *  sub-expressions (one level only, not recursively). */
 ex basic::map(map_function & f) const
 {
-	unsigned num = nops();
+	size_t num = nops();
 	if (num == 0)
 		return *this;
 
@@ -243,7 +243,7 @@ ex basic::map(map_function & f) const
 	copy->setflag(status_flags::dynallocated);
 	copy->clearflag(status_flags::hash_calculated | status_flags::expanded);
 	ex e(*copy);
-	for (unsigned i=0; i<num; i++)
+	for (size_t i=0; i<num; i++)
 		e.let_op(i) = f(e.op(i));
 	return e.eval();
 }
@@ -286,7 +286,7 @@ ex basic::collect(const ex & s, bool distributed) const
 		else if (distributed) {
 
 			// Get lower/upper degree of all symbols in list
-			int num = s.nops();
+			size_t num = s.nops();
 			struct sym_info {
 				ex sym;
 				int ldeg, deg;
@@ -295,7 +295,7 @@ ex basic::collect(const ex & s, bool distributed) const
 			};
 			sym_info *si = new sym_info[num];
 			ex c = *this;
-			for (int i=0; i<num; i++) {
+			for (size_t i=0; i<num; i++) {
 				si[i].sym = s.op(i);
 				si[i].ldeg = si[i].cnt = this->ldegree(si[i].sym);
 				si[i].deg = this->degree(si[i].sym);
@@ -306,14 +306,14 @@ ex basic::collect(const ex & s, bool distributed) const
 
 				// Calculate coeff*x1^c1*...*xn^cn
 				ex y = _ex1;
-				for (int i=0; i<num; i++) {
+				for (size_t i=0; i<num; i++) {
 					int cnt = si[i].cnt;
 					y *= power(si[i].sym, cnt);
 				}
 				x += y * si[num - 1].coeff;
 
 				// Increment counters
-				int n = num - 1;
+				size_t n = num - 1;
 				while (true) {
 					++si[n].cnt;
 					if (si[n].cnt <= si[n].deg) {
@@ -323,7 +323,7 @@ ex basic::collect(const ex & s, bool distributed) const
 							c = *this;
 						else
 							c = si[n - 1].coeff;
-						for (int i=n; i<num; i++)
+						for (size_t i=n; i<num; i++)
 							c = si[i].coeff = c.coeff(si[i].sym, si[i].cnt);
 						break;
 					}
@@ -340,8 +340,13 @@ done:		delete[] si;
 
 			// Recursive form
 			x = *this;
-			for (int n=s.nops()-1; n>=0; n--)
+			size_t n = s.nops() - 1;
+			while (true) {
 				x = x.collect(s[n]);
+				if (n == 0)
+					break;
+				n--;
+			}
 		}
 
 	} else {
@@ -479,7 +484,7 @@ bool basic::match(const ex & pattern, lst & repl_lst) const
 		// Wildcard matches anything, but check whether we already have found
 		// a match for that wildcard first (if so, it the earlier match must
 		// be the same expression)
-		for (unsigned i=0; i<repl_lst.nops(); i++) {
+		for (size_t i=0; i<repl_lst.nops(); i++) {
 			if (repl_lst.op(i).op(0).is_equal(pattern))
 				return is_equal(ex_to<basic>(repl_lst.op(i).op(1)));
 		}
@@ -506,7 +511,7 @@ bool basic::match(const ex & pattern, lst & repl_lst) const
 			return false;
 
 		// Otherwise the subexpressions must match one-to-one
-		for (unsigned i=0; i<nops(); i++)
+		for (size_t i=0; i<nops(); i++)
 			if (!op(i).match(pattern.op(i), repl_lst))
 				return false;
 
@@ -522,12 +527,12 @@ ex basic::subs(const lst & ls, const lst & lr, unsigned options) const
 	GINAC_ASSERT(ls.nops() == lr.nops());
 
 	if (options & subs_options::subs_no_pattern) {
-		for (unsigned i=0; i<ls.nops(); i++) {
+		for (size_t i=0; i<ls.nops(); i++) {
 			if (is_equal(ex_to<basic>(ls.op(i))))
 				return lr.op(i);
 		}
 	} else {
-		for (unsigned i=0; i<ls.nops(); i++) {
+		for (size_t i=0; i<ls.nops(); i++) {
 			lst repl_lst;
 			if (match(ex_to<basic>(ls.op(i)), repl_lst))
 				return lr.op(i).subs(repl_lst, options | subs_options::subs_no_pattern); // avoid infinite recursion when re-substituting the wildcards
@@ -653,7 +658,7 @@ unsigned basic::return_type_tinfo(void) const
 unsigned basic::calchash(void) const
 {
 	unsigned v = golden_ratio_hash(tinfo());
-	for (unsigned i=0; i<nops(); i++) {
+	for (size_t i=0; i<nops(); i++) {
 		v = rotate_left(v);
 		v ^= (const_cast<basic *>(this))->op(i).gethash();
 	}
@@ -708,7 +713,7 @@ ex basic::subs(const ex & e, unsigned options) const
 	}
 	lst ls;
 	lst lr;
-	for (unsigned i=0; i<e.nops(); i++) {
+	for (size_t i=0; i<e.nops(); i++) {
 		ex r = e.op(i);
 		if (!r.info(info_flags::relation_equal)) {
 			throw(std::invalid_argument("basic::subs(ex): argument must be a list of equations"));
diff --git a/ginac/basic.h b/ginac/basic.h
index 57cf5734..c68234f3 100644
--- a/ginac/basic.h
+++ b/ginac/basic.h
@@ -23,6 +23,7 @@
 #ifndef __GINAC_BASIC_H__
 #define __GINAC_BASIC_H__
 
+#include <cstddef> // for size_t
 #include <vector>
 // CINT needs <algorithm> to work properly with <vector>
 #include <algorithm>
@@ -96,13 +97,13 @@ public: // only const functions please (may break reference counting)
 	virtual void dbgprinttree(void) const;
 	virtual unsigned precedence(void) const;
 	virtual bool info(unsigned inf) const;
-	virtual unsigned nops() const;
-	virtual ex op(int i) const;
+	virtual size_t nops() const;
+	virtual ex op(size_t i) const;
 	virtual ex operator[](const ex & index) const;
-	virtual ex operator[](int i) const;
-	virtual ex & let_op(int i);
+	virtual ex operator[](size_t i) const;
+	virtual ex & let_op(size_t i);
 	virtual ex & operator[](const ex & index);
-	virtual ex & operator[](int i);
+	virtual ex & operator[](size_t i);
 	virtual ex expand(unsigned options = 0) const;
 	virtual bool has(const ex & other) const;
 	virtual ex map(map_function & f) const;
@@ -140,7 +141,7 @@ protected: // functions that should be called from class ex only
 	// non-virtual functions in this class
 public:
 	ex subs(const ex & e, unsigned options = 0) const;
-	ex diff(const symbol & s, unsigned nth=1) const;
+	ex diff(const symbol & s, unsigned nth = 1) const;
 	int compare(const basic & other) const;
 	bool is_equal(const basic & other) const;
 	const basic & hold(void) const;
diff --git a/ginac/clifford.cpp b/ginac/clifford.cpp
index 0918fe42..9a1330e4 100644
--- a/ginac/clifford.cpp
+++ b/ginac/clifford.cpp
@@ -525,7 +525,7 @@ static bool is_clifford_tinfo(unsigned ti)
 
 /** Take trace of a string of an even number of Dirac gammas given a vector
  *  of indices. */
-static ex trace_string(exvector::const_iterator ix, unsigned num)
+static ex trace_string(exvector::const_iterator ix, size_t num)
 {
 	// Tr gamma.mu gamma.nu = 4 g.mu.nu
 	if (num == 2)
@@ -547,8 +547,8 @@ static ex trace_string(exvector::const_iterator ix, unsigned num)
 	exvector v(num - 2);
 	int sign = 1;
 	ex result;
-	for (unsigned i=1; i<num; i++) {
-		for (unsigned n=1, j=0; n<num; n++) {
+	for (size_t i=1; i<num; i++) {
+		for (size_t n=1, j=0; n<num; n++) {
 			if (n == i)
 				continue;
 			v[j++] = ix[n];
@@ -577,7 +577,7 @@ ex dirac_trace(const ex & e, unsigned char rl, const ex & trONE)
 
 		// Trace of product: pull out non-clifford factors
 		ex prod = _ex1;
-		for (unsigned i=0; i<e.nops(); i++) {
+		for (size_t i=0; i<e.nops(); i++) {
 			const ex &o = e.op(i);
 			if (is_clifford_tinfo(o.return_type_tinfo(), rl))
 				prod *= dirac_trace(o, rl, trONE);
@@ -601,7 +601,7 @@ ex dirac_trace(const ex & e, unsigned char rl, const ex & trONE)
 
 		// gamma5 gets moved to the front so this check is enough
 		bool has_gamma5 = is_a<diracgamma5>(e.op(0).op(0));
-		unsigned num = e.nops();
+		size_t num = e.nops();
 
 		if (has_gamma5) {
 
@@ -625,23 +625,23 @@ ex dirac_trace(const ex & e, unsigned char rl, const ex & trONE)
 			//   I/4! * epsilon0123.mu1.mu2.mu3.mu4 * Tr gamma.mu1 gamma.mu2 gamma.mu3 gamma.mu4 S_2k
 			// (the epsilon is always 4-dimensional)
 			exvector ix(num-1), bv(num-1);
-			for (unsigned i=1; i<num; i++)
+			for (size_t i=1; i<num; i++)
 				base_and_index(e.op(i), bv[i-1], ix[i-1]);
 			num--;
 			int *iv = new int[num];
 			ex result;
-			for (unsigned i=0; i<num-3; i++) {
+			for (size_t i=0; i<num-3; i++) {
 				ex idx1 = ix[i];
-				for (unsigned j=i+1; j<num-2; j++) {
+				for (size_t j=i+1; j<num-2; j++) {
 					ex idx2 = ix[j];
-					for (unsigned k=j+1; k<num-1; k++) {
+					for (size_t k=j+1; k<num-1; k++) {
 						ex idx3 = ix[k];
-						for (unsigned l=k+1; l<num; l++) {
+						for (size_t l=k+1; l<num; l++) {
 							ex idx4 = ix[l];
 							iv[0] = i; iv[1] = j; iv[2] = k; iv[3] = l;
 							exvector v;
 							v.reserve(num - 4);
-							for (unsigned n=0, t=4; n<num; n++) {
+							for (size_t n=0, t=4; n<num; n++) {
 								if (n == i || n == j || n == k || n == l)
 									continue;
 								iv[t++] = n;
@@ -672,7 +672,7 @@ ex dirac_trace(const ex & e, unsigned char rl, const ex & trONE)
 			}
 
 			exvector iv(num), bv(num);
-			for (unsigned i=0; i<num; i++)
+			for (size_t i=0; i<num; i++)
 				base_and_index(e.op(i), bv[i], iv[i]);
 
 			return trONE * (trace_string(iv.begin(), num) * mul(bv)).simplify_indexed();
@@ -693,7 +693,7 @@ ex canonicalize_clifford(const ex & e)
 	// Scan for any ncmul objects
 	lst srl;
 	ex aux = e.to_rational(srl);
-	for (unsigned i=0; i<srl.nops(); i++) {
+	for (size_t i=0; i<srl.nops(); i++) {
 
 		ex lhs = srl.op(i).lhs();
 		ex rhs = srl.op(i).rhs();
@@ -713,7 +713,7 @@ ex canonicalize_clifford(const ex & e)
 
 			exvector v;
 			v.reserve(rhs.nops());
-			for (unsigned j=0; j<rhs.nops(); j++)
+			for (size_t j=0; j<rhs.nops(); j++)
 				v.push_back(rhs.op(j));
 
 			// Stupid recursive bubble sort because we only want to swap adjacent gammas
diff --git a/ginac/color.cpp b/ginac/color.cpp
index 9a11b4f7..94f1871b 100644
--- a/ginac/color.cpp
+++ b/ginac/color.cpp
@@ -530,7 +530,7 @@ ex color_trace(const ex & e, unsigned char rl)
 
 		// Trace of product: pull out non-color factors
 		ex prod = _ex1;
-		for (unsigned i=0; i<e.nops(); i++) {
+		for (size_t i=0; i<e.nops(); i++) {
 			const ex &o = e.op(i);
 			if (is_color_tinfo(o.return_type_tinfo(), rl))
 				prod *= color_trace(o, rl);
@@ -549,7 +549,7 @@ ex color_trace(const ex & e, unsigned char rl)
 		if (!is_a<ncmul>(e_expanded))
 			return color_trace(e_expanded, rl);
 
-		unsigned num = e.nops();
+		size_t num = e.nops();
 
 		if (num == 2) {
 
@@ -573,7 +573,7 @@ ex color_trace(const ex & e, unsigned char rl)
 
 			exvector v1;
 			v1.reserve(num - 2);
-			for (unsigned i=0; i<num-2; i++)
+			for (size_t i=0; i<num-2; i++)
 				v1.push_back(e.op(i));
 
 			exvector v2 = v1;
diff --git a/ginac/container.pl b/ginac/container.pl
index 84c0a49d..945622fe 100755
--- a/ginac/container.pl
+++ b/ginac/container.pl
@@ -107,30 +107,28 @@ END_OF_SORT_IMPLEMENTATION
 
 if ($let_op) {
 	$LET_OP_DEFINITION=<<END_OF_LET_OP_DEFINITION;
-	ex & let_op(int i);
+	ex & let_op(size_t i);
 END_OF_LET_OP_DEFINITION
 
 	$LET_OP_IMPLEMENTATION=<<END_OF_LET_OP_IMPLEMENTATION;
-ex ${CONTAINER}::op(int i) const
+ex ${CONTAINER}::op(size_t i) const
 {
-	GINAC_ASSERT(i>=0);
 	GINAC_ASSERT(i<nops());
 
 	${STLT}::const_iterator it=seq.begin();
-	for (int j=0; j<i; j++) {
+	for (size_t j=0; j<i; j++) {
 		++it;
 	}
 	return *it;
 }
 
-ex & ${CONTAINER}::let_op(int i)
+ex & ${CONTAINER}::let_op(size_t i)
 {
-	GINAC_ASSERT(i>=0);
 	GINAC_ASSERT(i<nops());
 
 	ensure_if_modifiable();
 	${STLT}::iterator it=seq.begin();
-	for (int j=0; j<i; j++) {
+	for (size_t j=0; j<i; j++) {
 		++it;
 	}
 	return *it;
@@ -255,8 +253,8 @@ public:
 	void print(const print_context & c, unsigned level = 0) const;
 	unsigned precedence(void) const {return 10;}
 	bool info(unsigned inf) const;
-	unsigned nops() const;
-	ex op(int i) const;
+	size_t nops() const;
+	ex op(size_t i) const;
 ${LET_OP_DEFINITION}
 	ex map(map_function & f) const;
 	ex eval(int level=0) const;
@@ -462,7 +460,7 @@ void ${CONTAINER}::print(const print_context & c, unsigned level) const
 
 // ${CONTAINER}::info() will be implemented by user elsewhere";
 
-unsigned ${CONTAINER}::nops() const
+size_t ${CONTAINER}::nops() const
 {
 	return seq.size();
 }
diff --git a/ginac/ex.cpp b/ginac/ex.cpp
index d4c37a0a..2c31212c 100644
--- a/ginac/ex.cpp
+++ b/ginac/ex.cpp
@@ -124,14 +124,14 @@ bool ex::find(const ex & pattern, lst & found) const
 		return true;
 	}
 	bool any_found = false;
-	for (unsigned i=0; i<nops(); i++)
+	for (size_t i=0; i<nops(); i++)
 		if (op(i).find(pattern, found))
 			any_found = true;
 	return any_found;
 }
 
 /** Return modifyable operand/member at position i. */
-ex & ex::let_op(int i)
+ex & ex::let_op(size_t i)
 {
 	makewriteable();
 	GINAC_ASSERT(bp!=0);
@@ -145,7 +145,7 @@ ex & ex::operator[](const ex & index)
 	return (*bp)[index];
 }
 
-ex & ex::operator[](int i)
+ex & ex::operator[](size_t i)
 {
 	makewriteable();
 	GINAC_ASSERT(bp!=0);
diff --git a/ginac/ex.h b/ginac/ex.h
index 3d79c7ac..5e8ac0de 100644
--- a/ginac/ex.h
+++ b/ginac/ex.h
@@ -115,7 +115,7 @@ public:
 	void dbgprint(void) const;
 	void dbgprinttree(void) const;
 	bool info(unsigned inf) const { return bp->info(inf); }
-	unsigned nops() const { return bp->nops(); }
+	size_t nops() const { return bp->nops(); }
 	ex expand(unsigned options=0) const;
 	bool has(const ex & pattern) const { return bp->has(pattern); }
 	ex map(map_function & f) const { return bp->map(f); }
@@ -159,12 +159,12 @@ public:
 	ex symmetrize_cyclic(void) const;
 	ex symmetrize_cyclic(const lst & l) const;
 	ex eval_ncmul(const exvector & v) const { return bp->eval_ncmul(v); }
-	ex op(int i) const { return bp->op(i); }
+	ex op(size_t i) const { return bp->op(i); }
 	ex operator[](const ex & index) const { return (*bp)[index]; }
-	ex operator[](int i) const { return (*bp)[i]; }
-	ex & let_op(int i);
+	ex operator[](size_t i) const { return (*bp)[i]; }
+	ex & let_op(size_t i);
 	ex & operator[](const ex & index);
-	ex & operator[](int i);
+	ex & operator[](size_t i);
 	ex lhs(void) const;
 	ex rhs(void) const;
 	int compare(const ex & other) const;
@@ -371,7 +371,7 @@ inline bool are_ex_trivially_equal(const ex &e1, const ex &e2)
 }
 
 // wrapper functions around member functions
-inline unsigned nops(const ex & thisex)
+inline size_t nops(const ex & thisex)
 { return thisex.nops(); }
 
 inline ex expand(const ex & thisex, unsigned options = 0)
@@ -461,7 +461,7 @@ inline ex symmetrize_cyclic(const ex & thisex)
 inline ex symmetrize_cyclic(const ex & thisex, const lst & l)
 { return thisex.symmetrize_cyclic(l); }
 
-inline ex op(const ex & thisex, int i)
+inline ex op(const ex & thisex, size_t i)
 { return thisex.op(i); }
 
 inline ex lhs(const ex & thisex)
diff --git a/ginac/expairseq.cpp b/ginac/expairseq.cpp
index 5c3c2fe7..0fbbf543 100644
--- a/ginac/expairseq.cpp
+++ b/ginac/expairseq.cpp
@@ -197,8 +197,8 @@ void expairseq::print(const print_context &c, unsigned level) const
 		    << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 		    << ", nops=" << nops()
 		    << std::endl;
-		unsigned num = seq.size();
-		for (unsigned i=0; i<num; ++i) {
+		size_t num = seq.size();
+		for (size_t i=0; i<num; ++i) {
 			seq[i].rest.print(c, level + delta_indent);
 			seq[i].coeff.print(c, level + delta_indent);
 			if (i != num - 1)
@@ -276,7 +276,7 @@ bool expairseq::info(unsigned inf) const
 	return inherited::info(inf);
 }
 
-unsigned expairseq::nops() const
+size_t expairseq::nops() const
 {
 	if (overall_coeff.is_equal(default_overall_coeff()))
 		return seq.size();
@@ -284,9 +284,9 @@ unsigned expairseq::nops() const
 		return seq.size()+1;
 }
 
-ex expairseq::op(int i) const
+ex expairseq::op(size_t i) const
 {
-	if (unsigned(i)<seq.size())
+	if (i < seq.size())
 		return recombine_pair_to_ex(seq[i]);
 	GINAC_ASSERT(!overall_coeff.is_equal(default_overall_coeff()));
 	return overall_coeff;
@@ -333,7 +333,7 @@ bool expairseq::match(const ex & pattern, lst & repl_lst) const
 		// expression", like "*" above) is present
 		bool has_global_wildcard = false;
 		ex global_wildcard;
-		for (unsigned int i=0; i<pattern.nops(); i++) {
+		for (size_t i=0; i<pattern.nops(); i++) {
 			if (is_exactly_a<wildcard>(pattern.op(i))) {
 				has_global_wildcard = true;
 				global_wildcard = pattern.op(i);
@@ -347,12 +347,12 @@ bool expairseq::match(const ex & pattern, lst & repl_lst) const
 		// Chop into terms
 		exvector ops;
 		ops.reserve(nops());
-		for (unsigned i=0; i<nops(); i++)
+		for (size_t i=0; i<nops(); i++)
 			ops.push_back(op(i));
 
 		// Now, for every term of the pattern, look for a matching term in
 		// the expression and remove the match
-		for (unsigned i=0; i<pattern.nops(); i++) {
+		for (size_t i=0; i<pattern.nops(); i++) {
 			ex p = pattern.op(i);
 			if (has_global_wildcard && p.is_equal(global_wildcard))
 				continue;
@@ -373,13 +373,13 @@ found:		;
 			// Assign all the remaining terms to the global wildcard (unless
 			// it has already been matched before, in which case the matches
 			// must be equal)
-			unsigned num = ops.size();
+			size_t num = ops.size();
 			epvector *vp = new epvector();
 			vp->reserve(num);
-			for (unsigned i=0; i<num; i++)
+			for (size_t i=0; i<num; i++)
 				vp->push_back(split_ex_to_pair(ops[i]));
 			ex rest = thisexpairseq(vp, default_overall_coeff());
-			for (unsigned i=0; i<repl_lst.nops(); i++) {
+			for (size_t i=0; i<repl_lst.nops(); i++) {
 				if (repl_lst.op(i).op(0).is_equal(global_wildcard))
 					return rest.is_equal(repl_lst.op(i).op(1));
 			}
@@ -1291,7 +1291,7 @@ void expairseq::drop_coeff_0_terms(epvector::iterator &first_numeric,
 	// move terms with coeff 0 to end and remove them from hashtab
 	// check only those elements which have been touched
 	epp current = seq.begin();
-	unsigned i = 0;
+	size_t i = 0;
 	while (current!=first_numeric) {
 		if (!touched[i]) {
 			++current;
@@ -1387,7 +1387,7 @@ void expairseq::combine_same_terms(void)
 	epvector::iterator first_numeric = seq.end();
 	epvector::iterator last_non_zero = seq.end()-1;
 	
-	unsigned num = seq.size();
+	size_t num = seq.size();
 	std::vector<bool> touched(num);
 	
 	unsigned number_of_zeroes = 0;
@@ -1565,7 +1565,7 @@ epvector * expairseq::subschildren(const lst &ls, const lst &lr, unsigned option
 	// is a product or power. In this case we have to recombine the pairs
 	// because the numeric coefficients may be part of the search pattern.
 	bool complex_subs = false;
-	for (unsigned i=0; i<ls.nops(); ++i) {
+	for (size_t i=0; i<ls.nops(); ++i) {
 		if (is_exactly_a<mul>(ls.op(i)) || is_exactly_a<power>(ls.op(i))) {
 			complex_subs = true;
 			break;
diff --git a/ginac/expairseq.h b/ginac/expairseq.h
index add8a8f9..8f4b6a89 100644
--- a/ginac/expairseq.h
+++ b/ginac/expairseq.h
@@ -89,8 +89,8 @@ public:
 	void print(const print_context & c, unsigned level = 0) const;
 	unsigned precedence(void) const {return 10;}
 	bool info(unsigned inf) const;
-	unsigned nops() const;
-	ex op(int i) const;
+	size_t nops() const;
+	ex op(size_t i) const;
 	ex map(map_function & f) const;
 	ex eval(int level=0) const;
 	ex to_rational(lst &repl_lst) const;
diff --git a/ginac/exprseq_suppl.cpp b/ginac/exprseq_suppl.cpp
index e348853c..3653b582 100644
--- a/ginac/exprseq_suppl.cpp
+++ b/ginac/exprseq_suppl.cpp
@@ -33,11 +33,9 @@ bool exprseq::info(unsigned inf) const
 	return basic::info(inf);
 }
 
-ex exprseq::op(int i) const
+ex exprseq::op(size_t i) const
 {
-	GINAC_ASSERT(i>=0);
 	GINAC_ASSERT(i<nops());
-
 	return seq[i];
 }
 
diff --git a/ginac/fderivative.cpp b/ginac/fderivative.cpp
index d668858e..8216179d 100644
--- a/ginac/fderivative.cpp
+++ b/ginac/fderivative.cpp
@@ -117,7 +117,7 @@ void fderivative::print(const print_context & c, unsigned level) const
 			c.s << *i++ << ",";
 		c.s << *i << std::endl;
 		unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
-		for (unsigned i=0; i<seq.size(); ++i)
+		for (size_t i=0; i<seq.size(); ++i)
 			seq[i].print(c, level + delta_indent);
 		c.s << std::string(level + delta_indent, ' ') << "=====" << std::endl;
 
@@ -180,7 +180,7 @@ ex fderivative::thisexprseq(exvector * vp) const
 ex fderivative::derivative(const symbol & s) const
 {
 	ex result;
-	for (unsigned i=0; i!=seq.size(); i++) {
+	for (size_t i=0; i<seq.size(); i++) {
 		ex arg_diff = seq[i].diff(s);
 		if (!arg_diff.is_zero()) {
 			paramset ps = parameter_set;
diff --git a/ginac/function.pl b/ginac/function.pl
index fd0bc190..fba34e20 100755
--- a/ginac/function.pl
+++ b/ginac/function.pl
@@ -681,7 +681,7 @@ void function::print(const print_context & c, unsigned level) const
 		    << ", nops=" << nops()
 		    << std::endl;
 		unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
-		for (unsigned i=0; i<seq.size(); ++i)
+		for (size_t i=0; i<seq.size(); ++i)
 			seq[i].print(c, level + delta_indent);
 		c.s << std::string(level + delta_indent, ' ') << "=====" << std::endl;
 
@@ -689,8 +689,8 @@ void function::print(const print_context & c, unsigned level) const
 
 		// Print function name in lowercase
 		std::string lname = registered_functions()[serial].name;
-		unsigned num = lname.size();
-		for (unsigned i=0; i<num; i++)
+		size_t num = lname.size();
+		for (size_t i=0; i<num; i++)
 			lname[i] = tolower(lname[i]);
 		c.s << lname << "(";
 
@@ -811,7 +811,7 @@ ${evalf_switch_statement}
 unsigned function::calchash(void) const
 {
 	unsigned v = golden_ratio_hash(golden_ratio_hash(tinfo()) ^ serial);
-	for (unsigned i=0; i<nops(); i++) {
+	for (size_t i=0; i<nops(); i++) {
 		v = rotate_left(v);
 		v ^= this->op(i).gethash();
 	}
@@ -875,8 +875,8 @@ ex function::derivative(const symbol & s) const
 	} else {
 		// Chain rule
 		ex arg_diff;
-		unsigned num = seq.size();
-		for (unsigned i=0; i<num; i++) {
+		size_t num = seq.size();
+		for (size_t i=0; i<num; i++) {
 			arg_diff = seq[i].diff(s);
 			// We apply the chain rule only when it makes sense.  This is not
 			// just for performance reasons but also to allow functions to
@@ -1005,8 +1005,8 @@ void function::store_remember_table(ex const & result) const
 
 unsigned function::register_new(function_options const & opt)
 {
-	unsigned same_name = 0;
-	for (unsigned i=0; i<registered_functions().size(); ++i) {
+	size_t same_name = 0;
+	for (size_t i=0; i<registered_functions().size(); ++i) {
 		if (registered_functions()[i].name==opt.name) {
 			++same_name;
 		}
diff --git a/ginac/idx.cpp b/ginac/idx.cpp
index 2323c29c..277bc316 100644
--- a/ginac/idx.cpp
+++ b/ginac/idx.cpp
@@ -269,13 +269,13 @@ bool idx::info(unsigned inf) const
 	return inherited::info(inf);
 }
 
-unsigned idx::nops() const
+size_t idx::nops() const
 {
 	// don't count the dimension as that is not really a sub-expression
 	return 1;
 }
 
-ex idx::op(int i) const
+ex idx::op(size_t i) const
 {
 	GINAC_ASSERT(i == 0);
 	return value;
@@ -376,7 +376,7 @@ ex idx::subs(const lst & ls, const lst & lr, unsigned options) const
 	GINAC_ASSERT(ls.nops() == lr.nops());
 
 	// First look for index substitutions
-	for (unsigned i=0; i<ls.nops(); i++) {
+	for (size_t i=0; i<ls.nops(); i++) {
 		if (is_equal(ex_to<basic>(ls.op(i)))) {
 
 			// Substitution index->index
diff --git a/ginac/idx.h b/ginac/idx.h
index 7fc743b5..e9587e42 100644
--- a/ginac/idx.h
+++ b/ginac/idx.h
@@ -50,8 +50,8 @@ public:
 public:
 	void print(const print_context & c, unsigned level = 0) const;
 	bool info(unsigned inf) const;
-	unsigned nops() const;
-	ex op(int i) const;
+	size_t nops() const;
+	ex op(size_t i) const;
 	ex map(map_function & f) const;
 	ex evalf(int level = 0) const;
 	ex subs(const lst & ls, const lst & lr, unsigned options = 0) const;
@@ -249,7 +249,7 @@ inline void find_dummy_indices(const exvector & v, exvector & out_dummy)
 }
 
 /** Count the number of dummy index pairs in an index vector. */
-inline unsigned count_dummy_indices(const exvector & v)
+inline size_t count_dummy_indices(const exvector & v)
 {
 	exvector free_indices, dummy_indices;
 	find_free_and_dummy(v.begin(), v.end(), free_indices, dummy_indices);
@@ -257,7 +257,7 @@ inline unsigned count_dummy_indices(const exvector & v)
 }
 
 /** Count the number of dummy index pairs in an index vector. */
-inline unsigned count_free_indices(const exvector & v)
+inline size_t count_free_indices(const exvector & v)
 {
 	exvector free_indices, dummy_indices;
 	find_free_and_dummy(v.begin(), v.end(), free_indices, dummy_indices);
diff --git a/ginac/indexed.cpp b/ginac/indexed.cpp
index 6e53a589..5799642a 100644
--- a/ginac/indexed.cpp
+++ b/ginac/indexed.cpp
@@ -294,7 +294,7 @@ ex indexed::expand(unsigned options) const
 		// expand_indexed expands (a+b).i -> a.i + b.i
 		const ex & base = seq[0];
 		ex sum = _ex0;
-		for (unsigned i=0; i<base.nops(); i++) {
+		for (size_t i=0; i<base.nops(); i++) {
 			exvector s = seq;
 			s[0] = base.op(i);
 			sum += thisexprseq(s).expand();
@@ -459,7 +459,7 @@ exvector indexed::get_free_indices(void) const
 exvector add::get_free_indices(void) const
 {
 	exvector free_indices;
-	for (unsigned i=0; i<nops(); i++) {
+	for (size_t i=0; i<nops(); i++) {
 		if (i == 0)
 			free_indices = op(i).get_free_indices();
 		else {
@@ -475,7 +475,7 @@ exvector mul::get_free_indices(void) const
 {
 	// Concatenate free indices of all factors
 	exvector un;
-	for (unsigned i=0; i<nops(); i++) {
+	for (size_t i=0; i<nops(); i++) {
 		exvector free_indices_of_factor = op(i).get_free_indices();
 		un.insert(un.end(), free_indices_of_factor.begin(), free_indices_of_factor.end());
 	}
@@ -490,7 +490,7 @@ exvector ncmul::get_free_indices(void) const
 {
 	// Concatenate free indices of all factors
 	exvector un;
-	for (unsigned i=0; i<nops(); i++) {
+	for (size_t i=0; i<nops(); i++) {
 		exvector free_indices_of_factor = op(i).get_free_indices();
 		un.insert(un.end(), free_indices_of_factor.begin(), free_indices_of_factor.end());
 	}
@@ -517,8 +517,8 @@ exvector power::get_free_indices(void) const
  *    by the function */
 static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, exvector & local_dummy_indices)
 {
-	unsigned global_size = global_dummy_indices.size(),
-	         local_size = local_dummy_indices.size();
+	size_t global_size = global_dummy_indices.size(),
+	       local_size = local_dummy_indices.size();
 
 	// Any local dummy indices at all?
 	if (local_size == 0)
@@ -528,7 +528,7 @@ static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, ex
 
 		// More local indices than we encountered before, add the new ones
 		// to the global set
-		int old_global_size = global_size;
+		size_t old_global_size = global_size;
 		int remaining = local_size - global_size;
 		exvector::const_iterator it = local_dummy_indices.begin(), itend = local_dummy_indices.end();
 		while (it != itend && remaining > 0) {
@@ -548,10 +548,10 @@ static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, ex
 
 	// Construct lists of index symbols
 	exlist local_syms, global_syms;
-	for (unsigned i=0; i<local_size; i++)
+	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 (unsigned i=0; i<local_size; i++) // don't use more global symbols than necessary
+	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());
 
@@ -662,7 +662,7 @@ ex simplify_indexed_product(const ex & e, exvector & free_indices, exvector & du
 		v.push_back(e.op(0));
 		v.push_back(e.op(0));
 	} else {
-		for (unsigned i=0; i<e.nops(); i++) {
+		for (size_t i=0; i<e.nops(); i++) {
 			ex f = e.op(i);
 			if (is_exactly_a<power>(f) && f.op(1).is_equal(_ex2)) {
 				v.push_back(f.op(0));
@@ -670,7 +670,7 @@ ex simplify_indexed_product(const ex & e, exvector & free_indices, exvector & du
 			} else if (is_exactly_a<ncmul>(f)) {
 				// Noncommutative factor found, split it as well
 				non_commutative = true; // everything becomes noncommutative, ncmul will sort out the commutative factors later
-				for (unsigned j=0; j<f.nops(); j++)
+				for (size_t j=0; j<f.nops(); j++)
 					v.push_back(f.op(j));
 			} else
 				v.push_back(f);
@@ -711,7 +711,7 @@ try_again:
 			// Check whether the two factors share dummy indices
 			exvector free, dummy;
 			find_free_and_dummy(un, free, dummy);
-			unsigned num_dummies = dummy.size();
+			size_t num_dummies = dummy.size();
 			if (num_dummies == 0)
 				continue;
 
@@ -870,7 +870,7 @@ class symminfo {
 public:
 	symminfo() : num(0) {}
 
-	symminfo(const ex & symmterm_, const ex & orig_, unsigned num_) : orig(orig_), num(num_)
+	symminfo(const ex & symmterm_, const ex & orig_, size_t num_) : orig(orig_), num(num_)
 	{
 		if (is_exactly_a<mul>(symmterm_) && is_exactly_a<numeric>(symmterm_.op(symmterm_.nops()-1))) {
 			coeff = symmterm_.op(symmterm_.nops()-1);
@@ -884,7 +884,7 @@ public:
 	ex symmterm;  /**< symmetrized term */
 	ex coeff;     /**< coefficient of symmetrized term */
 	ex orig;      /**< original term */
-	unsigned num; /**< how many symmetrized terms resulted from the original term */
+	size_t num; /**< how many symmetrized terms resulted from the original term */
 };
 
 class symminfo_is_less_by_symmterm {
@@ -941,7 +941,7 @@ ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indi
 		ex sum;
 		free_indices.clear();
 
-		for (unsigned i=0; i<e_expanded.nops(); i++) {
+		for (size_t i=0; i<e_expanded.nops(); i++) {
 			exvector free_indices_of_term;
 			ex term = simplify_indexed(e_expanded.op(i), free_indices_of_term, dummy_indices, sp);
 			if (!term.is_zero()) {
@@ -971,7 +971,7 @@ ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indi
 		}
 
 		// More than one term and more than one dummy index?
-		int num_terms_orig = (is_exactly_a<add>(sum) ? sum.nops() : 1);
+		size_t num_terms_orig = (is_exactly_a<add>(sum) ? sum.nops() : 1);
 		if (num_terms_orig < 2 || dummy_indices.size() < 2)
 			return sum;
 
@@ -983,7 +983,7 @@ ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indi
 		// Chop the sum into terms and symmetrize each one over the dummy
 		// indices
 		std::vector<terminfo> terms;
-		for (unsigned i=0; i<sum.nops(); i++) {
+		for (size_t i=0; i<sum.nops(); i++) {
 			const ex & term = sum.op(i);
 			ex term_symm = term.symmetrize(dummy_syms);
 			if (term_symm.is_zero())
@@ -997,7 +997,7 @@ ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indi
 		// Combine equal symmetrized terms
 		std::vector<terminfo> terms_pass2;
 		for (std::vector<terminfo>::const_iterator i=terms.begin(); i!=terms.end(); ) {
-			unsigned num = 1;
+			size_t num = 1;
 			std::vector<terminfo>::const_iterator j = i + 1;
 			while (j != terms.end() && j->symm == i->symm) {
 				num++;
@@ -1015,8 +1015,8 @@ ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indi
 		std::vector<symminfo> sy;
 		for (std::vector<terminfo>::const_iterator i=terms_pass2.begin(); i!=terms_pass2.end(); ++i) {
 			if (is_exactly_a<add>(i->symm)) {
-				unsigned num = i->symm.nops();
-				for (unsigned j=0; j<num; j++)
+				size_t num = i->symm.nops();
+				for (size_t j=0; j<num; j++)
 					sy.push_back(symminfo(i->symm.op(j), i->orig, num));
 			} else
 				sy.push_back(symminfo(i->symm, i->orig, 1));
@@ -1063,7 +1063,7 @@ ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indi
 			for (std::vector<symminfo>::const_iterator i=sy_pass2.begin(); i!=sy_pass2.end(); ) {
 
 				// How many symmetrized terms of this original term are left?
-				unsigned num = 1;
+				size_t num = 1;
 				std::vector<symminfo>::const_iterator j = i + 1;
 				while (j != sy_pass2.end() && j->orig == i->orig) {
 					num++;
@@ -1224,10 +1224,10 @@ void scalar_products::add(const ex & v1, const ex & v2, const ex & dim, const ex
 void scalar_products::add_vectors(const lst & l, const ex & dim)
 {
 	// Add all possible pairs of products
-	unsigned num = l.nops();
-	for (unsigned i=0; i<num; i++) {
+	size_t num = l.nops();
+	for (size_t i=0; i<num; i++) {
 		ex a = l.op(i);
-		for (unsigned j=0; j<num; j++) {
+		for (size_t j=0; j<num; j++) {
 			ex b = l.op(j);
 			add(a, b, dim, a*b);
 		}
diff --git a/ginac/inifcns.cpp b/ginac/inifcns.cpp
index ffcb229f..73f9c221 100644
--- a/ginac/inifcns.cpp
+++ b/ginac/inifcns.cpp
@@ -324,7 +324,7 @@ static ex Li2_series(const ex &x, const relational &rel, int order, unsigned opt
 			seq.push_back(expair(Li2(x_pt), _ex0));
 			// compute the intermediate terms:
 			ex replarg = series(Li2(x), s==foo, order);
-			for (unsigned i=1; i<replarg.nops()-1; ++i)
+			for (size_t i=1; i<replarg.nops()-1; ++i)
 				seq.push_back(expair((replarg.op(i)/power(s-foo,i)).series(foo==point,1,options).op(0).subs(foo==s),i));
 			// append an order term:
 			seq.push_back(expair(Order(_ex1), replarg.nops()-1));
@@ -454,7 +454,7 @@ ex lsolve(const ex &eqns, const ex &symbols, unsigned options)
 	if (!eqns.info(info_flags::list)) {
 		throw(std::invalid_argument("lsolve(): 1st argument must be a list"));
 	}
-	for (unsigned i=0; i<eqns.nops(); i++) {
+	for (size_t i=0; i<eqns.nops(); i++) {
 		if (!eqns.op(i).info(info_flags::relation_equal)) {
 			throw(std::invalid_argument("lsolve(): 1st argument must be a list of equations"));
 		}
@@ -462,7 +462,7 @@ ex lsolve(const ex &eqns, const ex &symbols, unsigned options)
 	if (!symbols.info(info_flags::list)) {
 		throw(std::invalid_argument("lsolve(): 2nd argument must be a list"));
 	}
-	for (unsigned i=0; i<symbols.nops(); i++) {
+	for (size_t i=0; i<symbols.nops(); i++) {
 		if (!symbols.op(i).info(info_flags::symbol)) {
 			throw(std::invalid_argument("lsolve(): 2nd argument must be a list of symbols"));
 		}
@@ -473,10 +473,10 @@ ex lsolve(const ex &eqns, const ex &symbols, unsigned options)
 	matrix rhs(eqns.nops(),1);
 	matrix vars(symbols.nops(),1);
 	
-	for (unsigned r=0; r<eqns.nops(); r++) {
+	for (size_t r=0; r<eqns.nops(); r++) {
 		const ex eq = eqns.op(r).op(0)-eqns.op(r).op(1); // lhs-rhs==0
 		ex linpart = eq;
-		for (unsigned c=0; c<symbols.nops(); c++) {
+		for (size_t c=0; c<symbols.nops(); c++) {
 			const ex co = eq.coeff(ex_to<symbol>(symbols.op(c)),1);
 			linpart -= co*symbols.op(c);
 			sys(r,c) = co;
@@ -486,7 +486,7 @@ ex lsolve(const ex &eqns, const ex &symbols, unsigned options)
 	}
 	
 	// test if system is linear and fill vars matrix
-	for (unsigned i=0; i<symbols.nops(); i++) {
+	for (size_t i=0; i<symbols.nops(); i++) {
 		vars(i,0) = symbols.op(i);
 		if (sys.has(symbols.op(i)))
 			throw(std::logic_error("lsolve: system is not linear"));
@@ -507,7 +507,7 @@ ex lsolve(const ex &eqns, const ex &symbols, unsigned options)
 	
 	// return list of equations of the form lst(var1==sol1,var2==sol2,...)
 	lst sollist;
-	for (unsigned i=0; i<symbols.nops(); i++)
+	for (size_t i=0; i<symbols.nops(); i++)
 		sollist.append(symbols.op(i)==solution(i,0));
 	
 	return sollist;
diff --git a/ginac/matrix.cpp b/ginac/matrix.cpp
index ecf876c0..fe24e09e 100644
--- a/ginac/matrix.cpp
+++ b/ginac/matrix.cpp
@@ -95,9 +95,9 @@ matrix::matrix(unsigned r, unsigned c, const lst & l)
 {
 	m.resize(r*c, _ex0);
 
-	for (unsigned i=0; i<l.nops(); i++) {
-		unsigned x = i % c;
-		unsigned y = i / c;
+	for (size_t i=0; i<l.nops(); i++) {
+		size_t x = i % c;
+		size_t y = i / c;
 		if (y >= r)
 			break; // matrix smaller than list: throw away excessive elements
 		m[y*c+x] = l.op(i);
@@ -193,24 +193,22 @@ void matrix::print(const print_context & c, unsigned level) const
 }
 
 /** nops is defined to be rows x columns. */
-unsigned matrix::nops() const
+size_t matrix::nops() const
 {
-	return row*col;
+	return static_cast<size_t>(row) * static_cast<size_t>(col);
 }
 
 /** returns matrix entry at position (i/col, i%col). */
-ex matrix::op(int i) const
+ex matrix::op(size_t i) const
 {
-	GINAC_ASSERT(i>=0);
 	GINAC_ASSERT(i<nops());
 	
 	return m[i];
 }
 
 /** returns writable matrix entry at position (i/col, i%col). */
-ex & matrix::let_op(int i)
+ex & matrix::let_op(size_t i)
 {
-	GINAC_ASSERT(i>=0);
 	GINAC_ASSERT(i<nops());
 	
 	ensure_if_modifiable();
@@ -1455,7 +1453,7 @@ int matrix::pivot(unsigned ro, unsigned co, bool symbolic)
 ex lst_to_matrix(const lst & l)
 {
 	// Find number of rows and columns
-	unsigned rows = l.nops(), cols = 0, i, j;
+	size_t rows = l.nops(), cols = 0, i, j;
 	for (i=0; i<rows; i++)
 		if (l.op(i).nops() > cols)
 			cols = l.op(i).nops();
@@ -1474,11 +1472,11 @@ ex lst_to_matrix(const lst & l)
 
 ex diag_matrix(const lst & l)
 {
-	unsigned dim = l.nops();
+	size_t dim = l.nops();
 
 	matrix &m = *new matrix(dim, dim);
 	m.setflag(status_flags::dynallocated);
-	for (unsigned i=0; i<dim; i++)
+	for (size_t i=0; i<dim; i++)
 		m(i, i) = l.op(i);
 
 	return m;
diff --git a/ginac/matrix.h b/ginac/matrix.h
index 7fabe709..8176f882 100644
--- a/ginac/matrix.h
+++ b/ginac/matrix.h
@@ -44,9 +44,9 @@ public:
 	// functions overriding virtual functions from base classes
 public:
 	void print(const print_context & c, unsigned level = 0) const;
-	unsigned nops() const;
-	ex op(int i) const;
-	ex & let_op(int i);
+	size_t nops() const;
+	ex op(size_t i) const;
+	ex & let_op(size_t i);
 	ex eval(int level=0) const;
 	ex evalm(void) const {return *this;}
 	ex subs(const lst & ls, const lst & lr, unsigned options = 0) const;
@@ -98,7 +98,7 @@ protected:
 
 // wrapper functions around member functions
 
-inline unsigned nops(const matrix & m)
+inline size_t nops(const matrix & m)
 { return m.nops(); }
 
 inline ex expand(const matrix & m, unsigned options = 0)
diff --git a/ginac/mul.cpp b/ginac/mul.cpp
index d665db90..fac84d8b 100644
--- a/ginac/mul.cpp
+++ b/ginac/mul.cpp
@@ -181,7 +181,7 @@ void mul::print(const print_context & c, unsigned level) const
 	} else if (is_a<print_python_repr>(c)) {
 		c.s << class_name() << '(';
 		op(0).print(c);
-		for (unsigned i=1; i<nops(); ++i) {
+		for (size_t i=1; i<nops(); ++i) {
 			c.s << ',';
 			op(i).print(c);
 		}
@@ -566,7 +566,7 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con
 	std::vector<bool> subsed(seq.size(), false);
 	exvector subsresult(seq.size());
 
-	for (unsigned i=0; i<ls.nops(); i++) {
+	for (size_t i=0; i<ls.nops(); i++) {
 
 		if (is_exactly_a<mul>(ls.op(i))) {
 
@@ -575,9 +575,9 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con
 			bool succeed = true;
 			lst repls;
 
-			for (unsigned j=0; j<ls.op(i).nops(); j++) {
+			for (size_t j=0; j<ls.op(i).nops(); j++) {
 				bool found=false;
-				for (unsigned k=0; k<nops(); k++) {
+				for (size_t k=0; k<nops(); k++) {
 					if (currsubsed[k] || subsed[k])
 						continue;
 					if (tryfactsubs(op(k), ls.op(i).op(j), nummatches, repls)) {
@@ -595,7 +595,7 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con
 				continue;
 
 			bool foundfirstsubsedfactor = false;
-			for (unsigned j=0; j<subsed.size(); j++) {
+			for (size_t j=0; j<subsed.size(); j++) {
 				if (currsubsed[j]) {
 					if (foundfirstsubsedfactor)
 						subsresult[j] = op(j);
@@ -612,7 +612,7 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con
 			unsigned nummatches = std::numeric_limits<unsigned>::max();
 			lst repls;
 
-			for (unsigned j=0; j<this->nops(); j++) {
+			for (size_t j=0; j<this->nops(); j++) {
 				if (!subsed[j] && tryfactsubs(op(j), ls.op(i), nummatches, repls)) {
 					subsed[j] = true;
 					subsresult[j] = op(j) * power(lr.op(i).subs(ex(repls), subs_options::subs_no_pattern) / ls.op(i).subs(ex(repls), subs_options::subs_no_pattern), nummatches);
@@ -622,7 +622,7 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con
 	}
 
 	bool subsfound = false;
-	for (unsigned i=0; i<subsed.size(); i++) {
+	for (size_t i=0; i<subsed.size(); i++) {
 		if (subsed[i]) {
 			subsfound = true;
 			break;
@@ -632,7 +632,7 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con
 		return basic::subs(ls, lr, options | subs_options::subs_algebraic);
 
 	exvector ev; ev.reserve(nops());
-	for (unsigned i=0; i<nops(); i++) {
+	for (size_t i=0; i<nops(); i++) {
 		if (subsed[i])
 			ev.push_back(subsresult[i]);
 		else
@@ -648,7 +648,7 @@ ex mul::algebraic_subs_mul(const lst & ls, const lst & lr, unsigned options) con
  *  @see ex::diff */
 ex mul::derivative(const symbol & s) const
 {
-	unsigned num = seq.size();
+	size_t num = seq.size();
 	exvector addseq;
 	addseq.reserve(num);
 	
@@ -846,26 +846,7 @@ ex mul::expand(unsigned options) const
 			(cit->coeff.is_equal(_ex1))) {
 			++number_of_adds;
 			if (is_exactly_a<add>(last_expanded)) {
-#if 0
-				// Expand a product of two sums, simple and robust version.
-				const add & add1 = ex_to<add>(last_expanded);
-				const add & add2 = ex_to<add>(cit->rest);
-				const int n1 = add1.nops();
-				const int n2 = add2.nops();
-				ex tmp_accu;
-				exvector distrseq;
-				distrseq.reserve(n2);
-				for (int i1=0; i1<n1; ++i1) {
-					distrseq.clear();
-					// cache the first operand (for efficiency):
-					const ex op1 = add1.op(i1);
-					for (int i2=0; i2<n2; ++i2)
-						distrseq.push_back(op1 * add2.op(i2));
-					tmp_accu += (new add(distrseq))->
-					             setflag(status_flags::dynallocated);
-				}
-				last_expanded = tmp_accu;
-#else
+
 				// Expand a product of two sums, aggressive version.
 				// Caring for the overall coefficients in separate loops can
 				// sometimes give a performance gain of up to 15%!
@@ -917,7 +898,7 @@ ex mul::expand(unsigned options) const
 					tmp_accu += (new add(distrseq, oc))->setflag(status_flags::dynallocated);
 				}
 				last_expanded = tmp_accu;
-#endif
+
 			} else {
 				non_adds.push_back(split_ex_to_pair(last_expanded));
 				last_expanded = cit->rest;
@@ -935,9 +916,9 @@ ex mul::expand(unsigned options) const
 	if (is_exactly_a<add>(last_expanded)) {
 		const add & finaladd = ex_to<add>(last_expanded);
 		exvector distrseq;
-		int n = finaladd.nops();
+		size_t n = finaladd.nops();
 		distrseq.reserve(n);
-		for (int i=0; i<n; ++i) {
+		for (size_t i=0; i<n; ++i) {
 			epvector factors = non_adds;
 			factors.push_back(split_ex_to_pair(finaladd.op(i)));
 			distrseq.push_back((new mul(factors, overall_coeff))->
diff --git a/ginac/ncmul.cpp b/ginac/ncmul.cpp
index d2a2f21f..dc328818 100644
--- a/ginac/ncmul.cpp
+++ b/ginac/ncmul.cpp
@@ -144,15 +144,15 @@ ex ncmul::expand(unsigned options) const
 	intvector positions_of_adds(expanded_seq.size());
 	intvector number_of_add_operands(expanded_seq.size());
 
-	int number_of_adds = 0;
-	int number_of_expanded_terms = 1;
+	size_t number_of_adds = 0;
+	size_t number_of_expanded_terms = 1;
 
-	unsigned current_position = 0;
+	size_t current_position = 0;
 	exvector::const_iterator last = expanded_seq.end();
 	for (exvector::const_iterator cit=expanded_seq.begin(); cit!=last; ++cit) {
 		if (is_exactly_a<add>(*cit)) {
 			positions_of_adds[number_of_adds] = current_position;
-			unsigned num_ops = cit->nops();
+			size_t num_ops = cit->nops();
 			number_of_add_operands[number_of_adds] = num_ops;
 			number_of_expanded_terms *= num_ops;
 			number_of_adds++;
@@ -174,7 +174,7 @@ ex ncmul::expand(unsigned options) const
 
 	while (true) {
 		exvector term = expanded_seq;
-		for (int i=0; i<number_of_adds; i++)
+		for (size_t i=0; i<number_of_adds; i++)
 			term[positions_of_adds[i]] = expanded_seq[positions_of_adds[i]].op(k[i]);
 		distrseq.push_back((new ncmul(term, true))->
 		                    setflag(status_flags::dynallocated | (options == 0 ? status_flags::expanded : 0)));
@@ -251,12 +251,12 @@ ex ncmul::coeff(const ex & s, int n) const
 	return _ex0;
 }
 
-unsigned ncmul::count_factors(const ex & e) const
+size_t ncmul::count_factors(const ex & e) const
 {
 	if ((is_exactly_a<mul>(e)&&(e.return_type()!=return_types::commutative))||
 		(is_exactly_a<ncmul>(e))) {
-		unsigned factors=0;
-		for (unsigned i=0; i<e.nops(); i++)
+		size_t factors=0;
+		for (size_t i=0; i<e.nops(); i++)
 			factors += count_factors(e.op(i));
 		
 		return factors;
@@ -268,7 +268,7 @@ void ncmul::append_factors(exvector & v, const ex & e) const
 {
 	if ((is_exactly_a<mul>(e)&&(e.return_type()!=return_types::commutative))||
 		(is_exactly_a<ncmul>(e))) {
-		for (unsigned i=0; i<e.nops(); i++)
+		for (size_t i=0; i<e.nops(); i++)
 			append_factors(v,e.op(i));
 	} else 
 		v.push_back(e);
@@ -305,7 +305,7 @@ ex ncmul::eval(int level) const
 
 	// ncmul(...,*(x1,x2),...,ncmul(x3,x4),...) ->
 	//     ncmul(...,x1,x2,...,x3,x4,...)  (associativity)
-	unsigned factors = 0;
+	size_t factors = 0;
 	exvector::const_iterator cit = evaledseq.begin(), citend = evaledseq.end();
 	while (cit != citend)
 		factors += count_factors(*cit++);
@@ -325,10 +325,10 @@ ex ncmul::eval(int level) const
 	// determine return types
 	unsignedvector rettypes;
 	rettypes.reserve(assocseq.size());
-	unsigned i = 0;
-	unsigned count_commutative=0;
-	unsigned count_noncommutative=0;
-	unsigned count_noncommutative_composite=0;
+	size_t i = 0;
+	size_t count_commutative=0;
+	size_t count_noncommutative=0;
+	size_t count_noncommutative_composite=0;
 	cit = assocseq.begin(); citend = assocseq.end();
 	while (cit != citend) {
 		switch (rettypes[i] = cit->return_type()) {
@@ -355,8 +355,8 @@ ex ncmul::eval(int level) const
 		commutativeseq.reserve(count_commutative+1);
 		exvector noncommutativeseq;
 		noncommutativeseq.reserve(assocseq.size()-count_commutative);
-		unsigned num = assocseq.size();
-		for (unsigned i=0; i<num; ++i) {
+		size_t num = assocseq.size();
+		for (size_t i=0; i<num; ++i) {
 			if (rettypes[i]==return_types::commutative)
 				commutativeseq.push_back(assocseq[i]);
 			else
@@ -374,7 +374,7 @@ ex ncmul::eval(int level) const
 		// elements in assocseq
 		GINAC_ASSERT(count_commutative==0);
 
-		unsigned assoc_num = assocseq.size();
+		size_t assoc_num = assocseq.size();
 		exvectorvector evv;
 		unsignedvector rttinfos;
 		evv.reserve(assoc_num);
@@ -383,7 +383,7 @@ ex ncmul::eval(int level) const
 		cit = assocseq.begin(), citend = assocseq.end();
 		while (cit != citend) {
 			unsigned ti = cit->return_type_tinfo();
-			unsigned rtt_num = rttinfos.size();
+			size_t rtt_num = rttinfos.size();
 			// search type in vector of known types
 			for (i=0; i<rtt_num; ++i) {
 				if (ti == rttinfos[i]) {
@@ -401,11 +401,11 @@ ex ncmul::eval(int level) const
 			++cit;
 		}
 
-		unsigned evv_num = evv.size();
+		size_t evv_num = evv.size();
 #ifdef DO_GINAC_ASSERT
 		GINAC_ASSERT(evv_num == rttinfos.size());
 		GINAC_ASSERT(evv_num > 0);
-		unsigned s=0;
+		size_t s=0;
 		for (i=0; i<evv_num; ++i)
 			s += evv[i].size();
 		GINAC_ASSERT(s == assoc_num);
@@ -474,13 +474,13 @@ ex ncmul::thisexprseq(exvector * vp) const
  *  @see ex::diff */
 ex ncmul::derivative(const symbol & s) const
 {
-	unsigned num = seq.size();
+	size_t num = seq.size();
 	exvector addseq;
 	addseq.reserve(num);
 	
 	// D(a*b*c) = D(a)*b*c + a*D(b)*c + a*b*D(c)
 	exvector ncmulseq = seq;
-	for (unsigned i=0; i<num; ++i) {
+	for (size_t i=0; i<num; ++i) {
 		ex e = seq[i].diff(s);
 		e.swap(ncmulseq[i]);
 		addseq.push_back((new ncmul(ncmulseq))->setflag(status_flags::dynallocated));
diff --git a/ginac/ncmul.h b/ginac/ncmul.h
index 78a45ed6..b2ea9ca2 100644
--- a/ginac/ncmul.h
+++ b/ginac/ncmul.h
@@ -75,7 +75,7 @@ protected:
 
 	// non-virtual functions in this class
 protected:
-	unsigned count_factors(const ex & e) const;
+	size_t count_factors(const ex & e) const;
 	void append_factors(exvector & v, const ex & e) const;
 	exvector expandchildren(unsigned options) const;
 public:
diff --git a/ginac/normal.cpp b/ginac/normal.cpp
index c9e9352f..297a052b 100644
--- a/ginac/normal.cpp
+++ b/ginac/normal.cpp
@@ -97,7 +97,7 @@ static bool get_first_symbol(const ex &e, const symbol *&x)
 		x = &ex_to<symbol>(e);
 		return true;
 	} else if (is_exactly_a<add>(e) || is_exactly_a<mul>(e)) {
-		for (unsigned i=0; i<e.nops(); i++)
+		for (size_t i=0; i<e.nops(); i++)
 			if (get_first_symbol(e.op(i), x))
 				return true;
 	} else if (is_exactly_a<power>(e)) {
@@ -138,7 +138,7 @@ struct sym_desc {
 	int max_deg;
 
 	/** Maximum number of terms of leading coefficient of symbol in both polynomials */
-	int max_lcnops;
+	size_t max_lcnops;
 
 	/** Commparison operator for sorting */
 	bool operator<(const sym_desc &x) const
@@ -173,7 +173,7 @@ static void collect_symbols(const ex &e, sym_desc_vec &v)
 	if (is_a<symbol>(e)) {
 		add_symbol(&ex_to<symbol>(e), v);
 	} else if (is_exactly_a<add>(e) || is_exactly_a<mul>(e)) {
-		for (unsigned i=0; i<e.nops(); i++)
+		for (size_t i=0; i<e.nops(); i++)
 			collect_symbols(e.op(i), v);
 	} else if (is_exactly_a<power>(e)) {
 		collect_symbols(e.op(0), v);
@@ -234,12 +234,12 @@ static numeric lcmcoeff(const ex &e, const numeric &l)
 		return lcm(ex_to<numeric>(e).denom(), l);
 	else if (is_exactly_a<add>(e)) {
 		numeric c = _num1;
-		for (unsigned i=0; i<e.nops(); i++)
+		for (size_t i=0; i<e.nops(); i++)
 			c = lcmcoeff(e.op(i), c);
 		return lcm(c, l);
 	} else if (is_exactly_a<mul>(e)) {
 		numeric c = _num1;
-		for (unsigned i=0; i<e.nops(); i++)
+		for (size_t i=0; i<e.nops(); i++)
 			c *= lcmcoeff(e.op(i), _num1);
 		return lcm(c, l);
 	} else if (is_exactly_a<power>(e)) {
@@ -271,10 +271,10 @@ static numeric lcm_of_coefficients_denominators(const ex &e)
 static ex multiply_lcm(const ex &e, const numeric &lcm)
 {
 	if (is_exactly_a<mul>(e)) {
-		unsigned num = e.nops();
+		size_t num = e.nops();
 		exvector v; v.reserve(num + 1);
 		numeric lcm_accum = _num1;
-		for (unsigned i=0; i<e.nops(); i++) {
+		for (size_t i=0; i<num; i++) {
 			numeric op_lcm = lcmcoeff(e.op(i), _num1);
 			v.push_back(multiply_lcm(e.op(i), op_lcm));
 			lcm_accum *= op_lcm;
@@ -282,9 +282,9 @@ static ex multiply_lcm(const ex &e, const numeric &lcm)
 		v.push_back(lcm / lcm_accum);
 		return (new mul(v))->setflag(status_flags::dynallocated);
 	} else if (is_exactly_a<add>(e)) {
-		unsigned num = e.nops();
+		size_t num = e.nops();
 		exvector v; v.reserve(num);
-		for (unsigned i=0; i<num; i++)
+		for (size_t i=0; i<num; i++)
 			v.push_back(multiply_lcm(e.op(i), lcm));
 		return (new add(v))->setflag(status_flags::dynallocated);
 	} else if (is_exactly_a<power>(e)) {
@@ -1260,11 +1260,11 @@ ex gcd(const ex &a, const ex &b, ex *ca, ex *cb, bool check_args)
 		if (is_exactly_a<mul>(b) && b.nops() > a.nops())
 			goto factored_b;
 factored_a:
-		unsigned num = a.nops();
+		size_t num = a.nops();
 		exvector g; g.reserve(num);
 		exvector acc_ca; acc_ca.reserve(num);
 		ex part_b = b;
-		for (unsigned i=0; i<num; i++) {
+		for (size_t i=0; i<num; i++) {
 			ex part_ca, part_cb;
 			g.push_back(gcd(a.op(i), part_b, &part_ca, &part_cb, check_args));
 			acc_ca.push_back(part_ca);
@@ -1279,11 +1279,11 @@ factored_a:
 		if (is_exactly_a<mul>(a) && a.nops() > b.nops())
 			goto factored_a;
 factored_b:
-		unsigned num = b.nops();
+		size_t num = b.nops();
 		exvector g; g.reserve(num);
 		exvector acc_cb; acc_cb.reserve(num);
 		ex part_a = a;
-		for (unsigned i=0; i<num; i++) {
+		for (size_t i=0; i<num; i++) {
 			ex part_ca, part_cb;
 			g.push_back(gcd(part_a, b.op(i), &part_ca, &part_cb, check_args));
 			acc_cb.push_back(part_cb);
@@ -1620,15 +1620,15 @@ ex sqrfree_parfrac(const ex & a, const symbol & x)
 	// Factorize denominator and compute cofactors
 	exvector yun = sqrfree_yun(denom, x);
 //clog << "yun factors: " << exprseq(yun) << endl;
-	unsigned num_yun = yun.size();
+	size_t num_yun = yun.size();
 	exvector factor; factor.reserve(num_yun);
 	exvector cofac; cofac.reserve(num_yun);
-	for (unsigned i=0; i<num_yun; i++) {
+	for (size_t i=0; i<num_yun; i++) {
 		if (!yun[i].is_equal(_ex1)) {
-			for (unsigned j=0; j<=i; j++) {
+			for (size_t j=0; j<=i; j++) {
 				factor.push_back(pow(yun[i], j+1));
 				ex prod = _ex1;
-				for (unsigned k=0; k<num_yun; k++) {
+				for (size_t k=0; k<num_yun; k++) {
 					if (k == i)
 						prod *= pow(yun[k], i-j);
 					else
@@ -1638,7 +1638,7 @@ ex sqrfree_parfrac(const ex & a, const symbol & x)
 			}
 		}
 	}
-	unsigned num_factors = factor.size();
+	size_t num_factors = factor.size();
 //clog << "factors  : " << exprseq(factor) << endl;
 //clog << "cofactors: " << exprseq(cofac) << endl;
 
@@ -1647,7 +1647,7 @@ ex sqrfree_parfrac(const ex & a, const symbol & x)
 	matrix sys(max_denom_deg + 1, num_factors);
 	matrix rhs(max_denom_deg + 1, 1);
 	for (int i=0; i<=max_denom_deg; i++) {
-		for (unsigned j=0; j<num_factors; j++)
+		for (size_t j=0; j<num_factors; j++)
 			sys(i, j) = cofac[j].coeff(x, i);
 		rhs(i, 0) = red_numer.coeff(x, i);
 	}
@@ -1656,13 +1656,13 @@ ex sqrfree_parfrac(const ex & a, const symbol & x)
 
 	// Solve resulting linear system
 	matrix vars(num_factors, 1);
-	for (unsigned i=0; i<num_factors; i++)
+	for (size_t i=0; i<num_factors; i++)
 		vars(i, 0) = symbol();
 	matrix sol = sys.solve(vars, rhs);
 
 	// Sum up decomposed fractions
 	ex sum = 0;
-	for (unsigned i=0; i<num_factors; i++)
+	for (size_t i=0; i<num_factors; i++)
 		sum += sol(i, 0) / factor[i];
 
 	return red_poly + sum;
@@ -1689,7 +1689,7 @@ ex sqrfree_parfrac(const ex & a, const symbol & x)
 static ex replace_with_symbol(const ex &e, lst &sym_lst, lst &repl_lst)
 {
 	// Expression already in repl_lst? Then return the assigned symbol
-	for (unsigned i=0; i<repl_lst.nops(); i++)
+	for (size_t i=0; i<repl_lst.nops(); i++)
 		if (repl_lst.op(i).is_equal(e))
 			return sym_lst.op(i);
 	
@@ -1712,7 +1712,7 @@ static ex replace_with_symbol(const ex &e, lst &sym_lst, lst &repl_lst)
 static ex replace_with_symbol(const ex &e, lst &repl_lst)
 {
 	// Expression already in repl_lst? Then return the assigned symbol
-	for (unsigned i=0; i<repl_lst.nops(); i++)
+	for (size_t i=0; i<repl_lst.nops(); i++)
 		if (repl_lst.op(i).op(1).is_equal(e))
 			return repl_lst.op(i).op(0);
 	
@@ -2247,12 +2247,12 @@ static ex find_common_factor(const ex & e, ex & factor, lst & repl)
 {
 	if (is_exactly_a<add>(e)) {
 
-		unsigned num = e.nops();
+		size_t num = e.nops();
 		exvector terms; terms.reserve(num);
 		ex gc;
 
 		// Find the common GCD
-		for (unsigned i=0; i<num; i++) {
+		for (size_t i=0; i<num; i++) {
 			ex x = e.op(i).to_polynomial(repl);
 
 			if (is_exactly_a<add>(x) || is_exactly_a<mul>(x)) {
@@ -2276,16 +2276,16 @@ static ex find_common_factor(const ex & e, ex & factor, lst & repl)
 		factor *= gc;
 
 		// Now divide all terms by the GCD
-		for (unsigned i=0; i<num; i++) {
+		for (size_t i=0; i<num; i++) {
 			ex x;
 
 			// Try to avoid divide() because it expands the polynomial
 			ex &t = terms[i];
 			if (is_exactly_a<mul>(t)) {
-				for (unsigned j=0; j<t.nops(); j++) {
+				for (size_t j=0; j<t.nops(); j++) {
 					if (t.op(j).is_equal(gc)) {
 						exvector v; v.reserve(t.nops());
-						for (unsigned k=0; k<t.nops(); k++) {
+						for (size_t k=0; k<t.nops(); k++) {
 							if (k == j)
 								v.push_back(_ex1);
 							else
@@ -2305,10 +2305,10 @@ term_done:	;
 
 	} else if (is_exactly_a<mul>(e)) {
 
-		unsigned num = e.nops();
+		size_t num = e.nops();
 		exvector v; v.reserve(num);
 
-		for (unsigned i=0; i<num; i++)
+		for (size_t i=0; i<num; i++)
 			v.push_back(find_common_factor(e.op(i), factor, repl));
 
 		return (new mul(v))->setflag(status_flags::dynallocated);
diff --git a/ginac/power.cpp b/ginac/power.cpp
index 7312a692..f2fbe703 100644
--- a/ginac/power.cpp
+++ b/ginac/power.cpp
@@ -228,14 +228,13 @@ bool power::info(unsigned inf) const
 	return inherited::info(inf);
 }
 
-unsigned power::nops() const
+size_t power::nops() const
 {
 	return 2;
 }
 
-ex power::op(int i) const
+ex power::op(size_t i) const
 {
-	GINAC_ASSERT(i>=0);
 	GINAC_ASSERT(i<2);
 
 	return i==0 ? basis : exponent;
@@ -530,7 +529,7 @@ extern bool tryfactsubs(const ex &, const ex &, unsigned &, lst &);
 ex power::subs(const lst & ls, const lst & lr, unsigned options) const
 {
 	if (options & subs_options::subs_algebraic) {
-		for (unsigned i=0; i<ls.nops(); i++) {
+		for (size_t i=0; i<ls.nops(); i++) {
 			unsigned nummatches = std::numeric_limits<unsigned>::max();
 			lst repls;
 			if (tryfactsubs(*this, ls.op(i), nummatches, repls))
@@ -678,7 +677,7 @@ ex power::expand_add(const add & a, int n) const
 	if (n==2)
 		return expand_add_2(a);
 
-	const int m = a.nops();
+	const size_t m = a.nops();
 	exvector result;
 	// The number of terms will be the number of combinatorial compositions,
 	// i.e. the number of unordered arrangement of m nonnegative integers
@@ -690,7 +689,7 @@ ex power::expand_add(const add & a, int n) const
 	intvector upper_limit(m-1);
 	int l;
 
-	for (int l=0; l<m-1; ++l) {
+	for (size_t l=0; l<m-1; ++l) {
 		k[l] = 0;
 		k_cum[l] = 0;
 		upper_limit[l] = n;
@@ -746,10 +745,10 @@ ex power::expand_add(const add & a, int n) const
 		// recalc k_cum[] and upper_limit[]
 		k_cum[l] = (l==0 ? k[0] : k_cum[l-1]+k[l]);
 
-		for (int i=l+1; i<m-1; ++i)
+		for (size_t i=l+1; i<m-1; ++i)
 			k_cum[i] = k_cum[i-1]+k[i];
 
-		for (int i=l+1; i<m-1; ++i)
+		for (size_t i=l+1; i<m-1; ++i)
 			upper_limit[i] = n-k_cum[i-1];
 	}
 
@@ -763,7 +762,7 @@ ex power::expand_add(const add & a, int n) const
 ex power::expand_add_2(const add & a) const
 {
 	epvector sum;
-	unsigned a_nops = a.nops();
+	size_t a_nops = a.nops();
 	sum.reserve((a_nops*(a_nops+1))/2);
 	epvector::const_iterator last = a.seq.end();
 
diff --git a/ginac/power.h b/ginac/power.h
index 7a192c41..fa67350e 100644
--- a/ginac/power.h
+++ b/ginac/power.h
@@ -51,8 +51,8 @@ public:
 	void print(const print_context & c, unsigned level = 0) const;
 	unsigned precedence(void) const {return 60;}
 	bool info(unsigned inf) const;
-	unsigned nops() const;
-	ex op(int i) const;
+	size_t nops() const;
+	ex op(size_t i) const;
 	ex map(map_function & f) const;
 	int degree(const ex & s) const;
 	int ldegree(const ex & s) const;
diff --git a/ginac/pseries.cpp b/ginac/pseries.cpp
index 33b9cfdc..0f1f08d9 100644
--- a/ginac/pseries.cpp
+++ b/ginac/pseries.cpp
@@ -126,8 +126,8 @@ void pseries::print(const print_context & c, unsigned level) const
 		    << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 		    << std::endl;
 		unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
-		unsigned num = seq.size();
-		for (unsigned i=0; i<num; ++i) {
+		size_t num = seq.size();
+		for (size_t i=0; i<num; ++i) {
 			seq[i].rest.print(c, level + delta_indent);
 			seq[i].coeff.print(c, level + delta_indent);
 			c.s << std::string(level + delta_indent, ' ') << "-----" << std::endl;
@@ -141,8 +141,8 @@ void pseries::print(const print_context & c, unsigned level) const
 		c.s << ',';
 		point.print(c);
 		c.s << "),[";
-		unsigned num = seq.size();
-		for (unsigned i=0; i<num; ++i) {
+		size_t num = seq.size();
+		for (size_t i=0; i<num; ++i) {
 			if (i)
 				c.s << ',';
 			c.s << '(';
@@ -254,15 +254,15 @@ int pseries::compare_same_type(const basic & other) const
 }
 
 /** Return the number of operands including a possible order term. */
-unsigned pseries::nops(void) const
+size_t pseries::nops(void) const
 {
 	return seq.size();
 }
 
 /** Return the ith term in the series when represented as a sum. */
-ex pseries::op(int i) const
+ex pseries::op(size_t i) const
 {
-	if (i < 0 || unsigned(i) >= seq.size())
+	if (i >= seq.size())
 		throw (std::out_of_range("op() out of range"));
 
 	return seq[i].rest * power(var - point, seq[i].coeff);
diff --git a/ginac/pseries.h b/ginac/pseries.h
index 17ae934a..308d0457 100644
--- a/ginac/pseries.h
+++ b/ginac/pseries.h
@@ -44,8 +44,8 @@ public:
 public:
 	void print(const print_context & c, unsigned level = 0) const;
 	unsigned precedence(void) const {return 38;} // for clarity just below add::precedence
-	unsigned nops(void) const;
-	ex op(int i) const;
+	size_t nops(void) const;
+	ex op(size_t i) const;
 	int degree(const ex &s) const;
 	int ldegree(const ex &s) const;
 	ex coeff(const ex &s, int n = 1) const;
diff --git a/ginac/relational.cpp b/ginac/relational.cpp
index 7c395eea..70dc5bd0 100644
--- a/ginac/relational.cpp
+++ b/ginac/relational.cpp
@@ -160,14 +160,13 @@ bool relational::info(unsigned inf) const
 	return 0;
 }
 
-unsigned relational::nops() const
+size_t relational::nops() const
 {
 	return 2;
 }
 
-ex relational::op(int i) const
+ex relational::op(size_t i) const
 {
-	GINAC_ASSERT(i>=0);
 	GINAC_ASSERT(i<2);
 
 	return i==0 ? lh : rh;
diff --git a/ginac/relational.h b/ginac/relational.h
index 133a25f3..c8044f58 100644
--- a/ginac/relational.h
+++ b/ginac/relational.h
@@ -54,8 +54,8 @@ public:
 	void print(const print_context & c, unsigned level = 0) const;
 	unsigned precedence(void) const {return 20;}
 	bool info(unsigned inf) const;
-	unsigned nops() const;
-	ex op(int i) const;
+	size_t nops() const;
+	ex op(size_t i) const;
 	ex map(map_function & f) const;
 	ex eval(int level=0) const;
 
diff --git a/ginac/remember.cpp b/ginac/remember.cpp
index cdd346f7..cc8a6b88 100644
--- a/ginac/remember.cpp
+++ b/ginac/remember.cpp
@@ -44,8 +44,8 @@ bool remember_table_entry::is_equal(function const & f) const
 {
 	GINAC_ASSERT(f.seq.size()==seq.size());
 	if (f.gethash()!=hashvalue) return false;
-	unsigned num = seq.size();
-	for (unsigned i=0; i<num; ++i)
+	size_t num = seq.size();
+	for (size_t i=0; i<num; ++i)
 		if (!seq[i].is_equal(f.seq[i])) return false;
 	++last_access = access_counter;
 	++successful_hits;
diff --git a/ginac/symbol.cpp b/ginac/symbol.cpp
index 177e00ab..af98215a 100644
--- a/ginac/symbol.cpp
+++ b/ginac/symbol.cpp
@@ -110,7 +110,7 @@ ex symbol::unarchive(const archive_node &n, const lst &sym_lst)
 	ex s = (new symbol(n, sym_lst))->setflag(status_flags::dynallocated);
 	
 	// If symbol is in sym_lst, return the existing symbol
-	for (unsigned i=0; i<sym_lst.nops(); i++) {
+	for (size_t i=0; i<sym_lst.nops(); i++) {
 		if (is_a<symbol>(sym_lst.op(i)) && (ex_to<symbol>(sym_lst.op(i)).name == ex_to<symbol>(s).name))
 			return sym_lst.op(i);
 	}
diff --git a/ginac/symmetry.cpp b/ginac/symmetry.cpp
index f5623608..473bd5b9 100644
--- a/ginac/symmetry.cpp
+++ b/ginac/symmetry.cpp
@@ -205,8 +205,8 @@ void symmetry::print(const print_context & c, unsigned level) const
 				default: c.s << '?'; break;
 			}
 			c.s << '(';
-			unsigned num = children.size();
-			for (unsigned i=0; i<num; i++) {
+			size_t num = children.size();
+			for (size_t i=0; i<num; i++) {
 				children[i].print(c);
 				if (i != num - 1)
 					c.s << ",";
@@ -428,7 +428,7 @@ ex ex::symmetrize(const lst & l) const
 {
 	exvector v;
 	v.reserve(l.nops());
-	for (unsigned i=0; i<l.nops(); i++)
+	for (size_t i=0; i<l.nops(); i++)
 		v.push_back(l.op(i));
 	return symm(*this, v.begin(), v.end(), false);
 }
@@ -438,7 +438,7 @@ ex ex::antisymmetrize(const lst & l) const
 {
 	exvector v;
 	v.reserve(l.nops());
-	for (unsigned i=0; i<l.nops(); i++)
+	for (size_t i=0; i<l.nops(); i++)
 		v.push_back(l.op(i));
 	return symm(*this, v.begin(), v.end(), true);
 }
@@ -449,7 +449,7 @@ ex ex::symmetrize_cyclic(const lst & l) const
 {
 	exvector v;
 	v.reserve(l.nops());
-	for (unsigned i=0; i<l.nops(); i++)
+	for (size_t i=0; i<l.nops(); i++)
 		v.push_back(l.op(i));
 	return GiNaC::symmetrize_cyclic(*this, v.begin(), v.end());
 }
diff --git a/ginac/tensor.cpp b/ginac/tensor.cpp
index 08d368dd..eef89a1e 100644
--- a/ginac/tensor.cpp
+++ b/ginac/tensor.cpp
@@ -313,13 +313,13 @@ ex tensepsilon::eval_indexed(const basic & i) const
 		// a canonic order but we can't assume what exactly that order is)
 		std::vector<int> v;
 		v.reserve(i.nops() - 1);
-		for (unsigned j=1; j<i.nops(); j++)
+		for (size_t j=1; j<i.nops(); j++)
 			v.push_back(ex_to<numeric>(ex_to<idx>(i.op(j)).get_value()).to_int());
 		int sign = permutation_sign(v.begin(), v.end());
 
 		// In a Minkowski space, check for covariant indices
 		if (minkowski) {
-			for (unsigned j=1; j<i.nops(); j++) {
+			for (size_t j=1; j<i.nops(); j++) {
 				const ex & x = i.op(j);
 				if (!is_a<varidx>(x))
 					throw(std::runtime_error("indices of epsilon tensor in Minkowski space must be of type varidx"));
@@ -347,7 +347,7 @@ bool tensor::replace_contr_index(exvector::iterator self, exvector::iterator oth
 
 again:
 	if (self_idx->is_symbolic()) {
-		for (unsigned i=1; i<other->nops(); i++) {
+		for (size_t i=1; i<other->nops(); i++) {
 			const idx &other_idx = ex_to<idx>(other->op(i));
 			if (is_dummy_pair(*self_idx, other_idx)) {
 
@@ -462,7 +462,7 @@ bool spinmetric::contract_with(exvector::iterator self, exvector::iterator other
 
 again:
 	if (self_idx->is_symbolic()) {
-		for (unsigned i=1; i<other->nops(); i++) {
+		for (size_t i=1; i<other->nops(); i++) {
 			const idx &other_idx = ex_to<idx>(other->op(i));
 			if (is_dummy_pair(*self_idx, other_idx)) {
 
@@ -495,15 +495,15 @@ bool tensepsilon::contract_with(exvector::iterator self, exvector::iterator othe
 	GINAC_ASSERT(is_a<indexed>(*self));
 	GINAC_ASSERT(is_a<indexed>(*other));
 	GINAC_ASSERT(is_a<tensepsilon>(self->op(0)));
-	unsigned num = self->nops() - 1;
+	size_t num = self->nops() - 1;
 
 	if (is_exactly_a<tensepsilon>(other->op(0)) && num+1 == other->nops()) {
 
 		// Contraction of two epsilon tensors is a determinant
 		bool variance = is_a<varidx>(self->op(1));
 		matrix M(num, num);
-		for (unsigned i=0; i<num; i++) {
-			for (unsigned j=0; j<num; j++) {
+		for (size_t i=0; i<num; i++) {
+			for (size_t j=0; j<num; j++) {
 				if (minkowski)
 					M(i, j) = lorentz_g(self->op(i+1), other->op(j+1), pos_sig);
 				else if (variance)
diff --git a/ginsh/ginsh_parser.yy b/ginsh/ginsh_parser.yy
index d1485943..02318cbd 100644
--- a/ginsh/ginsh_parser.yy
+++ b/ginsh/ginsh_parser.yy
@@ -352,9 +352,9 @@ static ex f_determinant(const exprseq &e)
 
 static ex f_diag(const exprseq &e)
 {
-	unsigned dim = e.nops();
+	size_t dim = e.nops();
 	matrix &m = *new matrix(dim, dim);
-	for (unsigned i=0; i<dim; i++)
+	for (size_t i=0; i<dim; i++)
 		m.set(i, i, e.op(i));
 	return m;
 }
-- 
2.47.0