USACO DEC07 Problem 'mud' Analysis

by Rob Kolstad

This is a classic 'Flood fill' algorithm. The idea is that one seeds a 'flood fill' with an initial point followed by massive recursion to find the points around that point.

Using a small 5 x 5 matrix as an example:

 4 . . . . .
 3 . . . . .
 2 . . * . .
 1 . . . . .
 0 . . . . .
   0 1 2 3 4
one can easily imagine starting at point 2,2 and then examining points around it after marking those points with the 'distance' from the initial point:
 4 . . . . .
 3 . . 1 . .
 2 . 1 * 1 .
 1 . . 1 . .
 0 . . . . .
   0 1 2 3 4
and then examining points around those 1's:
 4 . . 2 . .
 3 . 2 1 2 .
 2 2 1 * 1 2
 1 . 2 1 2 .
 0 . . 2 . .
   0 1 2 3 4
and so on until all distances are determined:
 4|4 3 2 3 4
 3|3 2 1 2 3
 2|2 1 * 1 2
 1|3 2 1 2 3
 0|4 3 2 3 4
  +---------
   0 1 2 3 4
Adding in obstacles just makes some of the distances a bit different since certain cells can not be used in the recursive search. Here is the flooded version of the sample input data:
   4| 6  7  8  9 10 11 10  9  
   3| 5  M  9 10 11 10  9  8 
Y  2| 4  5  M  B  M 11  M  7 
   1| 3  M  1  M  3  M  5  6 
   0| 2  1  *  1  2  3  4  5 
  -1| 3  2  1  2  3  4  5  6 
    +-----------------------
     -2 -1  0  1  2  3  4  5 
                 X
The final '11' has not been written on to the 'B' square on which Bessie sits, but it's clear that 11 is the answer.

The obvious way to implement this algorithm is so simple (yet wrong):

check(x,y,dist) {
    if (dist > matrix[x,y]) return;
    matrix[x,y] = dist;
    if (validsubscript (x+1, y  ) check (x+1, y  , dist+1);
    if (validsubscript (x-1, y  ) check (x-1, y  , dist+1);
    if (validsubscript (x  , y+1) check (x  , y+1, dist+1);
    if (validsubscript (x  , y-1) check (x  , y-1, dist+1);
}
The pleasing intuitive brevity, though, is nothing but a trap. Why? Because, as shown, the subprogram will quickly recurse to check(x+1, ...), check (x+2, ...) and so on down the x axis instead of the geometrically nice circular pattern exhibited a few paragraphs back. Will the algorithm still get the right answer? Yes! Is there a way to implement the algorithm that will take longer to get to the right answer? Probably not.

Instead, one needs to implement a FIFO queue structure so that all four of the 'next elements to search' are queued before any of them is examined. Bonus performance for knowing if an element is already queued (and thus avoiding requeueing it). Note that one should requeue elements whose value has changed -- but it's hardly ever an advantage to have the same element on the queue more than once.

Aside: queues can be nicely implemented with linked lists. Such lists are slightly tricky to get 'right the first time' and sometimes are implemented to allocate a tiny bit of storage for each element added to the list. In contest programming, one might as well just make a big array and circle through it. If you need more space than the array allows, you would have needed more linked-list space -- and failed very late in the program's execution. Better to expose the potential failure up front!

The final challenge in this task was the coordinate scheme: it starts at some large negative number, so all subscripts need to be shifted into a valid range. The first program below does this with #define's; the second explicitly codes the shift into each subscript.

My solution is intended to be straightforward with the tiny trick of using a bit more storage to set up a 'border' around the edges of the matrix. This avoids having to validate every possible enqueueing to ensure the subscripts will not go out of bounds for the array:

#include <stdio.h>

#define RANGE 500
#define FIELD(i,j) field[(i)+RANGE+2][(j)+RANGE+2]
#define DONE(i,j) done[(i)+RANGE+2][(j)+RANGE+2]
#define QUEUED(i,j) queued[(i)+RANGE+2][(j)+RANGE+2]

short field[2*RANGE+5][2*RANGE+5];
char done[2*RANGE+5][2*RANGE+5], queued[2*RANGE+5][2*RANGE+5];	

struct q { int x, y; } queue[800000];
int nq = 0, head = 0; tail = 0;

tryenq (x,y,thisdist) {
    if (FIELD(x,y) != 0 && thisdist >= FIELD(x,y)) return;
    if (QUEUED(x,y)) return;
    QUEUED(x,y) = 1;
    FIELD(x,y) = thisdist;
    if (++nq >= 790000) exit (9);
    queue[tail].x = x;
    queue[tail].y = y;
    if (++tail >= 800000) tail = 0;
}

main () {
    FILE *fin = fopen ("mud.in", "r");
    FILE *fout = fopen ("mud.out", "w");
    int X, Y, N, i, j, A, B;

    fscanf (fin, "%d %d %d", &X, &Y, &N);
    for (i = 0; i < N; i++) {
	fscanf (fin, "%d %d", &A, &B);
	DONE(A,B) = -1;
    }
	/* set up a border so we don't look outside box: */
    for (i = -RANGE; i <= RANGE; i++)
	FIELD(-RANGE-1,i) = FIELD(RANGE+1,i) = FIELD(i,RANGE+1) =
					FIELD(i,-RANGE-1) = -1;
    tryenq (0,0,0);
    while (nq) {
	int x=queue[head].x; 		/* dequeue */
	int y=queue[head].y;
	QUEUED(x,y) = 0;
	nq--;
	if (++head >= 800000) head = 0;

	/* enqueue if needed: */
	if (DONE(x,y)) continue;
	DONE(x,y) = 1;
	tryenq(x+1, y  , FIELD(x,y)+1);
	tryenq(x-1, y  , FIELD(x,y)+1);
	tryenq(x  , y+1, FIELD(x,y)+1);
	tryenq(x  , y-1, FIELD(x,y)+1);
    }

    fprintf (fout, "%d\n", FIELD(X,Y));
    exit (0);
}

Richard Peng's solution is a bit shorter (as usual). Among the techniques for compaction:

Richard's program is over 1.5x faster than the program above.
#include <cstdio>
#include <cstring>

short vis[1020][1020], r, c, shif, q[1100000][2];
int a, b, x, y, n, tai;
int dir[4][2]={-1, 0, 1, 0, 0, -1, 0, 1};

void add (int r1, int c1, int d) {
    if (r1>=0 && r1<r && c1>=0 && c1<c && vis[r1][c1]==-1) {
        vis[r1][c1]=d;
        q[tai][0]=r1;
        q[tai][1]=c1;
        tai++;
    }
}

int main () {
    int i, j;
    freopen ("mud.in", "r", stdin);
    freopen ("mud.out", "w", stdout);
    r = 1010;  c = 1010;  shif = 505;
    tai = 0;
    scanf ("%d%d%d", &x, &y, &n);
    for (i = 0; i < r; i++)
        for (j = 0; j < c; j++)
            vis[i][j] = -1;
    for (i = 0; i < n; i++){
        scanf ("%d%d", &a, &b);
        vis[a+shif][b+shif] = -2;
    }
    tai = 0;
    add (shif, shif, 0);
    for (i = 0; i < tai; i++)
        for (j = 0; j < 4; j++)
            add (q[i][0]+dir[j][0], q[i][1]+dir[j][1], vis[q[i][0]][q[i][1]]+1);
    printf ("%d\n", vis[x+shif][y+shif]);
    return 0;
}