packages feed

posit 2022.0.1.4 → 2022.1.0.0

raw patch · 10 files changed

+300/−48 lines, 10 filesdep ~posit

Dependency ranges changed: posit

Files

ChangeLog.md view
@@ -1,8 +1,19 @@ # Changelog for Posit Numbers +# posit-2022.1.0.0++  * Optimized the `Show` instance to properly handle the tapered precision+  * Changes to round better+  * Improved accuracy of `exp`+  * Changed titles of accuracy charts to be more consistant+  * Added test to confirm `read.show == id` to perhaps optimise the textual representation, and verify that the implementation is sufficent+  * Test round trip with command to run: `stack test posit:test-posit-readShowId`+  * Added catagory to .cabal file+ # posit-2022.0.1.4    * `atan 1 :: Posit256` is very slow, added a few precice constants to speed it up+  * Improved accuracy of `sin` and `cos` at the extreme positive and negative values  # posit-2022.0.1.3
README.md view
@@ -1,4 +1,4 @@-# posit 2022.0.1.4+# posit 2022.1.0.0  The [Posit Standard 2022](https://posithub.org/docs/posit_standard-2.pdf), and [Posit Standard 3.2](https://posithub.org/docs/posit_standard.pdf),
posit.cabal view
@@ -1,7 +1,8 @@ cabal-version: 1.12  name:           posit-version:        2022.0.1.4+version:        2022.1.0.0+category:       Numeric, Math description:    The Posit Number format attempting to conform to the Posit Standard Versions 3.2 and 2022.  Where Real numbers are approximated by `Maybe Rational` and sampled in a similar way to the projective real line. homepage:       https://github.com/waivio/posit#readme bug-reports:    https://github.com/waivio/posit/issues@@ -15,7 +16,7 @@                      GHC == 8.10.7,                      GHC == 9.0.2,                      GHC == 9.2.8,-                     GHC == 9.4.5,+                     GHC == 9.4.6,                      GHC == 9.6.2 synopsis:       Posit Numbers extra-source-files:@@ -113,8 +114,19 @@   main-is: WeighPosit.hs   ghc-options: -Wall -O2   build-depends:-    base >=4.7 && <5,+    base >= 4.7 && < 5,     posit,     vector,     weigh+  default-language: Haskell2010++-- Test `(show.read == id)`+benchmark test-posit-readShowId+  type: exitcode-stdio-1.0+  hs-source-dirs: test+  main-is: TestReadShow.hs+  ghc-options: -Wall -O2+  build-depends:+    base >= 4.7 && < 5,+    posit   default-language: Haskell2010
src/Posit.hs view
@@ -201,7 +201,7 @@ -- instance PositC es => Show (Posit es) where   show NaR = "NaR"-  show (R r) = formatScientific Generic (Just $ decimalPrec @es) (fst.fromRationalRepetendUnlimited $ r)+  show p@(Posit int) = formatScientific Generic (Just $ decimalPrec @es int) (fst.fromRationalRepetendUnlimited $ toRational p) -- #endif @@ -710,11 +710,14 @@   approx_exp :: forall es. PositC es => Posit es -> Posit es     -- Comment by Abigale Emily:  xcddfffff-approx_exp (R bx) = fromIntegral (uSeed @es)^^m * (taylor_approx_expm1 (R c * log_USeed) + 1)+approx_exp NaR = NaR+approx_exp (R x) = 2^^k * funExpTaylor (R r)   where-    (m,c) = properFraction $ bx / luS-    R luS = log_USeed @es+    k = floor (x / ln2 + 0.5)  -- should be Integer or Int+    r = x - fromInteger k * ln2  -- should be Rational +ln2 :: Rational+ln2 = 0.6931471805599453094172321214581765680755001343602552541206800094933936219696947156058633269964186875420014810205706857336855202  log_USeed :: forall es. PositC es => Posit es log_USeed = approx_log $ fromIntegral (uSeed @es)@@ -870,7 +873,6 @@   | otherwise = tuma_approx_cos $ 2*approx_pi * x --  - -- Using the CORDIC domain reduction and some approximation function of log funLogDomainReduction :: forall es. PositC es => (Posit es -> Posit es) -> Posit es -> Posit es funLogDomainReduction _ NaR = NaR@@ -881,9 +883,9 @@     where       (ex, sig) = (int * fromIntegral (2^(exponentSize @es)) + fromIntegral nat + 1, fromRational rat / 2) -- move significand range from 1,2 to 0.5,1       (_,int,nat,rat) = (posit2TupPosit @es).toRational $ x -- sign should always be positive-     -  ++ -- natural log with log phi acurate to 9 ULP funLogTaylor :: forall es. PositC es => Posit es -> Posit es funLogTaylor NaR = NaR@@ -925,9 +927,25 @@ taylor_approx_expm1 x = go 0 [x^n / fromIntegral (fac n) | n <- [1..]]   where     go :: Posit es -> [Posit es] -> Posit es+    go !acc [] = acc     go !acc (h:t) | acc == acc + h = acc                   | otherwise = go (acc + h) t  -- need to use a Taylor Series, the `tanh` formulation doesn't work because it requires something that depends on `exp`++--+-- calculate exp, its most accurate near zero+-- sum k=0 to k=inf of the terms, iterate until a fixed point is reached+funExpTaylor :: forall es. PositC es => Posit es -> Posit es+funExpTaylor NaR = NaR+funExpTaylor 0 = 1+funExpTaylor z = go 0 0+  where+    go :: Natural -> Posit es -> Posit es+    go !k !acc+      | acc == (acc + term k) = acc  -- if x == x + dx then terminate and return x+      | otherwise = go (k+1) (acc + term k)+    term :: Natural -> Posit es+    term k = (z^k) / (fromIntegral.fac $ k)   -- =====================================================================
src/Posit/Internal/PositC.hs view
@@ -25,6 +25,7 @@ {-# LANGUAGE DerivingVia #-}  -- To Derive instances for newtypes to eliminate Orphan Instances {-# LANGUAGE UndecidableInstances #-}  -- For deriving DoubleWord {-# LANGUAGE CPP #-} -- To remove Storable instances to remove noise when performing analysis of Core+{-# LANGUAGE MultiWayIf #-} {-# OPTIONS_GHC -Wno-unticked-promoted-constructors #-}  -- Turn off noise {-# OPTIONS_GHC -Wno-type-defaults #-}  -- Turn off noise @@ -32,7 +33,6 @@ --  |Posit Class, implementing: -- --   * PositC---   * Orphan Instances of Storable for Word128, Int128, Int256 -- ----  module Posit.Internal.PositC@@ -61,7 +61,7 @@ import Data.Int (Int8,Int16,Int32,Int64)  -- Import standard Int sizes import Data.DoubleWord (Word128,Int128,Int256,fromHiAndLo,hiWord,loWord,DoubleWord,BinaryWord) -- Import large Int sizes import Data.Word (Word64)-import Data.Bits (Bits(..), shiftL, shift, testBit, (.&.), shiftR,FiniteBits)+import Data.Bits (Bits(..), shiftL, shiftR, testBit, (.&.), FiniteBits)   -- Import Naturals and Rationals@@ -291,7 +291,7 @@   formRegime :: Integer -> (IntN es, Integer)   formRegime power     | 0 <= power =-      let offset = (fromIntegral (nBits @es - 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 @es - 1) - abs power - 1)@@ -300,16 +300,23 @@   formExponent :: Natural -> Integer -> (IntN es, Integer)   formExponent power offset =     let offset' = offset - fromIntegral (exponentSize @es)-    in (fromIntegral power `shift` fromInteger offset', offset')+        result = fromIntegral power `shiftR` negate (fromInteger offset')+    in if offset' >= 0+       then (fromIntegral power `shiftL` fromInteger offset', offset')+       else if testBit (fromIntegral power :: Int) (pred.fromIntegral.negate $ offset')  -- What is "Twilight Zone" Posit Number?+            then (succ result, offset')+            else (result, offset')  -- rounding case      formFraction :: Rational -> Integer -> IntN es   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+    in if | numFractionBits >= 1 -> fromInteger normFraction+          | numFractionBits == 0 -> if 2 * (r - 1) > 1 -- rounding case -- Previously was 0, this should be "banker's rounding" per Gustafson's ch3+                                    then 1+                                    else 0+          | otherwise -> 0      tupPosit2Posit :: (Bool,Integer,Natural,Rational) -> Maybe Rational   tupPosit2Posit (sgn,regime,exponent,rat) = -- s = isNeg posit == True@@ -379,9 +386,17 @@                else "0" ++ go (idx - 1)      -- decimal Precision-  decimalPrec :: Int-  decimalPrec = fromIntegral $ 2 * (nBytes @es) + 1-  +  decimalPrec :: IntN es -> Int+  -- decimalPrec _ = fromIntegral $ 2 * (nBytes @es) + 1  -- The Olden way, 1.37% fail to round trip+  --+  -- `1 + ` Gustafson's Decimal Accuracy (from coorispondance with Alessandro)  -- Sucessfuly round triped for both Posit32 and P32, and below+  decimalPrec posit =+    let regimeFormat = findRegimeFormat @es posit+        regimeCount = countRegimeBits @es regimeFormat posit+        fractionSize = max (fromIntegral (nBits @es) - fromIntegral (signBitSize @es) - regimeCount - fromIntegral (exponentSize @es)) 0  -- fractionSize is at least zero+    in succ.ceiling $ (fromIntegral fractionSize + 1) * log10Of2 + halflog10Of2+  --+   {-# MINIMAL exponentSize | nBytes #-}  @@ -462,6 +477,14 @@  xnor :: Bool -> Bool -> Bool xnor a b = not $ (a || b) && not (b && a)+++log10Of2 :: Double+log10Of2 = 0.3010299956639811952137388947+++halflog10Of2 :: Double+halflog10Of2 = 0.1505149978319905976068694473   #ifndef O_NO_STORABLE_RANDOM
stack.yaml view
@@ -2,7 +2,7 @@ # that coorispond to a specific GHC or Stackage version  # resolver: nightly-2023-07-06  # ghc-9.6.2-resolver: lts-21.1 # ghc-9.4.5+resolver: lts-21.9 # ghc-9.4.6 # resolver: lts-20.26 # ghc-9.2.8 # resolver: lts-19.33 # ghc-9.0.2 # resolver: lts-18.28 # ghc-8.10.7
test/Test/Algorithms.hs view
@@ -126,11 +126,13 @@       z = 1 / 2^n  -----+{-+approx_exp :: forall es. PositC es => Posit es -> Posit es+approx_exp (R bx) = fromIntegral (uSeed @es)^^m * (taylor_approx_expm1 (R c * log_USeed) + 1)+  where+    (m,c) = properFraction $ bx / luS+    R luS = log_USeed @es+-}   @@ -301,8 +303,8 @@     go !prevA !prevS !n !a !s       | prevA == a = a       | prevS == s = a-      | abs (prevA - a) <= 2*eps = a  -- P256 or Posit128, will not reach a fixed point where `prevA == a` it sort of oscelates until divergence occurs, if we test for less than 2*eps it can stop early-      | abs (prevS - s) <= 2*eps = a+      | abs (prevA - a) <= 3*eps = a  -- P256 or Posit128, will not reach a fixed point where `prevA == a` it sort of oscelates until divergence occurs, if we test for less than 3*eps it can stop early+      | abs (prevS - s) <= 3*eps = a       | otherwise =         let x = 5 / s - 1             y = (x - 1)^2 + 7
test/TestElementaryFunctions.hs view
@@ -32,91 +32,91 @@ expPlotP16Posit16 = toFile def "./test/Results/Bits Accuracy of exp with P16 and Posit16.png" $ do     let expP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (exp p) exp) | p <- enumFrom (NaR :: P16)]         expPosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (exp p) exp) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of  P16 and Posit16 exp"+    layout_title .= "Number of Accurate Bits of  P16 and Posit16 'exp'"     plot (line "exp P16 error" [expP16])     plot (line "exp Posit16 error" [expPosit16])  logPlotP16Posit16 = toFile def "./test/Results/Bits Accuracy of log with P16 and Posit16.png" $ do     let lnP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (log p) log) | p <- enumFrom (0 :: P16)]         lnPosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (log p) log) | p <- enumFrom (0 :: Posit16)]-    layout_title .= "Number of Accurate Bits of  P16 and Posit16 Log"+    layout_title .= "Number of Accurate Bits of  P16 and Posit16 'log'"     plot (line "log P16 error" [lnP16])     plot (line "log Posit16 error" [lnPosit16])  sqrtPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of sqrt with P16 and Posit16.png" $ do     let sqrtP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (sqrt p) sqrt) | p <- enumFrom (0 :: P16)]         sqrtPosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (sqrt p) sqrt) | p <- enumFrom (0 :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Square Root"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'sqrt'"     plot (line "sqrt P16 error" [sqrtP16])     plot (line "sqrt Posit16 error" [sqrtPosit16])  sinPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of sin with P16 and Posit16.png" $ do     let sineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (sin p) sin) | p <- enumFrom (NaR :: P16)]         sinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (sin p) sin) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Sine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'sin'"     plot (line "sin P16 error" [sineP16])     plot (line "sin Posit16 error" [sinePosit16])  cosPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of cos with P16 and Posit16.png" $ do     let cosineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (cos p) cos) | p <- enumFrom (NaR :: P16)]         cosinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (cos p) cos) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Cosine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'cos'"     plot (line "cos P16 error" [cosineP16])     plot (line "cos Posit16 error" [cosinePosit16])  asinPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of asin with P16 and Posit16.png" $ do     let arcsineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (asin p) sin) | p <- enumFromTo (-1 :: P16) 1]         arcsinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (asin p) sin) | p <- enumFromTo (-1 :: Posit16) 1]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 arcSine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'asin'"     plot (line "asin P16 error" [arcsineP16])     plot (line "asin Posit16 error" [arcsinePosit16])  acosPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of acos with P16 and Posit16.png" $ do     let arccosineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (acos p) acos) | p <- enumFromTo (-1 :: P16) 1]         arccosinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (acos p) acos) | p <- enumFromTo (-1 :: Posit16) 1]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 arcCosine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'acos'"     plot (line "acos P16 error" [arccosineP16])     plot (line "acos Posit16 error" [arccosinePosit16])  atanPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of atan with P16 and Posit16.png" $ do     let arctangentP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (atan p) atan) | p <- enumFrom (NaR :: P16)]         arctangentPosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (atan p) atan) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 arcTangent"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'atan'"     plot (line "atan P16 error" [arctangentP16])     plot (line "atan Posit16 error" [arctangentPosit16])  sinhPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of sinh with P16 and Posit16.png" $ do     let hypsineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (sinh p) sinh) | p <- enumFrom (NaR :: P16)]         hypsinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (sinh p) sinh) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Hyperbolic Sine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'sinh'"     plot (line "sinh P16 error" [hypsineP16])     plot (line "sinh Posit16 error" [hypsinePosit16])  coshPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of cosh with P16 and Posit16.png" $ do     let hypcosineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (cosh p) cosh) | p <- enumFrom (NaR :: P16)]         hypcosinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (cosh p) cosh) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Hyperbolic Cosine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'cosh'"     plot (line "cosh P16 error" [hypcosineP16])     plot (line "cosh Posit16 error" [hypcosinePosit16])  asinhPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of asinh with P16 and Posit16.png" $ do     let archypsineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (asinh p) asinh) | p <- enumFrom (NaR :: P16)]         archypsinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (asinh p) asinh) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Hyperbolic Sine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'asinh'"     plot (line "asinh P16 error" [archypsineP16])     plot (line "asinh Posit16 error" [archypsinePosit16])  acoshPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of acosh with P16 and Posit16.png" $ do     let archypcosineP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (acosh p) acosh) | p <- enumFrom (NaR :: P16)]         archypcosinePosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (acosh p) acosh) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Inv Hyperbolic Cosine"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'acosh'"     plot (line "acosh P16 error" [archypcosineP16])     plot (line "acosh Posit16 error" [archypcosinePosit16])  atanhPlotsP16Posit16 = toFile def "./test/Results/Bits Accuracy of atanh with P16 and Posit16.png" $ do     let archyptangentP16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (atanh p) atanh) | p <- enumFrom (NaR :: P16)]         archyptangentPosit16 = filter (\(_,d) -> not $ nanOrInf d) [(read (displayIntegral p) :: Double, err p (atanh p) atanh) | p <- enumFrom (NaR :: Posit16)]-    layout_title .= "Number of Accurate Bits of P16 and Posit16 Inv Hyperbolic Tangent"+    layout_title .= "Number of Accurate Bits of P16 and Posit16 'atanh'"     plot (line "atanh P16 error" [archyptangentP16])     plot (line "atanh Posit16 error" [archyptangentPosit16])
test/TestPosit.hs view
@@ -75,13 +75,21 @@   print $ "Machine epsilon P64 ~1.0: " ++ show (eps :: P64) --    print $ "Machine epsilon P128 ~1.0: " ++ show (eps :: P128) --    print $ "Machine epsilon P256 ~1.0: " ++ show (eps :: P256) -- -  -- | Taylor vs. Tuma   print $ "Does (1 - 1) == 0 ?: " ++ show ((1 - 1) == (0 :: Posit256)) -- [(1 - 1) == zero | zero = 0 :: Posit es, es <- Z .. V]-  let sqrtTaylor = (funLogDomainReduction funLogTaylor).(/2).(funExp2 funExpTaylor).(/log 2)-  print $ "sqrt phi using a Taylor algorithm: " ++ show (sqrtTaylor phi)-  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))+  print $ "truth    :: Posit256: " ++ show (1.2720196495140689642524224617374914917156080418400962486166403825392975755360680118303842149884602585385141476367280265057103381 :: Posit256)+  print $ "sqrt phi :: Posit256: " ++ show (sqrt phi :: Posit256)+  print $ "sqrt phi :: Posit128: " ++ show (sqrt phi :: Posit128)+  print $ "sqrt phi :: Posit64: " ++ show (sqrt phi :: Posit64)+  print $ "sqrt phi :: Posit32: " ++ show (sqrt phi :: Posit32)+  print $ "sqrt phi :: Posit16: " ++ show (sqrt phi :: Posit16)+  print $ "sqrt phi :: Posit8: " ++ show (sqrt phi :: Posit8)+  print $ "truth    :: P256: " ++ show (1.2720196495140689642524224617374914917156080418400962486166403825392975755360680118303842149884602585385141476367280265057103381 :: P256)+  print $ "sqrt phi :: P256: " ++ show (sqrt phi :: P256)+  print $ "sqrt phi :: P128: " ++ show (sqrt phi :: P128)+  print $ "sqrt phi :: P64: " ++ show (sqrt phi :: P64)+  print $ "sqrt phi :: P32: " ++ show (sqrt phi :: P32)+  print $ "sqrt phi :: P16: " ++ show (sqrt phi :: P16)+  print $ "sqrt phi :: P8: " ++ show (sqrt phi :: P8)   {-   let truthPosit256 = 0.8956731517052878608869612167009786079379812529831641161347143256  :: Posit256  -- 0.89566032673209158354178209470474131001971567786620187475744721557  :: Posit256   -- 0.8956731517052878608869612167009786079379812529831641161347143256836782657295966290940929214799036260987761959338755143914935872 :: Posit256   let truthP256 = 0.8956731517052878608869612167009786079379812529831641161347143256 :: P256 --  0.89566032673209158354178209470474131001971567786620187475744721557 :: P256    -- 0.8956731517052878608869612167009786079379812529831641161347143256836782657295966290940929214799036260987761959338755143914935872 :: P256
+ test/TestReadShow.hs view
@@ -0,0 +1,178 @@++--------------------------------------------------------------------------------------------+-- | Posit Numbers+--   Copyright   :  (C) 2022-2023 Nathan Waivio+--   License     :  BSD3+--   Maintainer  :  Nathan Waivio <nathan.waivio@gmail.com>+--   Stability   :  Stable+--   Portability :  Portable+--+--   Test Suite for a Library implementing standard Posit Numbers+-- +---------------------------------------------------------------------------------------------++-- Hmm... phi is the most irrational number?++{-# LANGUAGE ScopedTypeVariables #-} --   To reduce some code duplication, this is important+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeApplications #-}++import Posit+import Posit.Internal.PositC++import Data.List (partition)++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++main :: IO ()+main = do+  let (pPosit8, fPosit8) = partition testCase [ x | x <- enumFrom (NaR :: Posit8)]+      lpPosit8 = length pPosit8 :: Int -- Int+      lfPosit8 = length fPosit8 :: Int -- Int+  print $ "Fraction Error (read.show) :: Posit8 : " ++ show (lfPosit8 % (lpPosit8 + lfPosit8))+  print $ "Failure case: " ++ show fPosit8+  let (pP8, fP8) = partition testCase [ x | x <- enumFrom (NaR :: P8)]+      lpP8 = length pP8 :: Int+      lfP8 = length fP8 :: Int+  print $ "Fraction Error (read.show) :: P8 : " ++ show (lfP8 % (lpP8 + lfP8))+  print $ "Failure case: " ++ show fP8+  let (pPosit16, fPosit16) = partition testCase [ x | x <- enumFrom (NaR :: Posit16)]+      lpPosit16 = length pPosit16 :: Int -- Int+      lfPosit16 = length fPosit16 :: Int -- Int+  print $ "Fraction Error (read.show) :: Posit16 : " ++ show (lfPosit16 % (lpPosit16 + lfPosit16))+  print $ "Failure case: " ++ show fPosit16+  let (pP16, fP16) = partition testCase [ x | x <- enumFrom (NaR :: P16)]+      lpP16 = length pP16 :: Int+      lfP16 = length fP16 :: Int+  print $ "Fraction Error (read.show) :: P16 : " ++ show (lfP16 % (lpP16 + lfP16))+  print $ "Failure case: " ++ show fP16+  -- Larger size can only test percentage+  --+  let lfPosit32 = length [ x | x <- enumFrom (NaR :: Posit32), not (testCase x)]+  print $ "Fraction Error (read.show) :: Posit32 : " ++ show (lfPosit32 % 2^32)+  let lfP32 = length [ x | x <- enumFrom (NaR :: P32), not (testCase x)]+  print $ "Fraction Error (read.show) :: P32 : " ++ show (lfP32 % 2^32)+  --+++testCase :: forall es. PositC es => Posit es -> Bool+testCase x = x == ((read @(Posit es)).show $ x)+++-- * Gustafson's coorispondance with Alessandro+-- +-- Decimal accuracy is not the same as the minimal number of decimals+-- needed to assure conversion back to the original posit. The latter is+-- larger. For example, a 32-bit posit has 28 bits of significance in+-- its highest-accuracy region, which is about:+-- > 28 * log10(2) + 0.5 log10(2) ≈ 8.6 +-- decimals of minimum relative accuracy, but if you only use 9 decimals+-- to express a posit and you convert those 9 decimals back into a+-- posit, occasionally you will get a different posit from what you+-- started with.++-- * And then responding again to Alessandro++-- First of all, there are two regime bits (01 or 10), not one. The+-- termination (opposite) bit is part of the regime. That leaves 27 bits+-- for the fraction, which is 28 bits of significand including the+-- hidden value left of the radix point.++-- Yes, your formula should then work, but remember that you have+-- wobbling relative accuracy (wobble of ±½ log10(2) ≈ ±0.15 decimals)+-- for both posits and floats within any given binade. You should+-- probably be pessimistic and use the lower bound on the sawtooth shape+-- if you are trying to prove something about the accuracy of a+-- computation.+++-- * During test of "round tripping" a Decimal text string, a rounding+-- issue was found.++-- After resolving a rounding error that added many failure cases when +-- there should have been rounding for when `2 * remainingSignificand - 1 > 1`+-- for Posits with no fraction.+-- >  ...+-- >      else if 2 * (r - 1) > 1 -- Previously was 0, this should be "banker's rounding" as per Gustafson (upcoming Ch3), where if `2*(r - 1) > 1` then 1 else 0+-- >           then 1+-- >           else 0+++-- The baseline results using `decimalPrec = fromIntegral $ 2 * (nBytes @es) + 1` -- [3,5,9,17,33,65]+-- Benchmark test-posit-readShowId: RUNNING...+-- "Fraction Error (read.show) :: Posit8 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: P8 : 1 % 128"+-- "Failure case: [-2.441e-4,2.441e-4]"+-- "Fraction Error (read.show) :: Posit16 : 1 % 8192"+-- "Failure case: [-2.68435e8,-3.35544e7,-5.96046e-8,-2.98023e-8,2.98023e-8,5.96046e-8,3.35544e7,2.68435e8]"+-- "Fraction Error (read.show) :: P16 : 3 % 16384"+-- "Failure case: [-1.40737e14,-7.03687e13,-5.68434e-14,-2.84217e-14,-3.55271e-15,-8.88178e-16,8.88178e-16,3.55271e-15,2.84217e-14,5.68434e-14,7.03687e13,1.40737e14]"+-- "Fraction Error (read.show) :: Posit32 : 37904825 % 2147483648"+-- "Fraction Error (read.show) :: P32 : 37904825 % 2147483648"+-- Benchmark test-posit-readShowId: FINISH+--+-- After bug fix:+-- Benchmark test-posit-readShowId: RUNNING...+-- "Fraction Error (read.show) :: Posit8 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: P8 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: Posit16 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: P16 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: Posit32 : 37904819 % 2147483648"+-- "Fraction Error (read.show) :: P32 : 37904819 % 2147483648"+-- Benchmark test-posit-readShowId: FINISH+++-- Gustafson's Decimal Accuracy: -- []+-- >   decimalPrec posit =+-- >    let regimeFormat = findRegimeFormat @es posit+-- >        regimeCount = countRegimeBits @es regimeFormat posit+-- >        fractionSize = max (fromIntegral (nBits @es) - fromIntegral (signBitSize @es) - regimeCount - fromIntegral (exponentSize @es)) 0  -- fractionSize is at least zero+-- >    in ceiling $ (fromIntegral fractionSize + 1) * log10Of2 + halflog10Of2+-- Benchmark test-posit-readShowId: RUNNING...+-- "Fraction Error (read.show) :: Posit8 : 1 % 128"+-- "Failure case: [-3.1e-2,3.12e-2]"+-- "Fraction Error (read.show) :: P8 : 7 % 128"+-- "Failure case: [-0.12,-2.4e-4,-1.2e-4,-6.1e-5,-1.5e-5,-3.8e-6,-9.5e-7,9.5e-7,3.8e-6,1.5e-5,6.1e-5,1.2e-4,2.44e-4,0.12]"+-- "Fraction Error (read.show) :: Posit16 : 3 % 65536"+-- "Failure case: [-6.7e7,5.96e-8,6.7e7]"+-- "Fraction Error (read.show) :: P16 : 23 % 65536"+-- "Failure case: [-7.2e16,-4.5e15,-1.1e15,-2.8e14,-1.4e14,-7.0e13,-3.5e13,-2.8e-14,-1.4e-14,-7.1e-15,-2.2e-16,2.2e-16,7.1e-15,1.4e-14,2.8e-14,5.68e-14,3.5e13,7.0e13,1.4e14,2.8e14,1.1e15,4.5e15,7.2e16]"+-- "Fraction Error (read.show) :: Posit32 : 75809653 % 4294967296"+-- "Fraction Error (read.show) :: P32 : 75809653 % 4294967296"+-- Benchmark test-posit-readShowId: FINISH+--+-- After bug fix:+-- Benchmark test-posit-readShowId: RUNNING...+-- "Fraction Error (read.show) :: Posit8 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: P8 : 1 % 128"+-- "Failure case: [-0.12,0.12]"+-- "Fraction Error (read.show) :: Posit16 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: P16 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: Posit32 : 37904819 % 2147483648"+-- "Fraction Error (read.show) :: P32 : 37904819 % 2147483648"+-- Benchmark test-posit-readShowId: FINISH+--+-- After bug fix and incrementing by 1 decimal digit+-- Benchmark test-posit-readShowId: RUNNING...+-- "Fraction Error (read.show) :: Posit8 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: P8 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: Posit16 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: P16 : 0 % 1"+-- "Failure case: []" +-- "Fraction Error (read.show) :: Posit32 : 0 % 1"+-- "Fraction Error (read.show) :: P32 : 0 % 1"+-- Benchmark test-posit-readShowId: FINISH+++