You are currently browsing the category archive for the ‘computing’ category.

**The Van Eck Sequence**

This is a nice sequence as discussed in the Numberphile video above. There are only 2 rules:

- If you have not seen the number in the sequence before, add a 0 to the sequence.
- If you have seen the number in the sequence before, count how long since you last saw it.

You start with a 0.

0

You have never seen a 0 before, so the next number is 0.

00

You have seen a 0 before, and it was 1 step ago, so the next number is 1.

001

You have never seen a 1 before, so the next number is 0.

0010

You have seen a 0 before, it was 2 steps ago, so the next number is 2.

00102.

etc.

I can run a quick Python program (adapted from the entry in the Online Encyclopedia of Integer Sequences here) to find the first 100 terms.

```
```A181391 = [0, 0]

for n in range(1, 10**2):

for m in range(n-1, -1, -1):

if A181391[m] == A181391[n]:

A181391.append(n-m)

break

else:

A181391.append(0)

print(A181391)

This returns:

[0, 0, 1, 0, 2, 0, 2, 2, 1, 6, 0, 5, 0, 2, 6, 5, 4, 0, 5, 3, 0, 3, 2, 9, 0, 4, 9, 3, 6, 14, 0, 6, 3, 5, 15, 0, 5, 3, 5, 2, 17, 0, 6, 11, 0, 3, 8, 0, 3, 3, 1, 42, 0, 5, 15, 20, 0, 4, 32, 0, 3, 11, 18, 0, 4, 7, 0, 3, 7, 3, 2, 31, 0, 6, 31, 3, 6, 3, 2, 8, 33, 0, 9, 56, 0, 3, 8, 7, 19, 0, 5, 37, 0, 3, 8, 8, 1, 46, 0, 6, 23]

I then assigned each term an x coordinate value, i.e.:

0 , 0

1 , 0

2 , 1

3 , 0

4 , 2

5 , 0

6 , 2

7 , 2

8 , 1

9 , 6

10 , 0

11 , 5

12 , 0

13 , 2

14 , 6

15 , 5

16 , 4

17 , 0

18 , 5

19 , 3

20 , 0

etc.

This means that you can then plot the sequence as a line graph, with the y values corresponding to the sequence terms. As you can see, every time we hit a new peak the following value is 0, leading to the peaks and troughs seen below:

Let’s extend the sequence to the first 1000 terms:

We can see that the line y = x provides a reasonably good upper bound for this data:

But it is not known if every number would actually appear in the sequence somewhere – so this bound may not hold for larger values.

**Length of steps before new numbers appear.**

We can also investigate how long we have to wait to see each number for the first time by running the following Python code:

```
```A181391 = [0, 0]

for n in range(1, 10**3):

for m in range(n-1, -1, -1):

if A181391[m] == A181391[n]:

A181391.append(n-m)

break

else:

A181391.append(0)

for m in range(1,50):

if A181391[n]==m:

print(m, ",", n+1)

break

This returns the following data:

1 , 3

2 , 5

6 , 10

5 , 12

4 , 17

3 , 20

9 , 24

14 , 30

15 , 35

17 , 41

11 , 44

8 , 47

42 , 52

20 , 56

32 , 59

18 , 63

7 , 66

31 , 72

33 , 81

19 , 89

etc.

The first coordinate tells us the number we are interested in, and the second number tells us how long we have to wait in the sequence until it appears. So (1 , 3) means that we have to wait until 3 terms in the sequence to see the number 1 for the first time.

Plotting this for numbers 1-50 on a graph returns the following:

So, we can see (for example that we wait 66 terms to first see a 7, and 173 terms to first see a 12. There seems to be a general trend that as the numbers get larger we have to wait longer to see them. Testing this with a linear regression we can see a weak to moderate correlation:

Checking for the numbers up to 300 we get the following:

For example this shows that we have to wait 9700 terms until we see the number 254 for the first time. Testing this with a linear correlation we have a weaker positive correlation than previously.

So, a nice and quick investigation using a combination of sequences, coding, graphing and regression, with lots of areas this could be developed further.

Computers can brute force a lot of simple mathematical problems, so I thought I’d try and write some code to solve some of them. In nearly all these cases there’s probably a more elegant way of coding the problem – but these all do the job! You can run all of these with a Python editor such as Repl.it. Just copy and paste the below code and see what happens.

1) **Happy Numbers.**

Happy numbers are defined by the rule that you start with any positive integer, square each of the digits then add them together. Now do the same with the new number. Happy numbers will eventually spiral down to a number of 1. Numbers that don’t eventually reach 1 are called unhappy numbers.

