The previous post looked at how many digits are in the reduced fraction for the nth harmonic number. I was curious about how long the cycle of digits in a harmonic number might be.

I wrote about the period length for the digits of fractions almost a decade ago. This post includes code so I can apply it to harmonic denominators.

from sympy import lcm, factorint, n_order

def period(n):
    factors = factorint(n)
    exp2 = factors.get(2, 0)
    exp5 = factors.get(5, 0)
    r = max(exp2, exp5)

    d = n // (2**exp2 * 5**exp5)
    s = 1 if d == 1 else n_order(10, d)
    return (r, s)

This function returns two numbers: r is the number of non-repeating digits at the beginning and s is the length of the repeating part.

The following code

from functools import reduce

def lcm_range(n):
    return reduce(lcm, range(1, n + 1))

print( period( lcm_range(50) ) )

prints (5, 1275120) meaning that 1/lcm(1, 2, 3, …, 49, 50) has five non-repeating digits following by 1,275,120 digits that repeat ad infinitum. And so the decimals in the expansion of H₅₀ have a cycle length of 1,275,120.

The previous post looked at writing the harmonic numbers as reduced fractions and estimating the number of digits in the numerator and denominator based on asymptotics. This is a follow up post with plots.

We’ll choose our base b to be 2. And we’ll look at the total number of bits in both the numerator and denominator, which we will use as the height of the fractions.

First, let’s look at the actual and estimated heights, using the estimates from the previous post.

Next let’s look at the difference between the actual and estimated heights.

In the previous post I looked at n = 50, which was kind of a lucky choice, the error being smaller than usual. I had also looked at, but didn’t publish, n = 100, which would be an unlucky choice.

Finally, let’s look at the relative error in the estimates, and plot over a larger range of n.

The error goes to zero, as predicted by the asymptotic estimates. And it goes noisily, which you’d expect since the heights are related to the distribution of primes.

The nth harmonic number is the sum of the reciprocals of the first n positive integers.

H_n = 1 + 1/2 + 1/3 + 1/4 + … + 1/n

The product of all the denominators is n!, so you could write H_n as a fraction

H_n = p/q

where p = n! H_n is an integer and q = n!.

While p/q is a way to write H_n as a fraction, it’s not the most efficient because p and n! will have common factors.

If we write H_n as a reduced fraction, the denominator will be the least common multiple of the integers 1 through n. That number is asymptotically exp(n). That estimate follows from the prime number theorem.

So for large n the denominator will be roughly exp(n), and in base b it would have around

n/log(b)

digits.

The numerator will be exp(n) H_n, and since H_n is asymptotically log(n) + γ, the numerator for large n will be roughly

exp(n) (log(n) + γ)

and will have around

(n + log log(n) ) / log(b)

digits.

Let’s see how well our asymptotic estimates work for n = 50. The 50th harmonic number is

H₅₀ = 13943237577224054960759 / 3099044504245996706400.

This fraction has 23 digits in the numerator and 22 in the denominator. We would have predicted around

(50 + log(log(50)))/log(10) = 22.3

digits in the numerator and

50/log(10) = 21.7

digits in the denominator.

Let’s try a larger example, looking at the 1000th harmonic number in binary. We’ll use the following Python code.

from fractions import Fraction

def bits(n):
    H = sum(Fraction(i, i+1) for i in range(1, n+1))
    p, q = H.numerator, H.denominator
    # subtract 2 because bin returns a string starting with 0b.
    return len(bin(p)) - 2, len(bin(q)) - 2

print(bits(1000))

This returns 1448 and 1438. We would have estimated

(1000 + log(log(1000)))/log(2) = 1445.4

bits in the numerator and

1000/log(2) = 1442.7

bits in the denominator.

Update: See the next post for plots as a function of n.

Related posts

• Computing harmonic numbers
• Sums of consecutive reciprocals
• AM over GM

In this post we find all Pythagorean triples that contain consecutive numbers, all Pythagorean triples (a, b, c) such that a + 1 = b or b + 1 = c.

a + 1 = b

George Osborne wrote a paper [1] addressing the question of when the squares of two consecutive numbers is also a square. Geometrically this is asking for primitive Pythagorean triples for which the legs are consecutive integers.

He proved that the sequence shorter legs satisfies the recurrence relation

