The idea of using computers to assist with proofs is not new. The first major computer-assisted proof was published in 1976, the proof of the four color theorem by Kenneth Appel and Wolfgang Haken. The authors reduced the proof of the four color theorem to verifying calculations on 1,834 configurations, each checked by a computer program.

The proof was simplified over the years, and formalized in Coq in 2005. Everyone is satisfied that the theorem is true, but there has never been a satisfying proof, one that a human could read and say “I see now why any map can be colored using only four colors.” And there may never be one, but see this post for a contrary prediction.

The IBM mainframe that ran the calculations completing the proof of the four color theorem did not generate the proof. It simply executed the FORTRAN program that Haken and Appel (and Koch [1]) gave it.

I don’t see the recent proof of finite-time blowup for solutions to the Navier-Stokes equations as entirely different. Computers did higher-level tasks for the OpenAI team than the mainframe did for Haken and Appel, and these tasks were not as directly programmed as the tasks that were given to the mainframe, but still machines do what they are told to do.

Related posts

• The seven color map theorem
• Detecting errors with the four color theorem

[1] John A. Koch was a programmer who worked on the four color proof with Haken and Appel. I don’t know how much credit he deserves, but I suspect it may be more than he was given.

I like AI-generated ASCII diagrams. Because nobody would ask AI to generate ASCII diagrams, and so, it’s congruous. I like incongruity [1].

Aside from the incongruity of using a gazillion-parameter neural network to make 1970’s style ASCII art, ASCII diagrams have some uses. They’re absolutely tiny compared to image files. But more importantly they can be inserted into plain text files, such as source code or markdown. A diagram embedded directly into a source file cannot become separated from the code.

ASCII diagrams are tedious to create, though there are tools to mitigate the tedium. But if an AI can generate the diagram, the tedium goes away.

I was curious how well Claude could create ASCII diagrams, so I tried a few examples. I hope these render well in whatever format you’re reading this post. They look fine for me previewing the post in a browser. I expect they might not turn out so well in an RSS reader.

I asked it to reproduce the graphs from my recently post on the graph imbalance theorem and the first diagram turned out nicely.

                 +-------+                                +-------+
                 |   A   |--------------------------------|   B   |
                 +-------+                                +-------+
                     |                                        |
                     |                                        |
   ------------------|------------------             ---------|---------
   |        |        |        |        |             |        |        |
   |        |        |        |        |             |        |        |
+-----+  +-----+  +-----+  +-----+  +-----+       +-----+  +-----+  +-----+
| a0  |  | a1  |  | a2  |  | a3  |  | a4  |       | b0  |  | b1  |  | b2  |
+-----+  +-----+  +-----+  +-----+  +-----+       +-----+  +-----+  +-----+

The second network is more complicated and so the corresponding ASCII diagram is hard to read.

   +----------------------------------------------------------------+
   |                                                                |
   |+-----------------------------------------------+               |
   ||                                               |               |
   ||+-------------------------------+              |               |
  +-------+       +-------+       +-------+       +-------+       +-------+
  |  R1   |-------|  R2   |-------|  R3   |-------|  R4   |-------|  R5   |
  +-------+       +-------+       +-------+       +-------+       +-------+
      |             | | |          |   |           |   |           |  |  |
     ++             | | |          |   |           |   |           |  |  |
     | +---------------------------+   |           |   |           |  |  |
     | |            ++| |              |           |   |           |  |  |
     | |             || +--------------+           |   |           |  |  |
     | |             || |  +---------------------------------------+  |  |
     | |             |+-|--|-------------+         |   |              |  |
     | |             |  |  |             |  +------+   |              |  |
     | |             |  |  |             |  |  +----------------------+  |
     | |             |  +--|-------------|--|--|-------------+           |
     | |             |  |  |             |  |  |       +-----|--+        |
     | |             |  |  |             |  |  |             |  |  +-----+
     | |             |  |  |             |  |  |             |  |  |
  +-------+         +-------+           +-------+           +-------+
  |  G1   |         |  B1   |           |  B2   |           |  B3   |
  +-------+         +-------+           +-------+           +-------+

For a third example, here is a fairly complicated diagram that nevertheless lends itself to a readable ASCII diagram. It’s a Feistel network diagram for DES encryption.

   +-------------+                    +-------------+
   |   L(i-1)    |                    |   R(i-1)    |--------
   +-------------+                    +-------------+       |
          |                                  |              |
          |                                  |              |
          |                      +-----------------------+  |
          |                      |   E (expand 32->48)   |  |
          |                      +-----------------------+  |
          |                                  |              |
          |                      +-----------------------+  |
          |                      |     XOR with K(i)     |  |
          |                      +-----------------------+  |
          |                                  |              |
          |                      +-----------------------+  |
          |                      |    S-boxes S1..S8     |  |
          |                      +-----------------------+  |
          |                                  |              |
          |                      +-----------------------+  |
          |                      |    P (permutation)    |  |
          |                      +-----------------------+  |
          |                                  |              |
          |                                  |              |
          |            +-------+             |              |
          +------------|  XOR  |-------------+              |
                       +-------+                            |
                           |                                |
          +----------------|--------------------------------+
          |                +------------------+
          |                                   |
   +-------------+                    +-------------+
   |    L(i)     |                    |    R(i)     |
   +-------------+                    +-------------+

Related posts

• The worst tool for the job
• ASCII armor
• Memorizing the ASCII table

[1] See Christian Wolff’s discussion of dogs playing poker in The Accountant (2016).

At the end of my post on how not to calculate cosine I said that the first step in calculating cosine, particularly cosine of a large number, would be to do range reduction. This post will present a simple range reduction method by Cody and Waite that is adequate for moderately large arguments.

