easytensor 0.3.0.0 → 0.3.1.0
raw patch · 44 files changed
+5088/−432 lines, 44 filesPVP: major bump suggested
API removals or changes: PVP suggests a major version bump
API changes (from Hackage documentation)
- Numeric.Commons: type DoubleBytes a = (PrimBytes a, ElemRep a ~ DoubleRep, ElemPrim a ~ Double#)
- Numeric.Commons: type FloatBytes a = (PrimBytes a, ElemRep a ~ FloatRep, ElemPrim a ~ Float#)
- Numeric.Commons: type IntBytes a = (PrimBytes a, ElemRep a ~ IntRep, ElemPrim a ~ Int#)
- Numeric.Commons: type WordBytes a = (PrimBytes a, ElemRep a ~ WordRep, ElemPrim a ~ Word#)
+ Numeric.DataFrame: AIDoubleX2 :: ArrayInstance t
+ Numeric.DataFrame: AIDoubleX3 :: ArrayInstance t
+ Numeric.DataFrame: AIDoubleX4 :: ArrayInstance t
+ Numeric.DataFrame.IO: SomeIODataFrame :: (IODataFrame t ns) -> SomeIODataFrame t
+ Numeric.DataFrame.IO: data SomeIODataFrame t (xns :: [XNat])
+ Numeric.DataFrame.ST: SomeSTDataFrame :: (STDataFrame s t ns) -> SomeSTDataFrame s t
+ Numeric.DataFrame.ST: data SomeSTDataFrame s t (xns :: [XNat])
+ Numeric.Matrix: class HomTransform4 t
+ Numeric.Matrix: fromHom :: HomTransform4 t => Vector t 4 -> Vector t 3
+ Numeric.Matrix: lookAt :: HomTransform4 t => Vector t 3 -> Vector t 3 -> Vector t 3 -> Matrix t 4 4
+ Numeric.Matrix: orthogonal :: HomTransform4 t => t -> t -> t -> t -> Matrix t 4 4
+ Numeric.Matrix: perspective :: HomTransform4 t => t -> t -> t -> t -> Matrix t 4 4
+ Numeric.Matrix: rotate :: HomTransform4 t => Vector t 3 -> t -> Matrix t 4 4
+ Numeric.Matrix: rotateEuler :: HomTransform4 t => t -> t -> t -> Matrix t 4 4
+ Numeric.Matrix: rotateX :: HomTransform4 t => t -> Matrix t 4 4
+ Numeric.Matrix: rotateY :: HomTransform4 t => t -> Matrix t 4 4
+ Numeric.Matrix: rotateZ :: HomTransform4 t => t -> Matrix t 4 4
+ Numeric.Matrix: toHomPoint :: HomTransform4 t => Vector t 3 -> Vector t 4
+ Numeric.Matrix: toHomVector :: HomTransform4 t => Vector t 3 -> Vector t 4
+ Numeric.Matrix: translate3 :: HomTransform4 t => Vector t 3 -> Matrix t 4 4
+ Numeric.Matrix: translate4 :: HomTransform4 t => Vector t 4 -> Matrix t 4 4
+ Numeric.Matrix: type Mat22d = Matrix Double 2 2
+ Numeric.Matrix: type Mat23d = Matrix Double 2 3
+ Numeric.Matrix: type Mat24d = Matrix Double 2 4
+ Numeric.Matrix: type Mat32d = Matrix Double 3 2
+ Numeric.Matrix: type Mat33d = Matrix Double 3 3
+ Numeric.Matrix: type Mat34d = Matrix Double 3 4
+ Numeric.Matrix: type Mat42d = Matrix Double 4 2
+ Numeric.Matrix: type Mat43d = Matrix Double 4 3
+ Numeric.Matrix: type Mat44d = Matrix Double 4 4
+ Numeric.Quaternion: axisRotation :: Quaternion t => Vector t 3 -> t -> Quater t
+ Numeric.Quaternion: class Quaternion t where data Quater t where {
+ Numeric.Quaternion: conjugate :: Quaternion t => Quater t -> Quater t
+ Numeric.Quaternion: data family Quater t;
+ Numeric.Quaternion: fromMatrix33 :: Quaternion t => Matrix t 3 3 -> Quater t
+ Numeric.Quaternion: fromMatrix44 :: Quaternion t => Matrix t 4 4 -> Quater t
+ Numeric.Quaternion: fromVec4 :: Quaternion t => Vector t 4 -> Quater t
+ Numeric.Quaternion: fromVecNum :: Quaternion t => Vector t 3 -> t -> Quater t
+ Numeric.Quaternion: getRotScale :: Quaternion t => Vector t 3 -> Vector t 3 -> Quater t
+ Numeric.Quaternion: im :: Quaternion t => Quater t -> Quater t
+ Numeric.Quaternion: imVec :: Quaternion t => Quater t -> Vector t 3
+ Numeric.Quaternion: packQ :: Quaternion t => t -> t -> t -> t -> Quater t
+ Numeric.Quaternion: qArg :: Quaternion t => Quater t -> t
+ Numeric.Quaternion: re :: Quaternion t => Quater t -> Quater t
+ Numeric.Quaternion: rotScale :: Quaternion t => Quater t -> Vector t 3 -> Vector t 3
+ Numeric.Quaternion: square :: Quaternion t => Quater t -> t
+ Numeric.Quaternion: takei :: Quaternion t => Quater t -> t
+ Numeric.Quaternion: takej :: Quaternion t => Quater t -> t
+ Numeric.Quaternion: takek :: Quaternion t => Quater t -> t
+ Numeric.Quaternion: taker :: Quaternion t => Quater t -> t
+ Numeric.Quaternion: toMatrix33 :: Quaternion t => Quater t -> Matrix t 3 3
+ Numeric.Quaternion: toMatrix44 :: Quaternion t => Quater t -> Matrix t 4 4
+ Numeric.Quaternion: toVec4 :: Quaternion t => Quater t -> Vector t 4
+ Numeric.Quaternion: type QDouble = Quater Double
+ Numeric.Quaternion: type QFloat = Quater Float
+ Numeric.Quaternion: unpackQ :: Quaternion t => Quater t -> (t, t, t, t)
+ Numeric.Quaternion: }
+ Numeric.Vector: normalized :: (Floating t, Fractional (Vector t n), ElementWise (Idx '[n]) t (Vector t n)) => Vector t n -> Vector t n
+ Numeric.Vector: unpackV2 :: ElementWise (Idx '[2]) t (Vector t 2) => Vector t 2 -> (t, t)
+ Numeric.Vector: unpackV3 :: ElementWise (Idx '[3]) t (Vector t 3) => Vector t 3 -> (t, t, t)
+ Numeric.Vector: unpackV4 :: ElementWise (Idx '[4]) t (Vector t 4) => Vector t 4 -> (t, t, t, t)
- Numeric.DataFrame.IO: copyDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]). (PrimBytes (DataFrame t as), ConcatList as bs asbs, Dimensions bs) => DataFrame t as -> Idx bs -> IODataFrame t asbs -> IO ()
+ Numeric.DataFrame.IO: copyDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]). (ConcatList as (b :+ bs) asbs, Dimensions (b :+ bs), PrimBytes (DataFrame t (as +: b'))) => DataFrame t (as +: b') -> Idx (b :+ bs) -> IODataFrame t asbs -> IO ()
- Numeric.DataFrame.IO: copyMutableDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]). (PrimBytes t, ConcatList as bs asbs, Dimensions bs) => IODataFrame t as -> Idx bs -> IODataFrame t asbs -> IO ()
+ Numeric.DataFrame.IO: copyMutableDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]). (PrimBytes t, ConcatList as (b :+ bs) asbs, Dimensions (b :+ bs)) => IODataFrame t (as +: b') -> Idx (b :+ bs) -> IODataFrame t asbs -> IO ()
- Numeric.DataFrame.ST: copyDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) s. (PrimBytes (DataFrame t as), ConcatList as bs asbs, Dimensions bs) => DataFrame t as -> Idx bs -> STDataFrame s t asbs -> ST s ()
+ Numeric.DataFrame.ST: copyDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s. (ConcatList as (b :+ bs) asbs, Dimensions (b :+ bs), PrimBytes (DataFrame t (as +: b'))) => DataFrame t (as +: b') -> Idx (b :+ bs) -> STDataFrame s t asbs -> ST s ()
- Numeric.DataFrame.ST: copyMutableDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) s. (PrimBytes t, ConcatList as bs asbs, Dimensions bs) => STDataFrame s t as -> Idx bs -> STDataFrame s t asbs -> ST s ()
+ Numeric.DataFrame.ST: copyMutableDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s. (PrimBytes t, ConcatList as (b :+ bs) asbs, Dimensions (b :+ bs)) => STDataFrame s t (as +: b') -> Idx (b :+ bs) -> STDataFrame s t asbs -> ST s ()
Files
- bench/misc.hs +20/−0
- easytensor.cabal +18/−6
- src-base/Numeric/Array.hs +4/−0
- src-base/Numeric/Array/Family.hs +62/−7
- src-base/Numeric/Array/Family/ArrayD.hs +1/−1
- src-base/Numeric/Array/Family/ArrayF.hs +1/−1
- src-base/Numeric/Array/Family/ArrayI64.hs +43/−12
- src-base/Numeric/Array/Family/ArrayW64.hs +46/−8
- src-base/Numeric/Array/Family/DoubleX2.hs +262/−0
- src-base/Numeric/Array/Family/DoubleX3.hs +299/−0
- src-base/Numeric/Array/Family/DoubleX4.hs +335/−0
- src-base/Numeric/DataFrame/Contraction.hs +83/−16
- src-base/Numeric/DataFrame/Inference.hs +12/−2
- src-base/Numeric/DataFrame/Mutable.hs +224/−0
- src-base/Numeric/Matrix/Mat44d.hs +37/−0
- src-base/Numeric/Matrix/Mat44f.hs +37/−0
- src-base/Numeric/Quaternion/QDouble.hs +542/−0
- src-base/Numeric/Quaternion/QFloat.hs +562/−0
- src-ghcjs/Numeric/Array/Family.hs +7/−3
- src-ghcjs/Numeric/Array/Family/ArrayT.hs +5/−2
- src-ghcjs/Numeric/Array/Family/ArrayT.js +64/−20
- src-ghcjs/Numeric/DataFrame/Contraction.hs +20/−5
- src-ghcjs/Numeric/DataFrame/Mutable.hs +296/−0
- src-ghcjs/Numeric/Matrix/Mat44.js +103/−0
- src-ghcjs/Numeric/Matrix/Mat44d.hs +69/−0
- src-ghcjs/Numeric/Matrix/Mat44f.hs +69/−0
- src-ghcjs/Numeric/Quaternion/QDouble.hs +285/−0
- src-ghcjs/Numeric/Quaternion/QFloat.hs +295/−0
- src-ghcjs/Numeric/Quaternion/Quaternion.js +205/−0
- src/Numeric/Commons.hs +106/−11
- src/Numeric/DataFrame.hs +2/−0
- src/Numeric/DataFrame/IO.hs +241/−17
- src/Numeric/DataFrame/Mutable.hs +0/−223
- src/Numeric/DataFrame/ST.hs +238/−14
- src/Numeric/DataFrame/SubSpace.hs +19/−9
- src/Numeric/DataFrame/Type.hs +30/−3
- src/Numeric/Matrix.hs +34/−22
- src/Numeric/Matrix/Class.hs +114/−0
- src/Numeric/Matrix/Type.hs +0/−49
- src/Numeric/Quaternion.hs +35/−0
- src/Numeric/Quaternion/Class.hs +85/−0
- src/Numeric/Vector.hs +30/−0
- test/Numeric/QuaternionTest.hs +124/−0
- test/Spec.hs +24/−1
bench/misc.hs view
@@ -9,8 +9,12 @@ import Numeric.DataFrame import Numeric.Dimensions +import qualified Control.Monad.ST as ST+import qualified Numeric.DataFrame.ST as ST+-- import qualified Numeric.Dimensions.Traverse.ST as ST + main :: IO () main = do putStrLn "Hello world!"@@ -75,3 +79,19 @@ -- Updating existing frames print $ update (2:!Z) (scalar 777) rVec print $ update (2:!3:!Z) (vec2 999 999) x++ let matX = iwgen (scalar . fromEnum) :: DataFrame Int '[2,5,4]+ matY = iwgen (scalar . fromEnum) :: DataFrame Int '[5,4]+ putStrLn "Check carefully that this returns no garbage"+ print matX+ print (ewmap (<+:> scalar 111) matX :: DataFrame Int '[3,5,4])+ print matY+ print (ewmap fromScalar matY :: DataFrame Int '[3,5,4])++ -- Working with mutable frames+ print $ ST.runST $ do+ sdf <- ST.thawDataFrame matY+ ST.writeDataFrame sdf (1:!1:!Z) 900101+ ST.writeDataFrame sdf (3:!3:!Z) 900303+ ST.writeDataFrame sdf (5:!3:!Z) 900503+ ST.unsafeFreezeDataFrame sdf
easytensor.cabal view
@@ -1,8 +1,8 @@ name: easytensor-version: 0.3.0.0+version: 0.3.1.0 cabal-version: >=1.20 build-type: Simple-license: MIT+license: BSD3 license-file: LICENSE copyright: (c) Artem Chirkin maintainer: chirkin@arch.ethz.ch@@ -42,6 +42,7 @@ Numeric.Matrix Numeric.Vector Numeric.Scalar+ Numeric.Quaternion build-depends: base >=4.9 && <5, ghc-prim >=0.5,@@ -65,7 +66,12 @@ Numeric.DataFrame.Type Numeric.DataFrame.Inference Numeric.DataFrame.Shape- Numeric.Matrix.Type+ Numeric.Matrix.Class+ Numeric.Matrix.Mat44d+ Numeric.Matrix.Mat44f+ Numeric.Quaternion.Class+ Numeric.Quaternion.QDouble+ Numeric.Quaternion.QFloat if impl(ghcjs) other-modules: Numeric.Array.Family.ArrayT@@ -86,18 +92,24 @@ Numeric.Array.Family.FloatX2 Numeric.Array.Family.FloatX3 Numeric.Array.Family.FloatX4+ Numeric.Array.Family.DoubleX2+ Numeric.Array.Family.DoubleX3+ Numeric.Array.Family.DoubleX4 js-sources: src-ghcjs/Numeric/Array/Family/ArrayT.js- ghc-options: -Wall -fwarn-tabs -O2+ src-ghcjs/Numeric/Quaternion/Quaternion.js+ src-ghcjs/Numeric/Matrix/Mat44.js+ ghc-options: -Wall -fwarn-tabs -fwarn-unused-do-bind -fwarn-monomorphism-restriction -O2 test-suite et-test - type: detailed-0.9- test-module: Spec+ type: exitcode-stdio-1.0+ main-is: Spec.hs other-modules: Numeric.DataFrame.Arbitraries Numeric.DataFrame.SubSpaceTest Numeric.DataFrame.BasicTest+ Numeric.QuaternionTest build-depends: base -any, Cabal >=1.20,
src-base/Numeric/Array.hs view
@@ -31,3 +31,7 @@ import Numeric.Array.Family.FloatX2 () import Numeric.Array.Family.FloatX3 () import Numeric.Array.Family.FloatX4 ()++import Numeric.Array.Family.DoubleX2 ()+import Numeric.Array.Family.DoubleX3 ()+import Numeric.Array.Family.DoubleX4 ()
src-base/Numeric/Array/Family.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE ConstraintKinds #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE ExistentialQuantification #-}@@ -31,17 +32,23 @@ , ArrayW (..), ArrayW8 (..), ArrayW16 (..), ArrayW32 (..), ArrayW64 (..) , Scalar (..) , FloatX2 (..), FloatX3 (..), FloatX4 (..)+ , DoubleX2 (..), DoubleX3 (..), DoubleX4 (..) , ArrayInstanceInference, ElemType (..), ArraySize (..) , ElemTypeInference (..), ArraySizeInference (..), ArrayInstanceEvidence , getArrayInstance, ArrayInstance (..), inferArrayInstance ) where +#include "MachDeps.h" import Data.Int (Int16, Int32, Int64, Int8) import Data.Type.Equality ((:~:) (..)) import Data.Word (Word16, Word32, Word64, Word8)-import GHC.Prim (ByteArray#, Double#, Float#, Int#,- Word#, unsafeCoerce#)+import GHC.Prim ( ByteArray#, Double#, Float#+#if WORD_SIZE_IN_BITS < 64+ , Int64#, Word64#+#endif+ , Int#, Word#, unsafeCoerce#)+import GHC.Exts (RuntimeRep(..)) import Numeric.Array.ElementWise import Numeric.Commons@@ -54,6 +61,9 @@ Array Float '[2] = FloatX2 Array Float '[3] = FloatX3 Array Float '[4] = FloatX4+ Array Double '[2] = DoubleX2+ Array Double '[3] = DoubleX3+ Array Double '[4] = DoubleX4 Array Float (d ': ds) = ArrayF (d ': ds) Array Double (d ': ds) = ArrayD (d ': ds) Array Int (d ': ds) = ArrayI (d ': ds)@@ -75,19 +85,47 @@ instance Show t => Show (Scalar t) where show (Scalar t) = "{ " ++ show t ++ " }" -type instance ElemRep (Scalar t) = ElemRep t+type instance ElemRep (Scalar Float ) = 'FloatRep+type instance ElemRep (Scalar Double) = 'DoubleRep+type instance ElemRep (Scalar Int ) = 'IntRep+type instance ElemRep (Scalar Int8 ) = 'IntRep+type instance ElemRep (Scalar Int16 ) = 'IntRep+type instance ElemRep (Scalar Int32 ) = 'IntRep+#if WORD_SIZE_IN_BITS < 64+type instance ElemRep (Scalar Int64 ) = 'Int64Rep+#else+type instance ElemRep (Scalar Int64 ) = 'IntRep+#endif+type instance ElemRep (Scalar Word ) = 'WordRep+type instance ElemRep (Scalar Word8 ) = 'WordRep+type instance ElemRep (Scalar Word16) = 'WordRep+type instance ElemRep (Scalar Word32) = 'WordRep+#if WORD_SIZE_IN_BITS < 64+type instance ElemRep (Scalar Word64) = 'Word64Rep+#else+type instance ElemRep (Scalar Word64) = 'WordRep+#endif+ type instance ElemPrim (Scalar Float ) = Float# type instance ElemPrim (Scalar Double) = Double# type instance ElemPrim (Scalar Int ) = Int# type instance ElemPrim (Scalar Int8 ) = Int# type instance ElemPrim (Scalar Int16 ) = Int# type instance ElemPrim (Scalar Int32 ) = Int#+#if WORD_SIZE_IN_BITS < 64+type instance ElemPrim (Scalar Int64 ) = Int64#+#else type instance ElemPrim (Scalar Int64 ) = Int#+#endif type instance ElemPrim (Scalar Word ) = Word# type instance ElemPrim (Scalar Word8 ) = Word# type instance ElemPrim (Scalar Word16) = Word# type instance ElemPrim (Scalar Word32) = Word#+#if WORD_SIZE_IN_BITS < 64+type instance ElemPrim (Scalar Word64) = Word64#+#else type instance ElemPrim (Scalar Word64) = Word#+#endif deriving instance PrimBytes (Scalar Float) deriving instance PrimBytes (Scalar Double)@@ -146,8 +184,13 @@ | FromScalarI16# Int# data ArrayI32 (ds :: [Nat]) = ArrayI32# Int# Int# ByteArray# | FromScalarI32# Int#+#if WORD_SIZE_IN_BITS < 64 data ArrayI64 (ds :: [Nat]) = ArrayI64# Int# Int# ByteArray#+ | FromScalarI64# Int64#+#else+data ArrayI64 (ds :: [Nat]) = ArrayI64# Int# Int# ByteArray# | FromScalarI64# Int#+#endif data ArrayW (ds :: [Nat]) = ArrayW# Int# Int# ByteArray# | FromScalarW# Word# data ArrayW8 (ds :: [Nat]) = ArrayW8# Int# Int# ByteArray#@@ -156,8 +199,13 @@ | FromScalarW16# Word# data ArrayW32 (ds :: [Nat]) = ArrayW32# Int# Int# ByteArray# | FromScalarW32# Word#+#if WORD_SIZE_IN_BITS < 64 data ArrayW64 (ds :: [Nat]) = ArrayW64# Int# Int# ByteArray#+ | FromScalarW64# Word64#+#else+data ArrayW64 (ds :: [Nat]) = ArrayW64# Int# Int# ByteArray# | FromScalarW64# Word#+#endif -- * Specialized types -- More efficient data types for small fixed-size tensors@@ -165,6 +213,10 @@ data FloatX3 = FloatX3# Float# Float# Float# data FloatX4 = FloatX4# Float# Float# Float# Float# +data DoubleX2 = DoubleX2# Double# Double#+data DoubleX3 = DoubleX3# Double# Double# Double#+data DoubleX4 = DoubleX4# Double# Double# Double# Double#+ -- * Recovering type instances at runtime -- A combination of `ElemType t` and `ArraySize ds` should -- define an instance of `Array t ds` unambiguously.@@ -218,6 +270,9 @@ | ( Array t ds ~ FloatX2, ds ~ '[2], t ~ Float) => AIFloatX2 | ( Array t ds ~ FloatX3, ds ~ '[3], t ~ Float) => AIFloatX3 | ( Array t ds ~ FloatX4, ds ~ '[4], t ~ Float) => AIFloatX4+ | ( Array t ds ~ DoubleX2, ds ~ '[2], t ~ Double) => AIDoubleX2+ | ( Array t ds ~ DoubleX3, ds ~ '[3], t ~ Double) => AIDoubleX3+ | ( Array t ds ~ DoubleX4, ds ~ '[4], t ~ Double) => AIDoubleX4 -- | A singleton type used to prove that the given Array family instance -- has a known instance@@ -346,7 +401,7 @@ (ETWord64 , ASScalar) -> AIScalar (ETFloat , ASX2) -> AIFloatX2- (ETDouble , ASX2) -> AIArrayD+ (ETDouble , ASX2) -> AIDoubleX2 (ETInt , ASX2) -> AIArrayI (ETInt8 , ASX2) -> AIArrayI8 (ETInt16 , ASX2) -> AIArrayI16@@ -359,7 +414,7 @@ (ETWord64 , ASX2) -> AIArrayW64 (ETFloat , ASX3) -> AIFloatX3- (ETDouble , ASX3) -> AIArrayD+ (ETDouble , ASX3) -> AIDoubleX3 (ETInt , ASX3) -> AIArrayI (ETInt8 , ASX3) -> AIArrayI8 (ETInt16 , ASX3) -> AIArrayI16@@ -372,7 +427,7 @@ (ETWord64 , ASX3) -> AIArrayW64 (ETFloat , ASX4) -> AIFloatX4- (ETDouble , ASX4) -> AIArrayD+ (ETDouble , ASX4) -> AIDoubleX4 (ETInt , ASX4) -> AIArrayI (ETInt8 , ASX4) -> AIArrayI8 (ETInt16 , ASX4) -> AIArrayI16@@ -385,7 +440,7 @@ (ETWord64 , ASX4) -> AIArrayW64 (ETFloat , ASXN) -> unsafeCoerce# (AIArrayF :: ArrayInstance Float '[5])- (ETDouble , ASXN) -> AIArrayD+ (ETDouble , ASXN) -> unsafeCoerce# (AIArrayD :: ArrayInstance Double '[5]) (ETInt , ASXN) -> AIArrayI (ETInt8 , ASXN) -> AIArrayI8 (ETInt16 , ASXN) -> AIArrayI16
src-base/Numeric/Array/Family/ArrayD.hs view
@@ -40,7 +40,7 @@ import Numeric.Dimensions import Numeric.Dimensions.Traverse import Numeric.TypeLits-import Numeric.Matrix.Type+import Numeric.Matrix.Class #include "MachDeps.h"
src-base/Numeric/Array/Family/ArrayF.hs view
@@ -42,7 +42,7 @@ import Numeric.Dimensions import Numeric.Dimensions.Traverse import Numeric.TypeLits-import Numeric.Matrix.Type+import Numeric.Matrix.Class #include "MachDeps.h" #define ARR_TYPE ArrayF
src-base/Numeric/Array/Family/ArrayI64.hs view
@@ -27,10 +27,14 @@ module Numeric.Array.Family.ArrayI64 () where +#include "MachDeps.h" import GHC.Base (runRW#) import GHC.Prim import GHC.Types (Int (..), RuntimeRep (..), isTrue#) import GHC.Int (Int64 (..))+#if WORD_SIZE_IN_BITS < 64+import GHC.IntWord64+#endif import Numeric.Array.ElementWise import Numeric.Array.Family@@ -39,11 +43,37 @@ import Numeric.Dimensions.Traverse -#include "MachDeps.h"++#if SIZEOF_HSWORD < 8 #define ARR_TYPE ArrayI64 #define ARR_FROMSCALAR FromScalarI64# #define ARR_CONSTR ArrayI64# #define EL_TYPE_BOXED Int64+#define EL_TYPE_PRIM Int64#+#define EL_RUNTIME_REP 'Int64Rep+#define EL_CONSTR I64#+#define EL_SIZE SIZEOF_INT64#+#define EL_ALIGNMENT ALIGNMENT_INT64#+#define EL_ZERO (intToInt64# 0#)+#define EL_ONE (intToInt64# 1#)+#define EL_MINUS_ONE (intToInt64# (-1#))+#define INDEX_ARRAY indexInt64Array#+#define WRITE_ARRAY writeInt64Array#+#define OP_EQ (eqInt64#)+#define OP_NE (neInt64#)+#define OP_GT (gtInt64#)+#define OP_GE (geInt64#)+#define OP_LT (ltInt64#)+#define OP_LE (leInt64#)+#define OP_PLUS (plusInt64#)+#define OP_MINUS (minusInt64#)+#define OP_TIMES (timesInt64#)+#define OP_NEGATE negateInt64#+#else+#define ARR_TYPE ArrayI64+#define ARR_FROMSCALAR FromScalarI64#+#define ARR_CONSTR ArrayI64#+#define EL_TYPE_BOXED Int64 #define EL_TYPE_PRIM Int# #define EL_RUNTIME_REP 'IntRep #define EL_CONSTR I64#@@ -64,28 +94,29 @@ #define OP_MINUS (-#) #define OP_TIMES (*#) #define OP_NEGATE negateInt#+#endif #include "Array.h" instance Num (ArrayI64 ds) where- (+) = zipV (+#)+ (+) = zipV OP_PLUS {-# INLINE (+) #-}- (-) = zipV (-#)+ (-) = zipV OP_MINUS {-# INLINE (-) #-}- (*) = zipV (*#)+ (*) = zipV OP_TIMES {-# INLINE (*) #-}- negate = mapV negateInt#+ negate = mapV OP_NEGATE {-# INLINE negate #-}- abs = mapV (\x -> if isTrue# (x >=# 0#)+ abs = mapV (\x -> if isTrue# (OP_GE x EL_ZERO) then x- else negateInt# x+ else OP_NEGATE x ) {-# INLINE abs #-}- signum = mapV (\x -> if isTrue# (x ># 0#)- then 1#- else if isTrue# (x <# 0#)- then -1#- else 0#+ signum = mapV (\x -> if isTrue# (OP_GT x EL_ZERO)+ then EL_ONE+ else if isTrue# (OP_LT x EL_ZERO)+ then EL_MINUS_ONE+ else EL_ZERO ) {-# INLINE signum #-} fromInteger = broadcastArray . fromInteger
src-base/Numeric/Array/Family/ArrayW64.hs view
@@ -27,10 +27,14 @@ module Numeric.Array.Family.ArrayW64 () where +#include "MachDeps.h" import GHC.Base (runRW#) import GHC.Prim import GHC.Types (Int (..), RuntimeRep (..), isTrue#) import GHC.Word (Word64 (..))+#if WORD_SIZE_IN_BITS < 64+import GHC.IntWord64+#endif import Numeric.Array.ElementWise import Numeric.Array.Family@@ -38,12 +42,42 @@ import Numeric.Dimensions import Numeric.Dimensions.Traverse +#if SIZEOF_HSWORD < 8 -#include "MachDeps.h"+plusWord64#, minusWord64#, timesWord64#+ :: Word64# -> Word64# -> Word64#+plusWord64# x y = int64ToWord64# (word64ToInt64# x `plusInt64#` word64ToInt64# y)+minusWord64# x y = int64ToWord64# (word64ToInt64# x `minusInt64#` word64ToInt64# y)+timesWord64# x y = int64ToWord64# (word64ToInt64# x `timesInt64#` word64ToInt64# y)++ #define ARR_TYPE ArrayW64 #define ARR_FROMSCALAR FromScalarW64# #define ARR_CONSTR ArrayW64# #define EL_TYPE_BOXED Word64+#define EL_TYPE_PRIM Word64#+#define EL_RUNTIME_REP 'Word64Rep+#define EL_CONSTR W64#+#define EL_SIZE SIZEOF_WORD64#+#define EL_ALIGNMENT ALIGNMENT_WORD64#+#define EL_ZERO (wordToWord64# 0##)+#define EL_ONE (wordToWord64# 1##)+#define INDEX_ARRAY indexWord64Array#+#define WRITE_ARRAY writeWord64Array#+#define OP_EQ (eqWord64#)+#define OP_NE (neWord64#)+#define OP_GT (gtWord64#)+#define OP_GE (geWord64#)+#define OP_LT (ltWord64#)+#define OP_LE (leWord64#)+#define OP_PLUS (plusWord64#)+#define OP_MINUS (minusWord64#)+#define OP_TIMES (timesWord64#)+#else+#define ARR_TYPE ArrayW64+#define ARR_FROMSCALAR FromScalarW64#+#define ARR_CONSTR ArrayW64#+#define EL_TYPE_BOXED Word64 #define EL_TYPE_PRIM Word# #define EL_RUNTIME_REP 'WordRep #define EL_CONSTR W64#@@ -51,7 +85,6 @@ #define EL_ALIGNMENT ALIGNMENT_WORD64# #define EL_ZERO 0## #define EL_ONE 1##-#define EL_MINUS_ONE -1# #define INDEX_ARRAY indexWord64Array# #define WRITE_ARRAY writeWord64Array# #define OP_EQ eqWord#@@ -63,22 +96,27 @@ #define OP_PLUS plusWord# #define OP_MINUS minusWord# #define OP_TIMES timesWord#+#endif #include "Array.h" instance Num (ArrayW64 ds) where- (+) = zipV plusWord#+ (+) = zipV OP_PLUS {-# INLINE (+) #-}- (-) = zipV minusWord#+ (-) = zipV OP_MINUS {-# INLINE (-) #-}- (*) = zipV timesWord#+ (*) = zipV OP_TIMES {-# INLINE (*) #-}+#if SIZEOF_HSWORD < 8+ negate = mapV (\x -> int64ToWord64# (negateInt64# (word64ToInt64# x)))+#else negate = mapV (\x -> int2Word# (negateInt# (word2Int# x)))+#endif {-# INLINE negate #-} abs = id {-# INLINE abs #-}- signum = mapV (\x -> if isTrue# (gtWord# x 0##)- then 1##- else 0##+ signum = mapV (\x -> if isTrue# (OP_GT x EL_ZERO)+ then EL_ONE+ else EL_ZERO ) {-# INLINE signum #-} fromInteger = broadcastArray . fromInteger
+ src-base/Numeric/Array/Family/DoubleX2.hs view
@@ -0,0 +1,262 @@+{-# LANGUAGE CPP #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UnboxedTuples #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+-----------------------------------------------------------------------------+-- |+-- Module : Numeric.Array.Family.DoubleX2+-- Copyright : (c) Artem Chirkin+-- License : BSD3+--+-- Maintainer : chirkin@arch.ethz.ch+--+--+-----------------------------------------------------------------------------++module Numeric.Array.Family.DoubleX2 () where+++#include "MachDeps.h"++import GHC.Base (runRW#)+import GHC.Prim+import GHC.Types (Double (..), RuntimeRep (..),+ isTrue#)++import Numeric.Array.ElementWise+import Numeric.Array.Family+import Numeric.Commons+import Numeric.Dimensions++++++instance Show DoubleX2 where+ show (DoubleX2# a1 a2) = "{ " ++ show (D# a1)+ ++ ", " ++ show (D# a2)+ ++ " }"++++instance Eq DoubleX2 where+ DoubleX2# a1 a2 == DoubleX2# b1 b2 = isTrue# ( (a1 ==## b1)+ `andI#` (a2 ==## b2)+ )+ {-# INLINE (==) #-}+ DoubleX2# a1 a2 /= DoubleX2# b1 b2 = isTrue# ( (a1 /=## b1)+ `orI#` (a2 /=## b2)+ )+ {-# INLINE (/=) #-}++++-- | Implement partial ordering for `>`, `<`, `>=`, `<=`+-- and lexicographical ordering for `compare`+instance Ord DoubleX2 where+ DoubleX2# a1 a2 > DoubleX2# b1 b2 = isTrue# ( (a1 >## b1)+ `andI#` (a2 >## b2)+ )+ {-# INLINE (>) #-}+ DoubleX2# a1 a2 < DoubleX2# b1 b2 = isTrue# ( (a1 <## b1)+ `andI#` (a2 <## b2)+ )+ {-# INLINE (<) #-}+ DoubleX2# a1 a2 >= DoubleX2# b1 b2 = isTrue# ( (a1 >=## b1)+ `andI#` (a2 >=## b2)+ )+ {-# INLINE (>=) #-}+ DoubleX2# a1 a2 <= DoubleX2# b1 b2 = isTrue# ( (a1 <=## b1)+ `andI#` (a2 <=## b2)+ )+ {-# INLINE (<=) #-}+ -- | Compare lexicographically+ compare (DoubleX2# a1 a2) (DoubleX2# b1 b2)+ | isTrue# (a1 >## b1) = GT+ | isTrue# (a1 <## b1) = LT+ | isTrue# (a2 >## b2) = GT+ | isTrue# (a2 <## b2) = LT+ | otherwise = EQ+ {-# INLINE compare #-}+ -- | Element-wise minimum+ min (DoubleX2# a1 a2) (DoubleX2# b1 b2) =+ DoubleX2# (if isTrue# (a1 >## b1) then b1 else a1)+ (if isTrue# (a2 >## b2) then b2 else a2)+ {-# INLINE min #-}+ -- | Element-wise maximum+ max (DoubleX2# a1 a2) (DoubleX2# b1 b2) =+ DoubleX2# (if isTrue# (a1 >## b1) then a1 else b1)+ (if isTrue# (a2 >## b2) then a2 else b2)+ {-# INLINE max #-}++++-- | element-wise operations for vectors+instance Num DoubleX2 where+ DoubleX2# a1 a2 + DoubleX2# b1 b2+ = DoubleX2# ((+##) a1 b1) ((+##) a2 b2)+ {-# INLINE (+) #-}+ DoubleX2# a1 a2 - DoubleX2# b1 b2+ = DoubleX2# ((-##) a1 b1) ((-##) a2 b2)+ {-# INLINE (-) #-}+ DoubleX2# a1 a2 * DoubleX2# b1 b2+ = DoubleX2# ((*##) a1 b1) ((*##) a2 b2)+ {-# INLINE (*) #-}+ negate (DoubleX2# a1 a2)+ = DoubleX2# (negateDouble# a1) (negateDouble# a2)+ {-# INLINE negate #-}+ abs (DoubleX2# a1 a2)+ = DoubleX2# (if isTrue# (a1 >=## 0.0##) then a1 else negateDouble# a1)+ (if isTrue# (a2 >=## 0.0##) then a2 else negateDouble# a2)+ {-# INLINE abs #-}+ signum (DoubleX2# a1 a2)+ = DoubleX2# (if isTrue# (a1 >## 0.0##)+ then 1.0##+ else if isTrue# (a1 <## 0.0##) then -1.0## else 0.0## )+ (if isTrue# (a2 >## 0.0##)+ then 1.0##+ else if isTrue# (a2 <## 0.0##) then -1.0## else 0.0## )+ {-# INLINE signum #-}+ fromInteger n = case fromInteger n of D# x -> DoubleX2# x x+ {-# INLINE fromInteger #-}++++instance Fractional DoubleX2 where+ DoubleX2# a1 a2 / DoubleX2# b1 b2 = DoubleX2# ((/##) a1 b1)+ ((/##) a2 b2)+ {-# INLINE (/) #-}+ recip (DoubleX2# a1 a2) = DoubleX2# ((/##) 1.0## a1)+ ((/##) 1.0## a2)+ {-# INLINE recip #-}+ fromRational r = case fromRational r of D# x -> DoubleX2# x x+ {-# INLINE fromRational #-}++++instance Floating DoubleX2 where+ pi = DoubleX2# 3.141592653589793238## 3.141592653589793238##+ {-# INLINE pi #-}+ exp (DoubleX2# a1 a2) = DoubleX2# (expDouble# a1)+ (expDouble# a2)+ {-# INLINE exp #-}+ log (DoubleX2# a1 a2) = DoubleX2# (logDouble# a1)+ (logDouble# a2)+ {-# INLINE log #-}+ sqrt (DoubleX2# a1 a2) = DoubleX2# (sqrtDouble# a1)+ (sqrtDouble# a2)+ {-# INLINE sqrt #-}+ sin (DoubleX2# a1 a2) = DoubleX2# (sinDouble# a1)+ (sinDouble# a2)+ {-# INLINE sin #-}+ cos (DoubleX2# a1 a2) = DoubleX2# (cosDouble# a1)+ (cosDouble# a2)+ {-# INLINE cos #-}+ tan (DoubleX2# a1 a2) = DoubleX2# (tanDouble# a1)+ (tanDouble# a2)+ {-# INLINE tan #-}+ asin (DoubleX2# a1 a2) = DoubleX2# (asinDouble# a1)+ (asinDouble# a2)+ {-# INLINE asin #-}+ acos (DoubleX2# a1 a2) = DoubleX2# (acosDouble# a1)+ (acosDouble# a2)+ {-# INLINE acos #-}+ atan (DoubleX2# a1 a2) = DoubleX2# (atanDouble# a1)+ (atanDouble# a2)+ {-# INLINE atan #-}+ sinh (DoubleX2# a1 a2) = DoubleX2# (sinDouble# a1)+ (sinDouble# a2)+ {-# INLINE sinh #-}+ cosh (DoubleX2# a1 a2) = DoubleX2# (coshDouble# a1)+ (coshDouble# a2)+ {-# INLINE cosh #-}+ tanh (DoubleX2# a1 a2) = DoubleX2# (tanhDouble# a1)+ (tanhDouble# a2)+ {-# INLINE tanh #-}+ DoubleX2# a1 a2 ** DoubleX2# b1 b2 = DoubleX2# ((**##) a1 b1)+ ((**##) a2 b2)+ {-# INLINE (**) #-}++ logBase x y = log y / log x+ {-# INLINE logBase #-}+ asinh x = log (x + sqrt (1.0+x*x))+ {-# INLINE asinh #-}+ acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))+ {-# INLINE acosh #-}+ atanh x = 0.5 * log ((1.0+x) / (1.0-x))+ {-# INLINE atanh #-}++++type instance ElemRep DoubleX2 = 'DoubleRep+type instance ElemPrim DoubleX2 = Double#+instance PrimBytes DoubleX2 where+ toBytes (DoubleX2# a1 a2) = case runRW#+ ( \s0 -> case newByteArray# (SIZEOF_HSFLOAT# *# 2#) s0 of+ (# s1, marr #) -> case writeDoubleArray# marr 0# a1 s1 of+ s2 -> case writeDoubleArray# marr 1# a2 s2 of+ s3 -> unsafeFreezeByteArray# marr s3+ ) of (# _, a #) -> (# 0#, 2#, a #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = DoubleX2#+ (indexDoubleArray# arr off)+ (indexDoubleArray# arr (off +# 1#))+ {-# INLINE fromBytes #-}+ byteSize _ = SIZEOF_HSFLOAT# *# 2#+ {-# INLINE byteSize #-}+ byteAlign _ = ALIGNMENT_HSFLOAT#+ {-# INLINE byteAlign #-}+ elementByteSize _ = SIZEOF_HSFLOAT#+ {-# INLINE elementByteSize #-}+ ix 0# (DoubleX2# a1 _) = a1+ ix 1# (DoubleX2# _ a2) = a2+ ix _ _ = undefined+ {-# INLINE ix #-}+++instance ElementWise (Idx '[2]) Double DoubleX2 where+ indexOffset# (DoubleX2# a1 _) 0# = D# a1+ indexOffset# (DoubleX2# _ a2) 1# = D# a2+ indexOffset# _ _ = undefined+ {-# INLINE indexOffset# #-}++ (!) (DoubleX2# a1 _) ( 1 :! Z) = D# a1+ (!) (DoubleX2# _ a2) ( 2 :! Z) = D# a2+ (!) _ ( _ :! Z) = undefined+ {-# INLINE (!) #-}++ broadcast (D# x) = DoubleX2# x x+ {-# INLINE broadcast #-}++ ewmap f (DoubleX2# x y) = case (f (1:!Z) (D# x), f (2:!Z) (D# y)) of+ (D# r1, D# r2) -> DoubleX2# r1 r2+ {-# INLINE ewmap #-}++ ewgen f = case (f (1:!Z), f (2:!Z)) of (D# r1, D# r2) -> DoubleX2# r1 r2+ {-# INLINE ewgen #-}++ ewgenA f = (\(D# r1) (D# r2) -> DoubleX2# r1 r2) <$> f (1:!Z) <*> f (2:!Z)+ {-# INLINE ewgenA #-}++ ewfoldl f x0 (DoubleX2# x y) = f (2:!Z) (f (1:!Z) x0 (D# x)) (D# y)+ {-# INLINE ewfoldl #-}++ ewfoldr f x0 (DoubleX2# x y) = f (1:!Z) (D# x) (f (2:!Z) (D# y) x0)+ {-# INLINE ewfoldr #-}++ elementWise f (DoubleX2# x y) = (\(D# a) (D# b) -> DoubleX2# a b)+ <$> f (D# x) <*> f (D# y)+ {-# INLINE elementWise #-}++ indexWise f (DoubleX2# x y) = (\(D# a) (D# b) -> DoubleX2# a b)+ <$> f (1:!Z) (D# x) <*> f (2:!Z) (D# y)+ {-# INLINE indexWise #-}++ update (1 :! Z) (D# q) (DoubleX2# _ y) = DoubleX2# q y+ update (2 :! Z) (D# q) (DoubleX2# x _) = DoubleX2# x q+ update (_ :! Z) _ x = x+ {-# INLINE update #-}
+ src-base/Numeric/Array/Family/DoubleX3.hs view
@@ -0,0 +1,299 @@+{-# LANGUAGE CPP #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UnboxedTuples #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+-----------------------------------------------------------------------------+-- |+-- Module : Numeric.Array.Family.DoubleX3+-- Copyright : (c) Artem Chirkin+-- License : BSD3+--+-- Maintainer : chirkin@arch.ethz.ch+--+--+-----------------------------------------------------------------------------++module Numeric.Array.Family.DoubleX3 () where+++#include "MachDeps.h"++import GHC.Base (runRW#)+import GHC.Prim+import GHC.Types (Double (..), RuntimeRep (..),+ isTrue#)++import Numeric.Array.ElementWise+import Numeric.Array.Family+import Numeric.Commons+import Numeric.Dimensions++++++instance Show DoubleX3 where+ show (DoubleX3# a1 a2 a3) = "{ " ++ show (D# a1)+ ++ ", " ++ show (D# a2)+ ++ ", " ++ show (D# a3)+ ++ " }"++++instance Eq DoubleX3 where+ DoubleX3# a1 a2 a3 == DoubleX3# b1 b2 b3 = isTrue# ( (a1 ==## b1)+ `andI#` (a2 ==## b2)+ `andI#` (a3 ==## b3)+ )+ {-# INLINE (==) #-}+ DoubleX3# a1 a2 a3 /= DoubleX3# b1 b2 b3 = isTrue# ( (a1 /=## b1)+ `orI#` (a2 /=## b2)+ `orI#` (a3 /=## b3)+ )+ {-# INLINE (/=) #-}++++-- | Implement partial ordering for `>`, `<`, `>=`, `<=`+-- and lexicographical ordering for `compare`+instance Ord DoubleX3 where+ DoubleX3# a1 a2 a3 > DoubleX3# b1 b2 b3 = isTrue# ( (a1 >## b1)+ `andI#` (a2 >## b2)+ `andI#` (a3 >## b3)+ )+ {-# INLINE (>) #-}+ DoubleX3# a1 a2 a3 < DoubleX3# b1 b2 b3 = isTrue# ( (a1 <## b1)+ `andI#` (a2 <## b2)+ `andI#` (a3 <## b3)+ )+ {-# INLINE (<) #-}+ DoubleX3# a1 a2 a3 >= DoubleX3# b1 b2 b3 = isTrue# ( (a1 >=## b1)+ `andI#` (a2 >=## b2)+ `andI#` (a3 >=## b3)+ )+ {-# INLINE (>=) #-}+ DoubleX3# a1 a2 a3 <= DoubleX3# b1 b2 b3 = isTrue# ( (a1 <=## b1)+ `andI#` (a2 <=## b2)+ `andI#` (a3 <=## b3)+ )+ {-# INLINE (<=) #-}+ -- | Compare lexicographically+ compare (DoubleX3# a1 a2 a3) (DoubleX3# b1 b2 b3)+ | isTrue# (a1 >## b1) = GT+ | isTrue# (a1 <## b1) = LT+ | isTrue# (a2 >## b2) = GT+ | isTrue# (a2 <## b2) = LT+ | isTrue# (a3 >## b3) = GT+ | isTrue# (a3 <## b3) = LT+ | otherwise = EQ+ {-# INLINE compare #-}+ -- | Element-wise minimum+ min (DoubleX3# a1 a2 a3) (DoubleX3# b1 b2 b3) =+ DoubleX3# (if isTrue# (a1 >## b1) then b1 else a1)+ (if isTrue# (a2 >## b2) then b2 else a2)+ (if isTrue# (a3 >## b3) then b3 else a3)+ {-# INLINE min #-}+ -- | Element-wise maximum+ max (DoubleX3# a1 a2 a3) (DoubleX3# b1 b2 b3) =+ DoubleX3# (if isTrue# (a1 >## b1) then a1 else b1)+ (if isTrue# (a2 >## b2) then a2 else b2)+ (if isTrue# (a3 >## b3) then a3 else b3)+ {-# INLINE max #-}++++-- | element-wise operations for vectors+instance Num DoubleX3 where+ DoubleX3# a1 a2 a3 + DoubleX3# b1 b2 b3+ = DoubleX3# ((+##) a1 b1) ((+##) a2 b2) ((+##) a3 b3)+ {-# INLINE (+) #-}+ DoubleX3# a1 a2 a3 - DoubleX3# b1 b2 b3+ = DoubleX3# ((-##) a1 b1) ((-##) a2 b2) ((-##) a3 b3)+ {-# INLINE (-) #-}+ DoubleX3# a1 a2 a3 * DoubleX3# b1 b2 b3+ = DoubleX3# ((*##) a1 b1) ((*##) a2 b2) ((*##) a3 b3)+ {-# INLINE (*) #-}+ negate (DoubleX3# a1 a2 a3)+ = DoubleX3# (negateDouble# a1) (negateDouble# a2) (negateDouble# a3)+ {-# INLINE negate #-}+ abs (DoubleX3# a1 a2 a3)+ = DoubleX3# (if isTrue# (a1 >=## 0.0##) then a1 else negateDouble# a1)+ (if isTrue# (a2 >=## 0.0##) then a2 else negateDouble# a2)+ (if isTrue# (a3 >=## 0.0##) then a3 else negateDouble# a3)+ {-# INLINE abs #-}+ signum (DoubleX3# a1 a2 a3)+ = DoubleX3# (if isTrue# (a1 >## 0.0##)+ then 1.0##+ else if isTrue# (a1 <## 0.0##) then -1.0## else 0.0## )+ (if isTrue# (a2 >## 0.0##)+ then 1.0##+ else if isTrue# (a2 <## 0.0##) then -1.0## else 0.0## )+ (if isTrue# (a3 >## 0.0##)+ then 1.0##+ else if isTrue# (a3 <## 0.0##) then -1.0## else 0.0## )+ {-# INLINE signum #-}+ fromInteger n = case fromInteger n of D# x -> DoubleX3# x x x+ {-# INLINE fromInteger #-}++++instance Fractional DoubleX3 where+ DoubleX3# a1 a2 a3 / DoubleX3# b1 b2 b3 = DoubleX3# ((/##) a1 b1)+ ((/##) a2 b2)+ ((/##) a3 b3)+ {-# INLINE (/) #-}+ recip (DoubleX3# a1 a2 a3) = DoubleX3# ((/##) 1.0## a1)+ ((/##) 1.0## a2)+ ((/##) 1.0## a3)+ {-# INLINE recip #-}+ fromRational r = case fromRational r of D# x -> DoubleX3# x x x+ {-# INLINE fromRational #-}++++instance Floating DoubleX3 where+ pi = DoubleX3# 3.141592653589793238## 3.141592653589793238## 3.141592653589793238##+ {-# INLINE pi #-}+ exp (DoubleX3# a1 a2 a3) = DoubleX3# (expDouble# a1)+ (expDouble# a2)+ (expDouble# a3)+ {-# INLINE exp #-}+ log (DoubleX3# a1 a2 a3) = DoubleX3# (logDouble# a1)+ (logDouble# a2)+ (logDouble# a3)+ {-# INLINE log #-}+ sqrt (DoubleX3# a1 a2 a3) = DoubleX3# (sqrtDouble# a1)+ (sqrtDouble# a2)+ (sqrtDouble# a3)+ {-# INLINE sqrt #-}+ sin (DoubleX3# a1 a2 a3) = DoubleX3# (sinDouble# a1)+ (sinDouble# a2)+ (sinDouble# a3)+ {-# INLINE sin #-}+ cos (DoubleX3# a1 a2 a3) = DoubleX3# (cosDouble# a1)+ (cosDouble# a2)+ (cosDouble# a3)+ {-# INLINE cos #-}+ tan (DoubleX3# a1 a2 a3) = DoubleX3# (tanDouble# a1)+ (tanDouble# a2)+ (tanDouble# a3)+ {-# INLINE tan #-}+ asin (DoubleX3# a1 a2 a3) = DoubleX3# (asinDouble# a1)+ (asinDouble# a2)+ (asinDouble# a3)+ {-# INLINE asin #-}+ acos (DoubleX3# a1 a2 a3) = DoubleX3# (acosDouble# a1)+ (acosDouble# a2)+ (acosDouble# a3)+ {-# INLINE acos #-}+ atan (DoubleX3# a1 a2 a3) = DoubleX3# (atanDouble# a1)+ (atanDouble# a2)+ (atanDouble# a3)+ {-# INLINE atan #-}+ sinh (DoubleX3# a1 a2 a3) = DoubleX3# (sinDouble# a1)+ (sinDouble# a2)+ (sinDouble# a3)+ {-# INLINE sinh #-}+ cosh (DoubleX3# a1 a2 a3) = DoubleX3# (coshDouble# a1)+ (coshDouble# a2)+ (coshDouble# a3)+ {-# INLINE cosh #-}+ tanh (DoubleX3# a1 a2 a3) = DoubleX3# (tanhDouble# a1)+ (tanhDouble# a2)+ (tanhDouble# a3)+ {-# INLINE tanh #-}+ DoubleX3# a1 a2 a3 ** DoubleX3# b1 b2 b3 = DoubleX3# ((**##) a1 b1)+ ((**##) a2 b2)+ ((**##) a3 b3)+ {-# INLINE (**) #-}++ logBase x y = log y / log x+ {-# INLINE logBase #-}+ asinh x = log (x + sqrt (1.0+x*x))+ {-# INLINE asinh #-}+ acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))+ {-# INLINE acosh #-}+ atanh x = 0.5 * log ((1.0+x) / (1.0-x))+ {-# INLINE atanh #-}++++type instance ElemRep DoubleX3 = 'DoubleRep+type instance ElemPrim DoubleX3 = Double#+instance PrimBytes DoubleX3 where+ toBytes (DoubleX3# a1 a2 a3) = case runRW#+ ( \s0 -> case newByteArray# (SIZEOF_HSFLOAT# *# 3#) s0 of+ (# s1, marr #) -> case writeDoubleArray# marr 0# a1 s1 of+ s2 -> case writeDoubleArray# marr 1# a2 s2 of+ s3 -> case writeDoubleArray# marr 2# a3 s3 of+ s4 -> unsafeFreezeByteArray# marr s4+ ) of (# _, a #) -> (# 0#, 3#, a #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = DoubleX3#+ (indexDoubleArray# arr off)+ (indexDoubleArray# arr (off +# 1#))+ (indexDoubleArray# arr (off +# 2#))+ {-# INLINE fromBytes #-}+ byteSize _ = SIZEOF_HSFLOAT# *# 3#+ {-# INLINE byteSize #-}+ byteAlign _ = ALIGNMENT_HSFLOAT#+ {-# INLINE byteAlign #-}+ elementByteSize _ = SIZEOF_HSFLOAT#+ {-# INLINE elementByteSize #-}+ ix 0# (DoubleX3# a1 _ _) = a1+ ix 1# (DoubleX3# _ a2 _) = a2+ ix 2# (DoubleX3# _ _ a3) = a3+ ix _ _ = undefined+ {-# INLINE ix #-}+++instance ElementWise (Idx '[3]) Double DoubleX3 where+ indexOffset# (DoubleX3# a1 _ _) 0# = D# a1+ indexOffset# (DoubleX3# _ a2 _) 1# = D# a2+ indexOffset# (DoubleX3# _ _ a3) 2# = D# a3+ indexOffset# _ _ = undefined+ {-# INLINE indexOffset# #-}++ (!) (DoubleX3# a1 _ _) ( 1 :! Z) = D# a1+ (!) (DoubleX3# _ a2 _) ( 2 :! Z) = D# a2+ (!) (DoubleX3# _ _ a3) ( 3 :! Z) = D# a3+ (!) _ ( _ :! Z) = undefined+ {-# INLINE (!) #-}++ broadcast (D# x) = DoubleX3# x x x+ {-# INLINE broadcast #-}++ ewmap f (DoubleX3# x y z) = case (f (1:!Z) (D# x), f (2:!Z) (D# y), f (3:!Z) (D# z)) of+ (D# r1, D# r2, D# r3) -> DoubleX3# r1 r2 r3+ {-# INLINE ewmap #-}++ ewgen f = case (f (1:!Z), f (2:!Z), f (3:!Z)) of (D# r1, D# r2, D# r3) -> DoubleX3# r1 r2 r3+ {-# INLINE ewgen #-}++ ewgenA f = (\(D# r1) (D# r2) (D# r3) -> DoubleX3# r1 r2 r3)+ <$> f (1:!Z) <*> f (2:!Z) <*> f (3:!Z)+ {-# INLINE ewgenA #-}++ ewfoldl f x0 (DoubleX3# x y z) = f (3:!Z) (f (2:!Z) (f (1:!Z) x0 (D# x)) (D# y)) (D# z)+ {-# INLINE ewfoldl #-}++ ewfoldr f x0 (DoubleX3# x y z) = f (1:!Z) (D# x) (f (2:!Z) (D# y) (f (3:!Z) (D# z) x0))+ {-# INLINE ewfoldr #-}++ elementWise f (DoubleX3# x y z) = (\(D# a) (D# b) (D# c) -> DoubleX3# a b c)+ <$> f (D# x) <*> f (D# y) <*> f (D# z)+ {-# INLINE elementWise #-}++ indexWise f (DoubleX3# x y z) = (\(D# a) (D# b) (D# c) -> DoubleX3# a b c)+ <$> f (1:!Z) (D# x) <*> f (2:!Z) (D# y) <*> f (3:!Z) (D# z)+ {-# INLINE indexWise #-}++ update (1 :! Z) (D# q) (DoubleX3# _ y z) = DoubleX3# q y z+ update (2 :! Z) (D# q) (DoubleX3# x _ z) = DoubleX3# x q z+ update (3 :! Z) (D# q) (DoubleX3# x y _) = DoubleX3# x y q+ update (_ :! Z) _ x = x+ {-# INLINE update #-}
+ src-base/Numeric/Array/Family/DoubleX4.hs view
@@ -0,0 +1,335 @@+{-# LANGUAGE CPP #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UnboxedTuples #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+-----------------------------------------------------------------------------+-- |+-- Module : Numeric.Array.Family.DoubleX4+-- Copyright : (c) Artem Chirkin+-- License : BSD3+--+-- Maintainer : chirkin@arch.ethz.ch+--+--+-----------------------------------------------------------------------------++module Numeric.Array.Family.DoubleX4 () where+++#include "MachDeps.h"++import GHC.Base (runRW#)+import GHC.Prim+import GHC.Types (Double (..), RuntimeRep (..),+ isTrue#)++import Numeric.Array.ElementWise+import Numeric.Array.Family+import Numeric.Commons+import Numeric.Dimensions++++++instance Show DoubleX4 where+ show (DoubleX4# a1 a2 a3 a4) = "{ " ++ show (D# a1)+ ++ ", " ++ show (D# a2)+ ++ ", " ++ show (D# a3)+ ++ ", " ++ show (D# a4)+ ++ " }"++++instance Eq DoubleX4 where+ DoubleX4# a1 a2 a3 a4 == DoubleX4# b1 b2 b3 b4 = isTrue# ( (a1 ==## b1)+ `andI#` (a2 ==## b2)+ `andI#` (a3 ==## b3)+ `andI#` (a4 ==## b4)+ )+ {-# INLINE (==) #-}+ DoubleX4# a1 a2 a3 a4 /= DoubleX4# b1 b2 b3 b4 = isTrue# ( (a1 /=## b1)+ `orI#` (a2 /=## b2)+ `orI#` (a3 /=## b3)+ `orI#` (a4 /=## b4)+ )+ {-# INLINE (/=) #-}++++-- | Implement partial ordering for `>`, `<`, `>=`, `<=`+-- and lexicographical ordering for `compare`+instance Ord DoubleX4 where+ DoubleX4# a1 a2 a3 a4 > DoubleX4# b1 b2 b3 b4 = isTrue# ( (a1 >## b1)+ `andI#` (a2 >## b2)+ `andI#` (a3 >## b3)+ `andI#` (a4 >## b4)+ )+ {-# INLINE (>) #-}+ DoubleX4# a1 a2 a3 a4 < DoubleX4# b1 b2 b3 b4 = isTrue# ( (a1 <## b1)+ `andI#` (a2 <## b2)+ `andI#` (a3 <## b3)+ `andI#` (a4 <## b4)+ )+ {-# INLINE (<) #-}+ DoubleX4# a1 a2 a3 a4 >= DoubleX4# b1 b2 b3 b4 = isTrue# ( (a1 >=## b1)+ `andI#` (a2 >=## b2)+ `andI#` (a3 >=## b3)+ `andI#` (a4 >=## b4)+ )+ {-# INLINE (>=) #-}+ DoubleX4# a1 a2 a3 a4 <= DoubleX4# b1 b2 b3 b4 = isTrue# ( (a1 <=## b1)+ `andI#` (a2 <=## b2)+ `andI#` (a3 <=## b3)+ `andI#` (a4 <=## b4)+ )+ {-# INLINE (<=) #-}+ -- | Compare lexicographically+ compare (DoubleX4# a1 a2 a3 a4) (DoubleX4# b1 b2 b3 b4)+ | isTrue# (a1 >## b1) = GT+ | isTrue# (a1 <## b1) = LT+ | isTrue# (a2 >## b2) = GT+ | isTrue# (a2 <## b2) = LT+ | isTrue# (a3 >## b3) = GT+ | isTrue# (a3 <## b3) = LT+ | isTrue# (a4 >## b4) = GT+ | isTrue# (a4 <## b4) = LT+ | otherwise = EQ+ {-# INLINE compare #-}+ -- | Element-wise minimum+ min (DoubleX4# a1 a2 a3 a4) (DoubleX4# b1 b2 b3 b4) =+ DoubleX4# (if isTrue# (a1 >## b1) then b1 else a1)+ (if isTrue# (a2 >## b2) then b2 else a2)+ (if isTrue# (a3 >## b3) then b3 else a3)+ (if isTrue# (a4 >## b4) then b4 else a4)+ {-# INLINE min #-}+ -- | Element-wise maximum+ max (DoubleX4# a1 a2 a3 a4) (DoubleX4# b1 b2 b3 b4) =+ DoubleX4# (if isTrue# (a1 >## b1) then a1 else b1)+ (if isTrue# (a2 >## b2) then a2 else b2)+ (if isTrue# (a3 >## b3) then a3 else b3)+ (if isTrue# (a4 >## b4) then a4 else b4)+ {-# INLINE max #-}++++-- | element-wise operations for vectors+instance Num DoubleX4 where+ DoubleX4# a1 a2 a3 a4 + DoubleX4# b1 b2 b3 b4+ = DoubleX4# ((+##) a1 b1) ((+##) a2 b2) ((+##) a3 b3) ((+##) a4 b4)+ {-# INLINE (+) #-}+ DoubleX4# a1 a2 a3 a4 - DoubleX4# b1 b2 b3 b4+ = DoubleX4# ((-##) a1 b1) ((-##) a2 b2) ((-##) a3 b3) ((-##) a4 b4)+ {-# INLINE (-) #-}+ DoubleX4# a1 a2 a3 a4 * DoubleX4# b1 b2 b3 b4+ = DoubleX4# ((*##) a1 b1) ((*##) a2 b2) ((*##) a3 b3) ((*##) a4 b4)+ {-# INLINE (*) #-}+ negate (DoubleX4# a1 a2 a3 a4)+ = DoubleX4# (negateDouble# a1) (negateDouble# a2) (negateDouble# a3) (negateDouble# a4)+ {-# INLINE negate #-}+ abs (DoubleX4# a1 a2 a3 a4)+ = DoubleX4# (if isTrue# (a1 >=## 0.0##) then a1 else negateDouble# a1)+ (if isTrue# (a2 >=## 0.0##) then a2 else negateDouble# a2)+ (if isTrue# (a3 >=## 0.0##) then a3 else negateDouble# a3)+ (if isTrue# (a4 >=## 0.0##) then a4 else negateDouble# a4)+ {-# INLINE abs #-}+ signum (DoubleX4# a1 a2 a3 a4)+ = DoubleX4# (if isTrue# (a1 >## 0.0##)+ then 1.0##+ else if isTrue# (a1 <## 0.0##) then -1.0## else 0.0## )+ (if isTrue# (a2 >## 0.0##)+ then 1.0##+ else if isTrue# (a2 <## 0.0##) then -1.0## else 0.0## )+ (if isTrue# (a3 >## 0.0##)+ then 1.0##+ else if isTrue# (a3 <## 0.0##) then -1.0## else 0.0## )+ (if isTrue# (a4 >## 0.0##)+ then 1.0##+ else if isTrue# (a4 <## 0.0##) then -1.0## else 0.0## )+ {-# INLINE signum #-}+ fromInteger n = case fromInteger n of D# x -> DoubleX4# x x x x+ {-# INLINE fromInteger #-}++++instance Fractional DoubleX4 where+ DoubleX4# a1 a2 a3 a4 / DoubleX4# b1 b2 b3 b4 = DoubleX4# ((/##) a1 b1)+ ((/##) a2 b2)+ ((/##) a3 b3)+ ((/##) a4 b4)+ {-# INLINE (/) #-}+ recip (DoubleX4# a1 a2 a3 a4) = DoubleX4# ((/##) 1.0## a1)+ ((/##) 1.0## a2)+ ((/##) 1.0## a3)+ ((/##) 1.0## a4)+ {-# INLINE recip #-}+ fromRational r = case fromRational r of D# x -> DoubleX4# x x x x+ {-# INLINE fromRational #-}++++instance Floating DoubleX4 where+ pi = DoubleX4# 3.141592653589793238## 3.141592653589793238## 3.141592653589793238## 3.141592653589793238##+ {-# INLINE pi #-}+ exp (DoubleX4# a1 a2 a3 a4) = DoubleX4# (expDouble# a1)+ (expDouble# a2)+ (expDouble# a3)+ (expDouble# a4)+ {-# INLINE exp #-}+ log (DoubleX4# a1 a2 a3 a4) = DoubleX4# (logDouble# a1)+ (logDouble# a2)+ (logDouble# a3)+ (logDouble# a4)+ {-# INLINE log #-}+ sqrt (DoubleX4# a1 a2 a3 a4) = DoubleX4# (sqrtDouble# a1)+ (sqrtDouble# a2)+ (sqrtDouble# a3)+ (sqrtDouble# a4)+ {-# INLINE sqrt #-}+ sin (DoubleX4# a1 a2 a3 a4) = DoubleX4# (sinDouble# a1)+ (sinDouble# a2)+ (sinDouble# a3)+ (sinDouble# a4)+ {-# INLINE sin #-}+ cos (DoubleX4# a1 a2 a3 a4) = DoubleX4# (cosDouble# a1)+ (cosDouble# a2)+ (cosDouble# a3)+ (cosDouble# a4)+ {-# INLINE cos #-}+ tan (DoubleX4# a1 a2 a3 a4) = DoubleX4# (tanDouble# a1)+ (tanDouble# a2)+ (tanDouble# a3)+ (tanDouble# a4)+ {-# INLINE tan #-}+ asin (DoubleX4# a1 a2 a3 a4) = DoubleX4# (asinDouble# a1)+ (asinDouble# a2)+ (asinDouble# a3)+ (asinDouble# a4)+ {-# INLINE asin #-}+ acos (DoubleX4# a1 a2 a3 a4) = DoubleX4# (acosDouble# a1)+ (acosDouble# a2)+ (acosDouble# a3)+ (acosDouble# a4)+ {-# INLINE acos #-}+ atan (DoubleX4# a1 a2 a3 a4) = DoubleX4# (atanDouble# a1)+ (atanDouble# a2)+ (atanDouble# a3)+ (atanDouble# a4)+ {-# INLINE atan #-}+ sinh (DoubleX4# a1 a2 a3 a4) = DoubleX4# (sinDouble# a1)+ (sinDouble# a2)+ (sinDouble# a3)+ (sinDouble# a4)+ {-# INLINE sinh #-}+ cosh (DoubleX4# a1 a2 a3 a4) = DoubleX4# (coshDouble# a1)+ (coshDouble# a2)+ (coshDouble# a3)+ (coshDouble# a4)+ {-# INLINE cosh #-}+ tanh (DoubleX4# a1 a2 a3 a4) = DoubleX4# (tanhDouble# a1)+ (tanhDouble# a2)+ (tanhDouble# a3)+ (tanhDouble# a4)+ {-# INLINE tanh #-}+ DoubleX4# a1 a2 a3 a4 ** DoubleX4# b1 b2 b3 b4 = DoubleX4# ((**##) a1 b1)+ ((**##) a2 b2)+ ((**##) a3 b3)+ ((**##) a4 b4)+ {-# INLINE (**) #-}++ logBase x y = log y / log x+ {-# INLINE logBase #-}+ asinh x = log (x + sqrt (1.0+x*x))+ {-# INLINE asinh #-}+ acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))+ {-# INLINE acosh #-}+ atanh x = 0.5 * log ((1.0+x) / (1.0-x))+ {-# INLINE atanh #-}++++type instance ElemRep DoubleX4 = 'DoubleRep+type instance ElemPrim DoubleX4 = Double#+instance PrimBytes DoubleX4 where+ toBytes (DoubleX4# a1 a2 a3 a4) = case runRW#+ ( \s0 -> case newByteArray# (SIZEOF_HSFLOAT# *# 3#) s0 of+ (# s1, marr #) -> case writeDoubleArray# marr 0# a1 s1 of+ s2 -> case writeDoubleArray# marr 1# a2 s2 of+ s3 -> case writeDoubleArray# marr 2# a3 s3 of+ s4 -> case writeDoubleArray# marr 3# a4 s4 of+ s5 -> unsafeFreezeByteArray# marr s5+ ) of (# _, a #) -> (# 0#, 4#, a #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = DoubleX4#+ (indexDoubleArray# arr off)+ (indexDoubleArray# arr (off +# 1#))+ (indexDoubleArray# arr (off +# 2#))+ (indexDoubleArray# arr (off +# 3#))+ {-# INLINE fromBytes #-}+ byteSize _ = SIZEOF_HSFLOAT# *# 4#+ {-# INLINE byteSize #-}+ byteAlign _ = ALIGNMENT_HSFLOAT#+ {-# INLINE byteAlign #-}+ elementByteSize _ = SIZEOF_HSFLOAT#+ {-# INLINE elementByteSize #-}+ ix 0# (DoubleX4# a1 _ _ _) = a1+ ix 1# (DoubleX4# _ a2 _ _) = a2+ ix 2# (DoubleX4# _ _ a3 _) = a3+ ix 3# (DoubleX4# _ _ _ a4) = a4+ ix _ _ = undefined+ {-# INLINE ix #-}+++instance ElementWise (Idx '[4]) Double DoubleX4 where+ indexOffset# (DoubleX4# a1 _ _ _) 0# = D# a1+ indexOffset# (DoubleX4# _ a2 _ _) 1# = D# a2+ indexOffset# (DoubleX4# _ _ a3 _) 2# = D# a3+ indexOffset# (DoubleX4# _ _ _ a4) 3# = D# a4+ indexOffset# _ _ = undefined+ {-# INLINE indexOffset# #-}++ (!) (DoubleX4# a1 _ _ _) ( 1 :! Z) = D# a1+ (!) (DoubleX4# _ a2 _ _) ( 2 :! Z) = D# a2+ (!) (DoubleX4# _ _ a3 _) ( 3 :! Z) = D# a3+ (!) (DoubleX4# _ _ _ a4) ( 4 :! Z) = D# a4+ (!) _ ( _ :! Z) = undefined+ {-# INLINE (!) #-}++ broadcast (D# x) = DoubleX4# x x x x+ {-# INLINE broadcast #-}++ ewmap f (DoubleX4# x y z w) = case (f (1:!Z) (D# x), f (2:!Z) (D# y), f (3:!Z) (D# z), f (3:!Z) (D# w)) of+ (D# r1, D# r2, D# r3, D# r4) -> DoubleX4# r1 r2 r3 r4+ {-# INLINE ewmap #-}++ ewgen f = case (f (1:!Z), f (2:!Z), f (3:!Z), f (4:!Z)) of (D# r1, D# r2, D# r3, D# r4) -> DoubleX4# r1 r2 r3 r4+ {-# INLINE ewgen #-}++ ewgenA f = (\(D# a) (D# b) (D# c) (D# d) -> DoubleX4# a b c d)+ <$> f (1:!Z) <*> f (2:!Z) <*> f (3:!Z) <*> f (4:!Z)+ {-# INLINE ewgenA #-}++ ewfoldl f x0 (DoubleX4# x y z w) = f (4:!Z) (f (3:!Z) (f (2:!Z) (f (1:!Z) x0 (D# x)) (D# y)) (D# z)) (D# w)+ {-# INLINE ewfoldl #-}++ ewfoldr f x0 (DoubleX4# x y z w) = f (1:!Z) (D# x) (f (2:!Z) (D# y) (f (3:!Z) (D# z) (f (4:!Z) (D# w) x0)))+ {-# INLINE ewfoldr #-}++ elementWise f (DoubleX4# x y z w) = (\(D# a) (D# b) (D# c) (D# d) -> DoubleX4# a b c d)+ <$> f (D# x) <*> f (D# y) <*> f (D# z) <*> f (D# w)+ {-# INLINE elementWise #-}++ indexWise f (DoubleX4# x y z w) = (\(D# a) (D# b) (D# c) (D# d) -> DoubleX4# a b c d)+ <$> f (1:!Z) (D# x) <*> f (2:!Z) (D# y) <*> f (3:!Z) (D# z) <*> f (4:!Z) (D# w)+ {-# INLINE indexWise #-}++ update (1 :! Z) (D# q) (DoubleX4# _ y z w) = DoubleX4# q y z w+ update (2 :! Z) (D# q) (DoubleX4# x _ z w) = DoubleX4# x q z w+ update (3 :! Z) (D# q) (DoubleX4# x y _ w) = DoubleX4# x y q w+ update (4 :! Z) (D# q) (DoubleX4# x y z _) = DoubleX4# x y z q+ update (_ :! Z) _ x = x+ {-# INLINE update #-}
src-base/Numeric/DataFrame/Contraction.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE FlexibleContexts #-} {-# LANGUAGE FlexibleInstances #-}@@ -30,6 +31,7 @@ ( Contraction (..), (%*) ) where +#include "MachDeps.h" import Data.Int (Int16, Int32, Int64, Int8) import Data.Word (Word16, Word32, Word64, Word8)@@ -198,10 +200,17 @@ , I# m <- intNatVal pm , I# n <- totalDim (Proxy @as) , I# k <- totalDim (Proxy @bs)+#if WORD_SIZE_IN_BITS < 64+ , Refl <- unsafeCoerce Refl :: ElemRep (Array Int64 (m : bs) ) :~: 'Int64Rep+ , Refl <- unsafeCoerce Refl :: ElemPrim (Array Int64 (m : bs) ) :~: Int64#+ , Refl <- unsafeCoerce Refl :: ElemRep (Array Int64 (as +: m)) :~: 'Int64Rep+ , Refl <- unsafeCoerce Refl :: ElemPrim (Array Int64 (as +: m)) :~: Int64#+#else , Refl <- unsafeCoerce Refl :: ElemRep (Array Int64 (m : bs) ) :~: 'IntRep , Refl <- unsafeCoerce Refl :: ElemPrim (Array Int64 (m : bs) ) :~: Int# , Refl <- unsafeCoerce Refl :: ElemRep (Array Int64 (as +: m)) :~: 'IntRep , Refl <- unsafeCoerce Refl :: ElemPrim (Array Int64 (as +: m)) :~: Int#+#endif = prodI64 n m k x y where getM :: forall m p . p (m ': bs) -> Proxy m@@ -290,18 +299,27 @@ , I# m <- intNatVal pm , I# n <- totalDim (Proxy @as) , I# k <- totalDim (Proxy @bs)+#if WORD_SIZE_IN_BITS < 64+ , Refl <- unsafeCoerce Refl :: ElemRep (Array Word64 (m : bs) ) :~: 'Word64Rep+ , Refl <- unsafeCoerce Refl :: ElemPrim (Array Word64 (m : bs) ) :~: Word64#+ , Refl <- unsafeCoerce Refl :: ElemRep (Array Word64 (as +: m)) :~: 'Word64Rep+ , Refl <- unsafeCoerce Refl :: ElemPrim (Array Word64 (as +: m)) :~: Word64#+#else , Refl <- unsafeCoerce Refl :: ElemRep (Array Word64 (m : bs) ) :~: 'WordRep , Refl <- unsafeCoerce Refl :: ElemPrim (Array Word64 (m : bs) ) :~: Word# , Refl <- unsafeCoerce Refl :: ElemRep (Array Word64 (as +: m)) :~: 'WordRep , Refl <- unsafeCoerce Refl :: ElemPrim (Array Word64 (as +: m)) :~: Word#+#endif = prodW64 n m k x y where getM :: forall m p . p (m ': bs) -> Proxy m getM _ = Proxy --prodF :: (FloatBytes a, FloatBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodF :: (PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Float#, ElemRep a ~ 'FloatRep+ , ElemPrim b ~ Float#, ElemRep b ~ 'FloatRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodF n m k x y = case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -317,7 +335,10 @@ bs = n *# k *# elementByteSize x {-# INLINE prodF #-} -prodD :: (DoubleBytes a, DoubleBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodD :: (PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Double#, ElemRep a ~ 'DoubleRep+ , ElemPrim b ~ Double#, ElemRep b ~ 'DoubleRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodD n m k x y= case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -333,7 +354,10 @@ bs = n *# k *# elementByteSize x {-# INLINE prodD #-} -prodI :: (IntBytes a, IntBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodI :: (PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Int#, ElemRep a ~ 'IntRep+ , ElemPrim b ~ Int#, ElemRep b ~ 'IntRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodI n m k x y= case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -349,7 +373,10 @@ bs = n *# k *# elementByteSize x {-# INLINE prodI #-} -prodI8 :: (IntBytes a, IntBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodI8 :: (PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Int#, ElemRep a ~ 'IntRep+ , ElemPrim b ~ Int#, ElemRep b ~ 'IntRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodI8 n m k x y= case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -366,7 +393,10 @@ {-# INLINE prodI8 #-} -prodI16 :: (IntBytes a, IntBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodI16 :: ( PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Int#, ElemRep a ~ 'IntRep+ , ElemPrim b ~ Int#, ElemRep b ~ 'IntRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodI16 n m k x y= case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -383,7 +413,10 @@ {-# INLINE prodI16 #-} -prodI32 :: (IntBytes a, IntBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodI32 :: ( PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Int#, ElemRep a ~ 'IntRep+ , ElemPrim b ~ Int#, ElemRep b ~ 'IntRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodI32 n m k x y= case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -399,7 +432,19 @@ bs = n *# k *# elementByteSize x {-# INLINE prodI32 #-} -prodI64 :: (IntBytes a, IntBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c++prodI64 :: ( PrimBytes a, PrimBytes b, PrimBytes c+#if WORD_SIZE_IN_BITS < 64+ , ElemPrim a ~ Int64#, ElemRep a ~ 'Int64Rep+ , ElemPrim b ~ Int64#, ElemRep b ~ 'Int64Rep+#else+ , ElemPrim a ~ Int#, ElemRep a ~ 'IntRep+ , ElemPrim b ~ Int#, ElemRep b ~ 'IntRep+#endif+ ) => Int# -> Int# -> Int# -> a -> b -> c+#if WORD_SIZE_IN_BITS < 64+prodI64 = undefined+#else prodI64 n m k x y= case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -414,9 +459,12 @@ where bs = n *# k *# elementByteSize x {-# INLINE prodI64 #-}-+#endif -prodW :: (WordBytes a, WordBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodW :: ( PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Word#, ElemRep a ~ 'WordRep+ , ElemPrim b ~ Word#, ElemRep b ~ 'WordRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodW n m k x y = case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -432,7 +480,10 @@ bs = n *# k *# elementByteSize x {-# INLINE prodW #-} -prodW8 :: (WordBytes a, WordBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodW8 :: ( PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Word#, ElemRep a ~ 'WordRep+ , ElemPrim b ~ Word#, ElemRep b ~ 'WordRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodW8 n m k x y = case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -449,7 +500,10 @@ {-# INLINE prodW8 #-} -prodW16 :: (WordBytes a, WordBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodW16 :: ( PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Word#, ElemRep a ~ 'WordRep+ , ElemPrim b ~ Word#, ElemRep b ~ 'WordRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodW16 n m k x y = case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -465,7 +519,10 @@ bs = n *# k *# elementByteSize x {-# INLINE prodW16 #-} -prodW32 :: (WordBytes a, WordBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodW32 :: ( PrimBytes a, PrimBytes b, PrimBytes c+ , ElemPrim a ~ Word#, ElemRep a ~ 'WordRep+ , ElemPrim b ~ Word#, ElemRep b ~ 'WordRep+ ) => Int# -> Int# -> Int# -> a -> b -> c prodW32 n m k x y = case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -481,7 +538,18 @@ bs = n *# k *# elementByteSize x {-# INLINE prodW32 #-} -prodW64 :: (WordBytes a, WordBytes b, PrimBytes c) => Int# -> Int# -> Int# -> a -> b -> c+prodW64 :: ( PrimBytes a, PrimBytes b, PrimBytes c+#if WORD_SIZE_IN_BITS < 64+ , ElemPrim a ~ Word64#, ElemRep a ~ 'Word64Rep+ , ElemPrim b ~ Word64#, ElemRep b ~ 'Word64Rep+#else+ , ElemPrim a ~ Word#, ElemRep a ~ 'WordRep+ , ElemPrim b ~ Word#, ElemRep b ~ 'WordRep+#endif+ ) => Int# -> Int# -> Int# -> a -> b -> c+#if WORD_SIZE_IN_BITS < 64+prodW64 = undefined+#else prodW64 n m k x y = case runRW# ( \s0 -> case newByteArray# bs s0 of (# s1, marr #) ->@@ -496,6 +564,7 @@ where bs = n *# k *# elementByteSize x {-# INLINE prodW64 #-}+#endif -- | Do something in a loop for int i from 0 to n-1 and j from 0 to m-1 loop2# :: Int# -> Int# -> (Int# -> Int#-> State# s -> State# s) -> State# s -> State# s@@ -505,5 +574,3 @@ | isTrue# (i ==# n) = loop' 0# (j +# 1#) s | otherwise = case f i j s of s1 -> loop' (i +# 1#) j s1 {-# INLINE loop2# #-}--
src-base/Numeric/DataFrame/Inference.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE ExistentialQuantification #-} {-# LANGUAGE FlexibleContexts #-}@@ -61,7 +62,7 @@ ETWord8 -> Evidence ETWord16 -> Evidence ETWord32 -> Evidence- ETWord64 -> Evidence+ ETWord64 -> Evidence AIArrayF -> Evidence AIArrayD -> Evidence AIArrayI -> Evidence@@ -77,6 +78,9 @@ AIFloatX2 -> Evidence AIFloatX3 -> Evidence AIFloatX4 -> Evidence+ AIDoubleX2 -> Evidence+ AIDoubleX3 -> Evidence+ AIDoubleX4 -> Evidence inferElementWise :: forall t (ds :: [Nat]) . ( ArrayInstanceInference t ds@@ -100,6 +104,9 @@ AIFloatX2 -> Evidence AIFloatX3 -> Evidence AIFloatX4 -> Evidence+ AIDoubleX2 -> Evidence+ AIDoubleX3 -> Evidence+ AIDoubleX4 -> Evidence inferNumericFrame :: forall t (ds :: [Nat])@@ -113,6 +120,9 @@ AIFloatX2 -> Evidence AIFloatX3 -> Evidence AIFloatX4 -> Evidence+ AIDoubleX2 -> Evidence+ AIDoubleX3 -> Evidence+ AIDoubleX4 -> Evidence AIScalar -> case elemTypeInstance @t of ETFloat -> Evidence ETDouble -> Evidence@@ -125,7 +135,7 @@ ETWord8 -> Evidence ETWord16 -> Evidence ETWord32 -> Evidence- ETWord64 -> Evidence+ ETWord64 -> Evidence AIArrayF -> Evidence AIArrayD -> Evidence AIArrayI -> Evidence
+ src-base/Numeric/DataFrame/Mutable.hs view
@@ -0,0 +1,224 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE ExistentialQuantification #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE KindSignatures #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UnboxedTuples #-}+{-# LANGUAGE TypeOperators #-}+-----------------------------------------------------------------------------+-- |+-- Module : Numeric.DataFrame.Mutable+-- Copyright : (c) Artem Chirkin+-- License : BSD3+--+-- Maintainer : chirkin@arch.ethz.ch+--+-- Interfrace to perform primitive stateful operations on mutable frames.+--+-----------------------------------------------------------------------------++module Numeric.DataFrame.Mutable+ ( MutableFrame (..), MDataFrame ()+ , newDataFrame#, copyDataFrame#, copyMDataFrame#, unsafeFreezeDataFrame#+ , freezeDataFrame#, thawDataFrame#+ , writeDataFrame#, readDataFrame#+ ) where+++import GHC.Int (Int16 (..), Int32 (..), Int64 (..),+ Int8 (..))+import GHC.Prim+import GHC.Types (Double (..), Float (..), Int (..),+ Word (..))+import GHC.Word (Word16 (..), Word32 (..), Word64 (..),+ Word8 (..))++import Numeric.Commons+import Numeric.DataFrame.Type+import Numeric.Dimensions+import Numeric.TypeLits+++-- | Mutable DataFrame type+data MDataFrame s t (ns :: [Nat]) = MDataFrame# Int# Int# (MutableByteArray# s)++-- | Create a new mutable DataFrame.+newDataFrame# :: forall t (ns :: [Nat]) s+ . ( PrimBytes t, Dimensions ns)+ => State# s -> (# State# s, MDataFrame s t ns #)+newDataFrame# s0+ | elS <- elementByteSize (undefined :: t)+ , I# n <- totalDim (Proxy @ns)+ , (# s1, mba #) <- newByteArray# (n *# elS) s0+ = (# s1, MDataFrame# 0# n mba #)+{-# INLINE newDataFrame# #-}++-- | Copy one DataFrame into another mutable DataFrame at specified position.+copyDataFrame# :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s+ . ( PrimBytes (DataFrame t (as +: b'))+ , ConcatList as (b :+ bs) asbs+ , Dimensions (b :+ bs)+ )+ => DataFrame t (as +: b') -> Idx (b :+ bs) -> MDataFrame s t asbs -> State# s -> (# State# s, () #)+copyDataFrame# df ei (MDataFrame# offM _ arrM) s+ | (# offA, lenA, arrA #) <- toBytes df+ , elS <- elementByteSize df+ , I# i <- fromEnum ei+ = (# copyByteArray# arrA (offA *# elS) arrM ((offM +# i) *# elS) (lenA *# elS) s, () #)+{-# INLINE copyDataFrame# #-}++-- | Copy one mutable DataFrame into another mutable DataFrame at specified position.+copyMDataFrame# :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s+ . ( PrimBytes t+ , ConcatList as (b :+ bs) asbs+ , Dimensions (b :+ bs)+ )+ => MDataFrame s t (as +: b') -> Idx (b :+ bs) -> MDataFrame s t asbs -> State# s -> (# State# s, () #)+copyMDataFrame# (MDataFrame# offA lenA arrA) ei (MDataFrame# offM _ arrM) s+ | elS <- elementByteSize (undefined :: t)+ , I# i <- fromEnum ei+ = (# copyMutableByteArray# arrA (offA *# elS) arrM ((offM +# i) *# elS) (lenA *# elS) s, () #)+{-# INLINE copyMDataFrame# #-}++-- | Make a mutable DataFrame immutable, without copying.+unsafeFreezeDataFrame# :: forall t (ns :: [Nat]) s+ . PrimBytes (DataFrame t ns)+ => MDataFrame s t ns -> State# s -> (# State# s, DataFrame t ns #)+unsafeFreezeDataFrame# (MDataFrame# offM lenM arrM) s1+ | (# s2, arrA #) <- unsafeFreezeByteArray# arrM s1+ = (# s2, fromBytes (# offM, lenM, arrA #) #)+{-# INLINE unsafeFreezeDataFrame# #-}++-- | Copy content of a mutable DataFrame into a new immutable DataFrame.+freezeDataFrame# :: forall t (ns :: [Nat]) s+ . PrimBytes (DataFrame t ns)+ => MDataFrame s t ns -> State# s -> (# State# s, DataFrame t ns #)+freezeDataFrame# (MDataFrame# offM n arrM) s0+ | elS <- elementByteSize (undefined :: DataFrame t ns)+ , (# s1, mba #) <- newByteArray# (n *# elS) s0+ , s2 <- copyMutableByteArray# arrM (offM *# elS) mba 0# (n *# elS) s1+ , (# s3, arrA #) <- unsafeFreezeByteArray# mba s2+ = (# s3, fromBytes (# 0#, n, arrA #) #)+{-# INLINE freezeDataFrame# #-}++-- | Create a new mutable DataFrame and copy content of immutable one in there.+thawDataFrame# :: forall t (ns :: [Nat]) s+ . PrimBytes (DataFrame t ns)+ => DataFrame t ns -> State# s -> (# State# s, MDataFrame s t ns #)+thawDataFrame# df s0+ | elS <- elementByteSize (undefined :: DataFrame t ns)+ , (# offA, n, arrA #) <- toBytes df+ , (# s1, arrM #) <- newByteArray# (n *# elS) s0+ , s2 <- copyByteArray# arrA (offA *# elS) arrM 0# (n *# elS) s1+ = (# s2, MDataFrame# 0# n arrM #)+{-# INLINE thawDataFrame# #-}++-- | Write a single element at the specified index+writeDataFrame# :: forall t (ns :: [Nat]) s+ . ( MutableFrame t ns, Dimensions ns )+ => MDataFrame s t ns -> Idx ns -> t -> State# s -> (# State# s, () #)+writeDataFrame# mdf ei x s | I# i <- fromEnum ei = (# writeDataFrameOff# mdf i x s, () #)+{-# INLINE writeDataFrame# #-}++-- | Read a single element at the specified index+readDataFrame# :: forall t (ns :: [Nat]) s+ . ( MutableFrame t ns, Dimensions ns )+ => MDataFrame s t ns -> Idx ns -> State# s -> (# State# s, t #)+readDataFrame# mdf ei | I# i <- fromEnum ei = readDataFrameOff# mdf i+{-# INLINE readDataFrame# #-}++class MutableFrame t (ns :: [Nat]) where+ -- | Write a single element at the specified element offset+ writeDataFrameOff# :: MDataFrame s t ns -> Int# -> t -> State# s -> State# s+ -- | Read a single element at the specified element offset+ readDataFrameOff# :: MDataFrame s t ns -> Int# -> State# s -> (# State# s, t #)++instance MutableFrame Float (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (F# x) = writeFloatArray# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readFloatArray# arr (off +# i) s0 = (# s1, F# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Double (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (D# x) = writeDoubleArray# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readDoubleArray# arr (off +# i) s0 = (# s1, D# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (I# x) = writeIntArray# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readIntArray# arr (off +# i) s0 = (# s1, I# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int8 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (I8# x) = writeInt8Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readInt8Array# arr (off +# i) s0 = (# s1, I8# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int16 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (I16# x) = writeInt16Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readInt16Array# arr (off +# i) s0 = (# s1, I16# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int32 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (I32# x) = writeInt32Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readInt32Array# arr (off +# i) s0 = (# s1, I32# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int64 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (I64# x) = writeInt64Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readInt64Array# arr (off +# i) s0 = (# s1, I64# r #)+ {-# INLINE readDataFrameOff# #-}+++instance MutableFrame Word (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (W# x) = writeWordArray# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readWordArray# arr (off +# i) s0 = (# s1, W# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word8 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (W8# x) = writeWord8Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readWord8Array# arr (off +# i) s0 = (# s1, W8# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word16 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (W16# x) = writeWord16Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readWord16Array# arr (off +# i) s0 = (# s1, W16# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word32 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (W32# x) = writeWord32Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readWord32Array# arr (off +# i) s0 = (# s1, W32# r #)+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word64 (ns :: [Nat]) where+ writeDataFrameOff# (MDataFrame# off _ arr) i (W64# x) = writeWord64Array# arr (off +# i) x+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# (MDataFrame# off _ arr) i s0+ | (# s1, r #) <- readWord64Array# arr (off +# i) s0 = (# s1, W64# r #)+ {-# INLINE readDataFrameOff# #-}
+ src-base/Numeric/Matrix/Mat44d.hs view
@@ -0,0 +1,37 @@+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Matrix.Mat44d () where+++import Numeric.Matrix.Class++notYet :: a+notYet = error "Sorry, this function is not implemented for current platform yet."++instance HomTransform4 Double where+ translate4 = notYet+ {-# INLINE translate4 #-}+ translate3 = notYet+ {-# INLINE translate3 #-}+ rotateX = notYet+ {-# INLINE rotateX #-}+ rotateY = notYet+ {-# INLINE rotateY #-}+ rotateZ = notYet+ {-# INLINE rotateZ #-}+ rotate = notYet+ {-# INLINE rotate #-}+ rotateEuler = notYet+ {-# INLINE rotateEuler #-}+ lookAt = notYet+ {-# INLINE lookAt #-}+ perspective = notYet+ {-# INLINE perspective #-}+ orthogonal = notYet+ {-# INLINE orthogonal #-}+ toHomPoint = notYet+ {-# INLINE toHomPoint #-}+ toHomVector = notYet+ {-# INLINE toHomVector #-}+ fromHom = notYet+ {-# INLINE fromHom #-}+
+ src-base/Numeric/Matrix/Mat44f.hs view
@@ -0,0 +1,37 @@+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Matrix.Mat44f () where+++import Numeric.Matrix.Class++notYet :: a+notYet = error "Sorry, this function is not implemented for current platform yet."++instance HomTransform4 Float where+ translate4 = notYet+ {-# INLINE translate4 #-}+ translate3 = notYet+ {-# INLINE translate3 #-}+ rotateX = notYet+ {-# INLINE rotateX #-}+ rotateY = notYet+ {-# INLINE rotateY #-}+ rotateZ = notYet+ {-# INLINE rotateZ #-}+ rotate = notYet+ {-# INLINE rotate #-}+ rotateEuler = notYet+ {-# INLINE rotateEuler #-}+ lookAt = notYet+ {-# INLINE lookAt #-}+ perspective = notYet+ {-# INLINE perspective #-}+ orthogonal = notYet+ {-# INLINE orthogonal #-}+ toHomPoint = notYet+ {-# INLINE toHomPoint #-}+ toHomVector = notYet+ {-# INLINE toHomVector #-}+ fromHom = notYet+ {-# INLINE fromHom #-}+
+ src-base/Numeric/Quaternion/QDouble.hs view
@@ -0,0 +1,542 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeSynonymInstances #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE UnboxedTuples #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Quaternion.QDouble+ ( QDouble, Quater (..)+ ) where++import GHC.Exts+import Data.Coerce (coerce)++import Numeric.Array+import Numeric.DataFrame.Type+import Numeric.Commons+import Numeric.Dimensions+import Numeric.Scalar+import Numeric.Vector+import Numeric.Matrix+import qualified Numeric.DataFrame.ST as ST+import qualified Numeric.Dimensions.Traverse.ST as ST+import qualified Control.Monad.ST as ST++import Numeric.Quaternion.Class+++type QDouble = Quater Double++instance Quaternion Double where+ newtype Quater Double = QDouble DoubleX4+ {-# INLINE packQ #-}+ packQ (D# x) (D# y) (D# z) (D# w) = QDouble (DoubleX4# x y z w)+ {-# INLINE unpackQ #-}+ unpackQ (QDouble (DoubleX4# x y z w)) = (D# x, D# y, D# z, D# w)+ {-# INLINE fromVecNum #-}+ fromVecNum (KnownDataFrame (DoubleX3# x y z)) (D# w) = QDouble (DoubleX4# x y z w)+ {-# INLINE fromVec4 #-}+ fromVec4 = coerce+ {-# INLINE toVec4 #-}+ toVec4 = coerce+ {-# INLINE square #-}+ square q = D# (qdot q)+ {-# INLINE im #-}+ im (QDouble (DoubleX4# x y z _)) = QDouble (DoubleX4# x y z 0.0##)+ {-# INLINE re #-}+ re (QDouble (DoubleX4# _ _ _ w)) = QDouble (DoubleX4# 0.0## 0.0## 0.0## w)+ {-# INLINE imVec #-}+ imVec (QDouble (DoubleX4# x y z _)) = KnownDataFrame (DoubleX3# x y z)+ {-# INLINE taker #-}+ taker (QDouble (DoubleX4# _ _ _ w)) = D# w+ {-# INLINE takei #-}+ takei (QDouble (DoubleX4# x _ _ _)) = D# x+ {-# INLINE takej #-}+ takej (QDouble (DoubleX4# _ y _ _)) = D# y+ {-# INLINE takek #-}+ takek (QDouble (DoubleX4# _ _ z _)) = D# z+ {-# INLINE conjugate #-}+ conjugate (QDouble (DoubleX4# x y z w)) = QDouble (DoubleX4#+ (negateDouble# x)+ (negateDouble# y)+ (negateDouble# z) w)+ {-# INLINE rotScale #-}+ rotScale (QDouble (DoubleX4# i j k t))+ (KnownDataFrame (DoubleX3# x y z))+ = let l = t*##t -## i*##i -## j*##j -## k*##k+ d = 2.0## *## ( i*##x +## j*##y +## k*##z)+ t2 = t *## 2.0##+ in KnownDataFrame+ ( DoubleX3#+ (l*##x +## d*##i +## t2 *## (z*##j -## y*##k))+ (l*##y +## d*##j +## t2 *## (x*##k -## z*##i))+ (l*##z +## d*##k +## t2 *## (y*##i -## x*##j))+ )+ {-# INLINE getRotScale #-}+ getRotScale _ (KnownDataFrame (DoubleX3# 0.0## 0.0## 0.0##))+ = QDouble (DoubleX4# 0.0## 0.0## 0.0## 0.0##)+ getRotScale (KnownDataFrame (DoubleX3# 0.0## 0.0## 0.0##)) _+ = case infty of D# x -> QDouble (DoubleX4# x x x x)+ getRotScale a@(KnownDataFrame (DoubleX3# a1 a2 a3))+ b@(KnownDataFrame (DoubleX3# b1 b2 b3))+ = let ma = sqrtDouble# (a1*##a1 +## a2*##a2 +## a3*##a3)+ mb = sqrtDouble# (b1*##b1 +## b2*##b2 +## b3*##b3)+ d = a1*##b1 +## a2*##b2 +## a3*##b3+ c = sqrtDouble# (ma*##mb +## d)+ ma2 = ma *## sqrtDouble# 2.0##+ r = 1.0## /## (ma2 *## c)+ in case cross a b of+ KnownDataFrame (DoubleX3# 0.0## 0.0## 0.0##) ->+ if isTrue# (d >## 0.0##)+ then QDouble (DoubleX4# 0.0## 0.0## 0.0## (sqrtDouble# (mb /## ma)))+ -- Shall we move result from k to i component?+ else QDouble (DoubleX4# 0.0## 0.0## (sqrtDouble# (mb /## ma)) 0.0##)+ KnownDataFrame (DoubleX3# t1 t2 t3) -> QDouble+ ( DoubleX4#+ (t1 *## r)+ (t2 *## r)+ (t3 *## r)+ (c /## ma2)+ )+ {-# INLINE axisRotation #-}+ axisRotation (KnownDataFrame (DoubleX3# 0.0## 0.0## 0.0##)) _+ = QDouble (DoubleX4# 0.0## 0.0## 0.0## 1.0##)+ axisRotation (KnownDataFrame (DoubleX3# x y z)) (D# a)+ = let c = cosDouble# (a *## 0.5##)+ s = sinDouble# (a *## 0.5##)+ /## sqrtDouble# (x*##x +## y*##y +## z*##z)+ in QDouble+ ( DoubleX4#+ (x *## s)+ (y *## s)+ (z *## s)+ c+ )+ {-# INLINE qArg #-}+ qArg (QDouble (DoubleX4# x y z w))+ = case atan2 (D# (sqrtDouble# (x*##x +## y*##y +## z*##z)))+ (D# w) of+ D# a -> D# (a *## 2.0##)+ {-# INLINE fromMatrix33 #-}+ fromMatrix33 m+ = let d =+ ( ix 0# m *## ( ix 4# m *## ix 8# m -## ix 5# m *## ix 7# m )+ -## ix 1# m *## ( ix 3# m *## ix 8# m -## ix 5# m *## ix 6# m )+ +## ix 2# m *## ( ix 3# m *## ix 7# m -## ix 4# m *## ix 6# m )+ ) **## 0.33333333333333333333333333333333##+ in QDouble+ ( DoubleX4#+ (sqrtDouble# (max# 0.0## (d +## ix 0# m -## ix 4# m -## ix 8# m )) *## sign# (ix 5# m -## ix 7# m) *## 0.5##)+ (sqrtDouble# (max# 0.0## (d -## ix 0# m +## ix 4# m -## ix 8# m )) *## sign# (ix 6# m -## ix 2# m) *## 0.5##)+ (sqrtDouble# (max# 0.0## (d -## ix 0# m -## ix 4# m +## ix 8# m )) *## sign# (ix 1# m -## ix 3# m) *## 0.5##)+ (sqrtDouble# (max# 0.0## (d +## ix 0# m +## ix 4# m +## ix 8# m )) *## 0.5##)+ )+ {-# INLINE fromMatrix44 #-}+ fromMatrix44 m+ = let d =+ ( ix 0# m *## ( ix 5# m *## ix 10# m -## ix 6# m *## ix 9# m )+ -## ix 1# m *## ( ix 4# m *## ix 10# m -## ix 6# m *## ix 8# m )+ +## ix 2# m *## ( ix 4# m *## ix 9# m -## ix 5# m *## ix 8# m )+ ) **## 0.33333333333333333333333333333333##+ c = 0.5## /## ix 15# m+ in QDouble+ ( DoubleX4#+ (sqrtDouble# (max# 0.0## (d +## ix 0# m -## ix 5# m -## ix 10# m )) *## sign# (ix 6# m -## ix 9# m) *## c)+ (sqrtDouble# (max# 0.0## (d -## ix 0# m +## ix 5# m -## ix 10# m )) *## sign# (ix 8# m -## ix 2# m) *## c)+ (sqrtDouble# (max# 0.0## (d -## ix 0# m -## ix 5# m +## ix 10# m )) *## sign# (ix 1# m -## ix 4# m) *## c)+ (sqrtDouble# (max# 0.0## (d +## ix 0# m +## ix 5# m +## ix 10# m )) *## c)+ )+ {-# INLINE toMatrix33 #-}+ toMatrix33 (QDouble (DoubleX4# 0.0## 0.0## 0.0## w)) = diag (scalar (D# (w *## w)))+ toMatrix33 (QDouble (DoubleX4# x' y' z' w')) =+ let x = scalar (D# x')+ y = scalar (D# y')+ z = scalar (D# z')+ w = scalar (D# w')+ x2 = x * x+ y2 = y * y+ z2 = z * z+ w2 = w * w+ l2 = x2 + y2 + z2 + w2+ in ST.runST $ do+ df <- ST.newDataFrame+ ST.writeDataFrameOff df 0 $ l2 - 2*(z2 + y2)+ ST.writeDataFrameOff df 1 $ 2*(x*y + z*w)+ ST.writeDataFrameOff df 2 $ 2*(x*z - y*w)+ ST.writeDataFrameOff df 3 $ 2*(x*y - z*w)+ ST.writeDataFrameOff df 4 $ l2 - 2*(z2 + x2)+ ST.writeDataFrameOff df 5 $ 2*(y*z + x*w)+ ST.writeDataFrameOff df 6 $ 2*(x*z + y*w)+ ST.writeDataFrameOff df 7 $ 2*(y*z - x*w)+ ST.writeDataFrameOff df 8 $ l2 - 2*(y2 + x2)+ ST.unsafeFreezeDataFrame df+ {-# INLINE toMatrix44 #-}+ toMatrix44 (QDouble (DoubleX4# 0.0## 0.0## 0.0## w)) = ST.runST $ do+ df <- ST.newDataFrame+ ST.overDimOff_ (dim :: Dim '[4,4]) (\i -> ST.writeDataFrameOff df (I# i) 0) 0# 1#+ let w2 = scalar (D# (w *## w))+ ST.writeDataFrameOff df 0 w2+ ST.writeDataFrameOff df 5 w2+ ST.writeDataFrameOff df 10 w2+ ST.writeDataFrameOff df 15 1+ ST.unsafeFreezeDataFrame df+ toMatrix44 (QDouble (DoubleX4# x' y' z' w')) =+ let x = scalar (D# x')+ y = scalar (D# y')+ z = scalar (D# z')+ w = scalar (D# w')+ x2 = x * x+ y2 = y * y+ z2 = z * z+ w2 = w * w+ l2 = x2 + y2 + z2 + w2+ in ST.runST $ do+ df <- ST.newDataFrame+ ST.writeDataFrameOff df 0 $ l2 - 2*(z2 + y2)+ ST.writeDataFrameOff df 1 $ 2*(x*y + z*w)+ ST.writeDataFrameOff df 2 $ 2*(x*z - y*w)+ ST.writeDataFrameOff df 3 0+ ST.writeDataFrameOff df 4 $ 2*(x*y - z*w)+ ST.writeDataFrameOff df 5 $ l2 - 2*(z2 + x2)+ ST.writeDataFrameOff df 6 $ 2*(y*z + x*w)+ ST.writeDataFrameOff df 7 0+ ST.writeDataFrameOff df 8 $ 2*(x*z + y*w)+ ST.writeDataFrameOff df 9 $ 2*(y*z - x*w)+ ST.writeDataFrameOff df 10 $ l2 - 2*(y2 + x2)+ ST.writeDataFrameOff df 11 0+ ST.writeDataFrameOff df 12 0+ ST.writeDataFrameOff df 13 0+ ST.writeDataFrameOff df 14 0+ ST.writeDataFrameOff df 15 1+ ST.unsafeFreezeDataFrame df++qdot :: QDouble -> Double#+qdot (QDouble (DoubleX4# x y z w)) = (x *## x) +##+ (y *## y) +##+ (z *## z) +##+ (w *## w)+{-# INLINE qdot #-}++infty :: Double+infty = read "Infinity"++max# :: Double# -> Double# -> Double#+max# a b | isTrue# (a >## b) = a+ | otherwise = b+{-# INLINE max# #-}++sign# :: Double# -> Double#+sign# a | isTrue# (a >## 0.0##) = 1.0##+ | isTrue# (a <## 0.0##) = negateDouble# 1.0##+ | otherwise = 0.0##+{-# INLINE sign# #-}++--------------------------------------------------------------------------+-- Num+--------------------------------------------------------------------------++instance Num QDouble where+ QDouble a + QDouble b+ = QDouble (a + b)+ {-# INLINE (+) #-}+ QDouble a - QDouble b+ = QDouble (a - b)+ {-# INLINE (-) #-}+ QDouble (DoubleX4# a1 a2 a3 a4) * QDouble (DoubleX4# b1 b2 b3 b4)+ = QDouble+ ( DoubleX4#+ ((a4 *## b1) +##+ (a1 *## b4) +##+ (a2 *## b3) -##+ (a3 *## b2)+ )+ ((a4 *## b2) -##+ (a1 *## b3) +##+ (a2 *## b4) +##+ (a3 *## b1)+ )+ ((a4 *## b3) +##+ (a1 *## b2) -##+ (a2 *## b1) +##+ (a3 *## b4)+ )+ ((a4 *## b4) -##+ (a1 *## b1) -##+ (a2 *## b2) -##+ (a3 *## b3)+ )+ )+ {-# INLINE (*) #-}+ negate (QDouble a) = QDouble (negate a)+ {-# INLINE negate #-}+ abs q = QDouble (DoubleX4# 0.0## 0.0## 0.0## (sqrtDouble# (qdot q)))+ {-# INLINE abs #-}+ signum q@(QDouble (DoubleX4# x y z w))+ = case qdot q of+ 0.0## -> QDouble (DoubleX4# 0.0## 0.0## 0.0## 0.0##)+ qd -> case 1.0## /## sqrtDouble# qd of+ s -> QDouble+ ( DoubleX4#+ (x *## s)+ (y *## s)+ (z *## s)+ (w *## s)+ )+ {-# INLINE signum #-}+ fromInteger n = case fromInteger n of+ D# x -> QDouble (DoubleX4# 0.0## 0.0## 0.0## x)+ {-# INLINE fromInteger #-}++++--------------------------------------------------------------------------+-- Fractional+--------------------------------------------------------------------------++instance Fractional QDouble where+ {-# INLINE recip #-}+ recip q@(QDouble (DoubleX4# x y z w)) = case -1.0## /## qdot q of+ c -> QDouble+ ( DoubleX4#+ (x *## c)+ (y *## c)+ (z *## c)+ (negateDouble# (w *## c))+ )+ {-# INLINE (/) #-}+ a / b = a * recip b+ {-# INLINE fromRational #-}+ fromRational q = case fromRational q of+ D# x -> QDouble (DoubleX4# 0.0## 0.0## 0.0## x)++--------------------------------------------------------------------------+-- Doubleing+--------------------------------------------------------------------------++instance Floating QDouble where+ {-# INLINE pi #-}+ pi = QDouble (DoubleX4# 0.0## 0.0## 0.0##+ 3.141592653589793##)+ {-# INLINE exp #-}+ exp (QDouble (DoubleX4# x y z w))+ = case (# (x *## x) +##+ (y *## y) +##+ (z *## z)+ , expDouble# w+ #) of+ (# 0.0##, et #) -> QDouble (DoubleX4# 0.0## 0.0## 0.0## et)+ (# mv2, et #) -> case sqrtDouble# mv2 of+ mv -> case et *## sinDouble# mv+ /## mv of+ l -> QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (et *## cosDouble# mv)+ )+ {-# INLINE log #-}+ log (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> if isTrue# (w >=## 0.0##)+ then QDouble (DoubleX4# 0.0## 0.0## 0.0## (logDouble# w))+ else QDouble (DoubleX4# 3.141592653589793## 0.0## 0.0##+ (logDouble# (negateDouble# w)))+ mv2 -> case (# sqrtDouble# (mv2 +## (w *## w))+ , sqrtDouble# mv2+ #) of+ (# mq, mv #) -> case atan2 (D# mv) (D# w) / D# mv of+ D# l -> QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (logDouble# mq)+ )+ {-# INLINE sqrt #-}+ sqrt (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> if isTrue# (w >=## 0.0##)+ then QDouble (DoubleX4# 0.0## 0.0## 0.0## (sqrtDouble# w))+ else QDouble (DoubleX4# (sqrtDouble# (negateDouble# w)) 0.0## 0.0## 0.0##)+ mv2 ->+ let mq = sqrtDouble# (mv2 +## w *## w)+ l2 = sqrtDouble# mq+ tq = w /## (mq *## 2.0##)+ sina = sqrtDouble# (0.5## -## tq) *## l2 /## sqrtDouble# mv2+ in QDouble+ ( DoubleX4#+ (x *## sina)+ (y *## sina)+ (z *## sina)+ (sqrtDouble# (0.5## +## tq) *## l2)+ )+ {-# INLINE sin #-}+ sin (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> QDouble (DoubleX4# 0.0## 0.0## 0.0## (sinDouble# w))+ mv2 -> case sqrtDouble# mv2 of+ mv -> case cosDouble# w *## sinhDouble# mv+ /## mv of+ l -> QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (sinDouble# w *## coshDouble# mv)+ )+ {-# INLINE cos #-}+ cos (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> QDouble (DoubleX4# 0.0## 0.0## 0.0## (cosDouble# w))+ mv2 -> case sqrtDouble# mv2 of+ mv -> case sinDouble# w *## sinhDouble# mv+ /## negateDouble# mv of+ l -> QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (cosDouble# w *## coshDouble# mv)+ )+ {-# INLINE tan #-}+ tan (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> QDouble (DoubleX4# 0.0## 0.0## 0.0## (tanDouble# w))+ mv2 ->+ let mv = sqrtDouble# mv2+ chv = coshDouble# mv+ shv = sinhDouble# mv+ ct = cosDouble# w+ st = sinDouble# w+ cq = 1.0## /##+ ( (ct *## ct *## chv *## chv)+ +##+ (st *## st *## shv *## shv)+ )+ l = chv *## shv *## cq+ /## mv+ in QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (ct *## st *## cq)+ )+ {-# INLINE sinh #-}+ sinh (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> QDouble (DoubleX4# 0.0## 0.0## 0.0## (sinhDouble# w))+ mv2 -> case sqrtDouble# mv2 of+ mv -> case coshDouble# w *## sinDouble# mv+ /## mv of+ l -> QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (sinhDouble# w *## cosDouble# mv)+ )+ {-# INLINE cosh #-}+ cosh (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> QDouble (DoubleX4# 0.0## 0.0## 0.0## (coshDouble# w))+ mv2 -> case sqrtDouble# mv2 of+ mv -> case sinhDouble# w *## sinDouble# mv+ /## mv of+ l -> QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (coshDouble# w *## cosDouble# mv)+ )+ {-# INLINE tanh #-}+ tanh (QDouble (DoubleX4# x y z w))+ = case (x *## x) +##+ (y *## y) +##+ (z *## z) of+ 0.0## -> QDouble (DoubleX4# 0.0## 0.0## 0.0## (tanhDouble# w))+ mv2 ->+ let mv = sqrtDouble# mv2+ cv = cosDouble# mv+ sv = sinDouble# mv+ cht = coshDouble# w+ sht = sinhDouble# w+ cq = 1.0## /##+ ( (cht *## cht *## cv *## cv)+ +##+ (sht *## sht *## sv *## sv)+ )+ l = cv *## sv *## cq+ /## mv+ in QDouble+ ( DoubleX4#+ (x *## l)+ (y *## l)+ (z *## l)+ (cht *## sht *## cq)+ )+ {-# INLINE asin #-}+ asin q = -i * log (i*q + sqrt (1 - q*q))+ where+ i = case signum . im $ q of+ 0 -> QDouble (DoubleX4# 1.0## 0.0## 0.0## 0.0##)+ i' -> i'+ {-# INLINE acos #-}+ acos q = pi/2 - asin q+ {-# INLINE atan #-}+ atan q@(QDouble (DoubleX4# _ _ _ w))+ = if square imq == 0+ then QDouble (DoubleX4# 0.0## 0.0## 0.0## (atanDouble# w))+ else i / 2 * log ( (i + q) / (i - q) )+ where+ i = signum imq+ imq = im q+ {-# INLINE asinh #-}+ asinh q = log (q + sqrt (q*q + 1))+ {-# INLINE acosh #-}+ acosh q = log (q + sqrt (q*q - 1))+ {-# INLINE atanh #-}+ atanh q = 0.5 * log ((1+q)/(1-q))++--------------------------------------------------------------------------+-- Eq+--------------------------------------------------------------------------++instance Eq QDouble where+ {-# INLINE (==) #-}+ QDouble a == QDouble b = a == b+ {-# INLINE (/=) #-}+ QDouble a /= QDouble b = a /= b++++--------------------------------------------------------------------------+-- Show+--------------------------------------------------------------------------++instance Show QDouble where+ show (QDouble (DoubleX4# x y z w)) =+ show (D# w) ++ ss x ++ "i"+ ++ ss y ++ "j"+ ++ ss z ++ "k"+ where+ ss a# = case D# a# of+ a -> if a >= 0 then " + " ++ show a+ else " - " ++ show (negate a)
+ src-base/Numeric/Quaternion/QFloat.hs view
@@ -0,0 +1,562 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeSynonymInstances #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE UnboxedTuples #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Quaternion.QFloat+ ( QFloat, Quater (..)+ ) where++import GHC.Exts+import Data.Coerce (coerce)++import Numeric.Array+import Numeric.DataFrame.Type+import Numeric.Commons+import Numeric.Dimensions+import Numeric.Scalar+import Numeric.Vector+import Numeric.Matrix+import qualified Numeric.DataFrame.ST as ST+import qualified Numeric.Dimensions.Traverse.ST as ST+import qualified Control.Monad.ST as ST++import Numeric.Quaternion.Class+++type QFloat = Quater Float++instance Quaternion Float where+ newtype Quater Float = QFloat FloatX4+ {-# INLINE packQ #-}+ packQ (F# x) (F# y) (F# z) (F# w) = QFloat (FloatX4# x y z w)+ {-# INLINE unpackQ #-}+ unpackQ (QFloat (FloatX4# x y z w)) = (F# x, F# y, F# z, F# w)+ {-# INLINE fromVecNum #-}+ fromVecNum (KnownDataFrame (FloatX3# x y z)) (F# w) = QFloat (FloatX4# x y z w)+ {-# INLINE fromVec4 #-}+ fromVec4 = coerce+ {-# INLINE toVec4 #-}+ toVec4 = coerce+ {-# INLINE square #-}+ square q = F# (qdot q)+ {-# INLINE im #-}+ im (QFloat (FloatX4# x y z _)) = QFloat (FloatX4# x y z 0.0#)+ {-# INLINE re #-}+ re (QFloat (FloatX4# _ _ _ w)) = QFloat (FloatX4# 0.0# 0.0# 0.0# w)+ {-# INLINE imVec #-}+ imVec (QFloat (FloatX4# x y z _)) = KnownDataFrame (FloatX3# x y z)+ {-# INLINE taker #-}+ taker (QFloat (FloatX4# _ _ _ w)) = F# w+ {-# INLINE takei #-}+ takei (QFloat (FloatX4# x _ _ _)) = F# x+ {-# INLINE takej #-}+ takej (QFloat (FloatX4# _ y _ _)) = F# y+ {-# INLINE takek #-}+ takek (QFloat (FloatX4# _ _ z _)) = F# z+ {-# INLINE conjugate #-}+ conjugate (QFloat (FloatX4# x y z w)) = QFloat (FloatX4#+ (negateFloat# x)+ (negateFloat# y)+ (negateFloat# z) w)+ {-# INLINE rotScale #-}+ rotScale (QFloat (FloatX4# i j k t))+ (KnownDataFrame (FloatX3# x y z))+ = let l = t*%t -% i*%i -% j*%j -% k*%k+ d = 2.0# *% ( i*%x +% j*%y +% k*%z)+ t2 = t *% 2.0#+ in KnownDataFrame+ ( FloatX3#+ (l*%x +% d*%i +% t2 *% (z*%j -% y*%k))+ (l*%y +% d*%j +% t2 *% (x*%k -% z*%i))+ (l*%z +% d*%k +% t2 *% (y*%i -% x*%j))+ )+ {-# INLINE getRotScale #-}+ getRotScale _ (KnownDataFrame (FloatX3# 0.0# 0.0# 0.0#))+ = QFloat (FloatX4# 0.0# 0.0# 0.0# 0.0#)+ getRotScale (KnownDataFrame (FloatX3# 0.0# 0.0# 0.0#)) _+ = case infty of F# x -> QFloat (FloatX4# x x x x)+ getRotScale a@(KnownDataFrame (FloatX3# a1 a2 a3))+ b@(KnownDataFrame (FloatX3# b1 b2 b3))+ = let ma = sqrtFloat# (a1*%a1 +% a2*%a2 +% a3*%a3)+ mb = sqrtFloat# (b1*%b1 +% b2*%b2 +% b3*%b3)+ d = a1*%b1 +% a2*%b2 +% a3*%b3+ c = sqrtFloat# (ma*%mb +% d)+ ma2 = ma *% sqrtFloat# 2.0#+ r = 1.0# /% (ma2 *% c)+ in case cross a b of+ KnownDataFrame (FloatX3# 0.0# 0.0# 0.0#) ->+ if isTrue# (gtFloat# d 0.0#)+ then QFloat (FloatX4# 0.0# 0.0# 0.0# (sqrtFloat# (mb /% ma)))+ -- Shall we move result from k to i component?+ else QFloat (FloatX4# 0.0# 0.0# (sqrtFloat# (mb /% ma)) 0.0#)+ KnownDataFrame (FloatX3# t1 t2 t3) -> QFloat+ ( FloatX4#+ (t1 *% r)+ (t2 *% r)+ (t3 *% r)+ (c /% ma2)+ )+ {-# INLINE axisRotation #-}+ axisRotation (KnownDataFrame (FloatX3# 0.0# 0.0# 0.0#)) _+ = QFloat (FloatX4# 0.0# 0.0# 0.0# 1.0#)+ axisRotation (KnownDataFrame (FloatX3# x y z)) (F# a)+ = let c = cosFloat# (a *% 0.5#)+ s = sinFloat# (a *% 0.5#)+ /% sqrtFloat# (x*%x +% y*%y +% z*%z)+ in QFloat+ ( FloatX4#+ (x *% s)+ (y *% s)+ (z *% s)+ c+ )+ {-# INLINE qArg #-}+ qArg (QFloat (FloatX4# x y z w))+ = case atan2 (F# (sqrtFloat# (x*%x +% y*%y +% z*%z)))+ (F# w) of+ F# a -> F# (a *% 2.0#)+ {-# INLINE fromMatrix33 #-}+ fromMatrix33 m+ = let d = powerFloat#+ ( ix 0# m *% ( ix 4# m *% ix 8# m -% ix 5# m *% ix 7# m )+ -% ix 1# m *% ( ix 3# m *% ix 8# m -% ix 5# m *% ix 6# m )+ +% ix 2# m *% ( ix 3# m *% ix 7# m -% ix 4# m *% ix 6# m )+ ) 0.33333333333333333333333333333333#+ in QFloat+ ( FloatX4#+ (sqrtFloat# (max# 0.0# (d +% ix 0# m -% ix 4# m -% ix 8# m )) *% sign# (ix 5# m -% ix 7# m) *% 0.5#)+ (sqrtFloat# (max# 0.0# (d -% ix 0# m +% ix 4# m -% ix 8# m )) *% sign# (ix 6# m -% ix 2# m) *% 0.5#)+ (sqrtFloat# (max# 0.0# (d -% ix 0# m -% ix 4# m +% ix 8# m )) *% sign# (ix 1# m -% ix 3# m) *% 0.5#)+ (sqrtFloat# (max# 0.0# (d +% ix 0# m +% ix 4# m +% ix 8# m )) *% 0.5#)+ )+ {-# INLINE fromMatrix44 #-}+ fromMatrix44 m+ = let d = powerFloat#+ ( ix 0# m *% ( ix 5# m *% ix 10# m -% ix 6# m *% ix 9# m )+ -% ix 1# m *% ( ix 4# m *% ix 10# m -% ix 6# m *% ix 8# m )+ +% ix 2# m *% ( ix 4# m *% ix 9# m -% ix 5# m *% ix 8# m )+ ) 0.33333333333333333333333333333333#+ c = 0.5# /% ix 15# m+ in QFloat+ ( FloatX4#+ (sqrtFloat# (max# 0.0# (d +% ix 0# m -% ix 5# m -% ix 10# m )) *% sign# (ix 6# m -% ix 9# m) *% c)+ (sqrtFloat# (max# 0.0# (d -% ix 0# m +% ix 5# m -% ix 10# m )) *% sign# (ix 8# m -% ix 2# m) *% c)+ (sqrtFloat# (max# 0.0# (d -% ix 0# m -% ix 5# m +% ix 10# m )) *% sign# (ix 1# m -% ix 4# m) *% c)+ (sqrtFloat# (max# 0.0# (d +% ix 0# m +% ix 5# m +% ix 10# m )) *% c)+ )+ {-# INLINE toMatrix33 #-}+ toMatrix33 (QFloat (FloatX4# 0.0# 0.0# 0.0# w)) = diag (scalar (F# (w *% w)))+ toMatrix33 (QFloat (FloatX4# x' y' z' w')) =+ let x = scalar (F# x')+ y = scalar (F# y')+ z = scalar (F# z')+ w = scalar (F# w')+ x2 = x * x+ y2 = y * y+ z2 = z * z+ w2 = w * w+ l2 = x2 + y2 + z2 + w2+ in ST.runST $ do+ df <- ST.newDataFrame+ ST.writeDataFrameOff df 0 $ l2 - 2*(z2 + y2)+ ST.writeDataFrameOff df 1 $ 2*(x*y + z*w)+ ST.writeDataFrameOff df 2 $ 2*(x*z - y*w)+ ST.writeDataFrameOff df 3 $ 2*(x*y - z*w)+ ST.writeDataFrameOff df 4 $ l2 - 2*(z2 + x2)+ ST.writeDataFrameOff df 5 $ 2*(y*z + x*w)+ ST.writeDataFrameOff df 6 $ 2*(x*z + y*w)+ ST.writeDataFrameOff df 7 $ 2*(y*z - x*w)+ ST.writeDataFrameOff df 8 $ l2 - 2*(y2 + x2)+ ST.unsafeFreezeDataFrame df+ {-# INLINE toMatrix44 #-}+ toMatrix44 (QFloat (FloatX4# 0.0# 0.0# 0.0# w)) = ST.runST $ do+ df <- ST.newDataFrame+ ST.overDimOff_ (dim :: Dim '[4,4]) (\i -> ST.writeDataFrameOff df (I# i) 0) 0# 1#+ let w2 = scalar (F# (w *% w))+ ST.writeDataFrameOff df 0 w2+ ST.writeDataFrameOff df 5 w2+ ST.writeDataFrameOff df 10 w2+ ST.writeDataFrameOff df 15 1+ ST.unsafeFreezeDataFrame df+ toMatrix44 (QFloat (FloatX4# x' y' z' w')) =+ let x = scalar (F# x')+ y = scalar (F# y')+ z = scalar (F# z')+ w = scalar (F# w')+ x2 = x * x+ y2 = y * y+ z2 = z * z+ w2 = w * w+ l2 = x2 + y2 + z2 + w2+ in ST.runST $ do+ df <- ST.newDataFrame+ ST.writeDataFrameOff df 0 $ l2 - 2*(z2 + y2)+ ST.writeDataFrameOff df 1 $ 2*(x*y + z*w)+ ST.writeDataFrameOff df 2 $ 2*(x*z - y*w)+ ST.writeDataFrameOff df 3 0+ ST.writeDataFrameOff df 4 $ 2*(x*y - z*w)+ ST.writeDataFrameOff df 5 $ l2 - 2*(z2 + x2)+ ST.writeDataFrameOff df 6 $ 2*(y*z + x*w)+ ST.writeDataFrameOff df 7 0+ ST.writeDataFrameOff df 8 $ 2*(x*z + y*w)+ ST.writeDataFrameOff df 9 $ 2*(y*z - x*w)+ ST.writeDataFrameOff df 10 $ l2 - 2*(y2 + x2)+ ST.writeDataFrameOff df 11 0+ ST.writeDataFrameOff df 12 0+ ST.writeDataFrameOff df 13 0+ ST.writeDataFrameOff df 14 0+ ST.writeDataFrameOff df 15 1+ ST.unsafeFreezeDataFrame df++qdot :: QFloat -> Float#+qdot (QFloat (FloatX4# x y z w)) = (x *% x) +%+ (y *% y) +%+ (z *% z) +%+ (w *% w)+{-# INLINE qdot #-}++(*%) :: Float# -> Float# -> Float#+(*%) = timesFloat#+{-# INLINE (*%) #-}+infixl 7 *%++(-%) :: Float# -> Float# -> Float#+(-%) = minusFloat#+{-# INLINE (-%) #-}+infixl 6 -%++(+%) :: Float# -> Float# -> Float#+(+%) = plusFloat#+{-# INLINE (+%) #-}+infixl 6 +%++(/%) :: Float# -> Float# -> Float#+(/%) = divideFloat#+{-# INLINE (/%) #-}+infixl 7 /%++infty :: Float+infty = read "Infinity"++max# :: Float# -> Float# -> Float#+max# a b | isTrue# (gtFloat# a b) = a+ | otherwise = b+{-# INLINE max# #-}++sign# :: Float# -> Float#+sign# a | isTrue# (gtFloat# a 0.0#) = 1.0#+ | isTrue# (ltFloat# a 0.0#) = negateFloat# 1.0#+ | otherwise = 0.0#+{-# INLINE sign# #-}++--------------------------------------------------------------------------+-- Num+--------------------------------------------------------------------------++instance Num QFloat where+ QFloat a + QFloat b+ = QFloat (a + b)+ {-# INLINE (+) #-}+ QFloat a - QFloat b+ = QFloat (a - b)+ {-# INLINE (-) #-}+ QFloat (FloatX4# a1 a2 a3 a4) * QFloat (FloatX4# b1 b2 b3 b4)+ = QFloat+ ( FloatX4#+ ((a4 *% b1) +%+ (a1 *% b4) +%+ (a2 *% b3) -%+ (a3 *% b2)+ )+ ((a4 *% b2) -%+ (a1 *% b3) +%+ (a2 *% b4) +%+ (a3 *% b1)+ )+ ((a4 *% b3) +%+ (a1 *% b2) -%+ (a2 *% b1) +%+ (a3 *% b4)+ )+ ((a4 *% b4) -%+ (a1 *% b1) -%+ (a2 *% b2) -%+ (a3 *% b3)+ )+ )+ {-# INLINE (*) #-}+ negate (QFloat a) = QFloat (negate a)+ {-# INLINE negate #-}+ abs q = QFloat (FloatX4# 0.0# 0.0# 0.0# (sqrtFloat# (qdot q)))+ {-# INLINE abs #-}+ signum q@(QFloat (FloatX4# x y z w))+ = case qdot q of+ 0.0# -> QFloat (FloatX4# 0.0# 0.0# 0.0# 0.0#)+ qd -> case 1.0# /% sqrtFloat# qd of+ s -> QFloat+ ( FloatX4#+ (x *% s)+ (y *% s)+ (z *% s)+ (w *% s)+ )+ {-# INLINE signum #-}+ fromInteger n = case fromInteger n of+ F# x -> QFloat (FloatX4# 0.0# 0.0# 0.0# x)+ {-# INLINE fromInteger #-}++++--------------------------------------------------------------------------+-- Fractional+--------------------------------------------------------------------------++instance Fractional QFloat where+ {-# INLINE recip #-}+ recip q@(QFloat (FloatX4# x y z w)) = case -1.0# /% qdot q of+ c -> QFloat+ ( FloatX4#+ (x *% c)+ (y *% c)+ (z *% c)+ (negateFloat# (w *% c))+ )+ {-# INLINE (/) #-}+ a / b = a * recip b+ {-# INLINE fromRational #-}+ fromRational q = case fromRational q of+ F# x -> QFloat (FloatX4# 0.0# 0.0# 0.0# x)++--------------------------------------------------------------------------+-- Floating+--------------------------------------------------------------------------++instance Floating QFloat where+ {-# INLINE pi #-}+ pi = QFloat (FloatX4# 0.0# 0.0# 0.0#+ 3.141592653589793#)+ {-# INLINE exp #-}+ exp (QFloat (FloatX4# x y z w))+ = case (# (x *% x) +%+ (y *% y) +%+ (z *% z)+ , expFloat# w+ #) of+ (# 0.0#, et #) -> QFloat (FloatX4# 0.0# 0.0# 0.0# et)+ (# mv2, et #) -> case sqrtFloat# mv2 of+ mv -> case et *% sinFloat# mv+ /% mv of+ l -> QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (et *% cosFloat# mv)+ )+ {-# INLINE log #-}+ log (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> if isTrue# (w `geFloat#` 0.0#)+ then QFloat (FloatX4# 0.0# 0.0# 0.0# (logFloat# w))+ else QFloat (FloatX4# 3.141592653589793# 0.0# 0.0#+ (logFloat# (negateFloat# w)))+ mv2 -> case (# sqrtFloat# (mv2 +% (w *% w))+ , sqrtFloat# mv2+ #) of+ (# mq, mv #) -> case atan2 (F# mv) (F# w) / F# mv of+ F# l -> QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (logFloat# mq)+ )+ {-# INLINE sqrt #-}+ sqrt (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> if isTrue# (w `geFloat#` 0.0#)+ then QFloat (FloatX4# 0.0# 0.0# 0.0# (sqrtFloat# w))+ else QFloat (FloatX4# (sqrtFloat# (negateFloat# w)) 0.0# 0.0# 0.0#)+ mv2 ->+ let mq = sqrtFloat# (mv2 +% w *% w)+ l2 = sqrtFloat# mq+ tq = w /% (mq *% 2.0#)+ sina = sqrtFloat# (0.5# -% tq) *% l2 /% sqrtFloat# mv2+ in QFloat+ ( FloatX4#+ (x *% sina)+ (y *% sina)+ (z *% sina)+ (sqrtFloat# (0.5# +% tq) *% l2)+ )+ {-# INLINE sin #-}+ sin (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> QFloat (FloatX4# 0.0# 0.0# 0.0# (sinFloat# w))+ mv2 -> case sqrtFloat# mv2 of+ mv -> case cosFloat# w *% sinhFloat# mv+ /% mv of+ l -> QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (sinFloat# w *% coshFloat# mv)+ )+ {-# INLINE cos #-}+ cos (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> QFloat (FloatX4# 0.0# 0.0# 0.0# (cosFloat# w))+ mv2 -> case sqrtFloat# mv2 of+ mv -> case sinFloat# w *% sinhFloat# mv+ /% negateFloat# mv of+ l -> QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (cosFloat# w *% coshFloat# mv)+ )+ {-# INLINE tan #-}+ tan (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> QFloat (FloatX4# 0.0# 0.0# 0.0# (tanFloat# w))+ mv2 ->+ let mv = sqrtFloat# mv2+ chv = coshFloat# mv+ shv = sinhFloat# mv+ ct = cosFloat# w+ st = sinFloat# w+ cq = 1.0# /%+ ( (ct *% ct *% chv *% chv)+ +%+ (st *% st *% shv *% shv)+ )+ l = chv *% shv *% cq+ /% mv+ in QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (ct *% st *% cq)+ )+ {-# INLINE sinh #-}+ sinh (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> QFloat (FloatX4# 0.0# 0.0# 0.0# (sinhFloat# w))+ mv2 -> case sqrtFloat# mv2 of+ mv -> case coshFloat# w *% sinFloat# mv+ /% mv of+ l -> QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (sinhFloat# w *% cosFloat# mv)+ )+ {-# INLINE cosh #-}+ cosh (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> QFloat (FloatX4# 0.0# 0.0# 0.0# (coshFloat# w))+ mv2 -> case sqrtFloat# mv2 of+ mv -> case sinhFloat# w *% sinFloat# mv+ /% mv of+ l -> QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (coshFloat# w *% cosFloat# mv)+ )+ {-# INLINE tanh #-}+ tanh (QFloat (FloatX4# x y z w))+ = case (x *% x) +%+ (y *% y) +%+ (z *% z) of+ 0.0# -> QFloat (FloatX4# 0.0# 0.0# 0.0# (tanhFloat# w))+ mv2 ->+ let mv = sqrtFloat# mv2+ cv = cosFloat# mv+ sv = sinFloat# mv+ cht = coshFloat# w+ sht = sinhFloat# w+ cq = 1.0# /%+ ( (cht *% cht *% cv *% cv)+ +%+ (sht *% sht *% sv *% sv)+ )+ l = cv *% sv *% cq+ /% mv+ in QFloat+ ( FloatX4#+ (x *% l)+ (y *% l)+ (z *% l)+ (cht *% sht *% cq)+ )+ {-# INLINE asin #-}+ asin q = -i * log (i*q + sqrt (1 - q*q))+ where+ i = case signum . im $ q of+ 0 -> QFloat (FloatX4# 1.0# 0.0# 0.0# 0.0#)+ i' -> i'+ {-# INLINE acos #-}+ acos q = pi/2 - asin q+ {-# INLINE atan #-}+ atan q@(QFloat (FloatX4# _ _ _ w))+ = if square imq == 0+ then QFloat (FloatX4# 0.0# 0.0# 0.0# (atanFloat# w))+ else i / 2 * log ( (i + q) / (i - q) )+ where+ i = signum imq+ imq = im q+ {-# INLINE asinh #-}+ asinh q = log (q + sqrt (q*q + 1))+ {-# INLINE acosh #-}+ acosh q = log (q + sqrt (q*q - 1))+ {-# INLINE atanh #-}+ atanh q = 0.5 * log ((1+q)/(1-q))++--------------------------------------------------------------------------+-- Eq+--------------------------------------------------------------------------++instance Eq QFloat where+ {-# INLINE (==) #-}+ QFloat a == QFloat b = a == b+ {-# INLINE (/=) #-}+ QFloat a /= QFloat b = a /= b++++--------------------------------------------------------------------------+-- Show+--------------------------------------------------------------------------++instance Show QFloat where+ show (QFloat (FloatX4# x y z w)) =+ show (F# w) ++ ss x ++ "i"+ ++ ss y ++ "j"+ ++ ss z ++ "k"+ where+ ss a# = case F# a# of+ a -> if a >= 0 then " + " ++ show a+ else " - " ++ show (negate a)
src-ghcjs/Numeric/Array/Family.hs view
@@ -14,6 +14,8 @@ {-# LANGUAGE UnboxedTuples #-} {-# LANGUAGE StandaloneDeriving #-} {-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE UnliftedFFITypes #-} ----------------------------------------------------------------------------- -- | -- Module : Numeric.Array.Family@@ -38,7 +40,7 @@ import Data.Type.Equality ((:~:) (..)) import Data.Word (Word16, Word32, Word8) import GHC.Prim (Double#, Float#, Int#,- Word#, unsafeCoerce#)+ Word#, unsafeCoerce#, ByteArray#) import GHC.Types (Int (..)) import GHCJS.Types @@ -74,9 +76,9 @@ type instance ElemRep Word8Clamped = ElemRep Int type instance ElemPrim Word8Clamped = Int# instance PrimBytes Word8Clamped where- toBytes (Clamped i) = toBytes (fromIntegral (min 0 (max 255 i)) :: Word8)+ toBytes v = (# 0#, 1#, js_wrapWord8Clamped v #) {-# INLINE toBytes #-}- fromBytes bs = fromIntegral (fromBytes bs :: Word8)+ fromBytes (# off, _, arr #) = js_unwrapWord8Clamped arr off {-# INLINE fromBytes #-} byteSize _ = 1# {-# INLINE byteSize #-}@@ -86,6 +88,8 @@ {-# INLINE elementByteSize #-} ix _ (Clamped (I# x)) = x {-# INLINE ix #-}+foreign import javascript unsafe "h$wrapBuffer((new Uint8ClampedArray([$1])).buffer)" js_wrapWord8Clamped :: Word8Clamped -> ByteArray#+foreign import javascript unsafe "($1.uc || new Uint8ClampedArray($1.buf))[$2]" js_unwrapWord8Clamped :: ByteArray# -> Int# -> Word8Clamped instance ElementWise (Idx ('[] :: [Nat])) Word8Clamped Word8Clamped where indexOffset# x _ = x
src-ghcjs/Numeric/Array/Family/ArrayT.hs view
@@ -15,7 +15,6 @@ {-# LANGUAGE JavaScriptFFI #-} {-# LANGUAGE GHCForeignImportPrim #-} {-# LANGUAGE UnliftedFFITypes #-}-{-# LANGUAGE GeneralizedNewtypeDeriving #-} {-# LANGUAGE Strict #-} {-# OPTIONS_GHC -fno-warn-orphans #-} module Numeric.Array.Family.ArrayT () where@@ -35,7 +34,7 @@ import Numeric.Dimensions import Numeric.Dimensions.Traverse import Numeric.TypeLits-import Numeric.Matrix.Type+import Numeric.Matrix.Class type instance ElemRep (ArrayT t ds) = ElemRep t@@ -255,11 +254,15 @@ 1 -> GT 0 -> EQ _ -> LT+ max = js_arrayMax+ min = js_arrayMin foreign import javascript unsafe "$1.every(function (e, i) { return e < $2[i]; })" js_arrayTLT :: ArrayT t ds -> ArrayT t ds -> Bool foreign import javascript unsafe "$1.every(function (e, i) { return e <= $2[i]; })" js_arrayTLE :: ArrayT t ds -> ArrayT t ds -> Bool foreign import javascript unsafe "$1.every(function (e, i) { return e > $2[i]; })" js_arrayTGT :: ArrayT t ds -> ArrayT t ds -> Bool foreign import javascript unsafe "$1.every(function (e, i) { return e >= $2[i]; })" js_arrayTGE :: ArrayT t ds -> ArrayT t ds -> Bool foreign import javascript unsafe "$1.reduce(function (r, e, i) { return r === 0 ? (e > $2[i] ? 1 : (e < $2[i] ? -1 : 0)) : r;}, 0)" js_arrayTCmp :: ArrayT t ds -> ArrayT t ds -> Int+foreign import javascript unsafe "$1.map(function (e, i) { return Math.max(e,$2[i]); })" js_arrayMax :: ArrayT t ds -> ArrayT t ds -> ArrayT t ds+foreign import javascript unsafe "$1.map(function (e, i) { return Math.min(e,$2[i]); })" js_arrayMin :: ArrayT t ds -> ArrayT t ds -> ArrayT t ds instance Dimensions ds => Num (ArrayT Float ds) where
src-ghcjs/Numeric/Array/Family/ArrayT.js view
@@ -170,6 +170,7 @@ function h$easytensor_transpose(n, mat) {+ 'use strict'; var nmat = new mat.constructor(mat.length), m = Math.round(mat.length / n); for(var i = 0; i < n; i++) {@@ -181,99 +182,119 @@ } function h$easytensor_eyeFloat32(n) {- var mat = new Float32Array(n*n).fill(0);+ 'use strict';+ var mat = new Float32Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeFloat64(n) {- var mat = new Float64Array(n*n).fill(0);+ 'use strict';+ var mat = new Float64Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeInt8(n) {- var mat = new Int8Array(n*n).fill(0);+ 'use strict';+ var mat = new Int8Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeInt16(n) {- var mat = new Int16Array(n*n).fill(0);+ 'use strict';+ var mat = new Int16Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeInt32(n) {- var mat = new Int32Array(n*n).fill(0);+ 'use strict';+ var mat = new Int32Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeUint8(n) {- var mat = new Uint8Array(n*n).fill(0);+ 'use strict';+ var mat = new Uint8Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeUint8Clamped(n) {- var mat = new Uint8ClampedArray(n*n).fill(0);+ 'use strict';+ var mat = new Uint8ClampedArray(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeUint16(n) {- var mat = new Uint16Array(n*n).fill(0);+ 'use strict';+ var mat = new Uint16Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_eyeUint32(n) {- var mat = new Uint32Array(n*n).fill(0);+ 'use strict';+ var mat = new Uint32Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=1;} return mat; } function h$easytensor_diagFloat32(n,x) {- var mat = new Float32Array(n*n).fill(0);+ 'use strict';+ 'use strict';+ var mat = new Float32Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagFloat64(n,x) {- var mat = new Float64Array(n*n).fill(0);+ 'use strict';+ var mat = new Float64Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagInt8(n,x) {- var mat = new Int8Array(n*n).fill(0);+ 'use strict';+ var mat = new Int8Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagInt16(n,x) {- var mat = new Int16Array(n*n).fill(0);+ 'use strict';+ var mat = new Int16Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagInt32(n,x) {- var mat = new Int32Array(n*n).fill(0);+ 'use strict';+ var mat = new Int32Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagUint8(n,x) {- var mat = new Uint8Array(n*n).fill(0);+ 'use strict';+ var mat = new Uint8Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagUint8Clamped(n,x) {- var mat = new Uint8ClampedArray(n*n).fill(0);+ 'use strict';+ var mat = new Uint8ClampedArray(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagUint16(n,x) {- var mat = new Uint16Array(n*n).fill(0);+ 'use strict';+ var mat = new Uint16Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_diagUint32(n,x) {- var mat = new Uint32Array(n*n).fill(0);+ 'use strict';+ var mat = new Uint32Array(n*n); for(var i = 0; i < n*n; i += n + 1){mat[i]=x;} return mat; } function h$easytensor_trace(mat, n) {+ 'use strict'; var r = 0; for(var i = 0; i < n*n; i += n + 1){r+=mat[i];} return r;@@ -281,6 +302,7 @@ function h$easytensor_det(mat, n) {+ 'use strict'; switch (n) { case 1: return mat[0];@@ -296,10 +318,12 @@ } function h$easytensor_detJSMat2(mat) {+ 'use strict'; return (mat[0]*mat[3] - mat[1]*mat[2]); } function h$easytensor_detJSMat3(mat) {+ 'use strict'; return ( mat[0]*(mat[4]*mat[8]-mat[5]*mat[7]) - mat[1]*(mat[3]*mat[8]-mat[5]*mat[6])@@ -308,6 +332,7 @@ } function h$easytensor_detJSMat4(mat) {+ 'use strict'; var n11 = mat[ 0 ], n12 = mat[ 4 ], n13 = mat[ 8 ], n14 = mat[ 12 ]; var n21 = mat[ 1 ], n22 = mat[ 5 ], n23 = mat[ 9 ], n24 = mat[ 13 ]; var n31 = mat[ 2 ], n32 = mat[ 6 ], n33 = mat[ 10 ], n34 = mat[ 14 ];@@ -353,6 +378,7 @@ function h$easytensor_inverse(mat, n) {+ 'use strict'; switch (n) { case 1: return 1 / mat[0];@@ -369,6 +395,7 @@ function h$easytensor_inverseJSM4(mat) {+ 'use strict'; var rez = new mat.constructor(16); rez[0] = mat[13]*(mat[ 6]*mat[11]-mat[10]*mat[ 7])+mat[ 9]*(mat[14]*mat[ 7]-mat[ 6]*mat[15])+mat[ 5]*(mat[10]*mat[15]-mat[14]*mat[11]); rez[4] = mat[12]*(mat[10]*mat[ 7]-mat[ 6]*mat[11])+mat[ 8]*(mat[ 6]*mat[15]-mat[14]*mat[ 7])+mat[ 4]*(mat[14]*mat[11]-mat[10]*mat[15]); @@ -390,12 +417,13 @@ if (det === 0) { return undefined; } else {- for(var i = 0; i < 16; i++) {rez[i] !== det;}+ for(var i = 0; i < 16; i++) {rez[i] /= det;} return rez; } } function h$easytensor_inverseJSM3(mat) {+ 'use strict'; var rez = new mat.constructor(9); rez[0] = mat[4]*mat[8] - mat[7]*mat[5]; rez[3] = mat[6]*mat[5] - mat[3]*mat[8]; @@ -410,12 +438,13 @@ if (det === 0) { return undefined; } else {- for(var i = 0; i < 9; i++) {rez[i] !== det;}+ for(var i = 0; i < 9; i++) {rez[i] /= det;} return rez; } } function h$easytensor_inverseJSM2(mat) {+ 'use strict'; var det = mat[0]*mat[3] - mat[1]*mat[2]; if (det === 0) { return undefined;@@ -430,6 +459,7 @@ function h$easytensor_contract(n,m,k,lhs,rhs) {+ 'use strict'; var t, rez = new lhs.constructor(n*k); for(var i = 0; i < n; i++) { for(var j = 0; j < k; j++) {@@ -441,4 +471,18 @@ } } return rez;+}+++function h$easytensor_dot(lhs, rhs) {+ 'use strict';+ return lhs.reduce(function (r, e, i) { return r + e*rhs[i];}, 0);+}++function h$easytensor_cross(a, b) {+ 'use strict';+ return [ a[1]*b[2]-a[2]*b[1]+ , a[2]*b[0]-a[0]*b[2]+ , a[0]*b[1]-a[1]*b[0]+ ]; }
src-ghcjs/Numeric/DataFrame/Contraction.hs view
@@ -11,6 +11,9 @@ {-# LANGUAGE UnboxedTuples #-} {-# LANGUAGE UndecidableInstances #-} {-# LANGUAGE InstanceSigs #-}+{-# LANGUAGE UnliftedFFITypes #-}+{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE GHCForeignImportPrim #-} ----------------------------------------------------------------------------- -- | -- Module : Numeric.DataFrame.Contraction@@ -62,17 +65,29 @@ {-# INLINE (%*) #-} infixl 7 %* -instance ( ConcatList as bs asbs+instance {-# OVERLAPPABLE #-}+ ( ConcatList as bs asbs , Dimensions as , Dimensions bs+ , asbs ~ (a' ': sbs') ) => Contraction t as bs asbs where contract :: forall m . KnownDim m => DataFrame t (as +: m) -> DataFrame t (m :+ bs) -> DataFrame t asbs contract dx dy- | Refl <- unsafeCoerce Refl :: Array t asbs :~: ArrayT t asbs- , Refl <- unsafeCoerce Refl :: Array t (as +: m) :~: ArrayT t (as +: m)+ | Refl <- unsafeCoerce Refl :: Array t (as +: m) :~: ArrayT t (as +: m) , Evidence <- inferConcatDimensions @as @bs- = KnownDataFrame $ js_conctract @t @as @m @bs (dimVal (dim @as)) (dimVal' @m) (dimVal (dim @bs)) (coerce dx) (coerce dy)+ = KnownDataFrame $ js_contract @t @as @m @bs (dimVal (dim @as)) (dimVal' @m) (dimVal (dim @bs)) (coerce dx) (coerce dy) foreign import javascript unsafe "h$easytensor_contract($1,$2,$3,$4,$5)"- js_conctract :: forall t as m bs . Int -> Int -> Int -> ArrayT t (as +: m) -> ArrayT t (m :+ bs) -> ArrayT t (as ++ bs)+ js_contract :: forall t as m bs . Int -> Int -> Int -> ArrayT t (as +: m) -> ArrayT t (m :+ bs) -> ArrayT t (as ++ bs)++instance {-# OVERLAPPING #-}+ Contraction t '[] '[] '[] where+ contract :: forall m . KnownDim m => DataFrame t '[m] -> DataFrame t '[m] -> DataFrame t ('[] :: [Nat])+ contract dx dy+ = KnownDataFrame $ unsafeCoerce (js_contract0 (coerce dx) (coerce dy))++foreign import javascript unsafe "$1.reduce(function (r, e, i) { return e*$2[i] + r;}, 0)"+ js_contract0 :: ArrayT t '[m] -> ArrayT t '[m] -> Any++
+ src-ghcjs/Numeric/DataFrame/Mutable.hs view
@@ -0,0 +1,296 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE ExistentialQuantification #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE KindSignatures #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UnboxedTuples #-}+{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE GHCForeignImportPrim #-}+{-# LANGUAGE UnliftedFFITypes #-}+{-# LANGUAGE TypeOperators #-}+-----------------------------------------------------------------------------+-- |+-- Module : Numeric.DataFrame.Mutable+-- Copyright : (c) Artem Chirkin+-- License : BSD3+--+-- Maintainer : chirkin@arch.ethz.ch+--+-- Interfrace to perform primitive stateful operations on mutable frames.+--+-----------------------------------------------------------------------------++module Numeric.DataFrame.Mutable+ ( MutableFrame (..), MDataFrame (..)+ , newDataFrame#, copyDataFrame#, copyMDataFrame#, unsafeFreezeDataFrame#+ , freezeDataFrame#, thawDataFrame#+ , writeDataFrame#, readDataFrame#+ , newArrayBuffer#, arrayBuffer#, viewFloatArray#, viewDoubleArray#+ , viewIntArray#, viewInt32Array#, viewInt16Array#, viewInt8Array#+ , viewWordArray#, viewWord32Array#, viewWord16Array#, viewWord8Array#, viewWord8ClampedArray#+ ) where++++import GHCJS.Types (IsJSVal(), JSVal)+import GHC.Int (Int16 (..), Int32 (..),Int8 (..))+import GHC.Prim+import GHC.Types (Double (..), Float (..), Int (..), Word (..))+import GHC.Word (Word16 (..), Word32 (..), Word8 (..))+import Unsafe.Coerce (unsafeCoerce)++import Numeric.DataFrame.Type+import Numeric.Array.Family+import Numeric.Dimensions+++-- | Mutable DataFrame type+newtype MDataFrame s t (ns :: [Nat]) = MDataFrame (MutableArrayT s t ns)+instance IsJSVal (MDataFrame s t ds)+++-- | Create a new mutable DataFrame.+newDataFrame# :: forall t (ns :: [Nat]) s+ . (ElemTypeInference t, Dimensions ns)+ => State# s -> (# State# s, MDataFrame s t ns #)+newDataFrame# = case elemTypeInstance @t of+ ETFloat -> js_createFloatArray n+ ETDouble -> js_createDoubleArray n+ ETInt -> js_createIntArray n+ ETInt8 -> js_createInt8Array n+ ETInt16 -> js_createInt16Array n+ ETInt32 -> js_createInt32Array n+ ETWord -> js_createWordArray n+ ETWord8 -> js_createWord8Array n+ ETWord16 -> js_createWord16Array n+ ETWord32 -> js_createWord32Array n+ ETWord8C -> js_createWord8ClampedArray n+ where+ n = dimVal (dim @ns)+{-# INLINE newDataFrame# #-}+++-- | Copy one DataFrame into another mutable DataFrame at specified position.+copyDataFrame# :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s+ . ( ArraySizeInference (as +: b')+ , ConcatList as (b :+ bs) asbs+ , Dimensions as+ , Dimensions (b :+ bs)+ )+ => DataFrame t (as +: b') -> Idx (b :+ bs) -> MDataFrame s t asbs -> State# s -> (# State# s, () #)+copyDataFrame# df i mdf s0 = case arraySizeInstance @(as +: b') of+ ASScalar -> df `seq` (# js_writeArrayOffsetJSVal# mdf (fromEnum i) (unsafeCoerce df) s0, () #)+ ASArray -> js_copyDataFrame (coerce df) (fromEnum i * dimVal (dim @as)) mdf s0+{-# INLINE copyDataFrame# #-}++++-- | Copy one mutable DataFrame into another mutable DataFrame at specified position.+copyMDataFrame# :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s+ . ( ConcatList as (b :+ bs) asbs+ , Dimensions as+ , Dimensions (b :+ bs)+ )+ => MDataFrame s t (as +: b') -> Idx (b :+ bs) -> MDataFrame s t asbs -> State# s -> (# State# s, () #)+copyMDataFrame# d i = js_copyMDataFrame d (fromEnum i * dimVal (dim @as))+{-# INLINE copyMDataFrame# #-}+++-- | Make a mutable DataFrame immutable, without copying.+unsafeFreezeDataFrame# :: forall t (ns :: [Nat]) s+ . (MutableFrame t ns, ArraySizeInference ns)+ => MDataFrame s t ns -> State# s -> (# State# s, DataFrame t ns #)+unsafeFreezeDataFrame# a s = case arraySizeInstance @ns of+ ASScalar -> case readDataFrameOff# a 0# s of+ (# s1, v #) -> (# s1, coerce v #)+ ASArray -> (# s, coerce a #)+{-# INLINE unsafeFreezeDataFrame# #-}+++--unsafeThawArrayT# :: ArrayT t ds -> State# s -> (#State# s, MutableArrayT s t ds #)+--unsafeThawArrayT# a s = (# s, coerce a #)+--{-# INLINE unsafeThawArrayT# #-}+++-- | Copy content of a mutable DataFrame into a new immutable DataFrame.+freezeDataFrame# :: forall t (ns :: [Nat]) s+ . (MutableFrame t ns, ArraySizeInference ns)+ => MDataFrame s t ns -> State# s -> (# State# s, DataFrame t ns #)+freezeDataFrame# a s = case arraySizeInstance @ns of+ ASScalar -> case readDataFrameOff# a 0# s of+ (# s1, v #) -> (# s1, coerce v #)+ ASArray -> case js_freeze a s of+ (# s1, v #) -> (# s1, coerce v #)+{-# INLINE freezeDataFrame# #-}++++-- | Create a new mutable DataFrame and copy content of immutable one in there.+thawDataFrame# :: forall t (ns :: [Nat]) s+ . (MutableFrame t ns, ArrayInstanceInference t ns)+ => DataFrame t ns -> State# s -> (# State# s, MDataFrame s t ns #)+thawDataFrame# a s = case arraySizeInstance @ns of+ ASScalar -> case newDataFrame# @t @'[] s of+ (# s1, df #) -> (# writeDataFrameOff# df 0# (coerce a) s1, df #)+ ASArray -> js_thaw (coerce a) s+{-# INLINE thawDataFrame# #-}++-- | Write a single element at the specified index+writeDataFrame# :: forall t (ns :: [Nat]) s+ . ( MutableFrame t ns, Dimensions ns )+ => MDataFrame s t ns -> Idx ns -> t -> State# s -> (# State# s, () #)+writeDataFrame# mdf ei x s | I# i <- fromEnum ei = (# writeDataFrameOff# mdf i x s, () #)+{-# INLINE writeDataFrame# #-}++-- | Read a single element at the specified index+readDataFrame# :: forall t (ns :: [Nat]) s+ . ( MutableFrame t ns, Dimensions ns )+ => MDataFrame s t ns -> Idx ns -> State# s -> (# State# s, t #)+readDataFrame# mdf ei | I# i <- fromEnum ei = readDataFrameOff# mdf i+{-# INLINE readDataFrame# #-}++class MutableFrame t (ns :: [Nat]) where+ -- | Write a single element at the specified element offset+ writeDataFrameOff# :: MDataFrame s t ns -> Int# -> t -> State# s -> State# s+ -- | Read a single element at the specified element offset+ readDataFrameOff# :: MDataFrame s t ns -> Int# -> State# s -> (# State# s, t #)++instance MutableFrame Float (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetFloat#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetFloat#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Double (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetDouble#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetDouble#+ {-# INLINE readDataFrameOff# #-}+++instance MutableFrame Int (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetInt#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetInt#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int8 (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetInt8#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetInt8#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int16 (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetInt16#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetInt16#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Int32 (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetInt32#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetInt32#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetWord#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetWord#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word8 (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetWord8#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetWord8#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word16 (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetWord16#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetWord16#+ {-# INLINE readDataFrameOff# #-}++instance MutableFrame Word32 (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetWord32#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# = js_readArrayOffsetWord32#+ {-# INLINE readDataFrameOff# #-}+++instance MutableFrame Word8Clamped (ns :: [Nat]) where+ writeDataFrameOff# = js_writeArrayOffsetWord8Clamped#+ {-# INLINE writeDataFrameOff# #-}+ readDataFrameOff# m o s = case js_readArrayOffsetWord8Clamped# m o s of+ (# s1, y #) -> (# s1, coerce y #)+ {-# INLINE readDataFrameOff# #-}++++foreign import javascript unsafe "new Float32Array($1)" js_createFloatArray :: Int -> State# s -> (# State# s, MDataFrame s Float ds #)+foreign import javascript unsafe "new Float64Array($1)" js_createDoubleArray :: Int -> State# s -> (# State# s, MDataFrame s Double ds #)+foreign import javascript unsafe "new Int32Array($1)" js_createIntArray :: Int -> State# s -> (# State# s, MDataFrame s Int ds #)+foreign import javascript unsafe "new Int32Array($1)" js_createInt32Array :: Int -> State# s -> (# State# s, MDataFrame s Int32 ds #)+foreign import javascript unsafe "new Int16Array($1)" js_createInt16Array :: Int -> State# s -> (# State# s, MDataFrame s Int16 ds #)+foreign import javascript unsafe "new Int8Array($1)" js_createInt8Array :: Int -> State# s -> (# State# s, MDataFrame s Int8 ds #)+foreign import javascript unsafe "new Uint32Array($1)" js_createWordArray :: Int -> State# s -> (# State# s, MDataFrame s Word ds #)+foreign import javascript unsafe "new Uint32Array($1)" js_createWord32Array :: Int -> State# s -> (# State# s, MDataFrame s Word32 ds #)+foreign import javascript unsafe "new Uint16Array($1)" js_createWord16Array :: Int -> State# s -> (# State# s, MDataFrame s Word16 ds #)+foreign import javascript unsafe "new Uint8Array($1)" js_createWord8Array :: Int -> State# s -> (# State# s, MDataFrame s Word8 ds #)+foreign import javascript unsafe "new Uint8ClampedArray($1)" js_createWord8ClampedArray :: Int -> State# s -> (# State# s, MDataFrame s Word8Clamped ds #)+++foreign import javascript unsafe "$1[$2]" js_readArrayOffsetFloat# :: MDataFrame s Float ds -> Int# -> State# s -> (# State# s, Float #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetDouble# :: MDataFrame s Double ds -> Int# -> State# s -> (# State# s, Double #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetInt# :: MDataFrame s Int ds -> Int# -> State# s -> (# State# s, Int #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetInt8# :: MDataFrame s Int8 ds -> Int# -> State# s -> (# State# s, Int8 #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetInt16# :: MDataFrame s Int16 ds -> Int# -> State# s -> (# State# s, Int16 #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetInt32# :: MDataFrame s Int32 ds -> Int# -> State# s -> (# State# s, Int32 #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetWord# :: MDataFrame s Word ds -> Int# -> State# s -> (# State# s, Word #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetWord8# :: MDataFrame s Word8 ds -> Int# -> State# s -> (# State# s, Word8 #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetWord8Clamped# :: MDataFrame s Word8Clamped ds -> Int# -> State# s -> (# State# s, Int #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetWord16# :: MDataFrame s Word16 ds -> Int# -> State# s -> (# State# s, Word16 #)+foreign import javascript unsafe "$1[$2]" js_readArrayOffsetWord32# :: MDataFrame s Word32 ds -> Int# -> State# s -> (# State# s, Word32 #)++++foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetFloat# :: MDataFrame s Float ds -> Int# -> Float -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetDouble# :: MDataFrame s Double ds -> Int# -> Double -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetInt# :: MDataFrame s Int ds -> Int# -> Int -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetInt8# :: MDataFrame s Int8 ds -> Int# -> Int8 -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetInt16# :: MDataFrame s Int16 ds -> Int# -> Int16 -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetInt32# :: MDataFrame s Int32 ds -> Int# -> Int32 -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetWord# :: MDataFrame s Word ds -> Int# -> Word -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetWord8# :: MDataFrame s Word8 ds -> Int# -> Word8 -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetWord8Clamped# :: MDataFrame s Word8Clamped ds -> Int# -> Word8Clamped -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetWord16# :: MDataFrame s Word16 ds -> Int# -> Word16 -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetWord32# :: MDataFrame s Word32 ds -> Int# -> Word32 -> State# s -> State# s+foreign import javascript unsafe "$1[$2] = $3;" js_writeArrayOffsetJSVal# :: MDataFrame s t ds -> Int -> JSVal -> State# s -> State# s+++foreign import javascript unsafe "$3.set($1, $2);" js_copyDataFrame :: ArrayT t as -> Int -> MDataFrame s t asbs -> State# s -> (# State# s, () #)+foreign import javascript unsafe "$3.set($1, $2);" js_copyMDataFrame :: MDataFrame s t as -> Int -> MDataFrame s t asbs -> State# s -> (# State# s, () #)+++foreign import javascript unsafe "$1.slice()" js_freeze :: MDataFrame s t as -> State# s -> (# State# s, ArrayT t ds #)+foreign import javascript unsafe "$1.slice()" js_thaw :: ArrayT t as -> State# s -> (# State# s, MDataFrame s t ds #)++++foreign import javascript unsafe "new ArrayBuffer($1)" newArrayBuffer# :: Int -> State# s -> (# State# s, JSVal #)+foreign import javascript unsafe "new Float32Array($1)" viewFloatArray# :: JSVal -> State# s -> (# State# s, MDataFrame s Float ds #)+foreign import javascript unsafe "new Float64Array($1)" viewDoubleArray# :: JSVal -> State# s -> (# State# s, MDataFrame s Double ds #)+foreign import javascript unsafe "new Int32Array($1)" viewIntArray# :: JSVal -> State# s -> (# State# s, MDataFrame s Int ds #)+foreign import javascript unsafe "new Int32Array($1)" viewInt32Array# :: JSVal -> State# s -> (# State# s, MDataFrame s Int32 ds #)+foreign import javascript unsafe "new Int16Array($1)" viewInt16Array# :: JSVal -> State# s -> (# State# s, MDataFrame s Int16 ds #)+foreign import javascript unsafe "new Int8Array($1)" viewInt8Array# :: JSVal -> State# s -> (# State# s, MDataFrame s Int8 ds #)+foreign import javascript unsafe "new Uint32Array($1)" viewWordArray# :: JSVal -> State# s -> (# State# s, MDataFrame s Word ds #)+foreign import javascript unsafe "new Uint32Array($1)" viewWord32Array# :: JSVal -> State# s -> (# State# s, MDataFrame s Word32 ds #)+foreign import javascript unsafe "new Uint16Array($1)" viewWord16Array# :: JSVal -> State# s -> (# State# s, MDataFrame s Word16 ds #)+foreign import javascript unsafe "new Uint8Array($1)" viewWord8Array# :: JSVal -> State# s -> (# State# s, MDataFrame s Word8 ds #)+foreign import javascript unsafe "new Uint8ClampedArray($1)" viewWord8ClampedArray# :: JSVal -> State# s -> (# State# s, MDataFrame s Word8Clamped ds #)+foreign import javascript unsafe "$1.buffer" arrayBuffer# :: MDataFrame s t ds -> State# s -> (# State# s, JSVal #)
+ src-ghcjs/Numeric/Matrix/Mat44.js view
@@ -0,0 +1,103 @@+function h$easytensor_m4fromHom(v) {+ 'use strict';+ var r = v.slice(0,3), t = v[3];+ if (t !== 0) {+ r[0] /= t; r[1] /= t; r[2] /= t;+ }+ return r; +}+++function h$easytensor_m4translate(v) {+ 'use strict';+ var m = new v.constructor(16);+ m.set(v, 12);+ m[0] = 1;+ m[5] = 1;+ m[10] = 1;+ m[15] = 1;+ return m;+} + +function h$easytensor_m4rotateX(a) {+ 'use strict';+ var c = Math.cos(a), s = Math.sin(a);+ return [ 1, 0, 0, 0+ , 0, c, s, 0+ , 0,-s, c, 0+ , 0, 0, 0, 1]; +}+ +function h$easytensor_m4rotateY(a) {+ 'use strict';+ var c = Math.cos(a), s = Math.sin(a);+ return [ c, 0,-s, 0+ , 0, 1, 0, 0+ , s, 0, c, 0+ , 0, 0, 0, 1]; +}++function h$easytensor_m4rotateZ(a) {+ 'use strict';+ var c = Math.cos(a), s = Math.sin(a);+ return [ c, s, 0, 0+ ,-s, c, 0, 0+ , 0, 0, 1, 0+ , 0, 0, 0, 1]; +}++function h$easytensor_m4rotate(vec, a) {+ 'use strict';+ var c = Math.cos(a);+ var s = Math.sin(a);+ var c1 = 1 - c;+ var x = vec[0], y = vec[1], z = vec[2];+ return [ c + c1*x*x, c1*x*y + s*z, c1*x*z - s*y, 0+ , c1*x*y - s*z, c + c1*y*y, c1*y*z + s*x, 0+ , c1*x*z + s*y, c1*y*z - s*x, c + c1*z*z, 0+ , 0, 0, 0, 1]; +}++function h$easytensor_m4rotateEuler(x, y, z) {+ 'use strict';+ var cx = Math.cos(x), sx = Math.sin(x), cy = Math.cos(y), sy = Math.sin(y), cz = Math.cos(z), sz = Math.sin(z);+ return [ cy*cz, -cy*sz, sy, 0+ , cx*sz + sx*sy*cz, cx*cz - sx*sy*sz, -sx*cy, 0+ , sx*sz - cx*sy*cz, sx*cz + cx*sy*sz, cx*cy, 0+ , 0, 0, 0, 1]; +}++function h$easytensor_m4lookAt(up,camera,point) {+ 'use strict';+ var zDir = [ camera[0] - point[0], camera[1] - point[1], camera[2] - point[2] ];+ var t = Math.hypot.apply(null,zDir);+ zDir = zDir.map(function (e){return e / t;});+ var xDir = h$easytensor_cross(up,zDir);+ t = Math.hypot.apply(null,xDir);+ xDir = xDir.map(function (e){return e / t;});+ var yDir = h$easytensor_cross(zDir,xDir);+ return [ xDir[0], yDir[0], zDir[0], 0+ , xDir[1], yDir[1], zDir[1], 0+ , xDir[2], yDir[2], zDir[2], 0+ , - h$easytensor_dot(xDir,camera), - h$easytensor_dot(yDir,camera), - h$easytensor_dot(zDir,camera), 1+ ]; +}++function h$easytensor_m4perspective(n, f, fovy, aspect) {+ 'use strict';+ var h2 = n*Math.tan(fovy/2); + var w2 = aspect*h2;+ return [ n/w2, 0, 0, 0+ , 0, n/h2, 0, 0+ , 0, 0, (n+f)/(n-f),-1+ , 0, 0, 2*n*f/(n-f), 0 ]; +}++function h$easytensor_m4orthogonal(n, f, w, h) {+ 'use strict';+ return [ 2/w, 0, 0, 0+ , 0, 2/h, 0, 0+ , 0, 0, 2/(n-f), 0+ , 0, 0, (n+f)/(n-f), 1 ]; +}+
+ src-ghcjs/Numeric/Matrix/Mat44d.hs view
@@ -0,0 +1,69 @@+{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE DataKinds #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Matrix.Mat44d () where+++import Numeric.DataFrame.Type+import Numeric.Vector+import Numeric.Matrix.Class+++instance HomTransform4 Double where+ translate4 = js_translate+ {-# INLINE translate4 #-}+ translate3 = js_translate+ {-# INLINE translate3 #-}+ rotateX = js_rotateX+ {-# INLINE rotateX #-}+ rotateY = js_rotateY+ {-# INLINE rotateY #-}+ rotateZ = js_rotateZ+ {-# INLINE rotateZ #-}+ rotate = js_rotate+ {-# INLINE rotate #-}+ rotateEuler = js_rotateEuler+ {-# INLINE rotateEuler #-}+ lookAt = js_lookAt+ {-# INLINE lookAt #-}+ perspective = js_perspective+ {-# INLINE perspective #-}+ orthogonal = js_orthogonal+ {-# INLINE orthogonal #-}+ toHomPoint = js_toHomPoint+ {-# INLINE toHomPoint #-}+ toHomVector = js_toHomVector+ {-# INLINE toHomVector #-}+ fromHom = js_fromHom+ {-# INLINE fromHom #-}++++foreign import javascript unsafe "h$easytensor_m4translate($1)" + js_translate :: Vector Double n -> Matrix Double 4 4 +foreign import javascript unsafe "new Float64Array(h$easytensor_m4rotateX($1))" + js_rotateX :: Double -> Matrix Double 4 4 +foreign import javascript unsafe "new Float64Array(h$easytensor_m4rotateY($1))" + js_rotateY :: Double -> Matrix Double 4 4+foreign import javascript unsafe "new Float64Array(h$easytensor_m4rotateZ($1))" + js_rotateZ :: Double -> Matrix Double 4 4+foreign import javascript unsafe "new Float64Array(h$easytensor_m4rotate($1, $2))" + js_rotate :: Vector Double 3 -> Double -> Matrix Double 4 4+foreign import javascript unsafe "new Float64Array(h$easytensor_m4rotateEuler($1, $2, $3))" + js_rotateEuler :: Double -> Double -> Double -> Matrix Double 4 4+foreign import javascript unsafe "new Float64Array(h$easytensor_m4lookAt($1,$2,$3))" + js_lookAt :: Vector Double 3 -> Vector Double 3 -> Vector Double 3 -> Matrix Double 4 4+foreign import javascript unsafe "new Float64Array(h$easytensor_m4perspective($1, $2, $3, $4))" + js_perspective :: Double -> Double -> Double -> Double -> Matrix Double 4 4+foreign import javascript unsafe "new Float64Array(h$easytensor_m4orthogonal($1, $2, $3, $4))" + js_orthogonal :: Double -> Double -> Double -> Double -> Matrix Double 4 4++foreign import javascript unsafe "$r = new $1.constructor(4); $r.set($1); $r[3] = 1;"+ js_toHomPoint :: Vector Double 3 -> Vector Double 4+foreign import javascript unsafe "$r = new $1.constructor(4); $r.set($1); $r[3] = 0;"+ js_toHomVector :: Vector Double 3 -> Vector Double 4+foreign import javascript unsafe "h$easytensor_m4fromHom($1)"+ js_fromHom :: Vector Double 4 -> Vector Double 3+++
+ src-ghcjs/Numeric/Matrix/Mat44f.hs view
@@ -0,0 +1,69 @@+{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE DataKinds #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Matrix.Mat44f () where+++import Numeric.DataFrame.Type+import Numeric.Vector+import Numeric.Matrix.Class+++instance HomTransform4 Float where+ translate4 = js_translate+ {-# INLINE translate4 #-}+ translate3 = js_translate+ {-# INLINE translate3 #-}+ rotateX = js_rotateX+ {-# INLINE rotateX #-}+ rotateY = js_rotateY+ {-# INLINE rotateY #-}+ rotateZ = js_rotateZ+ {-# INLINE rotateZ #-}+ rotate = js_rotate+ {-# INLINE rotate #-}+ rotateEuler = js_rotateEuler+ {-# INLINE rotateEuler #-}+ lookAt = js_lookAt+ {-# INLINE lookAt #-}+ perspective = js_perspective+ {-# INLINE perspective #-}+ orthogonal = js_orthogonal+ {-# INLINE orthogonal #-}+ toHomPoint = js_toHomPoint+ {-# INLINE toHomPoint #-}+ toHomVector = js_toHomVector+ {-# INLINE toHomVector #-}+ fromHom = js_fromHom+ {-# INLINE fromHom #-}++++foreign import javascript unsafe "h$easytensor_m4translate($1)" + js_translate :: Vector Float n -> Matrix Float 4 4 +foreign import javascript unsafe "new Float32Array(h$easytensor_m4rotateX($1))" + js_rotateX :: Float -> Matrix Float 4 4 +foreign import javascript unsafe "new Float32Array(h$easytensor_m4rotateY($1))" + js_rotateY :: Float -> Matrix Float 4 4+foreign import javascript unsafe "new Float32Array(h$easytensor_m4rotateZ($1))" + js_rotateZ :: Float -> Matrix Float 4 4+foreign import javascript unsafe "new Float32Array(h$easytensor_m4rotate($1, $2))" + js_rotate :: Vector Float 3 -> Float -> Matrix Float 4 4+foreign import javascript unsafe "new Float32Array(h$easytensor_m4rotateEuler($1, $2, $3))" + js_rotateEuler :: Float -> Float -> Float -> Matrix Float 4 4+foreign import javascript unsafe "new Float32Array(h$easytensor_m4lookAt($1,$2,$3))" + js_lookAt :: Vector Float 3 -> Vector Float 3 -> Vector Float 3 -> Matrix Float 4 4+foreign import javascript unsafe "new Float32Array(h$easytensor_m4perspective($1, $2, $3, $4))" + js_perspective :: Float -> Float -> Float -> Float -> Matrix Float 4 4+foreign import javascript unsafe "new Float32Array(h$easytensor_m4orthogonal($1, $2, $3, $4))" + js_orthogonal :: Float -> Float -> Float -> Float -> Matrix Float 4 4++foreign import javascript unsafe "$r = new $1.constructor(4); $r.set($1); $r[3] = 1;"+ js_toHomPoint :: Vector Float 3 -> Vector Float 4+foreign import javascript unsafe "$r = new $1.constructor(4); $r.set($1); $r[3] = 0;"+ js_toHomVector :: Vector Float 3 -> Vector Float 4+foreign import javascript unsafe "h$easytensor_m4fromHom($1)"+ js_fromHom :: Vector Float 4 -> Vector Float 3+++
+ src-ghcjs/Numeric/Quaternion/QDouble.hs view
@@ -0,0 +1,285 @@+{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeSynonymInstances #-}+{-# LANGUAGE FlexibleInstances #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Quaternion.QDouble+ ( QDouble, Quater (..)+ ) where++import Data.JSString (JSString, unpack')+import Data.Coerce (coerce)++import Numeric.Array+import Numeric.DataFrame.Type+import Numeric.Vector+import Numeric.Matrix++import Numeric.Quaternion.Class+++type QDouble = Quater Double++instance Quaternion Double where+ newtype Quater Double = QDouble (ArrayT Double '[4])+ {-# INLINE packQ #-}+ packQ = js_packQ+ {-# INLINE fromVecNum #-}+ fromVecNum = js_fromVecNum+ {-# INLINE fromVec4 #-}+ fromVec4 = coerce+ {-# INLINE toVec4 #-}+ toVec4 = coerce+ {-# INLINE unpackQ #-}+ unpackQ = js_unpackQ+ {-# INLINE square #-}+ square = js_square+ {-# INLINE im #-}+ im = js_im+ {-# INLINE re #-}+ re = js_re+ {-# INLINE imVec #-}+ imVec = js_imVec+ {-# INLINE taker #-}+ taker = js_taker+ {-# INLINE takei #-}+ takei = js_takei+ {-# INLINE takej #-}+ takej = js_takej+ {-# INLINE takek #-}+ takek = js_takek+ {-# INLINE conjugate #-}+ conjugate = js_conjugate+ {-# INLINE rotScale #-}+ rotScale = js_rotScale+ {-# INLINE getRotScale #-}+ getRotScale = js_getRotScale+ {-# INLINE axisRotation #-}+ axisRotation = js_axisRotation+ {-# INLINE qArg #-}+ qArg = js_qArg+ {-# INLINE fromMatrix33 #-}+ fromMatrix33 = js_fromMatrix33+ {-# INLINE fromMatrix44 #-}+ fromMatrix44 = js_fromMatrix44+ {-# INLINE toMatrix33 #-}+ toMatrix33 = js_toMatrix33+ {-# INLINE toMatrix44 #-}+ toMatrix44 = js_toMatrix44+++foreign import javascript unsafe "new Float64Array([$1,$2,$3,$4])"+ js_packQ :: Double -> Double -> Double -> Double -> QDouble++foreign import javascript unsafe "new Float64Array([$1[0],$1[1],$1[2],$2])"+ js_fromVecNum :: Vector Double 3 -> Double -> QDouble++foreign import javascript unsafe "$r1 = $1[0];$r2 = $1[1];$r3 = $1[2];$r4 = $1[3];"+ js_unpackQ :: QDouble -> (Double,Double,Double,Double)+++foreign import javascript unsafe "$1[0]*$1[0] + $1[1]*$1[1] + $1[2]*$1[2] + $1[3]*$1[3]"+ js_square :: QDouble -> Double++foreign import javascript unsafe "new Float64Array([$1[0],$1[1],$1[2],0])"+ js_im :: QDouble -> QDouble++foreign import javascript unsafe "$1.slice(0,3)"+ js_imVec :: QDouble -> Vector Double 3++foreign import javascript unsafe "$r = new Float64Array(4); $r[3] = $1[3];"+ js_re :: QDouble -> QDouble++foreign import javascript unsafe "$1[3]"+ js_taker :: QDouble -> Double++foreign import javascript unsafe "$1[0]"+ js_takei :: QDouble -> Double++foreign import javascript unsafe "$1[1]"+ js_takej :: QDouble -> Double++foreign import javascript unsafe "$1[2]"+ js_takek :: QDouble -> Double++foreign import javascript unsafe "new Float64Array([-$1[0],-$1[1],-$1[2],$1[3]])"+ js_conjugate :: QDouble -> QDouble+++foreign import javascript unsafe "new Float64Array(h$easytensor_rotScale($1,$2))"+ js_rotScale :: QDouble -> Vector Double 3 -> Vector Double 3++foreign import javascript unsafe "new Float64Array(h$easytensor_getRotScale($1,$2))"+ js_getRotScale :: Vector Double 3 -> Vector Double 3 -> QDouble++foreign import javascript unsafe "new Float64Array(h$easytensor_axisRotation($1,$2))"+ js_axisRotation :: Vector Double 3 -> Double -> QDouble++foreign import javascript unsafe "h$easytensor_qArg($1)"+ js_qArg :: QDouble -> Double++foreign import javascript unsafe "new Float64Array(h$easytensor_qfromMatrix33($1))"+ js_fromMatrix33 :: Matrix Double 3 3 -> QDouble++foreign import javascript unsafe "new Float64Array(h$easytensor_qfromMatrix44($1))"+ js_fromMatrix44 :: Matrix Double 4 4 -> QDouble++foreign import javascript unsafe "new Float64Array(h$easytensor_qtoMatrix33($1))"+ js_toMatrix33 :: QDouble -> Matrix Double 3 3++foreign import javascript unsafe "new Float64Array(h$easytensor_qtoMatrix44($1))"+ js_toMatrix44 :: QDouble -> Matrix Double 4 4+++--------------------------------------------------------------------------+-- Num+--------------------------------------------------------------------------++instance Num QDouble where+ {-# INLINE (+) #-}+ (+) = js_plus+ {-# INLINE (-) #-}+ (-) = js_minus+ {-# INLINE (*) #-}+ (*) = js_times+ {-# INLINE abs #-}+ abs = js_abs+ {-# INLINE signum #-}+ signum = js_signum+ {-# INLINE negate #-}+ negate = fromVec4 . negate . toVec4+ {-# INLINE fromInteger #-}+ fromInteger = js_toQuaternion . fromInteger++foreign import javascript unsafe "$1.map(function (e, i) {return e + $2[i];})"+ js_plus :: QDouble -> QDouble -> QDouble+++foreign import javascript unsafe "$1.map(function (e, i) {return e - $2[i];})"+ js_minus :: QDouble -> QDouble -> QDouble++foreign import javascript unsafe "new Float64Array(\+ \[ $1[3]*$2[0] + $1[0]*$2[3] + $1[1]*$2[2] - $1[2]*$2[1]\+ \, $1[3]*$2[1] - $1[0]*$2[2] + $1[1]*$2[3] + $1[2]*$2[0]\+ \, $1[3]*$2[2] + $1[0]*$2[1] - $1[1]*$2[0] + $1[2]*$2[3]\+ \, $1[3]*$2[3] - $1[0]*$2[0] - $1[1]*$2[1] - $1[2]*$2[2] ])"+ js_times :: QDouble -> QDouble -> QDouble++foreign import javascript unsafe "new Float64Array([0,0,0,Math.hypot($1[0],$1[1],$1[2],$1[3])])"+ js_abs :: QDouble -> QDouble++foreign import javascript unsafe "var l = Math.hypot($1[0],$1[1],$1[2],$1[3]); $r = (l == 0) ? (new Float64Array(4)) : $1.map(function (e) {return e/l;})"+ js_signum :: QDouble -> QDouble++foreign import javascript unsafe "$r = new Float64Array(4); $r[3] = $1;"+ js_toQuaternion :: Double -> QDouble++{-# RULES+"realToFrac/DoubleQDouble" realToFrac = js_toQuaternion+"realToFrac/aQDouble" realToFrac = js_toQuaternion . realToFrac+"fromIntegral/aQDouble" fromIntegral = js_toQuaternion . fromIntegral+ #-}+++--------------------------------------------------------------------------+-- Fractional+--------------------------------------------------------------------------++instance Fractional QDouble where+ {-# INLINE recip #-}+ recip = js_recip+ {-# INLINE (/) #-}+ (/) p = js_times p . js_recip+ {-# INLINE fromRational #-}+ fromRational = js_toQuaternion . fromRational++foreign import javascript unsafe "new Float64Array(h$easytensor_qrecip($1))"+ js_recip :: QDouble -> QDouble+++--------------------------------------------------------------------------+-- Floating+--------------------------------------------------------------------------++instance Floating QDouble where+ {-# INLINE pi #-}+ pi = js_toQuaternion pi+ {-# INLINE exp #-}+ exp = js_exp+ {-# INLINE log #-}+ log = js_log+ {-# INLINE sqrt #-}+ sqrt = js_sqrt+ {-# INLINE sin #-}+ sin = js_sin+ {-# INLINE cos #-}+ cos = js_cos+ {-# INLINE tan #-}+ tan = js_tan+ {-# INLINE sinh #-}+ sinh = js_sinh+ {-# INLINE cosh #-}+ cosh = js_cosh+ {-# INLINE tanh #-}+ tanh = js_tanh+ {-# INLINE asin #-}+ asin q = -i * log (i*q + sqrt (1 - q*q))+ where i = signum . im $ q+ {-# INLINE acos #-}+ acos q = pi/2 - asin q+ {-# INLINE atan #-}+ atan q = if square imq == 0+ then js_toQuaternion (atan $ taker q)+ else i / 2 * log ( (i + q) / (i - q) )+ where i = signum imq+ imq = im q+ {-# INLINE asinh #-}+ asinh q = log (q + sqrt (q*q + 1))+ {-# INLINE acosh #-}+ acosh q = log (q + sqrt (q*q - 1))+ {-# INLINE atanh #-}+ atanh q = 0.5 * log ((1+q)/(1-q))++foreign import javascript unsafe "new Float64Array(h$easytensor_qexp($1))" js_exp :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qlog($1))" js_log :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qsqrt($1))" js_sqrt :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qsin($1))" js_sin :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qcos($1))" js_cos :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qtan($1))" js_tan :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qsinh($1))" js_sinh :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qcosh($1))" js_cosh :: QDouble -> QDouble+foreign import javascript unsafe "new Float64Array(h$easytensor_qtanh($1))" js_tanh :: QDouble -> QDouble+++--------------------------------------------------------------------------+-- Eq+--------------------------------------------------------------------------++instance Eq QDouble where+ {-# INLINE (==) #-}+ (==) = js_eq+ {-# INLINE (/=) #-}+ (/=) = js_neq++++foreign import javascript unsafe "$1[0] === $2[0] && $1[1] === $2[1] && $1[2] === $2[2] && $1[3] === $2[3]"+ js_eq :: QDouble -> QDouble -> Bool+foreign import javascript unsafe "$1[0] !== $2[0] || $1[1] !== $2[1] || $1[2] !== $2[2] || $1[3] !== $2[3]"+ js_neq :: QDouble -> QDouble -> Bool++++--------------------------------------------------------------------------+-- Show+--------------------------------------------------------------------------++instance Show QDouble where+ show = unpack' . js_show++foreign import javascript unsafe "$1[3].toPrecision(8)\+ \ + ($1[0] >= 0 ? ' + ' : ' - ') + Math.abs($1[0]).toPrecision(8) + 'i'\+ \ + ($1[1] >= 0 ? ' + ' : ' - ') + Math.abs($1[1]).toPrecision(8) + 'j'\+ \ + ($1[2] >= 0 ? ' + ' : ' - ') + Math.abs($1[2]).toPrecision(8) + 'k'"+ js_show:: QDouble -> JSString
+ src-ghcjs/Numeric/Quaternion/QFloat.hs view
@@ -0,0 +1,295 @@+-- Note,+--+-- The whole module is made by copying Numeric.Quaternion.QDouble and replacing:+-- Double -> Float+-- Float64 -> Float32+--+-- If we are doing any refactoring of one of these modules, just do the same operation.++{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeSynonymInstances #-}+{-# LANGUAGE FlexibleInstances #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.Quaternion.QFloat+ ( QFloat, Quater (..)+ ) where++import Data.JSString (JSString, unpack')+import Data.Coerce (coerce)++import Numeric.Array+import Numeric.DataFrame.Type+import Numeric.Vector+import Numeric.Matrix++import Numeric.Quaternion.Class+++type QFloat = Quater Float++instance Quaternion Float where+ newtype Quater Float = QFloat (ArrayT Float '[4])+ {-# INLINE packQ #-}+ packQ = js_packQ+ {-# INLINE fromVecNum #-}+ fromVecNum = js_fromVecNum+ {-# INLINE fromVec4 #-}+ fromVec4 = coerce+ {-# INLINE toVec4 #-}+ toVec4 = coerce+ {-# INLINE unpackQ #-}+ unpackQ = js_unpackQ+ {-# INLINE square #-}+ square = js_square+ {-# INLINE im #-}+ im = js_im+ {-# INLINE re #-}+ re = js_re+ {-# INLINE imVec #-}+ imVec = js_imVec+ {-# INLINE taker #-}+ taker = js_taker+ {-# INLINE takei #-}+ takei = js_takei+ {-# INLINE takej #-}+ takej = js_takej+ {-# INLINE takek #-}+ takek = js_takek+ {-# INLINE conjugate #-}+ conjugate = js_conjugate+ {-# INLINE rotScale #-}+ rotScale = js_rotScale+ {-# INLINE getRotScale #-}+ getRotScale = js_getRotScale+ {-# INLINE axisRotation #-}+ axisRotation = js_axisRotation+ {-# INLINE qArg #-}+ qArg = js_qArg+ {-# INLINE fromMatrix33 #-}+ fromMatrix33 = js_fromMatrix33+ {-# INLINE fromMatrix44 #-}+ fromMatrix44 = js_fromMatrix44+ {-# INLINE toMatrix33 #-}+ toMatrix33 = js_toMatrix33+ {-# INLINE toMatrix44 #-}+ toMatrix44 = js_toMatrix44+++foreign import javascript unsafe "new Float32Array([$1,$2,$3,$4])"+ js_packQ :: Float -> Float -> Float -> Float -> QFloat++foreign import javascript unsafe "new Float32Array([$1[0],$1[1],$1[2],$2])"+ js_fromVecNum :: Vector Float 3 -> Float -> QFloat++foreign import javascript unsafe "$r1 = $1[0];$r2 = $1[1];$r3 = $1[2];$r4 = $1[3];"+ js_unpackQ :: QFloat -> (Float,Float,Float,Float)+++foreign import javascript unsafe "$1[0]*$1[0] + $1[1]*$1[1] + $1[2]*$1[2] + $1[3]*$1[3]"+ js_square :: QFloat -> Float++foreign import javascript unsafe "new Float32Array([$1[0],$1[1],$1[2],0])"+ js_im :: QFloat -> QFloat++foreign import javascript unsafe "$1.slice(0,3)"+ js_imVec :: QFloat -> Vector Float 3++foreign import javascript unsafe "$r = new Float32Array(4); $r[3] = $1[3];"+ js_re :: QFloat -> QFloat++foreign import javascript unsafe "$1[3]"+ js_taker :: QFloat -> Float++foreign import javascript unsafe "$1[0]"+ js_takei :: QFloat -> Float++foreign import javascript unsafe "$1[1]"+ js_takej :: QFloat -> Float++foreign import javascript unsafe "$1[2]"+ js_takek :: QFloat -> Float++foreign import javascript unsafe "new Float32Array([-$1[0],-$1[1],-$1[2],$1[3]])"+ js_conjugate :: QFloat -> QFloat+++foreign import javascript unsafe "new Float32Array(h$easytensor_rotScale($1,$2))"+ js_rotScale :: QFloat -> Vector Float 3 -> Vector Float 3++foreign import javascript unsafe "new Float32Array(h$easytensor_getRotScale($1,$2))"+ js_getRotScale :: Vector Float 3 -> Vector Float 3 -> QFloat++foreign import javascript unsafe "new Float32Array(h$easytensor_axisRotation($1,$2))"+ js_axisRotation :: Vector Float 3 -> Float -> QFloat++foreign import javascript unsafe "h$easytensor_qArg($1)"+ js_qArg :: QFloat -> Float++foreign import javascript unsafe "new Float32Array(h$easytensor_qfromMatrix33($1))"+ js_fromMatrix33 :: Matrix Float 3 3 -> QFloat++foreign import javascript unsafe "new Float32Array(h$easytensor_qfromMatrix44($1))"+ js_fromMatrix44 :: Matrix Float 4 4 -> QFloat++foreign import javascript unsafe "new Float32Array(h$easytensor_qtoMatrix33($1))"+ js_toMatrix33 :: QFloat -> Matrix Float 3 3++foreign import javascript unsafe "new Float32Array(h$easytensor_qtoMatrix44($1))"+ js_toMatrix44 :: QFloat -> Matrix Float 4 4+++--------------------------------------------------------------------------+-- Num+--------------------------------------------------------------------------++instance Num QFloat where+ {-# INLINE (+) #-}+ (+) = js_plus+ {-# INLINE (-) #-}+ (-) = js_minus+ {-# INLINE (*) #-}+ (*) = js_times+ {-# INLINE abs #-}+ abs = js_abs+ {-# INLINE signum #-}+ signum = js_signum+ {-# INLINE negate #-}+ negate = fromVec4 . negate . toVec4+ {-# INLINE fromInteger #-}+ fromInteger = js_toQuaternion . fromInteger++foreign import javascript unsafe "$1.map(function (e, i) {return e + $2[i];})"+ js_plus :: QFloat -> QFloat -> QFloat+++foreign import javascript unsafe "$1.map(function (e, i) {return e - $2[i];})"+ js_minus :: QFloat -> QFloat -> QFloat++foreign import javascript unsafe "new Float32Array(\+ \[ $1[3]*$2[0] + $1[0]*$2[3] + $1[1]*$2[2] - $1[2]*$2[1]\+ \, $1[3]*$2[1] - $1[0]*$2[2] + $1[1]*$2[3] + $1[2]*$2[0]\+ \, $1[3]*$2[2] + $1[0]*$2[1] - $1[1]*$2[0] + $1[2]*$2[3]\+ \, $1[3]*$2[3] - $1[0]*$2[0] - $1[1]*$2[1] - $1[2]*$2[2] ])"+ js_times :: QFloat -> QFloat -> QFloat++foreign import javascript unsafe "new Float32Array([0,0,0,Math.hypot($1[0],$1[1],$1[2],$1[3])])"+ js_abs :: QFloat -> QFloat++foreign import javascript unsafe "var l = Math.hypot($1[0],$1[1],$1[2],$1[3]); $r = (l == 0) ? (new Float32Array(4)) : $1.map(function (e) {return e/l;})"+ js_signum :: QFloat -> QFloat++foreign import javascript unsafe "$r = new Float32Array(4); $r[3] = $1;"+ js_toQuaternion :: Float -> QFloat++{-# RULES+"realToFrac/FloatQFloat" realToFrac = js_toQuaternion+"realToFrac/aQFloat" realToFrac = js_toQuaternion . realToFrac+"fromIntegral/aQFloat" fromIntegral = js_toQuaternion . fromIntegral+ #-}+++--------------------------------------------------------------------------+-- Fractional+--------------------------------------------------------------------------++instance Fractional QFloat where+ {-# INLINE recip #-}+ recip = js_recip+ {-# INLINE (/) #-}+ (/) p = js_times p . js_recip+ {-# INLINE fromRational #-}+ fromRational = js_toQuaternion . fromRational++foreign import javascript unsafe "new Float32Array(h$easytensor_qrecip($1))"+ js_recip :: QFloat -> QFloat+++--------------------------------------------------------------------------+-- Floating+--------------------------------------------------------------------------++instance Floating QFloat where+ {-# INLINE pi #-}+ pi = js_toQuaternion pi+ {-# INLINE exp #-}+ exp = js_exp+ {-# INLINE log #-}+ log = js_log+ {-# INLINE sqrt #-}+ sqrt = js_sqrt+ {-# INLINE sin #-}+ sin = js_sin+ {-# INLINE cos #-}+ cos = js_cos+ {-# INLINE tan #-}+ tan = js_tan+ {-# INLINE sinh #-}+ sinh = js_sinh+ {-# INLINE cosh #-}+ cosh = js_cosh+ {-# INLINE tanh #-}+ tanh = js_tanh+ {-# INLINE asin #-}+ asin q = -i * log (i*q + sqrt (1 - q*q))+ where i = signum . im $ q+ {-# INLINE acos #-}+ acos q = pi/2 - asin q+ {-# INLINE atan #-}+ atan q = if square imq == 0+ then js_toQuaternion (atan $ taker q)+ else i / 2 * log ( (i + q) / (i - q) )+ where i = signum imq+ imq = im q+ {-# INLINE asinh #-}+ asinh q = log (q + sqrt (q*q + 1))+ {-# INLINE acosh #-}+ acosh q = log (q + sqrt (q*q - 1))+ {-# INLINE atanh #-}+ atanh q = 0.5 * log ((1+q)/(1-q))++foreign import javascript unsafe "new Float32Array(h$easytensor_qexp($1))" js_exp :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qlog($1))" js_log :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qsqrt($1))" js_sqrt :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qsin($1))" js_sin :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qcos($1))" js_cos :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qtan($1))" js_tan :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qsinh($1))" js_sinh :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qcosh($1))" js_cosh :: QFloat -> QFloat+foreign import javascript unsafe "new Float32Array(h$easytensor_qtanh($1))" js_tanh :: QFloat -> QFloat+++--------------------------------------------------------------------------+-- Eq+--------------------------------------------------------------------------++instance Eq QFloat where+ {-# INLINE (==) #-}+ (==) = js_eq+ {-# INLINE (/=) #-}+ (/=) = js_neq++++foreign import javascript unsafe "$1[0] === $2[0] && $1[1] === $2[1] && $1[2] === $2[2] && $1[3] === $2[3]"+ js_eq :: QFloat -> QFloat -> Bool+foreign import javascript unsafe "$1[0] !== $2[0] || $1[1] !== $2[1] || $1[2] !== $2[2] || $1[3] !== $2[3]"+ js_neq :: QFloat -> QFloat -> Bool++++--------------------------------------------------------------------------+-- Show+--------------------------------------------------------------------------++instance Show QFloat where+ show = unpack' . js_show++foreign import javascript unsafe "$1[3].toPrecision(8)\+ \ + ($1[0] >= 0 ? ' + ' : ' - ') + Math.abs($1[0]).toPrecision(8) + 'i'\+ \ + ($1[1] >= 0 ? ' + ' : ' - ') + Math.abs($1[1]).toPrecision(8) + 'j'\+ \ + ($1[2] >= 0 ? ' + ' : ' - ') + Math.abs($1[2]).toPrecision(8) + 'k'"+ js_show:: QFloat -> JSString++
+ src-ghcjs/Numeric/Quaternion/Quaternion.js view
@@ -0,0 +1,205 @@+function h$easytensor_rotScale(quat, vec) {+ 'use strict';+ var i = quat[0], j = quat[1], k = quat[2], t = quat[3];+ var x = vec[0], y = vec[1], z = vec[2];+ var l = t*t - i*i - j*j - k*k;+ var d = 2*(i*x + j*y + k*z);+ t *= 2;+ return [ l*x + d*i + t*(z*j - y*k)+ , l*y + d*j + t*(x*k - z*i)+ , l*z + d*k + t*(y*i - x*j)+ ];+}++function h$easytensor_qArg(quat) {+ 'use strict';+ return Math.atan2( Math.hypot(quat[0],quat[1],quat[2]) , quat[3] ) * 2 ;+}++function h$easytensor_getRotScale(a, b) {+ 'use strict';+ if (b[0] === 0 && b[1] === 0 && b[2] === 0) { return [0,0,0,0];}+ if (a[0] === 0 && a[1] === 0 && a[2] === 0) { return [Infinity,Infinity,Infinity,Infinity];}+ var t = h$easytensor_cross(a, b);+ var ma = Math.hypot(a[0],a[1],a[2]);+ var mb = Math.hypot(b[0],b[1],b[2]);+ var dot = a[0]*b[0]+a[1]*b[1]+a[2]*b[2];+ if (t[0] === 0 && t[1] === 0 && t[2] === 0) {+ if (dot > 0) {return [0,0,0,Math.sqrt(mb/ma)];}+ else {return [0,0,Math.sqrt(mb/ma),0];}+ }+ var c = Math.sqrt(ma*mb + dot);+ ma *= Math.SQRT2;+ return [ t[0]/(ma*c)+ , t[1]/(ma*c)+ , t[2]/(ma*c)+ , c/ma+ ];+}++function h$easytensor_axisRotation(axis, a) {+ 'use strict';+ if (axis[0] === 0 && axis[1] === 0 && axis[2] === 0) { return [0,0,0,1];}+ var c = Math.cos(a*0.5), s = Math.sin(a*0.5) / Math.hypot(axis[0],axis[1],axis[2]);+ return [ axis[0]*s, axis[1]*s, axis[2]*s, c];+}++function h$easytensor_qfromMatrix33(m) {+ 'use strict';+ var d = Math.cbrt(+ m[0]*(m[4]*m[8]-m[5]*m[7])+ - m[1]*(m[3]*m[8]-m[5]*m[6])+ + m[2]*(m[3]*m[7]-m[4]*m[6]));+ return [ Math.sqrt(Math.max( 0, d + m[0] - m[4] - m[8] )) * Math.sign(m[5] - m[7]) * 0.5+ , Math.sqrt(Math.max( 0, d - m[0] + m[4] - m[8] )) * Math.sign(m[6] - m[2]) * 0.5+ , Math.sqrt(Math.max( 0, d - m[0] - m[4] + m[8] )) * Math.sign(m[1] - m[3]) * 0.5+ , Math.sqrt(Math.max( 0, d + m[0] + m[4] + m[8] )) * 0.5 ];+}++function h$easytensor_qfromMatrix44(m) {+ 'use strict';+ var d = Math.cbrt(+ m[0]*(m[5]*m[10]-m[6]*m[9])+ - m[1]*(m[4]*m[10]-m[6]*m[8])+ + m[2]*(m[4]*m[ 9]-m[5]*m[8]));+ return [ Math.sqrt(Math.max( 0, d + m[0] - m[5] - m[10] )) * Math.sign(m[6] - m[9]) * 0.5 / m[15]+ , Math.sqrt(Math.max( 0, d - m[0] + m[5] - m[10] )) * Math.sign(m[8] - m[2]) * 0.5 / m[15]+ , Math.sqrt(Math.max( 0, d - m[0] - m[5] + m[10] )) * Math.sign(m[1] - m[4]) * 0.5 / m[15]+ , Math.sqrt(Math.max( 0, d + m[0] + m[5] + m[10] )) * 0.5 / m[15] ];+}++function h$easytensor_qtoMatrix33(quat) {+ 'use strict';+ var x = quat[0], y = quat[1], z = quat[2], w = quat[3];+ var w2 = w*w;+ if (x === 0 && y === 0 && z === 0) {+ return [w2,0,0+ ,0,w2,0+ ,0,0,w2];+ }+ var x2 = x*x, y2 = y*y, z2 = z*z;+ var l2 = x2+y2+z2+w2;+ return [ l2 - 2*(z2 + y2), 2*(x*y + z*w), 2*(x*z - y*w)+ , 2*(x*y - z*w), l2 - 2*(z2 + x2), 2*(y*z + x*w)+ , 2*(x*z + y*w), 2*(y*z - x*w), l2 - 2*(y2 + x2) ]; +}++function h$easytensor_qtoMatrix44(quat) {+ 'use strict';+ var x = quat[0], y = quat[1], z = quat[2], w = quat[3];+ var w2 = w*w;+ if (x === 0 && y === 0 && z === 0) {+ return [w2,0,0,0+ ,0,w2,0,0+ ,0,0,w2,0+ ,0,0,0,1];+ }+ var x2 = x*x, y2 = y*y, z2 = z*z;+ var l2 = x2+y2+z2+w2;+ return [ l2 - 2*(z2 + y2), 2*(x*y + z*w), 2*(x*z - y*w), 0+ , 2*(x*y - z*w), l2 - 2*(z2 + x2), 2*(y*z + x*w), 0+ , 2*(x*z + y*w), 2*(y*z - x*w), l2 - 2*(y2 + x2), 0+ , 0, 0, 0, 1]; +}++function h$easytensor_qrecip(q) {+ 'use strict';+ var c = -1 / (q[0]*q[0] + q[1]*q[1] + q[2]*q[2] + q[3]*q[3]);+ return [q[0]*c,q[1]*c,q[2]*c,-q[3]*c];+}++function h$easytensor_qexp(q) {+ 'use strict';+ var mv = Math.hypot(q[0],q[1],q[2]), et = Math.exp(q[3]);+ if(mv === 0) {return [0,0,0,et];}+ var l = et * Math.sin(mv) / mv;+ return [q[0]*l,q[1]*l,q[2]*l,et*Math.cos(mv)];+}+function h$easytensor_qlog(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {+ if(q[3] >= 0){+ return [0,0,0,Math.log(q[3])];+ } else {+ return [Math.PI,0,0,Math.log(-q[3])];+ }+ }+ var mq = Math.sqrt(mv + q[3]*q[3]);+ mv = Math.sqrt(mv);+ var l = Math.atan2( mv, q[3] ) / mv;+ return [q[0]*l,q[1]*l,q[2]*l,Math.log(mq)];+}+function h$easytensor_qsqrt(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {+ if(q[3] >= 0){+ return [0,0,0,Math.sqrt(q[3])];+ } else {+ return [Math.sqrt(-q[3]),0,0,0];+ }+ }+ var l = Math.sqrt(mv + q[3]*q[3]);+ var l2 = Math.sqrt(l);+ var tq = q[3] / (l * 2);+ var sina = Math.sqrt(0.5 - tq) * l2 / Math.sqrt(mv);+ return [q[0]*sina,q[1]*sina,q[2]*sina,Math.sqrt(0.5 + tq) * l2];+}++// A good tutorial on complex number trigonometric functions is available here+// http://www.milefoot.com/math/complex/functionsofi.htm+// I extend it to complex numbers by replacing complex i with quaternion vector ijk++function h$easytensor_qsin(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {return [0,0,0,Math.sin(q[3])];}+ mv = Math.sqrt(mv);+ var l = Math.cos(q[3]) * Math.sinh(mv) / mv;+ return [q[0]*l,q[1]*l,q[2]*l, Math.sin(q[3])*Math.cosh(mv)];+}+function h$easytensor_qcos(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {return [0,0,0,Math.cos(q[3])];}+ mv = Math.sqrt(mv);+ var l = - Math.sin(q[3]) * Math.sinh(mv) / mv;+ return [q[0]*l,q[1]*l,q[2]*l, Math.cos(q[3])*Math.cosh(mv)];+}+function h$easytensor_qtan(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {return [0,0,0,Math.tan(q[3])];}+ mv = Math.sqrt(mv);+ var chv = Math.cosh(mv), shv = Math.sinh(mv), ct = Math.cos(q[3]), st = Math.sin(q[3]);+ var cq = 1 / (ct*ct*chv*chv + st*st*shv*shv);+ var l = chv * shv * cq / mv;+ return [q[0]*l,q[1]*l,q[2]*l, ct * st * cq];+}+function h$easytensor_qsinh(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {return [0,0,0,Math.sinh(q[3])];}+ mv = Math.sqrt(mv);+ var l = Math.cosh(q[3]) * Math.sin(mv) / mv;+ return [q[0]*l,q[1]*l,q[2]*l, Math.sinh(q[3])*Math.cos(mv)];+}+function h$easytensor_qcosh(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {return [0,0,0,Math.cosh(q[3])];}+ mv = Math.sqrt(mv);+ var l = Math.sinh(q[3]) * Math.sin(mv) / mv;+ return [q[0]*l,q[1]*l,q[2]*l, Math.cosh(q[3])*Math.cos(mv)];+}+function h$easytensor_qtanh(q) {+ 'use strict';+ var mv = q[0]*q[0] + q[1]*q[1] + q[2]*q[2];+ if(mv === 0) {return [0,0,0,Math.tanh(q[3])];}+ mv = Math.sqrt(mv);+ var cv = Math.cos(mv), sv = Math.sin(mv), cht = Math.cosh(q[3]), sht = Math.sinh(q[3]);+ var cq = 1 / (cht*cht*cv*cv + sht*sht*sv*sv);+ var l = cv * sv * cq / mv;+ return [q[0]*l,q[1]*l,q[2]*l, cht * sht * cq];+}
src/Numeric/Commons.hs view
@@ -12,6 +12,10 @@ {-# LANGUAGE TypeInType #-} {-# LANGUAGE UnboxedTuples #-} {-# LANGUAGE UndecidableInstances #-}+#ifdef ghcjs_HOST_OS+{-# LANGUAGE JavaScriptFFI #-}+{-# LANGUAGE UnliftedFFITypes #-}+#endif ----------------------------------------------------------------------------- -- | -- Module : Numeric.Commons@@ -25,13 +29,14 @@ module Numeric.Commons ( ElemRep, ElemPrim- , PrimBytes (..), FloatBytes, DoubleBytes, IntBytes, WordBytes+ , PrimBytes (..) ) where #include "MachDeps.h"-#include "HsBaseConfig.h" +#ifndef ghcjs_HOST_OS import GHC.Base (runRW#)+#endif import GHC.Int (Int16 (..), Int32 (..), Int64 (..), Int8 (..)) import GHC.Prim import GHC.Types (Double (..), Float (..), Int (..), RuntimeRep (..),@@ -47,7 +52,7 @@ type instance ElemRep Int8 = 'IntRep type instance ElemRep Int16 = 'IntRep type instance ElemRep Int32 = 'IntRep-#if SIZEOF_HSWORD < 8+#if WORD_SIZE_IN_BITS < 64 type instance ElemRep Int64 = 'Int64Rep #else type instance ElemRep Int64 = 'IntRep@@ -56,7 +61,7 @@ type instance ElemRep Word8 = 'WordRep type instance ElemRep Word16 = 'WordRep type instance ElemRep Word32 = 'WordRep-#if SIZEOF_HSWORD < 8+#if WORD_SIZE_IN_BITS < 64 type instance ElemRep Word64 = 'Word64Rep #else type instance ElemRep Word64 = 'WordRep@@ -69,7 +74,7 @@ type instance ElemPrim Int8 = Int# type instance ElemPrim Int16 = Int# type instance ElemPrim Int32 = Int#-#if SIZEOF_HSWORD < 8+#if WORD_SIZE_IN_BITS < 64 type instance ElemPrim Int64 = Int64# #else type instance ElemPrim Int64 = Int#@@ -78,18 +83,13 @@ type instance ElemPrim Word8 = Word# type instance ElemPrim Word16 = Word# type instance ElemPrim Word32 = Word#-#if SIZEOF_HSWORD < 8+#if WORD_SIZE_IN_BITS < 64 type instance ElemPrim Word64 = Word64# #else type instance ElemPrim Word64 = Word# #endif -type FloatBytes a = (PrimBytes a, ElemRep a ~ 'FloatRep , ElemPrim a ~ Float#)-type DoubleBytes a = (PrimBytes a, ElemRep a ~ 'DoubleRep, ElemPrim a ~ Double#)-type IntBytes a = (PrimBytes a, ElemRep a ~ 'IntRep , ElemPrim a ~ Int#)-type WordBytes a = (PrimBytes a, ElemRep a ~ 'WordRep , ElemPrim a ~ Word#)- -- | Facilities to convert to and from raw byte array. -- Warning! offsets and sizes are in elements, not in bytes! -- Therefore one must be really carefull if having a crazy idea of@@ -111,6 +111,12 @@ ix :: Int# -> a -> (ElemPrim a :: TYPE (ElemRep a)) instance PrimBytes Float where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapFloat v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapFloat arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(F# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeFloatArray# marr 0# x s1 of@@ -119,6 +125,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = F# (indexFloatArray# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_HSFLOAT# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_HSFLOAT#@@ -129,6 +136,12 @@ {-# INLINE ix #-} instance PrimBytes Double where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapDouble v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapDouble arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(D# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeDoubleArray# marr 0# x s1 of@@ -137,6 +150,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = D# (indexDoubleArray# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_HSDOUBLE# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_HSDOUBLE#@@ -147,6 +161,12 @@ {-# INLINE ix #-} instance PrimBytes Int where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapInt v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapInt arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(I# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeIntArray# marr 0# x s1 of@@ -155,6 +175,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = I# (indexIntArray# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_HSINT# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_HSINT#@@ -165,6 +186,12 @@ {-# INLINE ix #-} instance PrimBytes Int8 where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapInt8 v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapInt8 arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(I8# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeInt8Array# marr 0# x s1 of@@ -173,6 +200,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = I8# (indexInt8Array# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_INT8# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_INT8#@@ -183,6 +211,12 @@ {-# INLINE ix #-} instance PrimBytes Int16 where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapInt16 v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapInt16 arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(I16# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeInt16Array# marr 0# x s1 of@@ -191,6 +225,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = I16# (indexInt16Array# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_INT16# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_INT16#@@ -201,6 +236,12 @@ {-# INLINE ix #-} instance PrimBytes Int32 where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapInt32 v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapInt32 arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(I32# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeInt32Array# marr 0# x s1 of@@ -209,6 +250,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = I32# (indexInt32Array# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_INT32# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_INT32#@@ -239,6 +281,12 @@ #endif instance PrimBytes Word where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapWord v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapWord arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(W# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeWordArray# marr 0# x s1 of@@ -247,6 +295,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = W# (indexWordArray# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_HSWORD# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_HSWORD#@@ -257,6 +306,12 @@ {-# INLINE ix #-} instance PrimBytes Word8 where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapWord8 v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapWord8 arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(W8# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeWord8Array# marr 0# x s1 of@@ -265,6 +320,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = W8# (indexWord8Array# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_WORD8# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_WORD8#@@ -275,6 +331,12 @@ {-# INLINE ix #-} instance PrimBytes Word16 where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapWord16 v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapWord16 arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(W16# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeWord16Array# marr 0# x s1 of@@ -283,6 +345,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = W16# (indexWord16Array# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_WORD16# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_WORD16#@@ -293,6 +356,12 @@ {-# INLINE ix #-} instance PrimBytes Word32 where+#ifdef ghcjs_HOST_OS+ toBytes v = (# 0#, 1#, js_wrapWord32 v #)+ {-# INLINE toBytes #-}+ fromBytes (# off, _, arr #) = js_unwrapWord32 arr off+ {-# INLINE fromBytes #-}+#else toBytes v@(W32# x) = case runRW# ( \s0 -> case newByteArray# (byteSize v) s0 of (# s1, marr #) -> case writeWord32Array# marr 0# x s1 of@@ -301,6 +370,7 @@ {-# INLINE toBytes #-} fromBytes (# off, _, arr #) = W32# (indexWord32Array# arr off) {-# INLINE fromBytes #-}+#endif byteSize _ = SIZEOF_WORD32# {-# INLINE byteSize #-} byteAlign _ = ALIGNMENT_WORD32#@@ -331,3 +401,28 @@ {-# INLINE ix #-} #endif +#ifdef ghcjs_HOST_OS+foreign import javascript unsafe "h$wrapBuffer((new Float32Array([$1])).buffer)" js_wrapFloat :: Float -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Float64Array([$1])).buffer)" js_wrapDouble :: Double -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Int32Array([$1])).buffer)" js_wrapInt :: Int -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Int32Array([$1])).buffer)" js_wrapInt32 :: Int32 -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Int16Array([$1])).buffer)" js_wrapInt16 :: Int16 -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Int8Array([$1])).buffer)" js_wrapInt8 :: Int8 -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Uint32Array([$1])).buffer)" js_wrapWord :: Word -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Uint32Array([$1])).buffer)" js_wrapWord32 :: Word32 -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Uint16Array([$1])).buffer)" js_wrapWord16 :: Word16 -> ByteArray#+foreign import javascript unsafe "h$wrapBuffer((new Uint8Array([$1])).buffer)" js_wrapWord8 :: Word8 -> ByteArray#++++foreign import javascript unsafe "($1.f3 || new Float32Array($1.buf))[$2]" js_unwrapFloat :: ByteArray# -> Int# -> Float+foreign import javascript unsafe "($1.f6 || new Float64Array($1.buf))[$2]" js_unwrapDouble :: ByteArray# -> Int# -> Double+foreign import javascript unsafe "($1.i3 || new Int32Array($1.buf))[$2]" js_unwrapInt :: ByteArray# -> Int# -> Int+foreign import javascript unsafe "($1.i3 || new Int32Array($1.buf))[$2]" js_unwrapInt32 :: ByteArray# -> Int# -> Int32+foreign import javascript unsafe "($1.i1 || new Int16Array($1.buf))[$2]" js_unwrapInt16 :: ByteArray# -> Int# -> Int16+foreign import javascript unsafe "($1.i8 || new Int8Array($1.buf))[$2]" js_unwrapInt8 :: ByteArray# -> Int# -> Int8+foreign import javascript unsafe "($1.u3 || new Uint32Array($1.buf))[$2]" js_unwrapWord :: ByteArray# -> Int# -> Word+foreign import javascript unsafe "($1.u3 || new Uint32Array($1.buf))[$2]" js_unwrapWord32 :: ByteArray# -> Int# -> Word32+foreign import javascript unsafe "($1.u1 || new Uint16Array($1.buf))[$2]" js_unwrapWord16 :: ByteArray# -> Int# -> Word16+foreign import javascript unsafe "($1.u8 || new Uint8Array($1.buf))[$2]" js_unwrapWord8 :: ByteArray# -> Int# -> Word8+#endif
src/Numeric/DataFrame.hs view
@@ -41,3 +41,5 @@ import Numeric.Matrix import Numeric.Scalar import Numeric.Vector++
src/Numeric/DataFrame/IO.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE ExistentialQuantification #-} {-# LANGUAGE FlexibleContexts #-}@@ -9,6 +10,7 @@ {-# LANGUAGE TypeApplications #-} {-# LANGUAGE TypeFamilies #-} {-# LANGUAGE UnboxedTuples #-}+{-# LANGUAGE TypeOperators #-} ----------------------------------------------------------------------------- -- | -- Module : Numeric.DataFrame.IO@@ -17,22 +19,46 @@ -- -- Maintainer : chirkin@arch.ethz.ch ----- Mutable DataFrames living in IO.+---- Mutable DataFrames living in IO. -- ----------------------------------------------------------------------------- module Numeric.DataFrame.IO- ( MutableFrame (), IODataFrame ()+ (+#ifdef ghcjs_HOST_OS+ MutableFrame (), IODataFrame (..), MDataFrame (..), MutableArrayT (..)+#else+ MutableFrame (), IODataFrame ()+#endif+ , SomeIODataFrame (..) , newDataFrame, copyDataFrame, copyMutableDataFrame , unsafeFreezeDataFrame , freezeDataFrame, thawDataFrame , writeDataFrame, readDataFrame , writeDataFrameOff, readDataFrameOff+#ifdef ghcjs_HOST_OS+ -- * JavaScript-specific functions+ , MutableArrayBuffer+ , newArrayBuffer, arrayBuffer, viewFloatArray, viewDoubleArray+ , viewIntArray, viewInt32Array, viewInt16Array, viewInt8Array+ , viewWordArray, viewWord32Array, viewWord16Array, viewWord8Array, viewWord8ClampedArray+#endif ) where import GHC.Prim (RealWorld) import GHC.Types (Int (..), IO (..)) ++#ifdef ghcjs_HOST_OS+import Numeric.Array.Family hiding (Scalar)+import JavaScript.TypedArray.ArrayBuffer+import GHC.Prim+import Data.Int+import Data.Word+import Data.Maybe+import GHCJS.Types+import Numeric.DataFrame.Inference+#endif import Numeric.Commons import Numeric.DataFrame.Type import Numeric.DataFrame.Mutable@@ -42,32 +68,53 @@ -- | Mutable DataFrame that lives in IO. -- Internal representation is always a ByteArray. newtype IODataFrame t (ns :: [Nat]) = IODataFrame (MDataFrame RealWorld t (ns :: [Nat]))-+#ifdef ghcjs_HOST_OS+instance IsJSVal (IODataFrame t ds)+#endif+-- | Mutable DataFrame of unknown dimensionality+data SomeIODataFrame t (xns :: [XNat])+ = forall (ns :: [Nat])+ . ( FixedDim xns ns ~ ns+ , FixedXDim xns ns ~ xns+ , NumericFrame t ns+ )+ => SomeIODataFrame (IODataFrame t ns) -- | Create a new mutable DataFrame. newDataFrame :: forall t (ns :: [Nat])+#ifdef ghcjs_HOST_OS+ . ( ElemTypeInference t, Dimensions ns)+#else . ( PrimBytes t, Dimensions ns)+#endif => IO (IODataFrame t ns) newDataFrame = IODataFrame <$> IO (newDataFrame# @t @ns) {-# INLINE newDataFrame #-} -- | Copy one DataFrame into another mutable DataFrame at specified position.-copyDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat])- . ( PrimBytes (DataFrame t as)- , ConcatList as bs asbs- , Dimensions bs- )- => DataFrame t as -> Idx bs -> IODataFrame t asbs -> IO ()+copyDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat])+ . ( ConcatList as (b :+ bs) asbs, Dimensions (b :+ bs)+#ifdef ghcjs_HOST_OS+ , ArraySizeInference (as +: b'), Dimensions as+#else+ , PrimBytes (DataFrame t (as +: b'))+#endif+ )+ => DataFrame t (as +: b') -> Idx (b :+ bs) -> IODataFrame t asbs -> IO () copyDataFrame df ei (IODataFrame mdf) = IO (copyDataFrame# df ei mdf) {-# INLINE copyDataFrame #-} + -- | Copy one mutable DataFrame into another mutable DataFrame at specified position.-copyMutableDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat])+copyMutableDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) . ( PrimBytes t- , ConcatList as bs asbs- , Dimensions bs+ , ConcatList as (b :+ bs) asbs+ , Dimensions (b :+ bs)+#ifdef ghcjs_HOST_OS+ , Dimensions as+#endif )- => IODataFrame t as -> Idx bs -> IODataFrame t asbs -> IO ()+ => IODataFrame t (as +: b') -> Idx (b :+ bs) -> IODataFrame t asbs -> IO () copyMutableDataFrame (IODataFrame mdfA) ei (IODataFrame mdfB) = IO (copyMDataFrame# mdfA ei mdfB) {-# INLINE copyMutableDataFrame #-}@@ -75,7 +122,11 @@ -- | Make a mutable DataFrame immutable, without copying. unsafeFreezeDataFrame :: forall t (ns :: [Nat])- . PrimBytes (DataFrame t ns)+#ifdef ghcjs_HOST_OS+ . (MutableFrame t ns, ArraySizeInference ns)+#else+ . PrimBytes (DataFrame t ns)+#endif => IODataFrame t ns -> IO (DataFrame t ns) unsafeFreezeDataFrame (IODataFrame mdf) = IO (unsafeFreezeDataFrame# mdf) {-# INLINE unsafeFreezeDataFrame #-}@@ -83,14 +134,22 @@ -- | Copy content of a mutable DataFrame into a new immutable DataFrame. freezeDataFrame :: forall t (ns :: [Nat])- . PrimBytes (DataFrame t ns)- => IODataFrame t ns -> IO (DataFrame t ns)+#ifdef ghcjs_HOST_OS+ . (MutableFrame t ns, ArraySizeInference ns)+#else+ . PrimBytes (DataFrame t ns)+#endif+ => IODataFrame t ns -> IO (DataFrame t ns) freezeDataFrame (IODataFrame mdf) = IO (freezeDataFrame# mdf) {-# INLINE freezeDataFrame #-} -- | Create a new mutable DataFrame and copy content of immutable one in there. thawDataFrame :: forall t (ns :: [Nat])- . PrimBytes (DataFrame t ns)+#ifdef ghcjs_HOST_OS+ . (MutableFrame t ns, ArrayInstanceInference t ns)+#else+ . PrimBytes (DataFrame t ns)+#endif => DataFrame t ns -> IO (IODataFrame t ns) thawDataFrame df = IODataFrame <$> IO (thawDataFrame# df) {-# INLINE thawDataFrame #-}@@ -126,3 +185,168 @@ => IODataFrame t ns -> Int -> IO (Scalar t) readDataFrameOff (IODataFrame mdf) (I# i) = scalar <$> IO (readDataFrameOff# mdf i) {-# INLINE readDataFrameOff #-}+++#ifdef ghcjs_HOST_OS+newArrayBuffer :: Int -> IO MutableArrayBuffer+newArrayBuffer n = unsafeCoerce# <$> IO (newArrayBuffer# n)++viewFloatArray :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Float (AsXDims ds +: XN 0))+viewFloatArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Float))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewFloatArray# (jsval ab) :: IO (IODataFrame Float (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Float ds) pn+ of+ Evidence -> case inferNumericFrame @Float @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewDoubleArray :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Double (AsXDims ds +: XN 0))+viewDoubleArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Double))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewDoubleArray# (jsval ab) :: IO (IODataFrame Double (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Double ds) pn+ of+ Evidence -> case inferNumericFrame @Double @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewIntArray :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Int (AsXDims ds +: XN 0))+viewIntArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewIntArray# (jsval ab) :: IO (IODataFrame Int (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int ds) pn+ of+ Evidence -> case inferNumericFrame @Int @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewInt32Array :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Int32 (AsXDims ds +: XN 0))+viewInt32Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int32))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewInt32Array# (jsval ab) :: IO (IODataFrame Int32 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int32 ds) pn+ of+ Evidence -> case inferNumericFrame @Int32 @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewInt16Array :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Int16 (AsXDims ds +: XN 0))+viewInt16Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int16))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewInt16Array# (jsval ab) :: IO (IODataFrame Int16 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int16 ds) pn+ of+ Evidence -> case inferNumericFrame @Int16 @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewInt8Array :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Int8 (AsXDims ds +: XN 0))+viewInt8Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int8))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewInt8Array# (jsval ab) :: IO (IODataFrame Int8 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int8 ds) pn+ of+ Evidence -> case inferNumericFrame @Int8 @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewWordArray :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Word (AsXDims ds +: XN 0))+viewWordArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewWordArray# (jsval ab) :: IO (IODataFrame Word (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word ds) pn+ of+ Evidence -> case inferNumericFrame @Word @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewWord32Array :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Word32 (AsXDims ds +: XN 0))+viewWord32Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word32))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewWord32Array# (jsval ab) :: IO (IODataFrame Word32 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word32 ds) pn+ of+ Evidence -> case inferNumericFrame @Word32 @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewWord16Array :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Word16 (AsXDims ds +: XN 0))+viewWord16Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word16))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewWord16Array# (jsval ab) :: IO (IODataFrame Word16 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word16 ds) pn+ of+ Evidence -> case inferNumericFrame @Word16 @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewWord8Array :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Word8 (AsXDims ds +: XN 0))+viewWord8Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word8))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewWord8Array# (jsval ab) :: IO (IODataFrame Word8 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word8 ds) pn+ of+ Evidence -> case inferNumericFrame @Word8 @(ds +: n) of+ Evidence -> SomeIODataFrame df++viewWord8ClampedArray :: forall ds+ . ( Dimensions ds, ArraySizeInference ds)+ => MutableArrayBuffer -> IO (SomeIODataFrame Word8Clamped (AsXDims ds +: XN 0))+viewWord8ClampedArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word8Clamped))) (dim @ds) ab+ df <- fmap IODataFrame . IO $ viewWord8ClampedArray# (jsval ab) :: IO (IODataFrame Word8Clamped (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word8Clamped ds) pn+ of+ Evidence -> case inferNumericFrame @Word8Clamped @(ds +: n) of+ Evidence -> SomeIODataFrame df++arrayBuffer :: IODataFrame t ds -> IO MutableArrayBuffer+arrayBuffer (IODataFrame x) = unsafeCoerce# <$> IO (arrayBuffer# x)+++foreign import javascript unsafe "$1.length" js_abLength :: MutableArrayBuffer -> IO Int++abDim :: Int -> Dim (ds :: [Nat]) -> MutableArrayBuffer -> IO SomeDim+abDim elS d ab = fromMaybe (SomeDim (Dn :: Dim 0)) . someDimVal . (`quot` (elS * dimVal d)) <$> js_abLength ab++unsafeForceFixedDims :: forall ds n+ . Evidence ( FixedDim (AsXDims ds +: XN 0) (ds +: n) ~ (ds +: n)+ , FixedXDim (AsXDims ds +: XN 0) (ds +: n) ~ (AsXDims ds +: XN 0)+ )+unsafeForceFixedDims = unsafeCoerce# (Evidence :: Evidence ( (ds +: n) ~ (ds +: n) , (ds +: n) ~ (ds +: n) ))++#endif
− src/Numeric/DataFrame/Mutable.hs
@@ -1,223 +0,0 @@-{-# LANGUAGE DataKinds #-}-{-# LANGUAGE ExistentialQuantification #-}-{-# LANGUAGE FlexibleContexts #-}-{-# LANGUAGE FlexibleInstances #-}-{-# LANGUAGE KindSignatures #-}-{-# LANGUAGE MagicHash #-}-{-# LANGUAGE MultiParamTypeClasses #-}-{-# LANGUAGE ScopedTypeVariables #-}-{-# LANGUAGE TypeApplications #-}-{-# LANGUAGE TypeFamilies #-}-{-# LANGUAGE UnboxedTuples #-}--------------------------------------------------------------------------------- |--- Module : Numeric.DataFrame.Mutable--- Copyright : (c) Artem Chirkin--- License : BSD3------ Maintainer : chirkin@arch.ethz.ch------ Interfrace to perform primitive stateful operations on mutable frames.-----------------------------------------------------------------------------------module Numeric.DataFrame.Mutable- ( MutableFrame (..), MDataFrame ()- , newDataFrame#, copyDataFrame#, copyMDataFrame#, unsafeFreezeDataFrame#- , freezeDataFrame#, thawDataFrame#- , writeDataFrame#, readDataFrame#- ) where---import GHC.Int (Int16 (..), Int32 (..), Int64 (..),- Int8 (..))-import GHC.Prim-import GHC.Types (Double (..), Float (..), Int (..),- Word (..))-import GHC.Word (Word16 (..), Word32 (..), Word64 (..),- Word8 (..))--import Numeric.Commons-import Numeric.DataFrame.Type-import Numeric.Dimensions-import Numeric.TypeLits----- | Mutable DataFrame type-data MDataFrame s t (ns :: [Nat]) = MDataFrame# Int# Int# (MutableByteArray# s)---- | Create a new mutable DataFrame.-newDataFrame# :: forall t (ns :: [Nat]) s- . ( PrimBytes t, Dimensions ns)- => State# s -> (# State# s, MDataFrame s t ns #)-newDataFrame# s0- | elS <- elementByteSize (undefined :: t)- , I# n <- totalDim (Proxy @ns)- , (# s1, mba #) <- newByteArray# (n *# elS) s0- = (# s1, MDataFrame# 0# n mba #)-{-# INLINE newDataFrame# #-}---- | Copy one DataFrame into another mutable DataFrame at specified position.-copyDataFrame# :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) s- . ( PrimBytes (DataFrame t as)- , ConcatList as bs asbs- , Dimensions bs- )- => DataFrame t as -> Idx bs -> MDataFrame s t asbs -> State# s -> (# State# s, () #)-copyDataFrame# df ei (MDataFrame# offM _ arrM) s- | (# offA, lenA, arrA #) <- toBytes df- , elS <- elementByteSize df- , I# i <- fromEnum ei- = (# copyByteArray# arrA (offA *# elS) arrM ((offM +# i) *# elS) (lenA *# elS) s, () #)-{-# INLINE copyDataFrame# #-}---- | Copy one mutable DataFrame into another mutable DataFrame at specified position.-copyMDataFrame# :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) s- . ( PrimBytes t- , ConcatList as bs asbs- , Dimensions bs- )- => MDataFrame s t as -> Idx bs -> MDataFrame s t asbs -> State# s -> (# State# s, () #)-copyMDataFrame# (MDataFrame# offA lenA arrA) ei (MDataFrame# offM _ arrM) s- | elS <- elementByteSize (undefined :: t)- , I# i <- fromEnum ei- = (# copyMutableByteArray# arrA (offA *# elS) arrM ((offM +# i) *# elS) (lenA *# elS) s, () #)-{-# INLINE copyMDataFrame# #-}---- | Make a mutable DataFrame immutable, without copying.-unsafeFreezeDataFrame# :: forall t (ns :: [Nat]) s- . PrimBytes (DataFrame t ns)- => MDataFrame s t ns -> State# s -> (# State# s, DataFrame t ns #)-unsafeFreezeDataFrame# (MDataFrame# offM lenM arrM) s1- | (# s2, arrA #) <- unsafeFreezeByteArray# arrM s1- = (# s2, fromBytes (# offM, lenM, arrA #) #)-{-# INLINE unsafeFreezeDataFrame# #-}---- | Copy content of a mutable DataFrame into a new immutable DataFrame.-freezeDataFrame# :: forall t (ns :: [Nat]) s- . PrimBytes (DataFrame t ns)- => MDataFrame s t ns -> State# s -> (# State# s, DataFrame t ns #)-freezeDataFrame# (MDataFrame# offM n arrM) s0- | elS <- elementByteSize (undefined :: DataFrame t ns)- , (# s1, mba #) <- newByteArray# (n *# elS) s0- , s2 <- copyMutableByteArray# arrM (offM *# elS) mba 0# (n *# elS) s1- , (# s3, arrA #) <- unsafeFreezeByteArray# mba s2- = (# s3, fromBytes (# 0#, n, arrA #) #)-{-# INLINE freezeDataFrame# #-}---- | Create a new mutable DataFrame and copy content of immutable one in there.-thawDataFrame# :: forall t (ns :: [Nat]) s- . PrimBytes (DataFrame t ns)- => DataFrame t ns -> State# s -> (# State# s, MDataFrame s t ns #)-thawDataFrame# df s0- | elS <- elementByteSize (undefined :: DataFrame t ns)- , (# offA, n, arrA #) <- toBytes df- , (# s1, arrM #) <- newByteArray# (n *# elS) s0- , s2 <- copyByteArray# arrA (offA *# elS) arrM 0# (n *# elS) s1- = (# s2, MDataFrame# 0# n arrM #)-{-# INLINE thawDataFrame# #-}---- | Write a single element at the specified index-writeDataFrame# :: forall t (ns :: [Nat]) s- . ( MutableFrame t ns, Dimensions ns )- => MDataFrame s t ns -> Idx ns -> t -> State# s -> (# State# s, () #)-writeDataFrame# mdf ei x s | I# i <- fromEnum ei = (# writeDataFrameOff# mdf i x s, () #)-{-# INLINE writeDataFrame# #-}---- | Read a single element at the specified index-readDataFrame# :: forall t (ns :: [Nat]) s- . ( MutableFrame t ns, Dimensions ns )- => MDataFrame s t ns -> Idx ns -> State# s -> (# State# s, t #)-readDataFrame# mdf ei | I# i <- fromEnum ei = readDataFrameOff# mdf i-{-# INLINE readDataFrame# #-}--class MutableFrame t (ns :: [Nat]) where- -- | Write a single element at the specified element offset- writeDataFrameOff# :: MDataFrame s t ns -> Int# -> t -> State# s -> State# s- -- | Read a single element at the specified element offset- readDataFrameOff# :: MDataFrame s t ns -> Int# -> State# s -> (# State# s, t #)--instance MutableFrame Float (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (F# x) = writeFloatArray# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readFloatArray# arr (off +# i) s0 = (# s1, F# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Double (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (D# x) = writeDoubleArray# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readDoubleArray# arr (off +# i) s0 = (# s1, D# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Int (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (I# x) = writeIntArray# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readIntArray# arr (off +# i) s0 = (# s1, I# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Int8 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (I8# x) = writeInt8Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readInt8Array# arr (off +# i) s0 = (# s1, I8# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Int16 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (I16# x) = writeInt16Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readInt16Array# arr (off +# i) s0 = (# s1, I16# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Int32 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (I32# x) = writeInt32Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readInt32Array# arr (off +# i) s0 = (# s1, I32# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Int64 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (I64# x) = writeInt64Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readInt64Array# arr (off +# i) s0 = (# s1, I64# r #)- {-# INLINE readDataFrameOff# #-}---instance MutableFrame Word (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (W# x) = writeWordArray# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readWordArray# arr (off +# i) s0 = (# s1, W# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Word8 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (W8# x) = writeWord8Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readWord8Array# arr (off +# i) s0 = (# s1, W8# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Word16 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (W16# x) = writeWord16Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readWord16Array# arr (off +# i) s0 = (# s1, W16# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Word32 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (W32# x) = writeWord32Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readWord32Array# arr (off +# i) s0 = (# s1, W32# r #)- {-# INLINE readDataFrameOff# #-}--instance MutableFrame Word64 (ns :: [Nat]) where- writeDataFrameOff# (MDataFrame# off _ arr) i (W64# x) = writeWord64Array# arr (off +# i) x- {-# INLINE writeDataFrameOff# #-}- readDataFrameOff# (MDataFrame# off _ arr) i s0- | (# s1, r #) <- readWord64Array# arr (off +# i) s0 = (# s1, W64# r #)- {-# INLINE readDataFrameOff# #-}
src/Numeric/DataFrame/ST.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE ExistentialQuantification #-} {-# LANGUAGE FlexibleContexts #-}@@ -9,6 +10,7 @@ {-# LANGUAGE TypeApplications #-} {-# LANGUAGE TypeFamilies #-} {-# LANGUAGE UnboxedTuples #-}+{-# LANGUAGE TypeOperators #-} ----------------------------------------------------------------------------- -- | -- Module : Numeric.DataFrame.ST@@ -22,18 +24,41 @@ ----------------------------------------------------------------------------- module Numeric.DataFrame.ST- ( MutableFrame (), STDataFrame ()+ (+#ifdef ghcjs_HOST_OS+ MutableFrame (), STDataFrame (..), MDataFrame (..), MutableArrayT (..)+#else+ MutableFrame (), STDataFrame ()+#endif+ , SomeSTDataFrame (..) , newDataFrame, copyDataFrame, copyMutableDataFrame , unsafeFreezeDataFrame , freezeDataFrame, thawDataFrame , writeDataFrame, readDataFrame , writeDataFrameOff, readDataFrameOff+#ifdef ghcjs_HOST_OS+ -- * JavaScript-specific functions+ , STArrayBuffer+ , newArrayBuffer, arrayBuffer, viewFloatArray, viewDoubleArray+ , viewIntArray, viewInt32Array, viewInt16Array, viewInt8Array+ , viewWordArray, viewWord32Array, viewWord16Array, viewWord8Array, viewWord8ClampedArray+#endif ) where import GHC.Types (Int (..)) import GHC.ST (ST(..)) +#ifdef ghcjs_HOST_OS+import Numeric.Array.Family hiding (Scalar)+import JavaScript.TypedArray.ArrayBuffer.ST+import GHC.Prim+import Data.Int+import Data.Word+import Data.Maybe+import GHCJS.Types+import Numeric.DataFrame.Inference+#endif import Numeric.Commons import Numeric.DataFrame.Type import Numeric.DataFrame.Mutable@@ -44,32 +69,52 @@ -- | Mutable DataFrame that lives in ST. -- Internal representation is always a ByteArray. newtype STDataFrame s t (ns :: [Nat]) = STDataFrame (MDataFrame s t (ns :: [Nat]))-+#ifdef ghcjs_HOST_OS+instance IsJSVal (STDataFrame s t ds)+#endif+-- | Mutable DataFrame of unknown dimensionality+data SomeSTDataFrame s t (xns :: [XNat])+ = forall (ns :: [Nat])+ . ( FixedDim xns ns ~ ns+ , FixedXDim xns ns ~ xns+ , NumericFrame t ns+ )+ => SomeSTDataFrame (STDataFrame s t ns) -- | Create a new mutable DataFrame. newDataFrame :: forall t (ns :: [Nat]) s+#ifdef ghcjs_HOST_OS+ . ( ElemTypeInference t, Dimensions ns)+#else . ( PrimBytes t, Dimensions ns)+#endif => ST s (STDataFrame s t ns) newDataFrame = STDataFrame <$> ST (newDataFrame# @t @ns) {-# INLINE newDataFrame #-} -- | Copy one DataFrame into another mutable DataFrame at specified position.-copyDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) s- . ( PrimBytes (DataFrame t as)- , ConcatList as bs asbs- , Dimensions bs- )- => DataFrame t as -> Idx bs -> STDataFrame s t asbs -> ST s ()+copyDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s+ . ( ConcatList as (b :+ bs) asbs, Dimensions (b :+ bs)+#ifdef ghcjs_HOST_OS+ , ArraySizeInference (as +: b'), Dimensions as+#else+ , PrimBytes (DataFrame t (as +: b'))+#endif+ )+ => DataFrame t (as +: b') -> Idx (b :+ bs) -> STDataFrame s t asbs -> ST s () copyDataFrame df ei (STDataFrame mdf) = ST (copyDataFrame# df ei mdf) {-# INLINE copyDataFrame #-} -- | Copy one mutable DataFrame into another mutable DataFrame at specified position.-copyMutableDataFrame :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) s+copyMutableDataFrame :: forall t (as :: [Nat]) (b' :: Nat) (b :: Nat) (bs :: [Nat]) (asbs :: [Nat]) s . ( PrimBytes t- , ConcatList as bs asbs- , Dimensions bs+ , ConcatList as (b :+ bs) asbs+ , Dimensions (b :+ bs)+#ifdef ghcjs_HOST_OS+ , Dimensions as+#endif )- => STDataFrame s t as -> Idx bs -> STDataFrame s t asbs -> ST s ()+ => STDataFrame s t (as +: b') -> Idx (b :+ bs) -> STDataFrame s t asbs -> ST s () copyMutableDataFrame (STDataFrame mdfA) ei (STDataFrame mdfB) = ST (copyMDataFrame# mdfA ei mdfB) {-# INLINE copyMutableDataFrame #-}@@ -77,7 +122,11 @@ -- | Make a mutable DataFrame immutable, without copying. unsafeFreezeDataFrame :: forall t (ns :: [Nat]) s- . PrimBytes (DataFrame t ns)+#ifdef ghcjs_HOST_OS+ . (MutableFrame t ns, ArraySizeInference ns)+#else+ . PrimBytes (DataFrame t ns)+#endif => STDataFrame s t ns -> ST s (DataFrame t ns) unsafeFreezeDataFrame (STDataFrame mdf) = ST (unsafeFreezeDataFrame# mdf) {-# INLINE unsafeFreezeDataFrame #-}@@ -85,14 +134,22 @@ -- | Copy content of a mutable DataFrame into a new immutable DataFrame. freezeDataFrame :: forall t (ns :: [Nat]) s+#ifdef ghcjs_HOST_OS+ . (MutableFrame t ns, ArraySizeInference ns)+#else . PrimBytes (DataFrame t ns)+#endif => STDataFrame s t ns -> ST s (DataFrame t ns) freezeDataFrame (STDataFrame mdf) = ST (freezeDataFrame# mdf) {-# INLINE freezeDataFrame #-} -- | Create a new mutable DataFrame and copy content of immutable one in there. thawDataFrame :: forall t (ns :: [Nat]) s- . PrimBytes (DataFrame t ns)+#ifdef ghcjs_HOST_OS+ . (MutableFrame t ns, ArrayInstanceInference t ns)+#else+ . PrimBytes (DataFrame t ns)+#endif => DataFrame t ns -> ST s (STDataFrame s t ns) thawDataFrame df = STDataFrame <$> ST (thawDataFrame# df) {-# INLINE thawDataFrame #-}@@ -128,3 +185,170 @@ => STDataFrame s t ns -> Int -> ST s (Scalar t) readDataFrameOff (STDataFrame mdf) (I# i) = scalar <$> ST (readDataFrameOff# mdf i) {-# INLINE readDataFrameOff #-}+++#ifdef ghcjs_HOST_OS++newArrayBuffer :: Int -> ST s (STArrayBuffer s)+newArrayBuffer n = unsafeCoerce# <$> ST (newArrayBuffer# n)++viewFloatArray :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Float (AsXDims ds +: XN 0))+viewFloatArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Float))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewFloatArray# (jsval ab) :: ST s (STDataFrame s Float (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Float ds) pn+ of+ Evidence -> case inferNumericFrame @Float @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewDoubleArray :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Double (AsXDims ds +: XN 0))+viewDoubleArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Double))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewDoubleArray# (jsval ab) :: ST s (STDataFrame s Double (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Double ds) pn+ of+ Evidence -> case inferNumericFrame @Double @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewIntArray :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Int (AsXDims ds +: XN 0))+viewIntArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewIntArray# (jsval ab) :: ST s (STDataFrame s Int (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int ds) pn+ of+ Evidence -> case inferNumericFrame @Int @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewInt32Array :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Int32 (AsXDims ds +: XN 0))+viewInt32Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int32))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewInt32Array# (jsval ab) :: ST s (STDataFrame s Int32 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int32 ds) pn+ of+ Evidence -> case inferNumericFrame @Int32 @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewInt16Array :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Int16 (AsXDims ds +: XN 0))+viewInt16Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int16))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewInt16Array# (jsval ab) :: ST s (STDataFrame s Int16 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int16 ds) pn+ of+ Evidence -> case inferNumericFrame @Int16 @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewInt8Array :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Int8 (AsXDims ds +: XN 0))+viewInt8Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Int8))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewInt8Array# (jsval ab) :: ST s (STDataFrame s Int8 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Int8 ds) pn+ of+ Evidence -> case inferNumericFrame @Int8 @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewWordArray :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Word (AsXDims ds +: XN 0))+viewWordArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewWordArray# (jsval ab) :: ST s (STDataFrame s Word (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word ds) pn+ of+ Evidence -> case inferNumericFrame @Word @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewWord32Array :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Word32 (AsXDims ds +: XN 0))+viewWord32Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word32))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewWord32Array# (jsval ab) :: ST s (STDataFrame s Word32 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word32 ds) pn+ of+ Evidence -> case inferNumericFrame @Word32 @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewWord16Array :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Word16 (AsXDims ds +: XN 0))+viewWord16Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word16))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewWord16Array# (jsval ab) :: ST s (STDataFrame s Word16 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word16 ds) pn+ of+ Evidence -> case inferNumericFrame @Word16 @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewWord8Array :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Word8 (AsXDims ds +: XN 0))+viewWord8Array ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word8))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewWord8Array# (jsval ab) :: ST s (STDataFrame s Word8 (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word8 ds) pn+ of+ Evidence -> case inferNumericFrame @Word8 @(ds +: n) of+ Evidence -> SomeSTDataFrame df++viewWord8ClampedArray :: forall s ds+ . ( Dimensions ds, ArraySizeInference ds)+ => STArrayBuffer s -> ST s (SomeSTDataFrame s Word8Clamped (AsXDims ds +: XN 0))+viewWord8ClampedArray ab = do+ SomeDim (pn@Dn :: Dim (n :: Nat)) <- abDim (I# (byteSize (undefined :: Word8Clamped))) (dim @ds) ab+ df <- fmap STDataFrame . ST $ viewWord8ClampedArray# (jsval ab) :: ST s (STDataFrame s Word8Clamped (ds +: n))+ return $ case unsafeForceFixedDims @ds @n+ `sumEvs` inferSnocDimensions @ds @n+ `sumEvs` inferSnocArrayInstance (undefined :: DataFrame Word8Clamped ds) pn+ of+ Evidence -> case inferNumericFrame @Word8Clamped @(ds +: n) of+ Evidence -> SomeSTDataFrame df++arrayBuffer :: STDataFrame s t ds -> ST s (STArrayBuffer s)+arrayBuffer (STDataFrame x) = unsafeCoerce# <$> ST (arrayBuffer# x)+++foreign import javascript unsafe "$1.length" js_abLength :: STArrayBuffer s -> Int++abDim :: Int -> Dim (ds :: [Nat]) -> STArrayBuffer s -> ST s SomeDim+abDim elS d ab = fromMaybe (SomeDim (Dn :: Dim 0)) . someDimVal . (`quot` (elS * dimVal d)) <$> pure (js_abLength ab)+{-# NOINLINE abDim #-}++unsafeForceFixedDims :: forall ds n+ . Evidence ( FixedDim (AsXDims ds +: XN 0) (ds +: n) ~ (ds +: n)+ , FixedXDim (AsXDims ds +: XN 0) (ds +: n) ~ (AsXDims ds +: XN 0)+ )+unsafeForceFixedDims = unsafeCoerce# (Evidence :: Evidence ( (ds +: n) ~ (ds +: n) , (ds +: n) ~ (ds +: n) ))++#endif
src/Numeric/DataFrame/SubSpace.hs view
@@ -69,7 +69,7 @@ | asbs as -> bs, asbs bs -> as, as bs -> asbs where -- | Unsafely get a sub-dataframe by its primitive element subset. -- The offset is not checked to be aligned to the space structure or for bounds.- -- Arguments are zero-based element offset and element size (aka `totalDim` of sub dataframe)+ -- Arguments are zero-based primitive element offset and subset ("as" element) size (aka `totalDim` of sub dataframe) -- -- Normal indexing can be expressed in terms of `indexOffset#`: --@@ -152,13 +152,13 @@ ewfoldMap :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) m . (Monoid m, SubSpace t as bs asbs) => (DataFrame t as -> m) -> DataFrame t asbs -> m-ewfoldMap f = ewfoldl (\b -> mappend b . f) mempty+ewfoldMap f = ewfoldl (\m b -> m `seq` (mappend m $! f b)) mempty {-# INLINE ewfoldMap #-} iwfoldMap :: forall t (as :: [Nat]) (bs :: [Nat]) (asbs :: [Nat]) m . ( Monoid m, SubSpace t as bs asbs) => (Idx bs -> DataFrame t as -> m) -> DataFrame t asbs -> m-iwfoldMap f = iwfoldl (\i b -> mappend b . f i) mempty+iwfoldMap f = iwfoldl (\i m b -> m `seq` (mappend m $! f i b)) mempty {-# INLINE iwfoldMap #-} @@ -218,34 +218,44 @@ #endif {-# INLINE indexOffset# #-} + ewmap :: forall s (as' :: [Nat]) (asbs' :: [Nat])+ . SubSpace s as' bs asbs'+ => (DataFrame s as' -> DataFrame t as)+ -> DataFrame s asbs' -> DataFrame t asbs ewmap f df | elS <- elementByteSize (undefined :: DataFrame t asbs) , I# lenBS <- totalDim (Proxy @bs) , I# lenAS <- totalDim (Proxy @as)+ , I# lenAS' <- totalDim (Proxy @as') , lenASB <- lenAS *# elS = case runRW# ( \s0 -> case newByteArray# (lenAS *# lenBS *# elS) s0 of (# s1, marr #) -> case overDimOff_# (dim @bs)- ( \pos s -> case toBytes $ f (indexOffset# pos lenAS df) of- (# offX, _, arrX #) -> copyByteArray# arrX (offX *# elS) marr (pos *# elS) lenASB s- ) 0# lenAS s1 of+ ( \pos s -> case toBytes $ f (indexOffset# (pos *# lenAS') lenAS' df) of+ (# offX, _, arrX #) -> copyByteArray# arrX (offX *# elS) marr (pos *# lenASB) lenASB s+ ) 0# 1# s1 of s2 -> unsafeFreezeByteArray# marr s2 ) of (# _, r #) -> fromBytes (# 0#, lenAS *# lenBS, r #) {-# INLINE ewmap #-} + iwmap :: forall s (as' :: [Nat]) (asbs' :: [Nat])+ . SubSpace s as' bs asbs'+ => (Idx bs -> DataFrame s as' -> DataFrame t as)+ -> DataFrame s asbs' -> DataFrame t asbs iwmap f df | elS <- elementByteSize (undefined :: DataFrame t asbs) , I# lenBS <- totalDim (Proxy @bs) , I# lenAS <- totalDim (Proxy @as)+ , I# lenAS' <- totalDim (Proxy @as') , lenASB <- lenAS *# elS = case runRW# ( \s0 -> case newByteArray# (lenAS *# lenBS *# elS) s0 of (# s1, marr #) -> case overDim_# (dim @bs)- ( \i pos s -> case toBytes $ f i (indexOffset# pos lenAS df) of- (# offX, _, arrX #) -> copyByteArray# arrX (offX *# elS) marr (pos *# elS) lenASB s- ) 0# lenAS s1 of+ ( \i pos s -> case toBytes $ f i (indexOffset# (pos *# lenAS') lenAS' df) of+ (# offX, _, arrX #) -> copyByteArray# arrX (offX *# elS) marr (pos *# lenASB) lenASB s+ ) 0# 1# s1 of s2 -> unsafeFreezeByteArray# marr s2 ) of (# _, r #) -> fromBytes (# 0#, lenAS *# lenBS, r #)
src/Numeric/DataFrame/Type.hs view
@@ -38,10 +38,17 @@ , FPFRame, IntegralFrame, NumericVariantFrame, CommonOpFrame ) where +#include "MachDeps.h" import Data.Int (Int16, Int32, Int64, Int8) import Data.Word (Word16, Word32, Word64, Word8) import Foreign.Storable (Storable (..))-import GHC.Exts (Int (..), Ptr (..), Float#, Double#, Int#, Word#)+import GHC.Exts (Int (..), Ptr (..), Float#, Double#, Int#, Word#+#ifndef ghcjs_HOST_OS+#if WORD_SIZE_IN_BITS < 64+ , Int64#, Word64#+#endif+#endif+ ) import GHC.Prim (copyAddrToByteArray#, copyByteArrayToAddr#, newByteArray#, plusAddr#, quotInt#,@@ -210,8 +217,12 @@ type instance ElemPrim (DataFrame Int16 ds) = Int# type instance ElemPrim (DataFrame Int32 ds) = Int# #ifndef ghcjs_HOST_OS+#if WORD_SIZE_IN_BITS < 64+type instance ElemPrim (DataFrame Int64 ds) = Int64#+#else type instance ElemPrim (DataFrame Int64 ds) = Int# #endif+#endif type instance ElemPrim (DataFrame Word ds) = Word# type instance ElemPrim (DataFrame Word8 ds) = Word# type instance ElemPrim (DataFrame Word16 ds) = Word#@@ -219,8 +230,12 @@ #ifdef ghcjs_HOST_OS type instance ElemPrim (DataFrame Word8Clamped ds) = Int# #else+#if WORD_SIZE_IN_BITS < 64+type instance ElemPrim (DataFrame Word64 ds) = Word64#+#else type instance ElemPrim (DataFrame Word64 ds) = Word# #endif+#endif deriving instance ( PrimBytes (Array Float ds) , ElemPrim (Array Float ds) ~ Float# , ElemRep (Array Float ds) ~ 'FloatRep) => PrimBytes (DataFrame Float ds)@@ -241,9 +256,15 @@ , ElemRep (Array Int32 ds) ~ 'IntRep) => PrimBytes (DataFrame Int32 ds) #ifndef ghcjs_HOST_OS deriving instance ( PrimBytes (Array Int64 ds)+#if WORD_SIZE_IN_BITS < 64+ , ElemPrim (Array Int64 ds) ~ Int64#+ , ElemRep (Array Int64 ds) ~ 'Int64Rep+#else , ElemPrim (Array Int64 ds) ~ Int#- , ElemRep (Array Int64 ds) ~ 'IntRep) => PrimBytes (DataFrame Int64 ds)+ , ElemRep (Array Int64 ds) ~ 'IntRep #endif+ ) => PrimBytes (DataFrame Int64 ds)+#endif deriving instance ( PrimBytes (Array Word ds) , ElemPrim (Array Word ds) ~ Word# , ElemRep (Array Word ds) ~ 'WordRep) => PrimBytes (DataFrame Word ds)@@ -262,8 +283,14 @@ , ElemRep (Array Word8Clamped ds) ~ 'IntRep) => PrimBytes (DataFrame Word8Clamped ds) #else deriving instance ( PrimBytes (Array Word64 ds)+#if WORD_SIZE_IN_BITS < 64+ , ElemPrim (Array Word64 ds) ~ Word64#+ , ElemRep (Array Word64 ds) ~ 'Word64Rep+#else , ElemPrim (Array Word64 ds) ~ Word#- , ElemRep (Array Word64 ds) ~ 'WordRep) => PrimBytes (DataFrame Word64 ds)+ , ElemRep (Array Word64 ds) ~ 'WordRep+#endif+ ) => PrimBytes (DataFrame Word64 ds) #endif
src/Numeric/Matrix.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE ExistentialQuantification #-} {-# LANGUAGE FlexibleContexts #-}@@ -20,39 +21,40 @@ ( MatrixCalculus (..) , SquareMatrixCalculus (..) , MatrixInverse (..)+ , HomTransform4 (..) , Matrix , Mat22f, Mat23f, Mat24f , Mat32f, Mat33f, Mat34f , Mat42f, Mat43f, Mat44f+ , Mat22d, Mat23d, Mat24d+ , Mat32d, Mat33d, Mat34d+ , Mat42d, Mat43d, Mat44d , mat22, mat33, mat44 , (%*) ) where +++#ifdef ghcjs_HOST_OS+import Numeric.Array.Family (ElemTypeInference)+#endif+ import GHC.Types (Type) import Numeric.Commons import Numeric.DataFrame.Contraction ((%*)) import Numeric.DataFrame.Shape-import Numeric.Dimensions (Nat)-import Numeric.Matrix.Type+import Numeric.Dimensions (Nat, Idx (..))+import Numeric.Matrix.Class+import Numeric.Matrix.Mat44d ()+import Numeric.Matrix.Mat44f () import Numeric.Vector import Control.Monad.ST import Numeric.DataFrame.ST --- Type abbreviations -type Mat22f = Matrix Float 2 2-type Mat32f = Matrix Float 3 2-type Mat42f = Matrix Float 4 2-type Mat23f = Matrix Float 2 3-type Mat33f = Matrix Float 3 3-type Mat43f = Matrix Float 4 3-type Mat24f = Matrix Float 2 4-type Mat34f = Matrix Float 3 4-type Mat44f = Matrix Float 4 4 - -- | Compose a 2x2D matrix mat22 :: ( PrimBytes (Vector t 2) , PrimBytes (Matrix t 2 2)@@ -61,30 +63,40 @@ mat22 = (<::>) -- | Compose a 3x3D matrix-mat33 :: ( PrimBytes t+mat33 :: (+#ifdef ghcjs_HOST_OS+ ElemTypeInference t, MutableFrame t '[3,3]+#else+ PrimBytes t , PrimBytes (Vector t 3) , PrimBytes (Matrix t 3 3)+#endif ) => Vector t 3 -> Vector t 3 -> Vector t 3 -> Matrix t 3 3 mat33 a b c = runST $ do mmat <- newDataFrame- copyDataFrame a 1 mmat- copyDataFrame b 2 mmat- copyDataFrame c 3 mmat+ copyDataFrame a (1:!1:!Z) mmat+ copyDataFrame b (1:!2:!Z) mmat+ copyDataFrame c (1:!3:!Z) mmat unsafeFreezeDataFrame mmat -- | Compose a 4x4D matrix mat44 :: forall (t :: Type)- . ( PrimBytes t+ . (+#ifdef ghcjs_HOST_OS+ ElemTypeInference t, MutableFrame t '[4,4]+#else+ PrimBytes t , PrimBytes (Vector t (4 :: Nat)) , PrimBytes (Matrix t (4 :: Nat) (4 :: Nat))+#endif ) => Vector t (4 :: Nat) -> Vector t (4 :: Nat) -> Vector t (4 :: Nat) -> Vector t (4 :: Nat) -> Matrix t (4 :: Nat) (4 :: Nat) mat44 a b c d = runST $ do mmat <- newDataFrame- copyDataFrame a 1 mmat- copyDataFrame b 2 mmat- copyDataFrame c 3 mmat- copyDataFrame d 4 mmat+ copyDataFrame a (1:!1:!Z) mmat+ copyDataFrame b (1:!2:!Z) mmat+ copyDataFrame c (1:!3:!Z) mmat+ copyDataFrame d (1:!4:!Z) mmat unsafeFreezeDataFrame mmat
+ src/Numeric/Matrix/Class.hs view
@@ -0,0 +1,114 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE KindSignatures #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE TypeFamilies #-}+module Numeric.Matrix.Class+ ( MatrixCalculus (..)+ , SquareMatrixCalculus (..)+ , MatrixInverse (..)+ , Matrix+ , HomTransform4 (..)+ , Mat22f, Mat23f, Mat24f+ , Mat32f, Mat33f, Mat34f+ , Mat42f, Mat43f, Mat44f+ , Mat22d, Mat23d, Mat24d+ , Mat32d, Mat33d, Mat34d+ , Mat42d, Mat43d, Mat44d+ ) where++import Numeric.Commons+import Numeric.DataFrame.Type+import Numeric.Dimensions (Nat)+import Numeric.Scalar+import Numeric.Vector++-- | Alias for DataFrames of rank 2+type Matrix t (n :: Nat) (m :: Nat) = DataFrame t '[n,m]++class MatrixCalculus t (n :: Nat) (m :: Nat) where+ -- | Transpose Mat+ transpose :: (MatrixCalculus t m n, PrimBytes (Matrix t m n)) => Matrix t n m -> Matrix t m n+++class SquareMatrixCalculus t (n :: Nat) where+ -- | Mat with 1 on diagonal and 0 elsewhere+ eye :: Matrix t n n+ -- | Put the same value on the Mat diagonal, 0 otherwise+ diag :: Scalar t -> Matrix t n n+ -- | Determinant of Mat+ det :: Matrix t n n -> Scalar t+ -- | Sum of diagonal elements+ trace :: Matrix t n n -> Scalar t++class MatrixInverse t (n :: Nat) where+ -- | Matrix inverse+ inverse :: DataFrame t '[n,n] -> DataFrame t '[n,n]+++-- | Operations on 4x4 transformation matrices and vectors in homogeneous coordinates.+-- All angles are specified in radians.+class HomTransform4 t where+ -- | Create a translation matrix from a vector+ translate4 :: Vector t 4 -> Matrix t 4 4+ -- | Create a translation matrix from a vector+ translate3 :: Vector t 3 -> Matrix t 4 4+ -- | Rotation matrix for a rotation around the X axis, angle is given in radians.+ rotateX :: t -> Matrix t 4 4+ -- | Rotation matrix for a rotation around the Y axis, angle is given in radians.+ rotateY :: t -> Matrix t 4 4+ -- | Rotation matrix for a rotation around the Z axis, angle is given in radians.+ rotateZ :: t -> Matrix t 4 4+ -- | Rotation matrix for a rotation around an arbitrary normalized vector+ rotate :: Vector t 3 -> t -> Matrix t 4 4+ -- | Rotation matrix from the Euler angles yaw pitch and roll+ rotateEuler :: t -> t -> t -> Matrix t 4 4+ -- | Create a transform matrix using up direction, camera position and a point to look at. + -- Just the same as GluLookAt.+ lookAt :: Vector t 3 -- ^ The up direction, not necessary unit length or perpendicular to the view vector+ -> Vector t 3 -- ^ The viewers position+ -> Vector t 3 -- ^ The point to look at+ -> Matrix t 4 4+ -- | A perspective symmetric projection matrix. Right-handed coordinate system. (@x@ - right, @y@ - top) + -- http://en.wikibooks.org/wiki/GLSL_Programming/Vertex_Transformations+ perspective :: t -- ^ Near plane clipping distance (always positive) + -> t -- ^ Far plane clipping distance (always positive) + -> t -- ^ Field of view of the y axis, in radians + -> t -- ^ Aspect ratio, i.e. screen's width\/height + -> Matrix t 4 4+ -- | An orthogonal symmetric projection matrix. Right-handed coordinate system. (@x@ - right, @y@ - top) + -- http://en.wikibooks.org/wiki/GLSL_Programming/Vertex_Transformations+ orthogonal :: t -- ^ Near plane clipping distance + -> t -- ^ Far plane clipping distance + -> t -- ^ width + -> t -- ^ height + -> Matrix t 4 4+ -- | Add one more dimension and set it to 1.+ toHomPoint :: Vector t 3 -> Vector t 4+ -- | Add one more dimension and set it to 0.+ toHomVector :: Vector t 3 -> Vector t 4+ -- | Transform a homogenous vector or point into a normal 3D vector.+ -- If the last coordinate is not zero, divide the rest by it.+ fromHom :: Vector t 4 -> Vector t 3++-- Type abbreviations++type Mat22f = Matrix Float 2 2+type Mat32f = Matrix Float 3 2+type Mat42f = Matrix Float 4 2+type Mat23f = Matrix Float 2 3+type Mat33f = Matrix Float 3 3+type Mat43f = Matrix Float 4 3+type Mat24f = Matrix Float 2 4+type Mat34f = Matrix Float 3 4+type Mat44f = Matrix Float 4 4++type Mat22d = Matrix Double 2 2+type Mat32d = Matrix Double 3 2+type Mat42d = Matrix Double 4 2+type Mat23d = Matrix Double 2 3+type Mat33d = Matrix Double 3 3+type Mat43d = Matrix Double 4 3+type Mat24d = Matrix Double 2 4+type Mat34d = Matrix Double 3 4+type Mat44d = Matrix Double 4 4
− src/Numeric/Matrix/Type.hs
@@ -1,49 +0,0 @@-{-# LANGUAGE DataKinds #-}-{-# LANGUAGE FlexibleContexts #-}-{-# LANGUAGE KindSignatures #-}-{-# LANGUAGE MultiParamTypeClasses #-}-{-# LANGUAGE TypeFamilies #-}--------------------------------------------------------------------------------- |--- Module : Numeric.Matrix.Type--- Copyright : (c) Artem Chirkin--- License : BSD3------ Maintainer : chirkin@arch.ethz.ch--------------------------------------------------------------------------------------module Numeric.Matrix.Type- ( MatrixCalculus (..)- , SquareMatrixCalculus (..)- , MatrixInverse (..)- , Matrix- ) where--import Numeric.Commons-import Numeric.DataFrame.Type-import Numeric.Dimensions (Nat)-import Numeric.Scalar---- | Alias for DataFrames of rank 2-type Matrix t (n :: Nat) (m :: Nat) = DataFrame t '[n,m]--class MatrixCalculus t (n :: Nat) (m :: Nat) where- -- | Transpose Mat- transpose :: (MatrixCalculus t m n, PrimBytes (Matrix t m n)) => Matrix t n m -> Matrix t m n---class SquareMatrixCalculus t (n :: Nat) where- -- | Mat with 1 on diagonal and 0 elsewhere- eye :: Matrix t n n- -- | Put the same value on the Mat diagonal, 0 otherwise- diag :: Scalar t -> Matrix t n n- -- | Determinant of Mat- det :: Matrix t n n -> Scalar t- -- | Sum of diagonal elements- trace :: Matrix t n n -> Scalar t--class MatrixInverse t (n :: Nat) where- -- | Matrix inverse- inverse :: DataFrame t '[n,n] -> DataFrame t '[n,n]
+ src/Numeric/Quaternion.hs view
@@ -0,0 +1,35 @@+-----------------------------------------------------------------------------+-- |+-- Module : Numeric.Quaternion+-- Copyright : (c) Artem Chirkin+-- License : BSD3+--+-- Maintainer : chirkin@arch.ethz.ch+--+-- Quaternion operations implemented for Floats and Doubles.+--+-- The types `QDouble` and `QFloat` have the same representation as corresponding `Vector t 4`.+-- This means, you can do a cheap conversion between the types.+--+-- However, arithmetic instances, such as Num and Floating are implemented in substentially different ways.+-- For example, fromInteger fills a vector fully but sets only real part to a quaternion:+--+-- >>> 1 = vec4 1 1 1 1+-- >>> 1 = packQ 0 0 0 1+--+-- All other numeric operations for vectors are element-wise, but for quaternions I have implemented+-- the actual quaternion math.+-- Some of the operations (such as trigonometric operations) are ambiguous for quaternions;+-- the general rules I follow:+--+-- 1. Preserve imaginary vector axis same if possible;+-- 2. If both @+q@ and @-q@ are possible, prefer real value positive (@re q >= 0@).+-----------------------------------------------------------------------------+module Numeric.Quaternion+ ( Quaternion (..)+ , QDouble, QFloat+ ) where++import Numeric.Quaternion.Class+import Numeric.Quaternion.QFloat+import Numeric.Quaternion.QDouble
+ src/Numeric/Quaternion/Class.hs view
@@ -0,0 +1,85 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+module Numeric.Quaternion.Class+ ( Quaternion (..)+ ) where++import Numeric.Matrix (Matrix)+import Numeric.Vector (Vector)++-- | Quaternion operations+class Quaternion t where+ -- | Quaternion data type. The ordering of coordinates is @x y z w@,+ -- where @w@ is an argument, and @x y z@ are components of a 3D vector+ data Quater t+ -- | Set the quaternion in format (x,y,z,w)+ packQ :: t -> t -> t -> t -> Quater t+ -- | Get the values of the quaternion in format (x,y,z,w)+ unpackQ :: Quater t -> (t,t,t,t)+ -- | Set the quaternion from 3D axis vector and argument+ fromVecNum :: Vector t 3 -> t -> Quater t+ -- | Set the quaternion from 4D vector in format (x,y,z,w)+ fromVec4 :: Vector t 4 -> Quater t+ -- | Transform the quaternion to 4D vector in format (x,y,z,w)+ toVec4 :: Quater t -> Vector t 4+ -- | Get scalar square of the quaternion.+ --+ -- >> realToFrac (square q) == q * conjugate q+ square :: Quater t -> t+ -- | Imagine part of quaternion (orientation vector)+ im :: Quater t -> Quater t+ -- | Real part of the quaternion+ re :: Quater t -> Quater t+ -- | Get imagenery 3D vector of the quaternion+ imVec :: Quater t -> Vector t 3+ -- | Real part of the quaternion into number+ taker :: Quater t -> t+ -- | i-th component+ takei :: Quater t -> t+ -- | j-th component+ takej :: Quater t -> t+ -- | k-th component+ takek :: Quater t -> t+ -- | Conjugate quaternion (negate imaginary part)+ conjugate :: Quater t -> Quater t+ -- | Rotates and scales vector in 3D using quaternion.+ -- Let @q = (cos (a\/2), sin (a\/2) * v)@; then the rotation angle is @a@, and the axis of rotation is @v@.+ -- Scaling is proportional to @|v|^2@.+ --+ -- >>> rotScale q x == q * x * (conjugate q)+ rotScale :: Quater t -> Vector t 3 -> Vector t 3+ -- | Creates a quaternion @q@ from two vectors @a@ and @b@.+ -- @rotScale q a == b@+ getRotScale :: Vector t 3 -> Vector t 3 -> Quater t+ -- | Creates a rotation versor from an axis vector and an angle in radians.+ -- Result is always a unit quaternion (versor).+ -- If the argument vector is zero, then result is a real unit quaternion.+ axisRotation :: Vector t 3 -> t -> Quater t+ -- | Quaternion rotation angle+ --+ -- >>> q /= 0 ==> axisRotation (imVec q) (qArg q) == signum q+ qArg :: Quater t -> t+ -- | Create a quaternion from a rotation matrix.+ -- Note, that rotations of `q` and `-q` are equivalent, there result of this function may be+ -- ambiguious. I decided to force its real part be positive:+ --+ -- >>> taker (fromMatrix33 m) >= 0+ fromMatrix33 :: Matrix t 3 3 -> Quater t+ -- | Create a quaternion from a homogenious coordinates trasform matrix.+ -- Ignores matrix translation transform.+ -- Note, that rotations of `q` and `-q` are equivalent, there result of this function may be+ -- ambiguious. I decided to force its real part be positive:+ --+ -- >>> taker (fromMatrix44 m) >= 0+ fromMatrix44 :: Matrix t 4 4 -> Quater t+ -- | Create a rotation matrix from a quaternion.+ -- Note, that rotations of `q` and `-q` are equivalent, so the following property holds:+ --+ -- >>> toMatrix33 q == toMatrix33 (-q)+ toMatrix33 :: Quater t -> Matrix t 3 3+ -- | Create a homogenious coordinates trasform matrix from a quaternion.+ -- Translation of the output matrix is zero.+ -- Note, that rotations of `q` and `-q` are equivalent, so the following property holds:+ --+ -- >>> toMatrix44 q == toMatrix44 (-q)+ toMatrix44 :: Quater t -> Matrix t 4 4
src/Numeric/Vector.hs view
@@ -20,8 +20,10 @@ -- * Common operations , (.*.), dot, (·) , normL1, normL2, normLPInf, normLNInf, normLP+ , normalized , vec2, vec3, vec4 , det2, cross, (×)+ , unpackV2, unpackV3, unpackV4 ) where import Numeric.Array.ElementWise@@ -86,6 +88,16 @@ => Vector t n -> Scalar t normL2 = scalar . sqrt . ewfoldr (const (\a -> (a*a +))) 0 +-- | Normalize vector w.r.t. Euclidean metric (L2).+normalized :: ( Floating t+ , Fractional (Vector t n)+ , ElementWise (Idx '[n]) t (Vector t n)+ )+ => Vector t n -> Vector t n+normalized v = v / n+ where+ n = broadcast . sqrt $ ewfoldr (const (\a -> (a*a +))) 0 v+ -- | Maximum of absolute values normLPInf :: ( Ord t, Num t , ElementWise (Idx '[n]) t (Vector t n)@@ -169,3 +181,21 @@ f (2 :! Z) = b f (3 :! Z) = c f _ = d+++unpackV2 :: ElementWise (Idx '[2]) t (Vector t 2)+ => Vector t 2 -> (t, t)+unpackV2 v = (v ! 1, v ! 2)+{-# INLINE unpackV2 #-}+++unpackV3 :: ElementWise (Idx '[3]) t (Vector t 3)+ => Vector t 3 -> (t, t, t)+unpackV3 v = (v ! 1, v ! 2, v ! 3)+{-# INLINE unpackV3 #-}+++unpackV4 :: ElementWise (Idx '[4]) t (Vector t 4)+ => Vector t 4 -> (t, t, t, t)+unpackV4 v = (v ! 1, v ! 2, v ! 3, v ! 4)+{-# INLINE unpackV4 #-}
+ test/Numeric/QuaternionTest.hs view
@@ -0,0 +1,124 @@+{-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE TypeSynonymInstances #-}+{-# LANGUAGE FlexibleInstances #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}+module Numeric.QuaternionTest (runTests) where++import Test.QuickCheck+import Numeric.Quaternion+import Numeric.Vector+++instance (Quaternion t, Arbitrary t) => Arbitrary (Quater t) where+ arbitrary = packQ <$> arbitrary <*> arbitrary <*> arbitrary <*> arbitrary+ shrink q | (x,y,z,t) <- unpackQ q = packQ <$> shrink x <*> shrink y <*> shrink z <*> shrink t++instance Arbitrary Vec3d where+ arbitrary = vec3 <$> arbitrary <*> arbitrary <*> arbitrary+ shrink v | (x,y,z) <- unpackV3 v = vec3 <$> shrink x <*> shrink y <*> shrink z+++(=~=) :: (Quaternion a, Num a, Ord a, Num (Quater a), Fractional a) => Quater a -> Quater a -> Bool+(=~=) a b = taker (abs (a - b)) <= eps+ where+ s = max 1e-6 $ max (taker (abs a)) (taker (abs b))+ eps = s * 1e-6+infix 4 =~=+++prop_Eq :: QDouble -> Bool+prop_Eq q = and+ [ (\(x,y,z,t) -> packQ x y z t) (unpackQ q) == q+ , packQ (takei q) (takej q) (takek q) (taker q) == q+ , fromVecNum (imVec q) (taker q) == q+ , im q + re q == q+ ]+++prop_DoubleConjugate :: QDouble -> Bool+prop_DoubleConjugate q = conjugate (conjugate q) == q++prop_Square :: QDouble -> Bool+prop_Square q = q * conjugate q =~= realToFrac (square q)+++prop_RotScale :: QDouble -> Vec3d -> Bool+prop_RotScale q v = fromVecNum (rotScale q v) 0 =~= q * fromVecNum v 0 * conjugate q++prop_GetRotScale :: Vec3d -> Vec3d -> Bool+prop_GetRotScale a b = fromVecNum (rotScale q a) 0 =~= fromVecNum b 0+ where+ q = getRotScale a b+++prop_InverseRotScale :: QDouble -> Vec3d -> Bool+prop_InverseRotScale q v | q /= 0 = fromVecNum (rotScale (1/q) (rotScale q v)) 0 =~= fromVecNum v 0+ | otherwise = True++prop_FromToMatrix33 :: QDouble -> Bool+prop_FromToMatrix33 q | q == 0 = True+ | otherwise = fromMatrix33 (toMatrix33 (signum q)) =~= signum q * signum (re q)++prop_FromToMatrix44 :: QDouble -> Bool+prop_FromToMatrix44 q | q == 0 = True+ | otherwise = fromMatrix44 (toMatrix44 (signum q)) =~= signum q * signum (re q)++prop_RotationArg :: QDouble -> Bool+prop_RotationArg q | q == 0 = axisRotation (imVec q) (qArg q) =~= 1+ | otherwise = axisRotation (imVec q) (qArg q) =~= signum q+++prop_UnitQ :: QDouble -> Bool+prop_UnitQ q | q == 0 = True+ | otherwise = realToFrac (square (q / q)) =~= (1 :: QDouble)+++prop_ExpLog :: QDouble -> Bool+prop_ExpLog q | q == 0 = log (exp q) == q+ | otherwise = exp (log q) =~= q++prop_SinAsin :: QDouble -> Bool+prop_SinAsin q = sin (asin q') =~= q'+ where+ q' = signum q++prop_CosAcos :: QDouble -> Bool+prop_CosAcos q = cos (acos q') =~= q'+ where+ q' = signum q++prop_TanAtan :: QDouble -> Bool+prop_TanAtan q = tan (atan q') =~= q'+ where -- protect agains big numbers and rounding errors+ q' = if square q >= 1e6 then signum q else q+++prop_SinhAsinh :: QDouble -> Bool+prop_SinhAsinh q = sinh (asinh q') =~= q'+ where -- protect agains big numbers and rounding errors+ q' = if square q >= 1e6 then signum q else q++prop_CoshAcosh :: QDouble -> Bool+prop_CoshAcosh q = cosh (acosh q') =~= q'+ where -- protect agains big numbers and rounding errors+ q' = if square q >= 1e6 then signum q else q++prop_TanhAtanh :: QDouble -> Bool+prop_TanhAtanh q = tanh (atanh q') =~= q'+ where -- protect agains big numbers and rounding errors+ q' = if square q >= 1e6 then signum q else q++prop_SinCos :: QDouble -> Bool+prop_SinCos q = sin q' * sin q' + cos q' * cos q' =~= 1+ where+ q' = signum q++prop_SinhCosh :: QDouble -> Bool+prop_SinhCosh q = cosh q' * cosh q' - sinh q' * sinh q' =~= 1+ where+ q' = signum q++return []+runTests :: IO Bool+runTests = $quickCheckAll
test/Spec.hs view
@@ -1,17 +1,40 @@-module Spec (tests) where+module Main (tests, main) where +import System.Exit import Distribution.TestSuite import qualified Numeric.DataFrame.BasicTest import qualified Numeric.DataFrame.SubSpaceTest+import qualified Numeric.QuaternionTest ++-- | Collection of tests in detailed-0.9 format tests :: IO [Test] tests = return [ test "DataFrame.Basic" Numeric.DataFrame.BasicTest.runTests , test "DataFrame.SubSpace" Numeric.DataFrame.SubSpaceTest.runTests+ , test "Quaternion" Numeric.QuaternionTest.runTests ] +++-- | Run tests as exitcode-stdio-1.0+main :: IO ()+main = do+ ts <- tests+ trs <- mapM (\(Test ti) ->(,) (name ti) <$> run ti) ts+ case filter (not . isGood) trs of+ [] -> exitSuccess+ xs -> do+ putStrLn $ "Failed tests: " ++ unwords (fmap fst xs)+ exitFailure+ where+ isGood (_, Finished Pass) = True+ isGood _ = False+++-- | Convert QuickCheck props into Cabal tests test :: String -> IO Bool -> Test test tName propOp = Test testI where