geo3d NORMAL and NORMAL_SQUARED point distances are wrong on WGS84 because Vector.normalDistance assumes a unit-length vector

Vector.normalDistanceSquared computes the perpendicular distance from a vector v to the line through the origin and the receiver p as |p * (p . v) - v|, without dividing the projection by |p|^2. That is right only when |p| = 1, which holds for every surface point on PlanetModel.SPHERE but not on WGS84, CLARKE_1866 or other ellipsoids, where surface points are 0.998 to 1.0011 units long. So on WGS84 p.normalDistance(p) is |p| * ||p|^2 - 1| instead of 0: 0.0016507737062451833 planet units (10.5 km) at latitude 0.3 rad, 14.3 km at the equator, 28.4 km at the poles, and 0 only near 35 degrees (geocentric latitude), where |p| = 1. DistanceStyle.NORMAL and NORMAL_SQUARED point distances are exactly these calls, and circles and paths measure them from their centre or path points, so a WGS84 circle's NORMAL distance for a point 50 m from its centre is 0.00165 instead of 7.85e-6, and a path-centre comparison ranks a point on the path behind one 1.63 km off it. The error is a floor that combines in quadrature with the true distance, so near distances and small thresholds are hit hardest. Nothing fails loudly. The formula comes from the sphere-only first geo3d code; the 2015 change that added WGS84 (5.3.0) made arc distances magnitude-aware but kept this method's sphere-era formula. Lucene's own queries and sorts use ARC, so only applications that call the geom distance API with NORMAL or NORMAL_SQUARED on an ellipsoid are hit. Dividing t by |p|^2 fixes it.

1. The formula. Vector.normalDistanceSquared(Vector v) (Vector.java:424-430) computes double t = dotProduct(v); (:425), then this * t - v and its squared length; the (x, y, z) overload is the same (:441-447, t at :442). normalDistance is the square root (:456-458, :469-471). There is no division by |this|^2.
2. What it computes. Write p = |p| u with u of length 1 and split v into its part along u and the perpendicular rest v_perp. Then p * (p . v) - v = (|p|^2 - 1)(u . v) u - v_perp, so, apart from rounding, the computed square is the true perpendicular distance squared plus ((|p|^2 - 1)(u . v))^2. Hence the result is never below the true distance, nor below ||p|^2 - 1| * |v|, and for v = p it is |p| * ||p|^2 - 1|. The error is a floor combined in quadrature, not a constant offset: near p the result is about sqrt(true^2 + floor^2) with floor = |p| * ||p|^2 - 1|, so the excess shrinks with distance (at latitude 0.3: 10.5 km at 0 and 50 m, 551 m at 100 km, 54 m at 1,000 km). It is zero only when |p| = 1 or p . v = 0. From a receiver at latitude 0.3 no WGS84 surface point is closer than 0.00165 (10.5 km).
3. Who calls it. NormalDistance.computeDistance(point1, x2, y2, z2) returns point1.normalDistance(x2, y2, z2) (NormalDistance.java:33-41), and NormalSquaredDistance the squared form (NormalSquaredDistance.java:33-41); these are the only production callers. Circles call distanceStyle.computeDistance(this.center, x, y, z) (GeoStandardCircle.java:158, GeoExactCircle.java:235), and paths use their segment end points and path points as receivers (GeoStandardPath.java:914-948, :1593, :1688, :1745, :1778-1780; GeoDegeneratePath.java:504-536, :663, :713). On paths, computeDistance, computeNearestDistance, computePathCenterDistance and computeDeltaDistance reach these point-to-point calls; computeOutsideDistance is mixed: segment planes (exact) plus segment corners and end points (affected) (GeoStandardPath.java:1784-1838, :946-949; GeoDegeneratePath.java:866-874, :534-537). Path NORMAL_SQUARED values also carry two separate path bugs, LFV-271 (GeoStandardPath) and LFV-741 (GeoDegeneratePath); the path values below are NORMAL, which those do not touch.
4. Not affected. GeoPoint.arcDistance divides by both magnitudes (GeoPoint.java:230-244); LINEAR is a plain difference; the plane-based NORMAL distances (Plane.java:645-701) do not use this method. On SPHERE |p| = 1 and the result is exact up to rounding.

