[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: [Help-smalltalk] BUG gst2.2: Float asExactFraction asFloat inexact
From: |
Paolo Bonzini |
Subject: |
Re: [Help-smalltalk] BUG gst2.2: Float asExactFraction asFloat inexact |
Date: |
Sat, 27 Jan 2007 16:51:58 +0100 |
User-agent: |
Thunderbird 1.5.0.9 (Macintosh/20061207) |
Unfortunately, current implementation of both Largeinteger>>asFloat and
Fraction>>asFloat will cause successive round off errors...
With 2.3.1 and your fix to asExactFraction, this works. But there are
more bugs in general, and I tried to rewrite the methods to use
timesTwoPower: as much as possible, and to avoid problems with denormals.
More over, denormalized number (gradual underflow) will convert back to
zero...
| x |
x := 1.1d-310.
^x asExactFraction asFloat = 0
This is exposing a bug in the lexer.
st> 1.1d-310 = 0!
true
st> 1.1d-306 / 10000 = 0!
false
However, there is another problem, in that (1 << 1072) is not a
representable floating point number, but its reciprocal is a (denormal)
floating point number. Therefore, the problem is that for denormals
asExactFraction returns a denominator that is not representable, and the
bug is actually in Fraction>>#asFloat.
The attached patch fixes all these bugs, except the one in the lexer. I
would like a double-check in Fraction>>#asFloat, maybe the rounding
could be improved if the numerator or denominator do not fit in a Float.
If you can provide some regression tests for rounding errors in
#asExactFraction and #asFloat, please do. Just write the code for the
borderline cases you can think of, I'll check if the results are right
since they may still be wrong with your outdated image.
Paolo
2007-01-27 Paolo Bonzini <address@hidden>
Nicolas Cellier <address@hidden>
* kernel/Float.st: Rewrite #truncated and #asExactFraction.
* kernel/Fraction.st: Use #asFloatD in #hash. Be careful about
denormals in #asFloat:
* kernel/LargeInt.st: Rewrite #asFloat:.
* kernel/MappedColl.st: Use map size as collection size.
2007-01-27 Paolo Bonzini <address@hidden>
* libgst/prims.def: Use trunc{,f,l} and lrint{,f,l} to truncate
floats.
--- orig/kernel/Float.st
+++ mod/kernel/Float.st
@@ -194,27 +194,30 @@ truncated
ifFalse: [ float := self negated ].
exponent := float exponent.
- bytes := ByteArray new: exponent // 8 + 2.
- float := float timesTwoPower: (exponent bitClear: 7) negated.
+ bytes := LargePositiveInteger new: (self class precision + 7) // 8 + 1.
+ float := float timesTwoPower: float class precision - exponent - 8.
- bytes size - 1 to: 1 by: -1 do: [ :i |
- bytes at: i put: float truncated.
- float := float fractionPart timesTwoPower: 8
+ 1 to: bytes size - 1 do: [ :i |
+ bytes digitAt: i put: (float fractionPart timesTwoPower: 8) truncated.
+ float := float integerPart timesTwoPower: -8.
].
- ^positive
- ifTrue: [ (LargeInteger from: bytes) ]
- ifFalse: [ (LargeInteger from: bytes) negated ]
+ bytes := bytes bitShift: (exponent - float class precision).
+ positive ifFalse: [ bytes := bytes negated ].
+ ^bytes
!
asExactFraction
"Convert the receiver into a fraction with optimal approximation,
but with usually huge terms."
+ | shift mantissa |
self checkCoercion.
-
- ^(self timesTwoPower: self exponent + self class precision) truncated /
- (1 bitShift: self exponent + self class precision)
+ shift := self exponent negated + self class precision.
+ mantissa := (self timesTwoPower: shift) truncated.
+ ^shift negative
+ ifTrue: [mantissa * (1 bitShift: shift negated)]
+ ifFalse: [mantissa / (1 bitShift: shift)]
!
asFraction
--- orig/kernel/Fraction.st
+++ mod/kernel/Fraction.st
@@ -288,7 +288,7 @@ truncated
hash
"Answer an hash value for the receiver"
denominator = 1 ifTrue: [ ^numerator hash ].
- ^(numerator asFloatD / denominator asFloatD) hash
+ ^self asFloatD hash
! !
@@ -385,7 +385,7 @@ storeOn: aStream
asFloat: characterization
"Answer the receiver converted to a Float"
- | n d shift sign |
+ | n d nShift dShift sign x |
n := numerator.
d := denominator.
sign := n sign * d sign.
@@ -394,11 +394,15 @@ asFloat: characterization
1e1800 asFloat / (1e1799 + 1) asFloat = NaN, but
(1e1800 / (1e1799 + 1)) asFloat must be 10."
- shift := (characterization emax - numerator highBit).
- shift := shift min: (characterization emax - denominator highBit).
- shift < 0 ifTrue: [
- "Lose some more precision, but we MUST avoid infinites and NaNs!"
- shift := shift - 10. n := n bitShift: shift. d := d bitShift: shift ].
+ nShift := (characterization emax - numerator highBit).
+ nShift < 0
+ ifTrue: [ n := n bitShift: nShift ]
+ ifFalse: [ nShift := 0 ].
+
+ dShift := (characterization emax - denominator highBit).
+ dShift < 0
+ ifTrue: [ d := d bitShift: dShift ]
+ ifFalse: [ dShift := 0 ].
d = 0 ifTrue: [
^sign > 0
@@ -409,8 +413,8 @@ asFloat: characterization
ifTrue: [ characterization coerce: 0 ]
ifFalse: [ characterization negativeInfinity reciprocal ] ].
- ^(characterization coerce: n)
- / (characterization coerce: d)
+ x := (characterization coerce: n) / (characterization coerce: d).
+ ^x timesTwoPower: dShift - nShift
!
reduce
--- orig/kernel/LargeInt.st
+++ mod/kernel/LargeInt.st
@@ -977,14 +988,16 @@ sign
asFloat: characterization
"Answer the receiver converted to a Float"
- | adjust result byte |
- byte := characterization coerce: 256.
- adjust := byte raisedToInteger: self size - 1.
+ | result byte shift |
+ shift := 8 * self size.
result := 0.
self size to: 1 by: -1 do: [ :index |
- result := (self at: index) * adjust + result.
- adjust := adjust / byte.
+ byte := self at: index.
+ byte = 0 ifFalse: [
+ byte := characterization coerce: byte.
+ result := (byte timesTwoPower: shift) + result ].
+ shift := shift - 8
].
^result
!
--- orig/kernel/MappedColl.st
+++ mod/kernel/MappedColl.st
@@ -89,7 +89,7 @@ at: key put: value
size
"Answer the receiver's size"
- ^domain size
+ ^map size
!
add: anObject
--- orig/libgst/prims.def
+++ mod/libgst/prims.def
@@ -1531,14 +1531,9 @@ primitive VMpr_FloatD_truncated [succeed
if COMMON (RECEIVER_IS_CLASS (oop1, _gst_floatd_class))
{
double oopValue = FLOATD_OOP_VALUE (oop1);
- if COMMON (oopValue >= 0.0 && oopValue <= MAX_ST_INT)
+ if COMMON (oopValue >= MIN_ST_INT && oopValue <= MAX_ST_INT)
{
- PUSH_INT ((intptr_t) (oopValue + 0.000000000000005));
- PRIM_SUCCEEDED;
- }
- else if COMMON (oopValue < 0.0 && oopValue >= MIN_ST_INT)
- {
- PUSH_INT ((intptr_t) (oopValue - 0.000000000000005));
+ PUSH_INT (lrint (trunc (oopValue)));
PRIM_SUCCEEDED;
}
}
@@ -1732,14 +1727,9 @@ primitive VMpr_FloatE_truncated [succeed
if COMMON (RECEIVER_IS_CLASS (oop1, _gst_floate_class))
{
double oopValue = FLOATE_OOP_VALUE (oop1);
- if COMMON (oopValue >= 0.0 && oopValue <= MAX_ST_INT)
+ if COMMON (oopValue >= MIN_ST_INT && oopValue <= MAX_ST_INT)
{
- PUSH_INT ((intptr_t) (oopValue + 0.000000000000005));
- PRIM_SUCCEEDED;
- }
- else if COMMON (oopValue < 0.0 && oopValue >= MIN_ST_INT)
- {
- PUSH_INT ((intptr_t) (oopValue - 0.000000000000005));
+ PUSH_INT (lrintf (truncf (oopValue)));
PRIM_SUCCEEDED;
}
}
@@ -1933,14 +1923,9 @@ primitive VMpr_FloatQ_truncated [succeed
if COMMON (RECEIVER_IS_CLASS (oop1, _gst_floatq_class))
{
long double oopValue = FLOATQ_OOP_VALUE (oop1);
- if COMMON (oopValue >= 0.0 && oopValue <= MAX_ST_INT)
- {
- PUSH_INT ((intptr_t) (oopValue + 0.000000000000005));
- PRIM_SUCCEEDED;
- }
- else if COMMON (oopValue < 0.0 && oopValue >= MIN_ST_INT)
+ if COMMON (oopValue >= MIN_ST_INT && oopValue <= MAX_ST_INT)
{
- PUSH_INT ((intptr_t) (oopValue - 0.000000000000005));
+ PUSH_INT (lrintl (truncl (oopValue)));
PRIM_SUCCEEDED;
}
}