# geo.distance The great-circle distance between two latitude/longitude points, in metres, by the haversine formula, rounded to the millimetre. ## The model The Earth is taken as a sphere of radius 6,371,008.8 m, the IUGG mean radius R1 = (2a + b) / 3 of the GRS 80 ellipsoid (H. Moritz, "Geodetic Reference System 1980", Journal of Geodesy 74 (2000) 128-133; also Bulletin Géodésique 54 (1980)). Many snippets use 6,371,000 m or the equatorial 6,378,137 m; those give different answers, which is why the radius is part of the contract. A sphere is not the Earth. Against the WGS 84 ellipsoid (Vincenty, Karney) the haversine distance is off by up to about 0.5%. Use this for "how far is the nearest depot", delivery radii and sorting by distance, not for surveying. The millimetre rounding is about reproducibility, not accuracy. Longitudes do not need wrapping: 179.5 to -179.5 is 1 degree across the antimeridian, because the formula only uses sin² of half the difference. Coordinates outside -90..90 and -180..180 are errors rather than being wrapped, since they usually mean latitude and longitude were swapped. The result uses atan2(√a, √(1−a)) rather than asin(√a). For points closer than a metre or so it keeps millimetre accuracy, where the spherical law of cosines (acos) loses it. Near-antipodal points are ill-conditioned for any haversine in double precision (see Accuracy). ## What changed in 2.0.0 - **Rounding is on the exact value of the double.** 1.x multiplied by 1000 and rounded that product, so a distance stored just below a millimetre tie could round up: the double 14704080.96649999916... (Tokyo Haneda, 35.5494, 139.7798, to -33.9399, 18.786034 near Stellenbosch) times 1000 rounds to exactly 14704080966.5, and 1.x said 14704080.967. 2.0.0 rounds with `math.round-float`, half away from zero on the value actually held, and says 14704080.966, the way 2.675 to two places is now 2.67 rather than 2.68. The high-precision reference agrees: the true distance for those inputs is 14704080.96649999826... m, below the tie as well. - **Sine, cosine and arctangent come from `math.sin-cos` and `math.atan`** instead of private copies. `math.sin-cos` reduces the angle with a two-part π/2 and `math.atan` is fdlibm's arctangent, both from + − × ÷ only, so the unrounded distance can differ from 1.x in the last bits (by up to about 1e-9 m on the vectors) while staying the same double in every language. The one piece left here is a four-line atan2 for y, x ≥ 0 on top of `math.atan` (atan(y/x), and π/2 when x is 0, which only exactly antipodal points reach), since the registry has no atan2. - **Vectors.** Every 1.x vector keeps its expected answer; each was re-derived from an independent 60-digit reference (own series for π, sin, cos and atan, decimal square root), and the rounded true distance equals it. Two vectors were added that show the rounding change: Tokyo Haneda to near Stellenbosch (1.x .967, now .966) and London to 33.8688°S 151.714482°E in the Tasman Sea (17024346.35449999943..., 1.x .355, now .354). The vectors now compare exactly (`floats exact`): the three languages return the same double. Callers of 1.x who want the new rounding should move to ^2; 1.x is unchanged. ## Why the answers are identical in every language JavaScript's Math.sin, Python's math.sin and Rust's f64::sin come from different maths libraries, which are allowed to differ in the last bit. Rounding to millimetres hides that almost always, but not always: a value that lands on a rounding boundary can go either way. So this capability does not call them. Its trigonometry comes from `math.sin-cos` and `math.atan`, built only from operations IEEE 754 requires to be correctly rounded, in the same order in every language; the square roots here are IEEE's correctly rounded sqrt. No language fuses a multiply and an add behind your back (JavaScript and Python cannot; Rust does not without an explicit mul_add). The unrounded distance is therefore the same double in TypeScript, Python and Rust, and `math.round-float` rounds it the same way in all three. ## Accuracy Measured against the 60-digit reference, with the inputs taken as the exact doubles given, on 5,500 uniformly random pairs plus 1,500 of each special kind (coordinates to 6 decimals): - anywhere on the globe: the unrounded distance is within 4e-8 m of the true spherical distance (the largest errors are on the longest distances; at 20,000 km that is about ten units in the last place); - points less than about 150 m apart: within 5e-14 m; - 1° to 7° from the antipode: within 3e-7 m; within 1° of the antipode: within 3e-6 m, because the haversine itself is ill-conditioned there. 1.x had the same error bounds; they come from converting degrees to radians and from the formula's conditioning, not from the sine or arctangent. On every sampled pair the rounded answer equalled the true distance rounded half away from zero (1.x missed one, near the antipode). It cannot be guaranteed: when the true distance is within a few nanometres of a millimetre tie, the double can sit on the other side of it (35.5494, 139.7798 to -33.9399, 18.933865 is truly 14691235.9925000015... m, but its double is 14691235.99249999970... m, so the answer is 14691235.992). The rounding is exact on the double; the double is within the bounds above of the truth. ## Rounding Half away from zero to 3 decimals, by `math.round-float`, on the exact value of the double: the nearest multiple of 0.001 to the value held, an exact tie going away from zero. Not Python's round(), which rounds halves to even, and not floor(d × 1000 + 0.5), which rounds the product first. A result is never -0.