Programmer’s Challenge

By Bob Boonstra, Westford, Massachusetts

Mutant Life

Time for a little nostalgia this month. Most of you probably remember John Conway’s exploration of cellular automata known as the game of Life. The game is played on a grid of square cells. A cell has one of two states - it can be occupied (“alive”) or empty (“dead”). Time proceeds in discrete increments, or generations, and the state of a cell at time N+1 is determined by its state and that of its eight neighbors at time N. In the simplest variations of the game, a “birth” occurs in an empty cell if exactly three of its neighbors were alive in the previous generation. A “death” occurs in an occupied cell surrounded by four or more living cells, or by fewer than two living cells.

This month, the challenge is to write code that will compute the state of a Life-like world some number of generations into the future. The prototype for the code you should write is:

pascal long PropagateLife(
 BitMap cells,   /* the boundaries and population of your automata */
 long numGenerations,/* number of generations to propagate */
 short birthRules, /* defines when cells become alive */
 short deathRules/* defines when cells die */
);

Your automata live in a world defined by the rectangle cells.bounds (with top and left coordinates guaranteed to be 0). Their world is actually a torus instead of a rectangle: the cells.bounds.right-1 column of cells is adjacent to column 0, and the cells.bounds.bottom-1 row of cells is adjacent to row 0. The rules for birth and death are generalized from those in the first paragraph and defined by birthRules and deathRules. An empty cell with X occupied neighbors becomes alive in the next generation if the bit (birthRules & (1<<X)) is set. An occupied cell with Y occupied neighbors dies in the next generation if the bit (deathRules & (1<<Y)) is set. Any other cell retains its previous state (occupied or empty) from one generation to the next. As an example, the version of the game described in the first paragraph would have birthRules=0x0008 and deathRules=0x01F3.

The initial population of automata is pointed to by cells.baseAddr, one bit per cell, when PropagateLife is called. An occupied cell has the value 1, and an empty cell has the value 0. The cells BitMap is defined in the usual way, with row R found starting at *(cells.baseAddr + R*cells.rowBytes). You are to use birthRules and deathRules to propagate this population ahead for numGenerations generations, stopping only in the event that the population of generation N is identical to that of the immediately preceeding generation. Your code must return the number of generations processed (which will be numGenerations unless a static population was reached). When you return, the memory pointed to by cells.baseAddr must contain the propagated population.

You may allocate a reasonable amount of auxiliary storage if that is helpful, provided (as always) that you deallocate any memory before returning, as I will be calling your code many times.

This month, we continue the language experiment that permits your solution to the Challenge to be coded in C, C++, or Pascal, using your choice among the MPW, Metrowerks, or Symantec compilers for these languages. The environment you choose must support linking your solution with test code written in C. Along with your solution, you should provide a project file or make file that will generate a stand-alone application that calls your solution from C test code.

This will be a native PowerPC Challenge. Now, start propagating

Two Months Ago Winner

Congratulations to Ernst Munter (Kanata, Ontario) for submitting the fastest entry to the Intersecting Rectangles Challenge. Of the eighteen contestants who submitted entries, sixteen provided correct solutions. Recall that the Challenge was to provide code that would return a set of output rectangles containing all points inside in an odd number (or an even number, depending on an input parameter) of input rectangles.

A number of solutions scanned the list of input rectangles and created a list of rectangles formed by the intersections, keeping track of whether the resulting subrectangles were inside an odd or an even number of input rectangles. Other solutions used a bitmap approach, calculating the exclusive OR of the input rectangles (for the odd parity case). The bitmap technique tended to suffer when the rectangles spanned a large x/y space.

The winning solution combines these techniques in an interesting way. Ernst first scans the input rectangles to collect and sort the unique x and y vertex coordinates. He then forms a reduced-scale bitmap using these virtual pixels (dubbed “vixels”), applying the XOR technique to compute the odd or even parity intersections of the input rectangles. Finally, Ernst scans the “vixelMap” to form output rectangles of the appropriate parity. An innovative technique that was not only fast but also space-efficient compared with many of the other entries.

