A number-theoretic tangent for your Saturday afternoon

You may know that 1001 has a simple prime factorization, namely 1001 = (7)(11)(13).

And of course 1000 has a simple prime factorization, namely 1000 = (2³)(5³).

The natural question to ask at this point, to me, is "how often does this happen?" Of course one has to define "this". One question is "how often are there two consecutive numbers that have all their prime factors less than or equal to 13?" There seem to be finitely many, and after I compiled a list I googled a few numbers from it and discovered that David Rusin had also done so. In case you're wondering why I conjecture there are finitely many: the number of integers less than n with all prime factors at most 13 is proportional to (log n)⁶. (The exponent 6 comes from the fact that there are six primes which are 13 or less.) So the "probability" (in the usual heuristic number-theoretic sense) that a given number n is "13-smooth" is the derivative of this, and thus of the order of (log n)⁵/n. The "probability" that two consecutive numbers are "13-smooth" is on the order of the square of this "probability", i. e. (log n)¹⁰n^-2, and the integral of that from, say, 2 to infinity converges. Nothing is special about 13 here; there should be finitely many pairs of consecutive "p-smooth numbers" where p is any prime.

(Rusin's list, by the way, was something he compiled to illustrate how one might calculate logarithms; since 1000 ~ 1001, we can take common logs and get 3 ~ log(7) + log(11) + log(13). Given enough relations like this one can solve for the logarithms themselves.)

Alternatively, we can define the "roughness" of n as (log p)/(log n) where p is the largest prime factor of n. I call it "roughness", not "smoothness", because if it's smaller than the corresponding number is smoother, i. e. has comparatively small prime factors. This naturally compensates for the size of the number; if I recall correctly the proportion of numbers less than n with roughness less than r approaches some limit (a function of r) as n goes to infinity. (I only dabble in analytic number theory so I can't provide a source off the top of my head.)

It seems reasonable then that there should be a similar limit for pairs of consecutive numbers. Then empirically it seems that the probability that n and n+1 both have roughness at most (log 13)/(log 1001) is about 1 in 250. The pairs of numbers that are "at least as round" as (1000, 1001) in this sense are

(80, 81), (224, 225), (1000, 1001), (1715, 1716), (2079, 2080), (2400, 2401), ...

and it seems like about one in 230 pairs of consecutive numbers have this property (4291 in the first million). As is often the case, the interesting thing is that the limiting probability appears to exist and be neither zero nor one.


Postscript: I know some readers will be offended by my use of "probability" in this number-theoretic sense, because the prime factorization of a number is hardly a random object! But I'm a probabilist; this is how I think.

Odd periods in continued fractions

Here's a question. Why is the period of the quotients in the continued fraction of N^1/2 "usually" even? For example, if N runs over the ninety non-squares less than 100, then only 20 times does the continued fraction expansion of N^1/2 have an odd period. Of the 992 non-squares less than 1024, 157 have an odd period. Of the 9900 squares less than 10⁴, 1322 have an odd period. This is a sign that something is going on under the hood -- naively you'd expect half the periods to be odd.

Arnold has observed this, but only empirically; I first observed it from this problem from Project Euler.

The period of the continued fraction of N^1/2 is odd if and only if x² - Ny² = -1 has solutions in integers. All such integers, it turns out, have no prime factors congruent to 3 mod 4, which is pretty rare for large numbers. (The number of positive integers less than N with no prime factors congruent to 3 mod 4 is about N(log N)^-1/2.) For integers having no prime factors congruent to 3 mod 4, though, a paper of Etienne Fouvry and Jurgen Kluners shows that asymptotically at least 52% of such numbers have odd period, and at most two-thirds do.

Heuristic derivation of the Prime Number Theorem

One way of stating the Prime Number Theorem is that the "probability" that a large number near x is prime is 1/log(x). (Here, as always, all logs are natural.)

A heuristic derivation of the prime number theorem, Frank Morgan, via Andrew Gelman, with some embellishment by me: let P(x) be the "probability" that a large integer x is prime. Then how do P(x+1) and P(x) compare? Say x is prime, which it is with "probability" P(x); it "divides" x+1 (and all larger numbers) with probability 1/x. (Of course, this doesn't actually make sense; the idea is that we're modeling the numbers divisible by x as a random set with density 1/x, which should have the same large-scale properties.) Other than this, x and x+1 are equally "likely" to be prime. So the function P satisfies

(*) P(x+1) = P(x) [P(x) (1 - 1/x)] + (1 - P(x)) P(x)

since x+1 is "prime" with probability P(x) (1-1/x) if x is "prime", and P(x) otherwise. Divide through by P(x) to get

P(x+1)/P(x) = P(x) (1-1/x) + (1-P(x))
P(x+1)/P(x) = 1 - P(x)/x.

Now, P(x+1) can be approximated by P(x) + P'(x), so we have

1 + P'(x)/P(x) = 1 - P(x)/x
P'(x)/P(x) = -P(x)/x

which has the general solution P(x) = 1/(C + log x).

