harpie-numhask (empty) → 0.1.0.0

raw patch · 5 files changed

+397/−0 lines, 5 filesdep +adjunctionsdep +basedep +doctest-parallel

Dependencies added: adjunctions, base, doctest-parallel, first-class-families, harpie, harpie-numhask, numhask

Files

ChangeLog.md +6/−0
LICENSE +30/−0
harpie-numhask.cabal +69/−0
src/Harpie/NumHask.hs +284/−0
test/doctests.hs +8/−0

+ ChangeLog.md view

@@ -0,0 +1,6 @@++0.1+===++- Harpie.Fixed orphans +

+ LICENSE view

@@ -0,0 +1,30 @@+Copyright (c) 2016, Tony Day++All rights reserved.++Redistribution and use in source and binary forms, with or without+modification, are permitted provided that the following conditions are met:++    * Redistributions of source code must retain the above copyright+      notice, this list of conditions and the following disclaimer.++    * Redistributions in binary form must reproduce the above+      copyright notice, this list of conditions and the following+      disclaimer in the documentation and/or other materials provided+      with the distribution.++    * Neither the name of Tony Day nor the names of other+      contributors may be used to endorse or promote products derived+      from this software without specific prior written permission.++THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS+"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT+LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR+A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT+OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,+DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY+THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT+(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

+ harpie-numhask.cabal view

@@ -0,0 +1,69 @@+cabal-version: 3.0+name: harpie-numhask+version: 0.1.0.0+license: BSD-3-Clause+license-file: LICENSE+copyright: Tony Day (c) 2024+category: math+author: Tony Day+maintainer: tonyday567@gmail.com+homepage: https://github.com/tonyday567/harpie-numhask#readme+bug-reports: https://github.com/tonyday567/harpie-numhask/issues+synopsis: numhask shim for harpie.+description:+    This package provides orphans and interface into the [numhask](https://hackage.haskell.org/package/numhask) API for [harpie](https://hackage.haskell.org/package/harpie).++    == Usage++    >>> {-# LANGUAGE RebindableSyntax #-}+    >>> import NumHask.Prelude+    >>> import Harpie.Fixed as F+    >>> import Harpie.Numhask qualified as N++build-type: Simple+tested-with:+    , GHC == 9.10.1+    , GHC == 9.6.5+    , GHC == 9.8.2+extra-doc-files: ChangeLog.md++source-repository head+    type: git+    location: https://github.com/tonyday567/harpie-numhask++common ghc-options-stanza+    ghc-options:+        -Wall+        -Wcompat+        -Widentities+        -Wincomplete-record-updates+        -Wincomplete-uni-patterns+        -Wpartial-fields+        -Wredundant-constraints++common ghc2021-stanza+    default-language: GHC2021++library+    import: ghc-options-stanza+    import: ghc2021-stanza+    hs-source-dirs: src+    build-depends:+        , adjunctions          >=4.4.2 && <4.5+        , base                 >=4.14 && <5+        , first-class-families >=0.8 && <0.9+        , harpie               >=0.1.0 && <0.2+        , numhask              >=0.12 && <0.13+    exposed-modules: Harpie.NumHask++test-suite doctests+    import: ghc2021-stanza+    main-is: doctests.hs+    hs-source-dirs: test+    build-depends:+        , base             >=4.14 && <5+        , doctest-parallel >=0.3 && <0.4+        , harpie-numhask+    default-extensions: RebindableSyntax+    ghc-options: -threaded+    type: exitcode-stdio-1.0

+ src/Harpie/NumHask.hs view