The table below summarizes the results for entries that worked correctly. It shows the total time required for 60 test cases of up to 250 input rectangles per test case, the number of output rectangles produced, and the total code/data size of each entry. (The limit of 250 input rectangles resulted from the large memory requirements of some of the solutions.) Numbers in parentheses after a person’s name indicate that person’s cumulative point total for all previous Challenges, not including this one.

Top 20 Contestants of All Time

Here are the Top Contestants for the Programmer’s Challenges to date, including everyone who has accumulated more than 20 points. The numbers below include points awarded for this month’s entrants.

1. [Name deleted] 176 11. Mallett, Jeff 44

2. Munter, Ernst 132 12. Kasparian, Raffi 42

3. Gregg, Xan 92 13. Vineyard, Jeremy 42

4. Larsson, Gustav 87 14. Lengyel, Eric 40

5. Karsh, Bill 80 15. Darrah, Dave 31

6. Stenger, Allen 65 16. Brown, Jorg 30

7. Riha, Stepan 51 17. Landry, Larry 29

8. Cutts, Kevin 50 18. Elwertowski, Tom 24

9. Goebel, James 49 19. Lee, Johnny 22

10. Nepsund, Ronald 47 20. Noll, Robert 22

There are three ways to earn points: (1) scoring in the top 5 of any Challenge, (2) being the first person to find a bug in a published winning solution or, (3) being the first person to suggest a Challenge that I use. The points you can win are:

1st place 20 points 5th place 2 points

2nd place 10 points finding bug 2 points

3rd place 7 points suggesting Challenge 2 points

4th place 4 points

Xan Gregg earns two points this month for being the first to point out an error in the winning Find Again and Again solution by Gustav Larsson published in the February issue. The error occurs because the routines BMH_Search() and SimpleSearch() use signed declarations char * when they ought to use unsigned char *. As a result, processing is not correct in some cases when the textToSearch contains characters >= 0x80. There was confusion on this point in a number of the entries, and I did not penalize any of the solutions for making this error.

Here is Ernst’s winning Intersecting Rectangles solution:

IntersectRects.c

Copyright 1996, Ernst Munter, Kanata, ON, Canada

/*
  The Problem
  -----------
  Given a bunch of overlapping rectangles, compute a set
  of rectangles which covers the area of either an odd or
  an even number of overlaps.  The output rects should only
  use edges from the repertoire of edges contained in the
  input set of rects.

  General Strategy
  ----------------
  We create a virtual raster with a (variable) resolution,
  where each x or y coordinate value corresponds to an
  edge of at least one input rectangle.  Depending on the
  number of input rects, and their coincidence of edges,
  this raster may be very small, or fairly large, but never
  larger than the screen it represents.

  We then paint rectangles into the raster, each raster
  point being represented by 1 bit, regardless how many
  pixels are within the corresponding edges on the real
  screen.  I call these bits “virtual pixels” or “vixels”.

  After all vixels are painted, the bit map is scanned
  to identify rectangular areas of set bits.

  The vertical extent of each output rect is at least equal
  to the distance between the two neighboring input edges.

  We then follow the slice down over as many slices as
  possible to maximize the height of the rectangle.

  Memory Use
  ----------
  The maximum amount of memory allocated dynamically is
  determined by the number of input rects.  The actual
  amount will be less if some input rects share edge
  coordinate values.

  Approximate size of the index heap:
  (16 * numRectsIn)  bytes
     plus a few overhead bytes

  Approximate combined size of the two vixel maps:
  (numRectsIn * numRectsIn) bytes
     plus a few overhead bytes,
     minus gain from elimination of duplicate values

  A double size vixel map is always allocated although
  only the even parity case needs both.

  For example, total dynamic memory for 100 rectangles will
  be about 16K.  1000 rectangles might need 1MB, but on
  any reasonable size screen, 1000 rectangles will share
  a very large number of edges, and will have considerably
  less memory allocated.


  Other assumptions (these are not checked)
  -----------------------------------------
  There is at least one input rect.

  All input rects are legal and not empty, that is:
      top<bottom, and left<right.
*/

#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#define MAXLONG 0x7fffffff

void RectangleIntersections(
     const Rect inputRects[],
     const long numRectsIn,
     Rect  outputRects[],
     long  *numRectsOut,
     const Boolean oddParity);

