Lines Matching +full:2 +full:x
12 /* j1(x), y1(x)
14 * Method -- j1(x):
15 * 1. For tiny x, we use j1(x) = x/2 - x^3/16 + x^5/384 - ...
16 * 2. Reduce x to |x| since j1(x)=-j1(-x), and
17 * for x in (0,2)
18 * j1(x) = x/2 + x*z*R0/S0, where z = x*x;
19 * (precision: |j1/x - 1/2 - R0/S0 |<2**-61.51 )
20 * for x in (2,inf)
21 * j1(x) = sqrt(2/(pi*x))*(p1(x)*cos(x1)-q1(x)*sin(x1))
22 * y1(x) = sqrt(2/(pi*x))*(p1(x)*sin(x1)+q1(x)*cos(x1))
23 * where x1 = x-3*pi/4. It is better to compute sin(x1),cos(x1)
25 * cos(x1) = cos(x)cos(3pi/4)+sin(x)sin(3pi/4)
26 * = 1/sqrt(2) * (sin(x) - cos(x))
27 * sin(x1) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4)
28 * = -1/sqrt(2) * (sin(x) + cos(x))
30 * sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
38 * Method -- y1(x):
39 * 1. screen out x<=0 cases: y1(0)=-inf, y1(x<0)=NaN
40 * 2. For x<2.
42 * y1(x) = 2/pi*(j1(x)*(ln(x/2)+Euler)-1/x-x/2+5/64*x^3-...)
43 * therefore y1(x)-2/pi*j1(x)*ln(x)-1/x is an odd function.
45 * y1(x) = x*U(z)/V(z) + (2/pi)*(j1(x)*ln(x)-1/x), z= x^2
46 * where for x in [0,2] (abs err less than 2**-65.89)
49 * Note: For tiny x, 1/x dominate y1 and hence
50 * y1(tiny) = -2/pi/tiny, (choose tiny<2**-54)
51 * 3. For x>=2.
52 * y1(x) = sqrt(2/(pi*x))*(p1(x)*sin(x1)+q1(x)*cos(x1))
53 * where x1 = x-3*pi/4. It is better to compute sin(x1),cos(x1)
69 /* R0/S0 on [0,2] */
83 j1(double x) in j1() argument
88 GET_HIGH_WORD(hx,x); in j1()
90 if(ix>=0x7ff00000) return one/x; in j1()
91 y = fabs(x); in j1()
92 if(ix >= 0x40000000) { /* |x| >= 2.0 */ in j1()
102 * j1(x) = 1/sqrt(pi) * (P(1,x)*cc - Q(1,x)*ss) / sqrt(x) in j1()
103 * y1(x) = 1/sqrt(pi) * (P(1,x)*ss + Q(1,x)*cc) / sqrt(x) in j1()
113 if(ix<0x3e400000) { /* |x|<2**-27 */ in j1()
114 if(huge+x>one) return 0.5*x;/* inexact if x!=0 necessary */ in j1()
116 z = x*x; in j1()
119 r *= x; in j1()
120 return(x*0.5+r/s); in j1()
139 y1(double x) in y1() argument
144 EXTRACT_WORDS(hx,lx,x); in y1()
151 if(ix>=0x7ff00000) return vone/(x+x*x); in y1()
154 /* y1(x<0) = NaN and raise invalid exception. */ in y1()
156 if(ix >= 0x40000000) { /* |x| >= 2.0 */ in y1()
157 sincos(x, &s, &c); in y1()
160 if(ix<0x7fe00000) { /* make sure x+x not overflow */ in y1()
161 z = cos(x+x); in y1()
165 /* y1(x) = sqrt(2/(pi*x))*(p1(x)*sin(x0)+q1(x)*cos(x0)) in y1()
166 * where x0 = x-3pi/4 in y1()
168 * cos(x0) = cos(x)cos(3pi/4)+sin(x)sin(3pi/4) in y1()
169 * = 1/sqrt(2) * (sin(x) - cos(x)) in y1()
170 * sin(x0) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4) in y1()
171 * = -1/sqrt(2) * (cos(x) + sin(x)) in y1()
173 * sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x)) in y1()
176 if(ix>0x48000000) z = (invsqrtpi*ss)/sqrt(x); in y1()
178 u = pone(x); v = qone(x); in y1()
179 z = invsqrtpi*(u*ss+v*cc)/sqrt(x); in y1()
183 if(ix<=0x3c900000) { /* x < 2**-54 */ in y1()
184 return(-tpi/x); in y1()
186 z = x*x; in y1()
187 u = U0[0]+z*(U0[1]+z*(U0[2]+z*(U0[3]+z*U0[4]))); in y1()
188 v = one+z*(V0[0]+z*(V0[1]+z*(V0[2]+z*(V0[3]+z*V0[4])))); in y1()
189 return(x*(u/v) + tpi*(j1(x)*log(x)-one/x)); in y1()
192 /* For x >= 8, the asymptotic expansions of pone is
193 * 1 + 15/128 s^2 - 4725/2^15 s^4 - ..., where s = 1/x.
195 * pone(x) = 1 + (R/S)
196 * where R = pr0 + pr1*s^2 + pr2*s^4 + ... + pr5*s^10
197 * S = 1 + ps0*s^2 + ... + ps4*s^10
199 * | pone(x)-1-R/S | <= 2 ** ( -60.06)
202 static const double pr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
218 static const double pr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
250 static const double pr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
267 pone(double x) in pone() argument
272 GET_HIGH_WORD(ix,x); in pone()
278 z = one/(x*x); in pone()
279 r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); in pone()
280 s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*q[4])))); in pone()
285 /* For x >= 8, the asymptotic expansions of qone is
286 * 3/8 s - 105/1024 s^3 - ..., where s = 1/x.
288 * qone(x) = s*(0.375 + (R/S))
289 * where R = qr1*s^2 + qr2*s^4 + ... + qr5*s^10
290 * S = 1 + qs1*s^2 + ... + qs6*s^12
292 * | qone(x)/s -0.375-R/S | <= 2 ** ( -61.13)
295 static const double qr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
312 static const double qr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
346 static const double qr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
364 qone(double x) in qone() argument
369 GET_HIGH_WORD(ix,x); in qone()
375 z = one/(x*x); in qone()
376 r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); in qone()
377 s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*(q[4]+z*q[5]))))); in qone()
378 return (.375 + r/s)/x; in qone()