@@ -0,0 +1,284 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE RebindableSyntax #-}+{-# LANGUAGE TypeFamilies #-}+{-# OPTIONS_GHC -Wno-incomplete-uni-patterns #-}+{-# OPTIONS_GHC -Wno-orphans #-}+{-# OPTIONS_GHC -Wno-redundant-constraints #-}++-- | numhask orphan instances and shim for harpie.+module Harpie.NumHask+  ( -- * Usage+    -- $usage+    ident,+    undiag,+    mult,+    invtri,+    inverse,+    chol,+  )+where++import Data.Functor.Rep+import Fcf qualified+import GHC.TypeNats+import Harpie.Fixed as F hiding (chol, ident, inverse, invtri, mult, undiag)+import Harpie.Shape (DeleteDims, Eval, GetDims, KnownNats, Rank, type (++))+import Harpie.Shape qualified as S+import NumHask.Prelude as P hiding (Min, cycle, diff, drop, empty, find, length, repeat, sequence, take, zipWith)++-- $setup+--+-- >>> :m -Prelude+-- >>> :set -XDataKinds+-- >>> :set -XRebindableSyntax+-- >>> import NumHask.Prelude hiding (cycle, repeat, empty, diff, take, drop, zipWith, find)+-- >>> import Harpie.Fixed qualified as F+-- >>> import Harpie.Fixed (Array, array, range, shape, toDynamic)+-- >>> import Harpie.NumHask+-- >>> import Prettyprinter hiding (dot,fill)+--+-- >>> s = 1 :: Array '[] Int+-- >>> s+-- [1]+-- >>> shape s+-- []+-- >>> pretty s+-- 1+-- >>> let v = range @'[3]+-- >>> pretty v+-- [0,1,2]+-- >>> let m = range @[2,3]+-- >>> pretty m+-- [[0,1,2],+--  [3,4,5]]+-- >>> a = range @[2,3,4]+-- >>> a+-- [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23]+-- >>> pretty a+-- [[[0,1,2,3],+--   [4,5,6,7],+--   [8,9,10,11]],+--  [[12,13,14,15],+--   [16,17,18,19],+--   [20,21,22,23]]]+-- >>> toDynamic a+-- UnsafeArray [2,3,4] [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23]++-- * NumHask heirarchy++instance+  ( Additive a,+    KnownNats s+  ) =>+  Additive (Array s a)+  where+  (+) = liftR2 (+)++  zero = pureRep zero++instance+  ( Subtractive a,+    KnownNats s+  ) =>+  Subtractive (Array s a)+  where+  negate = fmapRep negate++instance+  (Multiplicative a) =>+  MultiplicativeAction (Array s a)+  where+  type Scalar (Array s a) = a+  (|*) r s = fmap (s *) r++instance (Additive a) => AdditiveAction (Array s a) where+  type AdditiveScalar (Array s a) = a+  (|+) r s = fmap (s +) r++instance+  (Subtractive a) =>+  SubtractiveAction (Array s a)+  where+  (|-) r s = fmap (\x -> x - s) r++instance+  (Divisive a) =>+  DivisiveAction (Array s a)+  where+  (|/) r s = fmap (/ s) r++instance (KnownNats s, JoinSemiLattice a) => JoinSemiLattice (Array s a) where+  (\/) = liftR2 (\/)++instance (KnownNats s, MeetSemiLattice a) => MeetSemiLattice (Array s a) where+  (/\) = liftR2 (/\)++instance (KnownNats s, Subtractive a, Epsilon a) => Epsilon (Array s a) where+  epsilon = konst epsilon++instance (FromInteger a) => FromInteger (Array ('[] :: [Nat]) a) where+  fromInteger x = toScalar (fromInteger x)++instance (FromRational a) => FromRational (Array ('[] :: [Nat]) a) where+  fromRational x = toScalar (fromRational x)++-- | The identity array.+--+-- >>> pretty $ ident @[3,3]+-- [[1,0,0],+--  [0,1,0],+--  [0,0,1]]+ident :: (KnownNats s, Additive a, Multiplicative a) => Array s a+ident = tabulate (bool zero one . S.isDiag . S.fromFins)++-- | Expand the array to form a diagonal array+--+-- >>> pretty $ undiag (range @'[3])+-- [[0,0,0],+--  [0,1,0],+--  [0,0,2]]+undiag ::+  forall s' a s.+  ( KnownNats s,+    KnownNats s',+    s' ~ Eval ((++) s s),+    Additive a+  ) =>+  Array s a ->+  Array s' a+undiag a = tabulate (\xs -> bool zero (index a (S.UnsafeFins $ pure $ S.getDim 0 (S.fromFins xs))) (S.isDiag (S.fromFins xs)))++-- | Array multiplication.+--+-- matrix multiplication+--+-- >>> pretty $ mult m (F.transpose m)+-- [[5,14],+--  [14,50]]+--+-- inner product+--+-- >>> pretty $ mult v v+-- 5+--+-- matrix-vector multiplication+--+-- >>> pretty $ mult v (F.transpose m)+-- [5,14]+--+-- >>> pretty $ mult m v+-- [5,14]+mult ::+  forall a ds0 ds1 s0 s1 so0 so1 st si.+  ( Ring a,+    KnownNats s0,+    KnownNats s1,+    KnownNats ds0,+    KnownNats ds1,+    KnownNats so0,+    KnownNats so1,+    KnownNats st,+    KnownNats si,+    so0 ~ Eval (DeleteDims ds0 s0),+    so1 ~ Eval (DeleteDims ds1 s1),+    si ~ Eval (GetDims ds0 s0),+    si ~ Eval (GetDims ds1 s1),+    st ~ Eval ((++) so0 so1),+    ds0 ~ '[Eval ((Fcf.-) (Eval (Rank s0)) 1)],+    ds1 ~ '[0]+  ) =>+  Array s0 a ->+  Array s1 a ->+  Array st a+mult = dot sum (*)++instance+  ( Multiplicative a,+    P.Distributive a,+    Subtractive a,+    KnownNat m+  ) =>+  Multiplicative (Matrix m m a)+  where+  (*) = mult++  one = ident++instance+  ( Multiplicative a,+    P.Distributive a,+    Subtractive a,+    Eq a,+    ExpField a,+    KnownNat m+  ) =>+  Divisive (Matrix m m a)+  where+  recip a = invtri (transpose (chol a)) * invtri (chol a)++-- | Inverse of a square matrix.+--+-- > mult (inverse a) a == a+--+-- >>> e = array @[3,3] @Double [4,12,-16,12,37,-43,-16,-43,98]+-- >>> pretty (inverse e)+-- [[49.36111111111111,-13.555555555555554,2.1111111111111107],+--  [-13.555555555555554,3.7777777777777772,-0.5555555555555555],+--  [2.1111111111111107,-0.5555555555555555,0.1111111111111111]]+inverse :: (Eq a, ExpField a, KnownNat m) => Matrix m m a -> Matrix m m a+inverse a = mult (invtri (transpose (chol a))) (invtri (chol a))++-- | [Inversion of a Triangular Matrix](https://math.stackexchange.com/questions/1003801/inverse-of-an-invertible-upper-triangular-matrix-of-order-3)+--+-- >>> t = array @[3,3] @Double [1,0,1,0,1,2,0,0,1]+-- >>> pretty (invtri t)+-- [[1.0,0.0,-1.0],+--  [0.0,1.0,-2.0],+--  [0.0,0.0,1.0]]+--+-- > ident == mult t (invtri t)+-- True+invtri :: forall a n. (KnownNat n, ExpField a, Eq a) => Matrix n n a -> Matrix n n a+invtri a = sum (fmap (l ^) (iota @n)) * ti+  where+    ti = undiag (fmap recip (diag a))+    tl = a - undiag (diag a)+    l = negate (ti * tl)++-- | cholesky decomposition+--+-- Uses the <https://en.wikipedia.org/wiki/Cholesky_decomposition#The_Cholesky_algorithm Cholesky-Crout> algorithm.+--+-- >>> e = array @[3,3] @Double [4,12,-16,12,37,-43,-16,-43,98]+-- >>> pretty (chol e)+-- [[2.0,0.0,0.0],+--  [6.0,1.0,0.0],+--  [-8.0,5.0,3.0]]+-- >>> mult (chol e) (F.transpose (chol e)) == e+-- True+chol :: (KnownNat m, ExpField a) => Matrix m m a -> Matrix m m a+chol a =+  let l =+        unsafeTabulate+          ( \[i, j] ->+              bool+                ( one+                    / unsafeIndex l [j, j]+                    * ( unsafeIndex a [i, j]+                          - sum+                            ( (\k -> unsafeIndex l [i, k] * unsafeIndex l [j, k])+                                <$> ([0 .. (j - 1)] :: [Int])+                            )+                      )+                )+                ( sqrt+                    ( unsafeIndex a [i, i]+                        - sum+                          ( (\k -> unsafeIndex l [j, k] ^ (2 :: Int))+                              <$> ([0 .. (j - 1)] :: [Int])+                          )+                    )+                )+                (i == j)+          )+   in l

+ test/doctests.hs view

@@ -0,0 +1,8 @@+module Main where++import System.Environment (getArgs)+import Test.DocTest (mainFromCabal)+import Prelude (IO, (=<<))++main :: IO ()+main = mainFromCabal "harpie-numhask" =<< getArgs