/src/gmpy_mpfr.c
C | 2919 lines | 2354 code | 492 blank | 73 comment | 356 complexity | dc9c4e536aaf2cff3e22b6ad4dd1a020 MD5 | raw file
Large files files are truncated, but you can click here to view the full file
1/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 2 * gmpy_mpfr.c * 3 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 4 * Python interface to the GMP or MPIR, MPFR, and MPC multiple precision * 5 * libraries. * 6 * * 7 * Copyright 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, * 8 * 2008, 2009 Alex Martelli * 9 * * 10 * Copyright 2008, 2009, 2010, 2011, 2012, 2013 Case Van Horsen * 11 * * 12 * This file is part of GMPY2. * 13 * * 14 * GMPY2 is free software: you can redistribute it and/or modify it under * 15 * the terms of the GNU Lesser General Public License as published by the * 16 * Free Software Foundation, either version 3 of the License, or (at your * 17 * option) any later version. * 18 * * 19 * GMPY2 is distributed in the hope that it will be useful, but WITHOUT * 20 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or * 21 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public * 22 * License for more details. * 23 * * 24 * You should have received a copy of the GNU Lesser General Public * 25 * License along with GMPY2; if not, see <http://www.gnu.org/licenses/> * 26 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ 27 28PyDoc_STRVAR(doc_g_mpfr_f2q, 29"f2q(x,[err]) -> mpq\n\n" 30"Return the 'best' mpq approximating x to within relative error 'err'.\n" 31"Default is the precision of x. Uses Stern-Brocot tree to find the\n" 32"'best' approximation. An 'mpz' is returned if the the denominator\n" 33"is 1. If 'err'<0, error sought is 2.0 ** err."); 34 35static PyObject * 36Pympfr_f2q(PyObject *self, PyObject *args) 37{ 38 PympfrObject *err = 0; 39 PyObject *result; 40 41 if (!PyArg_ParseTuple(args, "O&|O&", Pympfr_convert_arg, &self, 42 Pympfr_convert_arg, &err)) { 43 TYPE_ERROR("f2q() requires 'mpfr', ['mpfr'] arguments"); 44 return NULL; 45 } 46 47 result = (PyObject*)stern_brocot((PympfrObject*)self, err, 0, 1); 48 Py_DECREF(self); 49 Py_XDECREF((PyObject*)err); 50 return result; 51} 52 53PyDoc_STRVAR(doc_mpfr, 54"mpfr() -> mpfr(0.0)\n\n" 55" If no argument is given, return mpfr(0.0).\n\n" 56"mpfr(n[, precison=0]) -> mpfr\n\n" 57" Return an 'mpfr' object after converting a numeric value. If\n" 58" no precision, or a precision of 0, is specified; the precison\n" 59" is taken from the current context.\n\n" 60"mpfr(s[, precision=0[, [base=0]]) -> mpfr\n\n" 61" Return 'mpfr' object after converting a string 's' made up of\n" 62" digits in the given base, possibly with fraction-part (with\n" 63" period as a separator) and/or exponent-part (with exponent\n" 64" marker 'e' for base<=10, else '@'). If no precision, or a\n" 65" precision of 0, is specified; the precison is taken from the\n" 66" current context. The base of the string representation must\n" 67" be 0 or in the interval 2 ... 62. If the base is 0, the leading\n" 68" digits of the string are used to identify the base: 0b implies\n" 69" base=2, 0x implies base=16, otherwise base=10 is assumed.\n"); 70 71static PyObject * 72Pygmpy_mpfr(PyObject *self, PyObject *args, PyObject *keywds) 73{ 74 PympfrObject *result = 0; 75 PyObject *arg0; 76 int base = 0; 77 Py_ssize_t argc; 78 /* Assumes mpfr_prec_t is the same as a long. */ 79 mpfr_prec_t bits = 0; 80 static char *kwlist_s[] = {"s", "precision", "base", NULL}; 81 static char *kwlist_n[] = {"n", "precision", NULL}; 82 83 argc = PyTuple_Size(args); 84 if ((argc < 0) || (argc > 3)) { 85 TYPE_ERROR("mpfr() requires 0 to 3 arguments"); 86 return NULL; 87 } 88 89 if (argc == 0) { 90 if ((result = (PympfrObject*)Pympfr_new(0))) { 91 mpfr_set_ui(result->f, 0, context->ctx.mpfr_round); 92 } 93 return (PyObject*)result; 94 } 95 arg0 = PyTuple_GetItem(args, 0); 96 if (PyStrOrUnicode_Check(arg0)) { 97 /* Can have both precision and/or base as keyword arguments. */ 98 if (PyArg_ParseTupleAndKeywords(args, keywds, "O|li", kwlist_s, 99 &arg0, &bits, &base)) { 100 if ((base!=0) && ((base<2)||(base>62))) { 101 VALUE_ERROR("base for mpfr() must be 0 or in the " 102 "interval 2 ... 62"); 103 } 104 else if (bits < 0) { 105 VALUE_ERROR("precision for mpfr() must be >= 0"); 106 } 107 else { 108 result = Pympfr_From_PyStr(arg0, base, bits); 109 } 110 } 111 SUBNORMALIZE(result); 112 return (PyObject*)result; 113 } 114 115 /* Optimize the common case */ 116 if (isReal(arg0) && argc == 1 && !keywds) { 117 result = Pympfr_From_Real(arg0, bits); 118 SUBNORMALIZE(result); 119 return (PyObject*)result; 120 } 121 122 /* Can only have precision as keyword argument. */ 123 if (PyArg_ParseTupleAndKeywords(args, keywds, "O|l", kwlist_n, &arg0, &bits)) { 124 if (bits < 0) { 125 VALUE_ERROR("precision for mpfr() must be >= 0"); 126 } 127 else { 128 result = Pympfr_From_Real(arg0, bits); 129 if (!result) 130 TYPE_ERROR("mpfr() requires numeric or string argument"); 131 } 132 } 133 return (PyObject*)result; 134} 135 136/* Implement the .precision attribute of an mpfr. */ 137 138static PyObject * 139Pympfr_getprec_attrib(PympfrObject *self, void *closure) 140{ 141 return PyIntOrLong_FromSsize_t((Py_ssize_t)mpfr_get_prec(self->f)); 142} 143 144/* Implement the .rc attribute of an mpfr. */ 145 146static PyObject * 147Pympfr_getrc_attrib(PympfrObject *self, void *closure) 148{ 149 return PyIntOrLong_FromLong((long)self->rc); 150} 151 152/* Implement the .imag attribute of an mpfr. */ 153 154static PyObject * 155Pympfr_getimag_attrib(PympfrObject *self, void *closure) 156{ 157 PympfrObject *result; 158 159 if ((result = (PympfrObject*)Pympfr_new(0))) 160 mpfr_set_zero(result->f, 1); 161 return (PyObject*)result; 162} 163 164/* Implement the .real attribute of an mpfr. */ 165 166static PyObject * 167Pympfr_getreal_attrib(PympfrObject *self, void *closure) 168{ 169 return (PyObject*)Pympfr_From_Pympfr((PyObject*)self, 0); 170} 171 172/* Implement the nb_bool slot. */ 173 174static int 175Pympfr_nonzero(PympfrObject *self) 176{ 177 return !mpfr_zero_p(self->f); 178} 179 180/* Implement the conjugate() method. */ 181 182PyDoc_STRVAR(doc_mpfr_conjugate, 183"x.conjugate() -> mpfr\n\n" 184"Return the conjugate of x (which is just a copy of x since x is\n" 185"not a complex number)."); 186 187static PyObject * 188Pympfr_conjugate(PyObject *self, PyObject *args) 189{ 190 return (PyObject*)Pympfr_From_Pympfr(self, 0); 191} 192 193/* Implement the nb_positive slot. */ 194 195/* TODO: can probably just call Pympfr_From_Real. */ 196 197static PyObject * 198Pympfr_pos(PympfrObject *self) 199{ 200 PympfrObject *result; 201 202 if (!(result = (PympfrObject*)Pympfr_new(mpfr_get_prec(self->f)))) 203 return NULL; 204 205 mpfr_clear_flags(); 206 207 /* Since result has the same precision as self, no rounding occurs. */ 208 mpfr_set(result->f, self->f, context->ctx.mpfr_round); 209 result->round_mode = self->round_mode; 210 result->rc = self->rc; 211 /* Force the exponents to be valid. */ 212 result->rc = mpfr_check_range(result->f, result->rc, result->round_mode); 213 /* Now round result to the current precision. */ 214 result->rc = mpfr_prec_round(result->f, context->ctx.mpfr_prec, 215 context->ctx.mpfr_round); 216 217 SUBNORMALIZE(result); 218 MERGE_FLAGS; 219 CHECK_FLAGS("__pos__"); 220 done: 221 if (PyErr_Occurred()) { 222 Py_XDECREF((PyObject*)result); 223 result = NULL; 224 } 225 return (PyObject*)result; 226} 227 228PyDoc_STRVAR(doc_g_mpfr_get_emin_min, 229"get_emin_min() -> integer\n\n" 230"Return the minimum possible exponent that can be set for 'mpfr'."); 231 232static PyObject * 233Pympfr_get_emin_min(PyObject *self, PyObject *args) 234{ 235 return PyIntOrLong_FromSsize_t((Py_ssize_t)mpfr_get_emin_min()); 236} 237 238PyDoc_STRVAR(doc_g_mpfr_get_emax_max, 239"get_emax_max() -> integer\n\n" 240"Return the maximum possible exponent that can be set for 'mpfr'."); 241 242static PyObject * 243Pympfr_get_emax_max(PyObject *self, PyObject *args) 244{ 245 return PyIntOrLong_FromSsize_t((Py_ssize_t)mpfr_get_emax_max()); 246} 247 248PyDoc_STRVAR(doc_g_mpfr_get_max_precision, 249"get_max_precision() -> integer\n\n" 250"Return the maximum bits of precision that can be used for calculations.\n" 251"Note: to allow extra precision for intermediate calculations, avoid\n" 252"setting precision close the maximum precision."); 253 254static PyObject * 255Pympfr_get_max_precision(PyObject *self, PyObject *args) 256{ 257 return PyIntOrLong_FromSsize_t((Py_ssize_t)MPFR_PREC_MAX); 258} 259 260PyDoc_STRVAR(doc_g_mpfr_get_exp, 261"get_exp(mpfr) -> integer\n\n" 262"Return the exponent of an mpfr. Returns 0 for NaN or Infinity and\n" 263"sets the erange flag and will raise an exception if trap_erange\n" 264"is set."); 265 266static PyObject * 267Pympfr_get_exp(PyObject *self, PyObject *other) 268{ 269 PyObject *result = 0; 270 Py_ssize_t exp; 271 272 PARSE_ONE_MPFR_OTHER("get_exp() requires 'mpfr' argument"); 273 274 if (mpfr_regular_p(Pympfr_AS_MPFR(self))) { 275 exp = (Py_ssize_t)mpfr_get_exp(Pympfr_AS_MPFR(self)); 276 result = PyIntOrLong_FromSsize_t((Py_ssize_t)exp); 277 } 278 else if (mpfr_zero_p(Pympfr_AS_MPFR(self))) { 279 Py_DECREF(self); 280 result = PyIntOrLong_FromSsize_t(0); 281 } 282 else { 283 context->ctx.erange = 1; 284 if (context->ctx.trap_erange) { 285 GMPY_ERANGE("Can not get exponent from NaN or Infinity."); 286 } 287 else { 288 result = PyIntOrLong_FromSsize_t(0); 289 } 290 } 291 Py_DECREF(self); 292 return result; 293} 294 295PyDoc_STRVAR(doc_g_mpfr_set_exp, 296"set_exp(mpfr, n) -> mpfr\n\n" 297"Set the exponent of an mpfr to n. If n is outside the range of\n" 298"valid exponents, set_exp() will set the erange flag and either\n" 299"return the original value or raise an exception if trap_erange\n" 300"is set."); 301 302static PyObject * 303Pympfr_set_exp(PyObject *self, PyObject *args) 304{ 305 PympfrObject *result = 0; 306 long exp = 0; 307 308 if (!PyArg_ParseTuple(args, "O&l", Pympfr_convert_arg, &self, &exp)) { 309 TYPE_ERROR("set_exp() requires 'mpfr', 'integer' arguments"); 310 return NULL; 311 } 312 313 if (!(result = Pympfr_From_Pympfr(self, 0))) 314 return NULL; 315 Py_DECREF(self); 316 317 result->rc = mpfr_set_exp(Pympfr_AS_MPFR(result), exp); 318 319 if (result->rc) { 320 context->ctx.erange = 1; 321 if (context->ctx.trap_erange) { 322 GMPY_ERANGE("New exponent is out-of-bounds."); 323 Py_DECREF(result); 324 return NULL; 325 } 326 } 327 328 return (PyObject*)result; 329} 330 331PyDoc_STRVAR(doc_g_mpfr_set_sign, 332"set_sign(mpfr, bool) -> mpfr\n\n" 333"If 'bool' is True, then return an 'mpfr' with the sign bit set."); 334 335static PyObject * 336Pympfr_set_sign(PyObject *self, PyObject *args) 337{ 338 PympfrObject *result = 0; 339 PyObject *boolean = 0; 340 int s; 341 342 if (!PyArg_ParseTuple(args, "O&O", Pympfr_convert_arg, &self, &boolean)) { 343 TYPE_ERROR("set_sign() requires 'mpfr', 'boolean' arguments"); 344 return NULL; 345 } 346 347 if (!(result = (PympfrObject*)Pympfr_new(0))) 348 return NULL; 349 350 s = PyObject_IsTrue(boolean); 351 if (s == -1) { 352 TYPE_ERROR("set_sign() requires 'mpfr', 'boolean' arguments"); 353 Py_DECREF(self); 354 Py_DECREF(boolean); 355 Py_DECREF(result); 356 return NULL; 357 } 358 359 result->rc = mpfr_setsign(Pympfr_AS_MPFR(result), Pympfr_AS_MPFR(self), 360 s, context->ctx.mpfr_round); 361 362 Py_DECREF(self); 363 Py_DECREF(boolean); 364 return (PyObject*)result; 365} 366 367PyDoc_STRVAR(doc_g_mpfr_copy_sign, 368"copy_sign(mpfr, mpfr) -> mpfr\n\n" 369"Return an 'mpfr' composed of the first argument with the sign of the\n" 370"second argument."); 371 372static PyObject * 373Pympfr_copy_sign(PyObject *self, PyObject *args) 374{ 375 PympfrObject *result = 0; 376 PyObject *other = 0; 377 378 if (!PyArg_ParseTuple(args, "O&O&", Pympfr_convert_arg, &self, 379 Pympfr_convert_arg, &other)) { 380 TYPE_ERROR("copy_sign() requires 'mpfr', 'mpfr' arguments"); 381 return NULL; 382 } 383 384 if (!(result = (PympfrObject*)Pympfr_new(0))) 385 return NULL; 386 387 result->rc = mpfr_copysign(Pympfr_AS_MPFR(result), Pympfr_AS_MPFR(self), 388 Pympfr_AS_MPFR(other), context->ctx.mpfr_round); 389 390 Py_DECREF(self); 391 Py_DECREF(other); 392 return (PyObject*)result; 393} 394 395static PyObject * 396Pympfr_div_2exp(PyObject *self, PyObject *args) 397{ 398 PympfrObject *result = 0; 399 unsigned long exp = 0; 400 401 if (!PyArg_ParseTuple(args, "O&k", Pympfr_convert_arg, &self, &exp)) { 402 TYPE_ERROR("div_2exp() requires 'mpfr', 'integer' arguments"); 403 return NULL; 404 } 405 406 if (!(result = (PympfrObject*)Pympfr_new(0))) 407 return NULL; 408 409 mpfr_clear_flags(); 410 411 result->rc = mpfr_div_2ui(Pympfr_AS_MPFR(result), Pympfr_AS_MPFR(self), 412 exp, context->ctx.mpfr_round); 413 414 MPFR_CLEANUP_SELF("div_2exp()"); 415} 416 417static PyObject * 418Pympfr_mul_2exp(PyObject *self, PyObject *args) 419{ 420 PympfrObject *result = 0; 421 unsigned long exp = 0; 422 423 if (!PyArg_ParseTuple(args, "O&k", Pympfr_convert_arg, &self, &exp)) { 424 TYPE_ERROR("mul_2exp() requires 'mpfr', 'integer' arguments"); 425 return NULL; 426 } 427 428 if (!(result = (PympfrObject*)Pympfr_new(0))) 429 return NULL; 430 431 mpfr_clear_flags(); 432 433 result->rc = mpfr_mul_2ui(Pympfr_AS_MPFR(result), Pympfr_AS_MPFR(self), 434 exp, context->ctx.mpfr_round); 435 436 MPFR_CLEANUP_SELF("mul_2exp()"); 437} 438 439PyDoc_STRVAR(doc_g_mpfr_set_nan, 440"nan() -> mpfr\n\n" 441"Return an 'mpfr' initialized to NaN (Not-A-Number)."); 442 443static PyObject * 444Pympfr_set_nan(PyObject *self, PyObject *other) 445{ 446 PympfrObject *result; 447 448 if ((result = (PympfrObject*)Pympfr_new(0))) 449 mpfr_set_nan(result->f); 450 return (PyObject*)result; 451} 452 453PyDoc_STRVAR(doc_g_mpfr_set_inf, 454"inf(n) -> mpfr\n\n" 455"Return an 'mpfr' initialized to Infinity with the same sign as n.\n" 456"If n is not given, +Infinity is returned."); 457 458static PyObject * 459Pympfr_set_inf(PyObject *self, PyObject *args) 460{ 461 PympfrObject *result; 462 long s = 1; 463 464 if (PyTuple_Size(args) == 1) { 465 s = clong_From_Integer(PyTuple_GET_ITEM(args, 0)); 466 if (s == -1 && PyErr_Occurred()) { 467 TYPE_ERROR("inf() requires 'int' argument"); 468 return NULL; 469 } 470 } 471 472 if ((result = (PympfrObject*)Pympfr_new(0))) 473 mpfr_set_inf(result->f, s<0?-1:1); 474 return (PyObject*)result; 475} 476 477PyDoc_STRVAR(doc_g_mpfr_set_zero, 478"zero(n) -> mpfr\n\n" 479"Return an 'mpfr' inialized to 0.0 with the same sign as n.\n" 480"If n is not given, +0.0 is returned."); 481 482static PyObject * 483Pympfr_set_zero(PyObject *self, PyObject *args) 484{ 485 PympfrObject *result; 486 long s = 1; 487 488 if (PyTuple_Size(args) == 1) { 489 s = clong_From_Integer(PyTuple_GET_ITEM(args, 0)); 490 if (s == -1 && PyErr_Occurred()) { 491 TYPE_ERROR("zero() requires 'int' argument"); 492 return NULL; 493 } 494 } 495 496 if ((result = (PympfrObject*)Pympfr_new(0))) 497 mpfr_set_zero(result->f, s<0?-1:1); 498 return (PyObject*)result; 499} 500 501PyDoc_STRVAR(doc_g_mpfr_is_signed, 502"is_signed(x) -> boolean\n\n" 503"Return True if the sign bit of x is set."); 504 505static PyObject * 506Pympfr_is_signed(PyObject *self, PyObject *other) 507{ 508 int res; 509 if(self && Pympfr_Check(self)) { 510 Py_INCREF(self); 511 } 512 else if(Pympfr_Check(other)) { 513 self = other; 514 Py_INCREF((PyObject*)self); 515 } 516 else if (!(self = (PyObject*)Pympfr_From_Real(other, 0))) { 517 TYPE_ERROR("is_signed() requires 'mpfr' argument"); 518 return NULL; 519 } 520 res = mpfr_signbit(Pympfr_AS_MPFR(self)); 521 Py_DECREF(self); 522 if (res) 523 Py_RETURN_TRUE; 524 else 525 Py_RETURN_FALSE; 526} 527 528#define MPFR_TEST_OTHER(NAME, msg) \ 529static PyObject * \ 530Pympfr_is_##NAME(PyObject *self, PyObject *other)\ 531{\ 532 int res;\ 533 if(self && Pympfr_Check(self)) {\ 534 Py_INCREF(self);\ 535 }\ 536 else if(Pympfr_Check(other)) {\ 537 self = other;\ 538 Py_INCREF((PyObject*)self);\ 539 }\ 540 else if (!(self = (PyObject*)Pympfr_From_Real(other, 0))) {\ 541 PyErr_SetString(PyExc_TypeError, msg);\ 542 return NULL;\ 543 }\ 544 res = mpfr_##NAME##_p(Pympfr_AS_MPFR(self));\ 545 Py_DECREF(self);\ 546 if (res)\ 547 Py_RETURN_TRUE;\ 548 else\ 549 Py_RETURN_FALSE;\ 550} 551 552MPFR_TEST_OTHER(nan, "is_nan() requires 'mpfr' argument"); 553 554MPFR_TEST_OTHER(inf, "is_infinite() requires 'mpfr' argument"); 555 556PyDoc_STRVAR(doc_g_mpfr_is_number, 557"is_number(x) -> boolean\n\n" 558"Return True if x is an actual number (i.e. not NaN or Infinity);\n" 559"False otherwise.\n" 560"Note: is_number() is deprecated; please use is_finite()."); 561 562MPFR_TEST_OTHER(number, "is_finite() requires 'mpfr' argument"); 563 564MPFR_TEST_OTHER(zero, "is_zero() requires 'mpfr' argument"); 565 566PyDoc_STRVAR(doc_g_mpfr_is_regular, 567"is_regular(x) -> boolean\n\n" 568"Return True if x is not zero, NaN, or Infinity; False otherwise."); 569 570MPFR_TEST_OTHER(regular, "is_regular() requires 'mpfr' argument"); 571 572PyDoc_STRVAR(doc_mpfr_is_integer, 573"x.is_integer() -> boolean\n\n" 574"Return True if x is an integer; False otherwise."); 575 576PyDoc_STRVAR(doc_g_mpfr_is_integer, 577"is_integer(x) -> boolean\n\n" 578"Return True if x is an integer; False otherwise."); 579 580MPFR_TEST_OTHER(integer, "is_integer() requires 'mpfr' argument"); 581 582/* produce string for an mpfr with requested/defaulted parameters */ 583 584PyDoc_STRVAR(doc_mpfr_digits, 585"x.digits([base=10[, prec=0]]) -> (mantissa, exponent, bits)\n\n" 586"Returns up to 'prec' digits in the given base. If 'prec' is 0, as many\n" 587"digits that are available are returned. No more digits than available\n" 588"given x's precision are returned. 'base' must be between 2 and 62,\n" 589"inclusive. The result is a three element tuple containing the mantissa,\n" 590"the exponent, and the number of bits of precision."); 591 592/* TODO: support keyword arguments. */ 593 594static PyObject * 595Pympfr_digits(PyObject *self, PyObject *args) 596{ 597 int base = 10; 598 int prec = 0; 599 PyObject *result; 600 601 if (self && Pympfr_Check(self)) { 602 if (!PyArg_ParseTuple(args, "|ii", &base, &prec)) 603 return NULL; 604 Py_INCREF(self); 605 } 606 else { 607 if(!PyArg_ParseTuple(args, "O&|ii", Pympfr_convert_arg, &self, 608 &base, &prec)) 609 return NULL; 610 } 611 result = Pympfr_To_PyStr((PympfrObject*)self, base, prec); 612 Py_DECREF(self); 613 return result; 614} 615 616PyDoc_STRVAR(doc_mpfr_integer_ratio, 617"x.as_integer_ratio() -> (num, den)\n\n" 618"Return the exact rational equivalent of an mpfr. Value is a tuple\n" 619"for compatibility with Python's float.as_integer_ratio()."); 620 621static PyObject * 622Pympfr_integer_ratio(PyObject *self, PyObject *args) 623{ 624 PympzObject *num = 0, *den = 0; 625 mpfr_exp_t temp, twocount; 626 PyObject *result; 627 628 if (mpfr_nan_p(Pympfr_AS_MPFR(self))) { 629 VALUE_ERROR("Cannot pass NaN to mpfr.as_integer_ratio."); 630 return NULL; 631 } 632 633 if (mpfr_inf_p(Pympfr_AS_MPFR(self))) { 634 OVERFLOW_ERROR("Cannot pass Infinity to mpfr.as_integer_ratio."); 635 return NULL; 636 } 637 638 num = (PympzObject*)Pympz_new(); 639 den = (PympzObject*)Pympz_new(); 640 if (!num || !den) { 641 Py_XDECREF((PyObject*)num); 642 Py_XDECREF((PyObject*)den); 643 return NULL; 644 } 645 646 if (mpfr_zero_p(Pympfr_AS_MPFR(self))) { 647 mpz_set_ui(num->z, 0); 648 mpz_set_ui(den->z, 1); 649 } 650 else { 651 temp = mpfr_get_z_2exp(num->z, Pympfr_AS_MPFR(self)); 652 twocount = (mpfr_exp_t)mpz_scan1(num->z, 0); 653 if (twocount) { 654 temp += twocount; 655 mpz_div_2exp(num->z, num->z, twocount); 656 } 657 mpz_set_ui(den->z, 1); 658 if (temp > 0) 659 mpz_mul_2exp(num->z, num->z, temp); 660 else if (temp < 0) 661 mpz_mul_2exp(den->z, den->z, -temp); 662 } 663 result = Py_BuildValue("(NN)", (PyObject*)num, (PyObject*)den); 664 if (!result) { 665 Py_DECREF((PyObject*)num); 666 Py_DECREF((PyObject*)den); 667 } 668 return result; 669} 670 671PyDoc_STRVAR(doc_mpfr_mantissa_exp, 672"x.as_mantissa_exp() -> (mantissa,exponent)\n\n" 673"Return the mantissa and exponent of an mpfr."); 674 675static PyObject * 676Pympfr_mantissa_exp(PyObject *self, PyObject *args) 677{ 678 PympzObject *mantissa = 0, *exponent = 0; 679 mpfr_exp_t temp; 680 PyObject *result; 681 682 if (mpfr_nan_p(Pympfr_AS_MPFR(self))) { 683 VALUE_ERROR("Cannot pass NaN to mpfr.as_mantissa_exp."); 684 return NULL; 685 } 686 687 if (mpfr_inf_p(Pympfr_AS_MPFR(self))) { 688 OVERFLOW_ERROR("Cannot pass Infinity to mpfr.as_mantissa_exp."); 689 return NULL; 690 } 691 692 mantissa = (PympzObject*)Pympz_new(); 693 exponent = (PympzObject*)Pympz_new(); 694 if (!mantissa || !exponent) { 695 Py_XDECREF((PyObject*)mantissa); 696 Py_XDECREF((PyObject*)exponent); 697 return NULL; 698 } 699 700 if (mpfr_zero_p(Pympfr_AS_MPFR(self))) { 701 mpz_set_ui(mantissa->z, 0); 702 mpz_set_ui(exponent->z, 1); 703 } 704 else { 705 temp = mpfr_get_z_2exp(mantissa->z, Pympfr_AS_MPFR(self)); 706 mpz_set_si(exponent->z, temp); 707 } 708 result = Py_BuildValue("(NN)", (PyObject*)mantissa, (PyObject*)exponent); 709 if (!result) { 710 Py_DECREF((PyObject*)mantissa); 711 Py_DECREF((PyObject*)exponent); 712 } 713 return result; 714} 715 716PyDoc_STRVAR(doc_mpfr_simple_fraction, 717"x.as_simple_fraction([precision=0]) -> mpq\n\n" 718"Return a simple rational approximation to x. The result will be\n" 719"accurate to 'precision' bits. If 'precision' is 0, the precision\n" 720"of 'x' will be used."); 721 722static PyObject * 723Pympfr_simple_fraction(PyObject *self, PyObject *args, PyObject *keywds) 724{ 725 mpfr_prec_t prec = 0; 726 static char *kwlist[] = {"precision", NULL}; 727 728 if (!PyArg_ParseTupleAndKeywords(args, keywds, "|l", kwlist, &prec)) 729 return NULL; 730 731 return (PyObject*)stern_brocot((PympfrObject*)self, 0, prec, 0); 732} 733 734static Py_hash_t 735_mpfr_hash(mpfr_t f) 736{ 737#ifdef _PyHASH_MODULUS 738 Py_uhash_t hash = 0; 739 Py_ssize_t exp; 740 size_t msize; 741 int sign; 742 743 /* Handle special cases first */ 744 if (!mpfr_number_p(f)) { 745 if (mpfr_inf_p(f)) 746 if (mpfr_sgn(f) > 0) 747 return _PyHASH_INF; 748 else 749 return -_PyHASH_INF; 750 else 751 return _PyHASH_NAN; 752 } 753 754 /* Calculate the number of limbs in the mantissa. */ 755 msize = (f->_mpfr_prec + mp_bits_per_limb - 1) / mp_bits_per_limb; 756 757 /* Calculate the hash of the mantissa. */ 758 if (mpfr_sgn(f) > 0) { 759 hash = mpn_mod_1(f->_mpfr_d, msize, _PyHASH_MODULUS); 760 sign = 1; 761 } 762 else if (mpfr_sgn(f) < 0) { 763 hash = mpn_mod_1(f->_mpfr_d, msize, _PyHASH_MODULUS); 764 sign = -1; 765 } 766 else { 767 return 0; 768 } 769 770 /* Calculate the final hash. */ 771 exp = f->_mpfr_exp - (msize * mp_bits_per_limb); 772 exp = exp >= 0 ? exp % _PyHASH_BITS : _PyHASH_BITS-1-((-1-exp) % _PyHASH_BITS); 773 hash = ((hash << exp) & _PyHASH_MODULUS) | hash >> (_PyHASH_BITS - exp); 774 775 hash *= sign; 776 if (hash == (Py_uhash_t)-1) 777 hash = (Py_uhash_t)-2; 778 return (Py_hash_t)hash; 779#else 780 double temp; 781 temp = mpfr_get_d(f, context->ctx.mpfr_round); 782 return _Py_HashDouble(temp); 783#endif 784} 785 786static Py_hash_t 787Pympfr_hash(PympfrObject *self) 788{ 789 if (self->hash_cache == -1) 790 self->hash_cache = _mpfr_hash(self->f); 791 return self->hash_cache; 792} 793 794/* This function is used in gmpy_mpany. */ 795 796static PyObject * 797Pympfr_pow(PyObject *base, PyObject *exp, PyObject *m) 798{ 799 PympfrObject *tempb, *tempe, *result; 800#ifdef WITHMPC 801 PympcObject *mpc_result; 802#endif 803 804 if (m != Py_None) { 805 TYPE_ERROR("pow() 3rd argument not allowed unless all arguments are integers"); 806 return NULL; 807 } 808 809 tempb = Pympfr_From_Real(base, 0); 810 tempe = Pympfr_From_Real(exp, 0); 811 812 if (!tempe || !tempb) { 813 Py_XDECREF((PyObject*)tempe); 814 Py_XDECREF((PyObject*)tempb); 815 Py_RETURN_NOTIMPLEMENTED; 816 } 817 818 result = (PympfrObject*)Pympfr_new(0); 819 820 if (!result) { 821 Py_DECREF((PyObject*)tempe); 822 Py_DECREF((PyObject*)tempb); 823 return NULL; 824 } 825 826 if (mpfr_zero_p(tempb->f) && (mpfr_sgn(tempe->f) < 0)) { 827 context->ctx.divzero = 1; 828 if (context->ctx.trap_divzero) { 829 GMPY_DIVZERO("zero cannot be raised to a negative power"); 830 goto done; 831 } 832 } 833 834 mpfr_clear_flags(); 835 result->rc = mpfr_pow(result->f, tempb->f, 836 tempe->f, context->ctx.mpfr_round); 837#ifdef WITHMPC 838 if (result && mpfr_nanflag_p() && context->ctx.allow_complex) { 839 /* If we don't get a valid result, or the result is a nan, then just 840 * return the original mpfr value. */ 841 if (!(mpc_result = (PympcObject*)Pympc_pow(base, exp, m)) || 842 MPC_IS_NAN_P(mpc_result)) { 843 844 Py_XDECREF((PyObject*)mpc_result); 845 context->ctx.invalid = 1; 846 GMPY_INVALID("invalid operation in 'mpfr' pow()"); 847 goto done; 848 } 849 /* return a valid complex result */ 850 Py_DECREF(result); 851 result = (PympfrObject*)mpc_result; 852 goto done; 853 } 854#endif 855 856 SUBNORMALIZE(result) 857 MERGE_FLAGS 858 CHECK_FLAGS("pow()") 859 done: 860 Py_DECREF((PyObject*)tempe); 861 Py_DECREF((PyObject*)tempb); 862 if (PyErr_Occurred()) { 863 Py_XDECREF((PyObject*)result); 864 result = NULL; 865 } 866 return (PyObject*)result; 867} 868 869#define MPFR_CONST(NAME) \ 870static PyObject * \ 871Pympfr_##NAME(PyObject *self, PyObject *args, PyObject *keywds) \ 872{ \ 873 PympfrObject *result; \ 874 mpfr_prec_t bits = 0; \ 875 static char *kwlist[] = {"precision", NULL}; \ 876 if (!PyArg_ParseTupleAndKeywords(args, keywds, "|l", kwlist, &bits)) return NULL; \ 877 if ((result = (PympfrObject*)Pympfr_new(bits))) { \ 878 mpfr_clear_flags(); \ 879 result->rc = mpfr_##NAME(result->f, context->ctx.mpfr_round); \ 880 MERGE_FLAGS \ 881 CHECK_FLAGS(#NAME "()") \ 882 } \ 883 done: \ 884 return (PyObject*)result; \ 885} 886 887PyDoc_STRVAR(doc_mpfr_const_pi, 888"const_pi([precision=0]) -> mpfr\n\n" 889"Return the constant pi using the specified precision. If no\n" 890"precision is specified, the default precision is used."); 891 892MPFR_CONST(const_pi) 893 894PyDoc_STRVAR(doc_mpfr_const_euler, 895"const_euler([precision=0]) -> mpfr\n\n" 896"Return the euler constant using the specified precision. If no\n" 897"precision is specified, the default precision is used."); 898 899MPFR_CONST(const_euler) 900 901PyDoc_STRVAR(doc_mpfr_const_log2, 902"const_log2([precision=0]) -> mpfr\n\n" 903"Return the log2 constant using the specified precision. If no\n" 904"precision is specified, the default precision is used."); 905 906MPFR_CONST(const_log2) 907 908PyDoc_STRVAR(doc_mpfr_const_catalan, 909"const_catalan([precision=0]) -> mpfr\n\n" 910"Return the catalan constant using the specified precision. If no\n" 911"precision is specified, the default precision is used."); 912 913MPFR_CONST(const_catalan) 914 915static PyObject * 916Pympfr_sqrt(PyObject *self, PyObject *other) 917{ 918 PympfrObject *result; 919 920 PARSE_ONE_MPFR_OTHER("sqrt() requires 'mpfr' argument"); 921 922#ifdef WITHMPC 923 if (mpfr_sgn(Pympfr_AS_MPFR(self)) < 0 && context->ctx.allow_complex) { 924 Py_DECREF(self); 925 return Pympc_sqrt(self, other); 926 } 927#endif 928 929 if (!(result = (PympfrObject*)Pympfr_new(0))) { 930 Py_DECREF(self); 931 return NULL; 932 } 933 934 mpfr_clear_flags(); 935 result->rc = mpfr_sqrt(result->f, Pympfr_AS_MPFR(self), 936 context->ctx.mpfr_round); 937 938 MPFR_CLEANUP_SELF("sqrt()"); 939} 940 941PyDoc_STRVAR(doc_g_mpfr_rec_sqrt, 942"rec_sqrt(x) -> mpfr\n\n" 943"Return the reciprocal of the square root of x."); 944 945static PyObject * 946Pympfr_rec_sqrt(PyObject *self, PyObject *other) 947{ 948 PympfrObject *result; 949 950 PARSE_ONE_MPFR_OTHER("rec_sqrt() requires 'mpfr' argument"); 951 952 if (!(result = (PympfrObject*)Pympfr_new(0))) 953 goto done; 954 955 mpfr_clear_flags(); 956 result->rc = mpfr_rec_sqrt(result->f, Pympfr_AS_MPFR(self), 957 context->ctx.mpfr_round); 958 959 MPFR_CLEANUP_SELF("rec_sqrt()"); 960} 961 962PyDoc_STRVAR(doc_mpfr_root, 963"root(x, n) -> mpfr\n\n" 964"Return n-th root of x. The result always an 'mpfr'."); 965 966static PyObject * 967Pympfr_root(PyObject *self, PyObject *args) 968{ 969 long n; 970 PympfrObject *result; 971 972 PARSE_ONE_MPFR_REQ_CLONG(&n, "root() requires 'mpfr','int' arguments"); 973 974 if (!(result = (PympfrObject*)Pympfr_new(0))) 975 goto done; 976 977 if (n <= 0) { 978 VALUE_ERROR("n must be > 0"); 979 goto done; 980 } 981 982 mpfr_clear_flags(); 983 result->rc = mpfr_root(result->f, Pympfr_AS_MPFR(self), n, 984 context->ctx.mpfr_round); 985 986 MPFR_CLEANUP_SELF("root()"); 987} 988 989PyDoc_STRVAR(doc_g_mpfr_round2, 990"round2(x[, n]) -> mpfr\n\n" 991"Return x rounded to n bits. Uses default precision if n is not\n" 992"specified. See round_away() to access the mpfr_round() function."); 993 994static PyObject * 995Pympfr_round2(PyObject *self, PyObject *args) 996{ 997 mpfr_prec_t prec = context->ctx.mpfr_prec; 998 PympfrObject *result = 0; 999 1000 PARSE_ONE_MPFR_OPT_CLONG(&prec, 1001 "round2() requires 'mpfr',['int'] arguments"); 1002 1003 if (prec < MPFR_PREC_MIN || prec > MPFR_PREC_MAX) { 1004 VALUE_ERROR("invalid precision"); 1005 goto done; 1006 } 1007 1008 if (!(result = (PympfrObject*)Pympfr_new(mpfr_get_prec(Pympfr_AS_MPFR(self))))) { 1009 goto done; 1010 } 1011 1012 mpfr_clear_flags(); 1013 /* Duplicate the code from Pympfr_pos. */ 1014 mpfr_set(result->f, Pympfr_AS_MPFR(self), context->ctx.mpfr_round); 1015 result->round_mode = ((PympfrObject*)self)->round_mode; 1016 result->rc = ((PympfrObject*)self)->rc; 1017 result->rc = mpfr_check_range(result->f, result->rc, result->round_mode); 1018 result->rc = mpfr_prec_round(result->f, prec, context->ctx.mpfr_round); 1019 1020 MPFR_CLEANUP_SELF("round2()"); 1021} 1022 1023PyDoc_STRVAR(doc_g_mpfr_round10, 1024"__round__(x[, n = 0]) -> mpfr\n\n" 1025"Return x rounded to n decimal digits before (n < 0) or after (n > 0)\n" 1026"the decimal point. Rounds to an integer if n is not specified."); 1027 1028static PyObject * 1029Pympfr_round10(PyObject *self, PyObject *args) 1030{ 1031 Py_ssize_t digits = 0; 1032 mpz_t temp; 1033 PympfrObject *resultf = 0; 1034 PympzObject *resultz; 1035 1036 /* If the size of args is 0, we just return an mpz. */ 1037 1038 if (PyTuple_GET_SIZE(args) == 0) { 1039 if ((resultz = (PympzObject*)Pympz_new())) { 1040 if (mpfr_nan_p(Pympfr_AS_MPFR(self))) { 1041 Py_DECREF((PyObject*)resultz); 1042 VALUE_ERROR("'mpz' does not support NaN"); 1043 return NULL; 1044 } 1045 if (mpfr_inf_p(Pympfr_AS_MPFR(self))) { 1046 Py_DECREF((PyObject*)resultz); 1047 OVERFLOW_ERROR("'mpz' does not support Infinity"); 1048 return NULL; 1049 } 1050 /* return code is ignored */ 1051 mpfr_get_z(resultz->z, Pympfr_AS_MPFR(self), MPFR_RNDN); 1052 } 1053 return (PyObject*)resultz; 1054 } 1055 1056 /* Now we need to return an mpfr, so handle the simple cases. */ 1057 1058 if (!mpfr_regular_p(Pympfr_AS_MPFR(self))) { 1059 Py_INCREF(self); 1060 return self; 1061 } 1062 1063 if (PyTuple_GET_SIZE(args) > 1) { 1064 TYPE_ERROR("Too many arguments for __round__()."); 1065 return NULL; 1066 } 1067 1068 if (PyTuple_GET_SIZE(args) == 1) { 1069 digits = ssize_t_From_Integer(PyTuple_GET_ITEM(args, 0)); 1070 if (digits == -1 && PyErr_Occurred()) { 1071 TYPE_ERROR("__round__() requires 'int' argument"); 1072 return NULL; 1073 } 1074 } 1075 1076 /* TODO: better error analysis, or else convert the mpfr to an exact 1077 * fraction, round the fraction, and then convert back to an mpfr. 1078 */ 1079 1080 resultf = (PympfrObject*)Pympfr_new(mpfr_get_prec(Pympfr_AS_MPFR(self))+100); 1081 if (!resultf) 1082 return NULL; 1083 1084 mpz_inoc(temp); 1085 mpz_ui_pow_ui(temp, 10, digits > 0 ? digits : -digits); 1086 if (digits >= 0) { 1087 mpfr_mul_z(resultf->f, Pympfr_AS_MPFR(self), temp, MPFR_RNDN); 1088 } 1089 else { 1090 mpfr_div_z(resultf->f, Pympfr_AS_MPFR(self), temp, MPFR_RNDN); 1091 } 1092 1093 mpfr_rint(resultf->f, resultf->f, MPFR_RNDN); 1094 1095 if (digits >= 0) { 1096 mpfr_div_z(resultf->f, resultf->f, temp, MPFR_RNDN); 1097 } 1098 else { 1099 mpfr_mul_z(resultf->f, resultf->f, temp, MPFR_RNDN); 1100 } 1101 mpfr_prec_round(resultf->f, mpfr_get_prec(Pympfr_AS_MPFR(self)), MPFR_RNDN); 1102 1103 mpz_cloc(temp); 1104 return((PyObject*)resultf); 1105} 1106 1107PyDoc_STRVAR(doc_g_mpfr_reldiff, 1108"reldiff(x, y) -> mpfr\n\n" 1109"Return the relative difference between x and y. Result is equal to\n" 1110"abs(x-y)/x."); 1111 1112static PyObject * 1113Pympfr_reldiff(PyObject *self, PyObject *args) 1114{ 1115 PympfrObject *result; 1116 PyObject *other; 1117 1118 PARSE_TWO_MPFR_ARGS(other, "reldiff() requires 'mpfr','mpfr' arguments"); 1119 1120 if (!(result = (PympfrObject*)Pympfr_new(0))) { 1121 Py_DECREF(self); 1122 Py_DECREF(other); 1123 return NULL; 1124 } 1125 1126 /* mpfr_reldiff doesn't guarantee correct rounding and doesn't appear 1127 * to set any exceptions. 1128 */ 1129 mpfr_reldiff(result->f, Pympfr_AS_MPFR(self), Pympfr_AS_MPFR(other), 1130 context->ctx.mpfr_round); 1131 result->rc = 0; 1132 Py_DECREF(self); 1133 Py_DECREF(other); 1134 return (PyObject*)result; 1135} 1136 1137static PyObject * 1138Pympfr_sign(PyObject *self, PyObject *other) 1139{ 1140 long sign; 1141 1142 PARSE_ONE_MPFR_OTHER("sign() requires 'mpfr' argument"); 1143 1144 mpfr_clear_flags(); 1145 sign = mpfr_sgn(Pympfr_AS_MPFR(self)); 1146 1147 MERGE_FLAGS; 1148 CHECK_ERANGE("range error in 'mpfr' sign(), NaN argument"); 1149 1150 done: 1151 Py_DECREF((PyObject*)self); 1152 if (PyErr_Occurred()) 1153 return NULL; 1154 else 1155 return PyIntOrLong_FromLong(sign); 1156} 1157 1158#define MPFR_MONOP(NAME) \ 1159static PyObject * \ 1160Py##NAME(PympfrObject *x) \ 1161{ \ 1162 PympfrObject *r; \ 1163 if (!(r = (PympfrObject*)Pympfr_new(0))) \ 1164 return NULL; \ 1165 if (Pympfr_CheckAndExp(x)) { \ 1166 r->rc = NAME(r->f, x->f, context->ctx.mpfr_round); \ 1167 } \ 1168 else { \ 1169 mpfr_set(r->f, x->f, context->ctx.mpfr_round); \ 1170 r->round_mode = x->round_mode; \ 1171 r->rc = x->rc; \ 1172 mpfr_clear_flags(); \ 1173 mpfr_check_range(r->f, r->rc, r->round_mode); \ 1174 r->rc = NAME(r->f, r->f, context->ctx.mpfr_round); \ 1175 MERGE_FLAGS; \ 1176 CHECK_FLAGS(#NAME "()"); \ 1177 } \ 1178 done: \ 1179 return (PyObject *) r; \ 1180} 1181 1182MPFR_MONOP(mpfr_abs) 1183MPFR_MONOP(mpfr_neg) 1184 1185#define MPFR_UNIOP_NOROUND(NAME) \ 1186static PyObject * \ 1187Pympfr_##NAME(PyObject* self, PyObject *other) \ 1188{ \ 1189 PympfrObject *result; \ 1190 PARSE_ONE_MPFR_OTHER(#NAME "() requires 'mpfr' argument"); \ 1191 if (!(result = (PympfrObject*)Pympfr_new(0))) goto done; \ 1192 mpfr_clear_flags(); \ 1193 result->rc = mpfr_##NAME(result->f, Pympfr_AS_MPFR(self)); \ 1194 MPFR_CLEANUP_SELF(#NAME "()"); \ 1195} 1196 1197PyDoc_STRVAR(doc_mpfr_ceil, 1198"x.__ceil__() -> mpfr\n\n" 1199"Return an 'mpfr' that is the smallest integer >= x."); 1200 1201PyDoc_STRVAR(doc_g_mpfr_ceil, 1202"ceil(x) ->mpfr\n\n" 1203"Return an 'mpfr' that is the smallest integer >= x."); 1204 1205MPFR_UNIOP_NOROUND(ceil) 1206 1207PyDoc_STRVAR(doc_mpfr_floor, 1208"x.__floor__() -> mpfr\n\n" 1209"Return an 'mpfr' that is the smallest integer <= x."); 1210 1211PyDoc_STRVAR(doc_g_mpfr_floor, 1212"floor(x) -> mpfr\n\n" 1213"Return an 'mpfr' that is the smallest integer <= x."); 1214 1215MPFR_UNIOP_NOROUND(floor); 1216 1217PyDoc_STRVAR(doc_mpfr_trunc, 1218"x.__trunc__() -> mpfr\n\n" 1219"Return an 'mpfr' that is truncated towards 0. Same as\n" 1220"x.floor() if x>=0 or x.ceil() if x<0."); 1221 1222PyDoc_STRVAR(doc_g_mpfr_trunc, 1223"trunc(x) -> mpfr\n\n" 1224"Return an 'mpfr' that is x truncated towards 0. Same as\n" 1225"x.floor() if x>=0 or x.ceil() if x<0."); 1226 1227MPFR_UNIOP_NOROUND(trunc) 1228 1229PyDoc_STRVAR(doc_g_mpfr_round_away, 1230"round_away(x) -> mpfr\n\n" 1231"Return an 'mpfr' that is x rounded to the nearest integer,\n" 1232"with ties rounded away from 0."); 1233 1234static PyObject * 1235Pympfr_round_away(PyObject* self, PyObject *other) 1236{ 1237 PympfrObject *result; 1238 1239 PARSE_ONE_MPFR_OTHER("round_away() requires 'mpfr' argument"); 1240 if (!(result = (PympfrObject*)Pympfr_new(0))) 1241 goto done; 1242 mpfr_clear_flags(); 1243 result->rc = mpfr_round(result->f, Pympfr_AS_MPFR(self)); 1244 MPFR_CLEANUP_SELF("round_away()"); 1245} 1246 1247#define MPFR_UNIOP(NAME) \ 1248static PyObject * \ 1249Pympfr_##NAME(PyObject* self, PyObject *other) \ 1250{ \ 1251 PympfrObject *result; \ 1252 PARSE_ONE_MPFR_OTHER(#NAME "() requires 'mpfr' argument"); \ 1253 if (!(result = (PympfrObject*)Pympfr_new(0))) goto done; \ 1254 mpfr_clear_flags(); \ 1255 result->rc = mpfr_##NAME(result->f, Pympfr_AS_MPFR(self), context->ctx.mpfr_round); \ 1256 MPFR_CLEANUP_SELF(#NAME "()"); \ 1257} 1258 1259PyDoc_STRVAR(doc_g_mpfr_rint, 1260"rint(x) -> mpfr\n\n" 1261"Return x rounded to the nearest integer using the current rounding\n" 1262"mode."); 1263 1264MPFR_UNIOP(rint) 1265 1266PyDoc_STRVAR(doc_g_mpfr_rint_ceil, 1267"rint_ceil(x) -> mpfr\n\n" 1268"Return x rounded to the nearest integer by first rounding to the\n" 1269"next higher or equal integer and then, if needed, using the current\n" 1270"rounding mode."); 1271 1272MPFR_UNIOP(rint_ceil) 1273 1274PyDoc_STRVAR(doc_g_mpfr_rint_floor, 1275"rint_floor(x) -> mpfr\n\n" 1276"Return x rounded to the nearest integer by first rounding to the\n" 1277"next lower or equal integer and then, if needed, using the current\n" 1278"rounding mode."); 1279 1280MPFR_UNIOP(rint_floor) 1281 1282PyDoc_STRVAR(doc_g_mpfr_rint_round, 1283"rint_round(x) -> mpfr\n\n" 1284"Return x rounded to the nearest integer by first rounding to the\n" 1285"nearest integer (ties away from 0) and then, if needed, using\n" 1286"the current rounding mode."); 1287 1288MPFR_UNIOP(rint_round) 1289 1290PyDoc_STRVAR(doc_g_mpfr_rint_trunc, 1291"rint_trunc(x) -> mpfr\n\n" 1292"Return x rounded to the nearest integer by first rounding towards\n" 1293"zero and then, if needed, using the current rounding mode."); 1294 1295MPFR_UNIOP(rint_trunc) 1296 1297PyDoc_STRVAR(doc_g_mpfr_frac, 1298"frac(x) -> mpfr\n\n" 1299"Return fractional part of x."); 1300 1301MPFR_UNIOP(frac) 1302 1303PyDoc_STRVAR(doc_g_mpfr_modf, 1304"modf(x) -> (mpfr, mpfr)\n\n" 1305"Return a tuple containing the integer and fractional portions\n" 1306"of x."); 1307 1308static PyObject * 1309Pympfr_modf(PyObject *self, PyObject *other) 1310{ 1311 PympfrObject *s, *c; 1312 PyObject *result; 1313 int code; 1314 1315 PARSE_ONE_MPFR_OTHER("modf() requires 'mpfr' argument"); 1316 1317 s = (PympfrObject*)Pympfr_new(0); 1318 c = (PympfrObject*)Pympfr_new(0); 1319 result = PyTuple_New(2); 1320 if (!s || !c || !result) 1321 goto done; 1322 1323 mpfr_clear_flags(); 1324 code = mpfr_modf(s->f, c->f, Pympfr_AS_MPFR(self), 1325 context->ctx.mpfr_round); 1326 s->rc = code & 0x03; 1327 c->rc = code >> 2; 1328 if (s->rc == 2) s->rc = -1; 1329 if (c->rc == 2) c->rc = -1; 1330 SUBNORMALIZE(s); 1331 SUBNORMALIZE(c); 1332 MERGE_FLAGS; 1333 CHECK_FLAGS("modf()"); 1334 1335 done: 1336 Py_DECREF(self); 1337 if (PyErr_Occurred()) { 1338 Py_XDECREF((PyObject*)s); 1339 Py_XDECREF((PyObject*)c); 1340 Py_XDECREF(result); 1341 result = NULL; 1342 } 1343 else { 1344 PyTuple_SET_ITEM(result, 0, (PyObject*)s); 1345 PyTuple_SET_ITEM(result, 1, (PyObject*)c); 1346 } 1347 return result; 1348} 1349 1350/* Needed for square() in mpz_mpany.c */ 1351 1352MPFR_UNIOP(sqr) 1353 1354PyDoc_STRVAR(doc_g_mpfr_cbrt, 1355"cbrt(x) -> mpfr\n\n" 1356"Return the cube root of x."); 1357 1358MPFR_UNIOP(cbrt) 1359 1360/* Called via gmpy_mpany so doc-string is there. */ 1361 1362MPFR_UNIOP(log) 1363 1364PyDoc_STRVAR(doc_g_mpfr_log2, 1365"log2(x) -> mpfr\n\n" 1366"Return base-2 logarithm of x."); 1367 1368MPFR_UNIOP(log2) 1369 1370/* Called via gmpy_mpany so doc-string is there. */ 1371 1372MPFR_UNIOP(log10) 1373 1374/* Called via gmpy_mpany so doc-string is there. */ 1375 1376MPFR_UNIOP(exp) 1377 1378PyDoc_STRVAR(doc_g_mpfr_exp2, 1379"exp2(x) -> mpfr\n\n" 1380"Return 2**x."); 1381 1382MPFR_UNIOP(exp2) 1383 1384PyDoc_STRVAR(doc_g_mpfr_exp10, 1385"exp10(x) -> mpfr\n\n" 1386"Return 10**x."); 1387 1388MPFR_UNIOP(exp10) 1389 1390MPFR_UNIOP(sin) 1391 1392MPFR_UNIOP(cos) 1393 1394MPFR_UNIOP(tan) 1395 1396PyDoc_STRVAR(doc_g_mpfr_sec, 1397"sec(x) -> mpfr\n\n" 1398"Return secant of x; x in radians."); 1399 1400MPFR_UNIOP(sec) 1401 1402PyDoc_STRVAR(doc_g_mpfr_csc, 1403"csc(x) -> mpfr\n\n" 1404"Return cosecant of x; x in radians."); 1405 1406MPFR_UNIOP(csc) 1407 1408PyDoc_STRVAR(doc_g_mpfr_cot, 1409"cot(x) -> mpfr\n\n" 1410"Return cotangent of x; x in radians."); 1411 1412MPFR_UNIOP(cot) 1413 1414static PyObject * 1415Pympfr_acos(PyObject* self, PyObject *other) 1416{ 1417 PympfrObject *result; 1418 1419 PARSE_ONE_MPFR_OTHER("acos() requires 'mpfr' argument"); 1420 1421#ifdef WITHMPC 1422 if (!mpfr_nan_p(Pympfr_AS_MPFR(self)) && 1423 (mpfr_cmp_si(Pympfr_AS_MPFR(self), 1) > 0 || 1424 mpfr_cmp_si(Pympfr_AS_MPFR(self), -1) < 0) && 1425 context->ctx.allow_complex) { 1426 Py_DECREF(self); 1427 return Pympc_acos(self, other); 1428 } 1429#endif 1430 1431 if (!(result = (PympfrObject*)Pympfr_new(0))) { 1432 Py_DECREF(self); 1433 return NULL; 1434 } 1435 mpfr_clear_flags(); 1436 result->rc = mpfr_acos(result->f, Pympfr_AS_MPFR(self), 1437 context->ctx.mpfr_round); 1438 MPFR_CLEANUP_SELF("acos()"); 1439} 1440 1441static PyObject * 1442Pympfr_asin(PyObject* self, PyObject *other) 1443{ 1444 PympfrObject *result; 1445 1446 PARSE_ONE_MPFR_OTHER("asin() requires 'mpfr' argument"); 1447 1448#ifdef WITHMPC 1449 if (!mpfr_nan_p(Pympfr_AS_MPFR(self)) && 1450 (mpfr_cmp_si(Pympfr_AS_MPFR(self), 1) > 0 || 1451 mpfr_cmp_si(Pympfr_AS_MPFR(self), -1) < 0) && 1452 context->ctx.allow_complex) { 1453 Py_DECREF(self); 1454 return Pympc_asin(self, other); 1455 } 1456#endif 1457 1458 if (!(result = (PympfrObject*)Pympfr_new(0))) { 1459 Py_DECREF(self); 1460 return NULL; 1461 } 1462 mpfr_clear_flags(); 1463 result->rc = mpfr_asin(result->f, Pympfr_AS_MPFR(self), 1464 context->ctx.mpfr_round); 1465 MPFR_CLEANUP_SELF("asin()"); 1466} 1467 1468MPFR_UNIOP(atan) 1469 1470MPFR_UNIOP(cosh) 1471 1472MPFR_UNIOP(sinh) 1473 1474MPFR_UNIOP(tanh) 1475 1476PyDoc_STRVAR(doc_g_mpfr_sech, 1477"sech(x) -> mpfr\n\n" 1478"Returns hyperbolic secant of x."); 1479 1480MPFR_UNIOP(sech) 1481 1482PyDoc_STRVAR(doc_g_mpfr_csch, 1483"csch(x) -> mpfr\n\n" 1484"Return hyperbolic cosecant of x."); 1485 1486MPFR_UNIOP(csch) 1487 1488PyDoc_STRVAR(doc_g_mpfr_coth, 1489"coth(x) -> mpfr\n\n" 1490"Return hyperbolic cotangent of x."); 1491 1492MPFR_UNIOP(coth) 1493 1494MPFR_UNIOP(acosh) 1495 1496MPFR_UNIOP(asinh) 1497 1498static PyObject * 1499Pympfr_atanh(PyObject* self, PyObject *other) 1500{ 1501 PympfrObject *result; 1502 1503 PARSE_ONE_MPFR_OTHER("atanh() requires 'mpfr' argument"); 1504 1505#ifdef WITHMPC 1506 if (!mpfr_nan_p(Pympfr_AS_MPFR(self)) && 1507 (mpfr_cmp_si(Pympfr_AS_MPFR(self), 1) > 0 || 1508 mpfr_cmp_si(Pympfr_AS_MPFR(self), -1) < 0) && 1509 context->ctx.allow_complex) { 1510 Py_DECREF(self); 1511 return Pympc_atanh(self, other); 1512 } 1513#endif 1514 1515 if (!(result = (PympfrObject*)Pympfr_new(0))) { 1516 Py_DECREF(self); 1517 return NULL; 1518 } 1519 mpfr_clear_flags(); 1520 result->rc = mpfr_asin(result->f, Pympfr_AS_MPFR(self), 1521 context->ctx.mpfr_round); 1522 MPFR_CLEANUP_SELF("atanh()"); 1523} 1524 1525 1526PyDoc_STRVAR(doc_g_mpfr_log1p, 1527"log1p(x) -> mpfr\n\n" 1528"Return logarithm of (1+x)."); 1529 1530MPFR_UNIOP(log1p) 1531 1532PyDoc_STRVAR(doc_g_mpfr_expm1, 1533"expm1(x) -> mpfr\n\n" 1534"Return exponential(x) - 1."); 1535 1536MPFR_UNIOP(expm1) 1537 1538PyDoc_STRVAR(doc_g_mpfr_eint, 1539"eint(x) -> mpfr\n\n" 1540"Return exponential integral of x."); 1541 1542MPFR_UNIOP(eint) 1543 1544PyDoc_STRVAR(doc_g_mpfr_li2, 1545"li2(x) -> mpfr\n\n" 1546"Return real part of dilogarithm of x."); 1547 1548MPFR_UNIOP(li2) 1549 1550PyDoc_STRVAR(doc_g_mpfr_gamma, 1551"gamma(x) -> mpfr\n\n" 1552"Return gamma of x."); 1553 1554MPFR_UNIOP(gamma) 1555 1556PyDoc_STRVAR(doc_g_mpfr_lngamma, 1557"lngamma(x) -> mpfr\n\n" 1558"Return logarithm of gamma(x)."); 1559 1560MPFR_UNIOP(lngamma) 1561 1562PyDoc_STRVAR(doc_g_mpfr_lgamma, 1563"lgamma(x) -> (mpfr, int)\n\n" 1564"Return a tuple containing the logarithm of the absolute value of\n" 1565"gamma(x) and the sign of gamma(x)"); 1566 1567static PyObject * 1568Pympfr_lgamma(PyObject* self, PyObject *other) 1569{ 1570 PyObject *result; 1571 PympfrObject *value; 1572 int signp = 0; 1573 1574 PARSE_ONE_MPFR_OTHER("lgamma() requires 'mpfr' argument"); 1575 1576 value = (PympfrObject*)Pympfr_new(0); 1577 result = PyTuple_New(2); 1578 if (!value || !result) 1579 goto done; 1580 1581 mpfr_clear_flags(); 1582 value->rc = mpfr_lgamma(value->f, &signp, Pympfr_AS_MPFR(self), 1583 context->ctx.mpfr_round); 1584 SUBNORMALIZE(value); 1585 MERGE_FLAGS; 1586 CHECK_FLAGS("lgamma()"); 1587 1588 done: 1589 Py_DECREF(self); 1590 if (PyErr_Occurred()) { 1591 Py_XDECREF(result); 1592 Py_XDECREF((PyObject*)value); 1593 result = NULL; 1594 } 1595 else { 1596 PyTuple_SET_ITEM(result, 0, (PyObject*)value); 1597 PyTuple_SET_ITEM(result, 1, PyIntOrLong_FromLong((long)signp)); 1598 } 1599 return result; 1600} 1601 1602PyDoc_STRVAR(doc_g_mpfr_digamma, 1603"digamma(x) -> mpfr\n\n" 1604"Return digamma of x."); 1605 1606MPFR_UNIOP(digamma) 1607 1608PyDoc_STRVAR(doc_g_mpfr_zeta, 1609"zeta(x) -> mpfr\n\n" 1610"Return Riemann zeta of x."); 1611 1612MPFR_UNIOP(zeta) 1613 1614PyDoc_STRVAR(doc_g_mpfr_erf, 1615"erf(x) -> mpfr\n\n" 1616"Return error function of x."); 1617 1618MPFR_UNIOP(erf) 1619 1620PyDoc_STRVAR(doc_g_mpfr_erfc, 1621"erfc(x) -> mpfr\n\n" 1622"Return complementary error function of x."); 1623 1624MPFR_UNIOP(erfc) 1625 1626PyDoc_STRVAR(doc_g_mpfr_j0, 1627"j0(x) -> mpfr\n\n" 1628"Return first kind Bessel function of order 0 of x."); 1629 1630MPFR_UNIOP(j0) 1631 1632PyDoc_STRVAR(doc_g_mpfr_j1, 1633"j1(x) -> mpfr\n\n" 1634"Return first kind Bessel function of order 1 of x."); 1635 1636MPFR_UNIOP(j1) 1637 1638PyDoc_STRVAR(doc_g_mpfr_jn, 1639"jn(x,n) -> mpfr\n\n" 1640"Return the first kind Bessel function of order n of x."); 1641 1642static PyObject * 1643Pympfr_jn(PyObject *self, PyObject *args) 1644{ 1645 PympfrObject *result; 1646 long n = 0; 1647 1648 PARSE_ONE_MPFR_REQ_CLONG(&n, "jn() requires 'mpfr','int' arguments"); 1649 1650 if (!(result = (PympfrObject*)Pympfr_new(0))) 1651 goto done; 1652 1653 mpfr_clear_flags(); 1654 result->rc = mpfr_jn(result->f, n, Pympfr_AS_MPFR(self), 1655 context->ctx.mpfr_round); 1656 MPFR_CLEANUP_SELF("jn()"); 1657} 1658 1659PyDoc_STRVAR(doc_g_mpfr_y0, 1660"y0(x) -> mpfr\n\n" 1661"Return second kind Bessel function of order 0 of x."); 1662 1663MPFR_UNIOP(y0) 1664 1665PyDoc_STRVAR(doc_g_mpfr_y1, 1666"y1(x) -> mpfr\n\n" 1667"Return second kind Bessel function of order 1 of x."); 1668 1669MPFR_UNIOP(y1) 1670 1671PyDoc_STRVAR(doc_g_mpfr_yn, 1672"yn(x,n) -> mpfr\n\n" 1673"Return the second kind Bessel function of order n of x."); 1674 1675static PyObject * 1676Pympfr_yn(PyObject *self, PyObject *args) 1677{ 1678 PympfrObject *result; 1679 long n = 0; 1680 1681 PARSE_ONE_MPFR_REQ_CLONG(&n, "yn() requires 'mpfr','int' arguments"); 1682 1683 if (!(result = (PympfrObject*)Pympfr_new(0))) 1684 goto done; 1685 1686 mpfr_clear_flags(); 1687 result->rc = mpfr_yn(result->f, n, Pympfr_AS_MPFR(self), 1688 context->ctx.mpfr_round); 1689 MPFR_CLEANUP_SELF("yn()"); 1690} 1691 1692PyDoc_STRVAR(doc_g_mpfr_ai, 1693"ai(x) -> mpfr\n\n" 1694"Return Airy function of x."); 1695 1696MPFR_UNIOP(ai) 1697 1698static PyObject * 1699Pympfr_add_fast(PyObject *x, PyObject *y) 1700{ 1701 PympfrObject *result; 1702 1703 if (Pympfr_CheckAndExp(x) && Pympfr_CheckAndExp(y)) { 1704 if (!(result = (PympfrObject*)Pympfr_new(0))) { 1705 return NULL; 1706 } 1707 mpfr_clear_flags(); 1708 result->rc = mpfr_add(result->f, 1709 Pympfr_AS_MPFR(x), 1710 Pympfr_AS_MPFR(y), 1711 context->ctx.mpfr_round); 1712 MPFR_CLEANUP_RESULT("addition"); 1713 return (PyObject*)result; 1714 } 1715 else { 1716 return Pybasic_add(x, y); 1717 } 1718} 1719 1720static PyObject * 1721Pympfr_add(PyObject *self, PyObject *args) 1722{ 1723 PympfrObject *result; 1724 PyObject *other; 1725 1726 PARSE_TWO_MPFR_ARGS(other, "add() requires 'mpfr','mpfr' arguments"); 1727 1728 if (!(result = (PympfrObject*)Pympfr_new(0))) 1729 goto done; 1730 1731 mpfr_clear_flags(); 1732 result->rc = mpfr_add(result->f, Pympfr_AS_MPFR(self), 1733 Pympfr_AS_MPFR(other), context->ctx.mpfr_round); 1734 MPFR_CLEANUP_SELF_OTHER("add()"); 1735} 1736 1737static PyObject * 1738Pympfr_sub_fast(PyObject *x, PyObject *y) 1739{ 1740 PympfrObject *result; 1741 1742 if (Pympfr_CheckAndExp(x) && Pympfr_CheckAndExp(y)) { 1743 if (!(result = (PympfrObject*)Pympfr_new(0))) { 1744 return NULL; 1745 } 1746 mpfr_clear_flags(); 1747 result->rc = mpfr_sub(result->f, 1748 Pympfr_AS_MPFR(x), 1749 Pympfr_AS_MPFR(y), 1750 context->ctx.mpfr_round); 1751 MPFR_CLEANUP_RESULT("subtraction"); 1752 return (PyObject*)result; 1753 } 1754 else { 1755 return Pybasic_sub(x, y); 1756 } 1757} 1758 1759static PyObject * 1760Pympfr_sub(PyObject *self, PyObject *args) 1761{ 1762 PympfrObject *result; 1763 PyObject *other; 1764 1765 PARSE_TWO_MPFR_ARGS(other, "sub() requires 'mpfr','mpfr' arguments"); 1766 1767 if (!(result = (PympfrObject*)Pympfr_new(0))) 1768 goto done; 1769 1770 mpfr_clear_flags(); 1771 result->rc = mpfr_sub(result->f, Pympfr_AS_MPFR(self), 1772 Pympfr_AS_MPFR(other), context->ctx.mpfr_round); 1773 MPFR_CLEANUP_SELF_OTHER("sub()"); 1774} 1775 1776static PyObject * 1777Pympfr_mul_fast(PyObject *x, PyObject *y) 1778{ 1779 PympfrObject *result; 1780 1781 if (Pympfr_CheckAndExp(x) && Pympfr_CheckAndExp(y)) { 1782 if (!(result = (PympfrObject*)Pympfr_new(0))) { 1783 return NULL; 1784 } 1785 mpfr_clear_flags(); 1786 result->rc = mpfr_mul(result->f, 1787 Pympfr_AS_MPFR(x), 1788 Pympfr_AS_MPFR(y), 1789 context->ctx.mpfr_round); 1790 MPFR_CLEANUP_RESULT("multiplication"); 1791 return (PyObject*)result; 1792 } 1793 else { 1794 return Pybasic_mul(x, y); 1795 } 1796} 1797 1798static PyObject * 1799Pympfr_mul(PyObject *self, PyObject *args) 1800{ 1801 PympfrObject *result; 1802 PyObject *other; 1803 1804 PARSE_TWO_MPFR_ARGS(other, "mul() requires 'mpfr','mpfr' arguments"); 1805 1806 if (!(result = (PympfrObject*)Pympfr_new(0))) 1807 goto done; 1808 1809 mpfr_clear_flags(); 1810 result->rc = mpfr_mul(result->f, Pympfr_AS_MPFR(self), 1811 Pympfr_AS_MPFR(other), context->ctx.mpfr_round); 1812 MPFR_CLEANUP_SELF_OTHER("mul()"); 1813} 1814 1815static PyObject * 1816Pympfr_truediv_fast(PyObject *x, PyObject *y) 1817{ 1818 PympfrObject *result; 1819 1820 if (Pympfr_Check…
Large files files are truncated, but you can click here to view the full file