cyclotomic 0.2 → 0.3

raw patch · 5 files changed

+142/−27 lines, 5 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

• Examples.hs +17/−0
• Properties.hs +49/−0
• UnitTests.hs +28/−0
• cyclotomic.cabal +7/−5
• src/Data/Complex/Cyclotomic.hs +41/−22

@@ -0,0 +1,17 @@+{-# OPTIONS_GHC -Wall #-}++module Examples where++import Data.Complex.Cyclotomic+import Data.Ratio++-- | A Discrete Fourier Transform+dft :: [Cyclotomic] -> [Cyclotomic]+dft cs = [sum $ zipWith (*) [conj (e m^(k*n)) | n <- [0..]] cs | k <- [0..m-1]]+          where m = fromIntegral $ length cs++-- | Inverse Discrete Fourier Transform+dftInv :: [Cyclotomic] -> [Cyclotomic]+dftInv cs = [minv * sum (zipWith (*) [e m^(k*n) | n <- [0..]] cs) | k <- [0..m-1]]+          where m = fromIntegral $ length cs+                minv = fromRational (1 % m)

@@ -0,0 +1,49 @@+{-# OPTIONS_GHC -Wall #-}++-- Properties.hs is modified from Brent Yorgey's file of the same name+-- for his package Data.List.Split++module Properties where++import Data.Complex.Cyclotomic+import Test.QuickCheck+import Test.SmallCheck++import Control.Monad+import Text.Printf++import Examples++main :: IO ()+main = do+    results <- mapM (\(s,t) -> printf "%-40s" s >> t) tests+    when (not . all isSuccess $ results) $ fail "Not all tests passed!"+ where+    isSuccess (Success{}) = True+    isSuccess _ = False+    qc si su x = quickCheckWithResult (stdArgs { maxSize = si, maxSuccess = su }) x+    tests = [ ("square of sqrtRat",             qc  15  15 prop_square_sqrtRat)+            , ("dftInv . dft",                  qc  30  15 prop_dftInv_dft)+            , ("dft . dftInv",                  qc  30  15 prop_dft_dftInv)+            ]++two :: Integer+two = 2++prop_square_sqrtRat :: Int -> Bool+prop_square_sqrtRat n = sqrtRat (fromIntegral n) ^ two == fromIntegral n++prop_dftInv_dft :: [Rational] -> Bool+prop_dftInv_dft rs = dftInv (dft cs) == cs+    where cs = map fromRational rs++prop_dft_dftInv :: [Rational] -> Bool+prop_dft_dftInv rs = dft (dftInv cs) == cs+    where cs = map fromRational rs++smallCheckTests :: IO ()+smallCheckTests+    = smallCheck 10 prop_square_sqrtRat+--      smallCheck 10 prop_dftInv_dft     >>+--      smallCheck 10 prop_dft_dftInv+

@@ -0,0 +1,28 @@+{-# OPTIONS_GHC -Wall #-}++module UnitTests where++import Data.Complex.Cyclotomic+import Test.HUnit+--import Examples+import Data.Ratio++rationals :: [Rational]+rationals = 0 % 1 : [sign * k % j | n <- [0..], m <- [0..n-1], sign <- [1,-1]+                    , let k = m + 1, let j = n - m, gcd k j == 1]++test1 :: Test+test1 = TestLabel "polarRat" $ TestList [polarRat 1 (p % q) ~=? e q^p | p <- [0..10], q <- [1..10]]++test2 :: Test+test2 = TestList [sqrtRat r ^ (2::Int) ~=? fromRational r | r <- take 100 rationals]++test3 :: Test+test3 = TestList [sqrtRat (r*r) ~=? fromRational (abs r) | r <- take 100 rationals]++test4 :: Test+test4 = TestList [z * (1 / z) ~=? 1 | n <- [1..10], m <- [1..10], let z = e n + e m, z /= 0]++main :: IO Counts+main = runTestTT $ TestList [test1,test2,test3,test4]+

@@ -1,15 +1,17 @@ Name:                cyclotomic-Version:             0.2-Synopsis:            A subfield of the complex numbers for exact calculation+Version:             0.3+Synopsis:            A subfield of the complex numbers for exact calculation. Description:         The cyclotomic numbers are a subset of the-                     complex numbers with a number of nice properties.-                     They are represented exactly, enabling exact+                     complex numbers that are represented exactly, enabling exact                      computations and equality comparisons.  They                      contain the Gaussian rationals (complex numbers-                     of the form p + q i with p and q rational).  The+                     of the form p + q i with p and q rational), as well+                     as all complex roots of unity.  The                      cyclotomic numbers contain the square roots of                      all rational numbers.  They contain the sine and                      cosine of all rational multiples of pi.+                     The cyclotomic numbers form a field, being closed under+                     addition, subtraction, mutiplication, and division. License:             GPL-3 License-file:        LICENSE Author:              Scott N. Walck

@@ -69,7 +69,6 @@     ,conj     ,real     ,imag-    ,modSq     ,isReal     ,isRat     ,isGaussianRat@@ -90,16 +89,17 @@                              , coeffs :: M.Map Integer Rational                              } deriving (Eq) --- | @signum c@ is the complex number with magnitude 1 that has the same argument as c;+-- | @abs@ and @signum@ are partial functions.+--   A cyclotomic number is not guaranteed to have a cyclotomic absolute value.+--   When defined, @signum c@ is the complex number with magnitude 1 that has the same argument as c; --   @signum c = c / abs c@. instance Num Cyclotomic where     (+) = sumCyc     (*) = prodCyc     (-) c1 c2 = sumCyc c1 (aInvCyc c2)     negate = aInvCyc-    abs = sqrtRat . modSq-    signum 0 = zeroCyc-    signum c = c / abs c+    abs = absVal+    signum = sigNum     fromInteger 0 = zeroCyc     fromInteger n = Cyclotomic 1 (M.singleton 0 (fromIntegral n)) @@ -132,7 +132,7 @@ leadingTerm r _ 0 = showRat r leadingTerm r n ex     | r == 1     = t-    | r == (-1)  = "-" ++ t+    | r == (-1)  = '-' : t     | r > 0      = showRat r ++ "*" ++ t     | r < 0      = "-" ++ showRat (abs r) ++ "*" ++ t     | otherwise  = ""@@ -226,12 +226,18 @@ imag :: Cyclotomic -> Cyclotomic imag z = (z - conj z) / (2*i) --- | Modulus squared.-modSq :: Cyclotomic -> Rational-modSq z = case toRat (z * conj z) of-            Just msq -> msq-            Nothing  -> error $ "modSq:  tried z = " ++ show z+absVal :: Cyclotomic -> Cyclotomic+absVal c = let modsq = c * conj c+           in case toRat modsq of+                Just msq -> sqrtRat msq+                Nothing  -> error "abs not available for this number" +sigNum :: Cyclotomic -> Cyclotomic+sigNum 0 = zeroCyc+sigNum c = let modsq = c * conj c+           in if isRat modsq then c / absVal c+              else error "signum not available for this number"+ convertToBase :: Integer -> M.Map Integer Rational -> M.Map Integer Rational convertToBase n mp = foldr (\(p,r) m -> replace n p r m) mp (extraneousPowers n) @@ -250,7 +256,7 @@ gcdReduce :: Cyclotomic -> Cyclotomic gcdReduce cyc@(Cyclotomic n mp) = case gcdCyc cyc of                                     1 -> cyc-                                    d -> Cyclotomic (n `div` d) (M.mapKeys (\k -> k `div` d) mp)+                                    d -> Cyclotomic (n `div` d) (M.mapKeys (`div` d) mp)  gcdCyc :: Cyclotomic -> Integer gcdCyc (Cyclotomic n mp) = gcdList (n:M.keys mp)@@ -274,16 +280,14 @@     where       factors = factorise n       phi = foldr (\p n' -> n' `div` p * (p-1)) n [p | (p,_) <- factors]-      nrp = length (factors)+      nrp = length factors       sqfree = all (<=1) [m | (_,m) <- factors]  equalCoefficients :: Cyclotomic -> Maybe Rational equalCoefficients (Cyclotomic _ mp)     = case ts of         []    -> Nothing-        (x:_) -> case equal ts of-                   True  -> Just x-                   False -> Nothing+        (x:_) -> if equal ts then Just x else Nothing       where         ts = M.elems mp @@ -299,7 +303,7 @@  reduceByPrime :: Integer -> Cyclotomic -> Cyclotomic reduceByPrime p c@(Cyclotomic n _)-    = case sequence $ map (\r -> equalReplacements p r c) [0,p..n-p] of+    = case mapM (\r -> equalReplacements p r c) [0,p..n-p] of         Just cfs -> Cyclotomic (n `div` p) $ removeZeros $ M.fromList $ zip [0..(n `div` p)-1] (map negate cfs)         Nothing  -> c @@ -323,8 +327,8 @@ includeMods n q start = [start] ++ takeWhile (>= 0) [start-q,start-2*q..] ++ takeWhile (< n) [start+q,start+2*q..]  removeExps :: Integer -> Integer -> Integer -> [Integer]-removeExps n 2 q = concat $ map (includeMods n q) $ map ((n `div` q) *) [q `div` 2..q-1]-removeExps n p q = concat $ map (includeMods n q) $ map ((n `div` q) *) [-m..m]+removeExps n 2 q = concatMap (includeMods n q) $ map ((n `div` q) *) [q `div` 2..q-1]+removeExps n p q = concatMap (includeMods n q) $ map ((n `div` q) *) [-m..m]     where m = (q `div` p - 1) `div` 2  pqPairs :: Integer -> [(Integer,Integer)]@@ -333,7 +337,7 @@ extraneousPowers :: Integer -> [(Integer,Integer)] extraneousPowers n     | n < 1      = error "extraneousPowers needs a postive integer"-    | otherwise  = nub $ concat $ [[(p,r) | r <- removeExps n p q] | (p,q) <- pqPairs n]+    | otherwise  = nub $ concat [[(p,r) | r <- removeExps n p q] | (p,q) <- pqPairs n]  -- | Sum of two cyclotomic numbers. sumCyc :: Cyclotomic -> Cyclotomic -> Cyclotomic@@ -361,15 +365,30 @@  -- | Additive identity. zeroCyc :: Cyclotomic-zeroCyc = Cyclotomic 1 (M.empty)+zeroCyc = Cyclotomic 1 M.empty  -- | Additive inverse. aInvCyc :: Cyclotomic -> Cyclotomic aInvCyc = prodRatCyc (-1) +multiplyExponents :: Integer -> Cyclotomic -> Cyclotomic+multiplyExponents j (Cyclotomic n mp)+    | gcd j n /= 1  = error "multiplyExponents needs gcd == 1"+    | otherwise     = mkCyclotomic n (M.mapKeys (\k -> j*k `mod` n) mp)++productOfGaloisConjugates :: Cyclotomic -> Cyclotomic+productOfGaloisConjugates c+    = product [multiplyExponents j c | j <- [2..ord], gcd j ord == 1]+      where+        ord = order c+ -- | Multiplicative inverse. invCyc :: Cyclotomic -> Cyclotomic-invCyc z = prodRatCyc (1 / modSq z) (conj z)+invCyc z = case toRat (z * prod) of+             Just r  -> prodRatCyc (1 / r) prod+             Nothing -> error "invCyc:  product of Galois conjugates not rational; this is a bug, please inform package maintainer"+    where+      prod = productOfGaloisConjugates z  -- | Is the cyclotomic a real number? isReal :: Cyclotomic -> Bool