cpython-lru_cache-weakref: Modules/

annotate Modules/_math.c @ 87690:66c1c9f32567 default tip

Add use_weakref option to functools.lru_cache The use_weakref option makes lru_cache create weak references to the function arguments when creating the cache key, which prevents them from staying in the memory forever even after all other references to them are gone. This is a proof of concept for http://bugs.python.org/issue19859

author	Radomir Dopieralski <devel@sheep.art.pl>
date	Tue, 03 Dec 2013 16:02:48 +0100
parents	08418369da7b
children

rev	line source
57758 b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	1 /* Definitions of some C99 math library functions, for those platforms
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	2 that don't implement these functions already. */
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	3
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	4 #include "Python.h"
57758 b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	5 #include <float.h>
57869 1a25b420ec00 Merged revisions 76982 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57866 diff changeset	6 #include "_math.h"
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	7
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	8 /* The following copyright notice applies to the original
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	9 implementations of acosh, asinh and atanh. */
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	10
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	11 /*
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	12 * ====================================================
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	13 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	14 *
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	15 * Developed at SunPro, a Sun Microsystems, Inc. business.
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	16 * Permission to use, copy, modify, and distribute this
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	17 * software is freely granted, provided that this notice
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	18 * is preserved.
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	19 * ====================================================
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	20 */
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	21
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	22 static const double ln2 = 6.93147180559945286227E-01;
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	23 static const double two_pow_m28 = 3.7252902984619141E-09; /* 2*-28 /
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	24 static const double two_pow_p28 = 268435456.0; /* 2*28 /
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	25 static const double zero = 0.0;
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	26
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	27 /* acosh(x)
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	28 * Method :
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	29 * Based on
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	30 * acosh(x) = log [ x + sqrt(x*x-1) ]
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	31 * we have
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	32 * acosh(x) := log(x)+ln2, if x is large; else
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	33 * acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	34 * acosh(x) := log1p(t+sqrt(2.0t+tt)); where t=x-1.
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	35 *
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	36 * Special cases:
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	37 * acosh(x) is NaN with signal if x<1.
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	38 * acosh(NaN) is NaN without signal.
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	39 */
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	40
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	41 double
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	42 _Py_acosh(double x)
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	43 {
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	44 if (Py_IS_NAN(x)) {
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	45 return x+x;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	46 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	47 if (x < 1.) { /* x < 1; return a signaling NaN */
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	48 errno = EDOM;
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	49 #ifdef Py_NAN
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	50 return Py_NAN;
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	51 #else
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	52 return (x-x)/(x-x);
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	53 #endif
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	54 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	55 else if (x >= two_pow_p28) { /* x > 2*28 /
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	56 if (Py_IS_INFINITY(x)) {
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	57 return x+x;
62830 2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	58 }
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	59 else {
62820 a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	60 return log(x)+ln2; /* acosh(huge)=log(2x) */
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	61 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	62 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	63 else if (x == 1.) {
62820 a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	64 return 0.0; /* acosh(1) = 0 */
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	65 }
62820 a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	66 else if (x > 2.) { /* 2 < x < 2*28 /
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	67 double t = x*x;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	68 return log(2.0*x - 1.0 / (x + sqrt(t - 1.0)));
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	69 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	70 else { /* 1 < x <= 2 */
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	71 double t = x - 1.0;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	72 return m_log1p(t + sqrt(2.0t + tt));
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	73 }
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	74 }
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	75
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	76
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	77 /* asinh(x)
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	78 * Method :
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	79 * Based on
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	80 * asinh(x) = sign(x) * log [ \|x\| + sqrt(x*x+1) ]
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	81 * we have
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	82 * asinh(x) := x if 1+x*x=1,
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	83 * := sign(x)*(log(x)+ln2)) for large \|x\|, else
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	84 * := sign(x)log(2\|x\|+1/(\|x\|+sqrt(xx+1))) if\|x\|>2, else
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	85 * := sign(x)*log1p(\|x\| + x^2/(1 + sqrt(1+x^2)))
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	86 */
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	87
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	88 double
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	89 _Py_asinh(double x)
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	90 {
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	91 double w;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	92 double absx = fabs(x);
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	93
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	94 if (Py_IS_NAN(x) \|\| Py_IS_INFINITY(x)) {
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	95 return x+x;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	96 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	97 if (absx < two_pow_m28) { /* \|x\| < 2*-28 /
62820 a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	98 return x; /* return x inexact except 0 */
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	99 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	100 if (absx > two_pow_p28) { /* \|x\| > 2*28 /
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	101 w = log(absx)+ln2;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	102 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	103 else if (absx > 2.0) { /* 2 < \|x\| < 2*28 /
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	104 w = log(2.0absx + 1.0 / (sqrt(xx + 1.0) + absx));
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	105 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	106 else { /* 2*-28 <= \|x\| < 2= /
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	107 double t = x*x;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	108 w = m_log1p(absx + t / (1.0 + sqrt(1.0 + t)));
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	109 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	110 return copysign(w, x);
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	111
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	112 }
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	113
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	114 /* atanh(x)
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	115 * Method :
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	116 * 1.Reduced x to positive by atanh(-x) = -atanh(x)
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	117 * 2.For x>=0.5
62820 a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	118 * 1 2x x
a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	119 * atanh(x) = --- * log(1 + -------) = 0.5 * log1p(2 * -------)
a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	120 * 2 1 - x 1 - x
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	121 *
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	122 * For x<0.5
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	123 * atanh(x) = 0.5log1p(2x+2xx/(1-x))
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	124 *
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	125 * Special cases:
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	126 * atanh(x) is NaN if \|x\| >= 1 with signal;
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	127 * atanh(NaN) is that NaN with no signal;
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	128 *
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	129 */
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	130
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	131 double
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	132 _Py_atanh(double x)
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	133 {
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	134 double absx;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	135 double t;
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	136
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	137 if (Py_IS_NAN(x)) {
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	138 return x+x;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	139 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	140 absx = fabs(x);
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	141 if (absx >= 1.) { /* \|x\| >= 1 */
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	142 errno = EDOM;
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	143 #ifdef Py_NAN
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	144 return Py_NAN;
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	145 #else
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	146 return x/zero;
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	147 #endif
61463 fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	148 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	149 if (absx < two_pow_m28) { /* \|x\| < 2*-28 /
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	150 return x;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	151 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	152 if (absx < 0.5) { /* \|x\| < 0.5 */
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	153 t = absx+absx;
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	154 t = 0.5 * m_log1p(t + t*absx / (1.0 - absx));
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	155 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	156 else { /* 0.5 <= \|x\| <= 1.0 */
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	157 t = 0.5 * m_log1p((absx + absx) / (1.0 - absx));
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	158 }
fd25030223b0 Untabify C files. Will watch buildbots. Antoine Pitrou <solipsis@pitrou.net> parents: 57869 diff changeset	159 return copysign(t, x);
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	160 }
57758 b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	161
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	162 /* Mathematically, expm1(x) = exp(x) - 1. The expm1 function is designed
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	163 to avoid the significant loss of precision that arises from direct
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	164 evaluation of the expression exp(x) - 1, for x near 0. */
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	165
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	166 double
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	167 _Py_expm1(double x)
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	168 {
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	169 /* For abs(x) >= log(2), it's safe to evaluate exp(x) - 1 directly; this
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	170 also works fine for infinities and nans.
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	171
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	172 For smaller x, we can use a method due to Kahan that achieves close to
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	173 full accuracy.
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	174 */
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	175
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	176 if (fabs(x) < 0.7) {
62830 2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	177 double u;
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	178 u = exp(x);
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	179 if (u == 1.0)
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	180 return x;
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	181 else
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	182 return (u - 1.0) * x / log(u);
57758 b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	183 }
b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	184 else
62830 2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	185 return exp(x) - 1.0;
57758 b45ca11b9a3e Merged revisions 76861 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: diff changeset	186 }
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	187
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	188 /* log1p(x) = log(1+x). The log1p function is designed to avoid the
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	189 significant loss of precision that arises from direct evaluation when x is
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	190 small. */
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	191
78640 08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	192 #ifdef HAVE_LOG1P
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	193
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	194 double
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	195 _Py_log1p(double x)
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	196 {
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	197 /* Some platforms supply a log1p function but don't respect the sign of
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	198 zero: log1p(-0.0) gives 0.0 instead of the correct result of -0.0.
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	199
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	200 To save fiddling with configure tests and platform checks, we handle the
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	201 special case of zero input directly on all platforms.
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	202 */
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	203 if (x == 0.0) {
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	204 return x;
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	205 }
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	206 else {
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	207 return log1p(x);
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	208 }
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	209 }
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	210
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	211 #else
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	212
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	213 double
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	214 _Py_log1p(double x)
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	215 {
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	216 /* For x small, we use the following approach. Let y be the nearest float
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	217 to 1+x, then
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	218
62820 a07ef5eab90b Merged revisions 82595 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 61463 diff changeset	219 1+x = y * (1 - (y-1-x)/y)
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	220
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	221 so log(1+x) = log(y) + log(1-(y-1-x)/y). Since (y-1-x)/y is tiny, the
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	222 second term is well approximated by (y-1-x)/y. If abs(x) >=
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	223 DBL_EPSILON/2 or the rounding-mode is some form of round-to-nearest
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	224 then y-1-x will be exactly representable, and is computed exactly by
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	225 (y-1)-x.
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	226
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	227 If abs(x) < DBL_EPSILON/2 and the rounding mode is not known to be
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	228 round-to-nearest then this method is slightly dangerous: 1+x could be
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	229 rounded up to 1+DBL_EPSILON instead of down to 1, and in that case
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	230 y-1-x will not be exactly representable any more and the result can be
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	231 off by many ulps. But this is easily fixed: for a floating-point
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	232 number \|x\| < DBL_EPSILON/2., the closest floating-point number to
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	233 log(1+x) is exactly x.
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	234 */
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	235
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	236 double y;
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	237 if (fabs(x) < DBL_EPSILON/2.) {
62830 2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	238 return x;
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	239 }
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	240 else if (-0.5 <= x && x <= 1.) {
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	241 /* WARNING: it's possible than an overeager compiler
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	242 will incorrectly optimize the following two lines
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	243 to the equivalent of "return log(1.+x)". If this
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	244 happens, then results from log1p will be inaccurate
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	245 for small x. */
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	246 y = 1.+x;
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	247 return log(y)-((y-1.)-x)/y;
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	248 }
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	249 else {
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	250 /* NaNs and infinities should end up here */
2399e0451d0d Merged revisions 82607 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 62820 diff changeset	251 return log(1.+x);
57866 9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	252 }
9751ddcc1882 Merged revisions 76978 via svnmerge from Mark Dickinson <dickinsm@gmail.com> parents: 57758 diff changeset	253 }
78640 08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	254
08418369da7b Issue #15477: Add workaround for log1p(-0.0) on platforms where it's broken. Mark Dickinson <mdickinson@enthought.com> parents: 62830 diff changeset	255 #endif /* ifdef HAVE_LOG1P */

Mercurial > cpython-lru_cache-weakref

annotate Modules/_math.c @ 87690:66c1c9f32567 default tip