John D. Cook

Van der Pol’s differential equation is The equation describes a system with nonlinear damping, the degree of nonlinearity given by μ. If μ = 0 the system is linear and undamped, but as μ increases the strength of the nonlinearity increases. We will plot the phase portrait for the solution to Van der Pol’s equation […]

Today’s exponential sum looks sorta like a Christmas wreath with candles. Let me back up and say what “today’s exponential sum” means. I have a page that plots a graph each day, where the function being graphed takes parameters from the date. It plots the partial sums of where m is the month, d is […]

Consider a random walk on the integers mod p. At each step, you either take a step forward or a step back, with equal probability. After just a few steps, you’ll have to be near where you started. After a few more steps, you could be far from your starting point, but you’re probably not. […]

These have been the most popular privacy-related posts here this year. Hashing PII does not protect privacy Rare and strange IDC-10 codes Font fingerprinting Computed IDs and privacy Supercookies

A continuous linear system is stable if and only if all its eigenvalues have negative real parts. There are various other stability criteria, but they boil down to the statement above. For example, there is the Routh stability criterion and the Hurwitz stability criterion. There’s also a continued fraction criterion. But these criteria are just […]

Here’s a simple probability problem with an answer you may find surprising. (The statement of the problem and solution are simple. The proof is not as simple.) Suppose you draw 1,000 serial numbers at random from a set of 1,000,000. Then you make another random sample of 1,000. How likely is it that no numbers […]

Since you can describe a point in the plane with two numbers, why would you choose to use three numbers? Why would you ever want to use a coordinate system with more coordinates than necessary? Barycentric coordinates One way to indicate the location of a point inside a triangle is to give the distance to […]

My daughter had a homework problem the other day that gave the frequencies of several Fourier components and asked her to find the fundamental frequency. The numbers were nice enough that brute force worked, and I’m sure that’s what students were expected to do. But this could easily be a much more sophisticated problem. If […]

Top blog posts this year about command line tools. The hard part in becoming a command line wizard Computational survivalist Computing π with bc Set theory at the command line Working with wide text files Random sampling from a file

If you take a square matrix M, subtract x from the elements on the diagonal, and take the determinant, you get a polynomial in x called the characteristic polynomial of M. For example, let Then The characteristic equation is the equation that sets the characteristic polynomial to zero. The roots of this polynomial are eigenvalues […]

Suppose you’re looking for instances of -42 in a file foo.txt. The command grep -42 foo.txt won’t work. Instead you’ll get a warning message like the following. Usage: grep [OPTION]... PATTERN [FILE]... Try 'grep --help' for more information. Putting single or double quotes around -42 won’t help. The problem is that grep interprets 42 as […]

If you learned regular expressions by using a programming language like Perl or Python, you may be surprised when tools like grep seem broken. That’s because what you think of as simply regular expressions, these tools consider extended regular expressions. Tell them to search on extended regular expressions and some of your frustration will go […]

Toward the end of each year I write a post or two listing the most popular posts by category. This year the categories will be a little different. I’ll start by listing my most popular posts about cryptography this year. Reversing an MD5 hash Naming elliptic curves Inside the AES S-box What is an elliptic […]

How hard is it to factor large numbers? And how secure are encryption methods based on the difficulty of factoring large numbers? The RSA factoring challenges were set up to address these questions. Last year RSA-230 was factored, and this week RSA-240 was factored. This is a 240 digit (795 bit) number, the product of […]

Last night I was helping my daughter with calculus homework. I told her that a common mistake was to forget what the original problem was after getting absorbed in sub-problems that have to be solved. I saw this over and over when I taught college. Then a few minutes later, we both did exactly what […]

I recently got a review copy of Statistics, Data Mining, and Machine Learning in Astronomy. I’m sure the book is especially useful to astronomers, but those of us who are not astronomers could use it as a survey of data analysis techniques, especially using Python tools, where all the examples happen to come from astronomy. […]

Alex Kontorovich had an interesting tweet about the Riemann hypothesis a few days ago. Why is the Riemann Hypothesis hard? Just one reason (of very many): it’s not an analytic question. Here are the values of zeta on the 1/2-line (where at least 40% of the zeros are, all should be) and the 4/5-line, where […]

