1 /**
2  * Data structure for optimized "sparse" N-dimensional matrices
3  *
4  * License:
5  *   This Source Code Form is subject to the terms of
6  *   the Mozilla Public License, v. 2.0. If a copy of
7  *   the MPL was not distributed with this file, You
8  *   can obtain one at http://mozilla.org/MPL/2.0/.
9  *
10  * Authors:
11  *   Vladimir Panteleev <ae@cy.md>
12  */
13 
14 module ae.utils.mapset;
15 
16 static if (__VERSION__ >= 2083):
17 
18 import std.algorithm.iteration;
19 import std.algorithm.searching;
20 import std.algorithm.sorting;
21 import std.array;
22 import std.exception;
23 import std.typecons : tuple;
24 
25 import ae.utils.aa;
26 import ae.utils.array : amap;
27 
28 /**
29    Data structure for holding optimized "sparse" N-dimensional matrices.
30    The number of dimensions is arbitrary and can be varied at runtime.
31 
32    Unlike classical sparse arrays/matrices, this data structure is
33    optimized not just for arrays which are mostly zero (or any given
34    value), but for any matrices where data is likely to repeat across
35    dimensions.
36 
37    This is done by storing the data as a tree, with each depth
38    representing one dimension. A sub-tree is, thus, a sub-matrix
39    (slice of the matrix represented by the top-level tree) along the
40    top-level dimension. Each sub-tree is individually immutable, and
41    can therefore be shared within the same tree or even by several
42    instances of this type.
43 
44    Dimension order need not be the same for all sub-trees. Even the
45    number of dimensions in siblings' sub-trees need not be the same;
46    i.e, the data structure is itself sparse with regards to the
47    dimensions.
48 
49    As there is no explicit "value" type (associated with every set of
50    coordinates), this data structure can be seen as a representation
51    of a set of points, but it can be simulated by adding a "value"
52    dimension (in the same way how a gray-scale image can be
53    represented as a set of 3D points like a height-map).
54 
55    An alternative way to look at this data structure is as a set of
56    maps (i.e. associative arrays). Thus, each layer of the tree stores
57    one possible map key, and under it, all the maps that have this key
58    with a specific value. Note that for this use case, DimValue needs
59    to have a value (configured using nullValue) which indicates the
60    absence of a certain key (dimension) in a map.
61 
62    Params:
63      DimName = The type used to indicate the name (or index) of a dimension.
64      DimValue = The type for indicating a point on a dimension's axis.
65      nullValue = A value for DimValue indicating that a dimension is unset
66  */
67 struct MapSet(DimName, DimValue, DimValue nullValue = DimValue.init)
68 {
69 	/// Logically, each MapSet node has a map of values to a subset.
70 	/// However, it is faster to represent that map as an array of key-value pairs
71 	/// rather than a D associative array, so we do that here.
72 	struct Pair
73 	{
74 		DimValue value; ///
75 		MapSet set; ///
76 
77 		int opCmp(ref const typeof(this) other) const
78 		{
79 			static if (is(typeof(value.opCmp(other.value)) : int))
80 				return value.opCmp(other.value);
81 			else
82 				return value < other.value ? -1 : value > other.value ? 1 : 0;
83 		} ///
84 	}
85 
86 	struct Node
87 	{
88 		DimName dim; ///
89 		Pair[] children; ///
90 
91 		immutable this(DimName dim, immutable Pair[] children)
92 		{
93 			// Zero children doesn't make sense.
94 			// It would be the equivalent of an empty set,
95 			// but then we should just use MapSet.emptySet instead.
96 			assert(!children.empty, "Node with zero children");
97 
98 			// Nodes should be in their canonical form for
99 			// memoization and deduplication to be effective.
100 			assert(children.isSorted, "Children are not sorted");
101 
102 			this.dim = dim;
103 			this.children = children;
104 
105 			// Because each subset is immutable, we can
106 			// pre-calculate the hash during construction.
107 			hash = hashOf(dim) ^ hashOf(children);
108 
109 			size_t totalMembers = 0;
110 			foreach (i, ref pair; children)
111 			{
112 				if (i)
113 					assert(pair.value != children[i-1].value, "Duplicate value");
114 
115 				pair.set.assertDeduplicated();
116 				// Same as "Node with zero children"
117 				assert(pair.set !is emptySet, "Empty set as subset");
118 
119 				totalMembers += pair.set.count;
120 			}
121 			this.totalMembers = totalMembers;
122 		} ///
123 
124 		immutable this(DimName dim, immutable MapSet[DimValue] children)
125 		{
126 			auto childrenList = children.byKeyValue.map!(kv => Pair(kv.key, kv.value)).array;
127 			childrenList.sort();
128 			this(dim, cast(immutable) childrenList);
129 		} ///
130 
131 		void toString(scope void delegate(const(char)[]) sink) const
132 		{
133 			import std.format : formattedWrite;
134 			sink.formattedWrite("{ %(%s%), [ ", (&dim)[0..1]);
135 			bool first = true;
136 			foreach (ref pair; children)
137 			{
138 				if (first)
139 					first = false;
140 				else
141 					sink(", ");
142 				sink.formattedWrite("%s : %s", pair.value, pair.set);
143 			}
144 			sink(" ] }");
145 		} ///
146 
147 		private hash_t hash;
148 		private size_t totalMembers;
149 
150 		hash_t toHash() const
151 		{
152 			return hash;
153 		} ///
154 
155 		bool opEquals(ref const typeof(this) s) const
156 		{
157 			return hash == s.hash && dim == s.dim && children == s.children;
158 		} ///
159 	} ///
160 
161 	/// Indicates the empty set.
162 	/// May only occur at the top level (never as a subset).
163 	private enum emptySetRoot = cast(immutable(Node)*)1;
164 
165 	/// If emptySetRoot, zero values.
166 	/// If null, one value with zero dimensions.
167 	/// Otherwise, pointer to node describing the next dimension.
168 	immutable(Node)* root = emptySetRoot;
169 
170 	/// A set containing a single nil-dimensional element.
171 	/// Holds exactly one value (a point with all dimensions at
172 	/// nullValue).
173 	static immutable unitSet = MapSet(null);
174 
175 	/// The empty set. Represents a set which holds zero values.
176 	static immutable emptySet = MapSet(emptySetRoot);
177 
178 	/// Return the total number of items in this set.
179 	size_t count() const
180 	{
181 		if (this is emptySet)
182 			return 0;
183 		if (this is unitSet)
184 			return 1;
185 		return root.totalMembers;
186 	}
187 
188 	private struct SetSetOp { MapSet a, b; }
189 	private struct SetDimOp { MapSet set; DimName dim; }
190 	private struct SetIdxOp { MapSet set; size_t index; }
191 	private struct Cache
192 	{
193 		/// For deduplication - key is value
194 		MapSet[MapSet] instances;
195 		/// Operations - things that operate recursively on subtrees
196 		/// should be memoized here
197 		MapSet[SetSetOp] merge, subtract, cartesianProduct, reorderUsing;
198 		MapSet[SetDimOp] remove, bringToFront;
199 		MapSet[SetIdxOp] swapDepth;
200 		MapSet[MapSet] optimize, completeSuperset;
201 		size_t[MapSet] uniqueNodes, maxDepth;
202 	}
203 
204 	/// Because subtrees can be reused within the tree, a way of
205 	/// memoizing operations across the entire tree (instead of just
206 	/// across children of a single node, or siblings) is crucial for
207 	/// performance.
208 	static Cache cache;
209 
210 	/// Clear the global operations cache.
211 	static void clearCache()
212 	{
213 		cache = Cache.init;
214 	}
215 
216 	/// Because MapSet operations benefit greatly from memoization,
217 	/// maintaining a set of interned sets benefits performance. After
218 	/// calling `clearCache`, call this method to re-register extant
219 	/// live instances to the instance cache.
220 	void addToCache() const
221 	{
222 		cache.instances.updateVoid(
223 			this,
224 			{
225 				if (this !is emptySet && this !is unitSet)
226 					foreach (ref pair; root.children)
227 						pair.set.addToCache();
228 				return this;
229 			},
230 			(ref MapSet set)
231 			{
232 				assert(set is this);
233 			}
234 		);
235 	}
236 
237 	/// Intern and deduplicate this MapSet.
238 	/// Needs to be called only after constructing a MapSet manually
239 	/// (by allocating, populating, and setting `root`).
240 	MapSet deduplicate() const
241 	{
242 		MapSet deduplicated;
243 		cache.instances.updateVoid(
244 			this,
245 			{
246 				debug if (this !is emptySet && this !is unitSet)
247 					foreach (ref pair; root.children)
248 						pair.set.assertDeduplicated();
249 				deduplicated = this;
250 				return this;
251 			},
252 			(ref MapSet set)
253 			{
254 				deduplicated = set;
255 			}
256 		);
257 		return deduplicated;
258 	}
259 
260 	private void assertDeduplicated() const { debug assert(this is emptySet || this is unitSet || cache.instances[this] is this); }
261 
262 	/// Count and return the total number of unique nodes in this MapSet.
263 	size_t uniqueNodes() const
264 	{
265 		return cache.uniqueNodes.require(this, {
266 			HashSet!MapSet seen;
267 			void visit(MapSet set)
268 			{
269 				if (set is emptySet || set is unitSet || set in seen)
270 					return;
271 				seen.add(set);
272 				foreach (ref pair; set.root.children)
273 					visit(pair.set);
274 			}
275 			visit(this);
276 			return seen.length;
277 		}());
278 	}
279 
280 	/// Collect the names of all dimensions occurring in this tree.
281 	DimName[] getDims() const
282 	{
283 		HashSet!DimName dims;
284 		HashSet!MapSet seen;
285 		void visit(MapSet set)
286 		{
287 			if (set is emptySet || set is unitSet || set in seen)
288 				return;
289 			seen.add(set);
290 			dims.add(set.root.dim);
291 			foreach (ref pair; set.root.children)
292 				visit(pair.set);
293 		}
294 		visit(this);
295 		return dims.keys;
296 	}
297 
298 	/// Return all values for all dimensions occurring in this set.
299 	/// The Cartesian product of these values would thus be a superset
300 	/// of this set.
301 	/// This is equivalent to, but faster than, calling `getDims` and
302 	/// then `all` for each dim.
303 	HashSet!(immutable(DimValue))[DimName] getDimsAndValues(this This)()
304 	{
305 		if (this is emptySet) return null;
306 
307 		// Be careful to count the implicit nullValues on branches
308 		// where a dim doesn't occur.
309 		MapSet set = this.completeSuperset;
310 		HashSet!(immutable(DimValue))[DimName] result;
311 		while (set !is unitSet)
312 		{
313 			DimName dim = set.root.dim;
314 			HashSet!(immutable(DimValue)) values = set.root.children.map!((ref child) => child.value).toSet;
315 			bool added = result.addNew(dim, values);
316 			assert(added, "Duplicate dimension");
317 			set = set.root.children[0].set;
318 		}
319 		return result;
320 	}
321 
322 	/// Combine two matrices together, returning their union.
323 	/// If `other` is a subset of `this`, return `this` unchanged.
324 	MapSet merge(MapSet other) const
325 	{
326 		if (this is emptySet) return other;
327 		if (other is emptySet) return this;
328 		if (this is other) return this;
329 		if (!root) return bringToFront(other.root.dim).merge(other);
330 
331 		this.assertDeduplicated();
332 		other.assertDeduplicated();
333 		return cache.merge.require(SetSetOp(this, other), {
334 			other = other.bringToFront(root.dim);
335 
336 			if (root.children.length == 1 && other.root.children.length == 1 &&
337 				root.children[0].value == other.root.children[0].value)
338 			{
339 				// Single-child optimization
340 				auto mergeResult = root.children[0].set.merge(other.root.children[0].set);
341 				if (mergeResult !is root.children[0].set)
342 					return MapSet(new immutable Node(root.dim, [Pair(root.children[0].value, mergeResult)])).deduplicate;
343 				else
344 					return this;
345 			}
346 
347 			MapSet[DimValue] newChildren;
348 			foreach (ref pair; root.children)
349 				newChildren[pair.value] = pair.set;
350 
351 			bool modified;
352 			foreach (ref pair; other.root.children)
353 				newChildren.updateVoid(pair.value,
354 					{
355 						modified = true;
356 						return pair.set;
357 					},
358 					(ref MapSet oldSubmatrix)
359 					{
360 						auto mergeResult = oldSubmatrix.merge(pair.set);
361 						if (oldSubmatrix !is mergeResult)
362 						{
363 							oldSubmatrix = mergeResult;
364 							modified = true;
365 						}
366 					}
367 				);
368 
369 			if (!modified)
370 				return this;
371 
372 			return MapSet(new immutable Node(root.dim, cast(immutable) newChildren)).deduplicate;
373 		}());
374 	}
375 
376 	/// Return the difference between `this` and the given set.
377 	/// If `other` does not intersect with `this`, return `this` unchanged.
378 	MapSet subtract(MapSet other) const
379 	{
380 		if (this is emptySet) return this;
381 		if (other is emptySet) return this;
382 		if (this is other) return emptySet;
383 		if (!root) return bringToFront(other.root.dim).subtract(other);
384 
385 		this.assertDeduplicated();
386 		other.assertDeduplicated();
387 		return cache.subtract.require(SetSetOp(this, other), {
388 			other = other.bringToFront(root.dim);
389 
390 			if (root.children.length == 1 && other.root.children.length == 1 &&
391 				root.children[0].value == other.root.children[0].value)
392 			{
393 				// Single-child optimization
394 				auto subtractResult = root.children[0].set.subtract(other.root.children[0].set);
395 				if (subtractResult is emptySet)
396 					return emptySet;
397 				else
398 				if (subtractResult !is root.children[0].set)
399 					return MapSet(new immutable Node(root.dim, [Pair(root.children[0].value, subtractResult)])).deduplicate;
400 				else
401 					return this;
402 			}
403 
404 			MapSet[DimValue] newChildren;
405 			foreach (ref pair; root.children)
406 				newChildren[pair.value] = pair.set;
407 
408 			bool modified;
409 			foreach (ref pair; other.root.children)
410 				if (auto poldSubmatrix = pair.value in newChildren)
411 				{
412 					auto subtractResult = poldSubmatrix.subtract(pair.set);
413 					if (*poldSubmatrix !is subtractResult)
414 					{
415 						*poldSubmatrix = subtractResult;
416 						if (subtractResult is emptySet)
417 							newChildren.remove(pair.value);
418 						modified = true;
419 					}
420 				}
421 
422 			if (!modified)
423 				return this;
424 			if (!newChildren.length)
425 				return emptySet;
426 
427 			return MapSet(new immutable Node(root.dim, cast(immutable) newChildren)).deduplicate;
428 		}());
429 	}
430 
431 	/// Apply `fn` over subsets, and return a new `MapSet`.
432 	/*private*/ MapSet lazyMap(scope MapSet delegate(MapSet) fn) const
433 	{
434 		// Defer allocation until the need to mutate
435 		foreach (i, ref pair; root.children) // Read-only scan
436 		{
437 			auto newSet = fn(pair.set);
438 			if (newSet !is pair.set)
439 			{
440 				auto newChildren = new Pair[root.children.length];
441 				// Known to not need mutation
442 				newChildren[0 .. i] = root.children[0 .. i];
443 				// Reuse already calculated result
444 				newChildren[i] = Pair(pair.value, newSet);
445 				// Continue scan with mutation
446 				foreach (j, ref pair2; root.children[i + 1 .. $])
447 					newChildren[i + 1 + j] = Pair(pair2.value, fn(pair2.set));
448 				return MapSet(new immutable Node(root.dim, cast(immutable) newChildren)).deduplicate;
449 			}
450 		}
451 		return this; // No mutation necessary
452 	}
453 
454 	/// Like `lazyMap`, but allow `fn` to return `emptySet`.
455 	/*private*/ MapSet lazyDeletingMap(scope MapSet delegate(MapSet) fn) const
456 	{
457 		// Defer allocation until the need to mutate
458 		foreach (i, ref pair; root.children) // Read-only scan
459 		{
460 			auto newSet = fn(pair.set);
461 			if (newSet !is pair.set)
462 			{
463 				auto newChildren = new Pair[root.children.length];
464 				// Known to not need mutation
465 				newChildren[0 .. i] = root.children[0 .. i];
466 				// Reuse already calculated result
467 				auto j = i;
468 				if (newSet != emptySet)
469 					newChildren[j++] = Pair(pair.value, newSet);
470 				// Continue scan with mutation
471 				foreach (ref pair2; root.children[i + 1 .. $])
472 				{
473 					newSet = fn(pair2.set);
474 					if (newSet != emptySet)
475 						newChildren[j++] = Pair(pair2.value, newSet);
476 				}
477 				newChildren = newChildren[0 .. j];
478 				if (newChildren.length)
479 					return MapSet(new immutable Node(root.dim, cast(immutable) newChildren)).deduplicate;
480 				else
481 					return emptySet;
482 			}
483 		}
484 		return this; // No mutation necessary
485 	}
486 
487 	/// "Unset" a given dimension, removing it from the matrix.
488 	/// The result is the union of all sub-matrices for all values of `dim`.
489 	MapSet remove(DimName dim) const
490 	{
491 		if (this is emptySet || this is unitSet) return this;
492 		this.assertDeduplicated();
493 		return cache.remove.require(SetDimOp(this, dim), {
494 			if (root.dim == dim)
495 			{
496 				MapSet result;
497 				foreach (ref pair; root.children)
498 					result = result.merge(pair.set);
499 				return result;
500 			}
501 			return lazyMap(set => set.remove(dim));
502 		}());
503 	}
504 
505 	/// Unset dimensions according to a predicate.
506 	/// This is faster than removing dimensions individually, however,
507 	/// unlike the `DimName` overload, this one does not benefit from global memoization.
508 	MapSet remove(bool delegate(DimName) pred) const
509 	{
510 		if (this is emptySet) return this;
511 		this.assertDeduplicated();
512 
513 		MapSet[MapSet] cache;
514 		MapSet visit(MapSet set)
515 		{
516 			if (set is unitSet)
517 				return set;
518 			return cache.require(set, {
519 				if (pred(set.root.dim))
520 				{
521 					MapSet result;
522 					foreach (ref pair; set.root.children)
523 						result = result.merge(visit(pair.set));
524 					return result;
525 				}
526 				return set.lazyMap(set => visit(set));
527 			}());
528 		}
529 		return visit(this);
530 	}
531 
532 	/// Set the given dimension to always have the given value,
533 	/// collapsing (returning the union of) all sub-matrices
534 	/// for previously different values of `dim`.
535 	MapSet set(DimName dim, DimValue value) const
536 	{
537 		return this.remove(dim).addDim(dim, value);
538 	}
539 
540 	/// Adds a new dimension with the given value.
541 	/// The dimension must not have previously existed in `this`.
542 	private MapSet addDim(DimName dim, DimValue value) const
543 	{
544 		if (this is emptySet) return this;
545 		this.assertDeduplicated();
546 		assert(this is this.remove(dim), "Duplicate dimension");
547 		if (value == nullValue)
548 			return this;
549 		return MapSet(new immutable Node(dim, [Pair(value, this)])).deduplicate;
550 	}
551 
552 	/// Return a sub-matrix for all points where the given dimension has this value.
553 	/// The dimension itself is not included in the result.
554 	/// Note: if `dim` doesn't occur (e.g. because `this` is the unit set) and `value` is `nullValue`,
555 	/// this returns the unit set (not the empty set).
556 	MapSet slice(DimName dim, DimValue value) const
557 	{
558 		if (this is emptySet) return emptySet;
559 		this.assertDeduplicated();
560 		foreach (ref pair; bringToFront(dim).root.children)
561 			if (pair.value == value)
562 				return pair.set;
563 		return emptySet;
564 	}
565 
566 	/// Return a subset of this set for all points where the given dimension has this value.
567 	/// Unlike `slice`, the dimension itself is included in the result (with the given value).
568 	MapSet get(DimName dim, DimValue value) const
569 	{
570 		if (this is emptySet) return emptySet;
571 		this.assertDeduplicated();
572 		enum reorder = false;
573 		static if (reorder)
574 		{
575 			foreach (ref pair; bringToFront(dim).root.children)
576 				if (pair.value == value)
577 					return MapSet(new immutable Node(dim, [Pair(value, pair.set)])).deduplicate;
578 			return emptySet;
579 		}
580 		else
581 		{
582 			MapSet[MapSet] cache;
583 			MapSet visit(MapSet set)
584 			{
585 				return cache.require(set, {
586 					if (set == unitSet)
587 					{
588 						if (value == nullValue)
589 							return set;
590 						else
591 							return emptySet;
592 					}
593 					if (set.root.dim == dim)
594 					{
595 						foreach (ref pair; set.root.children)
596 							if (pair.value == value)
597 								return MapSet(new immutable Node(dim, [Pair(value, pair.set)])).deduplicate;
598 						return emptySet;
599 					}
600 					return set.lazyDeletingMap(&visit);
601 				}());
602 			}
603 			return visit(this);
604 		}
605 	}
606 
607 	/// Return all unique values occurring for a given dimension.
608 	/// Unless this is the empty set, the return value is always non-empty.
609 	/// If `dim` doesn't occur, it will be `[nullValue]`.
610 	const(DimValue)[] all(DimName dim) const
611 	{
612 		// return bringToFront(dim).root.children.keys;
613 		if (this is emptySet) return null;
614 		if (this is unitSet) return [nullValue];
615 		if (root.dim == dim) return root.children.amap!(child => child.value);
616 		this.assertDeduplicated();
617 
618 		HashSet!DimValue allValues;
619 		HashSet!MapSet seen;
620 		void visit(MapSet set)
621 		{
622 			if (set is unitSet)
623 			{
624 				allValues.add(nullValue);
625 				return;
626 			}
627 			if (set in seen)
628 				return;
629 			seen.add(set);
630 
631 			if (set.root.dim == dim)
632 				foreach (ref pair; set.root.children)
633 					allValues.add(pair.value);
634 			else
635 				foreach (ref pair; set.root.children)
636 					visit(pair.set);
637 		}
638 		visit(this);
639 		return allValues.keys;
640 	}
641 
642 	/// Return a set which represents the Cartesian product between
643 	/// this set and the given `values` across the specified
644 	/// dimension.
645 	MapSet cartesianProduct(DimName dim, DimValue[] values) const
646 	{
647 		if (this is emptySet) return emptySet;
648 		if (values.length == 0) return emptySet;
649 		this.assertDeduplicated();
650 		auto unset = this.remove(dim);
651 		auto children = values.map!(value => Pair(value, unset)).array;
652 		children.sort();
653 		return MapSet(new immutable Node(dim, cast(immutable) children)).deduplicate;
654 	}
655 
656 	/// Assumes that `dim` does not occur in `this` set.
657 	private MapSet uncheckedCartesianProduct(DimName dim, DimValue[] values) const
658 	{
659 		if (this is emptySet) return emptySet;
660 		if (values.length == 0) return emptySet;
661 		this.assertDeduplicated();
662 		auto children = values.map!(value => Pair(value, this)).array;
663 		children.sort();
664 		return MapSet(new immutable Node(dim, cast(immutable) children)).deduplicate;
665 	}
666 
667 	/// Return a set which represents the Cartesian product between
668 	/// this and the given set.
669 	/// Duplicate dimensions are first removed from `this` set.
670 	/// For best performance, call `big.cartesianProduct(small)`
671 	MapSet cartesianProduct(MapSet other) const
672 	{
673 		if (this is emptySet || other is emptySet) return emptySet;
674 		if (this is unitSet) return other;
675 		if (other is unitSet) return this;
676 
677 		MapSet unset = this;
678 		foreach (dim; other.getDims())
679 			unset = unset.remove(dim);
680 
681 		return other.uncheckedCartesianProduct(unset);
682 	}
683 
684 	// Assumes that the dimensions in `this` and `other` are disjoint.
685 	private MapSet uncheckedCartesianProduct(MapSet other) const
686 	{
687 		if (this is unitSet) return other;
688 		if (other is unitSet) return this;
689 
690 		this.assertDeduplicated();
691 		other.assertDeduplicated();
692 
693 		return cache.cartesianProduct.require(SetSetOp(this, other), {
694 			return lazyMap(set => set.uncheckedCartesianProduct(other));
695 		}());
696 	}
697 
698 	/// Return a superset of this set, consisting of a Cartesian
699 	/// product of every value of every dimension. The total number of
700 	/// unique nodes is thus equal to the number of dimensions.
701 	MapSet completeSuperset() const
702 	{
703 		if (this is emptySet || this is unitSet) return this;
704 		this.assertDeduplicated();
705 
706 		return cache.completeSuperset.require(this, {
707 			MapSet child = MapSet.emptySet;
708 			foreach (ref pair; root.children)
709 				child = child.merge(pair.set);
710 			child = child.completeSuperset();
711 			auto newChildren = root.children.dup;
712 			foreach (ref pair; newChildren)
713 				pair.set = child;
714 			return MapSet(new immutable Node(root.dim, cast(immutable) newChildren)).deduplicate;
715 		}());
716 	}
717 
718 	/// Refactor this matrix into one with the same data,
719 	/// but putting the given dimension in front.
720 	/// This will speed up access to values with the given dimension.
721 	/// If the dimension does not yet occur in the set (or any subset),
722 	/// it is instantiated with a single `nullValue` value.
723 	/// The set must be non-empty.
724 	MapSet bringToFront(DimName dim) const
725 	{
726 		assert(this !is emptySet, "Empty sets may not have dimensions");
727 
728 		if (this is unitSet)
729 		{
730 			// We reached the bottom, and did not find `dim` along the way.
731 			// Create it now.
732 			return MapSet(new immutable Node(dim, [Pair(nullValue, unitSet)])).deduplicate;
733 		}
734 
735 		if (dim == root.dim)
736 		{
737 			// Already at the front.
738 			return this;
739 		}
740 
741 		// 1. Recurse.
742 		// 2. After recursion, all children should have `dim` at the front.
743 		// So, just swap this layer with the next one.
744 
745 		this.assertDeduplicated();
746 		return cache.bringToFront.require(SetDimOp(this, dim), {
747 			if (root.children.length == 1)
748 			{
749 				auto newSubmatrix = root.children[0].set.bringToFront(dim);
750 				auto newChildren = newSubmatrix.root.children.dup;
751 				foreach (ref pair; newChildren)
752 					pair.set = MapSet(new immutable Node(root.dim, [Pair(root.children[0].value, pair.set)])).deduplicate;
753 				return MapSet(new immutable Node(dim, cast(immutable) newChildren)).deduplicate;
754 			}
755 
756 			Pair[][DimValue] subsets;
757 			foreach (ref pair; root.children)
758 			{
759 				auto newSubmatrix = pair.set.bringToFront(dim);
760 				assert(newSubmatrix.root.dim == dim);
761 				foreach (ref pair2; newSubmatrix.root.children)
762 					subsets[pair2.value] ~= Pair(pair.value, pair2.set);
763 			}
764 			Pair[] newChildren;
765 			foreach (value, children; subsets)
766 			{
767 				children.sort();
768 				newChildren ~= Pair(value, MapSet(new immutable Node(root.dim, cast(immutable) children)).deduplicate);
769 			}
770 			newChildren.sort();
771 			return MapSet(new immutable Node(dim, cast(immutable) newChildren)).deduplicate;
772 		}());
773 	}
774 
775 	/// Refactor this matrix into one with the same data,
776 	/// but attempting to lower the total number of nodes.
777 	MapSet optimize() const
778 	{
779 		if (this is emptySet || this is unitSet) return this;
780 		this.assertDeduplicated();
781 
782 		return cache.optimize.require(this, {
783 			struct Result { bool done; MapSet[MapSet] map; }
784 			static Result optimizeLayer(HashSet!MapSet sets0)
785 			{
786 				// - At the bottom?
787 				//   - Yes:
788 				//     - return failure
789 				//   - No:
790 				//     - Try to swap this layer. Success?
791 				//       - Yes:
792 				//         - return success
793 				//       - No:
794 				//         - Recurse and try to swap next layer. Success?
795 				//           - Yes: Retry this layer
796 				//           - No: return failure (bottom reached)
797 
798 				assert(!sets0.empty);
799 				if (sets0.byKey.front is unitSet)
800 					return Result(true, null); // at the bottom
801 				auto dim0 = sets0.byKey.front.root.dim;
802 				assert(sets0.byKey.all!(set => set !is unitSet), "Leaf/non-leaf nodes mismatch");
803 				assert(sets0.byKey.all!(set => set.root.dim == dim0), "Dim mismatch");
804 
805 				auto sets1 = sets0.byKey.map!(set => set.root.children.map!(function MapSet (ref child) => child.set)).joiner.toSet;
806 				assert(!sets1.empty);
807 				if (sets1.byKey.front is unitSet)
808 					return Result(true, null); // one layer away from the bottom, nothing to swap with
809 				auto dim1 = sets1.byKey.front.root.dim;
810 				assert(sets1.byKey.all!(set => set !is unitSet), "Leaf/non-leaf nodes mismatch");
811 				assert(sets1.byKey.all!(set => set.root.dim == dim1), "Dim mismatch");
812 
813 				auto currentNodes = sets0.length + sets1.length;
814 
815 				MapSet[MapSet] swappedSets;
816 				HashSet!MapSet sets0new, sets1new;
817 
818 				foreach (set0; sets0.byKey)
819 				{
820 					Pair[][DimValue] subsets;
821 					foreach (ref pair0; set0.root.children)
822 						foreach (ref pair1; pair0.set.root.children)
823 							subsets[pair1.value] ~= Pair(pair0.value, pair1.set);
824 
825 					Pair[] newChildren;
826 					foreach (value, children; subsets)
827 					{
828 						children.sort();
829 						auto set1new = MapSet(new immutable Node(dim0, cast(immutable) children)).deduplicate;
830 						sets1new.add(set1new);
831 						newChildren ~= Pair(value, set1new);
832 					}
833 					newChildren.sort();
834 					auto set0new = MapSet(new immutable Node(dim1, cast(immutable) newChildren)).deduplicate;
835 					sets0new.add(set0new);
836 					swappedSets[set0] = set0new;
837 				}
838 
839 				auto newNodes = sets0new.length + sets1new.length;
840 
841 				if (newNodes < currentNodes)
842 					return Result(false, swappedSets); // Success, retry above layer
843 
844 				// Failure, descend
845 
846 				auto result1 = optimizeLayer(sets1);
847 				if (!result1.map)
848 				{
849 					assert(result1.done);
850 					return Result(true, null); // Done, bottom reached
851 				}
852 
853 				// Apply result
854 				sets0new.clear();
855 				foreach (set0; sets0.byKey)
856 				{
857 					set0.assertDeduplicated();
858 					auto newChildren = set0.root.children.dup;
859 					foreach (ref pair0; newChildren)
860 						pair0.set = result1.map[pair0.set];
861 					auto set0new = MapSet(new immutable Node(dim0, cast(immutable) newChildren)).deduplicate;
862 					sets0new.add(set0new);
863 					swappedSets[set0] = set0new;
864 				}
865 
866 				if (result1.done)
867 					return Result(true, swappedSets);
868 
869 				// Retry this layer
870 				auto result0 = optimizeLayer(sets0new);
871 				if (!result0.map)
872 				{
873 					assert(result0.done);
874 					return Result(true, swappedSets); // Bottom was reached upon retry, just return our results unchanged
875 				}
876 
877 				MapSet[MapSet] compoundedResult;
878 				foreach (set0; sets0.byKey)
879 					compoundedResult[set0] = result0.map[swappedSets[set0]];
880 				return Result(result0.done, compoundedResult);
881 			}
882 
883 			MapSet[1] root = [this.normalize];
884 			while (true)
885 			{
886 				auto result = optimizeLayer(root[].toSet());
887 				if (result.map)
888 					root[0] = result.map[root[0]];
889 				if (result.done)
890 					return root[0];
891 			}
892 		}());
893 	}
894 
895 	/// Return a set equivalent to a unit set (all dimensions
896 	/// explicitly set to `nullValue`), with all dimensions in `this`,
897 	/// in an order approximately following that of the dimensions in
898 	/// `this`.  Implicitly normalized.
899 	private MapSet dimOrderReference() const
900 	{
901 		if (this is unitSet || this is emptySet) return this;
902 
903 		DimName[] dims;
904 		HashSet!MapSet seen;
905 		void visit(MapSet set, size_t pos)
906 		{
907 			if (set is unitSet) return;
908 			if (set in seen) return;
909 			seen.add(set);
910 			if ((pos == dims.length || dims[pos] != set.root.dim) && !dims.canFind(set.root.dim))
911 			{
912 				dims.insertInPlace(pos, set.root.dim);
913 				pos++;
914 			}
915 			foreach (ref pair; set.root.children)
916 				visit(pair.set, pos);
917 		}
918 		visit(this, 0);
919 
920 		MapSet result = unitSet;
921 		foreach_reverse (dim; dims)
922 			result = MapSet(new immutable Node(dim, [Pair(nullValue, result)])).deduplicate();
923 
924 		return result;
925 	}
926 
927 	/// Refactor this matrix into one in which dimensions always occur
928 	/// in the same order, no matter what path is taken.
929 	MapSet normalize() const
930 	{
931 		if (this is unitSet || this is emptySet) return this;
932 
933 		return reorderUsing(dimOrderReference());
934 	}
935 
936 	private size_t maxDepth() const
937 	{
938 		import std.algorithm.comparison : max;
939 
940 		if (this is emptySet || this is unitSet)
941 			return 0;
942 		this.assertDeduplicated();
943 		return cache.maxDepth.require(this, {
944 			size_t maxDepth = 0;
945 			foreach (ref pair; root.children)
946 				maxDepth = max(maxDepth, pair.set.maxDepth());
947 			return 1 + maxDepth;
948 		}());
949 	}
950 
951 	/// Swap the two adjacent dimensions which are `depth` dimensions away.
952 	/*private*/ MapSet swapDepth(size_t depth) const
953 	{
954 		if (this is emptySet || this is unitSet) return this;
955 		this.assertDeduplicated();
956 
957 		return cache.swapDepth.require(SetIdxOp(this, depth), {
958 			if (depth == 0)
959 			{
960 				foreach (ref pair; root.children)
961 					if (pair.set !is unitSet)
962 						return bringToFront(pair.set.root.dim);
963 				return this;
964 			}
965 			else
966 			{
967 				auto newChildren = root.children.dup;
968 				foreach (ref pair; newChildren)
969 					pair.set = pair.set.swapDepth(depth - 1);
970 				return MapSet(new immutable Node(root.dim, cast(immutable) newChildren)).deduplicate;
971 			}
972 		}());
973 	}
974 
975 	/// Refactor this matrix into one with the same data, but in which
976 	/// the dimensions always occur as in `reference` (which is
977 	/// assumed to be normalized).
978 	MapSet reorderUsing(MapSet reference) const
979 	{
980 		if (this is emptySet || reference is emptySet || reference is unitSet) return this;
981 		this.assertDeduplicated();
982 		reference.assertDeduplicated();
983 
984 		return cache.reorderUsing.require(SetSetOp(this, reference), {
985 			return bringToFront(reference.root.dim).lazyMap(set => set.reorderUsing(reference.root.children[0].set));
986 		}());
987 	}
988 
989 	void toString(scope void delegate(const(char)[]) sink) const
990 	{
991 		import std.format : formattedWrite;
992 		if (this is emptySet)
993 			sink("{}");
994 		else
995 		if (this is unitSet)
996 			sink("{[]}");
997 		else
998 			sink.formattedWrite!"%s"(*root);
999 	} ///
1000 
1001 	hash_t toHash() const
1002 	{
1003 		return
1004 			this is emptySet ? 0 :
1005 			this is unitSet ? 1 :
1006 			root.toHash();
1007 	} ///
1008 
1009 	bool opEquals(const typeof(this) s) const
1010 	{
1011 		if (root is s.root)
1012 			return true;
1013 		if (this is emptySet || this is unitSet || s is emptySet || s is unitSet)
1014 			return this is s;
1015 		return *root == *s.root;
1016 	} ///
1017 }
1018 
1019 unittest
1020 {
1021 	import std.algorithm.sorting : sort;
1022 
1023 	alias M = MapSet!(string, int);
1024 	M m, n;
1025 	m = m.merge(M.unitSet.set("x", 1).set("y", 5));
1026 	m = m.merge(M.unitSet.set("x", 1).set("y", 6));
1027 	assert(m.all("x") == [1]);
1028 	assert(m.all("y").dup.sort.release == [5, 6]);
1029 
1030 	m = m.merge(M.unitSet.set("x", 2).set("y", 6));
1031 	assert(m.get("x", 1).all("y").dup.sort.release == [5, 6]);
1032 	assert(m.get("y", 6).all("x").dup.sort.release == [1, 2]);
1033 
1034 	m = m.subtract(M.unitSet.set("x", 1).set("y", 6));
1035 	assert(m.all("x").dup.sort.release == [1, 2]);
1036 	assert(m.all("y").dup.sort.release == [5, 6]);
1037 	assert(m.get("x", 1).all("y") == [5]);
1038 	assert(m.get("y", 6).all("x") == [2]);
1039 
1040 	m = M.emptySet;
1041 	m = m.merge(M.unitSet.set("x", 10).set("y", 20));
1042 	m = m.merge(M.unitSet.set("x", 10).set("y", 21));
1043 	m = m.merge(M.unitSet.set("x", 11).set("y", 21));
1044 	m = m.merge(M.unitSet.set("x", 11).set("y", 22));
1045 	auto o = m.optimize();
1046 	assert(o.uniqueNodes < m.uniqueNodes);
1047 
1048 	m = M.emptySet;
1049 	assert(m.all("x") == []);
1050 	m = M.unitSet;
1051 	assert(m.all("x") == [0]);
1052 	m = m.merge(M.unitSet.set("x", 1));
1053 	assert(m.all("x").dup.sort.release == [0, 1]);
1054 
1055 	m = M.unitSet;
1056 	assert(m.set("x", 1).set("x", 1).all("x") == [1]);
1057 
1058 	m = M.unitSet;
1059 	m = m.cartesianProduct("x", [1, 2, 3]);
1060 	m = m.cartesianProduct("y", [1, 2, 3]);
1061 	m = m.cartesianProduct("z", [1, 2, 3]);
1062 	assert(m.count == 3 * 3 * 3);
1063 	assert(m            .all("x").dup.sort.release == [1, 2, 3]);
1064 	assert(m.set("z", 1).all("x").dup.sort.release == [1, 2, 3]);
1065 	assert(m.set("x", 1).all("z").dup.sort.release == [1, 2, 3]);
1066 
1067 	m = M.unitSet;
1068 	m = m.cartesianProduct("a", [1, 2, 3]);
1069 	m = m.cartesianProduct("b", [1, 2, 3]);
1070 	n = M.unitSet;
1071 	n = n.cartesianProduct("c", [1, 2, 3]);
1072 	n = n.cartesianProduct("d", [1, 2, 3]);
1073 	m = m.cartesianProduct(n);
1074 	assert(m.count == 3 * 3 * 3 * 3);
1075 
1076 	assert(M.unitSet != M.unitSet.set("x", 1));
1077 
1078 	m = M.emptySet;
1079 	m = m.merge(M.unitSet.set("x", 1).set("y", 11));
1080 	m = m.merge(M.unitSet.set("x", 2).set("y", 12).set("z", 22));
1081 	m = m.completeSuperset();
1082 	assert(m.uniqueNodes == 3);
1083 	assert(m.count == 8); // 2 ^^ 3
1084 	assert(m.all("x").dup.sort.release == [1, 2]);
1085 	assert(m.all("y").dup.sort.release == [11, 12]);
1086 	assert(m.all("z").dup.sort.release == [0, 22]);
1087 }
1088 
1089 
1090 /// Allows executing a deterministic algorithm over all states in a given MapSet.
1091 /// If a variable is not queried by the algorithm, states for all
1092 /// variations of that variable are processed in one iteration.
1093 struct MapSetVisitor(A, V, V nullValue = V.init)
1094 {
1095 	/// Underlying `MapSet`.
1096 	alias Set = MapSet!(A, V, nullValue);
1097 	Set set; /// ditto
1098 
1099 	/// Internal state.
1100 	/*private*/ public
1101 	{
1102 		// Iteration state for resolved values
1103 		struct Var
1104 		{
1105 			A name; ///
1106 			const(V)[] values; ///
1107 			size_t pos; ///
1108 		}
1109 		Var[] stack;
1110 		size_t stackPos;
1111 
1112 		private enum Maybe { no, maybe, yes }
1113 
1114 		struct VarState
1115 		{
1116 			// Variable holds this one value if `haveValue` is true.
1117 			V value = nullValue;
1118 
1119 			// Optimization.
1120 			// For variables which have been resolved, or have been
1121 			// set to some specific value, remember that value here
1122 			// (in the `value` field).
1123 			// Faster than checking workingSet.all(name)[0].
1124 			// If this is set, it has a concrete value because
1125 			// - we are iterating over it (it's in the stack)
1126 			// - due to a `put` call
1127 			// - it was never in the set (so it's implicitly at nullValue).
1128 			bool haveValue = true;
1129 
1130 			// Optimization.
1131 			// Accumulate MapSet.set calls, and flush then in bulk.
1132 			// The value to flush is stored in the `value` field.
1133 			// If `dirty` is true, `haveValue` must be true.
1134 			bool dirty;
1135 
1136 			// Optimization.
1137 			// Remember whether this variable is in the set or not.
1138 			Maybe inSet = Maybe.no;
1139 
1140 			void setInSet(bool inSet)
1141 			{
1142 				if (inSet)
1143 					this.inSet = Maybe.yes;
1144 				else
1145 				{
1146 					this.inSet = Maybe.no;
1147 					assert(!this.dirty, "TODO");
1148 					if (this.haveValue)
1149 						assert(this.value == nullValue);
1150 					else
1151 					{
1152 						this.haveValue = true;
1153 						this.value = nullValue;
1154 					}
1155 				}
1156 			}
1157 		}
1158 		VarState[A] varState, initialVarState;
1159 
1160 		// The version of the set for the current iteration.
1161 		private Set workingSet;
1162 	}
1163 
1164 	this(Set set)
1165 	{
1166 		this.set = set;
1167 		foreach (dim, values; set.getDimsAndValues())
1168 		{
1169 			auto pstate = &initialVarState.require(dim);
1170 			pstate.inSet = Maybe.yes;
1171 			if (values.length == 1)
1172 				pstate.value = values.byKey.front;
1173 			else
1174 				pstate.haveValue = false;
1175 		}
1176 	} ///
1177 
1178 	@disable this(this);
1179 
1180 	MapSetVisitor dup()
1181 	{
1182 		MapSetVisitor r;
1183 		r.set = set;
1184 		r.stack = stack.dup;
1185 		r.stackPos = stackPos;
1186 		r.varState = varState.dup;
1187 		r.workingSet = workingSet;
1188 		return r;
1189 	}
1190 
1191 	/// Resets iteration to the beginning.
1192 	/// Equivalent to but faster than constructing a new MapSetVisitor
1193 	/// instance (`visitor = MapSetVisitor(visitor.set)`).
1194 	void reset()
1195 	{
1196 		workingSet = Set.emptySet;
1197 		stack = null;
1198 	}
1199 
1200 	/// Returns true if there are more states to iterate over,
1201 	/// otherwise returns false
1202 	bool next()
1203 	{
1204 		if (set is Set.emptySet)
1205 			return false;
1206 		if (workingSet is Set.emptySet)
1207 		{
1208 			// first iteration
1209 		}
1210 		else
1211 			while (true)
1212 			{
1213 				if (!stack.length)
1214 					return false; // All possibilities exhausted
1215 				auto last = &stack[$-1];
1216 				last.pos++;
1217 				if (last.pos == last.values.length)
1218 				{
1219 					stack = stack[0 .. $ - 1];
1220 					continue;
1221 				}
1222 				break;
1223 			}
1224 
1225 		workingSet = set;
1226 		varState = initialVarState.dup;
1227 		stackPos = 0;
1228 		return true;
1229 	}
1230 
1231 	private void flush()
1232 	{
1233 		auto toRemove = varState.byKeyValue
1234 			.filter!(pair => pair.value.dirty && pair.value.inSet >= Maybe.maybe)
1235 			.map!(pair => pair.key)
1236 			.toSet;
1237 		workingSet = workingSet.remove((A name) => name in toRemove);
1238 		foreach (name, ref state; varState)
1239 			if (state.dirty)
1240 			{
1241 				workingSet = workingSet.addDim(name, state.value);
1242 				state.dirty = false;
1243 				state.setInSet(state.value != nullValue); // addDim is a no-op with nullValue
1244 			}
1245 	}
1246 
1247 	private void flush(A name)
1248 	{
1249 		if (auto pstate = name in varState)
1250 			if (pstate.dirty)
1251 			{
1252 				pstate.dirty = false;
1253 				if (pstate.inSet >= Maybe.maybe)
1254 				{
1255 					if (pstate.inSet == Maybe.yes)
1256 					{
1257 						auto oldSet = workingSet;
1258 						auto newSet = workingSet.remove(name);
1259 						assert(oldSet != newSet, "Actually wasn't in the set");
1260 						workingSet = newSet;
1261 					}
1262 					else
1263 						workingSet = workingSet.remove(name);
1264 				}
1265 				workingSet = workingSet.addDim(name, pstate.value);
1266 				pstate.setInSet(pstate.value != nullValue); // addDim is a no-op with nullValue
1267 			}
1268 	}
1269 
1270 	/// Peek at the subset the algorithm is currently working with.
1271 	@property Set currentSubset()
1272 	{
1273 		assert(workingSet !is Set.emptySet, "Not iterating");
1274 		flush();
1275 		return workingSet;
1276 	}
1277 
1278 	/// Get all possible values for this variable at this point.
1279 	/// Should be used mainly for debugging.
1280 	/*private*/ const(V)[] getAll(A name)
1281 	{
1282 		assert(workingSet !is Set.emptySet, "Not iterating");
1283 		if (auto pstate = name in varState)
1284 			if (pstate.haveValue)
1285 				return (&pstate.value)[0 .. 1];
1286 
1287 		return workingSet.all(name);
1288 	}
1289 
1290 	/// Algorithm interface - get a value by name
1291 	V get(A name)
1292 	{
1293 		assert(workingSet !is Set.emptySet, "Not iterating");
1294 		auto pstate = &varState.require(name);
1295 		if (pstate.haveValue)
1296 			return pstate.value;
1297 
1298 		// We are going to narrow the workingSet - update inSet appropriately
1299 		foreach (varName, ref state; varState)
1300 			if (varName == name)
1301 				state.inSet = Maybe.maybe;
1302 			else
1303 			if (state.inSet == Maybe.yes)
1304 				state.inSet = Maybe.maybe;
1305 
1306 		if (stackPos == stack.length)
1307 		{
1308 			// Expand new variable
1309 			auto values = workingSet.all(name);
1310 			auto value = values[0];
1311 			// auto pstate = varState
1312 			pstate.value = value;
1313 			pstate.haveValue = true;
1314 			stack ~= Var(name, values, 0);
1315 			stackPos++;
1316 			if (values.length > 1)
1317 				workingSet = workingSet.get(name, value);
1318 			return value;
1319 		}
1320 
1321 		// Iterate over known variable
1322 		auto var = &stack[stackPos];
1323 		assert(var.name == name, "Mismatching get order");
1324 		auto value = var.values[var.pos];
1325 		workingSet = workingSet.get(var.name, value);
1326 		assert(workingSet !is Set.emptySet, "Empty set after restoring");
1327 		pstate.value = value;
1328 		pstate.haveValue = true;
1329 		stackPos++;
1330 		return value;
1331 	}
1332 
1333 	/// Algorithm interface - set a value by name
1334 	void put(A name, V value)
1335 	{
1336 		assert(workingSet !is Set.emptySet, "Not iterating");
1337 
1338 		auto pstate = &varState.require(name);
1339 
1340 		if (pstate.haveValue && pstate.value == value)
1341 			return; // All OK
1342 
1343 		pstate.value = value;
1344 		pstate.haveValue = pstate.dirty = true;
1345 
1346 		// Flush null values ASAP, to avoid non-null values
1347 		// accumulating in the set and increasing the set size.
1348 		if (value == nullValue)
1349 			flush(name);
1350 	}
1351 
1352 	/// Prepare an unresolved variable for overwriting (with more than
1353 	/// one value).
1354 	private void destroy(A name)
1355 	{
1356 		auto pstate = &varState.require(name);
1357 		pstate.haveValue = pstate.dirty = false;
1358 		if (pstate.inSet >= Maybe.maybe)
1359 		{
1360 			workingSet = workingSet.remove(name);
1361 			pstate.inSet = Maybe.no;
1362 		}
1363 	}
1364 
1365 	/// A smarter `workingSet = workingSet.bringToFront(name)`, which
1366 	/// checks if `name` is in the set first.
1367 	private void bringToFront(A name)
1368 	{
1369 		auto pState = &varState.require(name);
1370 		if (pState.inSet == Maybe.no)
1371 			workingSet = Set(new immutable Set.Node(name, [Set.Pair(nullValue, workingSet)])).deduplicate;
1372 		else
1373 			workingSet = workingSet.bringToFront(name);
1374 		pState.inSet = Maybe.yes;
1375 	}
1376 
1377 	/// Algorithm interface - copy a value target another name,
1378 	/// without resolving it (unless it's already resolved).
1379 	void copy(bool reorder = false)(A source, A target)
1380 	{
1381 		if (source == target)
1382 			return;
1383 
1384 		auto pSourceState = &varState.require(source);
1385 		auto pTargetState = &varState.require(target);
1386 		if (pSourceState.haveValue)
1387 		{
1388 			put(target, pSourceState.value);
1389 			return;
1390 		}
1391 
1392 		assert(workingSet !is Set.emptySet, "Not iterating");
1393 
1394 		static if (reorder)
1395 		{
1396 			destroy(target);
1397 
1398 			bringToFront(source);
1399 			auto newChildren = workingSet.root.children.dup;
1400 			foreach (ref pair; newChildren)
1401 			{
1402 				auto set = Set(new immutable Set.Node(source, [Set.Pair(pair.value, pair.set)])).deduplicate;
1403 				pair = Set.Pair(pair.value, set);
1404 			}
1405 			workingSet = Set(new immutable Set.Node(target, cast(immutable) newChildren)).deduplicate;
1406 			pSourceState.inSet = Maybe.yes;
1407 			pTargetState.inSet = Maybe.yes;
1408 		}
1409 		else
1410 		{
1411 			targetTransform!true(source, target,
1412 				(ref const V inputValue, out V outputValue)
1413 				{
1414 					outputValue = inputValue;
1415 				}
1416 			);
1417 		}
1418 	}
1419 
1420 	/// Apply a function over every possible value of the given
1421 	/// variable, without resolving it (unless it's already resolved).
1422 	void transform(A name, scope void delegate(ref V value) fun)
1423 	{
1424 		assert(workingSet !is Set.emptySet, "Not iterating");
1425 		auto pState = &varState.require(name);
1426 		if (pState.haveValue)
1427 		{
1428 			pState.dirty = true;
1429 			fun(pState.value);
1430 			return;
1431 		}
1432 
1433 		bringToFront(name);
1434 		Set[V] newChildren;
1435 		foreach (ref child; workingSet.root.children)
1436 		{
1437 			V value = child.value;
1438 			fun(value);
1439 			newChildren.updateVoid(value,
1440 				() => child.set,
1441 				(ref Set set)
1442 				{
1443 					set = set.merge(child.set);
1444 				});
1445 		}
1446 		workingSet = Set(new immutable Set.Node(name, cast(immutable) newChildren)).deduplicate;
1447 		pState.inSet = Maybe.yes;
1448 	}
1449 
1450 	/// Apply a function over every possible value of the given
1451 	/// variable, without resolving it (unless it's already resolved).
1452 	/// The function is assumed to be injective (does not produce
1453 	/// duplicate outputs for distinct inputs).
1454 	void injectiveTransform(A name, scope void delegate(ref V value) fun)
1455 	{
1456 		assert(workingSet !is Set.emptySet, "Not iterating");
1457 		auto pState = &varState.require(name);
1458 		if (pState.haveValue)
1459 		{
1460 			pState.dirty = true;
1461 			fun(pState.value);
1462 			return;
1463 		}
1464 
1465 		bringToFront(name);
1466 		auto newChildren = workingSet.root.children.dup;
1467 		foreach (ref child; newChildren)
1468 			fun(child.value);
1469 		newChildren.sort();
1470 		workingSet = Set(new immutable Set.Node(name, cast(immutable) newChildren)).deduplicate;
1471 		pState.inSet = Maybe.yes;
1472 	}
1473 
1474 	/// Perform a transformation with one input and one output.
1475 	/// Does not reorder the MapSet.
1476 	void targetTransform(bool injective = false)(A input, A output, scope void delegate(ref const V inputValue, out V outputValue) fun)
1477 	{
1478 		assert(input != output, "Input is the same as output - use transform instead");
1479 
1480 		flush(input);
1481 		destroy(output);
1482 
1483 		auto pInputState = &varState.require(input);
1484 		auto pOutputState = &varState.require(output);
1485 
1486 		if (pInputState.haveValue)
1487 		{
1488 			V outputValue;
1489 			fun(pInputState.value, outputValue);
1490 			put(output, outputValue);
1491 			return;
1492 		}
1493 
1494 		bool sawInput, addedOutput;
1495 		Set[Set] cache;
1496 		Set visit(Set set)
1497 		{
1498 			return cache.require(set, {
1499 				if (set == Set.unitSet)
1500 				{
1501 					V inputValue = nullValue;
1502 					V outputValue;
1503 					fun(inputValue, outputValue);
1504 					if (outputValue != nullValue)
1505 						addedOutput = true;
1506 					return set.addDim(output, outputValue);
1507 				}
1508 				if (set.root.dim == input)
1509 				{
1510 					sawInput = true;
1511 					static if (injective)
1512 					{
1513 						auto outputChildren = new Set.Pair[set.root.children.length];
1514 						foreach (i, ref outputPair; outputChildren)
1515 						{
1516 							auto inputPair = &set.root.children[i];
1517 							fun(inputPair.value, outputPair.value);
1518 							outputPair.set = Set(new immutable Set.Node(input, inputPair[0..1])).deduplicate;
1519 						}
1520 						outputChildren.sort();
1521 						addedOutput = true;
1522 						return Set(new immutable Set.Node(output, cast(immutable) outputChildren)).deduplicate;
1523 					}
1524 					else
1525 					{
1526 						auto inputChildren = set.root.children;
1527 						set = Set.emptySet;
1528 						foreach (i, ref inputPair; inputChildren)
1529 						{
1530 							V outputValue;
1531 							fun(inputPair.value, outputValue);
1532 							auto inputSet = Set(new immutable Set.Node(input, (&inputPair)[0..1])).deduplicate;
1533 							auto outputSet = inputSet.addDim(output, outputValue);
1534 							set = set.merge(outputSet);
1535 							if (outputValue != nullValue)
1536 								addedOutput = true;
1537 						}
1538 						return set;
1539 					}
1540 				}
1541 				else
1542 					return set.lazyMap(&visit);
1543 			}());
1544 		}
1545 		workingSet = visit(workingSet);
1546 		pInputState.setInSet(sawInput);
1547 		pOutputState.setInSet(addedOutput);
1548 	}
1549 
1550 	/// Perform a transformation with multiple inputs and outputs.
1551 	/// Inputs and outputs must not overlap.
1552 	/// Can be used to perform binary operations, copy-transforms, and more.
1553 	void multiTransform(scope A[] inputs, scope A[] outputs, scope void delegate(scope V[] inputValues, scope V[] outputValues) fun)
1554 	{
1555 		assert(inputs.length > 0 && outputs.length > 0, "");
1556 		foreach (output; outputs)
1557 			assert(!inputs.canFind(output), "Inputs and outputs overlap");
1558 
1559 		foreach (input; inputs)
1560 			flush(input);
1561 		foreach (output; outputs)
1562 			destroy(output);
1563 		foreach_reverse (input; inputs)
1564 			bringToFront(input);
1565 
1566 		Set resultSet = Set.emptySet;
1567 		auto inputValues = new V[inputs.length];
1568 		auto outputValues = new V[outputs.length];
1569 		auto addedInput = new bool[inputs.length];
1570 		auto addedOutput = new bool[outputs.length];
1571 
1572 		void visit(Set set, size_t depth)
1573 		{
1574 			if (depth == inputs.length)
1575 			{
1576 				fun(inputValues, outputValues);
1577 				foreach_reverse (i, input; inputs)
1578 				{
1579 					set = set.addDim(input, inputValues[i]);
1580 					if (inputValues[i] != nullValue)
1581 						addedInput[i] = true;
1582 				}
1583 				foreach_reverse (i, output; outputs)
1584 				{
1585 					set = set.addDim(output, outputValues[i]);
1586 					if (outputValues[i] != nullValue)
1587 						addedOutput[i] = true;
1588 				}
1589 				resultSet = resultSet.merge(set);
1590 			}
1591 			else
1592 			{
1593 				assert(set.root.dim == inputs[depth]);
1594 				foreach (ref pair; set.root.children)
1595 				{
1596 					inputValues[depth] = pair.value;
1597 					visit(pair.set, depth + 1);
1598 				}
1599 			}
1600 		}
1601 		visit(workingSet, 0);
1602 		workingSet = resultSet;
1603 		foreach (i, input; inputs)
1604 			varState.require(input).setInSet(addedInput[i]);
1605 		foreach (i, output; outputs)
1606 			varState.require(output).setInSet(addedOutput[i]);
1607 	}
1608 
1609 	/// Inject a variable and values to iterate over.
1610 	/// The variable must not have been resolved yet.
1611 	void inject(A name, V[] values)
1612 	{
1613 		assert(values.length > 0, "Injecting zero values would result in an empty set");
1614 		destroy(name);
1615 		workingSet = workingSet.uncheckedCartesianProduct(name, values);
1616 		varState.require(name).inSet = Maybe.yes;
1617 	}
1618 }
1619 
1620 /// An algorithm which divides two numbers.
1621 /// When the divisor is zero, we don't even query the dividend,
1622 /// therefore processing all dividends in one iteration.
1623 unittest
1624 {
1625 	alias M = MapSet!(string, int);
1626 	M m = M.unitSet
1627 		.cartesianProduct("divisor" , [0, 1, 2])
1628 		.cartesianProduct("dividend", [0, 1, 2]);
1629 	assert(m.count == 9);
1630 
1631 	auto v = MapSetVisitor!(string, int)(m);
1632 	M results;
1633 	int iterations;
1634 	while (v.next())
1635 	{
1636 		iterations++;
1637 		auto divisor = v.get("divisor");
1638 		if (divisor == 0)
1639 			continue;
1640 		auto dividend = v.get("dividend");
1641 		v.put("quotient", dividend / divisor);
1642 		results = results.merge(v.currentSubset);
1643 	}
1644 
1645 	assert(iterations == 7); // 1 for division by zero + 3 for division by one + 3 for division by two
1646 	assert(results.get("divisor", 2).get("dividend", 2).all("quotient") == [1]);
1647 	assert(results.get("divisor", 0).count == 0);
1648 }
1649 
1650 unittest
1651 {
1652 	import std.algorithm.sorting : sort;
1653 
1654 	alias M = MapSet!(string, int);
1655 	M m = M.unitSet.cartesianProduct("x", [1, 2, 3]);
1656 	auto v = MapSetVisitor!(string, int)(m);
1657 	v.next();
1658 	v.transform("x", (ref int v) { v *= 2; });
1659 	assert(v.currentSubset.all("x").dup.sort.release == [2, 4, 6]);
1660 }
1661 
1662 unittest
1663 {
1664 	alias M = MapSet!(string, int);
1665 	M m = M.unitSet.cartesianProduct("x", [1, 2, 3]);
1666 	auto v = MapSetVisitor!(string, int)(m);
1667 	while (v.next())
1668 	{
1669 		v.transform("x", (ref int v) { v *= 2; });
1670 		v.put("y", v.get("x"));
1671 	}
1672 }
1673 
1674 // Test that initialVarState does not interfere with flushing
1675 unittest
1676 {
1677 	alias M = MapSet!(string, int);
1678 	M m = M.unitSet.cartesianProduct("x", [1]);
1679 	auto v = MapSetVisitor!(string, int)(m);
1680 	assert(v.initialVarState["x"].haveValue);
1681 	v.next();
1682 	v.put("x", 2);
1683 	v.currentSubset;
1684 	assert(v.get("x") == 2);
1685 }
1686 
1687 // Test resolving the same variable several times
1688 unittest
1689 {
1690 	alias M = MapSet!(string, int);
1691 	M m = M.unitSet.cartesianProduct("x", [10, 20, 30]);
1692 	auto v = MapSetVisitor!(string, int)(m);
1693 	int[] result;
1694 	while (v.next())
1695 	{
1696 		auto a = v.get("x"); // First resolve
1697 		v.inject("x", [1, 2, 3]);
1698 		auto b = v.get("x"); // Second resolve
1699 		result ~= a + b;
1700 	}
1701 	assert(result == [11, 12, 13, 21, 22, 23, 31, 32, 33]);
1702 }
1703 
1704 // targetTransform
1705 unittest
1706 {
1707 	import std.algorithm.sorting : sort;
1708 
1709 	alias M = MapSet!(string, int);
1710 	M m = M.unitSet.cartesianProduct("x", [1, 2, 3, 4, 5]);
1711 	auto v = MapSetVisitor!(string, int)(m);
1712 	v.next();
1713 	v.targetTransform("x", "y", (ref const int input, out int output) { output = input + 1; });
1714 	assert(v.currentSubset.all("y").dup.sort.release == [2, 3, 4, 5, 6]);
1715 	assert(!v.next());
1716 }
1717 
1718 // multiTransform
1719 unittest
1720 {
1721 	alias M = MapSet!(string, int);
1722 	M m = M.unitSet.cartesianProduct("x", [1, 2, 3, 4, 5]);
1723 	auto v = MapSetVisitor!(string, int)(m);
1724 	v.next();
1725 	v.copy("x", "y");
1726 	v.transform("y", (ref int v) { v = 6 - v; });
1727 	v.multiTransform(["x", "y"], ["z"], (int[] inputs, int[] outputs) { outputs[0] = inputs[0] + inputs[1]; });
1728 	assert(v.currentSubset.all("z") == [6]);
1729 	assert(!v.next());
1730 }