diff --git a/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistance.java b/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistance.java
index 8f9295bef0..d1f82a97c0 100644
--- a/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistance.java
+++ b/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistance.java
@@ -31,9 +31,16 @@ public double compute(double[] a, double[] b) {
MathArrays.checkEqualLength(a, b);
double sum = 0;
for (int i = 0; i < a.length; i++) {
+ final double absA = JdkMath.abs(a[i]);
+ final double absB = JdkMath.abs(b[i]);
final double num = JdkMath.abs(a[i] - b[i]);
- final double denom = JdkMath.abs(a[i]) + JdkMath.abs(b[i]);
- sum += num == 0.0 && denom == 0.0 ? 0.0 : num / denom;
+ final double denom = absA + absB;
+ if (Double.isInfinite(denom) && Double.isFinite(a[i]) && Double.isFinite(b[i])) {
+ // Rescale large finite coordinates to avoid intermediate overflow.
+ sum += JdkMath.abs(a[i] * 0.5 - b[i] * 0.5) / (absA * 0.5 + absB * 0.5);
+ } else {
+ sum += denom == 0.0 ? 0.0 : num / denom;
+ }
}
return sum;
}
diff --git a/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistanceTest.java b/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistanceTest.java
index 85384ba116..8b38df937b 100644
--- a/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistanceTest.java
+++ b/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/ml/distance/CanberraDistanceTest.java
@@ -16,6 +16,9 @@
*/
package org.apache.commons.math4.legacy.ml.distance;
+import java.math.BigDecimal;
+import java.math.MathContext;
+
import org.junit.Assert;
import org.junit.Test;
@@ -45,4 +48,81 @@ public void test() {
Assert.assertEquals(expected, distance.compute(a, b), 0d);
Assert.assertEquals(expected, distance.compute(b, a), 0d);
}
+
+ @Test
+ public void testExtremeOppositeSigns() {
+ // Both |MAX - (-MAX)| and |MAX| + |-MAX| overflow to infinity: the true ratio is 1.
+ final double max = Double.MAX_VALUE;
+ Assert.assertEquals(1, distance.compute(new double[] {max}, new double[] {-max}), 0d);
+ Assert.assertEquals(1, distance.compute(new double[] {-max}, new double[] {max}), 0d);
+ Assert.assertEquals(1, distance.compute(new double[] {max}, new double[] {-max / 2}), 0d);
+ Assert.assertEquals(1, distance.compute(new double[] {-max / 2}, new double[] {max}), 0d);
+ }
+
+ @Test
+ public void testExtremeSameSign() {
+ // |MAX - MAX/2| / (|MAX| + |MAX/2|) == 1/3; only the denominator overflows.
+ final double max = Double.MAX_VALUE;
+ // 1/3 is not representable, so allow one ulp instead of an exact match.
+ final double expected = 1.0 / 3.0;
+ final double tol = Math.ulp(expected);
+ Assert.assertEquals(expected, distance.compute(new double[] {max}, new double[] {max / 2}), tol);
+ Assert.assertEquals(expected, distance.compute(new double[] {max / 2}, new double[] {max}), tol);
+ Assert.assertEquals(expected, distance.compute(new double[] {-max}, new double[] {-max / 2}), tol);
+ Assert.assertEquals(expected, distance.compute(new double[] {-max / 2}, new double[] {-max}), tol);
+ // Identical extreme values have distance 0.
+ Assert.assertEquals(0, distance.compute(new double[] {max}, new double[] {max}), 0d);
+ Assert.assertEquals(0, distance.compute(new double[] {-max}, new double[] {-max}), 0d);
+ }
+
+ @Test
+ public void testExtremeVector() {
+ final double max = Double.MAX_VALUE;
+ final double[] a = { max, max / 2, max };
+ final double[] b = { -max, max, max };
+ // 1 + 1/3 + 0
+ Assert.assertEquals(4.0 / 3.0, distance.compute(a, b), 1e-15);
+ Assert.assertEquals(4.0 / 3.0, distance.compute(b, a), 1e-15);
+ }
+
+ @Test
+ public void testExtremeFiniteValuesAgainstBigDecimalReference() {
+ final double max = Double.MAX_VALUE;
+ final double[] a = {max, max, 0.75 * max, max, -max, max};
+ final double[] b = {max / 2, Math.nextDown(max), Math.nextDown(0.75 * max), -max, -max / 2, max / 16};
+ for (int i = 0; i < a.length; i++) {
+ final double expected = referenceCanberra(a[i], b[i]);
+ assertWithinOneUlp(expected, distance.compute(new double[] {a[i]}, new double[] {b[i]}));
+ assertWithinOneUlp(expected, distance.compute(new double[] {b[i]}, new double[] {a[i]}));
+ }
+ }
+
+ /**
+ * Computes a high-precision reference from the exact decimal representations of the input doubles.
+ */
+ private static double referenceCanberra(double a, double b) {
+ final BigDecimal bdA = new BigDecimal(a);
+ final BigDecimal bdB = new BigDecimal(b);
+ final BigDecimal num = bdA.subtract(bdB).abs();
+ final BigDecimal den = bdA.abs().add(bdB.abs());
+ return num.divide(den, MathContext.DECIMAL128).doubleValue();
+ }
+
+ private static void assertWithinOneUlp(double expected, double actual) {
+ Assert.assertEquals("expected " + expected + " but was " + actual,
+ expected, actual, Math.ulp(expected));
+ }
+
+ @Test
+ public void testSpecialValuesPreserved() {
+ final double min = Double.MIN_VALUE;
+ // Subnormal ratios are unaffected by the extreme-value handling.
+ Assert.assertEquals(1.0 / 3.0, distance.compute(new double[] {min}, new double[] {2 * min}), 0d);
+ Assert.assertEquals(1, distance.compute(new double[] {min}, new double[] {-min}), 0d);
+ // Signed zero: the zero/zero convention contributes 0.
+ Assert.assertEquals(0, distance.compute(new double[] {0.0}, new double[] {-0.0}), 0d);
+ // Infinite and NaN inputs keep their historical NaN result.
+ Assert.assertTrue(Double.isNaN(distance.compute(new double[] {Double.POSITIVE_INFINITY}, new double[] {1})));
+ Assert.assertTrue(Double.isNaN(distance.compute(new double[] {Double.NaN}, new double[] {0})));
+ }
}
diff --git a/src/changes/changes.xml b/src/changes/changes.xml
index 077104d34d..b111a20cc7 100644
--- a/src/changes/changes.xml
+++ b/src/changes/changes.xml
@@ -101,6 +101,9 @@ Caveat:
Instances not cannot be constructed that would silently overwrite entries
due to key duplication from interger overflow.
+
+ "CanberraDistance": Avoid overflow for large finite coordinates.
+
"SynchronizedDescriptiveStatistics": Synchronize all mutators
of the internal state. Updates previously unsynchronized methods