// Local function prototypes:
void  PaintOdd(long* vm,long H,long L,long R,long mapWidth);
void  PaintEven(long* vm,long H,long L,long R,long mapWidth);
void  PackMap(long* vm,long mapSize);
void  Insert(long* h,long size,long x);
long* Sort(long* h,long size);
long  GetIndex(long size,long* index,long z);

//Some shorthand macros:
#define IRT (inputRects[i].top)
#define IRL (inputRects[i].left)
#define IRB (inputRects[i].bottom)
#define IRR (inputRects[i].right)

#define ORT (ORptr->top)
#define ORL (ORptr->left)
#define ORB (ORptr->bottom)
#define ORR (ORptr->right)

/* Masks needed to process the edges of vixel blocks
   which are not necessarily aligned with bitmap words.
*/
long leftMask[32] =
   {0xFFFFFFFF, 0x7FFFFFFF, 0x3FFFFFFF, 0x1FFFFFFF,
    0x0FFFFFFF, 0x07FFFFFF, 0x03FFFFFF, 0x01FFFFFF,
    0x00FFFFFF, 0x007FFFFF, 0x003FFFFF, 0x001FFFFF,
    0x000FFFFF, 0x0007FFFF, 0x0003FFFF, 0x0001FFFF,
    0x0000FFFF, 0x00007FFF, 0x00003FFF, 0x00001FFF,
    0x00000FFF, 0x000007FF, 0x000003FF, 0x000001FF,
    0x000000FF, 0x0000007F, 0x0000003F, 0x0000001F,
    0x0000000F, 0x00000007, 0x00000003, 0x00000001};

long rightMask[32] =
   {0x80000000, 0xC0000000, 0xE0000000, 0xF0000000,
    0xF8000000, 0xFC000000, 0xFE000000, 0xFF000000,
    0xFF800000, 0xFFC00000, 0xFFE00000, 0xFFF00000,
    0xFFF80000, 0xFFFC0000, 0xFFFE0000, 0xFFFF0000,
    0xFFFF8000, 0xFFFFC000, 0xFFFFE000, 0xFFFFF000,
    0xFFFFF800, 0xFFFFFC00, 0xFFFFFE00, 0xFFFFFF00,
    0xFFFFFF80, 0xFFFFFFC0, 0xFFFFFFE0, 0xFFFFFFF0,
    0xFFFFFFF8, 0xFFFFFFFC, 0xFFFFFFFE, 0xFFFFFFFF};


