NavList:

Tony wrote:
"I was quite surprised to discover this site a ~year ago. I did not occur to me before that calculators could fail"

I believe this is worth some discussion (Tony, please do not worry that this note is directed to you personally -- it is a note intended for general consumption).

All computational methods and tools fail on floating point operations. It's the nature of the math. It has very little to do with issues in navvigation computations. Also, the peculiar "forensic" calculation, which Tony borrowed from this website, is not something that anyone should worry about and certainly not something that you should use --unless you're interested in the purpose they intended for it, which is "fingerprinting" certain chipsets found inside commercially available calculators without tearing them apart. And for any readers who are hyper-fascinated by calculators, then that may be a useful goal. But it's completely irrelevant to issues in navigation. 

Let's do a little forensics in a different direction...

The web page on "calculator forensics" that Tony linked hasn't been updated in a dozen years, and even at a glance we can tell that the original is much older than that with a bit of "website forensics". For example, the page uses a background bitmap texture, like "wallpaper". This was common and fashionable about twenty years ago, but it's now quite old-fashioned, at least as a pure page background. For the most recent evidence of updates, we can view its source (right-click or two-finger tap in most web browsers will reveal a "view source" option), and discover that it was created and edited using a product called Microsoft Frontpage 12.0. Long defunct, that was still available about fifteen years ago. Does it matter? Is the page obsolete? No. It's archived information that reflects some useful details about calculators primarily manufactured from 20-40 years ago. But just as we will learn much that is wrong by digging around in the NavList archives from fifteen years ago, we should also be wary of fossilized webpages like this one.

Do we learn something about the infamous haversine by looking at the keystroke sequence on that webpage? No. That's really wrong. So let's do it right...

How can you test a calculator to see if it has a problem with small angles? You could try taking the sine of an angle, and then the inverse sine. But if you do, you'll discover that there is no evidence of error even down to absurdly small angles like 10^(-50) degrees. How can our devices and apps be so accurate? There are two factors at work here. First, there is the simple fact that the series expansion of the sine of an angle (expressed as a pure ration) is given by:
  sin(x) = x - x³/6 + ...,
which for small enough angles is equivalent to sin(x) = x. So apart from the process of conversion from degrees to "radians", taking the sine of a number gives the same number back. The second factor involves the representation of real numbers ("floating point" numbers) on nearly all modern computing devices. Unlike integers, which are represented "perfectly", real numbers are represented in a fashion similar to "scientific notation" with a  mantissa greater than or equal to 1 and less than 10, accompanied by an exponent indicating the powers of ten to multiply the mantissa. When the accuracy of a result is judged by significant digits, as is normal, it's not hard to see that floating point numbers in computers and calculators are just as good when they're astronomical as when they're sub-microscopic. If I take the sine of 10^(-50) degrees, it's accuracy will be no worse than the sine of 0.01 degrees.

Unlike the sine of an angle, the cosine of an angle is limited by actual calculating accuracy of the device. As with the sine, there are two reasons for this. The series expansion of the cosine for small angles is given by:
  cos(x) = 1 - x²/2 + ...,
which for small enough angles is the same as cos(x) = 1. The catch is that this implies the same unit value for the cosine for a range of small angles. If you take the cos of a small angle and then the inverse cos, you'll quickly discover that there is some small error and then below a certain limit, for every angle x, acos(cos(x)) yields 0. So how do we deal with that? First we need to know under what conditions this occurs. Then we need to know if it ever matters given the data we're dealing with. And finally we need to ask what the appropriate algorithm solution would be.

What is the limiting small angle x for which acos(cos(x)) yields zero on common devices? You can try this out on your own devices, apps, and tools. There are three broad categories here: single precision devices, modern calculators, and double precision devices. 

The earliest home computers used a type of floating point math called "single precision". Any non-integer math typically was accurate to eight significant digits (though with specialized software it was possible to go well beyond those limits). This was the standard of computing that many of us first experimented with in the early to mid 1980s, and it was normal for about a decade. Home computers soon included "math coprocessors" as optional, expensive components, and within a few years these chipsets were integrated on the main processing chip as the FPU or floating point unit. Starting around 1987, these new math coprocessors were common but still quite expensive. By 1990 this capability was integrated in the CPU. By 1993 this capability was cheap enough to be standard on almost all new home computers. After that date, single precision floating point math was fundamentally obsolete, and nearly every computation could be expected to have sixteen significant digits of accuracy. Finally modern, ultra-cheap calculators typically manage 10-12 significant digits for floating point operations, and it's interesting to note that inexpensive calculators often had better accuracy in the mid-1980s than typical home computers. 

