|
@@ -676,17 +676,17 @@ type
|
|
|
* to produce the hexadecimal values shown.
|
|
|
*
|
|
|
}
|
|
|
- function fpc_exp_real(x: Double):Double;compilerproc;
|
|
|
+ function fpc_exp_real(d: Double):Double;compilerproc;
|
|
|
const
|
|
|
- one = 1.0,
|
|
|
+ one = 1.0;
|
|
|
halF : array[0..1] of double = (0.5,-0.5);
|
|
|
huge = 1.0e+300;
|
|
|
twom1000 = 9.33263618503218878990e-302; { 2**-1000=0x01700000,0}
|
|
|
o_threshold = 7.09782712893383973096e+02; { 0x40862E42, 0xFEFA39EF }
|
|
|
u_threshold = -7.45133219101941108420e+02; { 0xc0874910, 0xD52D3051 }
|
|
|
- ln2HI : array[0..1] of double = ( 6.93147180369123816490e-01: { 0x3fe62e42, 0xfee00000 }
|
|
|
+ ln2HI : array[0..1] of double = ( 6.93147180369123816490e-01, { 0x3fe62e42, 0xfee00000 }
|
|
|
-6.93147180369123816490e-01); { 0xbfe62e42, 0xfee00000 }
|
|
|
- ln2LO : array[0..1] of double = (1.90821492927058770002e-10; { 0x3dea39ef, 0x35793c76 }
|
|
|
+ ln2LO : array[0..1] of double = (1.90821492927058770002e-10, { 0x3dea39ef, 0x35793c76 }
|
|
|
-1.90821492927058770002e-10); { 0xbdea39ef, 0x35793c76 }
|
|
|
invln2 = 1.44269504088896338700e+00; { 0x3ff71547, 0x652b82fe }
|
|
|
P1 = 1.66666666666666019037e-01; { 0x3FC55555, 0x5555553E }
|
|
@@ -695,44 +695,43 @@ type
|
|
|
P4 = -1.65339022054652515390e-06; { 0xBEBBBD41, 0xC5D26BF1 }
|
|
|
P5 = 4.13813679705723846039e-08; { 0x3E663769, 0x72BEA4D0 }
|
|
|
var
|
|
|
- double hi : double = 0.0;
|
|
|
- double lo = 0.0;
|
|
|
- c : double;
|
|
|
- t : double;
|
|
|
- int32_t k=0;
|
|
|
- xsb : longint;
|
|
|
+ c,hi,lo,t,y : double;
|
|
|
+ k,xsb : longint;
|
|
|
hx,hy,lx : dword;
|
|
|
begin
|
|
|
- hx:=float64(x).high;
|
|
|
- xsb := (hx shr 31) and 1; { sign bit of x }
|
|
|
- hx := hx and $7fffffff; { high word of |x| }
|
|
|
+ hi:=0.0;
|
|
|
+ lo:=0.0;
|
|
|
+ k:=0;
|
|
|
+ hx:=float64(d).high;
|
|
|
+ xsb := (hx shr 31) and 1; { sign bit of d }
|
|
|
+ hx := hx and $7fffffff; { high word of |d| }
|
|
|
|
|
|
{ filter out non-finite argument }
|
|
|
if hx >= $40862E42 then
|
|
|
- begin { if |x|>=709.78... }
|
|
|
- if hx >= $7ff00000
|
|
|
+ begin { if |d|>=709.78... }
|
|
|
+ if hx >= $7ff00000 then
|
|
|
begin
|
|
|
- lx:=float64(x).low;
|
|
|
+ lx:=float64(d).low;
|
|
|
if ((hx and $fffff) or lx)<>0 then
|
|
|
begin
|
|
|
- result:=x+x; { NaN }
|
|
|
+ result:=d+d; { NaN }
|
|
|
exit;
|
|
|
- else
|
|
|
+ end
|
|
|
else
|
|
|
begin
|
|
|
if xsb=0 then
|
|
|
- result:=x
|
|
|
+ result:=d
|
|
|
else
|
|
|
result:=0.0; { exp(+-inf)=begininf,0end }
|
|
|
exit;
|
|
|
end;
|
|
|
end;
|
|
|
- if x > o_threshold then
|
|
|
+ if d > o_threshold then
|
|
|
begin
|
|
|
result:=huge*huge; { overflow }
|
|
|
exit;
|
|
|
end;
|
|
|
- if x < u_threshold then
|
|
|
+ if d < u_threshold then
|
|
|
begin
|
|
|
result:=twom1000*twom1000; { underflow }
|
|
|
exit;
|
|
@@ -741,45 +740,45 @@ type
|
|
|
|
|
|
{ argument reduction }
|
|
|
if hx > $3fd62e42 then
|
|
|
- begin { if |x| > 0.5 ln2 }
|
|
|
- if hx < $3FF0A2B2 then { and |x| < 1.5 ln2 }
|
|
|
+ begin { if |d| > 0.5 ln2 }
|
|
|
+ if hx < $3FF0A2B2 then { and |d| < 1.5 ln2 }
|
|
|
begin
|
|
|
- hi := x-ln2HI[xsb];
|
|
|
+ hi := d-ln2HI[xsb];
|
|
|
lo:=ln2LO[xsb];
|
|
|
k := 1-xsb-xsb;
|
|
|
end
|
|
|
else
|
|
|
begin
|
|
|
- k := invln2*x+halF[xsb];
|
|
|
+ k := round(invln2*d+halF[xsb]);
|
|
|
t := k;
|
|
|
- hi := x - t*ln2HI[0]; { t*ln2HI is exact here }
|
|
|
+ hi := d - t*ln2HI[0]; { t*ln2HI is exact here }
|
|
|
lo := t*ln2LO[0];
|
|
|
end;
|
|
|
- x := hi - lo;
|
|
|
+ d := hi - lo;
|
|
|
end
|
|
|
else if hx < $3e300000 then
|
|
|
- begin { when |x|<2**-28 }
|
|
|
- if huge+x>one then
|
|
|
+ begin { when |d|<2**-28 }
|
|
|
+ if huge+d>one then
|
|
|
begin
|
|
|
- result:=one+x;{ trigger inexact }
|
|
|
+ result:=one+d;{ trigger inexact }
|
|
|
exit;
|
|
|
end;
|
|
|
end
|
|
|
else
|
|
|
k := 0;
|
|
|
|
|
|
- { x is now in primary range }
|
|
|
- t:=x*x;
|
|
|
- c:=x - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
|
|
|
+ { d is now in primary range }
|
|
|
+ t:=d*d;
|
|
|
+ c:=d - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
|
|
|
if k=0 then
|
|
|
begin
|
|
|
- result:=one-((x*c)/(c-2.0)-x);
|
|
|
+ result:=one-((d*c)/(c-2.0)-d);
|
|
|
exit;
|
|
|
end
|
|
|
else
|
|
|
- y := one-((lo-(x*c)/(2.0-c))-hi);
|
|
|
+ y := one-((lo-(d*c)/(2.0-c))-hi);
|
|
|
|
|
|
- if k >= -1021
|
|
|
+ if k >= -1021 then
|
|
|
begin
|
|
|
hy:=float64(y).high;
|
|
|
float64(y).high:=hy+(k shl 20); { add k to y's exponent }
|