#include "Python.h"

#include <ctype.h>

#include <float.h>

#ifndef PyFloat_MAXFREELIST

#define PyFloat_MAXFREELIST 100

static int numfree = 0;

static PyFloatObject *free_list = NULL;

PyFloat_GetMax(void)

{

return DBL_MAX;

}

PyFloat_GetMin(void)

{

return DBL_MIN;

}

static PyTypeObject FloatInfoType;

PyDoc_STRVAR(floatinfo__doc__,

your system's :file:`float.h` for more information.");

static PyStructSequence_Field floatinfo_fields[] = {

{"max", "DBL_MAX -- maximum representable finite float"},

{"max_exp", "DBL_MAX_EXP -- maximum int e such that radix**(e-1) "

"is representable"},

{"max_10_exp", "DBL_MAX_10_EXP -- maximum int e such that 10**e "

"is representable"},

{"min", "DBL_MIN -- Minimum positive normalizer float"},

{"min_exp", "DBL_MIN_EXP -- minimum int e such that radix**(e-1) "

"is a normalized float"},

{"min_10_exp", "DBL_MIN_10_EXP -- minimum int e such that 10**e is "

"a normalized"},

{"dig", "DBL_DIG -- digits"},

{"mant_dig", "DBL_MANT_DIG -- mantissa digits"},

{"epsilon", "DBL_EPSILON -- Difference between 1 and the next "

"representable float"},

{"radix", "FLT_RADIX -- radix of exponent"},

{"rounds", "FLT_ROUNDS -- addition rounds"},

{0}

};

static PyStructSequence_Desc floatinfo_desc = {

"sys.float_info", /* name */

floatinfo__doc__, /* doc */

floatinfo_fields, /* fields */

11

};

PyObject *

PyFloat_GetInfo(void)

{

PyObject* floatinfo;

int pos = 0;

floatinfo = PyStructSequence_New(&FloatInfoType);

if (floatinfo == NULL) {

return NULL;

}

#define SetIntFlag(flag) \

PyStructSequence_SET_ITEM(floatinfo, pos++, PyLong_FromLong(flag))

#define SetDblFlag(flag) \

PyStructSequence_SET_ITEM(floatinfo, pos++, PyFloat_FromDouble(flag))

SetDblFlag(DBL_MAX);

SetIntFlag(DBL_MAX_EXP);

SetIntFlag(DBL_MAX_10_EXP);

SetDblFlag(DBL_MIN);

SetIntFlag(DBL_MIN_EXP);

SetIntFlag(DBL_MIN_10_EXP);

SetIntFlag(DBL_DIG);

SetIntFlag(DBL_MANT_DIG);

SetDblFlag(DBL_EPSILON);

SetIntFlag(FLT_RADIX);

SetIntFlag(FLT_ROUNDS);

#undef SetIntFlag

#undef SetDblFlag

if (PyErr_Occurred()) {

Py_CLEAR(floatinfo);

return NULL;

}

return floatinfo;

}

PyObject *

PyFloat_FromDouble(double fval)

{

PyFloatObject *op = free_list;

if (op != NULL) {

free_list = (PyFloatObject *) Py_TYPE(op);

numfree--;

} else {

op = (PyFloatObject*) PyObject_MALLOC(sizeof(PyFloatObject));

if (!op)

return PyErr_NoMemory();

}

/* Inline PyObject_New */

PyObject_INIT(op, &PyFloat_Type);

op->ob_fval = fval;

return (PyObject *) op;

}

PyObject *

PyFloat_FromString(PyObject *v)

{

const char *s, *last, *end;

double x;

PyObject *s_buffer = NULL;

Py_ssize_t len;

PyObject *result = NULL;

if (PyUnicode_Check(v)) {

s_buffer = _PyUnicode_TransformDecimalAndSpaceToASCII(v);

if (s_buffer == NULL)

return NULL;

s = PyUnicode_AsUTF8AndSize(s_buffer, &len);

if (s == NULL) {

Py_DECREF(s_buffer);

return NULL;

}

}

else if (PyObject_AsCharBuffer(v, &s, &len)) {

PyErr_SetString(PyExc_TypeError,

"float() argument must be a string or a number");

return NULL;

}

last = s + len;

/* strip space */

while (s < last && Py_ISSPACE(*s))

s++;

while (s < last - 1 && Py_ISSPACE(last[-1]))

last--;

/* We don't care about overflow or underflow. If the platform

x = PyOS_string_to_double(s, (char **)&end, NULL);

if (end != last) {

PyErr_Format(PyExc_ValueError,

"could not convert string to float: "

"%R", v);

result = NULL;

}

else if (x == -1.0 && PyErr_Occurred())

result = NULL;

else

result = PyFloat_FromDouble(x);

Py_XDECREF(s_buffer);

return result;

}

static void

float_dealloc(PyFloatObject *op)

{

if (PyFloat_CheckExact(op)) {

if (numfree >= PyFloat_MAXFREELIST) {

PyObject_FREE(op);

return;

}

numfree++;

Py_TYPE(op) = (struct _typeobject *)free_list;

free_list = op;

}

else

Py_TYPE(op)->tp_free((PyObject *)op);

}

PyFloat_AsDouble(PyObject *op)

{

PyNumberMethods *nb;

PyFloatObject *fo;

double val;

if (op && PyFloat_Check(op))

return PyFloat_AS_DOUBLE((PyFloatObject*) op);

if (op == NULL) {

PyErr_BadArgument();

return -1;

}

if ((nb = Py_TYPE(op)->tp_as_number) == NULL || nb->nb_float == NULL) {

PyErr_SetString(PyExc_TypeError, "a float is required");

return -1;

}

fo = (PyFloatObject*) (*nb->nb_float) (op);

if (fo == NULL)

return -1;

if (!PyFloat_Check(fo)) {

PyErr_SetString(PyExc_TypeError,

"nb_float should return float object");

return -1;

}

val = PyFloat_AS_DOUBLE(fo);

Py_DECREF(fo);

return val;

}

#define CONVERT_TO_DOUBLE(obj, dbl) \

if (PyFloat_Check(obj)) \

dbl = PyFloat_AS_DOUBLE(obj); \

else if (convert_to_double(&(obj), &(dbl)) < 0) \

return obj;

static int

convert_to_double(PyObject **v, double *dbl)

{

PyObject *obj = *v;

if (PyLong_Check(obj)) {

*dbl = PyLong_AsDouble(obj);

if (*dbl == -1.0 && PyErr_Occurred()) {

*v = NULL;

return -1;

}

}

else {

Py_INCREF(Py_NotImplemented);

*v = Py_NotImplemented;

return -1;

}

return 0;

}

static PyObject *

float_repr(PyFloatObject *v)

{

PyObject *result;

char *buf;

buf = PyOS_double_to_string(PyFloat_AS_DOUBLE(v),

'r', 0,

Py_DTSF_ADD_DOT_0,

NULL);

if (!buf)

return PyErr_NoMemory();

result = _PyUnicode_FromASCII(buf, strlen(buf));

PyMem_Free(buf);

return result;

}

static PyObject*

float_richcompare(PyObject *v, PyObject *w, int op)

{

double i, j;

int r = 0;

assert(PyFloat_Check(v));

i = PyFloat_AS_DOUBLE(v);

/* Switch on the type of w. Set i and j to doubles to be compared,

if (PyFloat_Check(w))

j = PyFloat_AS_DOUBLE(w);

else if (!Py_IS_FINITE(i)) {

if (PyLong_Check(w))

/* If i is an infinity, its magnitude exceeds any

j = 0.0;

else

goto Unimplemented;

}

else if (PyLong_Check(w)) {

int vsign = i == 0.0 ? 0 : i < 0.0 ? -1 : 1;

int wsign = _PyLong_Sign(w);

size_t nbits;

int exponent;

if (vsign != wsign) {

/* Magnitudes are irrelevant -- the signs alone

i = (double)vsign;

j = (double)wsign;

goto Compare;

}

/* The signs are the same. */

/* Convert w to a double if it fits. In particular, 0 fits. */

nbits = _PyLong_NumBits(w);

if (nbits == (size_t)-1 && PyErr_Occurred()) {

/* This long is so large that size_t isn't big enough

PyErr_Clear();

i = (double)vsign;

assert(wsign != 0);

j = wsign * 2.0;

goto Compare;

}

if (nbits <= 48) {

j = PyLong_AsDouble(w);

/* It's impossible that <= 48 bits overflowed. */

assert(j != -1.0 || ! PyErr_Occurred());

goto Compare;

}

assert(wsign != 0); /* else nbits was 0 */

assert(vsign != 0); /* if vsign were 0, then since wsign is

/* We want to work with non-negative numbers. */

if (vsign < 0) {

/* "Multiply both sides" by -1; this also swaps the

i = -i;

op = _Py_SwappedOp[op];

}

assert(i > 0.0);

(void) frexp(i, &exponent);

/* exponent is the # of bits in v before the radix point;

if (exponent < 0 || (size_t)exponent < nbits) {

i = 1.0;

j = 2.0;

goto Compare;

}

if ((size_t)exponent > nbits) {

i = 2.0;

j = 1.0;

goto Compare;

}

/* v and w have the same number of bits before the radix

{

double fracpart;

double intpart;

PyObject *result = NULL;

PyObject *one = NULL;

PyObject *vv = NULL;

PyObject *ww = w;

if (wsign < 0) {

ww = PyNumber_Negative(w);

if (ww == NULL)

goto Error;

}

else

Py_INCREF(ww);

fracpart = modf(i, &intpart);

vv = PyLong_FromDouble(intpart);

if (vv == NULL)

goto Error;

if (fracpart != 0.0) {

/* Shift left, and or a 1 bit into vv

PyObject *temp;

one = PyLong_FromLong(1);

if (one == NULL)

goto Error;

temp = PyNumber_Lshift(ww, one);

if (temp == NULL)

goto Error;

Py_DECREF(ww);

ww = temp;

temp = PyNumber_Lshift(vv, one);

if (temp == NULL)

goto Error;

Py_DECREF(vv);

vv = temp;

temp = PyNumber_Or(vv, one);

if (temp == NULL)

goto Error;

Py_DECREF(vv);

vv = temp;

}

r = PyObject_RichCompareBool(vv, ww, op);

if (r < 0)

goto Error;

result = PyBool_FromLong(r);

Error:

Py_XDECREF(vv);

Py_XDECREF(ww);

Py_XDECREF(one);

return result;

}

} /* else if (PyLong_Check(w)) */

else /* w isn't float or int */

goto Unimplemented;

Compare:

PyFPE_START_PROTECT("richcompare", return NULL)

switch (op) {

case Py_EQ:

r = i == j;

break;

case Py_NE:

r = i != j;

break;

case Py_LE:

r = i <= j;

break;

case Py_GE:

r = i >= j;

break;

case Py_LT:

r = i < j;

break;

case Py_GT:

r = i > j;

break;

}

PyFPE_END_PROTECT(r)

return PyBool_FromLong(r);

Unimplemented:

Py_RETURN_NOTIMPLEMENTED;

}

static Py_hash_t

float_hash(PyFloatObject *v)

{

return _Py_HashDouble(v->ob_fval);

}

static PyObject *

float_add(PyObject *v, PyObject *w)

{

double a,b;

CONVERT_TO_DOUBLE(v, a);

CONVERT_TO_DOUBLE(w, b);

PyFPE_START_PROTECT("add", return 0)

a = a + b;

PyFPE_END_PROTECT(a)

return PyFloat_FromDouble(a);

}

static PyObject *

float_sub(PyObject *v, PyObject *w)

{

double a,b;

CONVERT_TO_DOUBLE(v, a);

CONVERT_TO_DOUBLE(w, b);

PyFPE_START_PROTECT("subtract", return 0)

a = a - b;

PyFPE_END_PROTECT(a)

return PyFloat_FromDouble(a);

}

static PyObject *

float_mul(PyObject *v, PyObject *w)

{

double a,b;

CONVERT_TO_DOUBLE(v, a);

CONVERT_TO_DOUBLE(w, b);

PyFPE_START_PROTECT("multiply", return 0)

a = a * b;

PyFPE_END_PROTECT(a)

return PyFloat_FromDouble(a);

}

static PyObject *

float_div(PyObject *v, PyObject *w)

{

double a,b;

CONVERT_TO_DOUBLE(v, a);

CONVERT_TO_DOUBLE(w, b);

if (b == 0.0) {

PyErr_SetString(PyExc_ZeroDivisionError,

"float division by zero");

return NULL;

}

PyFPE_START_PROTECT("divide", return 0)

a = a / b;

PyFPE_END_PROTECT(a)

return PyFloat_FromDouble(a);

}

static PyObject *

float_rem(PyObject *v, PyObject *w)

{

double vx, wx;

double mod;

CONVERT_TO_DOUBLE(v, vx);

CONVERT_TO_DOUBLE(w, wx);

if (wx == 0.0) {

PyErr_SetString(PyExc_ZeroDivisionError,

"float modulo");

return NULL;

}

PyFPE_START_PROTECT("modulo", return 0)

mod = fmod(vx, wx);

if (mod) {

/* ensure the remainder has the same sign as the denominator */

if ((wx < 0) != (mod < 0)) {

mod += wx;

}

}

else {

/* the remainder is zero, and in the presence of signed zeroes

mod = copysign(0.0, wx);

}

PyFPE_END_PROTECT(mod)

return PyFloat_FromDouble(mod);

}

static PyObject *

float_divmod(PyObject *v, PyObject *w)

{

double vx, wx;

double div, mod, floordiv;

CONVERT_TO_DOUBLE(v, vx);

CONVERT_TO_DOUBLE(w, wx);

if (wx == 0.0) {

PyErr_SetString(PyExc_ZeroDivisionError, "float divmod()");

return NULL;

}

PyFPE_START_PROTECT("divmod", return 0)

mod = fmod(vx, wx);

/* fmod is typically exact, so vx-mod is *mathematically* an

div = (vx - mod) / wx;

if (mod) {

/* ensure the remainder has the same sign as the denominator */

if ((wx < 0) != (mod < 0)) {

mod += wx;

div -= 1.0;

}

}

else {

/* the remainder is zero, and in the presence of signed zeroes

mod = copysign(0.0, wx);

}

/* snap quotient to nearest integral value */

if (div) {

floordiv = floor(div);

if (div - floordiv > 0.5)

floordiv += 1.0;

}

else {

/* div is zero - get the same sign as the true quotient */

floordiv = copysign(0.0, vx / wx); /* zero w/ sign of vx/wx */

}

PyFPE_END_PROTECT(floordiv)

return Py_BuildValue("(dd)", floordiv, mod);

}

static PyObject *

float_floor_div(PyObject *v, PyObject *w)

{

PyObject *t, *r;

t = float_divmod(v, w);

if (t == NULL || t == Py_NotImplemented)

return t;

assert(PyTuple_CheckExact(t));

r = PyTuple_GET_ITEM(t, 0);

Py_INCREF(r);

Py_DECREF(t);

return r;

}

#define DOUBLE_IS_ODD_INTEGER(x) (fmod(fabs(x), 2.0) == 1.0)

static PyObject *

float_pow(PyObject *v, PyObject *w, PyObject *z)

{

double iv, iw, ix;

int negate_result = 0;

if ((PyObject *)z != Py_None) {

PyErr_SetString(PyExc_TypeError, "pow() 3rd argument not "

"allowed unless all arguments are integers");

return NULL;

}

CONVERT_TO_DOUBLE(v, iv);

CONVERT_TO_DOUBLE(w, iw);

/* Sort out special cases here instead of relying on pow() */

if (iw == 0) { /* v**0 is 1, even 0**0 */

return PyFloat_FromDouble(1.0);

}

if (Py_IS_NAN(iv)) { /* nan**w = nan, unless w == 0 */

return PyFloat_FromDouble(iv);

}

if (Py_IS_NAN(iw)) { /* v**nan = nan, unless v == 1; 1**nan = 1 */

return PyFloat_FromDouble(iv == 1.0 ? 1.0 : iw);

}

if (Py_IS_INFINITY(iw)) {

/* v**inf is: 0.0 if abs(v) < 1; 1.0 if abs(v) == 1; inf if

iv = fabs(iv);

if (iv == 1.0)

return PyFloat_FromDouble(1.0);

else if ((iw > 0.0) == (iv > 1.0))

return PyFloat_FromDouble(fabs(iw)); /* return inf */

else

return PyFloat_FromDouble(0.0);

}

if (Py_IS_INFINITY(iv)) {

/* (+-inf)**w is: inf for w positive, 0 for w negative; in

int iw_is_odd = DOUBLE_IS_ODD_INTEGER(iw);

if (iw > 0.0)

return PyFloat_FromDouble(iw_is_odd ? iv : fabs(iv));

else

return PyFloat_FromDouble(iw_is_odd ?

copysign(0.0, iv) : 0.0);

}

if (iv == 0.0) { /* 0**w is: 0 for w positive, 1 for w zero

int iw_is_odd = DOUBLE_IS_ODD_INTEGER(iw);

if (iw < 0.0) {

PyErr_SetString(PyExc_ZeroDivisionError,

"0.0 cannot be raised to a "

"negative power");

return NULL;

}

/* use correct sign if iw is odd */

return PyFloat_FromDouble(iw_is_odd ? iv : 0.0);

}

if (iv < 0.0) {

/* Whether this is an error is a mess, and bumps into libm

if (iw != floor(iw)) {

/* Negative numbers raised to fractional powers

return PyComplex_Type.tp_as_number->nb_power(v, w, z);

}

/* iw is an exact integer, albeit perhaps a very large

iv = -iv;

negate_result = DOUBLE_IS_ODD_INTEGER(iw);

}

if (iv == 1.0) { /* 1**w is 1, even 1**inf and 1**nan */

/* (-1) ** large_integer also ends up here. Here's an

return PyFloat_FromDouble(negate_result ? -1.0 : 1.0);

}

/* Now iv and iw are finite, iw is nonzero, and iv is

errno = 0;

PyFPE_START_PROTECT("pow", return NULL)

ix = pow(iv, iw);

PyFPE_END_PROTECT(ix)

Py_ADJUST_ERANGE1(ix);

if (negate_result)

ix = -ix;

if (errno != 0) {

/* We don't expect any errno value other than ERANGE, but

PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError :

PyExc_ValueError);

return NULL;

}

return PyFloat_FromDouble(ix);

}

#undef DOUBLE_IS_ODD_INTEGER

static PyObject *

float_neg(PyFloatObject *v)

{

return PyFloat_FromDouble(-v->ob_fval);

}

static PyObject *

float_abs(PyFloatObject *v)

{

return PyFloat_FromDouble(fabs(v->ob_fval));

}

static int

float_bool(PyFloatObject *v)

{

return v->ob_fval != 0.0;

}

static PyObject *

float_is_integer(PyObject *v)

{

double x = PyFloat_AsDouble(v);

PyObject *o;

if (x == -1.0 && PyErr_Occurred())

return NULL;

if (!Py_IS_FINITE(x))

Py_RETURN_FALSE;

errno = 0;

PyFPE_START_PROTECT("is_integer", return NULL)

o = (floor(x) == x) ? Py_True : Py_False;

PyFPE_END_PROTECT(x)

if (errno != 0) {

PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError :

PyExc_ValueError);

return NULL;

}

Py_INCREF(o);

return o;

}

#if 0

static PyObject *

float_is_inf(PyObject *v)

{

double x = PyFloat_AsDouble(v);

if (x == -1.0 && PyErr_Occurred())

return NULL;

return PyBool_FromLong((long)Py_IS_INFINITY(x));

}

static PyObject *

float_is_nan(PyObject *v)

{

double x = PyFloat_AsDouble(v);

if (x == -1.0 && PyErr_Occurred())

return NULL;

return PyBool_FromLong((long)Py_IS_NAN(x));

}

static PyObject *

float_is_finite(PyObject *v)

{

double x = PyFloat_AsDouble(v);

if (x == -1.0 && PyErr_Occurred())

return NULL;

return PyBool_FromLong((long)Py_IS_FINITE(x));

}

static PyObject *

float_trunc(PyObject *v)

{

double x = PyFloat_AsDouble(v);

double wholepart; /* integral portion of x, rounded toward 0 */

(void)modf(x, &wholepart);

/* Try to get out cheap if this fits in a Python int. The attempt

if (LONG_MIN < wholepart && wholepart < LONG_MAX) {

const long aslong = (long)wholepart;

return PyLong_FromLong(aslong);

}

return PyLong_FromDouble(wholepart);

}

#ifndef PY_NO_SHORT_FLOAT_REPR

static PyObject *

double_round(double x, int ndigits) {

double rounded;

Py_ssize_t buflen, mybuflen=100;

char *buf, *buf_end, shortbuf[100], *mybuf=shortbuf;

int decpt, sign;

PyObject *result = NULL;

_Py_SET_53BIT_PRECISION_HEADER;

/* round to a decimal string */

_Py_SET_53BIT_PRECISION_START;

buf = _Py_dg_dtoa(x, 3, ndigits, &decpt, &sign, &buf_end);

_Py_SET_53BIT_PRECISION_END;

if (buf == NULL) {

PyErr_NoMemory();

return NULL;

}

/* Get new buffer if shortbuf is too small. Space needed <= buf_end -

buflen = buf_end - buf;

if (buflen + 8 > mybuflen) {

mybuflen = buflen+8;

mybuf = (char *)PyMem_Malloc(mybuflen);

if (mybuf == NULL) {

PyErr_NoMemory();

goto exit;

}

}

/* copy buf to mybuf, adding exponent, sign and leading 0 */

PyOS_snprintf(mybuf, mybuflen, "%s0%se%d", (sign ? "-" : ""),

buf, decpt - (int)buflen);

/* and convert the resulting string back to a double */

errno = 0;

_Py_SET_53BIT_PRECISION_START;

rounded = _Py_dg_strtod(mybuf, NULL);

_Py_SET_53BIT_PRECISION_END;

if (errno == ERANGE && fabs(rounded) >= 1.)

PyErr_SetString(PyExc_OverflowError,

"rounded value too large to represent");

else

result = PyFloat_FromDouble(rounded);

/* done computing value; now clean up */

if (mybuf != shortbuf)

PyMem_Free(mybuf);

exit:

_Py_dg_freedtoa(buf);

return result;

}

#else /* PY_NO_SHORT_FLOAT_REPR */

static PyObject *

double_round(double x, int ndigits) {

double pow1, pow2, y, z;

if (ndigits >= 0) {

if (ndigits > 22) {

/* pow1 and pow2 are each safe from overflow, but

pow1 = pow(10.0, (double)(ndigits-22));

pow2 = 1e22;

}

else {

pow1 = pow(10.0, (double)ndigits);

pow2 = 1.0;

}

y = (x*pow1)*pow2;

/* if y overflows, then rounded value is exactly x */

if (!Py_IS_FINITE(y))

return PyFloat_FromDouble(x);

}

else {

pow1 = pow(10.0, (double)-ndigits);

pow2 = 1.0; /* unused; silences a gcc compiler warning */

y = x / pow1;

}

z = round(y);

if (fabs(y-z) == 0.5)

/* halfway between two integers; use round-half-even */

z = 2.0*round(y/2.0);

if (ndigits >= 0)

z = (z / pow2) / pow1;

else

z *= pow1;

/* if computation resulted in overflow, raise OverflowError */

if (!Py_IS_FINITE(z)) {

PyErr_SetString(PyExc_OverflowError,

"overflow occurred during round");

return NULL;

}

return PyFloat_FromDouble(z);

}

#endif /* PY_NO_SHORT_FLOAT_REPR */

static PyObject *

float_round(PyObject *v, PyObject *args)

{

double x, rounded;

PyObject *o_ndigits = NULL;

Py_ssize_t ndigits;

x = PyFloat_AsDouble(v);

if (!PyArg_ParseTuple(args, "|O", &o_ndigits))

return NULL;

if (o_ndigits == NULL) {

/* single-argument round: round to nearest integer */

rounded = round(x);

if (fabs(x-rounded) == 0.5)

/* halfway case: round to even */

rounded = 2.0*round(x/2.0);

return PyLong_FromDouble(rounded);

}

/* interpret second argument as a Py_ssize_t; clips on overflow */

ndigits = PyNumber_AsSsize_t(o_ndigits, NULL);

if (ndigits == -1 && PyErr_Occurred())

return NULL;

/* nans and infinities round to themselves */

if (!Py_IS_FINITE(x))

return PyFloat_FromDouble(x);

/* Deal with extreme values for ndigits. For ndigits > NDIGITS_MAX, x