The Sierpinski triangle has come up several times on this blog. Here’s yet another way to produce a Sierpinski triangle, this time by taking the bitwise-and of two integers. The ith bit of x&y is 1 if and only if the ith bit of x and the ith bit of y are both 1. The […]

Suppose ages in some database are reported in decades: 0, 10, 20, etc. You need to add a 27 year old woman to the data set. How do you record her age? A reasonable approach would to round-to-nearest. In this case, 27 would be rounded up to 30. Another approach would be stochastic rounding. In […]

Often you have two lists, and you want to know what items belong to both lists, or which things belong to one list but not the other. The command line utility comm was written for this task. Given two files, A and B, comm returns its output in three columns: A – B B – […]

I’ve been reading Avi Wigdenson’s new book Mathematics and Computation. It brings a lot of interesting ideas together in compact but readable form. A couple days ago I quoted Widgerson’s comment about problems in class P often having small exponents. Here I wanted to point out another interesting comment [1]. The chapter on quantum computing […]

A crinkle crankle wall, also called a serpentine wall, is a wavy wall that may seem to sacrifice some efficiency for aesthetics. The curves add visual interest, but they use more material than a straight way. Except they don’t! They use more bricks than a straight wall of the same thickness but they don’t have […]

An m by n rook graph is formed by associating a node with each square of an m by n chess board and connecting two nodes with an edge if a rook can legally move from one to the other. That is, two nodes are connected if they are either in the same row or […]

The P vs NP conjecture has always seemed a little odd to me. Or rather the interpretation of the conjecture that any problem in P is tractable. How reassuring is to to know a problem can be solved in time less than some polynomial function of the size of the input if that polynomial has […]

Programming languages differ in the names they use for inverse trig functions, and in some cases differ in their return values for the same inputs. Choice of range If sin(θ) = x, then θ is the inverse sine of x, right? Maybe. Depends on how you define the inverse sine. If -1 ≤ x ≤ […]

I recently had to work with the function f(x; k) = arctan( k tan(x) ) The semicolon between x and k is meant to imply that we’re thinking of x as the argument and k as a parameter. This function is an example of the coming full circle theme I’ve written about several times. Here’s […]

Phys.com published an article a couple days ago NASA finds Neptune moons locked in ‘dance of avoidance’. The article is based on the scholarly paper Orbits and resonances of the regular moons of Neptune. The two moons closest to Neptune, named Naiad and Thalassa, orbit at nearly the same distance, 48,224 km for Naiad and […]

Suppose you want to compute the natural logarithms of every floating point number, correctly truncated to a floating point result. Here by floating point number we mean an IEEE standard 64-bit float, what C calls a double. Which logarithm is hardest to compute? We’ll get to the hardest logarithm shortly, but we’ll first start with […]

The Sketchviz site lets you enter GraphViz code and get out something that looks hand-drawn. I decided to try it on some of the GraphViz diagrams I’ve made for this blog. Here’s a diagram from my post on Rock, Paper, Scissors, Lizard, Spock. The Sketchviz site defaults to Tinos font, which does not look hand-drawn. […]

If you have two sequences, and one is a permutation of the other, how do you measure how far apart the two sequences are? This post was motivated by the previous post on anagram statistics. For a certain dictionary, the code used in that post found that the longest anagram pairs were the following: certifications, […]

An anagram of a word is another word formed by rearranging its letters. For example, “restful” and “fluster” are anagrams. How common are anagrams? What is the largest set of words that are anagrams of each other? What are the longest words which have anagrams? You’ll get different answers to these questions depending on what […]

Many of you have been following my blog for years. Thank you for your time and your feedback. But there are always newcomers, and so periodically I like to write a sort of orientation post. You can subscribe to this blog by RSS or email. I also have a monthly newsletter where I highlight the […]

The previous post looked at an example of strange floating point behavior taking from book End of Error. This post looks at another example. This example, by Siegfried Rump, asks us to evaluate 333.75 y6 + x2 (11 x2 y2 – y6 – 121 y4 – 2) + 5.5 y8 + x/(2y) at x = […]

Here is an example fiendishly designed to point out what could go wrong with floating point arithmetic, even with high precision. I found the following example by Jean-Michel Muller in John Gustafson’s book End of Error. The task is to evaluate the following function at 15, 16, 17, and 9999. Here e(0) is defined by […]

