# Duncan’s blog

## June 19, 2011

### Project Euler problem 26

Filed under: Project Euler,Python — duncan @ 9:27 pm

A unit fraction contains 1 in the numerator. The decimal representation of the unit fractions with denominators 2 to 10 are given:

1/2 = 0.5
1/3 = 0.(3)
1/4 = 0.25
1/5 = 0.2
1/6 = 0.1(6)
1/7 = 0.(142857)
1/8 = 0.125
1/9 = 0.(1)
1/10 = 0.1

Where 0.1(6) means 0.166666…, and has a 1-digit recurring cycle. It can be seen that 1/7 has a 6-digit recurring cycle.

Find the value of d < 1000 for which 1/d contains the longest recurring cycle in its decimal fraction part.

Another one I initially attempted with ColdFusion, but it didn’t give decimals to enough precision (at least without using the underlying Java), so turned to Python.

To me this seemed like a good place to use a regular expression, where you can match a group and any repetitions of it, using back-references. Some problems around which parts should be greedy or not. For instance if you have 0.01010101 is the repeating pattern 0101 twice, or 01 four times?

Here’s the regular expression I used:

“^0\.[0-9]*([0-9]{7,}?)\1+[0-9]*?\$”

And here’s an explanation of the various parts of it:

^0\.[0-9]*([0-9]{7,}?)\1+[0-9]*?\$

The ^ matches the beginning of the string, and the \$ matches the end of the string. They’re not always strictly necessary, but it can be good practice for examples like this.

^0\.[0-9]*([0-9]{7,}?)\1+[0-9]*?\$

The value is 1/d. As d starts at 2 and increases to 1000, the value starts at 0.5 and gets closer and closer to zero. So the value is always going to start with a zero. Normally a . in a regular expression indicates ‘any character’. The \. just tells the regular expression engine to ‘escape’ it and treat it as a string literal, not a regex special character.

^0\.[0-9]*([0-9]{7,}?)\1+[0-9]*?\$

[0-9] just means match any digit between 0 and 9. The * means 0 or many. Notice this is before the part in the ( ) that I’m really interested in. So in other words, there may be some digits before the pattern. e.g. for 1/6, which is 0.166666… I don’t care about that first digit 1.

^0\.[0-9]*([0-9]{7,}?)\1+[0-9]*?\$

The ( ) is a way of grouping parts of a regex, in this case for the benefit of back-referencing.
The fact we’re given a 6 digit example means we can assume the answer has to be at least 7 digits. So I used the {7,} syntax to mean the group has to be a minimum of 7 digits long.
The ? means 0 or 1. In other words, only give me up to 1 matching group that is a minimum of 7 characters long, and don’t be greedy!

^0\.[0-9]*([0-9]{7,}?)\1+[0-9]*?\$

The \1 is a back reference to the earlier grouped part. If there were multiple parts each with their own ( ) you could back-reference them with \2, \3 etc.
And the + means 1 or more, so there has to be at least 1 repeating group.

^0\.[0-9]*([0-9]{7,}?)\1+[0-9]*?\$

The final [0-9] is for any trailing digits that don’t fit the pattern. Even though once the pattern starts in theory it continues infinitely, in reality you might get a value that is rounded, e.g. where you have 0.16666667 instead of 0.16666666…
The * means zero or many of these trailing digits, and the ? stops it from being greedy. Interestingly, you’d assume that just the ? on its own would be more efficient, as that extra digit is just going to be one optional digit and no more. But that reduces the performance somewhat.

And here’s the complete code.

```import re
from decimal import *

def displaymatch(match):
if match is None:
return 0

return len(match.group(1))

length = 0
maxlen = 0
d = 0
getcontext().prec = 2000

for i in range(1,1000):
fraction = 1/Decimal(i)
pattern = re.search(r"^[0-9]\.[0-9]*([0-9]{7,}?)(\1+)[0-9]*?\$", str(fraction))

length = displaymatch(pattern)

if length > maxlen:
maxlen = length
d = i

print(d)```