If you want to compute the sine or cosine of an angle x you could start by reducing x mod 2π since that would not change the result. However, accurately reducing a number mod 2π is not trivial; that’s why range reduction is an area of algorithm development.

Range reduction mod π/2

Even better would be to reduce x mod π/2. Reducing to a smaller range means that power series method, and other methods such as rational approximation, will be more efficient.

So suppose you can find an integer k such that

x − k π/2 = y

where 0 ≤ y ≤ π/2. Then sin(x) is ±sin(y) or ±cos(y), depending on k mod 4 equals 0, 1, 2, or 3.

from math import *

def reduced_sin(x, k):
     match k % 4:
        case 0: return sin(x)
        case 1: return cos(x)
        case 2: return -sin(x)
        case 3: return -cos(x)

Naive range reduction

Now let’s set x = 500. Then k = 318 because that’s the multiple of π/2 we need to subtract to bring x into range, and the sine of x should be the negative of the sine of the reduced value y because 318 = 2 mod 4.

The following code computes sin(x) with naive range reduction

def naive_sin(x):
    k = floor(x / (pi/2))
    y = x % (pi/2)
    return reduced_sin(y, k)

and when x = 500 the error is on the order of 1.7 × 10⁻¹⁴.

Better range reduction

The value of k above is fine, but we’d like to calculate y more accurately. The following code is much better.

def Cody_Waite_sin(x):
    C1 = 1686629713 / 2**30
    C2 = 4701928774853425 / 2**86

    k = floor(x / (pi/2))
    y = (x - k*C1) - k*C2
    return reduced_sin(y, k)

This will compute sin(500) to full machine precision. What kind of magic is this?

The trick is that the exact value of C1 + C2 equals π/2 to more precision than is possible in a single float [1]. You can confirm, with bc or some other extended precision software, that the difference between C1 + C2 and π/2 is roughly 2⁻⁸⁸, while the limit of float precision is 2⁻⁵².

If we compute

y = x - k*(C1 + C2)

then we’re doing the same calculation as naive_sin and will get the same error. But if we compute

y = (x - k*C1) - k*C2

we will get a more accurate result, provided x isn’t too large.

You can use the following code to play around and see how large x can be before errors start to creep in. For small enough x, like 500, the Cody and Waite sine returns full precision. For larger x it’s better than naive sine but does not return full precision. And for large enough x it completely breaks down.

def compare(x):
    y0 = naive_sin(x) 
    y1 = Cody_Waite_sin(x)
    y2 = sin(x)
    print("Naive error:     ", y2 - y0)
    print("Cody Waite error:", y2 - y1)

Now this may seem circular since we’re using math.sin as our gold standard. However, this function is calling the sine function on your CPU, which is using sophisticated range reduction to compute its result accurately down to the last bit, assuming you run the code on a computer that’s less than 40 years old.

The Cody and Waite algorithm is inadequate for large x, but it’s a good place to begin studying range reduction. It shows there are clever ways of squeezing out more precision than seems possible.

[1] The numerator n₁ of C1 is ⌊2³⁰ π/2⌋. The numerator n₂ of C2 is the solution to

2⁸⁶⁻³⁰ n₁ + n₂ = ⌊2⁸⁶ π/2⌋.

In a footnote to the previous post, I said that Python’s math library can calculate the logarithm of extremely large numbers but not the cosine. This post will expand on that comment.

In this post I’ll use n = 200! as my example rather than 1000! because this value of N is larger than the largest representable floating point number but small enough to be more convenient to work with.

Suppose someone calculates 200! for you:

78865786736479050355236321393218506229513597768717326329474253324435\
94499634033429203042840119846239041772121389196388302576427902426371\
05061926624952829931113462857270763317237396988943922445621451664240\
25403329186413122742829485327752424240757390324032125740557956866022\
60319041703240623517008587961789222227896237038973747200000000000000\
00000000000000000000000000000000000

You could now calculate log(n) using

n = 7.886578673647905 × 10³⁷⁴

and so

log(n) = log(7.886578673647905 × 10³⁷⁴)
= log(7.886578673647905) + 374 log(10) = 863.2319871924055.

The key thing that makes this possible is that the least significant digits of n only affect the least significant digits of log(n). In the calculation above I kept the first 16 digits of n. Python couldn’t make use of any more digits, and had no need of any more digits, in order to produce the logarithm to machine precision.

Cosine doesn’t work that way. The cosine of n depends on the remainder when n is divided by 2π, and that remainder depends on every single digit of n. I’ll illustrate that below.

Using bc -l and setting the scale to 400, I can calculated n then calculate

cos(n + 10ⁱ)

for i running from 0 to 374, tweaking each digit one at a time. (Except when a digit is a 9 and the addition results in a carry.)

    n = 1
    for (i = 1; i <= 200; i++) n *= i
    scale = 400
    for (i = 1; i <= 374; i++) {
        x = c(n+10^i)
        scale = 16
        print x/1, "\n"
        scale = 400
    }

Here’s what a plot of the results look like.

The value of cos(n) is about −0.985, but the values above are all over the map. We can look at the range by projecting all the points over to the left edge then rotating a quarter turn:

The remarkable thing about this image is that there are a few gaps, i.e. a few values the cosine does not take on.

Here’s a more sophisticated way to look at it. The sequence 10ⁱ mod 2π is dense in [0, 2π], and so by going far enough out in the sequence, we can find a value that shifts the phase of n by any desired amount within any given tolerance.

Every digit in n matters, and changing any digit can change the value of cosine to be essentially any value. You cannot calculate the cosine of an enormous number without using some kind of extended precision arithmetic. There are clever range reduction algorithms that minimize the amount of extended arithmetic necessary, but extended arithmetic cannot be completely eliminated.