As an example, say we start with the number 23. Next we do 2²+3² = 13. Now, 1²+3² = 10. Now 1²+0² = 1. 23 is therefore a happy number.

k= int(input("type a 2 digit number "))

a = int(k%10)

c = int(k//100)

b = int(k//10 -10*c)

print (a**2+b**2+c**2)

```
```for k in range (1,20):

` k = a**2+b**2 + c**2`

a = int(k%10)

c = int(k//100)

b = int(k//10 -10*c)

print (a**2+b**2+c**2)

2) **Sum of 3 cubes**

Most (though not all) numbers can be written as the sum of 3 cubes. For example:

1^{3} + 2^{3} + 2^{3} = 17. Therefore 17 can be written as the sum of 3 cubes.

This program allows you to see all the combinations possible when using the integers -10 to 10 and trying to make all the numbers up to 29.

for k in range(1,30):

```
```

` for a in range(-10, 10):`

for b in range(-10,10):

for c in range (-10, 10):

if a**3+b**3+c**3 == k :

print(k,a,b,c)

3) **Narcissistic Numbers**

A 3 digit narcissistic number is defined as one which the sum of the cubes of its digits equal the original number. This program allows you to see all 3 digit narcissistic numbers.

```
```for a in range (0,10):

for b in range(0, 10):

for c in range(0,10):

if a**3 + b**3 + c**3 ==100*a + 10*b + c:

print(int(100*a + 10*b + c))

4) **Pythagorean triples**

Pythagorean triples are integer solutions to Pythagoras’ Theorem. For example:

3^{2} + 4^{2} = 5^{2} is an integer solution to Pythagoras’ Theorem.

This code allows you to find all integer solutions to Pythagoras’ Theorem for the numbers in the range you specify.

```
```k = 100

`for a in range(1, k):`

for b in range(1,k):

for c in range (1, 2*k):

if a**2+b**2==c**2:

print(a,b,c)

5) **Perfect Numbers**

Perfect numbers are numbers whose proper factors (factors excluding the number itself) add to the number. This is easier to see with an example.

6 is a perfect number because its proper factors are 1,2,3 and 1+2+3 = 6

8 is not a perfect number because its proper factors are 1,2,4 and 1+2+4 = 7

Perfect numbers have been known about for about 2000 years – however they are exceptionally rare. The first 4 perfect numbers are 6, 28, 496, 8128. These were all known to the Greeks. The next perfect number wasn’t discovered until around 1500 years later – and not surprisingly as it’s 33,550,336.

The code below will find all the perfect numbers less than 10,000.

```
```for n in range(1,10000):

list = []

for i in range (1,n):

if n%i ==0:

list.append(i)

if sum(list)==n:

print(n)

**Friendly Numbers**

Friendly numbers are numbers which share a relationship with other numbers. They require the use of σ(a) which is called the divisor function and means the addition of all the factors of a. For example σ(7) = 1 + 7 = 8 and σ(10) = 1 +2 +5 + 10 = 18.

Friendly numbers therefore satisfy:

σ(a)/a = σ(b)/b

As an example:

σ(6) / 6 = (1+2+3+6) / 6 = 2,

σ(28) / 28 = (1+2+4+7+14+28) / 28 = 2

σ(496)/496 = (1+2+4+8+16+31+62+124+248+496)/496 = 2

Therefore 28 and 6 are friendly numbers because they share a common relationship.

This code will help find some Friendly numbers (though these are very difficult to find, as we need to check against every other integer until we find a relationship).

The code below will find some Friendly numbers less than 200, and their friendly pair less than 5000:

for n in range(1,5000):

list = []

```
``` for i in range (1,n+1):

if n%i ==0:

list.append(i)

Result1 = sum(list)

for m in range(1,200):

list2 = []

for j in range (1,m+1):

if m%j ==0:

list2.append(j)

Result2 = sum(list2)

` if Result2/m ==Result1/n:`

if n != m:

print(n,m)

**Hailstone numbers**

Hailstone numbers are created by the following rules:

if n is even: divide by 2

if n is odd: times by 3 and add 1

We can then generate a sequence from any starting number. For example, starting with 10:

10, 5, 16, 8, 4, 2, 1, 4, 2, 1…

we can see that this sequence loops into an infinitely repeating 4,2,1 sequence. Trying another number, say 58:

58, 29, 88, 44, 22, 11, 34, 17, 52, 26, 13, 40, 20, 10, 5, 16, 8, 4, 2, 1, 4, 2, 1…

