@@ -14,86 +14,109 @@ module Data.Aeson.Encode.Double
1414 double
1515 ) where
1616
17−import GHC.Float
18−
17+import Blaze.ByteString.Builder (Builder, fromByteString)
18+import Blaze.ByteString.Builder.Char8 (fromChar)
19+import Data.Aeson.Encode.Int (digit, int, minus)
1920import Data.ByteString.Char8 ()
20−import Data.Monoid
21−import Data.Aeson.Encode.Int
22−import Blaze.ByteString.Builder
23−import Blaze.ByteString.Builder.Char8
21+import Data.Monoid (mappend, mconcat, mempty)
22+import qualified Data.Vector as V
2423
25−double :: Double -> Builder
26−double = showpGFloat Nothing
24+-- The code below is originally from GHC.Float, but has been optimised
25+-- in quite a few ways.
2726
28−-- | Show a signed RealFloat value using decimal notation when the
29−-- absolute value lies between 0.1 and 9,999,999, and scientific
30−-- notation otherwise. The optional integer can be used to specify
31−-- precision.
32−showpGFloat :: RealFloat a => Maybe Int -> a -> Builder
33−showpGFloat = putFormattedFloat FFGeneric
27+data T = T [Int] {-# UNPACK #-} !Int
3428
35−-- | Show a signed RealFloat value using decimal notation. The optional
36−-- integer can be used to specify precision.
37−showpFFloat :: RealFloat a => Maybe Int -> a -> Builder
38−showpFFloat = putFormattedFloat FFFixed
29+double :: Double -> Builder
30+double f
31+ | isNaN f = fromByteString "NaN"
32+ | isInfinite f = fromByteString $
33+ if f < 0 then "-Infinity" else "Infinity"
34+ | f < 0 || isNegativeZero f = minus `mappend` goGeneric (floatToDigits (-f))
35+ | otherwise = goGeneric (floatToDigits f)
36+ where
37+ goGeneric p@(T _ e)
38+ | e < 0 || e > 7 = goExponent p
39+ | otherwise = goFixed p
40+ goExponent (T is e) =
41+ case is of
42+ [] -> error "putFormattedFloat"
43+ [0] -> fromByteString "0.0e0"
44+ [d] -> digit d `mappend` fromByteString ".0e" `mappend` int (e-1)
45+ (d:ds) -> digit d `mappend` fromChar '.' `mappend` digits ds `mappend`
46+ fromChar 'e' `mappend` int (e-1)
47+ goFixed (T is e)
48+ | e <= 0 = fromChar '0' `mappend` fromChar '.' `mappend`
49+ mconcat (replicate (-e) (fromChar '0')) `mappend`
50+ digits is
51+ | otherwise = let g 0 rs = fromChar '.' `mappend` mk0 rs
52+ g n [] = fromChar '0' `mappend` g (n-1) []
53+ g n (r:rs) = digit r `mappend` g (n-1) rs
54+ in g e is
55+ mk0 [] = fromChar '0'
56+ mk0 rs = digits rs
3957
40−-- | Show a signed RealFloat value using scientific (exponential) notation.
41−-- The optional integer can be used to specify precision.
42−showpEFloat :: RealFloat a => Maybe Int -> a -> Builder
43−showpEFloat = putFormattedFloat FFExponent
58+digits :: [Int] -> Builder
59+digits (d:ds) = digit d `mappend` digits ds
60+digits _ = mempty
61+{-# INLINE digits #-}
4462
45−putFormattedFloat :: RealFloat a => FFFormat -> Maybe Int -> a -> Builder
46−putFormattedFloat fmt decs f
47− | isNaN f = fromChar 'N' `mappend` fromChar 'a' `mappend` fromChar 'N'
48− | isInfinite f = fromByteString (if f < 0 then "-Infinity" else "Infinity")
49− | f < 0 || isNegativeZero f = fromChar '-' `mappend` go fmt (floatToDigits (toInteger base) (-f))
50− | otherwise = go fmt (floatToDigits (toInteger base) f)
63+floatToDigits :: Double -> T
64+floatToDigits 0 = T [0] 0
65+floatToDigits x = T (reverse rds) k
5166 where
52− base = 10
53−
54− go FFGeneric p@(_,e)
55− | e < 0 || e > 7 = go FFExponent p
56− | otherwise = go FFFixed p
57− go FFExponent (is, e) =
58− case decs of
59− Nothing -> case is of
60− [] -> error "putFormattedFloat"
61− [0] -> fromByteString "0.0e0"
62− [d] -> digit d `mappend` fromByteString ".0e" `mappend` int (e-1)
63− (d:ds) -> digit d `mappend` fromChar '.' `mappend` mconcat (map digit ds)
64− `mappend` fromChar 'e' `mappend` int (e-1)
65− Just dec ->
66− let dec' = max dec 1 in
67− case is of
68− [0] -> fromChar '0' `mappend` fromChar '.' `mappend` mconcat (replicate dec' (fromChar '0'))
69− `mappend` fromChar 'e' `mappend` fromChar '0'
70− _ ->
71− let (ei, is') = roundTo base (dec'+1) is
72− (d:ds) = if ei > 0 then init is' else is'
73− in digit d `mappend` fromChar '.' `mappend` mconcat (map digit ds)
74− `mappend` fromChar 'e' `mappend` int (e - 1 + ei)
75− go FFFixed (is, e) = case decs of
76− Nothing
77− | e <= 0 -> fromChar '0' `mappend` fromChar '.' `mappend` mconcat (replicate (-e) (fromChar '0'))
78− `mappend` mconcat (map digit is)
79− | otherwise -> let g 0 rs = fromChar '.' `mappend` mk0 rs
80− g n [] = fromChar '0' `mappend` g (n-1) []
81− g n (r:rs) = digit r `mappend` g (n-1) rs
82− in g e is
83− Just dec ->
84− let dec' = max dec 0 in
85− if e >= 0 then
86− let (ei, is') = roundTo base (dec' + e) is
87− (ls,rs) = splitAt (e+ei) is'
88− in if null ls
89− then mk0 ls
90− else mk0 ls `mappend` (fromChar '.' `mappend` mconcat (map digit rs))
91− else
92− let (ei, is') = roundTo base dec' (replicate (-e) 0 ++ is)
93− d:ds = if ei > 0 then is' else 0:is'
94− in if null ds
95− then digit d
96− else digit d `mappend` (fromChar '.' `mappend` mconcat (map digit ds))
97−
98− mk0 [] = fromChar '0'
99− mk0 rs = mconcat (map digit rs)
67+ (f0, e0) = decodeFloat x
68+ (minExp0, _) = floatRange (undefined::Double)
69+ p = floatDigits x
70+ b = floatRadix x
71+ minExp = minExp0 - p -- the real minimum exponent
72+ -- Haskell requires that f be adjusted so denormalized numbers
73+ -- will have an impossibly low exponent. Adjust for this.
74+ (# f, e #) =
75+ let n = minExp - e0 in
76+ if n > 0 then (# f0 `div` (b^n), e0+n #) else (# f0, e0 #)
77+ (# r, s, mUp, mDn #) =
78+ if e >= 0
79+ then let be = b^ e
80+ in if f == b^(p-1)
81+ then (# f*be*b*2, 2*b, be*b, b #)
82+ else (# f*be*2, 2, be, be #)
83+ else if e > minExp && f == b^(p-1)
84+ then (# f*b*2, b^(-e+1)*2, b, 1 #)
85+ else (# f*2, b^(-e)*2, 1, 1 #)
86+ k = fixup k0
87+ where
88+ k0 | b == 2 = (p - 1 + e0) * 3 `div` 10
89+ -- logBase 10 2 is slightly bigger than 3/10 so the following
90+ -- will err on the low side. Ignoring the fraction will make
91+ -- it err even more. Haskell promises that p-1 <= logBase b f
92+ -- < p.
93+ | otherwise = ceiling ((log (fromInteger (f+1) :: Double) +
94+ fromIntegral e * log (fromInteger b)) / log 10)
95+ fixup n
96+ | n >= 0 = if r + mUp <= exp10 n * s then n else fixup (n+1)
97+ | otherwise = if exp10 (-n) * (r + mUp) <= s then n else fixup (n+1)
98+
99+ gen ds !rn !sN !mUpN !mDnN =
100+ let (dn0, rn') = (rn * 10) `divMod` sN
101+ mUpN' = mUpN * 10
102+ mDnN' = mDnN * 10
103+ !dn = fromInteger dn0
104+ !dn' = dn + 1
105+ in case (# rn' < mDnN', rn' + mUpN' > sN #) of
106+ (# True, False #) -> dn : ds
107+ (# False, True #) -> dn' : ds
108+ (# True, True #) -> if rn' * 2 < sN then dn : ds else dn' : ds
109+ (# False, False #) -> gen (dn:ds) rn' sN mUpN' mDnN'
110+
111+ rds | k >= 0 = gen [] r (s * exp10 k) mUp mDn
112+ | otherwise = gen [] (r * bk) s (mUp * bk) (mDn * bk)
113+ where bk = exp10 (-k)
114+
115+exp10 :: Int -> Integer
116+exp10 n
117+ | n >= 0 && n < maxExpt = V.unsafeIndex expts n
118+ | otherwise = 10 ^ n
119+ where expts = V.generate maxExpt (10^)
120+ {-# NOINLINE expts #-}
121+ maxExpt = 17
122+{-# INLINE exp10 #-}