About Me
Posts tagged duke
One of my greatest achievements at Duke was getting an A+ in COMPSCI 100E (I think it is 201, nowadays?). It's a bit sad how little I have to be proud of, but my only other A+ came from an intro physics course unless you count study abroad. I'd also like to make a shoutout to my fellow math major Robert Won, who took COMPSCI 100E with me.
I did all the extra credit and every single Algorithmic Problem-solving Testing (APT) for full credit. I remember the hardest one being APT CountPaths. I remember neglecting all my other coursework and spending several days on this problem, which apparently isn't suppose to be that hard. I recently encountered a variation of this problem again, here, so I figured that I would spend some time writing up a solution. My original solution that I wrote back at Duke is a whopping 257 lines of Java. My current solution is less than 100 lines in C++, which includes comments.
Beginner Version
In the most basic version of the problem, we're trying to get from one corner to the other corner.
Going from point A to B, only down and right moves are permitted. If the grid has $R$ rows and $C$ columns. We need to make $R-1$ down moves and $C-1$ right moves. The order doesn't matter, and any down move is indistinguishable from another down move. Thus, the number of paths is $$\boxed{\frac{(R + C - 2)!}{(R-1)!(C-1)!} = {{R + C - 2}\choose{R-1}}}.$$
Intermediate Version
In another version, there are some squares that we are not allowed to visit, and the grid is not too big, say less than a billion squares.
Here, we are not allowed to visit black squares, $x_1$, $x_2$, $x_3$, and $x_4$. This version admits a dynamic programming solution that is $O(RC)$ in time complexity and $O\left(\max(R,C)\right)$ in memory since we only need to keep track of two rows or columns at a time. Let $P(r, c)$ be the number of paths to square $(r,c)$. We have that $$ P(r,c) = \begin{cases} 1, &\text{if $r = 0$ and $c=0$} \\ 0, &\text{if $(r,c)$ is a black square} \\ P(r-1,c) + P(r,c-1), &\text{if $0 \leq r < R$ and $0 \leq c < C$ and $(r,c) \neq (0,0)$} \\ 0, &\text{otherwise}. \end{cases} $$
We simply loop through rows and columns in order with our loop invariant being that at square $(r,c)$, we know the number of paths to squares $(r^\prime, c^\prime)$, where either $r^\prime < r$ or $r^\prime = r$ and $c^\prime < c$. See the code in C++ below with some variations since the problem asks for number of paths modulo 1,000,000,007.
#include <algorithm>
#include <iostream>
#include <vector>
using namespace std;
const int MOD = 1000000007;
int main(int argc, char *argv[]) {
int h, w, n;
cin >> h >> w >> n;
vector<vector<int>> grid(2, vector<int>(w, 0));
vector<pair<int, int>> blackSquares;
for (int i = 0; i < n; ++i) {
int r, c;
cin >> r >> c;
blackSquares.emplace_back(r - 1, c - 1);
}
sort(blackSquares.begin(), blackSquares.end());
grid.front()[0] = 1; // guaranteed that 0,0 is white
auto blackIt = blackSquares.begin();
for (int c = 1; c < w; ++c) {
if (blackIt != blackSquares.end() && blackIt -> first == 0 && blackIt -> second == c) {
grid.front()[c] = -1;
++blackIt;
} else if (grid.front()[c-1] != -1) {
grid.front()[c] = 1;
}
}
for (int r = 1; r < h; ++r) {
if (blackIt != blackSquares.end() && blackIt -> first == r && blackIt -> second == 0) {
grid.back()[0] = -1;
++blackIt;
} else if (grid.front()[0] != -1) {
grid.back()[0] = 1;
} else {
grid.back()[0] = 0;
}
for (int c = 1; c < w; ++c) {
if (blackIt != blackSquares.end() && blackIt -> first == r && blackIt -> second == c) {
grid.back()[c] = -1;
++blackIt;
} else {
grid.back()[c] = 0;
if (grid.front()[c] != -1) grid.back()[c] += grid.front()[c];
if (grid.back()[c-1] != -1) grid.back()[c] += grid.back()[c-1];
grid.back()[c] %= MOD;
}
}
grid.front().swap(grid.back());
}
cout << grid.front()[w-1] << endl;
return 0;
}
Advanced Version
Now, if we have a very large grid, the above algorithm will take too long. If $1 \leq R,C \leq 10^5$, we could have as many as 10 billion squares. Assume that number of black squares $N \ll RC$ ($N$ much less than $RC$). Then, we can write an algorithm that runs in $O\left(N^2\log\left(\max(R,C)\right)\right)$.
First, sort the black squares. Given $x = (r_1, c_1)$ and $y = (r_2, c_2)$, let $x < y$ if either $r_1 < r_2$ or $r_1 = r_2$ and $c_1 < c_2$. Let $x_1,x_2,\ldots,x_N$ be the black squares, where $x_i = (r_i,c_i)$. For each $x_i$, associate the set $$S_i = \left\{x_j = (r_j, y_j) : r_j \leq r_i~\text{and}~c_j \leq c_i~\text{and}~x_j \neq x_i\right\}.$$
Recalling our previous diagram,
$S_1 = \emptyset$, $S_2 = S_3 = \{ x_1 \}$, $S_4 = \{x_1,x_3\}$, and $S_B = \{x_1, x_2, x_3, x_4\}$. Let let us loop through these points in order, and calculate the number of paths to $x_i$ avoiding all the squares in $S_i$. It's clear that for all $x_j \in S_i$, $j < i$, so we have our loop invariant: at $i$, for all $j < i$, we know the number of paths to $x_j$ that do not pass through any square in $S_j$. Assuming this, let us calculate the number of paths to $x_i$.
If $S_i$ is empty, the problem reduces to the beginner version and $$P(x_i) = P(r_i,c_i) = {{r_i + c_i - 2}\choose{r_i-1}}.$$ If $S_i$ is not empty, we need to subtract any paths that go through any $x_j \in S_i$, so if $i = 4$, we want to subtract paths like these red ones:
Thus, there are 3 types of paths we need to subtract. Paths that go through $x_1$ and not $x_3$, paths that go through $x_3$ and not $x_1$, and paths that go through both $x_1$ and $x_3$.
It might be helpful to think about this as a Venn diagram:
We see that the set of these paths is the same as any path that goes through $x_1$ and paths that go through $x_3$ but not $x_1$. The number of paths that go through $x_1$ is $$P(x_1)\cdot{(r_4 - r_1) + (c_4 - c_1)\choose{r_4 - r_1}},$$ and since $x_1 \in S_3$, the number of paths that go through $x_3$ but not $x_1$ is $$P(x_3)\cdot{(r_4 - r_3) + (c_4 - c_3)\choose{r_4 - r_3}}.$$
Generalizing this, we have that $$P(x_i) = {{r_i + c_i - 2}\choose{r_i-1}} - \sum_{x_j \in S_i}P(x_j)\cdot{{(r_i-r_j) + (c_i - c_j)}\choose{r_i-r_j}},$$ which leads to the $N^2$ factor in the time complexity. Treating B as a black square, we can calculate the number of paths to B that avoid $x_1, x_2,\ldots,x_N$.
Now, the problem is further complicated by the fact that we need to return the number of paths modulo 1,000,000,007. Computing the binomial coefficient involves division, so we need to be able to find the modular multiplicative inverse. We can do this with the Euclidean algorithm, which gives us the $\log\left(\max(R, C)\right)$ factor in the time complexity.
Since $m = 10^9 + 7$ is prime, given any $n < m$, $\gcd(n,m) = 1$, so there exists $w,x \in \mathbb{Z}$ such that $wn + xm = 1$. We can find $w$ and $x$ by setting up a linear system, $$ \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix}\begin{pmatrix} n \\ m \end{pmatrix} = \begin{pmatrix} n \\ m \end{pmatrix}, $$ which can be rewritten as an augmented matrix, $$ \left(\begin{array}{cc|c} 1 & 0 & n\\ 0 & 1 & m \end{array}\right). $$ You can perform the Euclidean algorithm by doing elementary row operations until a $1$ appears in the rightmost column. The first and second column of that row is $w$ and $x$, respectively. Since the system remains consistent with $(n,m)^\intercal$ as a solution, we have $$wn + xm = 1 \Rightarrow wn = 1 - xm \equiv 1 \bmod m,$$ so the modular multiplicative inverse of $n$ is $w$.
If you memoize this step and precalculate all the modular multiplicative inverses up to $\max(R,C)$, you can actually achieve a $O(N^2)$ algorithm, but this is unnecessary for this problem. Here's the code below.
#include <algorithm>
#include <iostream>
#include <vector>
using namespace std;
const long long MOD = 1000000007LL;
const int FACTORIALS_LENGTH = 199999;
long long factorials[FACTORIALS_LENGTH];
long long invertMod(int n, int m) {
/*
* find the modular inverse with
* extended euclidean algorithm
* [w x][a] = [a]
* [y z][b] = [b]
*/
n %= m;
long w = 1;
long x = 0;
long y = 0;
long z = 1;
long a = n;
long b = m;
while (a != 0 && b != 0) {
if (a <= b) {
int q = b / a;
int r = b % a;
b = r;
z -= q*x; z %= m;
y -= q*w; y %= m;
} else if (a > b) {
int q = a / b;
int r = a % b;
a = r;
x -= q*z; x %= m;
w -= q*y; w %= m;
}
}
int inverseMod = a != 0 ? w : y;
if (inverseMod < 0) inverseMod += m;
return inverseMod;
}
long long choose(int n, int r) {
long long res = (factorials[n]*invertMod(factorials[r], MOD)) % MOD;
res *= invertMod(factorials[n-r], MOD);
return res % MOD;
}
int main(int argc, char *argv[]) {
factorials[0] = 1;
for (int n = 1; n < FACTORIALS_LENGTH; ++n) factorials[n] = (factorials[n-1]*n) % MOD;
int h, w, n;
cin >> h >> w >> n;
vector<pair<int, int>> blackSquares;
vector<long long> blackSquarePaths;
for (int i = 0; i < n; ++i) {
int r, c;
cin >> r >> c;
blackSquares.emplace_back(r - 1, c - 1);
blackSquarePaths.push_back(0);
}
// just treat the end point as a black square
blackSquares.emplace_back(h - 1, w - 1);
blackSquarePaths.push_back(0);
sort(blackSquares.begin(), blackSquares.end());
for (int i = 0; i < blackSquares.size(); ++i) {
// start by considering all paths
blackSquarePaths[i] = choose(blackSquares[i].first + blackSquares[i].second, blackSquares[i].first);
for (int j = 0; j < i; ++j) {
// subtract unallowed paths, we already know row is smaller because we sorted already
if (blackSquares[j].second <= blackSquares[i].second) {
/*
* A
* B
* C
* if we're trying to get to C, we need to subtract all the paths that go through, A, B, or both
* we hit A first and subtract paths that go through A and both A and B
* we hit B and subtract all paths that go through B but not A
*
* loop invariant is that we know how many paths to get to a
* black square by avoiding other black squares
*/
int rDiff = blackSquares[i].first - blackSquares[j].first;
int cDiff = blackSquares[i].second - blackSquares[j].second;
blackSquarePaths[i] -= (blackSquarePaths[j]*choose(rDiff + cDiff, rDiff)) % MOD;
blackSquarePaths[i] %= MOD;
if (blackSquarePaths[i] < 0) blackSquarePaths[i] += MOD;
}
}
}
cout << blackSquarePaths.back() << endl;
return 0;
}
Duke APT CountPath Version
Now the version that I did at Duke has a few key differences. First, you start in the lower left and end up in the upper right. More importantly, you need to calculate the paths that go through 1 black square, 2 black squares, 3 black squares, etc. We need to do this up to paths that go through all $N$ black squares. Also, the grid is much smaller, and the black squares must be visited in the order given.
Since the grid is much smaller, a dynamic programming solution works where the state is four variables $(r,c,k,l)$, where $k$ is the index of last black square visited and $l$ is the number of black squares on that path. In general, if $T$ is the set of indicies of black squares that are to the lower-left of $k$ and appear earlier in the order given, we have that $$ P(r,c,k,l) = \begin{cases} P(r-1,c,k,l) + P(r,c-1,k,l), &\text{if $(r,c)$ is not a black square} \\ \sum_{k^\prime \in T} P(r-1,c,k^\prime,l-1) + P(r,c-1,k^\prime,l-1), &\text{if $(r,c)$ is a black square.} \end{cases} $$ This is $O(RCN^2)$ in both time and space. Here is the Java code.
import java.util.Arrays;
public class CountPaths {
private static final int MOD = 1000007;
public int[] difPaths(int r, int c, int[] fieldrow, int[] fieldcol) {
// row, col, last special field, num of previous special fields,
// need to add 1 to third index, in case there are 0 special fields
int[][][][] memo = new int[r+1][c+1][fieldrow.length+1][fieldrow.length+1];
int[][] isSpecial = new int[r+1][c+1];
for (int i = 0; i <= r; ++i) { Arrays.fill(isSpecial[i], -1); };
for (int i = 0; i < fieldrow.length; ++i) { isSpecial[fieldrow[i]][fieldcol[i]] = i; }
// 00 is special, marks no special fields yet
memo[1][0][0][0] = 1; // so that memo[1][1] gets 1 path
for (int i = 1; i <= r; ++i) {
for (int j = 1; j <= c; ++j) {
int kLimit; int lStart;
// handle 00 case where no special fields have been reached yet
if (isSpecial[i][j] == -1) {
// l = 0 case, no special field came before it
memo[i][j][0][0] += memo[i-1][j][0][0] + memo[i][j-1][0][0];
memo[i][j][0][0] %= MOD;
kLimit = fieldrow.length; lStart = 1;
} else {
// since this field is special min length is 1
memo[i][j][isSpecial[i][j]][1] += memo[i-1][j][0][0] + memo[i][j-1][0][0];
memo[i][j][isSpecial[i][j]][1] %= MOD;
kLimit = isSpecial[i][j]; lStart = 2;
}
for (int l = lStart; l <= fieldrow.length; ++l) {
// only check special fields that come before
for (int k = 0; k < kLimit; ++k) {
// this point is only reachable from special field if it is up and to the right
if (i >= fieldrow[k] && j >= fieldcol[k]) {
if (isSpecial[i][j] == -1) {
// when the field is not special, length is the same
memo[i][j][k][l] += memo[i-1][j][k][l] + memo[i][j-1][k][l];
memo[i][j][k][l] %= MOD;
} else if (isSpecial[i][j] > -1) {
// when the field is special add 1 to the length
memo[i][j][kLimit][l] += memo[i-1][j][k][l-1] + memo[i][j-1][k][l-1];
memo[i][j][kLimit][l] %= MOD;
}
}
}
}
}
}
int[] paths = new int[fieldrow.length+1];
for (int k = 0; k < fieldrow.length; ++k) {
for (int l = 1; l < fieldrow.length + 1; ++l) {
paths[l] += memo[r][c][k][l];
paths[l] %= MOD;
}
}
paths[0] += memo[r][c][0][0]; paths[0] %= MOD;
return paths;
}
}
Using ideas from the advanced version, we can significantly reduce the memory needed. Basically, we represent each black square as a node and construct a directed graph, represented by an adjacency matrix, where the edges are weighted by the number of paths from one black square to another that avoid all other black squares. This takes $O(N^4)$ time and $O(N^2)$ space. Now, let the nodes be numbered $1,2,\ldots,N$. Now, we can calculate the paths that go through any number of black squares with some dynamic programming. Let $G(j,k)$ be the number of paths from black square $j$ and black square $k$ that don't go through any other black squares. Let $Q(k, l)$ be the number paths of length $l$ where $k$ is the last black square. We have that $$Q(k,l) = \sum_{j < k} Q(j, l-1)\cdot G(j,k), $$ keeping in mind that the nodes are sorted, and $G(j,k) = 0$ if the black squares are in the wrong order. Here's the C++ code below.
#include <algorithm>
#include <cassert>
#include <iostream>
#include <sstream>
#include <vector>
using namespace std;
class CountPaths {
private:
const int MOD = 1000007;
vector<vector<int>> choose;
vector<vector<int>> initializeChoose(int N) {
choose = vector<vector<int>>();
choose.push_back(vector<int>{1}); // choose[0]
for (int n = 1; n < N; ++n) {
choose.push_back(vector<int>(n+1));
choose[n][n] = choose[n][0] = 1;
for (int r = 1; r < n; ++r) {
// out of the previous n - 1 things choose r things
// and out of the previous n - 1 things choose r - 1 things and choose the nth thing
choose[n][r] = choose[n-1][r] + choose[n-1][r-1];
choose[n][r] %= MOD;
}
}
return choose;
}
int computePaths(int R, int C, vector<pair<int, int>> specialFields) {
// specialFields should already be sorted in this case
specialFields.emplace_back(R, C);
int N = specialFields.size();
vector<int> paths = vector<int>(N);
for (int i = 0; i < N; ++i) {
paths[i] = choose[specialFields[i].first + specialFields[i].second][specialFields[i].first];
for (int j = 0; j < i; ++j) {
int rDiff = specialFields[i].first - specialFields[j].first;
int cDiff = specialFields[i].second - specialFields[j].second;
assert(rDiff >= 0);
if (cDiff >= 0) {
paths[i] -= (((long long) paths[j])*choose[rDiff + cDiff][rDiff]) % MOD;
paths[i] %= MOD;
if (paths[i] < 0) paths[i] += MOD;
}
}
}
return paths.back();
}
vector<tuple<int, int, int>> makeSpecialFields(int r, int c, vector<int> fieldrow, vector<int> fieldcol) {
vector<tuple<int, int, int>> specialFields;
bool originIsSpecial = false;
bool endIsSpecial = false;
int N = fieldrow.size();
for (int i = 0; i < N; ++i) {
if (fieldrow[i] == 1 && fieldcol[i] == 1) originIsSpecial = true;
if (fieldrow[i] == r && fieldcol[i] == c) endIsSpecial = true;
specialFields.emplace_back(fieldrow[i] - 1, fieldcol[i] - 1, i);
}
if (!originIsSpecial) specialFields.emplace_back(0, 0, -1);
if (!endIsSpecial) specialFields.emplace_back(r - 1, c - 1, N);
sort(specialFields.begin(), specialFields.end(),
[](tuple<int, int, int> a, tuple<int, int, int> b) {
int ra = get<0>(a);
int rb = get<0>(b);
int ca = get<1>(a);
int cb = get<1>(b);
return ra < rb || (ra == rb && ca < cb) || (ra == rb && ca == cb && get<2>(a) < get<2>(b));
});
return specialFields;
}
public:
vector<int> difPaths(int r, int c, vector<int> fieldrow, vector<int> fieldcol) {
initializeChoose(r + c - 1);
// make a directed graph of paths between without going through any other fields
// tuple is (row, col, idx)
vector<tuple<int, int, int>> specialFields = makeSpecialFields(r, c, fieldrow, fieldcol);
int N = specialFields.size();
// graph[i][j] is the number of paths going from i to j without going through any other special fields
vector<vector<int>> graph(N, vector<int>(N, 0));
for (int i = 0; i < N - 1; ++i) {
for (int j = i + 1; j < N; ++j) {
int rDiff = get<0>(specialFields[j]) - get<0>(specialFields[i]);
int cDiff = get<1>(specialFields[j]) - get<1>(specialFields[i]);
vector<pair<int, int>> intermediateFields;
assert(rDiff >= 0); // since fields are sorted
if (cDiff >= 0 && get<2>(specialFields[i]) < get<2>(specialFields[j])) {
for (int k = i + 1; k < j; ++k) {
assert(get<0>(specialFields[i]) <= get<0>(specialFields[k]) && get<0>(specialFields[k]) <= get<0>(specialFields[j])); // since fields are sorted
if (get<1>(specialFields[i]) <= get<1>(specialFields[k]) && get<1>(specialFields[k]) <= get<1>(specialFields[j])) {
intermediateFields.emplace_back(get<0>(specialFields[k]) - get<0>(specialFields[i]),
get<1>(specialFields[k]) - get<1>(specialFields[i]));
}
}
graph[i][j] = computePaths(rDiff, cDiff, intermediateFields);
}
}
}
// first index refers to last special index (index in specialFields)
// second index refers to length of path
vector<vector<int>> pathsMemo(N, vector<int>(N + 1, 0));
pathsMemo[0][1] = 1; // counting the origin as a special field, we have 1 path going through the origin that goes through 1 special field (the origin)
for (int l = 2; l <= N; ++l) {
for (int k = 0; k < N; ++k) {
for (int j = 0; j < k; ++j) {
pathsMemo[k][l] += ((long long) pathsMemo[j][l-1]*graph[j][k]) % MOD;
pathsMemo[k][l] %= MOD;
}
}
}
return vector<int>(pathsMemo[N-1].end() - fieldrow.size() - 1, pathsMemo[N-1].end());
}
};