and we see the same loop of 4,2,1.

The question is, does every number end in this loop? Well, we don’t know. Every number mathematicians have checked do indeed lead to this loop, but that is not a proof. Perhaps there is a counter-example, we just haven’t found it yet.

Run the code below, and by changing the value of n you can see how quickly the number enters the 4,2,1 loop.

n = 300

for k in range(1,40):

```
```

` if n%2 ==0:`

print(n/2)

n =n/2

elif n%2 !=0:

print(3*n+1)

n =3*n+1

**Generating the Golden ratio**

The Golden ratio can be approximated by dividing any 2 successive terms of the Fibonacci sequence. As we divide ever larger successive terms we get a better approximation for the Golden ratio. This code returns successive terms of the Fibonacci sequence and the corresponding approximation for the Golden ratio.

a = 0

b = 1

print(a)

print(b)

for k in range(1,30):

```
``` a = a+b

b = a+b

` print(a,b, b/a)`

**Partial sums**

We can use programs to see if sums to infinity converge. For example with the sequence 1/n, if I add the terms together I get: 1/1 + 1/2 + 1/3 + 1/4…In this case the series (surprisingly) diverges. The code below shows that the sum of the sequence 1/n^{2} converges to a number (pi^{2}/6).

```
```

`list = []`

for n in range(1,100):

n = 1/(n**2)

list.append(n)

print(sum(list))

**Returning to 6174**

This is a nice number trick. You take any 4 digit number, then rearrange the digits so that you make the largest number possible and also the smallest number possible. You then take away the smallest number from the largest number, and then start again. For example with the number 6785, the largest number we can make is 8765 and the smallest is 5678. So we do 8765 – 5678 = 3087. We then carry on with the same method. Eventually we will arrive at the number 6174!

k= int(input("type a 4 digit number "))

a = int(k%10)

d = int(k//1000)

c = int(k//100 - 10*d)

b = int(k//10 -10*c-100*d)

```
```for n in range(1,10):

list = []

list = [a,b,c,d]

list.sort()

a = list[0]

d = list[3]

c = list[2]

b = list[1]

print(1000*d+100*c+10*b+a -1000*a-100*b-10*c-d)

k = int(1000*d+100*c+10*b+a -1000*a-100*b-10*c-d)

a = int(k%10)

d = int(k//1000)

c = int(k//100 - 10*d)

b = int(k//10 -10*c-100*d)

list = []

list = [a,b,c,d]

list.sort()

a = list[0]

d = list[3]

c = list[2]

b = list[1]

` print(1000*d+100*c+10*b+a -1000*a-100*b-10*c-d)`

**Maximising the volume of a cuboid**

If we take a cuboid of length n, and cut squares of size x from the corner, what value of x will give the maximum volume? This code will look at initial squares of size 10×10 up to 90×90 and find the value of x for each which give the maximum volume.

def compute():

```
``` list1=[]

k=6

z = int(0.5*a*10**k)

for x in range(1,z):

list1.append((10*a-2*x/10**(k-1))*(10*a-2*x/10**(k-1))*(x/10**(k-1)))

print("length of original side is, ", 10*a)

y= max(list1)

print("maximum volume is, ", max(list1))

q = list1.index(y)

print("length of square removed from corner is, ", (q+1)/10**(k-1))

`for a in range(1,10):`

print(compute())

**Stacking cannonballs – solving maths with code**

Numberphile have recently done a video looking at the maths behind stacking cannonballs – so in this post I’ll look at the code needed to solve this problem.

**Triangular based pyramid.**

A triangular based pyramid would have:

1 ball on the top layer

1 + 3 balls on the second layer

1 + 3 + 6 balls on the third layer

1 + 3 + 6 + 10 balls on the fourth layer.

Therefore a triangular based pyramid is based on the sum of the first n triangular numbers.

The formula for the triangular numbers is:

and the formula for the sum of the first n triangular numbers is:

We can simplify this by using the identity for the sum of the first n square numbers and also the identity for the sum of the first n natural numbers:

Therefore:

and the question we want to find out is whether there is triangular based pyramid with a certain number of cannonballs which can be rearranged into a triangular number i.e.:

here n and m can be any natural number. For example if we choose n = 3 and m = 4 we see that we have the following:

Therefore we can have a triangular pyramid of height 3, which has 10 cannonballs. There 10 cannonballs can then be rearranged into a triangular number.

**Square based pyramids and above.**

For a square based pyramid we would have:

1 ball on the top layer

1 + 4 balls on the second layer

1 + 4 + 9 balls on the third layer

1 + 4 + 9 + 16 balls on the fourth layer.

This is the sum of the first n square numbers. So the formula for the square numbers is:

and the sum of the first n square numbers is:

**For a pentagonal based pyramid we have:**

1 ball on the top layer

1 + 5 balls on the second layer

1 + 5 + 12 balls on the third layer

1 + 5 + 12 + 22 balls on the fourth layer.

This is the sum of the first n pentagonal numbers. So the formula for the pentagonal numbers is:

and the formula for the first n pentagonal numbers is:

**For a hexagonal based pyramid we have:**

The formula for the first n hexagonal numbers:

and the formula for the sum of the first n hexagonal numbers:

For a **k-agon based pyramid we have**

and the formula for the sum of the first n k-agon numbers:

Therefore the general case is to ask if a k-agonal pyramid can be rearranged into a k-agon number i.e:

**Computers to the rescue**

We can then use some coding to brute force some solutions by running through large numbers of integers and seeing if any values give a solution. Here is the Python code. Type it (taking care with the spacing) into a Python editor and you can run it yourself.

You can then change the k range to check larger k-agons and also change the range for a and b. Running this we can find the following. (The first number is the value of k, the second the height of a k-agonal pyramid, the third number the k-agon number and the last number the number of cannonballs used).

**Solutions:**

3 , 3 , 4 , 10

3 , 8 , 15 , 120

3 , 20 , 55 , 1540

3 , 34 , 119 , 7140

4 , 24 , 70 , 4900

6 , 11 , 22 , 946

8 , 10 , 19 , 1045

8 , 18 , 45 , 5985

10 , 5 , 7 , 175

11 , 25 , 73 , 23725

14 , 6 , 9 , 441

14 , 46 , 181 , 195661

17 , 73 , 361 , 975061

20 , 106 , 631 , 3578401

23 , 145 , 1009 , 10680265

26 , 190 , 1513 , 27453385

29 , 241 , 2161 , 63016921

30 , 17 , 41 , 23001

32 , 298 , 2971 , 132361021

35 , 361 , 3961 , 258815701

38 , 430 , 5149 , 477132085

41 , 204 , 1683 , 55202400

41 , 505 , 6553 , 837244045

43 , 33 , 110 , 245905

44 , 586 , 8191 , 1408778281

50 , 34 , 115 , 314755

88 , 15 , 34 , 48280

145, 162, 1191, 101337426

276, 26, 77, 801801)

322, 28, 86, 1169686

823, 113, 694, 197427385

2378, 103, 604, 432684460

31265, 259, 2407, 90525801730

For example we can see a graphical representation of this. When k is 6, we have a hexagonal pyramid with height 11 or the 22nd hexagonal number – both of which give a solution of 946. These are all the solutions I can find – can you find any others? Leave a comment below if you do find any others and I’ll add them to the list!

**What’s so special about 277777788888899?**

Numberphile have just done a nice video which combines mathematics and computer programing. The challenge is to choose any number (say 347)

Then we do 3x4x7 = 84

next we do 8×4 = 32

next we do 3×2 = 6.

And when we get to a single digit number then we have finished. It took 3 steps to get from 347 to a single digit number, therefore 347 has a *persistence* of 3. The challenge is to find a number with as big a persistence as possible. The current world record is 277777788888899 which is the smallest number with a persistence of 11. No numbers with a persistence of greater than 11 have ever been found. In the video Matt writes a Python program to check this, though I tried to make my own version below. It’s not very sophisticated, but it gets the job done (with a small glitch of returning a 0 followed by 1s when it should just return 0s!)

The full code should be available to run here, or download here. If you run the program above in an online Python site like repl.it you can choose any number you like as see what its persistence is.

If you find any number that hasn’t gone to a single digit after 11 rounds, you’ve found a new world record persistence!

To very briefly explain the code used above:

We start by defining “result” as 1. We then have some add any integer number on the screen (let’s use 347). We then do 347 mod 10 (number % 10) which gives 7, and do result (which is 1) multiplied by 7. We then do 347 divided by 10 ignoring remainders (number//10). This gives 34.

We then start the process again. 34 mod 10 = 4. So now we have 1 x 7 x 4. Next we do 34 divided by 10 ignoring remainders which gives 3. Last we do 3 mod 10 = 3. So we have 1 x 7 x 4 x 3. If we carried on the loop we would next have 3/10 = 0 ignoring remainders, therefore our loop would stop.

The program then defines the *new* starting number as 7x4x3 = 84 and then starts again. So, a nice use of mathematics and computing – see what levels of persistence you can find!