A killer app is an program so useful that people will buy a larger system just to use it. For example, the VisiCalc spreadsheet was a killer app for the Apple II, maybe the first program to be called a “killer app.” People would buy an Apple II just so they could run VisiCalc. Microsoft […]

Robert Bosch has a new book coming out next week titled OPT ART: From Mathematical Optimization to Visual Design. This post will look at just one example from the book: creating images of Frankenstein’s monster [1] using dominoes. The book includes two images, a low-resolution image made from 3 sets of dominoes and a high […]

A couple weeks ago I wrote about how De Bruijn sequences can be used to attack locks where there is no “enter” key, i.e. the lock will open once the right symbols have been entered. Here’s a variation on this theme: what about locks that let you press more than one button at a time? […]

Adding up an array of numbers is simple: just loop over the array, adding each element to an accumulator. Usually that’s good enough, but sometimes you need more precision, either because your application need a high-precision answer, or because your problem is poorly conditioned and you need a moderate-precision answer [1]. There are many algorithms […]

Any positive number can be found at the beginning of a factorial. That is, for every positive positive integer n, there is an integer m such that the leading digits of m! are the digits of n. There’s a tradition in math to use the current year when you need an arbitrary numbers; you’ll see […]

Having access to arbitrary precision arithmetic does not necessarily make numerical calculation difficulties go away. You still need to know what you’re doing. Before you extend precision, you have to know how far to extend it. If you don’t extend it enough, your answers won’t be accurate enough. If you extend it too far, you […]

Thirty years ago, a lot of US states thought it would be a good idea to compute someone’s drivers license number (DLN) from their personal information [1]. In 1991, fifteen states simply used your Social Security Number as your DLN. Eleven other states computed DLNs by applying a hash function to personal information such as […]

Yesterday I wrote about Householder’s higher-order generalizations of Newton’s root finding method. For n at least 2, define and iterate Hn to find a root of f(x). When n = 2, this is Newton’s method. In yesterday’s post I used Mathematica to find expressions for H3 and H4, then used Mathematica’s FortranForm[] function to export […]

Are there infinitely many positive integers n such that tan(n) > n? David P. Bellamy and Felix Lazebnik [1] asked in 1998 whether there were infinitely many solutions to |tan(n)| > n and tan(n) > n/4. In both cases the answer is yes. But at least as recently as 2014 the problem at the top […]

There are variations on Newton’s root finding method that use higher derivatives and converge faster. Alston Householder developed a sequence of such methods of the form where the superscripts in parentheses indicate derivatives. When n = 2, Householder’s method reduces to Newton’s method. Once Newton’s method is close enough to a root, the error is […]

The utility numfmt, part of Gnu Coreutils, formats numbers. The main uses are grouping digits and converting to and from unit suffixes like k for kilo and M for mega. This is somewhat useful for individual invocations, but like most command line utilities, the real value is using it as part of pipeline. The --grouping […]

I wanted to stress test the bc calculator a little and so I calculated π to 10,000 digits a couple different ways. First I ran time bc -l <<< "scale=10000;4*a(1)" which calculates π as 4 arctan(1). This took 2 minutes and 38 seconds. I imagine bc is using some sort of power series to compute […]

The previous post looked at an algorithm for generating De Bruijn sequences B(k, n) where k is a prime number. These are optimal sequences that contain every possible consecutive sequence of n symbols from an alphabet of size k. As noted near the end of the post, the case k = 2 is especially important […]

Last week I wrote about a way to use De Bruijn cycles. Today I’ll write about a way to generate De Bruijn cycles. De Bruijn cycles Start with an alphabet of k symbols. B(k, n) is the set of cycles of that contain every sequence of n symbols, and is as short as possible. Since […]

Google announced today that it has demonstrated “quantum supremacy,” i.e. that they have solved a problem on a quantum computer that could not be solved on a classical computer. Google says Our machine performed the target computation in 200 seconds, and from measurements in our experiment we determined that it would take the world’s fastest […]

As I mentioned in my computational survivalist post, I’m working on a project where I have a dedicated computer with little more than basic Unix tools, ported to Windows. It’s given me new appreciation for how the standard Unix tools fit together; I’ve had to rely on them for tasks I’d usually do a different […]