A rational triangle is a triangle whose sides have rational length and whose area is rational. Can any two rational numbers be sizes of a rational triangle? Surprisingly no. You can always find a third side of rational length, but it might not be possible to do so while keeping the area rational.

The following theorem comes from [1].

Let q be a prime number not congruent to 1 or −1 mod 8. And let ω be an odd integer such that no prime factor of q^2ω − 1 is congruent to 1 mod 4. Then a = q^ω + 1 and b = q^ω − 1 are not the sides of any rational triangle.

To make the theorem explicitly clear, here’s a Python implementation.

from sympy import isprime, primefactors

def test(q, omega):
    if not isprime(q):
        return False
    if q % 8 == 1 or q % 8 == 7:
        return False
    if omega % 2 == 0:
        return False
    factors = primefactors(q**(2*omega) - 1)
    for f in factors:
        if f % 4 == 1:
            return False
    return True

We can see that q = 2 and ω = 1 passes the test, so there is no rational triangle with sides a = 2¹ + 1 = 3 and b = 2¹ − 1 = 1.

You could use the code above to search for (a, b) pairs systematically.

from sympy import primerange

for q in primerange(20):
    for omega in range(1, 20, 2):
        if test(q, omega):
            a, b = q**omega + 1, q**omega - 1
            print(a, b)

This outputs the following:

3 1
9 7
33 31
129 127
8193 8191
32769 32767
131073 131071
524289 524287
4 2
6 4
126 124
14 12

Note that one of the outputs is (4, 2). Since multiplying the sides of a rational triangle by a rational number gives another rational triangle, there cannot be a rational triangle in which one side is twice as long as another.

As the author in [1] points out, “There are undoubtedly many pairs not covered by [the theorem cited above]. It would be of interest to characterize all pairs.” The paper was published in 1976. Maybe there has been additional work since then.

Related posts

• Rational right triangles
• The congruent number problem
• Do perimeter and area determine a triangle?
• Approximation with Pythagorean triangles

[1] N. J. Fine. On Rational Triangles. The American Mathematical Monthly. Vol. 83. No. 7, pp. 517–521.

I recently ran across an article by Katherine Stange [1] on Lehmer’s factoring stencils [2]. These stencils were the basis of an effective method for factoring moderately large numbers before the advent of electronic computers.

This post will describe the stencils and simulate their use with a Python script.

Misconceptions

When I started reading [1] I had three misconceptions that slowed down my understanding of the article. I’ll address these first in case any of you come to this article with the same unhelpful expectations.

First, when I saw something about using stencils for factoring I thought they would be something like a Sieve of Eratosthenes, stencils with multiples of 2, 3, 5, etc. cut out. That is not what’s going on; Lehmer’s method is more sophisticated than that.

Second, related to the first, I thought that the factoring method would involve placing stencils on top of a large grid of numbers that included the number you wanted to factor. That is also not the case. The stencils can be used to factor nine-digit numbers, and a sheet of over a billion numbers would not be practical.

Third, I thought there might be some geometric significance to the position of the holes in the circular stencils. There isn’t. The stencils were circular for convenience.

How the stencils worked

Lehmer made a set of 295 stencils. Katherine Stange has put a set of SVG files for such stencils on her Github repo. The image below is the stencil corresponding to the number 42.

I want focus on how the stencils work, not why they work.

To factor a number N, you would first find a handful of numbers X₁, X₂, X₃, … X_k, related to N, numbers that are fairly easy to compute by hand. Then you would stack the stencils corresponding to the X_i and hold the stack up to the light.

If y is a quadratic residue mod N, then y is a quadratic residue mod each of the prime factors of N. Each stencil cuts the number of potential prime factors by about a half.

Once you found a probable prime factor, you would test it by long division, and then you know with certainty whether the number is a factor.

What is each stencil?

So what are these X_i and what are the holes in them?

The X_i are quadratic residues mod N, i.e. numbers such that

X_i = x² mod N

has a solution. This would have required some hand calculation, but not nearly as much calculation as trying to factor N unaided. A couple algorithms for finding quadratic residues are described in this video and in this PDF by Katherine Stange.

NB: The hard part isn’t finding quadratic residues mod N; you could do that by squaring random integers mod N. The hard part is finding quadratic residues that correspond to your set of available stencils.

The holes in the stencil correspond to the prime numbers up to 48593. Disk X has a hole for prime p if X is a quadratic residue mod X. If X is a quadratic residue mod p and mod N, there’s a 50-50 chance p divides N.

The number 48593 is the 4999th prime, so there are locations on each disk corresponding to 4999 primes. Lehmer arranged these locations in a spiral, but they could have been laid out in a grid or any other pattern.

Python simulation

See the Github repo referenced above for Python code to make the SVG files for the physical stencils. The Python code below simulates the operation of the stencils, not their geometry. I wrote the code to be illustrative, not to be practical efficient code for factoring integers.

First, here’s code to make our stencils. Here a 1 corresponds to a hole.

from sympy import primerange, legendre_symbol
from numpy import array

primes = array([p for p in primerange(3, 48594)])

def stencil(x):
    # stencil has a 1 in position k if x is a quadratic residue
    # modulo the kth odd prime and a 0 otherwise
    return array([1 if legendre_symbol(x, p) >= 0 else 0 for p in primes])

