You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by lu...@apache.org on 2010/09/05 00:59:21 UTC

svn commit: r992697 [1/4] - in /commons/proper/math/trunk/src: main/java/org/apache/commons/math/dfp/ site/xdoc/ test/java/org/apache/commons/math/dfp/

Author: luc
Date: Sat Sep  4 22:59:21 2010
New Revision: 992697

URL: http://svn.apache.org/viewvc?rev=992697&view=rev
Log:
added the dfp library
The Dfp class represent the high precision numbers, it implements our existing Field interface and hence each instance is associated with a DfpField that provides the constants at the required precision as well as factory methods. This allowed to remove the compile-time constraint in the library. Users can even use at the same time a field for 20 decimals digits precision and another field for 100 digits precision. Dfp instances with different precision CANNOT be mixed in the same computation (doing so creates a NaN). A few utility methods have been added, like constructors from integral types, isInfinite and isNaN methods, equal has been renames to equals and its signature changed to match the general Object method (a hashcode method has been added too).
JIRA: MATH-412

Added:
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/Dfp.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/DfpDec.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/DfpField.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/DfpMath.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/package.html   (with props)
    commons/proper/math/trunk/src/test/java/org/apache/commons/math/dfp/
    commons/proper/math/trunk/src/test/java/org/apache/commons/math/dfp/Decimal10.java   (with props)
    commons/proper/math/trunk/src/test/java/org/apache/commons/math/dfp/DfpDecTest.java   (with props)
    commons/proper/math/trunk/src/test/java/org/apache/commons/math/dfp/DfpMathTest.java   (with props)
    commons/proper/math/trunk/src/test/java/org/apache/commons/math/dfp/DfpTest.java   (with props)
Modified:
    commons/proper/math/trunk/src/site/xdoc/changes.xml

Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/Dfp.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/Dfp.java?rev=992697&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/Dfp.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/Dfp.java Sat Sep  4 22:59:21 2010
@@ -0,0 +1,2400 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *      http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math.dfp;
+
+import java.util.Arrays;
+
+import org.apache.commons.math.FieldElement;
+
+/**
+ *  Decimal floating point library for Java
+ *
+ *  <p>Another floating point class.  This one is built using radix 10000
+ *  which is 10<sup>4</sup>, so its almost decimal.</p>
+ *
+ *  <p>The design goals here are:
+ *  <ol>
+ *    <li>Decimal math, or close to it</li>
+ *    <li>Settable precision (but no mix between numbers using different settings)</li>
+ *    <li>Portability.  Code should be keep as portable as possible.</li>
+ *    <li>Performance</li>
+ *    <li>Accuracy  - Results should always be +/- 1 ULP for basic
+ *         algebraic operation</li>
+ *    <li>Comply with IEEE 854-1987 as much as possible.
+ *         (See IEEE 854-1987 notes below)</li>
+ *  </ol></p>
+ *
+ *  <p>Trade offs:
+ *  <ol>
+ *    <li>Memory foot print.  I'm using more memory than necessary to
+ *         represent numbers to get better performance.</li>
+ *    <li>Digits are bigger, so rounding is a greater loss.  So, if you
+ *         really need 12 decimal digits, better use 4 base 10000 digits
+ *         there can be one partially filled.</li>
+ *  </ol></p>
+ *
+ *  <p>Numbers are represented  in the following form:
+ *  <pre>
+ *  n  =  sign &times; mant &times; (radix)<sup>exp</sup>;</p>
+ *  </pre>
+ *  where sign is &plusmn;1, mantissa represents a fractional number between
+ *  zero and one.  mant[0] is the least significant digit.
+ *  exp is in the range of -32767 to 32768</p>
+ *
+ *  <p>IEEE 854-1987  Notes and differences</p>
+ *
+ *  <p>IEEE 854 requires the radix to be either 2 or 10.  The radix here is
+ *  10000, so that requirement is not met, but  it is possible that a
+ *  subclassed can be made to make it behave as a radix 10
+ *  number.  It is my opinion that if it looks and behaves as a radix
+ *  10 number then it is one and that requirement would be met.</p>
+ *
+ *  <p>The radix of 10000 was chosen because it should be faster to operate
+ *  on 4 decimal digits at once instead of one at a time.  Radix 10 behavior
+ *  can be realized by add an additional rounding step to ensure that
+ *  the number of decimal digits represented is constant.</p>
+ *
+ *  <p>The IEEE standard specifically leaves out internal data encoding,
+ *  so it is reasonable to conclude that such a subclass of this radix
+ *  10000 system is merely an encoding of a radix 10 system.</p>
+ *
+ *  <p>IEEE 854 also specifies the existence of "sub-normal" numbers.  This
+ *  class does not contain any such entities.  The most significant radix
+ *  10000 digit is always non-zero.  Instead, we support "gradual underflow"
+ *  by raising the underflow flag for numbers less with exponent less than
+ *  expMin, but don't flush to zero until the exponent reaches MIN_EXP-digits.
+ *  Thus the smallest number we can represent would be:
+ *  1E(-(MIN_EXP-digits-1)*4),  eg, for digits=5, MIN_EXP=-32767, that would
+ *  be 1e-131092.</p>
+ *
+ *  <p>IEEE 854 defines that the implied radix point lies just to the right
+ *  of the most significant digit and to the left of the remaining digits.
+ *  This implementation puts the implied radix point to the left of all
+ *  digits including the most significant one.  The most significant digit
+ *  here is the one just to the right of the radix point.  This is a fine
+ *  detail and is really only a matter of definition.  Any side effects of
+ *  this can be rendered invisible by a subclass.</p>
+ * @see DfpField
+ * @version $Revision$ $Date$
+ * @since 2.2
+ */
+public class Dfp implements FieldElement<Dfp> {
+
+    /** The radix, or base of this system.  Set to 10000 */
+    public static final int RADIX = 10000;
+
+    /** The minimum exponent before underflow is signaled.  Flush to zero
+     *  occurs at minExp-DIGITS */
+    public static final int MIN_EXP = -32767;
+
+    /** The maximum exponent before overflow is signaled and results flushed
+     *  to infinity */
+    public static final int MAX_EXP =  32768;
+
+    /** The amount under/overflows are scaled by before going to trap handler */
+    public static final int ERR_SCALE = 32760;
+
+    /** Indicator value for normal finite numbers. */
+    public static final byte FINITE = 0;
+
+    /** Indicator value for Infinity. */
+    public static final byte INFINITE = 1;
+
+    /** Indicator value for signaling NaN. */
+    public static final byte SNAN = 2;
+
+    /** Indicator value for quiet NaN. */
+    public static final byte QNAN = 3;
+
+    /** String for NaN representation. */
+    private static final String NAN_STRING = "NaN";
+
+    /** String for positive infinity representation. */
+    private static final String POS_INFINITY_STRING = "Infinity";
+
+    /** String for negative infinity representation. */
+    private static final String NEG_INFINITY_STRING = "-Infinity";
+
+    /** Name for traps triggered by addition. */
+    private static final String ADD_TRAP = "add";
+
+    /** Name for traps triggered by multiplication. */
+    private static final String MULTIPLY_TRAP = "multiply";
+
+    /** Name for traps triggered by division. */
+    private static final String DIVIDE_TRAP = "divide";
+
+    /** Name for traps triggered by square root. */
+    private static final String SQRT_TRAP = "sqrt";
+
+    /** Name for traps triggered by alignment. */
+    private static final String ALIGN_TRAP = "align";
+
+    /** Name for traps triggered by truncation. */
+    private static final String TRUNC_TRAP = "trunc";
+
+    /** Name for traps triggered by nextAfter. */
+    private static final String NEXT_AFTER_TRAP = "nextAfter";
+
+    /** Name for traps triggered by lessThan. */
+    private static final String LESS_THAN_TRAP = "lessThan";
+
+    /** Name for traps triggered by greaterThan. */
+    private static final String GREATER_THAN_TRAP = "greaterThan";
+
+    /** Name for traps triggered by newInstance. */
+    private static final String NEW_INSTANCE_TRAP = "newInstance";
+
+    /** Mantissa. */
+    protected int[] mant;
+
+    /** Sign bit: & for positive, -1 for negative. */
+    protected byte sign;
+
+    /** Exponent. */
+    protected int exp;
+
+    /** Indicator for non-finite / non-number values. */
+    protected byte nans;
+
+    /** Factory building similar Dfp's. */
+    private final DfpField field;
+
+    /** Makes an instance with a value of zero.
+     * @param field field to which this instance belongs
+     */
+    protected Dfp(final DfpField field) {
+        mant = new int[field.getRadixDigits()];
+        sign = 1;
+        exp = 0;
+        nans = FINITE;
+        this.field = field;
+    }
+
+    /** Create an instance from a byte value.
+     * @param field field to which this instance belongs
+     * @param x value to convert to an instance
+     */
+    protected Dfp(final DfpField field, byte x) {
+        this(field, (long) x);
+    }
+
+    /** Create an instance from an int value.
+     * @param field field to which this instance belongs
+     * @param x value to convert to an instance
+     */
+    protected Dfp(final DfpField field, int x) {
+        this(field, (long) x);
+    }
+
+    /** Create an instance from a long value.
+     * @param field field to which this instance belongs
+     * @param x value to convert to an instance
+     */
+    protected Dfp(final DfpField field, long x) {
+
+        // initialize as if 0
+        mant = new int[field.getRadixDigits()];
+        nans = FINITE;
+        this.field = field;
+
+        boolean isLongMin = false;
+        if (x == Long.MIN_VALUE) {
+            // special case for Long.MIN_VALUE (-9223372036854775808)
+            // we must shift it before taking its absolute value
+            isLongMin = true;
+            ++x;
+        }
+
+        // set the sign
+        if (x < 0) {
+            sign = -1;
+            x = -x;
+        } else {
+            sign = 1;
+        }
+
+        exp = 0;
+        while (x != 0) {
+            System.arraycopy(mant, mant.length - exp, mant, mant.length - 1 - exp, exp);
+            mant[mant.length - 1] = (int) (x % RADIX);
+            x /= RADIX;
+            exp++;
+        }
+
+        if (isLongMin) {
+            // remove the shift added for Long.MIN_VALUE
+            // we know in this case that fixing the last digit is sufficient
+            for (int i = 0; i < mant.length - 1; i++) {
+                if (mant[i] != 0) {
+                    mant[i]++;
+                    break;
+                }
+            }
+        }
+    }
+
+    /** Create an instance from a double value.
+     * @param field field to which this instance belongs
+     * @param x value to convert to an instance
+     */
+    protected Dfp(final DfpField field, double x) {
+
+        // initialize as if 0
+        mant = new int[field.getRadixDigits()];
+        sign = 1;
+        exp = 0;
+        nans = FINITE;
+        this.field = field;
+
+        long bits = Double.doubleToLongBits(x);
+        long mantissa = bits & 0x000fffffffffffffL;
+        int exponent = (int) ((bits & 0x7ff0000000000000L) >> 52) - 1023;
+
+        if (exponent == -1023) {
+            // Zero or sub-normal
+            if (x == 0) {
+                return;
+            }
+
+            exponent++;
+
+            // Normalize the subnormal number
+            while ( (mantissa & 0x0010000000000000L) == 0) {
+                exponent--;
+                mantissa <<= 1;
+            }
+            mantissa &= 0x000fffffffffffffL;
+        }
+
+        if (exponent == 1024) {
+            // infinity or NAN
+            if (x != x) {
+                sign = (byte) 1;
+                nans = QNAN;
+            } else if (x < 0) {
+                sign = (byte) -1;
+                nans = INFINITE;
+            } else {
+                sign = (byte) 1;
+                nans = INFINITE;
+            }
+            return;
+        }
+
+        Dfp xdfp = new Dfp(field, mantissa);
+        xdfp = xdfp.divide(new Dfp(field, 4503599627370496l)).add(field.getOne());  // Divide by 2^52, then add one
+        xdfp = xdfp.multiply(DfpMath.pow(field.getTwo(), exponent));
+
+        if ((bits & 0x8000000000000000L) != 0) {
+            xdfp = xdfp.negate();
+        }
+
+        System.arraycopy(xdfp.mant, 0, mant, 0, mant.length);
+        sign = xdfp.sign;
+        exp  = xdfp.exp;
+        nans = xdfp.nans;
+
+    }
+
+    /** Copy constructor.
+     * @param d instance to copy
+     */
+    public Dfp(final Dfp d) {
+        mant  = d.mant.clone();
+        sign  = d.sign;
+        exp   = d.exp;
+        nans  = d.nans;
+        field = d.field;
+    }
+
+    /** Create an instance from a String representation.
+     * @param field field to which this instance belongs
+     * @param s string representation of the instance
+     */
+    protected Dfp(final DfpField field, final String s) {
+
+        // initialize as if 0
+        mant = new int[field.getRadixDigits()];
+        sign = 1;
+        exp = 0;
+        nans = FINITE;
+        this.field = field;
+
+        boolean decimalFound = false;
+        final int rsize = 4;   // size of radix in decimal digits
+        final int offset = 4;  // Starting offset into Striped
+        final char[] striped = new char[getRadixDigits() * rsize + offset * 2];
+
+        // Check some special cases
+        if (s.equals(POS_INFINITY_STRING)) {
+            sign = (byte) 1;
+            nans = (byte) INFINITE;
+            return;
+        }
+
+        if (s.equals(NEG_INFINITY_STRING)) {
+            sign = (byte) -1;
+            nans = (byte) INFINITE;
+            return;
+        }
+
+        if (s.equals(NAN_STRING)) {
+            sign = (byte) 1;
+            nans = (byte) QNAN;
+            return;
+        }
+
+        // Check for scientific notation
+        int p = s.indexOf("e");
+        if (p == -1) { // try upper case?
+            p = s.indexOf("E");
+        }
+
+        final String fpdecimal;
+        int sciexp = 0;
+        if (p != -1) {
+            // scientific notation
+            fpdecimal = s.substring(0, p);
+            String fpexp = s.substring(p+1);
+            boolean negative = false;
+
+            for (int i=0; i<fpexp.length(); i++)
+            {
+                if (fpexp.charAt(i) == '-')
+                {
+                    negative = true;
+                    continue;
+                }
+                if (fpexp.charAt(i) >= '0' && fpexp.charAt(i) <= '9')
+                    sciexp = sciexp * 10 + fpexp.charAt(i) - '0';
+            }
+
+            if (negative) {
+                sciexp = -sciexp;
+            }
+        } else {
+            // normal case
+            fpdecimal = s;
+        }
+
+        // If there is a minus sign in the number then it is negative
+        if (fpdecimal.indexOf("-") !=  -1) {
+            sign = -1;
+        }
+
+        // First off, find all of the leading zeros, trailing zeros, and significant digits
+        p = 0;
+
+        // Move p to first significant digit
+        int decimalPos = 0;
+        for (;;) {
+            if (fpdecimal.charAt(p) >= '1' && fpdecimal.charAt(p) <= '9') {
+                break;
+            }
+
+            if (decimalFound && fpdecimal.charAt(p) == '0') {
+                decimalPos--;
+            }
+
+            if (fpdecimal.charAt(p) == '.') {
+                decimalFound = true;
+            }
+
+            p++;
+
+            if (p == fpdecimal.length()) {
+                break;
+            }
+        }
+
+        // Copy the string onto Stripped
+        int q = offset;
+        striped[0] = '0';
+        striped[1] = '0';
+        striped[2] = '0';
+        striped[3] = '0';
+        int significantDigits=0;
+        for(;;) {
+            if (p == (fpdecimal.length())) {
+                break;
+            }
+
+            // Don't want to run pass the end of the array
+            if (q == mant.length*rsize+offset+1) {
+                break;
+            }
+
+            if (fpdecimal.charAt(p) == '.') {
+                decimalFound = true;
+                decimalPos = significantDigits;
+                p++;
+                continue;
+            }
+
+            if (fpdecimal.charAt(p) < '0' || fpdecimal.charAt(p) > '9') {
+                p++;
+                continue;
+            }
+
+            striped[q] = fpdecimal.charAt(p);
+            q++;
+            p++;
+            significantDigits++;
+        }
+
+
+        // If the decimal point has been found then get rid of trailing zeros.
+        if (decimalFound && q != offset) {
+            for (;;) {
+                q--;
+                if (q == offset) {
+                    break;
+                }
+                if (striped[q] == '0') {
+                    significantDigits--;
+                } else {
+                    break;
+                }
+            }
+        }
+
+        // special case of numbers like "0.00000"
+        if (decimalFound && significantDigits == 0) {
+            decimalPos = 0;
+        }
+
+        // Implicit decimal point at end of number if not present
+        if (!decimalFound) {
+            decimalPos = q-offset;
+        }
+
+        // Find the number of significant trailing zeros
+        q = offset;  // set q to point to first sig digit
+        p = significantDigits-1+offset;
+
+        int trailingZeros = 0;
+        while (p > q) {
+            if (striped[p] != '0') {
+                break;
+            }
+            trailingZeros++;
+            p--;
+        }
+
+        // Make sure the decimal is on a mod 10000 boundary
+        int i = ((rsize * 100) - decimalPos - sciexp % rsize) % rsize;
+        q -= i;
+        decimalPos += i;
+
+        // Make the mantissa length right by adding zeros at the end if necessary
+        while ((p - q) < (mant.length * rsize)) {
+            for (i = 0; i < rsize; i++) {
+                striped[++p] = '0';
+            }
+        }
+
+        // Ok, now we know how many trailing zeros there are,
+        // and where the least significant digit is
+        for (i = mant.length - 1; i >= 0; i--) {
+            mant[i] = (striped[q]   - '0') * 1000 +
+                      (striped[q+1] - '0') * 100  +
+                      (striped[q+2] - '0') * 10   +
+                      (striped[q+3] - '0');
+            q += 4;
+        }
+
+
+        exp = (decimalPos+sciexp) / rsize;
+
+        if (q < striped.length) {
+            // Is there possible another digit?
+            round((striped[q] - '0')*1000);
+        }
+
+    }
+
+    /** Creates an instance with a non-finite value.
+     * @param field field to which this instance belongs
+     * @param sign sign of the Dfp to create
+     * @param nans code of the value, must be one of {@link #INFINITE},
+     * {@link #SNAN},  {@link #QNAN}
+     */
+    protected Dfp(final DfpField field, final byte sign, final byte nans) {
+        this.field = field;
+        this.mant    = new int[field.getRadixDigits()];
+        this.sign    = sign;
+        this.exp     = 0;
+        this.nans    = nans;
+    }
+
+    /** Create an instance with a value of 0.
+     * Use this internally in preference to constructors to facilitate subclasses
+     * @return a new instance with a value of 0
+     */
+    public Dfp newInstance() {
+        return new Dfp(getField());
+    }
+
+    /** Create an instance from a byte value.
+     * @param x value to convert to an instance
+     * @return a new instance with value x
+     */
+    public Dfp newInstance(final byte x) {
+        return new Dfp(getField(), x);
+    }
+
+    /** Create an instance from an int value.
+     * @param x value to convert to an instance
+     * @return a new instance with value x
+     */
+    public Dfp newInstance(final int x) {
+        return new Dfp(getField(), x);
+    }
+
+    /** Create an instance from a long value.
+     * @param x value to convert to an instance
+     * @return a new instance with value x
+     */
+    public Dfp newInstance(final long x) {
+        return new Dfp(getField(), x);
+    }
+
+    /** Create an instance from a double value.
+     * @param x value to convert to an instance
+     * @return a new instance with value x
+     */
+    public Dfp newInstance(final double x) {
+        return new Dfp(getField(), x);
+    }
+
+    /** Create an instance by copying an existing one.
+     * Use this internally in preference to constructors to facilitate subclasses.
+     * @param d instance to copy
+     * @return a new instance with the same value as d
+     */
+    public Dfp newInstance(final Dfp d) {
+
+        // make sure we don't mix number with different precision
+        if (field.getRadixDigits() != d.field.getRadixDigits()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            final Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            return dotrap(DfpField.FLAG_INVALID, NEW_INSTANCE_TRAP, d, result);
+        }
+
+        return new Dfp(d);
+
+    }
+
+    /** Create an instance from a String representation.
+     * Use this internally in preference to constructors to facilitate subclasses.
+     * @param s string representation of the instance
+     * @return a new instance parsed from specified string
+     */
+    public Dfp newInstance(final String s) {
+        return new Dfp(field, s);
+    }
+
+    /** Creates an instance with a non-finite value.
+     * @param sig sign of the Dfp to create
+     * @param code code of the value, must be one of {@link #INFINITE},
+     * {@link #SNAN},  {@link #QNAN}
+     * @return a new instance with a non-finite value
+     */
+    public Dfp newInstance(final byte sig, final byte code) {
+        return field.newDfp(sig, code);
+    }
+
+    /** Get the {@link org.apache.commons.math.Field Field} (really a {@link DfpField}) to which the instance belongs.
+     * <p>
+     * The field is linked to the number of digits and acts as a factory
+     * for {@link Dfp} instances.
+     * </p>
+     * @return {@link org.apache.commons.math.Field Field} (really a {@link DfpField}) to which the instance belongs
+     */
+    public DfpField getField() {
+        return field;
+    }
+
+    /** Get the number of radix digits of the instance.
+     * @return number of radix digits
+     */
+    public int getRadixDigits() {
+        return field.getRadixDigits();
+    }
+
+    /** Get the constant 0.
+     * @return a Dfp with value zero
+     */
+    public Dfp getZero() {
+        return field.getZero();
+    }
+
+    /** Get the constant 1.
+     * @return a Dfp with value one
+     */
+    public Dfp getOne() {
+        return field.getOne();
+    }
+
+    /** Get the constant 2.
+     * @return a Dfp with value two
+     */
+    public Dfp getTwo() {
+        return field.getTwo();
+    }
+
+    /** Shift the mantissa left, and adjust the exponent to compensate.
+     */
+    protected void shiftLeft() {
+        for (int i = mant.length - 1; i > 0; i--) {
+            mant[i] = mant[i-1];
+        }
+        mant[0] = 0;
+        exp--;
+    }
+
+    /* Note that shiftRight() does not call round() as that round() itself
+     uses shiftRight() */
+    /** Shift the mantissa right, and adjust the exponent to compensate.
+     */
+    protected void shiftRight() {
+        for (int i = 0; i < mant.length - 1; i++) {
+            mant[i] = mant[i+1];
+        }
+        mant[mant.length - 1] = 0;
+        exp++;
+    }
+
+    /** Make our exp equal to the supplied one, this may cause rounding.
+     *  Also causes de-normalized numbers.  These numbers are generally
+     *  dangerous because most routines assume normalized numbers.
+     *  Align doesn't round, so it will return the last digit destroyed
+     *  by shifting right.
+     *  @param e desired exponent
+     *  @return last digit destroyed by shifting right
+     */
+    protected int align(int e) {
+        int lostdigit = 0;
+        boolean inexact = false;
+
+        int diff = exp - e;
+
+        int adiff = diff;
+        if (adiff < 0) {
+            adiff = -adiff;
+        }
+
+        if (diff == 0) {
+            return 0;
+        }
+
+        if (adiff > (mant.length + 1)) {
+            // Special case
+            Arrays.fill(mant, 0);
+            exp = e;
+
+            field.setIEEEFlagsBits(DfpField.FLAG_INEXACT);
+            dotrap(DfpField.FLAG_INEXACT, ALIGN_TRAP, this, this);
+
+            return 0;
+        }
+
+        for (int i = 0; i < adiff; i++) {
+            if (diff < 0) {
+                /* Keep track of loss -- only signal inexact after losing 2 digits.
+                 * the first lost digit is returned to add() and may be incorporated
+                 * into the result.
+                 */
+                if (lostdigit != 0) {
+                    inexact = true;
+                }
+
+                lostdigit = mant[0];
+
+                shiftRight();
+            } else {
+                shiftLeft();
+            }
+        }
+
+        if (inexact) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INEXACT);
+            dotrap(DfpField.FLAG_INEXACT, ALIGN_TRAP, this, this);
+        }
+
+        return lostdigit;
+
+    }
+
+    /** Check if instance is less than x.
+     * @param x number to check instance against
+     * @return true if instance is less than x and neither are NaN, false otherwise
+     */
+    public boolean lessThan(final Dfp x) {
+
+        // make sure we don't mix number with different precision
+        if (field.getRadixDigits() != x.field.getRadixDigits()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            final Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            dotrap(DfpField.FLAG_INVALID, LESS_THAN_TRAP, x, result);
+            return false;
+        }
+
+        /* if a nan is involved, signal invalid and return false */
+        if (isNaN() || x.isNaN()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            dotrap(DfpField.FLAG_INVALID, LESS_THAN_TRAP, x, newInstance(getZero()));
+            return false;
+        }
+
+        return compare(this, x) < 0;
+    }
+
+    /** Check if instance is greater than x.
+     * @param x number to check instance against
+     * @return true if instance is greater than x and neither are NaN, false otherwise
+     */
+    public boolean greaterThan(final Dfp x) {
+
+        // make sure we don't mix number with different precision
+        if (field.getRadixDigits() != x.field.getRadixDigits()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            final Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            dotrap(DfpField.FLAG_INVALID, GREATER_THAN_TRAP, x, result);
+            return false;
+        }
+
+        /* if a nan is involved, signal invalid and return false */
+        if (isNaN() || x.isNaN()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            dotrap(DfpField.FLAG_INVALID, GREATER_THAN_TRAP, x, newInstance(getZero()));
+            return false;
+        }
+
+        return compare(this, x) > 0;
+    }
+
+    /** Check if instance is infinite.
+     * @return true if instance is infinite
+     */
+    public boolean isInfinite() {
+        return nans == INFINITE;
+    }
+
+    /** Check if instance is not a number.
+     * @return true if instance is not a number
+     */
+    public boolean isNaN() {
+        return (nans == QNAN) || (nans == SNAN);
+    }
+
+    /** Check if instance is equal to x.
+     * @param other object to check instance against
+     * @return true if instance is equal to x and neither are NaN, false otherwise
+     */
+    @Override
+    public boolean equals(final Object other) {
+
+        if (other instanceof Dfp) {
+            final Dfp x = (Dfp) other;
+            if (isNaN() || x.isNaN() || field.getRadixDigits() != x.field.getRadixDigits()) {
+                return false;
+            }
+
+            return compare(this, x) == 0;
+        }
+
+        return false;
+
+    }
+
+    /**
+     * Gets a hashCode for the instance.
+     * @return a hash code value for this object
+     */
+    @Override
+    public int hashCode() {
+        return 17 + (sign << 8) + (nans << 16) + exp + mant.hashCode();
+    }
+
+    /** Check if instance is not equal to x.
+     * @param x number to check instance against
+     * @return true if instance is not equal to x and neither are NaN, false otherwise
+     */
+    public boolean unequal(final Dfp x) {
+        if (isNaN() || x.isNaN() || field.getRadixDigits() != x.field.getRadixDigits()) {
+            return false;
+        }
+
+        return greaterThan(x) || lessThan(x);
+    }
+
+    /** Compare two instances.
+     * @param a first instance in comparison
+     * @param b second instance in comparison
+     * @return -1 if a<b, 1 if a>b and 0 if a==b
+     *  Note this method does not properly handle NaNs or numbers with different precision.
+     */
+    private static int compare(final Dfp a, final Dfp b) {
+        // Ignore the sign of zero
+        if (a.mant[a.mant.length - 1] == 0 && b.mant[b.mant.length - 1] == 0 &&
+            a.nans == FINITE && b.nans == FINITE) {
+            return 0;
+        }
+
+        if (a.sign != b.sign) {
+            if (a.sign == -1) {
+                return -1;
+            } else {
+                return 1;
+            }
+        }
+
+        // deal with the infinities
+        if (a.nans == INFINITE && b.nans == FINITE) {
+            return a.sign;
+        }
+
+        if (a.nans == FINITE && b.nans == INFINITE) {
+            return -b.sign;
+        }
+
+        if (a.nans == INFINITE && b.nans == INFINITE) {
+            return 0;
+        }
+
+        // Handle special case when a or b is zero, by ignoring the exponents
+        if (b.mant[b.mant.length-1] != 0 && a.mant[b.mant.length-1] != 0) {
+            if (a.exp < b.exp) {
+                return -a.sign;
+            }
+
+            if (a.exp > b.exp) {
+                return a.sign;
+            }
+        }
+
+        // compare the mantissas
+        for (int i = a.mant.length - 1; i >= 0; i--) {
+            if (a.mant[i] > b.mant[i]) {
+                return a.sign;
+            }
+
+            if (a.mant[i] < b.mant[i]) {
+                return -a.sign;
+            }
+        }
+
+        return 0;
+
+    }
+
+    /** Round to nearest integer using the round-half-even method.
+     *  That is round to nearest integer unless both are equidistant.
+     *  In which case round to the even one.
+     *  @return rounded value
+     */
+    public Dfp rint() {
+        return trunc(DfpField.RoundingMode.ROUND_HALF_EVEN);
+    }
+
+    /** Round to an integer using the round floor mode.
+     * That is, round toward -Infinity
+     *  @return rounded value
+     */
+    public Dfp floor() {
+        return trunc(DfpField.RoundingMode.ROUND_FLOOR);
+    }
+
+    /** Round to an integer using the round ceil mode.
+     * That is, round toward +Infinity
+     *  @return rounded value
+     */
+    public Dfp ceil() {
+        return trunc(DfpField.RoundingMode.ROUND_CEIL);
+    }
+
+    /** Returns the IEEE remainder.
+     * @param d divisor
+     * @return this less n &times; d, where n is the integer closest to this/d
+     */
+    public Dfp remainder(final Dfp d) {
+
+        final Dfp result = this.subtract(this.divide(d).rint().multiply(d));
+
+        // IEEE 854-1987 says that if the result is zero, then it carries the sign of this
+        if (result.mant[mant.length-1] == 0) {
+            result.sign = sign;
+        }
+
+        return result;
+
+    }
+
+    /** Does the integer conversions with the specified rounding.
+     * @param rmode rounding mode to use
+     * @return truncated value
+     */
+    protected Dfp trunc(final DfpField.RoundingMode rmode) {
+        boolean changed = false;
+
+        if (isNaN()) {
+            return newInstance(this);
+        }
+
+        if (nans == INFINITE) {
+            return newInstance(this);
+        }
+
+        if (mant[mant.length-1] == 0) {
+            // a is zero
+            return newInstance(this);
+        }
+
+        /* If the exponent is less than zero then we can certainly
+         * return zero */
+        if (exp < 0) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INEXACT);
+            Dfp result = newInstance(getZero());
+            result = dotrap(DfpField.FLAG_INEXACT, TRUNC_TRAP, this, result);
+            return result;
+        }
+
+        /* If the exponent is greater than or equal to digits, then it
+         * must already be an integer since there is no precision left
+         * for any fractional part */
+
+        if (exp >= mant.length) {
+            return newInstance(this);
+        }
+
+        /* General case:  create another dfp, result, that contains the
+         * a with the fractional part lopped off.  */
+
+        Dfp result = newInstance(this);
+        for (int i = 0; i < mant.length-result.exp; i++) {
+            changed |= result.mant[i] != 0;
+            result.mant[i] = 0;
+        }
+
+        if (changed) {
+            switch (rmode) {
+                case ROUND_FLOOR:
+                    if (result.sign == -1) {
+                        // then we must increment the mantissa by one
+                        result = result.add(newInstance(-1));
+                    }
+                    break;
+
+                case ROUND_CEIL:
+                    if (result.sign == 1) {
+                        // then we must increment the mantissa by one
+                        result = result.add(getOne());
+                    }
+                    break;
+
+                case ROUND_HALF_EVEN:
+                default:
+                    final Dfp half = newInstance("0.5");
+                    Dfp a = subtract(result);  // difference between this and result
+                    a.sign = 1;            // force positive (take abs)
+                    if (a.greaterThan(half)) {
+                        a = newInstance(getOne());
+                        a.sign = sign;
+                        result = result.add(a);
+                    }
+
+                    /** If exactly equal to 1/2 and odd then increment */
+                    if (a.equals(half) && result.exp > 0 && (result.mant[mant.length-result.exp]&1) != 0) {
+                        a = newInstance(getOne());
+                        a.sign = sign;
+                        result = result.add(a);
+                    }
+                    break;
+            }
+
+            field.setIEEEFlagsBits(DfpField.FLAG_INEXACT);  // signal inexact
+            result = dotrap(DfpField.FLAG_INEXACT, TRUNC_TRAP, this, result);
+            return result;
+        }
+
+        return result;
+    }
+
+    /** Convert this to an integer.
+     * If greater than 2147483647, it returns 2147483647. If less than -2147483648 it returns -2147483648.
+     * @return converted number
+     */
+    public int intValue() {
+        Dfp rounded;
+        int result = 0;
+
+        rounded = rint();
+
+        if (rounded.greaterThan(newInstance(2147483647))) {
+            return 2147483647;
+        }
+
+        if (rounded.lessThan(newInstance(-2147483648))) {
+            return -2147483648;
+        }
+
+        for (int i = mant.length - 1; i >= mant.length - rounded.exp; i--) {
+            result = result * RADIX + rounded.mant[i];
+        }
+
+        if (rounded.sign == -1) {
+            result = -result;
+        }
+
+        return result;
+    }
+
+    /** Get the exponent of the greatest power of 10000 that is
+     *  less than or equal to the absolute value of this.  I.E.  if
+     *  this is 10<sup>6</sup> then log10K would return 1.
+     *  @return integer base 10000 logarithm
+     */
+    public int log10K() {
+        return exp - 1;
+    }
+
+    /** Get the specified  power of 10000.
+     * @param e desired power
+     * @return 10000<sup>e</sup>
+     */
+    public Dfp power10K(final int e) {
+        Dfp d = newInstance(getOne());
+        d.exp = e + 1;
+        return d;
+    }
+
+    /** Get the exponent of the greatest power of 10 that is less than or equal to abs(this).
+     *  @return integer base 10 logarithm
+     */
+    public int log10()  {
+        if (mant[mant.length-1] > 1000) {
+            return exp * 4 - 1;
+        }
+        if (mant[mant.length-1] > 100) {
+            return exp * 4 - 2;
+        }
+        if (mant[mant.length-1] > 10) {
+            return exp * 4 - 3;
+        }
+        return exp * 4 - 4;
+    }
+
+    /** Return the specified  power of 10.
+     * @param e desired power
+     * @return 10<sup>e</sup>
+     */
+    public Dfp power10(final int e) {
+        Dfp d = newInstance(getOne());
+
+        if (e >= 0) {
+            d.exp = e / 4 + 1;
+        } else {
+            d.exp = (e + 1) / 4;
+        }
+
+        switch ((e % 4 + 4) % 4) {
+            case 0:
+                break;
+            case 1:
+                d = d.multiply(10);
+                break;
+            case 2:
+                d = d.multiply(100);
+                break;
+            default:
+                d = d.multiply(1000);
+        }
+
+        return d;
+    }
+
+    /** Negate the mantissa of this by computing the complement.
+     *  Leaves the sign bit unchanged, used internally by add.
+     *  Denormalized numbers are handled properly here.
+     *  @param extra ???
+     *  @return ???
+     */
+    protected int complement(int extra) {
+
+        extra = RADIX-extra;
+        for (int i = 0; i < mant.length; i++) {
+            mant[i] = RADIX-mant[i]-1;
+        }
+
+        int rh = extra / RADIX;
+        extra = extra - rh * RADIX;
+        for (int i = 0; i < mant.length; i++) {
+            final int r = mant[i] + rh;
+            rh = r / RADIX;
+            mant[i] = r - rh * RADIX;
+        }
+
+        return extra;
+    }
+
+    /** Add x to this.
+     * @param x number to add
+     * @return sum of this and x
+     */
+    public Dfp add(final Dfp x) {
+
+        // make sure we don't mix number with different precision
+        if (field.getRadixDigits() != x.field.getRadixDigits()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            final Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            return dotrap(DfpField.FLAG_INVALID, ADD_TRAP, x, result);
+        }
+
+        /* handle special cases */
+        if (nans != FINITE || x.nans != FINITE) {
+            if (isNaN()) {
+                return this;
+            }
+
+            if (x.isNaN()) {
+                return x;
+            }
+
+            if (nans == INFINITE && x.nans == FINITE) {
+                return this;
+            }
+
+            if (x.nans == INFINITE && nans == FINITE) {
+                return x;
+            }
+
+            if (x.nans == INFINITE && nans == INFINITE && sign == x.sign) {
+                return x;
+            }
+
+            if (x.nans == INFINITE && nans == INFINITE && sign != x.sign) {
+                field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+                Dfp result = newInstance(getZero());
+                result.nans = QNAN;
+                result = dotrap(DfpField.FLAG_INVALID, ADD_TRAP, x, result);
+                return result;
+            }
+        }
+
+        /* copy this and the arg */
+        Dfp a = newInstance(this);
+        Dfp b = newInstance(x);
+
+        /* initialize the result object */
+        Dfp result = newInstance(getZero());
+
+        /* Make all numbers positive, but remember their sign */
+        final byte asign = a.sign;
+        final byte bsign = b.sign;
+
+        a.sign = 1;
+        b.sign = 1;
+
+        /* The result will be signed like the arg with greatest magnitude */
+        byte rsign = bsign;
+        if (compare(a, b) > 0) {
+            rsign = asign;
+        }
+
+        /* Handle special case when a or b is zero, by setting the exponent
+       of the zero number equal to the other one.  This avoids an alignment
+       which would cause catastropic loss of precision */
+        if (b.mant[mant.length-1] == 0) {
+            b.exp = a.exp;
+        }
+
+        if (a.mant[mant.length-1] == 0) {
+            a.exp = b.exp;
+        }
+
+        /* align number with the smaller exponent */
+        int aextradigit = 0;
+        int bextradigit = 0;
+        if (a.exp < b.exp) {
+            aextradigit = a.align(b.exp);
+        } else {
+            bextradigit = b.align(a.exp);
+        }
+
+        /* complement the smaller of the two if the signs are different */
+        if (asign != bsign) {
+            if (asign == rsign) {
+                bextradigit = b.complement(bextradigit);
+            } else {
+                aextradigit = a.complement(aextradigit);
+            }
+        }
+
+        /* add the mantissas */
+        int rh = 0; /* acts as a carry */
+        for (int i = 0; i < mant.length; i++) {
+            final int r = a.mant[i]+b.mant[i]+rh;
+            rh = r / RADIX;
+            result.mant[i] = r - rh * RADIX;
+        }
+        result.exp = a.exp;
+        result.sign = rsign;
+
+        /* handle overflow -- note, when asign!=bsign an overflow is
+         * normal and should be ignored.  */
+
+        if (rh != 0 && (asign == bsign)) {
+            final int lostdigit = result.mant[0];
+            result.shiftRight();
+            result.mant[mant.length-1] = rh;
+            final int excp = result.round(lostdigit);
+            if (excp != 0) {
+                result = dotrap(excp, ADD_TRAP, x, result);
+            }
+        }
+
+        /* normalize the result */
+        for (int i = 0; i < mant.length; i++) {
+            if (result.mant[mant.length-1] != 0) {
+                break;
+            }
+            result.shiftLeft();
+            if (i == 0) {
+                result.mant[0] = aextradigit+bextradigit;
+                aextradigit = 0;
+                bextradigit = 0;
+            }
+        }
+
+        /* result is zero if after normalization the most sig. digit is zero */
+        if (result.mant[mant.length-1] == 0) {
+            result.exp = 0;
+
+            if (asign != bsign) {
+                // Unless adding 2 negative zeros, sign is positive
+                result.sign = 1;  // Per IEEE 854-1987 Section 6.3
+            }
+        }
+
+        /* Call round to test for over/under flows */
+        final int excp = result.round(aextradigit + bextradigit);
+        if (excp != 0) {
+            result = dotrap(excp, ADD_TRAP, x, result);
+        }
+
+        return result;
+    }
+
+    /** Returns a number that is this number with the sign bit reversed.
+     * @return the opposite of this
+     */
+    public Dfp negate() {
+        Dfp result = newInstance(this);
+        result.sign = (byte) - result.sign;
+        return result;
+    }
+
+    /** Subtract x from this.
+     * @param x number to subtract
+     * @return difference of this and a
+     */
+    public Dfp subtract(final Dfp x) {
+        return add(x.negate());
+    }
+
+    /** Round this given the next digit n using the current rounding mode.
+     * @param n ???
+     * @return the IEEE flag if an exception occurred
+     */
+    protected int round(int n) {
+        boolean inc = false;
+        switch (field.getRoundingMode()) {
+            case ROUND_DOWN:
+                inc = false;
+                break;
+
+            case ROUND_UP:
+                inc = n != 0;       // round up if n!=0
+                break;
+
+            case ROUND_HALF_UP:
+                inc = n >= 5000;  // round half up
+                break;
+
+            case ROUND_HALF_DOWN:
+                inc = n > 5000;  // round half down
+                break;
+
+            case ROUND_HALF_EVEN:
+                inc = n > 5000 || (n == 5000 && (mant[0] & 1) == 1);  // round half-even
+                break;
+
+            case ROUND_HALF_ODD:
+                inc = n > 5000 || (n == 5000 && (mant[0] & 1) == 0);  // round half-odd
+                break;
+
+            case ROUND_CEIL:
+                inc = sign == 1 && n != 0;  // round ceil
+                break;
+
+            case ROUND_FLOOR:
+            default:
+                inc = sign == -1 && n != 0;  // round floor
+                break;
+        }
+
+        if (inc) {
+            // increment if necessary
+            int rh = 1;
+            for (int i = 0; i < mant.length; i++) {
+                final int r = mant[i] + rh;
+                rh = r / RADIX;
+                mant[i] = r - rh * RADIX;
+            }
+
+            if (rh != 0) {
+                shiftRight();
+                mant[mant.length-1] = rh;
+            }
+        }
+
+        // check for exceptional cases and raise signals if necessary
+        if (exp < MIN_EXP) {
+            // Gradual Underflow
+            field.setIEEEFlagsBits(DfpField.FLAG_UNDERFLOW);
+            return DfpField.FLAG_UNDERFLOW;
+        }
+
+        if (exp > MAX_EXP) {
+            // Overflow
+            field.setIEEEFlagsBits(DfpField.FLAG_OVERFLOW);
+            return DfpField.FLAG_OVERFLOW;
+        }
+
+        if (n != 0) {
+            // Inexact
+            field.setIEEEFlagsBits(DfpField.FLAG_INEXACT);
+            return DfpField.FLAG_INEXACT;
+        }
+
+        return 0;
+
+    }
+
+    /** Multiply this by x.
+     * @param x multiplicand
+     * @return product of this and x
+     */
+    public Dfp multiply(final Dfp x) {
+
+        // make sure we don't mix number with different precision
+        if (field.getRadixDigits() != x.field.getRadixDigits()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            final Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            return dotrap(DfpField.FLAG_INVALID, MULTIPLY_TRAP, x, result);
+        }
+
+        Dfp result = newInstance(getZero());
+
+        /* handle special cases */
+        if (nans != FINITE || x.nans != FINITE) {
+            if (isNaN()) {
+                return this;
+            }
+
+            if (x.isNaN()) {
+                return x;
+            }
+
+            if (nans == INFINITE && x.nans == FINITE && x.mant[mant.length-1] != 0) {
+                result = newInstance(this);
+                result.sign = (byte) (sign * x.sign);
+                return result;
+            }
+
+            if (x.nans == INFINITE && nans == FINITE && mant[mant.length-1] != 0) {
+                result = newInstance(x);
+                result.sign = (byte) (sign * x.sign);
+                return result;
+            }
+
+            if (x.nans == INFINITE && nans == INFINITE) {
+                result = newInstance(this);
+                result.sign = (byte) (sign * x.sign);
+                return result;
+            }
+
+            if ( (x.nans == INFINITE && nans == FINITE && mant[mant.length-1] == 0) ||
+                    (nans == INFINITE && x.nans == FINITE && x.mant[mant.length-1] == 0) ) {
+                field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+                result = newInstance(getZero());
+                result.nans = QNAN;
+                result = dotrap(DfpField.FLAG_INVALID, MULTIPLY_TRAP, x, result);
+                return result;
+            }
+        }
+
+        int[] product = new int[mant.length*2];  // Big enough to hold even the largest result
+
+        for (int i = 0; i < mant.length; i++) {
+            int rh = 0;  // acts as a carry
+            for (int j=0; j<mant.length; j++) {
+                int r = mant[i] * x.mant[j];    // multiply the 2 digits
+                r = r + product[i+j] + rh;  // add to the product digit with carry in
+
+                rh = r / RADIX;
+                product[i+j] = r - rh * RADIX;
+            }
+            product[i+mant.length] = rh;
+        }
+
+        // Find the most sig digit
+        int md = mant.length * 2 - 1;  // default, in case result is zero
+        for (int i = mant.length * 2 - 1; i >= 0; i--) {
+            if (product[i] != 0) {
+                md = i;
+                break;
+            }
+        }
+
+        // Copy the digits into the result
+        for (int i = 0; i < mant.length; i++) {
+            result.mant[mant.length - i - 1] = product[md - i];
+        }
+
+        // Fixup the exponent.
+        result.exp = exp + x.exp + md - 2 * mant.length + 1;
+        result.sign = (byte)((sign == x.sign)?1:-1);
+
+        if (result.mant[mant.length-1] == 0) {
+            // if result is zero, set exp to zero
+            result.exp = 0;
+        }
+
+        final int excp;
+        if (md > (mant.length-1)) {
+            excp = result.round(product[md-mant.length]);
+        } else {
+            excp = result.round(0); // has no effect except to check status
+        }
+
+        if (excp != 0) {
+            result = dotrap(excp, MULTIPLY_TRAP, x, result);
+        }
+
+        return result;
+
+    }
+
+    /** Multiply this by a single digit 0&lt;=x&lt;radix.
+     * There are speed advantages in this special case
+     * @param x multiplicand
+     * @return product of this and x
+     */
+    public Dfp multiply(final int x) {
+        Dfp result = newInstance(this);
+
+        /* handle special cases */
+        if (nans != FINITE) {
+            if (isNaN()) {
+                return this;
+            }
+
+            if (nans == INFINITE && x != 0) {
+                result = newInstance(this);
+                return result;
+            }
+
+            if (nans == INFINITE && x == 0) {
+                field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+                result = newInstance(getZero());
+                result.nans = QNAN;
+                result = dotrap(DfpField.FLAG_INVALID, MULTIPLY_TRAP, newInstance(getZero()), result);
+                return result;
+            }
+        }
+
+        /* range check x */
+        if (x < 0 || x >= RADIX) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            result = newInstance(getZero());
+            result.nans = QNAN;
+            result = dotrap(DfpField.FLAG_INVALID, MULTIPLY_TRAP, result, result);
+            return result;
+        }
+
+        int rh = 0;
+        for (int i = 0; i < mant.length; i++) {
+            final int r = mant[i] * x + rh;
+            rh = r / RADIX;
+            result.mant[i] = r - rh * RADIX;
+        }
+
+        int lostdigit = 0;
+        if (rh != 0) {
+            lostdigit = result.mant[0];
+            result.shiftRight();
+            result.mant[mant.length-1] = rh;
+        }
+
+        if (result.mant[mant.length-1] == 0) { // if result is zero, set exp to zero
+            result.exp = 0;
+        }
+
+        final int excp = result.round(lostdigit);
+        if (excp != 0) {
+            result = dotrap(excp, MULTIPLY_TRAP, result, result);
+        }
+
+        return result;
+    }
+
+    /** Divide this by divisor.
+     * @param divisor divisor
+     * @return quotient of this by divisor
+     */
+    public Dfp divide(Dfp divisor) {
+        int dividend[]; // current status of the dividend
+        int quotient[]; // quotient
+        int remainder[];// remainder
+        int qd;         // current quotient digit we're working with
+        int nsqd;       // number of significant quotient digits we have
+        int trial=0;    // trial quotient digit
+        int minadj;     // minimum adjustment
+        boolean trialgood; // Flag to indicate a good trail digit
+        int md=0;       // most sig digit in result
+        int excp;       // exceptions
+
+        // make sure we don't mix number with different precision
+        if (field.getRadixDigits() != divisor.field.getRadixDigits()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            final Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            return dotrap(DfpField.FLAG_INVALID, DIVIDE_TRAP, divisor, result);
+        }
+
+        Dfp result = newInstance(getZero());
+
+        /* handle special cases */
+        if (nans != FINITE || divisor.nans != FINITE) {
+            if (isNaN()) {
+                return this;
+            }
+
+            if (divisor.isNaN()) {
+                return divisor;
+            }
+
+            if (nans == INFINITE && divisor.nans == FINITE) {
+                result = newInstance(this);
+                result.sign = (byte) (sign * divisor.sign);
+                return result;
+            }
+
+            if (divisor.nans == INFINITE && nans == FINITE) {
+                result = newInstance(getZero());
+                result.sign = (byte) (sign * divisor.sign);
+                return result;
+            }
+
+            if (divisor.nans == INFINITE && nans == INFINITE) {
+                field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+                result = newInstance(getZero());
+                result.nans = QNAN;
+                result = dotrap(DfpField.FLAG_INVALID, DIVIDE_TRAP, divisor, result);
+                return result;
+            }
+        }
+
+        /* Test for divide by zero */
+        if (divisor.mant[mant.length-1] == 0) {
+            field.setIEEEFlagsBits(DfpField.FLAG_DIV_ZERO);
+            result = newInstance(getZero());
+            result.sign = (byte) (sign * divisor.sign);
+            result.nans = INFINITE;
+            result = dotrap(DfpField.FLAG_DIV_ZERO, DIVIDE_TRAP, divisor, result);
+            return result;
+        }
+
+        dividend = new int[mant.length+1];  // one extra digit needed
+        quotient = new int[mant.length+2];  // two extra digits needed 1 for overflow, 1 for rounding
+        remainder = new int[mant.length+1]; // one extra digit needed
+
+        /* Initialize our most significant digits to zero */
+
+        dividend[mant.length] = 0;
+        quotient[mant.length] = 0;
+        quotient[mant.length+1] = 0;
+        remainder[mant.length] = 0;
+
+        /* copy our mantissa into the dividend, initialize the
+       quotient while we are at it */
+
+        for (int i = 0; i < mant.length; i++) {
+            dividend[i] = mant[i];
+            quotient[i] = 0;
+            remainder[i] = 0;
+        }
+
+        /* outer loop.  Once per quotient digit */
+        nsqd = 0;
+        for (qd = mant.length+1; qd >= 0; qd--) {
+            /* Determine outer limits of our quotient digit */
+
+            // r =  most sig 2 digits of dividend
+            final int divMsb = dividend[mant.length]*RADIX+dividend[mant.length-1];
+            int min = divMsb       / (divisor.mant[mant.length-1]+1);
+            int max = (divMsb + 1) / divisor.mant[mant.length-1];
+
+            trialgood = false;
+            while (!trialgood) {
+                // try the mean
+                trial = (min+max)/2;
+
+                /* Multiply by divisor and store as remainder */
+                int rh = 0;
+                for (int i = 0; i < mant.length + 1; i++) {
+                    int dm = (i<mant.length)?divisor.mant[i]:0;
+                    final int r = (dm * trial) + rh;
+                    rh = r / RADIX;
+                    remainder[i] = r - rh * RADIX;
+                }
+
+                /* subtract the remainder from the dividend */
+                rh = 1;  // carry in to aid the subtraction
+                for (int i = 0; i < mant.length + 1; i++) {
+                    final int r = ((RADIX-1) - remainder[i]) + dividend[i] + rh;
+                    rh = r / RADIX;
+                    remainder[i] = r - rh * RADIX;
+                }
+
+                /* Lets analyze what we have here */
+                if (rh == 0) {
+                    // trial is too big -- negative remainder
+                    max = trial-1;
+                    continue;
+                }
+
+                /* find out how far off the remainder is telling us we are */
+                minadj = (remainder[mant.length] * RADIX)+remainder[mant.length-1];
+                minadj = minadj / (divisor.mant[mant.length-1]+1);
+
+                if (minadj >= 2) {
+                    min = trial+minadj;  // update the minimum
+                    continue;
+                }
+
+                /* May have a good one here, check more thoroughly.  Basically
+           its a good one if it is less than the divisor */
+                trialgood = false;  // assume false
+                for (int i = mant.length - 1; i >= 0; i--) {
+                    if (divisor.mant[i] > remainder[i]) {
+                        trialgood = true;
+                    }
+                    if (divisor.mant[i] < remainder[i]) {
+                        break;
+                    }
+                }
+
+                if (remainder[mant.length] != 0) {
+                    trialgood = false;
+                }
+
+                if (trialgood == false) {
+                    min = trial+1;
+                }
+            }
+
+            /* Great we have a digit! */
+            quotient[qd] = trial;
+            if (trial != 0 || nsqd != 0) {
+                nsqd++;
+            }
+
+            if (field.getRoundingMode() == DfpField.RoundingMode.ROUND_DOWN && nsqd == mant.length) {
+                // We have enough for this mode
+                break;
+            }
+
+            if (nsqd > mant.length) {
+                // We have enough digits
+                break;
+            }
+
+            /* move the remainder into the dividend while left shifting */
+            dividend[0] = 0;
+            for (int i = 0; i < mant.length; i++) {
+                dividend[i + 1] = remainder[i];
+            }
+        }
+
+        /* Find the most sig digit */
+        md = mant.length;  // default
+        for (int i = mant.length + 1; i >= 0; i--) {
+            if (quotient[i] != 0) {
+                md = i;
+                break;
+            }
+        }
+
+        /* Copy the digits into the result */
+        for (int i=0; i<mant.length; i++) {
+            result.mant[mant.length-i-1] = quotient[md-i];
+        }
+
+        /* Fixup the exponent. */
+        result.exp = exp - divisor.exp + md - mant.length;
+        result.sign = (byte) ((sign == divisor.sign) ? 1 : -1);
+
+        if (result.mant[mant.length-1] == 0) { // if result is zero, set exp to zero
+            result.exp = 0;
+        }
+
+        if (md > (mant.length-1)) {
+            excp = result.round(quotient[md-mant.length]);
+        } else {
+            excp = result.round(0);
+        }
+
+        if (excp != 0) {
+            result = dotrap(excp, DIVIDE_TRAP, divisor, result);
+        }
+
+        return result;
+    }
+
+    /** Divide by a single digit less than radix.
+     *  Special case, so there are speed advantages. 0 &lt;= divisor &lt; radix
+     * @param divisor divisor
+     * @return quotient of this by divisor
+     */
+    public Dfp divide(int divisor) {
+
+        // Handle special cases
+        if (nans != FINITE) {
+            if (isNaN()) {
+                return this;
+            }
+
+            if (nans == INFINITE) {
+                return newInstance(this);
+            }
+        }
+
+        // Test for divide by zero
+        if (divisor == 0) {
+            field.setIEEEFlagsBits(DfpField.FLAG_DIV_ZERO);
+            Dfp result = newInstance(getZero());
+            result.sign = sign;
+            result.nans = INFINITE;
+            result = dotrap(DfpField.FLAG_DIV_ZERO, DIVIDE_TRAP, getZero(), result);
+            return result;
+        }
+
+        // range check divisor
+        if (divisor < 0 || divisor >= RADIX) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            result = dotrap(DfpField.FLAG_INVALID, DIVIDE_TRAP, result, result);
+            return result;
+        }
+
+        Dfp result = newInstance(this);
+
+        int rl = 0;
+        for (int i = mant.length-1; i >= 0; i--) {
+            final int r = rl*RADIX + result.mant[i];
+            final int rh = r / divisor;
+            rl = r - rh * divisor;
+            result.mant[i] = rh;
+        }
+
+        if (result.mant[mant.length-1] == 0) {
+            // normalize
+            result.shiftLeft();
+            final int r = rl * RADIX;        // compute the next digit and put it in
+            final int rh = r / divisor;
+            rl = r - rh * divisor;
+            result.mant[0] = rh;
+        }
+
+        final int excp = result.round(rl * RADIX / divisor);  // do the rounding
+        if (excp != 0) {
+            result = dotrap(excp, DIVIDE_TRAP, result, result);
+        }
+
+        return result;
+
+    }
+
+    /** Compute the square root.
+     * @return square root of the instance
+     */
+    public Dfp sqrt() {
+
+        // check for unusual cases
+        if (nans == FINITE && mant[mant.length-1] == 0) {
+            // if zero
+            return newInstance(this);
+        }
+
+        if (nans != FINITE) {
+            if (nans == INFINITE && sign == 1) {
+                // if positive infinity
+                return newInstance(this);
+            }
+
+            if (nans == QNAN) {
+                return newInstance(this);
+            }
+
+            if (nans == SNAN) {
+                Dfp result;
+
+                field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+                result = newInstance(this);
+                result = dotrap(DfpField.FLAG_INVALID, SQRT_TRAP, null, result);
+                return result;
+            }
+        }
+
+        if (sign == -1) {
+            // if negative
+            Dfp result;
+
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            result = newInstance(this);
+            result.nans = QNAN;
+            result = dotrap(DfpField.FLAG_INVALID, SQRT_TRAP, null, result);
+            return result;
+        }
+
+        Dfp x = newInstance(this);
+
+        /* Lets make a reasonable guess as to the size of the square root */
+        if (x.exp < -1 || x.exp > 1) {
+            x.exp = this.exp / 2;
+        }
+
+        /* Coarsely estimate the mantissa */
+        switch (x.mant[mant.length-1] / 2000) {
+            case 0:
+                x.mant[mant.length-1] = x.mant[mant.length-1]/2+1;
+                break;
+            case 2:
+                x.mant[mant.length-1] = 1500;
+                break;
+            case 3:
+                x.mant[mant.length-1] = 2200;
+                break;
+            default:
+                x.mant[mant.length-1] = 3000;
+        }
+
+        Dfp dx = newInstance(x);
+
+        /* Now that we have the first pass estimate, compute the rest
+       by the formula dx = (y - x*x) / (2x); */
+
+        Dfp px  = getZero();
+        Dfp ppx = getZero();
+        while (x.unequal(px)) {
+            dx = newInstance(x);
+            dx.sign = -1;
+            dx = dx.add(this.divide(x));
+            dx = dx.divide(2);
+            ppx = px;
+            px = x;
+            x = x.add(dx);
+
+            if (x.equals(ppx)) {
+                // alternating between two values
+                break;
+            }
+
+            // if dx is zero, break.  Note testing the most sig digit
+            // is a sufficient test since dx is normalized
+            if (dx.mant[mant.length-1] == 0) {
+                break;
+            }
+        }
+
+        return x;
+
+    }
+
+    /** Get a string representation of the instance.
+     * @return string representation of the instance
+     */
+    public String toString() {
+        if (nans != FINITE) {
+            // if non-finite exceptional cases
+            if (nans == INFINITE) {
+                return (sign < 0) ? NEG_INFINITY_STRING : POS_INFINITY_STRING;
+            } else {
+                return NAN_STRING;
+            }
+        }
+
+        if (exp > mant.length || exp < -1) {
+            return dfp2sci();
+        }
+
+        return dfp2string();
+
+    }
+
+    /** Convert an instance to a string using scientific notation.
+     * @param a instance to convert
+     * @return string representation of the instance in scientific notation
+     */
+    protected String dfp2sci() {
+        char rawdigits[]    = new char[mant.length * 4];
+        char outputbuffer[] = new char[mant.length * 4 + 20];
+        int p;
+        int q;
+        int e;
+        int ae;
+        int shf;
+
+        // Get all the digits
+        p = 0;
+        for (int i = mant.length - 1; i >= 0; i--) {
+            rawdigits[p++] = (char) ((mant[i] / 1000) + '0');
+            rawdigits[p++] = (char) (((mant[i] / 100) %10) + '0');
+            rawdigits[p++] = (char) (((mant[i] / 10) % 10) + '0');
+            rawdigits[p++] = (char) (((mant[i]) % 10) + '0');
+        }
+
+        // Find the first non-zero one
+        for (p = 0; p < rawdigits.length; p++) {
+            if (rawdigits[p] != '0') {
+                break;
+            }
+        }
+        shf = p;
+
+        // Now do the conversion
+        q = 0;
+        if (sign == -1) {
+            outputbuffer[q++] = '-';
+        }
+
+        if (p != rawdigits.length) {
+            // there are non zero digits...
+            outputbuffer[q++] = rawdigits[p++];
+            outputbuffer[q++] = '.';
+
+            while (p<rawdigits.length) {
+                outputbuffer[q++] = rawdigits[p++];
+            }
+        } else {
+            outputbuffer[q++] = '0';
+            outputbuffer[q++] = '.';
+            outputbuffer[q++] = '0';
+            outputbuffer[q++] = 'e';
+            outputbuffer[q++] = '0';
+            return new String(outputbuffer, 0, 5);
+        }
+
+        outputbuffer[q++] = 'e';
+
+        // Find the msd of the exponent
+
+        e = exp * 4 - shf - 1;
+        ae = e;
+        if (e < 0) {
+            ae = -e;
+        }
+
+        // Find the largest p such that p < e
+        for (p = 1000000000; p > ae; p /= 10) {
+            // nothing to do
+        }
+
+        if (e < 0) {
+            outputbuffer[q++] = '-';
+        }
+
+        while (p > 0) {
+            outputbuffer[q++] = (char)(ae / p + '0');
+            ae = ae % p;
+            p = p / 10;
+        }
+
+        return new String(outputbuffer, 0, q);
+
+    }
+
+    /** Convert an instance to a string using normal notation.
+     * @param a instance to convert
+     * @return string representation of the instance in normal notation
+     */
+    protected String dfp2string() {
+        char buffer[] = new char[mant.length*4 + 20];
+        int p = 1;
+        int q;
+        int e = exp;
+        boolean pointInserted = false;
+
+        buffer[0] = ' ';
+
+        if (e <= 0) {
+            buffer[p++] = '0';
+            buffer[p++] = '.';
+            pointInserted = true;
+        }
+
+        while (e < 0) {
+            buffer[p++] = '0';
+            buffer[p++] = '0';
+            buffer[p++] = '0';
+            buffer[p++] = '0';
+            e++;
+        }
+
+        for (int i = mant.length - 1; i >= 0; i--) {
+            buffer[p++] = (char) ((mant[i] / 1000) + '0');
+            buffer[p++] = (char) (((mant[i] / 100) % 10) + '0');
+            buffer[p++] = (char) (((mant[i] / 10) % 10) + '0');
+            buffer[p++] = (char) (((mant[i]) % 10) + '0');
+            if (--e == 0) {
+                buffer[p++] = '.';
+                pointInserted = true;
+            }
+        }
+
+        while (e > 0) {
+            buffer[p++] = '0';
+            buffer[p++] = '0';
+            buffer[p++] = '0';
+            buffer[p++] = '0';
+            e--;
+        }
+
+        if (!pointInserted) {
+            // Ensure we have a radix point!
+            buffer[p++] = '.';
+        }
+
+        // Suppress leading zeros
+        q = 1;
+        while (buffer[q] == '0') {
+            q++;
+        }
+        if (buffer[q] == '.') {
+            q--;
+        }
+
+        // Suppress trailing zeros
+        while (buffer[p-1] == '0') {
+            p--;
+        }
+
+        // Insert sign
+        if (sign < 0) {
+            buffer[--q] = '-';
+        }
+
+        return new String(buffer, q, p - q);
+
+    }
+
+    /** Raises a trap.  This does not set the corresponding flag however.
+     *  @param type the trap type
+     *  @param what - name of routine trap occurred in
+     *  @param oper - input operator to function
+     *  @param result - the result computed prior to the trap
+     *  @return The suggested return value from the trap handler
+     */
+    public Dfp dotrap(int type, String what, Dfp oper, Dfp result) {
+        Dfp def = result;
+
+        switch (type) {
+            case DfpField.FLAG_INVALID:
+                def = newInstance(getZero());
+                def.sign = result.sign;
+                def.nans = QNAN;
+                break;
+
+            case DfpField.FLAG_DIV_ZERO:
+                if (nans == FINITE && mant[mant.length-1] != 0) {
+                    // normal case, we are finite, non-zero
+                    def = newInstance(getZero());
+                    def.sign = (byte)(sign*oper.sign);
+                    def.nans = INFINITE;
+                }
+
+                if (nans == FINITE && mant[mant.length-1] == 0) {
+                    //  0/0
+                    def = newInstance(getZero());
+                    def.nans = QNAN;
+                }
+
+                if (nans == INFINITE || nans == QNAN) {
+                    def = newInstance(getZero());
+                    def.nans = QNAN;
+                }
+
+                if (nans == INFINITE || nans == SNAN) {
+                    def = newInstance(getZero());
+                    def.nans = QNAN;
+                }
+                break;
+
+            case DfpField.FLAG_UNDERFLOW:
+                if ( (result.exp+mant.length) < MIN_EXP) {
+                    def = newInstance(getZero());
+                    def.sign = result.sign;
+                } else {
+                    def = newInstance(result);  // gradual underflow
+                }
+                result.exp = result.exp + ERR_SCALE;
+                break;
+
+            case DfpField.FLAG_OVERFLOW:
+                result.exp = result.exp - ERR_SCALE;
+                def = newInstance(getZero());
+                def.sign = result.sign;
+                def.nans = INFINITE;
+                break;
+
+            default: def = result; break;
+        }
+
+        return trap(type, what, oper, def, result);
+
+    }
+
+    /** Trap handler.  Subclasses may override this to provide trap
+     *  functionality per IEEE 854-1987.
+     *
+     *  @param type  The exception type - e.g. FLAG_OVERFLOW
+     *  @param what  The name of the routine we were in e.g. divide()
+     *  @param oper  An operand to this function if any
+     *  @param def   The default return value if trap not enabled
+     *  @param result    The result that is specified to be delivered per
+     *                   IEEE 854, if any
+     *  @return the value that should be return by the operation triggering the trap
+     */
+    protected Dfp trap(int type, String what, Dfp oper, Dfp def, Dfp result) {
+        return def;
+    }
+
+    /** Returns the type - one of FINITE, INFINITE, SNAN, QNAN.
+     * @return type of the number
+     */
+    public int classify() {
+        return nans;
+    }
+
+    /** Creates an instance that is the same as x except that it has the sign of y.
+     * abs(x) = dfp.copysign(x, dfp.one)
+     * @param x number to get the value from
+     * @param y number to get the sign from
+     * @return a number with the value of x and the sign of y
+     */
+    public static Dfp copysign(final Dfp x, final Dfp y) {
+        Dfp result = x.newInstance(x);
+        result.sign = y.sign;
+        return result;
+    }
+
+    /** Returns the next number greater than this one in the direction of x.
+     * If this==x then simply returns this.
+     * @param x direction where to look at
+     * @return closest number next to instance in the direction of x
+     */
+    public Dfp nextAfter(final Dfp x) {
+
+        // make sure we don't mix number with different precision
+        if (field.getRadixDigits() != x.field.getRadixDigits()) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INVALID);
+            final Dfp result = newInstance(getZero());
+            result.nans = QNAN;
+            return dotrap(DfpField.FLAG_INVALID, NEXT_AFTER_TRAP, x, result);
+        }
+
+        // if this is greater than x
+        boolean up = false;
+        if (this.lessThan(x)) {
+            up = true;
+        }
+
+        if (compare(this, x) == 0) {
+            return newInstance(x);
+        }
+
+        if (lessThan(getZero())) {
+            up = !up;
+        }
+
+        final Dfp inc;
+        Dfp result;
+        if (up) {
+            inc = newInstance(getOne());
+            inc.exp = this.exp-mant.length+1;
+            inc.sign = this.sign;
+
+            if (this.equals(getZero())) {
+                inc.exp = MIN_EXP-mant.length;
+            }
+
+            result = add(inc);
+        } else {
+            inc = newInstance(getOne());
+            inc.exp = this.exp;
+            inc.sign = this.sign;
+
+            if (this.equals(inc)) {
+                inc.exp = this.exp-mant.length;
+            } else {
+                inc.exp = this.exp-mant.length+1;
+            }
+
+            if (this.equals(getZero())) {
+                inc.exp = MIN_EXP-mant.length;
+            }
+
+            result = this.subtract(inc);
+        }
+
+        if (result.classify() == INFINITE && this.classify() != INFINITE) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INEXACT);
+            result = dotrap(DfpField.FLAG_INEXACT, NEXT_AFTER_TRAP, x, result);
+        }
+
+        if (result.equals(getZero()) && this.equals(getZero()) == false) {
+            field.setIEEEFlagsBits(DfpField.FLAG_INEXACT);
+            result = dotrap(DfpField.FLAG_INEXACT, NEXT_AFTER_TRAP, x, result);
+        }
+
+        return result;
+
+    }
+
+    /** Convert the instance into a double.
+     * @return a double approximating the instance
+     * @see #toSplitDouble()
+     */
+    public double toDouble() {
+
+        if (isInfinite()) {
+            if (lessThan(getZero())) {
+                return Double.NEGATIVE_INFINITY;
+            } else {
+                return Double.POSITIVE_INFINITY;
+            }
+        }
+
+        if (isNaN()) {
+            return Double.NaN;
+        }
+
+        Dfp y = this;
+        boolean negate = false;
+        if (lessThan(getZero())) {
+            y = negate();
+            negate = true;
+        }
+
+        /* Find the exponent, first estimate by integer log10, then adjust.
+         Should be faster than doing a natural logarithm.  */
+        int exponent = (int)(y.log10() * 3.32);
+        if (exponent < 0) {
+            exponent--;
+        }
+
+        Dfp tempDfp = DfpMath.pow(getTwo(), exponent);
+        while (tempDfp.lessThan(y) || tempDfp.equals(y)) {
+            tempDfp = tempDfp.multiply(2);
+            exponent++;
+        }
+        exponent--;
+
+        /* We have the exponent, now work on the mantissa */
+
+        y = y.divide(DfpMath.pow(getTwo(), exponent));
+        if (exponent > -1023) {
+            y = y.subtract(getOne());
+        }
+
+        if (exponent < -1074) {
+            return 0;
+        }
+
+        if (exponent > 1023) {
+            return negate ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
+        }
+
+
+        y = y.multiply(newInstance(4503599627370496l)).rint();
+        String str = y.toString();
+        str = str.substring(0, str.length()-1);
+        long mantissa = Long.parseLong(str);
+
+        if (mantissa == 4503599627370496L) {
+            // Handle special case where we round up to next power of two
+            mantissa = 0;
+            exponent++;
+        }
+
+        /* Its going to be subnormal, so make adjustments */
+        if (exponent <= -1023) {
+            exponent--;
+        }
+
+        while (exponent < -1023) {
+            exponent++;
+            mantissa >>>= 1;
+        }
+
+        long bits = mantissa | ((exponent + 1023L) << 52);
+        double x = Double.longBitsToDouble(bits);
+
+        if (negate) {
+            x = -x;
+        }
+
+        return x;
+
+    }
+
+    /** Convert the instance into a split double.
+     * @return an array of two doubles which sum represent the instance
+     * @see #toDouble()
+     */
+    public double[] toSplitDouble() {
+        double split[] = new double[2];
+        long mask = 0xffffffffc0000000L;
+
+        split[0] = Double.longBitsToDouble(Double.doubleToLongBits(toDouble()) & mask);
+        split[1] = subtract(newInstance(split[0])).toDouble();
+
+        return split;
+    }
+
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/Dfp.java
------------------------------------------------------------------------------
    svn:eol-style = native

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math/dfp/Dfp.java
------------------------------------------------------------------------------
    svn:keywords = Author Date Id Revision