I’m using the Python regular expression and Decimal libraries. I also had to increase the precision of the decimal context to up to 2000, as lower values were giving me incorrect answers. i.e. it was truncating the fraction a bit, which meant my pattern matching wasn’t finding the correct answer. So by increasing the precision, I increased the length of the fraction being returned, and therefore improved the chances of finding the correct answer.

## December 13, 2009

### Project Euler problem 58

Filed under: Coldfusion,Project Euler — duncan @ 9:00 am
Tags: , , ,

Starting with 1 and spiralling anticlockwise in the following way, a square spiral with side length 7 is formed.

```37 36 35 34 33 32 31
38 17 16 15 14 13 30
39 18  5  4  3 12 29
40 19  6  1  2 11 28
41 20  7  8  9 10 27
42 21 22 23 24 25 26
43 44 45 46 47 48 49
```

It is interesting to note that the odd squares lie along the bottom right diagonal, but what is more interesting is that 8 out of the 13 numbers lying along both diagonals are prime; that is, a ratio of 8/13 ≈ 62%.

If one complete new layer is wrapped around the spiral above, a square spiral with side length 9 will be formed. If this process is continued, what is the side length of the square spiral for which the ratio of primes along both diagonals first falls below 10%?

This problem has some similarities to problem 28. In that case the spiral went clockwise; this one it goes anti-clockwise. For our purposes that’s a red herring.

We want to loop round, calculating the value for the corners of our spiral.
For each value, work out if it is prime (re-using the IsPrime function I first used back in problem 3).
Keep count of how many primes there are.
Also keep track of how long the side of our spiral is.
At each loop round the spiral, check if the ratio of primes / length is less than 10%.

```<cfset numPrimes = 0>
<cfset length = 1>
<cfset number = 1>
<cfset increment = 2>
<cfset ratio = 100>

<cfloop condition="ratio GT 10">
<cfset length += 2>

<!--- add our increment 4 times, calculating each time if it's prime --->
<cfloop index="i" from="1" to="4">
<cfset number += increment>

<!--- the 4th corner is never prime --->
<cfif i LT 4 AND IsPrime(number)>
<cfset numPrimes++>
</cfif>
</cfloop>

<cfset ratio = numPrimes / ((2*length)-1) * 100>

<cfset increment += 2>
</cfloop>

<cfoutput>#length#</cfoutput>

<cfscript>
function isPrime(x)
{
var isPrime = true;
var i = 0;

if (x LT 2)
{
return false;
}

if ((NOT x MOD 2) AND (x GT 2))
{	// a multiple of 2, but not 2 itself
return false;
}

for (i = 3; i LTE SQR(x); i = i + 2)
{
if (NOT x MOD i)
{	// found a factor of x
return false;
}
}

return isPrime;
}
</cfscript>```

Like in problem 28, we can see the corner values go:
3, 5, 7, 9, [increment: 2]
13, 17, 21, 25, [increment: 4]
31, 37, 43, 49 [increment: 6]
etc.

So we can work out the corner values by starting at 1, adding 2 to our increment each time round the spiral, then adding the increment four times to calculate the corners.

We don’t need to work out if the 4th corner is prime; it will always be a square number. e.g. on spiral with length 3, the 4th corner is 32 = 9. On length 5, the 4th corner is 52 = 25, and so on. So the if statement <cfif i LT 4 AND IsPrime(number)> will only execute the IsPrime() function if i LT 4 evaluates to true.

## December 12, 2009

### Project Euler problem 44

Filed under: Project Euler,Python — duncan @ 3:08 pm
Tags: , ,

Pentagonal numbers are generated by the formula, Pn=n(3n−1)/2. The first ten pentagonal numbers are:

1, 5, 12, 22, 35, 51, 70, 92, 117, 145, …

It can be seen that P4 + P7 = 22 + 70 = 92 = P8. However, their difference, 70 − 22 = 48, is not pentagonal.

Find the pair of pentagonal numbers, Pj and Pk, for which their sum and difference is pentagonal and D = |Pk − Pj| is minimised; what is the value of D?

I’d originally made a start on this in Python, then gave up on it. Came back to it much later and tried doing it in ColdFusion. My solution there took too long to run, so I went back to my Python version and finished it off.

Here’s the code:

