I’ve been tinkering with the α·max + β·min approximation to the distance formula. As a quick refresher, if you know the x and y distances to some point and want to estimate the straight-line distance, you can take the dot product of ⟨x, y⟩ with an appropriately chosen vector ⟨α, β⟩. Since the dot product equals the product of the magnitudes times the cosine of the angle between the vectors, then as long as the vectors are nearly parallel and the second vector is close to unit length, the result will be approximately the magnitude of the first vector.

The “max” and “min” parts comes in when you let u = max(|x|, |y|) and v = min(|x|, |y|), so ⟨u, v⟩ is essentially ⟨x, y⟩ reflected into the lower half of the first quadrant. In other words, ⟨u, v⟩ lies along the arc from 0° to 45°. For simplicity of notation, let us call that vector u = ⟨u, v⟩. By symmetry we should choose the vector α = ⟨α, β⟩ to have angle 22.5°. Then the angle between α and u will have cosine at least cos(π/8) = ½√(2+√2) ≈ 0.924, so our estimate will be pretty good.

But what if want the best possible estimate? Well then we need to figure out the optimal length for α = c·⟨cos(π/8), sin(π/8)⟩. And we need to decide what we’ll optimize for: maximum error, root-mean-squared error, average error?

• • •

Maximum error:

Spoiler:

The traditional choice is to minimize the maximum error, meaning to choose a length for α which minimizes max(|1 − α•u|) over all vectors u of unit length on the 45° arc. Clearly the maximum error occurs either at the midpoint (22.5°, where the vectors are colinear) or at the endpoints (0° and 45°, where the vectors are misaligned by 22.5°). The error at the midpoint is |1 − c| and at the ends it is |1 − c·cos(π/8)|. If one of these is larger than the other, then we can adjust c to shrink it while still keeping it larger than the other, so the optimum must occur when both are equal.

Moreover, the errors will have opposite signs, and we will have c ≥ 1, since otherwise the errors at the endpoints could not equal the error in the middle. Hence |1 — c| = c − 1 and |1 − c·cos(π/8)| = 1 − c·cos(π/8) and we want these to equal each other, so c − 1 = 1 − c·cos(π/8) and therefore c_{max} = 2 / (1 + cos(π/8)). Thus we have found the optimal vector to minimize the maximum error:

If we run the numbers we find that, for this choice of α we have:

3.9566% maximum error (minimum possible on 45° arc) 2.7001% RMS error 2.4083% average error

RMS error:

Spoiler:

If we want to minimize the root-mean-squared error we can find the c that minimizes ∫ (1 − c·cos(θ))^{2} dθ from 0 to π/8, which turns out to be c_{rms} = 2 / ((π/8)/sin(π/8) + cos(π/8)). The numeric values are:

α_{max} ≈ ⟨0.9475, 0.3925⟩

5.2456% maximum error 2.3325% RMS error (minimum possible on 45° arc) 2.0005% average error

Average error:

Spoiler:

And similarly the average error can be minimized by choosing c to minimize ∫ |1 − c·cos(θ)| dθ from 0 to π/8, which is somewhat trickier due to the absolute values but ends up being c_{avg} = 2 / √(3 + cos^{2}(π/8)). The corresponding numeric values are:

α_{avg} = ⟨0.9414, 0.3899⟩

5.8729% maximum error 2.4246% RMS error 1.9458% average error (minimum possible on 45° arc)

• • •

Now this is all well and good, but in practical use the purpose of the estimate is to speed up computation. If we have to do 2 multiplies, then we might as well just calculate the squared distance and work with that, at least for most purposes. No, the real benefit of this method comes in when we can replace the multiplies by a small number of bitshifts and adds. Specifically, we only want to multiply by powers of 2.

Numerical stuff:

Spoiler:

If we let α = ⟨15/16, 15/32⟩ = ⟨0.9375, 0.46875⟩ then we are are at least in the right ballpark for all our estimates, and now we can quickly calculate α•u = αu + βv = (15/16)u + (15/32)v by doing the following:

let w = u + v/2 return w − w/16

Since dividing by powers of 2 is computationally cheap (equivalent to a bitshift for integers, or subtracting exponents for floating point) now we can get a good approximation of the distance formula ridiculously fast. How good is it though? Well with this α we have:

6.2500% maximum error 3.4580% RMS error 3.0825% average error

To the best of my knowledge these are all optimal for only 4 arithmetic operations.

• • •

If we drop down to 2 operations then the best we can get is:

α = ⟨1, 1/4⟩ 11.6117% maximum error 4.1548% RMS error 3.2026% average error

But if we go up to 6 operations, then we can do substantially better. For the max error:

α = ⟨31/32, 97/256⟩ = ⟨0.96875, 0.37890625⟩ 4.7063% maximum error (minimum I know of with 6 operations) 2.7609% RMS error 2.4557% average error

Realized as: let w = u − v/8 let z = w − w/32 return z + v/2

For the RMS error:

α = ⟨61/64, 3/8⟩ = ⟨0.953125, 0.375⟩ 6.0874% maximum error 2.3727% RMS error (minimum I know of with 6 operations) 1.9969% average error

Realized as: let w = v/4 − u/32 let z = z + z/2 return u + z

And for the average error:

α = ⟨15/16, 13/32⟩ = ⟨0.9375, 0.40625⟩ 6.2500% maximum error 2.3938% RMS error 1.9693% average error (minimum I know of with 6 operations)

Realized as: let w = u + v/2 let z = w − w/16 return w − v/16

Greater numbers of operations of course allow ever-finer approximations, whereas fewer operations run faster.

• • •

Of course, we could do the same thing over other arcs besides 45°. For example, if we check whether u < 2*v then we split the 45° arc into two arcs, one from 0 to arctan(1/2) and another from arctan(1/2) to π/4. We can do similar calculations on those arcs, or any other. For an arc which subtends an angle δ, the optimal α vectors (facing the midpoint of the arc) have lengths:

Good approximations can be found for any particular arc of interest, though the mere act of splitting the circle into arcs entails additional comparisons. The cheapest are along the lines of (u < 2*v) or (u < 4*v).

• • •

I don’t really have a grand conclusion here, I just felt like sharing what I’ve been toying around with. If you’re looking for a challenge though, I’d appreciate having someone check my calculation for c_{avg} minimizing the average error.