fixed issue where some zeroes were not being removed causing matching errors

Author: Dominic Price

core/algorithms/young_reduce.cc

@@ -159,21 +159,21 @@ namespace cadabra {
 			return factor;
 		}
 
-		void ProjectedForm::combine(const ProjectedForm& other, mpq_class factor)
+		void ProjectedForm::combine(const ProjectedForm& other)
 		{
-			if (factor == 1) {
-				for (const auto& kv : other.data) {
-					data[kv.first] += kv.second;
-					if (data[kv.first] == 0)
-						data.erase(kv.first);
-				}
+			for (const auto& kv : other.data) {
+				data[kv.first] += kv.second;
+				if (data[kv.first] == 0)
+					data.erase(kv.first);
 			}
-			else {
-				for (const auto& kv : other.data) {
-					data[kv.first] += kv.second * factor;
-					if (data[kv.first] == 0)
-						data.erase(kv.first);
-				}
+		}
+
+		void ProjectedForm::combine(const ProjectedForm& other, mpq_class factor)
+		{
+			for (const auto& kv : other.data) {
+				data[kv.first] += kv.second * factor;
+				if (data[kv.first] == 0)
+					data.erase(kv.first);
 			}
 		}
 
@@ -190,7 +190,8 @@ namespace cadabra {
 
 		void ProjectedForm::insert(adjform_t adjform, mpq_class value)
 		{
-			data[adjform] = value;
+			if (value != 0)
+				data[adjform] = value;
 		}
 
 		void ProjectedForm::apply_young_symmetry(const adjform_t& indices, bool antisymmetric)
@@ -221,6 +222,8 @@ namespace cadabra {
 					}
 					cdebug << "\tMade term " << adjform_to_string(ret) << " * " << (parity * kv.second) << '\n';
 					data[ret] += parity * kv.second;
+					if (data[ret] == 0)
+						data.erase(ret);
 				} while (swaps = next_perm(perm));
 			}
 		}
