Fibonacci sums and unrepeatable digits

Fibonacci sums

Task 1 of week 149 asks us for the first 20 numbers, the sum of whose digits is a Fibonacci number.

It's easy to see that such numbers are quite common, so scanning up from 0 will work quickly enough. So only two functions are required:

sum_of_digits($j)

is_a_fibonacci_number($k)

The first of these doesn't need a sub as it can be done in one line of Perl:

$s = eval(join('+', split(//, $j)))

The 'split' splits the number into an array of digits, the 'join' joins them together with plus signs, and the 'eval' evaluates the string of plusses.

The Fibonacci series is easy to extend, so let's start with just (0, 1), and when we are testing a number for being in the series, we'll extend the series far enough that the last member is more than the number we're testing.

And there you have it.

Unrepeatable digits

Task 2 says given a number base, derive the largest perfect square with no repeated digits and return it as a string.

A number base b is made up from a pallete of b different digits - for example a base 4 number can only contain digits from the set 0, 1, 2 and 3.  So the largest number, base b, with no repeated digits can't be more than b digits long and it comprises the digits b - 1, b -2 ... 0.  For example, if b is 6, the largest number with no repeated digits is 543210. Generalising that, the number (call it z) is:

sum from j = 1 to b - 1 of j * (b**j)

z may not of course be a perfect square, but any number meeting the 'largest perfect square with no repeated digits' criterion must be no more than z.  So let's take the square root of z and truncate it to an integer, m.  We can then scan downwards from m to 1, checking each m**2 to see if it has no repeating digits, and the first one that meets that test is our answer.

This works fine up to b = 16.  Beyond that, Perl (or possibly my computer) cannot represent z as an integer and while my algorithm still holds, converting the initial m**2 to a string fails because it isn't held to sufficient precision. Fortunately Mohammad doesn't ask us to go that high.

Numbers in words and Cardano triplets

 Numbers not containing 'e' when spelt out in words

We are asked for all the numbers under 100 meeting the above condition.

Fortunately, I wrote a numbers-to-words function in Visual Basic about 20 years ago for the purpose of printing the words on cheques (or checks).  Even more fortunately, I was able to find it, and it took little time to translate it into Perl.

Thereafter the solution is trivial.

It's mildly interesting that there are 19 Eban numbers under 100, but the twentieth one is 2000.

Cardano triplets

We are asked to compute the first 5 Cardano triplets (a, b, c) specified by:

Exactly how you order them and thus determine which are the first 5 is not specified, but I took it to mean the first 5 if you order them by a + b + c.

I chose to do this by testing all possible triplets containing numbers up to 100, of which there are 1 million.  By iterating first over c, then b, then a I saved a little time by calculating sqrt(c) only 100 times, and b * sqrt(c) only 10 000 times, but in fact that made little difference: my little Raspberry Pi did it in under 5 seconds with or without my optimisation (maybe because Perl does it anyway).

If we write the definition of the triplets in the form t1 + t2 = 1 it is interesting to look at the pattern. In every case t1 is >1 and t2 is thus negative, for example:

triplet (32, 11, 85) yields 5.10977222864644 + -4.10977222864644 = 1

There are some sets of triplets where their the t1 of each and thus the t2 of each are the same, for example:

(47, 40, 20) :: 6.09016994374947 + -5.09016994374947 = 1

(47, 80, 5) ::   6.09016994374947 + -5.09016994374947 = 1

(47, 20, 80) :: 6.09016994374947 + -5.09016994374947 = 1

(47, 16, 125) :: 6.09016994374947 + -5.09016994374947 = 1

By inspection, it appears that there is only one t1, t2 pair for a given a: that is, any triplets where a is, for example, 47, will have the same t1 and t2.  That implies that they also share b * sqrt(c) and thus b**2 * c, and indeed from the 47 example:

b = 40, c= 20   so b**2 * c = 1600 * 20 = 32000
b = 80, c = 5    so b**2 *c   = 6400 * 5  = 32000

b = 20, c = 80  so b**2 *c   = 400 * 80  = 32000

b = 16, c = 125  so b**2 *c   = 256 * 125  = 32000

From that pattern we can deduce that the following are also Cardano triplets:

(47, 10, 320) and (47, 5, 1280)

(47, 8, 500), (47, 4, 2000), (47, 2, 8000) and (47, 1, 32000)

To get those I took the original set of 47s and successively divided b by 2 and multiplied c by 4 until b was indivisible by 2 (with an integral result), and similarly multiplied b by 2 and divided c by 4 until c was indivisible.

Can we then say that we have identified all the triplets where a = 47?  I postulate that this is the case. 

It must be the case that any factor of 32000 that is a square will generate at least one triplet.  The prime factors of 32000 are 2**8 * 5**3.  So the possible square factors of 32000, ie the possible b**2s, are:

1, 4, 16, 64, 256, 25, 100, 400, 1600, 6400

... and thus the possible values of b are:
1, 2, 4, 8, 16, 5, 10, 20, 40 and 80

... and all of these are in the original set (for which I only searched a, b, and c up to 100), or in the 6 others that I deduced.

But remember that I postulated that all triplets with a = 47 have the same t1 and t2. I can't prove that, although the evidence suggests that it is the case.

So could we do the same for, say a = 48?  At this point I can't see (a) how we could do that, or (b) prove whether in fact any triplets exist for a = 48, but I leave that - as one of my textbooks used to say - as an exercise for the reader.

Chop off their heads and conquer the Pentagon (week 147)

Chop off their heads!

Task 1 this week asks us to write a script to generate first 20 left-truncatable prime numbers in base 10.

A left-truncatable prime (LTP) is a prime number which, in a given base, contains no 0, and if the leading left digit is successively removed, then all resulting numbers are primes.

In principle this is easy enough: generate the primes using a sieve of our good friend Eratosthenes, discard any containing zeroes, and check the remainder for being left-truncatable. 

I chose to exclude single digit primes, ie 2, 3, 5, and 7 because I felt that they did not meet the criterion of 'all resulting numbers are primes', because there are no resulting numbers.

The next 16 LTPs are only 2 digits - they range from 11 to 97, and 4 more - 113, 131, 137 and 167 complete the required 20. The computing time is negligible.

However, had the task asked for the first 1000, or 10 000 LTPs it would be worth applying a little optimisation.  For example, the last digit of an LTP can only be 1, 3 or 7, since any number ending in 0, 2, 4, 5, 6 or 8 is not prime, any ending in 9 will fail the LTP test when only that digit remains.

So doing that I found 1000 LTPs in 3 seconds.

I then tried for 10 000.  It is interesting to note that 33333331 is an LTP, meaning that 33333331, 3333331, 333331, 33331, 3331, 331, 31 are all primes.  Could it be that any number of 3s followed by 1 is an LTP?

It took quite a while because they get steadily sparser, and my computer ran out of memory before getting to the 10 000th. The last one it found was 999 779 337.

Conquer the Pentagon

For task 2 we are asked to write a script to find the first pair of Pentagon Numbers (PNs) whose sum and difference are also a Pentagon Number. Pentagon Numbers are defined by P(n) = n(3n - 1)/2.

My first attempt at this worked, but was rather slow. I generated successive PNs, checked the difference between each new one and all the preceding ones, thus creating a queue of pairs which met the difference criterion.

I then continued generating PNs, and when I got to one that matched one of the queued sums, I had a result. 

Although we were only asked for the first pair, I wanted to find the second pair, so I deleted from the queue the one I had identified as a solution to the task, and also deleted any whose sums were less than that, but despite leaving this running for 15 minutes or so, I found no other qualifying pair. So either my method is faulty, or they are pretty rare, or maybe there is only one.

I mentioned that that was my first attempt.  My second attempt - the one I have submitted - was faster. I used the fact, gleaned from Wikipedia, that any PN gives an integer result in this formula:

(sqrt(24 * PN + 1) + 1) / 6

Using this avoids the need for my queue and the need to identify any PNs beyond the larger member of the pair.  

I am happy to say that both methods came up with the same answer, and despite its improved speed the second method also failed to find a second qualifying pair.

I did worry about applying the Wikipedia formula, because determining that a floating point operation results in an integer result is a little tricky.  For example, in some contexts, this:

$x = 1/3;
$y = $x + $x + $x;
say int($y);

results in 0 rather than the expected 1.  This is, roughly, because 0.333 + 0.333 + 0.333 = 0.999 which truncates to 0.  Happily Perl, on my system at least, does not make that error.

Large primes and curious fractions

The 10001st prime number 

Task 1 for week 146 is a more of a challenge than in recent weeks: we are asked to calculate the 10001st prime number.

The first question to ask is whether the desired number can be represented as an integer in Perl. The largest number which Perl treats as an integer is given by ~1 + 1 and on my computer that is 18446744073709551615 or about 2 * (10**19).  Primes get a bit thin on the ground as they get larger, but I think we can safely postulate that we will find 10**5 of them before we hit that limit.

My 'go to' algorithm for primes is the sieve of Eratosthenes.  Born in 276 BC, Eratosthenes was the first person to estimate the circumference of the earth - remarkably accurately - while managing the largest library in the world, in Alexandria.

His sieve works like this. Write down all the numbers from 2 to however large you want - say n. Take the leftmost number - 2 - and rub out all its multiples - 2, 4, 6, 8 and so on.  Now do it again: the leftmost number remaining is 3, because you rubbed out 2, so now you rub out 3, (not 6 because you already rubbed it out), 9, 15 and so on.  You repeat that until the leftmost number is more than the square root of n, and what remains is a list of all the prime numbers less than or equal to n.

So that gets us all the primes less than n, but what we want is the nth prime number. Eratosthenes unfortunately neglected to address that problem, but we can improvise.

After a bit of experimenting I came up with this:

• Make a sieve out of 1-1000
• Count the primes
• If you haven't reached 10001 of them then ...
• Make a sieve of 1-2000
• and so on

This works fine and gets the answer in 3 seconds on my computer. However if we want, say, the 1000001th prime it takes about 27 seconds, So I improved on it by extending the existing sieve by 1000 rather than making a new one. This makes the code somewhat harder to follow, but it reduced the time for the 10001th prime to under a second and the million and first one to under 7 seconds.

Curious fractions

Task 2 gives us a rather curious binary tree of fractions expressed as a/b, with 1/1 at the root, and asks for a script that can find the parent and grandparent of a given node. 

The relationship between parent and children turns out to be as follows:

• The left child of parent a/b is a/(a+b)
• The right child is (a+b)/b.
• It follows that if a member is < 1 then it is a left child 
• and if > 1 then a right child.

From that, the solution is fairly trivial and no doubt someone will come up with a one-liner.

But do all values for (positive integer) a and b occur in the tree? The answer is no: although it may be possible to calculate parent and grandparents for a given a/b, it does not follow that there is a path back to the root 1/1.  In fact, I think that a and b have to be mutually prime - ie share no prime factor - for a/b to appear in the tree.