Programming Praxis

In the study of probability, the Birthday Paradox states that in a group of 23 people, there is a 50% chance that two will have the same birthday; in a group of 57 people, the odds rise to 99%.

Your task is to simulate the birthday paradox over many trials and verify the odds. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

It’s been a while since we played word games. We have two today:

1) Find all the words in a dictionary that contain exactly five vowels (a, e, i, o and u) in ascending order.

2) Find all the words in a dictionary of length at least six that contain letters in strictly ascending alphabetical order.

These games are easy to play using regular expressions, so you should solve them without regular expressions, using only simple string manipulations. If your system doesn’t provide a dictionary, you can find one at http://icon.shef.ac.uk/Moby/mwords.html.

Your task is to play the two games. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

In the previous exercise we wrote some auxiliary functions needed for the implementation of the AKS primality prover. Today we will implement the AKS algorithm:

AKS Primality Prover: Given an integer n > 1, determine if it is prime or composite.
1. If n = a^b for integers a > 0 and b > 1, output COMPOSITE.
2. Find the smallest r such that ord_r(n) > (log₂ n)².
3. If 1 < gcd(a, n) < n for some a ≤ r, output COMPOSITE.
4. If n ≤ r, output PRIME.
5. For each a from 1 to floor √φ(r) · log₂ n), if (x + a)ⁿ &neq; xⁿ + a (mod x^r − 1, n), output COMPOSITE.
6. Output PRIME.

Here ord_r(n) is the multiplicative order of n modulo r, both logarithms are to the base 2, φ(r) is Euler’s totient function, and the polynomial modular exponentiation is done as in the previous exercise.

Your task is to write a program to prove primality using the AKS algorithm. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

In the previous exercise, we looked at Pocklington’s algorithm for proving the primality of a number n; the algorithm depends on at least a partial factorization of n−1. In today’s exercise we look at an algorithm of Manindra Agrawal, Neeraj Kayal, and Nitin Saxena of the Indian Institute of Technology in Kanpur, which does not rely on n to be of any special form, is purely deterministic, of provable polynomial time, not conditional on the Riemann Hypothesis, and unfortunately very slow.

The entire mathematical community was astonished in 2002 when Agrawal, Kayal and Saxena published their algorithm. In the first place, Kayal and Saxena were undergraduate students at the time; Agrawal was their professor. And in the second place, their proof is simple, not relying on any complicated math; the general reaction among mathematicians was “How could we have missed that? What else have we missed?” The three won numerous awards, and their algorithm has been thoroughly studied and analyzed.

We will look at the AKS algorithm in the next exercise. Today, we write two functions that will be needed to implement the AKS algorithm.

The first function computes the multiplicative order of n modulo r, written as ord_r(n), which is the smallest number k such that n^k ≡ 1 (mod r). The algorithm is straight forward: starting from k = 2, increment k until you find the answer. The multiplicative order exists only when n and r are co-prime.

The second function computes the modular exponentiation of a polynomial; it’s not hard to write, but requires some explanation. To compute the modular exponentiation of a polynomial is simple; it uses the same square-and-multiply algorithm that is used for the modular exponentiation of an integer, which is part of our Standard Prelude and is described in Algorithm 3.C of the essay Programming with Prime Numbers. All we need is the modular multiplication of two polynomials, which we will explain by example.

Consider the problem of squaring polynomial x³ + 4x² + 12x + 3 modulo (x⁵ − 1, 17). Polynomial multiplication is exactly the same as grade-school multiplication, except there are no carries, so the process looks like this:

                 1    4   12    3
                 1    4   12    3
               ---  ---  ---  ---
                 3   12   36    9
           12   48  144   36
       4   16   48   12
  1    4   12    3
---  ---  ---  ---  ---  ---  ---
  1    8   40  102  168   72    9

Thus, x³ + 4x² + 12x + 3 squared is x⁶ + 8x⁵ + 40x⁴ + 102x³ + 168x² + 72x + 9. To reduce the result modulo x⁵ − 1 we divide using the grade-school long-division algorithm and take the remainder, which gives x + 8 with a remainder of 40x⁴ + 102x³ + 168x² + 73x + 17:

                                          1    8
             + ---  ---  ---  ---  ---  ---  ---
