Volume Number: 16 (2000)
Issue Number: 5
Column Tag: Programmer's Challenge

Programmer's Challenge

by Bob Boonstra, Westford, MA

Back in September, 1995, we conducted an RSA Challenge that involved raising large integers to integral powers, modulo a third integer. The representation we used for those large integers was a BigNum type, where each digit of the large integer was stored in a byte. That representation and the operations on it were not particularly efficient, and this month we will belatedly recitfy that situation. Your Challenge is to implement a new BigNum type, of your own design, along with a number of arithmetic operations on these BigNums..

The prototype for the code you should write is:

typedef struct BigNum {
	long lengthInDigits;	/* length of the BigNum in digits */
	void *bigNumData;			/* pointer to BigNum data */
} BigNum;

BigNum NewBigNum (			/* create a BigNum */
	char sign,						/* +1 or -1 */
	char digits[],				/* digits to be made into a BigNum */
	long numDigits				/* number of digits */
);

void DisposeBigNum (		/* dispose of a BigNum */
	BigNum theBigNum			/* the BigNum to be disposed of */
);

BigNum AddBigNums (			/* sum two BigNums, returning a new one */
	BigNum bigNumA,				/* return the sum A+B */
	BigNum bigNumB
);

BigNum SubtractBigNums (	/* subtract two BigNums, returning a new one */
	BigNum bigNumA,				/* return the difference A-B */
	BigNum bigNumB
);

BigNum MultiplyBigNums (	/* multiply two BigNums, returning a new one */
	BigNum bigNumA,				/* return the product A*B */
	BigNum bigNumB
);

BigNum DivideBigNums (		/* divide two BigNums, returning a new one */
	BigNum bigNumA,				/* return the quotient A/B, discarding the remainder */
	BigNum bigNumB
);

BigNum ModBigNums (			/* divide two BigNums, returning a new one */
	BigNum bigNumA,				/* return the remainder A%B, discarding the quotient */
	BigNum bigNumB
);

BigNum PowerBigNums (		/* calculate one Bignum to the power of another, returning a new one */
	BigNum bigNumA,				/* return A raised to the power B, discarding the quotient */
	BigNum bigNumB
);

BigNum SqrtBigNum (			/* find the sqrt of a BigNum, returning a new one */
	BigNum bigNumA				/* return the square root of A */
);

long /* numDigits */ BigNumToDigits( /* convert a bigNum to decimal digits */
	BigNum bigNumA,				/* bigNum to be converted to decimal digits 0-9 */
	char *sign,						/* return +1 or -1 */
	char digits[]					/* decimal digits of bigNumA, preceeded by '-' if negative */
									/* storage for digits preallocated based on bigNumA.lengthInDigits */
);

The first thing you need to do is decide on an internal representation for BigNums. Then you need to write a NewBigNum routine that will create a BigNum from a sequence of numDigits digits and a sign value. Your NewBigNum code is responsible for allocating memory for the BigNumData. The DisposeBigNum routine is responsible for deallocating that memory. The caller of your code is responsible for pairing every NewBigNum call with a DisposeBigNum call, and the two routines should be implemented so as not to create any memory leaks. In addition to these allocation and deallocation routines, you need to write code to perform addition (AddBigNums), subtraction (SubtractBigNums), multiplication (MultiplyBigNums), division (DivideBigNums), remainders (ModBigNums), and exponentiation (PowerBigNums). Each of these routines takes two arguments, calculates the result, and returns the result in a new BigNum allocated by your code. Each of these returned BigNums will also be disposed of by a call to DisposeBigNum before the test is over, although they might be used for calculations in the interim.

Just to spice things up, you also need to provide a SqrtBigNum routine that calculates and returns the integer square root of a BigNum, the largest BigNum whose square is no larger than the original number.

And finally, to help me decipher your BigNums, you need to provide a BigNumToDigits conversion routine that converts your private BigNum data structure into a sequence of digits, along with a sign, and returns the number of digits in the decimal representation of the BigNum.

I'm not providing information on the distribution of calls to the various routines, except to say that the arithmetic routines will significantly outnumber the allocation and deallocation routines. The winner will be the solution that correctly completes a sequence of arithmetic operations on BigNums in the least amount of time. You are strongly encouraged to adequately comment the code in your submissions. Not only does that make your code more understandable if it is published as the winning solution, but it also helps me track down any minor problems that might occur.

