Tcl Source Code

Artifact [90233aa058]
Login

Artifact 90233aa058104811b65451fbe34bedfcc64f45f4:

Attachment "double.patch" to ticket [3120139fff] added by kennykb 2010-11-28 00:50:46.
? kaboodle
Index: generic/tclClock.c
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclClock.c,v
retrieving revision 1.76
diff -u -r1.76 tclClock.c
--- generic/tclClock.c	1 Oct 2010 12:52:49 -0000	1.76
+++ generic/tclClock.c	27 Nov 2010 17:40:52 -0000
@@ -1503,7 +1503,7 @@
 	fields->julianDay = JDAY_1_JAN_1_CE_JULIAN - 1
 		+ fields->dayOfMonth
 		+ daysInPriorMonths[year%4 == 0][month - 1]
-		+ (365 * ym1)
+		+ (ONE_YEAR * ym1)
 		+ ym1o4;
     }
 }
Index: generic/tclInt.decls
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclInt.decls,v
retrieving revision 1.150
diff -u -r1.150 tclInt.decls
--- generic/tclInt.decls	2 Oct 2010 00:23:44 -0000	1.150
+++ generic/tclInt.decls	27 Nov 2010 17:40:53 -0000
@@ -996,6 +996,11 @@
     int TclCopyChannel(Tcl_Interp *interp, Tcl_Channel inChan,
 	    Tcl_Channel outChan, Tcl_WideInt toRead, Tcl_Obj *cmdPtr)
 }
+
+declare 249 {
+    char* TclDoubleDigits(double dv, int ndigits, int flags,
+			  int* decpt, int* signum, char** endPtr)
+}
 
 ##############################################################################
 
Index: generic/tclInt.h
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclInt.h,v
retrieving revision 1.486
diff -u -r1.486 tclInt.h
--- generic/tclInt.h	15 Nov 2010 21:34:54 -0000	1.486
+++ generic/tclInt.h	27 Nov 2010 17:40:53 -0000
@@ -2795,6 +2795,32 @@
 				/* Procedure that unloads a loaded module */
 };
 
+/* Flags for conversion of doubles to digit strings */
+
+#define TCL_DD_SHORTEST 		0x4
+				/* Use the shortest possible string */
+#define TCL_DD_STEELE   		0x5
+				/* Use the original Steele&White algorithm */
+#define TCL_DD_E_FORMAT 		0x2
+				/* Use a fixed-length string of digits,
+				 * suitable for E format*/
+#define TCL_DD_F_FORMAT 		0x3
+				/* Use a fixed number of digits after the
+				 * decimal point, suitable for F format */
+
+#define TCL_DD_SHORTEN_FLAG 		0x4
+				/* Allow return of a shorter digit string
+				 * if it converts losslessly */
+#define TCL_DD_NO_QUICK 		0x8
+				/* Debug flag: forbid quick FP conversion */
+
+#define TCL_DD_CONVERSION_TYPE_MASK	0x3
+				/* Mask to isolate the conversion type */
+#define TCL_DD_STEELE0 			0x1
+				/* 'Steele&White' after masking */
+#define TCL_DD_SHORTEST0		0x0
+				/* 'Shortest possible' after masking */
+
 /*
  *----------------------------------------------------------------
  * Procedures shared among Tcl modules but not used by the outside world:
@@ -2843,7 +2869,6 @@
 MODULE_SCOPE ContLineLoc *TclContinuationsGet(Tcl_Obj *objPtr);
 MODULE_SCOPE void	TclContinuationsCopy(Tcl_Obj *objPtr,
 			    Tcl_Obj *originObjPtr);
-MODULE_SCOPE int	TclDoubleDigits(char *buf, double value, int *signum);
 MODULE_SCOPE void	TclDeleteNamespaceVars(Namespace *nsPtr);
 /* TIP #280 - Modified token based evulation, with line information. */
 MODULE_SCOPE int	TclEvalEx(Tcl_Interp *interp, const char *script,
Index: generic/tclIntDecls.h
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclIntDecls.h,v
retrieving revision 1.144
diff -u -r1.144 tclIntDecls.h
--- generic/tclIntDecls.h	2 Oct 2010 00:23:45 -0000	1.144
+++ generic/tclIntDecls.h	27 Nov 2010 17:40:53 -0000
@@ -596,6 +596,9 @@
 EXTERN int		TclCopyChannel(Tcl_Interp *interp,
 				Tcl_Channel inChan, Tcl_Channel outChan,
 				Tcl_WideInt toRead, Tcl_Obj *cmdPtr);
+/* 249 */
+EXTERN char*		TclDoubleDigits(double dv, int ndigits, int flags,
+				int*decpt, int*signum, char**endPtr);
 
 typedef struct TclIntStubs {
     int magic;
@@ -850,6 +853,7 @@
     int (*tclInitRewriteEnsemble) (Tcl_Interp *interp, int numRemoved, int numInserted, Tcl_Obj *const *objv); /* 246 */
     void (*tclResetRewriteEnsemble) (Tcl_Interp *interp, int isRootEnsemble); /* 247 */
     int (*tclCopyChannel) (Tcl_Interp *interp, Tcl_Channel inChan, Tcl_Channel outChan, Tcl_WideInt toRead, Tcl_Obj *cmdPtr); /* 248 */
+    char* (*tclDoubleDigits) (double dv, int ndigits, int flags, int*decpt, int*signum, char**endPtr); /* 249 */
 } TclIntStubs;
 
 #ifdef __cplusplus
@@ -1269,6 +1273,8 @@
 	(tclIntStubsPtr->tclResetRewriteEnsemble) /* 247 */
 #define TclCopyChannel \
 	(tclIntStubsPtr->tclCopyChannel) /* 248 */
+#define TclDoubleDigits \
+	(tclIntStubsPtr->tclDoubleDigits) /* 249 */
 
 #endif /* defined(USE_TCL_STUBS) */
 
Index: generic/tclStrToD.c
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclStrToD.c,v
retrieving revision 1.46
diff -u -r1.46 tclStrToD.c
--- generic/tclStrToD.c	21 May 2010 12:43:29 -0000	1.46
+++ generic/tclStrToD.c	27 Nov 2010 17:40:53 -0000
@@ -90,6 +90,66 @@
  * runtime).
  */
 
+/* Magic constants */
+
+#define LOG10_2 0.3010299956639812
+#define TWO_OVER_3LOG10 0.28952965460216784
+#define LOG10_3HALVES_PLUS_FUDGE 0.1760912590558
+
+/* Definitions of the parts of an IEEE754-format floating point number */
+
+#define SIGN_BIT 		0x80000000
+				/* Mask for the sign bit in the first
+				 * word of a double */
+#define EXP_MASK	     	0x7ff00000
+				/* Mask for the exponent field in the
+				 * first word of a double */
+#define EXP_SHIFT		20
+				/* Shift count to make the exponent an 
+				 * integer */
+#define HIDDEN_BIT		(((Tcl_WideUInt) 0x00100000) << 32)
+				/* Hidden 1 bit for the significand */
+#define HI_ORDER_SIG_MASK	0x000fffff
+				/* Mask for the high-order part of the
+				 * significand in the first word of a
+				 * double */
+#define SIG_MASK		(((Tcl_WideUInt) HI_ORDER_SIG_MASK << 32) \
+				 | 0xffffffff)
+				/* Mask for the 52-bit significand. */
+#define FP_PRECISION		53
+				/* Number of bits of significand plus the
+				 * hidden bit */
+#define EXPONENT_BIAS		0x3ff
+				/* Bias of the exponent 0 */
+
+/* Derived quantities */
+
+#define TEN_PMAX		22
+				/* floor(FP_PRECISION*log(2)/log(5)) */
+#define QUICK_MAX		14
+				/* floor((FP_PRECISION-1)*log(2)/log(10)) - 1 */
+#define BLETCH			0x10
+				/* Highest power of two that is greater than
+				 * DBL_MAX_10_EXP, divided by 16 */
+#define DIGIT_GROUP		8
+				/* floor(DIGIT_BIT*log(2)/log(10)) */
+
+/* Union used to dismantle floating point numbers. */
+
+typedef union Double {
+    struct {
+#ifdef WORDS_BIGENDIAN
+	int word0;
+	int word1;
+#else
+	int word1;
+	int word0;
+#endif
+    } w;
+    double d;
+    Tcl_WideUInt q;
+} Double;
+
 static int maxpow10_wide;	/* The powers of ten that can be represented
 				 * exactly as wide integers. */
 static Tcl_WideUInt *pow10_wide;
@@ -123,6 +183,7 @@
     1.0e+128,
     1.0e+256
 };
+
 static int n770_fp;		/* Flag is 1 on Nokia N770 floating point.
 				 * Nokia's floating point has the words
 				 * reversed: if big-endian is 7654 3210,
@@ -131,27 +192,161 @@
 				 * little-endian within the 32-bit words
 				 * but big-endian between them. */
 
+/* Table of powers of 5 that are small enough to fit in an mp_digit. */
+
+static const mp_digit dpow5[13] = {
+               1,              5,             25,            125,
+             625,           3125,          15625,          78125,
+          390625,        1953125,        9765625,       48828125,
+       244140625
+};
+
+/* Table of powers: pow5_13[n] = 5**(13*2**(n+1)) */
+static mp_int pow5_13[5];	/* Table of powers: 5**13, 5**26, 5**52,
+				 * 5**104, 5**208 */
+static const double tens[] = {
+    1e00, 1e01, 1e02, 1e03, 1e04, 1e05, 1e06, 1e07, 1e08, 1e09,
+    1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19,
+    1e20, 1e21, 1e22
+};
+
+static const int itens [] = {
+    1,
+    10,
+    100,
+    1000,
+    10000,
+    100000,
+    1000000,
+    10000000,
+    100000000
+};
+
+static const Tcl_WideUInt wtens[] = {
+    1, 10, 100, 1000, 10000, 100000, 1000000,
+    (Tcl_WideUInt) 1000000*10, 		(Tcl_WideUInt) 1000000*100,
+    (Tcl_WideUInt) 1000000*1000, 	(Tcl_WideUInt) 1000000*10000,
+    (Tcl_WideUInt) 1000000*100000, 	(Tcl_WideUInt) 1000000*1000000,
+    (Tcl_WideUInt) 1000000*1000000*10, 	(Tcl_WideUInt) 1000000*1000000*100,
+    (Tcl_WideUInt) 1000000*1000000*1000,(Tcl_WideUInt) 1000000*1000000*10000
+    
+};
+
+static const double bigtens[] = {
+    1e016, 1e032, 1e064, 1e128, 1e256
+};
+#define N_BIGTENS 5
+
+static const int log2pow5[27] = {
+    01,  3,  5,  7, 10, 12, 14, 17, 19, 21,
+    24, 26, 28, 31, 33, 35, 38, 40, 42, 45,
+    47, 49, 52, 54, 56, 59, 61
+};
+#define N_LOG2POW5 27
+
+static const Tcl_WideUInt wuipow5[27] = {
+    (Tcl_WideUInt) 1,		/* 5**0 */
+    (Tcl_WideUInt) 5,
+    (Tcl_WideUInt) 25,
+    (Tcl_WideUInt) 125,
+    (Tcl_WideUInt) 625,
+    (Tcl_WideUInt) 3125,	/* 5**5 */
+    (Tcl_WideUInt) 3125*5,
+    (Tcl_WideUInt) 3125*25,
+    (Tcl_WideUInt) 3125*125,
+    (Tcl_WideUInt) 3125*625,
+    (Tcl_WideUInt) 3125*3125,	/* 5**10 */
+    (Tcl_WideUInt) 3125*3125*5,
+    (Tcl_WideUInt) 3125*3125*25,
+    (Tcl_WideUInt) 3125*3125*125,
+    (Tcl_WideUInt) 3125*3125*625,
+    (Tcl_WideUInt) 3125*3125*3125, /* 5**15 */
+    (Tcl_WideUInt) 3125*3125*3125*5,
+    (Tcl_WideUInt) 3125*3125*3125*25,
+    (Tcl_WideUInt) 3125*3125*3125*125,
+    (Tcl_WideUInt) 3125*3125*3125*625,
+    (Tcl_WideUInt) 3125*3125*3125*3125,	/* 5**20 */
+    (Tcl_WideUInt) 3125*3125*3125*3125*5,
+    (Tcl_WideUInt) 3125*3125*3125*3125*25,
+    (Tcl_WideUInt) 3125*3125*3125*3125*125,
+    (Tcl_WideUInt) 3125*3125*3125*3125*625,
+    (Tcl_WideUInt) 3125*3125*3125*3125*3125,  /* 5**25 */
+    (Tcl_WideUInt) 3125*3125*3125*3125*3125*5 /* 5**26 */
+};
+
 /*
  * Static functions defined in this file.
  */
 
-static double		AbsoluteValue(double v, int *signum);
 static int		AccumulateDecimalDigit(unsigned, int, 
 			    Tcl_WideUInt *, mp_int *, int);
-static double		BignumToBiasedFrExp(const mp_int *big, int *machexp);
-static int		GetIntegerTimesPower(double v, mp_int *r, int *e);
 static double		MakeHighPrecisionDouble(int signum,
 			    mp_int *significand, int nSigDigs, int exponent);
 static double		MakeLowPrecisionDouble(int signum,
 			    Tcl_WideUInt significand, int nSigDigs,
 			    int exponent);
 static double		MakeNaN(int signum, Tcl_WideUInt tag);
-static Tcl_WideUInt	Nokia770Twiddle(Tcl_WideUInt w);
-static double		Pow10TimesFrExp(int exponent, double fraction,
-			    int *machexp);
 static double		RefineApproximation(double approx,
 			    mp_int *exactSignificand, int exponent);
+static void		MulPow5(mp_int*, unsigned, mp_int*);
+static int 		NormalizeRightward(Tcl_WideUInt*);
+static int		RequiredPrecision(Tcl_WideUInt);
+static void		DoubleToExpAndSig(double, Tcl_WideUInt*, int*, int*);
+static void		TakeAbsoluteValue(Double*, int*);
+static char*		FormatInfAndNaN(Double*, int*, char**);
+static char*		FormatZero(int*, char**);
+static int		ApproximateLog10(Tcl_WideUInt, int, int);
+static int		BetterLog10(double, int, int*);
+static void		ComputeScale(int, int, int*, int*, int*, int*);
+static void		SetPrecisionLimits(int, int, int*, int*, int*, int*);
+static char*		BumpUp(char*, char*, int*);
+static int		AdjustRange(double*, int);
+static char*		ShorteningQuickFormat(double, int, int, double, 
+			    char*, int*);
+static char*		StrictQuickFormat(double, int, int, double,
+			    char*, int*);
+static char*		QuickConversion(double, int, int, int, int, int, int,
+			    int*, char**);
+static void		CastOutPowersOf2(int*, int*, int*);
+static char*		ShorteningInt64Conversion(Double*, int, Tcl_WideUInt,
+			    int, int, int, int, int, int, int, int, int,
+			    int, int, int*, char**);
+static char*		StrictInt64Conversion(Double*, int, Tcl_WideUInt,
+			    int, int, int, int, int, int,
+			    int, int, int*, char**);
+static int		ShouldBankerRoundUpPowD(mp_int*, int, int);
+static int		ShouldBankerRoundUpToNextPowD(mp_int*, mp_int*,
+			    int, int, int, mp_int*);
+static char*		ShorteningBignumConversionPowD(Double* dPtr, 
+			    int convType, Tcl_WideUInt bw, int b2, int b5,
+			    int m2plus, int m2minus, int m5,
+			    int sd, int k, int len, 
+			    int ilim, int ilim1, int* decpt,
+			    char** endPtr);
+static char*		StrictBignumConversionPowD(Double* dPtr, int convType,
+			    Tcl_WideUInt bw, int b2, int b5,
+			    int sd, int k, int len, 
+			    int ilim, int ilim1, int* decpt,
+			    char** endPtr);
+static int		ShouldBankerRoundUp(mp_int*, mp_int*, int);
+static int		ShouldBankerRoundUpToNext(mp_int*, mp_int*, mp_int*,
+			    int, int, mp_int*);
+static char*		ShorteningBignumConversion(Double* dPtr, int convType,
+			    Tcl_WideUInt bw, int b2,
+			    int m2plus, int m2minus,
+			    int s2, int s5, int k, int len, 
+			    int ilim, int ilim1, int* decpt,
+			    char** endPtr);
+static char*		StrictBignumConversion(Double* dPtr, int convType,
+			    Tcl_WideUInt bw, int b2,
+			    int s2, int s5, int k, int len, 
+			    int ilim, int ilim1, int* decpt,
+			    char** endPtr);
+static double		BignumToBiasedFrExp(const mp_int *big, int *machexp);
+static double		Pow10TimesFrExp(int exponent, double fraction,
+			    int *machexp);
 static double		SafeLdExp(double fraction, int exponent);
