hafar 0.1.0.0 → 0.1.1.0

raw patch · 4 files changed

+163/−18 lines, 4 filesdep ~basePVP: major bump suggested

API removals or changes: PVP suggests a major version bump

Dependency ranges changed: base

API changes (from Hackage documentation)

Files

• hafar.cabal +6/−4
• src/Numeric/AffineForm.hs +2/−1
• src/Numeric/AffineForm/Internal.hs +30/−13
• src/Numeric/AffineForm/Subdivision.hs +125/−0

@@ -4,10 +4,10 @@ -- -- see: https://github.com/sol/hpack ----- hash: 11f7d969e8f00b5e392eda2b0897debe20363e72178cbd6a8a87a7ab29e2b819+-- hash: 6bf5e0c8d0a1cae72591bdd92749728d178e27c39856dbdc12dc983ea2f322f6  name:           hafar-version:        0.1.0.0+version:        0.1.1.0 synopsis:       Affine arithmetic library for Haskell description:    Hafar is an affine arithmetic library for Haskell. It is an efficient way to work with ranges of values or imprecise values. category:       Numeric@@ -31,13 +31,14 @@   exposed-modules:       Numeric.AffineForm       Numeric.AffineForm.ExplicitRounding+      Numeric.AffineForm.Subdivision   other-modules:       Numeric.AffineForm.Utils       Numeric.AffineForm.Internal   hs-source-dirs:       src   build-depends:-      base >=4.12 && <4.14+      base >=4.12 && <4.13     , intervals >=0.8 && <0.9     , mtl >=2.2 && <2.3   default-language: Haskell2010@@ -50,13 +51,14 @@       Numeric.AffineForm.Internal       Numeric.AffineForm.Utils       Numeric.AffineForm.ExplicitRounding+      Numeric.AffineForm.Subdivision   hs-source-dirs:       test       src   ghc-options: -threaded -rtsopts -with-rtsopts=-N   build-depends:       QuickCheck >=2.13 && <2.14-    , base >=4.12 && <4.14+    , base >=4.12 && <4.13     , hafar     , intervals >=0.8 && <0.9     , mtl >=2.2 && <2.3

@@ -6,12 +6,13 @@                            midpoint,                            inf, sup,                            interval,+                           afError,                            member,                            epscount_,                            setMidpoint,                            fix,                            addError,-                           (.+), (.*)+                           (.+), (.*),                           ) where  import Numeric.AffineForm.Internal