If we turn the clock back to 1984 and imagine working a problem like a great circle distance using a single precision device, like a typical "pc clone" such as the popular Radio Shack TRS-80 (affectionately or not, often called the "Trash 80"), we will immediately discover that cosines were a problem back then. At single precision accuracy, an angle as large as 0.01° is the approximate limit for the inverse cosine computation. This is about half a minute of arc. In celestial navigation, and for great circle distances on the globe, this problem would almost never arise, but in astronomical computations, this could be a serious issue for tasks like computing the occultation of a star by an asteroid. The standard great circle distance computation would yield exactly zero angular separation (implying a perfect occultation) whenever the actual angular distance was less than half a minute of arc. In other words, an algorithm like this would yield many false positives. It was this issue that inspired Roger Sinnott (now a NavList member) to write a little article on the "Virtues of Haversines" in the "Astronomical Computing" column, then just a few months old, in the August 1984 issue of Sky & Telescope (sample attached). It was an interesting issue at this time, but even in his article, Roger Sinnott points the way to the more effective solution: a variant of the Pythagorean theorem for special cases. That article should have become a historical curiosity within a decade.

After 1993 or so, more than a quarter century ago, double precision computation mostly consigned the whole issue to the dustbin of history. There's always some angle, x, at which the acos(cos(x)) computation will always yield zero, but with double precision floating point math, that limit is now about a millionth of a degree, less than 0.01 seconds of arc. On common scientific calculators, the limit is about 2 seconds of arc.

In celestial navigation or other related navigation problems, would we ever be concerned by such issues? Only if we're doing something exotic by brute force without thinking. For an exotic example, if you compute the distance from the bow of your sailboat to the stern using the standard cosine formula and great circle distances on a handheld calculator based on the latitude and longitude of those two points as inputs, then you will find that your boat has zero length (if it's less than about 200 feet in actual length). Otherwise? Never. It just isn't an issue. 

So why do we still hear about haversines in this context in various online articles? Mostly it's because the Internet works like medieval lore -- by copying. This is the origin of viral internet rumors, and it's also the origin of lingering incorrect or misleading information in technical topics. Celestial navigation itself, as a subject and the object of interest for a social community, has a long history of "lore" and is especially susceptible to this problem. It's the "Two people said it's true, so it must be true" problem. And issues of great circle distance computation were among the earliest tackled by the Internet starting even in the late 1990s. One tiny website, in particular, has had an out-sized influence on these issues, thanks to its generally good quality and the fact that it was originally created and referenced rather early in the history of Internet calculations. Created by Chris Veness, who also first learned professional computing in the early era of the pc clones, you can read his variant, based on Roger Sinnott's article, here: https://www.movable-type.co.uk/scripts/latlong.html. This author updated his article a few times in the early 21st century, and he's quite clear that this is now essentially irrelevant. Nonetheless, the internet is "lore-driven" and his account of the haversine as a panacea for great circle distances has become part of the culture. It's a shame, too, since in his updates he clearly makes the original article pointless, but it's been referenced so many times, and probably been of tremendous value to him professionally, that there's no incentive to take it down or merely replace it. 

So what's the right way to proceed? It depends on what sorts of problems you're looking at. Do you need to calculate angular distances between points that are exceedingly close together? Then use the standard Pythagorean theorem for angular distances, which can be proved by simple series expansion of the ordinary cosine law for angular distance in the case where the separation is small. If you must have the ability to cover very large angles and very small, too, or if you're simply made "nervous" by the thought of special cases, then write your algorithm to accommodate it:

Dealing mostly with large distances and some very short distances:
Use the usual cosine law:
D = acos(....)  
And then add:
If D<minD then D = Pyth(....)

Dealing with mostly small distances (like occultations):
D = Pyth(....)
if D>minD then ignore (or, more rarely, use the cosine law if needed for some other purpose).

Here minD is simply whatever appropriate minimum angle qualifies as the cutoff, e.g. one foot (0.3 meters), for problems on the Earth's surface, if you're using standard double precision math (any modern computer). The function Pyth() simply calculates the square root of the sum of the squares, as in the common Pythagorean Theorem for distance between two points: Pyth(dLat, dLon, Lat) = sqrt(dLat² + dLon²cos²Lat). Note, too, that this Pythagorean expression can easily incorporate differences in altitude when you're dealing with points on the Earth's surface, which the great circle distance cannot address in any way.

That's how it's done. That's how it should be done. No haversines needed.

I should conclude by saying that there are other cases where haversines come up in celestial navigation, primarily in certain popular methods from the early 20th century, but that usage is not significantly, necessarily connected to this issue of the cosine of small angles.

Frank Reed

PS: Of course, it's not just trigonometric functions that face the limits of floating point precision. Even something as simple as adding a pair of numbers will show issues if you know how to display results with all internal digits (usually hidden). I found a website with an adorable domain address that provides a nice set of examples here: https://0.30000000000000004.com/.

File: