File: factor_out.cc

package info (click to toggle)
cadabra2 2.4.3.2-2
  • links: PTS, VCS
  • area: main
  • in suites: sid
  • size: 78,732 kB
  • sloc: ansic: 133,450; cpp: 92,064; python: 1,530; javascript: 203; sh: 184; xml: 182; objc: 53; makefile: 51
file content (200 lines) | stat: -rw-r--r-- 5,890 bytes parent folder | download | duplicates (3)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200

#include "Cleanup.hh"
#include "Functional.hh"
#include "algorithms/factor_out.hh"
#include "algorithms/sort_product.hh"
#include <map>

using namespace cadabra;

factor_out::factor_out(const Kernel& k, Ex& e, Ex& args, bool right)
	: Algorithm(k, e), to_right(right)
	{
	cadabra::do_list(args, args.begin(), [&](Ex::iterator arg) {
		to_factor_out.push_back(Ex(arg));
		return true;
		}
	                );
	}

/// Check if the expression is a sum with more than one term
bool factor_out::can_apply(iterator st)
	{
	if(*st->name=="\\sum")  return true;
	return false;
	}

Algorithm::result_t factor_out::apply(iterator& it)
	{
	result_t result=result_t::l_no_action;

	// For every term in the sum, we look at the factors in the product
	// (or at the single object if there is no product). If this factor
	// needs to be factored out, we determine if it can be moved all the
	// way to the left of the expression. If so, move the object to
	// a 'factored_out' temporary, and take out of the tree. Rinse/repeat.
	// What's left at the end is two objects: the stuff factored out,
	// and the rest. Look up if we already have 'the stuff factored out'.
	// If not, create new. If so, add this term.

	Ex_comparator comparator(kernel.properties);

	typedef std::pair<Ex, std::vector<Ex> > new_term_t;
	std::vector<new_term_t> new_terms;

	auto term=tr.begin(it);
	while(term!=tr.end(it)) {
		auto next_term=term;
		++next_term;

		iterator prod=term;
		prod_wrap_single_term(prod);

		Ex collector("\\prod"); // collect all factors that we have taken out

		// Insert a dummy symbol at the very front or back.
		// FIXME: there is now a 'can_move_to_front', use that.
		iterator dummy;
		if(to_right) dummy = tr.append_child(prod, str_node("dummy"));
		else         dummy = tr.prepend_child(prod, str_node("dummy"));

		// Look at all factors in turn and determine if they should be taken out.
		if(to_right) {
			auto fac=tr.end(prod);
			auto next=fac;
			--next;
			do {
				fac=next;
				--next;
				for(size_t i=0; i<to_factor_out.size(); ++i) {
					auto match=comparator.equal_subtree(fac, to_factor_out[i].begin());
					if(match==Ex_comparator::match_t::subtree_match) {
						int sign=comparator.can_move_adjacent(prod, dummy, fac, false);
						if(sign!=0) {
							collector.append_child(collector.begin(), iterator(fac));
							multiply(prod->multiplier, sign);
							next=tr.erase(fac);
							result=result_t::l_applied;
							break;
							}
						}
					}
				}
			while(fac!=tr.begin(prod));
			}
		else {
			auto fac=tr.begin(prod);
			while(fac!=tr.end(prod)) {
				auto next=fac;
				++next;
				for(size_t i=0; i<to_factor_out.size(); ++i) {
					auto match=comparator.equal_subtree(fac, to_factor_out[i].begin());
					if(match==Ex_comparator::match_t::subtree_match) {
						int sign=comparator.can_move_adjacent(prod, dummy, fac, true);
						if(sign!=0) {
							collector.append_child(collector.begin(), iterator(fac));
							multiply(prod->multiplier, sign);
							next=tr.erase(fac);
							result=result_t::l_applied;
							break;
							}
						}
					}
				fac=next;
				}
			}

		tr.erase(dummy);
		if(tr.number_of_children(prod)==0)
			tr.append_child(prod, str_node("1"));

		// std::cerr << "product after factoring out " << Ex(prod) << std::endl;

		if(collector.number_of_children(collector.begin())!=0) {
			// The stuff factored out of this term is in 'collector'. See if we have
			// factored out that thing before. Because we may not always have collected
			// factors in the same order (the original expression may not have had
			// its product sorted), we first sort the collector product.

			sort_product sp(kernel, collector);
			sp.dont_cleanup(); // otherwise single-factor products will get stripped of the \prod wrapper.
			auto coltop=collector.begin();
			if(sp.can_apply(coltop)) {
				sp.apply(coltop);
				}
			multiply(prod->multiplier, *coltop->multiplier);
			one(coltop->multiplier);

			// Scan through the things factored out so far.
			bool found=false;
			for(auto& nt: new_terms) {
				if(nt.first==collector) { // have that factored out already, add the other factors
					nt.second.push_back(Ex(prod));
					found=true;
					break;
					}
				}
			// We hadn't factored this bit out before, make a new term.
			if(!found) {
				std::vector<Ex> v;
				v.push_back(Ex(prod));
				new_term_t nt(collector, v);
				new_terms.push_back(nt);
				}

			// All info is now in new_terms; can remove the original.
			tr.erase(prod);
			}
		else {
			prod_unwrap_single_term(prod);
			}
		term=next_term;
		}

	// Everything has been collected now into new_terms. Expand those out
	// into a proper sum of products.

	for(auto& nt: new_terms) {
		auto prod = tr.append_child(it, nt.first.begin());
		if(nt.second.size()==1) { // factored, but only one term found.
			auto top = nt.second[0].begin(); // prod node
			if(to_right) {
				auto ins = tr.end(top);
				--ins;
				while(tr.is_valid(ins)) {
					tr.prepend_child(prod, iterator(ins));
					--ins;
					}
				}
			else {
				auto ins = tr.begin(top);
				while(ins!=tr.end(top)) {
					tr.append_child(prod, iterator(ins));
					++ins;
					}
				}
			multiply(prod->multiplier, *(nt.second[0].begin()->multiplier));
			// FIXME: append_children has a BUG! Messes up the tree. But it is needed to
			// handle terms where the sub-factor is not a simple element.
			//			tr.append_children(prod, nt.second[0].begin(top), nt.second[0].end(top));

			cleanup_dispatch(kernel, tr, prod);
			}
		else {
			iterator sum;
			if(to_right)
				sum = tr.prepend_child(prod, str_node("\\sum"));
			else
				sum = tr.append_child(prod, str_node("\\sum"));
			for(auto& term: nt.second) {
				auto tmp = tr.append_child(sum, term.begin());
				cleanup_dispatch(kernel, tr, tmp);
				}
			}
		}

	// std::cerr << "end of factor_out: \n" << Ex(it) << std::endl;

	return result;
	}