online advertising

Monday, November 3, 2014

Problem 32 - Pandigital products

Problem: 

Please find the problem here.

Solution:

My current solution is kind of slow - it runs through all permutations and all possible splits to check for the product formula, and when it does accumulate the results, it is that simple.

It would be great if we could speed this up using some more constraints, but for now, this is good enough to produce and answer.

It would be an interesting detour to discuss the permutation generating algorithm. In some sense it is simple, pick the head, recursively compute the permutation of remaining elements, and concat with the head. Note that we could have been faster if we optimize the permutation method to reuse the same buffer instead of creating new list for each recursive result!

Code:

namespace Euler
{
    using System;
    using System.Collections.Generic;
    using System.Linq;
    using System.Numerics;

    internal static partial class Program
    {
        public static IEnumerable<IEnumerable<T>> Permutations<T>(IEnumerable<T> available)
        {
            List<T> availableList = available.ToList();
            return Permutations(availableList.Select(t => Pair<bool, T>.Create(true, t)).ToArray(), availableList.Count);
        }

        private static IEnumerable<IEnumerable<T>> Permutations<T>(Pair<bool, T>[] available, int count)
        {
            if (count == 0)
            {
                yield return new List<T> { };
            }

            for (int i = 0; i < available.Length; i++)
            {
                if (available[i].Item1)
                {
                    available[i].Item1 = false;
                    foreach (IEnumerable<T> recursiveResult in Permutations(available, count - 1))
                    {
                        yield return new List<T> { available[i].Item2 }.Concat(recursiveResult);
                    }

                    available[i].Item1 = true;
                }
            }
        }

        private class Pair<T, U>
        {
            public static Pair<T, U> Create(T item1, U item2)
            {
                return new Pair<T, U>(item1, item2);
            }

            private Pair(T item1, U item2)
            {
                this.Item1 = item1;
                this.Item2 = item2;
            }

            public T Item1 { get; set; }

            public U Item2 { get; set; }

        }

        public static void Problem032()
        {
            HashSet<int> results = new HashSet<int>();
            foreach (var permutation in Permutations(Enumerable.Range(1, 9).ToArray()))
            {
                // Break down the permutation into list of different partitions and check product formula
                for (int firstPart = 1; firstPart < 9; firstPart++)
                {
                    for (int secondPart = 1; secondPart < 9; secondPart++)
                    {
                        int thirdPart = 9 - firstPart - secondPart;
                        if (thirdPart > 0)
                        {
                            string permutationString = permutation.Aggregate("", (s, x) => s + x);
                            int firstPartValue = int.Parse(permutationString.Substring(0, firstPart));
                            int secondPartValue = int.Parse(permutationString.Substring(firstPart, secondPart));
                            int thirdPartValue = int.Parse(permutationString.Substring(9 - thirdPart));
                            if (firstPartValue == secondPartValue * thirdPartValue)
                            {
                                results.Add(firstPartValue);
                                //Console.WriteLine("{0} = {1} x {2}", firstPartValue, secondPartValue, thirdPartValue);
                            }
                        }
                    }
                }
            }
            Console.WriteLine(results.Aggregate((x, y) => x + y));
        }
    }
}

Problem 31 - Coin Sum

Problem:

Please find the problem here.

Solution:

It is easier to think when the number is smaller.

The number of ways to represent 15 using any number of coin is:

number of ways to represent 15 using coins (1, 2, 5) [by choosing using 0 10 coin]
+
number of ways to represent 5 using coins (1, 2, 5) [by choosing using 1 10 coin]

Similarly, the number of ways to represent 15 using coin (1, 2, 5) is:

number of ways to represent 15 using coins (1, 2, 5) [by choosing using 0 5 coin]
+
number of ways to represent 10 using coins (1, 2) [by choosing using 1 5 coin]
+
number of ways to represent 5 using coins (1, 2) [by choosing using 2 5 coin]
+
number of ways to represent 0 using coins (1, 2) [by choosing using 3 5 coin]

Of course the last term is 1, and now you get the idea how to use recursion to find the number of ways. The code do exactly that.

Code:

namespace Euler
{
    using System;
    using System.Linq;
    using System.Numerics;