1 0 0 0 0 -1 |   1    8   40  102  168   72    9
                 1    0    0    0    0   -1
               ---  ---  ---  ---  ---  ---
                      8   40  102  168   73    9
                      8    0    0    0    0   -8
                    ---  ---  ---  ---  ---  ---
                          40  102  168   73   17

We can confirm the calculation by multiplying and adding the remainder:

       1    0    0    0    0   -1
                           1    8
     ---  ---  ---  ---  ---  ---
       8    0    0    0    0   -8
  1    0    0    0    0   -1
---  ---  ---  ---  ---  ---  ---
  1    8    0    0    0   -1   -8
           40  102  168   73   17
---  ---  ---  ---  ---  ---  ---
  1    8   40  102  168   72    9

Then we simply reduce each of the coefficients of the remainder modulo 17, giving the result 6x⁴ + 0x³ + 15x² + 5x + 0. The whole computation can be organized as shown below. Note how the division and reduction modulo x⁵ − 1 is accomplished, eliminating the high-order coefficients and adding them back to the low-order coefficients; we are relying here on the simple form of the polynomial modulus, and would have to do more work if it was more complicated:

                 1    4   12    3   multiplicand
                 1    4   12    3   multiplier
               ---  ---  ---  ---
                 3   12   36    9    3 * 1 4 12 3 * x^0
           12   48  144   36        12 * 1 4 12 3 * x^1
       4   16   48   12              4 * 1 4 12 3 * x^2
  1    4   12    3                   1 * 1 4 12 3 * x^3
---  ---  ---  ---  ---  ---  ---
  1    8   40  102  168   72    9   1 4 12 3 * 1 4 12 3
 -1   -8                   1    8   reduce modulo x^5 - 1
---  ---  ---  ---  ---  ---  ---
           40  102  168   73   17   reduce modulo 17
            6    0   15    5    0   final result

Your task is to write functions that compute the multiplicative order of n modulo r and the modular exponentiation of polynomials. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution in the comments below.

Pages: 1 2

In previous exercises we have seen two methods for determining if a number is probably prime, the Miller-Rabin and Baillie-Wagstaff methods, and two methods for proving a number prime, trial division and the Lucas N-1 prover, which requires the complete factorization of n−1, where n is the number to be proved prime. In today’s exercise we examine a method of Pocklington in which a number can be proved prime with only a partial factorization of n−1.

We begin by expressing the number to be proved prime as n = f · r + 1, where the complete factorization of f is known and exceeds the square root of n. Pocklington proved that if, for some a, aⁿ⁻¹ ≡ 1 (mod n) and gcd(a^(n−1)/q − 1, n) = 1 for each prime q|f, then n is prime if f > √n.

As an example, let’s consider the primality of n = 2⁸⁹−1. A probable-prime test suggests that n is prime. Trial division by the primes less than a thousand tells us that n = 2 · 3 · 5 · 17 · 23 · 89 · 353 · 397 · 683 · 6194349127121 + 1 = 99924948842910 · 6194349127121 + 1, where the primality of 6194349127121 is indeterminate (it’s composite, but we don’t care). Since f = 99924948842910 is greater than the square root of n, we will be able to prove the primality of n by Pocklington’s theorem.

Consider first the situation a = 2; we can choose any a randomly from the range 1 < a < n − 1, but it is easiest to choose a from the sequence 2, 3, 4, …. Now 2ⁿ⁻¹ ≡ 1 (mod n), so the first criterion holds, but 2^(n−1)/2 ≡ 1 (mod n), so the gcd is n, and the second criterion fails. Next we consider a = 3. Now 3ⁿ⁻¹ ≡ 1 (mod n), so the first criterion holds, and the second criterion holds for each q ∈ {2, 3, 5, 17, 23, 89, 353, 397, 683}, so n = 2⁸⁹ − 1 is prime.

Brillhart, Lehmer and Selfridge later extended Pocklington’s theorem to work with f between the cube root and square root of n. The idea is to determine the constants c₁ and c₂ such that n = c₂ · f² + c₁ · f¹ + 1, which is easily done by extracting the “digits” of n to the base f. Then n is composite if c₁² − 4 c₂ is a perfect square and prime if it is not; Pocklington’s criteria must also hold.

