ad 0.40.1 → 0.44.0

raw patch · 13 files changed

+446/−265 lines, 13 filesdep −mlistPVP ok

version bump matches the API change (PVP)

Dependencies removed: mlist

API changes (from Hackage documentation)

Files

• Numeric/AD/Internal/Identity.hs +1/−0
• Numeric/AD/Internal/Reverse.hs +34/−0
• Numeric/AD/Internal/Sparse.hs +89/−4
• Numeric/AD/Internal/Stream.hs +1/−0
• Numeric/AD/Internal/Tensors.hs +0/−1
• Numeric/AD/Mode/Forward.hs +0/−14
• Numeric/AD/Mode/Mixed.hs +86/−84
• Numeric/AD/Mode/Reverse.hs +26/−84
• Numeric/AD/Mode/Sparse.hs +123/−0
• Numeric/AD/Mode/Tower.hs +0/−12
• Numeric/AD/Newton.hs +2/−60
• Numeric/AD/Types.hs +26/−0
• ad.cabal +58/−6

@@ -1,4 +1,5 @@ {-# LANGUAGE GeneralizedNewtypeDeriving, FlexibleContexts, FlexibleInstances, MultiParamTypeClasses, DeriveDataTypeable #-}+{-# OPTIONS_HADDOCK hide #-} ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.AD.Internal.Identity

@@ -32,6 +32,8 @@     , unbindMap     , unbindWith     , unbindMapWithDefault+    , vgrad, vgrad'+    , Grad(..)     ) where  import Prelude hiding (mapM)@@ -215,3 +217,35 @@  unbindMapWithDefault :: (Functor f, Var v, Num a) => b -> (a -> b -> c) -> f (v a) -> IntMap b -> f c unbindMapWithDefault z f xs ys = fmap (\v -> f (primal v) $ findWithDefault z (varId v) ys) xs++class Num a => Grad i o o' a | i -> a o o', o -> a i o', o' -> a i o where+    pack :: i -> [AD Reverse a] -> AD Reverse a+    unpack :: ([a] -> [a]) -> o+    unpack' :: ([a] -> (a, [a])) -> o'++instance Num a => Grad (AD Reverse a) [a] (a, [a]) a where+    pack i _ = i+    unpack f = f []+    unpack' f = f []++instance Grad i o o' a => Grad (AD Reverse a -> i) (a -> o) (a -> o') a where+    pack f (a:as) = pack (f a) as+    pack _ [] = error "Grad.pack: logic error"+    unpack f a = unpack (f . (a:))+    unpack' f a = unpack' (f . (a:))++vgrad :: Grad i o o' a => i -> o+vgrad i = unpack (unsafeGrad (pack i))+    where+        unsafeGrad f as = unbind vs (partialArray bds $ f vs)+            where+                (vs,bds) = bind as++vgrad' :: Grad i o o' a => i -> o'+vgrad' i = unpack' (unsafeGrad' (pack i))+    where+        unsafeGrad' f as = (primal r, unbind vs (partialArray bds r))+            where+                r = f vs+                (vs,bds) = bind as+

@@ -1,17 +1,21 @@-{-# LANGUAGE BangPatterns, TemplateHaskell, TypeFamilies, TypeOperators, FlexibleContexts, UndecidableInstances, DeriveDataTypeable #-}+{-# LANGUAGE BangPatterns, TemplateHaskell, TypeFamilies, TypeOperators, FlexibleContexts, UndecidableInstances, DeriveDataTypeable, MultiParamTypeClasses, FunctionalDependencies, FlexibleInstances #-} module Numeric.AD.Internal.Sparse      ( Index(..)     , emptyIndex     , addToIndex     , indices     , Sparse(..)+    , apply     , vars-    , d-    , d'-    , ds+    , d, d', ds     , skeleton     , spartial     , partial+    , vgrad+    , vgrad'+    , vgrads+    , Grad(..)+    , Grads(..)     ) where  import Prelude hiding (lookup)@@ -65,6 +69,10 @@         var !n a = (n + 1, AD $ Sparse a $ singleton n $ lift 1) {-# INLINE vars #-} +apply :: (Traversable f, Num a) => (f (AD Sparse a) -> b) -> f a -> b+apply f = f . vars +{-# INLINE apply #-}+ skeleton :: Traversable f => f a -> f Int skeleton = snd . mapAccumL (\ !n _ -> (n + 1, n)) 0 {-# INLINE skeleton #-}@@ -86,6 +94,32 @@             where ix' = addToIndex i ix {-# INLINE ds #-} +{-+vvars :: Num a => Vector a -> Vector (AD Sparse a)+vvars = Vector.imap (\n a -> AD $ Sparse a $ singleton n $ lift 1)+{-# INLINE vvars #-}++vapply :: Num a => (Vector (AD Sparse a) -> b) -> Vector a -> b+vapply f = f . vvars +{-# INLINE vapply #-}+++vd :: Num a => Int -> AD Sparse a -> Vector a+vd n (AD (Sparse _ da)) = Vector.generate n $ \i -> maybe 0 primal $ lookup i da+{-# INLINE vd #-}++vd' :: Num a => Int -> AD Sparse a -> (a, Vector a)+vd' n (AD (Sparse a da)) = (a , Vector.generate n $ \i -> maybe 0 primal $ lookup i da)+{-# INLINE vd' #-}++vds :: Num a => Int -> AD Sparse a -> Stream Vector a+vds n (AD as@(Sparse a _)) = a :< Vector.generate n (go emptyIndex)+    where+        go ix i = partial (indices ix') as :< Vector.generate n (go ix')+            where ix' = addToIndex i ix+{-# INLINE vds #-}+-}+ partial :: Num a => [Int] -> Sparse a -> a partial []     (Sparse a _)  = a partial (n:ns) (Sparse _ da) = partial ns $ findWithDefault (lift 0) n da@@ -98,6 +132,8 @@     spartial ns a' {-# INLINE spartial #-} ++ instance Primal Sparse where     primal (Sparse a _) = a @@ -134,3 +170,52 @@             (mapWithKey (times dadc) dc)  deriveLifted id $ conT ''Sparse+++class Num a => Grad i o o' a | i -> a o o', o -> a i o', o' -> a i o where+    pack :: i -> [AD Sparse a] -> AD Sparse a+    unpack :: ([a] -> [a]) -> o+    unpack' :: ([a] -> (a, [a])) -> o'++instance Num a => Grad (AD Sparse a) [a] (a, [a]) a where+    pack i _ = i+    unpack f = f []+    unpack' f = f []++instance Grad i o o' a => Grad (AD Sparse a -> i) (a -> o) (a -> o') a where+    pack f (a:as) = pack (f a) as+    pack _ [] = error "Grad.pack: logic error"+    unpack f a = unpack (f . (a:))+    unpack' f a = unpack' (f . (a:))++vgrad :: Grad i o o' a => i -> o+vgrad i = unpack (unsafeGrad (pack i))+    where+        unsafeGrad f as = d as $ apply f as+{-# INLINE vgrad #-}++vgrad' :: Grad i o o' a => i -> o'+vgrad' i = unpack' (unsafeGrad' (pack i))+    where+        unsafeGrad' f as = d' as $ apply f as+{-# INLINE vgrad' #-}++class Num a => Grads i o a | i -> a o, o -> a i where+    packs :: i -> [AD Sparse a] -> AD Sparse a+    unpacks :: ([a] -> Stream [] a) -> o++instance Num a => Grads (AD Sparse a) (Stream [] a) a where+    packs i _ = i+    unpacks f = f []++instance Grads i o a => Grads (AD Sparse a -> i) (a -> o) a where+    packs f (a:as) = packs (f a) as+    packs _ [] = error "Grad.pack: logic error"+    unpacks f a = unpacks (f . (a:))++vgrads :: Grads i o a => i -> o+vgrads i = unpacks (unsafeGrads (packs i))+    where+        unsafeGrads f as = ds as $ apply f as+{-# INLINE vgrads #-}+

@@ -59,6 +59,7 @@ tailS (_ :< as) = as {-# INLINE tailS #-} + unfoldS :: Functor f => (a -> (b, f a)) -> a -> Stream f b unfoldS f a = h :< unfoldS f <$> t      where

@@ -62,7 +62,6 @@         distribute :: Functor f => f (Tensors f a) -> Tensors f (f a)         distribute x = (headT <$> x) :- distribute (tailT <$> x) - instance Typeable1 f => Typeable1 (Tensors f) where     typeOf1 tfa = mkTyConApp tensorsTyCon [typeOf1 (undefined `asArgsType` tfa)]         where asArgsType :: f a -> t f a -> f a

@@ -40,9 +40,6 @@     , du'     , duF     , duF'-    -- * Monadic Combinators-    , diffM-    , diffM'     -- * Exposed Types     , UU, UF, FU, FF     , AD(..)@@ -51,7 +48,6 @@  import Data.Traversable (Traversable) import Control.Applicative-import Control.Monad (liftM) import Numeric.AD.Types import Numeric.AD.Internal.Classes import Numeric.AD.Internal.Composition@@ -97,16 +93,6 @@ diffF' :: (Functor f, Num a) => UF f a -> a -> f (a, a) diffF' f a = unbundle <$> apply f a {-# INLINE diffF' #-}---- | The 'dUM' function calculates the first derivative of scalar-to-scalar monadic function by F'orward' 'AD'-diffM :: (Monad m, Num a) => UF m a -> a -> m a-diffM f a = tangent `liftM` apply f a-{-# INLINE diffM #-}---- | The 'd'UM' function calculates the result and first derivative of a scalar-to-scalar monadic function by F'orward' 'AD'-diffM' :: (Monad m, Num a) => UF m a -> a -> m (a, a)-diffM' f a = unbundle `liftM` apply f a-{-# INLINE diffM' #-}  -- | A fast, simple transposed Jacobian computed with forward-mode AD. jacobianT :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> f (g a)

@@ -1,4 +1,4 @@-{-# LANGUAGE Rank2Types, TypeFamilies #-}+{-# LANGUAGE Rank2Types, TypeFamilies, PatternGuards #-} ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.AD.Mode.Mixed@@ -11,6 +11,31 @@ -- Mixed-Mode Automatic Differentiation. -- -- Each combinator exported from this module chooses an appropriate AD mode.+-- The following basic operations are supported, modified as appropriate by the suffixes below:+-- +-- * 'grad' computes the gradient (partial derivatives) of a function at a point+--+-- * 'jacobian' computes the Jacobian matrix of a function at a point+--+-- * 'diff' computes the derivative of a function at a point+--+-- * 'du' computes a directional derivative of a function at a point+-- +-- * 'hessian' compute the Hessian matrix (matrix of second partial derivatives) of a function at a point+-- +-- The suffixes have the following meanings:+-- +-- * @\'@ -- also return the answer+--+-- * @With@ lets the user supply a function to blend the input with the output+--+-- * @F@ is a version of the base function lifted to return a 'Traversable' (or 'Functor') result+--+-- * @s@ means the function returns all higher derivatives in a list or f-branching 'Stream'+--+-- * @T@ means the result is transposed with respect to the traditional formulation.+--+-- * @0@ means that the resulting derivative list is padded with 0s at the end. -----------------------------------------------------------------------------  module Numeric.AD.Mode.Mixed@@ -21,41 +46,35 @@     , gradWith     , gradWith' -    -- * Jacobians (Mixed Mode)+    -- * Higher Order Gradients (Sparse-on-Reverse)+    , grads++    -- * Jacobians (Sparse or Reverse)     , jacobian     , jacobian'     , jacobianWith     , jacobianWith' -    -- * Monadic Gradient/Jacobian (Reverse Mode)-    , gradM-    , gradM'-    , gradWithM-    , gradWithM'--    -- * Functorial Gradient/Jacobian (Reverse Mode)-    , gradF-    , gradF'-    , gradWithF-    , gradWithF'+    -- * Higher Order Jacobian (Sparse-on-Reverse)+    , jacobians      -- * Transposed Jacobians (Forward Mode)     , jacobianT     , jacobianWithT -    -- * Hessian (Forward-On-Reverse)+    -- * Hessian (Sparse-On-Reverse)     , hessian+    , hessian' -    -- * Hessian Tensors (Forward-On-Mixed)-    , hessianTensor+    -- * Hessian Tensors (Sparse or Sparse-On-Reverse)+    , hessianF+    -- * Hessian Tensors (Sparse)+    , hessianF'      -- * Hessian Vector Products (Forward-On-Reverse)     , hessianProduct     , hessianProduct' -    -- * Higher Order Gradients/Hessians (Sparse Forward)-    , gradients-     -- * Derivatives (Forward Mode)     , diff     , diffF@@ -90,13 +109,16 @@     , maclaurin     , maclaurin0 -    -- * Monadic Combinators (Forward Mode)-    , diffM-    , diffM'+    -- * Unsafe Variadic Grad+    , vgrad+    , vgrad'+    , vgrads      -- * Exposed Types     , module Numeric.AD.Types     , Mode(..)+    , Grad+    , Grads     ) where  import Data.Traversable (Traversable)@@ -104,88 +126,84 @@ import Control.Applicative  import Numeric.AD.Types-import Numeric.AD.Internal.Identity (probed, unprobe) import Numeric.AD.Internal.Composition import Numeric.AD.Classes (Mode(..)) -import qualified Numeric.AD.Mode.Forward as Forward import Numeric.AD.Mode.Forward      ( diff, diff', diffF, diffF'     , du, du', duF, duF'-    , diffM, diffM'-    , jacobianT, jacobianWithT-    ) +    , jacobianT, jacobianWithT )   import Numeric.AD.Mode.Tower      ( diffsF, diffs0F, diffs, diffs0     , taylor, taylor0, maclaurin, maclaurin0-    , dus, dus0, dusF, dus0F-    )+    , dus, dus0, dusF, dus0F )  import qualified Numeric.AD.Mode.Reverse as Reverse import Numeric.AD.Mode.Reverse -    ( grad, grad', gradWith, gradWith'-    , gradM, gradM', gradWithM, gradWithM'-    , gradF, gradF', gradWithF, gradWithF'-    )+    ( grad, grad', gradWith, gradWith', vgrad, vgrad', Grad)  -- temporary until we make a full sparse mode-import qualified Numeric.AD.Internal.Sparse as Sparse-+import qualified Numeric.AD.Mode.Sparse as Sparse+import Numeric.AD.Mode.Sparse+    ( grads, jacobians, hessian', hessianF', vgrads, Grads)+     -- | Calculate the Jacobian of a non-scalar-to-non-scalar function, automatically choosing between forward and reverse mode AD based on the number of inputs and outputs. ----- If you need to support functions where the output is only a 'Functor' or 'Monad', consider 'Numeric.AD.Reverse.jacobian' or 'Numeric.AD.Reverse.gradM' from "Numeric.AD.Reverse".-jacobian :: (Traversable f, Traversable g, Num a) => FF f g a -> f a -> g (f a)+-- If you know the relative number of inputs and outputs, consider 'Numeric.AD.Reverse.jacobian' or 'Nuneric.AD.Sparse.jacobian'.+jacobian :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f a) jacobian f bs = snd <$> jacobian' f bs {-# INLINE jacobian #-} --- | Calculate both the answer and Jacobian of a non-scalar-to-non-scalar function, automatically choosing between forward- and reverse- mode AD based on the relative, number of inputs and outputs.+data Nat = Z | S Nat deriving (Eq, Ord)++size :: Foldable f => f a -> Nat+size = foldr' (\_ b -> S b) Z ++big :: Nat -> Bool+big (S (S (S (S (S (S (S (S (S (S _)))))))))) = True+big _ = False++-- | Calculate both the answer and Jacobian of a non-scalar-to-non-scalar function, automatically choosing between forward- and reverse- mode AD based on the relative, based on the number of inputs ----- If you need to support functions where the output is only a 'Functor' or 'Monad', consider 'Numeric.AD.Reverse.jacobian'' or 'Numeric.AD.Reverse.gradM'' from "Numeric.AD.Reverse".-jacobian' :: (Traversable f, Traversable g, Num a) => FF f g a -> f a -> g (a, f a)-jacobian' f bs | n == 0    = fmap (\x -> (unprobe x, bs)) as-               | n > m     = Reverse.jacobian' f bs-               | otherwise = Forward.jacobian' f bs+-- If you know the relative number of inputs and outputs, consider 'Numeric.AD.Reverse.jacobian'' or 'Nuneric.AD.Sparse.jacobian''.+jacobian' :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (a, f a)+jacobian' f bs | Z <- n = fmap (\x -> (unprobe x, bs)) (f (probed bs))+               | big n  = Reverse.jacobian' f bs+               | otherwise = Sparse.jacobian' f bs     where-        as = f (probed bs)         n = size bs-        m = size as-        size :: Foldable f => f a -> Int-        size = foldr' (\_ b -> 1 + b) 0 {-# INLINE jacobian' #-}  -- | @'jacobianWith' g f@ calculates the Jacobian of a non-scalar-to-non-scalar function, automatically choosing between forward and reverse mode AD based on the number of inputs and outputs. -- -- The resulting Jacobian matrix is then recombined element-wise with the input using @g@. ----- If you need to support functions where the output is only a 'Functor' or 'Monad', consider 'Numeric.AD.Reverse.jacobianWith' or 'Numeric.AD.Reverse.gradWithM' from "Numeric.AD.Reverse".-jacobianWith :: (Traversable f, Traversable g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (f b)+-- If you know the relative number of inputs and outputs, consider 'Numeric.AD.Reverse.jacobianWith' or 'Nuneric.AD.Sparse.jacobianWith'.+jacobianWith :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (f b) jacobianWith g f bs = snd <$> jacobianWith' g f bs {-# INLINE jacobianWith #-} --- | @'jacobianWith'' g f@ calculates the answer and Jacobian of a non-scalar-to-non-scalar function, automatically choosing between forward and reverse mode AD based on the number of inputs and outputs.+-- | @'jacobianWith'' g f@ calculates the answer and Jacobian of a non-scalar-to-non-scalar function, automatically choosing between sparse and reverse mode AD based on the number of inputs and outputs. -- -- The resulting Jacobian matrix is then recombined element-wise with the input using @g@. ----- If you need to support functions where the output is only a 'Functor' or 'Monad', consider 'Numeric.AD.Reverse.jacobianWith'' or 'Numeric.AD.Reverse.gradWithM'' from "Numeric.AD.Reverse".-jacobianWith' :: (Traversable f, Traversable g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (a, f b)+-- If you know the relative number of inputs and outputs, consider 'Numeric.AD.Reverse.jacobianWith'' or 'Nuneric.AD.Sparse.jacobianWith''.+jacobianWith' :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (a, f b) jacobianWith' g f bs-    | n == 0    = fmap (\x -> (unprobe x, undefined <$> bs)) as-    | n > m     = Reverse.jacobianWith' g f bs-    | otherwise = Forward.jacobianWith' g f bs+    | Z <- n = fmap (\x -> (unprobe x, undefined <$> bs)) (f (probed bs))+    | big n  = Reverse.jacobianWith' g f bs+    | otherwise = Sparse.jacobianWith' g f bs     where-        as = f (probed bs)         n = size bs-        m = size as-        size :: Foldable f => f a -> Int-        size = foldr' (\_ b -> 1 + b) 0 {-# INLINE jacobianWith' #-}  -- | @'hessianProduct' f wv@ computes the product of the hessian @H@ of a non-scalar-to-scalar function @f@ at @w = 'fst' <$> wv@ with a vector @v = snd <$> wv@ using \"Pearlmutter\'s method\" from <http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.29.6143>, which states: -- -- > H v = (d/dr) grad_w (w + r v) | r = 0 -- --- Or in other words, we take the directional derivative of the gradient.+-- Or in other words, we take the directional derivative of the gradient. The gradient is calculated in reverse mode, then the directional derivative is calculated in forward mode.+-- hessianProduct :: (Traversable f, Num a) => FU f a -> f (a, a) -> f a hessianProduct f = duF (grad (decomposeMode . f . fmap composeMode)) @@ -193,32 +211,16 @@ -- -- > H v = (d/dr) grad_w (w + r v) | r = 0 -- --- Or in other words, we take the directional derivative of the gradient.--- +-- Or in other words, we return the gradient and the directional derivative of the gradient. The gradient is calculated in reverse mode, then the directional derivative is calculated in forward mode. hessianProduct' :: (Traversable f, Num a) => FU f a -> f (a, a) -> f (a, a) hessianProduct' f = duF' (grad (decomposeMode . f . fmap composeMode)) --- hessianProductWith' :: (Traversable f, Num a) => (a -> a -> a -> a -> b) -> (forall s. Mode s. f (AD s a) -> AD s a) -> f (a, a) -> f b---- | Compute the hessian via the jacobian of the gradient. gradient is computed in reverse mode and then the jacobian is computed in forward mode.+-- | Compute the Hessian via the Jacobian of the gradient. gradient is computed in reverse mode and then the Jacobian is computed in sparse (forward) mode. hessian :: (Traversable f, Num a) => FU f a -> f a -> f (f a)-hessian f = Forward.jacobian (grad (decomposeMode . f . fmap composeMode))---- | Compute the order 3 Hessian tensor on a non-scalar-to-non-scalar function via the forward-mode Jacobian of the mixed-mode Jacobian of the function.-hessianTensor :: (Traversable f, Traversable g, Num a) => FF f g a -> f a -> g (f (f a))-hessianTensor f = decomposeFunctor . Forward.jacobian (ComposeFunctor . jacobian (fmap decomposeMode . f . fmap composeMode))---- data f :> a = a :< f (f :> a)--- data f :- a = a :- (f :- f a) | Zero-{--flatten :: (f :> a) -> (f :- a)-grads :: (Traversable f, Num a) => (forall s. Mode s => f (AD s a) -> AD s a) -> f a -> (f :- a) -grads f b = a :- da :- d2a :- Zero-    (a, da) = grad2 f a-    dda = Forward.jacobian (grad (decomposeMode . f . fmap composeMode)-    ddda = Forward--}+hessian f = Sparse.jacobian (grad (decomposeMode . f . fmap composeMode)) -gradients :: (Traversable f, Num a) => FU f a -> f a -> Stream f a-gradients f as = Sparse.ds as $ f $ Sparse.vars as-    +-- | Compute the order 3 Hessian tensor on a non-scalar-to-non-scalar function using Sparse or Sparse-on-Reverse +hessianF :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f (f a))+hessianF f as +    | big (size as) = decomposeFunctor $ Sparse.jacobian (ComposeFunctor . Reverse.jacobian (fmap decomposeMode . f . fmap composeMode)) as+    | otherwise = Sparse.hessianF f as

@@ -1,4 +1,5 @@-{-# LANGUAGE Rank2Types, TemplateHaskell, BangPatterns #-}+-- {-# LANGUAGE Rank2Types, TemplateHaskell, BangPatterns #-}+{-# LANGUAGE Rank2Types, TemplateHaskell, BangPatterns, MultiParamTypeClasses, FunctionalDependencies, FlexibleInstances, UndecidableInstances, ScopedTypeVariables #-} ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.AD.Mode.Reverse@@ -24,6 +25,7 @@     , grad'     , gradWith     , gradWith'+     -- * Jacobian     , jacobian     , jacobian'@@ -31,34 +33,22 @@     , jacobianWith'     -- * Hessian     , hessian-    , hessianM-    , hessianTensor-    +    , hessianF     -- * Derivatives     , diff     , diff'     , diffF     , diffF'-    -- * Monadic Combinators-    , diffM-    , diffM'-    , gradM-    , gradM'-    , gradWithM-    , gradWithM'-    -- * Synonyms-    , gradF-    , gradF'-    , gradWithF-    , gradWithF'+    -- * Unsafe Variadic Gradient+    , vgrad, vgrad'     -- * Exposed Types     , UU, UF, FU, FF     , AD(..)     , Mode(..)+    , Grad     ) where -import Control.Monad (liftM)-import Control.Applicative (WrappedMonad(..),(<$>))+import Control.Applicative ((<$>)) import Data.Traversable (Traversable)  import Numeric.AD.Types@@ -99,23 +89,14 @@           r = f vs {-# INLINE gradWith' #-} --- | The 'gradF' function calculates the jacobian of a non-scalar-to-non-scalar function with reverse AD lazily in @m@ passes for @m@ outputs.-gradF :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f a)-gradF = jacobian-{-# INLINE gradF #-}---- | An alias for 'gradF'+-- | The 'jacobian' function calculates the jacobian of a non-scalar-to-non-scalar function with reverse AD lazily in @m@ passes for @m@ outputs. jacobian :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f a) jacobian f as = unbind vs . partialArray bds <$> f vs where     (vs, bds) = bind as {-# INLINE jacobian #-} --- | The 'gradF'' function calculates both the result and the Jacobian of a nonscalar-to-nonscalar function, using @m@ invocations of reverse AD,--- where @m@ is the output dimensionality. Applying @fmap snd@ to the result will recover the result of 'gradF'-gradF' :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (a, f a)-gradF' = jacobian' -{-# INLINE gradF' #-}-+-- | The 'jacobian'' function calculates both the result and the Jacobian of a nonscalar-to-nonscalar function, using @m@ invocations of reverse AD,+-- where @m@ is the output dimensionality. Applying @fmap snd@ to the result will recover the result of 'jacobian' -- | An alias for 'gradF'' jacobian' :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (a, f a) jacobian' f as = row <$> f vs where@@ -123,39 +104,29 @@     row a = (primal a, unbind vs (partialArray bds a)) {-# INLINE jacobian' #-} --- | 'gradWithF g f' calculates the Jacobian of a non-scalar-to-non-scalar function @f@ with reverse AD lazily in @m@ passes for @m@ outputs.+-- | 'jacobianWith g f' calculates the Jacobian of a non-scalar-to-non-scalar function @f@ with reverse AD lazily in @m@ passes for @m@ outputs. -- -- Instead of returning the Jacobian matrix, the elements of the matrix are combined with the input using the @g@. ----- > gradF == gradWithF (\_ dx -> dx)--- > gradWithF const == (\f x -> const x <$> f x)+-- > jacobian == jacobianWith (\_ dx -> dx)+-- > jacobianWith const == (\f x -> const x <$> f x) ---gradWithF :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (f b)-gradWithF g f as = unbindWith g vs . partialArray bds <$> f vs where-    (vs, bds) = bind as-{-# INLINE gradWithF #-}---- | An alias for 'gradWithF'. jacobianWith :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (f b)-jacobianWith = gradWithF +jacobianWith g f as = unbindWith g vs . partialArray bds <$> f vs where+    (vs, bds) = bind as {-# INLINE jacobianWith #-} --- | 'gradWithF' g f' calculates both the result and the Jacobian of a nonscalar-to-nonscalar function @f@, using @m@ invocations of reverse AD,--- where @m@ is the output dimensionality. Applying @fmap snd@ to the result will recover the result of 'gradWithF'+-- | 'jacobianWith' g f' calculates both the result and the Jacobian of a nonscalar-to-nonscalar function @f@, using @m@ invocations of reverse AD,+-- where @m@ is the output dimensionality. Applying @fmap snd@ to the result will recover the result of 'jacobianWith' -- -- Instead of returning the Jacobian matrix, the elements of the matrix are combined with the input using the @g@. ----- > jacobian' == gradWithF' (\_ dx -> dx)+-- > jacobian' == jacobianWith' (\_ dx -> dx) ---gradWithF' :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (a, f b)-gradWithF' g f as = row <$> f vs where+jacobianWith' :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (a, f b)+jacobianWith' g f as = row <$> f vs where     (vs, bds) = bind as     row a = (primal a, unbindWith g vs (partialArray bds a))-{-# INLINE gradWithF' #-}---- | An alias for 'gradWithF''-jacobianWith' :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (a, f b)-jacobianWith' = gradWithF' {-# INLINE jacobianWith' #-}  diff :: Num a => UU a -> a -> a@@ -176,44 +147,15 @@ diffF' f a = derivative' <$> f (var a 0) {-# INLINE diffF' #-} --- * Monadic Combinators--diffM :: (Monad m, Num a) => UF m a -> a -> m a-diffM f a = liftM derivative $ f (var a 0)-{-# INLINE diffM #-}--diffM' :: (Monad m, Num a) => UF m a -> a -> m (a, a)-diffM' f a = liftM derivative' $ f (var a 0)-{-# INLINE diffM' #-}--gradM :: (Traversable f, Monad m, Num a) => FF f m a -> f a -> m (f a)-gradM f = unwrapMonad . jacobian (WrapMonad . f)-{-# INLINE gradM #-}--gradM' :: (Traversable f, Monad m, Num a) => FF f m a -> f a -> m (a, f a)-gradM' f = unwrapMonad . jacobian' (WrapMonad . f)-{-# INLINE gradM' #-}--gradWithM :: (Traversable f, Monad m, Num a) => (a -> a -> b) -> FF f m a -> f a -> m (f b)-gradWithM g f = unwrapMonad . jacobianWith g (WrapMonad . f)--gradWithM' :: (Traversable f, Monad m, Num a) => (a -> a -> b) -> FF f m a -> f a -> m (a, f b)-gradWithM' g f = unwrapMonad . jacobianWith' g (WrapMonad . f)- -- | Compute the hessian via the jacobian of the gradient. gradient is computed in reverse mode and then the jacobian is computed in reverse mode. ----- However, since the @'grad f :: f a -> f a'@ is square this is not as fast as using the forward-mode Jacobian of a reverse mode gradient provided by 'Numeric.AD.hessian' in "Numeric.AD".+-- However, since the @'grad f :: f a -> f a'@ is square this is not as fast as using the forward-mode Jacobian of a reverse mode gradient provided by 'Numeric.AD.hessian'. hessian :: (Traversable f, Num a) => FU f a -> f a -> f (f a) hessian f = jacobian (grad (decomposeMode . f . fmap composeMode)) --- | Compute the order 3 Hessian tensor on a non-scalar-to-non-scalar function via the forward-mode Jacobian of the mixed-mode Jacobian of the function.+-- | Compute the order 3 Hessian tensor on a non-scalar-to-non-scalar function via the reverse-mode Jacobian of the reverse-mode Jacobian of the function. ----- While this is less efficient than 'Numeric.AD.hessianTensor' from "Numeric.AD" or 'Numeric.AD.Forward.hessianTensor' from "Numeric.AD.Forward", the type signature is more permissive with regards to the output non-scalar, and it may be more efficient if only a few coefficients of the result are consumed.-hessianTensor :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f (f a))-hessianTensor f = decomposeFunctor . jacobian (ComposeFunctor . jacobian (fmap decomposeMode . f . fmap composeMode))+-- Less efficient than 'Numeric.AD.Mode.Mixed.hessianF'.+hessianF :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f (f a))+hessianF f = decomposeFunctor . jacobian (ComposeFunctor . jacobian (fmap decomposeMode . f . fmap composeMode)) --- | Compute the hessian via the reverse-mode jacobian of the reverse-mode gradient of a non-scalar-to-scalar monadic action. ------ While this is less efficient than 'Numeric.AD.hessianTensor' from "Numeric.AD" or 'Numeric.AD.Forward.hessianTensor' from "Numeric.AD.Forward", the type signature is more permissive with regards to the output non-scalar, and it may be more efficient if only a few coefficients of the result are consumed.-hessianM :: (Traversable f, Monad m, Num a) => FF f m a -> f a -> m (f (f a))-hessianM f = unwrapMonad . hessianTensor (WrapMonad . f)

@@ -0,0 +1,123 @@+{-# LANGUAGE Rank2Types, BangPatterns #-}+-----------------------------------------------------------------------------+-- |+-- Module      : Numeric.AD.Mode.Sparse+-- Copyright   : (c) Edward Kmett 2010+-- License     : BSD3+-- Maintainer  : ekmett@gmail.com+-- Stability   : experimental+-- Portability : GHC only+--+-- Higher order derivatives via a \"dual number tower\".+--+-----------------------------------------------------------------------------++module Numeric.AD.Mode.Sparse+    (+    -- * Sparse Gradients+      grad+    , grad'+    , gradWith+    , gradWith'+    , grads+    +    -- * Sparse Jacobians (synonyms)+    , jacobian+    , jacobian'+    , jacobianWith+    , jacobianWith'+    , jacobians++    -- * Sparse Hessians+    , hessian+    , hessian'++    , hessianF+    , hessianF'++    -- * Unsafe gradients+    , vgrad+    , vgrads++    -- * Exposed Types+    , module Numeric.AD.Types+    , Mode(..)+    , Grad+    , Grads+    ) where++import Control.Applicative ((<$>))+import Data.Traversable+import Numeric.AD.Types+import Numeric.AD.Classes+import Numeric.AD.Internal.Sparse+import Numeric.AD.Internal.Combinators++second :: (a -> b) -> (c, a) -> (c, b)+second g (a,b) = (a, g b)+{-# INLINE second #-}++grad :: (Traversable f, Num a) => FU f a -> f a -> f a+grad f as = d as $ apply f as+{-# INLINE grad #-}++grad' :: (Traversable f, Num a) => FU f a -> f a -> (a, f a)+grad' f as = d' as $ apply f as+{-# INLINE grad' #-}++gradWith :: (Traversable f, Num a) => (a -> a -> b) -> FU f a -> f a -> f b+gradWith g f as = zipWithT g as $ grad f as+{-# INLINE gradWith #-}++gradWith' :: (Traversable f, Num a) => (a -> a -> b) -> FU f a -> f a -> (a, f b)+gradWith' g f as = second (zipWithT g as) $ grad' f as+{-# INLINE gradWith' #-}++jacobian :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f a)+jacobian f as = d as <$> apply f as+{-# INLINE jacobian #-}++jacobian' :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (a, f a)+jacobian' f as = d' as <$> apply f as+{-# INLINE jacobian' #-}++jacobianWith :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (f b)+jacobianWith g f as = zipWithT g as <$> jacobian f as+{-# INLINE jacobianWith #-}++jacobianWith' :: (Traversable f, Functor g, Num a) => (a -> a -> b) -> FF f g a -> f a -> g (a, f b)+jacobianWith' g f as = second (zipWithT g as) <$> jacobian' f as+{-# INLINE jacobianWith' #-}++grads :: (Traversable f, Num a) => FU f a -> f a -> Stream f a+grads f as = ds as $ apply f as+{-# INLINE grads #-}++jacobians :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (Stream f a)+jacobians f as = ds as <$> apply f as+{-# INLINE jacobians #-}++d2 :: Functor f => Stream f a -> f (f a)+d2 = headT . tailT . tailT . tensors +{-# INLINE d2 #-}++d2' :: Functor f => Stream f a -> (a, f (a, f a))+d2' (a :< as) = (a, fmap (\(da :< das) -> (da, fmap headS das)) as)+{-# INLINE d2' #-}++hessian :: (Traversable f, Num a) => FU f a -> f a -> f (f a)+hessian f as = d2 $ grads f as+{-# INLINE hessian #-}++hessian' :: (Traversable f, Num a) => FU f a -> f a -> (a, f (a, f a))+hessian' f as = d2' $ grads f as+{-# INLINE hessian' #-}++hessianF :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (f (f a))+hessianF f as = d2 <$> jacobians f as+{-# INLINE hessianF #-}++hessianF' :: (Traversable f, Functor g, Num a) => FF f g a -> f a -> g (a, f (a, f a))+hessianF' f as = d2' <$> jacobians f as+{-# INLINE hessianF' #-}+

@@ -36,16 +36,12 @@     , duF'    -- answer and directional derivative of (a -> f a)     , dusF    -- answer and all directional derivatives of (a -> f a)     , dus0F   -- answer and all zero padded directional derivatives of (a -> a)-    -- * Monadic Combinators-    , diffsM  -- answer and all derivatives of the monadic action (a -> m a)-    , diffs0M -- answer and all zero padded derivatives of (a -> m a)     -- * Exposed Types     , UU, UF, FU, FF     , Mode(..)     , AD(..)     ) where -import Control.Monad (liftM) import Control.Applicative ((<$>)) import Numeric.AD.Types import Numeric.AD.Classes@@ -66,14 +62,6 @@ diffs0F :: (Functor f, Num a) => UF f a -> a -> f [a] diffs0F f a = (zeroPad . getADTower) <$> apply f a {-# INLINE diffs0F #-}--diffsM :: (Monad m, Num a) => UF m a -> a -> m [a]-diffsM f a = getADTower `liftM` apply f a-{-# INLINE diffsM #-}--diffs0M :: (Monad m, Num a) => UF m a -> a -> m [a]-diffs0M f a = (zeroPad . getADTower) `liftM` apply f a-{-# INLINE diffs0M #-}  taylor :: Fractional a => UU a -> a -> a -> [a] taylor f x dx = go 1 1 (diffs f x)

@@ -14,18 +14,12 @@     (     -- * Newton's Method (Forward AD)       findZero-    , findZeroM     , inverse-    , inverseM     , fixedPoint-    , fixedPointM     , extremum-    , extremumM     -- * Gradient Ascent/Descent (Reverse AD)     , gradientDescent-    , gradientDescentM     , gradientAscent-    , gradientAscentM     -- * Exposed Types     , UU, UF, FU, FF     , AD(..)@@ -33,14 +27,12 @@     ) where  import Prelude hiding (all)-import Control.Monad (liftM)-import Data.MList import Data.Foldable (all) import Data.Traversable (Traversable) import Numeric.AD.Types import Numeric.AD.Classes-import Numeric.AD.Mode.Forward (diff, diff', diffM, diffM')-import Numeric.AD.Mode.Reverse (gradWith', gradWithM')+import Numeric.AD.Mode.Forward (diff, diff')+import Numeric.AD.Mode.Reverse (gradWith') import Numeric.AD.Internal.Composition  -- | The 'findZero' function finds a zero of a scalar function using@@ -62,16 +54,6 @@                 (y,y') = diff' f x {-# INLINE findZero #-} -findZeroM :: (Monad m, Fractional a) => UF m a -> a -> MList m a-findZeroM f x0 = MList (go x0)-    where-        go x = return $ -               MCons x $ -               MList $ do-                (y,y') <- diffM' f x-                go (x - y/y')-{-# INLINE findZeroM #-}- -- | The 'inverseNewton' function inverts a scalar function using -- Newton's method; its output is a stream of increasingly accurate -- results.  (Modulo the usual caveats.)@@ -84,10 +66,6 @@ inverse f x0 y = findZero (\x -> f x - lift y) x0 {-# INLINE inverse  #-} -inverseM :: (Monad m, Fractional a) => UF m a -> a -> a -> MList m a-inverseM f x0 y = findZeroM (\x -> subtract (lift y) `liftM` f x) x0-{-# INLINE inverseM  #-}- -- | The 'fixedPoint' function find a fixedpoint of a scalar -- function using Newton's method; its output is a stream of -- increasingly accurate results.  (Modulo the usual caveats.)@@ -97,10 +75,6 @@ fixedPoint f = findZero (\x -> f x - x) {-# INLINE fixedPoint #-} -fixedPointM :: (Monad m, Fractional a) => UF m a -> a -> MList m a-fixedPointM f = findZeroM (\x -> subtract x `liftM` f x)-{-# INLINE fixedPointM #-}- -- | The 'extremum' function finds an extremum of a scalar -- function using Newton's method; produces a stream of increasingly -- accurate results.  (Modulo the usual caveats.)@@ -110,10 +84,6 @@ extremum f = findZero (diff (decomposeMode . f . composeMode)) {-# INLINE extremum #-} -extremumM :: (Monad m, Fractional a) => UF m a -> a -> MList m a-extremumM f = findZeroM (diffM (liftM decomposeMode . f . composeMode))-{-# INLINE extremumM #-}- -- | The 'gradientDescent' function performs a multivariate -- optimization, based on the naive-gradient-descent in the file -- @stalingrad\/examples\/flow-tests\/pre-saddle-1a.vlad@ from the@@ -141,31 +111,3 @@ gradientAscent :: (Traversable f, Fractional a, Ord a) => FU f a -> f a -> [f a] gradientAscent f = gradientDescent (negate . f) {-# INLINE gradientAscent #-}---- monadic gradient descent-gradientDescentM :: (Traversable f, Monad m, Fractional a, Ord a) => FF f m a -> f a -> MList m (f a)-gradientDescentM f x0 = MList $ do-        (fx0, xgx0) <- gradWithM' (,) f x0-        go x0 fx0 xgx0 0.1 (0 :: Int)-    where-        go x fx xgx !eta !i-            | eta == 0  = return MNil -- step size is 0-            | otherwise = do-                (fx1, xgx1) <- gradWithM' (,) f x1-                case () of-                 _ | fx1 > fx     -> go x fx xgx (eta/2) 0 -- we stepped too far-                   | zeroGrad xgx -> return MNil -- gradient is 0-                   | otherwise    -> return $ -                                     MCons x1 $ -                                     MList $-                                     if i == 10-                                     then go x1 fx1 xgx1 (eta*2) 0-                                     else go x1 fx1 xgx1 eta (i+1)-            where-                x1 = fmap (\(xi,gxi) -> xi - eta * gxi) xgx-                zeroGrad = all (\(_,g) -> g == 0)-{-# INLINE gradientDescentM #-}--gradientAscentM :: (Traversable f, Monad m, Fractional a, Ord a) => FF f m a -> f a -> MList m (f a)-gradientAscentM f = gradientDescentM (liftM negate . f)-{-# INLINE gradientAscentM #-}

@@ -1,3 +1,4 @@+{-# LANGUAGE Rank2Types #-} ----------------------------------------------------------------------------- -- | -- Module      :  Numeric.AD.Types@@ -24,8 +25,33 @@     , headS     , tailS     , unfoldS+    -- * An Identity Mode. +    , Id(..)+    , probe, unprobe+    , probed, unprobed+    -- * Apply functions that use 'lift'+    , lowerUU, lowerUF, lowerFU, lowerFF     ) where +import Numeric.AD.Internal.Identity import Numeric.AD.Internal.Types import Numeric.AD.Internal.Stream import Numeric.AD.Internal.Tensors++-- these exploit the 'magic' that is probed to avoid the need for Functor, etc.++lowerUU :: UU a -> a -> a+lowerUU f = unprobe . f . probe+{-# INLINE lowerUU #-}++lowerUF :: UF f a -> a -> f a+lowerUF f = unprobed . f . probe+{-# INLINE lowerUF #-}++lowerFU :: FU f a -> f a -> a+lowerFU f = unprobe . f . probed+{-# INLINE lowerFU #-}++lowerFF :: FF f g a -> f a -> g a+lowerFF f = unprobed . f . probed+{-# INLINE lowerFF #-}

@@ -1,5 +1,5 @@ name:         ad-version:      0.40.1+version:      0.44.0 license:      BSD3 license-File: LICENSE copyright:    (c) Edward Kmett 2010,@@ -16,8 +16,61 @@     Type-level \"branding\" is used to both prevent the end user from confusing infinitesimals      and to limit unsafe access to the implementation details of each Mode.     .-    The combinators in "Numeric.AD" choose from a variety of automatic differentiation modes,-    based on the arity of their inputs and outputs.+    Each mode has a separate module full of combinators.+    .+    * @Numeric.AD.Mode.Forward@ provides basic forward-mode AD. It is good for computing simple derivatives.+    .+    * @Numeric.AD.Mode.Reverse@ uses benign side-effects to compute reverse-mode AD. It is good for computing gradients in one pass.+    .+    * @Numeric.AD.Mode.Sparse@ computes a sparse forward-mode AD tower. It is good for higher derivatives or large numbers of outputs.+    .+    * @Numeric.AD.Mode.Tower@ computes a dense forward-mode AD tower useful for higher derivatives of single input functions.+    .+    * @Numeric.AD.Mode.Mixed@ computes using whichever mode or combination thereof is suitable to each individual combinator. This mode is the default, re-exported by @Numeric.AD@+    .+    . +    While not every mode can provide all operations, the following basic operations are supported, modified as +    appropriate by the suffixes below:+    .+    * 'grad' computes the gradient (partial derivatives) of a function at a point.+    .+    * 'jacobian' computes the Jacobian matrix of a function at a point.+    .+    * 'diff' computes the derivative of a function at a point.+    .+    * 'du' computes a directional derivative of a function at a point.+    .+    * 'hessian' computes the Hessian matrix (matrix of second partial derivatives) of a function at a point.+    .+    The following suffixes alter the meanings of the functions above as follows:+    .+    * @\'@ -- also return the answer+    .+    * @With@ lets the user supply a function to blend the input with the output+    .+    * @F@ is a version of the base function lifted to return a 'Traversable' (or 'Functor') result+    .+    * @s@ means the function returns all higher derivatives in a list or f-branching 'Stream'+    .+    * @T@ means the result is transposed with respect to the traditional formulation.+    .+    * @0@ means that the resulting derivative list is padded with 0s at the end.+    .+    Changes since 0.40.0+    .+    * Bug fix in the derivative calculated for @'(/)' :: (Mode s, Fractional a) => AD s a@+    .+    * Improved documentation+    .+    * Regularized naming conventions+    .+    * Exposed 'Id', probe, and lower methods via @Numeric.AD.Types@+    .+    * Removed monadic combinators+    .+    * Retuned the 'Mixed' mode jacobian calculations to only require a 'Functor'-based result.+    .+    * Added unsafe variadic 'vgrad', 'vgrad'', and 'vgrads' combinators  build-type:   Simple build-depends:       @@ -25,7 +78,6 @@     data-reify >= 0.5 && < 0.6,      containers >= 0.2 && < 0.4,     template-haskell >= 2.4 && < 2.5,-    mlist >= 0.0.2 && <= 0.1,     array >= 0.2 && < 0.4  exposed-modules:@@ -43,20 +95,20 @@     Numeric.AD.Internal.Sparse     Numeric.AD.Internal.Dense     Numeric.AD.Internal.Composition-    Numeric.AD.Internal.Identity      Numeric.AD.Mode.Directed     Numeric.AD.Mode.Forward     Numeric.AD.Mode.Mixed     Numeric.AD.Mode.Reverse     Numeric.AD.Mode.Tower-+    Numeric.AD.Mode.Sparse  other-modules:     Numeric.AD.Internal.Types     Numeric.AD.Internal.Comonad     Numeric.AD.Internal.Stream     Numeric.AD.Internal.Tensors+    Numeric.AD.Internal.Identity  Extra-Source-Files: TODO GHC-Options: -Wall -fspec-constr -fdicts-cheap -O2