    internal static partial class Program
    {
        public static void Problem031()
        {
            // Count the number of using different coins
            int[] coinValues = new int[] { 1, 2, 5, 10, 20, 50, 100, 200 };
            Console.WriteLine(Way(200, coinValues, coinValues.Length - 1));
        }

        private static int Way(int value, int[] coinValues, int maxIndex)
        {
            // This is special because we know it is 1, the value is can always be represented as sum of 1 in only 1 way.
            if (maxIndex == 0)
            {
                return 1;
            }
            int maxValue = coinValues[maxIndex];
            return Enumerable.Range(0, value / maxValue + 1).Select(t => Way(value - maxValue * t, coinValues, maxIndex - 1)).Aggregate((x, y) => x + y);
        }
    }
}

Thursday, October 2, 2014

Problem 30 - Digit fifth powers

Problem:

Please find the problem here.

Solution:

The problem can be solved with brute force, just testing if the number matches the digit sum. The only problem is when do we stop the search.

Observe that the right hand side, despite raised to high power, are relatively small in magnitude since it got broken down into digits. For a six digit number, the maximum possible right hand side is
6x95 =354294, which is a 6 digit number, but 7x95 = 413343 is not a seven digit number. So all we need to do is to test the 6 digit numbers, in fact, only those less than 354294.

The rest is trivial. Caching the fifth powers of the digit can help quite a bit.

Code:

namespace Euler
{
    using System;
    using System.Linq;
    using System.Collections.Generic;

    internal static partial class Program
    {
        private class Solution030
        {
            public void Run()
            {
                Console.WriteLine(GetPowers().Aggregate((x, y) => x + y));
            }

            private IEnumerable<int> GetPowers()
            {
                int max = 9 * 9 * 9 * 9 * 9 * 6;
                int[] powers = Enumerable.Range(0, 10).Select(x => x * x * x * x * x).ToArray();
                for (int i = 1; i <= max; i++)
                {
                    IEnumerable<int> digits = i.ToString().Select(c => c - '0');
                    int digitPowerSum = digits.Select(d => powers[d]).Aggregate((x, y) => x + y);
                    if (digitPowerSum == i && i != 1)
                    {
                        yield return i;
                    }
                }
            }
        }

        public static void Problem030()
        {
            new Solution030().Run();
        }
    }
}

In fact, listing the found sums numbers is also fun:

4150 = 45 + 15 + 55 + 05
4151 = 45 + 15 + 55 + 15
54748 = 55 + 45 + 75 + 45 + 85
92727 = 95 + 25 + 75 + 25 + 75
93084 = 95 + 35 + 05 + 85 + 45
194979 = 15 + 95 + 45 + 95 + 75 + 95

Answer: 443839

Problem 29: Distinct Powers

Problem:

Please find the problem here.

Solution:

This is a counting problem. The elements can be big (e.g. 100100), but the number of elements is small, just at most 10,000 of them.

First, note that the set of powers can intersect only if they have the same prime factors with same ratios. For example, power of 2 will intersect with power of 4, but it will never intersect with power of 6. Similarly, power of 6 will never intersect with power of 12.

So we break the problem down to computing the count of union of some complicated disjoint sets.

{power of 2 union power of 4 union .. } union ... union { power of 6 union power of 36 union ... } union ...

As those sets are disjoint, the count of union is just the sum of counts. Now we focus on the compute the count of a particular term.

To compute the size of a union of a set of sets, we can use the inclusion-exclusion principle. In problem 1, we already applied the principle to compute the sum. This one is more complicated.

The hard part is to compute the size of an arbitrary intersection. For example, what is the size of {power of 2 union power of 4}?

Taking logarithm (to base 2 in this case) on the values will make the problem look simpler. Now we have the compute the size of {2, 3, ..., 100}, {4, 6, ..., 200}, so that is just {4, 6, ... 100}, which is 49, easy to compute in this case, but what is the general rule?

Abstractly, we are given sets of multiples and we wanted to find intersection. They are simply the common multiples. To compute the number of common multiples. We simply find the least common multiple, and then count how many multiple of that least common multiple can stay within the set.

Last but not least, we have to apply the inclusion-exclusion principle formula. The formula itself look daunting when applied with arbitrary number of sets, but fortunately, the formula can be implemented using gray code enumeration.

