Package gov.nih.mipav.util

Class DoubleDouble

java.lang.Object
gov.nih.mipav.util.DoubleDouble
All Implemented Interfaces:
Serializable, Cloneable, Comparable
public class DoubleDouble extends Object implements Serializable, Comparable, Cloneable
Immutable, extended-precision floating-point numbers which maintain 106 bits (approximately 30 decimal digits) of precision.

A DoubleDouble uses a representation containing two double-precision values. A number x is represented as a pair of doubles, x.hi and x.lo, such that the number represented by x is x.hi + x.lo, where 

    |x.lo| invalid input: '<'= 0.5*ulp(x.hi)
and ulp(y) means "unit in the last place of y". The basic arithmetic operations are implemented using convenient properties of IEEE-754 floating-point arithmetic.

The range of values which can be represented is the same as in IEEE-754. The precision of the representable numbers is twice as great as IEEE-754 double precision.

The correctness of the arithmetic algorithms relies on operations being performed with standard IEEE-754 double precision and rounding. This is the Java standard arithmetic model, but for performance reasons Java implementations are not constrained to using this standard by default. Some processors (notably the Intel Pentium architecure) perform floating point operations in (non-IEEE-754-standard) extended-precision. A JVM implementation may choose to use the non-standard extended-precision as its default arithmetic mode. To prevent this from happening, this code uses the Java strictfp modifier, which forces all operations to take place in the standard IEEE-754 rounding model.

The API provides a value-oriented interface. DoubleDouble values are immutable; operations on them return new objects carrying the result of the operation. This provides a much simpler semantics for writing DoubleDouble expressions, and Java memory management is efficient enough that this imposes very little performance penalty.

This implementation uses algorithms originally designed variously by Knuth, Kahan, Dekker, and Linnainmaa. Douglas Priest developed the first C implementation of these techniques. Other more recent C++ implementation are due to Keith M. Briggs and David Bailey et al.

References

  • Priest, D., Algorithms for Arbitrary Precision Floating Point Arithmetic, in P. Kornerup and D. Matula, Eds., Proc. 10th Symposium on Computer Arithmetic, IEEE Computer Society Press, Los Alamitos, Calif., 1991.
  • Yozo Hida, Xiaoye S. Li and David H. Bailey, Quad-Double Arithmetic: Algorithms, Implementation, and Application, manuscript, Oct 2000; Lawrence Berkeley National Laboratory Report BNL-46996.
  • David Bailey, High Precision Software Directory; http://crd.lbl.gov/~dhbailey/mpdist/index.html
This class is based on the DoubleDouble created by Martin Davis and distributed under the following license: 
 Copyright (c) 2008-2012 Martin Davis. All Rights Reserved.
 
 Redistribution and use in source and binary forms, with or without modification, are permitted provided that the
 following conditions are met:
 
 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following
 disclaimer.
 
 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the
 following disclaimer in the documentation and/or other materials provided with the distribution.
 
 3. The name of the author may not be used to endorse or promote products derived from this software without specific
 prior written permission.
 
 THIS SOFTWARE IS PROVIDED BY Martin Davis "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
 TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL
 THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
 BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
 INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 POSSIBILITY OF SUCH DAMAGE.