I'll close with a plug for the Challenge mailing list, where you can receive notice of the problems before the hard-copy magazine reaches your mailbox, and where any post-publication clarifications are distributed. Subscription instructions can be found at www.mactech.com/progchallenge/. This will be a native PowerPC Challenge, using the CodeWarrior Pro 5 environment. Solutions may be coded in C, C++, or Pascal.

The February Challenge required readers to calculate a minimal Latin square of a given order. Latin Squares are nxn arrays of integers, where each row and each column contains each integer from 1 to n exactly once. Congratulations to Willeke Rieken (The Netherlands) for coming up with the winning solution to the Latin Squares Challenge.

Eleven readers submitted entries to this Challenge, and their performance varied widely in efficiency. My test scenario was based on 28 test cases, consisting of the Latin Squares of orders 4, 5, 8, 9, 12, 13, 16, 17, 20, 21, 24, 25, 28, 29, 32, 33, 36, 37, 40, 41, 44, and 45. I selected those numbers because they formed a regular pattern that could be continued as far as the solutions would allow, and because they contained a mix of odd numbers, even numbers, perfect squares, prime numbers, and powers of two. My original intent was to test even larger numbers, but even the best solutions took too long to calculate some of the larger numbers.

Even limiting the tests to these cases, some of the solutions took a long time to execute, so I divided the tests into three sets. The first set consisted of the first ten test cases, and I ran all of the entries against that set. Three of the entries either did not complete all of the cases, or calculated a Latin Square that was larger than the squares calculated by other solutions. Three of the entries had fast execution times for the ten cases, and one more had an execution time within roughly two orders of magnitude of the best ones. So I ran the top four solutions against the next six test cases. Two of the entries completed those cases correctly, so I ran those cases against the final six test cases. The second place solution by Ernst Munter was by far the faster of the two, but unfortunately, it did not compute the minimal solution for the square of order 37. Where Ernst calculated a solution that included the following as the 28th row:

 28 27 26 25 32 31 30 35 36 37 33 34 19 20 17 21 7 8 9 5 6 4 ...

... Willeke's entry produced the following smaller value:

 28 27 26 25 32 31 30 35 36 37 33 34 19 20 17 21 7 8 9 5 6 3 ...

I decided not to disqualify solutions that produced suboptimal Latin Squares, or that failed to produce a result in a reasonable time. Instead, I ranked solutions by how many test cases they were able to complete, then how many they completed correctly, and then in order of increasing execution time. The problem statement called for the use of execution time only for correct solutions, but I felt that it was fairest to allow solutions that produced a suboptimal result to compete based on how well they did.

Willeke's algorithm takes advantage of the fact that squares whose size is a power of two can be generated with a systematic pattern of switching pairs of numbers in row n to create rows a power of 2 away from row n. He accomplishes this in his FillSquare2 routine. Squares of other sizes are filled by first filing the largest subsquare of size k (k a power of 2), filling the top right n-k square optimally, filling the diagonal, and then completing the square by trial and error. Ernst's entry makes more efficient use of information about which digits are forced into use before a particular column in a given row because the digit has already been used in subsequent columns. Ernst observes in his entry that execution time does not grow with problem size, and that problems of certain sizes (e.g., 41) take much longer to execute than one might expect based on the time required for squares of dimensions close in value.

The first table below lists, for each of the entries submitted, the final ranking based on all test cases completed, total execution time for the first ten cases, the number of test cases completed, the number completed incorrectly, and the code size, data size, and language parameters. As usual, the number in parentheses after the entrant's name is the total number of Challenge points earned in all Challenges prior to this one. The second and third tables provide the results for the remaining twelve test cases.

Note that while the top four positions in this Challenge were won by four of our top contestants in the points standing (the fifth did not compete), there are a number of new names in the list of contestants. Keep trying, folks, I know from personal experience that it takes a while to become good at this, but it is possible to knock the leaders from their perches.

Top Contestants

Listed here are the Top Contestants for the Programmer's Challenge, including everyone who has accumulated 10 or more points during the past two years. The numbers below include points awarded over the 24 most recent contests, including points earned by this month's entrants.

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:

Here is Willeke's winning Latin Squares solution:

LatinSquares.cp
Copyright © 2000
Welleke Rieken

