Home | History | Annotate | Download | only in Objects
      1 
      2 /* Float object implementation */
      3 
      4 /* XXX There should be overflow checks here, but it's hard to check
      5    for any kind of float exception without losing portability. */
      6 
      7 #include "Python.h"
      8 
      9 #include <ctype.h>
     10 #include <float.h>
     11 
     12 
     13 /* Special free list
     14    free_list is a singly-linked list of available PyFloatObjects, linked
     15    via abuse of their ob_type members.
     16 */
     17 
     18 #ifndef PyFloat_MAXFREELIST
     19 #define PyFloat_MAXFREELIST    100
     20 #endif
     21 static int numfree = 0;
     22 static PyFloatObject *free_list = NULL;
     23 
     24 double
     25 PyFloat_GetMax(void)
     26 {
     27     return DBL_MAX;
     28 }
     29 
     30 double
     31 PyFloat_GetMin(void)
     32 {
     33     return DBL_MIN;
     34 }
     35 
     36 static PyTypeObject FloatInfoType;
     37 
     38 PyDoc_STRVAR(floatinfo__doc__,
     39 "sys.float_info\n\
     40 \n\
     41 A structseq holding information about the float type. It contains low level\n\
     42 information about the precision and internal representation. Please study\n\
     43 your system's :file:`float.h` for more information.");
     44 
     45 static PyStructSequence_Field floatinfo_fields[] = {
     46     {"max",             "DBL_MAX -- maximum representable finite float"},
     47     {"max_exp",         "DBL_MAX_EXP -- maximum int e such that radix**(e-1) "
     48                     "is representable"},
     49     {"max_10_exp",      "DBL_MAX_10_EXP -- maximum int e such that 10**e "
     50                     "is representable"},
     51     {"min",             "DBL_MIN -- Minimum positive normalizer float"},
     52     {"min_exp",         "DBL_MIN_EXP -- minimum int e such that radix**(e-1) "
     53                     "is a normalized float"},
     54     {"min_10_exp",      "DBL_MIN_10_EXP -- minimum int e such that 10**e is "
     55                     "a normalized"},
     56     {"dig",             "DBL_DIG -- digits"},
     57     {"mant_dig",        "DBL_MANT_DIG -- mantissa digits"},
     58     {"epsilon",         "DBL_EPSILON -- Difference between 1 and the next "
     59                     "representable float"},
     60     {"radix",           "FLT_RADIX -- radix of exponent"},
     61     {"rounds",          "FLT_ROUNDS -- addition rounds"},
     62     {0}
     63 };
     64 
     65 static PyStructSequence_Desc floatinfo_desc = {
     66     "sys.float_info",           /* name */
     67     floatinfo__doc__,           /* doc */
     68     floatinfo_fields,           /* fields */
     69     11
     70 };
     71 
     72 PyObject *
     73 PyFloat_GetInfo(void)
     74 {
     75     PyObject* floatinfo;
     76     int pos = 0;
     77 
     78     floatinfo = PyStructSequence_New(&FloatInfoType);
     79     if (floatinfo == NULL) {
     80         return NULL;
     81     }
     82 
     83 #define SetIntFlag(flag) \
     84     PyStructSequence_SET_ITEM(floatinfo, pos++, PyLong_FromLong(flag))
     85 #define SetDblFlag(flag) \
     86     PyStructSequence_SET_ITEM(floatinfo, pos++, PyFloat_FromDouble(flag))
     87 
     88     SetDblFlag(DBL_MAX);
     89     SetIntFlag(DBL_MAX_EXP);
     90     SetIntFlag(DBL_MAX_10_EXP);
     91     SetDblFlag(DBL_MIN);
     92     SetIntFlag(DBL_MIN_EXP);
     93     SetIntFlag(DBL_MIN_10_EXP);
     94     SetIntFlag(DBL_DIG);
     95     SetIntFlag(DBL_MANT_DIG);
     96     SetDblFlag(DBL_EPSILON);
     97     SetIntFlag(FLT_RADIX);
     98     SetIntFlag(FLT_ROUNDS);
     99 #undef SetIntFlag
    100 #undef SetDblFlag
    101 
    102     if (PyErr_Occurred()) {
    103         Py_CLEAR(floatinfo);
    104         return NULL;
    105     }
    106     return floatinfo;
    107 }
    108 
    109 PyObject *
    110 PyFloat_FromDouble(double fval)
    111 {
    112     PyFloatObject *op = free_list;
    113     if (op != NULL) {
    114         free_list = (PyFloatObject *) Py_TYPE(op);
    115         numfree--;
    116     } else {
    117         op = (PyFloatObject*) PyObject_MALLOC(sizeof(PyFloatObject));
    118         if (!op)
    119             return PyErr_NoMemory();
    120     }
    121     /* Inline PyObject_New */
    122     (void)PyObject_INIT(op, &PyFloat_Type);
    123     op->ob_fval = fval;
    124     return (PyObject *) op;
    125 }
    126 
    127 static PyObject *
    128 float_from_string_inner(const char *s, Py_ssize_t len, void *obj)
    129 {
    130     double x;
    131     const char *end;
    132     const char *last = s + len;
    133     /* strip space */
    134     while (s < last && Py_ISSPACE(*s)) {
    135         s++;
    136     }
    137 
    138     while (s < last - 1 && Py_ISSPACE(last[-1])) {
    139         last--;
    140     }
    141 
    142     /* We don't care about overflow or underflow.  If the platform
    143      * supports them, infinities and signed zeroes (on underflow) are
    144      * fine. */
    145     x = PyOS_string_to_double(s, (char **)&end, NULL);
    146     if (end != last) {
    147         PyErr_Format(PyExc_ValueError,
    148                      "could not convert string to float: "
    149                      "%R", obj);
    150         return NULL;
    151     }
    152     else if (x == -1.0 && PyErr_Occurred()) {
    153         return NULL;
    154     }
    155     else {
    156         return PyFloat_FromDouble(x);
    157     }
    158 }
    159 
    160 PyObject *
    161 PyFloat_FromString(PyObject *v)
    162 {
    163     const char *s;
    164     PyObject *s_buffer = NULL;
    165     Py_ssize_t len;
    166     Py_buffer view = {NULL, NULL};
    167     PyObject *result = NULL;
    168 
    169     if (PyUnicode_Check(v)) {
    170         s_buffer = _PyUnicode_TransformDecimalAndSpaceToASCII(v);
    171         if (s_buffer == NULL)
    172             return NULL;
    173         s = PyUnicode_AsUTF8AndSize(s_buffer, &len);
    174         if (s == NULL) {
    175             Py_DECREF(s_buffer);
    176             return NULL;
    177         }
    178     }
    179     else if (PyBytes_Check(v)) {
    180         s = PyBytes_AS_STRING(v);
    181         len = PyBytes_GET_SIZE(v);
    182     }
    183     else if (PyByteArray_Check(v)) {
    184         s = PyByteArray_AS_STRING(v);
    185         len = PyByteArray_GET_SIZE(v);
    186     }
    187     else if (PyObject_GetBuffer(v, &view, PyBUF_SIMPLE) == 0) {
    188         s = (const char *)view.buf;
    189         len = view.len;
    190         /* Copy to NUL-terminated buffer. */
    191         s_buffer = PyBytes_FromStringAndSize(s, len);
    192         if (s_buffer == NULL) {
    193             PyBuffer_Release(&view);
    194             return NULL;
    195         }
    196         s = PyBytes_AS_STRING(s_buffer);
    197     }
    198     else {
    199         PyErr_Format(PyExc_TypeError,
    200             "float() argument must be a string or a number, not '%.200s'",
    201             Py_TYPE(v)->tp_name);
    202         return NULL;
    203     }
    204     result = _Py_string_to_number_with_underscores(s, len, "float", v, v,
    205                                                    float_from_string_inner);
    206     PyBuffer_Release(&view);
    207     Py_XDECREF(s_buffer);
    208     return result;
    209 }
    210 
    211 static void
    212 float_dealloc(PyFloatObject *op)
    213 {
    214     if (PyFloat_CheckExact(op)) {
    215         if (numfree >= PyFloat_MAXFREELIST)  {
    216             PyObject_FREE(op);
    217             return;
    218         }
    219         numfree++;
    220         Py_TYPE(op) = (struct _typeobject *)free_list;
    221         free_list = op;
    222     }
    223     else
    224         Py_TYPE(op)->tp_free((PyObject *)op);
    225 }
    226 
    227 double
    228 PyFloat_AsDouble(PyObject *op)
    229 {
    230     PyNumberMethods *nb;
    231     PyObject *res;
    232     double val;
    233 
    234     if (op == NULL) {
    235         PyErr_BadArgument();
    236         return -1;
    237     }
    238 
    239     if (PyFloat_Check(op)) {
    240         return PyFloat_AS_DOUBLE(op);
    241     }
    242 
    243     nb = Py_TYPE(op)->tp_as_number;
    244     if (nb == NULL || nb->nb_float == NULL) {
    245         PyErr_Format(PyExc_TypeError, "must be real number, not %.50s",
    246                      op->ob_type->tp_name);
    247         return -1;
    248     }
    249 
    250     res = (*nb->nb_float) (op);
    251     if (res == NULL) {
    252         return -1;
    253     }
    254     if (!PyFloat_CheckExact(res)) {
    255         if (!PyFloat_Check(res)) {
    256             PyErr_Format(PyExc_TypeError,
    257                          "%.50s.__float__ returned non-float (type %.50s)",
    258                          op->ob_type->tp_name, res->ob_type->tp_name);
    259             Py_DECREF(res);
    260             return -1;
    261         }
    262         if (PyErr_WarnFormat(PyExc_DeprecationWarning, 1,
    263                 "%.50s.__float__ returned non-float (type %.50s).  "
    264                 "The ability to return an instance of a strict subclass of float "
    265                 "is deprecated, and may be removed in a future version of Python.",
    266                 op->ob_type->tp_name, res->ob_type->tp_name)) {
    267             Py_DECREF(res);
    268             return -1;
    269         }
    270     }
    271 
    272     val = PyFloat_AS_DOUBLE(res);
    273     Py_DECREF(res);
    274     return val;
    275 }
    276 
    277 /* Macro and helper that convert PyObject obj to a C double and store
    278    the value in dbl.  If conversion to double raises an exception, obj is
    279    set to NULL, and the function invoking this macro returns NULL.  If
    280    obj is not of float or int type, Py_NotImplemented is incref'ed,
    281    stored in obj, and returned from the function invoking this macro.
    282 */
    283 #define CONVERT_TO_DOUBLE(obj, dbl)                     \
    284     if (PyFloat_Check(obj))                             \
    285         dbl = PyFloat_AS_DOUBLE(obj);                   \
    286     else if (convert_to_double(&(obj), &(dbl)) < 0)     \
    287         return obj;
    288 
    289 /* Methods */
    290 
    291 static int
    292 convert_to_double(PyObject **v, double *dbl)
    293 {
    294     PyObject *obj = *v;
    295 
    296     if (PyLong_Check(obj)) {
    297         *dbl = PyLong_AsDouble(obj);
    298         if (*dbl == -1.0 && PyErr_Occurred()) {
    299             *v = NULL;
    300             return -1;
    301         }
    302     }
    303     else {
    304         Py_INCREF(Py_NotImplemented);
    305         *v = Py_NotImplemented;
    306         return -1;
    307     }
    308     return 0;
    309 }
    310 
    311 static PyObject *
    312 float_repr(PyFloatObject *v)
    313 {
    314     PyObject *result;
    315     char *buf;
    316 
    317     buf = PyOS_double_to_string(PyFloat_AS_DOUBLE(v),
    318                                 'r', 0,
    319                                 Py_DTSF_ADD_DOT_0,
    320                                 NULL);
    321     if (!buf)
    322         return PyErr_NoMemory();
    323     result = _PyUnicode_FromASCII(buf, strlen(buf));
    324     PyMem_Free(buf);
    325     return result;
    326 }
    327 
    328 /* Comparison is pretty much a nightmare.  When comparing float to float,
    329  * we do it as straightforwardly (and long-windedly) as conceivable, so
    330  * that, e.g., Python x == y delivers the same result as the platform
    331  * C x == y when x and/or y is a NaN.
    332  * When mixing float with an integer type, there's no good *uniform* approach.
    333  * Converting the double to an integer obviously doesn't work, since we
    334  * may lose info from fractional bits.  Converting the integer to a double
    335  * also has two failure modes:  (1) an int may trigger overflow (too
    336  * large to fit in the dynamic range of a C double); (2) even a C long may have
    337  * more bits than fit in a C double (e.g., on a 64-bit box long may have
    338  * 63 bits of precision, but a C double probably has only 53), and then
    339  * we can falsely claim equality when low-order integer bits are lost by
    340  * coercion to double.  So this part is painful too.
    341  */
    342 
    343 static PyObject*
    344 float_richcompare(PyObject *v, PyObject *w, int op)
    345 {
    346     double i, j;
    347     int r = 0;
    348 
    349     assert(PyFloat_Check(v));
    350     i = PyFloat_AS_DOUBLE(v);
    351 
    352     /* Switch on the type of w.  Set i and j to doubles to be compared,
    353      * and op to the richcomp to use.
    354      */
    355     if (PyFloat_Check(w))
    356         j = PyFloat_AS_DOUBLE(w);
    357 
    358     else if (!Py_IS_FINITE(i)) {
    359         if (PyLong_Check(w))
    360             /* If i is an infinity, its magnitude exceeds any
    361              * finite integer, so it doesn't matter which int we
    362              * compare i with.  If i is a NaN, similarly.
    363              */
    364             j = 0.0;
    365         else
    366             goto Unimplemented;
    367     }
    368 
    369     else if (PyLong_Check(w)) {
    370         int vsign = i == 0.0 ? 0 : i < 0.0 ? -1 : 1;
    371         int wsign = _PyLong_Sign(w);
    372         size_t nbits;
    373         int exponent;
    374 
    375         if (vsign != wsign) {
    376             /* Magnitudes are irrelevant -- the signs alone
    377              * determine the outcome.
    378              */
    379             i = (double)vsign;
    380             j = (double)wsign;
    381             goto Compare;
    382         }
    383         /* The signs are the same. */
    384         /* Convert w to a double if it fits.  In particular, 0 fits. */
    385         nbits = _PyLong_NumBits(w);
    386         if (nbits == (size_t)-1 && PyErr_Occurred()) {
    387             /* This long is so large that size_t isn't big enough
    388              * to hold the # of bits.  Replace with little doubles
    389              * that give the same outcome -- w is so large that
    390              * its magnitude must exceed the magnitude of any
    391              * finite float.
    392              */
    393             PyErr_Clear();
    394             i = (double)vsign;
    395             assert(wsign != 0);
    396             j = wsign * 2.0;
    397             goto Compare;
    398         }
    399         if (nbits <= 48) {
    400             j = PyLong_AsDouble(w);
    401             /* It's impossible that <= 48 bits overflowed. */
    402             assert(j != -1.0 || ! PyErr_Occurred());
    403             goto Compare;
    404         }
    405         assert(wsign != 0); /* else nbits was 0 */
    406         assert(vsign != 0); /* if vsign were 0, then since wsign is
    407                              * not 0, we would have taken the
    408                              * vsign != wsign branch at the start */
    409         /* We want to work with non-negative numbers. */
    410         if (vsign < 0) {
    411             /* "Multiply both sides" by -1; this also swaps the
    412              * comparator.
    413              */
    414             i = -i;
    415             op = _Py_SwappedOp[op];
    416         }
    417         assert(i > 0.0);
    418         (void) frexp(i, &exponent);
    419         /* exponent is the # of bits in v before the radix point;
    420          * we know that nbits (the # of bits in w) > 48 at this point
    421          */
    422         if (exponent < 0 || (size_t)exponent < nbits) {
    423             i = 1.0;
    424             j = 2.0;
    425             goto Compare;
    426         }
    427         if ((size_t)exponent > nbits) {
    428             i = 2.0;
    429             j = 1.0;
    430             goto Compare;
    431         }
    432         /* v and w have the same number of bits before the radix
    433          * point.  Construct two ints that have the same comparison
    434          * outcome.
    435          */
    436         {
    437             double fracpart;
    438             double intpart;
    439             PyObject *result = NULL;
    440             PyObject *one = NULL;
    441             PyObject *vv = NULL;
    442             PyObject *ww = w;
    443 
    444             if (wsign < 0) {
    445                 ww = PyNumber_Negative(w);
    446                 if (ww == NULL)
    447                     goto Error;
    448             }
    449             else
    450                 Py_INCREF(ww);
    451 
    452             fracpart = modf(i, &intpart);
    453             vv = PyLong_FromDouble(intpart);
    454             if (vv == NULL)
    455                 goto Error;
    456 
    457             if (fracpart != 0.0) {
    458                 /* Shift left, and or a 1 bit into vv
    459                  * to represent the lost fraction.
    460                  */
    461                 PyObject *temp;
    462 
    463                 one = PyLong_FromLong(1);
    464                 if (one == NULL)
    465                     goto Error;
    466 
    467                 temp = PyNumber_Lshift(ww, one);
    468                 if (temp == NULL)
    469                     goto Error;
    470                 Py_DECREF(ww);
    471                 ww = temp;
    472 
    473                 temp = PyNumber_Lshift(vv, one);
    474                 if (temp == NULL)
    475                     goto Error;
    476                 Py_DECREF(vv);
    477                 vv = temp;
    478 
    479                 temp = PyNumber_Or(vv, one);
    480                 if (temp == NULL)
    481                     goto Error;
    482                 Py_DECREF(vv);
    483                 vv = temp;
    484             }
    485 
    486             r = PyObject_RichCompareBool(vv, ww, op);
    487             if (r < 0)
    488                 goto Error;
    489             result = PyBool_FromLong(r);
    490          Error:
    491             Py_XDECREF(vv);
    492             Py_XDECREF(ww);
    493             Py_XDECREF(one);
    494             return result;
    495         }
    496     } /* else if (PyLong_Check(w)) */
    497 
    498     else        /* w isn't float or int */
    499         goto Unimplemented;
    500 
    501  Compare:
    502     PyFPE_START_PROTECT("richcompare", return NULL)
    503     switch (op) {
    504     case Py_EQ:
    505         r = i == j;
    506         break;
    507     case Py_NE:
    508         r = i != j;
    509         break;
    510     case Py_LE:
    511         r = i <= j;
    512         break;
    513     case Py_GE:
    514         r = i >= j;
    515         break;
    516     case Py_LT:
    517         r = i < j;
    518         break;
    519     case Py_GT:
    520         r = i > j;
    521         break;
    522     }
    523     PyFPE_END_PROTECT(r)
    524     return PyBool_FromLong(r);
    525 
    526  Unimplemented:
    527     Py_RETURN_NOTIMPLEMENTED;
    528 }
    529 
    530 static Py_hash_t
    531 float_hash(PyFloatObject *v)
    532 {
    533     return _Py_HashDouble(v->ob_fval);
    534 }
    535 
    536 static PyObject *
    537 float_add(PyObject *v, PyObject *w)
    538 {
    539     double a,b;
    540     CONVERT_TO_DOUBLE(v, a);
    541     CONVERT_TO_DOUBLE(w, b);
    542     PyFPE_START_PROTECT("add", return 0)
    543     a = a + b;
    544     PyFPE_END_PROTECT(a)
    545     return PyFloat_FromDouble(a);
    546 }
    547 
    548 static PyObject *
    549 float_sub(PyObject *v, PyObject *w)
    550 {
    551     double a,b;
    552     CONVERT_TO_DOUBLE(v, a);
    553     CONVERT_TO_DOUBLE(w, b);
    554     PyFPE_START_PROTECT("subtract", return 0)
    555     a = a - b;
    556     PyFPE_END_PROTECT(a)
    557     return PyFloat_FromDouble(a);
    558 }
    559 
    560 static PyObject *
    561 float_mul(PyObject *v, PyObject *w)
    562 {
    563     double a,b;
    564     CONVERT_TO_DOUBLE(v, a);
    565     CONVERT_TO_DOUBLE(w, b);
    566     PyFPE_START_PROTECT("multiply", return 0)
    567     a = a * b;
    568     PyFPE_END_PROTECT(a)
    569     return PyFloat_FromDouble(a);
    570 }
    571 
    572 static PyObject *
    573 float_div(PyObject *v, PyObject *w)
    574 {
    575     double a,b;
    576     CONVERT_TO_DOUBLE(v, a);
    577     CONVERT_TO_DOUBLE(w, b);
    578     if (b == 0.0) {
    579         PyErr_SetString(PyExc_ZeroDivisionError,
    580                         "float division by zero");
    581         return NULL;
    582     }
    583     PyFPE_START_PROTECT("divide", return 0)
    584     a = a / b;
    585     PyFPE_END_PROTECT(a)
    586     return PyFloat_FromDouble(a);
    587 }
    588 
    589 static PyObject *
    590 float_rem(PyObject *v, PyObject *w)
    591 {
    592     double vx, wx;
    593     double mod;
    594     CONVERT_TO_DOUBLE(v, vx);
    595     CONVERT_TO_DOUBLE(w, wx);
    596     if (wx == 0.0) {
    597         PyErr_SetString(PyExc_ZeroDivisionError,
    598                         "float modulo");
    599         return NULL;
    600     }
    601     PyFPE_START_PROTECT("modulo", return 0)
    602     mod = fmod(vx, wx);
    603     if (mod) {
    604         /* ensure the remainder has the same sign as the denominator */
    605         if ((wx < 0) != (mod < 0)) {
    606             mod += wx;
    607         }
    608     }
    609     else {
    610         /* the remainder is zero, and in the presence of signed zeroes
    611            fmod returns different results across platforms; ensure
    612            it has the same sign as the denominator. */
    613         mod = copysign(0.0, wx);
    614     }
    615     PyFPE_END_PROTECT(mod)
    616     return PyFloat_FromDouble(mod);
    617 }
    618 
    619 static PyObject *
    620 float_divmod(PyObject *v, PyObject *w)
    621 {
    622     double vx, wx;
    623     double div, mod, floordiv;
    624     CONVERT_TO_DOUBLE(v, vx);
    625     CONVERT_TO_DOUBLE(w, wx);
    626     if (wx == 0.0) {
    627         PyErr_SetString(PyExc_ZeroDivisionError, "float divmod()");
    628         return NULL;
    629     }
    630     PyFPE_START_PROTECT("divmod", return 0)
    631     mod = fmod(vx, wx);
    632     /* fmod is typically exact, so vx-mod is *mathematically* an
    633        exact multiple of wx.  But this is fp arithmetic, and fp
    634        vx - mod is an approximation; the result is that div may
    635        not be an exact integral value after the division, although
    636        it will always be very close to one.
    637     */
    638     div = (vx - mod) / wx;
    639     if (mod) {
    640         /* ensure the remainder has the same sign as the denominator */
    641         if ((wx < 0) != (mod < 0)) {
    642             mod += wx;
    643             div -= 1.0;
    644         }
    645     }
    646     else {
    647         /* the remainder is zero, and in the presence of signed zeroes
    648            fmod returns different results across platforms; ensure
    649            it has the same sign as the denominator. */
    650         mod = copysign(0.0, wx);
    651     }
    652     /* snap quotient to nearest integral value */
    653     if (div) {
    654         floordiv = floor(div);
    655         if (div - floordiv > 0.5)
    656             floordiv += 1.0;
    657     }
    658     else {
    659         /* div is zero - get the same sign as the true quotient */
    660         floordiv = copysign(0.0, vx / wx); /* zero w/ sign of vx/wx */
    661     }
    662     PyFPE_END_PROTECT(floordiv)
    663     return Py_BuildValue("(dd)", floordiv, mod);
    664 }
    665 
    666 static PyObject *
    667 float_floor_div(PyObject *v, PyObject *w)
    668 {
    669     PyObject *t, *r;
    670 
    671     t = float_divmod(v, w);
    672     if (t == NULL || t == Py_NotImplemented)
    673         return t;
    674     assert(PyTuple_CheckExact(t));
    675     r = PyTuple_GET_ITEM(t, 0);
    676     Py_INCREF(r);
    677     Py_DECREF(t);
    678     return r;
    679 }
    680 
    681 /* determine whether x is an odd integer or not;  assumes that
    682    x is not an infinity or nan. */
    683 #define DOUBLE_IS_ODD_INTEGER(x) (fmod(fabs(x), 2.0) == 1.0)
    684 
    685 static PyObject *
    686 float_pow(PyObject *v, PyObject *w, PyObject *z)
    687 {
    688     double iv, iw, ix;
    689     int negate_result = 0;
    690 
    691     if ((PyObject *)z != Py_None) {
    692         PyErr_SetString(PyExc_TypeError, "pow() 3rd argument not "
    693             "allowed unless all arguments are integers");
    694         return NULL;
    695     }
    696 
    697     CONVERT_TO_DOUBLE(v, iv);
    698     CONVERT_TO_DOUBLE(w, iw);
    699 
    700     /* Sort out special cases here instead of relying on pow() */
    701     if (iw == 0) {              /* v**0 is 1, even 0**0 */
    702         return PyFloat_FromDouble(1.0);
    703     }
    704     if (Py_IS_NAN(iv)) {        /* nan**w = nan, unless w == 0 */
    705         return PyFloat_FromDouble(iv);
    706     }
    707     if (Py_IS_NAN(iw)) {        /* v**nan = nan, unless v == 1; 1**nan = 1 */
    708         return PyFloat_FromDouble(iv == 1.0 ? 1.0 : iw);
    709     }
    710     if (Py_IS_INFINITY(iw)) {
    711         /* v**inf is: 0.0 if abs(v) < 1; 1.0 if abs(v) == 1; inf if
    712          *     abs(v) > 1 (including case where v infinite)
    713          *
    714          * v**-inf is: inf if abs(v) < 1; 1.0 if abs(v) == 1; 0.0 if
    715          *     abs(v) > 1 (including case where v infinite)
    716          */
    717         iv = fabs(iv);
    718         if (iv == 1.0)
    719             return PyFloat_FromDouble(1.0);
    720         else if ((iw > 0.0) == (iv > 1.0))
    721             return PyFloat_FromDouble(fabs(iw)); /* return inf */
    722         else
    723             return PyFloat_FromDouble(0.0);
    724     }
    725     if (Py_IS_INFINITY(iv)) {
    726         /* (+-inf)**w is: inf for w positive, 0 for w negative; in
    727          *     both cases, we need to add the appropriate sign if w is
    728          *     an odd integer.
    729          */
    730         int iw_is_odd = DOUBLE_IS_ODD_INTEGER(iw);
    731         if (iw > 0.0)
    732             return PyFloat_FromDouble(iw_is_odd ? iv : fabs(iv));
    733         else
    734             return PyFloat_FromDouble(iw_is_odd ?
    735                                       copysign(0.0, iv) : 0.0);
    736     }
    737     if (iv == 0.0) {  /* 0**w is: 0 for w positive, 1 for w zero
    738                          (already dealt with above), and an error
    739                          if w is negative. */
    740         int iw_is_odd = DOUBLE_IS_ODD_INTEGER(iw);
    741         if (iw < 0.0) {
    742             PyErr_SetString(PyExc_ZeroDivisionError,
    743                             "0.0 cannot be raised to a "
    744                             "negative power");
    745             return NULL;
    746         }
    747         /* use correct sign if iw is odd */
    748         return PyFloat_FromDouble(iw_is_odd ? iv : 0.0);
    749     }
    750 
    751     if (iv < 0.0) {
    752         /* Whether this is an error is a mess, and bumps into libm
    753          * bugs so we have to figure it out ourselves.
    754          */
    755         if (iw != floor(iw)) {
    756             /* Negative numbers raised to fractional powers
    757              * become complex.
    758              */
    759             return PyComplex_Type.tp_as_number->nb_power(v, w, z);
    760         }
    761         /* iw is an exact integer, albeit perhaps a very large
    762          * one.  Replace iv by its absolute value and remember
    763          * to negate the pow result if iw is odd.
    764          */
    765         iv = -iv;
    766         negate_result = DOUBLE_IS_ODD_INTEGER(iw);
    767     }
    768 
    769     if (iv == 1.0) { /* 1**w is 1, even 1**inf and 1**nan */
    770         /* (-1) ** large_integer also ends up here.  Here's an
    771          * extract from the comments for the previous
    772          * implementation explaining why this special case is
    773          * necessary:
    774          *
    775          * -1 raised to an exact integer should never be exceptional.
    776          * Alas, some libms (chiefly glibc as of early 2003) return
    777          * NaN and set EDOM on pow(-1, large_int) if the int doesn't
    778          * happen to be representable in a *C* integer.  That's a
    779          * bug.
    780          */
    781         return PyFloat_FromDouble(negate_result ? -1.0 : 1.0);
    782     }
    783 
    784     /* Now iv and iw are finite, iw is nonzero, and iv is
    785      * positive and not equal to 1.0.  We finally allow
    786      * the platform pow to step in and do the rest.
    787      */
    788     errno = 0;
    789     PyFPE_START_PROTECT("pow", return NULL)
    790     ix = pow(iv, iw);
    791     PyFPE_END_PROTECT(ix)
    792     Py_ADJUST_ERANGE1(ix);
    793     if (negate_result)
    794         ix = -ix;
    795 
    796     if (errno != 0) {
    797         /* We don't expect any errno value other than ERANGE, but
    798          * the range of libm bugs appears unbounded.
    799          */
    800         PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError :
    801                              PyExc_ValueError);
    802         return NULL;
    803     }
    804     return PyFloat_FromDouble(ix);
    805 }
    806 
    807 #undef DOUBLE_IS_ODD_INTEGER
    808 
    809 static PyObject *
    810 float_neg(PyFloatObject *v)
    811 {
    812     return PyFloat_FromDouble(-v->ob_fval);
    813 }
    814 
    815 static PyObject *
    816 float_abs(PyFloatObject *v)
    817 {
    818     return PyFloat_FromDouble(fabs(v->ob_fval));
    819 }
    820 
    821 static int
    822 float_bool(PyFloatObject *v)
    823 {
    824     return v->ob_fval != 0.0;
    825 }
    826 
    827 static PyObject *
    828 float_is_integer(PyObject *v)
    829 {
    830     double x = PyFloat_AsDouble(v);
    831     PyObject *o;
    832 
    833     if (x == -1.0 && PyErr_Occurred())
    834         return NULL;
    835     if (!Py_IS_FINITE(x))
    836         Py_RETURN_FALSE;
    837     errno = 0;
    838     PyFPE_START_PROTECT("is_integer", return NULL)
    839     o = (floor(x) == x) ? Py_True : Py_False;
    840     PyFPE_END_PROTECT(x)
    841     if (errno != 0) {
    842         PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError :
    843                              PyExc_ValueError);
    844         return NULL;
    845     }
    846     Py_INCREF(o);
    847     return o;
    848 }
    849 
    850 #if 0
    851 static PyObject *
    852 float_is_inf(PyObject *v)
    853 {
    854     double x = PyFloat_AsDouble(v);
    855     if (x == -1.0 && PyErr_Occurred())
    856         return NULL;
    857     return PyBool_FromLong((long)Py_IS_INFINITY(x));
    858 }
    859 
    860 static PyObject *
    861 float_is_nan(PyObject *v)
    862 {
    863     double x = PyFloat_AsDouble(v);
    864     if (x == -1.0 && PyErr_Occurred())
    865         return NULL;
    866     return PyBool_FromLong((long)Py_IS_NAN(x));
    867 }
    868 
    869 static PyObject *
    870 float_is_finite(PyObject *v)
    871 {
    872     double x = PyFloat_AsDouble(v);
    873     if (x == -1.0 && PyErr_Occurred())
    874         return NULL;
    875     return PyBool_FromLong((long)Py_IS_FINITE(x));
    876 }
    877 #endif
    878 
    879 static PyObject *
    880 float_trunc(PyObject *v)
    881 {
    882     double x = PyFloat_AsDouble(v);
    883     double wholepart;           /* integral portion of x, rounded toward 0 */
    884 
    885     (void)modf(x, &wholepart);
    886     /* Try to get out cheap if this fits in a Python int.  The attempt
    887      * to cast to long must be protected, as C doesn't define what
    888      * happens if the double is too big to fit in a long.  Some rare
    889      * systems raise an exception then (RISCOS was mentioned as one,
    890      * and someone using a non-default option on Sun also bumped into
    891      * that).  Note that checking for >= and <= LONG_{MIN,MAX} would
    892      * still be vulnerable:  if a long has more bits of precision than
    893      * a double, casting MIN/MAX to double may yield an approximation,
    894      * and if that's rounded up, then, e.g., wholepart=LONG_MAX+1 would
    895      * yield true from the C expression wholepart<=LONG_MAX, despite
    896      * that wholepart is actually greater than LONG_MAX.
    897      */
    898     if (LONG_MIN < wholepart && wholepart < LONG_MAX) {
    899         const long aslong = (long)wholepart;
    900         return PyLong_FromLong(aslong);
    901     }
    902     return PyLong_FromDouble(wholepart);
    903 }
    904 
    905 /* double_round: rounds a finite double to the closest multiple of
    906    10**-ndigits; here ndigits is within reasonable bounds (typically, -308 <=
    907    ndigits <= 323).  Returns a Python float, or sets a Python error and
    908    returns NULL on failure (OverflowError and memory errors are possible). */
    909 
    910 #ifndef PY_NO_SHORT_FLOAT_REPR
    911 /* version of double_round that uses the correctly-rounded string<->double
    912    conversions from Python/dtoa.c */
    913 
    914 static PyObject *
    915 double_round(double x, int ndigits) {
    916 
    917     double rounded;
    918     Py_ssize_t buflen, mybuflen=100;
    919     char *buf, *buf_end, shortbuf[100], *mybuf=shortbuf;
    920     int decpt, sign;
    921     PyObject *result = NULL;
    922     _Py_SET_53BIT_PRECISION_HEADER;
    923 
    924     /* round to a decimal string */
    925     _Py_SET_53BIT_PRECISION_START;
    926     buf = _Py_dg_dtoa(x, 3, ndigits, &decpt, &sign, &buf_end);
    927     _Py_SET_53BIT_PRECISION_END;
    928     if (buf == NULL) {
    929         PyErr_NoMemory();
    930         return NULL;
    931     }
    932 
    933     /* Get new buffer if shortbuf is too small.  Space needed <= buf_end -
    934     buf + 8: (1 extra for '0', 1 for sign, 5 for exp, 1 for '\0').  */
    935     buflen = buf_end - buf;
    936     if (buflen + 8 > mybuflen) {
    937         mybuflen = buflen+8;
    938         mybuf = (char *)PyMem_Malloc(mybuflen);
    939         if (mybuf == NULL) {
    940             PyErr_NoMemory();
    941             goto exit;
    942         }
    943     }
    944     /* copy buf to mybuf, adding exponent, sign and leading 0 */
    945     PyOS_snprintf(mybuf, mybuflen, "%s0%se%d", (sign ? "-" : ""),
    946                   buf, decpt - (int)buflen);
    947 
    948     /* and convert the resulting string back to a double */
    949     errno = 0;
    950     _Py_SET_53BIT_PRECISION_START;
    951     rounded = _Py_dg_strtod(mybuf, NULL);
    952     _Py_SET_53BIT_PRECISION_END;
    953     if (errno == ERANGE && fabs(rounded) >= 1.)
    954         PyErr_SetString(PyExc_OverflowError,
    955                         "rounded value too large to represent");
    956     else
    957         result = PyFloat_FromDouble(rounded);
    958 
    959     /* done computing value;  now clean up */
    960     if (mybuf != shortbuf)
    961         PyMem_Free(mybuf);
    962   exit:
    963     _Py_dg_freedtoa(buf);
    964     return result;
    965 }
    966 
    967 #else /* PY_NO_SHORT_FLOAT_REPR */
    968 
    969 /* fallback version, to be used when correctly rounded binary<->decimal
    970    conversions aren't available */
    971 
    972 static PyObject *
    973 double_round(double x, int ndigits) {
    974     double pow1, pow2, y, z;
    975     if (ndigits >= 0) {
    976         if (ndigits > 22) {
    977             /* pow1 and pow2 are each safe from overflow, but
    978                pow1*pow2 ~= pow(10.0, ndigits) might overflow */
    979             pow1 = pow(10.0, (double)(ndigits-22));
    980             pow2 = 1e22;
    981         }
    982         else {
    983             pow1 = pow(10.0, (double)ndigits);
    984             pow2 = 1.0;
    985         }
    986         y = (x*pow1)*pow2;
    987         /* if y overflows, then rounded value is exactly x */
    988         if (!Py_IS_FINITE(y))
    989             return PyFloat_FromDouble(x);
    990     }
    991     else {
    992         pow1 = pow(10.0, (double)-ndigits);
    993         pow2 = 1.0; /* unused; silences a gcc compiler warning */
    994         y = x / pow1;
    995     }
    996 
    997     z = round(y);
    998     if (fabs(y-z) == 0.5)
    999         /* halfway between two integers; use round-half-even */
   1000         z = 2.0*round(y/2.0);
   1001 
   1002     if (ndigits >= 0)
   1003         z = (z / pow2) / pow1;
   1004     else
   1005         z *= pow1;
   1006 
   1007     /* if computation resulted in overflow, raise OverflowError */
   1008     if (!Py_IS_FINITE(z)) {
   1009         PyErr_SetString(PyExc_OverflowError,
   1010                         "overflow occurred during round");
   1011         return NULL;
   1012     }
   1013 
   1014     return PyFloat_FromDouble(z);
   1015 }
   1016 
   1017 #endif /* PY_NO_SHORT_FLOAT_REPR */
   1018 
   1019 /* round a Python float v to the closest multiple of 10**-ndigits */
   1020 
   1021 static PyObject *
   1022 float_round(PyObject *v, PyObject *args)
   1023 {
   1024     double x, rounded;
   1025     PyObject *o_ndigits = NULL;
   1026     Py_ssize_t ndigits;
   1027 
   1028     x = PyFloat_AsDouble(v);
   1029     if (!PyArg_ParseTuple(args, "|O", &o_ndigits))
   1030         return NULL;
   1031     if (o_ndigits == NULL || o_ndigits == Py_None) {
   1032         /* single-argument round or with None ndigits:
   1033          * round to nearest integer */
   1034         rounded = round(x);
   1035         if (fabs(x-rounded) == 0.5)
   1036             /* halfway case: round to even */
   1037             rounded = 2.0*round(x/2.0);
   1038         return PyLong_FromDouble(rounded);
   1039     }
   1040 
   1041     /* interpret second argument as a Py_ssize_t; clips on overflow */
   1042     ndigits = PyNumber_AsSsize_t(o_ndigits, NULL);
   1043     if (ndigits == -1 && PyErr_Occurred())
   1044         return NULL;
   1045 
   1046     /* nans and infinities round to themselves */
   1047     if (!Py_IS_FINITE(x))
   1048         return PyFloat_FromDouble(x);
   1049 
   1050     /* Deal with extreme values for ndigits. For ndigits > NDIGITS_MAX, x
   1051        always rounds to itself.  For ndigits < NDIGITS_MIN, x always
   1052        rounds to +-0.0.  Here 0.30103 is an upper bound for log10(2). */
   1053 #define NDIGITS_MAX ((int)((DBL_MANT_DIG-DBL_MIN_EXP) * 0.30103))
   1054 #define NDIGITS_MIN (-(int)((DBL_MAX_EXP + 1) * 0.30103))
   1055     if (ndigits > NDIGITS_MAX)
   1056         /* return x */
   1057         return PyFloat_FromDouble(x);
   1058     else if (ndigits < NDIGITS_MIN)
   1059         /* return 0.0, but with sign of x */
   1060         return PyFloat_FromDouble(0.0*x);
   1061     else
   1062         /* finite x, and ndigits is not unreasonably large */
   1063         return double_round(x, (int)ndigits);
   1064 #undef NDIGITS_MAX
   1065 #undef NDIGITS_MIN
   1066 }
   1067 
   1068 static PyObject *
   1069 float_float(PyObject *v)
   1070 {
   1071     if (PyFloat_CheckExact(v))
   1072         Py_INCREF(v);
   1073     else
   1074         v = PyFloat_FromDouble(((PyFloatObject *)v)->ob_fval);
   1075     return v;
   1076 }
   1077 
   1078 /* turn ASCII hex characters into integer values and vice versa */
   1079 
   1080 static char
   1081 char_from_hex(int x)
   1082 {
   1083     assert(0 <= x && x < 16);
   1084     return Py_hexdigits[x];
   1085 }
   1086 
   1087 static int
   1088 hex_from_char(char c) {
   1089     int x;
   1090     switch(c) {
   1091     case '0':
   1092         x = 0;
   1093         break;
   1094     case '1':
   1095         x = 1;
   1096         break;
   1097     case '2':
   1098         x = 2;
   1099         break;
   1100     case '3':
   1101         x = 3;
   1102         break;
   1103     case '4':
   1104         x = 4;
   1105         break;
   1106     case '5':
   1107         x = 5;
   1108         break;
   1109     case '6':
   1110         x = 6;
   1111         break;
   1112     case '7':
   1113         x = 7;
   1114         break;
   1115     case '8':
   1116         x = 8;
   1117         break;
   1118     case '9':
   1119         x = 9;
   1120         break;
   1121     case 'a':
   1122     case 'A':
   1123         x = 10;
   1124         break;
   1125     case 'b':
   1126     case 'B':
   1127         x = 11;
   1128         break;
   1129     case 'c':
   1130     case 'C':
   1131         x = 12;
   1132         break;
   1133     case 'd':
   1134     case 'D':
   1135         x = 13;
   1136         break;
   1137     case 'e':
   1138     case 'E':
   1139         x = 14;
   1140         break;
   1141     case 'f':
   1142     case 'F':
   1143         x = 15;
   1144         break;
   1145     default:
   1146         x = -1;
   1147         break;
   1148     }
   1149     return x;
   1150 }
   1151 
   1152 /* convert a float to a hexadecimal string */
   1153 
   1154 /* TOHEX_NBITS is DBL_MANT_DIG rounded up to the next integer
   1155    of the form 4k+1. */
   1156 #define TOHEX_NBITS DBL_MANT_DIG + 3 - (DBL_MANT_DIG+2)%4
   1157 
   1158 static PyObject *
   1159 float_hex(PyObject *v)
   1160 {
   1161     double x, m;
   1162     int e, shift, i, si, esign;
   1163     /* Space for 1+(TOHEX_NBITS-1)/4 digits, a decimal point, and the
   1164        trailing NUL byte. */
   1165     char s[(TOHEX_NBITS-1)/4+3];
   1166 
   1167     CONVERT_TO_DOUBLE(v, x);
   1168 
   1169     if (Py_IS_NAN(x) || Py_IS_INFINITY(x))
   1170         return float_repr((PyFloatObject *)v);
   1171 
   1172     if (x == 0.0) {
   1173         if (copysign(1.0, x) == -1.0)
   1174             return PyUnicode_FromString("-0x0.0p+0");
   1175         else
   1176             return PyUnicode_FromString("0x0.0p+0");
   1177     }
   1178 
   1179     m = frexp(fabs(x), &e);
   1180     shift = 1 - Py_MAX(DBL_MIN_EXP - e, 0);
   1181     m = ldexp(m, shift);
   1182     e -= shift;
   1183 
   1184     si = 0;
   1185     s[si] = char_from_hex((int)m);
   1186     si++;
   1187     m -= (int)m;
   1188     s[si] = '.';
   1189     si++;
   1190     for (i=0; i < (TOHEX_NBITS-1)/4; i++) {
   1191         m *= 16.0;
   1192         s[si] = char_from_hex((int)m);
   1193         si++;
   1194         m -= (int)m;
   1195     }
   1196     s[si] = '\0';
   1197 
   1198     if (e < 0) {
   1199         esign = (int)'-';
   1200         e = -e;
   1201     }
   1202     else
   1203         esign = (int)'+';
   1204 
   1205     if (x < 0.0)
   1206         return PyUnicode_FromFormat("-0x%sp%c%d", s, esign, e);
   1207     else
   1208         return PyUnicode_FromFormat("0x%sp%c%d", s, esign, e);
   1209 }
   1210 
   1211 PyDoc_STRVAR(float_hex_doc,
   1212 "float.hex() -> string\n\
   1213 \n\
   1214 Return a hexadecimal representation of a floating-point number.\n\
   1215 >>> (-0.1).hex()\n\
   1216 '-0x1.999999999999ap-4'\n\
   1217 >>> 3.14159.hex()\n\
   1218 '0x1.921f9f01b866ep+1'");
   1219 
   1220 /* Convert a hexadecimal string to a float. */
   1221 
   1222 static PyObject *
   1223 float_fromhex(PyObject *cls, PyObject *arg)
   1224 {
   1225     PyObject *result;
   1226     double x;
   1227     long exp, top_exp, lsb, key_digit;
   1228     char *s, *coeff_start, *s_store, *coeff_end, *exp_start, *s_end;
   1229     int half_eps, digit, round_up, negate=0;
   1230     Py_ssize_t length, ndigits, fdigits, i;
   1231 
   1232     /*
   1233      * For the sake of simplicity and correctness, we impose an artificial
   1234      * limit on ndigits, the total number of hex digits in the coefficient
   1235      * The limit is chosen to ensure that, writing exp for the exponent,
   1236      *
   1237      *   (1) if exp > LONG_MAX/2 then the value of the hex string is
   1238      *   guaranteed to overflow (provided it's nonzero)
   1239      *
   1240      *   (2) if exp < LONG_MIN/2 then the value of the hex string is
   1241      *   guaranteed to underflow to 0.
   1242      *
   1243      *   (3) if LONG_MIN/2 <= exp <= LONG_MAX/2 then there's no danger of
   1244      *   overflow in the calculation of exp and top_exp below.
   1245      *
   1246      * More specifically, ndigits is assumed to satisfy the following
   1247      * inequalities:
   1248      *
   1249      *   4*ndigits <= DBL_MIN_EXP - DBL_MANT_DIG - LONG_MIN/2
   1250      *   4*ndigits <= LONG_MAX/2 + 1 - DBL_MAX_EXP
   1251      *
   1252      * If either of these inequalities is not satisfied, a ValueError is
   1253      * raised.  Otherwise, write x for the value of the hex string, and
   1254      * assume x is nonzero.  Then
   1255      *
   1256      *   2**(exp-4*ndigits) <= |x| < 2**(exp+4*ndigits).
   1257      *
   1258      * Now if exp > LONG_MAX/2 then:
   1259      *
   1260      *   exp - 4*ndigits >= LONG_MAX/2 + 1 - (LONG_MAX/2 + 1 - DBL_MAX_EXP)
   1261      *                    = DBL_MAX_EXP
   1262      *
   1263      * so |x| >= 2**DBL_MAX_EXP, which is too large to be stored in C
   1264      * double, so overflows.  If exp < LONG_MIN/2, then
   1265      *
   1266      *   exp + 4*ndigits <= LONG_MIN/2 - 1 + (
   1267      *                      DBL_MIN_EXP - DBL_MANT_DIG - LONG_MIN/2)
   1268      *                    = DBL_MIN_EXP - DBL_MANT_DIG - 1
   1269      *
   1270      * and so |x| < 2**(DBL_MIN_EXP-DBL_MANT_DIG-1), hence underflows to 0
   1271      * when converted to a C double.
   1272      *
   1273      * It's easy to show that if LONG_MIN/2 <= exp <= LONG_MAX/2 then both
   1274      * exp+4*ndigits and exp-4*ndigits are within the range of a long.
   1275      */
   1276 
   1277     s = PyUnicode_AsUTF8AndSize(arg, &length);
   1278     if (s == NULL)
   1279         return NULL;
   1280     s_end = s + length;
   1281 
   1282     /********************
   1283      * Parse the string *
   1284      ********************/
   1285 
   1286     /* leading whitespace */
   1287     while (Py_ISSPACE(*s))
   1288         s++;
   1289 
   1290     /* infinities and nans */
   1291     x = _Py_parse_inf_or_nan(s, &coeff_end);
   1292     if (coeff_end != s) {
   1293         s = coeff_end;
   1294         goto finished;
   1295     }
   1296 
   1297     /* optional sign */
   1298     if (*s == '-') {
   1299         s++;
   1300         negate = 1;
   1301     }
   1302     else if (*s == '+')
   1303         s++;
   1304 
   1305     /* [0x] */
   1306     s_store = s;
   1307     if (*s == '0') {
   1308         s++;
   1309         if (*s == 'x' || *s == 'X')
   1310             s++;
   1311         else
   1312             s = s_store;
   1313     }
   1314 
   1315     /* coefficient: <integer> [. <fraction>] */
   1316     coeff_start = s;
   1317     while (hex_from_char(*s) >= 0)
   1318         s++;
   1319     s_store = s;
   1320     if (*s == '.') {
   1321         s++;
   1322         while (hex_from_char(*s) >= 0)
   1323             s++;
   1324         coeff_end = s-1;
   1325     }
   1326     else
   1327         coeff_end = s;
   1328 
   1329     /* ndigits = total # of hex digits; fdigits = # after point */
   1330     ndigits = coeff_end - coeff_start;
   1331     fdigits = coeff_end - s_store;
   1332     if (ndigits == 0)
   1333         goto parse_error;
   1334     if (ndigits > Py_MIN(DBL_MIN_EXP - DBL_MANT_DIG - LONG_MIN/2,
   1335                          LONG_MAX/2 + 1 - DBL_MAX_EXP)/4)
   1336         goto insane_length_error;
   1337 
   1338     /* [p <exponent>] */
   1339     if (*s == 'p' || *s == 'P') {
   1340         s++;
   1341         exp_start = s;
   1342         if (*s == '-' || *s == '+')
   1343             s++;
   1344         if (!('0' <= *s && *s <= '9'))
   1345             goto parse_error;
   1346         s++;
   1347         while ('0' <= *s && *s <= '9')
   1348             s++;
   1349         exp = strtol(exp_start, NULL, 10);
   1350     }
   1351     else
   1352         exp = 0;
   1353 
   1354 /* for 0 <= j < ndigits, HEX_DIGIT(j) gives the jth most significant digit */
   1355 #define HEX_DIGIT(j) hex_from_char(*((j) < fdigits ?            \
   1356                      coeff_end-(j) :                                    \
   1357                      coeff_end-1-(j)))
   1358 
   1359     /*******************************************
   1360      * Compute rounded value of the hex string *
   1361      *******************************************/
   1362 
   1363     /* Discard leading zeros, and catch extreme overflow and underflow */
   1364     while (ndigits > 0 && HEX_DIGIT(ndigits-1) == 0)
   1365         ndigits--;
   1366     if (ndigits == 0 || exp < LONG_MIN/2) {
   1367         x = 0.0;
   1368         goto finished;
   1369     }
   1370     if (exp > LONG_MAX/2)
   1371         goto overflow_error;
   1372 
   1373     /* Adjust exponent for fractional part. */
   1374     exp = exp - 4*((long)fdigits);
   1375 
   1376     /* top_exp = 1 more than exponent of most sig. bit of coefficient */
   1377     top_exp = exp + 4*((long)ndigits - 1);
   1378     for (digit = HEX_DIGIT(ndigits-1); digit != 0; digit /= 2)
   1379         top_exp++;
   1380 
   1381     /* catch almost all nonextreme cases of overflow and underflow here */
   1382     if (top_exp < DBL_MIN_EXP - DBL_MANT_DIG) {
   1383         x = 0.0;
   1384         goto finished;
   1385     }
   1386     if (top_exp > DBL_MAX_EXP)
   1387         goto overflow_error;
   1388 
   1389     /* lsb = exponent of least significant bit of the *rounded* value.
   1390        This is top_exp - DBL_MANT_DIG unless result is subnormal. */
   1391     lsb = Py_MAX(top_exp, (long)DBL_MIN_EXP) - DBL_MANT_DIG;
   1392 
   1393     x = 0.0;
   1394     if (exp >= lsb) {
   1395         /* no rounding required */
   1396         for (i = ndigits-1; i >= 0; i--)
   1397             x = 16.0*x + HEX_DIGIT(i);
   1398         x = ldexp(x, (int)(exp));
   1399         goto finished;
   1400     }
   1401     /* rounding required.  key_digit is the index of the hex digit
   1402        containing the first bit to be rounded away. */
   1403     half_eps = 1 << (int)((lsb - exp - 1) % 4);
   1404     key_digit = (lsb - exp - 1) / 4;
   1405     for (i = ndigits-1; i > key_digit; i--)
   1406         x = 16.0*x + HEX_DIGIT(i);
   1407     digit = HEX_DIGIT(key_digit);
   1408     x = 16.0*x + (double)(digit & (16-2*half_eps));
   1409 
   1410     /* round-half-even: round up if bit lsb-1 is 1 and at least one of
   1411        bits lsb, lsb-2, lsb-3, lsb-4, ... is 1. */
   1412     if ((digit & half_eps) != 0) {
   1413         round_up = 0;
   1414         if ((digit & (3*half_eps-1)) != 0 ||
   1415             (half_eps == 8 && (HEX_DIGIT(key_digit+1) & 1) != 0))
   1416             round_up = 1;
   1417         else
   1418             for (i = key_digit-1; i >= 0; i--)
   1419                 if (HEX_DIGIT(i) != 0) {
   1420                     round_up = 1;
   1421                     break;
   1422                 }
   1423         if (round_up) {
   1424             x += 2*half_eps;
   1425             if (top_exp == DBL_MAX_EXP &&
   1426                 x == ldexp((double)(2*half_eps), DBL_MANT_DIG))
   1427                 /* overflow corner case: pre-rounded value <
   1428                    2**DBL_MAX_EXP; rounded=2**DBL_MAX_EXP. */
   1429                 goto overflow_error;
   1430         }
   1431     }
   1432     x = ldexp(x, (int)(exp+4*key_digit));
   1433 
   1434   finished:
   1435     /* optional trailing whitespace leading to the end of the string */
   1436     while (Py_ISSPACE(*s))
   1437         s++;
   1438     if (s != s_end)
   1439         goto parse_error;
   1440     result = PyFloat_FromDouble(negate ? -x : x);
   1441     if (cls != (PyObject *)&PyFloat_Type && result != NULL) {
   1442         Py_SETREF(result, PyObject_CallFunctionObjArgs(cls, result, NULL));
   1443     }
   1444     return result;
   1445 
   1446   overflow_error:
   1447     PyErr_SetString(PyExc_OverflowError,
   1448                     "hexadecimal value too large to represent as a float");
   1449     return NULL;
   1450 
   1451   parse_error:
   1452     PyErr_SetString(PyExc_ValueError,
   1453                     "invalid hexadecimal floating-point string");
   1454     return NULL;
   1455 
   1456   insane_length_error:
   1457     PyErr_SetString(PyExc_ValueError,
   1458                     "hexadecimal string too long to convert");
   1459     return NULL;
   1460 }
   1461 
   1462 PyDoc_STRVAR(float_fromhex_doc,
   1463 "float.fromhex(string) -> float\n\
   1464 \n\
   1465 Create a floating-point number from a hexadecimal string.\n\
   1466 >>> float.fromhex('0x1.ffffp10')\n\
   1467 2047.984375\n\
   1468 >>> float.fromhex('-0x1p-1074')\n\
   1469 -5e-324");
   1470 
   1471 
   1472 static PyObject *
   1473 float_as_integer_ratio(PyObject *v, PyObject *unused)
   1474 {
   1475     double self;
   1476     double float_part;
   1477     int exponent;
   1478     int i;
   1479 
   1480     PyObject *py_exponent = NULL;
   1481     PyObject *numerator = NULL;
   1482     PyObject *denominator = NULL;
   1483     PyObject *result_pair = NULL;
   1484     PyNumberMethods *long_methods = PyLong_Type.tp_as_number;
   1485 
   1486     CONVERT_TO_DOUBLE(v, self);
   1487 
   1488     if (Py_IS_INFINITY(self)) {
   1489         PyErr_SetString(PyExc_OverflowError,
   1490                         "cannot convert Infinity to integer ratio");
   1491         return NULL;
   1492     }
   1493     if (Py_IS_NAN(self)) {
   1494         PyErr_SetString(PyExc_ValueError,
   1495                         "cannot convert NaN to integer ratio");
   1496         return NULL;
   1497     }
   1498 
   1499     PyFPE_START_PROTECT("as_integer_ratio", goto error);
   1500     float_part = frexp(self, &exponent);        /* self == float_part * 2**exponent exactly */
   1501     PyFPE_END_PROTECT(float_part);
   1502 
   1503     for (i=0; i<300 && float_part != floor(float_part) ; i++) {
   1504         float_part *= 2.0;
   1505         exponent--;
   1506     }
   1507     /* self == float_part * 2**exponent exactly and float_part is integral.
   1508        If FLT_RADIX != 2, the 300 steps may leave a tiny fractional part
   1509        to be truncated by PyLong_FromDouble(). */
   1510 
   1511     numerator = PyLong_FromDouble(float_part);
   1512     if (numerator == NULL)
   1513         goto error;
   1514     denominator = PyLong_FromLong(1);
   1515     if (denominator == NULL)
   1516         goto error;
   1517     py_exponent = PyLong_FromLong(Py_ABS(exponent));
   1518     if (py_exponent == NULL)
   1519         goto error;
   1520 
   1521     /* fold in 2**exponent */
   1522     if (exponent > 0) {
   1523         Py_SETREF(numerator,
   1524                   long_methods->nb_lshift(numerator, py_exponent));
   1525         if (numerator == NULL)
   1526             goto error;
   1527     }
   1528     else {
   1529         Py_SETREF(denominator,
   1530                   long_methods->nb_lshift(denominator, py_exponent));
   1531         if (denominator == NULL)
   1532             goto error;
   1533     }
   1534 
   1535     result_pair = PyTuple_Pack(2, numerator, denominator);
   1536 
   1537 error:
   1538     Py_XDECREF(py_exponent);
   1539     Py_XDECREF(denominator);
   1540     Py_XDECREF(numerator);
   1541     return result_pair;
   1542 }
   1543 
   1544 PyDoc_STRVAR(float_as_integer_ratio_doc,
   1545 "float.as_integer_ratio() -> (int, int)\n"
   1546 "\n"
   1547 "Return a pair of integers, whose ratio is exactly equal to the original\n"
   1548 "float and with a positive denominator.\n"
   1549 "Raise OverflowError on infinities and a ValueError on NaNs.\n"
   1550 "\n"
   1551 ">>> (10.0).as_integer_ratio()\n"
   1552 "(10, 1)\n"
   1553 ">>> (0.0).as_integer_ratio()\n"
   1554 "(0, 1)\n"
   1555 ">>> (-.25).as_integer_ratio()\n"
   1556 "(-1, 4)");
   1557 
   1558 
   1559 static PyObject *
   1560 float_subtype_new(PyTypeObject *type, PyObject *args, PyObject *kwds);
   1561 
   1562 static PyObject *
   1563 float_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
   1564 {
   1565     PyObject *x = Py_False; /* Integer zero */
   1566     static char *kwlist[] = {"x", 0};
   1567 
   1568     if (type != &PyFloat_Type)
   1569         return float_subtype_new(type, args, kwds); /* Wimp out */
   1570     if (!PyArg_ParseTupleAndKeywords(args, kwds, "|O:float", kwlist, &x))
   1571         return NULL;
   1572     /* If it's a string, but not a string subclass, use
   1573        PyFloat_FromString. */
   1574     if (PyUnicode_CheckExact(x))
   1575         return PyFloat_FromString(x);
   1576     return PyNumber_Float(x);
   1577 }
   1578 
   1579 /* Wimpy, slow approach to tp_new calls for subtypes of float:
   1580    first create a regular float from whatever arguments we got,
   1581    then allocate a subtype instance and initialize its ob_fval
   1582    from the regular float.  The regular float is then thrown away.
   1583 */
   1584 static PyObject *
   1585 float_subtype_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
   1586 {
   1587     PyObject *tmp, *newobj;
   1588 
   1589     assert(PyType_IsSubtype(type, &PyFloat_Type));
   1590     tmp = float_new(&PyFloat_Type, args, kwds);
   1591     if (tmp == NULL)
   1592         return NULL;
   1593     assert(PyFloat_Check(tmp));
   1594     newobj = type->tp_alloc(type, 0);
   1595     if (newobj == NULL) {
   1596         Py_DECREF(tmp);
   1597         return NULL;
   1598     }
   1599     ((PyFloatObject *)newobj)->ob_fval = ((PyFloatObject *)tmp)->ob_fval;
   1600     Py_DECREF(tmp);
   1601     return newobj;
   1602 }
   1603 
   1604 static PyObject *
   1605 float_getnewargs(PyFloatObject *v)
   1606 {
   1607     return Py_BuildValue("(d)", v->ob_fval);
   1608 }
   1609 
   1610 /* this is for the benefit of the pack/unpack routines below */
   1611 
   1612 typedef enum {
   1613     unknown_format, ieee_big_endian_format, ieee_little_endian_format
   1614 } float_format_type;
   1615 
   1616 static float_format_type double_format, float_format;
   1617 static float_format_type detected_double_format, detected_float_format;
   1618 
   1619 static PyObject *
   1620 float_getformat(PyTypeObject *v, PyObject* arg)
   1621 {
   1622     char* s;
   1623     float_format_type r;
   1624 
   1625     if (!PyUnicode_Check(arg)) {
   1626         PyErr_Format(PyExc_TypeError,
   1627          "__getformat__() argument must be string, not %.500s",
   1628                          Py_TYPE(arg)->tp_name);
   1629         return NULL;
   1630     }
   1631     s = PyUnicode_AsUTF8(arg);
   1632     if (s == NULL)
   1633         return NULL;
   1634     if (strcmp(s, "double") == 0) {
   1635         r = double_format;
   1636     }
   1637     else if (strcmp(s, "float") == 0) {
   1638         r = float_format;
   1639     }
   1640     else {
   1641         PyErr_SetString(PyExc_ValueError,
   1642                         "__getformat__() argument 1 must be "
   1643                         "'double' or 'float'");
   1644         return NULL;
   1645     }
   1646 
   1647     switch (r) {
   1648     case unknown_format:
   1649         return PyUnicode_FromString("unknown");
   1650     case ieee_little_endian_format:
   1651         return PyUnicode_FromString("IEEE, little-endian");
   1652     case ieee_big_endian_format:
   1653         return PyUnicode_FromString("IEEE, big-endian");
   1654     default:
   1655         Py_FatalError("insane float_format or double_format");
   1656         return NULL;
   1657     }
   1658 }
   1659 
   1660 PyDoc_STRVAR(float_getformat_doc,
   1661 "float.__getformat__(typestr) -> string\n"
   1662 "\n"
   1663 "You probably don't want to use this function.  It exists mainly to be\n"
   1664 "used in Python's test suite.\n"
   1665 "\n"
   1666 "typestr must be 'double' or 'float'.  This function returns whichever of\n"
   1667 "'unknown', 'IEEE, big-endian' or 'IEEE, little-endian' best describes the\n"
   1668 "format of floating point numbers used by the C type named by typestr.");
   1669 
   1670 static PyObject *
   1671 float_setformat(PyTypeObject *v, PyObject* args)
   1672 {
   1673     char* typestr;
   1674     char* format;
   1675     float_format_type f;
   1676     float_format_type detected;
   1677     float_format_type *p;
   1678 
   1679     if (!PyArg_ParseTuple(args, "ss:__setformat__", &typestr, &format))
   1680         return NULL;
   1681 
   1682     if (strcmp(typestr, "double") == 0) {
   1683         p = &double_format;
   1684         detected = detected_double_format;
   1685     }
   1686     else if (strcmp(typestr, "float") == 0) {
   1687         p = &float_format;
   1688         detected = detected_float_format;
   1689     }
   1690     else {
   1691         PyErr_SetString(PyExc_ValueError,
   1692                         "__setformat__() argument 1 must "
   1693                         "be 'double' or 'float'");
   1694         return NULL;
   1695     }
   1696 
   1697     if (strcmp(format, "unknown") == 0) {
   1698         f = unknown_format;
   1699     }
   1700     else if (strcmp(format, "IEEE, little-endian") == 0) {
   1701         f = ieee_little_endian_format;
   1702     }
   1703     else if (strcmp(format, "IEEE, big-endian") == 0) {
   1704         f = ieee_big_endian_format;
   1705     }
   1706     else {
   1707         PyErr_SetString(PyExc_ValueError,
   1708                         "__setformat__() argument 2 must be "
   1709                         "'unknown', 'IEEE, little-endian' or "
   1710                         "'IEEE, big-endian'");
   1711         return NULL;
   1712 
   1713     }
   1714 
   1715     if (f != unknown_format && f != detected) {
   1716         PyErr_Format(PyExc_ValueError,
   1717                      "can only set %s format to 'unknown' or the "
   1718                      "detected platform value", typestr);
   1719         return NULL;
   1720     }
   1721 
   1722     *p = f;
   1723     Py_RETURN_NONE;
   1724 }
   1725 
   1726 PyDoc_STRVAR(float_setformat_doc,
   1727 "float.__setformat__(typestr, fmt) -> None\n"
   1728 "\n"
   1729 "You probably don't want to use this function.  It exists mainly to be\n"
   1730 "used in Python's test suite.\n"
   1731 "\n"
   1732 "typestr must be 'double' or 'float'.  fmt must be one of 'unknown',\n"
   1733 "'IEEE, big-endian' or 'IEEE, little-endian', and in addition can only be\n"
   1734 "one of the latter two if it appears to match the underlying C reality.\n"
   1735 "\n"
   1736 "Override the automatic determination of C-level floating point type.\n"
   1737 "This affects how floats are converted to and from binary strings.");
   1738 
   1739 static PyObject *
   1740 float_getzero(PyObject *v, void *closure)
   1741 {
   1742     return PyFloat_FromDouble(0.0);
   1743 }
   1744 
   1745 static PyObject *
   1746 float__format__(PyObject *self, PyObject *args)
   1747 {
   1748     PyObject *format_spec;
   1749     _PyUnicodeWriter writer;
   1750     int ret;
   1751 
   1752     if (!PyArg_ParseTuple(args, "U:__format__", &format_spec))
   1753         return NULL;
   1754 
   1755     _PyUnicodeWriter_Init(&writer);
   1756     ret = _PyFloat_FormatAdvancedWriter(
   1757         &writer,
   1758         self,
   1759         format_spec, 0, PyUnicode_GET_LENGTH(format_spec));
   1760     if (ret == -1) {
   1761         _PyUnicodeWriter_Dealloc(&writer);
   1762         return NULL;
   1763     }
   1764     return _PyUnicodeWriter_Finish(&writer);
   1765 }
   1766 
   1767 PyDoc_STRVAR(float__format__doc,
   1768 "float.__format__(format_spec) -> string\n"
   1769 "\n"
   1770 "Formats the float according to format_spec.");
   1771 
   1772 
   1773 static PyMethodDef float_methods[] = {
   1774     {"conjugate",       (PyCFunction)float_float,       METH_NOARGS,
   1775      "Return self, the complex conjugate of any float."},
   1776     {"__trunc__",       (PyCFunction)float_trunc, METH_NOARGS,
   1777      "Return the Integral closest to x between 0 and x."},
   1778     {"__round__",       (PyCFunction)float_round, METH_VARARGS,
   1779      "Return the Integral closest to x, rounding half toward even.\n"
   1780      "When an argument is passed, work like built-in round(x, ndigits)."},
   1781     {"as_integer_ratio", (PyCFunction)float_as_integer_ratio, METH_NOARGS,
   1782      float_as_integer_ratio_doc},
   1783     {"fromhex", (PyCFunction)float_fromhex,
   1784      METH_O|METH_CLASS, float_fromhex_doc},
   1785     {"hex", (PyCFunction)float_hex,
   1786      METH_NOARGS, float_hex_doc},
   1787     {"is_integer",      (PyCFunction)float_is_integer,  METH_NOARGS,
   1788      "Return True if the float is an integer."},
   1789 #if 0
   1790     {"is_inf",          (PyCFunction)float_is_inf,      METH_NOARGS,
   1791      "Return True if the float is positive or negative infinite."},
   1792     {"is_finite",       (PyCFunction)float_is_finite,   METH_NOARGS,
   1793      "Return True if the float is finite, neither infinite nor NaN."},
   1794     {"is_nan",          (PyCFunction)float_is_nan,      METH_NOARGS,
   1795      "Return True if the float is not a number (NaN)."},
   1796 #endif
   1797     {"__getnewargs__",          (PyCFunction)float_getnewargs,  METH_NOARGS},
   1798     {"__getformat__",           (PyCFunction)float_getformat,
   1799      METH_O|METH_CLASS,                 float_getformat_doc},
   1800     {"__setformat__",           (PyCFunction)float_setformat,
   1801      METH_VARARGS|METH_CLASS,           float_setformat_doc},
   1802     {"__format__",          (PyCFunction)float__format__,
   1803      METH_VARARGS,                  float__format__doc},
   1804     {NULL,              NULL}           /* sentinel */
   1805 };
   1806 
   1807 static PyGetSetDef float_getset[] = {
   1808     {"real",
   1809      (getter)float_float, (setter)NULL,
   1810      "the real part of a complex number",
   1811      NULL},
   1812     {"imag",
   1813      (getter)float_getzero, (setter)NULL,
   1814      "the imaginary part of a complex number",
   1815      NULL},
   1816     {NULL}  /* Sentinel */
   1817 };
   1818 
   1819 PyDoc_STRVAR(float_doc,
   1820 "float(x) -> floating point number\n\
   1821 \n\
   1822 Convert a string or number to a floating point number, if possible.");
   1823 
   1824 
   1825 static PyNumberMethods float_as_number = {
   1826     float_add,          /*nb_add*/
   1827     float_sub,          /*nb_subtract*/
   1828     float_mul,          /*nb_multiply*/
   1829     float_rem,          /*nb_remainder*/
   1830     float_divmod,       /*nb_divmod*/
   1831     float_pow,          /*nb_power*/
   1832     (unaryfunc)float_neg, /*nb_negative*/
   1833     (unaryfunc)float_float, /*nb_positive*/
   1834     (unaryfunc)float_abs, /*nb_absolute*/
   1835     (inquiry)float_bool, /*nb_bool*/
   1836     0,                  /*nb_invert*/
   1837     0,                  /*nb_lshift*/
   1838     0,                  /*nb_rshift*/
   1839     0,                  /*nb_and*/
   1840     0,                  /*nb_xor*/
   1841     0,                  /*nb_or*/
   1842     float_trunc,        /*nb_int*/
   1843     0,                  /*nb_reserved*/
   1844     float_float,        /*nb_float*/
   1845     0,                  /* nb_inplace_add */
   1846     0,                  /* nb_inplace_subtract */
   1847     0,                  /* nb_inplace_multiply */
   1848     0,                  /* nb_inplace_remainder */
   1849     0,                  /* nb_inplace_power */
   1850     0,                  /* nb_inplace_lshift */
   1851     0,                  /* nb_inplace_rshift */
   1852     0,                  /* nb_inplace_and */
   1853     0,                  /* nb_inplace_xor */
   1854     0,                  /* nb_inplace_or */
   1855     float_floor_div, /* nb_floor_divide */
   1856     float_div,          /* nb_true_divide */
   1857     0,                  /* nb_inplace_floor_divide */
   1858     0,                  /* nb_inplace_true_divide */
   1859 };
   1860 
   1861 PyTypeObject PyFloat_Type = {
   1862     PyVarObject_HEAD_INIT(&PyType_Type, 0)
   1863     "float",
   1864     sizeof(PyFloatObject),
   1865     0,
   1866     (destructor)float_dealloc,                  /* tp_dealloc */
   1867     0,                                          /* tp_print */
   1868     0,                                          /* tp_getattr */
   1869     0,                                          /* tp_setattr */
   1870     0,                                          /* tp_reserved */
   1871     (reprfunc)float_repr,                       /* tp_repr */
   1872     &float_as_number,                           /* tp_as_number */
   1873     0,                                          /* tp_as_sequence */
   1874     0,                                          /* tp_as_mapping */
   1875     (hashfunc)float_hash,                       /* tp_hash */
   1876     0,                                          /* tp_call */
   1877     (reprfunc)float_repr,                       /* tp_str */
   1878     PyObject_GenericGetAttr,                    /* tp_getattro */
   1879     0,                                          /* tp_setattro */
   1880     0,                                          /* tp_as_buffer */
   1881     Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE, /* tp_flags */
   1882     float_doc,                                  /* tp_doc */
   1883     0,                                          /* tp_traverse */
   1884     0,                                          /* tp_clear */
   1885     float_richcompare,                          /* tp_richcompare */
   1886     0,                                          /* tp_weaklistoffset */
   1887     0,                                          /* tp_iter */
   1888     0,                                          /* tp_iternext */
   1889     float_methods,                              /* tp_methods */
   1890     0,                                          /* tp_members */
   1891     float_getset,                               /* tp_getset */
   1892     0,                                          /* tp_base */
   1893     0,                                          /* tp_dict */
   1894     0,                                          /* tp_descr_get */
   1895     0,                                          /* tp_descr_set */
   1896     0,                                          /* tp_dictoffset */
   1897     0,                                          /* tp_init */
   1898     0,                                          /* tp_alloc */
   1899     float_new,                                  /* tp_new */
   1900 };
   1901 
   1902 int
   1903 _PyFloat_Init(void)
   1904 {
   1905     /* We attempt to determine if this machine is using IEEE
   1906        floating point formats by peering at the bits of some
   1907        carefully chosen values.  If it looks like we are on an
   1908        IEEE platform, the float packing/unpacking routines can
   1909        just copy bits, if not they resort to arithmetic & shifts
   1910        and masks.  The shifts & masks approach works on all finite
   1911        values, but what happens to infinities, NaNs and signed
   1912        zeroes on packing is an accident, and attempting to unpack
   1913        a NaN or an infinity will raise an exception.
   1914 
   1915        Note that if we're on some whacked-out platform which uses
   1916        IEEE formats but isn't strictly little-endian or big-
   1917        endian, we will fall back to the portable shifts & masks
   1918        method. */
   1919 
   1920 #if SIZEOF_DOUBLE == 8
   1921     {
   1922         double x = 9006104071832581.0;
   1923         if (memcmp(&x, "\x43\x3f\xff\x01\x02\x03\x04\x05", 8) == 0)
   1924             detected_double_format = ieee_big_endian_format;
   1925         else if (memcmp(&x, "\x05\x04\x03\x02\x01\xff\x3f\x43", 8) == 0)
   1926             detected_double_format = ieee_little_endian_format;
   1927         else
   1928             detected_double_format = unknown_format;
   1929     }
   1930 #else
   1931     detected_double_format = unknown_format;
   1932 #endif
   1933 
   1934 #if SIZEOF_FLOAT == 4
   1935     {
   1936         float y = 16711938.0;
   1937         if (memcmp(&y, "\x4b\x7f\x01\x02", 4) == 0)
   1938             detected_float_format = ieee_big_endian_format;
   1939         else if (memcmp(&y, "\x02\x01\x7f\x4b", 4) == 0)
   1940             detected_float_format = ieee_little_endian_format;
   1941         else
   1942             detected_float_format = unknown_format;
   1943     }
   1944 #else
   1945     detected_float_format = unknown_format;
   1946 #endif
   1947 
   1948     double_format = detected_double_format;
   1949     float_format = detected_float_format;
   1950 
   1951     /* Init float info */
   1952     if (FloatInfoType.tp_name == NULL) {
   1953         if (PyStructSequence_InitType2(&FloatInfoType, &floatinfo_desc) < 0)
   1954             return 0;
   1955     }
   1956     return 1;
   1957 }
   1958 
   1959 int
   1960 PyFloat_ClearFreeList(void)
   1961 {
   1962     PyFloatObject *f = free_list, *next;
   1963     int i = numfree;
   1964     while (f) {
   1965         next = (PyFloatObject*) Py_TYPE(f);
   1966         PyObject_FREE(f);
   1967         f = next;
   1968     }
   1969     free_list = NULL;
   1970     numfree = 0;
   1971     return i;
   1972 }
   1973 
   1974 void
   1975 PyFloat_Fini(void)
   1976 {
   1977     (void)PyFloat_ClearFreeList();
   1978 }
   1979 
   1980 /* Print summary info about the state of the optimized allocator */
   1981 void
   1982 _PyFloat_DebugMallocStats(FILE *out)
   1983 {
   1984     _PyDebugAllocatorStats(out,
   1985                            "free PyFloatObject",
   1986                            numfree, sizeof(PyFloatObject));
   1987 }
   1988 
   1989 
   1990 /*----------------------------------------------------------------------------
   1991  * _PyFloat_{Pack,Unpack}{2,4,8}.  See floatobject.h.
   1992  * To match the NPY_HALF_ROUND_TIES_TO_EVEN behavior in:
   1993  * https://github.com/numpy/numpy/blob/master/numpy/core/src/npymath/halffloat.c
   1994  * We use:
   1995  *       bits = (unsigned short)f;    Note the truncation
   1996  *       if ((f - bits > 0.5) || (f - bits == 0.5 && bits % 2)) {
   1997  *           bits++;
   1998  *       }
   1999  */
   2000 
   2001 int
   2002 _PyFloat_Pack2(double x, unsigned char *p, int le)
   2003 {
   2004     unsigned char sign;
   2005     int e;
   2006     double f;
   2007     unsigned short bits;
   2008     int incr = 1;
   2009 
   2010     if (x == 0.0) {
   2011         sign = (copysign(1.0, x) == -1.0);
   2012         e = 0;
   2013         bits = 0;
   2014     }
   2015     else if (Py_IS_INFINITY(x)) {
   2016         sign = (x < 0.0);
   2017         e = 0x1f;
   2018         bits = 0;
   2019     }
   2020     else if (Py_IS_NAN(x)) {
   2021         /* There are 2046 distinct half-precision NaNs (1022 signaling and
   2022            1024 quiet), but there are only two quiet NaNs that don't arise by
   2023            quieting a signaling NaN; we get those by setting the topmost bit
   2024            of the fraction field and clearing all other fraction bits. We
   2025            choose the one with the appropriate sign. */
   2026         sign = (copysign(1.0, x) == -1.0);
   2027         e = 0x1f;
   2028         bits = 512;
   2029     }
   2030     else {
   2031         sign = (x < 0.0);
   2032         if (sign) {
   2033             x = -x;
   2034         }
   2035 
   2036         f = frexp(x, &e);
   2037         if (f < 0.5 || f >= 1.0) {
   2038             PyErr_SetString(PyExc_SystemError,
   2039                             "frexp() result out of range");
   2040             return -1;
   2041         }
   2042 
   2043         /* Normalize f to be in the range [1.0, 2.0) */
   2044         f *= 2.0;
   2045         e--;
   2046 
   2047         if (e >= 16) {
   2048             goto Overflow;
   2049         }
   2050         else if (e < -25) {
   2051             /* |x| < 2**-25. Underflow to zero. */
   2052             f = 0.0;
   2053             e = 0;
   2054         }
   2055         else if (e < -14) {
   2056             /* |x| < 2**-14. Gradual underflow */
   2057             f = ldexp(f, 14 + e);
   2058             e = 0;
   2059         }
   2060         else /* if (!(e == 0 && f == 0.0)) */ {
   2061             e += 15;
   2062             f -= 1.0; /* Get rid of leading 1 */
   2063         }
   2064 
   2065         f *= 1024.0; /* 2**10 */
   2066         /* Round to even */
   2067         bits = (unsigned short)f; /* Note the truncation */
   2068         assert(bits < 1024);
   2069         assert(e < 31);
   2070         if ((f - bits > 0.5) || ((f - bits == 0.5) && (bits % 2 == 1))) {
   2071             ++bits;
   2072             if (bits == 1024) {
   2073                 /* The carry propagated out of a string of 10 1 bits. */
   2074                 bits = 0;
   2075                 ++e;
   2076                 if (e == 31)
   2077                     goto Overflow;
   2078             }
   2079         }
   2080     }
   2081 
   2082     bits |= (e << 10) | (sign << 15);
   2083 
   2084     /* Write out result. */
   2085     if (le) {
   2086         p += 1;
   2087         incr = -1;
   2088     }
   2089 
   2090     /* First byte */
   2091     *p = (unsigned char)((bits >> 8) & 0xFF);
   2092     p += incr;
   2093 
   2094     /* Second byte */
   2095     *p = (unsigned char)(bits & 0xFF);
   2096 
   2097     return 0;
   2098 
   2099   Overflow:
   2100     PyErr_SetString(PyExc_OverflowError,
   2101                     "float too large to pack with e format");
   2102     return -1;
   2103 }
   2104 
   2105 int
   2106 _PyFloat_Pack4(double x, unsigned char *p, int le)
   2107 {
   2108     if (float_format == unknown_format) {
   2109         unsigned char sign;
   2110         int e;
   2111         double f;
   2112         unsigned int fbits;
   2113         int incr = 1;
   2114 
   2115         if (le) {
   2116             p += 3;
   2117             incr = -1;
   2118         }
   2119 
   2120         if (x < 0) {
   2121             sign = 1;
   2122             x = -x;
   2123         }
   2124         else
   2125             sign = 0;
   2126 
   2127         f = frexp(x, &e);
   2128 
   2129         /* Normalize f to be in the range [1.0, 2.0) */
   2130         if (0.5 <= f && f < 1.0) {
   2131             f *= 2.0;
   2132             e--;
   2133         }
   2134         else if (f == 0.0)
   2135             e = 0;
   2136         else {
   2137             PyErr_SetString(PyExc_SystemError,
   2138                             "frexp() result out of range");
   2139             return -1;
   2140         }
   2141 
   2142         if (e >= 128)
   2143             goto Overflow;
   2144         else if (e < -126) {
   2145             /* Gradual underflow */
   2146             f = ldexp(f, 126 + e);
   2147             e = 0;
   2148         }
   2149         else if (!(e == 0 && f == 0.0)) {
   2150             e += 127;
   2151             f -= 1.0; /* Get rid of leading 1 */
   2152         }
   2153 
   2154         f *= 8388608.0; /* 2**23 */
   2155         fbits = (unsigned int)(f + 0.5); /* Round */
   2156         assert(fbits <= 8388608);
   2157         if (fbits >> 23) {
   2158             /* The carry propagated out of a string of 23 1 bits. */
   2159             fbits = 0;
   2160             ++e;
   2161             if (e >= 255)
   2162                 goto Overflow;
   2163         }
   2164 
   2165         /* First byte */
   2166         *p = (sign << 7) | (e >> 1);
   2167         p += incr;
   2168 
   2169         /* Second byte */
   2170         *p = (char) (((e & 1) << 7) | (fbits >> 16));
   2171         p += incr;
   2172 
   2173         /* Third byte */
   2174         *p = (fbits >> 8) & 0xFF;
   2175         p += incr;
   2176 
   2177         /* Fourth byte */
   2178         *p = fbits & 0xFF;
   2179 
   2180         /* Done */
   2181         return 0;
   2182 
   2183     }
   2184     else {
   2185         float y = (float)x;
   2186         const unsigned char *s = (unsigned char*)&y;
   2187         int i, incr = 1;
   2188 
   2189         if (Py_IS_INFINITY(y) && !Py_IS_INFINITY(x))
   2190             goto Overflow;
   2191 
   2192         if ((float_format == ieee_little_endian_format && !le)
   2193             || (float_format == ieee_big_endian_format && le)) {
   2194             p += 3;
   2195             incr = -1;
   2196         }
   2197 
   2198         for (i = 0; i < 4; i++) {
   2199             *p = *s++;
   2200             p += incr;
   2201         }
   2202         return 0;
   2203     }
   2204   Overflow:
   2205     PyErr_SetString(PyExc_OverflowError,
   2206                     "float too large to pack with f format");
   2207     return -1;
   2208 }
   2209 
   2210 int
   2211 _PyFloat_Pack8(double x, unsigned char *p, int le)
   2212 {
   2213     if (double_format == unknown_format) {
   2214         unsigned char sign;
   2215         int e;
   2216         double f;
   2217         unsigned int fhi, flo;
   2218         int incr = 1;
   2219 
   2220         if (le) {
   2221             p += 7;
   2222             incr = -1;
   2223         }
   2224 
   2225         if (x < 0) {
   2226             sign = 1;
   2227             x = -x;
   2228         }
   2229         else
   2230             sign = 0;
   2231 
   2232         f = frexp(x, &e);
   2233 
   2234         /* Normalize f to be in the range [1.0, 2.0) */
   2235         if (0.5 <= f && f < 1.0) {
   2236             f *= 2.0;
   2237             e--;
   2238         }
   2239         else if (f == 0.0)
   2240             e = 0;
   2241         else {
   2242             PyErr_SetString(PyExc_SystemError,
   2243                             "frexp() result out of range");
   2244             return -1;
   2245         }
   2246 
   2247         if (e >= 1024)
   2248             goto Overflow;
   2249         else if (e < -1022) {
   2250             /* Gradual underflow */
   2251             f = ldexp(f, 1022 + e);
   2252             e = 0;
   2253         }
   2254         else if (!(e == 0 && f == 0.0)) {
   2255             e += 1023;
   2256             f -= 1.0; /* Get rid of leading 1 */
   2257         }
   2258 
   2259         /* fhi receives the high 28 bits; flo the low 24 bits (== 52 bits) */
   2260         f *= 268435456.0; /* 2**28 */
   2261         fhi = (unsigned int)f; /* Truncate */
   2262         assert(fhi < 268435456);
   2263 
   2264         f -= (double)fhi;
   2265         f *= 16777216.0; /* 2**24 */
   2266         flo = (unsigned int)(f + 0.5); /* Round */
   2267         assert(flo <= 16777216);
   2268         if (flo >> 24) {
   2269             /* The carry propagated out of a string of 24 1 bits. */
   2270             flo = 0;
   2271             ++fhi;
   2272             if (fhi >> 28) {
   2273                 /* And it also progagated out of the next 28 bits. */
   2274                 fhi = 0;
   2275                 ++e;
   2276                 if (e >= 2047)
   2277                     goto Overflow;
   2278             }
   2279         }
   2280 
   2281         /* First byte */
   2282         *p = (sign << 7) | (e >> 4);
   2283         p += incr;
   2284 
   2285         /* Second byte */
   2286         *p = (unsigned char) (((e & 0xF) << 4) | (fhi >> 24));
   2287         p += incr;
   2288 
   2289         /* Third byte */
   2290         *p = (fhi >> 16) & 0xFF;
   2291         p += incr;
   2292 
   2293         /* Fourth byte */
   2294         *p = (fhi >> 8) & 0xFF;
   2295         p += incr;
   2296 
   2297         /* Fifth byte */
   2298         *p = fhi & 0xFF;
   2299         p += incr;
   2300 
   2301         /* Sixth byte */
   2302         *p = (flo >> 16) & 0xFF;
   2303         p += incr;
   2304 
   2305         /* Seventh byte */
   2306         *p = (flo >> 8) & 0xFF;
   2307         p += incr;
   2308 
   2309         /* Eighth byte */
   2310         *p = flo & 0xFF;
   2311         /* p += incr; */
   2312 
   2313         /* Done */
   2314         return 0;
   2315 
   2316       Overflow:
   2317         PyErr_SetString(PyExc_OverflowError,
   2318                         "float too large to pack with d format");
   2319         return -1;
   2320     }
   2321     else {
   2322         const unsigned char *s = (unsigned char*)&x;
   2323         int i, incr = 1;
   2324 
   2325         if ((double_format == ieee_little_endian_format && !le)
   2326             || (double_format == ieee_big_endian_format && le)) {
   2327             p += 7;
   2328             incr = -1;
   2329         }
   2330 
   2331         for (i = 0; i < 8; i++) {
   2332             *p = *s++;
   2333             p += incr;
   2334         }
   2335         return 0;
   2336     }
   2337 }
   2338 
   2339 double
   2340 _PyFloat_Unpack2(const unsigned char *p, int le)
   2341 {
   2342     unsigned char sign;
   2343     int e;
   2344     unsigned int f;
   2345     double x;
   2346     int incr = 1;
   2347 
   2348     if (le) {
   2349         p += 1;
   2350         incr = -1;
   2351     }
   2352 
   2353     /* First byte */
   2354     sign = (*p >> 7) & 1;
   2355     e = (*p & 0x7C) >> 2;
   2356     f = (*p & 0x03) << 8;
   2357     p += incr;
   2358 
   2359     /* Second byte */
   2360     f |= *p;
   2361 
   2362     if (e == 0x1f) {
   2363 #ifdef PY_NO_SHORT_FLOAT_REPR
   2364         if (f == 0) {
   2365             /* Infinity */
   2366             return sign ? -Py_HUGE_VAL : Py_HUGE_VAL;
   2367         }
   2368         else {
   2369             /* NaN */
   2370 #ifdef Py_NAN
   2371             return sign ? -Py_NAN : Py_NAN;
   2372 #else
   2373             PyErr_SetString(
   2374                 PyExc_ValueError,
   2375                 "can't unpack IEEE 754 NaN "
   2376                 "on platform that does not support NaNs");
   2377             return -1;
   2378 #endif  /* #ifdef Py_NAN */
   2379         }
   2380 #else
   2381         if (f == 0) {
   2382             /* Infinity */
   2383             return _Py_dg_infinity(sign);
   2384         }
   2385         else {
   2386             /* NaN */
   2387             return _Py_dg_stdnan(sign);
   2388         }
   2389 #endif  /* #ifdef PY_NO_SHORT_FLOAT_REPR */
   2390     }
   2391 
   2392     x = (double)f / 1024.0;
   2393 
   2394     if (e == 0) {
   2395         e = -14;
   2396     }
   2397     else {
   2398         x += 1.0;
   2399         e -= 15;
   2400     }
   2401     x = ldexp(x, e);
   2402 
   2403     if (sign)
   2404         x = -x;
   2405 
   2406     return x;
   2407 }
   2408 
   2409 double
   2410 _PyFloat_Unpack4(const unsigned char *p, int le)
   2411 {
   2412     if (float_format == unknown_format) {
   2413         unsigned char sign;
   2414         int e;
   2415         unsigned int f;
   2416         double x;
   2417         int incr = 1;
   2418 
   2419         if (le) {
   2420             p += 3;
   2421             incr = -1;
   2422         }
   2423 
   2424         /* First byte */
   2425         sign = (*p >> 7) & 1;
   2426         e = (*p & 0x7F) << 1;
   2427         p += incr;
   2428 
   2429         /* Second byte */
   2430         e |= (*p >> 7) & 1;
   2431         f = (*p & 0x7F) << 16;
   2432         p += incr;
   2433 
   2434         if (e == 255) {
   2435             PyErr_SetString(
   2436                 PyExc_ValueError,
   2437                 "can't unpack IEEE 754 special value "
   2438                 "on non-IEEE platform");
   2439             return -1;
   2440         }
   2441 
   2442         /* Third byte */
   2443         f |= *p << 8;
   2444         p += incr;
   2445 
   2446         /* Fourth byte */
   2447         f |= *p;
   2448 
   2449         x = (double)f / 8388608.0;
   2450 
   2451         /* XXX This sadly ignores Inf/NaN issues */
   2452         if (e == 0)
   2453             e = -126;
   2454         else {
   2455             x += 1.0;
   2456             e -= 127;
   2457         }
   2458         x = ldexp(x, e);
   2459 
   2460         if (sign)
   2461             x = -x;
   2462 
   2463         return x;
   2464     }
   2465     else {
   2466         float x;
   2467 
   2468         if ((float_format == ieee_little_endian_format && !le)
   2469             || (float_format == ieee_big_endian_format && le)) {
   2470             char buf[4];
   2471             char *d = &buf[3];
   2472             int i;
   2473 
   2474             for (i = 0; i < 4; i++) {
   2475                 *d-- = *p++;
   2476             }
   2477             memcpy(&x, buf, 4);
   2478         }
   2479         else {
   2480             memcpy(&x, p, 4);
   2481         }
   2482 
   2483         return x;
   2484     }
   2485 }
   2486 
   2487 double
   2488 _PyFloat_Unpack8(const unsigned char *p, int le)
   2489 {
   2490     if (double_format == unknown_format) {
   2491         unsigned char sign;
   2492         int e;
   2493         unsigned int fhi, flo;
   2494         double x;
   2495         int incr = 1;
   2496 
   2497         if (le) {
   2498             p += 7;
   2499             incr = -1;
   2500         }
   2501 
   2502         /* First byte */
   2503         sign = (*p >> 7) & 1;
   2504         e = (*p & 0x7F) << 4;
   2505 
   2506         p += incr;
   2507 
   2508         /* Second byte */
   2509         e |= (*p >> 4) & 0xF;
   2510         fhi = (*p & 0xF) << 24;
   2511         p += incr;
   2512 
   2513         if (e == 2047) {
   2514             PyErr_SetString(
   2515                 PyExc_ValueError,
   2516                 "can't unpack IEEE 754 special value "
   2517                 "on non-IEEE platform");
   2518             return -1.0;
   2519         }
   2520 
   2521         /* Third byte */
   2522         fhi |= *p << 16;
   2523         p += incr;
   2524 
   2525         /* Fourth byte */
   2526         fhi |= *p  << 8;
   2527         p += incr;
   2528 
   2529         /* Fifth byte */
   2530         fhi |= *p;
   2531         p += incr;
   2532 
   2533         /* Sixth byte */
   2534         flo = *p << 16;
   2535         p += incr;
   2536 
   2537         /* Seventh byte */
   2538         flo |= *p << 8;
   2539         p += incr;
   2540 
   2541         /* Eighth byte */
   2542         flo |= *p;
   2543 
   2544         x = (double)fhi + (double)flo / 16777216.0; /* 2**24 */
   2545         x /= 268435456.0; /* 2**28 */
   2546 
   2547         if (e == 0)
   2548             e = -1022;
   2549         else {
   2550             x += 1.0;
   2551             e -= 1023;
   2552         }
   2553         x = ldexp(x, e);
   2554 
   2555         if (sign)
   2556             x = -x;
   2557 
   2558         return x;
   2559     }
   2560     else {
   2561         double x;
   2562 
   2563         if ((double_format == ieee_little_endian_format && !le)
   2564             || (double_format == ieee_big_endian_format && le)) {
   2565             char buf[8];
   2566             char *d = &buf[7];
   2567             int i;
   2568 
   2569             for (i = 0; i < 8; i++) {
   2570                 *d-- = *p++;
   2571             }
   2572             memcpy(&x, buf, 8);
   2573         }
   2574         else {
   2575             memcpy(&x, p, 8);
   2576         }
   2577 
   2578         return x;
   2579     }
   2580 }
   2581