001/**
002 * MathTools -- A collection of useful math utility routines.
003 * 
004 * Copyright (C) 1999-2023, Joseph A. Huwaldt. All rights reserved.
005 * 
006 * This library is free software; you can redistribute it and/or modify it under the terms
007 * of the GNU Lesser General Public License as published by the Free Software Foundation;
008 * either version 2 of the License, or (at your option) any later version.
009 * 
010 * This library is distributed in the hope that it will be useful, but WITHOUT ANY
011 * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
012 * PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
013 * 
014 * You should have received a copy of the GNU Lesser General Public License along with
015 * this program; if not, write to the Free Software Foundation, Inc., 59 Temple Place -
016 * Suite 330, Boston, MA 02111-1307, USA. Or visit: http://www.gnu.org/licenses/lgpl.html
017 */
018package jahuwaldt.tools.math;
019
020import java.util.BitSet;
021import java.math.BigInteger;
022
023/**
024 * A collection of useful static routines of a general mathematical nature. This file
025 * includes functions that accomplish all types of wondrous mathematical stuff.
026 * 
027 * <p> Modified by: Joseph A. Huwaldt </p>
028 * 
029 * @author Joseph A. Huwaldt, Date: September 29, 1997
030 * @version December 4, 2023
031 */
032public final class MathTools {
033
034    /**
035     * The natural logarithm of 10.
036     */
037    public static final double LOG10 = Math.log(10);
038
039    /**
040     * The natural logarithm of 2.
041     */
042    public static final double LOG2 = Math.log(2);
043
044    /**
045     * The natural logarithm of the maximum double value: log(MAX_VALUE).
046     */
047    public static final double MAX_LOG = Math.log(Double.MAX_VALUE);
048
049    /**
050     * The natural logarithm of the minimum double value: log(MIN_VALUE).
051     */
052    public static final double MIN_LOG = Math.log(Double.MIN_VALUE);
053
054    /**
055     * The machine epsilon (macheps) or unit roundoff for <code>double</code> in the Java
056     * environment. Machine epsilon gives an upper bound on the relative error due to
057     * rounding in floating point arithmetic. Machine epsilon is the smallest number such
058     * that (1.0 + EPS != 1.0).
059     */
060    public static final double EPS = epsilon(1.0);
061
062    /**
063     * Square-root of the machine epsilon for <code>double</code>.
064     */
065    public static final double SQRT_EPS = Math.sqrt(EPS);
066
067    /**
068     * The machine epsilon (macheps) or unit roundoff for <code>float</code> in the Java
069     * environment. Machine epsilon gives an upper bound on the relative error due to
070     * rounding in floating point arithmetic. Machine epsilon is the smallest number such
071     * that (1F + EPSF != 1F).
072     */
073    public static final float EPSF = epsilon(1F);
074
075    /**
076     * Square-root of the machine epsilon for <code>float</code>.
077     */
078    public static final float SQRT_EPSF = (float)Math.sqrt(EPSF);
079
080    /**
081     * Prevent anyone from instantiating this utility class.
082     */
083    private MathTools() { }
084
085    //-----------------------------------------------------------------------------------
086    /**
087     * Test to see if a given long integer is even.
088     *
089     * @param n Integer number to be tested.
090     * @return True if the number is even, false if it is odd.
091     */
092    public static boolean even(long n) {
093        return (n & 1) == 0;
094    }
095
096    /**
097     * Test to see if a given long integer is odd.
098     *
099     * @param n Integer number to be tested.
100     * @return True if the number is odd, false if it is even.
101     */
102    public static boolean odd(long n) {
103        return (n & 1) != 0;
104    }
105
106    /**
107     * Calculates the square (x^2) of the argument.
108     *
109     * @param x Argument to be squared.
110     * @return Returns the square (x^2) of the argument.
111     */
112    public static double sqr(double x) {
113        if (x == 0.)
114            return 0.;
115        else
116            return x * x;
117    }
118
119    /**
120     * Computes the cube root of the specified real number. If the argument is negative,
121     * then the cube root is negative.
122     *
123     * @param x Argument for which the cube root is to be found.
124     * @return The cube root of the argument is returned.
125     * @see Math#cbrt(double)
126     * @deprecated 
127     */
128    @Deprecated
129    public static double cubeRoot(double x) {
130//      double value;
131//      
132//      if ( x < 0. )
133//          value = -Math.exp( Math.log(-x) / 3. );
134//      else
135//          value = Math.exp( Math.log( x ) / 3. );
136//          
137//      return value;
138        return Math.cbrt(x);    //  Added with Java 1.5.
139    }
140
141    /**
142     * Returns a number "a" raised to the power "b". A "long" version of Math.pow(). This
143     * is much faster than using Math.pow() if the operands are integers.
144     *
145     * @param a Number to be raised to the power "b".
146     * @param b Power to raise number "a" to.
147     * @return A long integer "a" raised to the integer power "b".
148     * @throws ArithmeticException if "b" is negative.
149     */
150    public static long pow(long a, long b) throws ArithmeticException {
151        if (b < 0)
152            throw new ArithmeticException("Exponent must be positive.");
153
154        long r = 1;
155        while (b != 0) {
156            if (odd(b))
157                r *= a;
158
159            b >>>= 1;
160            a *= a;
161        }
162        return r;
163    }
164
165    /**
166     * Raises 2 to the small integer power indicated (eg: 2^3 = 8). This is MUCH faster
167     * than calling Math.pow(2, x).
168     *
169     * @param x Amount to raise 2 to the power of.
170     * @return Returns 2 raised to the power indicated.
171     */
172    public static long pow2(long x) {
173        long value = 1;
174        for (long i = 0; i < x; ++i)
175            value *= 2;
176        return value;
177    }
178
179    /**
180     * Raises 10 to the small integer power indicated (eg: 10^5 = 100000). This is faster
181     * than calling Math.pow(10, x).
182     *
183     * @param x Amount to raise 10 to the power of.
184     * @return Returns 10 raised to the power indicated.
185     */
186    public static double pow10(int x) {
187        double pow10 = 10.;
188
189        if (x != 0) {
190            boolean neg = false;
191            if (x < 0) {
192                x *= -1;
193                neg = true;
194            }
195
196            for (int i = 1; i < x; ++i)
197                pow10 *= 10.;
198
199            if (neg)
200                pow10 = 1. / pow10;
201
202        } else
203            pow10 = 1.;
204
205        return (pow10);
206    }
207
208    /**
209     * Find the base 10 logarithm of the given double.
210     *
211     * @param x Value to find the base 10 logarithm of.
212     * @return The base 10 logarithm of x.
213     * @see Math#log10(double) 
214     * @deprecated 
215     */
216    @Deprecated
217    public static double log10(double x) {
218        //return Math.log(x)/LOG10;
219        return Math.log10(x);   //  Added with Java 1.5.
220    }
221
222    /**
223     * Find the base 2 logarithm of the given double.
224     *
225     * @param x Value to find the base 2 logarithm of.
226     * @return The base 2 logarithm of x.
227     */
228    public static double log2(double x) {
229        return Math.log(x) / LOG2;
230    }
231
232    /**
233     * Rounds a floating point number to the desired decimal place. Example: 1346.4667
234     * rounded to the 2nd place = 1300.
235     *
236     * @param value The value to be rounded.
237     * @param place Number of decimal places to round value to. A place of 1 rounds to
238     *              10's place, 2 to 100's place, -2 to 1/100th place, et cetera.
239     * @return A floating point number rounded to the desired decimal place.
240     * @see #roundToSigFig(double, int) 
241     * @see #roundDownToPlace(double, int) 
242     * @see #roundUpToPlace(double, int) 
243     */
244    public static double roundToPlace(double value, int place) {
245
246        //  If the value is zero, just pass the number back out.
247        if (value != 0.) {
248
249            //  If the place is zero, round to the one's place.
250            if (place == 0)
251                value = Math.floor(value + 0.5);
252
253            else {
254                double pow10 = MathTools.pow10(place);  //  = 10 ^ place
255                double holdvalue = value / pow10;
256
257                value = Math.floor(holdvalue + 0.5);        // Round number to nearest integer
258                value *= pow10;
259            }
260        }
261
262        return value;
263    }
264
265    /**
266     * Rounds a floating point number up to the desired decimal place. Example: 1346.4667
267     * rounded up to the 2nd place = 1400.
268     *
269     * @param value The value to be rounded up.
270     * @param place Number of decimal places to round value to. A place of 1 rounds to
271     *              10's place, 2 to 100's place, -2 to 1/100th place, et cetera.
272     * @return A floating point number rounded up to the desired decimal place.
273     * @see #roundToPlace(double, int) 
274     * @see #roundDownToPlace(double, int) 
275     */
276    public static double roundUpToPlace(double value, int place) {
277
278        //  If the value is zero, just pass the number back out.
279        if (value != 0.) {
280
281            //  If the place is zero, round to the one's place.
282            if (place == 0)
283                value = Math.ceil(value);
284
285            else {
286                double pow10 = MathTools.pow10(place);  //  = 10 ^ place
287                double holdvalue = value / pow10;
288
289                value = Math.ceil(holdvalue);           // Round number up to nearest integer
290                value *= pow10;
291            }
292        }
293
294        return value;
295    }
296
297    /**
298     * Rounds a floating point number down to the desired decimal place. Example:
299     * 1346.4667 rounded down to the 1st place = 1340.
300     *
301     * @param value The value to be rounded down.
302     * @param place Number of decimal places to round value to. A place of 1 rounds to
303     *              10's place, 2 to 100's place, -2 to 1/100th place, et cetera.
304     * @return A floating point number rounded down to the desired decimal place.
305     * @see #roundToPlace(double, int) 
306     * @see #roundUpToPlace(double, int) 
307     */
308    public static double roundDownToPlace(double value, int place) {
309
310        //  If the value is zero, just pass the number back out.
311        if (value != 0.) {
312
313            //  If the place is zero, round to the one's place.
314            if (place == 0)
315                value = Math.floor(value);
316
317            else {
318                double pow10 = MathTools.pow10(place);  //  = 10 ^ place
319                double holdvalue = value / pow10;
320
321                value = Math.floor(holdvalue);          // Round number down to nearest integer
322                value *= pow10;
323            }
324        }
325
326        return value;
327    }
328
329    /**
330     * Rounds a floating point number to the desired number of significant digits.
331     * For example: 1346.4667 rounded to 3 significant digits is 1350.0.
332     * 
333     * @param value The value to be rounded.
334     * @param sigFig The number of significant digits/figures to retain.
335     * @return A floating point number rounded to the specified number of significant digits.
336     * @see #roundToPlace(double, int) 
337     */
338    public static double roundToSigFig(double value, int sigFig) {
339        int p = (int)Math.round(Math.log10(value) + 0.5);
340        p -= sigFig;
341        return MathTools.roundToPlace(value, p);
342    }
343    
344    /**
345     * Returns the factorial of n (usually written as n!) for any value of n. The
346     * factorial of n is the product of all integers up to and including n.
347     *
348     * @param n The number to calculate the factorial for.
349     * @return The factorial of n.
350     * @see #intFactorial(int)
351     */
352    public static BigInteger factorial(int n) {
353        if (n < 0)
354            throw new IllegalArgumentException("n must be positive");
355        if (n <= 12)
356            return BigInteger.valueOf(intFactorial(n));
357
358        BigInteger fact = BigInteger.ONE;
359        for (int i = 2; i <= n; ++i)
360            fact = fact.multiply(BigInteger.valueOf(i));
361        return fact;
362    }
363
364    /**
365     * Returns the factorial of n (usually written as n!) for values of n up to 12. The
366     * factorial of n is the product of all integers up to and including n. This method is
367     * more efficient than "factorial()" for values of n &le; 12. Values of n > 12 cause a
368     * numerical overflow with this method. Use "factorial()" for larger values of n (or
369     * any value of n where performance is not critical).
370     *
371     * @param n The number to calculate the factorial for (0 &le; n &le; 12).
372     * @return The factorial of n.
373     * @see #factorial(int) 
374    *
375     */
376    public static int intFactorial(int n) {
377        if (n < 0)
378            throw new IllegalArgumentException("n must be positive");
379        if (n > 12)
380            throw new IllegalArgumentException("Integer overflow in factorial()!");
381        int fact = 1;
382        for (int i = 2; i <= n; ++i)
383            fact *= i;
384        return fact;
385    }
386
387    /**
388     * Calculates the greatest common divisor between two input integers. The GCD is the
389     * largest number that can be divided into both input numbers.
390     *
391     * @param xval First integer
392     * @param yval Second integer
393     * @return The largest number that can be divided into both input values.
394     */
395    public static long greatestCommonDivisor(long xval, long yval) {
396        // Uses Euler's method.
397        long value = 0;
398        while (value != xval) {
399            if (xval < yval)
400                yval = yval - xval;
401
402            else {
403                if (xval > yval)
404                    xval = xval - yval;
405                else
406                    value = xval;
407            }
408        }
409        return (value);
410    }
411
412    /**
413     * Returns the fractional part of a floating point number (removes the integer part).
414     *
415     * @param x Argument for which the fractional part is to be returned.
416     * @return The fractional part of the argument is returned.
417     */
418    public static double frac(double x) {
419        x = x - (long)x;
420        if (x < 0.)
421            ++x;
422
423        return x;
424    }
425
426    /**
427     * Straight linear 1D interpolation between two points.
428     *
429     * @param x1 X-coordinate of the 1st point (the high point).
430     * @param y1 Y-coordinate of the 1st point (the high point).
431     * @param x2 X-coordinate of the 2nd point (the low point).
432     * @param y2 Y-coordinate of the 2nd point (the low point).
433     * @param x  The X coordinate of the point for which we want to interpolate to
434     *           determine a Y coordinate. Will extrapolate if X is outside of the bounds
435     *           of the point arguments.
436     * @return The linearly interpolated Y value corresponding to the input X value is
437     *         returned.
438     */
439    public static double lineInterp(double x1, double y1, double x2, double y2, double x) {
440        return ((y2 - y1) / (x2 - x1) * (x - x1) + y1);
441    }
442
443    /**
444     * Converts a positive decimal number to it's binary equivalent.
445     *
446     * @param decimal The positive decimal number to be encoded in binary.
447     * @param bits    The bitset to encode the number in.
448     */
449    public static void dec2bin(int decimal, BitSet bits) {
450        if (decimal < 0)
451            throw new IllegalArgumentException("Cannot convert a negative number to binary.");
452
453        int i = 0;
454        int value = decimal;
455        while (value > 0) {
456
457            if (value % 2 > 0)
458                bits.set(i);
459            else
460                bits.clear(i);
461
462            value /= 2;
463            ++i;
464        }
465
466        for (; i < bits.size(); ++i)
467            bits.clear(i);
468    }
469
470    /**
471     * Converts binary number to it's base 10 decimal equivalent.
472     *
473     * @param bits The bitset that encodes the number to be converted.
474     * @return Returns the decimal equivalent of the given binary number.
475     */
476    public static long bin2dec(BitSet bits) {
477        long value = 0;
478        int length = bits.size();
479
480        for (int i = 0; i < length; ++i) {
481            if (bits.get(i))
482                value += pow2(i);
483        }
484        return value;
485    }
486
487    /**
488     * Return the hyperbolic cosine of the specified argument.
489     *
490     * @param x The number whose hyperbolic cosine is to be returned.
491     * @return The hyperbolic cosine of x.
492     * @see Math#cosh(double) 
493     * @deprecated 
494     */
495    @Deprecated
496    public static double cosh(double x) {
497        return Math.cosh(x);
498    }
499
500    /**
501     * Return the hyperbolic sine of the specified argument.
502     *
503     * @param x The number whose hyperbolic sine is to be returned.
504     * @return The hyperbolic sine of x.
505     * @see Math#sinh(double) 
506     * @deprecated 
507     */
508    @Deprecated
509    public static double sinh(double x) {
510        return Math.sinh(x);
511    }
512
513    /**
514     * Returns the hyperbolic tangent of the specified argument.
515     *
516     * @param x The number whose hyperbolic tangent is to be returned.
517     * @return The hyperbolic tangent of x.
518     * @see Math#tanh(double) 
519     * @deprecated 
520     */
521    @Deprecated
522    public static double tanh(double x) {
523        return Math.tanh(x);
524    }
525
526    /**
527     * Returns the inverse hyperbolic cosine of the specified argument. The inverse
528     * hyperbolic cosine is defined as: <code>acosh(x) = log(x + sqrt((x-1)*(x+1)))</code>
529     *
530     * @param x Value to return inverse hyperbolic cosine of.
531     * @return The inverse hyperbolic cosine of x.
532     * @throws IllegalArgumentException if x is less than 1.0.
533     */
534    public static double acosh(double x) {
535        if (Double.isNaN(x))
536            return Double.NaN;
537        if (Double.isInfinite(x))
538            return x;
539        if (x < 1.0)
540            throw new IllegalArgumentException("x may not be less than 1.0");
541
542        double y;
543        if (x > 1.0E8) {
544            y = Math.log(x) + LOG2;
545
546        } else {
547            double a = Math.sqrt((x - 1.0) * (x + 1.0));
548            y = Math.log(x + a);
549        }
550
551        return y;
552    }
553
554    /**
555     * Returns the inverse hyperbolic sine of the specified argument. The inverse
556     * hyperbolic sine is defined as: <code>asinh(x) = log(x + sqrt(1 + x*x))</code>
557     *
558     * @param xx Value to return inverse hyperbolic cosine of.
559     * @return The inverse hyperbolic sine of x.
560     */
561    public static double asinh(double xx) {
562        if (Double.isNaN(xx))
563            return Double.NaN;
564        if (Double.isInfinite(xx))
565            return xx;
566        if (xx == 0)
567            return 0;
568
569        int sign = 1;
570        double x = xx;
571        if (xx < 0) {
572            sign = -1;
573            x = -xx;
574        }
575
576        double y;
577        if (x > 1.0E8) {
578            y = sign * (Math.log(x) + LOG2);
579
580        } else {
581            double a = Math.sqrt(x * x + 1.0);
582            y = sign * Math.log(x + a);
583        }
584
585        return y;
586    }
587
588    /**
589     * Returns the inverse hyperbolic tangent of the specified argument. The inverse
590     * hyperbolic tangent is defined as: <code>atanh(x) = 0.5*log((1 + x)/(1 - x))</code>
591     *
592     * @param x Value to return inverse hyperbolic cosine of.
593     * @return The inverse hyperbolic tangent of x.
594     * @throws IllegalArgumentException if x is outside the range -1, to +1.
595     */
596    public static double atanh(double x) {
597        if (Double.isNaN(x))
598            return Double.NaN;
599        if (x == 0)
600            return 0;
601
602        double z = Math.abs(x);
603        if (z >= 1.0) {
604            if (x == 1.0)
605                return Double.POSITIVE_INFINITY;
606            if (x == -1.0)
607                return Double.NEGATIVE_INFINITY;
608
609            throw new IllegalArgumentException("x outside of range -1 to +1");
610        }
611
612        if (z < 1.0E-7)
613            return x;
614
615        double y = 0.5 * Math.log((1.0 + x) / (1.0 - x));
616
617        return y;
618    }
619
620    /**
621     * Returns the absolute value of "a" times the sign of "b".
622     *
623     * @param a The value for which the magnitude is returned.
624     * @param b The value for which the sign is returned.
625     * @return The absolute value of "a" times the sign of "b".
626     */
627    public static double sign(double a, double b) {
628        return b >= 0 ? (a >= 0 ? a : -a) : (a >= 0 ? -a : a);
629    }
630
631    /**
632     * Returns the absolute value of "a" times the sign of "b".
633     *
634     * @param a The value for which the magnitude is returned.
635     * @param b The value for which the sign is returned.
636     * @return The absolute value of "a" times the sign of "b".
637     */
638    public static float sign(float a, double b) {
639        return b >= 0 ? (a >= 0 ? a : -a) : (a >= 0 ? -a : a);
640    }
641
642    /**
643     * Returns the absolute value of "a" times the sign of "b".
644     *
645     * @param a The value for which the magnitude is returned.
646     * @param b The value for which the sign is returned.
647     * @return The absolute value of "a" times the sign of "b".
648     */
649    public static long sign(long a, double b) {
650        return b >= 0 ? (a >= 0 ? a : -a) : (a >= 0 ? -a : a);
651    }
652
653    /**
654     * Returns the absolute value of "a" times the sign of "b".
655     *
656     * @param a The value for which the magnitude is returned.
657     * @param b The value for which the sign is returned.
658     * @return The absolute value of "a" times the sign of "b".
659     */
660    public static int sign(int a, double b) {
661        return b >= 0 ? (a >= 0 ? a : -a) : (a >= 0 ? -a : a);
662    }
663
664    //  Used by "lnGamma()".
665    private static final double[] gamCoef = {76.18009172947146, -86.50532032941677, 24.01409824083091,
666        -1.231739572450155, 0.1208650973866179e-2, -0.5395239384953e-5};
667
668    /**
669     * Returns the natural log of the Gamma Function defined by the integral:
670     * <code>Gamma(z) = integral from 0 to infinity of t^(z-1)*e^-t dt</code>.
671     * <p>
672     * It is better to implement ln(Gamma(x)) rather than Gamma(x) since the latter will
673     * overflow most computer floating point representations at quite modest values of x.
674     * </p>
675     *
676     * @param value The value to evaluate the log of the gamma function for.
677     * @return The natural log of the Gamma Function.
678     */
679    public static double lnGamma(double value) {
680        double x = value, y = value;
681        double tmp = x + 5.5;
682        tmp -= (x + 0.5) * Math.log(tmp);
683        double ser = 1.000000000190015;
684        for (int j = 0; j <= 5; ++j)
685            ser += gamCoef[j] / ++y;
686        return -tmp + Math.log(2.5066282746310005 * ser / x);
687    }
688
689    /**
690     * Returns the smallest roundoff in quantities of size x, EPSILON, such that
691     * <code>x + EPSILON > x</code>. This is the Units in the Last Place (ULP) and is
692     * different than the machine roundoff EPS (macheps).
693     *
694     * @param x The value to return the roundoff quantity for.
695     * @return The smallest roundoff in quantities of size x, EPSILON, such that x +
696     *         EPSILON > x
697     * @see #EPS
698     * @see Math#ulp(double) 
699     */
700    public static double epsilon(double x) {
701        return Math.ulp(x);
702    }
703
704    /**
705     * Returns the smallest roundoff in quantities of size x, EPSILON, such that
706     * <code>x + EPSILON > x</code>. This is the Units in the Last Place (ULP) and is
707     * different than the machine roundoff EPS (macheps).
708     *
709     * @param x The value to return the roundoff quantity for.
710     * @return The smallest roundoff in quantities of size x, EPSILON, such that x +
711     *         EPSILON > x
712     * @see #EPSF
713     * @see Math#ulp(float) 
714     */
715    public static float epsilon(float x) {
716        return Math.ulp(x);
717    }
718
719    /**
720     * Returns true if the two supplied numbers are approximately equal to within machine
721     * precision.
722     *
723     * @param a The 1st value to compare for approximate equality.
724     * @param b The 2nd value to compare for approximate equality.
725     * @return true if the two supplied numbers are approximately equal to within machine
726     *         precision
727     */
728    public static boolean isApproxEqual(double a, double b) {
729        double eps2 = epsilon(a);
730        double eps = (eps2 > EPS ? eps2 : EPS);
731        return Math.abs(a - b) <= eps;
732    }
733
734    /**
735     * Returns true if the two supplied numbers are approximately equal to within the
736     * specified tolerance.
737     *
738     * @param a   The 1st value to compare for approximate equality.
739     * @param b   The 2nd value to compare for approximate equality.
740     * @param tol The tolerance for equality.
741     * @return true if the two supplied numbers are approximately equal to within the
742     *         specified tolerance
743     */
744    public static boolean isApproxEqual(double a, double b, double tol) {
745        return Math.abs(a - b) <= tol;
746    }
747
748    /**
749     * Returns true if the two supplied numbers are approximately equal to within machine
750     * precision.
751     *
752     * @param a The 1st value to compare for approximate equality.
753     * @param b The 2nd value to compare for approximate equality.
754     * @return true if the two supplied numbers are approximately equal to within machine
755     *         precision
756     */
757    public static boolean isApproxEqual(float a, float b) {
758        float eps2 = epsilon(a);
759        float eps = (eps2 > EPSF ? eps2 : EPSF);
760        return Math.abs(a - b) <= eps;
761    }
762
763    /**
764     * Returns true if the two supplied numbers are approximately equal to within the
765     * specified tolerance.
766     *
767     * @param a   The 1st value to compare for approximate equality.
768     * @param b   The 2nd value to compare for approximate equality.
769     * @param tol The tolerance for equality.
770     * @return true if the two supplied numbers are approximately equal to within the
771     *         specified tolerance
772     */
773    public static boolean isApproxEqual(float a, float b, float tol) {
774        return Math.abs(a - b) <= tol;
775    }
776
777    /**
778     * Returns true if the supplied number is approximately zero to within machine
779     * precision.
780     *
781     * @param a The value to be compared with zero.
782     * @return true if the supplied number is approximately zero to within machine
783     *         precision
784     */
785    public static boolean isApproxZero(double a) {
786        return Math.abs(a) <= EPS;
787    }
788
789    /**
790     * Returns true if the supplied number is approximately zero to within the specified
791     * tolerance.
792     *
793     * @param a   The value to be compared with zero.
794     * @param tol The tolerance for equality.
795     * @return true if the supplied number is approximately zero to within the specified
796     *         tolerance
797     */
798    public static boolean isApproxZero(double a, double tol) {
799        return Math.abs(a) <= tol;
800    }
801
802    /**
803     * Returns true if the supplied number is approximately zero to within machine
804     * precision.
805     *
806     * @param a The value to be compared with zero.
807     * @return true if the supplied number is approximately zero to within machine
808     *         precision
809     */
810    public static boolean isApproxZero(float a) {
811        return Math.abs(a) <= EPSF;
812    }
813
814    /**
815     * Returns true if the supplied number is approximately zero to within the specified
816     * tolerance.
817     *
818     * @param a   The value to be compared with zero.
819     * @param tol The tolerance for equality.
820     * @return true if the supplied number is approximately zero to within the specified
821     *         tolerance
822     */
823    public static boolean isApproxZero(float a, float tol) {
824        return Math.abs(a) <= tol;
825    }
826
827    /**
828     * Used to test out the methods in this class.
829     *
830     * @param args Command line arguments (not used).
831     */
832    public static void main(String args[]) {
833
834        System.out.println();
835        System.out.println("Testing MathTools...");
836
837        System.out.println("  2 is an " + (even(2) ? "even" : "odd") + " number.");
838        System.out.println("  3 is an " + (odd(3) ? "odd" : "even") + " number.");
839        System.out.println("  The square of 3.8 is " + sqr(3.8) + ".");
840        System.out.println("  The integer 3^7 is " + pow(3, 7) + ".");
841        System.out.println("  The integer 2^8 is " + pow2(8) + ".");
842        System.out.println("  The double 10^-3 is " + pow10(-3) + ".");
843        System.out.println("  The base 2 logarithm of 8 is " + log2(8) + ".");
844        System.out.println("  1346.4667 rounded to the nearest 100 is "
845                + roundToPlace(1346.4667, 2) + ".");
846        System.out.println("  1346.4667 rounded up to the nearest 100 is "
847                + roundUpToPlace(1346.4667, 2) + ".");
848        System.out.println("  1346.4667 rounded down to the nearest 10 is "
849                + roundDownToPlace(1346.4667, 1) + ".");
850        System.out.println("  1346.4667 rounded to 3 significant digits is "
851                + roundToSigFig(1346.4667, 3) + ".");
852        System.out.println("  The GCD of 125 and 45 is " + greatestCommonDivisor(125, 45) + ".");
853        System.out.println("  The fractional part of 3.141593 is " + frac(3.141593) + ".");
854        double x = Math.cosh(5);
855        System.out.println("  The inv. hyperbolic cosine of " + (float)x + " = " + (float)acosh(x) + ".");
856        System.out.println("  The inv. hyperbolic sine of " + (float)x + " = " + (float)asinh(x) + ".");
857        x = Math.tanh(-0.25);
858        System.out.println("  The inv. hyperbolic tangent of " + (float)x + " = " + (float)atanh(x) + ".");
859
860        System.out.println("  4.56 with the sign of -6.33 is " + sign(4.56F, -6.33));
861        System.out.println("  epsilon(0) is " + epsilon(0));
862        System.out.println("  epsilon(2387.8) is " + epsilon(2387.8));
863
864    }
865
866}