Lines Matching +full:2 +full:x
12 /* j0(x), y0(x)
14 * Method -- j0(x):
15 * 1. For tiny x, we use j0(x) = 1 - x^2/4 + x^4/64 - ...
16 * 2. Reduce x to |x| since j0(x)=j0(-x), and
17 * for x in (0,2)
18 * j0(x) = 1-z/4+ z^2*R0/S0, where z = x*x;
19 * (precision: |j0-1+z/4-z^2R0/S0 |<2**-63.67 )
20 * for x in (2,inf)
21 * j0(x) = sqrt(2/(pi*x))*(p0(x)*cos(x0)-q0(x)*sin(x0))
22 * where x0 = x-pi/4. It is better to compute sin(x0),cos(x0)
24 * cos(x0) = cos(x)cos(pi/4)+sin(x)sin(pi/4)
25 * = 1/sqrt(2) * (cos(x) + sin(x))
26 * sin(x0) = sin(x)cos(pi/4)-cos(x)sin(pi/4)
27 * = 1/sqrt(2) * (sin(x) - cos(x))
29 * sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
37 * Method -- y0(x):
38 * 1. For x<2.
40 * y0(x) = 2/pi*(j0(x)*(ln(x/2)+Euler) + x^2/4 - ...)
41 * therefore y0(x)-2/pi*j0(x)*ln(x) is an even function.
43 * y0(x) = U(z)/V(z) + (2/pi)*(j0(x)*ln(x)), z= x^2
47 * with absolute approximation error bounded by 2**-72.
48 * Note: For tiny x, U/V = u0 and j0(x)~1, hence
49 * y0(tiny) = u0 + (2/pi)*ln(tiny), (choose tiny<2**-27)
50 * 2. For x>=2.
51 * y0(x) = sqrt(2/(pi*x))*(p0(x)*cos(x0)+q0(x)*sin(x0))
52 * where x0 = x-pi/4. It is better to compute sin(x0),cos(x0)
54 * 3. Special cases: y0(0)=-inf, y0(x<0)=NaN, y0(inf)=0.
82 j0(double x) in j0() argument
87 GET_HIGH_WORD(hx,x); in j0()
89 if(ix>=0x7ff00000) return one/(x*x); in j0()
90 x = fabs(x); in j0()
91 if(ix >= 0x40000000) { /* |x| >= 2.0 */ in j0()
92 sincos(x, &s, &c); in j0()
95 if(ix<0x7fe00000) { /* Make sure x+x does not overflow. */ in j0()
96 z = -cos(x+x); in j0()
101 * j0(x) = 1/sqrt(pi) * (P(0,x)*cc - Q(0,x)*ss) / sqrt(x) in j0()
102 * y0(x) = 1/sqrt(pi) * (P(0,x)*ss + Q(0,x)*cc) / sqrt(x) in j0()
104 if(ix>0x48000000) z = (invsqrtpi*cc)/sqrt(x); in j0()
106 u = pzero(x); v = qzero(x); in j0()
107 z = invsqrtpi*(u*cc-v*ss)/sqrt(x); in j0()
111 if(ix<0x3f200000) { /* |x| < 2**-13 */ in j0()
112 if(huge+x>one) { /* raise inexact if x != 0 */ in j0()
113 if(ix<0x3e400000) return one; /* |x|<2**-27 */ in j0()
114 else return one - x*x/4; in j0()
117 z = x*x; in j0()
120 if(ix < 0x3FF00000) { /* |x| < 1.00 */ in j0()
123 u = x/2; in j0()
142 y0(double x) in y0() argument
147 EXTRACT_WORDS(hx,lx,x); in y0()
154 if(ix>=0x7ff00000) return vone/(x+x*x); in y0()
157 /* y0(x<0) = NaN and raise invalid exception. */ in y0()
159 if(ix >= 0x40000000) { /* |x| >= 2.0 */ in y0()
160 /* y0(x) = sqrt(2/(pi*x))*(p0(x)*sin(x0)+q0(x)*cos(x0)) in y0()
161 * where x0 = x-pi/4 in y0()
163 * cos(x0) = cos(x)cos(pi/4)+sin(x)sin(pi/4) in y0()
164 * = 1/sqrt(2) * (sin(x) + cos(x)) in y0()
165 * sin(x0) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4) in y0()
166 * = 1/sqrt(2) * (sin(x) - cos(x)) in y0()
168 * sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x)) in y0()
171 sincos(x, &s, &c); in y0()
175 * j0(x) = 1/sqrt(pi) * (P(0,x)*cc - Q(0,x)*ss) / sqrt(x) in y0()
176 * y0(x) = 1/sqrt(pi) * (P(0,x)*ss + Q(0,x)*cc) / sqrt(x) in y0()
178 if(ix<0x7fe00000) { /* make sure x+x not overflow */ in y0()
179 z = -cos(x+x); in y0()
183 if(ix>0x48000000) z = (invsqrtpi*ss)/sqrt(x); in y0()
185 u = pzero(x); v = qzero(x); in y0()
186 z = invsqrtpi*(u*ss+v*cc)/sqrt(x); in y0()
190 if(ix<=0x3e400000) { /* x < 2**-27 */ in y0()
191 return(u00 + tpi*log(x)); in y0()
193 z = x*x; in y0()
196 return(u/v + tpi*(j0(x)*log(x))); in y0()
200 * 1 - 9/128 s^2 + 11025/98304 s^4 - ..., where s = 1/x.
201 * For x >= 2, We approximate pzero by
202 * pzero(x) = 1 + (R/S)
203 * where R = pR0 + pR1*s^2 + pR2*s^4 + ... + pR5*s^10
204 * S = 1 + pS0*s^2 + ... + pS4*s^10
206 * | pzero(x)-1-R/S | <= 2 ** ( -60.26)
208 static const double pR8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
224 static const double pR5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
240 static const double pR3[6] = {/* for x in [4.547,2.8571]=1/[0.2199,0.35001] */
256 static const double pR2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
273 pzero(double x) in pzero() argument
278 GET_HIGH_WORD(ix,x); in pzero()
284 z = one/(x*x); in pzero()
285 r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); in pzero()
286 s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*q[4])))); in pzero()
291 /* For x >= 8, the asymptotic expansions of qzero is
292 * -1/8 s + 75/1024 s^3 - ..., where s = 1/x.
294 * qzero(x) = s*(-1.25 + (R/S))
295 * where R = qR0 + qR1*s^2 + qR2*s^4 + ... + qR5*s^10
296 * S = 1 + qS0*s^2 + ... + qS5*s^12
298 * | qzero(x)/s +1.25-R/S | <= 2 ** ( -61.22)
300 static const double qR8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
317 static const double qR5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
334 static const double qR3[6] = {/* for x in [4.547,2.8571]=1/[0.2199,0.35001] */
351 static const double qR2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
369 qzero(double x) in qzero() argument
375 GET_HIGH_WORD(ix,x); in qzero()
381 z = one/(x*x); in qzero()
382 r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); in qzero()
383 s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*(q[4]+z*q[5]))))); in qzero()
384 return (r/s-eighth)/x; in qzero()