Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP

Comparing changes

Choose two branches to see what’s changed or to start a new pull request. If you need to, you can also compare across forks.

Open a pull request

Create a new pull request by comparing changes across two branches. If you need to, you can also compare across forks.
base fork: bjornbm/dimensional-experimental
base: 341abbcb9b
...
head fork: bjornbm/dimensional-experimental
compare: master
  • 17 commits
  • 7 files changed
  • 0 commit comments
  • 2 contributors
View
2  LICENSE
@@ -1,4 +1,4 @@
-Copyright (c) 2008, Bjorn Buckwalter.
+Copyright (c) 2008-2013, Bjorn Buckwalter.
All rights reserved.
Redistribution and use in source and binary forms, with or without
View
15 Numeric/Units/Dimensional/AD.hs
@@ -3,20 +3,23 @@
{-# LANGUAGE ExistentialQuantification #-}
{-# LANGUAGE Rank2Types #-}
-module Numeric.Units.Dimensional.AD (diff, Lift (lift)) where
+module Numeric.Units.Dimensional.AD (FAD, diff, Lift (lift)) where
import Numeric.Units.Dimensional (Dimensional (Dimensional), Quantity, Div)
-import Numeric.AD (AD, Mode)
-import qualified Numeric.AD (diff, lift)
+import Numeric.AD (AD, Mode, auto, Scalar)
+import qualified Numeric.AD (diff)
+import Numeric.AD.Mode.Forward (Forward)
-- | Unwrap a Dimensional's numeric representation.
undim :: Dimensional v d a -> a
undim (Dimensional a) = a
+type FAD tag a = AD tag (Forward a)
+
-- | @diff f x@ computes the derivative of the function @f(x)@ for the
-- given value of @x@.
diff :: (Num a, Div d2 d1 d3)
- => (forall tag. Mode tag => Quantity d1 (AD tag a) -> Quantity d2 (AD tag a))
+ => (forall tag. Quantity d1 (FAD tag a) -> Quantity d2 (FAD tag a))
-> Quantity d1 a -> Quantity d3 a
diff f = Dimensional . Numeric.AD.diff (undim . f . Dimensional) . undim
@@ -25,7 +28,7 @@ diff f = Dimensional . Numeric.AD.diff (undim . f . Dimensional) . undim
-- function).
class Lift w where
-- | Embed a constant data structure.
- lift :: (Num a, Mode t) => w a -> w (t a)
+ lift :: (Mode t) => w (Scalar t) -> w t
instance Lift (Dimensional v d)
- where lift = Dimensional . Numeric.AD.lift . undim
+ where lift = Dimensional . auto . undim
View
11 Numeric/Units/Dimensional/AEq.hs
@@ -1,4 +1,5 @@
{-# LANGUAGE TypeSynonymInstances #-}
+{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE GeneralizedNewtypeDeriving #-}
{-# LANGUAGE StandaloneDeriving #-}
@@ -7,12 +8,16 @@ module Numeric.Units.Dimensional.AEq where
import Numeric.Units.Dimensional (Dimensional (Dimensional), Quantity)
import Numeric.Units.Dimensional.LinearAlgebra.Vector (Vec (ListVec))
+import Numeric.Units.Dimensional.LinearAlgebra.Matrix (Mat (ListMat))
import Data.AEq
deriving instance AEq a => AEq (Quantity d a)
-instance (Floating a, AEq a) => AEq (Vec ds a) -- CPos et al
+instance (Floating a, AEq a) => AEq (Vec ds a)
where
- ListVec xs === ListVec ys = and $ zipWith (===) xs ys
- ListVec xs ~== ListVec ys = and $ zipWith (~==) xs ys
+ ListVec xs ~== ListVec ys = xs ~== ys
+
+instance (Floating a, AEq a) => AEq (Mat ds a)
+ where
+ ListMat xs ~== ListMat ys = xs ~== ys
View
87 Numeric/Units/Dimensional/LinearAlgebra/PosVel.lhs
@@ -7,15 +7,14 @@ The convention in this module is that a @C@ denotes cartesian coordinates and an
> module Numeric.Units.Dimensional.LinearAlgebra.PosVel where
> import qualified Prelude
+> import Numeric.NumType (Pos2)
> import Numeric.Units.Dimensional.Prelude
-> import Numeric.Units.Dimensional (Dimensional (Dimensional))
> import Numeric.Units.Dimensional.LinearAlgebra
> import Numeric.Units.Dimensional.LinearAlgebra.VectorAD (applyLinear)
Type synonyms for clearer documentation.
-> type DRadius = DLength; type Radius = Length
> -- | Angle from Z-axis.
> type DZenith = DPlaneAngle; type Zenith = Angle
> -- | Angle from X-axis towards Y-axis (positive about Z).
@@ -25,50 +24,48 @@ Type synonyms for clearer documentation.
Some type synonyms for convenience.
> type Vec3 d1 d2 d3 = Vec (d1:*:d2:*.d3)
-> type CPos = Vec3 DLength DLength DLength -- ^ x y z
-> type CVel = Vec3 DVelocity DVelocity DVelocity
-> type SPos = Vec3 DRadius DZenith DRightAscension
+> type Car d = Vec3 d d d -- ^ Cartesian vector (x,y,z).
+> type Sph d = Vec3 d DPlaneAngle DPlaneAngle -- ^ Spherical vector.
+> type CPos = Car DLength
+> type CVel = Car DVelocity
+> type SPos = Sph DRadius
> type SVel = Vec3 DVelocity DAngularVelocity DAngularVelocity
-Data type combining position and velocity into a state vector (minus epoch).
-
-> type CPosVel a = (CPos a, CVel a)
-> type SPosVel a = (SPos a, SVel a)
Querying
--------
-Cartesian position.
+Cartesian coordinates.
-> x :: CPos a -> Length a
+> x :: Car d a -> Quantity d a
> x = vElemAt zero
-> y :: CPos a -> Length a
+> y :: Car d a -> Quantity d a
> y = vElemAt pos1
-> z :: CPos a -> Length a
+> z :: Car d a -> Quantity d a
> z = vElemAt pos2
-Spherical position.
+Spherical coordinates.
-> radius :: SPos a -> Radius a
-> radius = vElemAt zero
+> magnitude :: Sph d a -> Quantity d a
+> magnitude = vElemAt zero
-> zenith :: SPos a -> Zenith a
+> zenith :: Sph d a -> Zenith a
> zenith = vElemAt pos1
-> colatitude :: SPos a -> Zenith a
+> colatitude :: Sph d a -> Zenith a
> colatitude = zenith
-> polarAngle :: SPos a -> Zenith a
+> polarAngle :: Sph d a -> Zenith a
> polarAngle = zenith
-> latitude :: Floating a => SPos a -> Angle a
-> latitude s = pi / _2 - colatitude s
-> declination :: Floating a => SPos a -> Angle a
+> latitude :: Floating a => Sph d a -> Angle a
+> latitude s = tau / _4 - colatitude s
+> declination :: Floating a => Sph d a -> Angle a
> declination = latitude
-> rightAscension :: SPos a -> RightAscension a
+> rightAscension :: Sph d a -> RightAscension a
> rightAscension = vElemAt pos2
-> longitude :: SPos a -> RightAscension a
+> longitude :: Sph d a -> RightAscension a
> longitude = rightAscension
-> hourAngle :: SPos a -> RightAscension a
+> hourAngle :: Sph d a -> RightAscension a
> hourAngle = rightAscension
@@ -76,20 +73,19 @@ Converting
----------
Converts a cartesian position vector into a spherical position vector.
-> c2s :: RealFloat a => CPos a -> SPos a
+> c2s :: (Div d d DOne, Mul d d d2, Root d2 Pos2 d, RealFloat a)
+> => Car d a -> Sph d a
> c2s c = fromTuple (r, zen, az)
> where
> (x, y, z) = toTuple c
> r = vNorm c -- sqrt (x^pos2 + y^pos2 + z^pos2)
-> zen = if r == 0 *~ meter then _0 else acos (z / r)
+> zen = if r == _0 then _0 else acos (z / r)
> az = atan2 y x
-> c2sEphem :: RealFloat a => CPosVel a -> SPosVel a
-> c2sEphem = applyLinear c2s -- unlinearize (c2s . linearize c :: RealFloat b => Time b -> SPos b)
Converts a spherical position vector into a cartesian position vector.
-> s2c :: Floating a => SPos a -> CPos a
+> s2c :: (Mul d DOne d) => Floating a => Sph d a -> Car d a
> s2c s = fromTuple (x, y, z)
> where
> (r, zen, az) = toTuple s
@@ -97,7 +93,34 @@ Converts a spherical position vector into a cartesian position vector.
> y = r * sin zen * sin az
> z = r * cos zen
-> s2cEphem :: RealFloat a => SPosVel a -> CPosVel a
-> s2cEphem = applyLinear s2c -- unlinearize (s2c . linearize s :: RealFloat b => Time b -> CPos b)
+Convert declination/latitude to zenith.
+
+> toZenith :: Floating a => Angle a -> Zenith a
+> toZenith decl = tau / _4 - decl
+
+Position vectors
+================
+> type DRadius = DLength; type Radius = Length
+
+> radius :: SPos a -> Radius a
+> radius = magnitude
+
+
+Position and velocity ephemeris
+===============================
+Data type combining position and velocity into a state vector (minus epoch).
+
+> type CPosVel a = (CPos a, CVel a)
+> type SPosVel a = (SPos a, SVel a)
+
+
+Conversions
+-----------
+
+> c2sEphem :: RealFloat a => CPosVel a -> SPosVel a
+> c2sEphem = applyLinear c2s -- unlinearize (c2s . linearize c :: RealFloat b => Time b -> SPos b)
+
+> s2cEphem :: RealFloat a => SPosVel a -> CPosVel a
+> s2cEphem = applyLinear s2c -- unlinearize (s2c . linearize s :: RealFloat b => Time b -> CPos b)
View
30 Numeric/Units/Dimensional/LinearAlgebra/Rotation.hs
@@ -28,16 +28,28 @@ type Homo33 d = Mat ((d:*:d:*.d) :*:
-- to rotating the coordinate system in opposite direction).
rotX :: Floating a => PlaneAngle a -> Homo33 DOne a
-rotX a = consRow (vCons _1 $ vCons _0 $ vSing _0)
- $ consRow (vCons _0 $ vCons (cos a) $ vSing (negate (sin a)))
- $ rowMatrix (vCons _0 $ vCons (sin a) $ vSing (cos a))
+rotX a = (_1 <: _0 <:. _0 )
+ |: (_0 <: cos a <:. negate (sin a))
+ |:. (_0 <: sin a <:. cos a )
rotY :: Floating a => PlaneAngle a -> Homo33 DOne a
-rotY a = consRow (vCons (cos a) $ vCons _0 $ vSing (sin a))
- $ consRow (vCons _0 $ vCons _1 $ vSing _0)
- $ rowMatrix (vCons (negate (sin a)) $ vCons _0 $ vSing (cos a))
+rotY a = ( cos a <: _0 <:. sin a)
+ |: ( _0 <: _1 <:. _0)
+ |:. (negate (sin a) <: _0 <:. cos a)
rotZ :: Floating a => PlaneAngle a -> Homo33 DOne a
-rotZ a = consRow (vCons (cos a) $ vCons (negate (sin a)) $ vSing _0)
- $ consRow (vCons (sin a) $ vCons (cos a) $ vSing _0)
- $ rowMatrix (vCons _0 $ vCons _0 $ vSing _1)
+rotZ a = (cos a <: negate (sin a) <:. _0)
+ |: (sin a <: cos a <:. _0)
+ |:. ( _0 <: _0 <:. _1)
+
+-- | @rotV v a@ generates the rotation matrix for a counterclockwise
+-- rotation through the angle @a@ about an axis in R^3, which is
+-- determined by the arbitrary unit vector @v@ that has its initial
+-- point at the origin (from Anton, Rorres, 1994, pages 190-191).
+rotV :: Floating a => Homo3 DOne a -> Angle a -> Homo33 DOne a
+rotV v t = (a*a*omct + cos t <: a*b*omct - c*sin t <:. a*c*omct + b*sin t)
+ |: (b*a*omct + c*sin t <: b*b*omct + cos t <:. b*c*omct - a*sin t)
+ |:. (c*a*omct - b*sin t <: c*b*omct + a*sin t <:. c*c*omct + cos t)
+ where
+ (a,b,c) = toTuple v
+ omct = _1 - cos t
View
15 Numeric/Units/Dimensional/LinearAlgebra/VectorAD.hs
@@ -8,8 +8,7 @@ import Numeric.Units.Dimensional.Prelude
import Numeric.Units.Dimensional (Dimensional (Dimensional))
import Numeric.Units.Dimensional.LinearAlgebra.Vector (Vec (ListVec), MulD, DivD, Homo, elemAdd, scaleVec)
import Numeric.Units.Dimensional.LinearAlgebra.HListExtras (HZipWith)
-import Numeric.AD (AD, diffF', Mode)
-import qualified Numeric.AD (lift)
+import Numeric.AD (auto, diffF')
import Numeric.Units.Dimensional.AD
@@ -17,7 +16,7 @@ import Numeric.Units.Dimensional.AD
-- @diff f@ is a function of the same type of quantity that returns
-- the first derivative of the result.
diffV :: (Num a, HMap (DivD,d) ds ds')
- => (forall tag. Mode tag => Quantity d (AD tag a) -> Vec ds (AD tag a)) -> Quantity d a -> Vec ds' a
+ => (forall tag. Quantity d (FAD tag a) -> Vec ds (FAD tag a)) -> Quantity d a -> Vec ds' a
diffV f = snd . diffV' f
@@ -25,7 +24,7 @@ unfzip as = (fmap fst as, fmap snd as)
-- | Like 'diffV' but returns a pair of the result and its first derivative.
diffV' :: (Num a, HMap (DivD,d) ds ds') -- Constraint could be changed to infer d instead (or also) if desired.
- => (forall tag. Mode tag => Quantity d (AD tag a) -> Vec ds (AD tag a))
+ => (forall tag. Quantity d (FAD tag a) -> Vec ds (FAD tag a))
-> Quantity d a -> (Vec ds a, Vec ds' a)
diffV' f (Dimensional x) = (ListVec ys, ListVec ys')
where
@@ -50,10 +49,10 @@ applyLinear :: forall a t ds ds' ds2 ds2' ts. (
HMap (MulD,t) ds' ds, -- Used in linearization.
HMap (DivD,t) ds2 ds2', -- Used in differentiation.
HZipWith DivD ds ds' ts, Homo ts t -- Necessary to infer t (the dimension w r t which we are differentiating).
- ) => (forall tag. Mode tag => Vec ds (AD tag a) -> Vec ds2 (AD tag a)) -> (Vec ds a, Vec ds' a) -> (Vec ds2 a, Vec ds2' a)
+ ) => (forall tag. Vec ds (FAD tag a) -> Vec ds2 (FAD tag a)) -> (Vec ds a, Vec ds' a) -> (Vec ds2 a, Vec ds2' a)
applyLinear f (p,v) = diffV' (\t -> f (lift p `elemAdd` scaleVec t (lift v))) t_0
where
- t_0 = Dimensional 0 :: Quantity t a
+ t_0 = _0 :: Quantity t a
-- 'applyLinearAt is analogous to 'applyLinear' but should be used when
-- the function is also dependent on the variable w r t which the vector
@@ -63,7 +62,7 @@ applyLinearAt :: forall a t ds ds' ds2 ds2' ts. (
HMap (MulD,t) ds' ds, -- Used in linearization.
HMap (DivD,t) ds2 ds2', -- Used in differentiation.
HZipWith DivD ds ds' ts, Homo ts t -- Necessary to infer t (the dimension w r t which we are differentiating).
- ) => (forall tag. Mode tag => Quantity t (AD tag a) -> Vec ds (AD tag a) -> Vec ds2 (AD tag a))
+ ) => (forall tag. Quantity t (FAD tag a) -> Vec ds (FAD tag a) -> Vec ds2 (FAD tag a))
-> Quantity t a -> (Vec ds a, Vec ds' a) -> (Vec ds2 a, Vec ds2' a)
applyLinearAt f t (p,v) = diffV' (\t' -> f t' (lift p `elemAdd` scaleVec (t' - lift t) (lift v))) t
@@ -72,4 +71,4 @@ applyLinearAt f t (p,v) = diffV' (\t' -> f t' (lift p `elemAdd` scaleVec (t' - l
-- -------
-- | Lift the elements of a vector to 'AD.AD's.
-instance Lift (Vec ds) where lift (ListVec xs) = ListVec (map Numeric.AD.lift xs)
+instance Lift (Vec ds) where lift (ListVec xs) = ListVec (map auto xs)
View
12 dimensional-experimental.cabal
@@ -1,8 +1,8 @@
Name: dimensional-experimental
-Version: 0.5
+Version: 0.10.0.0
License: BSD3
License-File: LICENSE
-Copyright: Bjorn Buckwalter 2008-2011
+Copyright: Bjorn Buckwalter 2008-2014
Author: Bjorn Buckwalter
Maintainer: bjorn.buckwalter@gmail.com
Stability: experimental
@@ -23,13 +23,13 @@ Build-Depends:
base,
time,
- dimensional,
+ dimensional >= 0.12,
dimensional-vectors >= 0.9,
numtype,
- HList >= 0.2.1,
+ HList-classic >= 0.2.1 && <1.1,
numeric-quest,
- ad >= 0.44.4,
- ieee754 >= 0.7.3
+ ad >= 4.2 && < 4.3,
+ ieee754 >= 0.7.4
Exposed-Modules:

No commit comments for this range

Something went wrong with that request. Please try again.