```import math

def pentagonal(x):
return x * (3 * x - 1) / 2

def isPentagonal(x):
n = (math.sqrt((24 * x) + 1) + 1) / 6

if n == int(n) and n > 0:
return 1
else:
return 0

i = 0
pentagonals = []
solution = 0

while 1:
i += 1
p = pentagonal(i)

for j in pentagonals:
diff = abs(p - j)

if isPentagonal(diff) == 1:
intSum = p + j

if isPentagonal(intSum) == 1:
print("solution:", int(diff))
solution = 1
break

if solution:
break

pentagonals.append(p)
```

The formula for calculating the pentagonal numbers is given to us. I already had a formula for the reverse, checking if a number is pentagonal, from problem 45.
To calculate the square root, I needed to import the math library.

Basically I loop until I find a solution. I calculate the pentagonal value of i. Then I subtract each pentagonal number smaller from it (but use the abs() function to get the absolute value).
If that difference is pentagonal, then I do a sum as well.
If that is also pentagonal, then our solution is the difference.
As I go along I append each pentagonal number into the array.
This solution runs in about 12 seconds for me.

## December 9, 2009

### Project Euler problem 55

Filed under: Project Euler,Python — duncan @ 12:00 am
Tags: , , , ,

Problem 55:

If we take 47, reverse and add, 47 + 74 = 121, which is palindromic.

Not all numbers produce palindromes so quickly. For example,

349 + 943 = 1292,
1292 + 2921 = 4213
4213 + 3124 = 7337

That is, 349 took three iterations to arrive at a palindrome.

Although no one has proved it yet, it is thought that some numbers, like 196, never produce a palindrome. A number that never forms a palindrome through the reverse and add process is called a Lychrel number. Due to the theoretical nature of these numbers, and for the purpose of this problem, we shall assume that a number is Lychrel until proven otherwise. In addition you are given that for every number below ten-thousand, it will either (i) become a palindrome in less than fifty iterations, or, (ii) no one, with all the computing power that exists, has managed so far to map it to a palindrome. In fact, 10677 is the first number to be shown to require over fifty iterations before producing a palindrome: 4668731596684224866951378664 (53 iterations, 28-digits).

Surprisingly, there are palindromic numbers that are themselves Lychrel numbers; the first example is 4994.

How many Lychrel numbers are there below ten-thousand?

Yet another problem which was easy to code in ColdFusion, but it couldn’t deal with the large numbers generated. So rewrote it in Python.

Just for interest, here’s what I did in CFML. NB: there’s a logic error here, which I fixed once I moved to Python:

```<cfscript>
function isPalindrome(x)
{
if (x EQ Reverse(x))
return true;
else
return false;
}
</cfscript>

<cfset Lychrel = ArrayNew(1)>

<cfloop index="i" from="1" to="9999">
<cfset k = i>

<cfloop index="j" from="1" to="50">

<cfset temp = k + Reverse(k)>

<cfif IsPalindrome(temp)>
<cfoutput>#k# + #Reverse(k)# = #temp#</cfoutput><br>
<cfset ArrayAppend(Lychrel, i)>
<cfbreak>
</cfif>

<cfset k = temp>
</cfloop>
</cfloop>

<cfoutput><strong>#ArrayLen(Lychrel)#</strong></cfoutput>
```

Running this it quickly throws an error:
‘The value 210+E11200002108.1 cannot be converted to a number.’

Python doesn’t seem to have a built-in Reverse function like ColdFusion, so you have to write your own.

```def Reverse(x):
num = str(x)
newnum = ''

for i in range(len(num)-1, 0, -1):
newnum += num[i]

newnum += num[0]

return newnum

Lychrel = []
limit = 10000

for i in range(1, limit):
k = i

for j in range(50):
temp = k + int(Reverse(k))

if str(temp) == str(Reverse(temp)):
break

k = temp

if j == 49:
Lychrel.append(i)

print(len(Lychrel))
```

This code is pretty fast, but it could be better. For instance I’m casting everything to a string or an int as required, but that could maybe be tidied up somehow.
Also in my Reverse function, I somehow couldn’t work out how to properly loop backwards through the string, so ended up having to add on the first character outside of the loop.

