hmatrix and ad don't (currently) mix. The problem is that hmatrix uses a packed structure that can't hold any of the AD mode variants we have as an Element. =(
I've been working with Alex Lang to explore in ad 4.0 ways that we can support monomorphic AD modes and still present largely the same API. We've seen a number of performance gains off of this refactoring already, but it doesn't go far enough to address what you need. A goal a bit farther out is to support AD on vector/matrix operations, but it is a much bigger refactoring than the one currently in the pipeline. =/ To support automatic differentiation on vector-based operations in a form that works with standard BLAS-like storage like the packed matrix rep used in hmatrix we need to convert from a 'matrix of AD variables' to an 'AD mode over of matrices'. This is similar to the difference between a matrix of complex numbers and a real matrix plus an imaginary matrix. This is a long term goal, but not one you're likely to see support for out of 'ad' in the short term. I can't build AD on hmatrix itself due in part to licensing restrictions and differing underlying storage requirements, so there are a lot of little issues in making that latter vision a reality. -Edward On Tue, Apr 9, 2013 at 10:46 AM, Dominic Steinitz <domi...@steinitz.org>wrote: > Hi Cafe, > > Suppose I want to find the grad of a function then it's easy I just > use http://hackage.haskell.org/package/ad-3.4: > > import Numeric.AD > import Data.Foldable (Foldable) > import Data.Traversable (Traversable) > > data MyMatrix a = MyMatrix (a, a) > deriving (Show, Functor, Foldable, Traversable) > > f :: Floating a => MyMatrix a -> a > f (MyMatrix (x, y)) = exp $ negate $ (x^2 + y^2) / 2.0 > > main :: IO () > main = do > putStrLn $ show $ f $ MyMatrix (0.0, 0.0) > putStrLn $ show $ grad f $ MyMatrix (0.0, 0.0) > > But now suppose I am doing some matrix calculations > http://hackage.haskell.org/package/hmatrix-0.14.1.0 and I want to find > the grad of a function of a matrix: > > import Numeric.AD > import Numeric.LinearAlgebra > import Data.Foldable (Foldable) > import Data.Traversable (Traversable) > > g :: (Element a, Floating a) => Matrix a -> a > g m = exp $ negate $ (x^2 + y^2) / 2.0 > where r = (toLists m)!!0 > x = r!!0 > y = r!!1 > > main :: IO () > main = do > putStrLn $ show $ g $ (1 >< 2) ([0.0, 0.0] :: [Double]) > putStrLn $ show $ grad g $ (1 >< 2) ([0.0, 0.0] :: [Double]) > > Then I am in trouble: > > /Users/dom/Dropbox/Private/Whales/MyAD.hs:24:21: > No instance for (Traversable Matrix) arising from a use of `grad' > Possible fix: add an instance declaration for (Traversable Matrix) > In the expression: grad g > In the second argument of `($)', namely > `grad g $ (1 >< 2) ([0.0, 0.0] :: [Double])' > In the second argument of `($)', namely > `show $ grad g $ (1 >< 2) ([0.0, 0.0] :: [Double])' > > /Users/dom/Dropbox/Private/Whales/MyAD.hs:24:26: > Could not deduce (Element > (ad-3.4:Numeric.AD.Internal.Types.AD s Double)) > arising from a use of `g' > from the context (Numeric.AD.Internal.Classes.Mode s) > bound by a type expected by the context: > Numeric.AD.Internal.Classes.Mode s => > Matrix (ad-3.4:Numeric.AD.Internal.Types.AD s Double) > -> ad-3.4:Numeric.AD.Internal.Types.AD s Double > at /Users/dom/Dropbox/Private/Whales/MyAD.hs:24:21-26 > Possible fix: > add an instance declaration for > (Element (ad-3.4:Numeric.AD.Internal.Types.AD s Double)) > In the first argument of `grad', namely `g' > In the expression: grad g > In the second argument of `($)', namely > `grad g $ (1 >< 2) ([0.0, 0.0] :: [Double])' > > What are my options here? Clearly I can convert my matrix into a list > (which is traversable), find the grad and convert it back into a > matrix but given I am doing numerical calculations and speed is an > important factor, this seems undesirable. > > I think I would have the same problem with: > > http://hackage.haskell.org/package/repa > http://hackage.haskell.org/package/yarr-1.3.1 > > although I haven'¯t checked. > > Thanks, Dominic. > > > _______________________________________________ > Haskell-Cafe mailing list > Haskell-Cafe@haskell.org > http://www.haskell.org/mailman/listinfo/haskell-cafe >
_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe