source: trunk/essentials/dev-lang/python/Objects/floatobject.c@ 3506

Last change on this file since 3506 was 3225, checked in by bird, 19 years ago

Python 2.5
File size: 40.4 KB
Line 
1
2/* Float object implementation */
3
4/* XXX There should be overflow checks here, but it's hard to check
5 for any kind of float exception without losing portability. */
6
7#include "Python.h"
8
9#include <ctype.h>
10
11#if !defined(__STDC__)
12extern double fmod(double, double);
13extern double pow(double, double);
14#endif
15
16/* Special free list -- see comments for same code in intobject.c. */
17#define BLOCK_SIZE 1000 /* 1K less typical malloc overhead */
18#define BHEAD_SIZE 8 /* Enough for a 64-bit pointer */
19#define N_FLOATOBJECTS ((BLOCK_SIZE - BHEAD_SIZE) / sizeof(PyFloatObject))
20
21struct _floatblock {
22 struct _floatblock *next;
23 PyFloatObject objects[N_FLOATOBJECTS];
24};
25
26typedef struct _floatblock PyFloatBlock;
27
28static PyFloatBlock *block_list = NULL;
29static PyFloatObject *free_list = NULL;
30
31static PyFloatObject *
32fill_free_list(void)
33{
34 PyFloatObject *p, *q;
35 /* XXX Float blocks escape the object heap. Use PyObject_MALLOC ??? */
36 p = (PyFloatObject *) PyMem_MALLOC(sizeof(PyFloatBlock));
37 if (p == NULL)
38 return (PyFloatObject *) PyErr_NoMemory();
39 ((PyFloatBlock *)p)->next = block_list;
40 block_list = (PyFloatBlock *)p;
41 p = &((PyFloatBlock *)p)->objects[0];
42 q = p + N_FLOATOBJECTS;
43 while (--q > p)
44 q->ob_type = (struct _typeobject *)(q-1);
45 q->ob_type = NULL;
46 return p + N_FLOATOBJECTS - 1;
47}
48
49PyObject *
50PyFloat_FromDouble(double fval)
51{
52 register PyFloatObject *op;
53 if (free_list == NULL) {
54 if ((free_list = fill_free_list()) == NULL)
55 return NULL;
56 }
57 /* Inline PyObject_New */
58 op = free_list;
59 free_list = (PyFloatObject *)op->ob_type;
60 PyObject_INIT(op, &PyFloat_Type);
61 op->ob_fval = fval;
62 return (PyObject *) op;
63}
64
65/**************************************************************************
66RED_FLAG 22-Sep-2000 tim
67PyFloat_FromString's pend argument is braindead. Prior to this RED_FLAG,
68
691. If v was a regular string, *pend was set to point to its terminating
70 null byte. That's useless (the caller can find that without any
71 help from this function!).
72
732. If v was a Unicode string, or an object convertible to a character
74 buffer, *pend was set to point into stack trash (the auto temp
75 vector holding the character buffer). That was downright dangerous.
76
77Since we can't change the interface of a public API function, pend is
78still supported but now *officially* useless: if pend is not NULL,
79*pend is set to NULL.
80**************************************************************************/
81PyObject *
82PyFloat_FromString(PyObject *v, char **pend)
83{
84 const char *s, *last, *end;
85 double x;
86 char buffer[256]; /* for errors */
87#ifdef Py_USING_UNICODE
88 char s_buffer[256]; /* for objects convertible to a char buffer */
89#endif
90 Py_ssize_t len;
91
92 if (pend)
93 *pend = NULL;
94 if (PyString_Check(v)) {
95 s = PyString_AS_STRING(v);
96 len = PyString_GET_SIZE(v);
97 }
98#ifdef Py_USING_UNICODE
99 else if (PyUnicode_Check(v)) {
100 if (PyUnicode_GET_SIZE(v) >= (Py_ssize_t)sizeof(s_buffer)) {
101 PyErr_SetString(PyExc_ValueError,
102 "Unicode float() literal too long to convert");
103 return NULL;
104 }
105 if (PyUnicode_EncodeDecimal(PyUnicode_AS_UNICODE(v),
106 PyUnicode_GET_SIZE(v),
107 s_buffer,
108 NULL))
109 return NULL;
110 s = s_buffer;
111 len = strlen(s);
112 }
113#endif
114 else if (PyObject_AsCharBuffer(v, &s, &len)) {
115 PyErr_SetString(PyExc_TypeError,
116 "float() argument must be a string or a number");
117 return NULL;
118 }
119
120 last = s + len;
121 while (*s && isspace(Py_CHARMASK(*s)))
122 s++;
123 if (*s == '\0') {
124 PyErr_SetString(PyExc_ValueError, "empty string for float()");
125 return NULL;
126 }
127 /* We don't care about overflow or underflow. If the platform supports
128 * them, infinities and signed zeroes (on underflow) are fine.
129 * However, strtod can return 0 for denormalized numbers, where atof
130 * does not. So (alas!) we special-case a zero result. Note that
131 * whether strtod sets errno on underflow is not defined, so we can't
132 * key off errno.
133 */
134 PyFPE_START_PROTECT("strtod", return NULL)
135 x = PyOS_ascii_strtod(s, (char **)&end);
136 PyFPE_END_PROTECT(x)
137 errno = 0;
138 /* Believe it or not, Solaris 2.6 can move end *beyond* the null
139 byte at the end of the string, when the input is inf(inity). */
140 if (end > last)
141 end = last;
142 if (end == s) {
143 PyOS_snprintf(buffer, sizeof(buffer),
144 "invalid literal for float(): %.200s", s);
145 PyErr_SetString(PyExc_ValueError, buffer);
146 return NULL;
147 }
148 /* Since end != s, the platform made *some* kind of sense out
149 of the input. Trust it. */
150 while (*end && isspace(Py_CHARMASK(*end)))
151 end++;
152 if (*end != '\0') {
153 PyOS_snprintf(buffer, sizeof(buffer),
154 "invalid literal for float(): %.200s", s);
155 PyErr_SetString(PyExc_ValueError, buffer);
156 return NULL;
157 }
158 else if (end != last) {
159 PyErr_SetString(PyExc_ValueError,
160 "null byte in argument for float()");
161 return NULL;
162 }
163 if (x == 0.0) {
164 /* See above -- may have been strtod being anal
165 about denorms. */
166 PyFPE_START_PROTECT("atof", return NULL)
167 x = PyOS_ascii_atof(s);
168 PyFPE_END_PROTECT(x)
169 errno = 0; /* whether atof ever set errno is undefined */
170 }
171 return PyFloat_FromDouble(x);
172}
173
174static void
175float_dealloc(PyFloatObject *op)
176{
177 if (PyFloat_CheckExact(op)) {
178 op->ob_type = (struct _typeobject *)free_list;
179 free_list = op;
180 }
181 else
182 op->ob_type->tp_free((PyObject *)op);
183}
184
185double
186PyFloat_AsDouble(PyObject *op)
187{
188 PyNumberMethods *nb;
189 PyFloatObject *fo;
190 double val;
191
192 if (op && PyFloat_Check(op))
193 return PyFloat_AS_DOUBLE((PyFloatObject*) op);
194
195 if (op == NULL) {
196 PyErr_BadArgument();
197 return -1;
198 }
199
200 if ((nb = op->ob_type->tp_as_number) == NULL || nb->nb_float == NULL) {
201 PyErr_SetString(PyExc_TypeError, "a float is required");
202 return -1;
203 }
204
205 fo = (PyFloatObject*) (*nb->nb_float) (op);
206 if (fo == NULL)
207 return -1;
208 if (!PyFloat_Check(fo)) {
209 PyErr_SetString(PyExc_TypeError,
210 "nb_float should return float object");
211 return -1;
212 }
213
214 val = PyFloat_AS_DOUBLE(fo);
215 Py_DECREF(fo);
216
217 return val;
218}
219
220/* Methods */
221
222static void
223format_float(char *buf, size_t buflen, PyFloatObject *v, int precision)
224{
225 register char *cp;
226 char format[32];
227 /* Subroutine for float_repr and float_print.
228 We want float numbers to be recognizable as such,
229 i.e., they should contain a decimal point or an exponent.
230 However, %g may print the number as an integer;
231 in such cases, we append ".0" to the string. */
232
233 assert(PyFloat_Check(v));
234 PyOS_snprintf(format, 32, "%%.%ig", precision);
235 PyOS_ascii_formatd(buf, buflen, format, v->ob_fval);
236 cp = buf;
237 if (*cp == '-')
238 cp++;
239 for (; *cp != '\0'; cp++) {
240 /* Any non-digit means it's not an integer;
241 this takes care of NAN and INF as well. */
242 if (!isdigit(Py_CHARMASK(*cp)))
243 break;
244 }
245 if (*cp == '\0') {
246 *cp++ = '.';
247 *cp++ = '0';
248 *cp++ = '\0';
249 }
250}
251
252/* XXX PyFloat_AsStringEx should not be a public API function (for one
253 XXX thing, its signature passes a buffer without a length; for another,
254 XXX it isn't useful outside this file).
255*/
256void
257PyFloat_AsStringEx(char *buf, PyFloatObject *v, int precision)
258{
259 format_float(buf, 100, v, precision);
260}
261
262/* Macro and helper that convert PyObject obj to a C double and store
263 the value in dbl; this replaces the functionality of the coercion
264 slot function. If conversion to double raises an exception, obj is
265 set to NULL, and the function invoking this macro returns NULL. If
266 obj is not of float, int or long type, Py_NotImplemented is incref'ed,
267 stored in obj, and returned from the function invoking this macro.
268*/
269#define CONVERT_TO_DOUBLE(obj, dbl) \
270 if (PyFloat_Check(obj)) \
271 dbl = PyFloat_AS_DOUBLE(obj); \
272 else if (convert_to_double(&(obj), &(dbl)) < 0) \
273 return obj;
274
275static int
276convert_to_double(PyObject **v, double *dbl)
277{
278 register PyObject *obj = *v;
279
280 if (PyInt_Check(obj)) {
281 *dbl = (double)PyInt_AS_LONG(obj);
282 }
283 else if (PyLong_Check(obj)) {
284 *dbl = PyLong_AsDouble(obj);
285 if (*dbl == -1.0 && PyErr_Occurred()) {
286 *v = NULL;
287 return -1;
288 }
289 }
290 else {
291 Py_INCREF(Py_NotImplemented);
292 *v = Py_NotImplemented;
293 return -1;
294 }
295 return 0;
296}
297
298/* Precisions used by repr() and str(), respectively.
299
300 The repr() precision (17 significant decimal digits) is the minimal number
301 that is guaranteed to have enough precision so that if the number is read
302 back in the exact same binary value is recreated. This is true for IEEE
303 floating point by design, and also happens to work for all other modern
304 hardware.
305
306 The str() precision is chosen so that in most cases, the rounding noise
307 created by various operations is suppressed, while giving plenty of
308 precision for practical use.
309
310*/
311
312#define PREC_REPR 17
313#define PREC_STR 12
314
315/* XXX PyFloat_AsString and PyFloat_AsReprString should be deprecated:
316 XXX they pass a char buffer without passing a length.
317*/
318void
319PyFloat_AsString(char *buf, PyFloatObject *v)
320{
321 format_float(buf, 100, v, PREC_STR);
322}
323
324void
325PyFloat_AsReprString(char *buf, PyFloatObject *v)
326{
327 format_float(buf, 100, v, PREC_REPR);
328}
329
330/* ARGSUSED */
331static int
332float_print(PyFloatObject *v, FILE *fp, int flags)
333{
334 char buf[100];
335 format_float(buf, sizeof(buf), v,
336 (flags & Py_PRINT_RAW) ? PREC_STR : PREC_REPR);
337 fputs(buf, fp);
338 return 0;
339}
340
341static PyObject *
342float_repr(PyFloatObject *v)
343{
344 char buf[100];
345 format_float(buf, sizeof(buf), v, PREC_REPR);
346 return PyString_FromString(buf);
347}
348
349static PyObject *
350float_str(PyFloatObject *v)
351{
352 char buf[100];
353 format_float(buf, sizeof(buf), v, PREC_STR);
354 return PyString_FromString(buf);
355}
356
357/* Comparison is pretty much a nightmare. When comparing float to float,
358 * we do it as straightforwardly (and long-windedly) as conceivable, so
359 * that, e.g., Python x == y delivers the same result as the platform
360 * C x == y when x and/or y is a NaN.
361 * When mixing float with an integer type, there's no good *uniform* approach.
362 * Converting the double to an integer obviously doesn't work, since we
363 * may lose info from fractional bits. Converting the integer to a double
364 * also has two failure modes: (1) a long int may trigger overflow (too
365 * large to fit in the dynamic range of a C double); (2) even a C long may have
366 * more bits than fit in a C double (e.g., on a a 64-bit box long may have
367 * 63 bits of precision, but a C double probably has only 53), and then
368 * we can falsely claim equality when low-order integer bits are lost by
369 * coercion to double. So this part is painful too.
370 */
371
372static PyObject*
373float_richcompare(PyObject *v, PyObject *w, int op)
374{
375 double i, j;
376 int r = 0;
377
378 assert(PyFloat_Check(v));
379 i = PyFloat_AS_DOUBLE(v);
380
381 /* Switch on the type of w. Set i and j to doubles to be compared,
382 * and op to the richcomp to use.
383 */
384 if (PyFloat_Check(w))
385 j = PyFloat_AS_DOUBLE(w);
386
387 else if (!Py_IS_FINITE(i)) {
388 if (PyInt_Check(w) || PyLong_Check(w))
389 /* If i is an infinity, its magnitude exceeds any
390 * finite integer, so it doesn't matter which int we
391 * compare i with. If i is a NaN, similarly.
392 */
393 j = 0.0;
394 else
395 goto Unimplemented;
396 }
397
398 else if (PyInt_Check(w)) {
399 long jj = PyInt_AS_LONG(w);
400 /* In the worst realistic case I can imagine, C double is a
401 * Cray single with 48 bits of precision, and long has 64
402 * bits.
403 */
404#if SIZEOF_LONG > 6
405 unsigned long abs = (unsigned long)(jj < 0 ? -jj : jj);
406 if (abs >> 48) {
407 /* Needs more than 48 bits. Make it take the
408 * PyLong path.
409 */
410 PyObject *result;
411 PyObject *ww = PyLong_FromLong(jj);
412
413 if (ww == NULL)
414 return NULL;
415 result = float_richcompare(v, ww, op);
416 Py_DECREF(ww);
417 return result;
418 }
419#endif
420 j = (double)jj;
421 assert((long)j == jj);
422 }
423
424 else if (PyLong_Check(w)) {
425 int vsign = i == 0.0 ? 0 : i < 0.0 ? -1 : 1;
426 int wsign = _PyLong_Sign(w);
427 size_t nbits;
428 int exponent;
429
430 if (vsign != wsign) {
431 /* Magnitudes are irrelevant -- the signs alone
432 * determine the outcome.
433 */
434 i = (double)vsign;
435 j = (double)wsign;
436 goto Compare;
437 }
438 /* The signs are the same. */
439 /* Convert w to a double if it fits. In particular, 0 fits. */
440 nbits = _PyLong_NumBits(w);
441 if (nbits == (size_t)-1 && PyErr_Occurred()) {
442 /* This long is so large that size_t isn't big enough
443 * to hold the # of bits. Replace with little doubles
444 * that give the same outcome -- w is so large that
445 * its magnitude must exceed the magnitude of any
446 * finite float.
447 */
448 PyErr_Clear();
449 i = (double)vsign;
450 assert(wsign != 0);
451 j = wsign * 2.0;
452 goto Compare;
453 }
454 if (nbits <= 48) {
455 j = PyLong_AsDouble(w);
456 /* It's impossible that <= 48 bits overflowed. */
457 assert(j != -1.0 || ! PyErr_Occurred());
458 goto Compare;
459 }
460 assert(wsign != 0); /* else nbits was 0 */
461 assert(vsign != 0); /* if vsign were 0, then since wsign is
462 * not 0, we would have taken the
463 * vsign != wsign branch at the start */
464 /* We want to work with non-negative numbers. */
465 if (vsign < 0) {
466 /* "Multiply both sides" by -1; this also swaps the
467 * comparator.
468 */
469 i = -i;
470 op = _Py_SwappedOp[op];
471 }
472 assert(i > 0.0);
473 (void) frexp(i, &exponent);
474 /* exponent is the # of bits in v before the radix point;
475 * we know that nbits (the # of bits in w) > 48 at this point
476 */
477 if (exponent < 0 || (size_t)exponent < nbits) {
478 i = 1.0;
479 j = 2.0;
480 goto Compare;
481 }
482 if ((size_t)exponent > nbits) {
483 i = 2.0;
484 j = 1.0;
485 goto Compare;
486 }
487 /* v and w have the same number of bits before the radix
488 * point. Construct two longs that have the same comparison
489 * outcome.
490 */
491 {
492 double fracpart;
493 double intpart;
494 PyObject *result = NULL;
495 PyObject *one = NULL;
496 PyObject *vv = NULL;
497 PyObject *ww = w;
498
499 if (wsign < 0) {
500 ww = PyNumber_Negative(w);
501 if (ww == NULL)
502 goto Error;
503 }
504 else
505 Py_INCREF(ww);
506
507 fracpart = modf(i, &intpart);
508 vv = PyLong_FromDouble(intpart);
509 if (vv == NULL)
510 goto Error;
511
512 if (fracpart != 0.0) {
513 /* Shift left, and or a 1 bit into vv
514 * to represent the lost fraction.
515 */
516 PyObject *temp;
517
518 one = PyInt_FromLong(1);
519 if (one == NULL)
520 goto Error;
521
522 temp = PyNumber_Lshift(ww, one);
523 if (temp == NULL)
524 goto Error;
525 Py_DECREF(ww);
526 ww = temp;
527
528 temp = PyNumber_Lshift(vv, one);
529 if (temp == NULL)
530 goto Error;
531 Py_DECREF(vv);
532 vv = temp;
533
534 temp = PyNumber_Or(vv, one);
535 if (temp == NULL)
536 goto Error;
537 Py_DECREF(vv);
538 vv = temp;
539 }
540
541 r = PyObject_RichCompareBool(vv, ww, op);
542 if (r < 0)
543 goto Error;
544 result = PyBool_FromLong(r);
545 Error:
546 Py_XDECREF(vv);
547 Py_XDECREF(ww);
548 Py_XDECREF(one);
549 return result;
550 }
551 } /* else if (PyLong_Check(w)) */
552
553 else /* w isn't float, int, or long */
554 goto Unimplemented;
555
556 Compare:
557 PyFPE_START_PROTECT("richcompare", return NULL)
558 switch (op) {
559 case Py_EQ:
560 r = i == j;
561 break;
562 case Py_NE:
563 r = i != j;
564 break;
565 case Py_LE:
566 r = i <= j;
567 break;
568 case Py_GE:
569 r = i >= j;
570 break;
571 case Py_LT:
572 r = i < j;
573 break;
574 case Py_GT:
575 r = i > j;
576 break;
577 }
578 PyFPE_END_PROTECT(r)
579 return PyBool_FromLong(r);
580
581 Unimplemented:
582 Py_INCREF(Py_NotImplemented);
583 return Py_NotImplemented;
584}
585
586static long
587float_hash(PyFloatObject *v)
588{
589 return _Py_HashDouble(v->ob_fval);
590}
591
592static PyObject *
593float_add(PyObject *v, PyObject *w)
594{
595 double a,b;
596 CONVERT_TO_DOUBLE(v, a);
597 CONVERT_TO_DOUBLE(w, b);
598 PyFPE_START_PROTECT("add", return 0)
599 a = a + b;
600 PyFPE_END_PROTECT(a)
601 return PyFloat_FromDouble(a);
602}
603
604static PyObject *
605float_sub(PyObject *v, PyObject *w)
606{
607 double a,b;
608 CONVERT_TO_DOUBLE(v, a);
609 CONVERT_TO_DOUBLE(w, b);
610 PyFPE_START_PROTECT("subtract", return 0)
611 a = a - b;
612 PyFPE_END_PROTECT(a)
613 return PyFloat_FromDouble(a);
614}
615
616static PyObject *
617float_mul(PyObject *v, PyObject *w)
618{
619 double a,b;
620 CONVERT_TO_DOUBLE(v, a);
621 CONVERT_TO_DOUBLE(w, b);
622 PyFPE_START_PROTECT("multiply", return 0)
623 a = a * b;
624 PyFPE_END_PROTECT(a)
625 return PyFloat_FromDouble(a);
626}
627
628static PyObject *
629float_div(PyObject *v, PyObject *w)
630{
631 double a,b;
632 CONVERT_TO_DOUBLE(v, a);
633 CONVERT_TO_DOUBLE(w, b);
634 if (b == 0.0) {
635 PyErr_SetString(PyExc_ZeroDivisionError, "float division");
636 return NULL;
637 }
638 PyFPE_START_PROTECT("divide", return 0)
639 a = a / b;
640 PyFPE_END_PROTECT(a)
641 return PyFloat_FromDouble(a);
642}
643
644static PyObject *
645float_classic_div(PyObject *v, PyObject *w)
646{
647 double a,b;
648 CONVERT_TO_DOUBLE(v, a);
649 CONVERT_TO_DOUBLE(w, b);
650 if (Py_DivisionWarningFlag >= 2 &&
651 PyErr_Warn(PyExc_DeprecationWarning, "classic float division") < 0)
652 return NULL;
653 if (b == 0.0) {
654 PyErr_SetString(PyExc_ZeroDivisionError, "float division");
655 return NULL;
656 }
657 PyFPE_START_PROTECT("divide", return 0)
658 a = a / b;
659 PyFPE_END_PROTECT(a)
660 return PyFloat_FromDouble(a);
661}
662
663static PyObject *
664float_rem(PyObject *v, PyObject *w)
665{
666 double vx, wx;
667 double mod;
668 CONVERT_TO_DOUBLE(v, vx);
669 CONVERT_TO_DOUBLE(w, wx);
670 if (wx == 0.0) {
671 PyErr_SetString(PyExc_ZeroDivisionError, "float modulo");
672 return NULL;
673 }
674 PyFPE_START_PROTECT("modulo", return 0)
675 mod = fmod(vx, wx);
676 /* note: checking mod*wx < 0 is incorrect -- underflows to
677 0 if wx < sqrt(smallest nonzero double) */
678 if (mod && ((wx < 0) != (mod < 0))) {
679 mod += wx;
680 }
681 PyFPE_END_PROTECT(mod)
682 return PyFloat_FromDouble(mod);
683}
684
685static PyObject *
686float_divmod(PyObject *v, PyObject *w)
687{
688 double vx, wx;
689 double div, mod, floordiv;
690 CONVERT_TO_DOUBLE(v, vx);
691 CONVERT_TO_DOUBLE(w, wx);
692 if (wx == 0.0) {
693 PyErr_SetString(PyExc_ZeroDivisionError, "float divmod()");
694 return NULL;
695 }
696 PyFPE_START_PROTECT("divmod", return 0)
697 mod = fmod(vx, wx);
698 /* fmod is typically exact, so vx-mod is *mathematically* an
699 exact multiple of wx. But this is fp arithmetic, and fp
700 vx - mod is an approximation; the result is that div may
701 not be an exact integral value after the division, although
702 it will always be very close to one.
703 */
704 div = (vx - mod) / wx;
705 if (mod) {
706 /* ensure the remainder has the same sign as the denominator */
707 if ((wx < 0) != (mod < 0)) {
708 mod += wx;
709 div -= 1.0;
710 }
711 }
712 else {
713 /* the remainder is zero, and in the presence of signed zeroes
714 fmod returns different results across platforms; ensure
715 it has the same sign as the denominator; we'd like to do
716 "mod = wx * 0.0", but that may get optimized away */
717 mod *= mod; /* hide "mod = +0" from optimizer */
718 if (wx < 0.0)
719 mod = -mod;
720 }
721 /* snap quotient to nearest integral value */
722 if (div) {
723 floordiv = floor(div);
724 if (div - floordiv > 0.5)
725 floordiv += 1.0;
726 }
727 else {
728 /* div is zero - get the same sign as the true quotient */
729 div *= div; /* hide "div = +0" from optimizers */
730 floordiv = div * vx / wx; /* zero w/ sign of vx/wx */
731 }
732 PyFPE_END_PROTECT(floordiv)
733 return Py_BuildValue("(dd)", floordiv, mod);
734}
735
736static PyObject *
737float_floor_div(PyObject *v, PyObject *w)
738{
739 PyObject *t, *r;
740
741 t = float_divmod(v, w);
742 if (t == NULL || t == Py_NotImplemented)
743 return t;
744 assert(PyTuple_CheckExact(t));
745 r = PyTuple_GET_ITEM(t, 0);
746 Py_INCREF(r);
747 Py_DECREF(t);
748 return r;
749}
750
751static PyObject *
752float_pow(PyObject *v, PyObject *w, PyObject *z)
753{
754 double iv, iw, ix;
755
756 if ((PyObject *)z != Py_None) {
757 PyErr_SetString(PyExc_TypeError, "pow() 3rd argument not "
758 "allowed unless all arguments are integers");
759 return NULL;
760 }
761
762 CONVERT_TO_DOUBLE(v, iv);
763 CONVERT_TO_DOUBLE(w, iw);
764
765 /* Sort out special cases here instead of relying on pow() */
766 if (iw == 0) { /* v**0 is 1, even 0**0 */
767 PyFPE_START_PROTECT("pow", return NULL)
768 if ((PyObject *)z != Py_None) {
769 double iz;
770 CONVERT_TO_DOUBLE(z, iz);
771 ix = fmod(1.0, iz);
772 if (ix != 0 && iz < 0)
773 ix += iz;
774 }
775 else
776 ix = 1.0;
777 PyFPE_END_PROTECT(ix)
778 return PyFloat_FromDouble(ix);
779 }
780 if (iv == 0.0) { /* 0**w is error if w<0, else 1 */
781 if (iw < 0.0) {
782 PyErr_SetString(PyExc_ZeroDivisionError,
783 "0.0 cannot be raised to a negative power");
784 return NULL;
785 }
786 return PyFloat_FromDouble(0.0);
787 }
788 if (iv < 0.0) {
789 /* Whether this is an error is a mess, and bumps into libm
790 * bugs so we have to figure it out ourselves.
791 */
792 if (iw != floor(iw)) {
793 PyErr_SetString(PyExc_ValueError, "negative number "
794 "cannot be raised to a fractional power");
795 return NULL;
796 }
797 /* iw is an exact integer, albeit perhaps a very large one.
798 * -1 raised to an exact integer should never be exceptional.
799 * Alas, some libms (chiefly glibc as of early 2003) return
800 * NaN and set EDOM on pow(-1, large_int) if the int doesn't
801 * happen to be representable in a *C* integer. That's a
802 * bug; we let that slide in math.pow() (which currently
803 * reflects all platform accidents), but not for Python's **.
804 */
805 if (iv == -1.0 && Py_IS_FINITE(iw)) {
806 /* Return 1 if iw is even, -1 if iw is odd; there's
807 * no guarantee that any C integral type is big
808 * enough to hold iw, so we have to check this
809 * indirectly.
810 */
811 ix = floor(iw * 0.5) * 2.0;
812 return PyFloat_FromDouble(ix == iw ? 1.0 : -1.0);
813 }
814 /* Else iv != -1.0, and overflow or underflow are possible.
815 * Unless we're to write pow() ourselves, we have to trust
816 * the platform to do this correctly.
817 */
818 }
819 errno = 0;
820 PyFPE_START_PROTECT("pow", return NULL)
821 ix = pow(iv, iw);
822 PyFPE_END_PROTECT(ix)
823 Py_ADJUST_ERANGE1(ix);
824 if (errno != 0) {
825 /* We don't expect any errno value other than ERANGE, but
826 * the range of libm bugs appears unbounded.
827 */
828 PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError :
829 PyExc_ValueError);
830 return NULL;
831 }
832 return PyFloat_FromDouble(ix);
833}
834
835static PyObject *
836float_neg(PyFloatObject *v)
837{
838 return PyFloat_FromDouble(-v->ob_fval);
839}
840
841static PyObject *
842float_pos(PyFloatObject *v)
843{
844 if (PyFloat_CheckExact(v)) {
845 Py_INCREF(v);
846 return (PyObject *)v;
847 }
848 else
849 return PyFloat_FromDouble(v->ob_fval);
850}
851
852static PyObject *
853float_abs(PyFloatObject *v)
854{
855 return PyFloat_FromDouble(fabs(v->ob_fval));
856}
857
858static int
859float_nonzero(PyFloatObject *v)
860{
861 return v->ob_fval != 0.0;
862}
863
864static int
865float_coerce(PyObject **pv, PyObject **pw)
866{
867 if (PyInt_Check(*pw)) {
868 long x = PyInt_AsLong(*pw);
869 *pw = PyFloat_FromDouble((double)x);
870 Py_INCREF(*pv);
871 return 0;
872 }
873 else if (PyLong_Check(*pw)) {
874 double x = PyLong_AsDouble(*pw);
875 if (x == -1.0 && PyErr_Occurred())
876 return -1;
877 *pw = PyFloat_FromDouble(x);
878 Py_INCREF(*pv);
879 return 0;
880 }
881 else if (PyFloat_Check(*pw)) {
882 Py_INCREF(*pv);
883 Py_INCREF(*pw);
884 return 0;
885 }
886 return 1; /* Can't do it */
887}
888
889static PyObject *
890float_long(PyObject *v)
891{
892 double x = PyFloat_AsDouble(v);
893 return PyLong_FromDouble(x);
894}
895
896static PyObject *
897float_int(PyObject *v)
898{
899 double x = PyFloat_AsDouble(v);
900 double wholepart; /* integral portion of x, rounded toward 0 */
901
902 (void)modf(x, &wholepart);
903 /* Try to get out cheap if this fits in a Python int. The attempt
904 * to cast to long must be protected, as C doesn't define what
905 * happens if the double is too big to fit in a long. Some rare
906 * systems raise an exception then (RISCOS was mentioned as one,
907 * and someone using a non-default option on Sun also bumped into
908 * that). Note that checking for >= and <= LONG_{MIN,MAX} would
909 * still be vulnerable: if a long has more bits of precision than
910 * a double, casting MIN/MAX to double may yield an approximation,
911 * and if that's rounded up, then, e.g., wholepart=LONG_MAX+1 would
912 * yield true from the C expression wholepart<=LONG_MAX, despite
913 * that wholepart is actually greater than LONG_MAX.
914 */
915 if (LONG_MIN < wholepart && wholepart < LONG_MAX) {
916 const long aslong = (long)wholepart;
917 return PyInt_FromLong(aslong);
918 }
919 return PyLong_FromDouble(wholepart);
920}
921
922static PyObject *
923float_float(PyObject *v)
924{
925 if (PyFloat_CheckExact(v))
926 Py_INCREF(v);
927 else
928 v = PyFloat_FromDouble(((PyFloatObject *)v)->ob_fval);
929 return v;
930}
931
932
933static PyObject *
934float_subtype_new(PyTypeObject *type, PyObject *args, PyObject *kwds);
935
936static PyObject *
937float_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
938{
939 PyObject *x = Py_False; /* Integer zero */
940 static char *kwlist[] = {"x", 0};
941
942 if (type != &PyFloat_Type)
943 return float_subtype_new(type, args, kwds); /* Wimp out */
944 if (!PyArg_ParseTupleAndKeywords(args, kwds, "|O:float", kwlist, &x))
945 return NULL;
946 if (PyString_Check(x))
947 return PyFloat_FromString(x, NULL);
948 return PyNumber_Float(x);
949}
950
951/* Wimpy, slow approach to tp_new calls for subtypes of float:
952 first create a regular float from whatever arguments we got,
953 then allocate a subtype instance and initialize its ob_fval
954 from the regular float. The regular float is then thrown away.
955*/
956static PyObject *
957float_subtype_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
958{
959 PyObject *tmp, *newobj;
960
961 assert(PyType_IsSubtype(type, &PyFloat_Type));
962 tmp = float_new(&PyFloat_Type, args, kwds);
963 if (tmp == NULL)
964 return NULL;
965 assert(PyFloat_CheckExact(tmp));
966 newobj = type->tp_alloc(type, 0);
967 if (newobj == NULL) {
968 Py_DECREF(tmp);
969 return NULL;
970 }
971 ((PyFloatObject *)newobj)->ob_fval = ((PyFloatObject *)tmp)->ob_fval;
972 Py_DECREF(tmp);
973 return newobj;
974}
975
976static PyObject *
977float_getnewargs(PyFloatObject *v)
978{
979 return Py_BuildValue("(d)", v->ob_fval);
980}
981
982/* this is for the benefit of the pack/unpack routines below */
983
984typedef enum {
985 unknown_format, ieee_big_endian_format, ieee_little_endian_format
986} float_format_type;
987
988static float_format_type double_format, float_format;
989static float_format_type detected_double_format, detected_float_format;
990
991static PyObject *
992float_getformat(PyTypeObject *v, PyObject* arg)
993{
994 char* s;
995 float_format_type r;
996
997 if (!PyString_Check(arg)) {
998 PyErr_Format(PyExc_TypeError,
999 "__getformat__() argument must be string, not %.500s",
1000 arg->ob_type->tp_name);
1001 return NULL;
1002 }
1003 s = PyString_AS_STRING(arg);
1004 if (strcmp(s, "double") == 0) {
1005 r = double_format;
1006 }
1007 else if (strcmp(s, "float") == 0) {
1008 r = float_format;
1009 }
1010 else {
1011 PyErr_SetString(PyExc_ValueError,
1012 "__getformat__() argument 1 must be "
1013 "'double' or 'float'");
1014 return NULL;
1015 }
1016
1017 switch (r) {
1018 case unknown_format:
1019 return PyString_FromString("unknown");
1020 case ieee_little_endian_format:
1021 return PyString_FromString("IEEE, little-endian");
1022 case ieee_big_endian_format:
1023 return PyString_FromString("IEEE, big-endian");
1024 default:
1025 Py_FatalError("insane float_format or double_format");
1026 return NULL;
1027 }
1028}
1029
1030PyDoc_STRVAR(float_getformat_doc,
1031"float.__getformat__(typestr) -> string\n"
1032"\n"
1033"You probably don't want to use this function. It exists mainly to be\n"
1034"used in Python's test suite.\n"
1035"\n"
1036"typestr must be 'double' or 'float'. This function returns whichever of\n"
1037"'unknown', 'IEEE, big-endian' or 'IEEE, little-endian' best describes the\n"
1038"format of floating point numbers used by the C type named by typestr.");
1039
1040static PyObject *
1041float_setformat(PyTypeObject *v, PyObject* args)
1042{
1043 char* typestr;
1044 char* format;
1045 float_format_type f;
1046 float_format_type detected;
1047 float_format_type *p;
1048
1049 if (!PyArg_ParseTuple(args, "ss:__setformat__", &typestr, &format))
1050 return NULL;
1051
1052 if (strcmp(typestr, "double") == 0) {
1053 p = &double_format;
1054 detected = detected_double_format;
1055 }
1056 else if (strcmp(typestr, "float") == 0) {
1057 p = &float_format;
1058 detected = detected_float_format;
1059 }
1060 else {
1061 PyErr_SetString(PyExc_ValueError,
1062 "__setformat__() argument 1 must "
1063 "be 'double' or 'float'");
1064 return NULL;
1065 }
1066
1067 if (strcmp(format, "unknown") == 0) {
1068 f = unknown_format;
1069 }
1070 else if (strcmp(format, "IEEE, little-endian") == 0) {
1071 f = ieee_little_endian_format;
1072 }
1073 else if (strcmp(format, "IEEE, big-endian") == 0) {
1074 f = ieee_big_endian_format;
1075 }
1076 else {
1077 PyErr_SetString(PyExc_ValueError,
1078 "__setformat__() argument 2 must be "
1079 "'unknown', 'IEEE, little-endian' or "
1080 "'IEEE, big-endian'");
1081 return NULL;
1082
1083 }
1084
1085 if (f != unknown_format && f != detected) {
1086 PyErr_Format(PyExc_ValueError,
1087 "can only set %s format to 'unknown' or the "
1088 "detected platform value", typestr);
1089 return NULL;
1090 }
1091
1092 *p = f;
1093 Py_RETURN_NONE;
1094}
1095
1096PyDoc_STRVAR(float_setformat_doc,
1097"float.__setformat__(typestr, fmt) -> None\n"
1098"\n"
1099"You probably don't want to use this function. It exists mainly to be\n"
1100"used in Python's test suite.\n"
1101"\n"
1102"typestr must be 'double' or 'float'. fmt must be one of 'unknown',\n"
1103"'IEEE, big-endian' or 'IEEE, little-endian', and in addition can only be\n"
1104"one of the latter two if it appears to match the underlying C reality.\n"
1105"\n"
1106"Overrides the automatic determination of C-level floating point type.\n"
1107"This affects how floats are converted to and from binary strings.");
1108
1109static PyMethodDef float_methods[] = {
1110 {"__getnewargs__", (PyCFunction)float_getnewargs, METH_NOARGS},
1111 {"__getformat__", (PyCFunction)float_getformat,
1112 METH_O|METH_CLASS, float_getformat_doc},
1113 {"__setformat__", (PyCFunction)float_setformat,
1114 METH_VARARGS|METH_CLASS, float_setformat_doc},
1115 {NULL, NULL} /* sentinel */
1116};
1117
1118PyDoc_STRVAR(float_doc,
1119"float(x) -> floating point number\n\
1120\n\
1121Convert a string or number to a floating point number, if possible.");
1122
1123
1124static PyNumberMethods float_as_number = {
1125 float_add, /*nb_add*/
1126 float_sub, /*nb_subtract*/
1127 float_mul, /*nb_multiply*/
1128 float_classic_div, /*nb_divide*/
1129 float_rem, /*nb_remainder*/
1130 float_divmod, /*nb_divmod*/
1131 float_pow, /*nb_power*/
1132 (unaryfunc)float_neg, /*nb_negative*/
1133 (unaryfunc)float_pos, /*nb_positive*/
1134 (unaryfunc)float_abs, /*nb_absolute*/
1135 (inquiry)float_nonzero, /*nb_nonzero*/
1136 0, /*nb_invert*/
1137 0, /*nb_lshift*/
1138 0, /*nb_rshift*/
1139 0, /*nb_and*/
1140 0, /*nb_xor*/
1141 0, /*nb_or*/
1142 float_coerce, /*nb_coerce*/
1143 float_int, /*nb_int*/
1144 float_long, /*nb_long*/
1145 float_float, /*nb_float*/
1146 0, /* nb_oct */
1147 0, /* nb_hex */
1148 0, /* nb_inplace_add */
1149 0, /* nb_inplace_subtract */
1150 0, /* nb_inplace_multiply */
1151 0, /* nb_inplace_divide */
1152 0, /* nb_inplace_remainder */
1153 0, /* nb_inplace_power */
1154 0, /* nb_inplace_lshift */
1155 0, /* nb_inplace_rshift */
1156 0, /* nb_inplace_and */
1157 0, /* nb_inplace_xor */
1158 0, /* nb_inplace_or */
1159 float_floor_div, /* nb_floor_divide */
1160 float_div, /* nb_true_divide */
1161 0, /* nb_inplace_floor_divide */
1162 0, /* nb_inplace_true_divide */
1163};
1164
1165PyTypeObject PyFloat_Type = {
1166 PyObject_HEAD_INIT(&PyType_Type)
1167 0,
1168 "float",
1169 sizeof(PyFloatObject),
1170 0,
1171 (destructor)float_dealloc, /* tp_dealloc */
1172 (printfunc)float_print, /* tp_print */
1173 0, /* tp_getattr */
1174 0, /* tp_setattr */
1175 0, /* tp_compare */
1176 (reprfunc)float_repr, /* tp_repr */
1177 &float_as_number, /* tp_as_number */
1178 0, /* tp_as_sequence */
1179 0, /* tp_as_mapping */
1180 (hashfunc)float_hash, /* tp_hash */
1181 0, /* tp_call */
1182 (reprfunc)float_str, /* tp_str */
1183 PyObject_GenericGetAttr, /* tp_getattro */
1184 0, /* tp_setattro */
1185 0, /* tp_as_buffer */
1186 Py_TPFLAGS_DEFAULT | Py_TPFLAGS_CHECKTYPES |
1187 Py_TPFLAGS_BASETYPE, /* tp_flags */
1188 float_doc, /* tp_doc */
1189 0, /* tp_traverse */
1190 0, /* tp_clear */
1191 float_richcompare, /* tp_richcompare */
1192 0, /* tp_weaklistoffset */
1193 0, /* tp_iter */
1194 0, /* tp_iternext */
1195 float_methods, /* tp_methods */
1196 0, /* tp_members */
1197 0, /* tp_getset */
1198 0, /* tp_base */
1199 0, /* tp_dict */
1200 0, /* tp_descr_get */
1201 0, /* tp_descr_set */
1202 0, /* tp_dictoffset */
1203 0, /* tp_init */
1204 0, /* tp_alloc */
1205 float_new, /* tp_new */
1206};
1207
1208void
1209_PyFloat_Init(void)
1210{
1211 /* We attempt to determine if this machine is using IEEE
1212 floating point formats by peering at the bits of some
1213 carefully chosen values. If it looks like we are on an
1214 IEEE platform, the float packing/unpacking routines can
1215 just copy bits, if not they resort to arithmetic & shifts
1216 and masks. The shifts & masks approach works on all finite
1217 values, but what happens to infinities, NaNs and signed
1218 zeroes on packing is an accident, and attempting to unpack
1219 a NaN or an infinity will raise an exception.
1220
1221 Note that if we're on some whacked-out platform which uses
1222 IEEE formats but isn't strictly little-endian or big-
1223 endian, we will fall back to the portable shifts & masks
1224 method. */
1225
1226#if SIZEOF_DOUBLE == 8
1227 {
1228 double x = 9006104071832581.0;
1229 if (memcmp(&x, "\x43\x3f\xff\x01\x02\x03\x04\x05", 8) == 0)
1230 detected_double_format = ieee_big_endian_format;
1231 else if (memcmp(&x, "\x05\x04\x03\x02\x01\xff\x3f\x43", 8) == 0)
1232 detected_double_format = ieee_little_endian_format;
1233 else
1234 detected_double_format = unknown_format;
1235 }
1236#else
1237 detected_double_format = unknown_format;
1238#endif
1239
1240#if SIZEOF_FLOAT == 4
1241 {
1242 float y = 16711938.0;
1243 if (memcmp(&y, "\x4b\x7f\x01\x02", 4) == 0)
1244 detected_float_format = ieee_big_endian_format;
1245 else if (memcmp(&y, "\x02\x01\x7f\x4b", 4) == 0)
1246 detected_float_format = ieee_little_endian_format;
1247 else
1248 detected_float_format = unknown_format;
1249 }
1250#else
1251 detected_float_format = unknown_format;
1252#endif
1253
1254 double_format = detected_double_format;
1255 float_format = detected_float_format;
1256}
1257
1258void
1259PyFloat_Fini(void)
1260{
1261 PyFloatObject *p;
1262 PyFloatBlock *list, *next;
1263 unsigned i;
1264 int bc, bf; /* block count, number of freed blocks */
1265 int frem, fsum; /* remaining unfreed floats per block, total */
1266
1267 bc = 0;
1268 bf = 0;
1269 fsum = 0;
1270 list = block_list;
1271 block_list = NULL;
1272 free_list = NULL;
1273 while (list != NULL) {
1274 bc++;
1275 frem = 0;
1276 for (i = 0, p = &list->objects[0];
1277 i < N_FLOATOBJECTS;
1278 i++, p++) {
1279 if (PyFloat_CheckExact(p) && p->ob_refcnt != 0)
1280 frem++;
1281 }
1282 next = list->next;
1283 if (frem) {
1284 list->next = block_list;
1285 block_list = list;
1286 for (i = 0, p = &list->objects[0];
1287 i < N_FLOATOBJECTS;
1288 i++, p++) {
1289 if (!PyFloat_CheckExact(p) ||
1290 p->ob_refcnt == 0) {
1291 p->ob_type = (struct _typeobject *)
1292 free_list;
1293 free_list = p;
1294 }
1295 }
1296 }
1297 else {
1298 PyMem_FREE(list); /* XXX PyObject_FREE ??? */
1299 bf++;
1300 }
1301 fsum += frem;
1302 list = next;
1303 }
1304 if (!Py_VerboseFlag)
1305 return;
1306 fprintf(stderr, "# cleanup floats");
1307 if (!fsum) {
1308 fprintf(stderr, "\n");
1309 }
1310 else {
1311 fprintf(stderr,
1312 ": %d unfreed float%s in %d out of %d block%s\n",
1313 fsum, fsum == 1 ? "" : "s",
1314 bc - bf, bc, bc == 1 ? "" : "s");
1315 }
1316 if (Py_VerboseFlag > 1) {
1317 list = block_list;
1318 while (list != NULL) {
1319 for (i = 0, p = &list->objects[0];
1320 i < N_FLOATOBJECTS;
1321 i++, p++) {
1322 if (PyFloat_CheckExact(p) &&
1323 p->ob_refcnt != 0) {
1324 char buf[100];
1325 PyFloat_AsString(buf, p);
1326 /* XXX(twouters) cast refcount to
1327 long until %zd is universally
1328 available
1329 */
1330 fprintf(stderr,
1331 "# <float at %p, refcnt=%ld, val=%s>\n",
1332 p, (long)p->ob_refcnt, buf);
1333 }
1334 }
1335 list = list->next;
1336 }
1337 }
1338}
1339
1340/*----------------------------------------------------------------------------
1341 * _PyFloat_{Pack,Unpack}{4,8}. See floatobject.h.
1342 *
1343 * TODO: On platforms that use the standard IEEE-754 single and double
1344 * formats natively, these routines could simply copy the bytes.
1345 */
1346int
1347_PyFloat_Pack4(double x, unsigned char *p, int le)
1348{
1349 if (float_format == unknown_format) {
1350 unsigned char sign;
1351 int e;
1352 double f;
1353 unsigned int fbits;
1354 int incr = 1;
1355
1356 if (le) {
1357 p += 3;
1358 incr = -1;
1359 }
1360
1361 if (x < 0) {
1362 sign = 1;
1363 x = -x;
1364 }
1365 else
1366 sign = 0;
1367
1368 f = frexp(x, &e);
1369
1370 /* Normalize f to be in the range [1.0, 2.0) */
1371 if (0.5 <= f && f < 1.0) {
1372 f *= 2.0;
1373 e--;
1374 }
1375 else if (f == 0.0)
1376 e = 0;
1377 else {
1378 PyErr_SetString(PyExc_SystemError,
1379 "frexp() result out of range");
1380 return -1;
1381 }
1382
1383 if (e >= 128)
1384 goto Overflow;
1385 else if (e < -126) {
1386 /* Gradual underflow */
1387 f = ldexp(f, 126 + e);
1388 e = 0;
1389 }
1390 else if (!(e == 0 && f == 0.0)) {
1391 e += 127;
1392 f -= 1.0; /* Get rid of leading 1 */
1393 }
1394
1395 f *= 8388608.0; /* 2**23 */
1396 fbits = (unsigned int)(f + 0.5); /* Round */
1397 assert(fbits <= 8388608);
1398 if (fbits >> 23) {
1399 /* The carry propagated out of a string of 23 1 bits. */
1400 fbits = 0;
1401 ++e;
1402 if (e >= 255)
1403 goto Overflow;
1404 }
1405
1406 /* First byte */
1407 *p = (sign << 7) | (e >> 1);
1408 p += incr;
1409
1410 /* Second byte */
1411 *p = (char) (((e & 1) << 7) | (fbits >> 16));
1412 p += incr;
1413
1414 /* Third byte */
1415 *p = (fbits >> 8) & 0xFF;
1416 p += incr;
1417
1418 /* Fourth byte */
1419 *p = fbits & 0xFF;
1420
1421 /* Done */
1422 return 0;
1423
1424 Overflow:
1425 PyErr_SetString(PyExc_OverflowError,
1426 "float too large to pack with f format");
1427 return -1;
1428 }
1429 else {
1430 float y = (float)x;
1431 const char *s = (char*)&y;
1432 int i, incr = 1;
1433
1434 if ((float_format == ieee_little_endian_format && !le)
1435 || (float_format == ieee_big_endian_format && le)) {
1436 p += 3;
1437 incr = -1;
1438 }
1439
1440 for (i = 0; i < 4; i++) {
1441 *p = *s++;
1442 p += incr;
1443 }
1444 return 0;
1445 }
1446}
1447
1448int
1449_PyFloat_Pack8(double x, unsigned char *p, int le)
1450{
1451 if (double_format == unknown_format) {
1452 unsigned char sign;
1453 int e;
1454 double f;
1455 unsigned int fhi, flo;
1456 int incr = 1;
1457
1458 if (le) {
1459 p += 7;
1460 incr = -1;
1461 }
1462
1463 if (x < 0) {
1464 sign = 1;
1465 x = -x;
1466 }
1467 else
1468 sign = 0;
1469
1470 f = frexp(x, &e);
1471
1472 /* Normalize f to be in the range [1.0, 2.0) */
1473 if (0.5 <= f && f < 1.0) {
1474 f *= 2.0;
1475 e--;
1476 }
1477 else if (f == 0.0)
1478 e = 0;
1479 else {
1480 PyErr_SetString(PyExc_SystemError,
1481 "frexp() result out of range");
1482 return -1;
1483 }
1484
1485 if (e >= 1024)
1486 goto Overflow;
1487 else if (e < -1022) {
1488 /* Gradual underflow */
1489 f = ldexp(f, 1022 + e);
1490 e = 0;
1491 }
1492 else if (!(e == 0 && f == 0.0)) {
1493 e += 1023;
1494 f -= 1.0; /* Get rid of leading 1 */
1495 }
1496
1497 /* fhi receives the high 28 bits; flo the low 24 bits (== 52 bits) */
1498 f *= 268435456.0; /* 2**28 */
1499 fhi = (unsigned int)f; /* Truncate */
1500 assert(fhi < 268435456);
1501
1502 f -= (double)fhi;
1503 f *= 16777216.0; /* 2**24 */
1504 flo = (unsigned int)(f + 0.5); /* Round */
1505 assert(flo <= 16777216);
1506 if (flo >> 24) {
1507 /* The carry propagated out of a string of 24 1 bits. */
1508 flo = 0;
1509 ++fhi;
1510 if (fhi >> 28) {
1511 /* And it also progagated out of the next 28 bits. */
1512 fhi = 0;
1513 ++e;
1514 if (e >= 2047)
1515 goto Overflow;
1516 }
1517 }
1518
1519 /* First byte */
1520 *p = (sign << 7) | (e >> 4);
1521 p += incr;
1522
1523 /* Second byte */
1524 *p = (unsigned char) (((e & 0xF) << 4) | (fhi >> 24));
1525 p += incr;
1526
1527 /* Third byte */
1528 *p = (fhi >> 16) & 0xFF;
1529 p += incr;
1530
1531 /* Fourth byte */
1532 *p = (fhi >> 8) & 0xFF;
1533 p += incr;
1534
1535 /* Fifth byte */
1536 *p = fhi & 0xFF;
1537 p += incr;
1538
1539 /* Sixth byte */
1540 *p = (flo >> 16) & 0xFF;
1541 p += incr;
1542
1543 /* Seventh byte */
1544 *p = (flo >> 8) & 0xFF;
1545 p += incr;
1546
1547 /* Eighth byte */
1548 *p = flo & 0xFF;
1549 p += incr;
1550
1551 /* Done */
1552 return 0;
1553
1554 Overflow:
1555 PyErr_SetString(PyExc_OverflowError,
1556 "float too large to pack with d format");
1557 return -1;
1558 }
1559 else {
1560 const char *s = (char*)&x;
1561 int i, incr = 1;
1562
1563 if ((double_format == ieee_little_endian_format && !le)
1564 || (double_format == ieee_big_endian_format && le)) {
1565 p += 7;
1566 incr = -1;
1567 }
1568
1569 for (i = 0; i < 8; i++) {
1570 *p = *s++;
1571 p += incr;
1572 }
1573 return 0;
1574 }
1575}
1576
1577double
1578_PyFloat_Unpack4(const unsigned char *p, int le)
1579{
1580 if (float_format == unknown_format) {
1581 unsigned char sign;
1582 int e;
1583 unsigned int f;
1584 double x;
1585 int incr = 1;
1586
1587 if (le) {
1588 p += 3;
1589 incr = -1;
1590 }
1591
1592 /* First byte */
1593 sign = (*p >> 7) & 1;
1594 e = (*p & 0x7F) << 1;
1595 p += incr;
1596
1597 /* Second byte */
1598 e |= (*p >> 7) & 1;
1599 f = (*p & 0x7F) << 16;
1600 p += incr;
1601
1602 if (e == 255) {
1603 PyErr_SetString(
1604 PyExc_ValueError,
1605 "can't unpack IEEE 754 special value "
1606 "on non-IEEE platform");
1607 return -1;
1608 }
1609
1610 /* Third byte */
1611 f |= *p << 8;
1612 p += incr;
1613
1614 /* Fourth byte */
1615 f |= *p;
1616
1617 x = (double)f / 8388608.0;
1618
1619 /* XXX This sadly ignores Inf/NaN issues */
1620 if (e == 0)
1621 e = -126;
1622 else {
1623 x += 1.0;
1624 e -= 127;
1625 }
1626 x = ldexp(x, e);
1627
1628 if (sign)
1629 x = -x;
1630
1631 return x;
1632 }
1633 else {
1634 float x;
1635
1636 if ((float_format == ieee_little_endian_format && !le)
1637 || (float_format == ieee_big_endian_format && le)) {
1638 char buf[4];
1639 char *d = &buf[3];
1640 int i;
1641
1642 for (i = 0; i < 4; i++) {
1643 *d-- = *p++;
1644 }
1645 memcpy(&x, buf, 4);
1646 }
1647 else {
1648 memcpy(&x, p, 4);
1649 }
1650
1651 return x;
1652 }
1653}
1654
1655double
1656_PyFloat_Unpack8(const unsigned char *p, int le)
1657{
1658 if (double_format == unknown_format) {
1659 unsigned char sign;
1660 int e;
1661 unsigned int fhi, flo;
1662 double x;
1663 int incr = 1;
1664
1665 if (le) {
1666 p += 7;
1667 incr = -1;
1668 }
1669
1670 /* First byte */
1671 sign = (*p >> 7) & 1;
1672 e = (*p & 0x7F) << 4;
1673
1674 p += incr;
1675
1676 /* Second byte */
1677 e |= (*p >> 4) & 0xF;
1678 fhi = (*p & 0xF) << 24;
1679 p += incr;
1680
1681 if (e == 2047) {
1682 PyErr_SetString(
1683 PyExc_ValueError,
1684 "can't unpack IEEE 754 special value "
1685 "on non-IEEE platform");
1686 return -1.0;
1687 }
1688
1689 /* Third byte */
1690 fhi |= *p << 16;
1691 p += incr;
1692
1693 /* Fourth byte */
1694 fhi |= *p << 8;
1695 p += incr;
1696
1697 /* Fifth byte */
1698 fhi |= *p;
1699 p += incr;
1700
1701 /* Sixth byte */
1702 flo = *p << 16;
1703 p += incr;
1704
1705 /* Seventh byte */
1706 flo |= *p << 8;
1707 p += incr;
1708
1709 /* Eighth byte */
1710 flo |= *p;
1711
1712 x = (double)fhi + (double)flo / 16777216.0; /* 2**24 */
1713 x /= 268435456.0; /* 2**28 */
1714
1715 if (e == 0)
1716 e = -1022;
1717 else {
1718 x += 1.0;
1719 e -= 1023;
1720 }
1721 x = ldexp(x, e);
1722
1723 if (sign)
1724 x = -x;
1725
1726 return x;
1727 }
1728 else {
1729 double x;
1730
1731 if ((double_format == ieee_little_endian_format && !le)
1732 || (double_format == ieee_big_endian_format && le)) {
1733 char buf[8];
1734 char *d = &buf[7];
1735 int i;
1736
1737 for (i = 0; i < 8; i++) {
1738 *d-- = *p++;
1739 }
1740 memcpy(&x, buf, 8);
1741 }
1742 else {
1743 memcpy(&x, p, 8);
1744 }
1745
1746 return x;
1747 }
1748}
Note: See TracBrowser for help on using the repository browser.