posit 3.2.0.0 → 3.2.0.1
raw patch · 6 files changed
+257/−747 lines, 6 files
Files
- ChangeLog.md +11/−0
- README.md +7/−6
- posit.cabal +3/−3
- src/Posit.hs +71/−37
- src/Posit/Internal/PositC.hs +152/−701
- test/TestPosit.hs +13/−0
ChangeLog.md view
@@ -1,5 +1,16 @@ # Changelog for Posit Numbers +# posit-3.2.0.1++ * Refactored `IntN` Type Family to be a closed type family instead of an associated type family+ * Refactored `IntN` constraints to use `ConstraintKinds` and made that to be a Super Class of `PositC` to improve the encapsulation the Constraints of the internal implementation+ * Refactored `PositC` to make use of `ConstrainedClassMethods` vastly reducing code duplication+ * Eliminated the `FlexableContexts` Language Extension from Posit.hs Interface, since the `InN` constraints no longer bleed into that file+ * Added test of Heegner numbers (almost integers)+ * Added test of various Gamma Function approximations+ * Improved function names in the Orphan Instance for `Storable` ( `Word128` )+ * Improved documentation+ ## posit-3.2.0.0 * Posit Standard 3.2 [Posit Standard] (https://posithub.org/docs/posit_standard.pdf)
README.md view
@@ -1,10 +1,11 @@-# posit 3.2.0.0+# posit 3.2.0.1 The [Posit Standard 3.2](https://posithub.org/docs/posit_standard.pdf),-where Real numbers are approximated by Maybe Rational. The Posit type-is mapped to a 2's complement integer type; smoothly and with tapering-precision, in a similar way to the projective real line. The 'posit'-library implements the following standard classes:+where Real numbers are approximated by Maybe Rational. The Posit +Numbers are a drop in replacement for `Float` or `Double` mapped to a +2's complement integer type; smoothly and with tapering precision, in a +similar way to the projective real line. The 'posit' library implements+the following standard classes: * Show * Eq@@ -24,7 +25,7 @@ * Floating -- Mathematical functions such as logarithm, exponential, trigonometric, and hyperbolic functions. Warning! May induce trance. The Posits are indexed by the type (es :: ES) where exponent size and-word size are related. In `posit-3.2.0.0` es is instantiated as Z, I,+word size are related. In `posit-3.2.0.1` es is instantiated as Z, I, II, III, IV, V. The word size (in bits) of the value is `= 8 * 2^es`, that is `2^es` bytes. The Types: 'Posit8', 'Posit16', 'Posit32', 'Posit64', 'Posit128', and 'Posit256' are implemented and include a
posit.cabal view
@@ -1,7 +1,7 @@ cabal-version: 1.12 name: posit-version: 3.2.0.0+version: 3.2.0.1 description: The Posit Number format. Please see the README on GitHub at <https://github.com/waivio/posit#readme> homepage: https://github.com/waivio/posit#readme bug-reports: https://github.com/waivio/posit/issues@@ -55,7 +55,7 @@ ghc-options: -Wall -O2 if flag(do-liquid)- ghc-options: -fplugin=LiquidHaskell -fplugin-opt=LiquidHaskell:--fast -fplugin-opt=LiquidHaskell:--max-case-expand=4 -fplugin-opt=LiquidHaskell:--no-termination -fplugin-opt=LiquidHaskell:--short-names+ ghc-options: -fplugin=LiquidHaskell -fplugin-opt=LiquidHaskell:--fast -fplugin-opt=LiquidHaskell:--no-termination -fplugin-opt=LiquidHaskell:--max-case-expand=4 -fplugin-opt=LiquidHaskell:--short-names if flag(do-no-storable) cpp-options: -DO_NO_STORABLE@@ -64,7 +64,7 @@ cpp-options: -DO_NO_ORPHANS if flag(do-liquid)- cpp-options: -DO_LIQUID -DO_NO_STORABLE -DO_NO_READ -DO_NO_SHOW+ cpp-options: -DO_LIQUID -DO_NO_STORABLE if flag(do-test) cpp-options: -DO_TEST
src/Posit.hs view
@@ -21,7 +21,6 @@ {-# LANGUAGE BangPatterns #-} -- Added Strictness for some fixed point algorithms {-# LANGUAGE PatternSynonyms #-} -- for a nice NaR interface {-# LANGUAGE FlexibleInstances #-} -- To make instances for each specific type [Posit8 .. Posit256]-{-# LANGUAGE FlexibleContexts #-} -- Allow non-type variables in the constraints {-# LANGUAGE TypeApplications #-} -- To apply types: @Type, it seems to select the specific class instance, when GHC is not able to reason about things, commenting this out shows an interesting interface {-# LANGUAGE MultiParamTypeClasses #-} -- To convert between Posit Types {-# LANGUAGE ScopedTypeVariables #-} -- To reduce some code duplication@@ -31,6 +30,7 @@ {-# OPTIONS_GHC -Wno-type-defaults #-} -- Turn off noise {-# OPTIONS_GHC -Wno-unused-top-binds #-} -- Turn off noise + -- ---- -- Posit numbers implementing: --@@ -54,7 +54,8 @@ -- ---- module Posit-(-- * Main Exported Types+(Posit(),+ -- * Main Exported Types Posit8, -- |An 8-bit Posit number with 'es' ~ 'Z' Posit16, -- |An 16-bit Posit number with 'es' ~ 'I' Posit32, -- |An 32-bit Posit number with 'es' ~ 'II'@@ -72,8 +73,10 @@ -- * Posits are Convertable between different Posit representations Convertible(..), +#ifndef O_NO_SHOW -- * Additional functions to show the Posit in different formats AltShow(..),+#endif -- * Additional Special Functions AltFloating(..),@@ -86,6 +89,7 @@ viaRational4, viaRational6, viaRational8,+ #ifdef O_TEST -- * Alternative algorithms for test purposes funExp,@@ -93,6 +97,13 @@ funExpTaylor, funLogTaylor, funExpTuma,+ funGammaSeriesFused,+ funGammaRamanujan,+ funGammaCalc,+ funGammaNemes,+ funGammaYang,+ funGammaChen,+ funGammaXminus1, funLogTuma, funLogDomainReduction, funPi1,@@ -103,14 +114,14 @@ funPsiSha2, funPsiSha3 #endif+ ) where import Prelude hiding (rem) -- Imports for Show and Read Instances-import Data.Scientific (Scientific- ,scientificP+import Data.Scientific (scientificP ,fromRationalRepetendUnlimited ,formatScientific ,FPFormat(Generic)) -- Used to print/show and read the rational value@@ -130,7 +141,7 @@ -- Imports for Storable Instance import Foreign.Storable (Storable, sizeOf, alignment, peek, poke) -- Used for Storable Instances of Posit-import Foreign.Ptr (Ptr, plusPtr, castPtr) -- Used for dealing with Pointers for the Posit Storable Instance+import Foreign.Ptr (Ptr, castPtr) -- Used for dealing with Pointers for the Posit Storable Instance -- would like to:@@ -138,7 +149,7 @@ -- Perhaps on the chopping block if we are moving to ElementaryFunctions -- Imports for implementing the Transcendental Functions import GHC.Natural (Natural) -- Import the Natural Numbers ℕ (u+2115) for some of the Transcendental Functions-import Data.Ratio (Rational, (%)) -- Import the Rational Numbers ℚ (u+211A), ℚ can get arbitrarily close to Real numbers ℝ (u+211D), used for some of the Transcendental Functions+import Data.Ratio ((%)) -- Import the Rational Numbers ℚ (u+211A), ℚ can get arbitrarily close to Real numbers ℝ (u+211D), used for some of the Transcendental Functions import Debug.Trace (trace) -- temporary for debug purposes @@ -148,7 +159,7 @@ -- ===================================================================== -- The machine implementation of the Posit encoding/decoding-import Posit.Internal.PositC (ES(..), PositC(..)) -- The main internal implementation details+import Posit.Internal.PositC -- The main internal implementation details -- |Base GADT rapper type, that uses the Exponent Size kind to index the various implementations@@ -156,7 +167,7 @@ Posit :: PositC es => !(IntN es) -> Posit es -- |Not a Real Number, the Posit is like a Maybe type, it's either a real number or not-pattern NaR :: (PositC es) => Posit es+pattern NaR :: PositC es => Posit es pattern NaR <- (Posit (decode -> Nothing)) where NaR = Posit unReal --@@ -182,7 +193,7 @@ #ifndef O_NO_SHOW -- Show ---instance forall es. (PositC es) => Show (Posit es) where+instance PositC es => Show (Posit es) where show NaR = "NaR" show (R r) = formatScientific Generic (Just $ decimalPrec @es) (fst.fromRationalRepetendUnlimited $ r) --@@ -193,7 +204,7 @@ -- Two Posit Numbers are Equal if their Finite Precision Integer representation is Equal -- -- All things equal I would rather write it like this:-instance forall es. (Eq (IntN es)) => Eq (Posit es) where+instance PositC es => Eq (Posit es) where (Posit int1) == (Posit int2) = int1 == int2 -- @@ -202,7 +213,7 @@ -- Two Posit Numbers are ordered by their Finite Precision Integer representation -- -- Ordinarily I would only like one instance to cover them all-instance forall es. (Ord (IntN es), PositC es) => Ord (Posit es) where+instance PositC es => Ord (Posit es) where compare (Posit int1) (Posit int2) = compare int1 int2 -- @@ -211,7 +222,7 @@ -- Num -- -- I'm num trying to get this definition:-instance forall es. (Num (IntN es), Ord (IntN es), PositC es) => Num (Posit es) where+instance PositC es => Num (Posit es) where -- Addition (+) = viaRational2 (+) -- Multiplication@@ -227,7 +238,7 @@ -- -- deriving via Integral Class, for the Integral representation of the posit-viaIntegral :: forall es. PositC es => (IntN es -> IntN es) -> Posit es -> Posit es+viaIntegral :: PositC es => (IntN es -> IntN es) -> Posit es -> Posit es viaIntegral f (Posit int) = Posit $ f int -- @@ -235,7 +246,7 @@ -- Enum-ish, A Posit has a Successor and Predecessor so its an ordinal number, as per Posit standard next, prior -- The Posit Standard requires 2's complement integer overflow to be ignored-instance forall es. (Num (IntN es), Enum (IntN es), Ord (IntN es), PositC es) => Enum (Posit es) where+instance PositC es => Enum (Posit es) where -- succ (Posit int) = Posit (int + 1) succ = viaIntegral (+1) -- succ = viaIntegral succ -- Non-compliant, runtime error pred NaR, and worse it is Int64 for types of greater precision, probably because of Preludes gross abomination of toEnum/fromEnum@@ -281,7 +292,7 @@ -- Fractional Instances; (Num => Fractional) -- -- How the Frac do I get this definition:-instance forall es. (Num (IntN es), Ord (IntN es), Eq (IntN es), PositC es) => Fractional (Posit es) where+instance PositC es => Fractional (Posit es) where fromRational = R recip 0 = NaR@@ -291,7 +302,7 @@ -- Rational Instances; Num & Ord Instanced => Real -- -- I for real want this definition:-instance forall es. (Num (IntN es), Ord (IntN es), PositC es) => Real (Posit es) where+instance PositC es => Real (Posit es) where toRational NaR = error "Your input is Not a Real or Rational (NaR) number, please try again!" toRational (R r) = r --@@ -299,25 +310,25 @@ -- Implementing instances via Rational Data Type's instance, -- The function checks for NaR, to protect against the runtime error 'toRational' would generate if called with a NaR value -- Unary::Arity NaR guarded pass through with wrapping and unwrapping use of a Rational function-viaRational :: (PositC es, Ord (IntN es), Num (IntN es)) => (Rational -> Rational) -> Posit es -> Posit es+viaRational :: PositC es => (Rational -> Rational) -> Posit es -> Posit es viaRational _ NaR = NaR viaRational f (R r) = fromRational $ f r -- Binary NaR guarded pass through with wrapping and unwrapping use of a Rational function-viaRational2 :: (PositC es, Ord (IntN es), Num (IntN es)) => (Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es+viaRational2 :: PositC es => (Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es viaRational2 _ NaR _ = NaR viaRational2 _ _ NaR = NaR viaRational2 f (R r1) (R r2) = R $ r1 `f` r2 -- Ternary NaR guarded pass through with wrapping and unwrapping use of a Rational function-viaRational3 :: (PositC es, Ord (IntN es), Num (IntN es)) => (Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es+viaRational3 :: PositC es => (Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es viaRational3 _ NaR _ _ = NaR viaRational3 _ _ NaR _ = NaR viaRational3 _ _ _ NaR = NaR viaRational3 f (R r1) (R r2) (R r3) = R $ f r1 r2 r3 -- Quaternary NaR guarded pass through with wrapping and unwrapping use of a Rational function-viaRational4 :: (PositC es, Ord (IntN es), Num (IntN es)) => (Rational -> Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es+viaRational4 :: PositC es => (Rational -> Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es viaRational4 _ NaR _ _ _ = NaR viaRational4 _ _ NaR _ _ = NaR viaRational4 _ _ _ NaR _ = NaR@@ -325,7 +336,7 @@ viaRational4 f (R r0) (R r1) (R r2) (R r3) = R $ f r0 r1 r2 r3 -- Senary NaR guarded pass through with wrapping and unwrapping use of a Rational function-viaRational6 :: (PositC es, Ord (IntN es), Num (IntN es)) => (Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es+viaRational6 :: PositC es => (Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es viaRational6 _ NaR _ _ _ _ _ = NaR viaRational6 _ _ NaR _ _ _ _ = NaR viaRational6 _ _ _ NaR _ _ _ = NaR@@ -335,7 +346,7 @@ viaRational6 f (R a1) (R a2) (R a3) (R b1) (R b2) (R b3) = R $ f a1 a2 a3 b1 b2 b3 -- Octonary NaR guarded pass through with wrapping and unwrapping use of a Rational function-viaRational8 :: (PositC es, Ord (IntN es), Num (IntN es)) => (Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es+viaRational8 :: PositC es => (Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational -> Rational) -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es -> Posit es viaRational8 _ NaR _ _ _ _ _ _ _ = NaR viaRational8 _ _ NaR _ _ _ _ _ _ = NaR viaRational8 _ _ _ NaR _ _ _ _ _ = NaR@@ -351,7 +362,7 @@ -- Bounded, bounded to what?!? To the ℝ! NaR is out of bounds!!! -- -- I'm bound to want this definition:-instance forall es. PositC es => Bounded (Posit es) where+instance PositC es => Bounded (Posit es) where -- 'minBound' the most negative number represented minBound = Posit mostNegVal -- 'maxBound' the most positive number represented@@ -411,7 +422,7 @@ -- ---instance forall es. (Num (IntN es), Ord (IntN es), PositC es) => FusedOps (Posit es) where+instance PositC es => FusedOps (Posit es) where -- Fused Subtract Multiply fsm = viaRational3 fsm -- Fuse Multiply Add@@ -460,7 +471,7 @@ class Convertible a b where convert :: a -> b -instance forall es1 es2. (PositC es1, PositC es2, Ord (IntN es1), Ord (IntN es2), Num (IntN es1), Num (IntN es2)) => Convertible (Posit es1) (Posit es2) where+instance (PositC es1, PositC es2) => Convertible (Posit es1) (Posit es2) where convert NaR = NaR convert (R r) = R r --@@ -484,7 +495,7 @@ -- ---instance forall es. (Show (IntN es), Ord (IntN es), Num (IntN es), PositC es) => AltShow (Posit es) where+instance PositC es => AltShow (Posit es) where displayBinary (Posit int) = displayBin int displayIntegral (Posit int) = show int@@ -494,7 +505,7 @@ displayDecimal = viaShowable (fst.fromRationalRepetendUnlimited) -- -viaShowable :: forall es a. (Show a, Ord (IntN es), Num (IntN es), PositC es) => (Rational -> a) -> Posit es -> String+viaShowable :: (Show a, PositC es) => (Rational -> a) -> Posit es -> String viaShowable _ NaR = "NaR" viaShowable f (R r) = show $ f r #endif@@ -505,7 +516,7 @@ -- ===================================================================== ---instance forall es. (PositC es) => Read (Posit es) where+instance PositC es => Read (Posit es) where readPrec = parens $ do x <- lexP@@ -528,7 +539,7 @@ -- #ifndef O_NO_STORABLE ---instance forall es. (Storable (IntN es), PositC es) => Storable (Posit es) where+instance PositC es => Storable (Posit es) where sizeOf _ = fromIntegral $ nBytes @es alignment _ = fromIntegral $ nBytes @es peek ptr = do@@ -545,12 +556,12 @@ -- ===================================================================== ---instance forall es. (Num (IntN es), Ord (IntN es), PositC es) => RealFrac (Posit es) where+instance PositC es => RealFrac (Posit es) where -- properFraction :: Integral b => a -> (b, a) properFraction = viaRationalErrTrunkation "NaR value is not a RealFrac" properFraction -- -viaRationalErrTrunkation :: forall es a. (Num (IntN es), (Ord (IntN es)), PositC es, Integral a) => String -> (Rational -> (a, Rational)) -> Posit es -> (a, Posit es)+viaRationalErrTrunkation :: PositC es => String -> (Rational -> (a, Rational)) -> Posit es -> (a, Posit es) viaRationalErrTrunkation err _ NaR = error err viaRationalErrTrunkation _ f (R r) = let (int, r') = f r@@ -560,7 +571,7 @@ -- === Real Float === -- ===================================================================== ---instance forall es. (Eq (IntN es), Ord (IntN es), Num (IntN es), Floating (Posit es), PositC es) => RealFloat (Posit es) where+instance (Floating (Posit es), PositC es) => RealFloat (Posit es) where isIEEE _ = False isDenormalized _ = False isNegativeZero _ = False@@ -773,6 +784,8 @@ | otherwise = funPhi (Posit x') where (Posit x') = (px^2 + 2*px) / (px^2 + 1)+ -- LiquidHaskell is telling me this is unsafe if px is imaginary+ -- lucky for us Posit256 is not imaginary -- calculate atan(1/2^n)@@ -1131,6 +1144,7 @@ in go next ((a' + b')^2 / (4 * t')) a' b' t' p' -- +#ifndef O_NO_SHOW -- Borwein's algorithm, with quintic convergence, -- gets to 7 ULP in 4 iterations, but really slow due to expensive function evaluations -- quite unstable and will not converge if sqrt is not accurate, which means log must be accurate@@ -1148,6 +1162,7 @@ s' = 25 / ((z + x/z + 1)^2 * s) in go a (n+1) (trace (show a') a') s' --+#endif -- Bailey–Borwein–Plouffe (BBP) formula, to 1-2 ULP, and blazing fast, converges in 60 iterations@@ -1253,7 +1268,9 @@ -- Interestingly enough, wikipedia defines two alternative solutions -- for the Shannon Wavelet, eventhough there are infinite solutions--- where the functions are equal, they are not equal+-- where the functions are equal, they are not equal. It a class of +-- functions with the charicteristic of being a band pass filter in the +-- frequency space. -- Shannon wavelet funPsiSha1 :: Posit256 -> Posit256 funPsiSha1 NaR = NaR@@ -1266,18 +1283,21 @@ funPsiSha2 t = funSinc (t/2) * cos (3*pi*t/2) -- --- Shannon wavelet+-- Shannon wavelet, same as funPsiSha1 but with a factor of pi, with the+-- Law: funPsiSha1.(pi*) === funPsiSha3+-- or : funPsiSha1 === funpsiSha3.(/pi)+-- Posit256 seems to hold to a few ULP funPsiSha3 :: Posit256 -> Posit256 funPsiSha3 NaR = NaR funPsiSha3 0 = 1 -- Why the hell not! funPsiSha3 t = let pit = pi * t- invpit = recip $ pit + invpit = recip pit in invpit * (sin (2 * pit) - sin pit) -- --- funPsiSha1.(pi*) === funPsiSha3 + -- Using the CORDIC domain reduction and some approximation function funLogDomainReduction :: (Posit256 -> Posit256) -> Posit256 -> Posit256 funLogDomainReduction _ NaR = NaR@@ -1346,6 +1366,17 @@ len = if lenA == lenB then lenA else error "Seiries Numerator and Denominator do not have the same length."++funGammaSeriesFused :: Posit256 -> Posit256+funGammaSeriesFused z = sqrt(2 * pi) * (z**(z - 0.5)) * exp (negate z) * (1 + series)+ where+ series :: Posit256+ series = fsumL $ zipWith (*) [fromRational (a % b) | (a,b) <- zip a001163 a001164] [recip $ z^n | n <- [1..len]] -- zipWith (\x y -> ) a001163 a001164+ lenA = length a001163+ lenB = length a001164+ len = if lenA == lenB+ then lenA+ else error "Seiries Numerator and Denominator do not have the same length." -- funGammaCalc :: Posit256 -> Posit256@@ -1365,4 +1396,7 @@ where x = z - 1 -+funGammaXminus1 :: Posit256 -> Posit256+funGammaXminus1 x = go (x - 1)+ where+ go z = sqrt (2 * pi) * exp z ** (negate z) * z ** (z + 0.5)
src/Posit/Internal/PositC.hs view
@@ -20,6 +20,8 @@ {-# LANGUAGE AllowAmbiguousTypes #-} -- The Haskell/GHC Type checker seems to have trouble things in the PositC class {-# LANGUAGE ScopedTypeVariables #-} -- To reduce some code duplication {-# LANGUAGE FlexibleContexts #-} -- To reduce some code duplication by claiming the type family provides some constraints, that GHC can't do without fully evaluating the type family+{-# LANGUAGE ConstrainedClassMethods #-} -- Allows constraints on class methods so default implementations of methods with Type Families can be implemented+{-# LANGUAGE ConstraintKinds #-} -- Simplify all of the constraints into a combinded constraint for the super class constraint {-# LANGUAGE CPP #-} -- To remove Storable instances to remove noise when performing analysis of Core {-# OPTIONS_GHC -Wno-unticked-promoted-constructors #-} -- Turn off noise {-# OPTIONS_GHC -Wno-type-defaults #-} -- Turn off noise@@ -33,7 +35,9 @@ module Posit.Internal.PositC (PositC(..),- ES(..)+ ES(..),+ IntN,+ FixedWidthInteger() ) where import Prelude hiding (exponent,significand)@@ -48,13 +52,14 @@ import Data.Int (Int8,Int16,Int32,Int64) -- Import standard Int sizes import Data.DoubleWord (Word128,Int128,Int256,fromHiAndLo,hiWord,loWord) -- Import large Int sizes import Data.Word (Word64)-import Data.Bits ((.|.), shiftL, shift, testBit, (.&.), shiftR)+import Data.Bits (Bits(..), (.|.), shiftL, shift, testBit, (.&.), shiftR) -- Import Naturals and Rationals {-@ embed Natural * as int @-} import GHC.Natural (Natural) -- Import the Natural Numbers ℕ (u+2115)-{-@ embed Ratio * as int @-}-import Data.Ratio (Rational, (%)) -- Import the Rational Numbers ℚ (u+211A), ℚ can get arbitrarily close to Real numbers ℝ (u+211D)+{-@ embed Ratio * as real @-}+{-@ embed Rational * as real @-}+import Data.Ratio ((%)) -- Import the Rational Numbers ℚ (u+211A), ℚ can get arbitrarily close to Real numbers ℝ (u+211D) -- | The Exponent Size 'ES' kind, the constructor for the Type is a Roman Numeral.@@ -65,43 +70,103 @@ | IV | V +-- | Type of the Finite Precision Representation, in our case Int8, +-- Int16, Int32, Int64, Int128, Int256. The 'es' of kind 'ES' will +-- determine a result of 'r' such that you can determine the 'es' by the+-- 'r'+{-@ embed IntN * as int @-}+type family IntN (es :: ES) = r | r -> es+ where+ IntN Z = Int8+ IntN I = Int16+ IntN II = Int32+ IntN III = Int64+ IntN IV = Int128+ IntN V = Int256 --- | The 'Posit' class is an approximation of ℝ, it is like a sampling on the Projective Real line ℙ(ℝ) with Maybe ℚ as the internal type.+-- | The 'FixedWidthInteger' is a Constraint Synonym that contains all+-- of the constraints provided by the 'IntN' Type Family. It is a super+-- class for the Posit Class.+type FixedWidthInteger a = + (Bits a+ ,Bounded a+ ,Enum a+ ,Integral a+ ,Eq a+ ,Ord a+ ,Num a+ ,Read a+ ,Show a+#ifndef O_NO_STORABLE+ ,Storable a+#endif+ )+++-- | The 'Posit' class is an approximation of ℝ, it is like a sampling +-- on the Projective Real line ℙ(ℝ) with Maybe ℚ as the internal type. -- The 'es' is an index that controlls the log2 word size of the Posit's -- fininte precision representation.-class PositC (es :: ES) where- -- | Type of the Finite Precision Representation, in our case Int8, Int16, Int32, Int64, Int128, Int256. The 'es' of kind 'ES' will determine a result of 'r' such that you can determine the 'es' by the 'r'- type IntN es = r | r -> es- - +class (FixedWidthInteger (IntN es)) => PositC (es :: ES) where+ -- | Transform to/from the Infinite Precision Representation encode :: Maybe Rational -> IntN es -- ^ Maybe you have some Rational Number and you want to encode it as some integer with a finite integer log2 word size.+ encode Nothing = unReal @es+ encode (Just 0) = 0+ encode (Just r)+ | r > maxPosRat @es = mostPosVal @es+ | r < minNegRat @es = mostNegVal @es+ | r > 0 && r < minPosRat @es = leastPosVal @es+ | r < 0 && r > maxNegRat @es = leastNegVal @es+ | otherwise = buildIntRep @es r+ decode :: IntN es -> Maybe Rational -- ^ You have an integer with a finite integer log2 word size decode it and Maybe it is Rational- + decode int+ | int == unReal @es = Nothing+ | int == 0 = Just 0+ | otherwise =+ let sgn = int < 0+ int' = if sgn+ then negate int+ else int+ (regime,nR) = regime2Integer @es int'+ exponent = exponent2Nat @es nR int' -- if no e or some bits missing, then they are considered zero+ rat = fraction2Posit @es nR int' -- if no fraction or some bits missing, then the missing bits are zero, making the significand p=1+ in tupPosit2Posit @es (sgn,regime,exponent,rat)+ + -- | Exponent Size based on the Posit Exponent kind ES exponentSize :: Natural -- ^ The exponent size, 'es' is a Natural number- + -- | Various other size definitions used in the Posit format with their default definitions nBytes :: Natural -- ^ 'nBytes' the number of bytes of the Posit Representation nBytes = 2^(exponentSize @es)- + nBits :: Natural -- ^ 'nBits' the number of bits of the Posit Representation nBits = 8 * (nBytes @es)- + signBitSize :: Natural -- ^ 'signBitSize' the size of the sign bit signBitSize = 1- + uSeed :: Natural -- ^ 'uSeed' scaling factor for the regime of the Posit Representation uSeed = 2^(nBytes @es)- + -- | Integer Representation of common bounds unReal :: IntN es -- ^ 'unReal' is something that is not Real, the integer value that is not a Real number- + unReal = minBound @(IntN es)+ mostPosVal :: IntN es+ mostPosVal = maxBound @(IntN es)+ leastPosVal :: IntN es+ leastPosVal = 1+ leastNegVal :: IntN es+ leastNegVal = -1+ mostNegVal :: IntN es- + mostNegVal = negate mostPosVal+ -- Rational Value of common bounds maxPosRat :: Rational maxPosRat = 2^((nBytes @es) * ((nBits @es) - 2)) % 1@@ -111,9 +176,9 @@ maxNegRat = negate (minPosRat @es) minNegRat :: Rational minNegRat = negate (maxPosRat @es)- + -- Functions to support encode and decode- + -- log base uSeed -- After calculating the regime the rational should be in the range [1,uSeed), it starts with (0,rational) log_uSeed :: (Integer, Rational) -> (Integer, Rational)@@ -121,91 +186,47 @@ | r < 1 = log_uSeed @es (regime-1,r * fromRational (toInteger (uSeed @es) % 1)) | r >= fromRational (toInteger (uSeed @es) % 1) = log_uSeed @es (regime+1,r * fromRational (1 % toInteger (uSeed @es))) | otherwise = (regime,r)- + getRegime :: Rational -> (Integer, Rational) getRegime r = log_uSeed @es (0,r)- + posit2TupPosit :: Rational -> (Bool, Integer, Natural, Rational) posit2TupPosit r = let (sgn,r') = getSign r -- returns the sign and a positive rational (regime,r'') = getRegime @es r' -- returns the regime and a rational between uSeed^-1 to uSeed^1 (exponent,significand) = getExponent r'' -- returns the exponent and a rational between [1,2), the significand in (sgn,regime,exponent,significand)- + buildIntRep :: Rational -> IntN es- mkIntRep :: Integer -> Natural -> Rational -> IntN es- formRegime :: Integer -> (IntN es, Integer)- formExponent :: Natural -> Integer -> (IntN es, Integer)- formFraction :: Rational -> Integer -> IntN es- - tupPosit2Posit :: (Bool,Integer,Natural,Rational) -> Maybe Rational- tupPosit2Posit (sgn,regime,exponent,rat) = -- s = isNeg posit == True- let pow2 = toRational (uSeed @es)^^regime * 2^exponent- scale = if sgn- then negate pow2- else pow2- in Just $ scale * rat- - regime2Integer :: IntN es -> (Integer, Int)- findRegimeFormat :: IntN es -> Bool- countRegimeBits :: Bool -> IntN es -> Int- exponent2Nat :: Int -> IntN es -> Natural- fraction2Posit :: Int -> IntN es -> Rational- - -- prints out the IntN es value in 0b... format- displayBin :: IntN es -> String- -- decimal Precision- decimalPrec :: Int- decimalPrec = fromIntegral $ 2 * (nBytes @es) + 1----instance PositC Z where- type IntN Z = Int8- exponentSize = 0- - -- Posit Integer Rep of various values- unReal = minBound @Int8- - mostPosVal = maxBound @Int8- leastPosVal = 1- leastNegVal = -1- mostNegVal = negate mostPosVal- - encode Nothing = unReal @Z- encode (Just 0) = 0- encode (Just r)- | r > maxPosRat @Z = mostPosVal @Z- | r < minNegRat @Z = mostNegVal @Z- | r > 0 && r < minPosRat @Z = leastPosVal @Z- | r < 0 && r > maxNegRat @Z = leastNegVal @Z- | otherwise = buildIntRep @Z r- buildIntRep r =- let (signBit,regime,exponent,significand) = posit2TupPosit @Z r- intRep = mkIntRep @Z regime exponent significand+ let (signBit,regime,exponent,significand) = posit2TupPosit @es r+ intRep = mkIntRep @es regime exponent significand in if signBit then negate intRep else intRep- + + mkIntRep :: Integer -> Natural -> Rational -> IntN es mkIntRep regime exponent significand =- let (regime', offset) = formRegime @Z regime -- offset is the number of binary digits remaining after the regime is formed- (exponent', offset') = formExponent @Z exponent offset -- offset' is the number of binary digits remaining after the exponent is formed- fraction = formFraction @Z significand offset'+ let (regime', offset) = formRegime @es regime -- offset is the number of binary digits remaining after the regime is formed+ (exponent', offset') = formExponent @es exponent offset -- offset' is the number of binary digits remaining after the exponent is formed+ fraction = formFraction @es significand offset' in regime' .|. exponent' .|. fraction- + + formRegime :: Integer -> (IntN es, Integer) formRegime power | 0 <= power =- let offset = (fromIntegral (nBits @Z - 1) - power - 1)+ let offset = (fromIntegral (nBits @es - 1) - power - 1) in (fromIntegral (2^(power + 1) - 1) `shiftL` fromInteger offset, offset - 1) | otherwise =- let offset = (fromIntegral (nBits @Z - 1) - abs power - 1)+ let offset = (fromIntegral (nBits @es - 1) - abs power - 1) in (1 `shiftL` fromInteger offset, offset)- + + formExponent :: Natural -> Integer -> (IntN es, Integer) formExponent power offset =- let offset' = offset - fromIntegral (exponentSize @Z)+ let offset' = offset - fromIntegral (exponentSize @es) in (fromIntegral power `shift` fromInteger offset', offset')- + + formFraction :: Rational -> Integer -> IntN es formFraction r offset = let numFractionBits = offset fractionSize = 2^numFractionBits@@ -213,63 +234,65 @@ in if numFractionBits >= 1 then fromInteger normFraction else 0- - decode int- | int == unReal @Z = Nothing- | int == 0 = Just 0- | otherwise =- let sgn = int < 0- int' = if sgn- then negate int- else int- (regime,nR) = regime2Integer @Z int'- exponent = exponent2Nat @Z nR int' -- if no e or some bits missing, then they are considered zero- rat = fraction2Posit @Z nR int' -- if no fraction or some bits missing, then the missing bits are zero, making the significand p=1- in tupPosit2Posit @Z (sgn,regime,exponent,rat)- + + tupPosit2Posit :: (Bool,Integer,Natural,Rational) -> Maybe Rational+ tupPosit2Posit (sgn,regime,exponent,rat) = -- s = isNeg posit == True+ let pow2 = toRational (uSeed @es)^^regime * 2^exponent+ scale = if sgn+ then negate pow2+ else pow2+ in Just $ scale * rat+ + regime2Integer :: IntN es -> (Integer, Int) regime2Integer posit =- let regimeFormat = findRegimeFormat @Z posit- regimeCount = countRegimeBits @Z regimeFormat posit+ let regimeFormat = findRegimeFormat @es posit+ regimeCount = countRegimeBits @es regimeFormat posit regime = calcRegimeInt regimeFormat regimeCount in (regime, regimeCount + 1) -- a rational representation of the regime and the regimeCount plus rBar which is the numBitsRegime- + -- will return the format of the regime, either HI or LO; it could get refactored in the future -- True means a 1 is the first bit in the regime- findRegimeFormat posit = testBit posit (fromIntegral (nBits @Z) - 1 - fromIntegral (signBitSize @Z))- - countRegimeBits format posit = go (fromIntegral (nBits @Z) - 1 - fromIntegral (signBitSize @Z)) 0+ findRegimeFormat :: IntN es -> Bool+ findRegimeFormat posit = testBit posit (fromIntegral (nBits @es) - 1 - fromIntegral (signBitSize @es))+ + countRegimeBits :: Bool -> IntN es -> Int+ countRegimeBits format posit = go (fromIntegral (nBits @es) - 1 - fromIntegral (signBitSize @es)) 0 where go (-1) acc = acc go index acc | xnor format (testBit posit index) = go (index - 1) (acc + 1) | otherwise = acc- + -- knowing the number of the regime bits, and the sign bit we can extract -- the exponent. We mask to the left of the exponent to remove the sign and regime, and -- then shift to the right to remove the fraction.+ exponent2Nat :: Int -> IntN es -> Natural exponent2Nat numBitsRegime posit =- let bitsRemaining = fromIntegral (nBits @Z) - numBitsRegime - fromIntegral (signBitSize @Z)+ let bitsRemaining = fromIntegral (nBits @es) - numBitsRegime - fromIntegral (signBitSize @es) signNRegimeMask = 2^bitsRemaining - 1 int = posit .&. signNRegimeMask- nBitsToTheRight = fromIntegral (nBits @Z) - numBitsRegime - fromIntegral (signBitSize @Z) - fromIntegral (exponentSize @Z)+ nBitsToTheRight = fromIntegral (nBits @es) - numBitsRegime - fromIntegral (signBitSize @es) - fromIntegral (exponentSize @es) in if bitsRemaining <=0 then 0 else if nBitsToTheRight < 0 then fromIntegral $ int `shiftL` negate nBitsToTheRight else fromIntegral $ int `shiftR` nBitsToTheRight- + -- knowing the number of the regime bits, sign bit, and the number of the -- exponent bits we can extract the fraction. We mask to the left of the fraction to -- remove the sign, regime, and exponent. If there is no fraction then the value is 1.+ fraction2Posit :: Int -> IntN es -> Rational fraction2Posit numBitsRegime posit =- let offset = fromIntegral $ (signBitSize @Z) + fromIntegral numBitsRegime + (exponentSize @Z)- fractionSize = fromIntegral (nBits @Z) - offset+ let offset = fromIntegral $ (signBitSize @es) + fromIntegral numBitsRegime + (exponentSize @es)+ fractionSize = fromIntegral (nBits @es) - offset fractionBits = posit .&. (2^fractionSize - 1) in if fractionSize >= 1 then (2^fractionSize + toInteger fractionBits) % 2^fractionSize else 1 % 1- - displayBin int = "0b" ++ go (fromIntegral (nBits @Z) - 1)+ + -- prints out the IntN es value in 0b... format+ displayBin :: IntN es -> String+ displayBin int = "0b" ++ go (fromIntegral (nBits @es) - 1) where go :: Int -> String go 0 = if testBit int 0@@ -277,614 +300,44 @@ else "0" go idx = if testBit int idx then "1" ++ go (idx - 1)- else "0" ++ go (idx -1)+ else "0" ++ go (idx - 1)+ + -- decimal Precision+ decimalPrec :: Int+ decimalPrec = fromIntegral $ 2 * (nBytes @es) + 1+ + {-# MINIMAL exponentSize #-} +-- =====================================================================+-- === PositC Instances ===+-- ===================================================================== +instance PositC Z where+ exponentSize = 0++ instance PositC I where- type IntN I = Int16 exponentSize = 1- - -- Posit Integer Rep of various values- unReal = minBound @Int16- - mostPosVal = maxBound @Int16- leastPosVal = 1- leastNegVal = -1- mostNegVal = negate mostPosVal- - encode Nothing = unReal @I- encode (Just 0) = 0- encode (Just r)- | r > maxPosRat @I = mostPosVal @I- | r < minNegRat @I = mostNegVal @I- | r > 0 && r < minPosRat @I = leastPosVal @I- | r < 0 && r > maxNegRat @I = leastNegVal @I- | otherwise = buildIntRep @I r- - buildIntRep r =- let (signBit,regime,exponent,significand) = posit2TupPosit @I r- intRep = mkIntRep @I regime exponent significand- in if signBit- then negate intRep- else intRep- - mkIntRep regime exponent significand =- let (regime', offset) = formRegime @I regime -- offset is the number of binary digits remaining after the regime is formed- (exponent', offset') = formExponent @I exponent offset -- offset' is the number of binary digits remaining after the exponent is formed- fraction = formFraction @I significand offset'- in regime' .|. exponent' .|. fraction- - formRegime power- | 0 <= power =- let offset = (fromIntegral (nBits @I - 1) - power - 1)- in (fromIntegral (2^(power + 1) - 1) `shiftL` fromInteger offset, offset - 1)- | otherwise =- let offset = (fromIntegral (nBits @I - 1) - abs power - 1)- in (1 `shiftL` fromInteger offset, offset)- - formExponent power offset =- let offset' = offset - fromIntegral (exponentSize @I)- in (fromIntegral power `shift` fromInteger offset', offset')- - formFraction r offset =- let numFractionBits = offset- fractionSize = 2^numFractionBits- normFraction = round $ (r - 1) * fractionSize -- "posit - 1" is changing it from the significand to the fraction: [1,2) -> [0,1)- in if numFractionBits >= 1- then fromInteger normFraction- else 0- - decode int- | int == unReal @I = Nothing- | int == 0 = Just 0- | otherwise =- let sgn = int < 0- int' = if sgn- then negate int- else int- (regime,nR) = regime2Integer @I int'- exponent = exponent2Nat @I nR int' -- if no e or some bits missing, then they are considered zero- rat = fraction2Posit @I nR int' -- if no fraction or some bits missing, then the missing bits are zero, making the significand p=1- in tupPosit2Posit @I (sgn,regime,exponent,rat)- - regime2Integer posit =- let regimeFormat = findRegimeFormat @I posit- regimeCount = countRegimeBits @I regimeFormat posit- regime = calcRegimeInt regimeFormat regimeCount- in (regime, regimeCount + 1) -- a rational representation of the regime and the regimeCount plus rBar which is the numBitsRegime- - -- will return the format of the regime, either HI or LO; it could get refactored in the future- -- True means a 1 is the first bit in the regime- findRegimeFormat posit = testBit posit (fromIntegral (nBits @I) - 1 - fromIntegral (signBitSize @I))- - countRegimeBits format posit = go (fromIntegral (nBits @I) - 1 - fromIntegral (signBitSize @I)) 0- where- go (-1) acc = acc- go index acc- | xnor format (testBit posit index) = go (index - 1) (acc + 1)- | otherwise = acc- - -- knowing the number of the regime bits, and the sign bit we can extract- -- the exponent. We mask to the left of the exponent to remove the sign and regime, and- -- then shift to the right to remove the fraction.- exponent2Nat numBitsRegime posit =- let bitsRemaining = fromIntegral (nBits @I) - numBitsRegime - fromIntegral (signBitSize @I)- signNRegimeMask = 2^bitsRemaining - 1- int = posit .&. signNRegimeMask- nBitsToTheRight = fromIntegral (nBits @I) - numBitsRegime - fromIntegral (signBitSize @I) - fromIntegral (exponentSize @I)- in if bitsRemaining <=0- then 0- else if nBitsToTheRight < 0- then fromIntegral $ int `shiftL` negate nBitsToTheRight- else fromIntegral $ int `shiftR` nBitsToTheRight- - -- knowing the number of the regime bits, sign bit, and the number of the- -- exponent bits we can extract the fraction. We mask to the left of the fraction to- -- remove the sign, regime, and exponent. If there is no fraction then the value is 1.- fraction2Posit numBitsRegime posit =- let offset = fromIntegral $ (signBitSize @I) + fromIntegral numBitsRegime + (exponentSize @I)- fractionSize = fromIntegral (nBits @I) - offset- fractionBits = posit .&. (2^fractionSize - 1)- in if fractionSize >= 1- then (2^fractionSize + toInteger fractionBits) % 2^fractionSize- else 1 % 1- - displayBin int = "0b" ++ go (fromIntegral (nBits @I) - 1)- where- go :: Int -> String- go 0 = if testBit int 0- then "1"- else "0"- go idx = if testBit int idx- then "1" ++ go (idx - 1)- else "0" ++ go (idx -1) - instance PositC II where- type IntN II = Int32 exponentSize = 2- - -- Posit Integer Rep of various values- unReal = minBound @Int32- - mostPosVal = maxBound @Int32- leastPosVal = 1- leastNegVal = -1- mostNegVal = negate mostPosVal- - encode Nothing = unReal @II- encode (Just 0) = 0- encode (Just r)- | r > maxPosRat @II = mostPosVal @II- | r < minNegRat @II = mostNegVal @II- | r > 0 && r < minPosRat @II = leastPosVal @II- | r < 0 && r > maxNegRat @II = leastNegVal @II- | otherwise = buildIntRep @II r- - buildIntRep r =- let (signBit,regime,exponent,significand) = posit2TupPosit @II r- intRep = mkIntRep @II regime exponent significand- in if signBit- then negate intRep- else intRep- - mkIntRep regime exponent significand =- let (regime', offset) = formRegime @II regime -- offset is the number of binary digits remaining after the regime is formed- (exponent', offset') = formExponent @II exponent offset -- offset' is the number of binary digits remaining after the exponent is formed- fraction = formFraction @II significand offset'- in regime' .|. exponent' .|. fraction- - formRegime power- | 0 <= power =- let offset = (fromIntegral (nBits @II - 1) - power - 1)- in (fromIntegral (2^(power + 1) - 1) `shiftL` fromInteger offset, offset - 1)- | otherwise =- let offset = (fromIntegral (nBits @II - 1) - abs power - 1)- in (1 `shiftL` fromInteger offset, offset)- - formExponent power offset =- let offset' = offset - fromIntegral (exponentSize @II)- in (fromIntegral power `shift` fromInteger offset', offset')- - formFraction r offset =- let numFractionBits = offset- fractionSize = 2^numFractionBits- normFraction = round $ (r - 1) * fractionSize -- "posit - 1" is changing it from the significand to the fraction: [1,2) -> [0,1)- in if numFractionBits >= 1- then fromInteger normFraction- else 0- - decode int- | int == unReal @II = Nothing- | int == 0 = Just 0- | otherwise =- let sgn = int < 0- int' = if sgn- then negate int- else int- (regime,nR) = regime2Integer @II int'- exponent = exponent2Nat @II nR int' -- if no e or some bits missing, then they are considered zero- rat = fraction2Posit @II nR int' -- if no fraction or some bits missing, then the missing bits are zero, making the significand p=1- in tupPosit2Posit @II (sgn,regime,exponent,rat)- - regime2Integer posit =- let regimeFormat = findRegimeFormat @II posit- regimeCount = countRegimeBits @II regimeFormat posit- regime = calcRegimeInt regimeFormat regimeCount- in (regime, regimeCount + 1) -- a rational representation of the regime and the regimeCount plus rBar which is the numBitsRegime- - -- will return the format of the regime, either HI or LO; it could get refactored in the future- -- True means a 1 is the first bit in the regime- findRegimeFormat posit = testBit posit (fromIntegral (nBits @II) - 1 - fromIntegral (signBitSize @II))- - countRegimeBits format posit = go (fromIntegral (nBits @II) - 1 - fromIntegral (signBitSize @II)) 0- where- go (-1) acc = acc- go index acc- | xnor format (testBit posit index) = go (index - 1) (acc + 1)- | otherwise = acc- - -- knowing the number of the regime bits, and the sign bit we can extract- -- the exponent. We mask to the left of the exponent to remove the sign and regime, and- -- then shift to the right to remove the fraction.- exponent2Nat numBitsRegime posit =- let bitsRemaining = fromIntegral (nBits @II) - numBitsRegime - fromIntegral (signBitSize @II)- signNRegimeMask = 2^bitsRemaining - 1- int = posit .&. signNRegimeMask- nBitsToTheRight = fromIntegral (nBits @II) - numBitsRegime - fromIntegral (signBitSize @II) - fromIntegral (exponentSize @II)- in if bitsRemaining <=0- then 0- else if nBitsToTheRight < 0- then fromIntegral $ int `shiftL` negate nBitsToTheRight- else fromIntegral $ int `shiftR` nBitsToTheRight- - -- knowing the number of the regime bits, sign bit, and the number of the- -- exponent bits we can extract the fraction. We mask to the left of the fraction to- -- remove the sign, regime, and exponent. If there is no fraction then the value is 1.- fraction2Posit numBitsRegime posit =- let offset = fromIntegral $ (signBitSize @II) + fromIntegral numBitsRegime + (exponentSize @II)- fractionSize = fromIntegral (nBits @II) - offset- fractionBits = posit .&. (2^fractionSize - 1)- in if fractionSize >= 1- then (2^fractionSize + toInteger fractionBits) % 2^fractionSize- else 1 % 1- - displayBin int = "0b" ++ go (fromIntegral (nBits @II) - 1)- where- go :: Int -> String- go 0 = if testBit int 0- then "1"- else "0"- go idx = if testBit int idx- then "1" ++ go (idx - 1)- else "0" ++ go (idx -1) - instance PositC III where- type IntN III = Int64 exponentSize = 3- - -- Posit Integer Rep of various values- unReal = minBound @Int64- - mostPosVal = maxBound @Int64- leastPosVal = 1- leastNegVal = -1- mostNegVal = negate mostPosVal- - encode Nothing = unReal @III- encode (Just 0) = 0- encode (Just r)- | r > maxPosRat @III = mostPosVal @III- | r < minNegRat @III = mostNegVal @III- | r > 0 && r < minPosRat @III = leastPosVal @III- | r < 0 && r > maxNegRat @III = leastNegVal @III- | otherwise = buildIntRep @III r- - buildIntRep r =- let (signBit,regime,exponent,significand) = posit2TupPosit @III r- intRep = mkIntRep @III regime exponent significand- in if signBit- then negate intRep- else intRep- - mkIntRep regime exponent significand =- let (regime', offset) = formRegime @III regime -- offset is the number of binary digits remaining after the regime is formed- (exponent', offset') = formExponent @III exponent offset -- offset' is the number of binary digits remaining after the exponent is formed- fraction = formFraction @III significand offset'- in regime' .|. exponent' .|. fraction- - formRegime power- | 0 <= power =- let offset = (fromIntegral (nBits @III - 1) - power - 1)- in (fromIntegral (2^(power + 1) - 1) `shiftL` fromInteger offset, offset - 1)- | otherwise =- let offset = (fromIntegral (nBits @III - 1) - abs power - 1)- in (1 `shiftL` fromInteger offset, offset)- - formExponent power offset =- let offset' = offset - fromIntegral (exponentSize @III)- in (fromIntegral power `shift` fromInteger offset', offset')- - formFraction r offset =- let numFractionBits = offset- fractionSize = 2^numFractionBits- normFraction = round $ (r - 1) * fractionSize -- "posit - 1" is changing it from the significand to the fraction: [1,2) -> [0,1)- in if numFractionBits >= 1- then fromInteger normFraction- else 0- - decode int- | int == unReal @III = Nothing- | int == 0 = Just 0- | otherwise =- let sgn = int < 0- int' = if sgn- then negate int- else int- (regime,nR) = regime2Integer @III int'- exponent = exponent2Nat @III nR int' -- if no e or some bits missing, then they are considered zero- rat = fraction2Posit @III nR int' -- if no fraction or some bits missing, then the missing bits are zero, making the significand p=1- in tupPosit2Posit @III (sgn,regime,exponent,rat)- - regime2Integer posit =- let regimeFormat = findRegimeFormat @III posit- regimeCount = countRegimeBits @III regimeFormat posit- regime = calcRegimeInt regimeFormat regimeCount- in (regime, regimeCount + 1) -- a rational representation of the regime and the regimeCount plus rBar which is the numBitsRegime- - -- will return the format of the regime, either HI or LO; it could get refactored in the future- -- True means a 1 is the first bit in the regime- findRegimeFormat posit = testBit posit (fromIntegral (nBits @III) - 1 - fromIntegral (signBitSize @III))- - countRegimeBits format posit = go (fromIntegral (nBits @III) - 1 - fromIntegral (signBitSize @III)) 0- where- go (-1) acc = acc- go index acc- | xnor format (testBit posit index) = go (index - 1) (acc + 1)- | otherwise = acc- - -- knowing the number of the regime bits, and the sign bit we can extract- -- the exponent. We mask to the left of the exponent to remove the sign and regime, and- -- then shift to the right to remove the fraction.- exponent2Nat numBitsRegime posit =- let bitsRemaining = fromIntegral (nBits @III) - numBitsRegime - fromIntegral (signBitSize @III)- signNRegimeMask = 2^bitsRemaining - 1- int = posit .&. signNRegimeMask- nBitsToTheRight = fromIntegral (nBits @III) - numBitsRegime - fromIntegral (signBitSize @III) - fromIntegral (exponentSize @III)- in if bitsRemaining <=0- then 0- else if nBitsToTheRight < 0- then fromIntegral $ int `shiftL` negate nBitsToTheRight- else fromIntegral $ int `shiftR` nBitsToTheRight- - -- knowing the number of the regime bits, sign bit, and the number of the- -- exponent bits we can extract the fraction. We mask to the left of the fraction to- -- remove the sign, regime, and exponent. If there is no fraction then the value is 1.- fraction2Posit numBitsRegime posit =- let offset = fromIntegral $ (signBitSize @III) + fromIntegral numBitsRegime + (exponentSize @III)- fractionSize = fromIntegral (nBits @III) - offset- fractionBits = posit .&. (2^fractionSize - 1)- in if fractionSize >= 1- then (2^fractionSize + toInteger fractionBits) % 2^fractionSize- else 1 % 1- - displayBin int = "0b" ++ go (fromIntegral (nBits @III) - 1)- where- go :: Int -> String- go 0 = if testBit int 0- then "1"- else "0"- go idx = if testBit int idx- then "1" ++ go (idx - 1)- else "0" ++ go (idx -1) - instance PositC IV where- type IntN IV = Int128 exponentSize = 4- - -- Posit Integer Rep of various values- unReal = minBound @Int128- - mostPosVal = maxBound @Int128- leastPosVal = 1- leastNegVal = -1- mostNegVal = negate mostPosVal- - encode Nothing = unReal @IV- encode (Just 0) = 0- encode (Just r)- | r > maxPosRat @IV = mostPosVal @IV- | r < minNegRat @IV = mostNegVal @IV- | r > 0 && r < minPosRat @IV = leastPosVal @IV- | r < 0 && r > maxNegRat @IV = leastNegVal @IV- | otherwise = buildIntRep @IV r- - buildIntRep r =- let (signBit,regime,exponent,significand) = posit2TupPosit @IV r- intRep = mkIntRep @IV regime exponent significand- in if signBit- then negate intRep- else intRep- - mkIntRep regime exponent significand =- let (regime', offset) = formRegime @IV regime -- offset is the number of binary digits remaining after the regime is formed- (exponent', offset') = formExponent @IV exponent offset -- offset' is the number of binary digits remaining after the exponent is formed- fraction = formFraction @IV significand offset'- in regime' .|. exponent' .|. fraction- - formRegime power- | 0 <= power =- let offset = (fromIntegral (nBits @IV - 1) - power - 1)- in (fromIntegral (2^(power + 1) - 1) `shiftL` fromInteger offset, offset - 1)- | otherwise =- let offset = (fromIntegral (nBits @IV - 1) - abs power - 1)- in (1 `shiftL` fromInteger offset, offset)- - formExponent power offset =- let offset' = offset - fromIntegral (exponentSize @IV)- in (fromIntegral power `shift` fromInteger offset', offset')- - formFraction r offset =- let numFractionBits = offset- fractionSize = 2^numFractionBits- normFraction = round $ (r - 1) * fractionSize -- "posit - 1" is changing it from the significand to the fraction: [1,2) -> [0,1)- in if numFractionBits >= 1- then fromInteger normFraction- else 0- - decode int- | int == unReal @IV = Nothing- | int == 0 = Just 0- | otherwise =- let sgn = int < 0- int' = if sgn- then negate int- else int- (regime,nR) = regime2Integer @IV int'- exponent = exponent2Nat @IV nR int' -- if no e or some bits missing, then they are considered zero- rat = fraction2Posit @IV nR int' -- if no fraction or some bits missing, then the missing bits are zero, making the significand p=1- in tupPosit2Posit @IV (sgn,regime,exponent,rat)- - regime2Integer posit =- let regimeFormat = findRegimeFormat @IV posit- regimeCount = countRegimeBits @IV regimeFormat posit- regime = calcRegimeInt regimeFormat regimeCount- in (regime, regimeCount + 1) -- a rational representation of the regime and the regimeCount plus rBar which is the numBitsRegime- - -- will return the format of the regime, either HI or LO; it could get refactored in the future- -- True means a 1 is the first bit in the regime- findRegimeFormat posit = testBit posit (fromIntegral (nBits @IV) - 1 - fromIntegral (signBitSize @IV))- - countRegimeBits format posit = go (fromIntegral (nBits @IV) - 1 - fromIntegral (signBitSize @IV)) 0- where- go (-1) acc = acc- go index acc- | xnor format (testBit posit index) = go (index - 1) (acc + 1)- | otherwise = acc- - -- knowing the number of the regime bits, and the sign bit we can extract- -- the exponent. We mask to the left of the exponent to remove the sign and regime, and- -- then shift to the right to remove the fraction.- exponent2Nat numBitsRegime posit =- let bitsRemaining = fromIntegral (nBits @IV) - numBitsRegime - fromIntegral (signBitSize @IV)- signNRegimeMask = 2^bitsRemaining - 1- int = posit .&. signNRegimeMask- nBitsToTheRight = fromIntegral (nBits @IV) - numBitsRegime - fromIntegral (signBitSize @IV) - fromIntegral (exponentSize @IV)- in if bitsRemaining <=0- then 0- else if nBitsToTheRight < 0- then fromIntegral $ int `shiftL` negate nBitsToTheRight- else fromIntegral $ int `shiftR` nBitsToTheRight- - -- knowing the number of the regime bits, sign bit, and the number of the- -- exponent bits we can extract the fraction. We mask to the left of the fraction to- -- remove the sign, regime, and exponent. If there is no fraction then the value is 1.- fraction2Posit numBitsRegime posit =- let offset = fromIntegral $ (signBitSize @IV) + fromIntegral numBitsRegime + (exponentSize @IV)- fractionSize = fromIntegral (nBits @IV) - offset- fractionBits = posit .&. (2^fractionSize - 1)- in if fractionSize >= 1- then (2^fractionSize + toInteger fractionBits) % 2^fractionSize- else 1 % 1- - displayBin int = "0b" ++ go (fromIntegral (nBits @IV) - 1)- where- go :: Int -> String- go 0 = if testBit int 0- then "1"- else "0"- go idx = if testBit int idx- then "1" ++ go (idx - 1)- else "0" ++ go (idx -1) - instance PositC V where- type IntN V = Int256 exponentSize = 5- - -- Posit Integer Rep of various values- unReal = minBound @Int256- - mostPosVal = maxBound @Int256- leastPosVal = 1- leastNegVal = -1- mostNegVal = negate mostPosVal- - encode Nothing = unReal @V- encode (Just 0) = 0- encode (Just r)- | r > maxPosRat @V = mostPosVal @V- | r < minNegRat @V = mostNegVal @V- | r > 0 && r < minPosRat @V = leastPosVal @V- | r < 0 && r > maxNegRat @V = leastNegVal @V- | otherwise = buildIntRep @V r- - buildIntRep r =- let (signBit,regime,exponent,significand) = posit2TupPosit @V r- intRep = mkIntRep @V regime exponent significand- in if signBit- then negate intRep- else intRep- - mkIntRep regime exponent significand =- let (regime', offset) = formRegime @V regime -- offset is the number of binary digits remaining after the regime is formed- (exponent', offset') = formExponent @V exponent offset -- offset' is the number of binary digits remaining after the exponent is formed- fraction = formFraction @V significand offset'- in regime' .|. exponent' .|. fraction- - formRegime power- | 0 <= power =- let offset = (fromIntegral (nBits @V - 1) - power - 1)- in (fromIntegral (2^(power + 1) - 1) `shiftL` fromInteger offset, offset - 1)- | otherwise =- let offset = (fromIntegral (nBits @V - 1) - abs power - 1)- in (1 `shiftL` fromInteger offset, offset)- - formExponent power offset =- let offset' = offset - fromIntegral (exponentSize @V)- in (fromIntegral power `shift` fromInteger offset', offset')- - formFraction r offset =- let numFractionBits = offset- fractionSize = 2^numFractionBits- normFraction = round $ (r - 1) * fractionSize -- "posit - 1" is changing it from the significand to the fraction: [1,2) -> [0,1)- in if numFractionBits >= 1- then fromInteger normFraction- else 0- - decode int- | int == unReal @V = Nothing- | int == 0 = Just 0- | otherwise =- let sgn = int < 0- int' = if sgn- then negate int- else int- (regime,nR) = regime2Integer @V int'- exponent = exponent2Nat @V nR int' -- if no e or some bits missing, then they are considered zero- rat = fraction2Posit @V nR int' -- if no fraction or some bits missing, then the missing bits are zero, making the significand p=1- in tupPosit2Posit @V (sgn,regime,exponent,rat)- - regime2Integer posit =- let regimeFormat = findRegimeFormat @V posit- regimeCount = countRegimeBits @V regimeFormat posit- regime = calcRegimeInt regimeFormat regimeCount- in (regime, regimeCount + 1) -- a rational representation of the regime and the regimeCount plus rBar which is the numBitsRegime- - -- will return the format of the regime, either HI or LO; it could get refactored in the future- -- True means a 1 is the first bit in the regime- findRegimeFormat posit = testBit posit (fromIntegral (nBits @V) - 1 - fromIntegral (signBitSize @V))- - countRegimeBits format posit = go (fromIntegral (nBits @V) - 1 - fromIntegral (signBitSize @V)) 0- where- go (-1) acc = acc- go index acc- | xnor format (testBit posit index) = go (index - 1) (acc + 1)- | otherwise = acc- - -- knowing the number of the regime bits, and the sign bit we can extract- -- the exponent. We mask to the left of the exponent to remove the sign and regime, and- -- then shift to the right to remove the fraction.- exponent2Nat numBitsRegime posit =- let bitsRemaining = fromIntegral (nBits @V) - numBitsRegime - fromIntegral (signBitSize @V)- signNRegimeMask = 2^bitsRemaining - 1- int = posit .&. signNRegimeMask- nBitsToTheRight = fromIntegral (nBits @V) - numBitsRegime - fromIntegral (signBitSize @V) - fromIntegral (exponentSize @V)- in if bitsRemaining <=0- then 0- else if nBitsToTheRight < 0- then fromIntegral $ int `shiftL` negate nBitsToTheRight- else fromIntegral $ int `shiftR` nBitsToTheRight- - -- knowing the number of the regime bits, sign bit, and the number of the- -- exponent bits we can extract the fraction. We mask to the left of the fraction to- -- remove the sign, regime, and exponent. If there is no fraction then the value is 1.- fraction2Posit numBitsRegime posit =- let offset = fromIntegral $ (signBitSize @V) + fromIntegral numBitsRegime + (exponentSize @V)- fractionSize = fromIntegral (nBits @V) - offset- fractionBits = posit .&. (2^fractionSize - 1)- in if fractionSize >= 1- then (2^fractionSize + toInteger fractionBits) % 2^fractionSize- else 1 % 1- - displayBin int = "0b" ++ go (fromIntegral (nBits @V) - 1)- where- go :: Int -> String- go 0 = if testBit int 0- then "1"- else "0"- go idx = if testBit int idx- then "1" ++ go (idx - 1)- else "0" ++ go (idx -1) + -- ===================================================================== -- === Encode and Decode Helpers === -- =====================================================================@@ -900,7 +353,7 @@ else r in (s,absPosit) -- pretty much the same as 'abs') --- Exponent should be an integer in the range of [0,uSeed), and also return the posit [1,2)+-- Exponent should be an integer in the range of [0,uSeed), and also return an exponent and a rational in the range of [1,2) getExponent :: Rational -> (Natural, Rational) getExponent r = log_2 (0,r) @@ -930,17 +383,15 @@ sizeOf _ = 16 alignment _ = 16 peek ptr = do- hi <- peek $ offsetInt 0+ hi <- peek $ offsetWord 0 lo <- peek $ offsetWord 1 return $ fromHiAndLo hi lo where- offsetInt i = (castPtr ptr :: Ptr Word64) `plusPtr` (i*8) offsetWord i = (castPtr ptr :: Ptr Word64) `plusPtr` (i*8) poke ptr int = do- poke (offsetInt 0) (hiWord int)+ poke (offsetWord 0) (hiWord int) poke (offsetWord 1) (loWord int) where- offsetInt i = (castPtr ptr :: Ptr Word64) `plusPtr` (i*8) offsetWord i = (castPtr ptr :: Ptr Word64) `plusPtr` (i*8) -- Orphan Instance for Int128 using the DoubleWord type class
test/TestPosit.hs view
@@ -19,6 +19,9 @@ main :: IO () main = do --+ print $ "exp(1)**(pi*sqrt 43): " ++ show (exp(1 :: Posit256) ** (pi * sqrt 43)) -- + print $ "exp(1)**(pi*sqrt 67): " ++ show (exp(1 :: Posit256) ** (pi * sqrt 67)) -- + print $ "exp(1)**(pi*sqrt 163): " ++ show (exp(1 :: Posit256) ** (pi * sqrt 163)) -- print $ "Machine Alpha Posit8 ~1.0: " ++ show (1.0 - succ (1.0 :: Posit8)) -- succ (Posit int) = Posit (succ int) print $ "Machine Alpha Posit16 ~1.0: " ++ show (1.0 - succ (1.0 :: Posit16)) -- print $ "Machine Alpha Posit32 ~1.0: " ++ show (1.0 - succ (1.0 :: Posit32)) -- @@ -31,6 +34,16 @@ let sqrtTuma = (funLogDomainReduction funLogTuma).(/2).(funExp2 funExpTuma).(/log 2) print $ "sqrt phi using a Tuma algorithm: " ++ show (sqrtTuma phi) print $ "Tuma is fasta: " ++ show (sqrtTaylor (1/1000000) - sqrtTuma (1/1000000))+ let truth = 0.8956731517052878608869612167009786079379812529831641161347143256836782657295966290940929214799036260987761959338755143914935872 :: Posit256+ eval "Standard: gamma(phi): " (gamma (phi)) truth+ eval "Fused Gamma: gamma(phi): " (funGammaSeriesFused (phi)) truth+ eval "Ramanujan Gamma: gamma(phi): " (funGammaRamanujan (phi)) truth+ eval "Calc Gamma: gamma(phi): " (funGammaCalc (phi)) truth+ eval "Nemes Gamma: gamma(phi): " (funGammaNemes (phi)) truth+ eval "Yang Gamma: gamma(phi): " (funGammaYang (phi)) truth+ eval "Chen Gamma: gamma(phi): " (funGammaChen (phi)) truth+ eval "Gamma (x - 1): gamma(phi): " (funGammaXminus1 (phi)) truth+ eval "Wolfram alpha: gamma(phi): " truth truth let truth = 5.0431656433600286513118821892854247103235901754138463603020001967777869609108929428415187821843384653305404495551887666992776792 :: Posit256 eval "Standard: exp(phi):" (exp (phi)) truth eval "Taylor: exp(phi):" (funExp2 funExpTaylor (phi / log 2)) truth