This is a bit of a lie. For one thing, the difference equation (*) actually seems to have solutions that differ from those of the differential equation by a constant factor, which seems to depend on the initial conditions. (This amounts to changing the base.) For another thing, the assumption that x+1 might be divisible by x is, um, stupid if we're actually talking about prime numbers. (It's probably possible to rephrase this heuristic derivation as a rigorous result about random sets, though.) Still, it gets the prime number theorem to within a constant factor, which isn't bad for such a simple argument.

A couple amusements from the Princeton Companion

From The Princeton Companion to Mathematics, specifically the article "Algebraic Geometry", by János Kollár:

Finally, if we marry a scheme to an orbifold, the outcome is a stack. The study of stacks is strongly recommended to people who would have been flagellants in earlier times.

I feel this way about most of algebraic geometry, but that's only because Penn has a high enough concentration of algebraic geometers that I get tired of hearing people walking around and talking about it all the time.

Also, I find myself screaming at the book quite often. But I do it in a good way; it's often some crucial insight that makes me think "why didn't anybody tell me that before?" (Of course, it's possible they did and I wasn't listening.) And sometimes it's "oh, damn, I thought I came up with that myself". For example, I just read the article on computational number theory, by Carl Pomerance, in which he explains a heuristic reason why Fermat's last theorem is true. I won't give it in full here, but it's basically the following. First, Euler showed it for n = 3. Second, consider all the positive integers which are nth powers for some n ≥ 4. The "probability" that a number m is in this set is about m^-3/4. So replace the set of fourth-or-higher powers with a random set S, which contains m with probability m^-3/4. Then the probability that a givennumber n can be written as a sum of two such elements of this set S is proportional to n^-1/2; independently, it has probability n^-3/4 of being in S. So the probability that n can be written as a sum of two elements of S and is also in S itself is proportional to n^-5/4, and so we expect finitely many examples. This isn't quite true, because the set of fourth-or-higher powers has some nontrivial structure. But it also took a couple hundred words, instead of a couple hundred pages like the real proof.

But I figured out something like this when I was in college, and I was so proud of myself! So it saddens me to learn I'm not the only one who thought of it. (It also makes me happy, though, because the idea is due to Erdos and Ulam, and there are worse people to be imitating.)

Playing around with sums of squares

An interesting random fact, learned from Arcadian Functor:

1² + 2² + 3² + ... + 24² = 70²

This doesn't seem quite so weird if you recall that the sum of the first n squares is n(n+1)(2n+1)/6; then the sum of the first 24 squares is (24)(25)(49)/6.

But are there any other n for which the sum of the first n squares is also a square? Not under a million. I tested. It seems kind of tricky to prove this directly, because knowing the factorization of any of n, n+1, or 2n+1 doesn't tell you much about the factorizations of the others.

Let f(n) = n(n+1)(2n+1)/6; for what proportion of values of n is f(n) squarefree? For random integers it's π²/6. Let g(n) be the number of integers [1,n] such that f(n) is squarefree; then g(1000) = 504, g(10000) = 5029, g(100000) = 50187.

