Tag Archives: prime

Problem Twenty Seven

18 Nov

Here’s problem twenty seven:

Question: Find the product of the coefficients, a and b, for the quadratic expression that produces the maximum number of primes for consecutive values of n, starting with n = 0.

I’m not exactly sure why, but this was a pretty fun problem to do. Basically I just had to brute force my way through this problem, and it actually ran much faster than I anticipated.

Take a look:

def isPrime(x): 
	divisor = 3 
	limit = int(x**0.5) + 1

	while divisor < limit:
		if x % divisor == 0:
			return False
		else:
			divisor += 2 

	return True

result = [0, 0, 0] # a, b, counter

for a in range(-999, 1000):
	for b in range(-999, 1000):
		n = 0
		while True:
			x = n*n + a*n + b
			if x > 0 and isPrime(x):
				n += 1	
			else:
				if n > result[2]:
					result = [a, b, n]
				break

print(result[0]*result[1])

As you can see, I reused the ever-so-popular isPrime() function for this problem. We use two for loops to test every combination of coefficients, and initialize n, the independent variable of the quadratic, to 0. We then continue to generate values of the quadratic until we reach one that isn’t prime. The length of the chain generated is tested to see if it’s greater then the previous longest chain, and if it is we overwrite the variable result.

Answer: -59231

Advertisements

Problem Twelve

7 Oct

Here’s problem twelve:

Question: What is the value of the first triangle number to have over five hundred divisors?

This problem has two parts: first is to find some way to generate triangle numbers, and second is to find the number of divisors of that number. As I soon realized, brute forcing the factors of a number takes a very long time, so a more efficient algorithm was definitely needed there. I also imported a module to help solve this problem, but its job was trivial in finding the solution, and I used it only because coding the function myself would involve more python list manipulation than  mathematics.

Here’s the code:

from collections import Counter

def primeFactors(n): # takes a number n and returns a list of prime factors
	factors = []
	d = 2

	while d  500:
		print(tn)
		break

	n += 1

Let’s look at the bottom first.

n = 1 # incrementer for sequence of triangle numbers
while True:
	tn = (n*(n+1))/2 # formula for the nth triangle number

	factors = countFactors(primeFactors(tn))

	if factors > 500:
		print(tn)
		break

	n += 1

This first part of this snippet of code generates triangle numbers. The formula for the nth triangle number is (n*(n+1))/2, as seen in the code. Then something strange happens. To find the number of factors of tn, we put it into a function inside a function. So let’s look at the innermost function first, primeFactors().

def primeFactors(n): # takes a number n and returns a list of prime factors
	factors = []
	d = 2

	while d 		if n % d == 0:
			factors.append(d)
			n /= d
			d = 2
		elif d == n:
			factors.append(d)
			break
		else:
			d += 1

	return factors

This code is actually taken from the solution for problem three, where we needed to find the largest prime factor of a number. It’s essentially brute forcing the prime factors of a number, which is much faster than brute forcing all the factors of a number. For a more detailed explanation of this function, go to the post for problem three.

At this point you may be wondering why we’re finding the prime factorization of a number, when the question requires the number of factors of a number. Well, brute forcing the number of factors of a number takes a considerably long time, so a more clever solution was needed. Using the prime factorization of a number, it is possible to find the number of factors of a number (this may or may not be called the divisor function). Basically, the prime factorization of a number is expressed as a product of its prime factors and their exponents; for example, 168 = 23 * 31 * 71. Once here we can discard the bases, as we only need to deal with the exponents. The number of factors is found by taking these exponents, adding one of each of them, and multiplying these resulting values. The product will be the number of factors of the initial number. So the number of factors of 168 would be (3+1)(1+1)(1+1) = 16. Who would’ve thought?

def countFactors(n): # takes a list of prime factors and returns the total number of factors for the product of the prime factors
	n = dict(Counter(n))
	n = list(n.values())

	factors = 1
	for i in n:
		i += 1
		factors *= i

	return factors

Now that we know this, we can look at the outer function, or countFactors(). Remember that primeFactors(tn) returns a list of the prime factors of tn, which is then put into countFactors(). Here’s where the imported Counter() function comes in. It takes a list, and returns a dictionary with every different value that appeared in that list, with the number of times it occurred. In other words, it returns Counter{value:occurrences}. dict() changes this counter object (dictionary subclass) into a dictionary, which we can take the values of to put into a list. Remember that values() returns a list of the values of a dictionary, which in our case is the occurrences, or the value of the exponents when the prime factorization is written out. Then we take the product of each of these ‘exponents’ plus one, in this case elements of the list n. This value is returned, which is the number of factors of the initial number tn.

