annotate Modules/_math.c @ 73323:dad02a080bbc

Improve byte-compilation in packaging to be independent of -O or -B. The code I fixed to comply with PEP 3147 still had one bug: When run under python -O, some paths for pyc files would be pyo, because I called imp.cache_from_source without explicit debug_override argument in some places, and under -O that would return .pyo (this is well explained in the imp docs). Now all code (util.byte_compile, build_py, install_lib) can create .pyo files according to options given by users, without interference from the calling Python’s own optimize mode. On a related topic, I also removed the code that prevented byte compilation under python -B. The rationale is that packaging gives control over the creation of pyc files to the user with its own explicit option, and the behavior should not be changed if the calling Python happens to run with -B for whatever reason. I will argue that this is a bug fix and ask to be allowed to backport this change to distutils. Finally, I moved one nugget of information about the --compile and --optimize options from the source into the doc. It clears up a misunderstanding that I (and maybe other people) had.
author Éric Araujo <merwok@netwok.org>
date Thu, 03 Nov 2011 05:08:28 +0100
parents 2399e0451d0d
children 08418369da7b
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
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.0*t+t*t)); 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.0*t + t*t));
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(x*x+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.0*absx + 1.0 / (sqrt(x*x + 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.5*log1p(2x+2x*x/(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
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
192 double
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
193 _Py_log1p(double x)
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
194 {
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
195 /* 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
196 to 1+x, then
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
197
62820
a07ef5eab90b Merged revisions 82595 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 61463
diff changeset
198 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
199
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
200 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
201 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
202 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
203 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
204 (y-1)-x.
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
205
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
206 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
207 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
208 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
209 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
210 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
211 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
212 log(1+x) is exactly x.
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
213 */
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
214
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
215 double y;
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
216 if (fabs(x) < DBL_EPSILON/2.) {
62830
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
217 return x;
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
218 }
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
219 else if (-0.5 <= x && x <= 1.) {
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
220 /* WARNING: it's possible than an overeager compiler
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
221 will incorrectly optimize the following two lines
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
222 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
223 happens, then results from log1p will be inaccurate
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
224 for small x. */
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
225 y = 1.+x;
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
226 return log(y)-((y-1.)-x)/y;
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
227 }
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
228 else {
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
229 /* NaNs and infinities should end up here */
2399e0451d0d Merged revisions 82607 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 62820
diff changeset
230 return log(1.+x);
57866
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
231 }
9751ddcc1882 Merged revisions 76978 via svnmerge from
Mark Dickinson <dickinsm@gmail.com>
parents: 57758
diff changeset
232 }