All Downloads are FREE. Search and download functionalities are using the official Maven repository.

org.apache.commons.math3.util.Precision Maven / Gradle / Ivy

There is a newer version: 5.0.70
Show newest version
/*
 * 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.math3.util;

import java.math.BigDecimal;

import org.apache.commons.math3.exception.MathArithmeticException;
import org.apache.commons.math3.exception.MathIllegalArgumentException;
import org.apache.commons.math3.exception.util.LocalizedFormats;

/**
 * Utilities for comparing numbers.
 *
 * @since 3.0
 */
public class Precision {
    /**
     * 

* Largest double-precision floating-point number such that * {@code 1 + EPSILON} is numerically equal to 1. This value is an upper * bound on the relative error due to rounding real numbers to double * precision floating-point numbers. *

*

* In IEEE 754 arithmetic, this is 2-53. *

* * @see Machine epsilon */ public static final double EPSILON; /** * Safe minimum, such that {@code 1 / SAFE_MIN} does not overflow. *
* In IEEE 754 arithmetic, this is also the smallest normalized * number 2-1022. */ public static final double SAFE_MIN; /** Exponent offset in IEEE754 representation. */ private static final long EXPONENT_OFFSET = 1023l; /** Offset to order signed double numbers lexicographically. */ private static final long SGN_MASK = 0x8000000000000000L; /** Offset to order signed double numbers lexicographically. */ private static final int SGN_MASK_FLOAT = 0x80000000; /** Positive zero. */ private static final double POSITIVE_ZERO = 0d; /** Positive zero bits. */ private static final long POSITIVE_ZERO_DOUBLE_BITS = Double.doubleToRawLongBits(+0.0); /** Negative zero bits. */ private static final long NEGATIVE_ZERO_DOUBLE_BITS = Double.doubleToRawLongBits(-0.0); /** Positive zero bits. */ private static final int POSITIVE_ZERO_FLOAT_BITS = Float.floatToRawIntBits(+0.0f); /** Negative zero bits. */ private static final int NEGATIVE_ZERO_FLOAT_BITS = Float.floatToRawIntBits(-0.0f); static { /* * This was previously expressed as = 0x1.0p-53; * However, OpenJDK (Sparc Solaris) cannot handle such small * constants: MATH-721 */ EPSILON = Double.longBitsToDouble((EXPONENT_OFFSET - 53l) << 52); /* * This was previously expressed as = 0x1.0p-1022; * However, OpenJDK (Sparc Solaris) cannot handle such small * constants: MATH-721 */ SAFE_MIN = Double.longBitsToDouble((EXPONENT_OFFSET - 1022l) << 52); } /** * Private constructor. */ private Precision() {} /** * Compares two numbers given some amount of allowed error. * * @param x the first number * @param y the second number * @param eps the amount of error to allow when checking for equality * @return
  • 0 if {@link #equals(double, double, double) equals(x, y, eps)}
  • *
  • < 0 if !{@link #equals(double, double, double) equals(x, y, eps)} && x < y
  • *
  • > 0 if !{@link #equals(double, double, double) equals(x, y, eps)} && x > y or * either argument is NaN
*/ public static int compareTo(double x, double y, double eps) { if (equals(x, y, eps)) { return 0; } else if (x < y) { return -1; } return 1; } /** * Compares two numbers given some amount of allowed error. * Two float numbers are considered equal if there are {@code (maxUlps - 1)} * (or fewer) floating point numbers between them, i.e. two adjacent floating * point numbers are considered equal. * Adapted from * Bruce Dawson. Returns {@code false} if either of the arguments is NaN. * * @param x first value * @param y second value * @param maxUlps {@code (maxUlps - 1)} is the number of floating point * values between {@code x} and {@code y}. * @return
  • 0 if {@link #equals(double, double, int) equals(x, y, maxUlps)}
  • *
  • < 0 if !{@link #equals(double, double, int) equals(x, y, maxUlps)} && x < y
  • *
  • > 0 if !{@link #equals(double, double, int) equals(x, y, maxUlps)} && x > y * or either argument is NaN
*/ public static int compareTo(final double x, final double y, final int maxUlps) { if (equals(x, y, maxUlps)) { return 0; } else if (x < y) { return -1; } return 1; } /** * Returns true iff they are equal as defined by * {@link #equals(float,float,int) equals(x, y, 1)}. * * @param x first value * @param y second value * @return {@code true} if the values are equal. */ public static boolean equals(float x, float y) { return equals(x, y, 1); } /** * Returns true if both arguments are NaN or they are * equal as defined by {@link #equals(float,float) equals(x, y, 1)}. * * @param x first value * @param y second value * @return {@code true} if the values are equal or both are NaN. * @since 2.2 */ public static boolean equalsIncludingNaN(float x, float y) { return (x != x || y != y) ? !(x != x ^ y != y) : equals(x, y, 1); } /** * Returns true if the arguments are equal or within the range of allowed * error (inclusive). Returns {@code false} if either of the arguments * is NaN. * * @param x first value * @param y second value * @param eps the amount of absolute error to allow. * @return {@code true} if the values are equal or within range of each other. * @since 2.2 */ public static boolean equals(float x, float y, float eps) { return equals(x, y, 1) || FastMath.abs(y - x) <= eps; } /** * Returns true if the arguments are both NaN, are equal, or are within the range * of allowed error (inclusive). * * @param x first value * @param y second value * @param eps the amount of absolute error to allow. * @return {@code true} if the values are equal or within range of each other, * or both are NaN. * @since 2.2 */ public static boolean equalsIncludingNaN(float x, float y, float eps) { return equalsIncludingNaN(x, y) || (FastMath.abs(y - x) <= eps); } /** * Returns true if the arguments are equal or within the range of allowed * error (inclusive). * Two float numbers are considered equal if there are {@code (maxUlps - 1)} * (or fewer) floating point numbers between them, i.e. two adjacent floating * point numbers are considered equal. * Adapted from * Bruce Dawson. Returns {@code false} if either of the arguments is NaN. * * @param x first value * @param y second value * @param maxUlps {@code (maxUlps - 1)} is the number of floating point * values between {@code x} and {@code y}. * @return {@code true} if there are fewer than {@code maxUlps} floating * point values between {@code x} and {@code y}. * @since 2.2 */ public static boolean equals(final float x, final float y, final int maxUlps) { final int xInt = Float.floatToRawIntBits(x); final int yInt = Float.floatToRawIntBits(y); final boolean isEqual; if (((xInt ^ yInt) & SGN_MASK_FLOAT) == 0) { // number have same sign, there is no risk of overflow isEqual = FastMath.abs(xInt - yInt) <= maxUlps; } else { // number have opposite signs, take care of overflow final int deltaPlus; final int deltaMinus; if (xInt < yInt) { deltaPlus = yInt - POSITIVE_ZERO_FLOAT_BITS; deltaMinus = xInt - NEGATIVE_ZERO_FLOAT_BITS; } else { deltaPlus = xInt - POSITIVE_ZERO_FLOAT_BITS; deltaMinus = yInt - NEGATIVE_ZERO_FLOAT_BITS; } if (deltaPlus > maxUlps) { isEqual = false; } else { isEqual = deltaMinus <= (maxUlps - deltaPlus); } } return isEqual && !Float.isNaN(x) && !Float.isNaN(y); } /** * Returns true if the arguments are both NaN or if they are equal as defined * by {@link #equals(float,float,int) equals(x, y, maxUlps)}. * * @param x first value * @param y second value * @param maxUlps {@code (maxUlps - 1)} is the number of floating point * values between {@code x} and {@code y}. * @return {@code true} if both arguments are NaN or if there are less than * {@code maxUlps} floating point values between {@code x} and {@code y}. * @since 2.2 */ public static boolean equalsIncludingNaN(float x, float y, int maxUlps) { return (x != x || y != y) ? !(x != x ^ y != y) : equals(x, y, maxUlps); } /** * Returns true iff they are equal as defined by * {@link #equals(double,double,int) equals(x, y, 1)}. * * @param x first value * @param y second value * @return {@code true} if the values are equal. */ public static boolean equals(double x, double y) { return equals(x, y, 1); } /** * Returns true if the arguments are both NaN or they are * equal as defined by {@link #equals(double,double) equals(x, y, 1)}. * * @param x first value * @param y second value * @return {@code true} if the values are equal or both are NaN. * @since 2.2 */ public static boolean equalsIncludingNaN(double x, double y) { return (x != x || y != y) ? !(x != x ^ y != y) : equals(x, y, 1); } /** * Returns {@code true} if there is no double value strictly between the * arguments or the difference between them is within the range of allowed * error (inclusive). Returns {@code false} if either of the arguments * is NaN. * * @param x First value. * @param y Second value. * @param eps Amount of allowed absolute error. * @return {@code true} if the values are two adjacent floating point * numbers or they are within range of each other. */ public static boolean equals(double x, double y, double eps) { return equals(x, y, 1) || FastMath.abs(y - x) <= eps; } /** * Returns {@code true} if there is no double value strictly between the * arguments or the relative difference between them is less than or equal * to the given tolerance. Returns {@code false} if either of the arguments * is NaN. * * @param x First value. * @param y Second value. * @param eps Amount of allowed relative error. * @return {@code true} if the values are two adjacent floating point * numbers or they are within range of each other. * @since 3.1 */ public static boolean equalsWithRelativeTolerance(double x, double y, double eps) { if (equals(x, y, 1)) { return true; } final double absoluteMax = FastMath.max(FastMath.abs(x), FastMath.abs(y)); final double relativeDifference = FastMath.abs((x - y) / absoluteMax); return relativeDifference <= eps; } /** * Returns true if the arguments are both NaN, are equal or are within the range * of allowed error (inclusive). * * @param x first value * @param y second value * @param eps the amount of absolute error to allow. * @return {@code true} if the values are equal or within range of each other, * or both are NaN. * @since 2.2 */ public static boolean equalsIncludingNaN(double x, double y, double eps) { return equalsIncludingNaN(x, y) || (FastMath.abs(y - x) <= eps); } /** * Returns true if the arguments are equal or within the range of allowed * error (inclusive). *

* Two float numbers are considered equal if there are {@code (maxUlps - 1)} * (or fewer) floating point numbers between them, i.e. two adjacent * floating point numbers are considered equal. *

*

* Adapted from * Bruce Dawson. Returns {@code false} if either of the arguments is NaN. *

* * @param x first value * @param y second value * @param maxUlps {@code (maxUlps - 1)} is the number of floating point * values between {@code x} and {@code y}. * @return {@code true} if there are fewer than {@code maxUlps} floating * point values between {@code x} and {@code y}. */ public static boolean equals(final double x, final double y, final int maxUlps) { final long xInt = Double.doubleToRawLongBits(x); final long yInt = Double.doubleToRawLongBits(y); final boolean isEqual; if (((xInt ^ yInt) & SGN_MASK) == 0l) { // number have same sign, there is no risk of overflow isEqual = FastMath.abs(xInt - yInt) <= maxUlps; } else { // number have opposite signs, take care of overflow final long deltaPlus; final long deltaMinus; if (xInt < yInt) { deltaPlus = yInt - POSITIVE_ZERO_DOUBLE_BITS; deltaMinus = xInt - NEGATIVE_ZERO_DOUBLE_BITS; } else { deltaPlus = xInt - POSITIVE_ZERO_DOUBLE_BITS; deltaMinus = yInt - NEGATIVE_ZERO_DOUBLE_BITS; } if (deltaPlus > maxUlps) { isEqual = false; } else { isEqual = deltaMinus <= (maxUlps - deltaPlus); } } return isEqual && !Double.isNaN(x) && !Double.isNaN(y); } /** * Returns true if both arguments are NaN or if they are equal as defined * by {@link #equals(double,double,int) equals(x, y, maxUlps)}. * * @param x first value * @param y second value * @param maxUlps {@code (maxUlps - 1)} is the number of floating point * values between {@code x} and {@code y}. * @return {@code true} if both arguments are NaN or if there are less than * {@code maxUlps} floating point values between {@code x} and {@code y}. * @since 2.2 */ public static boolean equalsIncludingNaN(double x, double y, int maxUlps) { return (x != x || y != y) ? !(x != x ^ y != y) : equals(x, y, maxUlps); } /** * Rounds the given value to the specified number of decimal places. * The value is rounded using the {@link BigDecimal#ROUND_HALF_UP} method. * * @param x Value to round. * @param scale Number of digits to the right of the decimal point. * @return the rounded value. * @since 1.1 (previously in {@code MathUtils}, moved as of version 3.0) */ public static double round(double x, int scale) { return round(x, scale, BigDecimal.ROUND_HALF_UP); } /** * Rounds the given value to the specified number of decimal places. * The value is rounded using the given method which is any method defined * in {@link BigDecimal}. * If {@code x} is infinite or {@code NaN}, then the value of {@code x} is * returned unchanged, regardless of the other parameters. * * @param x Value to round. * @param scale Number of digits to the right of the decimal point. * @param roundingMethod Rounding method as defined in {@link BigDecimal}. * @return the rounded value. * @throws ArithmeticException if {@code roundingMethod == ROUND_UNNECESSARY} * and the specified scaling operation would require rounding. * @throws IllegalArgumentException if {@code roundingMethod} does not * represent a valid rounding mode. * @since 1.1 (previously in {@code MathUtils}, moved as of version 3.0) */ public static double round(double x, int scale, int roundingMethod) { try { final double rounded = (new BigDecimal(Double.toString(x)) .setScale(scale, roundingMethod)) .doubleValue(); // MATH-1089: negative values rounded to zero should result in negative zero return rounded == POSITIVE_ZERO ? POSITIVE_ZERO * x : rounded; } catch (NumberFormatException ex) { if (Double.isInfinite(x)) { return x; } else { return Double.NaN; } } } /** * Rounds the given value to the specified number of decimal places. * The value is rounded using the {@link BigDecimal#ROUND_HALF_UP} method. * * @param x Value to round. * @param scale Number of digits to the right of the decimal point. * @return the rounded value. * @since 1.1 (previously in {@code MathUtils}, moved as of version 3.0) */ public static float round(float x, int scale) { return round(x, scale, BigDecimal.ROUND_HALF_UP); } /** * Rounds the given value to the specified number of decimal places. * The value is rounded using the given method which is any method defined * in {@link BigDecimal}. * * @param x Value to round. * @param scale Number of digits to the right of the decimal point. * @param roundingMethod Rounding method as defined in {@link BigDecimal}. * @return the rounded value. * @since 1.1 (previously in {@code MathUtils}, moved as of version 3.0) * @throws MathArithmeticException if an exact operation is required but result is not exact * @throws MathIllegalArgumentException if {@code roundingMethod} is not a valid rounding method. */ public static float round(float x, int scale, int roundingMethod) throws MathArithmeticException, MathIllegalArgumentException { final float sign = FastMath.copySign(1f, x); final float factor = (float) FastMath.pow(10.0f, scale) * sign; return (float) roundUnscaled(x * factor, sign, roundingMethod) / factor; } /** * Rounds the given non-negative value to the "nearest" integer. Nearest is * determined by the rounding method specified. Rounding methods are defined * in {@link BigDecimal}. * * @param unscaled Value to round. * @param sign Sign of the original, scaled value. * @param roundingMethod Rounding method, as defined in {@link BigDecimal}. * @return the rounded value. * @throws MathArithmeticException if an exact operation is required but result is not exact * @throws MathIllegalArgumentException if {@code roundingMethod} is not a valid rounding method. * @since 1.1 (previously in {@code MathUtils}, moved as of version 3.0) */ private static double roundUnscaled(double unscaled, double sign, int roundingMethod) throws MathArithmeticException, MathIllegalArgumentException { switch (roundingMethod) { case BigDecimal.ROUND_CEILING : if (sign == -1) { unscaled = FastMath.floor(FastMath.nextAfter(unscaled, Double.NEGATIVE_INFINITY)); } else { unscaled = FastMath.ceil(FastMath.nextAfter(unscaled, Double.POSITIVE_INFINITY)); } break; case BigDecimal.ROUND_DOWN : unscaled = FastMath.floor(FastMath.nextAfter(unscaled, Double.NEGATIVE_INFINITY)); break; case BigDecimal.ROUND_FLOOR : if (sign == -1) { unscaled = FastMath.ceil(FastMath.nextAfter(unscaled, Double.POSITIVE_INFINITY)); } else { unscaled = FastMath.floor(FastMath.nextAfter(unscaled, Double.NEGATIVE_INFINITY)); } break; case BigDecimal.ROUND_HALF_DOWN : { unscaled = FastMath.nextAfter(unscaled, Double.NEGATIVE_INFINITY); double fraction = unscaled - FastMath.floor(unscaled); if (fraction > 0.5) { unscaled = FastMath.ceil(unscaled); } else { unscaled = FastMath.floor(unscaled); } break; } case BigDecimal.ROUND_HALF_EVEN : { double fraction = unscaled - FastMath.floor(unscaled); if (fraction > 0.5) { unscaled = FastMath.ceil(unscaled); } else if (fraction < 0.5) { unscaled = FastMath.floor(unscaled); } else { // The following equality test is intentional and needed for rounding purposes if (FastMath.floor(unscaled) / 2.0 == FastMath.floor(FastMath.floor(unscaled) / 2.0)) { // even unscaled = FastMath.floor(unscaled); } else { // odd unscaled = FastMath.ceil(unscaled); } } break; } case BigDecimal.ROUND_HALF_UP : { unscaled = FastMath.nextAfter(unscaled, Double.POSITIVE_INFINITY); double fraction = unscaled - FastMath.floor(unscaled); if (fraction >= 0.5) { unscaled = FastMath.ceil(unscaled); } else { unscaled = FastMath.floor(unscaled); } break; } case BigDecimal.ROUND_UNNECESSARY : if (unscaled != FastMath.floor(unscaled)) { throw new MathArithmeticException(); } break; case BigDecimal.ROUND_UP : // do not round if the discarded fraction is equal to zero if (unscaled != FastMath.floor(unscaled)) { unscaled = FastMath.ceil(FastMath.nextAfter(unscaled, Double.POSITIVE_INFINITY)); } break; default : throw new MathIllegalArgumentException(LocalizedFormats.INVALID_ROUNDING_METHOD, roundingMethod, "ROUND_CEILING", BigDecimal.ROUND_CEILING, "ROUND_DOWN", BigDecimal.ROUND_DOWN, "ROUND_FLOOR", BigDecimal.ROUND_FLOOR, "ROUND_HALF_DOWN", BigDecimal.ROUND_HALF_DOWN, "ROUND_HALF_EVEN", BigDecimal.ROUND_HALF_EVEN, "ROUND_HALF_UP", BigDecimal.ROUND_HALF_UP, "ROUND_UNNECESSARY", BigDecimal.ROUND_UNNECESSARY, "ROUND_UP", BigDecimal.ROUND_UP); } return unscaled; } /** * Computes a number {@code delta} close to {@code originalDelta} with * the property that

     *   x + delta - x
     * 
* is exactly machine-representable. * This is useful when computing numerical derivatives, in order to reduce * roundoff errors. * * @param x Value. * @param originalDelta Offset value. * @return a number {@code delta} so that {@code x + delta} and {@code x} * differ by a representable floating number. */ public static double representableDelta(double x, double originalDelta) { return x + originalDelta - x; } }




© 2015 - 2024 Weber Informatics LLC | Privacy Policy