Now that we’ve examined the inner workings of the functions, we can go back to look at the run code one more time.

n = 1 # incrementer for sequence of triangle numbers
while True:
	tn = (n*(n+1))/2 # formula for the nth triangle number

	factors = countFactors(primeFactors(tn))

	if factors > 500:
		print(tn)
		break

	n += 1

So again, a new tn, or triangle number is generated until we decide to break the loop. The number of factors of tn are found using the relatively efficient technique of countFactors(primeFactors(tn)), and if this number is greater than 500, we print the number and break the loop.

Answer: 76576500

Problem Eleven

6 Oct

On to problem eleven.

Question: What is the greatest product of four adjacent numbers in any direction (up, down, left, right, or diagonally) in the 20×20 grid?

Ah, this question was actually pretty difficult. It took me a while to figure out a way to be able to find vertically and diagonally adjacent numbers in the grid, and my eventual solution was to use a list of lists. Essentially, each row in the grid is represented as a list of numbers, and 20 of these lists are put into another list in the order of the grid. After this was done, the rest entailed manipulating individual elements in the lists.

It might be clearer to see the code, so here it is:

nums = [
[ 8,  2, 22, 97, 38, 15,  0, 40,  0, 75,  4,  5,  7, 78, 52, 12, 50, 77, 91,  8], 
[49, 49, 99, 40, 17, 81, 18, 57, 60, 87, 17, 40, 98, 43, 69, 48,  4, 56, 62,  0], 
[81, 49, 31, 73, 55, 79, 14, 29, 93, 71, 40, 67, 53, 88, 30,  3, 49, 13, 36, 65], 
[52, 70, 95, 23,  4, 60, 11, 42, 69, 24, 68, 56,  1, 32, 56, 71, 37,  2, 36, 91], 
[22, 31, 16, 71, 51, 67, 63, 89, 41, 92, 36, 54, 22, 40, 40, 28, 66, 33, 13, 80], 
[24, 47, 32, 60, 99,  3, 45,  2, 44, 75, 33, 53, 78, 36, 84, 20, 35, 17, 12, 50], 
[32, 98, 81, 28, 64, 23, 67, 10, 26, 38, 40, 67, 59, 54, 70, 66, 18, 38, 64, 70], 
[67, 26, 20, 68,  2, 62, 12, 20, 95, 63, 94, 39, 63,  8, 40, 91, 66, 49, 94, 21], 
[24, 55, 58,  5, 66, 73, 99, 26, 97, 17, 78, 78, 96, 83, 14, 88, 34, 89, 63, 72], 
[21, 36, 23,  9, 75,  0, 76, 44, 20, 45, 35, 14,  0, 61, 33, 97, 34, 31, 33, 95], 
[78, 17, 53, 28, 22, 75, 31, 67, 15, 94,  3, 80,  4, 62, 16, 14,  9, 53, 56, 92], 
[16, 39,  5, 42, 96, 35, 31, 47, 55, 58, 88, 24,  0, 17, 54, 24, 36, 29, 85, 57], 
[86, 56,  0, 48, 35, 71, 89,  7,  5, 44, 44, 37, 44, 60, 21, 58, 51, 54, 17, 58], 
[19, 80, 81, 68,  5, 94, 47, 69, 28, 73, 92, 13, 86, 52, 17, 77,  4, 89, 55, 40], 
[ 4, 52,  8, 83, 97, 35, 99, 16,  7, 97, 57, 32, 16, 26, 26, 79, 33, 27, 98, 66], 
[88, 36, 68, 87, 57, 62, 20, 72,  3, 46, 33, 67, 46, 55, 12, 32, 63, 93, 53, 69], 
[ 4, 42, 16, 73, 38, 25, 39, 11, 24, 94, 72, 18,  8, 46, 29, 32, 40, 62, 76, 36], 
[20, 69, 36, 41, 72, 30, 23, 88, 34, 62, 99, 69, 82, 67, 59, 85, 74,  4, 36, 16], 
[20, 73, 35, 29, 78, 31, 90,  1, 74, 31, 49, 71, 48, 86, 81, 16, 23, 57,  5, 54], 
[ 1, 70, 54, 71, 83, 51, 54, 69, 16, 92, 33, 48, 61, 43, 52,  1, 89, 19, 67, 48]]