Next we pick an N and a few small quadratic residues mod N.

N = 19972009 
x = [2, 61, 79, 185, 238, 255, 313, 421, 491]

Multiplying the stencils together produces a 1 where every stencil has a hole.

mask = stencil(x[0])
for i in range(1, len(x)):
    mask *= stencil(x[i])
ps = mask*primes
print(ps[ps != 0])

This prints the following:

[ 97 103 1367 1999 15241 28351 35353 36847 44953]

We then test whether each of these primes is a factor of N, and in fact

19972009 = 97 × 103 × 1999.

Each stencil has 4999 primes, and if each of the 9 stencils eliminates about half the potential prime factors, we’d expect around 4999/2⁹ ≈ 10 candidate primes, which is very close to what we got.

Related posts

• Conway’s mental factoring method
• Factored random numbers
• Numbers don’t typically have many prime factors

[1] The Ingenious Physical Factoring Devices of D. N. Lehmer. Math Horizons. Volume 30 Issue 2. November 2022.

[2] D. N. Lehmer and his son D. H. Lehmer were both number theorists. I’ve referred to the later several times on the blog, most recently here. I also cited Katherine Stange a couple weeks ago.

The following prime number looks like a black-and-white image of an American flag.

The number mostly consists of the digits 1, 3, and 8, but there are a few 9’s. The following image colors the 8’s blue, the 3’s red, and the 1’s white. The background is gray so you can see the 1s.

I found this in [1]. The article includes a link to a text version of the number, but the link is broken, so I had to convert the image to text.

Unfortunately tesseract did a poor job, and so did Grok, so I ended up more or less doing the conversion by hand. I made a text version and posted it here.

See the next post for a calculation showing that the number is indeed prime.

[1] The United States of America Prime. Vadim Ponomarenko. Math Horizons, Vol. 28, No. 4 (April 2021)

In the process of writing the previous post, I ran across the Landau-Ramanujan theorem. It turned out to not be what I needed, but it’s an interesting result.

What portion of the numbers less than N can be written as the sum of the squares of two non-negative integers? Edmund Landau discovered in 1908, and Srinivasa Ramanujan independently rediscovered in 1913, that the proportion is asymptotically

c / (log N)^1/2

where c, the Landau-Ramanujan constant, equals 0.76422….

Let’s see how the number of squares less than 1000 compares to the estimate given by the theorem.

from math import sqrt, log

N = 1000
c = 0.76422
print("Predicted: ", c / sqrt(log(N)))

sumsq = N*[0]
for i in range(N):
    for j in range(N):
        n = i**2 + j**2
        if n < N:
            sumsq[n] = 1
            
print("Exact: ", sum(sumsq)/N)

The predicted proportion is 0.291 and the exact proportion is 0.330.

The script takes O(N²) time to run, so we'd need to do something more clever if we wanted to investigate very large N.

The sum of the first n odd numbers is n². For example,

1 + 3 + 5 + 7 + 9 = 25 = 5²

Here’s another way to say the same thing. Start with the list of positive integers, remove every second integer from the list, and take the cumulative sum. The resulting list is the list of squares.

We begin with

1, 2, 3, 4, 5, …

then have

1, 3, 5, 7, 9, …

and end up with

1, 4, 9, 16, 25, …

Cubes

Now start over with the list of positive integers. We’re going to remove every third integer, take the cumulative sum, then remove every second integer, and take the cumulative sum. The result is a list of cubes.

Here are the steps described above for the positive integers up to 20:

1, 2, 3, 4, 5, 6, 7, 8, 10, 11, …, 20
1, 2, 4, 5, 7, 8, 10, 11, 13, 14, 16, 17, 19, 20
1, 3, 7, 12, 19, 27, 37, 48, 61, 75, 91, 108, 127, 147
1, 7, 19, 37, 61, 91, 127
1, 8, 27, 64, 125, 216, 343

General case

In general, given an integer k we remove every kth element, take the cumulative sum. Then we reduce k by 1 and repeat. We keep doing this until k = 1, in which case we stop. The end result will be the list of kth powers.

The result at the top of the post, corresponding to k = 2, has been known since ancient times. The generalization to larger k wasn’t known until Alfred Moessner observed it in 1951. Oskar Perron proved Moessner’s conjecture later the same year. The result has been called Moesnner’s theorem or Moesnner’s magic.

Python implementation

A little Python code will make it easier to play with Moessner’s process for larger numbers.

from numpy import cumsum

def moessner(N, k):
    x = range(1, N+1)
    while k > 1:
        x = [x[n] for n in range(len(x)) if n % k != k - 1]
        # print(x)
        x = cumsum(x)
        # print(x)
        k = k - 1
    return x

To see the intermediate steps, uncomment the print statements.

If we run the following code

for k in [2, 3, 4, 5]:
    print( moessner(25, k) )

we get the following.

[1   4   9  16  25  36  49  64  81 100 121 144 169]
[1   8  27  64 125 216 343 512 729]
[1   16   81  256  625 1296 2401]
[1   32  243 1024 3125]

Related posts

• Easy to guess, hard to prove
• Recamán’s sequence