RectangleIntersections
void RectangleIntersections(
     const Rect inputRects[],
     const long numRectsIn,
     Rect  outputRects[],
     long  *numRectsOut,
     const Boolean oddParity) {

    long* xHeap;
    long* yHeap;
    long* vixelMap;

    long* xIndex;
    long* yIndex;
    long  xIndexMax;
    long  yIndexMax;
    long  xHeapSize;
    long  yHeapSize;

    long  i;
    long  mapWidth;
    long  mapSize;
    Rect* ORptr = outputRects;

//  First, we collect all X and Y coordinate values of
//  all input rectangles in a heap (priority queue), which
//  is then sorted into an index without duplicates for each
//  dimension, using a modified form of Heapsort.

    *numRectsOut=0;
    if (0==(yHeap=(long*)malloc((numRectsIn+3)*sizeof(long)*4)))
      return;
    xHeap=yHeap+(numRectsIn+3)*2;

    *xHeap=*yHeap=MAXLONG;
    xHeapSize=yHeapSize=1;

    for (i=0;i<numRectsIn;i++) {
      Insert(yHeap,yHeapSize,IRT);    yHeapSize++;
      Insert(yHeap,yHeapSize,IRB);    yHeapSize++;
      Insert(xHeap,xHeapSize,IRL);    xHeapSize++;
      Insert(xHeap,xHeapSize,IRR);    xHeapSize++;
    }

    xIndex=Sort(xHeap,xHeapSize);
    xIndexMax=xHeapSize-(xIndex-xHeap);
    yIndex=Sort(yHeap,yHeapSize);
    yIndexMax=yHeapSize-(yIndex-yHeap);

//note: IndexMax indexes to the last entry index[indexMax]
//      in each index list.  index[0] and index[indexMax]
//      are the edges of the enclosing rectangle.

//  Each block of real pixels that is defined by adjacent
//  X and Y edges (from any input rectangle) is considered
//  as a single virtual pixel (“vixel”).  The map of these
//  vixels will then be populated by the input rectangles.
//  Each vixel is represented by a bit in vixelMap.

//  We get memory for the vixel map and clear it to 0.
//  Vixels are stored as bitmaps in 32-bit words.

//  The vixel map is initially organized as either 1 word
//  per 32 vixels (odd parity) or 2 words (even parity).

    mapWidth=(32+xIndexMax) >> 5;
    mapSize=mapWidth*(yIndexMax+1);
    if (0==(vixelMap=
       (long*)malloc(2*mapSize*sizeof(long)))) return;
    if (oddParity) memset(vixelMap,0,mapSize*sizeof(long));
    else memset(vixelMap,0,2*mapSize*sizeof(long));

//  With odd parity, it is only necessary to XOR the vixels
//  of all input rects (PaintOdd).

//  With even parity, we also need to OR all vixels.  This
//  is done in the alternate words of vixelMap (PaintEven);
//  hence the vixelMap is stretched with even parity.

//  Accumulate the enclosed vixels of each input rect:

    for (i=0;i<numRectsIn;i++) {
      long  T,L,B,R,x,y;
      long* vm;
      T=GetIndex(yIndexMax,yIndex,IRT);
      L=GetIndex(xIndexMax,xIndex,IRL);
      B=GetIndex(yIndexMax,yIndex,IRB);
      R=GetIndex(xIndexMax,xIndex,IRR);
      if (oddParity) {
        vm=vixelMap+mapWidth*T+(L>>5);
        PaintOdd(vm,B-T,L,R-1,mapWidth);
      } else {
        vm=vixelMap+2*(mapWidth*T+(L>>5));
        PaintEven(vm,B-T,L,R-1,mapWidth);
      }
    }

//  For even parity, XOR all pairs of words in the vixelMap
//  to pack it into the same format as the odd parity
//  vixelMap.

    if (!oddParity) PackMap(vixelMap,mapSize);

//  Now the vixelMap (the bitmap of all vixels, that is
//  areas of the screen), is set to 1 for every vixel
//  meeting the criteria of either odd or even parity.

//  We scan the vixel map to find contiguous regions of
//  non-zero vixels in order to generate the output
//  rectangles.  For each row, we successively look for
//  blocks of set vixels. This will define one output rect.
//  The X/Y index arrays serve to convert the vixel
//  coordinates back to the real pixel coordinates which
//  define the output rectangles.

  { long word,bit,bb,c,u,L,B;
    long* vm=vixelMap;
    for (i=0;i<yIndexMax;i++) {
      bit=0;
      c=0;
      for (word=0;word<mapWidth;word++) {
        u=vm[word];c=0;
        if (u) {
          long* vmBelow=vm+word+mapWidth;
          bb=0;
          do {
            while (u>0) {bb++;u<<=1;}
            if (c==0) {
              L=bb;
              ORL=xIndex[bit+L]; c--;
            } else {
              long* vmx=vmBelow;
              long mask=~(leftMask[L] & rightMask[bb-1]);
              B=i+1;
//Default: the rectangle is 1 vixel high.
//We try to extend rectangle down as far as possible:
              while (-1==(mask | *vmx)) {
                B++;*vmx &= mask;vmx+=mapWidth;
              }
              ORB=yIndex[B];
              ORR=xIndex[bit+bb];
              ORT=yIndex[i];
              ORptr++;
              c=0;
            }
            if (0==(u=(~u) & rightMask[31-bb])) break;
          } while(bb<32);
          if (c) {
            long* vmx=vmBelow;
            long mask=~leftMask[L];
            B=i+1;
            while (-1==(mask | *vmx)) {
              B++;*vmx &= mask;vmx+=mapWidth;
            }
            ORB=yIndex[B];
            ORR=xIndex[bit+32];
            ORT=yIndex[i];
            ORptr++;
          }
        }
        bit+=32;
      }
      vm+=mapWidth;
    }

  }
  free(yHeap);                  //free allocated memory
  free(vixelMap);
  *numRectsOut=ORptr-outputRects;
}