@@ -258,8 +261,6 @@ namespace cadabra {
 			Ex::iterator l1 = lhs.begin(), l2 = lhs.end();
 			Ex::iterator r1 = rhs.begin(), r2 = rhs.end();
 
-			std::vector<Ex::iterator> l_indices, r_indices;
-
 			// Loop over all tree nodes using a depth first iterator. If the
 			// entry is an index ensure that it has the same parent_rel, if it
 			// is any other type of node check that the names match.
@@ -270,8 +271,6 @@ namespace cadabra {
 					if (l1->fl.parent_rel != r1->fl.parent_rel) {
 						return false;
 					}
-					l_indices.push_back(l1);
-					r_indices.push_back(r1);
 				}
 				else {
 					if (l1->name != r1->name || l1->multiplier != r1->multiplier) {
@@ -284,6 +283,18 @@ namespace cadabra {
 			return l1 == l2 && r1 == r2;
 		}
 
+		bool has_TableauBase(Ex::iterator it, const cadabra::Kernel& kernel)
+		{
+			if (*it->name == "\\prod" || *it->name == "\\sum") {
+				for (Ex::sibling_iterator beg = it.begin(), end = it.end(); beg != end; ++beg)
+					if (has_TableauBase(beg, kernel))
+						return true;
+			}
+			else {
+				return (kernel.properties.get_composite<cadabra::TableauBase>(it) != nullptr);
+			}
+		}
+
 		std::vector<Ex::iterator> split_ex(Ex::iterator it, const std::string& delim)
 		{
 			if (*it->name == delim) {
@@ -368,7 +379,12 @@ young_reduce::~young_reduce()
 
 bool young_reduce::can_apply(iterator it)
 {
-	return true;
+	if (pat == Ex::iterator()) {
+		// check for TableauBase
+		return has_TableauBase(it, kernel);
+	}
+	else
+		return true;
 }
 
 young_reduce::result_t young_reduce::apply(iterator& it)
@@ -513,17 +529,13 @@ bool young_reduce::set_pattern(Ex::iterator new_pat)
 		throw std::runtime_error("pat is empty");
 
 	cdebug << "Checking for TableauBase property...";
-	bool has_tableau_symmetry = false;
-	for (const auto& term: collect) {
-		if (kernel.properties.get_composite<TableauBase>(term) != nullptr) {
-			has_tableau_symmetry = true;
-			break;
-		}
-	}
-	cdebug << (has_tableau_symmetry ? "true!" : "false - exiting...") << '\n';
-
-	if (!has_tableau_symmetry) 
+	if (!has_TableauBase(new_pat, kernel)) {
+		cdebug << "false, returning...\n";
 		return false;
+	}
+	else {
+		cdebug << "true!\n";
+	}
 
 	pat = new_pat;
 	pat_sym = symmetrize(new_pat);	

core/algorithms/young_reduce.hh

@@ -22,7 +22,8 @@ namespace cadabra {
 			mpq_class compare(const ProjectedForm& other) const;
 
 			// Add all contributions from 'other' into 'this'
-			void combine(const ProjectedForm& other, mpq_class factor = 1);
+			void combine(const ProjectedForm& other);
+			void combine(const ProjectedForm& other, mpq_class factor);
 
 			// Multiply all terms by a constant factor
 			void multiply(mpq_class k);
@@ -47,17 +48,17 @@ namespace cadabra {
 		// Checks if 'lhs' and 'rhs' are identical up to index names
 		bool check_structure(Ex::iterator lhs, Ex::iterator rhs);
 
+		// Check for a TableauBase property (potentially deep)
+		bool has_TableauBase(Ex::iterator it, const Kernel& kernel);
+
 		// Returns a list of iterators of the children of 'it' if it->name
 		// is delim, otherwise returns a list containing 'it' as its only
 		// entry. If 'pat' is specified, any terms not matching 'pat' are
 		// ignored
 		std::vector<Ex::iterator> split_ex(Ex::iterator it, const std::string& delim);
 		std::vector<Ex::iterator> split_ex(Ex::iterator it, const std::string& delim, Ex::iterator pat);
 
-		// Rewrite adjform type indices as dummy indices
 		adjform_t collapse_dummy_indices(adjform_t adjform);
-
-		// Rewrite dummy indices as adjform type indices
 		adjform_t expand_dummy_indices(adjform_t adjform);
 	}
 

tests/youngreduce.cdb

@@ -6,7 +6,6 @@ def test01():
    S_{a b c}::Symmetric;
    ex:= S_{a b c}S_{b a c} + S_{b a c}S_{c a b}:
    young_reduce(ex, $S_{a b c}S_{a b c}$)
-   display(ex)
    assert(ex == $2S_{a b c}S_{a b c}$)
    print("Test 01 passed")
 
@@ -138,7 +137,8 @@ def test13():
    R_{a b c d}::RiemannTensor.
    A^{a b c}::AntiSymmetric.
    ex:= R_{a b c d} A^{a b c}:
-   young_reduce(ex, $R_{a b c d} A^{a b c}$)
+   young_reduce(ex)
+   #young_reduce(ex, $R_{a b c d} A^{a b c}$)
    assert ex==0
    print("Test 13 passed")
 
@@ -147,8 +147,7 @@ test13()
 def test14():
 	__cdbkernel__ = create_scope()
 	A_{a b c d e}::AntiSymmetric.
-	ex = young_reduce($9A_{b a c e d} - 3A_{e d c b a}$, $A_{b a c e d}$)
-	display(ex)
+	ex = young_reduce($9A_{b a c e d} - 3A_{e d c b a}$, $7A_{b a c e d}$)
 	assert ex == $6A_{b a c e d}$
 	print("Test 14 passed")
 
@@ -164,3 +163,12 @@ def test15():
 	print("Test 15 passed")
 
 test15()
+
+def test16():
+	__cdbkernel__ = create_scope()
+	S_{a b c}::Symmetric.
+	assert young_reduce($a + b + c$) == $a + b + c$
+	assert young_reduce($0$) == $0$
+	print("Test 16 passed")
+
+test16()