result = 0
factors = []

# finds highest in rows
for row in nums:
	for i in range(0, 17):
		product = row[i]*row[i+1]*row[i+2]*row[i+3] 

		if product > result:
			result = product
			factors = [row[i],row[i+1],row[i+2],row[i+3]] 
		i += 1
	
# finds highest in columns
for i in range(0, 17):
	for j in range(0, 20):
		product = nums[i][j]*nums[i+1][j]*nums[i+2][j]*nums[i+3][j]	

		if product > result:
			result = product
			factors = [nums[i][j],nums[i+1][j],nums[i+2][j],nums[i+3][j]]

# finds diagonal to the right
for i in range(0, 17):
	for j in range(0, 17):
		product = nums[i][j]*nums[i+1][j+1]*nums[i+2][j+2]*nums[i+3][j+3]

		if product > result:
			result = product
			factors = [nums[i][j],nums[i+1][j+1],nums[i+2][j+2],nums[i+3][j+3]]


# finds diagonal to the left
for i in range(16, -1, -1):
	for j in range(19, 2, -1):
		product = nums[i][j]*nums[i+1][j-1]*nums[i+2][j-2]*nums[i+3][j-3]

		if product > result:
			result = product
			factors = [nums[i][j],nums[i+1][j-1],nums[i+2][j-2],nums[i+3][j-3]]


print(result, factors)

The factors list is not necessary to solve the problem, but I prefer to keep it in there if only to see what exactly the computer is calculating. This code has four parts, which check the products of the adjacent horizontal, vertical, left-diagonal, and right-diagonal digits respectively.

# finds highest in rows
for row in nums:
	for i in range(0, 17):
		product = row[i]*row[i+1]*row[i+2]*row[i+3] 

		if product > result:
			result = product
			factors = [row[i],row[i+1],row[i+2],row[i+3]] 
		i += 1

Let’s look at line 26’s block first. The for loop finds each row in the grid, and puts it (as a list) into the variable row. For each row, i is the variable used to retrieve the ith element in the row (list). Since we’re looking for groups of four digits, the last four digits we need to test in each row are in positions 16, 17, 18, and 19 (remember, list positions start at 0). So i only has to go up to 16, since we find next three elements after i to multiply them. Every time four digits are multiplied, we check it against the previous largest product, and if it’s larger, the variable is overwritten with the that number. The four digits that are factors of the product are also stored.

# finds highest in columns
for i in range(0, 17):
	for j in range(0, 20):
		product = nums[i][j]*nums[i+1][j]*nums[i+2][j]*nums[i+3][j]	

		if product > result:
			result = product
			factors = [nums[i][j],nums[i+1][j],nums[i+2][j],nums[i+3][j]]

Next we tested vertically adjacent digits. In this case, i is the variable used to retrieve the ith element from nums, or in other words it retrieves the ith row in the grid. j is the variable used to find the jth element of that row. From a similar logic as before, the last group of rows we need to test are rows 16, 17, 18, and 19. Therefore, i only needs to go up to 16. However, j needs to go up to 19 because we need to go all the way to the end of each row and multiply the vertical column of digits below it. Four vertical digits are found by retrieving four adjacent rows (numbered i to i+3), then finding the same (jth) element of each. These four digits are then multiplied, and tested against product.

# finds diagonal to the right
for i in range(0, 17):
	for j in range(0, 17):
		product = nums[i][j]*nums[i+1][j+1]*nums[i+2][j+2]*nums[i+3][j+3]

		if product > result:
			result = product
			factors = [nums[i][j],nums[i+1][j+1],nums[i+2][j+2],nums[i+3][j+3]]

Diagonal to the right means that starting from any digit and going down diagonally, the last digit in that diagonal will be to the right of the initial digit. Again, i is used to find the ith row, and j is the jth element in that row. The ranges for i and j will become more obvious after explaining how the diagonal is found. So the first digit we find is in the ith row, at position j. The second digit is a row below that, at the i+1th row, and one position to the right, at j+1. This is repeated until the i+3rd row, where we will have a four digits that are diagonal to each other. Since the last digit is 3 down and 3 to the right of the initial digit, i (the row) and j (the element position) only have to go up to 16. Each product is checked against product.

