1 /* 2 * Licensed to the Apache Software Foundation (ASF) under one or more 3 * contributor license agreements. See the NOTICE file distributed with 4 * this work for additional information regarding copyright ownership. 5 * The ASF licenses this file to You under the Apache License, Version 2.0 6 * (the "License"); you may not use this file except in compliance with 7 * the License. You may obtain a copy of the License at 8 * 9 * https://www.apache.org/licenses/LICENSE-2.0 10 * 11 * Unless required by applicable law or agreed to in writing, software 12 * distributed under the License is distributed on an "AS IS" BASIS, 13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 14 * See the License for the specific language governing permissions and 15 * limitations under the License. 16 */ 17 18 /* 19 * This is not the original file distributed by the Apache Software Foundation 20 * It has been modified by the Hipparchus project 21 */ 22 23 package org.hipparchus.util; 24 25 import java.io.IOException; 26 import java.io.InputStream; 27 import java.util.Arrays; 28 import java.util.Properties; 29 30 import org.hipparchus.CalculusFieldElement; 31 import org.hipparchus.FieldElement; 32 import org.hipparchus.exception.Localizable; 33 import org.hipparchus.exception.LocalizedCoreFormats; 34 import org.hipparchus.exception.MathIllegalArgumentException; 35 import org.hipparchus.exception.MathRuntimeException; 36 import org.hipparchus.exception.NullArgumentException; 37 38 /** 39 * Miscellaneous utility functions. 40 * 41 * @see ArithmeticUtils 42 * @see Precision 43 * @see MathArrays 44 */ 45 public final class MathUtils { 46 /** \(2\pi\) */ 47 public static final double TWO_PI = 2 * FastMath.PI; 48 49 /** \(\pi^2\) */ 50 public static final double PI_SQUARED = FastMath.PI * FastMath.PI; 51 52 /** \(\pi/2\). */ 53 public static final double SEMI_PI = 0.5 * FastMath.PI; 54 55 /** 56 * Class contains only static methods. 57 */ 58 private MathUtils() {} 59 60 61 /** 62 * Returns an integer hash code representing the given double value. 63 * 64 * @param value the value to be hashed 65 * @return the hash code 66 */ 67 public static int hash(double value) { 68 return Double.hashCode(value); 69 } 70 71 /** 72 * Returns {@code true} if the values are equal according to semantics of 73 * {@link Double#equals(Object)}. 74 * 75 * @param x Value 76 * @param y Value 77 * @return {@code Double.valueOf(x).equals(Double.valueOf(y))} 78 */ 79 public static boolean equals(double x, double y) { 80 return Double.valueOf(x).equals(y); 81 } 82 83 /** 84 * Returns an integer hash code representing the given double array. 85 * 86 * @param value the value to be hashed (may be null) 87 * @return the hash code 88 */ 89 public static int hash(double[] value) { 90 return Arrays.hashCode(value); 91 } 92 93 /** 94 * Normalize an angle in a 2π wide interval around a center value. 95 * <p>This method has three main uses:</p> 96 * <ul> 97 * <li>normalize an angle between 0 and 2π:<br> 98 * {@code a = MathUtils.normalizeAngle(a, FastMath.PI);}</li> 99 * <li>normalize an angle between -π and +π<br> 100 * {@code a = MathUtils.normalizeAngle(a, 0.0);}</li> 101 * <li>compute the angle between two defining angular positions:<br> 102 * {@code angle = MathUtils.normalizeAngle(end, start) - start;}</li> 103 * </ul> 104 * <p>Note that due to numerical accuracy and since π cannot be represented 105 * exactly, the result interval is <em>closed</em>, it cannot be half-closed 106 * as would be more satisfactory in a purely mathematical view.</p> 107 * @param a angle to normalize 108 * @param center center of the desired 2π interval for the result 109 * @return a-2kπ with integer k and center-π <= a-2kπ <= center+π 110 */ 111 public static double normalizeAngle(double a, double center) { 112 return a - TWO_PI * FastMath.floor((a + FastMath.PI - center) / TWO_PI); 113 } 114 115 /** 116 * Normalize an angle in a 2π wide interval around a center value. 117 * <p>This method has three main uses:</p> 118 * <ul> 119 * <li>normalize an angle between 0 and 2π:<br> 120 * {@code a = MathUtils.normalizeAngle(a, FastMath.PI);}</li> 121 * <li>normalize an angle between -π and +π<br> 122 * {@code a = MathUtils.normalizeAngle(a, zero);}</li> 123 * <li>compute the angle between two defining angular positions:<br> 124 * {@code angle = MathUtils.normalizeAngle(end, start).subtract(start);}</li> 125 * </ul> 126 * <p>Note that due to numerical accuracy and since π cannot be represented 127 * exactly, the result interval is <em>closed</em>, it cannot be half-closed 128 * as would be more satisfactory in a purely mathematical view.</p> 129 * @param <T> the type of the field elements 130 * @param a angle to normalize 131 * @param center center of the desired 2π interval for the result 132 * @return a-2kπ with integer k and center-π <= a-2kπ <= center+π 133 */ 134 public static <T extends CalculusFieldElement<T>> T normalizeAngle(T a, T center) { 135 return a.subtract(FastMath.floor(a.add(FastMath.PI).subtract(center).divide(TWO_PI)).multiply(TWO_PI)); 136 } 137 138 /** Find the maximum of two field elements. 139 * @param <T> the type of the field elements 140 * @param e1 first element 141 * @param e2 second element 142 * @return max(a1, e2) 143 */ 144 public static <T extends CalculusFieldElement<T>> T max(final T e1, final T e2) { 145 return e1.subtract(e2).getReal() >= 0 ? e1 : e2; 146 } 147 148 /** Find the minimum of two field elements. 149 * @param <T> the type of the field elements 150 * @param e1 first element 151 * @param e2 second element 152 * @return min(a1, e2) 153 */ 154 public static <T extends CalculusFieldElement<T>> T min(final T e1, final T e2) { 155 return e1.subtract(e2).getReal() >= 0 ? e2 : e1; 156 } 157 158 /** 159 * <p>Reduce {@code |a - offset|} to the primary interval 160 * {@code [0, |period|)}.</p> 161 * 162 * <p>Specifically, the value returned is <br> 163 * {@code a - |period| * floor((a - offset) / |period|) - offset}.</p> 164 * 165 * <p>If any of the parameters are {@code NaN} or infinite, the result is 166 * {@code NaN}.</p> 167 * 168 * @param a Value to reduce. 169 * @param period Period. 170 * @param offset Value that will be mapped to {@code 0}. 171 * @return the value, within the interval {@code [0 |period|)}, 172 * that corresponds to {@code a}. 173 */ 174 public static double reduce(double a, 175 double period, 176 double offset) { 177 final double p = FastMath.abs(period); 178 return a - p * FastMath.floor((a - offset) / p) - offset; 179 } 180 181 /** 182 * Returns the first argument with the sign of the second argument. 183 * 184 * @param magnitude Magnitude of the returned value. 185 * @param sign Sign of the returned value. 186 * @return a value with magnitude equal to {@code magnitude} and with the 187 * same sign as the {@code sign} argument. 188 * @throws MathRuntimeException if {@code magnitude == Byte.MIN_VALUE} 189 * and {@code sign >= 0}. 190 */ 191 public static byte copySign(byte magnitude, byte sign) 192 throws MathRuntimeException { 193 if ((magnitude >= 0 && sign >= 0) || 194 (magnitude < 0 && sign < 0)) { // Sign is OK. 195 return magnitude; 196 } else if (sign >= 0 && 197 magnitude == Byte.MIN_VALUE) { 198 throw new MathRuntimeException(LocalizedCoreFormats.OVERFLOW); 199 } else { 200 return (byte) -magnitude; // Flip sign. 201 } 202 } 203 204 /** 205 * Returns the first argument with the sign of the second argument. 206 * 207 * @param magnitude Magnitude of the returned value. 208 * @param sign Sign of the returned value. 209 * @return a value with magnitude equal to {@code magnitude} and with the 210 * same sign as the {@code sign} argument. 211 * @throws MathRuntimeException if {@code magnitude == Short.MIN_VALUE} 212 * and {@code sign >= 0}. 213 */ 214 public static short copySign(short magnitude, short sign) 215 throws MathRuntimeException { 216 if ((magnitude >= 0 && sign >= 0) || 217 (magnitude < 0 && sign < 0)) { // Sign is OK. 218 return magnitude; 219 } else if (sign >= 0 && 220 magnitude == Short.MIN_VALUE) { 221 throw new MathRuntimeException(LocalizedCoreFormats.OVERFLOW); 222 } else { 223 return (short) -magnitude; // Flip sign. 224 } 225 } 226 227 /** 228 * Returns the first argument with the sign of the second argument. 229 * 230 * @param magnitude Magnitude of the returned value. 231 * @param sign Sign of the returned value. 232 * @return a value with magnitude equal to {@code magnitude} and with the 233 * same sign as the {@code sign} argument. 234 * @throws MathRuntimeException if {@code magnitude == Integer.MIN_VALUE} 235 * and {@code sign >= 0}. 236 */ 237 public static int copySign(int magnitude, int sign) 238 throws MathRuntimeException { 239 if ((magnitude >= 0 && sign >= 0) || 240 (magnitude < 0 && sign < 0)) { // Sign is OK. 241 return magnitude; 242 } else if (sign >= 0 && 243 magnitude == Integer.MIN_VALUE) { 244 throw new MathRuntimeException(LocalizedCoreFormats.OVERFLOW); 245 } else { 246 return -magnitude; // Flip sign. 247 } 248 } 249 250 /** 251 * Returns the first argument with the sign of the second argument. 252 * 253 * @param magnitude Magnitude of the returned value. 254 * @param sign Sign of the returned value. 255 * @return a value with magnitude equal to {@code magnitude} and with the 256 * same sign as the {@code sign} argument. 257 * @throws MathRuntimeException if {@code magnitude == Long.MIN_VALUE} 258 * and {@code sign >= 0}. 259 */ 260 public static long copySign(long magnitude, long sign) 261 throws MathRuntimeException { 262 if ((magnitude >= 0 && sign >= 0) || 263 (magnitude < 0 && sign < 0)) { // Sign is OK. 264 return magnitude; 265 } else if (sign >= 0 && 266 magnitude == Long.MIN_VALUE) { 267 throw new MathRuntimeException(LocalizedCoreFormats.OVERFLOW); 268 } else { 269 return -magnitude; // Flip sign. 270 } 271 } 272 /** 273 * Check that the argument is a real number. 274 * 275 * @param x Argument. 276 * @throws MathIllegalArgumentException if {@code x} is not a 277 * finite real number. 278 */ 279 public static void checkFinite(final double x) 280 throws MathIllegalArgumentException { 281 if (Double.isInfinite(x) || Double.isNaN(x)) { 282 throw new MathIllegalArgumentException(LocalizedCoreFormats.NOT_FINITE_NUMBER, x); 283 } 284 } 285 286 /** 287 * Check that all the elements are real numbers. 288 * 289 * @param val Arguments. 290 * @throws MathIllegalArgumentException if any values of the array is not a 291 * finite real number. 292 */ 293 public static void checkFinite(final double[] val) 294 throws MathIllegalArgumentException { 295 for (final double x : val) { 296 if (Double.isInfinite(x) || Double.isNaN(x)) { 297 throw new MathIllegalArgumentException(LocalizedCoreFormats.NOT_FINITE_NUMBER, x); 298 } 299 } 300 } 301 302 /** 303 * Checks that an object is not null. 304 * 305 * @param o Object to be checked. 306 * @param pattern Message pattern. 307 * @param args Arguments to replace the placeholders in {@code pattern}. 308 * @throws NullArgumentException if {@code o} is {@code null}. 309 */ 310 public static void checkNotNull(Object o, 311 Localizable pattern, 312 Object ... args) 313 throws NullArgumentException { 314 if (o == null) { 315 throw new NullArgumentException(pattern, args); 316 } 317 } 318 319 /** 320 * Checks that an object is not null. 321 * 322 * @param o Object to be checked. 323 * @throws NullArgumentException if {@code o} is {@code null}. 324 */ 325 public static void checkNotNull(Object o) 326 throws NullArgumentException { 327 if (o == null) { 328 throw new NullArgumentException(LocalizedCoreFormats.NULL_NOT_ALLOWED); 329 } 330 } 331 332 /** 333 * Checks that the given value is strictly within the range [lo, hi]. 334 * 335 * @param value value to be checked. 336 * @param lo the lower bound (inclusive). 337 * @param hi the upper bound (inclusive). 338 * @throws MathIllegalArgumentException if {@code value} is strictly outside [lo, hi]. 339 */ 340 public static void checkRangeInclusive(long value, long lo, long hi) { 341 if (value < lo || value > hi) { 342 throw new MathIllegalArgumentException(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, 343 value, lo, hi); 344 } 345 } 346 347 /** 348 * Checks that the given value is strictly within the range [lo, hi]. 349 * 350 * @param value value to be checked. 351 * @param lo the lower bound (inclusive). 352 * @param hi the upper bound (inclusive). 353 * @throws MathIllegalArgumentException if {@code value} is strictly outside [lo, hi]. 354 */ 355 public static void checkRangeInclusive(double value, double lo, double hi) { 356 if (value < lo || value > hi) { 357 throw new MathIllegalArgumentException(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, 358 value, lo, hi); 359 } 360 } 361 362 /** 363 * Checks that the given dimensions match. 364 * 365 * @param dimension the first dimension. 366 * @param otherDimension the second dimension. 367 * @throws MathIllegalArgumentException if length != otherLength. 368 */ 369 public static void checkDimension(int dimension, int otherDimension) { 370 if (dimension != otherDimension) { 371 throw new MathIllegalArgumentException(LocalizedCoreFormats.DIMENSIONS_MISMATCH, 372 dimension, otherDimension); 373 } 374 } 375 376 /** 377 * Get Hipparchus version. 378 * <p> 379 * The version is automatically retrieved from a properties file generated 380 * at maven compilation time. When using an IDE not configured to use 381 * maven, then a default value {@code "unknown"} will be returned. 382 * </p> 383 * @return hipparchus version 384 * @since 4.0 385 */ 386 public static String getHipparchusVersion() { 387 String version = "unknown"; 388 final Properties properties = new Properties(); 389 try (InputStream stream = MathUtils.class.getResourceAsStream("/assets/org/hipparchus/hipparchus.properties")) { 390 if (stream != null) { 391 properties.load(stream); 392 version = properties.getProperty("hipparchus.version", version); 393 } 394 } catch (IOException ioe) { // NOPMD 395 // ignored 396 } 397 return version; 398 } 399 400 /** 401 * Sums {@code a} and {@code b} using Møller's 2Sum algorithm. 402 * <p> 403 * References: 404 * <ul> 405 * <li>Møller, Ole. "Quasi double-precision in floating point addition." BIT 406 * 5, 37–50 (1965).</li> 407 * <li>Shewchuk, Richard J. "Adaptive Precision Floating-Point Arithmetic 408 * and Fast Robust Geometric Predicates." Discrete & Computational Geometry 409 * 18, 305–363 (1997).</li> 410 * <li><a href= 411 * "https://en.wikipedia.org/wiki/2Sum">https://en.wikipedia.org/wiki/2Sum</a></li> 412 * </ul> 413 * @param a first summand 414 * @param b second summand 415 * @return sum and residual error in the sum 416 */ 417 public static SumAndResidual twoSum(final double a, final double b) { 418 final double s = a + b; 419 final double aPrime = s - b; 420 final double bPrime = s - aPrime; 421 final double deltaA = a - aPrime; 422 final double deltaB = b - bPrime; 423 final double t = deltaA + deltaB; 424 return new SumAndResidual(s, t); 425 } 426 427 /** 428 * Sums {@code a} and {@code b} using Møller's 2Sum algorithm. 429 * <p> 430 * References: 431 * <ul> 432 * <li>Møller, Ole. "Quasi double-precision in floating point addition." BIT 433 * 5, 37–50 (1965).</li> 434 * <li>Shewchuk, Richard J. "Adaptive Precision Floating-Point Arithmetic 435 * and Fast Robust Geometric Predicates." Discrete & Computational Geometry 436 * 18, 305–363 (1997).</li> 437 * <li><a href= 438 * "https://en.wikipedia.org/wiki/2Sum">https://en.wikipedia.org/wiki/2Sum</a></li> 439 * </ul> 440 * @param <T> field element type 441 * @param a first summand 442 * @param b second summand 443 * @return sum and residual error in the sum 444 */ 445 public static <T extends FieldElement<T>> FieldSumAndResidual<T> twoSum(final T a, final T b) { 446 final T s = a.add(b); 447 final T aPrime = s.subtract(b); 448 final T bPrime = s.subtract(aPrime); 449 final T deltaA = a.subtract(aPrime); 450 final T deltaB = b.subtract(bPrime); 451 final T t = deltaA.add(deltaB); 452 return new FieldSumAndResidual<>(s, t); 453 } 454 455 /** 456 * Result class for {@link MathUtils#twoSum(double, double)} containing the 457 * sum and the residual error in the sum. 458 */ 459 public static final class SumAndResidual { 460 461 /** Sum. */ 462 private final double sum; 463 /** Residual error in the sum. */ 464 private final double residual; 465 466 /** 467 * Constructs a {@link SumAndResidual} instance. 468 * @param sum sum 469 * @param residual residual error in the sum 470 */ 471 private SumAndResidual(final double sum, final double residual) { 472 this.sum = sum; 473 this.residual = residual; 474 } 475 476 /** 477 * Returns the sum. 478 * @return sum 479 */ 480 public double getSum() { 481 return sum; 482 } 483 484 /** 485 * Returns the residual error in the sum. 486 * @return residual error in the sum 487 */ 488 public double getResidual() { 489 return residual; 490 } 491 492 } 493 494 /** 495 * Result class for 496 * {@link MathUtils#twoSum(FieldElement, FieldElement)} containing 497 * the sum and the residual error in the sum. 498 * @param <T> field element type 499 */ 500 public static final class FieldSumAndResidual<T extends FieldElement<T>> { 501 502 /** Sum. */ 503 private final T sum; 504 /** Residual error in the sum. */ 505 private final T residual; 506 507 /** 508 * Constructs a {@link FieldSumAndResidual} instance. 509 * @param sum sum 510 * @param residual residual error in the sum 511 */ 512 private FieldSumAndResidual(final T sum, final T residual) { 513 this.sum = sum; 514 this.residual = residual; 515 } 516 517 /** 518 * Returns the sum. 519 * @return sum 520 */ 521 public T getSum() { 522 return sum; 523 } 524 525 /** 526 * Returns the residual error in the sum. 527 * @return residual error in the sum 528 */ 529 public T getResidual() { 530 return residual; 531 } 532 533 } 534 535 }