Initially I was appending to my Lychrel array at the point where I break out of the loop. However on reading the question again I realised the Lychrel numbers are the numbers which aren’t solved within 50 iterations, not the other way around. I don’t even need an array, I could just have a counter that I increment. But having it in the array gives us the advantage of being able to output the array contents to make sure we’re on the right track.

## February 6, 2009

### Project Euler problem 56

Filed under: Project Euler,Python — duncan @ 12:00 am
Tags: , , ,

Problem 56:

A googol (10100) is a massive number: one followed by one-hundred zeros; 100100 is almost unimaginably large: one followed by two-hundred zeros. Despite their size, the sum of the digits in each number is only 1.

Considering natural numbers of the form, ab, where a, b < 100, what is the maximum digital sum?

So, let’s do nested loops, both going 1 to 100. Inside the inner loop, calculate a^b (or i^j in this case). Then loop through the digits of that value, adding them up. Keep track of which is the largest value this produces.

Running time about 3 seconds in Python:

```import time
tStart = time.time()

maxsum = 0

for i in range(100):
for j in range(100):
num = i ** j
total = 0

for k in str(num):
total += int(k)

if total > maxsum:
maxsum = total

print(maxsum)
print("time:" + str(time.time() - tStart))
```

## February 4, 2009

### Project Euler problem 53

Filed under: Project Euler,Python — duncan @ 12:00 am
Tags: , ,

Problem 53:

There are exactly ten ways of selecting three from five, 12345:

123, 124, 125, 134, 135, 145, 234, 235, 245, and 345

In combinatorics, we use the notation, 5C3 = 10.

In general,
nCr = n! / (r!(n-r)!)
where r <= n, n! = n*(n-1)*…*3*2*1, and 0! = 1.

It is not until n = 23, that a value exceeds one-million: 23C10 = 1144066.

How many, not necessarily distinct, values of nCr, for 1 <= n <= 100, are greater than one-million?

At first I looked at this and thought it was all about iterating through the combinations of digits, just like how it lists 123, 124, etc. Which sounds tricky. But it’s much simpler than that. In fact they give you the hardest part of the algorithm, the equation needed:
n! / (r!(n-r)!)

In Python this solution takes about 0.36 seconds. I’m going to start using the time module to calculate running times.

```import time
tStart = time.time()

values = 0

def factorial(x):
y = 1
for i in range(1,x+1):
y = y * i

return y

def calc(n, r):
return  factorial(n) / (factorial(r) * factorial(n-r))

for i in range(1,101):
for j in range(1, i):
if calc(i, j) > 1000000:
values += 1

print(values)
print("time:" + str(time.time() - tStart))
```

So, two functions, one calculates the factorial of x, the other just calculates our equation. The second function could instead have just been expressed as one line in our procedural code.

Then loop from 1 to 100. Within that, do an inner loop from 1 to whatever our current outer loop is on, passing the two loop counters into our function. We don’t really care which numbers produce values greater than a million; we only care how many numbers there are. So value is just a simple counter.

## February 3, 2009

### Project Euler problem 40

Filed under: Coldfusion,Project Euler,Python — duncan @ 12:00 am
Tags: , , , ,

Problem 40:

An irrational decimal fraction is created by concatenating the positive integers:

0.123456789101112131415161718192021…

It can be seen that the 12th digit of the fractional part is 1.

If dn represents the nth digit of the fractional part, find the value of the following expression.

d1 * d10 * d100 * d1000 * d10000 * d100000 * d1000000

Firstly I completed this in Coldfusion. I got the correct answer by brute-forcing it, but the solution was slow. This took about 11 minutes on a CF 5 server; it would undoubtedly be a lot quicker on any CFMX server. I then rewrote it in Python, which spits out the correct answer within a couple of seconds.

I thought it might be interesting to do a code comparison between the two languages.
CFML:

```<cfset fraction = "">
<cfset total = 1>

<cfloop index="i" from="1" to="1000000">
<!--- we're going to break out of the loop sooner than that --->
<cfset fraction = fraction & i>

<cfif Len(fraction) GTE 1000000>
<cfbreak>
</cfif>
</cfloop>

<cfset total = total *	Mid(fraction, 1, 1) *
Mid(fraction, 10, 1) *
Mid(fraction, 100, 1) *
Mid(fraction, 1000, 1) *
Mid(fraction, 10000, 1) *
Mid(fraction, 100000, 1) *
Mid(fraction, 1000000, 1)>

<cfoutput>#total#</cfoutput>
```