/*
	After generating several squares a pattern emerged.
	If n is even, every second row can be generated by
	switching pairs of numbers of the row above.
	If n can be divided by 4, every third and fourth
	row can be generated by switching squares of 2 by 2
	numbers of the 2 rows above.
	Example: n = 12 is generated by generating
	n = 3 and replacing every number by a square with
	n = 4.
	Other n's are generated by generating the biggest
	power of n that fits in the square and generating
	a square of n - 2^x at the top right. This square can
	be repeated to the bottom left till the first square ends.
	the numbers in the first column are in ascending order.
	the diagonal from top rigth to bottom left is filled with n.
	Example: n = 7
	1234567
	2143675
	3412756
	4567xxx
	5x7xxxx
	67xxxxx
	7xxxxxx
	The remaining numbers are generated by trial and error.
*/

#include "LatinSquares.h"

FillSquare2
static void FillSquare2(long n, short *latinSquare,
												long theDim,
												long theStartRow, long theStartCol,
												long theStartVal, long theNrOfRows)
// n is a power of 2. fill the first row with ascending numbers
// and switch them around to generate the other rows.
{
	short *aFrom1, *aTo1, *aFrom2, *aTo2;
	long 	aValue = theStartVal + 1, aRowsDone, aMultiple;
	short	*aStartSquare = latinSquare + (theStartRow * n) + 
			theStartCol;
	
	// fill first row
	aFrom1 = aStartSquare;
	for (long aCol = 0; aCol < theDim; aCol++)
	{
		*aFrom1  = aValue;
		aValue++;
		aFrom1++;
	}
	aRowsDone = 1;
	aMultiple = 1;
	while (aRowsDone < theNrOfRows)
	{
		for (long aRow = 0; aRow < aMultiple; aRow++)
		{
			if (aRow >= theNrOfRows)
				break;
			for (long anOffset = 0; anOffset < theDim; anOffset += 
							(aMultiple * 2))
			{
				aFrom2 = aStartSquare + (aRow * n) + anOffset;
				aFrom1 = aFrom2 + aMultiple;
		aTo1 = aStartSquare + ((aMultiple + aRow) * n) + anOffset;
				aTo2 = aTo1 + aMultiple;
				for (long aCol = 0; aCol < aMultiple; aCol++)
				{
					*aTo1 = *aFrom1;
					aFrom1++;
					aTo1++;
					*aTo2 = *aFrom2;
					aFrom2++;
					aTo2++;
				}
			}
		}
		aRowsDone += aMultiple;
		aMultiple <<= 1;
	}
}

CopySquare
static inline void CopySquare(long n, short *theFrom, short *theTo,
															long theDim)
{
// copy a square of size theDim from theFrom to theTo
	short *aFrom, *aTo;
	
	for (long aRow = 0; aRow < theDim; aRow++)
	{
		aFrom = theFrom + (aRow * n);
		aTo = theTo + (aRow * n);
		for (long aCol = 0; aCol < theDim; aCol++)
		{
			*aTo = *aFrom;
			aFrom++;
			aTo++;
		}
	}
}

CantFillRow
static short CantFillRow(long theDim, short *theValInCol,
													short *theValInRow, long theCol,
													long *theValue)
// check if there are numbers that can't be placed and if there
// are enough columns for the bigger numbers
{
	long	aGreaterPlacesNeeded = 0;
	short	aValOK = 0;
	for (long i = *theValue + 1; i < theDim; i++)
		if (!theValInRow[i])
		{
			aGreaterPlacesNeeded++;
			aValOK = 0;
			for (long j = theCol + 1; j < theDim; j++)
				if (!theValInCol[j * theDim + i])
				{
					aValOK = 1;
					break;
				}
			if (!aValOK)
			{
				*theValue = i - 1;
				return 1;
			}
		}
	for (long j = theCol + 1; j < theDim; j++)
	{
		aValOK = 0;
		for (long i = *theValue + 1; i < theDim; i++)
			if (!(theValInRow[i] || theValInCol[j * theDim + i]))
			{
				aValOK = 1;
				break;
			}
		if (aValOK)
			aGreaterPlacesNeeded-;
	}
	if (aGreaterPlacesNeeded > 0)
		return 1;
	return 0;	
}

CompleteSquare
static void CompleteSquare(long n, short *latinSquare,
						long theDim, long theSubDim,
						long theStartRow, long theStartCol,
						long theStartVal)
// fill remaining numbers by trial and error
{
	short	*aStartSquare = latinSquare +
												((theStartRow * n) << theSubDim) +
												(theStartCol << theSubDim);
	short	*aValInRow = new short[theDim];
	short	*aValInCol = new short[theDim * theDim];
	short	*aToBeFilled = new short[theDim * theDim];
	long	aRow, aCol, aValue, aSubDimvalue;
	short	*p, *q;

	aSubDimvalue = 1 << theSubDim;
	// fill left row and diagonal
	p = aStartSquare + ((n + theDim - 2) << theSubDim);
	q = aStartSquare + (n << theSubDim);
	for (aCol = 1; aCol < theDim; aCol++)
	{
		CopySquare(n, aStartSquare + ((theDim - 1) << theSubDim),
								p, aSubDimvalue);
		p += ((n - 1) << theSubDim);
		CopySquare(n, aStartSquare + (aCol << theSubDim),
								q, aSubDimvalue);
		q += (n << theSubDim);
	}
	// which numbers are used and which numbers have to be filled in
	for (aCol = 0; aCol < theDim * theDim; aCol++)
	{
		aValInCol[aCol] = 0;
		aToBeFilled[aCol] = 1;
	}
	for (aRow = 0; aRow < theDim; aRow++)
	{
		p = aStartSquare + ((aRow * n) << theSubDim);
		for (aCol = 0; aCol < theDim; aCol++)
		{
			if (*p)
			{
				aValue = aCol * theDim +
									(((*p - 1) >> theSubDim) - theStartVal);
				aValInCol[aValue] = 1;
				aToBeFilled[aValue] = 0;
			}
			p += (aSubDimvalue);
		}
	}

	// which numbers are in this row
	for (aValue = 0; aValue < theDim; aValue++)
		aValInRow[aValue] = 0;
	aValue = 0;
	aRow = 1;
	aCol = 0;
	p = aStartSquare + (n << theSubDim);
	while (1)
	{
		// find next place to ve filled
		while (*p)
		{
			aValInRow[((*p - 1) >> theSubDim) - theStartVal] = 1;
			aCol++;
			p += (aSubDimvalue);
			if (aCol >= theDim)
			{
				aCol = 0;
				aRow++;
				p = aStartSquare + ((aRow * n) << theSubDim);
				for (aValue = 0; aValue < theDim; aValue++)
					aValInRow[aValue] = 0;
				aValue = 0;
			}
		}
		// find next posible value
		while ((aValue < theDim) &&
						(aValInCol[aCol * theDim + aValue] ||
							aValInRow[aValue] ||
							CantFillRow(theDim, aValInCol, aValInRow,
													aCol, &aValue)))
			aValue++;
		if (aValue < theDim)
		{
			// place value
			aValInCol[aCol * theDim + aValue] = 1;
			aValInRow[aValue] = 1;
			CopySquare(n, aStartSquare + (aValue << theSubDim),
									p, aSubDimvalue);
			
			// next column
			aCol++;
			p += (aSubDimvalue);
			if (aCol >= theDim)
			{
				// next row
				aRow++;
				if (aRow < theDim)
				{
					p = aStartSquare + ((aRow * n) << theSubDim);
					for (aValue = 0; aValue < theDim; aValue++)
						aValInRow[aValue] = 0;
					for (aCol = 0; aCol < theDim; aCol++)
					{
						aValInRow[aValue] = 0;
						if (*p)
					aValInRow[((*p - 1) >> theSubDim) - theStartVal] = 1;
						p += (aSubDimvalue);
					}
					aCol = 0;
					p = aStartSquare + ((aRow * n) << theSubDim);
				}
				else
				{
					return;
				}
			}
			aValue = 0;
		}
		else
		{
			// undo
			aCol-;
			p -= (aSubDimvalue);
			aValue = ((*p - 1) >> theSubDim) - theStartVal;
			while (aCol >= 0 && !aToBeFilled[aCol * theDim + aValue])
			{
				aCol-;
				if (aCol >= 0)
				{
					p -= (aSubDimvalue);
					aValue = ((*p - 1) >> theSubDim) - theStartVal;
				}
			}
			if (aCol < 0)
			{
				aRow-;
				p = aStartSquare +
						(((aRow * n) + theDim - 1) << theSubDim);
				aCol = theDim - 1;
				for (aValue = 0; aValue < theDim; aValue++)
					aValInRow[aValue] = 1;
			}
			aValue = ((*p - 1) >> theSubDim) - theStartVal;
			*p = 0;
			aValInCol[aCol * theDim + aValue] = 0;
			aValInRow[aValue] = 0;
			aValue++;
		}
	}
	delete[] aValInCol;
	delete[] aValInRow;
	delete[] aToBeFilled;
}

