118 o_threshold = 7.09782712893383973096e+02,
119 ln2_hi = 6.93147180369123816490e-01,
120 ln2_lo = 1.90821492927058770002e-10,
121 invln2 = 1.44269504088896338700e+00,
123 Q1 = -3.33333333333331316428e-02,
124 Q2 = 1.58730158725481460165e-03,
125 Q3 = -7.93650757867487942473e-05,
126 Q4 = 4.00821782732936239552e-06,
127 Q5 = -2.01099218183624371326e-07;
130 double expm1(
double x)
136 double y,hi,lo,c=0,t,e,hxs,hfx,r1;
142 if(xsb==0) y=x;
else y= -x;
146 if(hx >= 0x4043687A) {
147 if(hx >= 0x40862E42) {
149 if(((hx&0xfffff)|
__LO(x))!=0)
151 else return (xsb==0)? x:-1.0;
153 if(x > o_threshold)
return huge*huge;
162 if(hx > 0x3fd62e42) {
163 if(hx < 0x3FF0A2B2) {
165 {hi = x - ln2_hi; lo = ln2_lo; k = 1;}
167 {hi = x + ln2_hi; lo = -ln2_lo; k = -1;}
169 k = invln2*x+((xsb==0)?0.5:-0.5);
177 else if(hx < 0x3c900000) {
179 return x - (t-(huge+x));
186 r1 =
one+hxs*(Q1+hxs*(Q2+hxs*(Q3+hxs*(Q4+hxs*Q5))));
188 e = hxs*((r1-t)/(6.0 - x*t));
189 if(k==0)
return x - (x*e-hxs);
193 if(k== -1)
return 0.5*(x-e)-0.5;
195 if(x < -0.25)
return -2.0*(e-(x+0.5));
196 else return one+2.0*(x-e);}
197 if (k <= -2 || k>56) {
204 __HI(t) = 0x3ff00000 - (0x200000>>k);
208 __HI(t) = ((0x3ff-k)<<20);