Skip to content

Commit

Permalink
v0.1.0.3
Browse files Browse the repository at this point in the history
  • Loading branch information
mstksg committed Mar 21, 2018
1 parent 52b9eca commit d4eca93
Show file tree
Hide file tree
Showing 3 changed files with 13 additions and 37 deletions.
3 changes: 2 additions & 1 deletion CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,8 @@ Version 0.1.0.3

* Compatibility with *base-4.11.0.0* and GHC 8.4
* Compatibility with *vector-sized-1.0.0.0*
* Internal conversion functions refactored using *hmatrix-vector-sized*.
* Internal conversion functions refactored using *hmatrix-vector-sized*,
*hessianF*.

Version 0.1.0.2
---------------
Expand Down
2 changes: 0 additions & 2 deletions package.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -35,8 +35,6 @@ library:
- Numeric.Hamilton
dependencies:
- ad
- comonad
- free
- hmatrix-gsl >=0.18
- hmatrix-vector-sized
- typelits-witnesses >=0.2.3
Expand Down
45 changes: 11 additions & 34 deletions src/Numeric/Hamilton.hs
Original file line number Diff line number Diff line change
Expand Up @@ -70,16 +70,14 @@ import Data.Foldable
import Data.Kind
import Data.Maybe
import Data.Proxy
import Data.Type.Equality hiding (sym)
import Data.Type.Equality
import GHC.Generics (Generic)
import GHC.TypeLits
import GHC.TypeLits.Compare
import Numeric.AD
import Numeric.GSL.ODE
import Numeric.LinearAlgebra.Static as H
import Numeric.LinearAlgebra.Static.Vector
import qualified Control.Comonad as C
import qualified Control.Comonad.Cofree as C
import qualified Data.Vector.Generic.Sized as VG
import qualified Data.Vector.Sized as V
import qualified Numeric.LinearAlgebra as LA
Expand Down Expand Up @@ -156,25 +154,12 @@ data System :: Nat -> Nat -> Type where
Sys :: { _sysInertia :: R m
, _sysCoords :: R n -> R m
, _sysJacobian :: R n -> L m n
, _sysJacobian2 :: R n -> V.Vector m (Sym n)
, _sysJacobian2 :: R n -> V.Vector m (L n n)
, _sysPotential :: R n -> Double
, _sysPotentialGrad :: R n -> R n
}
-> System m n

-- coordShift
-- :: (KnownNat m, KnownNat n, KnownNat o)
-- => (R o -> R n)
-- -> (R o -> L n o)
-- -> (R o -> V.Vector n (Sym o))
-- -> System m n
-- -> System m o
-- coordShift c j j2 = \case
-- Sys i c0 j0 j20 p g -> Sys i (c0 . c)
-- ((<>) <$> j0 . c <*> j)
-- ((\d -> fmap _) <$> j2 <*> j20 . c)
-- p g

-- | Converts the position of generalized coordinates of a system to the
-- coordinates of the system's underlying cartesian coordinate system.
-- Useful for plotting/drawing the system in cartesian space.
Expand Down Expand Up @@ -218,22 +203,14 @@ mkSystem
-- ^ The potential energy of the system as a function of
-- the generalized coordinate space's positions.
-> System m n
mkSystem m f u =
Sys m
(vecR . VG.convert . f . VG.convert . rVec)
(tr . vec2l . jacobianT f . VG.convert . rVec)
(fmap (sym . vec2l . j2 . C.hoistCofree VG.convert)
. VG.convert
. jacobians f
. VG.convert
. rVec
)
(u . VG.convert . rVec)
(vecR . VG.convert . grad u . VG.convert . rVec)
where
j2 :: C.Cofree (V.Vector n) Double
-> V.Vector n (V.Vector n Double)
j2 = fmap (fmap C.extract . C.unwrap) . C.unwrap
mkSystem m f u = Sys
{ _sysInertia = m
, _sysCoords = vecR . VG.convert . f . VG.convert . rVec
, _sysJacobian = tr . vec2l . jacobianT f . VG.convert . rVec
, _sysJacobian2 = fmap vec2l . hessianF f . VG.convert . rVec
, _sysPotential = u . VG.convert . rVec
, _sysPotentialGrad = vecR . VG.convert . grad u . VG.convert . rVec
}

-- | Convenience wrapper over 'mkSystem' that allows you to specify the
-- potential energy function in terms of the underlying cartesian
Expand Down Expand Up @@ -374,7 +351,7 @@ hamEqs Sys{..} Phs{..} = (dHdp, -dHdq)
mm = diag _sysInertia
j = _sysJacobian phsPositions
trj = tr j
j' = unSym <$> _sysJacobian2 phsPositions
j' = _sysJacobian2 phsPositions
jmj = trj H.<> mm H.<> j
ijmj = inv jmj
dTdq = vecR . VG.convert
Expand Down

0 comments on commit d4eca93

Please sign in to comment.