# finds diagonal to the left
for i in range(16, -1, -1):
	for j in range(19, 2, -1):
		product = nums[i][j]*nums[i+1][j-1]*nums[i+2][j-2]*nums[i+3][j-3]

		if product > result:
			result = product
			factors = [nums[i][j],nums[i+1][j-1],nums[i+2][j-2],nums[i+3][j-3]]

Finding the diagonal to the left is a bit strange. We start from the 4th row from the bottom (nums[16]), and work our way up to the top row (nums[0]). Therefore, i, or the ith row, starts at 16 and goes to 0. Following a similar logic as before, since the final digit in the left-diagonal is 3 down and 3 to the left from the initial digit, j, or the jth element must start from 19 (the end of the row), up to 3 (4 from the left of the row). We test the left-diagonals starting from the right, and moving to the left, or incrementing j by -1. Each product is checked against product.

So after we test the horizontal, vertical, and diagonal products, we will have exhausted all possibilities and the variable product will be storing the largest product it found, and so we print.

Answer: 70600674 [89, 94, 97, 87]

Problem Ten

8 Sep

Here’s problem ten:

Question: Find the sum of all the primes below two million.

This problem was fairly easy, at least when implementing my inelegant brute force solution. Ninety percent of the work was already done, in the form of a prime testing function that I created for problem seven. I just copied and pasted the function, and wrote a short while loop to find the solution.

Here’s the code:

def isPrime(x): 
	
	divisor = 3 
	limit = int(x**0.5) + 1

	while divisor < limit:
		if x % divisor == 0:
			return False
		else:
			divisor += 2 

	return True

i = 3
sum = 0

while i < 2000000:
	if isPrime(i):
		sum += i
	i += 2 # add two to skip evens

print(sum+2)

I’m not going to explain the isPrime() function again; if you’re curious about how it works, you can read my post about question seven. We initialize i to 3 and set sum to 0. While the current value of i is less than 2,000,000, we test to see if it’s prime. If true, we add it to sum. After prime testing we add two to i, which allows us to skip the even numbers. Finally we print sum, but since we skipped the prime number 2, we add it to sum before printing.

This was relatively inefficient, and took around 16 seconds to run. There’s an algorithm called the sieve of Eratosthenes, which finds primes remarkably fast. However, school just started, so in the interest of time (for doing homework D:), I didn’t research this method much. But if the ancient Greeks could come up with this method, how complicated could it be? Google it if you’re curious and see what you can come up with.

Answer: 142913828922

Problem Seven

30 Aug

On to problem seven.

Question: What is the 10,001st prime number?

This problem was more tricky that I had first expected. I came up with a solution after an hour of tinkering around, but it took around 110 seconds to run (much too long). I probably should’ve done some research for calculating primes beforehand though, since there was a simple trick that made the result 100 times faster (it finished less than a second!). I’ll discuss this theorem/postulate later.

This was also the first problem I had to write a new function for, in order to keep it somewhat readable.

Let’s think how to solve this problem. We’ll need some way to determine if a number is prime, and also a counter to keep track of when we reach the 10,001st prime. Since we don’t know how many times we need to iterate, we use a while true loop.

Here’s the code:

def isPrime(x): 
	
	divisor = 3
	limit = int(x**0.5) + 1 # the highest factor we need to test is the square root of the number

	while divisor < limit: 
		if x % divisor == 0:
			return False
		else:
			divisor += 2 

	return True 

counter = 1 # keeps track of number of primes, 3 is second
i = 3 # start testing at 3 

while counter != 10001: 
	if isPrime(i):
		counter += 1 

	i += 2 #iterates over even numbers

print(i - 2)

Let’s look from line 14-15 first. We start the counter at one rather than zero, since we don’t include two as a prime, so three will count as the second prime we find. We keep track of the current number we’re on with the variable i. Moving down, we see a while loop that runs as long as the counter is under 10001, meaning that we haven’t found the 10001st prime yet. Inside the loop, we check if i is a prime using the function isPrime(), and if it is, we add one to counter. If not, we add two to i and try again. The reason for adding two is that every other number is even, which will not be a prime, so we can completely skip calculating those.

When it’s done running, we print the curious line i – 2. Why? Because when we find the 10001st prime and counter is set to 10001, before the loop stops running, we add two to i, or the current number (which is the 10001st prime). Printing i – 2 is simpler than dealing with this using if statements.

