diff --git a/Data/Number/LogFloat.hs b/Data/Number/LogFloat.hs deleted file mode 100644 --- a/Data/Number/LogFloat.hs +++ /dev/null @@ -1,336 +0,0 @@ - --- Needed by our RealToFrac contexts -{-# LANGUAGE FlexibleContexts #-} - --- Removed -Wall because -fno-warn-orphans was removed in GHC 6.10 -{-# OPTIONS_GHC -fwarn-tabs #-} - --- Unfortunately we need -fglasgow-exts in order to actually pick --- up on the rules (see -ddump-rules). The -frewrite-rules flag --- doesn't do what you want. --- cf --- cf -{-# OPTIONS_GHC -O2 -fvia-C -optc-O3 -fexcess-precision -fglasgow-exts #-} - --- Version History --- (v0.11) Broke Data.Number.RealToFrac out --- (v0.10) Fixed bugs in Hugs for PartialOrd and Transfinite. --- Also added maxPO, minPO, comparingPO --- (v0.9.1) Fixed some PartialOrd stuff and sanitized documentation --- (v0.9.0) s/toFractional/realToFrac/g. --- Also moved realToFrac and log to Transfinite --- (v0.8.6) Removed buggy RULES --- (v0.8.5) Gave up and converted from lhs to hs so Hackage docs work --- (v0.8.4) Broke out Transfinite --- (v0.8.3) Documentation updates --- (v0.8.2) Announced release --- (v0.8) Did a bunch of tweaking. Things should be decent now --- (v0.7) Haddockified --- (v0.6) Fixed monomorphism. --- (v0.5) Added optimization rules. --- (v0.4) Translated to Haskell at revision 2007.12.20. --- (v0.3) Converted extensive comments to POD format. --- (v0.2) Did a bunch of profiling, optimizing, and debugging. --- (v0.1) Initial version created for hw5 for NLP with Jason Eisner. --- ----------------------------------------------------------------- --- ~ 2008.08.29 --- | --- Module : Data.Number.LogFloat --- Copyright : Copyright (c) 2007--2009 wren ng thornton --- License : BSD3 --- Maintainer : wren@community.haskell.org --- Stability : stable --- Portability : portable --- --- This module presents a type for storing numbers in the log-domain. --- The main reason for doing this is to prevent underflow when --- multiplying many small probabilities as is done in Hidden Markov --- Models and other statistical models often used for natural --- language processing. The log-domain also helps prevent overflow --- when multiplying many large numbers. In rare cases it can speed --- up numerical computation (since addition is faster than --- multiplication, though logarithms are exceptionally slow), but --- the primary goal is to improve accuracy of results. A secondary --- goal has been to maximize efficiency since these computations --- are frequently done within a /O(n^3)/ loop. --- --- The 'LogFloat' of this module is restricted to non-negative --- numbers for efficiency's sake, see the forthcoming --- "Data.Number.LogFloat.Signed" for doing signed log-domain --- calculations. ----------------------------------------------------------------- - -module Data.Number.LogFloat - ( - -- * Exceptional numeric values - module Data.Number.Transfinite - , module Data.Number.RealToFrac - - -- * @LogFloat@ data type and conversion functions - , LogFloat - , logFloat, logToLogFloat - , fromLogFloat, logFromLogFloat - ) where - -import Prelude hiding (log, realToFrac, isInfinite, isNaN) - -import Data.Number.RealToFrac -import Data.Number.Transfinite -import Data.Number.PartialOrd - ----------------------------------------------------------------- --- --- Try to add in some optimizations. Why these need to be --- down here and localized to the module, I don't know. We don't --- do anything foolish like this, but our clients might, or they --- might be generated by other code transformations. Note that due --- to the fuzz, these equations are not strictly true, even though --- they are mathematically correct. - -{-# RULES -"log/exp" forall x. log (exp x) = x -"log.exp" log . exp = id - -"exp/log" forall x. exp (log x) = x -"exp.log" exp . log = id - #-} - --- We'd like to be able to take advantage of general rule versions --- of our operators for 'LogFloat', with rules like @log x + log y --- = log (x * y)@ and @log x - log y = log (x / y)@. However the --- problem is that those equations could be driven in either direction --- depending on whether we think time performance or non-underflow --- performance is more important, and the answers may be different --- at every call site. --- --- Since we implore users to do normal-domain computations whenever --- it would not degenerate accuracy, we should not rewrite their --- decisions in any way. The log\/exp fusion strictly improves both --- time and accuracy, so those are safe. But the buck stops with --- them. - - --- These should only fire when it's type-safe --- This should already happen, but... --- TODO: Check the logs to see if it ever fires --- N.B. these are orphaned -{-# RULES -"toRational/fromRational" forall x. toRational (fromRational x) = x -"toRational.fromRational" toRational . fromRational = id - #-} - - ----------------------------------------------------------------- - --- | Reduce the number of constant string literals we need to store. -errorOutOfRange :: String -> a -errorOutOfRange fun = error $! "Data.Number.LogFloat."++fun - ++ ": argument out of range" - - --- | We need these guards in order to ensure some invariants. -guardNonNegative :: String -> Double -> Double -guardNonNegative fun x | x >= 0 = x - | otherwise = errorOutOfRange fun - - --- | It's unfortunate that 'notANumber' is not equal to itself, but --- we can hack around that. GHC gives NaN for the log of negatives --- and so we could ideally take advantage of @log . guardNonNegative --- fun = guardIsANumber fun . log@ to simplify things, but Hugs --- raises an error so that's non-portable. -guardIsANumber :: String -> Double -> Double -guardIsANumber fun x | isNaN x = errorOutOfRange fun - | otherwise = x - ----------------------------------------------------------------- --- --- | A @LogFloat@ is just a 'Double' with a special interpretation. --- The 'logFloat' function is presented instead of the constructor, --- in order to ensure semantic conversion. At present the 'Show' --- instance will convert back to the normal-domain, and so will --- underflow at that point. This behavior may change in the future. --- --- Performing operations in the log-domain is cheap, prevents --- underflow, and is otherwise very nice for dealing with miniscule --- probabilities. However, crossing into and out of the log-domain --- is expensive and should be avoided as much as possible. In --- particular, if you're doing a series of multiplications as in --- @lp * logFloat q * logFloat r@ it's faster to do @lp * logFloat --- (q * r)@ if you're reasonably sure the normal-domain multiplication --- won't underflow, because that way you enter the log-domain only --- once, instead of twice. --- --- Even more particularly, you should /avoid addition/ whenever --- possible. Addition is provided because it's necessary at times --- and the proper implementation is not immediately transparent. --- However, between two @LogFloat@s addition requires crossing the --- exp\/log boundary twice; with a @LogFloat@ and a regular number --- it's three times since the regular number needs to enter the --- log-domain first. This makes addition incredibly slow. Again, --- if you can parenthesize to do plain operations first, do it! - -newtype LogFloat = LogFloat Double - deriving (Eq, Ord) -- Should we really perpetuate the Ord lie? - -instance PartialOrd LogFloat where - cmp (LogFloat x) (LogFloat y) - | isNaN x || isNaN y = Nothing - | otherwise = Just $! x `compare` y - - --- | A constructor which does semantic conversion from normal-domain --- to log-domain. -logFloat :: (Real a, RealToFrac a Double) => a -> LogFloat -{-# SPECIALIZE logFloat :: Double -> LogFloat #-} -logFloat = LogFloat . log . guardNonNegative "logFloat" . realToFrac - - --- This is simply a polymorphic version of the 'LogFloat' data --- constructor. We present it mainly because we hide the constructor --- in order to make the type a bit more opaque. If the polymorphism --- turns out to be a performance liability because the rewrite rules --- can't remove it, then we need to rethink all four --- constructors\/destructors. --- --- | Constructor which assumes the argument is already in the --- log-domain. -logToLogFloat :: (Real a, RealToFrac a Double) => a -> LogFloat -{-# SPECIALIZE logToLogFloat :: Double -> LogFloat #-} -logToLogFloat = LogFloat . guardIsANumber "logToLogFloat" . realToFrac - - --- | Return our log-domain value back into normal-domain. Beware --- of overflow\/underflow. -fromLogFloat :: (Fractional a, Transfinite a, RealToFrac Double a) - => LogFloat -> a -{-# SPECIALIZE fromLogFloat :: LogFloat -> Double #-} -fromLogFloat (LogFloat x) = realToFrac (exp x) - - --- | Return the log-domain value itself without costly conversion -logFromLogFloat :: (Fractional a, Transfinite a, RealToFrac Double a) - => LogFloat -> a -{-# SPECIALIZE logFromLogFloat :: LogFloat -> Double #-} -logFromLogFloat (LogFloat x) = realToFrac x - - --- These are our module-specific versions of "log\/exp" and "exp\/log"; --- They do the same things but also have a @LogFloat@ in between --- the logarithm and exponentiation. --- --- In order to ensure these rules fire we may need to delay inlining --- of the four con-\/destructors, like we do for 'realToFrac'. --- Unfortunately, I'm not entirely sure whether they will be inlined --- already or not (and whether they are may be fragile) and I don't --- want to inline them excessively and lead to code bloat in the --- off chance that we could prune some of it away. --- TODO: thoroughly investigate this. - -{-# RULES --- Out of log-domain and back in -"log/fromLogFloat" forall x. log (fromLogFloat x) = logFromLogFloat x -"log.fromLogFloat" log . fromLogFloat = logFromLogFloat - -"logFloat/fromLogFloat" forall x. logFloat (fromLogFloat x) = x -"logFloat.fromLogFloat" logFloat . fromLogFloat = id - --- Into log-domain and back out -"fromLogFloat/logFloat" forall x. fromLogFloat (logFloat x) = x -"fromLogFloat.logFloat" fromLogFloat . logFloat = id - #-} - ----------------------------------------------------------------- --- To show it, we want to show the normal-domain value rather than --- the log-domain value. Also, if someone managed to break our --- invariants (e.g. by passing in a negative and noone's pulled on --- the thunk yet) then we want to crash before printing the --- constructor, rather than after. N.B. This means the show will --- underflow\/overflow in the same places as normal doubles since --- we underflow at the @exp@. Perhaps this means we should show the --- log-domain value instead. - -instance Show LogFloat where - show (LogFloat x) = let y = exp x - in y `seq` "LogFloat "++show y - - ----------------------------------------------------------------- --- These all work without causing underflow. However, do note that --- they tend to induce more of the floating-point fuzz than using --- regular floating numbers because @exp . log@ doesn't really equal --- @id@. In any case, our main aim is for preventing underflow when --- multiplying many small numbers (and preventing overflow for --- multiplying many large numbers) so we're not too worried about --- +\/- 4e-16. - -instance Num LogFloat where - -- BUG? In Hugs (Sept2006) the (>=) always returns True if - -- either isNaN. Only questionably a bug, since we try to - -- ensure that notANumber never occurs. Still... perhaps - -- we should use `ge` and other PartialOrd things in order - -- to play it safe. - -- TODO: benchmark and check core to see how much that hurts GHC. - - - (*) (LogFloat x) (LogFloat y) = LogFloat (x+y) - - (+) (LogFloat x) (LogFloat y) - | x >= y = LogFloat (x + log (1 + exp (y - x))) - | otherwise = LogFloat (y + log (1 + exp (x - y))) - - -- Without the guard this would return NaN instead of error - (-) (LogFloat x) (LogFloat y) - | x >= y = LogFloat (x + log (1 - exp (y - x))) - | otherwise = errorOutOfRange "(-)" - - signum (LogFloat x) - | x == negativeInfinity = 0 - | x > negativeInfinity = 1 - | otherwise = errorOutOfRange "signum" - -- The extra guard protects against NaN, in case someone - -- broke the invariant. That shouldn't be possible and - -- so noone else bothers to check, but we check here just - -- in case. - - negate _ = errorOutOfRange "negate" - - abs = id - - fromInteger = LogFloat . log - . guardNonNegative "fromInteger" . fromInteger - - -instance Fractional LogFloat where - -- n/0 is handled seamlessly for us; we must catch 0/0 though - (/) (LogFloat x) (LogFloat y) - | x == negativeInfinity - && y == negativeInfinity = errorOutOfRange "(/)" -- protect vs NaN - | otherwise = LogFloat (x-y) - - fromRational = LogFloat . log - . guardNonNegative "fromRational" . fromRational - - --- Just for fun. The more coersion functions the better. Though --- Rationals are very buggy when it comes to transfinite values -instance Real LogFloat where - toRational (LogFloat x) = toRational (exp x) - - -{- -- Commented out because I'm not sure about requiring MPTCs. Of course, those are already required by "Data.Number.Transfinite" so it's pretty moot... - --- LogFloat->LogFloat is already given via generic (a->a) --- No LogFloat->Rational since LogFloat can have 'infinity' --- Can't have LogFloat->a using fromLogFloat because Hugs dislikes incoherence. Adding an explicit LogFloat->LogFloat instance doesn't help like it does for GHC. - -instance RealToFrac LogFloat Double where - realToFrac = fromLogFloat - -instance RealToFrac LogFloat Float where - realToFrac = fromLogFloat --} - ----------------------------------------------------------------- ------------------------------------------------------------ fin. diff --git a/Data/Number/PartialOrd.hs b/Data/Number/PartialOrd.hs deleted file mode 100644 --- a/Data/Number/PartialOrd.hs +++ /dev/null @@ -1,148 +0,0 @@ -{-# LANGUAGE OverlappingInstances - , FlexibleInstances - , UndecidableInstances - #-} - -{-# OPTIONS_GHC -Wall -fwarn-tabs #-} - ----------------------------------------------------------------- --- ~ 2009.01.29 --- | --- Module : Data.Number.PartialOrd --- Copyright : Copyright (c) 2007--2009 wren ng thornton --- License : BSD3 --- Maintainer : wren@community.haskell.org --- Stability : experimental --- Portability : semi-portable (overlapping instances, etc) --- --- The Prelude's 'Ord' class for dealing with ordered types is often --- onerous to use because it requires 'Eq' as well as a total --- ordering. While such total orderings are common, partial orderings --- are moreso. This module presents a class for partially ordered --- types. ----------------------------------------------------------------- -module Data.Number.PartialOrd - ( - -- * Partial Ordering - PartialOrd(..) - -- * Functions - , comparingPO - ) where - --- Bugfix for Hugs (September 2006), see note below. -import Prelude hiding (isNaN) -import Hugs.RealFloat (isNaN) - ----------------------------------------------------------------- --- | This class defines a partially ordered type. The method names --- were chosen so as not to conflict with 'Ord' and 'Eq'. We use --- 'Maybe' instead of defining new types @PartialOrdering@ and --- @FuzzyBool@ because this way should make the class easier to --- use. --- --- Minimum complete definition: 'cmp' - -class PartialOrd a where - -- | like 'compare' - cmp :: a -> a -> Maybe Ordering - - -- | like ('>') - gt :: a -> a -> Maybe Bool - gt x y = case x `cmp` y of - Just GT -> Just True - Just _ -> Just False - Nothing -> Nothing - - -- | like ('>=') - ge :: a -> a -> Maybe Bool - ge x y = case x `cmp` y of - Just LT -> Just False - Just _ -> Just True - Nothing -> Nothing - - -- | like ('==') - eq :: a -> a -> Maybe Bool - eq x y = case x `cmp` y of - Just EQ -> Just True - Just _ -> Just False - Nothing -> Nothing - - -- | like ('/=') - ne :: a -> a -> Maybe Bool - ne x y = case x `cmp` y of - Just EQ -> Just False - Just _ -> Just True - Nothing -> Nothing - - -- | like ('<=') - le :: a -> a -> Maybe Bool - le x y = case x `cmp` y of - Just GT -> Just False - Just _ -> Just True - Nothing -> Nothing - - -- | like ('<') - lt :: a -> a -> Maybe Bool - lt x y = case x `cmp` y of - Just LT -> Just True - Just _ -> Just False - Nothing -> Nothing - - -- | like 'max'. The default instance returns the left argument - -- when they're equal. - maxPO :: a -> a -> Maybe a - maxPO x y = do o <- x `cmp` y - case o of - GT -> Just x - EQ -> Just x - LT -> Just y - - -- | like 'min'. The default instance returns the left argument - -- when they're equal. - minPO :: a -> a -> Maybe a - minPO x y = do o <- x `cmp` y - case o of - GT -> Just y - EQ -> Just x - LT -> Just x - -infix 4 `gt`, `ge`, `eq`, `ne`, `le`, `lt`, `maxPO`, `minPO` - -instance (Ord a) => PartialOrd a where - cmp x y = Just $! x `compare` y - gt x y = Just $! x > y - ge x y = Just $! x >= y - eq x y = Just $! x == y - ne x y = Just $! x /= y - le x y = Just $! x <= y - lt x y = Just $! x < y - maxPO x y = Just $! x `max` y - minPO x y = Just $! x `min` y - - --- N.B. Hugs (Sept 2006) has a buggy definition for 'isNaN' which --- always returns @False@. We use a fixed version, provided the CPP --- was run with the right arguments. See "Hugs.RealFloat". If 'cmp' --- returns @Just Eq@ for @notANumber@ then CPP was run wrongly. --- --- The instances inherited from Ord are wrong. So we'll fix them. -instance PartialOrd Float where - cmp x y | isNaN x || isNaN y = Nothing - | otherwise = Just $! x `compare` y - -instance PartialOrd Double where - cmp x y | isNaN x || isNaN y = Nothing - | otherwise = Just $! x `compare` y - ----------------------------------------------------------------- --- TODO? add maximumPO\/minimumPO via left or right fold? - --- BUG: Haddock doesn't link the `comparing` --- --- | Like @Data.Ord.comparing@. Helpful in conjunction with the --- @xxxBy@ family of functions from "Data.List" -comparingPO :: (PartialOrd b) => (a -> b) -> a -> a -> Maybe Ordering -comparingPO p x y = p x `cmp` p y - ----------------------------------------------------------------- ------------------------------------------------------------ fin. diff --git a/Data/Number/RealToFrac.hs b/Data/Number/RealToFrac.hs deleted file mode 100644 --- a/Data/Number/RealToFrac.hs +++ /dev/null @@ -1,116 +0,0 @@ --- Needed to ensure correctness, and because we can't guarantee rules fire -{-# LANGUAGE MultiParamTypeClasses - , OverlappingInstances - , CPP - #-} - --- Glasgow extensions needed to enable the # kind -{-# OPTIONS_GHC -fglasgow-exts #-} - -{-# OPTIONS_GHC -Wall -fwarn-tabs #-} - ----------------------------------------------------------------- --- ~ 2009.01.29 --- | --- Module : Data.Number.RealToFrac --- Copyright : Copyright (c) 2007--2009 wren ng thornton --- License : BSD3 --- Maintainer : wren@community.haskell.org --- Stability : experimental --- Portability : non-portable (CPP, MPTC, OverlappingInstances) --- --- This module presents a type class for generic conversion between --- numeric types, generalizing @realToFrac@ in order to overcome --- problems with pivoting through 'Rational' ----------------------------------------------------------------- -module Data.Number.RealToFrac (RealToFrac(..)) where - -import Prelude hiding (realToFrac, isInfinite, isNaN) -import qualified Prelude (realToFrac) - -import Data.Number.Transfinite - -#ifdef __GLASGOW_HASKELL__ -import GHC.Exts - ( Int(..), Float(..), Double(..) - , int2Double# - , int2Float# - , double2Float# - , float2Double# - ) -#endif - ----------------------------------------------------------------- --- | The 'Prelude.realToFrac' function is defined to pivot through --- a 'Rational' according to the haskell98 spec. This is non-portable --- and problematic as discussed in "Data.Number.Transfinite". Since --- there is resistance to breaking from the spec, this class defines --- a reasonable variant which deals with transfinite values --- appropriately. --- --- There is a generic instance from any Transfinite Real to any --- Transfinite Fractional, using checks to ensure correctness. GHC --- has specialized versions for some types which use primitive --- converters instead, for large performance gains. (These definitions --- are hidden from other compilers via CPP.) Due to a bug in Haddock --- the specialized instances are shown twice and the generic instance --- isn't shown at all. Since the instances are overlapped, you'll --- need to give type signatures if the arguments to 'realToFrac' --- are polymorphic. There's also a generic instance for any Real --- Fractional type to itself, thus if you write any generic instances --- beware of incoherence. --- --- If any of these restrictions (CPP, GHC-only optimizations, --- OverlappingInstances) are onerous to you, contact the maintainer --- (we like patches). Note that this /does/ work for Hugs with --- suitable options (e.g. @hugs -98 +o -F'cpp -P'@). However, Hugs --- doesn't allow @IncoherentInstances@ nor does it allow diamonds --- with @OverlappingInstances@, which restricts the ability to add --- additional generic instances. - -class (Real a, Fractional b) => RealToFrac a b where - realToFrac :: a -> b - -instance (Real a, Fractional a) => RealToFrac a a where - realToFrac = id - -instance (Real a, Transfinite a, Fractional b, Transfinite b) - => RealToFrac a b - where - realToFrac x - | isNaN x = notANumber - | isInfinite x = if x > 0 then infinity - else negativeInfinity - | otherwise = Prelude.realToFrac x - - -#ifdef __GLASGOW_HASKELL__ -instance RealToFrac Int Float where - {-# INLINE realToFrac #-} - realToFrac (I# i) = F# (int2Float# i) - -instance RealToFrac Int Double where - {-# INLINE realToFrac #-} - realToFrac (I# i) = D# (int2Double# i) - - -instance RealToFrac Integer Float where - -- TODO: is there a more primitive way? - realToFrac j = Prelude.realToFrac j - -instance RealToFrac Integer Double where - -- TODO: is there a more primitive way? - realToFrac j = Prelude.realToFrac j - - -instance RealToFrac Float Double where - {-# INLINE realToFrac #-} - realToFrac (F# f) = D# (float2Double# f) - -instance RealToFrac Double Float where - {-# INLINE realToFrac #-} - realToFrac (D# d) = F# (double2Float# d) -#endif - ----------------------------------------------------------------- ------------------------------------------------------------ fin. diff --git a/Data/Number/Transfinite.hs b/Data/Number/Transfinite.hs deleted file mode 100644 --- a/Data/Number/Transfinite.hs +++ /dev/null @@ -1,187 +0,0 @@ --- Glasgow extensions needed to enable the # kind -{-# OPTIONS_GHC -fglasgow-exts #-} - -{-# OPTIONS_GHC -Wall -fwarn-tabs #-} - ----------------------------------------------------------------- --- ~ 2009.01.29 --- | --- Module : Data.Number.Transfinite --- Copyright : Copyright (c) 2007--2009 wren ng thornton --- License : BSD3 --- Maintainer : wren@community.haskell.org --- Stability : experimental --- Portability : portable --- --- This module presents a type class for numbers which have --- representations for transfinite values. The idea originated from --- the IEEE-754 floating-point special values, used by --- "Data.Number.LogFloat". However not all 'Fractional' types --- necessarily support transfinite values. In particular, @Ratio@ --- types including 'Rational' do not have portable representations. --- --- For the Glasgow compiler (GHC 6.8.2), "GHC.Real" defines @1%0@ --- and @0%0@ as representations for 'infinity' and 'notANumber', --- but most operations on them will raise exceptions. If 'toRational' --- is used on an infinite floating value, the result is a rational --- with a numerator sufficiently large that it will overflow when --- converted back to a @Double@. If used on NaN, the result would --- buggily convert back as 'negativeInfinity'. For more discussion --- on why this approach is problematic, see: --- --- * --- --- * --- --- Hugs (September 2006) stays closer to the haskell98 spec and --- offers no way of constructing those values, raising arithmetic --- overflow errors if attempted. ----------------------------------------------------------------- -module Data.Number.Transfinite - ( Transfinite(..) - , log - ) where - -import Prelude hiding (log, isInfinite, isNaN) -import qualified Prelude (log) -import qualified Hugs.RealFloat as Prelude (isInfinite, isNaN) - -import Data.Number.PartialOrd - ----------------------------------------------------------------- --- | Many numbers are not 'Bounded' yet, even though they can --- represent arbitrarily large values, they are not necessarily --- able to represent transfinite values such as infinity itself. --- This class is for types which are capable of representing such --- values. Notably, this class does not require the type to be --- 'Fractional' nor 'Floating' since integral types could also have --- representations for transfinite values. By popular demand the --- 'Num' restriction has been lifted as well, due to complications --- of defining 'Show' or 'Eq' for some types. --- --- In particular, this class extends the ordered projection to have --- a maximum value 'infinity' and a minimum value 'negativeInfinity', --- as well as an exceptional value 'notANumber'. All the natural --- laws regarding @infinity@ and @negativeInfinity@ should pertain. --- (Some of these are discussed below.) --- --- Hugs (September 2006) has buggy Prelude definitions for --- 'Prelude.isNaN' and 'Prelude.isInfinite' on Float and Double. --- This module provides correct definitions, so long as "Hugs.RealFloat" --- is compiled correctly. - -class (PartialOrd a) => Transfinite a where - - -- | A transfinite value which is greater than all finite values. - -- Adding or subtracting any finite value is a no-op. As is - -- multiplying by any non-zero positive value (including - -- @infinity@), and dividing by any positive finite value. Also - -- obeys the law @negate infinity = negativeInfinity@ with all - -- appropriate ramifications. - - infinity :: a - - - -- | A transfinite value which is less than all finite values. - -- Obeys all the same laws as @infinity@ with the appropriate - -- changes for the sign difference. - - negativeInfinity :: a - - - -- | An exceptional transfinite value for dealing with undefined - -- results when manipulating infinite values. The following - -- operations must return @notANumber@, where @inf@ is any value - -- which @isInfinite@: - -- - -- * @inf + inf@ - -- - -- * @inf - inf@ - -- - -- * @inf * 0@ - -- - -- * @0 * inf@ - -- - -- * @inf \/ inf@ - -- - -- * @inf `div` inf@ - -- - -- * @0 \/ 0@ - -- - -- * @0 `div` 0@ - -- - -- Additionally, any mathematical operations on @notANumber@ - -- must also return @notANumber@, and any equality or ordering - -- comparison on @notANumber@ must return @False@ (violating - -- the law of the excluded middle, often assumed but not required - -- for 'Eq'; thus, 'eq' and 'ne' are preferred over ('==') and - -- ('/=')). Since it returns false for equality, there may be - -- more than one machine representation of this `value'. - - notANumber :: a - - - -- | Return true for both @infinity@ and @negativeInfinity@, - -- false for all other values. - isInfinite :: a -> Bool - - -- | Return true only for @notANumber@. - isNaN :: a -> Bool - - -instance Transfinite Double where - infinity = 1/0 - negativeInfinity = negate (1/0) - notANumber = 0/0 - isInfinite = Prelude.isInfinite - isNaN = Prelude.isNaN - - -instance Transfinite Float where - infinity = 1/0 - negativeInfinity = negate (1/0) - notANumber = 0/0 - isInfinite = Prelude.isInfinite - isNaN = Prelude.isNaN - - ----------------------------------------------------------------- --- | Since the normal 'Prelude.log' throws an error on zero, we --- have to redefine it in order for things to work right. Arguing --- from limits we can see that @log 0 == negativeInfinity@. Newer --- versions of GHC have this behavior already, but older versions --- and Hugs do not. --- --- This function will raise an error when taking the log of negative --- numbers, rather than returning 'notANumber' as the newer GHC --- implementation does. The reason being that typically this is a --- logical error, and @notANumber@ allows the error to propegate --- silently. --- --- In order to improve portability, the 'Transfinite' class is --- required to indicate that the 'Floating' type does in fact have --- a representation for negative infinity. Both native floating --- types ('Double' and 'Float') are supported. If you define your --- own instance of @Transfinite@, verify the above equation holds --- for your @0@ and @negativeInfinity@. If it doesn't, then you --- should avoid importing our @log@ and will probably want converters --- to handle the discrepancy. - -log :: (Floating a, Transfinite a) => a -> a -{-# SPECIALIZE log :: Double -> Double #-} -{-# SPECIALIZE log :: Float -> Float #-} -log x = case x `cmp` 0 of - Just GT -> Prelude.log x - Just EQ -> negativeInfinity - Just LT -> err "argument out of range" - Nothing -> err "argument not comparable to 0" - where - err e = error $! "Data.Number.Transfinite.log: "++e - --- Note, Floating ultimately requires Num, but not Ord. If PartialOrd --- proves to be an onerous requirement on Transfinite, we could --- hack our way around without using PartialOrd by using isNaN, (== --- 0), ((>0).signum) but that would be less efficient. - ----------------------------------------------------------------- ------------------------------------------------------------ fin. diff --git a/Hugs/RealFloat.hs b/Hugs/RealFloat.hs deleted file mode 100644 --- a/Hugs/RealFloat.hs +++ /dev/null @@ -1,63 +0,0 @@ - -{-# LANGUAGE CPP #-} - -{-# OPTIONS_GHC -Wall -fwarn-tabs #-} - ----------------------------------------------------------------- --- ~ 2009.01.29 --- | --- Module : Hugs.RealFloat --- Copyright : Copyright (c) 2007--2009 wren ng thornton --- License : BSD3 --- Maintainer : wren@community.haskell.org --- Stability : stable --- Portability : portable (with CPP) --- --- Hugs (September 2006) has buggy definitions for 'Prelude.isNaN' --- and 'Prelude.isInfinite' on Float and Double. If this module is --- run through CPP with the macro @__HUGS__@ set to a value no --- larger than 200609, then correct definitions are used. Otherwise --- the Prelude definitions are used (which should be correct for --- other compilers). For example, run Hugs with --- --- @hugs -F'cpp -P -D__HUGS__=200609' Hugs/RealFloat.hs@ --- --- N.B. The corrected definitions have only been tested to work for --- 'Float' and 'Double'. These definitions should probably not be --- used for other 'RealFloat' types. ----------------------------------------------------------------- -module Hugs.RealFloat - ( isInfinite - , isNaN - ) where - -import Prelude hiding (isInfinite, isNaN) -import qualified Prelude ----------------------------------------------------------------- - -isInfinite :: (RealFloat a) => a -> Bool -{-# INLINE isInfinite #-} -#if defined(__HUGS__) && (__HUGS__ <= 200609) -isInfinite x = (1/0) == abs x -#else -isInfinite = Prelude.isInfinite -#endif - - -isNaN :: (RealFloat a) => a -> Bool -{-# INLINE isNaN #-} -#if defined(__HUGS__) && (__HUGS__ <= 200609) -isNaN x = compareEQ x 0 && compareEQ x 1 - --- | In Hugs (September 2006), 'compare' always returns @EQ@ if one --- of the arguments is not a number. Thus, if a number is @compareEQ@ --- against multiple different numbers, then it must be @isNaN@. -compareEQ :: (Ord a) => a -> a -> Bool -compareEQ x y = case compare x y of - EQ -> True - _ -> False -#else -isNaN = Prelude.isNaN -#endif ----------------------------------------------------------------- ------------------------------------------------------------ fin. diff --git a/logfloat.cabal b/logfloat.cabal --- a/logfloat.cabal +++ b/logfloat.cabal @@ -3,7 +3,7 @@ ---------------------------------------------------------------- Name: logfloat -Version: 0.11.0.1 +Version: 0.11.1 Cabal-Version: >= 1.2 Build-Type: Simple Stability: experimental @@ -21,26 +21,25 @@ probabilities as is done in Hidden Markov Models. It is also helpful for preventing overflow. --- Doing it this way uncovers a bug in Cabal-1.2.3.0: --- "Setup.hs: 'parseField' called on a non-field. This is a bug." ---Flag hiddenPrim --- Description: Use GHC 6.10's newly hidden package for GHC.Prim --- if impl(ghc >= 6.10) --- Default: true --- else --- Default: false +Flag splitBase + Description: base-3.0 broke out array and other packages + Default: False + + Library + Hs-Source-Dirs: src Exposed-Modules: Data.Number.LogFloat , Data.Number.RealToFrac , Data.Number.Transfinite , Data.Number.PartialOrd , Hugs.RealFloat - Build-Depends: base --- No longer needed since we use GHC.Exts instead --- if flag(hiddenPrim) --- Build-Depends: ghc-prim - Hugs-Options: -98 +o -F'cpp -P -D__HUGS__=200609' + if flag(splitBase) + Build-depends: base >= 3.0, array + else + Build-depends: base < 3.0 + + Hugs-Options: -98 +o -F'cpp -P -traditional -D__HUGS__=200609' if impl(ghc < 6.10) GHC-Options: -fno-warn-orphans diff --git a/src/Data/Number/LogFloat.hs b/src/Data/Number/LogFloat.hs new file mode 100644 --- /dev/null +++ b/src/Data/Number/LogFloat.hs @@ -0,0 +1,437 @@ + +-- FlexibleContexts needed by our RealToFrac contexts +-- CPP needed for IArray UArray instance +{-# LANGUAGE FlexibleContexts + , CPP #-} + +-- Removed -Wall because -fno-warn-orphans was removed in GHC 6.10 +{-# OPTIONS_GHC -fwarn-tabs #-} + +-- Unfortunately we need -fglasgow-exts in order to actually pick +-- up on the rules (see -ddump-rules). The -frewrite-rules flag +-- doesn't do what you want. +-- cf +-- cf +{-# OPTIONS_GHC -O2 -fvia-C -optc-O3 -fexcess-precision -fglasgow-exts #-} + +-- Version History +-- (v0.11.1) Added IArray UArray instance +-- (v0.11) Broke Data.Number.RealToFrac out +-- (v0.10) Fixed bugs in Hugs for PartialOrd and Transfinite. +-- Also added maxPO, minPO, comparingPO +-- (v0.9.1) Fixed some PartialOrd stuff and sanitized documentation +-- (v0.9.0) s/toFractional/realToFrac/g. +-- Also moved realToFrac and log to Transfinite +-- (v0.8.6) Removed buggy RULES +-- (v0.8.5) Gave up and converted from lhs to hs so Hackage docs work +-- (v0.8.4) Broke out Transfinite +-- (v0.8.3) Documentation updates +-- (v0.8.2) Announced release +-- (v0.8) Did a bunch of tweaking. Things should be decent now +-- (v0.7) Haddockified +-- (v0.6) Fixed monomorphism. +-- (v0.5) Added optimization rules. +-- (v0.4) Translated to Haskell at revision 2007.12.20. +-- (v0.3) Converted extensive comments to POD format. +-- (v0.2) Did a bunch of profiling, optimizing, and debugging. +-- (v0.1) Initial version created for hw5 for NLP with Jason Eisner. +-- +---------------------------------------------------------------- +-- ~ 2009.03.07 +-- | +-- Module : Data.Number.LogFloat +-- Copyright : Copyright (c) 2007--2009 wren ng thornton +-- License : BSD3 +-- Maintainer : wren@community.haskell.org +-- Stability : stable +-- Portability : portable (with CPP) +-- +-- This module presents a type for storing numbers in the log-domain. +-- The main reason for doing this is to prevent underflow when +-- multiplying many small probabilities as is done in Hidden Markov +-- Models and other statistical models often used for natural +-- language processing. The log-domain also helps prevent overflow +-- when multiplying many large numbers. In rare cases it can speed +-- up numerical computation (since addition is faster than +-- multiplication, though logarithms are exceptionally slow), but +-- the primary goal is to improve accuracy of results. A secondary +-- goal has been to maximize efficiency since these computations +-- are frequently done within a /O(n^3)/ loop. +-- +-- The 'LogFloat' of this module is restricted to non-negative +-- numbers for efficiency's sake, see the forthcoming +-- "Data.Number.LogFloat.Signed" for doing signed log-domain +-- calculations. +---------------------------------------------------------------- + +module Data.Number.LogFloat + ( + -- * Exceptional numeric values + module Data.Number.Transfinite + , module Data.Number.RealToFrac + + -- * @LogFloat@ data type and conversion functions + , LogFloat + , logFloat, logToLogFloat + , fromLogFloat, logFromLogFloat + ) where + +import Prelude hiding (log, realToFrac, isInfinite, isNaN) + +import Data.Number.RealToFrac +import Data.Number.Transfinite +import Data.Number.PartialOrd + + +-- GHC can derive (IArray UArray LogFloat), but Hugs needs to coerce +-- TODO: see about nhc98/yhc, jhc/lhc +import Data.Array.Base (IArray(..)) +import Data.Array.Unboxed (UArray) + +-- Hugs (Sept 2006) doesn't use the generic wrapper in base:Unsafe.Coerce +-- so we'll just have to go back to the original source. +#ifdef __HUGS__ +import Hugs.IOExts (unsafeCoerce) +#elif __NHC__ +import NonStdUnsafeCoerce (unsafeCoerce) +#endif + +---------------------------------------------------------------- +-- +-- Try to add in some optimizations. Why these need to be +-- down here and localized to the module, I don't know. We don't +-- do anything foolish like this, but our clients might, or they +-- might be generated by other code transformations. Note that due +-- to the fuzz, these equations are not strictly true, even though +-- they are mathematically correct. + +{-# RULES +"log/exp" forall x. log (exp x) = x +"log.exp" log . exp = id + +"exp/log" forall x. exp (log x) = x +"exp.log" exp . log = id + #-} + +-- We'd like to be able to take advantage of general rule versions +-- of our operators for 'LogFloat', with rules like @log x + log y +-- = log (x * y)@ and @log x - log y = log (x / y)@. However the +-- problem is that those equations could be driven in either direction +-- depending on whether we think time performance or non-underflow +-- performance is more important, and the answers may be different +-- at every call site. +-- +-- Since we implore users to do normal-domain computations whenever +-- it would not degenerate accuracy, we should not rewrite their +-- decisions in any way. The log\/exp fusion strictly improves both +-- time and accuracy, so those are safe. But the buck stops with +-- them. + + +-- These should only fire when it's type-safe +-- This should already happen, but... +-- TODO: Check the logs to see if it ever fires +-- N.B. these are orphaned +{-# RULES +"toRational/fromRational" forall x. toRational (fromRational x) = x +"toRational.fromRational" toRational . fromRational = id + #-} + + +---------------------------------------------------------------- + +-- | Reduce the number of constant string literals we need to store. +errorOutOfRange :: String -> a +errorOutOfRange fun = error $! "Data.Number.LogFloat."++fun + ++ ": argument out of range" + + +-- | We need these guards in order to ensure some invariants. +guardNonNegative :: String -> Double -> Double +guardNonNegative fun x | x >= 0 = x + | otherwise = errorOutOfRange fun + + +-- TODO: since we're using Hugs.RealFloat instead of Prelude now, +-- is it still non-portable? +-- +-- | It's unfortunate that 'notANumber' is not equal to itself, but +-- we can hack around that. GHC gives NaN for the log of negatives +-- and so we could ideally take advantage of @log . guardNonNegative +-- fun = guardIsANumber fun . log@ to simplify things, but Hugs +-- raises an error so that's non-portable. +guardIsANumber :: String -> Double -> Double +guardIsANumber fun x | isNaN x = errorOutOfRange fun + | otherwise = x + +---------------------------------------------------------------- +-- +-- | A @LogFloat@ is just a 'Double' with a special interpretation. +-- The 'logFloat' function is presented instead of the constructor, +-- in order to ensure semantic conversion. At present the 'Show' +-- instance will convert back to the normal-domain, and so will +-- underflow at that point. This behavior may change in the future. +-- +-- Performing operations in the log-domain is cheap, prevents +-- underflow, and is otherwise very nice for dealing with miniscule +-- probabilities. However, crossing into and out of the log-domain +-- is expensive and should be avoided as much as possible. In +-- particular, if you're doing a series of multiplications as in +-- @lp * logFloat q * logFloat r@ it's faster to do @lp * logFloat +-- (q * r)@ if you're reasonably sure the normal-domain multiplication +-- won't underflow, because that way you enter the log-domain only +-- once, instead of twice. +-- +-- Even more particularly, you should /avoid addition/ whenever +-- possible. Addition is provided because it's necessary at times +-- and the proper implementation is not immediately transparent. +-- However, between two @LogFloat@s addition requires crossing the +-- exp\/log boundary twice; with a @LogFloat@ and a regular number +-- it's three times since the regular number needs to enter the +-- log-domain first. This makes addition incredibly slow. Again, +-- if you can parenthesize to do plain operations first, do it! + +newtype LogFloat = LogFloat Double + deriving + ( Eq + , Ord -- Should we really perpetuate the Ord lie? +#ifdef __GLASGOW_HASKELL__ + , IArray UArray + -- At least GHC 6.8.2 can derive IArray UArray (without + -- GeneralizedNewtypeDeriving). The H98 Report doesn't include + -- that among the options for automatic derivation though. +#endif + ) + + +#if __HUGS__ || __NHC__ + +-- These two operators make it much easier to read the instance. +-- Hopefully inlining everything will get rid of the eta overhead. +-- +{-# INLINE (~>) #-} +infixr 2 ~> +f ~> g = (. f) . (g .) + +{-# INLINE ($.) #-} +infixl 1 $. +($.) = flip ($) + + +{-# INLINE logFromLFAssocs #-} +logFromLFAssocs :: [(Int, LogFloat)] -> [(Int, Double)] +logFromLFAssocs = unsafeCoerce + +{-# INLINE logFromLFUArray #-} +logFromLFUArray :: UArray a LogFloat -> UArray a Double +logFromLFUArray = unsafeCoerce + +{-# INLINE logToLFUArray #-} +logToLFUArray :: UArray a Double -> UArray a LogFloat +logToLFUArray = unsafeCoerce + +{-# INLINE logToLFFunc #-} +logToLFFunc :: (LogFloat -> a -> LogFloat) -> (Double -> a -> Double) +logToLFFunc = ($. unsafeLogToLogFloat ~> id ~> logFromLogFloat) + +-- | Remove the extranious 'isNaN' test of 'logToLogFloat', when +-- we know we can. +{-# INLINE unsafeLogToLogFloat #-} +unsafeLogToLogFloat :: Double -> LogFloat +unsafeLogToLogFloat = LogFloat + + +instance IArray UArray LogFloat where + {-# INLINE bounds #-} + bounds = bounds . logFromLFUArray + +-- Apparently this method was added in base-2.0/GHC-6.6 but Hugs +-- (Sept 2006) doesn't have it. Not sure about NHC's base +#if __HUGS__ > 200609 + {-# INLINE numElements #-} + numElements = numElements . logFromLFUArray +#endif + + {-# INLINE unsafeArray #-} + unsafeArray = + unsafeArray $. id ~> logFromLFAssocs ~> logToLFUArray + + {-# INLINE unsafeAt #-} + unsafeAt = + unsafeAt $. logFromLFUArray ~> id ~> unsafeLogToLogFloat + + {-# INLINE unsafeReplace #-} + unsafeReplace = + unsafeReplace $. logFromLFUArray ~> logFromLFAssocs ~> logToLFUArray + + {-# INLINE unsafeAccum #-} + unsafeAccum = + unsafeAccum $. logToLFFunc ~> logFromLFUArray ~> id ~> logToLFUArray + + {-# INLINE unsafeAccumArray #-} + unsafeAccumArray = + unsafeAccumArray $. logToLFFunc ~> logFromLogFloat ~> id ~> id ~> logToLFUArray +#endif + + +instance PartialOrd LogFloat where + cmp (LogFloat x) (LogFloat y) + | isNaN x || isNaN y = Nothing + | otherwise = Just $! x `compare` y + + +---------------------------------------------------------------- +-- | A constructor which does semantic conversion from normal-domain +-- to log-domain. +logFloat :: (Real a, RealToFrac a Double) => a -> LogFloat +{-# SPECIALIZE logFloat :: Double -> LogFloat #-} +logFloat = LogFloat . log . guardNonNegative "logFloat" . realToFrac + + +-- This is simply a polymorphic version of the 'LogFloat' data +-- constructor. We present it mainly because we hide the constructor +-- in order to make the type a bit more opaque. If the polymorphism +-- turns out to be a performance liability because the rewrite rules +-- can't remove it, then we need to rethink all four +-- constructors\/destructors. +-- +-- | Constructor which assumes the argument is already in the +-- log-domain. +logToLogFloat :: (Real a, RealToFrac a Double) => a -> LogFloat +{-# SPECIALIZE logToLogFloat :: Double -> LogFloat #-} +logToLogFloat = LogFloat . guardIsANumber "logToLogFloat" . realToFrac + + +-- | Return our log-domain value back into normal-domain. Beware +-- of overflow\/underflow. +fromLogFloat :: (Fractional a, Transfinite a, RealToFrac Double a) + => LogFloat -> a +{-# SPECIALIZE fromLogFloat :: LogFloat -> Double #-} +fromLogFloat (LogFloat x) = realToFrac (exp x) + + +-- | Return the log-domain value itself without costly conversion +logFromLogFloat :: (Fractional a, Transfinite a, RealToFrac Double a) + => LogFloat -> a +{-# SPECIALIZE logFromLogFloat :: LogFloat -> Double #-} +logFromLogFloat (LogFloat x) = realToFrac x + + +-- These are our module-specific versions of "log\/exp" and "exp\/log"; +-- They do the same things but also have a @LogFloat@ in between +-- the logarithm and exponentiation. +-- +-- In order to ensure these rules fire we may need to delay inlining +-- of the four con-\/destructors, like we do for 'realToFrac'. +-- Unfortunately, I'm not entirely sure whether they will be inlined +-- already or not (and whether they are may be fragile) and I don't +-- want to inline them excessively and lead to code bloat in the +-- off chance that we could prune some of it away. +-- TODO: thoroughly investigate this. + +{-# RULES +-- Out of log-domain and back in +"log/fromLogFloat" forall x. log (fromLogFloat x) = logFromLogFloat x +"log.fromLogFloat" log . fromLogFloat = logFromLogFloat + +"logFloat/fromLogFloat" forall x. logFloat (fromLogFloat x) = x +"logFloat.fromLogFloat" logFloat . fromLogFloat = id + +-- Into log-domain and back out +"fromLogFloat/logFloat" forall x. fromLogFloat (logFloat x) = x +"fromLogFloat.logFloat" fromLogFloat . logFloat = id + #-} + +---------------------------------------------------------------- +-- To show it, we want to show the normal-domain value rather than +-- the log-domain value. Also, if someone managed to break our +-- invariants (e.g. by passing in a negative and noone's pulled on +-- the thunk yet) then we want to crash before printing the +-- constructor, rather than after. N.B. This means the show will +-- underflow\/overflow in the same places as normal doubles since +-- we underflow at the @exp@. Perhaps this means we should show the +-- log-domain value instead. + +instance Show LogFloat where + show (LogFloat x) = let y = exp x + in y `seq` "LogFloat "++show y + + +---------------------------------------------------------------- +-- These all work without causing underflow. However, do note that +-- they tend to induce more of the floating-point fuzz than using +-- regular floating numbers because @exp . log@ doesn't really equal +-- @id@. In any case, our main aim is for preventing underflow when +-- multiplying many small numbers (and preventing overflow for +-- multiplying many large numbers) so we're not too worried about +-- +\/- 4e-16. + +instance Num LogFloat where + -- BUG? In Hugs (Sept2006) the (>=) always returns True if + -- either isNaN. Only questionably a bug, since we try to + -- ensure that notANumber never occurs. Still... perhaps + -- we should use `ge` and other PartialOrd things in order + -- to play it safe. + -- TODO: benchmark and check core to see how much that hurts GHC. + + + (*) (LogFloat x) (LogFloat y) = LogFloat (x+y) + + (+) (LogFloat x) (LogFloat y) + | x >= y = LogFloat (x + log (1 + exp (y - x))) + | otherwise = LogFloat (y + log (1 + exp (x - y))) + + -- Without the guard this would return NaN instead of error + (-) (LogFloat x) (LogFloat y) + | x >= y = LogFloat (x + log (1 - exp (y - x))) + | otherwise = errorOutOfRange "(-)" + + signum (LogFloat x) + | x == negativeInfinity = 0 + | x > negativeInfinity = 1 + | otherwise = errorOutOfRange "signum" + -- The extra guard protects against NaN, in case someone + -- broke the invariant. That shouldn't be possible and + -- so noone else bothers to check, but we check here just + -- in case. + + negate _ = errorOutOfRange "negate" + + abs = id + + fromInteger = LogFloat . log + . guardNonNegative "fromInteger" . fromInteger + + +instance Fractional LogFloat where + -- n/0 is handled seamlessly for us; we must catch 0/0 though + (/) (LogFloat x) (LogFloat y) + | x == negativeInfinity + && y == negativeInfinity = errorOutOfRange "(/)" -- protect vs NaN + | otherwise = LogFloat (x-y) + + fromRational = LogFloat . log + . guardNonNegative "fromRational" . fromRational + + +-- Just for fun. The more coersion functions the better. Though +-- Rationals are very buggy when it comes to transfinite values +instance Real LogFloat where + toRational (LogFloat x) = toRational (exp x) + + +{- -- Commented out because I'm not sure about requiring MPTCs. Of course, those are already required by "Data.Number.Transfinite" so it's pretty moot... + +-- LogFloat->LogFloat is already given via generic (a->a) +-- No LogFloat->Rational since LogFloat can have 'infinity' +-- Can't have LogFloat->a using fromLogFloat because Hugs dislikes incoherence. Adding an explicit LogFloat->LogFloat instance doesn't help like it does for GHC. + +instance RealToFrac LogFloat Double where + realToFrac = fromLogFloat + +instance RealToFrac LogFloat Float where + realToFrac = fromLogFloat +-} + +---------------------------------------------------------------- +----------------------------------------------------------- fin. diff --git a/src/Data/Number/PartialOrd.hs b/src/Data/Number/PartialOrd.hs new file mode 100644 --- /dev/null +++ b/src/Data/Number/PartialOrd.hs @@ -0,0 +1,148 @@ +{-# LANGUAGE OverlappingInstances + , FlexibleInstances + , UndecidableInstances + #-} + +{-# OPTIONS_GHC -Wall -fwarn-tabs #-} + +---------------------------------------------------------------- +-- ~ 2009.01.29 +-- | +-- Module : Data.Number.PartialOrd +-- Copyright : Copyright (c) 2007--2009 wren ng thornton +-- License : BSD3 +-- Maintainer : wren@community.haskell.org +-- Stability : experimental +-- Portability : semi-portable (overlapping instances, etc) +-- +-- The Prelude's 'Ord' class for dealing with ordered types is often +-- onerous to use because it requires 'Eq' as well as a total +-- ordering. While such total orderings are common, partial orderings +-- are moreso. This module presents a class for partially ordered +-- types. +---------------------------------------------------------------- +module Data.Number.PartialOrd + ( + -- * Partial Ordering + PartialOrd(..) + -- * Functions + , comparingPO + ) where + +-- Bugfix for Hugs (September 2006), see note below. +import Prelude hiding (isNaN) +import Hugs.RealFloat (isNaN) + +---------------------------------------------------------------- +-- | This class defines a partially ordered type. The method names +-- were chosen so as not to conflict with 'Ord' and 'Eq'. We use +-- 'Maybe' instead of defining new types @PartialOrdering@ and +-- @FuzzyBool@ because this way should make the class easier to +-- use. +-- +-- Minimum complete definition: 'cmp' + +class PartialOrd a where + -- | like 'compare' + cmp :: a -> a -> Maybe Ordering + + -- | like ('>') + gt :: a -> a -> Maybe Bool + gt x y = case x `cmp` y of + Just GT -> Just True + Just _ -> Just False + Nothing -> Nothing + + -- | like ('>=') + ge :: a -> a -> Maybe Bool + ge x y = case x `cmp` y of + Just LT -> Just False + Just _ -> Just True + Nothing -> Nothing + + -- | like ('==') + eq :: a -> a -> Maybe Bool + eq x y = case x `cmp` y of + Just EQ -> Just True + Just _ -> Just False + Nothing -> Nothing + + -- | like ('/=') + ne :: a -> a -> Maybe Bool + ne x y = case x `cmp` y of + Just EQ -> Just False + Just _ -> Just True + Nothing -> Nothing + + -- | like ('<=') + le :: a -> a -> Maybe Bool + le x y = case x `cmp` y of + Just GT -> Just False + Just _ -> Just True + Nothing -> Nothing + + -- | like ('<') + lt :: a -> a -> Maybe Bool + lt x y = case x `cmp` y of + Just LT -> Just True + Just _ -> Just False + Nothing -> Nothing + + -- | like 'max'. The default instance returns the left argument + -- when they're equal. + maxPO :: a -> a -> Maybe a + maxPO x y = do o <- x `cmp` y + case o of + GT -> Just x + EQ -> Just x + LT -> Just y + + -- | like 'min'. The default instance returns the left argument + -- when they're equal. + minPO :: a -> a -> Maybe a + minPO x y = do o <- x `cmp` y + case o of + GT -> Just y + EQ -> Just x + LT -> Just x + +infix 4 `gt`, `ge`, `eq`, `ne`, `le`, `lt`, `maxPO`, `minPO` + +instance (Ord a) => PartialOrd a where + cmp x y = Just $! x `compare` y + gt x y = Just $! x > y + ge x y = Just $! x >= y + eq x y = Just $! x == y + ne x y = Just $! x /= y + le x y = Just $! x <= y + lt x y = Just $! x < y + maxPO x y = Just $! x `max` y + minPO x y = Just $! x `min` y + + +-- N.B. Hugs (Sept 2006) has a buggy definition for 'isNaN' which +-- always returns @False@. We use a fixed version, provided the CPP +-- was run with the right arguments. See "Hugs.RealFloat". If 'cmp' +-- returns @Just Eq@ for @notANumber@ then CPP was run wrongly. +-- +-- The instances inherited from Ord are wrong. So we'll fix them. +instance PartialOrd Float where + cmp x y | isNaN x || isNaN y = Nothing + | otherwise = Just $! x `compare` y + +instance PartialOrd Double where + cmp x y | isNaN x || isNaN y = Nothing + | otherwise = Just $! x `compare` y + +---------------------------------------------------------------- +-- TODO? add maximumPO\/minimumPO via left or right fold? + +-- BUG: Haddock doesn't link the `comparing` +-- +-- | Like @Data.Ord.comparing@. Helpful in conjunction with the +-- @xxxBy@ family of functions from "Data.List" +comparingPO :: (PartialOrd b) => (a -> b) -> a -> a -> Maybe Ordering +comparingPO p x y = p x `cmp` p y + +---------------------------------------------------------------- +----------------------------------------------------------- fin. diff --git a/src/Data/Number/RealToFrac.hs b/src/Data/Number/RealToFrac.hs new file mode 100644 --- /dev/null +++ b/src/Data/Number/RealToFrac.hs @@ -0,0 +1,117 @@ +-- Needed to ensure correctness, and because we can't guarantee rules fire +-- The MagicHash is for unboxed primitives (-fglasgow-exts also works) +-- We only need MagicHash if on GHC, but we can't hide it in an #ifdef +{-# LANGUAGE MultiParamTypeClasses + , OverlappingInstances + , FlexibleInstances + , CPP + , MagicHash + #-} + +{-# OPTIONS_GHC -Wall -fwarn-tabs #-} + +---------------------------------------------------------------- +-- ~ 2009.01.29 +-- | +-- Module : Data.Number.RealToFrac +-- Copyright : Copyright (c) 2007--2009 wren ng thornton +-- License : BSD3 +-- Maintainer : wren@community.haskell.org +-- Stability : experimental +-- Portability : non-portable (CPP, MPTC, OverlappingInstances) +-- +-- This module presents a type class for generic conversion between +-- numeric types, generalizing @realToFrac@ in order to overcome +-- problems with pivoting through 'Rational' +---------------------------------------------------------------- +module Data.Number.RealToFrac (RealToFrac(..)) where + +import Prelude hiding (realToFrac, isInfinite, isNaN) +import qualified Prelude (realToFrac) + +import Data.Number.Transfinite + +#ifdef __GLASGOW_HASKELL__ +import GHC.Exts + ( Int(..), Float(..), Double(..) + , int2Double# + , int2Float# + , double2Float# + , float2Double# + ) +#endif + +---------------------------------------------------------------- +-- | The 'Prelude.realToFrac' function is defined to pivot through +-- a 'Rational' according to the haskell98 spec. This is non-portable +-- and problematic as discussed in "Data.Number.Transfinite". Since +-- there is resistance to breaking from the spec, this class defines +-- a reasonable variant which deals with transfinite values +-- appropriately. +-- +-- There is a generic instance from any Transfinite Real to any +-- Transfinite Fractional, using checks to ensure correctness. GHC +-- has specialized versions for some types which use primitive +-- converters instead, for large performance gains. (These definitions +-- are hidden from other compilers via CPP.) Due to a bug in Haddock +-- the specialized instances are shown twice and the generic instance +-- isn't shown at all. Since the instances are overlapped, you'll +-- need to give type signatures if the arguments to 'realToFrac' +-- are polymorphic. There's also a generic instance for any Real +-- Fractional type to itself, thus if you write any generic instances +-- beware of incoherence. +-- +-- If any of these restrictions (CPP, GHC-only optimizations, +-- OverlappingInstances) are onerous to you, contact the maintainer +-- (we like patches). Note that this /does/ work for Hugs with +-- suitable options (e.g. @hugs -98 +o -F'cpp -P'@). However, Hugs +-- doesn't allow @IncoherentInstances@ nor does it allow diamonds +-- with @OverlappingInstances@, which restricts the ability to add +-- additional generic instances. + +class (Real a, Fractional b) => RealToFrac a b where + realToFrac :: a -> b + +instance (Real a, Fractional a) => RealToFrac a a where + realToFrac = id + +instance (Real a, Transfinite a, Fractional b, Transfinite b) + => RealToFrac a b + where + realToFrac x + | isNaN x = notANumber + | isInfinite x = if x > 0 then infinity + else negativeInfinity + | otherwise = Prelude.realToFrac x + + +#ifdef __GLASGOW_HASKELL__ +instance RealToFrac Int Float where + {-# INLINE realToFrac #-} + realToFrac (I# i) = F# (int2Float# i) + +instance RealToFrac Int Double where + {-# INLINE realToFrac #-} + realToFrac (I# i) = D# (int2Double# i) + + +instance RealToFrac Integer Float where + -- TODO: is there a more primitive way? + realToFrac j = Prelude.realToFrac j + +instance RealToFrac Integer Double where + -- TODO: is there a more primitive way? + realToFrac j = Prelude.realToFrac j + + +instance RealToFrac Float Double where + {-# INLINE realToFrac #-} + realToFrac (F# f) = D# (float2Double# f) + +instance RealToFrac Double Float where + {-# INLINE realToFrac #-} + realToFrac (D# d) = F# (double2Float# d) +#endif + +---------------------------------------------------------------- +----------------------------------------------------------- fin. diff --git a/src/Data/Number/Transfinite.hs b/src/Data/Number/Transfinite.hs new file mode 100644 --- /dev/null +++ b/src/Data/Number/Transfinite.hs @@ -0,0 +1,184 @@ +{-# OPTIONS_GHC -Wall -fwarn-tabs #-} + +---------------------------------------------------------------- +-- ~ 2009.01.29 +-- | +-- Module : Data.Number.Transfinite +-- Copyright : Copyright (c) 2007--2009 wren ng thornton +-- License : BSD3 +-- Maintainer : wren@community.haskell.org +-- Stability : experimental +-- Portability : portable +-- +-- This module presents a type class for numbers which have +-- representations for transfinite values. The idea originated from +-- the IEEE-754 floating-point special values, used by +-- "Data.Number.LogFloat". However not all 'Fractional' types +-- necessarily support transfinite values. In particular, @Ratio@ +-- types including 'Rational' do not have portable representations. +-- +-- For the Glasgow compiler (GHC 6.8.2), "GHC.Real" defines @1%0@ +-- and @0%0@ as representations for 'infinity' and 'notANumber', +-- but most operations on them will raise exceptions. If 'toRational' +-- is used on an infinite floating value, the result is a rational +-- with a numerator sufficiently large that it will overflow when +-- converted back to a @Double@. If used on NaN, the result would +-- buggily convert back as 'negativeInfinity'. For more discussion +-- on why this approach is problematic, see: +-- +-- * +-- +-- * +-- +-- Hugs (September 2006) stays closer to the haskell98 spec and +-- offers no way of constructing those values, raising arithmetic +-- overflow errors if attempted. +---------------------------------------------------------------- +module Data.Number.Transfinite + ( Transfinite(..) + , log + ) where + +import Prelude hiding (log, isInfinite, isNaN) +import qualified Prelude (log) +import qualified Hugs.RealFloat as Prelude (isInfinite, isNaN) + +import Data.Number.PartialOrd + +---------------------------------------------------------------- +-- | Many numbers are not 'Bounded' yet, even though they can +-- represent arbitrarily large values, they are not necessarily +-- able to represent transfinite values such as infinity itself. +-- This class is for types which are capable of representing such +-- values. Notably, this class does not require the type to be +-- 'Fractional' nor 'Floating' since integral types could also have +-- representations for transfinite values. By popular demand the +-- 'Num' restriction has been lifted as well, due to complications +-- of defining 'Show' or 'Eq' for some types. +-- +-- In particular, this class extends the ordered projection to have +-- a maximum value 'infinity' and a minimum value 'negativeInfinity', +-- as well as an exceptional value 'notANumber'. All the natural +-- laws regarding @infinity@ and @negativeInfinity@ should pertain. +-- (Some of these are discussed below.) +-- +-- Hugs (September 2006) has buggy Prelude definitions for +-- 'Prelude.isNaN' and 'Prelude.isInfinite' on Float and Double. +-- This module provides correct definitions, so long as "Hugs.RealFloat" +-- is compiled correctly. + +class (PartialOrd a) => Transfinite a where + + -- | A transfinite value which is greater than all finite values. + -- Adding or subtracting any finite value is a no-op. As is + -- multiplying by any non-zero positive value (including + -- @infinity@), and dividing by any positive finite value. Also + -- obeys the law @negate infinity = negativeInfinity@ with all + -- appropriate ramifications. + + infinity :: a + + + -- | A transfinite value which is less than all finite values. + -- Obeys all the same laws as @infinity@ with the appropriate + -- changes for the sign difference. + + negativeInfinity :: a + + + -- | An exceptional transfinite value for dealing with undefined + -- results when manipulating infinite values. The following + -- operations must return @notANumber@, where @inf@ is any value + -- which @isInfinite@: + -- + -- * @inf + inf@ + -- + -- * @inf - inf@ + -- + -- * @inf * 0@ + -- + -- * @0 * inf@ + -- + -- * @inf \/ inf@ + -- + -- * @inf `div` inf@ + -- + -- * @0 \/ 0@ + -- + -- * @0 `div` 0@ + -- + -- Additionally, any mathematical operations on @notANumber@ + -- must also return @notANumber@, and any equality or ordering + -- comparison on @notANumber@ must return @False@ (violating + -- the law of the excluded middle, often assumed but not required + -- for 'Eq'; thus, 'eq' and 'ne' are preferred over ('==') and + -- ('/=')). Since it returns false for equality, there may be + -- more than one machine representation of this `value'. + + notANumber :: a + + + -- | Return true for both @infinity@ and @negativeInfinity@, + -- false for all other values. + isInfinite :: a -> Bool + + -- | Return true only for @notANumber@. + isNaN :: a -> Bool + + +instance Transfinite Double where + infinity = 1/0 + negativeInfinity = negate (1/0) + notANumber = 0/0 + isInfinite = Prelude.isInfinite + isNaN = Prelude.isNaN + + +instance Transfinite Float where + infinity = 1/0 + negativeInfinity = negate (1/0) + notANumber = 0/0 + isInfinite = Prelude.isInfinite + isNaN = Prelude.isNaN + + +---------------------------------------------------------------- +-- | Since the normal 'Prelude.log' throws an error on zero, we +-- have to redefine it in order for things to work right. Arguing +-- from limits we can see that @log 0 == negativeInfinity@. Newer +-- versions of GHC have this behavior already, but older versions +-- and Hugs do not. +-- +-- This function will raise an error when taking the log of negative +-- numbers, rather than returning 'notANumber' as the newer GHC +-- implementation does. The reason being that typically this is a +-- logical error, and @notANumber@ allows the error to propegate +-- silently. +-- +-- In order to improve portability, the 'Transfinite' class is +-- required to indicate that the 'Floating' type does in fact have +-- a representation for negative infinity. Both native floating +-- types ('Double' and 'Float') are supported. If you define your +-- own instance of @Transfinite@, verify the above equation holds +-- for your @0@ and @negativeInfinity@. If it doesn't, then you +-- should avoid importing our @log@ and will probably want converters +-- to handle the discrepancy. + +log :: (Floating a, Transfinite a) => a -> a +{-# SPECIALIZE log :: Double -> Double #-} +{-# SPECIALIZE log :: Float -> Float #-} +log x = case x `cmp` 0 of + Just GT -> Prelude.log x + Just EQ -> negativeInfinity + Just LT -> err "argument out of range" + Nothing -> err "argument not comparable to 0" + where + err e = error $! "Data.Number.Transfinite.log: "++e + +-- Note, Floating ultimately requires Num, but not Ord. If PartialOrd +-- proves to be an onerous requirement on Transfinite, we could +-- hack our way around without using PartialOrd by using isNaN, (== +-- 0), ((>0).signum) but that would be less efficient. + +---------------------------------------------------------------- +----------------------------------------------------------- fin. diff --git a/src/Hugs/RealFloat.hs b/src/Hugs/RealFloat.hs new file mode 100644 --- /dev/null +++ b/src/Hugs/RealFloat.hs @@ -0,0 +1,63 @@ + +{-# LANGUAGE CPP #-} + +{-# OPTIONS_GHC -Wall -fwarn-tabs #-} + +---------------------------------------------------------------- +-- ~ 2009.01.29 +-- | +-- Module : Hugs.RealFloat +-- Copyright : Copyright (c) 2007--2009 wren ng thornton +-- License : BSD3 +-- Maintainer : wren@community.haskell.org +-- Stability : stable +-- Portability : portable (with CPP) +-- +-- Hugs (September 2006) has buggy definitions for 'Prelude.isNaN' +-- and 'Prelude.isInfinite' on Float and Double. If this module is +-- run through CPP with the macro @__HUGS__@ set to a value no +-- larger than 200609, then correct definitions are used. Otherwise +-- the Prelude definitions are used (which should be correct for +-- other compilers). For example, run Hugs with +-- +-- @hugs -F'cpp -P -D__HUGS__=200609' Hugs/RealFloat.hs@ +-- +-- N.B. The corrected definitions have only been tested to work for +-- 'Float' and 'Double'. These definitions should probably not be +-- used for other 'RealFloat' types. +---------------------------------------------------------------- +module Hugs.RealFloat + ( isInfinite + , isNaN + ) where + +import Prelude hiding (isInfinite, isNaN) +import qualified Prelude +---------------------------------------------------------------- + +isInfinite :: (RealFloat a) => a -> Bool +{-# INLINE isInfinite #-} +#if defined(__HUGS__) && (__HUGS__ <= 200609) +isInfinite x = (1/0) == abs x +#else +isInfinite = Prelude.isInfinite +#endif + + +isNaN :: (RealFloat a) => a -> Bool +{-# INLINE isNaN #-} +#if defined(__HUGS__) && (__HUGS__ <= 200609) +isNaN x = compareEQ x 0 && compareEQ x 1 + +-- | In Hugs (September 2006), 'compare' always returns @EQ@ if one +-- of the arguments is not a number. Thus, if a number is @compareEQ@ +-- against multiple different numbers, then it must be @isNaN@. +compareEQ :: (Ord a) => a -> a -> Bool +compareEQ x y = case compare x y of + EQ -> True + _ -> False +#else +isNaN = Prelude.isNaN +#endif +---------------------------------------------------------------- +----------------------------------------------------------- fin.