static pthread_t pthread_threads_counter;
-#ifdef NEED_SEPARATE_REGISTER_STACK
-/* Signal masks for the manager. These have to be global only when clone2
- is used since it's currently borken wrt signals in the child. */
-static sigset_t manager_mask; /* Manager normal signal mask */
-static sigset_t manager_mask_all; /* All bits set. */
-#endif
-
/* Forward declarations */
static int pthread_handle_create(pthread_t *thread, const pthread_attr_t *attr,
{
int reqfd = (int) (long int) arg;
struct pollfd ufd;
-#ifndef NEED_SEPARATE_REGISTER_STACK
sigset_t manager_mask;
-#endif
int n;
struct pthread_request request;
if (__pthread_threads_debug && __pthread_sig_debug > 0)
sigdelset(&manager_mask, __pthread_sig_debug);
sigprocmask(SIG_SETMASK, &manager_mask, NULL);
-#ifdef NEED_SEPARATE_REGISTER_STACK
- sigfillset(&manager_mask_all);
-#endif
/* Raise our priority to match that of main thread */
__pthread_manager_adjust_prio(__pthread_main_thread->p_priority);
/* Synchronize debugging of the thread manager */
And there is some argument for changing the __clone2
interface to pass sp and bsp instead, making it more IA64
specific, but allowing stacks to grow outward from each
- other, to get less paging and fewer mmaps. Clone2
- currently can't take signals in the child right after
- process creation. Mask them in the child. It resets the
- mask once it starts up. */
- sigprocmask(SIG_SETMASK, &manager_mask_all, NULL);
+ other, to get less paging and fewer mmaps. */
pid = __clone2(pthread_start_thread_event,
(void **)new_thread_bottom,
(char *)new_thread - new_thread_bottom,
CLONE_VM | CLONE_FS | CLONE_FILES | CLONE_SIGHAND |
__pthread_sig_cancel, new_thread);
- sigprocmask(SIG_SETMASK, &manager_mask, NULL);
#else
pid = __clone(pthread_start_thread_event, (void **) new_thread,
CLONE_VM | CLONE_FS | CLONE_FILES | CLONE_SIGHAND |
if (pid == 0)
{
#ifdef NEED_SEPARATE_REGISTER_STACK
- sigprocmask(SIG_SETMASK, &manager_mask_all, NULL);
pid = __clone2(pthread_start_thread,
(void **)new_thread_bottom,
(char *)new_thread - new_thread_bottom,
CLONE_VM | CLONE_FS | CLONE_FILES | CLONE_SIGHAND |
__pthread_sig_cancel, new_thread);
- sigprocmask(SIG_SETMASK, &manager_mask, NULL);
#else
pid = __clone(pthread_start_thread, (void **) new_thread,
CLONE_VM | CLONE_FS | CLONE_FILES | CLONE_SIGHAND |
ifloat: 1
double: 2
idouble: 2
+ldouble: 2
+ildouble: 2
Test "jn (0, 2.0) == 0.22389077914123566805":
float: 1
ifloat: 1
Test "jn (0, 8.0) == 0.17165080713755390609":
float: 1
ifloat: 1
+ldouble: 1
+ildouble: 1
Test "jn (1, 10.0) == 0.043472746168861436670":
float: 1
ifloat: 1
double: 2
idouble: 2
+ldouble: 2
+ildouble: 2
Test "jn (1, 2.0) == 0.57672480775687338720":
double: 1
idouble: 1
Test "jn (1, 8.0) == 0.23463634685391462438":
float: 1
ifloat: 1
+ldouble: 1
+ildouble: 1
Test "jn (10, -1.0) == 0.26306151236874532070e-9":
float: 1
ifloat: 1
+ldouble: 1
+ildouble: 1
Test "jn (10, 0.1) == 0.26905328954342155795e-19":
double: 4
float: 2
idouble: 4
ifloat: 2
+ldouble: 1
+ildouble: 1
Test "jn (10, 0.7) == 0.75175911502153953928e-11":
double: 4
float: 1
idouble: 4
ifloat: 1
+ldouble: 2
+ildouble: 2
Test "jn (10, 1.0) == 0.26306151236874532070e-9":
float: 1
ifloat: 1
+ldouble: 1
+ildouble: 1
Test "jn (10, 2.0) == 0.25153862827167367096e-6":
float: 1
ifloat: 1
double: 2
idouble: 2
+ldouble: 1
+ildouble: 1
Test "jn (10, 10.0) == 0.20748610663335885770":
float: 2
ifloat: 2
double: 4
idouble: 4
+ldouble: 2
+ildouble: 2
Test "jn (3, 0.1) == 0.000020820315754756261429":
double: 1
idouble: 1
+ldouble: 1
+ildouble: 1
Test "jn (3, 0.7) == 0.0069296548267508408077":
double: 2
idouble: 2
ifloat: 1
double: 3
idouble: 3
+ldouble: 1
+ildouble: 1
+Test "jn (3, -1.0) == -0.019563353982668405919":
+ldouble: 1
+ildouble: 1
+Test "jn (3, 1.0) == 0.019563353982668405919":
+ldouble: 1
+ildouble: 1
# lgamma
Test "lgamma (-0.5) == log(2*sqrt(pi))":
ifloat: 1
double: 3
idouble: 3
+ldouble: 2
+ildouble: 2
Test "yn (0, 1.0) == 0.088256964215676957983":
double: 2
float: 1
ifloat: 1
double: 2
idouble: 2
+ldouble: 1
+ildouble: 1
Test "yn (0, 8.0) == 0.22352148938756622053":
float: 1
ifloat: 1
double: 1
idouble: 1
+ldouble: 1
+ildouble: 1
Test "yn (1, 0.1) == -6.4589510947020269877":
double: 1
float: 1
idouble: 1
ifloat: 1
+ldouble: 1
+ildouble: 1
Test "yn (1, 0.7) == -1.1032498719076333697":
double: 1
idouble: 1
float: 2
idouble: 1
ifloat: 2
+ldouble: 1
+ildouble: 1
Test "yn (1, 8.0) == -0.15806046173124749426":
float: 2
ifloat: 2
+ldouble: 2
+ildouble: 2
Test "yn (10, 0.1) == -0.11831335132045197885e19":
double: 2
float: 1
idouble: 2
ifloat: 1
+ldouble: 2
+ildouble: 2
Test "yn (10, 0.7) == -0.42447194260703866924e10":
double: 6
float: 3
idouble: 6
ifloat: 3
+ldouble: 7
+ildouble: 7
Test "yn (10, 1.0) == -0.12161801427868918929e9":
double: 1
float: 1
ifloat: 1
double: 2
idouble: 2
+ldouble: 1
+ildouble: 1
Test "yn (3, 0.1) == -5099.3323786129048894":
double: 1
float: 1
idouble: 1
ifloat: 1
+ldouble: 2
+ildouble: 2
Test "yn (3, 0.7) == -15.819479052819633505":
double: 2
idouble: 2
+ldouble: 2
+ildouble: 2
Test "yn (3, 2.0) == -1.1277837768404277861":
float: 1
ifloat: 1
float: 2
idouble: 4
ifloat: 2
+ldouble: 2
+ildouble: 2
Function: "lgamma":
double: 1
float: 3
idouble: 6
ifloat: 3
+ldouble: 7
+ildouble: 7
# end of automatic generation
--- /dev/null
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/* Modifications for long double contributed by
+ Stephen L. Moshier <moshier@na-net.ornl.gov> */
+
+/*
+ * __ieee754_jn(n, x), __ieee754_yn(n, x)
+ * floating point Bessel's function of the 1st and 2nd kind
+ * of order n
+ *
+ * Special cases:
+ * y0(0)=y1(0)=yn(n,0) = -inf with division by zero signal;
+ * y0(-ve)=y1(-ve)=yn(n,-ve) are NaN with invalid signal.
+ * Note 2. About jn(n,x), yn(n,x)
+ * For n=0, j0(x) is called,
+ * for n=1, j1(x) is called,
+ * for n<x, forward recursion us used starting
+ * from values of j0(x) and j1(x).
+ * for n>x, a continued fraction approximation to
+ * j(n,x)/j(n-1,x) is evaluated and then backward
+ * recursion is used starting from a supposed value
+ * for j(n,x). The resulting value of j(0,x) is
+ * compared with the actual value to correct the
+ * supposed value of j(n,x).
+ *
+ * yn(n,x) is similar in all respects, except
+ * that forward recursion is used for all
+ * values of n>1.
+ *
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+#ifdef __STDC__
+static const long double
+#else
+static long double
+#endif
+ invsqrtpi = 5.64189583547756286948079e-1L, two = 2.0e0L, one = 1.0e0L;
+
+#ifdef __STDC__
+static const long double zero = 0.0L;
+#else
+static long double zero = 0.0L;
+#endif
+
+#ifdef __STDC__
+long double
+__ieee754_jnl (int n, long double x)
+#else
+long double
+__ieee754_jnl (n, x)
+ int n;
+ long double x;
+#endif
+{
+ u_int32_t se, i0, i1;
+ int32_t i, ix, sgn;
+ long double a, b, temp, di;
+ long double z, w;
+
+ /* J(-n,x) = (-1)^n * J(n, x), J(n, -x) = (-1)^n * J(n, x)
+ * Thus, J(-n,x) = J(n,-x)
+ */
+
+ GET_LDOUBLE_WORDS (se, i0, i1, x);
+ ix = se & 0x7fff;
+
+ /* if J(n,NaN) is NaN */
+ if ((ix == 0x7fff) && ((i0 & 0x7fffffff) != 0))
+ return x + x;
+ if (n < 0)
+ {
+ n = -n;
+ x = -x;
+ se ^= 0x8000;
+ }
+ if (n == 0)
+ return (__ieee754_j0l (x));
+ if (n == 1)
+ return (__ieee754_j1l (x));
+ sgn = (n & 1) & (se >> 15); /* even n -- 0, odd n -- sign(x) */
+ x = fabsl (x);
+ if ((ix | i0 | i1) == 0 || ix >= 0x7fff) /* if x is 0 or inf */
+ b = zero;
+ else if ((long double) n <= x)
+ {
+ /* Safe to use J(n+1,x)=2n/x *J(n,x)-J(n-1,x) */
+ if (ix >= 0x412D)
+ { /* x > 2**302 */
+
+ /* ??? This might be a futile gesture.
+ If x exceeds X_TLOSS anyway, the wrapper function
+ will set the result to zero. */
+
+ /* (x >> n**2)
+ * Jn(x) = cos(x-(2n+1)*pi/4)*sqrt(2/x*pi)
+ * Yn(x) = sin(x-(2n+1)*pi/4)*sqrt(2/x*pi)
+ * Let s=sin(x), c=cos(x),
+ * xn=x-(2n+1)*pi/4, sqt2 = sqrt(2),then
+ *
+ * n sin(xn)*sqt2 cos(xn)*sqt2
+ * ----------------------------------
+ * 0 s-c c+s
+ * 1 -s-c -c+s
+ * 2 -s+c -c-s
+ * 3 s+c c-s
+ */
+ long double s;
+ long double c;
+ __sincosl (x, &s, &c);
+ switch (n & 3)
+ {
+ case 0:
+ temp = c + s;
+ break;
+ case 1:
+ temp = -c + s;
+ break;
+ case 2:
+ temp = -c - s;
+ break;
+ case 3:
+ temp = c - s;
+ break;
+ }
+ b = invsqrtpi * temp / __ieee754_sqrtl (x);
+ }
+ else
+ {
+ a = __ieee754_j0l (x);
+ b = __ieee754_j1l (x);
+ for (i = 1; i < n; i++)
+ {
+ temp = b;
+ b = b * ((long double) (i + i) / x) - a; /* avoid underflow */
+ a = temp;
+ }
+ }
+ }
+ else
+ {
+ if (ix < 0x3fde)
+ { /* x < 2**-33 */
+ /* x is tiny, return the first Taylor expansion of J(n,x)
+ * J(n,x) = 1/n!*(x/2)^n - ...
+ */
+ if (n >= 400) /* underflow, result < 10^-4952 */
+ b = zero;
+ else
+ {
+ temp = x * 0.5;
+ b = temp;
+ for (a = one, i = 2; i <= n; i++)
+ {
+ a *= (long double) i; /* a = n! */
+ b *= temp; /* b = (x/2)^n */
+ }
+ b = b / a;
+ }
+ }
+ else
+ {
+ /* use backward recurrence */
+ /* x x^2 x^2
+ * J(n,x)/J(n-1,x) = ---- ------ ------ .....
+ * 2n - 2(n+1) - 2(n+2)
+ *
+ * 1 1 1
+ * (for large x) = ---- ------ ------ .....
+ * 2n 2(n+1) 2(n+2)
+ * -- - ------ - ------ -
+ * x x x
+ *
+ * Let w = 2n/x and h=2/x, then the above quotient
+ * is equal to the continued fraction:
+ * 1
+ * = -----------------------
+ * 1
+ * w - -----------------
+ * 1
+ * w+h - ---------
+ * w+2h - ...
+ *
+ * To determine how many terms needed, let
+ * Q(0) = w, Q(1) = w(w+h) - 1,
+ * Q(k) = (w+k*h)*Q(k-1) - Q(k-2),
+ * When Q(k) > 1e4 good for single
+ * When Q(k) > 1e9 good for double
+ * When Q(k) > 1e17 good for quadruple
+ */
+ /* determine k */
+ long double t, v;
+ long double q0, q1, h, tmp;
+ int32_t k, m;
+ w = (n + n) / (long double) x;
+ h = 2.0L / (long double) x;
+ q0 = w;
+ z = w + h;
+ q1 = w * z - 1.0L;
+ k = 1;
+ while (q1 < 1.0e11L)
+ {
+ k += 1;
+ z += h;
+ tmp = z * q1 - q0;
+ q0 = q1;
+ q1 = tmp;
+ }
+ m = n + n;
+ for (t = zero, i = 2 * (n + k); i >= m; i -= 2)
+ t = one / (i / x - t);
+ a = t;
+ b = one;
+ /* estimate log((2/x)^n*n!) = n*log(2/x)+n*ln(n)
+ * Hence, if n*(log(2n/x)) > ...
+ * single 8.8722839355e+01
+ * double 7.09782712893383973096e+02
+ * long double 1.1356523406294143949491931077970765006170e+04
+ * then recurrent value may overflow and the result is
+ * likely underflow to zero
+ */
+ tmp = n;
+ v = two / x;
+ tmp = tmp * __ieee754_logl (fabsl (v * tmp));
+
+ if (tmp < 1.1356523406294143949491931077970765006170e+04L)
+ {
+ for (i = n - 1, di = (long double) (i + i); i > 0; i--)
+ {
+ temp = b;
+ b *= di;
+ b = b / x - a;
+ a = temp;
+ di -= two;
+ }
+ }
+ else
+ {
+ for (i = n - 1, di = (long double) (i + i); i > 0; i--)
+ {
+ temp = b;
+ b *= di;
+ b = b / x - a;
+ a = temp;
+ di -= two;
+ /* scale b to avoid spurious overflow */
+ if (b > 1e100L)
+ {
+ a /= b;
+ t /= b;
+ b = one;
+ }
+ }
+ }
+ b = (t * __ieee754_j0l (x) / b);
+ }
+ }
+ if (sgn == 1)
+ return -b;
+ else
+ return b;
+}
+
+#ifdef __STDC__
+long double
+__ieee754_ynl (int n, long double x)
+#else
+long double
+__ieee754_ynl (n, x)
+ int n;
+ long double x;
+#endif
+{
+ u_int32_t se, i0, i1;
+ int32_t i, ix;
+ int32_t sign;
+ long double a, b, temp;
+
+
+ GET_LDOUBLE_WORDS (se, i0, i1, x);
+ ix = se & 0x7fff;
+ /* if Y(n,NaN) is NaN */
+ if ((ix == 0x7fff) && ((i0 & 0x7fffffff) != 0))
+ return x + x;
+ if ((ix | i0 | i1) == 0)
+ return -one / zero;
+ if (se & 0x8000)
+ return zero / zero;
+ sign = 1;
+ if (n < 0)
+ {
+ n = -n;
+ sign = 1 - ((n & 1) << 1);
+ }
+ if (n == 0)
+ return (__ieee754_y0l (x));
+ if (n == 1)
+ return (sign * __ieee754_y1l (x));
+ if (ix == 0x7fff)
+ return zero;
+ if (ix >= 0x412D)
+ { /* x > 2**302 */
+
+ /* ??? See comment above on the possible futility of this. */
+
+ /* (x >> n**2)
+ * Jn(x) = cos(x-(2n+1)*pi/4)*sqrt(2/x*pi)
+ * Yn(x) = sin(x-(2n+1)*pi/4)*sqrt(2/x*pi)
+ * Let s=sin(x), c=cos(x),
+ * xn=x-(2n+1)*pi/4, sqt2 = sqrt(2),then
+ *
+ * n sin(xn)*sqt2 cos(xn)*sqt2
+ * ----------------------------------
+ * 0 s-c c+s
+ * 1 -s-c -c+s
+ * 2 -s+c -c-s
+ * 3 s+c c-s
+ */
+ long double s;
+ long double c;
+ __sincosl (x, &s, &c);
+ switch (n & 3)
+ {
+ case 0:
+ temp = s - c;
+ break;
+ case 1:
+ temp = -s - c;
+ break;
+ case 2:
+ temp = -s + c;
+ break;
+ case 3:
+ temp = s + c;
+ break;
+ }
+ b = invsqrtpi * temp / __ieee754_sqrtl (x);
+ }
+ else
+ {
+ a = __ieee754_y0l (x);
+ b = __ieee754_y1l (x);
+ /* quit if b is -inf */
+ GET_LDOUBLE_WORDS (se, i0, i1, b);
+ for (i = 1; i < n && se != 0xffff; i++)
+ {
+ temp = b;
+ b = ((long double) (i + i) / x) * b - a;
+ GET_LDOUBLE_WORDS (se, i0, i1, b);
+ a = temp;
+ }
+ }
+ if (sign > 0)
+ return b;
+ else
+ return -b;
+}
git://sourceware.org / glibc.git / commitdiff