The Tax Historian
C++ using the BuDDy library
This seemed like a nice excuse to play with binary decision diagrams. The kingdom is converted into a big boolean formula where we have to count the number of ways in which it can be satisfied. That can (sometimes) be done more efficiently than it sounds.
The kingdom must be given as a program constant as a flat array and explicitely given dimensions. (Nice input is left as an execise for the reader :-)
Here is the embarrassingly simple code:
#include <iostream>
#include <bdd.h>
// describe the kingdom here:
constexpr int ROWS = 4;
constexpr int COLS = 10;
constexpr int a[] = {
1, 1, 0, 1, 0, 1, 0, 1, 1, 0,
1, 1, 0, 0, 0, 0, 1, 1, 1, 0,
1, 1, 0, 0, 0, 0, 0, 0, 0, 1,
0, 1, 0, 0, 0, 0, 1, 1, 0, 0,
};
// end of description
// check dimensions
static_assert(ROWS*COLS*sizeof(int)==sizeof(a),
"ROWS*COLS must be the number of entries of a");
// dimensions of previous generation
constexpr int R1 = ROWS+1;
constexpr int C1 = COLS+1;
// condition that exactly one is true
bdd one(bdd a, bdd b, bdd c, bdd d){
bdd q = a & !b & !c & !d;
q |= !a & b & !c & !d;
q |= !a & !b & c & !d;
q |= !a & !b & !c & d;
return q;
}
int main()
{
bdd_init(1000000, 10000); // tuneable, but not too important
bdd_setvarnum(R1*C1);
bdd q { bddtrue };
for(int j=COLS-1; j>=0; j--) // handle high vars first
for (int i=ROWS-1; i>=0; i--){
int x=i+R1*j;
bdd p=one(bdd_ithvar(x), bdd_ithvar(x+1),
bdd_ithvar(x+R1), bdd_ithvar(x+R1+1));
if (!a[COLS*i+j])
p = !p;
q &= p;
}
std::cout << "There are " << bdd_satcount(q) << " preimages\n";
bdd_done();
}
To compile with debian 8 (jessie), install libbdd-dev and do g++ -std=c++11 -o hist hist.cpp -lbdd. (Optimizing options make almost no difference because the real work is done in the library.)
Big examples can lead to messages about garbage collection. They could be suppressed, but I prefer to see them.
bdd_satcount returns a double, but that is good enough for the expected range of results. The same counting technique is possible with exact (big) integers.
The code is optimized for ROWS<COLS. If you have a lot more rows than columns, it might be a good idea to transpose the matrix.
Python 2.7
This is just a naïve first attempt. It's not particularly fast, but it is correct.
The first observation is that each cell is dependent on exactly four cells in the previous generation. We can represent those four cells as a 4-bit number (0-15). According to the rules, if exactly one neighboring cell in the previous generation is 1, then a given cell in the current generation will be 1, otherwise, it will be 0. Those correspond to the powers of two, namely, [1, 2, 4, 8]. When the four ancestors are represented as a 4-bit number, any other number will result in a 0 in the current generation. With this information, upon seeing a cell in the current generation, we can narrow down the possibilities of the neighborhood in the previous generation to one of four or one of twelve possibilities respectively.
I've chosen to represent the neighborhood as follows:
32
10
where 0 is the least-significant bit, and so on.
The second observation is that for two adjacent cells in the current generation, the two neighborhoods from the previous generation overlap:
32 32
10 10
or:
32
10
32
10
In the horizontal case, the 2 from the left neighborhood overlaps with the 3 from the right neighborhood, and similarly with the 0 on the left and the 1 on the right. In the vertical case, the 1 from the top neighborhood overlaps with the 3 from the bottom neighborhood, and similarly with the 0 on the top and the 2 on the bottom.
This overlap means that we can narrow down the possibilities for yet-unchosen neighborhoods based on what we've already chosen. The code works it's way from left to right, top to bottom, in a recursive depth-first search for possible preimages. The following diagram indicates which previous neighborhoods we have to consider when looking at the possible neighborhoods of a current cell:
f = free choice
h = only have to look at the neighborhood to the left
v = only have to look at the neighborhood to the top
b = have to look at both left and top neighborhoods
[f, h, h, h, h],
[v, b, b, b, b],
[v, b, b, b, b],
[v, b, b, b, b]
Here's the code:
def good_horizontal(left, right):
if (left & 4) >> 2 != (right & 8) >> 3:
return False
if left & 1 != (right & 2) >> 1:
return False
return True
def good_vertical(bottom, top):
if (bottom & 8) >> 3 != (top & 2) >> 1:
return False
if (bottom & 4) >> 2 != (top & 1):
return False
return True
ones = [1, 2, 4, 8]
zeros = [0, 3, 5, 6, 7, 9, 10, 11, 12, 13, 14, 15]
h = {}
v = {}
for i in range(16):
h[i] = [j for j in range(16) if good_horizontal(i, j)]
v[i] = [j for j in range(16) if good_vertical(i, j)]
def solve(arr):
height = len(arr)
width = len(arr[0])
if height == 1 and width == 1:
if arr[0][0] == 1:
return 4
else:
return 12
return solve_helper(arr)
def solve_helper(arr, i=0, j=0, partial=None):
height = len(arr)
width = len(arr[0])
if arr[i][j] == 1:
poss = ones
else:
poss = zeros
if i == height - 1 and j == width - 1: # We made it to the end of this chain
if height == 1:
return sum([1 for p in poss if p in h[partial[-1][-1]]])
else:
return sum([1 for p in poss if partial[-2][-1] in v[p] and p in h[partial[-1][-1]]])
if j == width - 1:
new_i, new_j = i + 1, 0
else:
new_i, new_j = i, j + 1
if i == 0:
if j == 0:
# first call
return sum([solve_helper(arr, new_i, new_j, [[p]]) for p in poss])
# still in the first row
return sum([solve_helper(arr, new_i, new_j, [partial[0] + [p]]) for p in poss if p in h[partial[0][-1]]])
if j == 0: # starting a new row
return sum([solve_helper(arr, new_i, new_j, [r for r in partial + [[p]]]) for p in poss if partial[i - 1][0] in v[p]])
return sum([solve_helper(arr, new_i, new_j, [r for r in partial[:-1] + ([partial[-1] + [p]])]) for p in poss if p in h[partial[i][-1]] and partial[i - 1][j] in v[p]])
To run it:
test3 = [
[1, 1, 0, 1, 0, 1, 0, 1, 1, 0],
[1, 1, 0, 0, 0, 0, 1, 1, 1, 0],
[1, 1, 0, 0, 0, 0, 0, 0, 0, 1],
[0, 1, 0, 0, 0, 0, 1, 1, 0, 0]
]
expected3 = 11567
assert(solve(test3) == expected3)