Author:
Martin Davis
See Also:

  • Field Details

    • PI

      public static final DoubleDouble PI
      The value nearest to the constant Pi.

    • TWO_PI

      public static final DoubleDouble TWO_PI
      The value nearest to the constant 2 * Pi.

    • PI_2

      public static final DoubleDouble PI_2
      The value nearest to the constant Pi / 2.

    • E

      public static final DoubleDouble E
      The value nearest to the constant e (the natural logarithm base).

    • NaN

      public static final DoubleDouble NaN
      A value representing the result of an operation which does not return a valid number.

    • POSITIVE_INFINITY

      public static final DoubleDouble POSITIVE_INFINITY

    • NEGATIVE_INFINITY

      public static final DoubleDouble NEGATIVE_INFINITY

    • EPS

      public static final double EPS
      The smallest representable relative difference between two {link @ DoubleDouble} values
      See Also:

    • SPLIT

      private static final double SPLIT
      The value to split a double-precision value on during multiplication
      See Also:

    • hi

      public double hi
      The high-order component of the double-double precision value.

    • lo

      public double lo
      The low-order component of the double-double precision value.

    • MAX_PRINT_DIGITS

      private static final int MAX_PRINT_DIGITS
      See Also:

    • TEN

      private static final DoubleDouble TEN

    • ONE

      private static final DoubleDouble ONE

    • SCI_NOT_EXPONENT_CHAR

      private static final String SCI_NOT_EXPONENT_CHAR
      See Also:

    • SCI_NOT_ZERO

      private static final String SCI_NOT_ZERO
      See Also:

  • Constructor Details

    • DoubleDouble

      public DoubleDouble()
      Creates a new DoubleDouble with value 0.0.

    • DoubleDouble

      public DoubleDouble(double x)
      Creates a new DoubleDouble with value x.
      Parameters:
      x - the value to initialize

    • DoubleDouble

      public DoubleDouble(double hi, double lo)
      Creates a new DoubleDouble with value (hi, lo).
      Parameters:
      hi - the high-order component
      lo - the high-order component

    • DoubleDouble

      public DoubleDouble(DoubleDouble dd)
      Creates a new DoubleDouble with value equal to the argument.
      Parameters:
      dd - the value to initialize

    • DoubleDouble

      public DoubleDouble(String str) throws NumberFormatException
      Creates a new DoubleDouble with value equal to the argument.
      Parameters:
      str - the value to initialize by
      Throws:
      NumberFormatException - if str is not a valid representation of a number

  • Method Details

    • valueOf

      public static DoubleDouble valueOf(String str) throws NumberFormatException
      Converts the string argument to a DoubleDouble number.
      Parameters:
      str - a string containing a representation of a numeric value
      Returns:
      the extended precision version of the value
      Throws:
      NumberFormatException - if s is not a valid representation of a number

    • valueOf

      public static DoubleDouble valueOf(double x)
      Converts the double argument to a DoubleDouble number.
      Parameters:
      x - a numeric value
      Returns:
      the extended precision version of the value

    • clone

      public Object clone()
      Creates and returns a copy of this value.
      Overrides:
      clone in class Object
      Returns:
      a copy of this value

    • init

      private void init(double x)

    • init

      private void init(double hi, double lo)

    • init

      private void init(DoubleDouble dd)

    • add

      public DoubleDouble add(DoubleDouble y)
      Returns a DoubleDouble whose value is (this + y).
      Parameters:
      y - the addend
      Returns:
      (this + y)

    • selfAdd

      private DoubleDouble selfAdd(DoubleDouble y)
      Adds the argument to the value of this. To prevent altering constants, this method must only be used on values known to be newly created.
      Parameters:
      y - the addend
      Returns:
      this, with its value incremented by y

    • subtract

      public DoubleDouble subtract(DoubleDouble y)
      Returns a DoubleDouble whose value is (this - y).
      Parameters:
      y - the subtrahend
      Returns:
      (this - y)

    • negate

      public DoubleDouble negate()
      Returns a DoubleDouble whose value is -this.
      Returns:
      -this

    • multiply

      public DoubleDouble multiply(DoubleDouble y)
      Returns a DoubleDouble whose value is (this * y).
      Parameters:
      y - the multiplicand
      Returns:
      (this * y)

    • selfMultiply

      private DoubleDouble selfMultiply(DoubleDouble y)
      Multiplies this by the argument, returning this. To prevent altering constants, this method must only be used on values known to be newly created.
      Parameters:
      y - a DoubleDouble value to multiply by
      Returns:
      this

    • divide

      public DoubleDouble divide(DoubleDouble y)
      Returns a DoubleDouble whose value is (this / y).
      Parameters:
      y - the divisor
      Returns:
      (this / y)

    • reciprocal

      public DoubleDouble reciprocal()
      Returns a DoubleDouble whose value is 1 / this.
      Returns:
      the reciprocal of this value

    • floor

      public DoubleDouble floor()
      Returns the largest (closest to positive infinity) value that is not greater than the argument and is equal to a mathematical integer. Special cases:
      • If this value is NaN, returns NaN.
      Returns:
      the largest (closest to positive infinity) value that is not greater than the argument and is equal to a mathematical integer.

    • ceil

      public DoubleDouble ceil()
      Returns the smallest (closest to negative infinity) value that is not less than the argument and is equal to a mathematical integer. Special cases:
      • If this value is NaN, returns NaN.
      Returns:
      the smallest (closest to negative infinity) value that is not less than the argument and is equal to a mathematical integer.

    • signum

      public int signum()
      Returns an integer indicating the sign of this value.
      • if this value is > 0, returns 1
      • if this value is invalid input: '<' 0, returns -1
      • if this value is = 0, returns 0
      • if this value is NaN, returns 0
      Returns:
      an integer indicating the sign of this value

    • rint

      public DoubleDouble rint()
      Rounds this value to the nearest integer. The value is rounded to an integer by adding 1/2 and taking the floor of the result. Special cases:
      • If this value is NaN, returns NaN.
      Returns:
      this value rounded to the nearest integer

    • trunc

      public DoubleDouble trunc()
      Returns the integer which is largest in absolute value and not further from zero than this value. Special cases:
      • If this value is NaN, returns NaN.
      Returns:
      the integer which is largest in absolute value and not further from zero than this value

    • abs

      public DoubleDouble abs()
      Returns the absolute value of this value. Special cases:
      • If this value is NaN, it is returned.
      Returns:
      the absolute value of this value

    • sqr

      public DoubleDouble sqr()
      Computes the square of this value.
      Returns:
      the square of this value.

    • sqrt

      public DoubleDouble sqrt()
      Computes the positive square root of this value. If the number is NaN or negative, NaN is returned.
      Returns:
      the positive square root of this number. If the argument is NaN or less than zero, the result is NaN.

    • exp

      public DoubleDouble exp()
      For all real, x exp(x) = 1 + x + x**2/2! + x**3/3! + x**4/4! + ...
      Returns:

    • log

      public DoubleDouble log()
      For x > 0, ratio = (x-1)/(x+1), log(x) = 2(ratio + ratio**3/3 + ratio**5/5 + ...)
      Returns:

    • sinh

      public DoubleDouble sinh()
      For all real x, sinh(x) = x + x**3/3! + x**5/5! + x**7/7! + ... + x**(2n+1)/(2n+1)! + ...
      Returns:

    • cosh

      public DoubleDouble cosh()
      For all real x, cosh(x) = 1 + x**2/2! + x**4/4! + x**6/6! + ... + x**(2*n)/((2*n)!) + ...
      Returns:

    • sin

      public DoubleDouble sin()
      For all real x, sin(x) = x - x**3/3! + x**5/5! - x**7/7! + ...
      Returns:

    • erf

      public DoubleDouble erf()

    • cos

      public DoubleDouble cos()
      For all real x, cos(x) = 1 - x**2/2! + x**4/4! - x**6/6! + ...
      Returns:

    • tan

      public DoubleDouble tan()
      For -PI/2 invalid input: '<' x invalid input: '<' PI/2, tan(x) = x + (x**3)/3 + 2*(x**5)/15 + 17*(x**7)/315 + 62*(x**9)/2835 + ... + (2**(2*n))*((2**(2*n)) - 1)*Bn*(x**(2*n-1))/((2*n)!) + ...
      Returns:

    • asin

      public DoubleDouble asin()
      For all -1 invalid input: '<' x invalid input: '<' 1, arcsin(x) = x + x**3/(2*3) + (1 * 3 * x**5)/(2 * 4 * 5) + (1 * 3 * 5 * x**7)/(2 * 4 * 6 * 7) + ...
      Returns:

    • acos

      public DoubleDouble acos()
      For all -1 invalid input: '<' x invalid input: '<' 1, arccos(x) = PI/2 - arcsin(x)
      Returns:

    • atan2

      public DoubleDouble atan2(DoubleDouble x)
      The atan of this(the imaginary y component) divided by the real component x Value will range from -PI to PI. Strays a bit from actual atan2 definition which recognizes positive zero and negative zero. Here both positive zero and negative zero are combined as zero.
      Parameters:
      x -
      Returns:

    • atan

      public DoubleDouble atan()
      For -1 invalid input: '<' x invalid input: '<' 1, arctan(x) = x - x**3/3 + x**5/5 - x**7/7 + ... For x > 1, arctan(x) = PI/2 - 1/x + 1/(3*x**3) - 1/(5*x**5) +1/(7*x**7) - ... * For x invalid input: '<' -1, arctan(x) = -PI/2 - 1/x + 1/(3*x**3) - 1/(5*x**5) +1/(7*x**7) - ...
      Returns:

    • BernoulliA

      public DoubleDouble BernoulliA(int n)

    • BernoulliB

      public DoubleDouble BernoulliB(int n)

    • Ci

      public DoubleDouble Ci()

    • Si

      public DoubleDouble Si()

    • cisia

      public void cisia(DoubleDouble x, DoubleDouble Ci, DoubleDouble Si)
      This is a port of subroutine CISIA which computes cosine and sine integrals from Computation of Special Functions by Shanjie Zhang and Jianming Jin. Waiting for Professor Jin's reply. Dear Professor Jianming Jin: There is an error in subroutine CISIA in Computation of Special Functions. Under ELSE IF (x .LE. 32.0D0) THEN in the DO 25 K =M,1,-1 loop the values of BJ(1) thru BJ(M) are set. You then have the loop DO 30 K=2,M,2 30 XS=XS+2.0D0*BJ(K+1) so for M even the value of BJ(M+1) will be used, but the value of BJ(M+1) has not been set. Sincerely, William Gandler

    • factorial

      public DoubleDouble factorial(int fac)
      Parameters:
      fac -
      Returns:

    • pow

      public DoubleDouble pow(int exp)
      Computes the value of this number raised to an integral power. Follows semantics of Java Math.pow as closely as possible.
      Parameters:
      exp - the integer exponent
      Returns:
      x raised to the integral power exp

    • pow

      public DoubleDouble pow(double x)
      Parameters:
      x - the double exponent
      Returns:
      a raised to the double power x For a > 0, base = x * log(a), a**x = 1 + base + base**2/2! + base**3/3! + ...

    • mod

      public DoubleDouble mod(DoubleDouble x)

    • pow

      public DoubleDouble pow(DoubleDouble x)
      Parameters:
      x - the DoubleDouble exponent
      Returns:
      a raised to the DoubleDouble power x For a > 0, base = x * log(a), a**x = 1 + base + base**2/2! + base**3/3! + ...

    • min

      public DoubleDouble min(DoubleDouble x)

    • max

      public DoubleDouble max(DoubleDouble x)

    • doubleValue

      public double doubleValue()
      Converts this value to the nearest double-precision number.
      Returns:
      the nearest double-precision number to this value

    • intValue

      public int intValue()
      Converts this value to the nearest integer.
      Returns:
      the nearest integer to this value

    • isPositiveFinite

      public boolean isPositiveFinite()

    • isNegativeFinite

      public boolean isNegativeFinite()

    • isPositiveInfinity

      public boolean isPositiveInfinity()

    • isNegativeInfinity

      public boolean isNegativeInfinity()

    • isZero

      public boolean isZero()
      Tests whether this value is equal to 0.
      Returns:
      true if this value is equal to 0

    • isNegative

      public boolean isNegative()
      Tests whether this value is less than 0.
      Returns:
      true if this value is less than 0

    • isPositive

      public boolean isPositive()
      Tests whether this value is greater than 0.
      Returns:
      true if this value is greater than 0

    • isNaN

      public boolean isNaN()
      Tests whether this value is NaN.
      Returns:
      true if this value is NaN

    • isInfinite

      public boolean isInfinite()

    • equals

      public boolean equals(DoubleDouble y)
      Tests whether this value is equal to another DoubleDouble value.
      Parameters:
      y - a DoubleDouble value
      Returns:
      true if this value = y

    • ne

      public boolean ne(DoubleDouble y)

    • gt

      public boolean gt(DoubleDouble y)
      Tests whether this value is greater than another DoubleDouble value.
      Parameters:
      y - a DoubleDouble value
      Returns:
      true if this value > y

    • ge

      public boolean ge(DoubleDouble y)
      Tests whether this value is greater than or equals to another DoubleDouble value.
      Parameters:
      y - a DoubleDouble value
      Returns:
      true if this value >= y

    • lt

      public boolean lt(DoubleDouble y)
      Tests whether this value is less than another DoubleDouble value.
      Parameters:
      y - a DoubleDouble value
      Returns:
      true if this value invalid input: '<' y

    • le

      public boolean le(DoubleDouble y)
      Tests whether this value is less than or equal to another DoubleDouble value.
      Parameters:
      y - a DoubleDouble value
      Returns:
      true if this value invalid input: '<'= y

    • compareTo

      public int compareTo(Object o)
      Compares two DoubleDouble objects numerically.
      Specified by:
      compareTo in interface Comparable
      Returns:
      -1,0 or 1 depending on whether this value is less than, equal to or greater than the value of o

    • dump

      public String dump()
      Dumps the components of this number to a string.
      Returns:
      a string showing the components of the number

    • toString

      public String toString()
      Returns a string representation of this number, in either standard or scientific notation. If the magnitude of the number is in the range [ 10-3, 108 ] standard notation will be used. Otherwise, scientific notation will be used.
      Overrides:
      toString in class Object
      Returns:
      a string representation of this number

    • toStandardNotation

      public String toStandardNotation()
      Returns the string representation of this value in standard notation.
      Returns:
      the string representation in standard notation

    • toSciNotation

      public String toSciNotation()
      Returns the string representation of this value in scientific notation.
      Returns:
      the string representation in scientific notation

    • extractSignificantDigits

      private String extractSignificantDigits(boolean insertDecimalPoint, int[] magnitude)
      Extracts the significant digits in the decimal representation of the argument. A decimal point may be optionally inserted in the string of digits (as long as its position lies within the extracted digits - if not, the caller must prepend or append the appropriate zeroes and decimal point).
      Parameters:
      y - the number to extract ( >= 0)
      decimalPointPos - the position in which to insert a decimal point
      Returns:
      the string containing the significant digits and possibly a decimal point

    • stringOfChar

      private static String stringOfChar(char ch, int len)
      Creates a string of a given length containing the given character
      Parameters:
      ch - the character to be repeated
      len - the len of the desired string
      Returns:
      the string

    • getSpecialNumberString

      private String getSpecialNumberString()
      Returns the string for this value if it has a known representation. (E.g. NaN or 0.0)
      Returns:
      the string for this special number

    • magnitude

      private static int magnitude(double x)
      Determines the decimal magnitude of a number. The magnitude is the exponent of the greatest power of 10 which is less than or equal to the number.
      Parameters:
      x - the number to find the magnitude of
      Returns:
      the decimal magnitude of x

    • parse

      public static DoubleDouble parse(String str) throws NumberFormatException
      Converts a string representation of a real number into a DoubleDouble value. The format accepted is similar to the standard Java real number syntax. It is defined by the following regular expression: 
       [+|-] {digit} [ . {digit} ] [ ( e | E ) [+|-] {digit}+
 
      

      

      Parameters:
      

      str - the string to parse
      

      Returns:
      

      the value of the parsed number
      

      Throws:
      

      NumberFormatException - if str is not a valid representation of a number
      

      

      

      • 

    

    

    •