@@ -15,6 +15,7 @@  import Numeric.AffineForm.Utils import Numeric.AffineForm.ExplicitRounding+import Numeric.AffineForm.Subdivision import qualified Numeric.Interval as IA import Numeric.Interval ((...)) import Data.Fixed (mod')@@ -69,32 +70,35 @@   acosh = minrange acosh (\x -> 1/((sqrt (x-1))*(sqrt (x+1)))) Concave   atanh = minrange atanh (\x -> 1/(1-x^2)) undefined +instance (Ord a, Fractional a, ExplicitRounding a) => Subdivisible (AF s a) where+  subdivide af n = mapToInterval af <$> subdivide (interval af) n+  combine2 l r   = fromInterval $ combine2 (interval l) (interval r)+ type AFIndex = Int  -- | AFM is a state monad that ensures that any new noise symbols have not been used by any previous affine form. -- All affine arithmetic calculations should be done inside the AFM monad. Affine forms do not make sense outside of their monad context.-newtype AFMT t s m a = AFMT {runAFMT :: s -> m (a, s)}-type AFM t a = AFMT t AFIndex Identity a+newtype AFMT t m a = AFMT {runAFMT :: AFIndex -> m (a, AFIndex)}+type AFM t a = AFMT t Identity a -instance (Monad m) => Functor (AFMT t s m) where+instance (Monad m) => Functor (AFMT t m) where   fmap = liftM -instance (Monad m) => Applicative (AFMT t s m) where+instance (Monad m) => Applicative (AFMT t m) where   pure = return   (<*>) = ap -instance (Monad m) => Monad (AFMT t s m) where+instance (Monad m) => Monad (AFMT t m) where   return a = AFMT $ \s -> return (a, s)   (AFMT x) >>= f = AFMT $ \s -> do     (v, s') <- x s     (AFMT x') <- return $ f v     x' s' -instance (Monad m) => MonadState s (AFMT t s m) where-  get   = AFMT $ \s -> return (s, s)-  put s = AFMT $ \_ -> return ((), s)+getIndex   = AFMT $ \s -> return (s, s)+putIndex s = AFMT $ \_ -> return ((), s) -instance MonadTrans (AFMT t s) where+instance MonadTrans (AFMT t) where   lift c = AFMT $ \s -> c >>= (\x -> return (x, s))  -- | This gives an affine form with midpoint 0 and radius 1.@@ -102,8 +106,8 @@ -- It can be used to instantiate new affine forms. newEps :: Num a => AFM t (AF t a) newEps = do-  idx <- get-  put $ idx + 1+  idx <- getIndex+  putIndex $ idx + 1   return $ AF 0 (replicate idx 0 ++ [1]) 0  -- | Creates a new affine form that covers the interval.@@ -119,12 +123,12 @@ singleton x = AF x [] 0  -- | Creates a new affine form that approximately represents some value.--- This function adds a small error to account for the 'wobble' in the computer representation of the value.+-- This function adds a small error to account for the `wobble` in the computer representation of the value. approxSingleton :: (ExplicitRounding a) => a -> AF s a approxSingleton x = AF x [] $ eps x  -- | Evaluates the AFM monad. It is not possible to get an AF out of an AFM monad.-evalAFM :: forall a b. (forall t. AFM t b) -> b+evalAFM :: (forall t. AFM t b) -> b evalAFM (AFMT x) = fst . runIdentity $ x 0  -- | Gives the radius of the affine form@@ -158,6 +162,10 @@ epscount_ :: AF s a -> Int epscount_ (AF _ xs _) = length xs +-- | Returns the value of the error term+afError :: AF s a -> a+afError (AF _ _ e) = e+ -- Affine arithmetic operations  -- | Sets the midpoint of the affine form@@ -169,6 +177,15 @@ addError (AF x xs xe) e   | e >= 0 = AF x xs (xe + e)   | otherwise = throw AddingNegativeError++mapToInterval :: (Fractional a, ExplicitRounding a) => AF s a -> IA.Interval a -> AF s a+mapToInterval af i = m .+ ((wi/waf) .* af)+  where waf = radius af+        wi  = IA.width i+        m   = IA.midpoint i++fromInterval :: (Fractional a) => IA.Interval a -> AF s a+fromInterval i = AF (IA.midpoint i) [] ((IA.width i)/2)  -- | Adds a scalar value to the affine form (.+) :: (Num a, ExplicitRounding a) => a -> AF s a -> AF s a

@@ -0,0 +1,125 @@+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE DeriveFunctor #-}++-- | Subdivision contains an implementation of the branch-and-bound algorithm.+-- This method divides interval-like values into smaller subdivisions and then applies a function to those values.+-- This results in a more accurate output at the cost of having to repeat the calculations multiple times+module Numeric.AffineForm.Subdivision+  ( Subdivisible(..)+  , SubdivisionEnvironment(..)+  , defaultEnvironment+  , branchAndBound+  ) where++import Control.Monad.Reader+import Control.Monad.State++import Data.List++import Numeric.Interval as IA++-- | The 'Subdivisible' class is used for datatypes that can be broken down into smaller pieces+-- and combined together+--+-- The 'subdivide' function subdivides the value into `n` smaller values and 'combine' joins+-- the subdivisions together into a bigger value+--+-- Subdividing and then combining a value is expected to give the initial value (sans rounding errors)+class Subdivisible a where+  subdivide    :: a -> Int -> [a]+  combine2     :: a -> a -> a+  combine      :: [a] -> a++  combine2 l r = combine [l,r]+  combine      = foldl1 combine2++instance (Ord a, Fractional a) => Subdivisible (Interval a) where+  subdivide i n = f . fromIntegral <$> [1..n]+    where f x = (((x-1)*w/n')IA....(x*w/n')) + (IA.singleton $ inf i)+          w   = width i+          n'  = fromIntegral n+  combine2 = hull++instance (Subdivisible a) => Subdivisible [a] where+  subdivide l n = sequenceA $ flip subdivide n <$> l+  combine2 l r  = uncurry combine2 <$> zip l r++-- | A data structure for configuring the 'branchAndBound' method+--+-- 'function' is the function that gets evaluated+-- 'errorFun' is the error measuring function of the result+-- 'maxError' specifies the maximum permitted error of an evaluation+-- 'maxDepth' specifies how deep the subdivision can go+-- 'subdivs'  specifies the number of subdivisions per value+data SubdivisionEnvironment a b e = SubdivisionEnvironment+  { function :: a -> b+  , errorFun :: b -> e+  , maxError :: e+  , maxDepth :: Int+  , subdivs  :: Int+  }++-- | This function creates a simple subdivision configuration.+-- It requires the evaluator function and an error measuring function as its parameters.+defaultEnvironment :: (Fractional e) => (a -> b) -> (b -> e) -> SubdivisionEnvironment a b e+defaultEnvironment f g = SubdivisionEnvironment+  { function = f+  , errorFun = g+  , maxError = 0.1+  , maxDepth = 3+  , subdivs  = 2+  }++type Subdivider a b e = Reader (SubdivisionEnvironment a b e) (SubdivisionTree a)++data SubdivisionTree a+  = Branch [SubdivisionTree a]+  | Node a+  deriving (Show, Functor)++deepen :: (Subdivisible a, Ord e) => SubdivisionTree a -> Int -> Subdivider a b e+deepen (Node x) depth =+  do+    env <- ask+    let res = function env $ x+        err = errorFun env $ res+        n   = subdivs env+    if err <= maxError env+       -- Accurate answer found+       then return $ Node x+       -- Subdivide and deepen+       else deepen (Branch $ Node <$> subdivide x n) $ depth+deepen (Branch l) depth =+  do+    env <- ask+    if depth >= maxDepth env+       then return $ Branch l+       else Branch <$> (sequence $ flip deepen (depth + 1) <$> l)++collapse :: (Subdivisible a) => SubdivisionTree a -> a+collapse (Node a)   = a+collapse (Branch l) = combine $ collapse <$> l++evalTree :: (Subdivisible b) => Subdivider a b e -> SubdivisionEnvironment a b e -> b+evalTree div cfg = flip runReader cfg $+  do+    val <- div+    env <- ask+    return . collapse $ function env <$> val++-- | This function iteratively subdivides a value to arrive at a more accurate result+branchAndBound :: (Subdivisible a, Subdivisible b, Ord e) => a -> SubdivisionEnvironment a b e -> b+branchAndBound x env = flip evalTree env $ deepen (Node x) 0++test :: Interval Double+test =+  flip evalTree env $ deepen (Node $ [3 IA.... 4, 1 IA....4]) 0+  where env = SubdivisionEnvironment+          { function = (\[x, y] -> x * y - x)+          , errorFun = width+          , maxError = 0.1+          , maxDepth = 5+          , subdivs  = 2+          }+