Each self-distance equals |p| * ||p|^2 - 1| to at least 11 significant digits. On WGS84 |p| = 1 at geocentric latitude 0.6143 rad (35.2 degrees), where the error vanishes. All latitudes are geo3d's geocentric latitudes.

• Who. Application code that calls the geom distance API (Vector.normalDistance, DistanceStyle.NORMAL.computeDistance, GeoCircle and GeoPath computeDistance, computeNearestDistance, computePathCenterDistance, ...) with DistanceStyle.NORMAL or NORMAL_SQUARED on a non-spherical planet model. The caller always picks the style, one of five public constants (DistanceStyle.java:30-46). The API is @lucene.experimental.
• Absolute values and thresholds. A point's distance to itself reads up to 28.4 km (at the poles), and a point 50 m from a circle's centre reads 10.5 km at latitude 0.3 rad. Because the computed value is never below the true one, a check such as "NORMAL distance below 1 km" drops points it should keep; from a centre at latitude 0.3 no surface point is closer than 10.5 km, so it is false for every point, the centre included. The excess shrinks with distance (551 m at 100 km), so near distances and small thresholds are hit hardest.
• Comparisons across receivers. computePathCenterDistance is "meant to allow comparisons between different paths to find the one that goes closest to a point" (GeoPath.java:60-61); on a WGS84 meridian path a point on the centre line ranks behind one 1.63 km off it. Distances measured from centres at different latitudes carry different floors.
• What stays right. Around one fixed circle centre, a verifier's 40-step ray control kept the ARC order (near the centre the floor is almost constant); order around one centre in general was not proved. SPHERE, ARC and LINEAR are correct up to rounding (ARC(p, p) reads up to 2.1e-8 from acos rounding); shape membership, Lucene's queries and sorts (Geo3DPointDistanceComparator and Geo3DPointOutsideDistanceComparator use ARC) and the index are unaffected.
• How likely. Plausible but rare: practically every opted-in call is wrong (all but those with |p| = 1 or orthogonal vectors), yet no Lucene code uses these styles, and the recorded searches found no downstream code using DistanceStyle.NORMAL and no report since 5.3.0 (2015). Registry priority low.

The class sits in org.apache.lucene.spatial3d.geom only to be next to TestGeoPoint. It compiles on 9.1 and later (org.apache.lucene.tests.util.LuceneTestCase); on 9.0 and older the base class is org.apache.lucene.util.LuceneTestCase (on 9.0.0 the test as written does not compile; on 9.1.0 it fails as the bug), and in 5.3 to 5.5 the package was org.apache.lucene.geo3d.

It shows the self-distance only. The circle, path, CLARKE_1866 and per-latitude values are the executed values in the table above, from a separate value probe, not extra tests.

Regression test: this same test. Its comment line is 101 columns, which google-java-format rejects (./gradlew tidy wraps it).

cp $REPRO/TestNormalDistanceOnEllipsoid.java $LUCENE/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/
cd $LUCENE && ./gradlew :lucene:spatial3d:test --tests org.apache.lucene.spatial3d.geom.TestNormalDistanceOnEllipsoid -Ptests.seed=DEADBEEF
# assertions off fails the same way: add -Ptests.asserts=false (the recorded runs added --offline --no-daemon)

# gradle, scratch worktree of 645d4fe14a20, seed DEADBEEF, assertions on
TestNormalDistanceOnEllipsoid > testNormalDistanceOfAPointToItself FAILED
    java.lang.AssertionError: expected:<0.0> but was:<0.0016507737062451833>
        at __randomizedtesting.SeedInfo.seed([DEADBEEF:A8ED8F369E317E22]:0)
        at junit@4.13.2/org.junit.Assert.fail(Assert.java:89)
        at junit@4.13.2/org.junit.Assert.failNotEquals(Assert.java:835)
        at junit@4.13.2/org.junit.Assert.assertEquals(Assert.java:555)
        at junit@4.13.2/org.junit.Assert.assertEquals(Assert.java:685)
        at org.apache.lucene.spatial3d.geom.TestNormalDistanceOnEllipsoid.testNormalDistanceOfAPointToItself(TestNormalDistanceOnEllipsoid.java:26)