FillSquare
static void FillSquare(long n, short *latinSquare,
						long theDim, long theSubDim,
						long theStartRow, long theStartCol,
						long theStartVal, long theNrOfRows)
// fill latin square
// if n can be divided by a power of 2,
// theSubDim is 2^x, theDim is n/(2^x)
{
	if (theDim == 1)	// n is a power of 2
		FillSquare2(n, latinSquare, 1 << theSubDim,
				theStartRow << theSubDim, theStartCol << theSubDim,
				theStartVal << theSubDim, theNrOfRows << theSubDim);
	else
	{
		long	aMaxPower2, aNrOfRows, aStartCol, aStartRow;
		short	*aStartSquare = latinSquare +
												((theStartRow * n) << theSubDim) +
															(theStartCol << theSubDim);
		aMaxPower2 = 1;
		while (aMaxPower2 <= theDim) aMaxPower2 <<= 1;
		aMaxPower2 >>= 1;
		// fill top left of the square with a square with n = 2^aMaxPower2
		FillSquare2(n, latinSquare, aMaxPower2 << theSubDim,
				theStartRow << theSubDim, theStartCol << theSubDim,
				theStartVal << theSubDim, aMaxPower2 << theSubDim);
		aNrOfRows = theDim - aMaxPower2;
		if (aNrOfRows > theNrOfRows) aNrOfRows = theNrOfRows;
		// fill top right of the square with a square with n = theDim - 2^aMaxPower2
		FillSquare(n, latinSquare, theDim - aMaxPower2, theSubDim,
								theStartRow, theStartCol + aMaxPower2,
								theStartVal + aMaxPower2, aNrOfRows);
		// copy the square from the top right along the diagonal to the bottom left
		aStartCol = aMaxPower2 - aNrOfRows;
		aStartRow = aNrOfRows;
		while (aStartCol >= 0 && aStartRow < theNrOfRows)
		{
			if (aStartRow + aNrOfRows > theNrOfRows)
				aNrOfRows = theNrOfRows - aStartRow;
			if (aNrOfRows > aStartCol && aStartCol > 0)
				aNrOfRows = aStartCol;
	for (long aRow = 0; aRow < (aNrOfRows << theSubDim); aRow++)
			{
				short	*aFrom = aStartSquare + (aRow * n) +
												(aMaxPower2 << theSubDim);
				short	*aTo = aStartSquare +
									(((aStartRow << theSubDim) + aRow) * n) +
											(aStartCol << theSubDim);
				for (long aCol = 0; aCol < ((theDim - aMaxPower2) << 
							theSubDim); aCol++)
				{
					*aTo = *aFrom;
					aFrom++;
					aTo++;
				}
			for (long aCol = ((aStartCol + (theDim - aMaxPower2)) <<
				theSubDim); aCol < (aMaxPower2 << theSubDim); aCol++)
				{
					*aTo = 0;
					aTo++;
				}
			}
			aStartCol -= (theDim - aMaxPower2);
			aStartRow += (theDim - aMaxPower2);
		}
		// generate the remaning numbers
		CompleteSquare(n, latinSquare, theDim, theSubDim,
										theStartRow, theStartCol, theStartVal);
	}
}

LatinSquares
void LatinSquares(
  short n, /* dimension of the latin square to be generated */
  short *latinSquare /* set latinSquare[c + r*n] to square value row r, col c */
) {
	short	*p = latinSquare;
	long	aSubDim = 0;
	// init
	for (long i = 0; i < n * n; i++, p++)
		*p = 0;
	// can n be divided by a power of 2
	while (!(n & (1 << aSubDim))) aSubDim++;
	FillSquare(n, latinSquare, n >> aSubDim, aSubDim,
							0, 0, 0, n >> aSubDim);
}