Hmm, it looks an awful lot like g(n) ~ n/2, doesn't it? I have no idea if this is true. It seems quite obvious that g(n) ~ cn for some constant c, but I'm not sure if c is exactly 1/2. (And I'd be surprised if it did come out to be 1/2; these sorts of probabilities are rarely rational.)

The density of abundant numbers

Let σ₁(n) be the sum of the divisors of n. So, for example, the divisors of 18 are 1, 2, 3, 6, 9, and 18; so σ₁(18) = 1 + 2 + 3 + 6 + 9 + 18 = 39.

If σ₁(n) > 2n, we call n abundant, so 18 is abundant. If σ₁(n) = 2n, we call n perfect; if σ₁(n) < 2n, we call it "deficient".

(If you're wondering about the subscript 1: σ_k(n), more generally, is the sum of kth powers of divisors of n. It turns out that σ₁(n)/n = σ_-1(n), so everything I'm about to say could be rephrased in that terminology.)

A basic result of analytic number theory (for example, see Apostol, Introduction to Analytic Number Theory, Theorem 3.4) is

. If we write

then the previous result says that Σ₁(n) is about n²(π²/12). But Σ₁ is in a sense an "integral" of σ₁; the "average value" of σ₁, for numbers near n, should just be Σ₁'(n), or about n(π²/6). That is, the sum of the divisors of a number near n is, on average, about 1.644n.

So, the logical question to ask is: what can we say about the distribution of σ(n)/n? It seems possible that almost all numbers could have σ(n)/n near π²/6; this would seem to mesh with the fact that large perfect numbers are very rare. But in the 1930's, Davenport and Erdos independently showed that this is in fact not the case, and σ(n)/n does seem to have a well-defined distribution; in particular, that

exists and is continuous. A(x) is just the proportion of numbers for which σ(n)/n is at least x, so A(2) is the proportion of numbers which are abundant. The best result I know of (although I just started thinking about this) is that of Deleglise, which is that A(2) is between 0.2474 and 0.2480. From numerical work one can get approximate A(x) for other x; for example, A(1.5) seems to be about 57%, and A(2.5) seems to be about 8.8%. I don't know the full distribution (although I haven't thoroughly read the last two papers I cite below, and I haven't read the first one at all) but I suspect that it is not one of the familiar ones. The Deleglise paper goes through a lot of work just to compute A(2); it seems like it would require equally much work to compute A(x) for each individual x, when what one would like is of course the entire distribution in one fell swoop.

Here are references for the papers I mentioned:
Davenport, H. "Über numeri abundantes." Sitzungsber. Preuss. Akad. Wiss., Phys.-Math. Kl., No. 6, 830-837, 1933.
Deléglise, M. "Bounds for the Density of Abundant Integers." Exp. Math. 7, 137-143, 1998.
Erdos, P. "On the Density of the Abundant Numbers." J. London Math. Soc. 9, 278-282, 1934.

Twin primes

Are There Infinitely Many Primes, by D. A. Goldston, via Casting Out Nines. This was based on a talk given to motivated high school students, so you don't need much background to understand it.

Since this is at least nominally a probability blog, I draw your attention to the following: Call p and p+2 a twin prime pair if they're both prime. Then it's conjectured that the number of twin prime pairs with smaller member less than or equal to x, denoted π₂(x), is asymptotically

and there is probabilistic reasoning that leads to this, which is a version of the Twin Prime Conjecture; basically, the "probability" that a number near t is prime is 1/(log t), and so the probability that two numbers near t are both prime should be the square of this. So the number of twin prime pairs under x should be just that integral, except that we have to take some divisiblity information into account, which is what the product out front does. If a number n is odd, then n+2 is twice as likely to be prime as it would ``otherwise" be; if n is not divisible by some odd prime p, then n+2 has ``probability" (p-2)/(p-1) of being not divisible by p (we must have that n isn't two less than a multiple of p), which is 1 - 1/(p-1)² times the "naive" probability that n+2 isn't divisible by p, namely (p-1)/p.

I rather like results of this form, where we can guess properties of the primes by assuming that the probability that a given integer is divisible by a prime p is 1/p and these events are independent for different primes p; I've talked about this before, when I sketched a heuristic argument for the Goldbach conjecture. Sure,.the primes have no business behaving randomly. But what reason is there that they should behave nonrandomly?

Probabilities of relative primality

From The Universe of Discourse -- the probability that two randomly-selected polynomials over Z₂ are relatively prime is 1/2. (More formally, the probability that two randomly-selected polynomials over Z₂ of degree less than n are relatively prime is 1/2 + 1/4ⁿ; the original result holds as n goes to infinity.) There's actually an almost one-to-one correspondence between pairs which aren't relatively prime and pairs which are. You can find this, and some more general results, in "The Probability of Relatively Prime Polynomials", from the June 2007 issue of Mathematics Magazine, (Arthur T. Benjamin and Curtis D Bennett, page 196).

The most general result is as follows: if a₁(x), ..., a_m(x) are randomly chosen polynomials of degree less than n in F[x], where F is the finite field with q elements, then the probability that they are relatively prime is 1 - 1/q^m-1 + (q-1)/q^mn; this is all shown with a counting argument that invokes the Euclidean algorithm. The case from The Universe of Discourse is the case m = 2, where n goes to infinity.

This reminds me of one of my favorite statements in mathematics: "the probability that two randomly chosen integers are relatively prime is 6/π²." This seems almost meaningless -- how can we choose integers "at random"? More formally, let A(n) be the number of relatively prime pairs of integers in {1, 2, ..., n}; then

The convergence is pretty quick; A(100)/100² = 0.6087, and A(1000)/1000² = 0.608383, while 6/π² = 0.6079...

Why is this true? Well, the probability that two integers are relatively prime is the probability that they're not both divisible by 2, times the probability that they're not both divisible by 3, times the probabiltiy that they're not both divisible by 5, and so on. The probability that two integers aren't both divisible by a given prime p is 1 - 1/p². So the probability that two "randomly chosen integers" aren't both divisible by the same prime p for all primes p is

and Euler showed that this was equal to 1/ζ(2), or 6/π².

Similarly, the probability that n randomly chosen integers are relatively prime (not pairwise relatively prime; I count a triple like (6, 10, 15), since there isn't a prime number dividing all three of those) is 1/ζ(n). We recall that
ζ(n) = 1 + (1/2)ⁿ + (1/3)ⁿ + ...
and if n is large we can throw away all but the first two terms; thus the probability that n randomly chosen integers are relatively prime is about 1 - (1/2)ⁿ. This makes sense; if there are n integers that aren't relatively prime, the reason they're not relatively prime is almost always because they're all even.

The logical next question is: what is the average of the greatest common divisor of randomly chosen integers? From some quick numerical work, it looks like

grows like log n, that is, the average of the GCDs of pairs of integers randomly chosen from {1, 2, ..., n} is about some constant times log n. (The constant, to one decimal place, looks like 0.6.) But I don't know the answer. If I get bored at some point in the next couple days I'll see if I can figure it out.