View Javadoc
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&pi;:<br>
98       *       {@code a = MathUtils.normalizeAngle(a, FastMath.PI);}</li>
99       *   <li>normalize an angle between -&pi; and +&pi;<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 &pi; 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&pi; interval for the result
109      * @return a-2k&pi; with integer k and center-&pi; &lt;= a-2k&pi; &lt;= center+&pi;
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&pi; 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&pi;:<br>
120       *       {@code a = MathUtils.normalizeAngle(a, FastMath.PI);}</li>
121       *   <li>normalize an angle between -&pi; and +&pi;<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 &pi; 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&pi; interval for the result
132       * @return a-2k&pi; with integer k and center-&pi; &lt;= a-2k&pi; &lt;= center+&pi;
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 &amp; 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 &amp; 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 }