Gray code is a specific sequence of binary strings such that adjacent element of the sequence differ by only one bit. A sample sequence of 3 element is shown below:

000
001
011
010
110
111
101
100

By definition, the gray code flip one bit every time, so the parity of the sequence alternates (number of 1 is even -> number of 1 is odd -> ...) This is exactly what we need in the inclusion exclusion sum! So I simply adopted the gray code generation algorithm from The Art of Computer Programming, and coded the formula.

As an aside - reverse engineering this piece of code (written two years ago) to the description above isn't easy at all. I should really have written the blog by the time I coded it.

Code:


namespace Euler
{
    using System;
    using System.Linq;
    using System.Collections.Generic;

    internal static partial class Program
    {

        class Solution029
        {
            public void Run()
            {
                int A = 100;
                int B = 100;
                long sum = 0;
                bool[] processed = new bool[A + 1];
                for (int a = 2; a <= A; a++)
                {
                    if (!processed[a])
                    {
                        int numberOfSets = (int)Math.Log(B, a);
                        int c = a;
                        for (int i = 1; i <= numberOfSets; i++)
                        {
                            processed[c] = true;
                            c *= a;
                        }
                        IEnumerable<long> allPowers = Enumerable.Range(1, numberOfSets).Select(t => (long)t);
                        foreach (var code in Gray(numberOfSets))
                        {
                            var joined = Enumerable.Zip(allPowers, code.Item1, (x, y) => Tuple.Create(x, y));
                            var powers = joined.Where(t => t.Item2).Select(t => t.Item1);
                            if (powers.Count() == 0)
                            {
                                // Gray code generate this sequence - but this should simply be discarded since it make
                                // no sense for having an intersection of no sets.
                                // One could argue this is simply empty set and contribute 0 to sum, that's fine too.
                                continue;
                            }
                            // The size of the intersection of a particular subset can be computed using this formula
                            // min(powers) * B / lcm
                            var lcm = powers.Aggregate((x, y) => CommonMultiple(x, y));
                            var min = powers.Min();
                            var intersectionSize = min * B / lcm;
                            if (powers.Contains(lcm))
                            {
                                intersectionSize--;
                            }
                            sum = sum + (code.Item2 ? -1 : 1) * intersectionSize;
                        }
                    }
                }
                Console.WriteLine(sum);
            }

            // From "The Art of Computer Programming, book 4"
            private IEnumerable<Tuple<List<bool>, bool>> Gray(int n)
            {
                int parity = 0;
                bool[] bits = new bool[n];
                int j = 0;
                while (true)
                {
                    // Iterating gray code and present the inclusion-exclusion sum
                    yield return Tuple.Create(bits.ToList(), parity == 0);
                    parity = 1 - parity;
                    if (parity == 1)
                    {
                        j = 0;
                    }
                    else
                    {
                        int k = 0;
                        while (true)
                        {
                            if (bits[k])
                            {
                                break;
                            }
                            else
                            {
                                k++;
                            }
                        }
                        // Now we know bits[k] is true, and the minimal such k, then 
                        j = k + 1;
                    }
                    if (j == n)
                    {
                        yield break;
                    }
                    else
                    {
                        bits[j] = !bits[j];
                    }
                }
            }
        }
        public static void Problem029()
        {
            new Solution029().Run();
        }
    }
}

Answer: 9183

Tuesday, September 23, 2014

Problem 28: Number spiral diagonals

Problem:

Please find the problem here.

Solution:

Since 10012 is small, I could have just fill the spiral and compute the sum. In fact, that is how I got the answer accepted for the first time. But I decide this is a chance for me to do some math, and so I do.



First notice the spiral (except the initial 1) can be decomposed into "Shells", each shell starts with the upper right corner and ends with a perfect square at the upper left corner. The values at the corners forms an arithmetic progression with common difference 2 * i, where i is the shell number, starting from 1.

With that observation, it is easy to translate that to the code below.

Apparently, the code can be further simplified, but I am not doing that in the code for simplicity.

Using Matlab, we can simplify this quickly using the commands below

clear
syms i
e = (2 * i + 1)^2
d = 2 * i
s = e - 3 * d
t = (s + e) * 4 / 2
expand(t)

Now t = 16i2+4i+4 is the term for summation, given the summation formula, we can even sum them up directly as follow

clear
syms n
q = n * (n + 1) * (2 * n + 1) / 6
l = n * (n + 1) / 2
o = 16 * q + 4 * l + 4 * n + 1
expand(o)

That gives the overall sum formula to be 16/3*n^3+10*n^2+26/3*n+1, an amazing simple cubic equation.

Code:

namespace Euler
{
    using System;

    internal static partial class Program
    {
        public static void Problem028()
        {
            long sum = 1;
            for (int i = 1; i < (1001 + 1) / 2; i++)
            {
                long end = (2 * i + 1) * (2 * i + 1);
                long diff = 2 * i;
                long start = end - 3 * diff;
                sum += (start + end) * 4 / 2;
            }
            Console.WriteLine(sum);
        }
    }
}

Answer: 669171001

Friday, September 5, 2014

Problem 27: Quadratic primes

Problem:

Please find the problem here.

Solution:

There isn't much we could do with this one, just try them all and find the best. The only optimization I have done is cache the is prime check, and it is quick enough for now.

namespace Euler
{
    using System;
    using System.Collections.Generic;
    using System.Linq;

    internal static partial class Program
    {
        private static Dictionary<long, bool> primeCache = new Dictionary<long, bool>();

        public static void Problem027()
        {
            Tuple<int, int> best = Tuple.Create(-1, -1);
            for (int a = -1000; a <= 1000; a++)
            {
                for (int b = -1000; b <= 1000; b++)
                {

                    long current = 0;
                    int length = 0;
                    while (true)
                    {
                        long value = current * current + a * current + b;
                        if (!IsPrime(value))
                        {
                            break;
                        }

                        length++;
                        current++;
                    }
                    if (length > best.Item1)
                    {
                        best = Tuple.Create(length, a * b);
                    }
                };
            };
            Console.WriteLine(best.Item1);
            Console.WriteLine(best.Item2);
        }

        private static bool IsPrime(long x)
        {
            if (x < 0)
            {
                x = -x;
            }
            bool result;
            if (primeCache.TryGetValue(x, out result))
            {
                return result;
            }
            for (long y = 2; y < x; y++)
            {
                if (x % y == 0)
                {
                    primeCache.Add(x, false);
                    return false;
                }
            }
            primeCache.Add(x, true);
            return true;
        }
    }
}

Answer: -59231

Given the restriction - the best quadratic formula can only output 71 primes.

Problem 26:Reciprocal cycles

Problem:

Please find the problem here.

Solution:

The key for this problem is to find the recurring cycle. For human, finding the recurring cycle involve long division and finding repeating dividend. We can code exactly that by dividing, multiply the remainder by 10 (much like adding zero in long division), and divide again, until we see the same dividend again.


namespace Euler
{
    using System;
    using System.Collections.Generic;
    using System.Linq;

    internal static partial class Program
    {
        public static void Problem026()
        {
            IEnumerable<Tuple<int, int>> recurringPartLength = GetRecurringPartLengths();
            int maxLength = recurringPartLength.Select(t => t.Item2).Max();
            Console.WriteLine(recurringPartLength.Single(t => t.Item2 == maxLength).Item1);
        }

        private static IEnumerable<Tuple<int, int>> GetRecurringPartLengths()
        {
            for (int d = 2; d < 1000; d++)
            {
                Tuple<List<int>, int> decimalRepresentation = GetDecimalRepresentation(d);
                if (decimalRepresentation != null)
                {
                    IEnumerable<int> recurringPart = decimalRepresentation.Item1.Skip(decimalRepresentation.Item2 - 1);
                    yield return Tuple.Create(d, recurringPart.Count());
                }
            }
        }

        private static Tuple<List<int>, int> GetDecimalRepresentation(int d)
        {
            int n = 1;
            List<int> quotients = new List<int>();
            // We need to stop when we see the same number to divide - we already know what would happen
            Dictionary<int, int> seen = new Dictionary<int, int>();
            int position = 0;
            while (true)
            {
                n = n * 10;
                position++;
                int previousIndex;

                if (seen.TryGetValue(n, out previousIndex))
                {
                    return Tuple.Create(quotients, previousIndex);
                }
                else
                {
                    seen.Add(n, position);
                    int quotient = n / d;
                    quotients.Add(quotient);
                }
                n = n % d;
                if (n == 0)
                {
                    return null;
                }
            }
        }
    }
}

Answer: 983

Thursday, September 4, 2014

Problem 25: 1000-digit Fibonacci number

Problem:

Please find the problem here.

Solution:

It is not hard to find online the Binet formula to compute Fibonacci number.
$ F_n = (\frac{\sqrt{5} + 1}{2})^n + (\frac{\sqrt{5} - 1}{2})^n $

Note that the second term is an exponential of a number less than 1, that goes to 0 quickly, so let's just ignore that and we get a good approximation

$ \begin{eqnarray*} F_n & \simeq & (\frac{\sqrt{5} + 1}{2})^n \\ \log(F_n) & \simeq & n \log (\frac{\sqrt{5} + 1}{2}) \end{eqnarray*} $

This gives us $ n \simeq 4785 $ when $ \log(F_n) \simeq 1000 $. That means if we use the trivial algorithm for computing Fibonacci number, it takes around 4785 additions and we will get the answer. However, since we know the answer will be approximately 4785, we have a much faster algorithm that allow us to compute $ F_{4785} $.

The algorithm start with observing the matrix

$ \begin{eqnarray*} \left(\begin{matrix}1 & 1 \\ 1 & 0\end{matrix}\right)^n = \left(\begin{matrix}F_{n+1} & F_{n} \\ F_{n} & 0\end{matrix}\right) \end{eqnarray*} $

With that, we can easily compute $ F_n $ by repeated squaring. Observing the matrix is symmetric we can reduce some arithmetic operations as well.

Once we get to the approximate solution, we need to walk to the solution. By some computation I know $ F_{4785} $ has 1,000 digits, so we need to walk backwards to see when it get started to become 999 digits, and that's when we can stop.

namespace Euler
{
    using System;
    using System.Numerics;

    internal static partial class Program
    {
        public static void Problem025()
        {
            new Solution025().Execute();
        }

        private class Solution025
        {
            public void Execute()
            {
                FibMatrix one = new FibMatrix
                {
                    a = 1,
                    b = 1,
                    d = 0
                };

                // By Binet formula approximation, we know the solution is around this
                int n = 4784;
                FibMatrix m = Power(one, n, one, Multiply);

                BigInteger big = m.a;
                BigInteger small = m.b;
                // small is the nth Fibonacci number

                // I already know small has 1000 digits by debugging, so I know I need to going backwards
                // in full generality I need to check if I need forward as well
                while (small.ToString().Length >= 1000)
                {
                    BigInteger prev = big - small;
                    big = small;
                    small = prev;
                    n--;
                    // small is still the nth Fibonacci number after the line
                }

                // Small has 999 digits, big has 1000 digit, so the answer is n + 1
                Console.WriteLine(n + 1);
            }

            private FibMatrix Multiply(FibMatrix p, FibMatrix q)
            {
                FibMatrix result = new FibMatrix();
                result.a = p.a * q.a + p.b * q.b;
                result.b = p.a * q.b + p.b * q.d;
                result.d = p.b * q.b + p.d * q.d;
                return result;
            }

            private class FibMatrix
            {
                public BigInteger a;
                public BigInteger b;
                public BigInteger d;
            }
        }
    }
}

Answer: 4782

Note how close is it for our analytic approximation.

Wednesday, September 3, 2014

Problem 24: Lexicographic permutations

Problem:

Please find the problem here.

Solution:

By dividing the set of permutations by the first digit into partitions, the set of all permutation of n digits can be divided into n partitions of size (n-1)!. The first thing to do is to determine which partition the kth permutation belongs to, it can easily computed to be floor((k - 1)/(n-1)!), that determine the first digit. Then, we know it is the (k-1) % (n-1)! permutation in the remaining digits, so we can just compute that recursively.

Once we have the analysis, the code came by really easily.

namespace Euler
{
    using System;
    using System.Collections.Generic;
    using System.Linq;

