Objects/floatobject.c

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