X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=blobdiff_plain;f=ginac%2Fmatrix.cpp;h=ca574cc32fd9f7355fe378dded9c9bde50cf9da8;hp=8d2e0b88b99e9fb7f8e5fb5e374b7c2dfcfd03e8;hb=1b5601580eba6925b0af9408c56d6b2027b56013;hpb=ed21ddd5e2bc0af018c10934342f526d0ae4b7a7

diff --git a/ginac/matrix.cpp b/ginac/matrix.cpp
index 8d2e0b88..ca574cc3 100644
--- a/ginac/matrix.cpp
+++ b/ginac/matrix.cpp
@@ -141,7 +141,7 @@ void matrix::archive(archive_node &n) const
 DEFAULT_UNARCHIVE(matrix)
 
 //////////
-// functions overriding virtual functions from bases classes
+// functions overriding virtual functions from base classes
 //////////
 
 // public
@@ -259,6 +259,16 @@ int matrix::compare_same_type(const basic & other) const
 	return 0;
 }
 
+bool matrix::match_same_type(const basic & other) const
+{
+	GINAC_ASSERT(is_exactly_of_type(other, matrix));
+	const matrix & o = static_cast<const matrix &>(other);
+	
+	// The number of rows and columns must be the same. This is necessary to
+	// prevent a 2x3 matrix from matching a 3x2 one.
+	return row == o.rows() && col == o.cols();
+}
+
 /** Automatic symbolic evaluation of an indexed matrix. */
 ex matrix::eval_indexed(const basic & i) const
 {
@@ -405,10 +415,8 @@ bool matrix::contract_with(exvector::iterator self, exvector::iterator other, ex
 	const matrix &other_matrix = ex_to<matrix>(other->op(0));
 
 	if (self->nops() == 2) {
-		unsigned self_dim = (self_matrix.col == 1) ? self_matrix.row : self_matrix.col;
 
 		if (other->nops() == 2) { // vector * vector (scalar product)
-			unsigned other_dim = (other_matrix.col == 1) ? other_matrix.row : other_matrix.col;
 
 			if (self_matrix.col == 1) {
 				if (other_matrix.col == 1) {
@@ -604,8 +612,11 @@ matrix matrix::pow(const ex & expn) const
 			for (unsigned r=0; r<row; ++r)
 				result(r,r) = _ex1();
 			numeric b(1);
-			// this loop computes the representation of k in base 2 and
-			// multiplies the factors whenever needed:
+			// This loop computes the representation of k in base 2 from right
+			// to left(!) and multiplies the factors whenever needed.  Note
+			// that this is not entirely optimal but close to optimal and
+			// "better" algorithms are much harder to implement.  (See Knuth,
+			// TAoCP2, section "Evaluation of Powers" for a good discussion.)
 			while (b.compare(k)<=0) {
 				b *= numeric(2);
 				numeric r(mod(k,b));