# assertions off (-Ptests.asserts=false): the same message
# standalone runner, seed DEADBEEF: 1 of 1 failed, same message, on main 645d4fe14a20 (-ea and -da), 10.6.0, 9.12.3
# with the recommended fix: the test passes (assertions on and off)

• The contract. Vector.normalDistance: "Compute the normal (perpendicular) distance to a vector described by a vector taken from the origin" (Vector.java:450-451, :461-462). Vector is "A 3d vector in space, not necessarily going through the origin" (:20), and GeoPoint "a point on the surface of a sphere or ellipsoid" (GeoPoint.java:24). No javadoc restricts NORMAL to SPHERE or to unit vectors.
• The author's definition. When the styles were added (LUCENE-6578, #7636, 2015, after WGS84 existed), the geo3d author Karl Wright (@DaddyWri) listed "NormalDistance (perpendicular to a plane/line through the center of the earth)" and said the squared styles "require no square roots, so if you need a fast distance estimate, that's what you'd use": our reading is that the thread ties "estimate" to skipping a square root and states no ellipsoid restriction for NORMAL.
• The tests. TestGeoCircle.testCircleDistance asserts a NORMAL distance of 0.0 at the centre and 0.049979 (sin 0.05) at arc 0.05 (TestGeoCircle.java:23-39, asserts :33 and :38), on SPHERE only; every upstream NORMAL assertion uses SPHERE, which is why no test caught it. The same centre assertion on WGS84 gives 0.002241464883617787.
• History. The formula is from the sphere-only first geo3d code (LUCENE-6196, bde056c7632, 2015-04-13). The WGS84 change (LUCENE-6487, 8c3e6b0ed32, 2015-06-02) made GeoPoint.arcDistance divide by both magnitudes but did not change Vector.normalDistance, although it edited the GeoCircle methods that return it (see the devil's advocate below); its author wrote that on the ellipsoid "distance measures continue to be computed as arc distances, normal distances, or direct (linear) distances".

The strongest case for "intended" (rated moderate by a cross-model review, weak to moderate by verifiers, weak by the other review): the geom API is @lucene.experimental; GeoDistance says its values "are meant to provide both actual distance values, as well as distance estimates that can be computed more cheaply" (GeoDistance.java:20-22), and the normalDistance javadoc promises only "Monotonically increasing for arc distances up to PI/2"; NORMAL was born on the unit sphere, and parameter docs still say "the point's unit x coordinate" (GeoDistance.java:51); and on WGS84 users should use ARC, since Geo3DPoint.newDistanceQuery already warns that the query "may have up to 0.5% error" (Geo3DPoint.java:112-113). The strongest contrary history evidence, in the refute-lens review's words: the WGS84 commit itself edited GeoCircle.computeNormalDistance and computeSquaredNormalDistance. It removed cutoffNormalDistance, replaced the normalDistance > cutoffNormalDistance checks with if (!isWithin(point)), kept return this.center.normalDistance(point);, and commented in the constructor "In an ellipsoidal world, cutoff distances make no sense, unfortunately. Only membership can be used to make in/out determination."

Why it fails: experimental status is about API stability, not permission to return wrong values. An estimate that gives 10.5 km for a point and itself estimates nothing: the result never drops below that floor, however close the points are. The javadoc's "monotonically increasing" note is naturally read for one receiver, and the path example does not contradict it; it shows that comparisons across receivers fail, which computePathCenterDistance's javadoc invites. The style was defined geometrically after WGS84 existed, and the same WGS84 change made arc distances magnitude-aware. The 0.5% note is about query shapes, not returned distances, and at latitude 0.3 the resulting error exceeds 0.5% of every true distance below about 105 km (the 10.5 km floor itself exceeds 0.5% of distances up to about 2,100 km). The 2015 comment is about the circle's membership cutoff, which the change replaced with isWithin and the circle plane; it says nothing about the values returned, and LUCENE-6578 (committed 2015-07-01) then defined NORMAL as a perpendicular distance. The independent reproducer, the cross-model reviews and the devil's-advocate verifier said real; the adjudication confirmed with high confidence.

The recommended fix was validated through gradle in a scratch worktree of 645d4fe14a20; the hunt validated the same division without the guard. The other two options were not built.

Project onto the receiver's direction: t = (this . v) / |this|^2. The guard keeps today's result for a zero-length receiver (t = 0, so the result is |v|^2, the squared distance to the origin) instead of NaN; surface GeoPoints are never zero length.

public double normalDistanceSquared(final Vector v) {
  // project onto the line through the origin and this vector, which need not be unit length
  final double magnitudeSquared = x * x + y * y + z * z;
  final double t = magnitudeSquared == 0.0 ? 0.0 : dotProduct(v) / magnitudeSquared;
  double deltaX = this.x * t - v.x;
  ...
// the (x, y, z) overload: the same, with this.x * this.x + ... and dotProduct(x, y, z)

Trade-offs

• Cost. Three multiplications, two additions, a comparison and one division per call; no square root, so the squared style stays root-free. Not benchmarked.
• Changed values. On SPHERE |p|^2 is 1 to within rounding, so values change at most in the last bits; on an ellipsoid NORMAL and NORMAL_SQUARED point distances change to the correct value (except where |p| = 1 or the vectors are orthogonal). It is a behaviour change of a public method: callers who tuned thresholds to the old numbers, or who call Vector.normalDistance on their own non-unit vectors, will see different results; a compatibility point, not evidence against the bug.
• Reach. Only NormalDistance and NormalSquaredDistance call this method in Lucene, so nothing else changes.

What the checks validated

• Gradle, scratch worktree of 645d4fe14a20, OpenJDK 25.0.3, seed DEADBEEF (2026-10-31). Stock: the repro fails with assertions on and off. With the fix: it passes with assertions on and off.
• :lucene:spatial3d:test, the whole module with the repro in place: 20 suites, 182 tests, 0 failures, 3 skipped (@Nightly tests, not enabled by default).
• :lucene:spatial3d checkGoogleJavaFormat, ecjLintMain, ecjLintTest, forbiddenApisMain, forbiddenApisTest: pass. The fix is format-clean; the repro's comment line is 101 columns and fails the format check, so it was wrapped for this run.
• The value probe gives the correct value for every row of the table above: every self-distance 0.0, the 50 m point 7.854520282177453E-6, the path A 0.0 < B 0.000256.

t = dotProduct(v) / (x * x + y * y + z * z) in both overloads, the finder's two-line sketch.

Trade-offs

Returns NaN for a zero-length receiver. Passed every hunt repro, 55 upstream tests in 7 spatial3d classes and the minimal repro on 645d4fe14a20.

GeoPoint caches its magnitude (GeoPoint.java:34-39, :281-288); an override could divide by magnitude() * magnitude().

Trade-offs

Fixes every style call, whose receivers are GeoPoints, but leaves Vector.normalDistance wrong for plain non-unit vectors, whose javadoc makes the same promise.

State that NORMAL and NORMAL_SQUARED are exact only on SPHERE, or for unit vectors.

Trade-offs

Leaves a 10 km self-distance on WGS84 and contradicts "perpendicular distance" in the javadoc.

• Translation validation. The hunt unit for geo3d distance functions ran the real spatial3d jar against an independent reference that computes the documented perpendicular distance, |q - p (p . q) / |p|^2|, with 50 significant digits, over circles, exact circles and paths on SPHERE and WGS84.
• What it found. On WGS84 all 4,050 inside NORMAL and NORMAL_SQUARED distances per circle type deviated, and more than 25,000 path nearest distances; SPHERE had none. With the fix shadowed, the WGS84 NORMAL deviations disappear except those of other, separately tracked causes.
• Root cause and test. Reading Vector.normalDistanceSquared showed the missing division; a public-API JUnit on the stock main jars failed with assertions on and off, and on 10.5.2 and 9.12.3.
• The search was bounded (a fixed grid of shapes and points), not a proof, so the label is semi-formal.