Python:

```fraction = ""
total = 1

for i in range(1,1000000):
fraction = fraction + str(i)

if len(fraction) >= 1000000:
break

total = total * \
int(fraction[0]) * \
int(fraction[9]) * \
int(fraction[99]) * \
int(fraction[999]) * \
int(fraction[9999]) * \
int(fraction[99999]) * \
int(fraction[999999])

print(total)
```

Almost identical code. I’m still not anywhere near writing pythonic code, but it’s early days.

Basically I loop from 1 to 1000000 (not really that far, because I break out of the loop once our string has length of 1000000). As I go I append each loop counter to the string.

Then I just multiply together the 1st, 10th, 100th etc values from the string. That part could also have been done in a loop to simplify the code slightly. The whole thing’s not particularly clever, but it works.

If you’re coming to this from a Coldfusion background, you’re maybe wondering why I’m doing a DIV between each of the values before multiplying them. Because whitespace is important in Python, the \ can be used to join two lines together where the line length gets too long.

## February 2, 2009

### Project Euler problem 34

Filed under: Coldfusion,Project Euler — duncan @ 12:00 am
Tags: , , , , ,

Problem 34:

145 is a curious number, as 1! + 4! + 5! = 1 + 24 + 120 = 145.

Find the sum of all numbers which are equal to the sum of the factorial of their digits.

Note: as 1! = 1 and 2! = 2 are not sums they are not included.

This one doesn’t seem too hard. However, the tricky part is determining what the upper bound should be. I wasn’t really sure how this one would progress, so decided to set a relatively low upper bound, 99999. I was then going to check what sort of values I was getting back, then increase the bound as necessary.

However running the following code, I saw there weren’t very many values at all. Just on the off-chance I’d got the final solution, I entered it into the Project Euler page, and I had indeed got it right, purely by luck.

I later discovered that the real upper bound that you would theoretically need to use is 9999999. Doing that with the following code would take much longer.

```<cfscript>
function factorial(x)
{
var i = 0;
var factorial = 1;

for (i = x; i GT 1; i = i - 1)
{
factorial = factorial * i;
}

return factorial;
}
</cfscript>

<cfset numbers = ArrayNew(1)>
<cfset total = 0>

<!--- how far do we need to loop?  Let's just start and see how it progresses... --->
<cfloop index="i" from="3" to="99999">
<cfset sum = 0>

<!--- get the factorial of each digit of i --->
<cfloop index="j" from="1" to="#Len(i)#">
<cfset sum = sum + factorial(Mid(i, j, 1))>
</cfloop>

<!--- does i = the sum of the factorials? --->
<cfif i EQ sum>
<cfset ArrayAppend(numbers, i)>
</cfif>
</cfloop>

<cfdump var="#numbers#">

<cfloop index="i" from="1" to="#ArrayLen(numbers)#">
<cfset total = total + numbers[i]>
</cfloop>

<cfoutput>#total#</cfoutput>
```

Factorials can get very large quickly. Fortunately I’m only ever going to be putting values from 0 to 9 into my factorial function.

Note that 0! = 1, not zero. This is the empty product. Fortunately my factorial function takes care of that.

Then I loop, ignoring 1 and 2, to my lucky choice of 99999 for the upper limit (the comment in my code gives an idea of how I was thinking!). Sum up the factorial of each digit. Check if that sum is the same as our number; if so add it to an array.

At the end I dumped out the array, just to see what sort of results I was getting. Then later I added in the final loop, adding up the values of the array to give the correct solution.

Judging by the forum discussion, there wasn’t any clear idea of what the upper limit should be, with most people seeming to take an educated guess.

Spoiler alert:
Interestingly, there are only four numbers that this works for, and they’re known as factorions. Two of those are 1 and 2, so the solution to this puzzle only requires you to work out two other numbers. And one of them, 145, has already been given to you in the problem!