+static Tcl_WideUInt	Nokia770Twiddle(Tcl_WideUInt w);
 
 /*
  *----------------------------------------------------------------------
@@ -1732,425 +1927,2373 @@
 }
 
 /*
- *----------------------------------------------------------------------
+ *-----------------------------------------------------------------------------
  *
- * TclDoubleDigits --
+ * MultPow5 --
+ *
+ *	Multiply a bignum by a power of 5.
  *
- *	Converts a double to a string of digits.
+ * Side effects:
+ *	Stores base*5**n in result
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static void
+MulPow5(mp_int* base, 		/* Number to multiply */
+	 unsigned n,		/* Power of 5 to multiply by */
+	 mp_int* result)	/* Place to store the result */
+{
+    mp_int* p = base;
+    int n13 = n / 13;
+    int r = n % 13;
+    if (r != 0) {
+	mp_mul_d(p, dpow5[r], result);
+	p = result;
+    }
+    r = 0;
+    while (n13 != 0) {
+	if (n13 & 1) {
+	    mp_mul(p, pow5_13+r, result);
+	    p = result;
+	}
+	n13 >>= 1;
+	++r;
+    }
+    if (p != result) {
+	mp_copy(p, result);
+    }
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * NormalizeRightward --
+ *
+ *	Shifts a number rightward until it is odd (that is, until the
+ *	least significant bit is nonzero.
  *
  * Results:
- *	Returns the position of the character in the string after which the
- *	decimal point should appear. Since the string contains only
- *	significant digits, the position may be less than zero or greater than
- *	the length of the string.
+ *	Returns the number of bit positions by which the number was shifted.
  *
  * Side effects:
- *	Stores the digits in the given buffer and sets 'signum' according to
- *	the sign of the number.
+ *	Shifts the number in place; *wPtr is replaced by the shifted number.
  *
- *----------------------------------------------------------------------
-
+ *-----------------------------------------------------------------------------
  */
 
-int
-TclDoubleDigits(
-    char *buffer,		/* Buffer in which to store the result, must
-				 * have at least 18 chars */
-    double v,			/* Number to convert. Must be finite, and not
-				 * NaN */
-    int *signum)		/* Output: 1 if the number is negative.
-				 * Should handle -0 correctly on the IEEE
-				 * architecture. */
-{
-    int e;			/* Power of FLT_RADIX that satisfies
-				 * v = f * FLT_RADIX**e */
-    int lowOK, highOK;
-    mp_int r;			/* Scaled significand. */
-    mp_int s;			/* Divisor such that v = r / s */
-    int smallestSig;		/* Flag == 1 iff v's significand is the
-				 * smallest that can be represented. */
-    mp_int mplus;		/* Scaled epsilon: (r + 2* mplus) == v(+)
-				 * where v(+) is the floating point successor
-				 * of v. */
-    mp_int mminus;		/* Scaled epsilon: (r - 2*mminus) == v(-)
-				 * where v(-) is the floating point
-				 * predecessor of v. */
-    mp_int temp;
-    int rfac2 = 0;		/* Powers of 2 and 5 by which large */
-    int rfac5 = 0;		/* integers should be scaled.	    */
-    int sfac2 = 0;
-    int sfac5 = 0;
-    int mplusfac2 = 0;
-    int mminusfac2 = 0;
-    char c;
-    int i, k, n;
-
-    /*
-     * Split the number into absolute value and signum.
-     */
-
-    v = AbsoluteValue(v, signum);
-
-    /*
-     * Handle zero specially.
-     */
-
-    if (v == 0.0) {
-	*buffer++ = '0';
-	*buffer++ = '\0';
-	return 1;
+inline static int
+NormalizeRightward(Tcl_WideUInt* wPtr)
+				/* INOUT: Number to shift */
+{
+    int rv = 0;
+    Tcl_WideUInt w = *wPtr;
+    if (!(w & (Tcl_WideUInt) 0xffffffff)) {
+	w >>= 32; rv += 32;
+    } 
+    if (!(w & (Tcl_WideUInt) 0xffff)) {
+	w >>= 16; rv += 16;
+    }
+    if (!(w & (Tcl_WideUInt) 0xff)) {
+	w >>= 8; rv += 8;
+    }
+    if (!(w & (Tcl_WideUInt) 0xf)) {
+	w >>= 4; rv += 4;
     }
-
-    /*
-     * Find a large integer r, and integer e, such that
-     *         v = r * FLT_RADIX**e
-     * and r is as small as possible. Also determine whether the significand
-     * is the smallest possible.
-     */
-
-    smallestSig = GetIntegerTimesPower(v, &r, &e);
-
-    lowOK = highOK = (mp_iseven(&r));
-
-    /*
-     * We are going to want to develop integers r, s, mplus, and mminus such
-     * that v = r / s, v(+)-v / 2 = mplus / s; v-v(-) / 2 = mminus / s and
-     * then scale either s or r, mplus, mminus by an appropriate power of ten.
-     *
-     * We actually do this by keeping track of the powers of 2 and 5 by which
-     * f is multiplied to yield v and by which 1 is multiplied to yield s,
-     * mplus, and mminus.
-     */
-
-    if (e >= 0) {
-	int bits = e * log2FLT_RADIX;
-
-	if (!smallestSig) {
-	    /*
-	     * Normal case, m+ and m- are both FLT_RADIX**e
-	     */
-
-	    rfac2 = bits + 1;
-	    sfac2 = 1;
-	    mplusfac2 = bits;
-	    mminusfac2 = bits;
-	} else {
-	    /*
-	     * If f is equal to the smallest significand, then we need another
-	     * factor of FLT_RADIX in s to cope with stepping to the next
-	     * smaller exponent when going to e's predecessor.
-	     */
-
-	    rfac2 = bits + log2FLT_RADIX + 1;
-	    sfac2 = 1 + log2FLT_RADIX;
-	    mplusfac2 = bits + log2FLT_RADIX;
-	    mminusfac2 = bits;
-	}
-    } else {
-	/*
-	 * v has digits after the binary point
-	 */
-
-	if (e <= DBL_MIN_EXP-DBL_MANT_DIG || !smallestSig) {
-	    /*
-	     * Either f isn't the smallest significand or e is the smallest
-	     * exponent. mplus and mminus will both be 1.
-	     */
-
-	    rfac2 = 1;
-	    sfac2 = 1 - e * log2FLT_RADIX;
-	    mplusfac2 = 0;
-	    mminusfac2 = 0;
-	} else {
-	    /*
-	     * f is the smallest significand, but e is not the smallest
-	     * exponent. We need to scale by FLT_RADIX again to cope with the
-	     * fact that v's predecessor has a smaller exponent.
-	     */
-
-	    rfac2 = 1 + log2FLT_RADIX;
-	    sfac2 = 1 + log2FLT_RADIX * (1 - e);
-	    mplusfac2 = FLT_RADIX;
-	    mminusfac2 = 0;
-	}
+    if (!(w & 0x3)) {
+	w >>= 2; rv += 2;
     }
+    if (!(w & 0x1)) {
+	w >>= 1; ++rv;
+    }
+    *wPtr = w;
+    return rv;
+}
+
+/*
+ *-----------------------------------------------------------------------------0
+ *
+ * RequiredPrecision --
+ *
+ *	Determines the number of bits needed to hold an intger.
+ *
+ * Results:
+ *	Returns the position of the most significant bit (0 - 63).
+ *	Returns 0 if the number is zero.
+ *
+ *----------------------------------------------------------------------------
+ */
 
-    /*
-     * Estimate the highest power of ten that will be needed to hold the
-     * result.
-     */
-
-    k = (int) ceil(log(v) / log(10.));
-    if (k >= 0) {
-	sfac2 += k;
-	sfac5 = k;
+static int
+RequiredPrecision(Tcl_WideUInt w)
+				/* Number to interrogate */
+{
+    int rv;
+    unsigned long wi;
+    if (w & ((Tcl_WideUInt) 0xffffffff << 32)) {
+	wi = w >> 32; rv = 32;
     } else {
-	rfac2 -= k;
-	mplusfac2 -= k;
-	mminusfac2 -= k;
-	rfac5 = -k;
+	wi = w; rv = 0;
     }
-
-    /*
-     * Scale r, s, mplus, mminus by the appropriate powers of 2 and 5.
-     */
-
-    mp_init_set(&mplus, 1);
-    for (i=0 ; i<=8 ; ++i) {
-	if (rfac5 & (1 << i)) {
-	    mp_mul(&mplus, pow5+i, &mplus);
-	}
+    if (wi & 0xffff0000) {
+	wi >>= 16; rv += 16;
     }
-    mp_mul(&r, &mplus, &r);
-    mp_mul_2d(&r, rfac2, &r);
-    mp_init_copy(&mminus, &mplus);
-    mp_mul_2d(&mplus, mplusfac2, &mplus);
-    mp_mul_2d(&mminus, mminusfac2, &mminus);
-    mp_init_set(&s, 1);
-    for (i=0 ; i<=8 ; ++i) {
-	if (sfac5 & (1 << i)) {
-	    mp_mul(&s, pow5+i, &s);
-	}
+    if (wi & 0xff00) {
+	wi >>= 8; rv += 8;
     }
-    mp_mul_2d(&s, sfac2, &s);
-
-    /*
-     * It is possible for k to be off by one because we used an inexact
-     * logarithm.
-     */
-
-    mp_init(&temp);
-    mp_add(&r, &mplus, &temp);
-    i = mp_cmp_mag(&temp, &s);
-    if (i>0 || (highOK && i==0)) {
-	mp_mul_d(&s, 10, &s);
-	k++;
-    } else {
-	mp_mul_d(&temp, 10, &temp);
-	i = mp_cmp_mag(&temp, &s);
-	if (i<0 || (highOK && i==0)) {
-	    mp_mul_d(&r, 10, &r);
-	    mp_mul_d(&mplus, 10, &mplus);
-	    mp_mul_d(&mminus, 10, &mminus);
-	    k--;
-	}
+    if (wi & 0xf0) {
+	wi >>= 4; rv += 4;
     }
-
-    /*
-     * At this point, k contains the power of ten by which we're scaling the
-     * result. r/s is at least 1/10 and strictly less than ten, and v = r/s *
-     * 10**k. mplus and mminus give the rounding limits.
-     */
-
-    for (;;) {
-	int tc1, tc2;
-
-	mp_mul_d(&r, 10, &r);
-	mp_div(&r, &s, &temp, &r);	/* temp = 10r / s; r = 10r mod s */
-	i = temp.dp[0];
-	mp_mul_d(&mplus, 10, &mplus);
-	mp_mul_d(&mminus, 10, &mminus);
-	tc1 = mp_cmp_mag(&r, &mminus);
-	if (lowOK) {
-	    tc1 = (tc1 <= 0);
-	} else {
-	    tc1 = (tc1 < 0);
-	}
-	mp_add(&r, &mplus, &temp);
-	tc2 = mp_cmp_mag(&temp, &s);
-	if (highOK) {
-	    tc2 = (tc2 >= 0);
-	} else {
-	    tc2 = (tc2 > 0);
-	}
-	if (!tc1) {
-	    if (!tc2) {
-		*buffer++ = '0' + i;
-	    } else {
-		c = (char) (i + '1');
-		break;
-	    }
-	} else {
-	    if (!tc2) {
-		c = (char) (i + '0');
-	    } else {
-		mp_mul_2d(&r, 1, &r);
-		n = mp_cmp_mag(&r, &s);
-		if (n < 0) {
-		    c = (char) (i + '0');
-		} else {
-		    c = (char) (i + '1');
-		}
-	    }
-	    break;
-	}
-    };
-    *buffer++ = c;
-    *buffer++ = '\0';
-
-    /*
-     * Free memory, and return.
-     */
-
-    mp_clear_multi(&r, &s, &mplus, &mminus, &temp, NULL);
-    return k;
+    if (wi & 0xc) {
+	wi >>= 2; rv += 2;
+    }
+    if (wi & 0x2) {
+	wi >>= 1; ++rv;
+    }
+    if (wi & 0x1) {
+	++rv;
+    }
+    return rv;
 }
 
 /*
- *----------------------------------------------------------------------
- *
- * AbsoluteValue --
+ *-----------------------------------------------------------------------------
  *
- *	Splits a 'double' into its absolute value and sign.
+ * DoubleToExpAndSig --
  *
- * Results:
- *	Returns the absolute value.
+ *	Separates a 'double' into exponent and significand.
  *
  * Side effects:
- *	Stores the signum in '*signum'.
+ *	Stores the significand in '*significand' and the exponent in
+ *	'*expon' so that dv == significand * 2.0**expon, and significand
+ *	is odd.  Also stores the position of the leftmost 1-bit in 'significand'
+ *	in 'bits'.
  *
- *----------------------------------------------------------------------
+ *-----------------------------------------------------------------------------
  */
 
-static double
-AbsoluteValue(
-    double v,			/* Number to split */
-    int *signum)		/* (Output) Sign of the number 1=-, 0=+ */
+inline static void
+DoubleToExpAndSig(double dv,	/* Number to convert */
+		  Tcl_WideUInt* significand,
+				/* OUTPUT: Significand of the number */
+		  int* expon,	/* OUTPUT: Exponent to multiply the number by */
+		  int* bits)	/* OUTPUT: Number of significant bits */
 {
-    /*
-     * Take the absolute value of the number, and report the number's sign.
-     * Take special steps to preserve signed zeroes in IEEE floating point.
-     * (We can't use fpclassify, because that's a C9x feature and we still
-     * have to build on C89 compilers.)
-     */
-
-#ifndef IEEE_FLOATING_POINT
-    if (v >= 0.0) {
-	*signum = 0;
+    Double d;			/* Number being converted */
+    Tcl_WideUInt z;		/* Significand under construction */
+    int de;			/* Exponent of the number */
+    int k;			/* Bit count */
+
+    d.d = dv;
+
+    /* Extract exponent and significand */
+
+    de = (d.w.word0 & EXP_MASK) >> EXP_SHIFT;
+    z = d.q & SIG_MASK;
+    if (de != 0) {
+	z |= HIDDEN_BIT;
+	k = NormalizeRightward(&z);
+	*bits = FP_PRECISION - k;
+	*expon = k + (de - EXPONENT_BIAS) - (FP_PRECISION-1);
     } else {
-	*signum = 1;
-	v = -v;
+	k = NormalizeRightward(&z);
+	*expon = k + (de - EXPONENT_BIAS) - (FP_PRECISION-1) + 1;
+	*bits = RequiredPrecision(z);
     }
-#else
-    union {
-	Tcl_WideUInt iv;
-	double dv;
-    } bitwhack;
-    bitwhack.dv = v;
-    if (n770_fp) {
-	bitwhack.iv = Nokia770Twiddle(bitwhack.iv);
+    *significand = z;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * TakeAbsoluteValue --
+ *
+ *	Takes the absolute value of a 'double' including 0, Inf and NaN
+ *
+ * Side effects:
+ *	The 'double' in *d is replaced with its absolute value. The
+ *	signum is stored in 'sign': 1 for negative, 0 for nonnegative.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static void
+TakeAbsoluteValue(Double* d,	/* Number to replace with absolute value */
+		  int* sign)	/* Place to put the signum */
+{
+    if (d->w.word0 & SIGN_BIT) {
+	*sign = 1;
+	d->w.word0 &= ~SIGN_BIT;
+    } else {
+	*sign = 0;
     }
-    if (bitwhack.iv & ((Tcl_WideUInt) 1 << 63)) {
-	*signum = 1;
-	bitwhack.iv &= ~((Tcl_WideUInt) 1 << 63);
-	if (n770_fp) {
-	    bitwhack.iv = Nokia770Twiddle(bitwhack.iv);
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * FormatInfAndNaN --
+ *
+ *	Bailout for formatting infinities and Not-A-Number.
+ *
+ * Results:
+ *	Returns one of the strings 'Infinity' and 'NaN'.
+ *
+ * Side effects:
+ *	Stores 9999 in *decpt, and sets '*endPtr' to designate the
+ *	terminating NUL byte of the string if 'endPtr' is not NULL. 
+ *
+ * The string returned must be freed by the caller using 'ckfree'.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+FormatInfAndNaN(Double* d,	/* Exceptional number to format */
+		int* decpt,	/* Decimal point to set to a bogus value */
+		char** endPtr)	/* Pointer to the end of the formatted data */
+{
+    char* retval;
+    *decpt = 9999;
+    if (!(d->w.word1) && !(d->w.word0 & HI_ORDER_SIG_MASK)) {
+	retval = ckalloc(9);
+	strcpy(retval, "Infinity");
+	if (endPtr) {
+	    *endPtr = retval + 8;
 	}
-	v = bitwhack.dv;
     } else {
-	*signum = 0;
+	retval = ckalloc(4);
+	strcpy(retval, "NaN");
+	if (endPtr) {
+	    *endPtr = retval + 3;
+	}
     }
-#endif
-    return v;
+    return retval;
 }
 
 /*
- *----------------------------------------------------------------------
+ *-----------------------------------------------------------------------------
  *
- * GetIntegerTimesPower --
+ * FormatZero --
  *
- *	Converts a floating point number to an exact integer times a power of
- *	the floating point radix.
+ *	Bailout to format a zero floating-point number.
  *
  * Results:
- *	Returns 1 if it converted the smallest significand, 0 otherwise.
+ *	Returns the constant string "0"
  *
  * Side effects:
- *	Initializes the integer value (does not just assign it), and stores
- *	the exponent.
+ *	Stores 1 in '*decpt' and puts a pointer to the NUL byte terminating
+ *	the string in '*endPtr' if 'endPtr' is not NULL.
  *
- *----------------------------------------------------------------------
+ *-----------------------------------------------------------------------------
  */
 
-static int
-GetIntegerTimesPower(
-    double v,			/* Value to convert */
-    mp_int *rPtr,		/* (Output) Integer value */
-    int *ePtr)			/* (Output) Power of FLT_RADIX by which r must
-				 * be multiplied to yield v*/
+inline static char*
+FormatZero(int* decpt,		/* Location of the decimal point */
+	   char** endPtr)	/* Pointer to the end of the formatted data */
 {
-    double a, f;
-    int e, i, n;
-
-    /*
-     * Develop f and e such that v = f * FLT_RADIX**e, with
-     * 1.0/FLT_RADIX <= f < 1.
-     */
-
-    f = frexp(v, &e);
-#if FLT_RADIX > 2
-    n = e % log2FLT_RADIX;
-    if (n > 0) {
-	n -= log2FLT_RADIX;
-	e += 1;
-	f *= ldexp(1.0, n);
+    char* retval = ckalloc(2);
+    strcpy(retval, "0");
+    if (endPtr) {
+	*endPtr = retval+1;
     }
-    e = (e - n) / log2FLT_RADIX;
-#endif
-    if (f == 1.0) {
-	f = 1.0 / FLT_RADIX;
-	e += 1;
+    *decpt = 0;
+    return retval;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ApproximateLog10 --
+ *
+ *	Computes a two-term Taylor series approximation to the common
+ *	log of a number, and computes the number's binary log.
+ *
+ * Results:
+ *	Return an approximation to floor(log10(bw*2**be)) that is either
+ *	exact or 1 too high.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static int
+ApproximateLog10(Tcl_WideUInt bw,
+				/* Integer significand of the number */
+		 int be,	/* Power of two to scale bw */
+		 int bbits)	/* Number of bits of precision in bw */
+{
+    int i;			/* Log base 2 of the number */
+    int k;			/* Floor(Log base 10 of the number) */
+    double ds;			/* Mantissa of the number */
+    Double d2;
+
+    /*
+     * Compute i and d2 such that d = d2*2**i, and 1 < d2 < 2.
+     * Compute an approximation to log10(d), 
+     *   log10(d) ~ log10(2) * i + log10(1.5) 
+     *            + (significand-1.5)/(1.5 * log(10))
+     */
+
+    d2.q = bw << (FP_PRECISION - bbits) & SIG_MASK;
+    d2.w.word0 |= (EXPONENT_BIAS) << EXP_SHIFT;
+    i = be + bbits - 1;
+    ds = (d2.d - 1.5) * TWO_OVER_3LOG10
+	+ LOG10_3HALVES_PLUS_FUDGE
+	+ LOG10_2 * i;
+    k = (int) ds;
+    if (k > ds) {
+	--k;
     }
+    return k;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * BetterLog10 --
+ *
+ *	Improves the result of ApproximateLog10 for numbers in the range
+ *	1 .. 10**(TEN_PMAX)-1
+ *
+ * Side effects:
+ *	Sets k_check to 0 if the new result is known to be exact, and to
+ *	1 if it may still be one too high.
+ *
+ * Results:
+ *	Returns the improved approximation to log10(d)
+ *
+ *-----------------------------------------------------------------------------
+ */
 
-    /*
-     * If the original number was denormalized, adjust e and f to be denormal
-     * as well.
+inline static int
+BetterLog10(double d,		/* Original number to format */
+	  int k,		/* Characteristic(Log base 10) of the number */
+	  int* k_check)		/* Flag == 1 if k is inexact */
+{
+    /* 
+     * Performance hack. If k is in the range 0..TEN_PMAX, then we can
+     * use a powers-of-ten table to check it.
      */
-
-    if (e < DBL_MIN_EXP) {
-	n = mantBits + (e - DBL_MIN_EXP)*log2FLT_RADIX;
-	f = ldexp(f, (e - DBL_MIN_EXP)*log2FLT_RADIX);
-	e = DBL_MIN_EXP;
-	n = (n + DIGIT_BIT - 1) / DIGIT_BIT;
+    if (k >= 0 && k <= TEN_PMAX) {
+	if (d < tens[k]) {
+	    k--;
+	}
+	*k_check = 0;
     } else {
-	n = mantDIGIT;
+	*k_check = 1;
     }
+    return k;
+}
+
+/* 
+ *-----------------------------------------------------------------------------
+ *
+ * ComputeScale --
+ *
+ *	Prepares to format a floating-point number as decimal.
+ *
+ * Parameters:
+ *	floor(log10*x) is k (or possibly k-1).  floor(log2(x) is i.
+ *	The significand of x requires bbits bits to represent.
+ *
+ * Results:
+ *	Determines integers b2, b5, s2, s5 so that sig*2**b2*5**b5/2**s2*2**s5
+ *	exactly represents the value of the x/10**k. This value will lie
+ *	in the range [1 .. 10), and allows for computing successive digits
+ *	by multiplying sig%10 by 10.
+ *
+ *-----------------------------------------------------------------------------
+ */
 
-    /*
-     * Now extract the base-2**DIGIT_BIT digits of f into a multi-precision
-     * integer r. Preserve the invariant v = r * 2**rfac2 * FLT_RADIX**e by
-     * adjusting e.
-     */
+inline static void
+ComputeScale(int be,		/* Exponent part of number: d = bw * 2**be */
+	     int k,		/* Characteristic of log10(number) */
+	     int* b2,		/* OUTPUT: Power of 2 in the numerator */
+	     int* b5,		/* OUTPUT: Power of 5 in the numerator */
+	     int* s2,		/* OUTPUT: Power of 2 in the denominator */
+	     int* s5)		/* OUTPUT: Power of 5 in the denominator */
+{
 
-    a = f;
-    n = mantDIGIT;
-    mp_init_size(rPtr, n);
-    rPtr->used = n;
-    rPtr->sign = MP_ZPOS;
-    i = (mantBits % DIGIT_BIT);
-    if (i == 0) {
-	i = DIGIT_BIT;
+    /* 
+     * Scale numerator and denominator powers of 2 so that the
+     * input binary number is the ratio of integers
+     */
+    if (be <= 0) {
+	*b2 = 0;
+	*s2 = -be;
+    } else {
+	*b2 = be;
+	*s2 = 0;
     }
-    while (n > 0) {
-	a *= ldexp(1.0, i);
-	i = DIGIT_BIT;
-	rPtr->dp[--n] = (mp_digit) a;
-	a -= (mp_digit) a;
+
+    /* 
+     * Scale numerator and denominator so that the output decimal number
+     * is the ratio of integers
+     */
+    if (k >= 0) {
+	*b5 = 0;
+	*s5 = k;
+	*s2 += k;
+    } else {
+	*b2 -= k;
+	*b5 = -k;
+	*s5 = 0;
     }
-    *ePtr = e - DBL_MANT_DIG;
-    return (f == 1.0 / FLT_RADIX);
 }
 
 /*
- *----------------------------------------------------------------------
- *
- * TclInitDoubleConversion --
+ *-----------------------------------------------------------------------------
  *
- *	Initializes constants that are needed for conversions to and from
- *	'double'
+ * SetPrecisionLimits --
  *
- * Results:
- *	None.
+ *	Determines how many digits of significance should be computed
+ *	(and, hence, how much memory need be allocated) for formatting a
+ *	floating point number.
+ *
+ * Given that 'k' is floor(log10(x)):
+ * if 'shortest' format is used, there will be at most 18 digits in the result.
+ * if 'F' format is used, there will be at most 'ndigits' + k + 1 digits
+ * if 'E' format is used, there will be exactly 'ndigits' digits.
+ *
+ * Side effects:
+ *	Adjusts '*ndigitsPtr' to have a valid value.
+ *	Stores the maximum memory allocation needed in *iPtr.
+ *	Sets '*iLimPtr' to the limiting number of digits to convert if k
+ *	has been guessed correctly, and '*iLim1Ptr' to the limiting number
+ *	of digits to convert if k has been guessed to be one too high.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static void
+SetPrecisionLimits(int convType,
+				/* Type of conversion:
+				 *   TCL_DD_SHORTEST
+				 *   TCL_DD_STEELE0
+				 *   TCL_DD_E_FMT
+				 *   TCL_DD_F_FMT */
+		   int k,	/* Floor(log10(number to convert)) */
+		   int* ndigitsPtr,
+				/* IN/OUT: Number of digits requested
+				 *         (Will be adjusted if needed) */
+		   int* iPtr,	/* OUT: Maximum number of digits
+				 *      to return */
+		   int *iLimPtr,/* OUT: Number of digits of significance
+				 *      if the bignum method is used.*/
+		   int *iLim1Ptr)
+				/* OUT: Number of digits of significance
+				 *      if the quick method is used. */
+{
+    switch(convType) {
+    case TCL_DD_SHORTEST0:
+    case TCL_DD_STEELE0:
+	*iLimPtr = *iLim1Ptr = -1;
+	*iPtr = 18;
+	*ndigitsPtr = 0;
+	break;
+    case TCL_DD_E_FORMAT:
+	if (*ndigitsPtr <= 0) {
+	    *ndigitsPtr = 1;
+	}
+	*iLimPtr = *iLim1Ptr = *iPtr = *ndigitsPtr;
+	break;
+    case TCL_DD_F_FORMAT:
+	*iPtr = *ndigitsPtr + k + 1;
+	*iLimPtr = *iPtr;
+	*iLim1Ptr = *iPtr - 1;
+	if (*iPtr <= 0) {
+	    *iPtr = 1;
+	}
+	break;
+    }
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * BumpUp --
+ *
+ *	Increases a string of digits ending in a series of nines to
+ *	designate the next higher number.  xxxxb9999... -> xxxx(b+1)0000...
+ *
+ * Results:
+ *	Returns a pointer to the end of the adjusted string.
+ *
+ * Side effects:
+ *	In the case that the string consists solely of '999999', sets it
+ *	to "1" and moves the decimal point (*kPtr) one place to the right.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+
+inline static char*
+BumpUp(char* s,		    	/* Cursor pointing one past the end of the
+				 * string */ 
+       char* retval,		/* Start of the string of digits */
+       int* kPtr)		/* Position of the decimal point */
+{
+    while (*--s == '9') {
+	if (s == retval) {
+	    ++(*kPtr);
+	    *s = '1';
+	    return s+1;
+	}
+    }
+    ++*s;
+    ++s;
+    return s;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * AdjustRange --
+ *
+ *	Rescales a 'double' in preparation for formatting it using the
+ *	'quick' double-to-string method.
+ *
+ * Results:
+ *	Returns the precision that has been lost in the prescaling as
+ *	a count of units in the least significant place.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static int
+AdjustRange(double* dPtr,	/* INOUT: Number to adjust */
+	    int k)		/* IN: floor(log10(d)) */
+{
+    int ieps;			/* Number of roundoff errors that have
+				 * accumulated */
+    double d = *dPtr;		/* Number to adjust */
+    double ds;
+    int i, j, j1;
+
+    ieps = 2;
+
+    if (k > 0) {
+	/*
+	 * The number must be reduced to bring it into range.
+	 */
+	ds = tens[k & 0xf];
+	j = k >> 4;
+	if (j & BLETCH) {
+	    j &= (BLETCH-1);
+	    d /= bigtens[N_BIGTENS - 1];
+	    ieps++;
+	}
+	i = 0;
+	for (; j != 0; j>>=1) {
+	    if (j & 1) {
+		ds *= bigtens[i];
+		++ieps;
+	    }
+	    ++i;
+	}
+	d /= ds;
+    } else if ((j1 = -k) != 0) {
+	/*
+	 * The number must be increased to bring it into range
+	 */
+	d *= tens[j1 & 0xf];
+	i = 0;
+	for (j = j1>>4; j; j>>=1) {
+	    if (j & 1) {
+		ieps++;
+		d *= bigtens[i];
+	    }
+	    ++i;
+	}
+    }
+
+    *dPtr = d;
+    return ieps;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShorteningQuickFormat --
+ *
+ *	Returns a 'quick' format of a double precision number to a string
+ *	of digits, preferring a shorter string of digits if the shorter
+ *	string is still within 1/2 ulp of the number. 
+ *
+ * Results:
+ *	Returns the string of digits. Returns NULL if the 'quick' method
+ *	fails and the bignum method must be used.
+ *
+ * Side effects:
+ *	Stores the position of the decimal point at '*kPtr'.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+ShorteningQuickFormat(double d,	/* Number to convert */
+		      int k,	/* floor(log10(d)) */
+		      int ilim,	/* Number of significant digits to return */
+		      double eps,
+				/* Estimated roundoff error */
+		      char* retval,
+				/* Buffer to receive the digit string */
+		      int* kPtr)
+				/* Pointer to stash the position of
+				 * the decimal point */
+{
+    char* s = retval;		/* Cursor in the return value */
+    int digit;			/* Current digit */
+    int i;
+
+    eps = 0.5 / tens[ilim-1] - eps;
+    i = 0;
+    for (;;) {
+	/* Convert a digit */
+
+	digit = d;
+	d -= digit;
+	*s++ = '0' + digit;
+
+	/* 
+	 * Truncate the conversion if the string of digits is within
+	 * 1/2 ulp of the actual value.
+	 */
+
+	if (d < eps) {
+	    *kPtr = k;
+	    return s;
+	}
+	if ((1. - d) < eps) {
+	    *kPtr = k;
+	    return BumpUp(s, retval, kPtr);
+	}
+
+	/*
+	 * Bail out if the conversion fails to converge to a sufficiently
+	 * precise value
+	 */
+
+	if (++i >= ilim) {
+	    return NULL;
+	}
+
+	/*
+	 * Bring the next digit to the integer part.
+	 */
+
+	eps *= 10;
+	d *= 10.0;
+    }
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * StrictQuickFormat --
+ *
+ *	Convert a double precision number of a string of a precise number
+ *	of digits, using the 'quick' double precision method.
+ *
+ * Results:
+ *	Returns the digit string, or NULL if the bignum method must be
+ *	used to do the formatting.
+ *
+ * Side effects:
+ *	Stores the position of the decimal point in '*kPtr'.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+StrictQuickFormat(double d,	/* Number to convert */
+		  int k,	/* floor(log10(d)) */
+		  int ilim,	/* Number of significant digits to return */
+		  double eps,	/* Estimated roundoff error */
+		  char* retval,	/* Start of the digit string */
+		  int* kPtr)	/* Pointer to stash the position of
+				 * the decimal point */
+{
+    char* s = retval;		/* Cursor in the return value */
+    int digit;			/* Current digit of the answer */
+    int i;
+
+    eps *= tens[ilim-1];
+    i = 1;
+    for (;;) {
+	/* Extract a digit */
+	digit = d;
+	d -= digit;
+	if (d == 0.0) {
+	    ilim = i;
+	}
+	*s++ = '0' + digit;
+
+	/* 
+	 * When the given digit count is reached, handle trailing strings
+	 * of 0 and 9.
+	 */
+	if (i == ilim) {
+	    if (d > 0.5 + eps) {
+		*kPtr = k;
+		return BumpUp(s, retval, kPtr);
+	    } else if (d < 0.5 - eps) {
+		while (*--s == '0') {
+		    /* do nothing */
+		}
+		s++;
+		*kPtr = k;
+		return s;
+	    } else {
+		return NULL;
+	    }
+	}
+
+	/* Advance to the next digit */
+	++i;
+	d *= 10.0;
+    }
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * QuickConversion --
+ *
+ *	Converts a floating point number the 'quick' way, when only a limited
+ *	number of digits is required and floating point arithmetic can
+ *	therefore be used for the intermediate results.
+ *
+ * Results:
+ *	Returns the converted string, or NULL if the bignum method must
+ *	be used.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+QuickConversion(double d,	/* Number to format */
+		int k,		/* floor(log10(d)), approximately */
+		int k_check,	/* 0 if k is exact, 1 if it may be too high */
+		int flags,	/* Flags passed to dtoa:
+				 *    TCL_DD_SHORTEN_FLAG */
+		int len,	/* Length of the return value */
+		int ilim,	/* Number of digits to store */
+		int ilim1,	/* Number of digits to store if we
+				 * musguessed k */
+		int* decpt,	/* OUTPUT: Location of the decimal point */
+		char** endPtr)	/* OUTPUT: Pointer to the terminal null byte */
+{
+    int ieps;			/* Number of 1-ulp roundoff errors that have
+				 * accumulated in the calculation*/
+    Double eps;			/* Estimated roundoff error */
+    char* retval;		/* Returned string */
+    char* end;			/* Pointer to the terminal null byte in the
+				 * returned string */
+
+    /*
+     * Bring d into the range [1 .. 10)
+     */
+    ieps = AdjustRange(&d, k);
+
+    /*
+     * If the guessed value of k didn't get d into range, adjust it
+     * by one. If that leaves us outside the range in which quick format
+     * is accurate, bail out.
+     */
+    if (k_check && d < 1. && ilim > 0) {
+	if (ilim1 < 0) {
+	    return NULL;
+	}
+	ilim = ilim1;
+	--k;
+	d *= 10.0;
+	++ieps;
+    }
+
+    /*
+     * Compute estimated roundoff error
+     */
+    eps.d = ieps * d + 7.;
+    eps.w.word0 -= (FP_PRECISION-1) << EXP_SHIFT;
+
+    /*
+     * Handle the peculiar case where the result has no significant
+     * digits.
+     */
+    retval = ckalloc(len + 1);
+    if (ilim == 0) {
+	d -= 5.;
+	if (d > eps.d) {
+	    *retval = '1';
+	    *decpt = k;
+	    return retval;
+	} else if (d < -eps.d) {
+	    *decpt = k;
+	    return retval;
+	} else {
+	    ckfree(retval);
+	    return NULL;
+	}
+    }
+
+    /* Format the digit string */
+
+    if (flags & TCL_DD_SHORTEN_FLAG) {
+	end = ShorteningQuickFormat(d, k, ilim, eps.d, retval, decpt);
+    } else {
+	end = StrictQuickFormat(d, k, ilim, eps.d, retval, decpt);
+    }
+    if (end == NULL) {
+	ckfree(retval);
+	return NULL;
+    }
+    *end = '\0';
+    if (endPtr != NULL) {
+	*endPtr = end;
+    }
+    return retval;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * CastOutPowersOf2 --
+ *
+ *	Adjust the factors 'b2', 'm2', and 's2' to cast out common powers
+ *	of 2 from numerator and denominator in preparation for the 'bignum'
+ *	method of floating point conversion.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static void
+CastOutPowersOf2(int* b2,	/* Power of 2 to multiply the significand  */
+		 int* m2,	/* Power of 2 to multiply 1/2 ulp  */
+		 int* s2)	/* Power of 2 to multiply the common
+				 * denominator */
+{
+    int i;
+    if (*m2 > 0 && *s2 > 0) {	/* Find the smallest power of 2 in the
+				 * numerator */
+        if (*m2 < *s2) {	/* Find the lowest common denominatorr */
+	    i = *m2;
+	} else {
+	    i = *s2;
+	}
+	*b2 -= i;		/* Reduce to lowest terms */
+	*m2 -= i;
+	*s2 -= i;
+    }
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShorteningInt64Conversion --
+ *
+ *	Converts a double-precision number to the shortest string of
+ *	digits that reconverts exactly to the given number, or to
+ *	'ilim' digits if that will yield a shorter result. The numerator and
+ *	denominator in David Gay's conversion algorithm are known to fit
+ *	in Tcl_WideUInt, giving considerably faster arithmetic than mp_int's.
+ *
+ * Results:
+ *	Returns the string of significant decimal digits, in newly
+ *	allocated memory
+ *
+ * Side effects:
+ *	Stores the location of the decimal point in '*decpt' and the
+ *      location of the terminal null byte in '*endPtr'.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+ShorteningInt64Conversion(Double* dPtr,
+				/* Original number to convert */
+			  int convType,
+				/* Type of conversion (shortest, Steele,
+				   E format, F format) */
+			  Tcl_WideUInt bw,
+				/* Integer significand */
+			  int b2, int b5,
+				/* Scale factor for the significand
+				 * in the numerator */
+			  int m2plus, int m2minus, int m5,
+			  	/* Scale factors for 1/2 ulp in
+				 * the numerator (will be different if
+				 * bw == 1 */
+			  int s2, int s5,
+				/* Scale factors for the denominator */
+			  int k,
+				/* Number of output digits before the decimal
+				 * point */
+			  int len,
+				/* Number of digits to allocate */
+			  int ilim,
+				/* Number of digits to convert if b >= s */
+			  int ilim1,
+				/* Number of digits to convert if b < s */
+			  int* decpt,
+				/* OUTPUT: Position of the decimal point */
+			  char** endPtr)
+				/* OUTPUT: Position of the terminal '\0'
+				 *         at the end of the returned string */
+{
+    
+    char* retval = ckalloc(len + 1);
+				/* Output buffer */
+    Tcl_WideUInt b = (bw * wuipow5[b5]) << b2;
+				/* Numerator of the fraction being converted */
+    Tcl_WideUInt S = wuipow5[s5] << s2;
+				/* Denominator of the fraction being 
+				 * converted */
+    Tcl_WideUInt mplus, mminus;	/* Ranges for testing whether the result
+				 * is within roundoff of being exact */
+    int digit;			/* Current output digit */
+    char* s = retval;		/* Cursor in the output buffer */
+    int i;			/* Current position in the output buffer */
+
+    /* Adjust if the logarithm was guessed wrong */
+
+    if (b < S) {
+	b = 10 * b;
+	++m2plus; ++m2minus; ++m5;
+	ilim = ilim1;
+	--k;
+    }
+
+    /* Compute roundoff ranges */
+
+    mplus = wuipow5[m5] << m2plus;
+    mminus = wuipow5[m5] << m2minus;
+
+    /* Loop through the digits */
+
+    i = 1;
+    for (;;) {
+	digit = (int)(b / S);
+	if (digit > 10) {
+	    Tcl_Panic("wrong digit!");
+	}
+	b = b % S;
+
+	/* 
+	 * Does the current digit put us on the low side of the exact value
+	 * but within within roundoff of being exact?
+	 */
+	if (b < mplus
+	    || (b == mplus
+		&& convType != TCL_DD_STEELE0
+		&& (dPtr->w.word1 & 1) == 0)) {
+	    /*
+	     * Make sure we shouldn't be rounding *up* instead,
+	     * in case the next number above is closer
+	     */
+	    if (2 * b > S
+		|| (2 * b == S
+		    && (digit & 1) != 0)) {
+		++digit;
+		if (digit == 10) {
+		    *s++ = '9';
+		    s = BumpUp(s, retval, &k);
+		    break;
+		}
+	    }
+
+	    /* Stash the current digit */
+
+	    *s++ = '0' + digit;
+	    break;
+	}
+
+	/*
+	 * Does one plus the current digit put us within roundoff of the
+	 * number?
+	 */
+	if (b > S - mminus
+	    || (b == S - mminus
+		&& convType != TCL_DD_STEELE0
+		&& (dPtr->w.word1 & 1) == 0)) {
+	    if (digit == 9) {
+		*s++ = '9';
+		s = BumpUp(s, retval, &k);
+		break;
+	    }
+	    ++digit;
+	    *s++ = '0' + digit;
+	    break;
+	}
+
+	/*
+	 * Have we converted all the requested digits?
+	 */
+	*s++ = '0' + digit;
+	if (i == ilim) {
+	    if (2*b > S
+		|| (2*b == S && (digit & 1) != 0)) {
+		s = BumpUp(s, retval, &k);
+	    }
+	    break;
+	}
+	
+	/* Advance to the next digit */
+	
+	b = 10 * b;
+	mplus = 10 * mplus;
+	mminus = 10 * mminus;
+	++i;
+    }
+
+    /* 
+     * Endgame - store the location of the decimal point and the end of the
+     * string.
+     */
+    *s = '\0';
+    *decpt = k;
+    if (endPtr) {
+	*endPtr = s;
+    }
+    return retval;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * StrictInt64Conversion --
+ *
+ *	Converts a double-precision number to a fixed-length string of
+ *	'ilim' digits that reconverts exactly to the given number.
+ *	('ilim' should be replaced with 'ilim1' in the case where
+ *	log10(d) has been overestimated).  The numerator and
+ *	denominator in David Gay's conversion algorithm are known to fit
+ *	in Tcl_WideUInt, giving considerably faster arithmetic than mp_int's.
+ *
+ * Results:
+ *	Returns the string of significant decimal digits, in newly
+ *	allocated memory
+ *
+ * Side effects:
+ *	Stores the location of the decimal point in '*decpt' and the
+ *      location of the terminal null byte in '*endPtr'.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+StrictInt64Conversion(Double* dPtr,
+				/* Original number to convert */
+		      int convType,
+				/* Type of conversion (shortest, Steele,
+				   E format, F format) */
+		      Tcl_WideUInt bw,
+				/* Integer significand */
+		      int b2, int b5,
+				/* Scale factor for the significand
+				 * in the numerator */
+		      int s2, int s5,
+				/* Scale factors for the denominator */
+		      int k,
+				/* Number of output digits before the decimal
+				 * point */
+		      int len,
+				/* Number of digits to allocate */
+		      int ilim,
+				/* Number of digits to convert if b >= s */
+		      int ilim1,
+				/* Number of digits to convert if b < s */
+		      int* decpt,
+				/* OUTPUT: Position of the decimal point */
+		      char** endPtr)
+				/* OUTPUT: Position of the terminal '\0'
+				 *         at the end of the returned string */
+{
+    
+    char* retval = ckalloc(len + 1);
+				/* Output buffer */
+    Tcl_WideUInt b = (bw * wuipow5[b5]) << b2;
+				/* Numerator of the fraction being converted */
+    Tcl_WideUInt S = wuipow5[s5] << s2;
+				/* Denominator of the fraction being 
+				 * converted */
+    int digit;			/* Current output digit */
+    char* s = retval;		/* Cursor in the output buffer */
+    int i;			/* Current position in the output buffer */
+
+    /* Adjust if the logarithm was guessed wrong */
+
+    if (b < S) {
+	b = 10 * b;
+	ilim = ilim1;
+	--k;
+    }
+
+    /* Loop through the digits */
+
+    i = 1;
+    for (;;) {
+	digit = (int)(b / S);
+	if (digit > 10) {
+	    Tcl_Panic("wrong digit!");
+	}
+	b = b % S;
+
+	/*
+	 * Have we converted all the requested digits?
+	 */
+	*s++ = '0' + digit;
+	if (i == ilim) {
+	    if (2*b > S
+		|| (2*b == S && (digit & 1) != 0)) {
+		s = BumpUp(s, retval, &k);
+	    }
+	    break;
+	}
+	
+	/* Advance to the next digit */
+	
+	b = 10 * b;
+	++i;
+    }
+
+    /* 
+     * Endgame - store the location of the decimal point and the end of the
+     * string.
+     */
+    *s = '\0';
+    *decpt = k;
+    if (endPtr) {
+	*endPtr = s;
+    }
+    return retval;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShouldBankerRoundUpPowD --
+ *
+ *	Test whether bankers' rounding should round a digit up. Assumption
+ *	is made that the denominator of the fraction being tested is
+ *	a power of 2**DIGIT_BIT.
+ *
+ * Results:
+ *	Returns 1 iff the fraction is more than 1/2, or if the fraction
+ *	is exactly 1/2 and the digit is odd.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static int
+ShouldBankerRoundUpPowD(mp_int* b,
+				/* Numerator of the fraction */
+			int sd,	/* Denominator is 2**(sd*DIGIT_BIT) */
+			int isodd)
+				/* 1 if the digit is odd, 0 if even */
+{
+    int i;
+    const static mp_digit topbit = (1<<(DIGIT_BIT-1));
+    if (b->used < sd || (b->dp[sd-1] & topbit) == 0) {
+	return 0;
+    }
+    if (b->dp[sd-1] != topbit) {
+	return 1;
+    }
+    for (i = sd-2; i >= 0; --i) {
+	if (b->dp[i] != 0) {
+	    return 1;
+	}
+    }
+    return isodd;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShouldBankerRoundUpToNextPowD --
+ *
+ *	Tests whether bankers' rounding will round down in the
+ *	"denominator is a power of 2**MP_DIGIT" case.
+ *
+ * Results:
+ *	Returns 1 if the rounding will be performed - which increases the
+ *	digit by one - and 0 otherwise.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static int
+ShouldBankerRoundUpToNextPowD(mp_int* b,
+				/* Numerator of the fraction */
+			      mp_int* m,
+				/* Numerator of the rounding tolerance */
+			      int sd, 
+				/* Common denominator is 2**(sd*DIGIT_BIT) */
+			      int convType,
+				/* Conversion type: STEELE defeats 
+				 * round-to-even (Not sure why one wants to
+				 * do this; I copied it from Gay) FIXME */
+			      int isodd,
+				/* 1 if the integer significand is odd */
+			      mp_int* temp)
+				/* Work area for the calculation */
+{
+    int i;
+
+    /* 
+     * Compare B and S-m -- which is the same as comparing B+m and S --
+     * which we do by computing b+m and doing a bitwhack compare against
+     * 2**(DIGIT_BIT*sd)
+     */
+    mp_add(b, m, temp);
+    if (temp->used <= sd) {	/* too few digits to be > S */
+	return 0;
+    }
+    if (temp->used > sd+1 || temp->dp[sd] > 1) {
+				/* >= 2s */
+	return 1;
+    }
+    for (i = sd-1; i >= 0; --i) {
+				/* check for ==s */
+	if (temp->dp[i] != 0) {	/* > s */
+	    return 1;
+	}
+    }
+    if (convType == TCL_DD_STEELE0) {
+				/* biased rounding */
+	return 0;
+    }
+    return isodd;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShorteningBignumConversionPowD --
+ *
+ *	Converts a double-precision number to the shortest string of
+ *	digits that reconverts exactly to the given number, or to
+ *	'ilim' digits if that will yield a shorter result. The denominator
+ *	in David Gay's conversion algorithm is known to be a power of
+ *	2**DIGIT_BIT, and hence the division in the main loop may be replaced
+ *	by a digit shift and mask.
+ *
+ * Results:
+ *	Returns the string of significant decimal digits, in newly
+ *	allocated memory
+ *
+ * Side effects:
+ *	Stores the location of the decimal point in '*decpt' and the
+ *      location of the terminal null byte in '*endPtr'.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+ShorteningBignumConversionPowD(Double* dPtr,
+				/* Original number to convert */
+			       int convType,
+				/* Type of conversion (shortest, Steele,
+				   E format, F format) */
+			       Tcl_WideUInt bw,
+				/* Integer significand */
+			       int b2, int b5,
+				/* Scale factor for the significand
+				 * in the numerator */
+			       int m2plus, int m2minus, int m5,
+			  	/* Scale factors for 1/2 ulp in
+				 * the numerator (will be different if
+				 * bw == 1 */
+			       int sd,
+				/* Scale factor for the denominator */
+			       int k,
+				/* Number of output digits before the decimal
+				 * point */
+			       int len,
+				/* Number of digits to allocate */
+			       int ilim,
+				/* Number of digits to convert if b >= s */
+			       int ilim1,
+				/* Number of digits to convert if b < s */
+			       int* decpt,
+				/* OUTPUT: Position of the decimal point */
+			       char** endPtr)
+				/* OUTPUT: Position of the terminal '\0'
+				 *         at the end of the returned string */
+{
+    
+    char* retval = ckalloc(len + 1);
+				/* Output buffer */
+    mp_int b;			/* Numerator of the fraction being converted */
+    mp_int mplus, mminus;	/* Bounds for roundoff */
+    mp_digit digit;		/* Current output digit */
+    char* s = retval;		/* Cursor in the output buffer */
+    int i;			/* Index in the output buffer */
+    mp_int temp;
+    int r1;
+
+    /* 
+     * b = bw * 2**b2 * 5**b5
+     * mminus = 5**m5
+     */
+
+    TclBNInitBignumFromWideUInt(&b, bw);
+    mp_init_set_int(&mminus, 1);
+    MulPow5(&b, b5, &b);
+    mp_mul_2d(&b, b2, &b);
+
+    /* Adjust if the logarithm was guessed wrong */
+
+    if (b.used <= sd) {
+	mp_mul_d(&b, 10, &b);
+	++m2plus; ++m2minus; ++m5;
+	ilim = ilim1;
+	--k;
+    }
+
+    /*
+     * mminus = 5**m5 * 2**m2minus
+     * mplus = 5**m5 * 2**m2plus
+     */
+
+    mp_mul_2d(&mminus, m2minus, &mminus);
+    MulPow5(&mminus, m5, &mminus);
+    if (m2plus > m2minus) {
+	mp_init_copy(&mplus, &mminus);
+	mp_mul_2d(&mplus, m2plus-m2minus, &mplus);
+    }
+    mp_init(&temp);
+
+    /* Loop through the digits. Do division and mod by s == 2**(sd*DIGIT_BIT)
+     * by mp_digit extraction */
+
+    i = 0;
+    for (;;) {
+	if (b.used <= sd) {
+	    digit = 0;
+	} else {
+	    digit = b.dp[sd];
+	    if (b.used > sd+1 || digit >= 10) {
+		Tcl_Panic("wrong digit!");
+	    }
+	    --b.used; mp_clamp(&b);
+	}
+
+	/* 
+	 * Does the current digit put us on the low side of the exact value
+	 * but within within roundoff of being exact?
+	 */
+	
+	r1 = mp_cmp_mag(&b, (m2plus > m2minus)? &mplus : &mminus);
+	if (r1 == MP_LT
+	    || (r1 == MP_EQ
+		&& convType != TCL_DD_STEELE0
+		&& (dPtr->w.word1 & 1) == 0)) {
+	    /*
+	     * Make sure we shouldn't be rounding *up* instead,
+	     * in case the next number above is closer
+	     */
+	    if (ShouldBankerRoundUpPowD(&b, sd, digit&1)) {
+		++digit;
+		if (digit == 10) {
+		    *s++ = '9';
+		    s = BumpUp(s, retval, &k);
+		    break;
+		}
+	    }
+
+	    /* Stash the last digit */
+
+	    *s++ = '0' + digit;
+	    break;
+	}
+
+	/*
+	 * Does one plus the current digit put us within roundoff of the
+	 * number?
+	 */
+	
+	if (ShouldBankerRoundUpToNextPowD(&b, &mminus, sd, 
+					   convType, dPtr->w.word1 & 1,
+					   &temp)) {
+	    if (digit == 9) {
+		*s++ = '9';
+		s = BumpUp(s, retval, &k);
+		break;
+	    }
+	    ++digit;
+	    *s++ = '0' + digit;
+	    break;
+	}
+
+	/*
+	 * Have we converted all the requested digits?
+	 */
+	*s++ = '0' + digit;
+	if (i == ilim) {
+	    if (ShouldBankerRoundUpPowD(&b, sd, digit&1)) {
+		s = BumpUp(s, retval, &k);
+	    }
+	    break;
+	}
+	
+	/* Advance to the next digit */
+	
+	mp_mul_d(&b, 10, &b);
+	mp_mul_d(&mminus, 10, &mminus);
+	if (m2plus > m2minus) {
+	    mp_mul_2d(&mminus, m2plus-m2minus, &mplus);
+	}
+	++i;
+    }
+
+    /* 
+     * Endgame - store the location of the decimal point and the end of the
+     * string.
+     */
+    if (m2plus > m2minus) {
+	mp_clear(&mplus);
+    }
+    mp_clear_multi(&b, &mminus, &temp, NULL);
+    *s = '\0';
+    *decpt = k;
+    if (endPtr) {
+	*endPtr = s;
+    }
+    return retval;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * StrictBignumConversionPowD --
+ *
+ *	Converts a double-precision number to a fixed-lengt string of
+ *	'ilim' digits (or 'ilim1' if log10(d) has been overestimated.)
+ *	The denominator in David Gay's conversion algorithm is known to
+ *	be a power of 2**DIGIT_BIT, and hence the division in the main 
+ *	loop may be replaced by a digit shift and mask.
+ *
+ * Results:
+ *	Returns the string of significant decimal digits, in newly
+ *	allocated memory.
+ *
+ * Side effects:
+ *	Stores the location of the decimal point in '*decpt' and the
+ *      location of the terminal null byte in '*endPtr'.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+StrictBignumConversionPowD(Double* dPtr,
+				/* Original number to convert */
+			   int convType,
+				/* Type of conversion (shortest, Steele,
+				   E format, F format) */
+			   Tcl_WideUInt bw,
+				/* Integer significand */
+			   int b2, int b5,
+				/* Scale factor for the significand
+				 * in the numerator */
+			   int sd,
+				/* Scale factor for the denominator */
+			   int k,
+				/* Number of output digits before the decimal
+				 * point */
+			   int len,
+				/* Number of digits to allocate */
+			   int ilim,
+				/* Number of digits to convert if b >= s */
+			   int ilim1,
+				/* Number of digits to convert if b < s */
+			   int* decpt,
+				/* OUTPUT: Position of the decimal point */
+			   char** endPtr)
+				/* OUTPUT: Position of the terminal '\0'
+				 *         at the end of the returned string */
+{
+    
+    char* retval = ckalloc(len + 1);
+				/* Output buffer */
+    mp_int b;			/* Numerator of the fraction being converted */
+    mp_digit digit;		/* Current output digit */
+    char* s = retval;		/* Cursor in the output buffer */
+    int i;			/* Index in the output buffer */
+    mp_int temp;
+
+    /* 
+     * b = bw * 2**b2 * 5**b5
+     */
+
+    TclBNInitBignumFromWideUInt(&b, bw);
+    MulPow5(&b, b5, &b);
+    mp_mul_2d(&b, b2, &b);
+
+    /* Adjust if the logarithm was guessed wrong */
+
+    if (b.used <= sd) {
+	mp_mul_d(&b, 10, &b);
+	ilim = ilim1;
+	--k;
+    }
+    mp_init(&temp);
+
+    /* 
+     * Loop through the digits. Do division and mod by s == 2**(sd*DIGIT_BIT)
+     * by mp_digit extraction 
+     */
+
+    i = 1;
+    for (;;) {
+	if (b.used <= sd) {
+	    digit = 0;
+	} else {
+	    digit = b.dp[sd];
+	    if (b.used > sd+1 || digit >= 10) {
+		Tcl_Panic("wrong digit!");
+	    }
+	    --b.used; mp_clamp(&b);
+	}
+
+	/*
+	 * Have we converted all the requested digits?
+	 */
+	*s++ = '0' + digit;
+	if (i == ilim) {
+	    if (ShouldBankerRoundUpPowD(&b, sd, digit&1)) {
+		s = BumpUp(s, retval, &k);
+	    }
+	    break;
+	}
+	
+	/* Advance to the next digit */
+	
+	mp_mul_d(&b, 10, &b);
+	++i;
+    }
+
+    /* 
+     * Endgame - store the location of the decimal point and the end of the
+     * string.
+     */
+    mp_clear_multi(&b, &temp, NULL);
+    *s = '\0';
+    *decpt = k;
+    if (endPtr) {
+	*endPtr = s;
+    }
+    return retval;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShouldBankerRoundUp --
+ *
+ *	Tests whether a digit should be rounded up or down when finishing
+ *	bignum-based floating point conversion.
+ *
+ * Results:
+ *	Returns 1 if the number needs to be rounded up, 0 otherwise.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static int
+ShouldBankerRoundUp(mp_int* twor,
+				/* 2x the remainder from thd division that
+				 * produced the last digit */
+		    mp_int* S,	/* Denominator */
+		    int isodd)	/* Flag == 1 if the last digit is odd */
+{
+    int r = mp_cmp_mag(twor, S);
+    switch (r) {
+    case MP_LT:
+	return 0;
+    case MP_EQ:
+	return isodd;
+    case MP_GT:
+	return 1;
+    }
+    Tcl_Panic("in ShouldBankerRoundUp, trichotomy fails!");
+    return 0;
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShouldBankerRoundUpToNext --
+ *
+ *	Tests whether the remainder is great enough to force rounding
+ *	to the next higher digit.
+ *
+ * Results:
+ *	Returns 1 if the number should be rounded up, 0 otherwise.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static int
+ShouldBankerRoundUpToNext(mp_int* b,
+				/* Remainder from the division that produced
+				 * the last digit. */
+			  mp_int* m,
+				/* Numerator of the rounding tolerance */
+			  mp_int* S,
+				/* Denominator */
+			  int convType,
+				/* Conversion type: STEELE0 defeats
+				 * round-to-even. (Not sure why one would
+				 * want this; I coped it from Gay. FIXME */
+			  int isodd,
+				/* 1 if the integer significand is odd */
+			  mp_int* temp)
+				/* Work area needed for the calculation */
+{
+    int r;
+    /* Compare b and S-m: this is the same as comparing B+m and S. */
+    mp_add(b, m, temp);
+    r = mp_cmp_mag(temp, S);
+    switch(r) {
+    case MP_LT:
+	return 0;
+    case MP_EQ:
+	if (convType == TCL_DD_STEELE0) {
+	    return 0;
+	} else {
+	    return isodd;
+	}
+    case MP_GT:
+	return 1;
+    }
+    Tcl_Panic("in ShouldBankerRoundUpToNext, trichotomy fails!");
+    return 0;
+}
+		 
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * ShorteningBignumConversion --
+ *
+ *	Convert a floating point number to a variable-length digit string
+ *	using the multiprecision method.
+ *
+ * Results:
+ *	Returns the string of digits.
+ *
+ * Side effects:
+ *	Stores the position of the decimal point in *decpt.
+ *	Stores a pointer to the end of the number in *endPtr.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+ShorteningBignumConversion(Double* dPtr,
+				/* Original number being converted */
+			   int convType,
+				/* Conversion type */
+			   Tcl_WideUInt bw,
+				/* Integer significand and exponent */
+			   int b2,
+				/* Scale factor for the significand */
+			   int m2plus, int m2minus,
+				/* Scale factors for 1/2 ulp in numerator */
+			   int s2, int s5,
+				/* Scale factors for denominator */
+			   int k,
+				/* Guessed position of the decimal point */
+			   int len,
+				/* Size of the digit buffer to allocate */
+			   int ilim,
+				/* Number of digits to convert if b >= s */
+			   int ilim1,
+				/* Number of digits to convert if b < s */
+			   int* decpt,
+				/* OUTPUT: Position of the decimal point */
+			   char** endPtr)
+				/* OUTPUT: Pointer to the end of the number */
+{
+    char* retval = ckalloc(len+1);
+				/* Buffer of digits to return */
+    char* s = retval;		/* Cursor in the return value */
+    mp_int b;			/* Numerator of the result */
+    mp_int mminus;		/* 1/2 ulp below the result */
+    mp_int mplus;		/* 1/2 ulp above the result */
+    mp_int S;			/* Denominator of the result */
+    mp_int dig;			/* Current digit of the result */
+    int digit;			/* Current digit of the result */
+    mp_int temp;		/* Work area */
+    int minit = 1;		/* Fudge factor for when we misguess k */
+    int i;
+    int r1;
+
+    /*
+     * b = bw * 2**b2 * 5**b5
+     * S = 2**s2 * 5*s5
+     */
+
+    TclBNInitBignumFromWideUInt(&b, bw);
+    mp_mul_2d(&b, b2, &b);
+    mp_init_set_int(&S, 1);
+    MulPow5(&S, s5, &S); mp_mul_2d(&S, s2, &S);
+
+    /*
+     * Handle the case where we guess the position of the decimal point 
+     * wrong. 
+     */
+    
+    if (mp_cmp_mag(&b, &S) == MP_LT) {
+	mp_mul_d(&b, 10, &b);
+	minit = 10;
+	ilim =ilim1;
+	--k;
+    }
+
+    /* mminus = 2**m2minus * 5**m5 */
+
+    mp_init_set_int(&mminus, minit);
+    mp_mul_2d(&mminus, m2minus, &mminus);
+    if (m2plus > m2minus) {
+	mp_init_copy(&mplus, &mminus);
+	mp_mul_2d(&mplus, m2plus-m2minus, &mplus);
+    }
+    mp_init(&temp);
+
+    /* Loop through the digits */
+
+    mp_init(&dig);
+    i = 1;
+    for (;;) {
+	mp_div(&b, &S, &dig, &b);
+	if (dig.used > 1 || dig.dp[0] >= 10) {
+	    Tcl_Panic("wrong digit!");
+	}
+	digit = dig.dp[0];
+
+	/* 
+	 * Does the current digit leave us with a remainder small enough to
+	 * round to it?
+	 */
+
+	r1 = mp_cmp_mag(&b, (m2plus > m2minus)? &mplus : &mminus);
+	if (r1 == MP_LT
+	    || (r1 == MP_EQ
+		&& convType != TCL_DD_STEELE0
+		&& (dPtr->w.word1 & 1) == 0)) {
+		mp_mul_2d(&b, 1, &b);
+	    if (ShouldBankerRoundUp(&b, &S, digit&1)) {
+		++digit;
+		if (digit == 10) {
+		    *s++ = '9';
+		    s = BumpUp(s, retval, &k);
+		    break;
+		}
+	    }
+	    *s++ = '0' + digit;
+	    break;
+	}
+
+	/*
+	 * Does the current digit leave us with a remainder large enough
+	 * to commit to rounding up to the next higher digit?
+	 */
+
+	if (ShouldBankerRoundUpToNext(&b, &mminus, &S, convType,
+				      dPtr->w.word1 & 1, &temp)) {
+	    ++digit;
+	    if (digit == 10) {
+		*s++ = '9';
+		s = BumpUp(s, retval, &k);
+		break;
+	    }
+	    *s++ = '0' + digit;
+	    break;
+	}
+
+	/* Have we converted all the requested digits? */
+
+	*s++ = '0' + digit;
+	if (i == ilim) {
+	    mp_mul_2d(&b, 1, &b);
+	    if (ShouldBankerRoundUp(&b, &S, digit&1)) {
+		s = BumpUp(s, retval, &k);  
+	    }
+	    break;
+	}
+
+	/* Advance to the next digit */
+
+	if (s5 > 0) {
+
+	    /* Can possibly shorten the denominator */
+	    mp_mul_2d(&b, 1, &b);
+	    mp_mul_2d(&mminus, 1, &mminus);
+	    if (m2plus > m2minus) {
+		mp_mul_2d(&mplus, 1, &mplus);
+	    }
+	    mp_div_d(&S, 5, &S, NULL);
+	    --s5;
+	    /* 
+	     * TODO: It might possibly be a win to fall back to
+	     *       int64 arithmetic here if S < 2**64/10. But it's
+	     *       a win only for a fairly narrow range of magnitudes
+	     *       so perhaps not worth bothering. We already know that
+	     *       we shorten the denominator by at least 1 mp_digit, perhaps
+	     *       2. as we do the conversion for 17 digits of significance.
+	     * Possible savings:
+	     * 10**26   1 trip through loop before fallback possible
+	     * 10**27   1 trip
+	     * 10**28   2 trips     
+	     * 10**29   3 trips
+	     * 10**30   4 trips
+	     * 10**31   5 trips
+	     * 10**32   6 trips
+	     * 10**33   7 trips
+	     * 10**34   8 trips
+	     * 10**35   9 trips
+	     * 10**36  10 trips
+	     * 10**37  11 trips
+	     * 10**38  12 trips
+	     * 10**39  13 trips
+	     * 10**40  14 trips
+	     * 10**41  15 trips
+	     * 10**42  16 trips
+	     * thereafter  no gain.
+	     */
+	} else {
+	    mp_mul_d(&b, 10, &b);
+	    mp_mul_d(&mminus, 10, &mminus);
+	    if (m2plus > m2minus) {
+		mp_mul_2d(&mplus, 10, &mplus);
+	    }
+	}
+
+	++i;
+    }
+
+
+    /* 
+     * Endgame - store the location of the decimal point and the end of the
+     * string.
+     */
+    if (m2plus > m2minus) {
+	mp_clear(&mplus);
+    }
+    mp_clear_multi(&b, &mminus, &temp, NULL);
+    *s = '\0';
+    *decpt = k;
+    if (endPtr) {
+	*endPtr = s;
+    }
+    return retval;
+
+}
+		 
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * StrictBignumConversion --
+ *
+ *	Convert a floating point number to a fixed-length digit string
+ *	using the multiprecision method.
+ *
+ * Results:
+ *	Returns the string of digits.
+ *
+ * Side effects:
+ *	Stores the position of the decimal point in *decpt.
+ *	Stores a pointer to the end of the number in *endPtr.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+inline static char*
+StrictBignumConversion(Double* dPtr,
+				/* Original number being converted */
+		       int convType,
+				/* Conversion type */
+		       Tcl_WideUInt bw,
+				/* Integer significand and exponent */
+		       int b2,	/* Scale factor for the significand */
+		       int s2, int s5,
+				/* Scale factors for denominator */
+		       int k,	/* Guessed position of the decimal point */
+		       int len,	/* Size of the digit buffer to allocate */
+		       int ilim,
+				/* Number of digits to convert if b >= s */
+		       int ilim1,
+				/* Number of digits to convert if b < s */
+		       int* decpt,
+				/* OUTPUT: Position of the decimal point */
+		       char** endPtr)
+				/* OUTPUT: Pointer to the end of the number */
+{
+    char* retval = ckalloc(len+1);
+				/* Buffer of digits to return */
+    char* s = retval;		/* Cursor in the return value */
+    mp_int b;			/* Numerator of the result */
+    mp_int S;			/* Denominator of the result */
+    mp_int dig;			/* Current digit of the result */
+    int digit;			/* Current digit of the result */
+    mp_int temp;		/* Work area */
+    int g;			/* Size of the current digit groun */
+    int i, j;
+    
+    /*
+     * b = bw * 2**b2 * 5**b5
+     * S = 2**s2 * 5*s5
+     */
+
+    TclBNInitBignumFromWideUInt(&b, bw);
+    mp_mul_2d(&b, b2, &b);
+    mp_init_set_int(&S, 1);
+    MulPow5(&S, s5, &S); mp_mul_2d(&S, s2, &S);
+
+    /*
+     * Handle the case where we guess the position of the decimal point 
+     * wrong. 
+     */
+    
+    if (mp_cmp_mag(&b, &S) == MP_LT) {
+	mp_mul_d(&b, 10, &b);
+	ilim =ilim1;
+	--k;
+    }
+    mp_init(&temp);
+
+    /* Convert the leading digit */
+
+    mp_init(&dig);
+    i = 0;
+    mp_div(&b, &S, &dig, &b);
+    if (dig.used > 1 || dig.dp[0] >= 10) {
+	    Tcl_Panic("wrong digit!");
+	}
+    digit = dig.dp[0];
+
+    /* Is a single digit all that was requested? */
+
+    *s++ = '0' + digit;
+    if (++i >= ilim) {
+	mp_mul_2d(&b, 1, &b);
+	if (ShouldBankerRoundUp(&b, &S, digit&1)) {
+	    s = BumpUp(s, retval, &k);  
+	}
+    } else {
+
+	for (;;) {
+
+	    /* Shift by a group of digits. */
+
+	    g = ilim - i;
+	    if (g > DIGIT_GROUP) {
+		g = DIGIT_GROUP;
+	    }
+	    if (s5 >= g) {
+		mp_div_d(&S, dpow5[g], &S, NULL);
+		s5 -= g;
+	    } else if (s5 > 0) {
+		mp_div_d(&S, dpow5[s5], &S, NULL);
+		mp_mul_d(&b, dpow5[g - s5], &b);
+		s5 = 0;
+	    } else {
+		mp_mul_d(&b, dpow5[g], &b);
+	    }
+	    mp_mul_2d(&b, g, &b);
+	    
+	    /*
+	     * As with the shortening bignum conversion, it's possible at
+	     * this point that we will have reduced the denominator to
+	     * less than 2**64/10, at which point it would be possible to
+	     * fall back to to int64 arithmetic. But the potential payoff
+	     * is tremendously less - unless we're working in F format -
+	     * because we know that three groups of digits will always
+	     * suffice for %#.17e, the longest format that doesn't introduce
+	     * empty precision.
+	     */
+
+	    /* Extract the next digit */
+	    
+	    mp_div(&b, &S, &dig, &b);
+	    if (dig.used > 1) {
+		Tcl_Panic("wrong digit!");
+	    }
+	    digit = dig.dp[0];
+	    for (j = g-1; j >= 0; --j) {
+		int t = itens[j];
+		*s++ = digit / t + '0';
+		digit %= t;
+	    }
+	    i += g;
+	    
+	    /* Have we converted all the requested digits? */
+	    
+	    if (i == ilim) {
+		mp_mul_2d(&b, 1, &b);
+		if (ShouldBankerRoundUp(&b, &S, digit&1)) {
+		    s = BumpUp(s, retval, &k);  
+		}
+		break;
+	    }
+	}
+    }
+    /* 
+     * Endgame - store the location of the decimal point and the end of the
+     * string.
+     */
+    mp_clear_multi(&b, &temp, NULL);
+    *s = '\0';
+    *decpt = k;
+    if (endPtr) {
+	*endPtr = s;
+    }
+    return retval;
+
+}
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * TclDoubleDigits --
+ *
+ *	Core of Tcl's conversion of double-precision floating point numbers
+ *	to decimal.
+ *
+ * Results:
+ *	Returns a newly-allocated string of digits.
+ *
+ * Side effects:
+ *	Sets *decpt to the index of the character in the string before the
+ *	place that the decimal point should go. If 'endPtr' is not NULL,
+ *	sets endPtr to point to the terminating '\0' byte of the string.
+ *	Sets *sign to 1 if a minus sign should be printed with the number,
+ *	or 0 if a plus sign (or no sign) should appear.
+ *
+ * This function is a service routine that produces the string of digits
+ * for floating-point-to-decimal conversion. It can do a number of things
+ * according to the 'flags' argument. Valid values for 'flags' include:
+ *	TCL_DD_SHORTEST - This is the default for floating point conversion
+ *		if ::tcl_precision is 0. It constructs the shortest string
+ *		of digits that will reconvert to the given number when scanned.
+ *		For floating point numbers that are exactly between two
+ *		decimal numbers, it resolves using the 'round to even' rule.
+ *		With this value, the 'ndigits' parameter is ignored.
+ *	TCL_DD_STEELE - This value is not recommended and may be removed
+ *		in the future. It follows the conversion algorithm outlined
+ *		in "How to Print Floating-Point Numbers Accurately" by
+ *		Guy L. Steele, Jr. and Jon L. White [Proc. ACM SIGPLAN '90, 
+ *		pp. 112-126]. This rule has the effect of rendering 1e23
+ *		as 9.9999999999999999e22 - which is a 'better' approximation
+ *		in the sense that it will reconvert correctly even if
+ *		a subsequent input conversion is 'round up' or 'round down'
+ *		rather than 'round to nearest', but is surprising otherwise.
+ *	TCL_DD_E_FORMAT - This value is used to prepare numbers for %e
+ *		format conversion (or for default floating->string if
+ *		tcl_precision is not 0). It constructs a string of at most
+ *		'ndigits' digits, choosing the one that is closest to the
+ *		given number (and resolving ties with 'round to even').
+ *		It is allowed to return fewer than 'ndigits' if the number
+ *		converts exactly; if the TCL_DD_E_FORMAT|TCL_DD_SHORTEN_FLAG 
+ *		is supplied instead, it is also allowed to return fewer digits
+ *		if the shorter string will still reconvert to the given
+ *		input number.
+ *	TCL_DD_F_FORMAT - This value is used to prepare numbers for %f
+ *		format conversion. It requests that conversion proceed until
+ *		'ndigits' digits after the decimal point have been converted.
+ *		It is possible for this format to result in a zero-length 
+ *		string if the number is sufficiently small. Again, it
+ *		is permissible for TCL_DD_F_FORMAT to return fewer digits
+ *		for a number that converts exactly, and changing the
+ *		argument to TCL_DD_F_FORMAT|TCL_DD_SHORTEN_FLAG will allow
+ *		the routine also to return fewer digits if the shorter string
+ *		will still reconvert without loss to the given input number.
+ *
+ *	To any of these flags may be OR'ed TCL_DD_NO_QUICK; this flag
+ *	requires all calculations to be done in exact arithmetic. Normally,
+ *	E and F format with fewer than about 14 digits will be done with
+ *	a quick floating point approximation and fall back on the exact
+ *	arithmetic only if the input number is close enough to the
+ *	midpoint between two decimal strings that more precision is needed
+ *	to resolve which string is correct.
+ *
+ * The value stored in the 'decpt' argument on return may be negative 
+ * (indicating that the decimal point falls to the left of the string) 
+ * or greater than the length of the string.  In addition, the value -9999
+ * is used as a sentinel to indicate that the string is one of the special
+ * values "Infinity" and "NaN", and that no decimal point should be inserted.
+ * 
+ *-----------------------------------------------------------------------------
+ */
+char*
+TclDoubleDigits(double dv,	/* Number to convert */
+		int ndigits,	/* Number of digits requested */
+		int flags,	/* Conversion flags */
+		int* decpt,	/* OUTPUT: Position of the decimal point */
+		int* sign,	/* OUTPUT: 1 if the result is negative */
+		char** endPtr)	/* OUTPUT: If not NULL, receives a pointer
+				 *         to one character beyond the end
+				 *         of the returned string */
+{
+    int convType = (flags & TCL_DD_CONVERSION_TYPE_MASK);
+				/* Type of conversion being performed
+				 * TCL_DD_SHORTEST0
+				 * TCL_DD_STEELE0
+				 * TCL_DD_E_FORMAT
+				 * TCL_DD_F_FORMAT */
+    Double d;			/* Union for deconstructing doubles */
+    Tcl_WideUInt bw;		/* Integer significand */
+    int be;			/* Power of 2 by which b must be multiplied */
+    int bbits;			/* Number of bits needed to represent b */
+    int denorm;			/* Flag == 1 iff the input number was
+				 * denormalized */
+    int k;			/* Estimate of floor(log10(d)) */
+    int k_check;		/* Flag == 1 if d is near enough to a 
+				 * power of ten that k must be checked */
+    int b2, b5, s2, s5;		/* Powers of 2 and 5 in the numerator and
+				 * denominator of intermediate results */
+    int ilim, ilim1;
+    char* retval;		/* Return value from this function */
+    int i;
+
+    /* Put the input number into a union for bit-whacking */
+
+    d.d = dv;
+
+    /* 
+     * Handle the cases of negative numbers (by taking the absolute value:
+     * this includes -Inf and -NaN!), infinity, Not a Number, and zero.
+     */
+
+    TakeAbsoluteValue(&d, sign);
+    if ((d.w.word0 & EXP_MASK) == EXP_MASK) {
+	return FormatInfAndNaN(&d, decpt, endPtr);
+    }
+    if (d.d == 0.0) {
+	return FormatZero(decpt, endPtr);
+    }
+
+    /* 
+     * Unpack the floating point into a wide integer and an exponent.
+     * Determine the number of bits that the big integer requires, and
+     * compute a quick approximation (which may be one too high) of
+     * ceil(log10(d.d)).
+     */
+    denorm = ((d.w.word0 & EXP_MASK) == 0);
+    DoubleToExpAndSig(d.d, &bw, &be, &bbits);
+    k = ApproximateLog10(bw, be, bbits);
+    k = BetterLog10(d.d, k, &k_check);
+
+    /* At this point, we have:
+     *	  d is the number to convert.
+     *    bw are significand and exponent: d == bw*2**be, 
+     *    bbits is the length of bw: 2**bbits-1 <= bw < 2**bbits
+     *	  k is either ceil(log10(d)) or ceil(log10(d))+1. k_check is 0
+     *      if we know that k is exactly ceil(log10(d)) and 1 if we need to
+     *	    check.
+     *    We want a rational number 
+     *      r = b * 10**(1-k) = bw * 2**b2 * 5**b5 / (2**s2 / 5**s5),
+     *    with b2, b5, s2, s5 >= 0.  Note that the most significant decimal
+     *    digit is floor(r) and that successive digits can be obtained
+     *    by setting r <- 10*floor(r) (or b <= 10 * (b % S)).
+     *    Find appropriate b2, b5, s2, s5.
+     */
+
+    ComputeScale(be, k, &b2, &b5, &s2, &s5);
+
+    /*
+     * Correct an incorrect caller-supplied 'ndigits'.
+     * Also determine:
+     *	i = The maximum number of decimal digits that will be returned in the
+     *      formatted string.  This is k + 1 + ndigits for F format, 18 for
+     *	    shortest and Steele, and ndigits for E format.
+     *  ilim = The number of significant digits to convert if
+     *         k has been guessed correctly. This is -1 for shortest and Steele
+     *         (which stop when all significance has been lost), 'ndigits'
+     *	       for E format, and 'k + 1 + ndigits' for F format.
+     *  ilim1 = The minimum number of significant digits to convert if
+     *	        k has been guessed 1 too high. This, too, is -1 for shortest
+     *	        and Steele, and 'ndigits' for E format, but it's 'ndigits-1'
+     *	        for F format.
+     */
+
+    SetPrecisionLimits(convType, k, &ndigits, &i, &ilim, &ilim1);
+
+    /* 
+     * Try to do low-precision conversion in floating point rather
+     * than resorting to expensive multiprecision arithmetic
+     */
+    if (ilim >= 0 && ilim <= QUICK_MAX && !(flags & TCL_DD_NO_QUICK)) {
+	if ((retval = QuickConversion(d.d, k, k_check, flags,
+				      i, ilim, ilim1,
+				      decpt, endPtr)) != NULL) {
+	    return retval;
+	}
+    }
+
+    /* 
+     * For shortening conversions, determine the upper and lower bounds
+     * for the remainder at which we can stop.
+     *   m+ = (2**m2plus * 5**m5) / (2**s2 * 5**s5) is the limit on the
+     *        high side, and
+     *   m- = (2**m2minus * 5**m5) / (2**s2 * 5**s5) is the limit on the
+     *        low side.
+     *   We may need to increase s2 to put m2plus, m2minus, b2 over a
+     *   common denominator.
+     */
+
+    if (flags & TCL_DD_SHORTEN_FLAG) {
+	int m2minus = b2;
+	int m2plus;
+	int m5 = b5;
+	int len = i;
+
+	/* 
+	 * Find the quantity i so that (2**i*5**b5)/(2**s2*5**s5)
+	 * is 1/2 unit in the least significant place of the floating 
+	 * point number.
+	 */
+	if (denorm) {
+	    i = be + EXPONENT_BIAS + (FP_PRECISION-1);
+	} else {
+	    i = 1 + FP_PRECISION - bbits;
+	}
+	b2 += i;
+	s2 += i;
+
+	/* 
+	 * Reduce the fractions to lowest terms, since the above calculation
+	 * may have left excess powers of 2 in numerator and denominator
+	 */
+	CastOutPowersOf2(&b2, &m2minus, &s2);
+
+	/*
+	 * In the special case where bw==1, the nearest floating point number
+	 * to it on the low side is 1/4 ulp below it. Adjust accordingly.
+	 */
+	m2plus = m2minus;
+	if (!denorm && bw == 1) {
+	    ++b2;
+	    ++s2;
+	    ++m2plus;
+	}
+
+	if (s5+1 < N_LOG2POW5
+	    && s2+1 + log2pow5[s5+1] <= 64) {
+	    /*
+	     * If 10*2**s2*5**s5 == 2**(s2+1)+5**(s5+1) fits in a 64-bit
+	     * word, then all our intermediate calculations can be done
+	     * using exact 64-bit arithmetic with no need for expensive
+	     * multiprecision operations. (This will be true for all numbers
+	     * in the range [1.0e-3 .. 1.0e+24]).
+	     */
+
+	    return ShorteningInt64Conversion(&d, convType, bw, b2, b5,
+					     m2plus, m2minus, m5,
+					     s2, s5, k, len, ilim, ilim1, 
+					     decpt, endPtr);
+	} else if (s5 == 0) {
+	    /*
+	     * The denominator is a power of 2, so we can replace division
+	     * by digit shifts. First we round up s2 to a multiple of
+	     * DIGIT_BIT, and adjust m2 and b2 accordingly. Then we launch
+	     * into a version of the comparison that's specialized for
+	     * the 'power of mp_digit in the denominator' case.
+	     */
+	    if (s2 % DIGIT_BIT != 0) {
+		int delta = DIGIT_BIT - (s2 % DIGIT_BIT);
+		b2 += delta;
+		m2plus += delta;
+		m2minus += delta;
+		s2 += delta;
+	    }
+	    return ShorteningBignumConversionPowD(&d, convType, bw, b2, b5,
+						  m2plus, m2minus, m5,
+						  s2/DIGIT_BIT, k, len, 
+						  ilim, ilim1, decpt, endPtr);
+	} else {
+
+	    /* 
+	     * Alas, there's no helpful special case; use full-up
+	     * bignum arithmetic for the conversion
+	     */
+
+	    return ShorteningBignumConversion(&d, convType, bw,
+					      b2, m2plus, m2minus,
+					      s2, s5, k, len,
+					      ilim, ilim1, decpt, endPtr);
+
+	}
+
+    } else {
+
+	/* Non-shortening conversion */
+
+	int len = i;
+
+	/* Reduce numerator and denominator to lowest terms */
+
+	if (b2 >= s2 && s2 > 0) {
+	    b2 -= s2; s2 = 0;
+	} else if (s2 >= b2 && b2 > 0) {
+	    s2 -= b2; b2 = 0;
+	}
+
+	if (s5+1 < N_LOG2POW5
+	    && s2+1 + log2pow5[s5+1] <= 64) {
+	    /*
+	     * If 10*2**s2*5**s5 == 2**(s2+1)+5**(s5+1) fits in a 64-bit
+	     * word, then all our intermediate calculations can be done
+	     * using exact 64-bit arithmetic with no need for expensive
+	     * multiprecision operations.
+	     */
+
+	    return StrictInt64Conversion(&d, convType, bw, b2, b5,
+					 s2, s5, k, len, ilim, ilim1, 
+					 decpt, endPtr);
+
+	} else if (s5 == 0) {
+	    /*
+	     * The denominator is a power of 2, so we can replace division
+	     * by digit shifts. First we round up s2 to a multiple of
+	     * DIGIT_BIT, and adjust m2 and b2 accordingly. Then we launch
+	     * into a version of the comparison that's specialized for
+	     * the 'power of mp_digit in the denominator' case.
+	     */
+	    if (s2 % DIGIT_BIT != 0) {
+		int delta = DIGIT_BIT - (s2 % DIGIT_BIT);
+		b2 += delta;
+		s2 += delta;
+	    }
+	    return StrictBignumConversionPowD(&d, convType, bw, b2, b5,
+					      s2/DIGIT_BIT, k, len, 
+					      ilim, ilim1, decpt, endPtr);
+	} else {
+	    /*
+	     * There are no helpful special cases, but at least we know
+	     * in advance how many digits we will convert. We can run the
+	     * conversion in steps of DIGIT_GROUP digits, so as to
+	     * have many fewer mp_int divisions.
+	     */
+	    return StrictBignumConversion(&d, convType, bw, b2, s2, s5,
+					  k, len, ilim, ilim1, decpt, endPtr);
+	}
+    } 
+}
+
+
+/*
+ *----------------------------------------------------------------------
+ *
+ * TclInitDoubleConversion --
+ *
+ *	Initializes constants that are needed for conversions to and from
+ *	'double'
+ *
+ * Results:
+ *	None.
  *
  * Side effects:
  *	The log base 2 of the floating point radix, the number of bits in a
@@ -2237,6 +4380,11 @@
     for (i=0; i<8; ++i) {
 	mp_sqr(pow5+i, pow5+i+1);
     }
+    mp_init_set_int(pow5_13, 1220703125);
+    for (i = 1; i < 5; ++i) {
+	mp_init(pow5_13 + i);
+	mp_sqr(pow5_13 + i - 1, pow5_13 + i);
+    }
 
     /*
      * Determine the number of decimal digits to the left and right of the
@@ -2293,7 +4441,7 @@
 {
     int i;
 
-    Tcl_Free((char *) pow10_wide);
+    ckfree((char *) pow10_wide);
     for (i=0; i<9; ++i) {
 	mp_clear(pow5 + i);
     }
@@ -2433,6 +4581,20 @@
 	return -r;
     }
 }
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * TclCeil --
+ *
+ *	Computes the smallest floating point number that is at least the
+ *	mp_int argument.
+ *
+ * Results:
+ *	Returns the floating point number.
+ *
+ *-----------------------------------------------------------------------------
+ */
 
 double
 TclCeil(
@@ -2476,6 +4638,20 @@
     mp_clear(&b);
     return r;
 }
+
+/*
+ *-----------------------------------------------------------------------------
+ *
+ * TclFloor --
+ *
+ *	Computes the largest floating point number less than or equal to
+ *	the mp_int argument.
+ *
+ * Results:
+ *	Returns the floating point value.
+ *
+ *-----------------------------------------------------------------------------
+ */
 
 double
 TclFloor(
Index: generic/tclStubInit.c
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclStubInit.c,v
retrieving revision 1.197
diff -u -r1.197 tclStubInit.c
--- generic/tclStubInit.c	16 Sep 2010 14:49:37 -0000	1.197
+++ generic/tclStubInit.c	27 Nov 2010 17:40:53 -0000
@@ -305,6 +305,7 @@
     TclInitRewriteEnsemble, /* 246 */
     TclResetRewriteEnsemble, /* 247 */
     TclCopyChannel, /* 248 */
+    TclDoubleDigits, /* 249 */
 };
 
 static const TclIntPlatStubs tclIntPlatStubs = {
@@ -460,6 +461,8 @@
     TclBN_s_mp_mul_digs, /* 58 */
     TclBN_s_mp_sqr, /* 59 */
     TclBN_s_mp_sub, /* 60 */
+    TclBN_mp_init_set_int, /* 61 */
+    TclBN_mp_set_int, /* 62 */
 };
 
 static const TclStubHooks tclStubHooks = {
Index: generic/tclTest.c
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclTest.c,v
retrieving revision 1.154
diff -u -r1.154 tclTest.c
--- generic/tclTest.c	27 Sep 2010 19:42:38 -0000	1.154
+++ generic/tclTest.c	27 Nov 2010 17:40:54 -0000
@@ -17,6 +17,8 @@
  * RCS: @(#) $Id: tclTest.c,v 1.154 2010/09/27 19:42:38 msofer Exp $
  */
 
+#include <math.h>
+
 #undef STATIC_BUILD
 #ifndef USE_TCL_STUBS
 #   define USE_TCL_STUBS
@@ -233,6 +235,9 @@
 			    Tcl_Interp *interp, int argc, const char **argv);
 static int		TestdelassocdataCmd(ClientData dummy,
 			    Tcl_Interp *interp, int argc, const char **argv);
+static int		TestdoubledigitsObjCmd(ClientData dummy,
+					       Tcl_Interp* interp,
+					       int objc, Tcl_Obj* const objv[]);
 static int		TestdstringCmd(ClientData dummy,
 			    Tcl_Interp *interp, int argc, const char **argv);
 static int		TestencodingObjCmd(ClientData dummy,
@@ -569,6 +574,8 @@
     Tcl_CreateCommand(interp, "testdel", TestdelCmd, NULL, NULL);
     Tcl_CreateCommand(interp, "testdelassocdata", TestdelassocdataCmd,
 	    NULL, NULL);
+    Tcl_CreateObjCommand(interp, "testdoubledigits", TestdoubledigitsObjCmd,
+			 NULL, NULL);
     Tcl_DStringInit(&dstring);
     Tcl_CreateCommand(interp, "testdstring", TestdstringCmd, NULL,
 	    NULL);
@@ -1594,6 +1601,102 @@
 }
 
 /*
+ *-----------------------------------------------------------------------------
+ *
+ * TestdoubledigitsCmd --
+ *
+ *	This procedure implements the 'testdoubledigits' command. It is
+ *	used to test the low-level floating-point formatting primitives
+ *	in Tcl.
+ *
+ * Usage:
+ *	testdoubledigits fpval ndigits type ?shorten"
+ *
+ * Parameters:
+ *	fpval - Floating-point value to format.
+ *	ndigits - Digit count to request from Tcl_DoubleDigits
+ *	type - One of 'shortest', 'Steele', 'e', 'f'
+ *	shorten - Indicates that the 'shorten' flag should be passed in.
+ *
+ *-----------------------------------------------------------------------------
+ */
+
+static int
+TestdoubledigitsObjCmd(ClientData unused,
+				/* NULL */
+		       Tcl_Interp* interp,
+				/* Tcl interpreter */
+		       int objc,
+				/* Parameter count */
+		       Tcl_Obj* const objv[])
+				/* Parameter vector */
+{
+    static const char* options[] = {
+	"shortest",
+	"Steele",
+	"e",
+	"f",
+	NULL
+    };
+    static const int types[] = {
+	TCL_DD_SHORTEST,
+	TCL_DD_STEELE,
+	TCL_DD_E_FORMAT,
+	TCL_DD_F_FORMAT
+    };
+
+    const Tcl_ObjType* doubleType;
+    double d;
+    int status;
+    int ndigits;
+    int type;
+    int decpt;
+    int signum;
+    char* str;
+    char* endPtr;
+    Tcl_Obj* strObj;
+    Tcl_Obj* retval;
+
+    if (objc < 4 || objc > 5) {
+	Tcl_WrongNumArgs(interp, 1, objv, "fpval ndigits type ?shorten?");
+	return TCL_ERROR;
+    }
+    status = Tcl_GetDoubleFromObj(interp, objv[1], &d);
+    if (status != TCL_OK) {
+	doubleType = Tcl_GetObjType("double");
+	if (objv[1]->typePtr == doubleType
+	    || TclIsNaN(objv[1]->internalRep.doubleValue)) {
+	    status = TCL_OK;
+	    memcpy(&d, &(objv[1]->internalRep.doubleValue), sizeof(double));
+	}
+    }
+    if (status != TCL_OK
+	|| Tcl_GetIntFromObj(interp, objv[2], &ndigits) != TCL_OK
+	|| Tcl_GetIndexFromObj(interp, objv[3], options, "conversion type",
+			       TCL_EXACT, &type) != TCL_OK) {
+	fprintf(stderr, "bad value? %g\n", d);
+	return TCL_ERROR;
+    }
+    type = types[type];
+    if (objc > 4) {
+	if (strcmp(Tcl_GetString(objv[4]), "shorten")) {
+	    Tcl_SetObjResult(interp, Tcl_NewStringObj("bad flag", -1));
+	    return TCL_ERROR;
+	}
+	type |= TCL_DD_SHORTEN_FLAG;
+    }
+    str = TclDoubleDigits(d, ndigits, type, &decpt, &signum, &endPtr);
+    strObj = Tcl_NewStringObj(str, endPtr-str);
+    ckfree(str);
+    retval = Tcl_NewListObj(1, &strObj);
+    Tcl_ListObjAppendElement(NULL, retval, Tcl_NewIntObj(decpt));
+    strObj = Tcl_NewStringObj(signum ? "-" : "+", 1);
+    Tcl_ListObjAppendElement(NULL, retval, strObj);
+    Tcl_SetObjResult(interp, retval);
+    return TCL_OK;
+}
+
+/*
  *----------------------------------------------------------------------
  *
  * TestdstringCmd --
Index: generic/tclTomMath.decls
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclTomMath.decls,v
retrieving revision 1.8
diff -u -r1.8 tclTomMath.decls
--- generic/tclTomMath.decls	15 Sep 2010 07:33:54 -0000	1.8
+++ generic/tclTomMath.decls	27 Nov 2010 17:40:54 -0000
@@ -214,3 +214,9 @@
 declare 60 {
     int TclBN_s_mp_sub(mp_int *a, mp_int *b, mp_int *c)
 }
+declare 61 {
+    int TclBN_mp_init_set_int(mp_int* a, unsigned long i)
+}
+declare 62 {
+    int TclBN_mp_set_int(mp_int* a, unsigned long i)
+}
\ No newline at end of file
Index: generic/tclTomMathDecls.h
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclTomMathDecls.h,v
retrieving revision 1.13
diff -u -r1.13 tclTomMathDecls.h
--- generic/tclTomMathDecls.h	19 Aug 2010 04:26:04 -0000	1.13
+++ generic/tclTomMathDecls.h	27 Nov 2010 17:40:54 -0000
@@ -79,6 +79,7 @@
 #define mp_init_copy TclBN_mp_init_copy
 #define mp_init_multi TclBN_mp_init_multi
 #define mp_init_set TclBN_mp_init_set
+#define mp_init_set_int TclBN_mp_init_set_int
 #define mp_init_size TclBN_mp_init_size
 #define mp_karatsuba_mul TclBN_mp_karatsuba_mul
 #define mp_karatsuba_sqr TclBN_mp_karatsuba_sqr
@@ -96,6 +97,7 @@
 #define mp_rshd TclBN_mp_rshd
 #define mp_s_rmap TclBNMpSRmap
 #define mp_set TclBN_mp_set
+#define mp_set_int TclBN_mp_set_int
 #define mp_shrink TclBN_mp_shrink
 #define mp_sqr TclBN_mp_sqr
 #define mp_sqrt TclBN_mp_sqrt
@@ -268,6 +270,10 @@
 EXTERN int		TclBN_s_mp_sqr(mp_int *a, mp_int *b);
 /* 60 */
 EXTERN int		TclBN_s_mp_sub(mp_int *a, mp_int *b, mp_int *c);
+/* 61 */
+EXTERN int		TclBN_mp_init_set_int(mp_int*a, unsigned long i);
+/* 62 */
+EXTERN int		TclBN_mp_set_int(mp_int*a, unsigned long i);
 
 typedef struct TclTomMathStubs {
     int magic;
@@ -334,6 +340,8 @@
     int (*tclBN_s_mp_mul_digs) (mp_int *a, mp_int *b, mp_int *c, int digs); /* 58 */
     int (*tclBN_s_mp_sqr) (mp_int *a, mp_int *b); /* 59 */
     int (*tclBN_s_mp_sub) (mp_int *a, mp_int *b, mp_int *c); /* 60 */
+    int (*tclBN_mp_init_set_int) (mp_int*a, unsigned long i); /* 61 */
+    int (*tclBN_mp_set_int) (mp_int*a, unsigned long i); /* 62 */
 } TclTomMathStubs;
 
 #ifdef __cplusplus
@@ -472,6 +480,10 @@
 	(tclTomMathStubsPtr->tclBN_s_mp_sqr) /* 59 */
 #define TclBN_s_mp_sub \
 	(tclTomMathStubsPtr->tclBN_s_mp_sub) /* 60 */
+#define TclBN_mp_init_set_int \
+	(tclTomMathStubsPtr->tclBN_mp_init_set_int) /* 61 */
+#define TclBN_mp_set_int \
+	(tclTomMathStubsPtr->tclBN_mp_set_int) /* 62 */
 
 #endif /* defined(USE_TCL_STUBS) */
 
Index: generic/tclUtil.c
===================================================================
RCS file: /cvsroot/tcl/tcl/generic/tclUtil.c,v
retrieving revision 1.118
diff -u -r1.118 tclUtil.c
--- generic/tclUtil.c	1 Oct 2010 12:52:50 -0000	1.118
+++ generic/tclUtil.c	27 Nov 2010 17:40:54 -0000
@@ -2234,129 +2234,100 @@
     char *p, c;
     int exponent;
     int signum;
-    char buffer[TCL_DOUBLE_SPACE];
+    char* digits;
+    char* end;
     Tcl_UniChar ch;
 
     int *precisionPtr = Tcl_GetThreadData(&precisionKey, (int)sizeof(int));
 
     /*
-     * If *precisionPtr == 0, then use TclDoubleDigits to develop a decimal
-     * significand and exponent, then format it in E or F format as
-     * appropriate. If *precisionPtr != 0, use the native sprintf and then add
-     * a trailing ".0" if there is no decimal point in the rep.
+     * Handle NaN.
      */
+    
+    if (TclIsNaN(value)) {
+	TclFormatNaN(value, dst);
+	return;
+    }
 
-    if (*precisionPtr == 0) {
+    /*
+     * Handle infinities.
+     */
+    
+    if (TclIsInfinite(value)) {
 	/*
-	 * Handle NaN.
+	 * Remember to copy the terminating NUL too.
 	 */
-
-	if (TclIsNaN(value)) {
-	    TclFormatNaN(value, dst);
-	    return;
+	
+	if (value < 0) {
+	    memcpy(dst, "-Inf", 5);
+	} else {
+	    memcpy(dst, "Inf", 4);
 	}
+	return;
+    }
 
+    /*
+     * Ordinary (normal and denormal) values.
+     */
+    
+    if (*precisionPtr == 0) {
+	digits = TclDoubleDigits(value, -1, TCL_DD_SHORTEST,
+				 &exponent, &signum, &end);
+    } else {
+	digits = TclDoubleDigits(value, *precisionPtr, TCL_DD_E_FORMAT, 
+				 &exponent, &signum, &end);
+    }
+    if (signum) {
+	*dst++ = '-';
+    }
+    p = digits;
+    if (exponent < -4 || exponent > 16) {
 	/*
-	 * Handle infinities.
+	 * E format for numbers < 1e-3 or >= 1e17.
 	 */
-
-	if (TclIsInfinite(value)) {
-	    /*
-	     * Remember to copy the terminating NUL too.
-	     */
-
-	    if (value < 0) {
-		memcpy(dst, "-Inf", 5);
-	    } else {
-		memcpy(dst, "Inf", 4);
+	
+	*dst++ = *p++;
+	c = *p;
+	if (c != '\0') {
+	    *dst++ = '.';
+	    while (c != '\0') {
+		*dst++ = c;
+		c = *++p;
 	    }
-	    return;
 	}
-
+	sprintf(dst, "e%+d", exponent);
+    } else {
 	/*
-	 * Ordinary (normal and denormal) values.
+	 * F format for others.
 	 */
-
-	exponent = TclDoubleDigits(buffer, value, &signum);
-	if (signum) {
-	    *dst++ = '-';
+	
+	if (exponent < 0) {
+	    *dst++ = '0';
 	}
-	p = buffer;
-	if (exponent < -3 || exponent > 17) {
-	    /*
-	     * E format for numbers < 1e-3 or >= 1e17.
-	     */
-
-	    *dst++ = *p++;
-	    c = *p;
+	c = *p;
+	while (exponent-- >= 0) {
 	    if (c != '\0') {
-		*dst++ = '.';
-		while (c != '\0') {
-		    *dst++ = c;
-		    c = *++p;
-		}
-	    }
-	    sprintf(dst, "e%+d", exponent-1);
-	} else {
-	    /*
-	     * F format for others.
-	     */
-
-	    if (exponent <= 0) {
-		*dst++ = '0';
-	    }
-	    c = *p;
-	    while (exponent-- > 0) {
-		if (c != '\0') {
-		    *dst++ = c;
-		    c = *++p;
-		} else {
-		    *dst++ = '0';
-		}
-	    }
-	    *dst++ = '.';
-	    if (c == '\0') {
-		*dst++ = '0';
+		*dst++ = c;
+		c = *++p;
 	    } else {
-		while (++exponent < 0) {
-		    *dst++ = '0';
-		}
-		while (c != '\0') {
-		    *dst++ = c;
-		    c = *++p;
-		}
+		*dst++ = '0';
 	    }
-	    *dst++ = '\0';
 	}
-    } else {
-	/*
-	 * tcl_precision is supplied, pass it to the native sprintf.
-	 */
-
-	sprintf(dst, "%.*g", *precisionPtr, value);
-
-	/*
-	 * If the ASCII result looks like an integer, add ".0" so that it
-	 * doesn't look like an integer anymore. This prevents floating-point
-	 * values from being converted to integers unintentionally. Check for
-	 * ASCII specifically to speed up the function.
-	 */
-
-	for (p = dst; *p != 0;) {
-	    if (UCHAR(*p) < 0x80) {
-		c = *p++;
-	    } else {
-		p += Tcl_UtfToUniChar(p, &ch);
-		c = UCHAR(ch);
+	*dst++ = '.';
+	if (c == '\0') {
+	    *dst++ = '0';
+	} else {
+	    while (++exponent < -1) {
+		*dst++ = '0';
 	    }
-	    if ((c == '.') || isalpha(UCHAR(c))) {	/* INTL: ISO only. */
-		return;
+	    while (c != '\0') {
+		*dst++ = c;
+		c = *++p;
 	    }
 	}
-	p[0] = '.';
-	p[1] = '0';
-	p[2] = 0;
+	*dst++ = '\0';
     }
+    ckfree(digits);
 }
 
 /*
Index: tests/util.test
===================================================================
RCS file: /cvsroot/tcl/tcl/tests/util.test,v
retrieving revision 1.20
diff -u -r1.20 util.test
--- tests/util.test	14 Oct 2008 16:35:44 -0000	1.20
+++ tests/util.test	27 Nov 2010 17:40:54 -0000
@@ -43,6 +43,10 @@
 		ieeeValues(+Infinity)
 	    binary scan \x00\x00\x00\x00\x00\x00\xf8\x7f d \
 		ieeeValues(NaN)
+	    binary scan \x00\x00\x00\x00\x00\x00\xf8\xff d \
+		ieeeValues(-NaN)
+	    binary scan \xef\xcd\xab\x89\x67\x45\xfb\xff d \
+		ieeeValues(-NaN(3456789abcdef))
 	    set ieeeValues(littleEndian) 1
 	    return 1
 	}
@@ -65,6 +69,10 @@
 		ieeeValues(+Infinity)
 	    binary scan \x7f\xf8\x00\x00\x00\x00\x00\x00 d \
 		ieeeValues(NaN)
+	    binary scan \xff\xf8\x00\x00\x00\x00\x00\x00 d \
+		ieeeValues(-NaN)
+	    binary scan \xff\xfb\x45\x67\x89\xab\xcd\xef d \
+		ieeeValues(-NaN(3456789abcdef))
 	    set ieeeValues(littleEndian) 0
 	    return 1
 	}
@@ -85,6 +93,30 @@
     return $result
 }
 
+proc verdonk_test {sig binexp shouldbe exp} {
+    regexp {([-+]?)([0-9a-f]+)} $sig -> signum sig
+    scan $sig %llx sig
+    if {$signum eq {-}} {
+	set signum [expr 1<<63]
+    } else {
+	set signum 0
+    }
+    regexp {E([-+]?[0-9]+)} $binexp -> binexp
+    set word [expr {$signum | (($binexp + 0x3ff)<<52)|($sig & ~(1<<52))}]
+    binary scan [binary format w $word] q double
+    regexp {([-+])(\d+)_(\d+)\&} $shouldbe -> signum digits1 digits2
+    regexp {E([-+]\d+)} $exp -> decexp
+    incr decexp [expr {[string length $digits1] - 1}]
+    lassign [testdoubledigits $double [string length $digits1] e] \
+	outdigits decpt outsign
+    if {[string index $digits2 0] >= 5} {
+	incr digits1
+    }
+    if {$outsign != $signum || $outdigits != $digits1 || $decpt != $decexp} {
+	return -code error "result is ${outsign}0.${outdigits}E$decpt\
+                            should be ${signum}0.${digits1}E$decexp"
+    }
+}
 
 test util-1.1 {TclFindElement procedure - binary element in middle of list} {
     lindex {0 foo\x00help 1} 1
@@ -1106,6 +1138,766 @@
     expr 1.1e17
 } {1.1e+17}
 
+test util-12.1 {Tcl_DoubleDigits - Inf} ieeeFloatingPoint {
+     testdoubledigits Inf -1 shortest
+} {Infinity 9999 +}
+test util-12.2 {Tcl_DoubleDigits - -Inf} ieeeFloatingPoint {
+     testdoubledigits -Inf -1 shortest
+} {Infinity 9999 -}
+test util-12.3 {Tcl_DoubleDigits - NaN} ieeeFloatingPoint {
+     testdoubledigits $ieeeValues(NaN) -1 shortest
+} {NaN 9999 +}
+test util-12.4 {Tcl_DoubleDigits - NaN} ieeeFloatingPoint {
+     testdoubledigits -NaN -1 shortest
+} {NaN 9999 -}
+test util-12.5 {Tcl_DoubleDigits - 0} {
+     testdoubledigits 0.0 -1 shortest
+} {0 0 +}
+test util-12.6 {Tcl_DoubleDigits - -0} {
+     testdoubledigits -0.0 -1 shortest
+} {0 0 -}
+
+# Verdonk test vectors
+
+test util-13.1 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1754e31cd072da E+1008 +4_000000000000000000& E+303
+    }
+    -result {}
+}
+test util-13.2 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1afcef51f0fb5f E+265 -1_000000000000000000& E+80
+    }
+    -result {}
+}
+test util-13.3 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1754e31cd072da E+1006 +1_000000000000000000& E+303
+    }
+    -result {}
+}
+test util-13.4 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1754e31cd072da E+1007 -2_000000000000000000& E+303
+    }
+    -result {}
+}
+test util-13.5 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1e07b27dd78b14 E-848 +1_00000000000000000& E-255
+    }
+    -result {}
+}
+test util-13.6 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1e29e9c56687fe E-709 -7_00000000000000000& E-214
+    }
+    -result {}
+}
+test util-13.7 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1be03d0bf225c7 E-137 +1_00000000000000000& E-41
+    }
+    -result {}
+}
+test util-13.8 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1a2fe76a3f9475 E-499 -1_00000000000000000& E-150
+    }
+    -result {}
+}
+test util-13.9 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 19a2028368022e E+1019 +8_999999999999999999& E+306
+    }
+    -result {}
+}
+test util-13.10 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1317e5ef3ab327 E+509 -1_999999999999999999& E+153
+    }
+    -result {}
+}
+test util-13.11 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1317e5ef3ab327 E+510 +3_99999999999999999& E+153
+    }
+    -result {}
+}
+test util-13.12 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1317e5ef3ab327 E+511 -7_99999999999999999& E+153
+    }
+    -result {}
+}
+test util-13.13 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1eb8e84fa0b278 E-1008 +6_999999999999999999& E-304
+    }
+    -result {}
+}
+test util-13.14 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -13339131c46f8b E-1004 -6_999999999999999999& E-303
+    }
+    -result {}
+}
+test util-13.15 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1c0f92a6276c9d E-162 +2_999999999999999999& E-49
+    }
+    -result {}
+}
+test util-13.16 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -15ce1f143d7ad2 E-443 -5_99999999999999999& E-134
+    }
+    -result {}
+}
+test util-13.17 {just over exact - 2 digits} {*}{
+    -body {
+        verdonk_test 1c0794d9d40e96 E-301 +43_000000000000000000& E-92
+    }
+    -result {}
+}
+test util-13.18 {just over exact - 2 digits} {*}{
+    -body {
+        verdonk_test -1c0794d9d40e96 E-300 -86_000000000000000000& E-92
+    }
+    -result {}
+}
+test util-13.19 {just over exact - 2 digits} {*}{
+    -body {
+        verdonk_test 1cd5bee57763e6 E-241 +51_000000000000000000& E-74
+    }
+    -result {}
+}
+test util-13.20 {just under exact - 2 digits} {*}{
+    -body {
+        verdonk_test 1d1c26db7d0dae E+651 +16_999999999999999999& E+195
+    }
+    -result {}
+}
+test util-13.21 {just under exact - 2 digits} {*}{
+    -body {
+        verdonk_test -13f7ced916872b E-5 -38_999999999999999999& E-3
+    }
+    -result {}
+}
+test util-13.22 {just over exact - 3 digits} {*}{
+    -body {
+        verdonk_test 17d93193f78fc6 E+588 +151_0000000000000000000& E+175
+    }
+    -result {}
+}
+test util-13.23 {just over exact - 3 digits} {*}{
+    -body {
+        verdonk_test -1a82a1631eeb30 E-625 -119_000000000000000000& E-190
+    }
+    -result {}
+}
+test util-13.24 {just under exact - 3 digits} {*}{
+    -body {
+        verdonk_test -16c309024bab4b E+290 -282_999999999999999999& E+85
+    }
+    -result {}
+}
+test util-13.25 {just over exact - 8 digits} {*}{
+    -body {
+        verdonk_test 1dbbac6f83a821 E-800 +27869147_0000000000000000000& E-248
+    }
+    -result {}
+}
+test util-13.26 {just under exact - 9 digits} {*}{
+    -body {
+        verdonk_test -1c569e968e0944 E+430 -491080653_9999999999999999999& E+121
+    }
+    -result {}
+}
+test util-13.27 {just under exact - 9 digits} {*}{
+    -body {
+        verdonk_test 1c569e968e0944 E+429 +245540326_9999999999999999999& E+121
+    }
+    -result {}
+}
+test util-13.28 {just over exact - 10 digits} {*}{
+    -body {
+        verdonk_test -1fc575867314ee E-330 -9078555839_0000000000000000000& E-109
+    }
+    -result {}
+}
+test util-13.29 {just under exact - 10 digits} {*}{
+    -body {
+        verdonk_test -1c569e968e0944 E+428 -1227701634_9999999999999999999& E+120
+    }
+    -result {}
+}
+test util-13.30 {just over exact - 11 digits} {*}{
+    -body {
+        verdonk_test 1fc575867314ee E-329 +18157111678_0000000000000000000& E-109
+    }
+    -result {}
+}
+test util-13.31 {just over exact - 14 digits} {*}{
+    -body {
+        verdonk_test -18bf7e7fa6f02a E-196 -15400733123779_0000000000000000000& E-72
+    }
+    -result {}
+}
+test util-13.32 {just over exact - 17 digits} {*}{
+    -body {
+        verdonk_test -13de005bd620df E+217 -26153245263757307_0000000000000000000& E+49
+    }
+    -result {}
+}
+test util-13.33 {just over exact - 18 digits} {*}{
+    -body {
+        verdonk_test 1f92bacb3cb40c E+718 +272104041512242479_0000000000000000000& E+199
+    }
+    -result {}
+}
+test util-13.34 {just over exact - 18 digits} {*}{
+    -body {
+        verdonk_test -1f92bacb3cb40c E+719 -544208083024484958_0000000000000000000& E+199
+    }
+    -result {}
+}
+test util-13.35 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 142dbf25096cf5 E+148 +4_500000000000000000& E+44
+    }
+    -result {}
+}
+test util-13.36 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1afcef51f0fb5f E+263 -2_500000000000000000& E+79
+    }
+    -result {}
+}
+test util-13.37 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 102498ea6df0c4 E+145 +4_500000000000000000& E+43
+    }
+    -result {}
+}
+test util-13.38 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1754e31cd072da E+1004 -2_500000000000000000& E+302
+    }
+    -result {}
+}
+test util-13.39 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 12deac01e2b4f7 E-557 +2_50000000000000000& E-168
+    }
+    -result {}
+}
+test util-13.40 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1b1df536c13eee E-307 -6_50000000000000000& E-93
+    }
+    -result {}
+}
+test util-13.41 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 10711fed5b19a4 E-154 +4_50000000000000000& E-47
+    }
+    -result {}
+}
+test util-13.42 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -148d67e8b1e00d E-151 -4_50000000000000000& E-46
+    }
+    -result {}
+}
+test util-13.43 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 1c8c574c0c6be7 E+187 +3_49999999999999999& E+56
+    }
+    -result {}
+}
+test util-13.44 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1756183c147514 E+206 -1_49999999999999999& E+62
+    }
+    -result {}
+}
+test util-13.45 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 12ab469676c410 E+203 +1_49999999999999999& E+61
+    }
+    -result {}
+}
+test util-13.46 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1539684e774b48 E+246 -1_49999999999999999& E+74
+    }
+    -result {}
+}
+test util-13.47 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 12e5f5dfa4fe9d E-286 +9_499999999999999999& E-87
+    }
+    -result {}
+}
+test util-13.48 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1bdc2417bf7787 E-838 -9_499999999999999999& E-253
+    }
+    -result {}
+}
+test util-13.49 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 1eb8e84fa0b278 E-1009 +3_499999999999999999& E-304
+    }
+    -result {}
+}
+test util-13.50 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1e3cbc9907fdc8 E-290 -9_499999999999999999& E-88
+    }
+    -result {}
+}
+test util-13.51 {just over half ulp - 2 digits} {*}{
+    -body {
+        verdonk_test 10ad836f269a17 E-324 +30_500000000000000000& E-99
+    }
+    -result {}
+}
+test util-13.52 {just over half ulp - 2 digits} {*}{
+    -body {
+        verdonk_test -1b39ae1909c31b E-687 -26_500000000000000000& E-208
+    }
+    -result {}
+}
+test util-13.53 {just over half ulp - 3 digits} {*}{
+    -body {
+        verdonk_test 1b2ab18615fcc6 E-576 +686_500000000000000000& E-176
+    }
+    -result {}
+}
+test util-13.54 {just over half ulp - 3 digits} {*}{
+    -body {
+        verdonk_test -13e1f90a573064 E-624 -178_500000000000000000& E-190
+    }
+    -result {}
+}
+test util-13.55 {just under half ulp - 3 digits} {*}{
+    -body {
+        verdonk_test 16c309024bab4b E+289 +141_499999999999999999& E+85
+    }
+    -result {}
+}
+test util-13.56 {just under half ulp - 4 digits} {*}{
+    -body {
+        verdonk_test -159bd3ad46e346 E+193 -1695_499999999999999999& E+55
+    }
+    -result {}
+}
+test util-13.57 {just under half ulp - 4 digits} {*}{
+    -body {
+        verdonk_test 1df4170f0fdecc E+124 +3981_499999999999999999& E+34
+    }
+    -result {}
+}
+test util-13.58 {just over half ulp - 6 digits} {*}{
+    -body {
+        verdonk_test 17e1e0f1c7a4ac E+415 +126300_5000000000000000000& E+120
+    }
+    -result {}
+}
+test util-13.59 {just over half ulp - 6 digits} {*}{
+    -body {
+        verdonk_test -1dda592e398dd7 E+418 -126300_5000000000000000000& E+121
+    }
+    -result {}
+}
+test util-13.60 {just under half ulp - 7 digits} {*}{
+    -body {
+        verdonk_test -1e597c0b94b7ae E+453 -4411845_499999999999999999& E+130
+    }
+    -result {}
+}
+test util-13.61 {just under half ulp - 9 digits} {*}{
+    -body {
+        verdonk_test 1c569e968e0944 E+427 +613850817_4999999999999999999& E+120
+    }
+    -result {}
+}
+test util-13.62 {just under half ulp - 9 digits} {*}{
+    -body {
+        verdonk_test -1c569e968e0944 E+428 -122770163_49999999999999999999& E+121
+    }
+    -result {}
+}
+test util-13.63 {just over half ulp - 18 digits} {*}{
+    -body {
+        verdonk_test 17ae0c186d8709 E+719 +408156062268363718_5000000000000000000& E+199
+    }
+    -result {}
+}
+test util-13.64 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 152d02c7e14af7 E+76 +1_0000000000000000& E+23
+    }
+    -result {}
+}
+test util-13.65 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -19d971e4fe8402 E+89 -1_0000000000000000& E+27
+    }
+    -result {}
+}
+test util-13.66 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 19d971e4fe8402 E+90 +2_0000000000000000& E+27
+    }
+    -result {}
+}
+test util-13.67 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -19d971e4fe8402 E+91 -4_0000000000000000& E+27
+    }
+    -result {}
+}
+test util-13.68 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 15798ee2308c3a E-27 +1_0000000000000000& E-8
+    }
+    -result {}
+}
+test util-13.69 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -15798ee2308c3a E-26 -2_0000000000000000& E-8
+    }
+    -result {}
+}
+test util-13.70 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test 15798ee2308c3a E-25 +4_0000000000000000& E-8
+    }
+    -result {}
+}
+test util-13.71 {just over exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1ef2d0f5da7dd9 E-84 -1_0000000000000000& E-25
+    }
+    -result {}
+}
+test util-13.72 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 1a784379d99db4 E+78 +4_9999999999999999& E+23
+    }
+    -result {}
+}
+test util-13.73 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1a784379d99db4 E+80 -1_9999999999999999& E+24
+    }
+    -result {}
+}
+test util-13.74 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 13da329b633647 E+81 +2_9999999999999999& E+24
+    }
+    -result {}
+}
+test util-13.75 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -1cf389cd46047d E+85 -6_9999999999999999& E+25
+    }
+    -result {}
+}
+test util-13.76 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 19999999999999 E-3 +1_99999999999999999& E-1
+    }
+    -result {}
+}
+test util-13.77 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -13333333333333 E-2 -2_99999999999999999& E-1
+    }
+    -result {}
+}
+test util-13.78 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test 16849b86a12b9b E-48 +4_99999999999999999& E-15
+    }
+    -result {}
+}
+test util-13.79 {just under exact - 1 digits} {*}{
+    -body {
+        verdonk_test -16849b86a12b9b E-46 -1_99999999999999999& E-14
+    }
+    -result {}
+}
+test util-13.80 {just over exact - 2 digits} {*}{
+    -body {
+        verdonk_test 17ccfc73126788 E-71 +63_00000000000000000& E-23
+    }
+    -result {}
+}
+test util-13.81 {just over exact - 2 digits} {*}{
+    -body {
+        verdonk_test -1dc03b8fd7016a E-68 -63_00000000000000000& E-22
+    }
+    -result {}
+}
+test util-13.82 {just under exact - 2 digits} {*}{
+    -body {
+        verdonk_test 13f7ced916872b E-5 +38_999999999999999999& E-3
+    }
+    -result {}
+}
+test util-13.83 {just over exact - 3 digits} {*}{
+    -body {
+        verdonk_test 1b297cad9f70b6 E+97 +269_000000000000000000& E+27
+    }
+    -result {}
+}
+test util-13.84 {just over exact - 3 digits} {*}{
+    -body {
+        verdonk_test -1b297cad9f70b6 E+98 -538_00000000000000000& E+27
+    }
+    -result {}
+}
+test util-13.85 {just over exact - 3 digits} {*}{
+    -body {
+        verdonk_test 1cdc06b20ef183 E-82 +373_00000000000000000& E-27
+    }
+    -result {}
+}
+test util-13.86 {just over exact - 4 digits} {*}{
+    -body {
+        verdonk_test 1b297cad9f70b6 E+96 +1345_00000000000000000& E+26
+    }
+    -result {}
+}
+# this one is not 4 digits, it is 3, and it is covered above.
+test util-13.87 {just over exact - 4 digits} {*}{
+    -constraints knownBadTest
+    -body {
+        verdonk_test -1b297cad9f70b6 E+97 -2690_00000000000000000& E+26
+    }
+    -result {}
+}
+test util-13.88 {just over exact - 5 digits} {*}{
+    -body {
+        verdonk_test -150a246ecd44f3 E-63 -14257_00000000000000000& E-23
+    }
+    -result {}
+}
+test util-13.89 {just under exact - 6 digits} {*}{
+    -body {
+        verdonk_test -119b96f36ec68b E-19 -209900_999999999999999999& E-11
+    }
+    -result {}
+}
+test util-13.90 {just over exact - 11 digits} {*}{
+    -body {
+        verdonk_test 1c06d366394441 E-35 +50980203373_000000000000000000& E-21
+    }
+    -result {}
+}
+test util-13.91 {just under exact - 12 digits} {*}{
+    -body {
+        verdonk_test -1f58ac4db68c90 E+122 -104166211810_99999999999999999& E+26
+    }
+    -result {}
+}
+test util-13.92 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 19d971e4fe8402 E+87 +2_5000000000000000& E+26
+    }
+    -result {}
+}
+test util-13.93 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1dc74be914d16b E+81 -4_500000000000000& E+24
+    }
+    -result {}
+}
+test util-13.94 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 14adf4b7320335 E+84 +2_500000000000000& E+25
+    }
+    -result {}
+}
+test util-13.95 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1ae22487c1042b E+85 -6_5000000000000000& E+25
+    }
+    -result {}
+}
+test util-13.96 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 187fe49aab41e0 E-54 +8_5000000000000000& E-17
+    }
+    -result {}
+}
+test util-13.97 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1f5c05e4b23fd7 E-61 -8_5000000000000000& E-19
+    }
+    -result {}
+}
+test util-13.98 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 1faa7ab552a552 E-42 +4_5000000000000000& E-13
+    }
+    -result {}
+}
+test util-13.99 {just over half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1b7cdfd9d7bdbb E-36 -2_5000000000000000& E-11
+    }
+    -result {}
+}
+test util-13.100 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 13da329b633647 E+80 +1_4999999999999999& E+24
+    }
+    -result {}
+}
+test util-13.101 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1cf389cd46047d E+84 -3_49999999999999999& E+25
+    }
+    -result {}
+}
+test util-13.102 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 1f04ef12cb04cf E+85 +7_4999999999999999& E+25
+    }
+    -result {}
+}
+test util-13.103 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -1f04ef12cb04cf E+86 -1_4999999999999999& E+26
+    }
+    -result {}
+}
+test util-13.104 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 13333333333333 E-3 +1_49999999999999999& E-1
+    }
+    -result {}
+}
+test util-13.105 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -107e1fe91b0b70 E-36 -1_49999999999999999& E-11
+    }
+    -result {}
+}
+test util-13.106 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test 149da7e361ce4c E-33 +1_49999999999999999& E-10
+    }
+    -result {}
+}
+test util-13.107 {just under half ulp - 1 digits} {*}{
+    -body {
+        verdonk_test -19c511dc3a41df E-30 -1_49999999999999999& E-9
+    }
+    -result {}
+}
+test util-13.108 {just over half ulp - 2 digits} {*}{
+    -body {
+        verdonk_test -1aa83d74267822 E+93 -16_5000000000000000& E+27
+    }
+    -result {}
+}
+test util-13.109 {just over half ulp - 2 digits} {*}{
+    -body {
+        verdonk_test 18f1d5969453de E+89 +96_5000000000000000& E+25
+    }
+    -result {}
+}
+test util-13.110 {just over half ulp - 2 digits} {*}{
+    -body {
+        verdonk_test 11d9bd564dcda6 E-70 +94_50000000000000000& E-23
+    }
+    -result {}
+}
+test util-13.111 {just over half ulp - 2 digits} {*}{
+    -body {
+        verdonk_test -1a58973ecbede6 E-48 -58_50000000000000000& E-16
+    }
+    -result {}
+}
+test util-13.112 {just over half ulp - 3 digits} {*}{
+    -body {
+        verdonk_test 1b297cad9f70b6 E+95 +672_50000000000000000& E+26
+    }
+    -result {}
+}
+test util-13.113 {just over half ulp - 3 digits} {*}{
+    -body {
+        verdonk_test -1b297cad9f70b6 E+96 -134_500000000000000000& E+27
+    }
+    -result {}
+}
+test util-13.114 {just over half ulp - 3 digits} {*}{
+    -body {
+        verdonk_test 1cdc06b20ef183 E-83 +186_50000000000000000& E-27
+    }
+    -result {}
+}
+test util-13.115 {just over half ulp - 3 digits} {*}{
+    -body {
+        verdonk_test -136071dcae4565 E-47 -860_50000000000000000& E-17
+    }
+    -result {}
+}
+test util-13.116 {just over half ulp - 6 digits} {*}{
+    -body {
+        verdonk_test 1cb968d297dde8 E+99 +113788_50000000000000000& E+25
+    }
+    -result {}
+}
+test util-13.117 {just over half ulp - 6 digits} {*}{
+    -body {
+        verdonk_test -11f3e1839eeab1 E+103 -113788_50000000000000000& E+26
+    }
+    -result {}
+}
+test util-13.118 {just under half ulp - 9 digits} {*}{
+    -body {
+        verdonk_test 1e9cec176c96f8 E+117 +317903333_49999999999999999& E+27
+    }
+    -result {}
+}
+test util-13.119 {just over half ulp - 11 digits} {*}{
+    -body {
+        verdonk_test 1c06d366394441 E-36 +25490101686_500000000000000000& E-21
+    }
+    -result {}
+}
+test util-13.120 {just under half ulp - 11 digits} {*}{
+    -body {
+        verdonk_test 1f58ac4db68c90 E+121 +52083105905_49999999999999999& E+26
+    }
+    -result {}
+}
+
+test util-14.1 {funky NaN} {*}{
+    -constraints ieeeFloatingPoint
+    -body {
+	set ieeeValues(-NaN)
+    }
+    -result -NaN
+}
+
+test util-14.2 {funky NaN} {*}{
+    -constraints ieeeFloatingPoint
+    -body {
+	set ieeeValues(-NaN(3456789abcdef))
+    }
+    -result -NaN(3456789abcdef)
+}
+
 # cleanup
 ::tcltest::cleanupTests
 return
Index: unix/Makefile.in
===================================================================
RCS file: /cvsroot/tcl/tcl/unix/Makefile.in,v
retrieving revision 1.310
diff -u -r1.310 Makefile.in
--- unix/Makefile.in	4 Nov 2010 15:55:45 -0000	1.310
+++ unix/Makefile.in	27 Nov 2010 17:40:55 -0000
@@ -322,12 +322,13 @@
 	bn_mp_div_2d.o bn_mp_div_3.o \
         bn_mp_exch.o bn_mp_expt_d.o bn_mp_grow.o bn_mp_init.o \
 	bn_mp_init_copy.o bn_mp_init_multi.o bn_mp_init_set.o \
-	bn_mp_init_size.o bn_mp_karatsuba_mul.o \
+	bn_mp_init_set_int.o bn_mp_init_size.o bn_mp_karatsuba_mul.o \
 	bn_mp_karatsuba_sqr.o \
         bn_mp_lshd.o bn_mp_mod.o bn_mp_mod_2d.o bn_mp_mul.o bn_mp_mul_2.o \
         bn_mp_mul_2d.o bn_mp_mul_d.o bn_mp_neg.o bn_mp_or.o \
 	bn_mp_radix_size.o bn_mp_radix_smap.o \
-        bn_mp_read_radix.o bn_mp_rshd.o bn_mp_set.o bn_mp_shrink.o \
+        bn_mp_read_radix.o bn_mp_rshd.o bn_mp_set.o bn_mp_set_int.o \
+	bn_mp_shrink.o \
 	bn_mp_sqr.o bn_mp_sqrt.o bn_mp_sub.o bn_mp_sub_d.o \
         bn_mp_to_unsigned_bin.o bn_mp_to_unsigned_bin_n.o \
 	bn_mp_toom_mul.o bn_mp_toom_sqr.o bn_mp_toradix_n.o \
@@ -496,6 +497,7 @@
 	$(TOMMATH_DIR)/bn_mp_init_copy.c \
 	$(TOMMATH_DIR)/bn_mp_init_multi.c \
 	$(TOMMATH_DIR)/bn_mp_init_set.c \
+	$(TOMMATH_DIR)/bn_mp_init_set_int.c \
 	$(TOMMATH_DIR)/bn_mp_init_size.c \
 	$(TOMMATH_DIR)/bn_mp_karatsuba_mul.c \
 	$(TOMMATH_DIR)/bn_mp_karatsuba_sqr.c \
@@ -513,6 +515,7 @@
 	$(TOMMATH_DIR)/bn_mp_read_radix.c \
 	$(TOMMATH_DIR)/bn_mp_rshd.c \
 	$(TOMMATH_DIR)/bn_mp_set.c \
+	$(TOMMATH_DIR)/bn_mp_set_int.c \
 	$(TOMMATH_DIR)/bn_mp_shrink.c \
 	$(TOMMATH_DIR)/bn_mp_sqr.c \
 	$(TOMMATH_DIR)/bn_mp_sqrt.c \
@@ -1380,6 +1383,9 @@
 bn_mp_init_set.o: $(TOMMATH_DIR)/bn_mp_init_set.c $(MATHHDRS)
 	$(CC) -c $(CC_SWITCHES) $(TOMMATH_DIR)/bn_mp_init_set.c
 
+bn_mp_init_set_int.o: $(TOMMATH_DIR)/bn_mp_init_set_int.c $(MATHHDRS)
+	$(CC) -c $(CC_SWITCHES) $(TOMMATH_DIR)/bn_mp_init_set_int.c
+
 bn_mp_init_size.o:$(TOMMATH_DIR)/bn_mp_init_size.c $(MATHHDRS)
 	$(CC) -c $(CC_SWITCHES) $(TOMMATH_DIR)/bn_mp_init_size.c
 
@@ -1431,6 +1437,9 @@
 bn_mp_set.o: $(TOMMATH_DIR)/bn_mp_set.c $(MATHHDRS)
 	$(CC) -c $(CC_SWITCHES) $(TOMMATH_DIR)/bn_mp_set.c
 
+bn_mp_set_int.o: $(TOMMATH_DIR)/bn_mp_set_int.c $(MATHHDRS)
+	$(CC) -c $(CC_SWITCHES) $(TOMMATH_DIR)/bn_mp_set_int.c
+
 bn_mp_shrink.o: $(TOMMATH_DIR)/bn_mp_shrink.c $(MATHHDRS)
 	$(CC) -c $(CC_SWITCHES) $(TOMMATH_DIR)/bn_mp_shrink.c
 
Index: win/Makefile.in
===================================================================
RCS file: /cvsroot/tcl/tcl/win/Makefile.in,v
retrieving revision 1.185
diff -u -r1.185 Makefile.in
--- win/Makefile.in	4 Nov 2010 21:48:23 -0000	1.185
+++ win/Makefile.in	27 Nov 2010 17:40:55 -0000
@@ -318,6 +318,7 @@
 	bn_mp_init_copy.${OBJEXT} \
 	bn_mp_init_multi.${OBJEXT} \
 	bn_mp_init_set.${OBJEXT} \
+	bn_mp_init_set_int.${OBJEXT} \
 	bn_mp_init_size.${OBJEXT} \
 	bn_mp_karatsuba_mul.${OBJEXT} \
 	bn_mp_karatsuba_sqr.$(OBJEXT) \
@@ -335,6 +336,7 @@
 	bn_mp_read_radix.${OBJEXT} \
 	bn_mp_rshd.${OBJEXT} \
 	bn_mp_set.${OBJEXT} \
+	bn_mp_set_int.${OBJEXT} \
 	bn_mp_shrink.${OBJEXT} \
 	bn_mp_sqr.${OBJEXT} \
 	bn_mp_sqrt.${OBJEXT} \
Index: win/makefile.vc
===================================================================
RCS file: /cvsroot/tcl/tcl/win/makefile.vc,v
retrieving revision 1.216
diff -u -r1.216 makefile.vc
--- win/makefile.vc	4 Nov 2010 21:48:23 -0000	1.216
+++ win/makefile.vc	27 Nov 2010 17:40:55 -0000
@@ -369,6 +369,7 @@
 	$(TMP_DIR)\bn_mp_init_copy.obj \
 	$(TMP_DIR)\bn_mp_init_multi.obj \
 	$(TMP_DIR)\bn_mp_init_set.obj \
+	$(TMP_DIR)\bn_mp_init_set_int.obj \
 	$(TMP_DIR)\bn_mp_init_size.obj \
 	$(TMP_DIR)\bn_mp_karatsuba_mul.obj \
 	$(TMP_DIR)\bn_mp_karatsuba_sqr.obj \
@@ -386,6 +387,7 @@
 	$(TMP_DIR)\bn_mp_read_radix.obj \
 	$(TMP_DIR)\bn_mp_rshd.obj \
 	$(TMP_DIR)\bn_mp_set.obj \
+	$(TMP_DIR)\bn_mp_set_int.obj \
 	$(TMP_DIR)\bn_mp_shrink.obj \
 	$(TMP_DIR)\bn_mp_sqr.obj \
 	$(TMP_DIR)\bn_mp_sqrt.obj \