with initial conditions u₀ = 0 and u₁ = 1. This is OEIS sequence A001652.

The method for solving recurrences like the one above is analogous to the method for solving linear differential equations. See a solution here. This gives us the following formula for the terms:

b + 1 = c

It’s also possible for the longer side and hypotenuse of a Pythagorean triangle to be consecutive numbers, as in the (5, 12, 13) triangle.

All primitive Pythagorean triples are given by Euclid’s formula

with integers m > n > 0. If b and c are consecutive numbers, then

and so m = n + 1. Therefore all possible values of b are given by 2n(n + 1) for n > 1.

[1] Geo. A. Osborne. A Problem in Number Theory. The American Mathematical Monthly, Vol. 21, No. 5 (May, 1914), pp. 148-150

The nth pentagonal number P_n is the number of dots in diagrams like those below with n concentric pentagons.

We have the formula

P_n = (3n² − n)/2

where n is a positive integer. If n is an integer but not positive, the equation above defines a generalized pentagonal number.

If you’re given an n, you can easily compute P_n. But suppose you’re given a large number x. How would you determine if it is a pentagonal number? And if it is a pentagonal number, how would you find n such that x = P_n?

Rejecting non-pentagonal numbers

If

x = P_n = (3n² − n)/2

then we can solve a quadratic equation for n:

n = (1 ± √(24x + 1))/6.

If 24x + 1 is not a perfect square, n is not an integer and x is not a pentagonal number, ordinary or generalized.

Small numbers

For example,

√(24 × 20260615 + 1)) = 22051.185…

and so 20260615 is not a pentagonal number nor a generalized pentagonal number.

Big numbers

Now suppose

x = 170141183460469231731687303715884105727.

Is this a pentagonal number? You can’t just compute √(24x + 1) in floating point arithmetic because the result is a 20-digit number, and floating point number have 15 digits of precision, so you can’t tell whether the result is an integer.

However, you can compute

⌊√(24x + 1)⌋

with only integer arithmetic using the sqrt_floor function from this post.

def sqrt_floor(n):
    a = n
    b = (n + 1) // 2
    while b < a:
        a = b
        b = (a*a + n) // (2*a)
    return a

The following prints a positive number,

x = 2**127 - 1
y = 24*x + 1
r = sqrt_floor(y)
print(y - r**2)

which tells us y is not a perfect square.

Finding the index

Now suppose y is a perfect square. Then the roots of

(1 ± √(24x + 1))/6

are rational, but are they integers? In fact one, and only one, of the roots will be an integer. If

24x + 1 = r²

then r is congruent to ±1 mod 6 because the left side is congruent to 1 mod 6. If r = 1 mod 6 then the smaller root is an integer, and if r = 5 mod 6 then the larger root is an integer.

So if 24x + 1 = r², then x is a pentagonal number if r = 5 mod 6 and a generalized pentagonal number otherwise.

The function pentagonal_index takes a number x and return n if x = P_n and None if no such n exists.

def pentagonal_index(x):
    y = 24*x + 1
    r = sqrt_floor(y)
    if r*r != y:
        return None
    if r % 6 == 5:
        return (1 + r) // 6
    else:
        return (1 - r) // 6

We can test this with the following code.

P = lambda n: (3*n**2 - n) // 2
for n in [2, 3, -4, -5, 10**200]:
    assert(pentagonal_index(P(n)) == n)

Integer division in Python

