packages feed

numeric-prelude 0.1.1 → 0.1.2

raw patch · 13 files changed

+714/−52 lines, 13 files

Files

docs/NOTES view
@@ -83,6 +83,17 @@  Juergen Bokowski <bokowski@mathematik.tu-darmstadt.de>    DMV-Nachrichten 2004/3 + Gerhard Navratil <Navratil@geoinfo.tuwien.ac.at>+   Partial derivatives, Haskell-Cafe 2006-05-08++ Paul Johnson <paul@cogito.org.uk>+   interval arithmetic+   http://sourceforge.net/projects/ranged-sets++ David Amos <polyomino@f2s.com>+   algebra+   http://polyomino.f2s.com/+ * RealFloat Defines the properties of a Floating type, thus should be named 'Floating'.
numeric-prelude.cabal view
@@ -1,5 +1,5 @@ Name:           numeric-prelude-Version:        0.1.1+Version:        0.1.2 License:        GPL License-File:   LICENSE Author:         Dylan Thurston <dpt@math.harvard.edu>, Henning Thielemann <numericprelude@henning-thielemann.de>, Mikael Johansson@@ -185,6 +185,7 @@     MathObj.DiscreteMap     MathObj.LaurentPolynomial     MathObj.Matrix+    MathObj.Monoid     MathObj.PartialFraction     MathObj.Permutation     MathObj.Permutation.CycleList@@ -204,6 +205,7 @@     Number.DimensionTerm.SI     Number.FixedPoint     Number.FixedPoint.Check+    Number.GaloisField2p32m5     Number.NonNegative     Number.NonNegativeChunky     Number.PartiallyTranscendental@@ -246,11 +248,14 @@   GHC-Options:    -Wall   Other-modules:     Test.NumericPrelude.Utility+    Test.Number.GaloisField2p32m5     Test.MathObj.PartialFraction+    Test.MathObj.Matrix     Test.MathObj.Polynomial     Test.MathObj.PowerSeries     Test.MathObj.Gaussian.Variance     Test.MathObj.Gaussian.Bell+    Test.MathObj.Gaussian.Polynomial   Main-Is: Test/Run.hs   If flag(buildTests)     Build-Depends: HUnit >=1 && <2
src/Algebra/Monoid.hs view
@@ -5,16 +5,59 @@ Stability    :   provisional Portability  : -Abstract concept of a Monoid. Will be used in order to generate-type classes for generic algebras. An algebra is a vector space-that also is a monoid.+Abstract concept of a Monoid.+Will be used in order to generate type classes for generic algebras.+An algebra is a vector space that also is a monoid.+Should we use the Monoid class from base library+despite its unfortunate method name @mappend@? -}  module Algebra.Monoid where -{- | We expect a monoid to adher to associativity and the identity-behaving decently. Nothing more, really. -}+import qualified Algebra.Additive as Additive+import qualified Algebra.Ring as Ring +import Data.Monoid as Mn++{- |+We expect a monoid to adher to associativity and+the identity behaving decently.+Nothing more, really.+-} class C a where   idt   :: a   (<*>) :: a -> a -> a++instance C Mn.All where+  idt = mempty+  (<*>) = mappend++instance C Any where+  idt = mempty+  (<*>) = mappend++instance C a => C (Dual a) where+  idt = Mn.Dual idt+  (Mn.Dual x) <*> (Mn.Dual y) = Mn.Dual (y <*> x)++instance C (Endo a) where+  idt = mempty+  (<*>) = mappend++instance C (First a) where+  idt = mempty+  (<*>) = mappend++instance C (Last a) where+  idt = mempty+  (<*>) = mappend+++instance Ring.C a => C (Product a) where+  idt = Mn.Product Ring.one+  (Mn.Product x) <*> (Mn.Product y) = Mn.Product (x Ring.* y)++instance Additive.C a => C (Sum a) where+  idt = Mn.Sum Additive.zero+  (Mn.Sum x) <*> (Mn.Sum y) = Mn.Sum (x Additive.+ y)+
src/Algebra/PrincipalIdealDomain.hs view
@@ -55,6 +55,8 @@ import Control.Monad (foldM, liftM) import Data.List (mapAccumL, mapAccumR, unfoldr) +import Data.Int  (Int,  Int8,  Int16,  Int32,  Int64,  )+ import PreludeBase import Prelude (Integer, Int) import qualified Prelude as P@@ -137,6 +139,11 @@    in  aux  -- could be implemented in a tail-recursive manner+{-+Unfortunately, with the normalization to the stdAssociate,+@gcd 0@ is no longer the identity function,+since @gcd 0 (-2) = 2@.+-} extendedEuclid :: (Units.C a, ZeroTestable.C a) =>    (a -> a -> (a,a)) -> a -> a -> (a,(a,a)) extendedEuclid genDivMod =@@ -288,6 +295,18 @@     gcd = euclid mod  instance C Int where+    gcd = euclid mod++instance C Int8 where+    gcd = euclid mod++instance C Int16 where+    gcd = euclid mod++instance C Int32 where+    gcd = euclid mod++instance C Int64 where     gcd = euclid mod
src/Algebra/Units.hs view
@@ -32,6 +32,8 @@ import Algebra.Additive       (negate) import Algebra.ZeroTestable   (isZero) +import Data.Int  (Int,  Int8,  Int16,  Int32,  Int64,  )+ import PreludeBase import Prelude (Integer, Int) import qualified Prelude as P@@ -96,6 +98,30 @@   stdUnitInv   = intStandardInverse  instance C Integer where+  isUnit       = intQuery+  stdAssociate = intAssociate+  stdUnit      = intStandard+  stdUnitInv   = intStandardInverse++instance C Int8 where+  isUnit       = intQuery+  stdAssociate = intAssociate+  stdUnit      = intStandard+  stdUnitInv   = intStandardInverse++instance C Int16 where+  isUnit       = intQuery+  stdAssociate = intAssociate+  stdUnit      = intStandard+  stdUnitInv   = intStandardInverse++instance C Int32 where+  isUnit       = intQuery+  stdAssociate = intAssociate+  stdUnit      = intStandard+  stdUnitInv   = intStandardInverse++instance C Int64 where   isUnit       = intQuery   stdAssociate = intAssociate   stdUnit      = intStandard
src/MathObj/Matrix.hs view
@@ -2,16 +2,45 @@ {-# LANGUAGE MultiParamTypeClasses #-} {-# LANGUAGE FlexibleInstances #-} {- |-Copyright    :   (c) Mikael Johansson 2006-Maintainer   :   mik@math.uni-jena.de+Copyright    :   (c) Henning Thielemann 2009, Mikael Johansson 2006+Maintainer   :   numericprelude@henning-thielemann.de Stability    :   provisional Portability  :   requires multi-parameter type classes  Routines and abstractions for Matrices and basic linear algebra over fields or rings.++We stick to simple Int indices.+Although advanced indices would be nice+e.g. for matrices with sub-matrices,+this is not easily implemented since arrays+do only support a lower and an upper bound+but no additional parameters.++ToDo:+ - Matrix inverse, determinant -} -module MathObj.Matrix where+module MathObj.Matrix (+   T, Dimension,+   format,+   transpose,+   rows,+   columns,+   fromRows,+   fromColumns,+   fromList,+   dimension,+   numRows,+   numColumns,+   zipWith,+   zero,+   one,+   diagonal,+   scale,+   random,+   randomR,+   ) where  import qualified Algebra.Module   as Module import qualified Algebra.Vector   as Vector@@ -20,32 +49,34 @@  import Algebra.Module((*>), ) import Algebra.Ring((*), fromInteger, scalarProduct, )-import Algebra.Additive((+), (-), zero, subtract, )+import Algebra.Additive((+), (-), subtract, ) -import Data.Array (Array, listArray, elems, bounds, (!), ixmap, range, )+import qualified System.Random as Rnd+import Data.Array (Array, array, listArray, accumArray, elems, bounds, (!), ixmap, range, ) import qualified Data.List as List  import Control.Monad (liftM2, ) import Control.Exception (assert, ) -import Data.Tuple.HT (swap, )+import Data.Tuple.HT (swap, mapFst, ) import Data.List.HT (outerProduct, )-import NumericPrelude (Integer, )++import NumericPrelude (Int, ) import PreludeBase hiding (zipWith, ) + {- |-A matrix is a twodimensional array of ring elements, indexed by integers.+A matrix is a twodimensional array, indexed by integers. -}+data T a =+   Cons (Array (Dimension, Dimension) a)+      deriving (Eq,Ord,Read) -data {-(Ring.C a) =>-}-       T a = Cons (Array (Integer, Integer) a) deriving (Eq,Ord,Read)+type Dimension = Int  {- |-Transposition of matrices is just transposition in the sense of-Data.List.+Transposition of matrices is just transposition in the sense of Data.List. -}-- transpose :: T a -> T a transpose (Cons m) =    let (lower,upper) = bounds m@@ -59,31 +90,63 @@ columns :: T a -> [[a]] columns (Cons m) =    let ((lr,lc), (ur,uc)) = bounds m-   in  outerProduct (curry(m!)) (range (lc,uc)) (range (lr,ur))+   in  outerProduct (flip(curry(m!))) (range (lc,uc)) (range (lr,ur)) -fromList :: Integer -> Integer -> [a] -> T a-fromList m n xs = Cons (listArray ((1,1),(m,n)) xs)+fromRows :: Dimension -> Dimension -> [[a]] -> T a+fromRows m n =+   Cons .+   array (indexBounds m n) .+   concat .+   List.zipWith (\r -> map (\(c,x) -> ((r,c),x))) allIndices .+   map (zip allIndices) -instance (Ring.C a, Show a) => Show (T a) where-  show m = List.unlines $ map (\r -> "(" ++ r ++ ")")-        $ map (unwords . map show) $ rows m+fromColumns :: Dimension -> Dimension -> [[a]] -> T a+fromColumns m n =+   Cons .+   array (indexBounds m n) .+   concat .+   List.zipWith (\r -> map (\(c,x) -> ((c,r),x))) allIndices .+   map (zip allIndices) +fromList :: Dimension -> Dimension -> [a] -> T a+fromList m n xs = Cons (listArray (indexBounds m n) xs) -dimension :: T a -> (Integer,Integer)+appPrec :: Int+appPrec = 10++instance (Show a) => Show (T a) where+   showsPrec p m =+      showParen (p >= appPrec)+         (showString "Matrix.fromRows " . showsPrec appPrec (rows m))++format :: (Show a) => T a -> String+format m = formatS m ""++formatS :: (Show a) => T a -> ShowS+formatS =+   concatS .+   map (\r -> showString "(" . concatS r . showString ")\n") .+   map (List.intersperse (' ':) . map (showsPrec 11)) .+   rows++concatS :: [ShowS] -> ShowS+concatS = flip (foldr ($))++dimension :: T a -> (Dimension,Dimension) dimension (Cons m) = uncurry subtract (bounds m) + (1,1) -numRows :: T a -> Integer+numRows :: T a -> Dimension numRows = fst . dimension -numColumns :: T a -> Integer+numColumns :: T a -> Dimension numColumns = snd . dimension  -- These implementations may benefit from a better exception than -- just assertions to validate dimensionalities instance (Additive.C a) => Additive.C (T a) where-  (+) = zipWith (+)-  (-) = zipWith (-)-  zero = zeroMatrix 1 1+   (+) = zipWith (+)+   (-) = zipWith (-)+   zero = zero 1 1  zipWith :: (a -> b -> c) -> T a -> T b -> T c zipWith op mM@(Cons m) nM@(Cons n) =@@ -93,30 +156,71 @@    in  assert (d == dimension nM) $          uncurry fromList d (List.zipWith op em en) -zeroMatrix :: (Additive.C a) => Integer -> Integer -> T a-zeroMatrix m n = fromList m n $-   List.repeat zero+zero :: (Additive.C a) => Dimension -> Dimension -> T a+zero m n =+   fromList m n $+   List.repeat Additive.zero --    List.replicate (fromInteger (m*n)) zero +one :: (Ring.C a) => Dimension -> T a+one n =+   Cons $+   accumArray (flip const) Additive.zero+      (indexBounds n n)+      (map (\i -> ((i,i), Ring.one)) (indexRange n))++diagonal :: (Additive.C a) => [a] -> T a+diagonal xs =+   let n = List.length xs+   in  Cons $+       accumArray (flip const) Additive.zero+          (indexBounds n n)+          (zip (map (\i -> (i,i)) allIndices) xs)++scale :: (Ring.C a) => a -> T a -> T a+scale s = Vector.functorScale s+ instance (Ring.C a) => Ring.C (T a) where-  mM * nM = assert (numRows mM == numColumns nM) $-        fromList (numColumns mM) (numRows nM)-           (liftM2 scalarProduct (rows mM) (columns nM))-  fromInteger n = fromList 1 1 [fromInteger n]+   mM * nM =+      assert (numColumns mM == numRows nM) $+      fromList (numRows mM) (numColumns nM)+         (liftM2 scalarProduct (rows mM) (columns nM))+   fromInteger n = fromList 1 1 [fromInteger n]  instance Functor T where    fmap f (Cons m) = Cons (fmap f m)  instance Vector.C T where-   zero  = zero+   zero  = Additive.zero    (<+>) = (+)-   (*>)  = Vector.functorScale+   (*>)  = scale  instance Module.C a b => Module.C a (T b) where    x *> m = fmap (x*>) m -{- |-What more do we need from our matrix class? We have addition,++random :: (Rnd.RandomGen g, Rnd.Random a) =>+   Dimension -> Dimension -> g -> (T a, g)+random =+   randomAux Rnd.random++randomR :: (Rnd.RandomGen g, Rnd.Random a) =>+   Dimension -> Dimension -> (a,a) -> g -> (T a, g)+randomR m n rng =+   randomAux (Rnd.randomR rng) m n++{-+could be made nicer with the State monad,+but I like to keep dependencies minimal+-}+randomAux :: (Rnd.RandomGen g, Rnd.Random a) =>+   (g -> (a,g)) -> Dimension -> Dimension -> g -> (T a, g)+randomAux rnd m n g0 =+   mapFst (fromList m n) $ swap $+   List.mapAccumL (\g _i -> swap $ rnd g) g0 (indexRange (m*n))++{-+What more do we need from our matrix type? We have addition, subtraction and multiplication, and thus composition of generic free-module-maps. We're going to want to solve linear equations with or without fields underneath, so we're going to want an implementation@@ -125,6 +229,7 @@ with the Gaussian algorithm or some other goodish method. -} +{- {- |  We'll want generic linear equation solving, returning one solution, any solution really, or nothing. Basically, this is asking for the@@ -144,7 +249,23 @@         (numRows a == numRows y &&     -- they match          numColumns y == 1)               -- and y is a column vector                 Nothing+-}  {- Cf. /usr/lib/hugs/demos/Matrix.hs -}+++-- these functions control whether we use 0 or 1 based indices++indexRange :: Dimension -> [Dimension]+indexRange n = [0..(n-1)]++indexBounds ::+   Dimension -> Dimension ->+   ((Dimension,Dimension), (Dimension,Dimension))+indexBounds m n =+   ((0,0), (m-1,n-1))++allIndices :: [Dimension]+allIndices = [0..]
+ src/MathObj/Monoid.hs view
@@ -0,0 +1,56 @@+{-# LANGUAGE NoImplicitPrelude #-}+module MathObj.Monoid where++import qualified Algebra.PrincipalIdealDomain as PID++import Algebra.PrincipalIdealDomain (gcd, lcm, )+import Algebra.Additive (zero, )+import Algebra.Monoid (C, idt, (<*>), )++import PreludeBase++{- |+It is only a monoid for non-negative numbers.++> idt <*> GCD (-2) = GCD 2++Thus, use this Monoid only for non-negative numbers!+-}+newtype GCD a = GCD {runGCD :: a}+   deriving (Show, Eq)++instance PID.C a => C (GCD a) where+   idt = GCD zero+   (GCD x) <*> (GCD y) = GCD (gcd x y)+++newtype LCM a = LCM {runLCM :: a}+   deriving (Show, Eq)++instance PID.C a => C (LCM a) where+   idt = LCM zero+   (LCM x) <*> (LCM y) = LCM (lcm x y)+++{- |+@Nothing@ is the largest element.+-}+newtype Min a = Min {runMin :: Maybe a}+   deriving (Show, Eq)++instance Ord a => C (Min a) where+   idt = Min Nothing+   (Min x) <*> (Min y) = Min $+      maybe y (\x' -> maybe x (Just . min x') y) x+++{- |+@Nothing@ is the smallest element.+-}+newtype Max a = Max {runMax :: Maybe a}+   deriving (Show, Eq)++instance Ord a => C (Max a) where+   idt = Max Nothing+   (Max x) <*> (Max y) = Max $+      maybe y (\x' -> maybe x (Just . max x') y) x
+ src/Number/GaloisField2p32m5.hs view
@@ -0,0 +1,93 @@+{-# LANGUAGE NoImplicitPrelude #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{- |+This number type is intended for tests of functions over fields,+where the field elements need constant space.+This way we can provide a Storable instance.+For 'Rational' this would not be possible.++However, be aware that sums of non-zero elements may yield zero.+Thus division is not always safe, where it is for rational numbers.+-}+module Number.GaloisField2p32m5 where++import qualified Number.ResidueClass as RC+import qualified Algebra.Module   as Module+import qualified Algebra.Field    as Field+import qualified Algebra.Ring     as Ring+import qualified Algebra.Additive as Additive++import Data.Int (Int64, )+import Data.Word (Word32, Word64, )++import qualified Foreign.Storable.Newtype as SN+import qualified Foreign.Storable as St++import Test.QuickCheck (Arbitrary(..), )++import PreludeBase+import NumericPrelude+++newtype T = Cons {decons :: Word32}+   deriving Eq++{-# INLINE appPrec #-}+appPrec :: Int+appPrec  = 10++instance Show T where+   showsPrec p (Cons x) =+      showsPrec p x+{-+      showParen (p >= appPrec)+         (showString "GF2p32m5.Cons " . shows x)+-}++instance Arbitrary T where+   arbitrary = fmap (Cons . fromInteger . flip mod base) arbitrary+   coarbitrary = undefined++instance St.Storable T where+   sizeOf = SN.sizeOf decons+   alignment = SN.alignment decons+   peek = SN.peek Cons+   poke = SN.poke decons+++base :: Ring.C a => a+base = 2^32-5+++{-# INLINE lift2 #-}+lift2 :: (Word64 -> Word64 -> Word64) -> (T -> T -> T)+lift2 f (Cons x) (Cons y) =+   Cons (fromIntegral (mod (f (fromIntegral x) (fromIntegral y)) base))++{-# INLINE lift2Integer #-}+lift2Integer :: (Int64 -> Int64 -> Int64) -> (T -> T -> T)+lift2Integer f (Cons x) (Cons y) =+   Cons (fromIntegral (mod (f (fromIntegral x) (fromIntegral y)) base))+++instance Additive.C T where+   zero = Cons zero+   (+) = lift2 (+)+--   (-) = lift2 (-)+   x-y = x + negate y+   negate n@(Cons x) =+      if x==0+        then n+        else Cons (base-x)++instance Ring.C T where+   one = Cons one+   (*) = lift2 (*)+   fromInteger =+      Cons . fromInteger . flip mod base++instance Field.C T where+   (/) = lift2Integer (RC.divide base)++instance Module.C T T where+   (*>) = (*)
+ test/Test/MathObj/Gaussian/Polynomial.hs view
@@ -0,0 +1,144 @@+{-# LANGUAGE NoImplicitPrelude #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE FlexibleInstances #-}+module Test.MathObj.Gaussian.Polynomial where++import qualified MathObj.Gaussian.Polynomial as G+import qualified MathObj.Gaussian.Bell as B++import qualified MathObj.Polynomial as Poly++-- import qualified Algebra.Ring           as Ring++import qualified Algebra.Laws as Laws++import qualified Number.Complex as Complex++import Test.NumericPrelude.Utility (testUnit)+import Test.QuickCheck (Testable, quickCheck, (==>))+import qualified Test.HUnit as HUnit++import qualified Number.NonNegative as NonNeg+import Data.Function.HT (nest, )+import Data.Tuple.HT (mapSnd, )++-- import Debug.Trace (trace, )++import PreludeBase as P+import NumericPrelude as NP+++simple ::+   (Testable t) =>+   (G.T Rational -> t) -> IO ()+simple f =+   quickCheck (\x -> f (x :: G.T Rational))++tests :: HUnit.Test+tests =+   HUnit.TestLabel "polynomial" $+   HUnit.TestList $+   map testUnit $+   testList++testList :: [(String, IO ())]+testList =+{-+      ("convolution, dirac",+          simple $ Laws.identity (+) zero) :+-}+      ("convolution, commutative",+          simple $ Laws.commutative G.convolve) :+--          simple $ \x -> Laws.commutative G.convolve (trace (show x) x)) :+      ("convolution, associative",+          simple $ Laws.associative G.convolve) :+      ("multiplication, one",+          simple $ Laws.identity G.multiply G.constant) :+      ("multiplication, commutative",+          simple $ Laws.commutative G.multiply) :+      ("multiplication, associative",+          simple $ Laws.associative G.multiply) :+      ("convolution, multplication, fourier",+          simple $ \x y ->+             G.fourier (G.convolve x y)+              == G.multiply (G.fourier x) (G.fourier y)) :+      ("fourier reverse",+          simple $ \x -> nest 2 G.fourier x == G.reverse x) :+      ("reverse identity",+          simple $ \x -> nest 2 G.reverse x == x) :+      ("fourier eigenfunction differential",+          quickCheck $ \m ->+             m <= 15 ==>+                let n = NonNeg.toNumber m+                    x = G.eigenfunctionDifferential n :: G.T Rational+                    k = Complex.conjugate Complex.imaginaryUnit ^ fromIntegral n+                in  G.fourier x  ==  G.scaleComplex k x) :+      ("fourier eigenfunction iterative",+          quickCheck $ \m ->+             m <= 15 ==>+                let n = NonNeg.toNumber m+                    x = G.eigenfunctionIterative n :: G.T Rational+                    k = Complex.conjugate Complex.imaginaryUnit ^ fromIntegral n+                in  G.fourier x  ==  G.scaleComplex k x) :+{- this does not hold, both functions compute different eigenbases+      ("fourier eigenfunction diff vs. iterative",+          quickCheck $ \n ->+             n <= 15 ==>+                G.eigenfunctionDifferential n ==+                (G.eigenfunctionIterative n :: G.T Rational)) :+-}+      ("translate additive",+          simple $ \x a b ->+             G.translate a (G.translate b x) == G.translate (a+b) x) :+      ("translateComplex additive",+          simple $ \x a b ->+             G.translateComplex a (G.translateComplex b x) == G.translateComplex (a+b) x) :+      ("translateComplex real",+          simple $ \x a ->+             G.translateComplex (Complex.fromReal a) x == G.translate a x) :+      ("modulate additive",+          simple $ \x a b ->+             G.modulate a (G.modulate b x) == G.modulate (a+b) x) :+      ("commute translate modulate",+          simple $ \x a b ->+             G.modulate b (G.translate a x)+              == G.turn (a*b) (G.translate a (G.modulate b x))) :+      ("fourier translate",+          simple $ \x a ->+             G.fourier (G.translate a x)+              == G.modulate a (G.fourier x)) :+      ("dilate multiplicative",+          simple $ \x a b -> a>0 && b>0 ==>+             G.dilate a (G.dilate b x) == G.dilate (a*b) x) :+      ("dilate by shrink",+          simple $ \x a -> a>0 ==>+             G.shrink a x == G.dilate (recip a) x) :+      ("fourier dilate",+          simple $ \x a -> a>0 ==>+             G.fourier (G.dilate a x) == G.amplify a (G.shrink a (G.fourier x))) :+      ("integrate differentiate",+          simple $ \x ->+             G.integrate (G.differentiate x) == (zero, x)) :+      ("fourier differentiate",+          simple $ \x ->+             G.fourier (G.differentiate x) ==+              let y = G.fourier x+              in  y{G.polynomial =+                      Poly.fromCoeffs [0, 0 Complex.+: 2] * G.polynomial y}) :+      ("approximate by bells, translate",+          simple $ \x unit d -> unit/=0 ==>+             G.approximateByBells unit (G.translateComplex d x) ==+             map (mapSnd (B.translateComplex d)) (G.approximateByBells unit x)) :+      ("approximate by bells, dilate",+          simple $ \x unit d -> unit/=0 && d/=0 ==>+             G.approximateByBells unit (G.dilate d x) ==+             map (mapSnd (B.dilate d)) (G.approximateByBells (unit/d) x)) :+      ("approximate by bells, shrink",+          simple $ \x unit d -> unit/=0 && d/=0 ==>+             G.approximateByBells unit (G.shrink d x) ==+             map (mapSnd (B.shrink d)) (G.approximateByBells (unit*d) x)) :+      ("approximate by bells, different implementations",+          quickCheck $ \unit d s p -> unit/=0 ==>+             G.approximateByBellsAtOnce unit d s (p::Poly.T (Complex.T Rational)) ==+             G.approximateByBellsByTranslation unit d s p) :+      []
+ test/Test/MathObj/Matrix.hs view
@@ -0,0 +1,96 @@+{-# LANGUAGE NoImplicitPrelude #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE FlexibleInstances #-}+module Test.MathObj.Matrix where++import qualified MathObj.Matrix as Matrix++import qualified Algebra.Ring           as Ring++import qualified Algebra.Laws as Laws++import qualified Number.NonNegative as NonNeg++import qualified System.Random as Random++import Test.NumericPrelude.Utility (testUnit, )+import Test.QuickCheck (quickCheck, )+import qualified Test.HUnit as HUnit+++import PreludeBase as P+import NumericPrelude as NP+++type Seed = Int+type Dimension = NonNeg.Int++random :: Dimension -> Dimension -> Seed -> Matrix.T Integer+random mn nn seed =+   fst $+   Matrix.random (NonNeg.toNumber mn) (NonNeg.toNumber nn) $+   Random.mkStdGen seed+++tests :: HUnit.Test+tests =+   HUnit.TestLabel "matrix" $+   HUnit.TestList $+   map testUnit $+      ("dimension",+          quickCheck (\m n a ->+             (NonNeg.toNumber m, NonNeg.toNumber n) == Matrix.dimension (random m n a))) :+      ("to and from rows",+          quickCheck (\m n a' ->+             let a = random m n a'+             in  a == Matrix.fromRows (NonNeg.toNumber m) (NonNeg.toNumber n) (Matrix.rows a))) :+      ("to and from columns",+          quickCheck (\m n a' ->+             let a = random m n a'+             in  a == Matrix.fromColumns (NonNeg.toNumber m) (NonNeg.toNumber n) (Matrix.columns a))) :+      ("transpose, rows, columns",+          quickCheck (\m n a' ->+             let a = random m n a'+             in  Matrix.rows a == Matrix.columns (Matrix.transpose a))) :+      ("transpose, columns, rows",+          quickCheck (\m n a' ->+             let a = random m n a'+             in  Matrix.columns a == Matrix.rows (Matrix.transpose a))) :+      ("addition, zero",+          quickCheck (\m n a ->+             Laws.identity (+) (Matrix.zero (NonNeg.toNumber m) (NonNeg.toNumber n)) (random m n a))) :+      ("addition, commutative",+          quickCheck (\m n a b ->+             Laws.commutative (+) (random m n a) (random m n b))) :+      ("addition, associative",+          quickCheck (\m n a b c ->+             Laws.associative (+) (random m n a) (random m n b) (random m n c))) :+      ("addition, transpose",+          quickCheck (\m n a b ->+             Laws.homomorphism Matrix.transpose (+) (+) (random m n a) (random m n b))) :+      ("one, diagonal",+          quickCheck (\n' ->+             let n = NonNeg.toNumber n'+             in Matrix.one n == (Matrix.diagonal $ replicate n Ring.one :: Matrix.T Integer))) :+      ("multiplication, one left",+          quickCheck (\m n a ->+             Laws.leftIdentity  (*) (Matrix.one (NonNeg.toNumber m)) (random m n a))) :+      ("multiplication, one right",+          quickCheck (\m n a ->+             Laws.rightIdentity (*) (Matrix.one (NonNeg.toNumber n)) (random m n a))) :+      ("multiplication, associative",+          quickCheck (\k l m n a b c ->+             Laws.associative (*) (random k l a) (random l m b) (random m n c))) :+      ("multiplication and addition, distributive left",+          quickCheck (\l m n a b c ->+             Laws.leftDistributive (*) (+) (random n l a) (random m n b) (random m n c))) :+      ("multiplication and addition, distributive right",+          quickCheck (\l m n a b c ->+             Laws.rightDistributive (*) (+) (random l m a) (random m n b) (random m n c))) :+      ("multiplication, transpose",+          quickCheck (\l m n a b ->+             Laws.homomorphism Matrix.transpose (*) (flip (*)) (random l m a) (random m n b))) :+{-+      ("division",       quickCheck (\x -> Integral.propInverse (x :: Poly.T Rational))) :+-}+      []
test/Test/MathObj/Polynomial.hs view
@@ -14,7 +14,7 @@ import qualified Data.List as List  import Test.NumericPrelude.Utility (testUnit)-import Test.QuickCheck (Property, quickCheck, (==>))+import Test.QuickCheck (Property, quickCheck, (==>), Testable, ) import qualified Test.HUnit as HUnit  @@ -32,6 +32,13 @@ mul xs ys  =  Poly.equal (Poly.mul xs ys) (Poly.mulShear xs ys)  +test :: Testable a => (Poly.T Integer -> a) -> IO ()+test = quickCheck++testRat :: Testable a => (Poly.T Rational -> a) -> IO ()+testRat = quickCheck++ tests :: HUnit.Test tests =    HUnit.TestLabel "polynomial" $@@ -39,12 +46,12 @@    map testUnit $       ("tensor product", quickCheck (tensorProductTranspose :: [Integer] -> [Integer] -> Property)) :       ("mul speed",      quickCheck (mul                    :: [Integer] -> [Integer] -> Bool)) :-      ("addition, zero",         quickCheck (\x -> Laws.identity (+) zero (x :: Poly.T Integer))) :-      ("addition, commutative",  quickCheck (\x -> Laws.commutative (+) (x :: Poly.T Integer))) :-      ("addition, associative",  quickCheck (\x -> Laws.associative (+) (x :: Poly.T Integer))) :-      ("multiplication, one",          quickCheck (\x -> Laws.identity (*) one (x :: Poly.T Integer))) :-      ("multiplication, commutative",  quickCheck (\x -> Laws.commutative (*) (x :: Poly.T Integer))) :-      ("multiplication, associative",  quickCheck (\x -> Laws.associative (*) (x :: Poly.T Integer))) :-      ("multiplication and addition, distributive",   quickCheck (\x -> Laws.leftDistributive (*) (+) (x :: Poly.T Integer))) :-      ("division",       quickCheck (\x -> Integral.propInverse (x :: Poly.T Rational))) :+      ("addition, zero",         test (Laws.identity (+) zero)) :+      ("addition, commutative",  test (Laws.commutative (+))) :+      ("addition, associative",  test (Laws.associative (+))) :+      ("multiplication, one",          test (Laws.identity (*) one)) :+      ("multiplication, commutative",  test (Laws.commutative (*))) :+      ("multiplication, associative",  test (Laws.associative (*))) :+      ("multiplication and addition, distributive",   test (Laws.leftDistributive (*) (+))) :+      ("division",       testRat (Integral.propInverse)) :       []
+ test/Test/Number/GaloisField2p32m5.hs view
@@ -0,0 +1,37 @@+{-# LANGUAGE NoImplicitPrelude #-}+module Test.Number.GaloisField2p32m5 where++import qualified Number.GaloisField2p32m5 as GF++import qualified Algebra.Laws as Laws++import Test.NumericPrelude.Utility (testUnit)+import Test.QuickCheck (Testable, quickCheck, (==>))+import qualified Test.HUnit as HUnit+++import PreludeBase as P+import NumericPrelude as NP+++test :: Testable a => (GF.T -> a) -> IO ()+test = quickCheck+++tests :: HUnit.Test+tests =+   HUnit.TestLabel "galois field 2^32-5" $+   HUnit.TestList $+   map testUnit $+      ("addition, zero",         test (Laws.identity (+) zero)) :+      ("addition, commutative",  test (Laws.commutative (+))) :+      ("addition, associative",  test (Laws.associative (+))) :+      ("addition, negate",       test (Laws.inverse (+) negate zero)) :+      ("addition, subtract",     test (\x -> Laws.inverse (+) (x-) x)) :+      ("multiplication, one",          test (Laws.identity (*) one)) :+      ("multiplication, commutative",  test (Laws.commutative (*))) :+      ("multiplication, associative",  test (Laws.associative (*))) :+      ("multiplication, recip",        test (\y -> y /= 0 ==> Laws.inverse (*) recip one y)) :+      ("multiplication, division",     test (\y x -> y /= 0 ==> Laws.inverse (*) (x/) x y)) :+      ("multiplication and addition, distributive",   test (Laws.leftDistributive (*) (+))) :+      []
test/Test/Run.hs view
@@ -4,8 +4,10 @@ import qualified Test.MathObj.Gaussian.Variance as GaussVariance import qualified Test.MathObj.Gaussian.Bell as GaussBell import qualified Test.MathObj.PartialFraction as PartialFraction+import qualified Test.MathObj.Matrix  as Matrix import qualified Test.MathObj.Polynomial  as Polynomial import qualified Test.MathObj.PowerSeries as PowerSeries+import qualified Test.Number.GaloisField2p32m5 as GF import qualified Test.HUnit.Text as HUnitText import qualified Test.HUnit as HUnit @@ -16,7 +18,9 @@          GaussBell.tests :          GaussPoly.tests :          PartialFraction.tests :+         Matrix.tests :          Polynomial.tests :          PowerSeries.tests :+         GF.tests :          [])       return ()