    internal static partial class Program
    {
        public static void Problem024()
        {
            int[] factorials = new int[11];
            factorials[0] = 1;
            for (int i = 1; i <= 10; i++)
            {
                factorials[i] = factorials[i - 1] * i;
            }

            // There are nine factorials number starts with 0
            // There are nine factorials number starts with 1, 
            // ...

            int current = 1000000 - 1;
            List<int> availableDigits = Enumerable.Range(0, 10).ToList();
            while (availableDigits.Count > 0)
            {
                int partitionSize = factorials[availableDigits.Count - 1];
                int quotient = current / partitionSize;
                int remainder = current % partitionSize;
                int partitionNumber = quotient + 1;
                current = remainder;
                Console.Write(availableDigits[partitionNumber - 1]);
                availableDigits.RemoveAt(partitionNumber - 1);
            }
            Console.WriteLine();
        }
    }
}

Answer: 2783915460

Problem 23: Non-abundant sums

Problem:

Please find the problem here.

Solution:

As we have already explained in Problem 21, it is easy to tell if a number is abundant by computing the sum of proper divisor function. Since we are given that every number larger than 28123 can be written as sum of two abundant numbers, so we don't have to consider any of the abundant numbers larger than that.

So we first creating a set of all abundant numbers less than 28124 (turn out there are around 7,000 of them), then we create the set of all abundant number sums. Finally, we evaluate all the number below 28124 and test if they're sum of abundant number, that's all.

namespace Euler
{
    using System;
    using System.Collections.Generic;
    using System.Linq;

    internal static partial class Program
    {
        public static void Problem023()
        {
            List<int> abundantNumbers = new List<int>();
            for (int i = 2; i <= 28124; i++)
            {
                if (IsAbundantNumbers(i))
                {
                    abundantNumbers.Add(i);
                }
            }

            HashSet<int> abundantSums = new HashSet<int>();
            foreach (int abundantNumber1 in abundantNumbers)
            {
                foreach (int abundantNumber2 in abundantNumbers)
                {
                    abundantSums.Add(abundantNumber1 + abundantNumber2);
                }
            }
            long answer = (1 + 23) * 23 / 2;
            for (int i = 24; i < 28124; i++)
            {
                if (!abundantSums.Contains(i))
                {
                    answer += i;
                }
            }
            Console.WriteLine(answer);
        }

        private static bool IsAbundantNumbers(int i)
        {
            return SumOfProperDivisor(i) > i;
        }
    }
}

Note that the code reused the SumOfProperDivisor function in Problem 21.

Answer: 4179871

Monday, September 1, 2014

Problem 22: Names scores

Problem:

Please find the problem here.

Solution:

Just compute all the scores and add them up. Nothing special, so I will just post the code.

namespace Euler
{
    using System;
    using System.Linq;

    internal static partial class Program
    {
        public static void Problem022()
        {
            string[] names = ReadResourceAsString("Euler.Problem022.txt").Split(new char[] { '"', ',' }, StringSplitOptions.RemoveEmptyEntries);
            Array.Sort(names);
            var nameValues = names.Select((s, i) => Tuple.Create(s, i + 1L)).Select(t => Tuple.Create(t.Item1, t.Item1.Select(c => c - 'A' + 1).Aggregate((x, y) => x + y), t.Item2));
            Console.WriteLine(nameValues.Select(t => t.Item2 * t.Item3).Aggregate((x, y) => x + y));
        }
    }
}

Answer: 871198282

Problem 21: Amicable numbers

Problem:

Please find the problem here.

Solution:

Finally, we got another problem that feel like more number theoretic. A straightforward implementation would be checking, for each number i under 10000, if the SumOfProperDivisor(SumOfProperDivisor(i)) == i, and if so, add it to the sum, and that's all.

Therefore the problem becomes computing sum of proper divisor, turn out it is easier to consider the sum of divisor first. The sum of divisor, like the divisor function we encountered in problem 12, is multiplicative. Just like in proving the divisor function is multiplicative, we use the one to one correspondence between the a divisor of pq and the divisor of p, divisor of q pair. Now group the pairs by the first component we can see the sum and be factorized and we can see the sum of divisor function is multiplicative. The sum of divisor function for a prime power is simply a geometric series and is easily evaluated. 

