Latest posts
- Mr. Bessel’s eponymous functionsJun 05, 2026John
Yesterday I wrote a post showing that the trapezoid rule evaluates the integral very efficiently. But how do we know what the exact integral is for comparison? If you ask Mathematica, it will tell you the integral equals −2π J1(1) where J1 is a Bessel function. This may seem like rabbit out of a hat, […] Mr. Bessel’s eponymous functions first appeared on John D. Cook.
- The Latin of LinuxJun 04, 2026John
One reason people study Latin is that it is the ancestor of many modern languages. English derives from West Germanic languages, not from Latin, but much of English vocabulary, perhaps as much as 60%, derives from Latin, either directly or indirectly through French. Knowing a bit of Latin makes sense of many things that would […] The Latin of Linux first appeared on John D. Cook.
- Integrating smooth periodic functionsJun 04, 2026John
Several posts lately have looked at the function f(x) = cos(sin(x) + x). This post will look at the function from a different angle. It’s a smooth function with period 2π. For reasons I wrote about here, this means that the trapezoid rule should find its integral very efficiently. In general, the error in the trapezoid […] Integrating smooth periodic functions first appeared on John D. Cook.
- Partitions over permutationsJun 04, 2026John
I was thinking more about the cosine approximation to the Gaussian exp(−z²) ≈ (1 + cos(sin(z) + z))/2 that I wrote about last week. The two expressions above are close along the real axis but not along the imaginary axis. If z = iy, the right side grows much faster than the left, behaving like exp(exp(y)). This led to […] Partitions over permutations first appeared on John D. Cook.
- Naively summing an alternating seriesJun 03, 2026John
Suppose you run across the power series for the exponential function and decide to code it up. Good idea: you’ll probably learn something, though maybe not what you expect. Maybe you decide a tolerance of 10−12 is good enough, and so you sum the terms until the next term to add is below the tolerance. […] Naively summing an alternating series first appeared on John D. Cook.
- It’s not just Taylor seriesJun 01, 2026John
There is still active discussion on X about the approximation exp(−x²) ≈ (1 + cos(sin(x) + x))/2 and some are saying this can just be explained by Taylor series: the series for the two sides differ for the first time at the x6 term, so that’s why you get a good approximation. As I wrote […] It’s not just Taylor series first appeared on John D. Cook.
- Subscribe by emailJun 01, 2026John
Readers have subscribed to this blog via email almost from its beginning in 2008, but how they have subscribed has changed several times. I’ve used several services to provide email subscription that have come and gone. For the past two years I’ve been using Substack to send out emails announcing new blog posts. That has […] Subscribe by email first appeared on John D. Cook.
- Another Gaussian approximationMay 31, 2026John
The function (1 + cos(x))/2 gives a fair approximation to the Gaussian density exp(−x²) You can make the approximation much better by raising it to a power. The function ((1 + cos(x))/2)4 gives a good lower bound and ((1 + cos(x))/2)3.5597 gives a good upper bound. More on that here. There are other ways of […] Another Gaussian approximation first appeared on John D. Cook.
- Spot checking polynomial identitiesMay 30, 2026John
If a polynomial identity holds at a few random points, it’s very like true. We’ll make this statement more precise, but first let’s look at some applications. You may want to test an identity that naturally presents itself as a statement that two polynomials are equal. Or you might use something like the binomial coefficient […] Spot checking polynomial identities first appeared on John D. Cook.
- Online (one-pass) algorithmsMay 29, 2026John
Canonical example The sample variance of a set of numbers is defined in terms of the sum of the squared distances from each point to the mean. So it would seem that you first need to calculate the mean, then go back and compute the squared differences from the mean. And yet sample variance can […] Online (one-pass) algorithms first appeared on John D. Cook.
- Turning K-L divergence into a metricMay 28, 2026John
Kullback-Leibler divergence Kullback-Leibler divergence is defined for two random variables X and Y by K-L divergence is non-negative, and it’s zero if and only if X and Y have the same distribution. But it is not a metric, for reasons explained here. For one thing, it’s not symmetric. Jeffreys divergence We can fix the symmetry problem by […] Turning K-L divergence into a metric first appeared on
- The Meta logo and fitting Besace curvesMay 27, 2026John
I saw a post yesterday saying that the Meta logo is a Besace curve. A Besace curve has the implicit form and the parametric form where t ranges over [0, 2π]. So given a Besace curve, such as the Meta logo, how do you find the parameters a and b to fit the curve? We can rewrite […] The Meta logo and fitting Besace curves first appeared on John D. Cook.
- Calculating the expected range of normal samplesMay 26, 2026John
The previous post looked at the expected IQ range in a jury of 12. This post will look more generally at computing the expected range of n samples from a N(0, 1) random variable. This will give the expected range in units of σ, i.e. multiply the results by σ if your σ isn’t 1. As mentioned […] Calculating the expected range of normal samples first appeared on John D. Cook.
- Expected IQ spread on a juryMay 26, 2026John
There’s been some discussion online lately about how a large difference in IQ makes it difficult for two people to communicate. There have been studies that confirm this effect. The difficulty is not insurmountable, but it takes deliberate effort to overcome. Someone dismissed this communication difficulty by pointing out that the expected difference in IQ […] Expected IQ spread on a jury first ap
- Hilbert transform as an infinite matrixMay 23, 2026John
The previous post linked to a post I wrote a few years ago about the Hilbert transform and Fourier series. That post says that if the Fourier series of a function is then the Fourier series of its Hilbert transform is When I looked back at that post I thought about how if you thought […] Hilbert transform as an infinite matrix first appeared on John D. Cook.
- Real and imaginary partsMay 23, 2026John
The previous post announced some notes I wrote up based on an article by Henry Baker implementing functions of a complex variable in terms of functions of a real variable. That is, it finds functions u(x, y) and v(x, y) such that f(x + iy) = u(x, y) + i v(x, y) where x, y, u, and v are […] Real and imaginary parts first appeared on John D. Cook.
- Building complex functions out of real partsMay 23, 2026John
A couple months ago I wrote about how to compute the sine and cosine of a complex number using only real functions of real variables using the equations You can do something analogous for all the elementary functions, though some of the equations are quite a bit more complicated than the ones above. See the […] Building complex functions out of real parts first appeared on John D. Cook.
- Couth and uncouth function pairsMay 21, 2026John
“You can’t always get what you want. But sometimes you get what you need.” — The Rolling Stones Circular functions and hyperbolic functions aren’t invertible, but we invert them anyway. These functions map many points in the domain to each point in the range, and we invert them by mapping a point in the range […] Couth and uncouth function pairs first appeared on John D. Cook.
- Circular and hyperbolic functions differ by rotationsMay 21, 2026John
The difference between a circular function and a hyperbolic function is a rotation or two. For example, cosh(z) = cos(iz). You can read that as saying that to find the hyperbolic cosine of z, first you rotate z a quarter turn to the left (i.e. multiply by i) and then take the cosine. For another example, […] Circular and hyperbolic functions differ by rotations first appeared on John D. Cook.
- Square root of x² − 1May 20, 2026John
How should we define √(z² − 1)? Well, you could square z, subtract 1, and take the square root. What else would you do?! The question turns out to be more subtle than it looks. When x is a non-negative real number, √x is defined to be the non-negative real number whose square is x. When x is […] Square root of x² − 1 first appeared on John D. Cook.
- Closer look at an identityMay 19, 2026John
The previous post derived the identity and said in a footnote that the identity holds at least for x > 1 and y > 1. That’s true, but let’s see why the footnote is necessary. Let’s have Mathematica plot The plot will be 0 where the identity above holds. The plot is indeed flat for x > 1 […] Closer look at an identity first appeared on John D. Cook.
- Approximating Markov’s equationMay 19, 2026John
Markov numbers are integer solutions to x² + y² + z² = 3xyz. The Wikipedia article on Markov numbers mentions that Don Zagier studied Markov numbers by looking the approximating equation x² + y² + z² = 3xyz + 4/9 which is equivalent to f(x) + f(y) = f(z) where f(t) is defined as arccosh(3t/2). It wasn’t clear to me why the […] Approximating Markov’s equation first appeared on John D. Cook.
- Recovering the state of xorshift128May 15, 2026John
I’ve written a couple posts lately about reverse engineering the internal state of a random number generator, first Mersenne Twister then lehmer64. This post will look at xorshift128, implemented below. import random # Seed the generator state a: int = random.getrandbits(32) b: int = random.getrandbits(32) c: int = random.getrandbits(32) d: int = random.getrandbits(32) MASK = […] Recovering the st
- Initialize and print 128-bit integers in CMay 12, 2026John
If you look very closely at my previous post, you’ll notice that I initialize a 128-bit integer with a 64-bit value. The 128-bit unsigned integer represents the internal state of a random number generator. Why not initialize it to a 128-bit value? I was trying to keep the code simple. A surprising feature of C […] Initialize and print 128-bit integers in C first appeared on John D. Cook.
- Hacking the lehmer64 RNGMay 12, 2026John
A couple days ago I wrote about hacking the Mersenne Twister. I explained how to recover the random number generator’s internal state from a stream of 640 outputs. This post will do something similar with the lehmer64 random number generator. This generator is very simple to implement. Daniel Lemire found it to be “the fastest […] Hacking the lehmer64 RNG first appeared on John D. Cook.
- Euler functionMay 12, 2026John
This morning I wrote a post about the probability that a random matrix over a finite field is invertible. If the field has q elements and the matrix has dimensions n × n then the probability is In that post I made observation that p(q, n) converges very quickly as a function of n [1]. […] Euler function first appeared on John D. Cook.
- Inverse shiftMay 11, 2026John
What is the inverse of shifting a sequence to the right? Shifting it to the left, obviously. But wait a minute. Suppose you have a sequence of eight bits abcdefgh and you shift it to the right. You get 0abcdefg. If you shift this sequence to the left you get abcdefg0 You can’t recover the […] Inverse shift first appeared on John D. Cook.
- Probability that a random binary matrix is invertibleMay 11, 2026John
The two latest posts have involved invertible matrices with 0 and 1 entries. If you fill an n × n matrix with 0s and 1s randomly, how likely is it to be invertible? What kind of inverse? There are a couple ways to find the probability that a binary matrix is invertible, depending on what […] Probability that a random binary matrix is invertible first appeared on John D. Cook.
- The linear algebra of bit twiddlingMay 10, 2026John
The previous post looked at the tempering step of the Mersenne Twister, formulating a sequence of bit operations as multiplication by a matrix mod 2. This post will look at the components more closely. The theorems of linear algebra generally hold independent of the field of scalars. Typically the field is ℝ or ℂ, but […] The linear algebra of bit twiddling first appeared on John D. Cook.
- Reverse engineering Mersenne Twister with Linear AlgebraMay 10, 2026John
The Mersenne Twister (MT) is a random number generator with good statistical properties but bad cryptographic properties. In buzzwords, it’s a PRNG but not a CSPRNG. This post will show how the internal state of a MT generator can be recovered from its output. We’ll do this using linear algebra. The bit twiddling approach is […] Reverse engineering Mersenne Twister with Linear Algebra first appear