From a1c1b0f46865bddd9a84dfae6c146a0e010a4496 Mon Sep 17 00:00:00 2001
From: Tim Peters <tim.peters@gmail.com>
Date: Tue, 4 Sep 2001 02:50:49 +0000
Subject: [PATCH] Introduce new private API function _PyLong_AsScaledDouble. 
 Not used yet, but will be the foundation for Good Things: + Speed
 PyLong_AsDouble. + Give PyLong_AsDouble the ability to detect overflow. +
 Make true division of long/long nearly as accurate as possible (no   spurious
 infinities or NaNs). + Return non-insane results from math.log and math.log10
 when passing a   long that can't be approximated by a double better than
 HUGE_VAL.

---
 Include/longobject.h |  9 ++++++++
 Objects/longobject.c | 52 ++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 61 insertions(+)

diff --git a/Include/longobject.h b/Include/longobject.h
index 5d1ea0a5174..e592891687a 100644
--- a/Include/longobject.h
+++ b/Include/longobject.h
@@ -18,6 +18,15 @@ extern DL_IMPORT(PyObject *) PyLong_FromUnsignedLong(unsigned long);
 extern DL_IMPORT(PyObject *) PyLong_FromDouble(double);
 extern DL_IMPORT(long) PyLong_AsLong(PyObject *);
 extern DL_IMPORT(unsigned long) PyLong_AsUnsignedLong(PyObject *);
+
+/* _PyLong_AsScaledDouble returns a double x and an exponent e such that
+   the true value is approximately equal to x * 2**(SHIFT*e).  e is >= 0.
+   x is 0.0 if and only if the input is 0 (in which case, e and x are both
+   zeroes).  Overflow is impossible.  Note that the exponent returned must
+   be multiplied by SHIFT!  There may not be enough room in an int to store
+   e*SHIFT directly. */
+extern DL_IMPORT(double) _PyLong_AsScaledDouble(PyObject *vv, int *e);
+
 extern DL_IMPORT(double) PyLong_AsDouble(PyObject *);
 extern DL_IMPORT(PyObject *) PyLong_FromVoidPtr(void *);
 extern DL_IMPORT(void *) PyLong_AsVoidPtr(PyObject *);
diff --git a/Objects/longobject.c b/Objects/longobject.c
index 4d4a9599bcf..b511928dbd4 100644
--- a/Objects/longobject.c
+++ b/Objects/longobject.c
@@ -474,6 +474,58 @@ Overflow:
 	
 }
 
+double
+_PyLong_AsScaledDouble(PyObject *vv, int *exponent)
+{
+/* NBITS_WANTED should be > the number of bits in a double's precision,
+   but small enough so that 2**NBITS_WANTED is within the normal double
+   range.  nbitsneeded is set to 1 less than that because the most-significant
+   Python digit contains at least 1 significant bit, but we don't want to
+   bother counting them (catering to the worst case cheaply).
+
+   57 is one more than VAX-D double precision; I (Tim) don't know of a double
+   format with more precision than that; it's 1 larger so that we add in at
+   least one round bit to stand in for the ignored least-significant bits.
+*/
+#define NBITS_WANTED 57
+	PyLongObject *v;
+	double x;
+	const double multiplier = (double)(1L << SHIFT);
+	int i, sign;
+	int nbitsneeded;
+
+	if (vv == NULL || !PyLong_Check(vv)) {
+		PyErr_BadInternalCall();
+		return -1;
+	}
+	v = (PyLongObject *)vv;
+	i = v->ob_size;
+	sign = 1;
+	if (i < 0) {
+		sign = -1;
+		i = -(i);
+	}
+	else if (i == 0) {
+		*exponent = 0;
+		return 0.0;
+	}
+	--i;
+	x = (double)v->ob_digit[i];
+	nbitsneeded = NBITS_WANTED - 1;
+	/* Invariant:  i Python digits remain unaccounted for. */
+	while (i > 0 && nbitsneeded > 0) {
+		--i;
+		x = x * multiplier + (double)v->ob_digit[i];
+		nbitsneeded -= SHIFT;
+	}
+	/* There are i digits we didn't shift in.  Pretending they're all
+	   zeroes, the true value is x * 2**(i*SHIFT). */
+	*exponent = i;
+	assert(x > 0.0);
+	return x * sign;
+#undef NBITS_WANTED
+}
+
 /* Get a C double from a long int object. */
 
 double