The Programmatic Side of Mathematica V: Segue into Primes

You might recall in the last post the expansion of a general binomial (a + b)^z:

binexpand[z_] := Expand[(a + b)^z]

You can also declare a loop in the way Mathematica does so that it gives a series of binomial expansions. Below is the output of the second statement from above, the first 6 expansions of a general binomial:

And if you want to mess around with the Map command, you can put square roots on everything:

You don’t always need to use Array[] to generate a list. You can also use Range[hi] which generates a list from 1 to “hi”:


Generates the first 100 primes:

{2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, \
73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, \
157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, \
239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, \
331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, \
421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, \
509, 521, 523, 541}

In earlier posts I had been curious about the distribution of prime numbers along the number line, and took my research to very large numbers, but I fear not large enough since I don’t have a supercomputer handy. My considered opinon first-hand, is that the primes are initially densely distributed, then gets gradually sparser in a way that is hard to predict. Mathematica claims it can easily access primes as high as on the order of 10^8 in magnitude, and likely does the higher primes more slowly due to the fact that it uses primes up to that magnitude to search for still higher primes. If a natural number n has no prime numbers as factors, then n itself is prime.

The Mersenne primes, those Mersenne primes which are still being sought after, numbers large enough to fill every page in a 200-page novel, are notoriously hard to prove indivisible, since I would fear that the primes up to 10^8 are not good enough. You need higher primes still, to many more orders of magnitude. Knowing the previous Mersenne prime is unreliable, since who knows how many primes lay between the current prime being checked and the previous one? There are likely to be hundreds, or even thousands of undiscovered primes which lay between 2^{p_{n-2}} - 1 and 2^{p_{n-1}} - 1, assuming 2^{p_{n-2}} and 2^{p_{n-1}} are prime (p, the exponent, is prime as well). The method of last resort, as I understand it, is to try every number up to \lfloor \sqrt{2^{p}} \rfloor past 10^8. While trying primes below 10^8 weeds out several numbers, I think would still take a very long time for certain numbers that remain unfactored past that, possibly days to months even on a fast computer processor — to weed out one number or to successfully prove it as prime.

All other prime numbers which lay between two consecutive Mersenne primes (especially the big numbers) are nearly impossible to find, and I am not aware of anyone who bothers trying to find them (that doesn’t mean no one is trying, though). To find the next prime up from a Mersenne prime means looking to the next consecutive odd number, and the next, and the next, until you find one. A lot of time — months of time — is possibly wasted on a lot of composite numbers. Mersenne numbers narrow the search down, but only to a certain kind of prime number.

Mathematica has, for what it’s worth, a function called Prime[n], which returns the nth prime. Prime[1] returns 2; Prime[2] returns 3, and so on. The 100th prime is 541:


And if we subtracted

Prime[100] - Prime[1]

we can say that 539 integers lay between the first and the 100th prime. We can now do this:

Prime[1000100] - Prime[1000000]

which tells us there is a difference of 1606 between the millionth prime and the 1,000,100th prime.

Prime[1000000100] - Prime[1000000000]

The air gets thinner, but only slightly, when counting between the billionth prime and the 1,000,000,100th prime. By the trillionth prime, we get a difference of 3546 by the 1,000,000,000,100th prime. Or you can say 100 primes per 3546 consecutive integers. The trillionth prime is 29,996,224,275,833, or just under 30 trillion. Mathematica’s Prime[] function chokes when getting to numbers on the order of a quadrillion, even though the documentation claims arguments to Prime[] are allowed to go into the quintillions: 1,152,921,504,606,846,976 or 260.

I could say that this still provides an incredible list of primes for the possible writing of a prime factorization function, however, Mathematica has beat me to it, wtih the function FactorInteger[]. Even the Mathematica website at Wolfram Math World has defined a function that modifies the output of this function further:


So, a call to FactorForm[] can be:


I don’t know what “Times@@” does, but I don’t seem to get anything other than the default terminal font. At any rate, the ordered pairs returned by FactorInteger[] are transformed into base integers and exponents.