So now let’s take a look at the bulk of the code, the isPrime() function.

We set the initial divisor to three, since anything will divide into one and we already sifted out the even numbers. We have an upper limit that’s a bit esoteric, which I will explain in a second. Inside the while loop, we keep testing to see if the number we passed in as a parameter is divisible by the divisor. We increment the divisor every time until it reaches a limit, and then it stops. If the number is divisible by any divisor under this limit, we know that it’s not prime and it returns False. If it’s not divisible by any of the divisors under the limit, then the while loop ends and we return True (it’s a prime).

Now for the limit variable. Initially I used no (defined) limit, and I would test all the numbers below the parameter to see if it was a factor. However, this was very slow and made the code run for more than 100 seconds. In the answers forum, I saw that someone had used a similar technique, but with the square root of the number as a limit. I googled this and found the mathematical reason. Basically, factors of a number come in pairs. One will be smaller, and the other will be larger. The only time when the two factors will be equal are when they are square roots of that number. So essentially, the largest factor of a number you need to worry about is the square root of that number. Any factor higher than that will have a corresponding lower factor, which has already been ruled out. That’s why the limit is set to the square root of the parameter (plus one).

Answer: 104743

Problem Three

25 Aug

Here’s problem three.

Question: What is the largest prime factor of the number 600851475143?

Ah, I must admit that I cheated a bit on this problem. My solution found the factors of that very large number, but not actually the prime factors. Simply by luck, the largest factor I found was a prime number, so it did work out in the end. However, a lot of people in the answers forum imported modules for primality testing, and one person even used a Mathematica function to list the prime factorization! So I guess to this end I’m cheating slightly less. 😛 (Correction to this statement at the bottom!)

So my solution to this problem was to test if a integer (incrementing every iteration) divided into the number cleanly. If true, I added that integer into a list called factors. I then divided the original number by that factor so that we could just put this code into a while loop.

Here’s the code: (d for divisor, facs for factors)

n = 600851475143
facs = [] # stores factors 
d = 2 # divisor must start at two

while d <= n:
	if n % d == 0: # if divisible, store factor in list and divide n by that factor
		facs.append(d)
		n /= d
		d = 2
	elif d == n: # last factor will be equal, so add to list and break
		facs.append(d)
		break
	else:
		d += 1 # if not divisible, increment divisor by one and test again

print(facs)

The number 600851475143 is stored in n, facs is an empty list, and d is set to 2 (division by zero if set to zero and infinite loop if set to one). Inside the while loop, it tests two things. If n can be divided by d (the current integer), it adds d to the list and sets n to the quotient of n and d. It also resets d (the divisor) to 2, since we are going to be testing a new number to divide in the next iteration. If it’s not divisible (and they are not equal), it increments d by one and repeats. When it reaches the last factor, d and n will be equal, so we add d the the list and break out of the loop. Finally we print the list.

So how would you solve this mathematically? I definitely would not like to, since the only way that I envision you could solve this is by brute force or some complex algorithm. I did google some mathematical algorithms initially to see if I could directly convert the math into python-speak, but they were almost all rooted in mathematical concepts beyond what I’ve learned, so I couldn’t really implement any of them. Some of these include the Pollard-Strassen method, the general number field sieve, and elliptic curve primality testing.

On the WolframMathworld website, I learned that this method is called direct search factorization, essentially brute forcing the solution by trial division. It states that this implementation is extremely inefficient and “simple-minded”. Ouch.

http://mathworld.wolfram.com/DirectSearchFactorization.html

Answer: 6857

Correction: I initially thought that my algorithm did not give prime factors, but on writing this post, I realized that it actually did. My solution was more elegant and complex that I had even percieved (ego boost)! Well, not really elegant (according to mr. wolfram), but it did do something I had not predicted. Guess my prescience could use a little work.

Thinking about it, the code tests for a factor by starting small and gradually getting larger. So by the time it reaches the large non-prime factors, it will have already found and eliminated the smaller primes (2, 3, 5, etc.) which are factors of these larger factors. My terminology is slightly arbitrary but I think you get the point. By the time it reaches a possible non-prime factor, it will have eliminated the smaller prime factors of that factor, so that particular number will not be counted as a factor. I think I might be just confusing you at this point, but Mathworld eloquently puts it as: “all possible factors (m) have had their cofactors already tested. …Therefore, m must be prime”.