The last twist in the problem is the SumOfProperDivisor(SumOfProperDivisor(i)) part. Turn out SumOfProperDivisor(i) can be larger than i (which are called abundant numbers), so we need to evaluate more SumOfProperDivisor(i) to its max to make sure we don't index out of bound.

Without further ado, here is the code:

namespace Euler
{
    using System;
    using System.Collections.Generic;
    using System.Linq;

    internal static partial class Program
    {
        public static void Problem021()
        {
            Dictionary<int, int> d = new Dictionary<int, int>();
            d.Add(1, 0);
            for (int i = 2; i <= 10000; i++)
            {
                d.Add(i, SumOfProperDivisor(i));
            }
            long maxSumOfDivisorWithin10000 = d.Values.Max();
            for (int i = 10001; i <= maxSumOfDivisorWithin10000; i++)
            {
                d.Add(i, SumOfProperDivisor(i));
            }
            Console.WriteLine(AmicableNumbers(d).Aggregate((x, y) => x + y));
        }

        private static IEnumerable<int> AmicableNumbers(Dictionary<int, int> d)
        {
            for (int i = 2; i <= 10000; i++)
            {
                if (d[d[i]] == i & i != d[i])
                {
                    yield return i;
                }
            }
        }

        private static int SumOfProperDivisor(int i)
        {
            return SumOfDivisor(i) - i;
        }

        private static int SumOfDivisor(int i)
        {
            var groupedFactors = BruteForceFactor(i).ToList();
            return groupedFactors.Select(kvp => (Power(kvp.Item1, kvp.Item2 + 1) - 1) / (kvp.Item1 - 1)).Aggregate((x, y) => x * y);
        }

        private static int Power(int x, int y)
        {
            if (y == 0)
            {
                return 1;
            }
            else if (y == 1)
            {
                return x;
            }
            else
            {
                int r = Power(x, y / 2);
                return y % 2 == 0 ? r * r : r * r * x;
            }

        }
    }
}

We could have optimized the factoring by eliminating trials using the Sieve_of_Eratosthenes. The idea is, instead of crossing out the number, we compute the power and store the information along. So at the end of the sieve, we got all factorizations we will need. It don't bother me to optimize right now as the performance of the previous code is fine.

Note that we reused the BruteForceFactor routine from problem 12.

Answer: 31626

Sunday, August 31, 2014

Problem 20: Factorial digit sum

Problem:

Please find the problem here.

Solution:

It is well known that we can use Strling's approximation to approximate factorial, let see what we'll get from it. First of all, we have
$ \begin{eqnarray*} \sqrt{2 \pi} n^{n+\frac{1}{2}} e^{-n} \leq n! \leq e n^{n+\frac{1}{2}} e^{-n} \end{eqnarray*} $

Only the term $n^{n+\frac{1}{2}}e^{-n}$ matters, the rest are small. That gives us 100! is approximately $ 1.01 \times 10^{158} $.

This number is not small, but it is manageable to just compute 100! I don't have other solution to compute the sum of the digits anyway. Here is the code

namespace Euler
{
    using System;
    using System.Linq;
    using System.Numerics;

    internal static partial class Program
    {
        public static void Problem020()
        {
            BigInteger prod = new BigInteger(1);
            for (int i = 2; i <= 100; i++)
            {
                prod = BigInteger.Multiply(prod, new BigInteger(i));
            }
            Console.WriteLine(prod.ToString().Select(i => i - '0').Aggregate((x, y) => x + y));
        }
    }
}

Answer: 648

The value of 100! is actually
93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000

It has 159 digits, so the Stirling approximation above is pretty close.

Problem 19: Counting Sundays

Problem:

Please find the problem here.

Solution:

We might be able to compute it by the description, but it is much easier to use the built-in calendar. Here is the code:

namespace Problem019
{
    using System;

    class Program
    {
        static void Main(string[] args)
        {
            int sum = 0;
            for (int year = 1901; year <= 2000; year++)
            {
                for (int month = 1; month <= 12; month++)
                {
                    if (new DateTime(year, month, 1).DayOfWeek == DayOfWeek.Sunday)
                    {
                        sum++;
                    }
                }
            }
            Console.WriteLine(sum);
            
        }
    }
}

Answer: 171