## February 1, 2009

### Project Euler problem 28

Filed under: Coldfusion,Project Euler — duncan @ 7:37 pm
Tags: , , ,

Problem 28:

Starting with the number 1 and moving to the right in a clockwise direction a 5 by 5 spiral is formed as follows:

```21 22 23 24 25
20  7  8  9 10
19  6  1  2 11
18  5  4  3 12
17 16 15 14 13```

It can be verified that the sum of both diagonals is 101.

What is the sum of both diagonals in a 1001 by 1001 spiral formed in the same way?

There’s various ways to look at this problem. One thing I noticed was that, in the example above, the numbers we’re adding go:
1,
3, 5, 7, 9,
13, 17, 21, 25

So, excluding the 1 in the centre, we have four numbers that increment by 2, then four numbers that increment by 4. If we made the spiral wider, the next set of four numbers would increment by 6, and so on. I used that observation as the basis for this CFML:

```<cfset total = 1>
<cfset runningtotal = 1>

<cfset length = 1001>

<cfloop index="i" from="2" to="#(length-1)#" step="2">
<cfloop index="j" from="1" to="4">
<cfset runningtotal = runningtotal + i>

<cfset total = total + runningtotal>
</cfloop>
</cfloop>

<cfoutput>#total#</cfoutput>```

I’m using two variables here, one is our grand total. The other is a running total. So on each iteration of the outer loop, our running total will have increased by (4 * 2), (4 * 4), (4 * 6), etc… The grand total will have all instances of running total added to it.

There’s neater ways this could have been expressed, but this was how I chose to do it; it’s simple and fast.

A word of warning: even though the problem states “What is the sum of both diagonals“, that’s not really the solution they’re looking for. If you were to add up the two diagonals in the example above, you’d have:
Diagonal 1: 21 + 7 + 1 + 3 + 13 = 45
Diagonal 2: 17 + 5 + 1 + 9 + 25 = 57
Adding those two values would give you 102, not 101. This is because you’d have added 1 twice. What the problem really wants here is the sum of all the unique values of both diagonals.

## January 31, 2009

### Project Euler problem 30

Filed under: Project Euler,Python — duncan @ 4:01 pm
Tags: , ,

Problem 30:

Surprisingly there are only three numbers that can be written as the sum of fourth powers of their digits:

1634 = 14 + 64 + 34 + 44
8208 = 84 + 24 + 04 + 84
9474 = 94 + 44 + 74 + 44

As 1 = 14 is not a sum it is not included.

The sum of these numbers is 1634 + 8208 + 9474 = 19316.

Find the sum of all the numbers that can be written as the sum of fifth powers of their digits.

Another problem where I’d worked out the logic in Coldfusion, but had to rewrite it in Python to get the solution.

We’re going to loop through numbers testing them out. We know we can skip 1. But how far do we need to loop up to? Let’s look at 4 digit numbers first; the lowest 4-digit number is 1000.  1^5 + 0^5 + 0^5 + 0^5 =1.  The highest 4-digit number is 9999.  9^5 + 9^5 + 9^5 + 9^5 = 236196.  Somewhere in between there is scope for the 5th powers of a 4 digit number to add up to a 4 digit number.

What about 5 digit numbers?  The highest sum of the 5th powers is 295,245. So all 5 digit numbers fall within that range.

And with 6 digit numbers?  The highest sum of the 5th powers is 354294.  So we can cover 6 digit numbers up to that far and no further.

With 7 digit numbers, the highest value is 413343, i.e. the highest sum of the 5th powers of a 7 digit number only adds up to a 6 digit number.  So we only need to loop from 2 to 354294.

```powers = []
grandtotal = 0

for i in range(2,354294):
total = 0
for j in str(i):
total += int(j) ** 5

if total == i:
powers.append(i)

for i in powers:
grandtotal += i

print("powers:", powers)
print("total:", grandtotal)
```

So we loop through our numbers. On each number, loop through the digits, adding up the values of each digit to the power of 5. If the total of those values = our number, then store that number in a list.
At the end of the loop, add up all the values in the list to give the solution to the problem.

Next Page »

Theme: Rubric. Blog at WordPress.com.