You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by er...@apache.org on 2021/06/12 22:50:20 UTC
[commons-numbers] branch master updated (9e0dc07 -> 19899e7)
This is an automated email from the ASF dual-hosted git repository.
erans pushed a change to branch master
in repository https://gitbox.apache.org/repos/asf/commons-numbers.git.
from 9e0dc07 Java 9 required for perftest
new a77a8eb NUMBERS-156: Norm API.
new 19899e7 Unit test.
The 2 revisions listed above as "new" are entirely new to this
repository and will be described in separate emails. The revisions
listed as "add" were already present in the repository and have only
been added to this reference.
Summary of changes:
.../org/apache/commons/numbers/angle/CosAngle.java | 4 +-
commons-numbers-arrays/pom.xml | 10 +
.../arrays/MultidimensionalCounterTest.java | 22 +
.../org/apache/commons/numbers/core/Norms.java | 448 ---------------
.../org/apache/commons/numbers/core/NormsTest.java | 625 ---------------------
5 files changed, 34 insertions(+), 1075 deletions(-)
delete mode 100644 commons-numbers-core/src/main/java/org/apache/commons/numbers/core/Norms.java
delete mode 100644 commons-numbers-core/src/test/java/org/apache/commons/numbers/core/NormsTest.java
[commons-numbers] 02/02: Unit test.
Posted by er...@apache.org.
This is an automated email from the ASF dual-hosted git repository.
erans pushed a commit to branch master
in repository https://gitbox.apache.org/repos/asf/commons-numbers.git
commit 19899e7124dcd182fb3ebe31f17405c15be6cd19
Author: Gilles Sadowski <gi...@gmail.com>
AuthorDate: Sun Jun 13 00:20:37 2021 +0200
Unit test.
---
commons-numbers-arrays/pom.xml | 10 ++++++++++
.../arrays/MultidimensionalCounterTest.java | 22 ++++++++++++++++++++++
2 files changed, 32 insertions(+)
diff --git a/commons-numbers-arrays/pom.xml b/commons-numbers-arrays/pom.xml
index 894b93c..756b76f 100644
--- a/commons-numbers-arrays/pom.xml
+++ b/commons-numbers-arrays/pom.xml
@@ -42,4 +42,14 @@
<numbers.parent.dir>${basedir}/..</numbers.parent.dir>
</properties>
+ <dependencies>
+
+ <dependency>
+ <groupId>org.apache.commons</groupId>
+ <artifactId>commons-math3</artifactId>
+ <scope>test</scope>
+ </dependency>
+
+ </dependencies>
+
</project>
diff --git a/commons-numbers-arrays/src/test/java/org/apache/commons/numbers/arrays/MultidimensionalCounterTest.java b/commons-numbers-arrays/src/test/java/org/apache/commons/numbers/arrays/MultidimensionalCounterTest.java
index 7307eb5..62afd23 100644
--- a/commons-numbers-arrays/src/test/java/org/apache/commons/numbers/arrays/MultidimensionalCounterTest.java
+++ b/commons-numbers-arrays/src/test/java/org/apache/commons/numbers/arrays/MultidimensionalCounterTest.java
@@ -117,4 +117,26 @@ class MultidimensionalCounterTest {
final MultidimensionalCounter c = MultidimensionalCounter.of(sizes);
Assertions.assertEquals(Arrays.toString(sizes), c.toString());
}
+
+ // Illustrates how to recover the iterator functionality that existed
+ // in Commons Math (v3.6.1) but was not ported to "Commons Numbers".
+ @Test
+ void testCommonsMathIterator() {
+ final int[] sizes = new int[] {3, 2, 5};
+ final org.apache.commons.math3.util.MultidimensionalCounter.Iterator cmIter =
+ new org.apache.commons.math3.util.MultidimensionalCounter(sizes).iterator();
+
+ final MultidimensionalCounter counter = MultidimensionalCounter.of(sizes);
+
+ Assertions.assertTrue(cmIter.hasNext());
+ Assertions.assertTrue(counter.getSize() > 0);
+
+ for (int i = 0; i < counter.getSize(); i++) {
+ cmIter.next();
+ Assertions.assertTrue(Arrays.equals(cmIter.getCounts(),
+ counter.toMulti(i)));
+ }
+
+ Assertions.assertFalse(cmIter.hasNext());
+ }
}
[commons-numbers] 01/02: NUMBERS-156: Norm API.
Posted by er...@apache.org.
This is an automated email from the ASF dual-hosted git repository.
erans pushed a commit to branch master
in repository https://gitbox.apache.org/repos/asf/commons-numbers.git
commit a77a8eb22b9527c50e50eac106a6821676a63fba
Author: Gilles Sadowski <gi...@gmail.com>
AuthorDate: Sun Jun 13 00:16:58 2021 +0200
NUMBERS-156: Norm API.
---
.../org/apache/commons/numbers/angle/CosAngle.java | 4 +-
.../org/apache/commons/numbers/core/Norms.java | 448 ---------------
.../org/apache/commons/numbers/core/NormsTest.java | 625 ---------------------
3 files changed, 2 insertions(+), 1075 deletions(-)
diff --git a/commons-numbers-angle/src/main/java/org/apache/commons/numbers/angle/CosAngle.java b/commons-numbers-angle/src/main/java/org/apache/commons/numbers/angle/CosAngle.java
index cfa569d..d9d7dd4 100644
--- a/commons-numbers-angle/src/main/java/org/apache/commons/numbers/angle/CosAngle.java
+++ b/commons-numbers-angle/src/main/java/org/apache/commons/numbers/angle/CosAngle.java
@@ -17,7 +17,7 @@
package org.apache.commons.numbers.angle;
import org.apache.commons.numbers.core.LinearCombination;
-import org.apache.commons.numbers.core.Norms;
+import org.apache.commons.numbers.core.Norm;
/**
* Computes the cosine of the angle between two vectors.
@@ -39,7 +39,7 @@ public final class CosAngle {
*/
public static double value(double[] v1,
double[] v2) {
- return LinearCombination.value(v1, v2) / Norms.euclidean(v1) / Norms.euclidean(v2);
+ return LinearCombination.value(v1, v2) / Norm.L2.of(v1) / Norm.L2.of(v2);
}
}
diff --git a/commons-numbers-core/src/main/java/org/apache/commons/numbers/core/Norms.java b/commons-numbers-core/src/main/java/org/apache/commons/numbers/core/Norms.java
deleted file mode 100644
index 2a357ca..0000000
--- a/commons-numbers-core/src/main/java/org/apache/commons/numbers/core/Norms.java
+++ /dev/null
@@ -1,448 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.numbers.core;
-
-/** Class providing methods to compute various norm values.
- *
- * <p>This class uses a variety of techniques to increase numerical accuracy
- * and reduce errors. A primary source for the included algorithms is the
- * 2005 paper <a href="https://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.2.1547">
- * Accurate Sum and Dot Product</a> by Takeshi Ogita, Siegfried M. Rump,
- * and Shin'ichi Oishi published in <em>SIAM J. Sci. Comput</em>.
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)">Norm</a>
- */
-public final class Norms {
-
- /** Threshold for scaling small numbers. This value is chosen such that doubles
- * set to this value can be squared without underflow. Values less than this must
- * be scaled up.
- */
- private static final double SMALL_THRESH = 0x1.0p-511;
-
- /** Threshold for scaling large numbers. This value is chosen such that 2^31 doubles
- * set to this value can be squared and added without overflow. Values greater than
- * this must be scaled down.
- */
- private static final double LARGE_THRESH = 0x1.0p+496;
-
- /** Threshold for scaling up a single value by {@link #SCALE_UP} without risking overflow
- * when the value is squared.
- */
- private static final double SAFE_SCALE_UP_THRESH = 0x1.0p-100;
-
- /** Value used to scale down large numbers. */
- private static final double SCALE_DOWN = 0x1.0p-600;
-
- /** Value used to scale up small numbers. */
- private static final double SCALE_UP = 0x1.0p+600;
-
- /** Threshold for the difference between the exponents of two Euclidean 2D input values
- * where the larger value dominates the calculation.
- */
- private static final int EXP_DIFF_THRESHOLD_2D = 54;
-
- /** Utility class; no instantiation. */
- private Norms() {}
-
- /** Compute the Manhattan norm (also known as the Taxicab norm or L1 norm) of the arguments.
- * The result is equal to \(|x| + |y|\), i.e., the sum of the absolute values of the arguments.
- *
- * <p>Special cases:
- * <ul>
- * <li>If either value is NaN, then the result is NaN.</li>
- * <li>If either value is infinite and the other value is not NaN, then the result is positive infinity.</li>
- * </ul>
- * @param x first input value
- * @param y second input value
- * @return Manhattan norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Taxicab_norm_or_Manhattan_norm">Manhattan norm</a>
- */
- public static double manhattan(final double x, final double y) {
- return Math.abs(x) + Math.abs(y);
- }
-
- /** Compute the Manhattan norm (also known as the Taxicab norm or L1 norm) of the arguments.
- * The result is equal to \(|x| + |y| + |z|\), i.e., the sum of the absolute values of the arguments.
- *
- * <p>Special cases:
- * <ul>
- * <li>If any value is NaN, then the result is NaN.</li>
- * <li>If any value is infinite and no value is NaN, then the result is positive infinity.</li>
- * </ul>
- * @param x first input value
- * @param y second input value
- * @param z third input value
- * @return Manhattan norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Taxicab_norm_or_Manhattan_norm">Manhattan norm</a>
- */
- public static double manhattan(final double x, final double y, final double z) {
- return Summation.value(
- Math.abs(x),
- Math.abs(y),
- Math.abs(z));
- }
-
- /** Compute the Manhattan norm (also known as the Taxicab norm or L1 norm) of the given values.
- * The result is equal to \(|v_0| + ... + |v_i|\), i.e., the sum of the absolute values of the input elements.
- *
- * <p>Special cases:
- * <ul>
- * <li>If any value is NaN, then the result is NaN.</li>
- * <li>If any value is infinite and no value is NaN, then the result is positive infinity.</li>
- * <li>If the array is empty, then the result is 0.</li>
- * </ul>
- * @param v input values
- * @return Manhattan norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Taxicab_norm_or_Manhattan_norm">Manhattan norm</a>
- */
- public static double manhattan(final double[] v) {
- double sum = 0d;
- double comp = 0d;
-
- for (int i = 0; i < v.length; ++i) {
- final double x = Math.abs(v[i]);
- final double sx = sum + x;
- comp += ExtendedPrecision.twoSumLow(sum, x, sx);
- sum = sx;
- }
-
- return Summation.summationResult(sum, comp);
- }
-
- /** Compute the Euclidean norm (also known as the L2 norm) of the arguments. The result is equal to
- * \(\sqrt{x^2 + y^2}\). This method correctly handles the possibility of overflow or underflow
- * during the computation.
- *
- * <p>Special cases:
- * <ul>
- * <li>If either value is NaN, then the result is NaN.</li>
- * <li>If either value is infinite and the other value is not NaN, then the result is positive infinity.</li>
- * </ul>
- *
- * <p><strong>Comparison with Math.hypot()</strong>
- * <p>While not guaranteed to return the same result, this method does provide similar error bounds to
- * the JDK's Math.hypot() method and may run faster on some JVMs.
- * @param x first input
- * @param y second input
- * @return Euclidean norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Euclidean_norm">Euclidean norm</a>
- */
- public static double euclidean(final double x, final double y) {
- final double xabs = Math.abs(x);
- final double yabs = Math.abs(y);
-
- final double max;
- final double min;
- // the compare method considers NaN greater than other values, meaning that our
- // check for if the max is finite later on will detect NaNs correctly
- if (Double.compare(xabs, yabs) > 0) {
- max = xabs;
- min = yabs;
- } else {
- max = yabs;
- min = xabs;
- }
-
- // if the max is not finite, then one of the inputs must not have
- // been finite
- if (!Double.isFinite(max)) {
- // let the standard multiply operation determine whether to return NaN or infinite
- return xabs * yabs;
- } else if (Math.getExponent(max) - Math.getExponent(min) > EXP_DIFF_THRESHOLD_2D) {
- // value is completely dominated by max; just return max
- return max;
- }
-
- // compute the scale and rescale values
- final double scale;
- final double rescale;
- if (max > LARGE_THRESH) {
- scale = SCALE_DOWN;
- rescale = SCALE_UP;
- } else if (max < SAFE_SCALE_UP_THRESH) {
- scale = SCALE_UP;
- rescale = SCALE_DOWN;
- } else {
- scale = 1d;
- rescale = 1d;
- }
-
- double sum = 0d;
- double comp = 0d;
-
- // add scaled x
- double sx = xabs * scale;
- final double px = sx * sx;
- comp += ExtendedPrecision.squareLowUnscaled(sx, px);
- final double sumPx = sum + px;
- comp += ExtendedPrecision.twoSumLow(sum, px, sumPx);
- sum = sumPx;
-
- // add scaled y
- double sy = yabs * scale;
- final double py = sy * sy;
- comp += ExtendedPrecision.squareLowUnscaled(sy, py);
- final double sumPy = sum + py;
- comp += ExtendedPrecision.twoSumLow(sum, py, sumPy);
- sum = sumPy;
-
- return Math.sqrt(sum + comp) * rescale;
- }
-
- /** Compute the Euclidean norm (also known as the L2 norm) of the arguments. The result is equal to
- * \(\sqrt{x^2 + y^2 + z^2}\). This method correctly handles the possibility of overflow or underflow
- * during the computation.
- *
- * <p>Special cases:
- * <ul>
- * <li>If any value is NaN, then the result is NaN.</li>
- * <li>If any value is infinite and no value is NaN, then the result is positive infinity.</li>
- * </ul>
- * @param x first input
- * @param y second input
- * @param z third input
- * @return Euclidean norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Euclidean_norm">Euclidean norm</a>
- */
- public static double euclidean(final double x, final double y, final double z) {
- final double xabs = Math.abs(x);
- final double yabs = Math.abs(y);
- final double zabs = Math.abs(z);
-
- final double max = Math.max(Math.max(xabs, yabs), zabs);
-
- // if the max is not finite, then one of the inputs must not have
- // been finite
- if (!Double.isFinite(max)) {
- // let the standard multiply operation determine whether to return NaN or infinite
- return xabs * yabs * zabs;
- }
-
- // compute the scale and rescale values
- final double scale;
- final double rescale;
- if (max > LARGE_THRESH) {
- scale = SCALE_DOWN;
- rescale = SCALE_UP;
- } else if (max < SAFE_SCALE_UP_THRESH) {
- scale = SCALE_UP;
- rescale = SCALE_DOWN;
- } else {
- scale = 1d;
- rescale = 1d;
- }
-
- double sum = 0d;
- double comp = 0d;
-
- // add scaled x
- double sx = xabs * scale;
- final double px = sx * sx;
- comp += ExtendedPrecision.squareLowUnscaled(sx, px);
- final double sumPx = sum + px;
- comp += ExtendedPrecision.twoSumLow(sum, px, sumPx);
- sum = sumPx;
-
- // add scaled y
- double sy = yabs * scale;
- final double py = sy * sy;
- comp += ExtendedPrecision.squareLowUnscaled(sy, py);
- final double sumPy = sum + py;
- comp += ExtendedPrecision.twoSumLow(sum, py, sumPy);
- sum = sumPy;
-
- // add scaled z
- final double sz = zabs * scale;
- final double pz = sz * sz;
- comp += ExtendedPrecision.squareLowUnscaled(sz, pz);
- final double sumPz = sum + pz;
- comp += ExtendedPrecision.twoSumLow(sum, pz, sumPz);
- sum = sumPz;
-
- return Math.sqrt(sum + comp) * rescale;
- }
-
- /** Compute the Euclidean norm (also known as the L2 norm) of the given values. The result is equal to
- * \(\sqrt{v_0^2 + ... + v_{n-1}^2}\). This method correctly handles the possibility of overflow or underflow
- * during the computation.
- *
- * <p>Special cases:
- * <ul>
- * <li>If any value is NaN, then the result is NaN.</li>
- * <li>If any value is infinite and no value is NaN, then the result is positive infinity.</li>
- * <li>If the array is empty, then the result is 0.</li>
- * </ul>
- * @param v input values
- * @return Euclidean norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Euclidean_norm">Euclidean norm</a>
- */
- public static double euclidean(final double[] v) {
- // sum of big, normal and small numbers
- double s1 = 0;
- double s2 = 0;
- double s3 = 0;
-
- // sum compensation values
- double c1 = 0;
- double c2 = 0;
- double c3 = 0;
-
- for (int i = 0; i < v.length; ++i) {
- final double x = Math.abs(v[i]);
- if (!Double.isFinite(x)) {
- // not finite; determine whether to return NaN or positive infinity
- return euclideanNormSpecial(v, i);
- } else if (x > LARGE_THRESH) {
- // scale down
- final double sx = x * SCALE_DOWN;
-
- // compute the product and product compensation
- final double p = sx * sx;
- final double cp = ExtendedPrecision.squareLowUnscaled(sx, p);
-
- // compute the running sum and sum compensation
- final double s = s1 + p;
- final double cs = ExtendedPrecision.twoSumLow(s1, p, s);
-
- // update running totals
- c1 += cp + cs;
- s1 = s;
- } else if (x < SMALL_THRESH) {
- // scale up
- final double sx = x * SCALE_UP;
-
- // compute the product and product compensation
- final double p = sx * sx;
- final double cp = ExtendedPrecision.squareLowUnscaled(sx, p);
-
- // compute the running sum and sum compensation
- final double s = s3 + p;
- final double cs = ExtendedPrecision.twoSumLow(s3, p, s);
-
- // update running totals
- c3 += cp + cs;
- s3 = s;
- } else {
- // no scaling
- // compute the product and product compensation
- final double p = x * x;
- final double cp = ExtendedPrecision.squareLowUnscaled(x, p);
-
- // compute the running sum and sum compensation
- final double s = s2 + p;
- final double cs = ExtendedPrecision.twoSumLow(s2, p, s);
-
- // update running totals
- c2 += cp + cs;
- s2 = s;
- }
- }
-
- // The highest sum is the significant component. Add the next significant.
- // Note that the "x * SCALE_DOWN * SCALE_DOWN" expressions must be executed
- // in the order given. If the two scale factors are multiplied together first,
- // they will underflow to zero.
- if (s1 != 0) {
- // add s1, s2, c1, c2
- final double s2Adj = s2 * SCALE_DOWN * SCALE_DOWN;
- final double sum = s1 + s2Adj;
- final double comp = ExtendedPrecision.twoSumLow(s1, s2Adj, sum) + c1 + (c2 * SCALE_DOWN * SCALE_DOWN);
- return Math.sqrt(sum + comp) * SCALE_UP;
- } else if (s2 != 0) {
- // add s2, s3, c2, c3
- final double s3Adj = s3 * SCALE_DOWN * SCALE_DOWN;
- final double sum = s2 + s3Adj;
- final double comp = ExtendedPrecision.twoSumLow(s2, s3Adj, sum) + c2 + (c3 * SCALE_DOWN * SCALE_DOWN);
- return Math.sqrt(sum + comp);
- }
- // add s3, c3
- return Math.sqrt(s3 + c3) * SCALE_DOWN;
- }
-
- /** Return a Euclidean norm value for special cases of non-finite input.
- * @param v input vector
- * @param start index to start examining the input vector from
- * @return Euclidean norm special value
- */
- private static double euclideanNormSpecial(final double[] v, final int start) {
- for (int i = start; i < v.length; ++i) {
- if (Double.isNaN(v[i])) {
- return Double.NaN;
- }
- }
- return Double.POSITIVE_INFINITY;
- }
-
- /** Compute the maximum norm (also known as the infinity norm or L<sub>inf</sub> norm) of the arguments.
- * The result is equal to \(\max{(|x|, |y|)}\), i.e., the maximum of the absolute values of the arguments.
- *
- * <p>Special cases:
- * <ul>
- * <li>If either value is NaN, then the result is NaN.</li>
- * <li>If either value is infinite and the other value is not NaN, then the result is positive infinity.</li>
- * </ul>
- * @param x first input
- * @param y second input
- * @return maximum norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Maximum_norm_(special_case_of:_infinity_norm,_uniform_norm,_or_supremum_norm)">Maximum norm</a>
- */
- public static double maximum(final double x, final double y) {
- return Math.max(Math.abs(x), Math.abs(y));
- }
-
- /** Compute the maximum norm (also known as the infinity norm or L<sub>inf</sub> norm) of the arguments.
- * The result is equal to \(\max{(|x|, |y|, |z|)}\), i.e., the maximum of the absolute values of the arguments.
- *
- * <p>Special cases:
- * <ul>
- * <li>If any value is NaN, then the result is NaN.</li>
- * <li>If any value is infinite and no value is NaN, then the result is positive infinity.</li>
- * </ul>
- * @param x first input
- * @param y second input
- * @param z third input
- * @return maximum norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Maximum_norm_(special_case_of:_infinity_norm,_uniform_norm,_or_supremum_norm)">Maximum norm</a>
- */
- public static double maximum(final double x, final double y, final double z) {
- return Math.max(
- Math.abs(x),
- Math.max(Math.abs(y), Math.abs(z)));
- }
-
- /** Compute the maximum norm (also known as the infinity norm or L<sub>inf</sub> norm) of the given values.
- * The result is equal to \(\max{(|v_0|, \ldots, |v_{n-1}|)}\), i.e., the maximum of the absolute values of the
- * input elements.
- *
- * <p>Special cases:
- * <ul>
- * <li>If any value is NaN, then the result is NaN.</li>
- * <li>If any value is infinite and no value is NaN, then the result is positive infinity.</li>
- * <li>If the array is empty, then the result is 0.</li>
- * </ul>
- * @param v input values
- * @return maximum norm
- * @see <a href="https://en.wikipedia.org/wiki/Norm_(mathematics)#Maximum_norm_(special_case_of:_infinity_norm,_uniform_norm,_or_supremum_norm)">Maximum norm</a>
- */
- public static double maximum(final double[] v) {
- double max = 0d;
- for (int i = 0; i < v.length; ++i) {
- max = Math.max(max, Math.abs(v[i]));
- }
- return max;
- }
-}
diff --git a/commons-numbers-core/src/test/java/org/apache/commons/numbers/core/NormsTest.java b/commons-numbers-core/src/test/java/org/apache/commons/numbers/core/NormsTest.java
deleted file mode 100644
index aa1496d..0000000
--- a/commons-numbers-core/src/test/java/org/apache/commons/numbers/core/NormsTest.java
+++ /dev/null
@@ -1,625 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.numbers.core;
-
-import java.math.BigDecimal;
-import java.math.MathContext;
-import java.util.Arrays;
-import java.util.function.ToDoubleFunction;
-
-import org.apache.commons.rng.UniformRandomProvider;
-import org.apache.commons.rng.simple.RandomSource;
-import org.junit.jupiter.api.Assertions;
-import org.junit.jupiter.api.Test;
-
-class NormsTest {
-
- private static final int SMALL_THRESH_EXP = -511;
-
- private static final int LARGE_THRESH_EXP = +496;
-
- private static final int RAND_VECTOR_CNT = 1_000;
-
- private static final int MAX_ULP_ERR = 1;
-
- private static final double HYPOT_COMPARE_EPS = 1e-2;
-
- private static final BigDecimal BD_MAX_VALUE = new BigDecimal(Double.MAX_VALUE);
- private static final BigDecimal BD_MIN_NORMAL = new BigDecimal(Double.MIN_NORMAL);
-
- /** The scale, used to scale the sqrt of the sum of squares. */
- private static final double SCALE = 0x1.0p200;
-
- /** The scale squared, used to scale the sum of squares. */
- private static final BigDecimal SCALE2 = new BigDecimal(SCALE * SCALE);
-
- @Test
- void testManhattan_2d() {
- // act/assert
- Assertions.assertEquals(0d, Norms.manhattan(0d, -0d));
- Assertions.assertEquals(3d, Norms.manhattan(-1d, 2d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.manhattan(Double.MAX_VALUE, Double.MAX_VALUE));
-
- Assertions.assertEquals(Double.NaN, Norms.manhattan(Double.NaN, 1d));
- Assertions.assertEquals(Double.NaN, Norms.manhattan(1d, Double.NaN));
- Assertions.assertEquals(Double.NaN, Norms.manhattan(Double.POSITIVE_INFINITY, Double.NaN));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(Double.POSITIVE_INFINITY, 0d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(0d, Double.POSITIVE_INFINITY));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY));
- }
-
- @Test
- void testManhattan_3d() {
- // act/assert
- Assertions.assertEquals(0d, Norms.manhattan(0d, -0d, 0d));
- Assertions.assertEquals(6d, Norms.manhattan(-1d, 2d, -3d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.manhattan(Double.MAX_VALUE, Double.MAX_VALUE, 0d));
-
- Assertions.assertEquals(Double.NaN, Norms.manhattan(Double.NaN, -2d, 1d));
- Assertions.assertEquals(Double.NaN, Norms.manhattan(-2d, Double.NaN, 1d));
- Assertions.assertEquals(Double.NaN, Norms.manhattan(-2d, 1d, Double.NaN));
- Assertions.assertEquals(Double.NaN, Norms.manhattan(-2d, Double.POSITIVE_INFINITY, Double.NaN));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(Double.POSITIVE_INFINITY, 2d, -4d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, -4d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY));
- }
-
- @Test
- void testManhattan_array() {
- // act/assert
- Assertions.assertEquals(0d, Norms.manhattan(new double[0]));
- Assertions.assertEquals(0d, Norms.manhattan(new double[] {0d, -0d}));
- Assertions.assertEquals(6d, Norms.manhattan(new double[] {-1d, 2d, -3d}));
- Assertions.assertEquals(10d, Norms.manhattan(new double[] {-1d, 2d, -3d, 4d}));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(new double[] {Double.MAX_VALUE, Double.MAX_VALUE}));
-
- Assertions.assertEquals(Double.NaN, Norms.manhattan(new double[] {-2d, Double.NaN, 1d}));
- Assertions.assertEquals(Double.NaN, Norms.manhattan(new double[] {Double.POSITIVE_INFINITY, Double.NaN, 1d}));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(new double[] {Double.POSITIVE_INFINITY, 0d}));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.manhattan(new double[] {Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY}));
- }
-
- @Test
- void testEuclidean_2d_simple() {
- // act/assert
- Assertions.assertEquals(0d, Norms.euclidean(0d, 0d));
- Assertions.assertEquals(1d, Norms.euclidean(1d, 0d));
- Assertions.assertEquals(1d, Norms.euclidean(0d, 1d));
- Assertions.assertEquals(5d, Norms.euclidean(-3d, 4d));
- Assertions.assertEquals(Double.MIN_VALUE, Norms.euclidean(0d, Double.MIN_VALUE));
- Assertions.assertEquals(Double.MAX_VALUE, Norms.euclidean(Double.MAX_VALUE, 0d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.euclidean(Double.MAX_VALUE, Double.MAX_VALUE));
-
- Assertions.assertEquals(Math.sqrt(2), Norms.euclidean(1d, -1d));
-
- Assertions.assertEquals(Double.NaN, Norms.euclidean(Double.NaN, -2d));
- Assertions.assertEquals(Double.NaN, Norms.euclidean(Double.NaN, Double.POSITIVE_INFINITY));
- Assertions.assertEquals(Double.NaN, Norms.euclidean(-2d, Double.NaN));
- Assertions.assertEquals(Double.NaN,
- Norms.euclidean(Double.NaN, Double.NEGATIVE_INFINITY));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(1d, Double.NEGATIVE_INFINITY));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(Double.POSITIVE_INFINITY, -1d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY));
- }
-
- @Test
- void testEuclidean_2d_scaled() {
- // arrange
- final double[] ones = new double[] {1, 1};
- final double[] multiplesOfTen = new double[] {1, 10};
- final ToDoubleFunction<double[]> fn = v -> Norms.euclidean(v[0], v[1]);
-
- // act/assert
- checkScaledEuclideanNorm(ones, 0, fn);
- checkScaledEuclideanNorm(ones, LARGE_THRESH_EXP, fn);
- checkScaledEuclideanNorm(ones, LARGE_THRESH_EXP + 1, fn);
- checkScaledEuclideanNorm(ones, -100, fn);
- checkScaledEuclideanNorm(ones, -101, fn);
- checkScaledEuclideanNorm(ones, SMALL_THRESH_EXP, fn);
- checkScaledEuclideanNorm(ones, SMALL_THRESH_EXP - 1, fn);
-
-
- checkScaledEuclideanNorm(multiplesOfTen, 0, fn);
- checkScaledEuclideanNorm(multiplesOfTen, -100, fn);
- checkScaledEuclideanNorm(multiplesOfTen, -101, fn);
- checkScaledEuclideanNorm(multiplesOfTen, LARGE_THRESH_EXP - 1, fn);
- checkScaledEuclideanNorm(multiplesOfTen, SMALL_THRESH_EXP, fn);
- }
-
- @Test
- void testEuclidean_2d_dominantValue() {
- // act/assert
- Assertions.assertEquals(Math.PI, Norms.euclidean(-Math.PI, 0x1.0p-55));
- Assertions.assertEquals(Math.PI, Norms.euclidean(0x1.0p-55, -Math.PI));
- }
-
- @Test
- void testEuclidean_2d_random() {
- // arrange
- final UniformRandomProvider rng = RandomSource.create(RandomSource.XO_RO_SHI_RO_1024_PP, 1L);
-
- // act/assert
- checkEuclideanRandom(2, rng, v -> Norms.euclidean(v[0], v[1]));
- }
-
- @Test
- void testEuclidean_2d_vsArray() {
- // arrange
- final double[][] inputs = {
- {-4.074598908124454E-9, 9.897869969944898E-28},
- {1.3472131556526359E-27, -9.064577177323565E9},
- {-3.9219339341360245E149, -7.132522817112096E148},
- {-1.4888098520466735E153, -2.9099184907796666E150},
- {-8.659395144898396E-152, -1.123275532302136E-150},
- {-3.660198254902351E-152, -6.656524053354807E-153}
- };
-
- // act/assert
- for (final double[] input : inputs) {
- Assertions.assertEquals(Norms.euclidean(input), Norms.euclidean(input[0], input[1]),
- () -> "Expected inline method result to equal array result for input " + Arrays.toString(input));
- }
- }
-
- @Test
- void testEuclidean_2d_vsHypot() {
- // arrange
- final int samples = 1000;
- final UniformRandomProvider rng = RandomSource.create(RandomSource.XO_RO_SHI_RO_1024_PP, 3L);
-
- // act/assert
- assertEuclidean2dVersusHypot(-10, +10, samples, rng);
- assertEuclidean2dVersusHypot(0, +20, samples, rng);
- assertEuclidean2dVersusHypot(-20, 0, samples, rng);
- assertEuclidean2dVersusHypot(-20, +20, samples, rng);
- assertEuclidean2dVersusHypot(-100, +100, samples, rng);
- assertEuclidean2dVersusHypot(LARGE_THRESH_EXP - 10, LARGE_THRESH_EXP + 10, samples, rng);
- assertEuclidean2dVersusHypot(SMALL_THRESH_EXP - 10, SMALL_THRESH_EXP + 10, samples, rng);
- assertEuclidean2dVersusHypot(-600, +600, samples, rng);
- }
-
- /** Assert that the Norms euclidean 2D computation produces similar error behavior to Math.hypot().
- * @param minExp minimum exponent for random inputs
- * @param maxExp maximum exponent for random inputs
- * @param samples sample count
- * @param rng random number generator
- */
- private static void assertEuclidean2dVersusHypot(final int minExp, final int maxExp, final int samples,
- final UniformRandomProvider rng) {
- // generate random inputs
- final double[][] inputs = new double[samples][];
- for (int i = 0; i < samples; ++i) {
- inputs[i] = DoubleTestUtils.randomArray(2, minExp, maxExp, rng);
- }
-
- // compute exact results
- final double[] exactResults = new double[samples];
- for (int i = 0; i < samples; ++i) {
- exactResults[i] = exactEuclideanNorm(inputs[i]);
- }
-
- // compute the std devs
- final UlpErrorStats hypotStats = computeUlpErrorStats(inputs, exactResults, v -> Math.hypot(v[0], v[1]));
- final UlpErrorStats normStats = computeUlpErrorStats(inputs, exactResults, v -> Norms.euclidean(v[0], v[1]));
-
- // ensure that we are within the ballpark of Math.hypot
- Assertions.assertTrue(normStats.getMean() <= (hypotStats.getMean() + HYPOT_COMPARE_EPS),
- () -> "Expected 2D norm result to have similar error mean to Math.hypot(): hypot error mean= " +
- hypotStats.getMean() + ", norm error mean= " + normStats.getMean());
-
- Assertions.assertTrue(normStats.getStdDev() <= (hypotStats.getStdDev() + HYPOT_COMPARE_EPS),
- () -> "Expected 2D norm result to have similar std deviation to Math.hypot(): hypot std dev= " +
- hypotStats.getStdDev() + ", norm std dev= " + normStats.getStdDev());
- }
-
- @Test
- void testEuclidean_3d_simple() {
- // act/assert
- Assertions.assertEquals(0d, Norms.euclidean(0d, 0d, 0d));
- Assertions.assertEquals(1d, Norms.euclidean(1d, 0d, 0d));
- Assertions.assertEquals(1d, Norms.euclidean(0d, 1d, 0d));
- Assertions.assertEquals(1d, Norms.euclidean(0d, 0d, 1d));
- Assertions.assertEquals(5 * Math.sqrt(2), Norms.euclidean(-3d, -4d, 5d));
- Assertions.assertEquals(Double.MIN_VALUE, Norms.euclidean(0d, 0d, Double.MIN_VALUE));
- Assertions.assertEquals(Double.MAX_VALUE, Norms.euclidean(Double.MAX_VALUE, 0d, 0d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(Double.MAX_VALUE, Double.MAX_VALUE, Double.MAX_VALUE));
-
- Assertions.assertEquals(Math.sqrt(3), Norms.euclidean(1d, -1d, 1d));
-
- Assertions.assertEquals(Double.NaN, Norms.euclidean(Double.NaN, -2d, 0d));
- Assertions.assertEquals(Double.NaN, Norms.euclidean(-2d, Double.NaN, 0d));
- Assertions.assertEquals(Double.NaN, Norms.euclidean(-2d, 0d, Double.NaN));
- Assertions.assertEquals(Double.NaN,
- Norms.euclidean(Double.POSITIVE_INFINITY, Double.NaN, 1d));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, 1d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY));
- }
-
- @Test
- void testEuclidean_3d_scaled() {
- // arrange
- final double[] ones = new double[] {1, 1, 1};
- final double[] multiplesOfTen = new double[] {1, 10, 100};
- final ToDoubleFunction<double[]> fn = v -> Norms.euclidean(v[0], v[1], v[2]);
-
- // act/assert
- checkScaledEuclideanNorm(ones, 0, fn);
- checkScaledEuclideanNorm(ones, LARGE_THRESH_EXP, fn);
- checkScaledEuclideanNorm(ones, LARGE_THRESH_EXP + 1, fn);
- checkScaledEuclideanNorm(ones, -100, fn);
- checkScaledEuclideanNorm(ones, -101, fn);
- checkScaledEuclideanNorm(ones, SMALL_THRESH_EXP, fn);
- checkScaledEuclideanNorm(ones, SMALL_THRESH_EXP - 1, fn);
-
- checkScaledEuclideanNorm(multiplesOfTen, 0, fn);
- checkScaledEuclideanNorm(multiplesOfTen, -100, fn);
- checkScaledEuclideanNorm(multiplesOfTen, -101, fn);
- checkScaledEuclideanNorm(multiplesOfTen, LARGE_THRESH_EXP - 1, fn);
- checkScaledEuclideanNorm(multiplesOfTen, SMALL_THRESH_EXP - 1, fn);
- }
-
- @Test
- void testEuclidean_3d_random() {
- // arrange
- final UniformRandomProvider rng = RandomSource.create(RandomSource.XO_RO_SHI_RO_1024_PP, 1L);
-
- // act/assert
- checkEuclideanRandom(3, rng, v -> Norms.euclidean(v[0], v[1], v[2]));
- }
-
- @Test
- void testEuclidean_3d_vsArray() {
- // arrange
- final double[][] inputs = {
- {-4.074598908124454E-9, 9.897869969944898E-28, 7.849935157082846E-14},
- {1.3472131556526359E-27, -9.064577177323565E9, 323771.526282239},
- {-3.9219339341360245E149, -7.132522817112096E148, -3.427334456813165E147},
- {-1.4888098520466735E153, -2.9099184907796666E150, 1.0144962310234785E152},
- {-8.659395144898396E-152, -1.123275532302136E-150, -2.151505326692001E-152},
- {-3.660198254902351E-152, -6.656524053354807E-153, -3.198606556986218E-154}
- };
-
- // act/assert
- for (final double[] input : inputs) {
- Assertions.assertEquals(Norms.euclidean(input), Norms.euclidean(input[0], input[1], input[2]),
- () -> "Expected inline method result to equal array result for input " + Arrays.toString(input));
- }
- }
-
- @Test
- void testEuclidean_array_simple() {
- // act/assert
- Assertions.assertEquals(0d, Norms.euclidean(new double[0]));
- Assertions.assertEquals(5d, Norms.euclidean(new double[] {-3d, 4d}));
-
- Assertions.assertEquals(Math.sqrt(2), Norms.euclidean(new double[] {1d, -1d}));
- Assertions.assertEquals(Math.sqrt(3), Norms.euclidean(new double[] {1d, -1d, 1d}));
- Assertions.assertEquals(2, Norms.euclidean(new double[] {1d, -1d, 1d, -1d}));
-
- final double[] longVec = new double[] {-0.9, 8.7, -6.5, -4.3, -2.1, 0, 1.2, 3.4, -5.6, 7.8, 9.0};
- Assertions.assertEquals(directEuclideanNorm(longVec), Norms.euclidean(longVec));
-
- Assertions.assertEquals(Double.MIN_VALUE, Norms.euclidean(new double[] {0d, Double.MIN_VALUE}));
- Assertions.assertEquals(Double.MAX_VALUE, Norms.euclidean(new double[] {Double.MAX_VALUE, 0d}));
-
- final double[] maxVec = new double[1000];
- Arrays.fill(maxVec, Double.MAX_VALUE);
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.euclidean(maxVec));
-
- final double[] largeThreshVec = new double[1000];
- Arrays.fill(largeThreshVec, 0x1.0p496);
- Assertions.assertEquals(Math.sqrt(largeThreshVec.length) * largeThreshVec[0], Norms.euclidean(largeThreshVec));
-
- Assertions.assertEquals(Double.NaN, Norms.euclidean(new double[] {-2d, Double.NaN, 1d}));
- Assertions.assertEquals(Double.NaN,
- Norms.euclidean(new double[] {Double.POSITIVE_INFINITY, Double.NaN}));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(new double[] {Double.POSITIVE_INFINITY, 1, 0}));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(new double[] {Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY}));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.euclidean(new double[] {Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY}));
- }
-
- @Test
- void testEuclidean_array_scaled() {
- // arrange
- final double[] ones = new double[] {1, 1, 1, 1};
- final double[] multiplesOfTen = new double[] {1, 10, 100, 1000};
-
- // act/assert
- checkScaledEuclideanNorm(ones, 0, Norms::euclidean);
- checkScaledEuclideanNorm(ones, LARGE_THRESH_EXP, Norms::euclidean);
- checkScaledEuclideanNorm(ones, LARGE_THRESH_EXP + 1, Norms::euclidean);
- checkScaledEuclideanNorm(ones, SMALL_THRESH_EXP, Norms::euclidean);
- checkScaledEuclideanNorm(ones, SMALL_THRESH_EXP - 1, Norms::euclidean);
-
- checkScaledEuclideanNorm(multiplesOfTen, 1, Norms::euclidean);
- checkScaledEuclideanNorm(multiplesOfTen, LARGE_THRESH_EXP - 1, Norms::euclidean);
- checkScaledEuclideanNorm(multiplesOfTen, SMALL_THRESH_EXP - 1, Norms::euclidean);
- }
-
- @Test
- void testEuclidean_array_random() {
- // arrange
- final UniformRandomProvider rng = RandomSource.create(RandomSource.XO_RO_SHI_RO_1024_PP, 1L);
-
- // act/assert
- checkEuclideanRandom(2, rng, Norms::euclidean);
- checkEuclideanRandom(3, rng, Norms::euclidean);
- checkEuclideanRandom(4, rng, Norms::euclidean);
- checkEuclideanRandom(10, rng, Norms::euclidean);
- checkEuclideanRandom(100, rng, Norms::euclidean);
- }
-
- @Test
- void testMaximum_2d() {
- // act/assert
- Assertions.assertEquals(0d, Norms.maximum(0d, -0d));
- Assertions.assertEquals(2d, Norms.maximum(1d, -2d));
- Assertions.assertEquals(3d, Norms.maximum(3d, 1d));
- Assertions.assertEquals(Double.MAX_VALUE, Norms.maximum(Double.MAX_VALUE, Double.MAX_VALUE));
-
- Assertions.assertEquals(Double.NaN, Norms.maximum(Double.NaN, 0d));
- Assertions.assertEquals(Double.NaN, Norms.maximum(0d, Double.NaN));
- Assertions.assertEquals(Double.NaN, Norms.maximum(Double.POSITIVE_INFINITY, Double.NaN));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.maximum(Double.POSITIVE_INFINITY, 0d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.maximum(Double.NEGATIVE_INFINITY, 0d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.maximum(0d, Double.NEGATIVE_INFINITY));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.maximum(Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY));
- }
-
- @Test
- void testMaximum_3d() {
- // act/assert
- Assertions.assertEquals(0d, Norms.maximum(0d, -0d, 0d));
- Assertions.assertEquals(3d, Norms.maximum(1d, -2d, 3d));
- Assertions.assertEquals(4d, Norms.maximum(-4d, -2d, 3d));
- Assertions.assertEquals(Double.MAX_VALUE, Norms.maximum(Double.MAX_VALUE, Double.MAX_VALUE, Double.MAX_VALUE));
-
- Assertions.assertEquals(Double.NaN, Norms.maximum(Double.NaN, 3d, 0d));
- Assertions.assertEquals(Double.NaN, Norms.maximum(3d, Double.NaN, 0d));
- Assertions.assertEquals(Double.NaN, Norms.maximum(3d, 0d, Double.NaN));
- Assertions.assertEquals(Double.NaN, Norms.maximum(Double.POSITIVE_INFINITY, 0d, Double.NaN));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.maximum(Double.POSITIVE_INFINITY, 0d, 1d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.maximum(0d, Double.POSITIVE_INFINITY, 1d));
- Assertions.assertEquals(Double.POSITIVE_INFINITY, Norms.maximum(0d, 1d, Double.NEGATIVE_INFINITY));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.maximum(Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY));
- }
-
- @Test
- void testMaximum_array() {
- // act/assert
- Assertions.assertEquals(0d, Norms.maximum(new double[0]));
- Assertions.assertEquals(0d, Norms.maximum(new double[] {0d, -0d}));
- Assertions.assertEquals(3d, Norms.maximum(new double[] {-1d, 2d, -3d}));
- Assertions.assertEquals(4d, Norms.maximum(new double[] {-1d, 2d, -3d, 4d}));
- Assertions.assertEquals(Double.MAX_VALUE, Norms.maximum(new double[] {Double.MAX_VALUE, Double.MAX_VALUE}));
-
- Assertions.assertEquals(Double.NaN, Norms.maximum(new double[] {-2d, Double.NaN, 1d}));
- Assertions.assertEquals(Double.NaN, Norms.maximum(new double[] {Double.POSITIVE_INFINITY, Double.NaN}));
-
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.maximum(new double[] {0d, Double.POSITIVE_INFINITY}));
- Assertions.assertEquals(Double.POSITIVE_INFINITY,
- Norms.maximum(new double[] {Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY}));
- }
-
- /** Check a number of random vectors of length {@code len} with various exponent
- * ranges.
- * @param len vector array length
- * @param rng random number generator
- * @param fn euclidean norm test function
- */
- private static void checkEuclideanRandom(final int len, final UniformRandomProvider rng,
- final ToDoubleFunction<double[]> fn) {
- checkEuclideanRandom(len, +600, +620, rng, fn);
- checkEuclideanRandom(len, LARGE_THRESH_EXP - 10, LARGE_THRESH_EXP + 10, rng, fn);
- checkEuclideanRandom(len, +400, +420, rng, fn);
- checkEuclideanRandom(len, +100, +120, rng, fn);
- checkEuclideanRandom(len, -10, +10, rng, fn);
- checkEuclideanRandom(len, -120, -100, rng, fn);
- checkEuclideanRandom(len, -420, -400, rng, fn);
- checkEuclideanRandom(len, SMALL_THRESH_EXP - 10, SMALL_THRESH_EXP + 10, rng, fn);
- checkEuclideanRandom(len, -620, -600, rng, fn);
-
- checkEuclideanRandom(len, -600, +600, rng, fn);
- }
-
- /** Check a number of random vectors of length {@code len} with elements containing
- * exponents in the range {@code [minExp, maxExp]}.
- * @param len vector array length
- * @param minExp min exponent
- * @param maxExp max exponent
- * @param rng random number generator
- * @param fn euclidean norm test function
- */
- private static void checkEuclideanRandom(final int len, final int minExp, final int maxExp,
- final UniformRandomProvider rng, final ToDoubleFunction<double[]> fn) {
- for (int i = 0; i < RAND_VECTOR_CNT; ++i) {
- // arrange
- final double[] v = DoubleTestUtils.randomArray(len, minExp, maxExp, rng);
-
- final double exact = exactEuclideanNorm(v);
- final double direct = directEuclideanNorm(v);
-
- // act
- final double actual = fn.applyAsDouble(v);
-
- // assert
- Assertions.assertTrue(Double.isFinite(actual), () ->
- "Computed norm was not finite; vector= " + Arrays.toString(v) + ", exact= " + exact +
- ", direct= " + direct + ", actual= " + actual);
-
- final int ulpError = Math.abs(DoubleTestUtils.computeUlpDifference(exact, actual));
-
- Assertions.assertTrue(ulpError <= MAX_ULP_ERR, () ->
- "Computed norm ulp error exceeds bounds; vector= " + Arrays.toString(v) +
- ", exact= " + exact + ", actual= " + actual + ", ulpError= " + ulpError);
- }
- }
-
- /** Assert that {@code directNorm(v) * 2^scaleExp = fn(v * 2^scaleExp)}.
- * @param v unscaled vector
- * @param scaleExp scale factor exponent
- * @param fn euclidean norm function
- */
- private static void checkScaledEuclideanNorm(final double[] v, final int scaleExp,
- final ToDoubleFunction<double[]> fn) {
-
- final double scale = Math.scalb(1d, scaleExp);
- final double[] scaledV = new double[v.length];
- for (int i = 0; i < v.length; ++i) {
- scaledV[i] = v[i] * scale;
- }
-
- final double norm = directEuclideanNorm(v);
- final double scaledNorm = fn.applyAsDouble(scaledV);
-
- Assertions.assertEquals(norm * scale, scaledNorm);
- }
-
- /** Direct euclidean norm computation.
- * @param v array
- * @return euclidean norm using direct summation.
- */
- private static double directEuclideanNorm(final double[] v) {
- double n = 0;
- for (int i = 0; i < v.length; i++) {
- n += v[i] * v[i];
- }
- return Math.sqrt(n);
- }
-
- /** Compute the exact double value of the vector norm using BigDecimals
- * with a math context of {@link MathContext#DECIMAL128}.
- * @param v array
- * @return euclidean norm using BigDecimal with MathContext.DECIMAL128
- */
- private static double exactEuclideanNorm(final double[] v) {
- final MathContext ctx = MathContext.DECIMAL128;
-
- BigDecimal sum = BigDecimal.ZERO;
- for (final double d : v) {
- sum = sum.add(new BigDecimal(d).pow(2), ctx);
- }
-
- // Java 9+:
- // sum.sqrt(ctx).doubleValue()
-
- // Require the sum to be in the range of a double for conversion before sqrt().
- // We scale by a power of 2. Rescaling uses the square root of this which is also
- // a power of 2 and can be accumulated for exact rescaling.
- double rescale = 1.0;
- if (sum.compareTo(BD_MIN_NORMAL) < 0) {
- while (sum.compareTo(BD_MIN_NORMAL) < 0) {
- sum = sum.multiply(SCALE2);
- rescale /= SCALE;
- }
- } else if (sum.compareTo(BD_MAX_VALUE) > 0) {
- while (sum.compareTo(BD_MAX_VALUE) > 0) {
- sum = sum.divide(SCALE2);
- rescale *= SCALE;
- }
- }
-
- return Math.sqrt(sum.doubleValue()) * rescale;
- }
-
- /** Compute statistics for the ulp error of {@code fn} for the given inputs and
- * array of exact results.
- * @param inputs sample inputs
- * @param exactResults array containing the exact expected results
- * @param fn function to perform the computation
- * @return ulp error statistics
- */
- private static UlpErrorStats computeUlpErrorStats(final double[][] inputs, final double[] exactResults,
- final ToDoubleFunction<double[]> fn) {
-
- // compute the ulp errors for each input
- final int[] ulpErrors = new int[inputs.length];
- int sum = 0;
- for (int i = 0; i < inputs.length; ++i) {
- final double exact = exactResults[i];
- final double actual = fn.applyAsDouble(inputs[i]);
-
- final int error = DoubleTestUtils.computeUlpDifference(exact, actual);
- ulpErrors[i] = error;
- sum += error;
- }
-
- // compute the mean
- final double mean = sum / (double) ulpErrors.length;
-
- // compute the std dev
- double diffSumSq = 0d;
- double diff;
- for (int ulpError : ulpErrors) {
- diff = ulpError - mean;
- diffSumSq += diff * diff;
- }
-
- final double stdDev = Math.sqrt(diffSumSq / (inputs.length - 1));
-
- return new UlpErrorStats(mean, stdDev);
- }
-
- /** Class containing ULP error statistics. */
- private static final class UlpErrorStats {
-
- private final double mean;
-
- private final double stdDev;
-
- UlpErrorStats(final double mean, final double stdDev) {
- this.mean = mean;
- this.stdDev = stdDev;
- }
-
- public double getMean() {
- return mean;
- }
-
- public double getStdDev() {
- return stdDev;
- }
- }
-}