## Google CodeJam: Numbers

The "Numbers" problem in Round 1A is considerably more devious than either Minimum Scalar Product, or Milkshakes problems. The problem is stated simply but the solution requires some innovative thinking to get around a precision problem. The fact that the exponent *n* can be up 2 x 10^{9} and √5 is irrational means that sufficient precision will be lost resulting in incorrect answers, or a stack overflow.

How do we know √5 is irrational? It can easily be proved using *Proof by Contradiction*.... Assume √5 is rational. In this case, by definition, it can be expressed as a quotient of 2 integers:

From number theory, we know that the **fundamental theorem of arithmetic** says that every positive composite integer has a unique prime factorization. This theorem also implies that you can't have both an even and odd number of prime factors. It is for this very reason that the theorem is also referred to as the **Unique Factorization Theorem**.

If *a* and *b* are in lowest terms (as supposed), their squares would each have an even number of prime factors. 5b² has one more prime factor than b², meaning it would have an odd number of prime factors. This contradiction forces the supposition wrong, so √5 cannot be rational. It is, therefore, irrational.

The irrationality of √5 is the heart of the precision problem. If it were an integer value then we could easily apply a modulus 1000 to interim calculations to truncate unnecessary precision and avoid overflow. So one strategy which helps remove the irrational component is the use of conjugates. Consider the following:

From (3) we note that ω = r - ω'. Now the problem becomes:

From this we can draw several observations. The first is that *r* is an integer, which is to be expected - that's why we introduced the conjugate in the first place. To see this you only need to examine the binomal expansion of the terms and you'll note that the even exponents of both √5 and -√5 produce *integers*, and the odd exponents of √5 and -√5 produce *irrational* components that cancel each other out since there is both a positive and negative value for each term.

The second observation is that (3 - √5) = (3 - 2.2306) < 1 which implies that ω' < 1. Since this term is less than 1 (and gets smaller as *n* gets larger) we can effectively ignore it since we are only interested in the first three digits in front of the decimal place. We have now reduced our problem to an integer calculation that won't run into any precision issues:

What is now needed is a way to calculate *r*:

We know the result is an integer but the irrational √5 component is still in the mix. The clue to solving this is given above in the binomial expansion of *r*. Looking at that it should be clear that the expansion of the first term in *r*, will produce a result in the following form:

Here the subscript *n* is used to denote the resulting co-efficients of the rational and irrational components when the original term, (3 + √5), is raised to the power *n*. When *n* equals 0 or 1 the values of the coefficients, *x* and *y*, are easy to determine. For n=0, the result must equal 1 so clearly x_{0}=1 and y_{0}=0. Since we know x_{0}, y_{0} it might be worthwhile deriving a recursive solution for higher values of *n*. To do this we can multiply the nth power by the original term and examine the results:

Equation (5) tells us that the co-efficients for the (n+1)-th power are a linear combination of the co-efficients for the n-th power. Since we know x_{0} and y_{0} we can solve for any value of n by using recursion. That is, repeatedly substituting values into (5) to get the desired x_{n} and y_{n} values. Here is the general equation in matrix form:

The astute reader will notice that we've only dealt with the first term in equation (4) above. However, the second term in (4) merely involves swapping the sign of the y_{n} component and in fact when the 2 terms are added together the y_{n} coefficients cancel out the irrational component leaving only the rational component (as shown earlier). Consequently,

So an elegant way to solve do this matrix multiplication problem quickly is to use **fast exponentiation** which requires far fewer (log n) operations than simple multiplication repeated n times. The difference between these two approaches is considerable when when *n* is very large.

Finally, here's some C# code to do this:

using System; using System.IO; namespace GoogleCodeJam.Numbers { class Program { static void Main() { string basePath = @""; StreamReader infile = new System.IO.StreamReader(basePath + "large.txt"); StreamWriter outfile = new System.IO.StreamWriter(basePath + "output.txt"); int testCases = Int32.Parse(infile.ReadLine()); for (int caseNo = 1; caseNo <= testCases; caseNo++) // note this is 1-based { // read in the input data uint n = UInt32.Parse(infile.ReadLine()); int[,] coeffs = new int[2, 2]; coeffs[0, 0] = 3; coeffs[0, 1] = 5; coeffs[1, 0] = 1; coeffs[1, 1] = 3; int[,] result = MatrixPow(coeffs, (int)n); ; double answer = (2 * result[0,0] + 999) % 1000; // write out the results outfile.WriteLine("Case #{0}: {1}", caseNo, answer.ToString("000")); } infile.Close(); outfile.Flush(); outfile.Close(); } // sqare matrix exponentiation (using fast exponentiation) static int[,] MatrixPow(int[,] a, int n) { if (n == 1) return a; if (n % 2 == 0) return MatrixPow(MatrixMultiply(2, a, a), n / 2); else return MatrixMultiply(2, MatrixPow(MatrixMultiply(2, a, a), n / 2), a); } // square matrix multiplication static int[,] MatrixMultiply(int size, int[,] m1, int[,] m2) { int[,] result = new int[size, size]; for (int i = 0; i < size; i++) { for (int j = 0; j < size; j++) { result[i, j] = 0; for (int k = 0; k < size; k++) { checked { result[i, j] += (m1[i, k] * m2[k, j]) % 1000; } } } } return result; } } }

*28 Jul 2008*
*Damien Wintour*

on 13 Jul 2009 at 2:42 am 1Some Things are Just Hard to Explain: The Taniyama-Shimura Conjecture | Necessary and Sufficient[...] In essence, this is proof by contradiction (a simple example of this can be found here). [...]

on 25 Feb 2010 at 9:16 am 2GastonGreat data, thanks for it!