diff --git a/brush-strokes/brush-strokes.cabal b/brush-strokes/brush-strokes.cabal index 08bf574..6362515 100644 --- a/brush-strokes/brush-strokes.cabal +++ b/brush-strokes/brush-strokes.cabal @@ -115,6 +115,7 @@ library , Math.Bezier.Stroke.EnvelopeEquation , Math.Differentiable , Math.Epsilon + , Math.Float.Utils , Math.Interval , Math.Linear , Math.Linear.Solve diff --git a/brush-strokes/src/lib/Math/Float/Utils.hs b/brush-strokes/src/lib/Math/Float/Utils.hs new file mode 100644 index 0000000..e99c9d6 --- /dev/null +++ b/brush-strokes/src/lib/Math/Float/Utils.hs @@ -0,0 +1,76 @@ +{-# LANGUAGE AllowAmbiguousTypes #-} +{-# LANGUAGE ScopedTypeVariables #-} + +module Math.Float.Utils + ( FPBits(..) + , nextAfter, succFP, prevFP + ) + where + +-- base +import Data.Bits + ( Bits((.&.), shiftL) ) +import Data.Word + ( Word32, Word64 ) +import GHC.Float + ( castFloatToWord32 , castWord32ToFloat + , castDoubleToWord64, castWord64ToDouble + ) + +-------------------------------------------------------------------------------- + +class ( RealFloat f, Num b, Bits b ) => FPBits f b | f -> b, b -> f where + toBits :: f -> b + fromBits :: b -> f + -- | Size in bytes. + sizeOf :: Int + +instance FPBits Float Word32 where + toBits = castFloatToWord32 + fromBits = castWord32ToFloat + sizeOf = 4 + +instance FPBits Double Word64 where + toBits = castDoubleToWord64 + fromBits = castWord64ToDouble + sizeOf = 8 + + +{-# SPECIALISE nextAfter :: Float -> Float -> Float #-} +{-# SPECIALISE nextAfter :: Double -> Double -> Double #-} +{-# INLINEABLE nextAfter #-} +-- | @nextAfter a b@ computes the next floating-point value after @a@ +-- in the direction of @b@. +nextAfter :: forall f b. FPBits f b => f -> f -> f +nextAfter a b + | isNaN a + = a + | isNaN b + = b + | a == b + = b + | otherwise + = let !res_bits + | a == 0 + , let !sgn_mask = 1 `shiftL` ( sizeOf @f * 8 - 1 ) + = ( toBits b .&. sgn_mask ) + 1 + | ( a < b ) == ( a > 0 ) + = toBits a + 1 + | otherwise + = toBits a - 1 + in fromBits res_bits + +{-# SPECIALISE succFP :: Float -> Float #-} +{-# SPECIALISE succFP :: Double -> Double #-} +{-# INLINEABLE succFP #-} +-- | The next floating-point number. +succFP :: forall f b. FPBits f b => f -> f +succFP x = nextAfter x (1/0) + + +{-# SPECIALISE prevFP :: Float -> Float #-} +{-# SPECIALISE prevFP :: Double -> Double #-} +{-# INLINEABLE prevFP #-} +-- | The previous floating-point number. +prevFP :: forall f b. FPBits f b => f -> f +prevFP x = nextAfter x (-1/0) diff --git a/brush-strokes/src/lib/Math/Interval/FMA.hs b/brush-strokes/src/lib/Math/Interval/FMA.hs index eea05f6..d5ca39d 100644 --- a/brush-strokes/src/lib/Math/Interval/FMA.hs +++ b/brush-strokes/src/lib/Math/Interval/FMA.hs @@ -1,77 +1,17 @@ {-# LANGUAGE AllowAmbiguousTypes #-} {-# LANGUAGE ScopedTypeVariables #-} -module Math.Interval.FMA ( addI, subI, prodI, divI, posPowI ) where +module Math.Interval.FMA + ( addI, subI, prodI, divI, posPowI ) + where -- base -import Data.Bits - ( Bits((.&.), shiftL) ) -import Data.Word - ( Word32, Word64 ) import GHC.Exts ( Double(D#), fmsubDouble#, fnmaddDouble# ) -import GHC.Float - ( castFloatToWord32 , castWord32ToFloat - , castDoubleToWord64, castWord64ToDouble - ) --------------------------------------------------------------------------------- - -class ( RealFloat f, Num b, Bits b ) => FPBits f b | f -> b, b -> f where - toBits :: f -> b - fromBits :: b -> f - -- | Size in bytes. - sizeOf :: Int - -instance FPBits Float Word32 where - toBits = castFloatToWord32 - fromBits = castWord32ToFloat - sizeOf = 4 - -instance FPBits Double Word64 where - toBits = castDoubleToWord64 - fromBits = castWord64ToDouble - sizeOf = 8 - - -{-# SPECIALISE nextAfter :: Float -> Float -> Float #-} -{-# SPECIALISE nextAfter :: Double -> Double -> Double #-} -{-# INLINEABLE nextAfter #-} --- | @nextAfter a b@ computes the next floating-point value after @a@ --- in the direction of @b@. -nextAfter :: forall f b. FPBits f b => f -> f -> f -nextAfter a b - | isNaN a - = a - | isNaN b - = b - | a == b - = b - | otherwise - = let !res_bits - | a == 0 - , let !sgn_mask = 1 `shiftL` ( sizeOf @f * 8 - 1 ) - = ( toBits b .&. sgn_mask ) + 1 - | ( a < b ) == ( a > 0 ) - = toBits a + 1 - | otherwise - = toBits a - 1 - in fromBits res_bits - -{-# SPECIALISE succFP :: Float -> Float #-} -{-# SPECIALISE succFP :: Double -> Double #-} -{-# INLINEABLE succFP #-} --- | The next floating-point number. -succFP :: forall f b. FPBits f b => f -> f -succFP x = nextAfter x (1/0) - - -{-# SPECIALISE prevFP :: Float -> Float #-} -{-# SPECIALISE prevFP :: Double -> Double #-} -{-# INLINEABLE prevFP #-} --- | The previous floating-point number. -prevFP :: forall f b. FPBits f b => f -> f -prevFP x = nextAfter x (-1/0) +-- brush-strokes +import Math.Float.Utils + ( prevFP, succFP ) --------------------------------------------------------------------------------