//////////////////////////////////////////////////////////////////
// Auxiliary functions called by RectangleIntersections: //
/////////////////////////////////////////////////////////////////

Insert
/* Insert grows a heap, that is a partially sorted balanced
   binary tree, where each node’s children must be less or
   equal, but not in any particular order.

   Each value x is inserted by appending it as the last node
   and then sifting it up (exchanging father and child
   nodes) until the heap property is restored.
*/
void Insert(long* h,long size,long x) {
  long i,j,z;
  i=size;
  do {
    j=i>>1;
    if (x<=(z=h[j])) break;
    h[i]=z;
    i=j;
  } while(1);
  h[i]=x;
}

Sort
/* The heap keeps the largest value at the root, at h[1].
   We sort as follows: each root value is removed and put
   at the end of the array;  then the last item in the heap
   is put into the root and sifted down until the heap
   property is restored.
   When we are done, the array is sorted.
   As we go along, we recognize duplicate values and remove
   them but do not put them back.  The result is that the
   start of the sorted list may be further up in the array.
*/
long* Sort(long* h,long size) {
  long  x,z,i,j;
  long* b=h+size+1;
  *b=MAXLONG;
  if (size>1) do {
    size--;
    i=1;
    j=2;
    if (*b != (z=h[1])) *(--b) = z;
    if (size<=1) break;
    x=h[size];
    h[size]=-MAXLONG;
    while (j<size) {
      long h0=h[j],h1=h[1+j];
      if (h0<h1) {j++; h0=h1;}
      if (x>=h0) break;
      h[i]=h0;
      i=j;
      j+=j;
    }
    h[i]=x;
  } while(1);
  return b;
}

GetIndex
/* GetIndex uses a binary search to locate a particular
   entry and returns its index.
*/
long GetIndex(long r,long* index,long z) {
  long l=0,m=r>>1,y;
  do {
    if (z>(y=index[m])) l=m+1;
    else if (z<y)       r=m-1;
    else return m;
    m=(l+r)>>1;
  } while (l<r);
  return r;
}


PaintOdd
/* The PaintOdd and PaintEven routines paint rectangles
   into the vixel map.
   PaintOdd only XORs a single bit map with a rectangle.
   PaintEven also ORs a second bit map with the same
   rectangle.  The 2 bit maps are word interleaved.
   It is hoped that this reduces cache misses by keeping
   to one area of memory for each row of a rectangle.
*/
void  PaintOdd(long* vm,long H,long L,long R,long mapWidth) {
long LM=leftMask[L & 31];
long RM=rightMask[R & 31];
long numMid=(>>R5)-(L>>5)-1;
long x,y,pad=mapWidth-numMid-2;

  if (numMid<0) {LM&=RM;RM=0;}

  for (y=0;y<H;y++) {
    *vm ^= LM; vm++;
    for (x=0;x<numMid;x++) {
      *vm ^= 0xFFFFFFFF; vm++;
    }
    if (RM) {
      *vm ^= RM; vm++;
    }
    vm+=pad;
  }
}

PaintEven 
void  PaintEven (long* vm,long H,long L,long R,long mapWidth) {
long LM=leftMask[L & 31];
long RM=rightMask[R & 31];
long numMid=(>>R5)-(L>>5)-1;
long x,y,pad=(mapWidth-numMid-2)<<1;

  if (numMid<0) {LM&=RM;RM=0;}

  for (y=0;y<H;y++) {
    *vm ^= LM; vm++;
    *vm |= LM; vm++;
    for (x=0;x<numMid;x++) {
      *vm ^= 0xFFFFFFFF; vm++;
      *vm |= 0xFFFFFFFF; vm++;
    }
    if (RM) {
      *vm ^= RM; vm++;
      *vm |= RM; vm++;
    }
    vm+=pad;
  }
}

PackMap
/* PackMap reduces the two interleaved bit maps used for
   the even parity case, into a single bit map.  Each
   pair of words, of the entire bitmap, is XORed together
   regardless of rectangle boundaries.
*/
void PackMap(long* vm,long mapSize) {
long* vmE=vm;
long* endOfMap=vm+mapSize;
  while (vm<endOfMap) {
    *vm++ = *vmE ^ vmE[1];
    vmE+=2;
  }
}