With this extension to Pocklington’s theorem we can prove the primality of 2⁸⁹ − 1 with only the factors less than 500. With a = 3, Pocklington’s first criterion holds, and Pocklington’s second criterion holds for each q ∈ {2, 3, 5, 17, 23, 89, 353, 397}, and f = 146302999770 is greater than the cube root of n, so we calculate n = 28917 f² + 96609474553 f + 1, then 96609474553² − 4 · 28917 = 9333390573406754434141 = 12611 · 95783 · 7726832165257 is not a perfect square, so n is prime.

Pocklington’s method fails if you can’t find sufficient factors of n − 1. We prefer trial division rather than some more powerful method because we have to prove that each of the factors is prime, and trial division does that for us automatically; if you can’t factor n−1 by trial division, but you can factor it by Pollard’s rho method or something even more powerful, you could still use Pocklington’s theorem, proving each of the factors prime by using Pocklington’s theorem on each of them recursively. For instance, if you are determined to prove the primality of n = 2⁵²¹ − 1, trial division by the primes less than a thousand finds the factors 2, 3, 5, 5, 11, 17, 31, 41, 53, 131, 157 and 521, and Pollard’s rho method quickly finds the factors 1613, 2731, 8191, 42641, 51481, 61681, 409891, and 858001, all of which can be proved prime by Pocklington’s method. Taken together, the product of those factors is greater than the cube root of n, so you can prove the primality of 2⁵²¹ − 1 even though you can’t factor n − 1 by trial division. Of course, it is still possible that you can’t find enough factors (as an example, try to prove the primality of 2⁶⁰⁷ − 1), in which case you need a more powerful method, but that’s a topic for a different exercise.

Your task is to write a program that proves the primality of an input number n, or shows that it is composite, using the method of Pocklington with extensions from Brillhart, Lehmer and Selfridge. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

Today is the autumn equinox, when the day and night are of equal lengths at the equator.

Your task is to write a program that calculates the length of the day and the night to show that they are of equal length. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

The ABC Conjecture has recently been in the news on math blogs because of the claim that it has been proved by Shinichi Mochizuki. Though the proof is being taken seriously, due to Mochizuki’s reputation, it is five hundred pages long, and confirmation will take several months. The conjecture states that given two positive integers a and b and their sum c with no common factors, the product of the distinct prime factors of abc is rarely much smaller than c.

The radical of a number n is the product of the distinct prime factors of n; for instance, rad(18) = 6 because 18 = 2 × 3 × 3 and, eliminating the duplicate occurrence of 3, 2 × 3 = 6. The quality of an (a,b,c) triple is given by q(a,b,c) = log(c) / log(rad(abc)). The precise statement of the ABC conjecture is that for every ε > 0 there are only finitely many triples (a,b,c) with a and b coprime positive integers and a + b = c such that rad(abc)^1+ε < c, or equivalently, such that q(a,b,c) > 1.

Your task is to write the functions rad and q and find all the triples with c < 1000 for which q(a,b,c) > 1. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

You will recall that fibonacci numbers are formed by a sequence starting with 0 and 1 where each succeeding number is the sum of the two preceeding numbers; that is, F[n] = F[n-1] + F[n-2] with F[0] = 0 and F[1] = 1. We studied fibonacci numbers in a previous exercise.

Tribonacci numbers are like fibonacci numbers except that the starting sequence is 0, 0 and 1 and each succeeding number is the sum of the three preceeding numbers; that is, T[n] = T[n-1] + T[n-2] + T[n-3] with T[-1] = 0, T[0] = 0 and T[1] = 1. The powering matrix for tribonacci numbers, used similarly to the powering matrix for fibonacci numbers, is:

The first ten terms of the tribonacci sequence, ignoring the two leading zeros, are 1, 1, 2, 4, 7, 13, 24, 44, 81 and 149 (A000073).

Your task is to write two functions that calculate the first n terms of the tribonacci sequence by iteration and the nth term by matrix powering; you should also calculate the tribonacci constant, which is the limit of the ratio between successive tribonacci numbers as n tends to infinity. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

The first ten prime numbers are 2, 3, 5, 7, 11, 13, 17, 19, 23, and 29. Their sum is 129.

Your task is to write a program that calculates the sum of the first billion primes. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2