Note that P(10**200) is too big to fit into a float, but the code works fine. This is because we use integer division (//) everywhere. If we had said

P = lambda n: (3*n**2 - n) / 2

the test above would pass for the small values of n but output

OverflowError: integer division result too large for a float

when it came to n = 10²⁰⁰.

Related posts

• Partitions and pentagons
• Certified Fibonacci numbers
• Truncated triangular numbers

There’s a lot of mathematics just off the beaten path. You can spend a career in math and yet not know all there is to know about even the most basic areas of math. For example, this post will demonstrate something you may not have seen about decimal forms of fractions.

Let p > 5 be a prime number and 0 < k < p. Then the digits in k/p might be the same for all k, varying only by cyclic permutations. This is the case, for example, when p = 7 or p = 17. More on these kinds of fractions here.

The digits in k/p repeat for every k, but different values of k might have sequences of digits that vary by more than cyclic permutations. For example, let’s look at the values of k/13.

>>> for i in range(1, 13):
...   print(i/13)
...
 1 0.0769230769230769
 2 0.1538461538461538
 3 0.2307692307692307
 4 0.3076923076923077
 5 0.3846153846153846
 6 0.4615384615384615
 7 0.5384615384615384
 8 0.6153846153846154
 9 0.6923076923076923
10 0.7692307692307693
11 0.8461538461538461
12 0.9230769230769231

One cycle goes through the digits 076923. You’ll see this when k = 1, 3, 4, 9, 10, or 11. The other cycle goes through 153846 for the rest of the values of k. The cycles 076923 and 153846 are called the distinct repeating sets of 13 in [1].

If we look at fractions with denominator 41, thee are six distinct repeating sets.

02439
04878
07317
09756
12195
14634
26829
36585

You could find these by modifying the Python code above. However, in general you’ll need more than default precision to see the full periods. You might want to shift over to bc, for example.

When you look at all the distinct repeating sets of a prime number, all digits appear almost the same number of times. Some digits may appear one more time than others, but that’s as uneven as you can get. A corollary in [1] states that if p = 10q + r, with 0 < r < 10, then 11 − r digits appear q times, and r − 1 digits appear q + 1 times.

Looking back at the example with p = 13, we have q = 1 and r = 3. The corollary says we should expect 8 digits to appear once and 2 digits to appear twice. And that’s what we see: in the sets 076923 and 153846 we have 3 and 6 repeated twice and the remaining 8 digits appear once.

In the example with p = 41, we have q = 4 and r = 1. So we expect all 10 digits to appear 4 times, which is the case.

Related posts

• Cyclic fractions
• Periods of fractions
• Calendars and continued fractions

[1] James K. Schiller. A Theorem in the Decimal Representation of Rationals. The American Mathematical Monthly
Vol. 66, No. 9 (Nov., 1959), pp. 797-798

I recently found out about Andrica’s conjecture: the square roots of consecutive primes are less than 1 apart.

In symbols, Andrica’s conjecture says that if p_n and p_n+1 are consecutive prime numbers, then

√p_n+1 − √p_n < 1.

This has been empirically verified for primes up to 2 × 10¹⁹.

If the conjecture is true, it puts an upper bound on how long you’d have to search to find the next prime:

p_n+1 < 1 + 2√p_n  + p_n,

which would be an improvement on the Bertrand-Chebyshev theorem that says

p_n+1 < 2p_n.

Pentagonal numbers are truncated triangular numbers. You can take the diagram that illustrates the nth pentagonal number and warp it into the base of the image that illustrates the (2n − 1)st triangular number. If you added a diagram for the (n − 1)st triangular number to the bottom of the image on the right, you’d have a diagram for the (2n − 1)st triangular number.

In short,

P_n = T_{2n − 1} − T_{n − 1}.

This is trivial to prove algebraically, though the visual proof above is more interesting.

The proof follows immediately from the definition of pentagonal numbers

P_n = (3n² − n)/2

and triangular numbers

T_n = (n² + n)/2.

Related posts

• Partitions and pentagons
• Tetrahedral numbers
• A pentagon, hexagon, and decagon walk into a bar …
  

If a number has a finite but nonzero fractional part, so do all its powers. I recently ran across a proof in [1] that is shorter than I expected.

Theorem: Suppose r is a real number that is not an integer, and the decimal part of r terminates. Then r^k is not an integer for any positive integer k.

Proof: The number r can be written as a reduced fraction a / 10^m for some positive m. If s = r^k were an integer, then

10^mk s = a^k.

Now the left side of this equation is divisible by 10 but the right side is not, and so s = r^k must not be an integer. QED.

The only thing special about base 10 is that we most easily think in terms of base 10, but you could replace 10 with any other base.

What about repeating decimals, like 1/7 = 0.142857142857…? They’re only repeating decimals in our chosen base. Pick the right base, i.e. 7 in this case, and they terminate. So the theorem above extends to repeating decimals.

[1] Eli Leher. √2 is Not 1.41421356237 or Anything of the Sort. The American Mathematical Monthly, Vol. 125, No. 4 (APRIL 2018), page 346.