padic (empty) → 0.1.0.0
raw patch · 17 files changed
+1936/−0 lines, 17 filesdep +QuickCheckdep +basedep +constraintssetup-changed
Dependencies added: QuickCheck, base, constraints, criterion, integer-gmp, mod, padic, tasty, tasty-expected-failure, tasty-hunit, tasty-quickcheck
Files
- ChangeLog.md +3/−0
- LICENSE +21/−0
- README.md +6/−0
- Setup.hs +2/−0
- bench/Bench.hs +115/−0
- padic.cabal +89/−0
- src/Math/NumberTheory/Padic.hs +162/−0
- src/Math/NumberTheory/Padic/Analysis.hs +400/−0
- src/Math/NumberTheory/Padic/Integer.hs +133/−0
- src/Math/NumberTheory/Padic/Rational.hs +171/−0
- src/Math/NumberTheory/Padic/Types.hs +280/−0
- test/Spec.hs +20/−0
- test/Test/Analysis.hs +50/−0
- test/Test/Base.hs +140/−0
- test/Test/Commons.hs +100/−0
- test/Test/Integer.hs +115/−0
- test/Test/Rational.hs +129/−0
+ ChangeLog.md view
@@ -0,0 +1,3 @@+# Changelog for padic++## Unreleased changes
+ LICENSE view
@@ -0,0 +1,21 @@+MIT License++Copyright (c) 2022 Sergey B. Samoylenko++Permission is hereby granted, free of charge, to any person obtaining a copy+of this software and associated documentation files (the "Software"), to deal+in the Software without restriction, including without limitation the rights+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell+copies of the Software, and to permit persons to whom the Software is+furnished to do so, subject to the following conditions:++The above copyright notice and this permission notice shall be included in all+copies or substantial portions of the Software.++THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE+SOFTWARE.
+ README.md view
@@ -0,0 +1,6 @@+# Math.NumberTheory.Padic++Module introduces p-adic integers and p-adic rational numbers of fixed and arbitratry precision, implementing basic arithmetic as well as some specific functions, i.e. detection of periodicity in digital sequence, rational reconstruction, square roots etc.++In order to gain efficiency the integer p-adic number with radix `p` is internally+represented in form `N mod p^k` as only one digit `N`, lifted to modulo `p^k`, where `k` is chosen so that within working precision numbers belogning to `Int` and `Ratio Int` types could be reconstructed by extended Euclidean algorithm. Canonical expansion is used for textual output only.
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ bench/Bench.hs view
@@ -0,0 +1,115 @@+{-# language TypeApplications #-}+{-# language DataKinds #-}+{-# language TypeFamilies #-}+{-# language TypeOperators #-}+{-# language FlexibleContexts #-}++import Criterion.Main+import qualified Math.NumberTheory.Padic.Fixed as F+import Math.NumberTheory.Padic+import Data.Maybe+import Data.Mod+import Data.List (transpose)++addBench :: (Show n, Num n) => n -> Int -> Int+addBench w n = length $ show (w + sum (take n (fib 0 1)))++mulBench :: (Show n, Num n) => n -> Integer -> Int+mulBench w n = let M ((x:_):_) = fibM n in length (show (w + x))++divBench :: (Show a, Num a) => a -> a -> Int -> String+divBench r w n = show (product $ take n $ logistic r w)++fib a b = a : fib b (a + b)+ +fibM 0 = I+fibM n = M [[1,1],[1,0]] <> fibM (n-1)++logistic r = iterate (\x -> r*x*(1-x))+ ++divMaybe a b | isInvertible b = Just (a `div` b)+ | otherwise = Nothing++fracMaybe a b | isInvertible b = Just (a / b)+ | otherwise = Nothing++data M a = I | M ![[a]]++dot a b = sum $! zipWith (*) a b ++instance Num a => Semigroup (M a) where+ I <> x = x+ x <> I = x+ M a <> M b = M $ [ [ dot x y | y <- transpose b ] | x <- a ]++instance Num a => Monoid (M a) where+ mempty = I++addN = 400+mulN = 400+divN = 100+ +suite :: [Benchmark]+suite =+ [ bgroup+ "add"+ [ bench "Integer" $ whnf (addBench (0 :: Integer)) addN+ , bench "Mod 2^20" $ whnf (addBench (0 :: Mod 2199023255552)) addN+ , bench "Z 2" $ whnf (addBench (0 :: Z' 2 20)) addN+ , bench "Z 13" $ whnf (addBench (0 :: Z' 13 20)) addN+ , bench "Z 251" $ whnf (addBench (0 :: Z' 251 20)) addN+ , bench "F.Z 2" $ whnf (addBench (0 :: F.Z' 2 20)) addN+ , bench "F.Z 2 100" $ whnf (addBench (0 :: F.Z' 2 1000)) addN+ , bench "F.Z 13" $ whnf (addBench (0 :: F.Z' 13 20)) addN+ , bench "F.Z 251" $ whnf (addBench (0 :: F.Z' 251 20)) addN+ , bench "Q 2" $ whnf (addBench (0 :: Q' 2 20)) addN+ , bench "Q 2 100" $ whnf (addBench (0 :: Q' 2 100)) addN+ , bench "Q 13" $ whnf (addBench (0 :: Q' 13 20)) addN+ , bench "Q 251" $ whnf (addBench (0 :: Q' 251 20)) addN+ , bench "F.Q 2" $ whnf (addBench (0 :: F.Q' 2 20)) addN+ , bench "F.Q 2 100" $ whnf (addBench (0 :: F.Q' 2 100)) addN+ , bench "F.Q 13" $ whnf (addBench (0 :: F.Q' 13 20)) addN+ , bench "F.Q 251" $ whnf (addBench (0 :: F.Q' 251 20)) addN+ ]+ , bgroup+ "mul"+ [ bench "Integer" $ whnf (mulBench (0 :: Integer)) mulN+ , bench "Mod 2^20" $ whnf (mulBench (0 :: Mod 2199023255552)) mulN+ , bench "Z 2" $ whnf (mulBench (0 :: Z' 2 20)) mulN+ , bench "Z 2 100" $ whnf (mulBench (0 :: Z' 2 100)) mulN+ , bench "Z 13" $ whnf (mulBench (0 :: Z' 13 20)) mulN+ , bench "Z 251" $ whnf (mulBench (0 :: Z' 251 20)) mulN+ , bench "F.Z 2" $ whnf (mulBench (0 :: F.Z' 2 20)) mulN+ , bench "F.Z 2 100" $ whnf (mulBench (0 :: F.Z' 2 100)) mulN+ , bench "F.Z 13" $ whnf (mulBench (0 :: F.Z' 13 20)) mulN+ , bench "F.Z 251" $ whnf (mulBench (0 :: F.Z' 251 20)) mulN+ , bench "Q 2" $ whnf (mulBench (0 :: Q' 2 20)) mulN+ , bench "Q 2 100" $ whnf (mulBench (0 :: Q' 2 100)) mulN+ , bench "Q 13" $ whnf (mulBench (0 :: Q' 13 20)) mulN+ , bench "Q 251" $ whnf (mulBench (0 :: Q' 251 20)) mulN+ , bench "F.Q 2" $ whnf (mulBench (0 :: F.Q' 2 20)) mulN+ , bench "F.Q 2 100" $ whnf (mulBench (0 :: F.Q' 2 100)) mulN+ , bench "F.Q 13" $ whnf (mulBench (0 :: F.Q' 13 20)) mulN+ , bench "F.Q 251" $ whnf (mulBench (0 :: F.Q' 251 20)) mulN+ ]+ , bgroup+ "div"+ [ bench "Double" $ whnf (divBench (13 / 5) (4 / 3 :: Double)) divN+ , bench "Z 2" $ whnf (divBench (13 `div` 5) (5 `div` 3 :: Z 2)) divN+ , bench "Z 17" $ whnf (divBench (13 `div` 4) (4 `div` 3 :: Z 17)) divN+ , bench "Z 251" $ whnf (divBench (13 `div` 4) (4 `div` 3 :: Z 251)) divN+ , bench "F.Z 2" $ whnf (divBench (13 `div` 5) (4 `div` 3 :: F.Z 2)) divN+ , bench "F.Z 13" $ whnf (divBench (13 `div` 4) (4 `div` 3 :: F.Z 13)) divN+ , bench "F.Z 251" $ whnf (divBench (13 `div` 4) (4 `div` 3 :: F.Z 251)) divN+ , bench "Q 2" $ whnf (divBench (13 / 4) (4 / 3 :: Q 2)) divN+ , bench "Q 13" $ whnf (divBench (13 / 4) (4 / 3 :: Q 13)) divN+ , bench "Q 251" $ whnf (divBench (13 / 4) (4 / 3 :: Q 251)) divN+ , bench "F.Q 2" $ whnf (divBench (13 / 4) (4 / 3 :: F.Q 2)) divN+ , bench "F.Q 13" $ whnf (divBench (13 / 4) (4 / 3 :: F.Q 13)) divN+ , bench "F.Q 251" $ whnf (divBench (13 / 4) (4 / 3 :: F.Q 251)) divN+ ]+ ]++main :: IO ()+main = defaultMain suite
+ padic.cabal view
@@ -0,0 +1,89 @@+cabal-version: 1.12++-- This file has been generated from package.yaml by hpack version 0.34.5.+--+-- see: https://github.com/sol/hpack++name: padic+version: 0.1.0.0+synopsis: Fast, type-safe p-adic arithmetic+description: Implementation of p-adic arithmetics on the base of fast modular arithmetics. Module introduces data types for p-adic integers and rationals with arbitrary precision as well as some specific functions (rational reconstruction, p-adic signum function, square roots etc.).+category: Math, Number Theory+homepage: https://github.com/samsergey/padic-0.1.0.0+bug-reports: https://github.com/samsergey/padic-0.1.0.0/issues+author: Sergey B. Samoylenko <samsergey@yandex.ru>+maintainer: samsergey@yandex.ru+copyright: 2022 Sergey B. Samoylenko+license: MIT+license-file: LICENSE+build-type: Simple+tested-with:+ GHC ==8.10.7+extra-source-files:+ README.md+ ChangeLog.md++source-repository head+ type: git+ location: https://github.com/samsergey/padic++library+ exposed-modules:+ Math.NumberTheory.Padic+ Math.NumberTheory.Padic.Analysis+ Math.NumberTheory.Padic.Integer+ Math.NumberTheory.Padic.Rational+ Math.NumberTheory.Padic.Types+ other-modules:+ Paths_padic+ hs-source-dirs:+ src+ build-depends:+ base >=4.14 && <4.17+ , constraints ==0.13.*+ , integer-gmp >=1.0.3 && <1.1+ , mod >=0.1.2.2 && <1.3+ default-language: Haskell2010++test-suite padic-test+ type: exitcode-stdio-1.0+ main-is: Spec.hs+ other-modules:+ Test.Analysis+ Test.Base+ Test.Commons+ Test.Integer+ Test.Rational+ Paths_padic+ hs-source-dirs:+ test+ ghc-options: -threaded -rtsopts -with-rtsopts=-N+ build-depends:+ QuickCheck >=2.14.2 && <2.15+ , base >=4.14+ , constraints ==0.13.*+ , integer-gmp >=1.0.3 && <1.1+ , mod >=0.1.2.2 && <1.3+ , padic+ , tasty >=1.4.2 && <1.5+ , tasty-expected-failure >=0.12 && <0.14+ , tasty-hunit >=0.10 && <0.12+ , tasty-quickcheck >=0.10 && <0.12+ default-language: Haskell2010++benchmark criterion-benchmarks+ type: exitcode-stdio-1.0+ main-is: Bench.hs+ other-modules:+ Paths_padic+ hs-source-dirs:+ bench+ ghc-options: -threaded -rtsopts -with-rtsopts=-N+ build-depends:+ base >=4.14 && <4.17+ , constraints ==0.13.*+ , criterion >=1.5.12 && <1.5.14+ , integer-gmp >=1.0.3 && <1.1+ , mod >=0.1.2.2 && <1.3+ , padic+ default-language: Haskell2010
+ src/Math/NumberTheory/Padic.hs view
@@ -0,0 +1,162 @@++{- |+Module : Math.NumberTheory.Padic.Fixed+Description : Representation and simple algebra for p-adic numbers.+Copyright : (c) Sergey Samoylenko, 2022+License : GPL-3+Maintainer : samsergey@yandex.ru+Stability : experimental+Portability : POSIX++Module introduces p-adic integers and rationals with basic p-adic arithmetics+and implments some specific functions (rational reconstruction, p-adic signum function, square roots etc.).++A truncated p-adic number \(x\) can be represented in three ways:++\[+\begin{align}+x &= p^v u & (1)\\+& = d_0 + d_1 p + d_2 p^2 + ... d_k p^k & (2)\\+&= N\ \mathrm{mod}\ p^k, & (3)+\end{align}+\]+where \(p > 1, k > 0, v \in \mathbb{Z},u \in \mathbb{Z_p},d_i \in \mathbb{Z}/p\mathbb{Z}, N \in \mathbb{Z}/p^k \mathbb{Z}\)++In order to gain efficiency the integer p-adic number with radix \(p\) is internally+represented in form \((3)\) as only one digit \(N\), lifted to modulo \(p^k\), where \(k\) is+chosen so that within working precision numbers belogning to @Int@ and @Ratio Int@ types could be+reconstructed by extended Euclidean method. Form \((2)\) is used for textual output only, and form \((1)\)+is used for transformations to and from rationals.++The documentation and the module bindings use following terminology:++ * `radix` -- modulus \(p\) of p-adic number,+ * `precision` -- maximal power \(k\) in p-adic expansion,+ * `unit` -- invertible muliplier \(u\) for prime \(p\),+ * `valuation` -- exponent \(v\),+ * `digits` -- list \(d_0,d_1,d_2,... d_k\) in the canonical p-adic expansion of a number,+ * `lifted` -- digit \(N\) lifted to modulo \(p^k\).++Rational p-adic number is represented as a unit (belonging to \(\mathbb{Z_p}\) ) and valuation, which may be negative.++The radix \(p\) of a p-adic number is specified at a type level via type-literals. In order to use them GHCi should be loaded with `-XDataKinds` extensions.++>>> :set -XDataKinds+>>> 45 :: Z 10+45+>>> 45 :: Q 5+140.0++Negative p-adic integers and rational p-adics have trailing periodic digit sequences, which are represented in parentheses.++>>> -45 :: Z 7+(6)04+>>> 1/7 :: Q 10+(285714)3.0++By default the precision of p-adics is computed so that it is possible to reconstruct integers and rationals using extended Euler's method. However precision could be specified explicitly via type-literal:++>>> sqrt 2 :: Q 7+…623164112011266421216213.0+>>> sqrt 2 :: Q' 7 5+…16213.0+>>> sqrt 2 :: Q' 7 50+…16244246442640361054365536623164112011266421216213.0+++Between types defined in the module there are bijective mappings as shown in the diagram:++@+ [Mod p]+ / \\+ digits / \\ digits+ fromDigits / \\ fromDigits+ / \\+ toInteger / fromUnit \\ fromRational+Z \<-----------\> Z p \<----------\> Q p \<------------\> Q+ fromInteger \\ unit / toRational+ \\ /+ lifted \\ / lifted+ mkUnit \\ / mkUnit+ \\ /+ Integer+@+++-}+------------------------------------------------------------+module Math.NumberTheory.Padic+( + -- * Data types+ -- ** p-Adic integers+ Z+ , Z'+ -- ** p-Adic rationals+ , Q+ , Q'+ , Padic+ , SufficientPrecision+ -- * Classes and functions+ -- ** Type synonyms and constraints+ , ValidRadix+ , KnownRadix+ , LiftedRadix+ , Radix+ -- ** p-adic numbers+ , PadicNum+ , Unit+ , Digit+ -- * Functions and utilities+ -- ** p-adic numbers and arithmetics+ , radix+ , precision+ , digits+ , firstDigit+ , reduce+ , fromDigits+ , lifted+ , mkUnit+ , splitUnit+ , fromUnit+ , unit+ , valuation+ , norm+ , normalize+ , inverse+ , isInvertible+ , isZero+ , getUnitZ+ , getUnitQ+ -- ** p-adic analysis+ , findSolutionMod+ , henselLifting+ , unityRoots+ , pSqrt+ , pPow+ , zPow+ , pExp+ , pLog+ , pSin+ , pCos+ , pSinh+ , pCosh+ , pTanh+ , pAsin+ , pAsinh+ , pAcosh+ , pAtanh+ -- ** Miscelleneos tools+ , fromRadix+ , toRadix+ , findCycle+ , sufficientPrecision+ ) where++import Math.NumberTheory.Padic.Types+import Math.NumberTheory.Padic.Integer+import Math.NumberTheory.Padic.Rational+import Math.NumberTheory.Padic.Analysis+import Data.Word+import Data.Ratio+import Data.Mod+
+ src/Math/NumberTheory/Padic/Analysis.hs view
@@ -0,0 +1,400 @@+{-# LANGUAGE LambdaCase #-}+{-# LANGUAGE FlexibleContexts #-}+{-# OPTIONS_HADDOCK hide, prune, ignore-exports #-}++{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE NoStarIsType #-}++module Math.NumberTheory.Padic.Analysis where++import Math.NumberTheory.Padic.Types+import Data.Mod+import Data.Ratio+import Data.List (unfoldr, genericLength, tails, inits,find)+import Data.Maybe+import GHC.TypeLits hiding (Mod)+import Control.Applicative ((<|>))++------------------------------------------------------------++-- | Unfolds a number to a list of digits (integers modulo @p@). +toRadix :: KnownRadix p => Integer -> [Mod p]+toRadix 0 = [0]+toRadix n = res+ where+ res = unfoldr go n+ p = fromIntegral $ natVal $ head $ 0 : res+ go 0 = Nothing+ go x =+ let (q, r) = quotRem x p+ in Just (fromIntegral r, q)+ +-- | Folds a list of digits (integers modulo @p@) to a number.+fromRadix :: KnownRadix p => [Mod p] -> Integer+fromRadix ds = foldr (\x r -> lifted x + r * p) 0 ds+ where+ p = fromIntegral $ natVal $ head $ 0 : ds++extEuclid :: Integral i => (Integer, Integer) -> Ratio i+extEuclid (n, m) = go (m, 0) (n, 1)+ where+ go (v1, v2) (w1, w2)+ | 2*w1*w1 > abs m =+ let q = v1 `div` w1+ in go (w1, w2) (v1 - q * w1, v2 - q * w2)+ | otherwise = fromRational (w1 % w2)++{- | Extracts p-adic unit from integer number. For radix \(p\) and integer \(n\) returns+pair \((u, k)\) such that \(n = u \cdot p^k\).++Examples:+ +>>> getUnitZ 10 120+(12,1)+>>> getUnitZ 2 120+(15,3)+>>> getUnitZ 3 120+(40,1)+-}+getUnitZ :: (Integral p, Integral n) => p -> n -> (p, Int)+getUnitZ _ 0 = (0, 0)+getUnitZ p x = (b, length v)+ where+ (v, b:_) = span (\n -> n `mod` p == 0) $ iterate (`div` p) $ fromIntegral x++{- | Extracts p-adic unit from a rational number. For radix \(p\) and rational number \(x\) returns+pair \((r/s, k)\) such that \(x = r/s \cdot p^k,\quad \gcd(r, s) = \gcd(s, p) = 1\) and \(p \nmid r\).+++Examples:++>>> getUnitQ 3 (75/157)+(25 % 157, 1)+>>> getUnitQ 5 (75/157)+(3 % 157, 2)+>>> getUnitQ 157 (75/157)+(75 % 1, -1)+>>> getUnitQ 10 (1/60)+(5 % 3, -2)+-}+getUnitQ :: Integral p => p -> Ratio p -> (Ratio p, Int)+getUnitQ _ 0 = (0, 0)+getUnitQ p x = ((n * n') % d, vn - vd)+ where+ go m (d, vd) = case gcd d p of+ 1 -> (m, d, vd)+ c -> let (d', vd') = getUnitZ c d+ in go (m * (p `div` c)^vd') (d', vd + vd')+ (n, vn) = getUnitZ p (numerator x)+ (n', d, vd) = go 1 $ getUnitZ p (denominator x)++-----------------------------------------------------------++{- | For a given list extracts prefix and a cycle, limiting length of prefix and cycle by @len@.+Uses the modified tortiose-and-hare method. -}+findCycle :: Eq a => Int -> [a] -> Maybe ([a], [a])+findCycle len lst =+ find test [ rollback (a, c)+ | (a, cs) <- tortoiseHare len lst+ , c <- take 1 [ c | c <- tail (inits cs)+ , and $ zipWith (==) cs (cycle c) ] ]+ where+ rollback (as, bs) = go (reverse as, reverse bs)+ where+ go =+ \case+ ([], ys) -> ([], reverse ys)+ (x:xs, y:ys)+ | x == y -> go (xs, ys ++ [x])+ (xs, ys) -> (reverse xs, reverse ys)+ test (_, []) = False+ test (pref, c) = and $ zipWith (==) (take len lst) (pref ++ cycle c)++tortoiseHare :: Eq a => Int -> [a] -> [([a], [a])]+tortoiseHare l x =+ map (fmap fst) $+ filter (\(_, (a, b)) -> concat (replicate 3 a) == b) $+ zip (inits x) $+ zipWith splitAt [1 .. l] $ zipWith take [4, 8 ..] $ tails x++ +{- | Returns p-adic solutions (if any) of the equation \(f(x) = 0\) using Hensel lifting method.+First, solutions of \(f(x) = 0\ \mathrm{mod}\ p\) are found, then by Newton's method this solution is get lifted to p-adic number (up to specified precision).++Examples:++>>> henselLifting (\x -> x*x - 2) (\x -> 2*x) :: [Z 7]+[…64112011266421216213,…02554655400245450454]+>>> henselLifting (\x -> x*x - x) (\x -> 2*x-1) :: [Q 10]+[0,1,…92256259918212890625,…07743740081787109376]+-}+henselLifting ::+ (Eq n, PadicNum n, KnownRadix p, Digit n ~ Mod p)+ => (n -> n) -- ^ Function to be vanished.+ -> (n -> n) -- ^ Derivative of the function.+ -> [n] -- ^ The result.+henselLifting f f' = res+ where+ pr = precision (head res)+ res = findSolutionMod f >>= newtonsMethod pr f f'++{- | Returns solution of the equation \(f(x) = 0\ \mathrm{mod}\ p\) in p-adics.+Used as a first step if `henselLifting` function and is usefull for introspection.++>>> findSolutionMod (\x -> x*x - 2) :: [Z 7]+[3,4]+>>> findSolutionMod (\x -> x*x - x) :: [Q 10]+[0.0,1.0,5.0,6.0]+-}+findSolutionMod :: (PadicNum n, KnownRadix p, Digit n ~ Mod p)+ => (n -> n) -> [n]+findSolutionMod f = [ fromMod d | d <- [0..], fm d == 0 ]+ where+ fm = firstDigit . f . fromMod+ fromMod x = fromDigits [x]++newtonsMethod+ :: PadicNum n => Int -> (n -> n) -> (n -> n) -> n -> [n]+newtonsMethod n f f' = iterateM n step+ where+ step x = do+ invf' <- maybeToList (inverse (f' x))+ return (x - f x * invf') ++iterateM :: (Eq a, Monad m) => Int -> (a -> m a) -> a -> m a+iterateM n f = go n+ where+ go 0 x = pure x+ go i x = do+ y <- f x+ if x == y then pure x else go (i - 1) y++{- | Returns a list of m-th roots of unity. -}+unityRoots :: (KnownRadix p, PadicNum n, Digit n ~ Mod p) => Integer -> [n]+unityRoots m = henselLifting f f'+ where+ f x = x^m - 1+ f' x = fromInteger m * x ^ (m - 1) ++pSignum :: (PadicNum n, KnownRadix p, Digit n ~ Mod p) => n -> n+pSignum n+ | d0 == 0 = 0+ | d0^p /= d0 = 1+ | otherwise = case res of+ [] -> 1+ x:_ -> x+ where+ d0 = firstDigit n+ p = radix n+ pr = precision n+ res = newtonsMethod pr (\x -> x^(p - 1) - 1) (\x -> fromInteger (p - 1)*x^(p-2)) (fromDigits [d0])+ +-------------------------------------------------------------++{- | Returns p-adic exponent function, calculated via Taylor series.+For given radix \(p\) converges for numbers which satisfy inequality:++\[|x|_p < p^\frac{1}{1-p}.\]++-}+pExp :: (Eq n, PadicNum n, Fractional n) => n -> Either String n+pExp x | fromRational (norm x) > p ** (-1/(p-1)) = Left "exp does not converge!"+ | otherwise = go (2 * precision x) 0 1 1+ where+ p = fromIntegral (radix x)+ go n s t i+ | n <= 0 = Left "exp failed to converge within precision!"+ | t == 0 = Right s+ | otherwise = go (n - 1) (s + t) (t*x/i) (i+1)++{- | Returns p-adic logarithm function, calculated via Taylor series.+For given radix \(p\) converges for numbers which satisfy inequality:++\[|x|_p < 1.\]++-}+pLog :: (Eq b, PadicNum b, Fractional b) => b -> Either String b+pLog x' | fromRational (norm (x' - 1)) >= 1 = Left "log does not converge!"+ | otherwise = f (x' - 1)+ where+ f x = go (2 * precision x) 0 x 1+ where+ nx = negate x+ go n s t i+ | n <= 0 = Left "log failed to converge within precision!"+ | t == 0 = Right s+ | otherwise = go (n - 1) (s + t/i) (nx*t) (i+1)++{- | Returns p-adic hyperbolic sine function, calculated via Taylor series.+For given radix \(p\) converges for numbers which satisfy inequality:++\[|x|_p < p^\frac{1}{1-p}.\]++-}+pSinh :: (PadicNum b, Fractional b) => b -> Either [Char] b+pSinh x+ | fromRational (norm x) > p ** (-1/(p-1)) = Left "sinh does not converge!"+ | otherwise = go (2 * precision x) 0 x 2+ where+ p = fromIntegral (radix x)+ x2 = x*x+ go n s t i+ | n <= 0 = Left "sinh failed to converge within precision!"+ | t == 0 = Right s+ | otherwise = go (n - 1) (s + t) (t*x2/(i*(i+1))) (i+2)++{- | Returns p-adic inverse hyperbolic sine function, calculated as++\[\mathrm{sinh}^{ -1} x = \log(x + \sqrt{x^2+1})\]++with convergence, corresponding to `pLog` and `pPow` functions.+-}+pAsinh :: (PadicNum b, Fractional b) => b -> Either String b+pAsinh x = do y <- pPow (x*x + 1) (1/2)+ pLog (x + y)++{- | Returns p-adic hyperbolic cosine function, calculated via Taylor series.+For given radix \(p\) converges for numbers which satisfy inequality:++\[|x|_p < p^\frac{1}{1-p}.\]++-}+pCosh :: (PadicNum b, Fractional b) => b -> Either [Char] b+pCosh x+ | fromRational (norm x) > p ** (-1/(p-1)) = Left "cosh does not converge!"+ | otherwise = go (2 * precision x) 0 1 1+ where+ p = fromIntegral (radix x)+ x2 = x*x+ go n s t i+ | n <= 0 = Left "cosh failed to converge within precision!"+ | t == 0 = Right s+ | otherwise = go (n - 1) (s + t) (t*x2/(i*(i+1))) (i+2)++{- | Returns p-adic inverse hyperbolic cosine function, calculated as++\[\mathrm{cosh}^{ -1}\ x = \log(x + \sqrt{x^2-1}),\]++with convergence, corresponding to `pLog` and `pPow` functions.++-}+pAcosh :: (PadicNum b, Fractional b) => b -> Either String b+pAcosh x = do y <- pPow (x*x - 1) (1/2)+ pLog (x + y)++{- | Returns p-adic hyperbolic tan function, calculated as++\[\mathrm{tanh}\ x = \frac{\mathrm{sinh}\ x}{\mathrm{cosh}\ x},\]++with convergence, corresponding to `pSinh` and `pCosh` functions.+-}+pTanh :: (Fractional b, PadicNum b) => b -> Either [Char] b+pTanh x = (/) <$> pSinh x <*> pCosh x++{- | Returns p-adic inverse hyperbolic tan function, calculated as++\[\mathrm{tanh}^{ -1 }\ x = \frac{1}{2} \log\left(\frac{x + 1}{x - 1}\right)\]++with convergence, corresponding to `pLog` function.+-}+pAtanh :: (PadicNum b, Fractional b) => b -> Either String b+pAtanh x = do y <- pLog ((x + 1) / (x - 1))+ return $ y / 2+++{- | Returns p-adic hyperbolic cosine function, calculated via Taylor series.+For given radix \(p\) converges for numbers which satisfy inequality:++\[|x|_p < p^\frac{1}{1-p}.\]++-}+pSin :: (PadicNum b, Fractional b) => b -> Either [Char] b+pSin x+ | fromRational (norm x) > p ** (-1/(p-1)) = Left "sin does not converge!"+ | otherwise = go (2 * precision x) 0 x 2+ where+ p = fromIntegral (radix x)+ x2 = negate x*x+ go n s t i+ | n <= 0 = Left "sin failed to converge within precision!"+ | t == 0 = Right s+ | otherwise = go (n - 1) (s + t) (t*x2/(i*(i+1))) (i+2)++{- | Returns p-adic cosine function, calculated via Taylor series.+For given radix \(p\) converges for numbers which satisfy inequality:++\[|x|_p < p^\frac{1}{1-p}.\]++-}+pCos :: (PadicNum b, Fractional b) => b -> Either [Char] b+pCos x+ | fromRational (norm x) > p ** (-1/(p-1)) = Left "cos does not converge!"+ | otherwise = go (2 * precision x) 0 1 1+ where+ p = fromIntegral (radix x)+ x2 = negate x*x+ go n s t i+ | n <= 0 = Left "cos failed to converge within precision!"+ | t == 0 = Right s+ | otherwise = go (n - 1) (s + t) (t*x2/(i*(i+1))) (i+2)++{- | Returns p-adic arcsine function, calculated via Taylor series.+For given radix \(p\) converges for numbers which satisfy inequality:++\[|x|_p < 1.\]++-}+pAsin x | norm x >= 1 = Left "asin does not converge!"+ | otherwise = Right $+ sum $ takeWhile (\t -> valuation t < pr) $+ take (2*pr) $ zipWith (*) xs cs+ where+ pr = precision x+ x2 = x*x+ xs = iterate (x2 *) x+ cs = zipWith (/) (zipWith (/) n2f nf2) [2*n+1 | n <- fromInteger <$> [0..]]+ n2f = scanl (*) 1 [n*(n+1) | n <- fromInteger <$> [1,3..]]+ nf2 = scanl (*) 1 [4*n^2 | n <- fromInteger <$> [1..]]++{- | Returns p-adic square root, calculated for odd radix via Hensel lifting,+and for \(p=2\) by recurrent product.+-}+pSqrt ::+ ( Fractional n+ , PadicNum n+ , KnownRadix p+ , Digit n ~ Mod p+ )+ => n -> [n]+pSqrt x+ | odd (radix x) && isSquareResidue x =+ henselLifting (\y -> y * y - x) (2 *)+ | lifted x `mod` 4 /= 3 && lifted x `mod` 8 == 1 =+ let r = pSqrt2 x in [r, -r]+ | otherwise = []++pSqrt2 :: (PadicNum a, Fractional a) => a -> a+pSqrt2 a = product $+ takeWhile (/= 1)+ $ take (2*precision a)+ $ go ((a - 1) / 8)+ where+ go x = (1 + 4*x) : go ((-2)*(x / (1 + 4*x))^2)++{- | Exponentiation for p-adic numbers, calculated as++\[ x^y = e^{y \log x}, \]++with convergence, corresponding to `pExp` and `pLog` functions.+-}+pPow :: (PadicNum p, Fractional p) => p -> p -> Either String p+pPow x y = case pLog x >>= \z -> pExp (z*y) of+ Right res -> Right res+ Left _ -> Left "exponentiation doesn't converge!"++{- | Returns @True@ for p-adics with square residue as a first digit.+-}+isSquareResidue :: (PadicNum n, KnownRadix p, Digit n ~ Mod p) => n -> Bool+isSquareResidue x = any (\m -> m*m == firstDigit x) [0..]+
+ src/Math/NumberTheory/Padic/Integer.hs view
@@ -0,0 +1,133 @@+{-# OPTIONS_HADDOCK hide, prune, ignore-exports #-}++{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE FlexibleInstances #-}++module Math.NumberTheory.Padic.Integer where++import Data.List (intercalate)+import Data.Mod+import Data.Word+import Data.Ratio+import GHC.TypeLits (Nat, natVal)+import GHC.Integer.GMP.Internals (recipModInteger)+import Math.NumberTheory.Padic.Types+import Math.NumberTheory.Padic.Analysis++------------------------------------------------------------++type instance Padic Integer p = Z' p (SufficientPrecision Word p)+type instance Padic Int p = Z' p (SufficientPrecision Int p)+type instance Padic Word8 p = Z' p (SufficientPrecision Word8 p)+type instance Padic Word16 p = Z' p (SufficientPrecision Word16 p)+type instance Padic Word32 p = Z' p (SufficientPrecision Word32 p)+type instance Padic Word64 p = Z' p (SufficientPrecision Word64 p)+type instance Padic Word p = Z' p (SufficientPrecision Word64 p)+++-- | Integer p-adic number (an element of \(\mathbb{Z}_p\)) with default precision.+type Z p = Z' p (SufficientPrecision Word32 p)++-- | Integer p-adic number with explicitly specified precision.+newtype Z' (p :: Nat) (prec :: Nat) = Z' (R prec p)+newtype R (prec :: Nat ) (p :: Nat) = R (Mod (LiftedRadix p prec))++instance Radix p prec => Eq (Z' p prec) where+ x@(Z' (R a)) == Z' (R b) = unMod a `mod` pk == unMod b `mod` pk+ where+ pk = radix x ^ precision x++instance Radix p prec => PadicNum (Z' p prec) where+ type Unit (Z' p prec) = Z' p prec+ type Digit (Z' p prec) = Mod p ++ {-# INLINE precision #-}+ precision = fromIntegral . natVal++ {-# INLINE radix #-}+ radix (Z' r) = fromIntegral $ natVal r+ + {-# INLINE fromDigits #-}+ fromDigits = mkUnit . fromRadix++ {-# INLINE digits #-}+ digits n = toRadix (lifted n)++ {-# INLINE lifted #-}+ lifted (Z' (R n)) = lifted n++ {-# INLINE mkUnit #-}+ mkUnit = Z' . R . fromInteger++ {-# INLINE fromUnit #-}+ fromUnit (u, v) = mkUnit $ radix u ^ fromIntegral v * lifted u++ splitUnit n = case getUnitZ (radix n) (lifted n) of+ (0, 0) -> (0, precision n)+ (u, v) -> (mkUnit u, v)+ + isInvertible n = (lifted n `mod` p) `gcd` p == 1+ where+ p = radix n+ + inverse (Z' (R n)) = Z' . R <$> invertMod n++instance Radix p prec => Show (Z' p prec) where+ show n = + case findCycle pr ds of+ Nothing | length ds > pr -> ell ++ toString (take pr ds)+ | otherwise -> toString ds+ Just ([],[0]) -> "0"+ Just (pref, [0]) -> toString pref+ Just (pref, cyc)+ | length pref + length cyc <= pr ->+ let sp = toString pref+ sc = toString cyc+ in "(" ++ sc ++ ")" ++ sep ++ sp+ | otherwise -> ell ++ toString (take pr $ pref ++ cyc)+ where+ pr = precision n+ ds = digits n+ showD = show . unMod+ toString = intercalate sep . map showD . reverse+ ell = "…" ++ sep + sep+ | radix n < 11 = ""+ | otherwise = " "++instance Radix p prec => Num (Z' p prec) where+ fromInteger = Z' . R . fromInteger+ Z' (R a) + Z' (R b) = Z' . R $ a + b+ Z' (R a) - Z' (R b) = Z' . R $ a - b+ Z' (R a) * Z' (R b) = Z' . R $ a * b+ negate (Z' (R a)) = Z' . R $ negate a+ abs x = if valuation x == 0 then 1 else 0+ signum = pSignum++instance Radix p prec => Enum (Z' p prec) where+ toEnum = fromIntegral+ fromEnum = fromIntegral . toInteger++instance Radix p prec => Real (Z' p prec) where+ toRational 0 = 0+ toRational n = extEuclid (lifted n, liftedRadix n)++instance Radix p prec => Integral (Z' p prec) where+ toInteger n = if denominator r == 1+ then numerator r+ else lifted n `mod` (radix n ^ precision n)+ where+ r = toRational n+ a `quotRem` b = case inverse b of+ Nothing -> error $ show b ++ " is not divisible modulo " ++ show (radix a) ++ "!" + Just r -> let q = a*r in (q, a - q * b)+++instance Radix p prec => Ord (Z' p prec) where+ compare = error "ordering is not defined for Z"++{-| Integer power function (analog of (^) operator ) -}+zPow :: Radix p prec => Z' p prec -> Z' p prec -> Z' p prec+zPow (Z' (R a)) (Z' (R b)) = Z' . R $ a ^% fromIntegral (unMod b)
+ src/Math/NumberTheory/Padic/Rational.hs view
@@ -0,0 +1,171 @@+{-# OPTIONS_HADDOCK hide, prune, ignore-exports #-}++{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UndecidableInstances #-}++module Math.NumberTheory.Padic.Rational where++import Data.List (intercalate)+import Data.Ratio+import Data.Mod+import Data.Word+import GHC.TypeLits (Nat, natVal)+import Math.NumberTheory.Padic.Types+import Math.NumberTheory.Padic.Analysis+import Math.NumberTheory.Padic.Integer++------------------------------------------------------------++type instance Padic Rational p = Q' p (SufficientPrecision Word p)+type instance Padic (Ratio Int) p = Q' p (SufficientPrecision Int p)+type instance Padic (Ratio Word8) p = Q' p (SufficientPrecision Word8 p)+type instance Padic (Ratio Word16) p = Q' p (SufficientPrecision Word16 p)+type instance Padic (Ratio Word32) p = Q' p (SufficientPrecision Word32 p)+type instance Padic (Ratio Word64) p = Q' p (SufficientPrecision Word64 p)+type instance Padic (Ratio Word) p = Q' p (SufficientPrecision Word64 p)++------------------------------------------------------------+-- | Rational p-adic number (an element of \(\mathbb{Q}_p\)) with default precision.+type Q p = Q' p (SufficientPrecision Word32 p)++-- | Rational p-adic number with explicitly specified precision.+data Q' (p :: Nat) (prec :: Nat) = Q' !(Z' p prec) !Int++instance Radix p prec => PadicNum (Q' p prec) where+ type Unit (Q' p prec) = Z' p prec+ type Digit (Q' p prec) = Digit (Z' p prec)++ {-# INLINE precision #-}+ precision = fromIntegral . natVal++ {-# INLINE radix #-}+ radix (Q' u _) = radix u++ {-# INLINE digits #-}+ digits (Q' u v) = replicate v 0 ++ toRadix (lifted u)++ {-# INLINE fromDigits #-}+ fromDigits ds = normalize $ Q' (fromDigits ds) 0++ {-# INLINE lifted #-}+ lifted (Q' u _) = lifted u++ {-# INLINE mkUnit #-}+ mkUnit ds = normalize $ Q' (mkUnit ds) 0++ {-# INLINE fromUnit #-}+ fromUnit (u, v) = Q' u v++ splitUnit n@(Q' u v) =+ let pr = precision n+ (u', v') = splitUnit u+ in if v + v' > pr then (0, pr) else (u', v + v') + + isInvertible = isInvertible . unit . normalize+ + inverse n = do r <- inverse (unit n)+ return $ fromUnit (r, - valuation n)++instance Radix p prec => Show (Q' p prec) where+ show n = si ++ sep ++ "." ++ sep ++ sf+ where+ (u, k) = splitUnit (normalize n)+ pr = precision n + ds = digits u+ (f, i) =+ case compare k 0 of+ EQ -> ([0], ds)+ GT -> ([0], replicate k 0 ++ ds)+ LT -> splitAt (-k) (ds ++ replicate (pr + k) 0)+ sf = intercalate sep $ showD <$> reverse f+ si =+ case findCycle pr i of+ Nothing+ | null i -> "0"+ | length i > pr -> ell ++ toString (take pr i)+ | otherwise -> toString i+ Just ([], [0]) -> "0"+ Just (pref, [0]) -> toString pref+ Just (pref, cyc)+ | length pref + length cyc <= pr ->+ let sp = toString pref+ sc = toString cyc+ in "(" ++ sc ++ ")" ++ sep ++ sp+ | otherwise -> ell ++ toString (take pr $ pref ++ cyc)+ showD = show . unMod+ toString = intercalate sep . map showD . reverse+ ell = "…" ++ sep+ sep+ | radix n < 11 = ""+ | otherwise = " "+ +instance Radix p prec => Eq (Q' p prec) where+ a' == b' =+ (isZero a && isZero b)+ || (valuation a == valuation b && unit a == unit b)+ where+ a = normalize a'+ b = normalize b'++instance Radix p prec => Ord (Q' p prec) where+ compare = error "Order is nor defined for p-adics."++instance Radix p prec => Num (Q' p prec) where+ fromInteger n = normalize $ Q' (fromInteger n) 0+ + x@(Q' (Z' (R a)) va) + Q' (Z' (R b)) vb =+ case compare va vb of+ LT -> Q' (Z' (R (a + p ^% (vb - va) * b))) va+ EQ -> Q' (Z' (R (a + b))) va+ GT -> Q' (Z' (R (p ^% (va - vb) * a + b))) vb+ where+ p = fromInteger (radix x)+ + Q' (Z' (R a)) va * Q' (Z' (R b)) vb = Q' (Z' (R (a * b))) (va + vb)+ + negate (Q' u v) = Q' (negate u) v+ abs = fromRational . norm+ signum = pSignum++newtype Delay prec p = Delay (Q' p prec)++instance Radix p prec => Fractional (Q' p prec) where+ fromRational 0 = 0+ fromRational x = res+ where+ res = Q' (n `div` d) v+ p = fromInteger $ natVal (Delay res)+ (q, v) = getUnitQ p x+ (n, d) = (mkUnit $ numerator q, mkUnit $ denominator q)+ a / b = Q' u (v + valuation a - valuation b')+ where+ b' = normalize b+ Q' u v = fromRational (lifted a % lifted b')++instance Radix p prec => Real (Q' p prec) where+ toRational n = toRational (unit n) / norm n++pUndefinedError s = error $ s ++ " is undifined for p-adics."++fromEither = either error id++instance Radix p prec => Floating (Q' p prec) where+ x ** y = fromEither $ pPow x y+ exp = fromEither . pExp+ log = fromEither . pLog+ sinh = fromEither . pSinh+ cosh = fromEither . pCosh+ sin = fromEither . pSin+ cos = fromEither . pCos+ asinh = fromEither . pAsinh+ acosh = fromEither . pCosh+ atanh = fromEither . pTanh+ asin = fromEither . pAsin+ sqrt x = case pSqrt x of+ res:_ -> res+ [] -> error $ "sqrt: digit " ++ show (firstDigit x) ++ " is not a square residue!"+ pi = pUndefinedError "pi"+ acos = pUndefinedError "acos"+ atan = pUndefinedError "atan"+
+ src/Math/NumberTheory/Padic/Types.hs view
@@ -0,0 +1,280 @@+{-# OPTIONS_HADDOCK hide, prune, ignore-exports #-}++{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE UndecidableSuperClasses #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE NoStarIsType #-}++module Math.NumberTheory.Padic.Types where++import Data.Ratio+import Data.Maybe (isJust, maybeToList)+import Data.Mod +import Data.Word +import GHC.TypeLits hiding (Mod)+import Data.Constraint (Constraint)+import GHC.Integer (smallInteger)+import GHC.Integer.Logarithms ( integerLogBase# )++------------------------------------------------------------+-- | Constraint for non-zero natural number which can be a radix.+type family ValidRadix (m :: Nat) :: Constraint where+ ValidRadix 0 = TypeError ('Text "Zero radix!")+ ValidRadix 1 = TypeError ('Text "Radix should be more then 1!")+ ValidRadix m = ()++-- | Constraint for valid radix of a number+type KnownRadix m = ( ValidRadix m , KnownNat m )+ +-- | Radix of the internal representation of integer p-adic number.+type family LiftedRadix p prec where+ LiftedRadix p prec = p ^ (2*prec + 1)++-- | Constraint for known valid radix of p-adic number as well as it's lifted radix.+type family Radix p prec :: Constraint where+ Radix p prec =+ ( KnownNat prec+ , KnownRadix p+ , KnownRadix (LiftedRadix p prec)+ )++{- | Precision sufficient for rational reconstruction of number belonging to a type @num@.+Used in a type declaration as follows:++>>> x = 1 `div` 1234567898765432123456789 :: Z 2 (Sufficientprecision Word32 2)+>>> toRational x+13822228938088947473 % 12702006275138148709+>>> x = 1 `div` 1234567898765432123456789 :: Z 2 (Sufficientprecision Int 2)+>>> toRational x+1 % 1234567898765432123456789++-} +type family SufficientPrecision num (p :: Nat) :: Nat where+ SufficientPrecision Word32 2 = 64+ SufficientPrecision Word32 3 = 43+ SufficientPrecision Word32 5 = 30+ SufficientPrecision Word32 6 = 26+ SufficientPrecision Word32 7 = 24+ SufficientPrecision Word8 p = Div (SufficientPrecision Word32 p) 4+ SufficientPrecision Word16 p = Div (SufficientPrecision Word32 p) 2+ SufficientPrecision Word64 p = 2 * (SufficientPrecision Word32 p) + 1+ SufficientPrecision Int p = 2 * SufficientPrecision Word32 p+ SufficientPrecision Word p = SufficientPrecision Word64 p+ SufficientPrecision (Ratio t) p = SufficientPrecision t p+ SufficientPrecision t p = Div (SufficientPrecision t 2) (Log2 p)++{- | Type family for p-adic numbers with precision defined by reconstructable number type.++>>>123456 :: Padic Int 7+1022634+>>> toInteger it+123456+>>> toRational (12345678987654321 :: Padic (Ratio Word16) 3)+537143292837 % 5612526479 -- insufficiend precision for proper reconstruction!!+>>> toRational (12345678987654321 :: Padic Rational 3)+12345678987654321 % 1++-}+type family Padic num (p :: Nat)++------------------------------------------------------------+{- | Typeclass for p-adic numbers.++-}+class (Eq n, Num n) => PadicNum n where+ -- | A type for p-adic unit.+ type Unit n+ -- | A type for digits of p-adic expansion.+ -- Associated type allows to assure that digits will agree with the radix @p@ of the number.+ type Digit n+ -- | Returns the precision of a number.+ --+ -- Examples:+ --+ -- >>> precision (123 :: Z 2)+ -- 20+ -- >>> precision (123 :: Z' 2 40)+ -- 40+ precision :: Integral i => n -> i+ -- | Returns the radix of a number+ --+ -- Examples:+ --+ -- >>> radix (5 :: Z 13)+ -- 13+ -- >>> radix (-5 :: Q' 3 40)+ -- 3+ radix :: Integral i => n -> i+ + -- | Constructor for a digital object from it's digits+ fromDigits :: [Digit n] -> n+ -- | Returns digits of a digital object+ --+ -- Examples:+ --+ -- >>> digits (123 :: Z 10)+ -- [(3 `modulo` 10),(2 `modulo` 10),(1 `modulo` 10),(0 `modulo` 10),(0 `modulo` 10)]+ -- >>> take 5 $ digits (-123 :: Z 2)+ -- [(1 `modulo` 2),(0 `modulo` 2),(1 `modulo` 2),(0 `modulo` 2),(0 `modulo` 2)]+ -- >>> take 5 $ digits (1/300 :: Q 10)+ -- [(7 `modulo` 10),(6 `modulo` 10),(6 `modulo` 10),(6 `modulo` 10),(6 `modulo` 10)]+ --+ digits :: n -> [Digit n]+ -- | Returns lifted digits+ --+ -- Examples:+ --+ -- >>> lifted (123 :: Z 10)+ -- 123+ -- >>> lifted (-123 :: Z 10)+ -- 9999999999999999999999999999999999999999877+ --+ lifted :: n -> Integer++ -- | Creates digital object from it's lifted digits.+ mkUnit :: Integer -> n++ -- | Creates p-adic number from given unit and valuation.+ --+ -- prop> fromUnit (u, v) = u * p^v+ fromUnit :: (Unit n, Int) -> n+ -- | Splits p-adic number into unit and valuation.+ --+ -- prop> splitUnit (u * p^v) = (u, v)+ splitUnit :: n -> (Unit n, Int)+ + -- | Returns @True@ for a p-adic number which is multiplicatively invertible.+ isInvertible :: n -> Bool++ -- | Partial multiplicative inverse of p-adic number (defined both for integer or rational p-adics).+ inverse :: n -> Maybe n+ ++{- | The least significant digit of a p-adic number.+--+-- >>> firstDigit (123 :: Z 10)+-- (3 `modulo` 10)+-- >>> firstDigit (123 :: Z 257)+-- (123 `modulo` 257)+-}+firstDigit :: PadicNum n => n -> Digit n+{-# INLINE firstDigit #-}+firstDigit = head . digits++{- | Returns p-adic number reduced modulo @p@++>>> reduce (123 :: Z 10) :: Mod 100+(23 `modulo` 100)+-}+reduce :: (KnownRadix p, PadicNum n) => n -> Mod p+reduce = fromIntegral . lifted+++{- | Returns the p-adic unit of a number++Examples:++>>> unit (120 :: Z 10)+12+>>> unit (75 :: Z 5)+3 -}+unit :: PadicNum n => n -> Unit n+{-# INLINE unit #-}+unit = fst . splitUnit++{- | Returns a p-adic valuation of a number++Examples:++>>> valuation (120 :: Z 10)+1+>>> valuation (75 :: Z 5)+2++Valuation of zero is equal to working precision++>>> valuation (0 :: Q 2)+64+>>> valuation (0 :: Q 10)+21 -}+valuation :: PadicNum n => n -> Int+{-# INLINE valuation #-}+valuation = snd . splitUnit++{- | Returns a rational p-adic norm of a number \(|x|_p\).++Examples:++>>> norm (120 :: Z 10)+0.1+>>> norm (75 :: Z 5)+4.0e-2+-}+norm :: (Integral i, PadicNum n) => n -> Ratio i+{-# INLINE norm #-}+norm n = (radix n % 1) ^^ (-valuation n)++{- | Adjusts unit and valuation of p-adic number, by removing trailing zeros from the right-side of the unit.++Examples:++>>> λ> x = 2313 + 1387 :: Q 10+>>> x+3700.0+>>> splitUnit x+(3700,0)+>>> splitUnit (normalize x)+(37,2) -}+normalize :: PadicNum n => n -> n+normalize = fromUnit . splitUnit++-- | Returns @True@ for a p-adic number which is equal to zero (within it's precision).+isZero :: PadicNum n => n -> Bool+{-# INLINE isZero #-}+isZero n = valuation n >= precision n++liftedRadix :: (PadicNum n, Integral a) => n -> a+{-# INLINE liftedRadix #-}+liftedRadix n = radix n ^ (2*precision n + 1)++{- | For given radix \(p\) and natural number \(m\) returns precision sufficient for rational+reconstruction of fractions with numerator and denominator not exceeding \(m\).++Examples:++>>> sufficientPrecision 2 (maxBound :: Int)+64+>>> sufficientPrecision 3 (maxBound :: Int)+41+>>> sufficientPrecision 10 (maxBound :: Int)+20+-}+sufficientPrecision :: Integral a => Integer -> a -> Integer+sufficientPrecision p m = ilog p (2 * fromIntegral m) + 2++ilog :: (Integral a, Integral b) => a -> a -> b+ilog a b =+ fromInteger $ smallInteger (integerLogBase# (fromIntegral a) (fromIntegral b))+++-----------------------------------------------------------++instance KnownRadix p => PadicNum (Mod p) where+ type Unit (Mod p) = Mod p+ type Digit (Mod p) = Mod p+ radix = fromIntegral . natVal+ precision _ = fromIntegral (maxBound :: Int)+ digits = pure+ fromDigits = head+ lifted = fromIntegral . unMod+ mkUnit = fromInteger+ inverse = invertMod+ isInvertible = isJust . invertMod + fromUnit (u, 0) = u+ fromUnit _ = 0+ splitUnit u = (u, 0)+
+ test/Spec.hs view
@@ -0,0 +1,20 @@+module Main where++import Test.Tasty+import qualified Test.Commons+import qualified Test.Integer+import qualified Test.Rational+import qualified Test.Analysis++-----------------------------------------------------------++testSuite :: TestTree+testSuite = testGroup "Padic module"+ [+ Test.Commons.testSuite+ , Test.Integer.testSuite+ , Test.Rational.testSuite + -- Test.Analysis.testSuite+ ]++main = defaultMain testSuite
+ test/Test/Analysis.hs view
@@ -0,0 +1,50 @@+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-}++module Test.Analysis (testSuite) where++import Math.NumberTheory.Padic+import Test.Base+import GHC.TypeLits hiding (Mod)+import Test.Tasty+import Test.Tasty.HUnit+import Test.Tasty.QuickCheck+import Test.QuickCheck+import Data.Mod+import Data.Ratio++a @/= b = assertBool "" (a /= b)++------------------------------------------------------------+mulSqrt :: (Eq n, Show n, PadicNum n, Fractional n, KnownRadix p, Digit n ~ Mod p) => n -> n -> Bool+mulSqrt w a = and $ do+ x <- pSqrt a+ return $ x*x == a++testSqrt = testGroup "pSqrt properties" $+ [+ testProperty "Q 2" $ mulSqrt (0 :: Q 2)+ , testProperty "Q 7" $ mulSqrt (0 :: Q 7)+ ]++mulExp :: (Eq n, Show n, PadicNum n, Fractional n) => n -> n -> n -> Bool+mulExp w a b = either (const True) id $ do+ x <- pExp a+ y <- pExp b+ z <- pExp (a + b)+ return $ x*y == z++testExp = testGroup "pExp properties" $+ [+ testProperty "multiplication" $ mulExp (0 :: Q 2)+ , testProperty "multiplication" $ mulExp (0 :: Q 7)+ ]++testSuite = testGroup "Analysis" $+ [+ testSqrt + ]
+ test/Test/Base.hs view
@@ -0,0 +1,140 @@+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-}++module Test.Base where++import Math.NumberTheory.Padic++import GHC.TypeLits hiding (Mod)+import Test.Tasty+import Test.Tasty.HUnit+import Test.Tasty.QuickCheck+import Test.Tasty.ExpectedFailure+import Test.QuickCheck+import Data.Mod+import Data.Word+import Data.Ratio+import Data.Maybe++instance KnownRadix m => Arbitrary (Mod m) where+ arbitrary = fromInteger <$> arbitrary++instance Radix p prec => Arbitrary (Z' p prec) where+ arbitrary = oneof [integerZ, rationalZ, arbitraryZ]+ where+ integerZ = fromInteger <$> arbitrary+ arbitraryZ = fromDigits . take 20 <$> infiniteList+ rationalZ = do+ a <- integerZ+ b <- suchThat integerZ isInvertible+ return $ a `div` b+ shrink _ = []++instance Radix p prec => Arbitrary (Q' p prec) where+ arbitrary = oneof [integerQ, rationalQ, arbitraryQ]+ where+ integerQ = fromInteger <$> arbitrary+ arbitraryQ = fromDigits . take 20 <$> infiniteList+ rationalQ = do+ SmallRational r <- arbitrary+ return $ fromRational r+ shrink _ = []++newtype SmallRational = SmallRational (Rational)+ deriving (Show, Eq, Num, Fractional)++instance Arbitrary SmallRational where+ arbitrary = do+ let m = fromIntegral (maxBound :: Word32)+ n <- chooseInteger (-m, m)+ d <- chooseInteger (1,m)+ return $ SmallRational (n % d)+ shrink (SmallRational r) = SmallRational <$> shrink r+ +a @/= b = assertBool "" (a /= b)++homo0 :: (Show a, Eq a) => (a -> t) -> (t -> a) -> t -> a -> Property+homo0 phi psi w a =+ let [x, _] = [phi a, w] in psi x === a++homo1 :: (Show t , Eq t) => (a -> t)+ -> (a -> a -> a)+ -> (t -> t -> t)+ -> t -> a -> a -> Property+homo1 phi opa opt w a b =+ let [x, y, _] = [phi a, phi b, w]+ in x `opt` y === phi (a `opa` b)++homo2 :: (Show a, Eq a) => (a -> t) -> (t -> a)+ -> (a -> a -> a)+ -> (t -> t -> t)+ -> t -> a -> a -> Property+homo2 phi psi opa opt w a b =+ let [x, y, _] = [phi a, phi b, w]+ in psi (x `opt` y) === a `opa` b++invOp :: (Show t, Eq t) => (a -> t) + -> (t -> t -> t) + -> (t -> t)+ -> (t -> Bool)+ -> t -> a -> a -> Property+invOp phi op inv p w a b =+ let [x, y, _] = [phi a, phi b, w]+ in p y ==> (x `op` y) `op` inv y === x + +addComm :: (Show a, Eq a, Num a) => a -> a -> a -> Bool+addComm t a b = a + b == b + a++addAssoc :: (Show a, Eq a, Num a) => a -> a -> a -> a -> Bool+addAssoc t a b c = a + (b + c) == (a + b) + c++negInvol :: (Show a, Eq a, Num a) => a -> a -> Bool+negInvol t a = - (- a) == a++negInvers :: (Eq a, Num a) => a -> a -> Bool+negInvers t a = a - a == 0++negScale :: (Eq a, Num a) => a -> a -> Bool+negScale t a = (-1) * a == - a++mulZero :: (Eq a, Num a) => a -> a -> Bool+mulZero t a = 0 * a == 0++mulOne :: (Eq a, Num a) => a -> a -> Bool+mulOne t a = 1 * a == a++mulComm :: (Eq a, Num a) => a -> a -> a -> Bool+mulComm t a b = a * b == b * a++mulAssoc :: (Eq a, Num a) => a -> a -> a -> a -> Bool+mulAssoc t a b c = a * (b * c) == (a * b) * c++mulDistr :: (Eq a, Num a) => a -> a -> a -> a -> Bool+mulDistr t a b c = a * (b + c) == a * b + a * c+ +divDistr :: (Eq a, Fractional a) => a -> a -> a -> a -> Property+divDistr t a b c = a /= 0 ==> (b + c) / a == b / a + c / a+ +divMul :: (Eq a, Fractional a) => a -> a -> a -> Property+divMul t a b = b /= 0 ==> (a / b) * b == a++mulSign :: (Eq a, Num a) => a -> a -> a -> Bool+mulSign t a b = and [a * (- b) == - (a * b), (- a) * (- b) == a * b]++ringLaws t = testGroup "Ring laws" $+ [ testProperty "Addition commutativity" $ addComm t+ , testProperty "Addition associativity" $ addAssoc t+ , testProperty "Negation involution" $ negInvol t+ , testProperty "Addition inverse" $ negInvers t+ , testProperty "Negative scaling" $ negScale t+ , testProperty "Multiplicative zero" $ mulZero t+ , testProperty "Multiplicative one" $ mulOne t+ , testProperty "Multiplication commutativity" $ mulComm t+ , testProperty "Multiplication associativity" $ mulAssoc t+ , testProperty "Multiplication distributivity" $ mulDistr t+ , testProperty "Multiplication signs" $ mulSign t+ ]
+ test/Test/Commons.hs view
@@ -0,0 +1,100 @@+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-}++module Test.Commons (testSuite) where++import Math.NumberTheory.Padic+import GHC.TypeLits hiding (Mod)+import Test.Tasty+import Test.Tasty.HUnit+import Test.Tasty.QuickCheck+import Test.QuickCheck+import Data.Mod+import Data.Ratio++instance KnownRadix m => Arbitrary (Mod m) where+ arbitrary = fromInteger <$> arbitrary++newtype AnyRadix = AnyRadix Integer+ deriving (Show, Eq, Num)+ +instance Arbitrary AnyRadix where+ arbitrary = AnyRadix <$> arbitrary `suchThat` (> 1)+ ++a @/= b = assertBool "" (a /= b)++------------------------------------------------------------++cycleTest :: TestTree+cycleTest = testGroup "findCycle tests"+ [ testCase "1" $ findCycle 10 [1..5] @?= Nothing+ , testCase "2" $ findCycle 10 [1] @?= Nothing+ , testCase "3" $ findCycle 10 (repeat 1) @?= Just ([],[1])+ , testCase "4" $ findCycle 10 ([1..5] ++ repeat 1) @?= Just ([1..5],[1])+ , testCase "5" $ findCycle 10 ([1..15] ++ repeat 1) @?= Nothing+ , testCase "6" $ findCycle 10 ([1,1,1] ++ repeat 1) @?= Just ([],[1])+ , testCase "7" $ findCycle 10 ([2,1,1] ++ repeat 1) @?= Just ([2],[1])+ , testCase "8" $ findCycle 10 ([1,2,3] ++ cycle [1,2,3]) @?= Just ([],[1,2,3])+ , testCase "9" $ findCycle 10 ([2,3] ++ cycle [1,2,3]) @?= Just ([],[2,3,1])+ , testCase "10" $ findCycle 10 ([0,1,2,3] ++ cycle [1,2,3]) @?= Just ([0],[1,2,3])+ , testCase "11" $ findCycle 10 ([0,2,3] ++ cycle [1,2,3]) @?= Just ([0],[2,3,1])+ , testCase "12" $ findCycle 200 ([1..99] ++ cycle [100..200]) @?= Just ([1..99],[100..200]) ]++------------------------------------------------------------++radixTest :: KnownRadix p => Mod p -> Positive Integer -> Bool+radixTest m (Positive n) =+ let ds = tail (m : toRadix n)+ n' = fromRadix ds + in n' == n && toRadix n' == ds++radixTests = testGroup "Conversion to and from digits"+ [ testProperty "2" $ radixTest (0 :: Mod 2)+ , testProperty "10" $ radixTest (0 :: Mod 10)+ , testProperty "65536" $ radixTest (0 :: Mod 65536)+ ]+ +------------------------------------------------------------+getUnitTests = testGroup "p-adic units."+ [ testGroup "Integers"+ [ testCase "2" $ getUnitZ 2 (4) @?= (1, 2)+ , testCase "7(2)" $ getUnitZ 2 (28) @?= (7, 2)+ , testCase "7(7)" $ getUnitZ 7 (28) @?= (4, 1)+ , testProperty "0" $ \(AnyRadix p) -> getUnitZ p 0 === (0, 0)+ , testProperty "1" $ \(AnyRadix p) -> getUnitZ p 1 === (1, 0)+ , testProperty "p" $+ \(AnyRadix p) r -> let (u, k) = getUnitZ p r+ in r === fromIntegral p ^ k * u+ ]+ , testGroup "Rationals"+ [ testCase "2" $ getUnitQ 2 (4%7) @?= (1 % 7, 2)+ , testCase "7" $ getUnitQ 7 (4%7) @?= (4 % 1, -1)+ , testCase "10" $ getUnitQ 10 (1%20) @?= (5 % 1, -2)+ , unitProperties+ ]+ ]++unitProperties = testGroup "Unit properties" $+ [ testProperty "reconstruction from unit" $+ \(AnyRadix p) r -> let (u, k) = getUnitQ p r+ in r === fromIntegral p^^k * u+ , testProperty "p does not divide numerator" $+ \(AnyRadix p) r -> let (u, k) = getUnitQ p r+ in u /= 0 ==> numerator u `mod` p =/= 0+ , testProperty "denominator is comprime with p" $+ \(AnyRadix p) r -> let (u, k) = getUnitQ p r+ in u /= 0 ==> gcd (denominator u) p == 1+ , testProperty "0" $ \(AnyRadix p) -> getUnitQ p 0 === (0, 0)+ , testProperty "1" $ \(AnyRadix p) -> getUnitQ p 1 === (1, 0)+ ]++------------------------------------------------------------+testSuite = testGroup "Commons"+ [ cycleTest+ , radixTests+ , getUnitTests ]
+ test/Test/Integer.hs view
@@ -0,0 +1,115 @@+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-}++module Test.Integer (testSuite) where++import Math.NumberTheory.Padic+import Test.Base+import GHC.TypeLits hiding (Mod)+import Test.Tasty+import Test.Tasty.HUnit+import Test.Tasty.QuickCheck+import Test.Tasty.ExpectedFailure+import Test.QuickCheck+import Data.Mod++------------------------------------------------------------+digitsTestZ :: (Show n, Eq n, PadicNum n) => n -> n -> Property+digitsTestZ t n = fromDigits (digits n) === n++digitsTests = testGroup "Conversion to and from digits"+ [ testProperty "Z 2" $ digitsTestZ (0 :: Z 2)+ , testProperty "Z 10" $ digitsTestZ (0 :: Z 10)+ , testProperty "Z 257" $ digitsTestZ (0 :: Z 257)+ ]+ +------------------------------------------------------------+equivTest :: TestTree+equivTest = testGroup "Equivalence tests"+ [ testCase "1" $ (0 :: Z' 10 5) @?= 432100000+ , testCase "2" $ (0 :: Z' 10 5) @/= 543210000+ , testCase "3" $ (87654321 :: Z' 10 5) @?= 87054321+ , testCase "4" $ (87654321 :: Z' 10 5) @/= 87604321+ ]++------------------------------------------------------------+showTests = testGroup "String representation"+ [ testCase "0" $ show (0 :: Z 3) @?= "0"+ , testCase "3" $ show (3 :: Z 3) @?= "10"+ , testCase "-3" $ show (-3 :: Z 3) @?= "(2)0"+ , testCase "123" $ show (123 :: Z 10) @?= "123"+ , testCase "123" $ show (123 :: Z 2) @?= "1111011"+ , testCase "123456789" $ show (123456789 :: Z' 10 5) @?= "…56789"+ , testCase "-123" $ show (-123 :: Z 10) @?= "(9)877"+ , testCase "1/23" $ show (1 `div` 23 :: Z 10) @?= "…565217391304347826087"+ , testCase "1/23" $ show (1 `div` 23 :: Z' 17 5) @?= "… 8 14 13 5 3"+ , testCase "123456" $ show (123456 :: Z' 257 4) @?= "1 223 96"+ , testCase "123456" $ show (-12345 :: Z 257) @?= "(256) 208 248"+ , testCase "123456" $ show (-123456 :: Z 257) @?= "… 256 256 256 256 256 255 33 161"+ ]++------------------------------------------------------------+ringIsoZ ::+ ( Integral n+ , PadicNum n+ , KnownRadix p+ , Digit n ~ Mod p+ , Arbitrary n+ , Show n+ )+ => TestName -> n -> TestTree+ringIsoZ s t = testGroup s + [ testProperty "Z <-> Zp" $ homo0 fromInteger toInteger t+ , testProperty "add Z -> Zp" $ homo1 fromInteger (+) (+) t+ , testProperty "add Zp -> Z" $ homo2 fromInteger toInteger (+) (+) t+ , testProperty "mul Z -> Zp" $ homo1 fromInteger (*) (*) t+ , testProperty "mul Zp -> Z" $ homo2 fromInteger toInteger (*) (*) t+ , testProperty "negation Zp" $ invOp fromInteger (+) negate (const True) t+ , testProperty "inversion Zp" $ invOp fromInteger (*) (div 1) isInvertible t+ , ringLaws t+ , testProperty "Division in the ring" $ divMulZ t+ ]++ringIsoZTests = testGroup "Ring isomorphism"+ [ ringIsoZ "Z 2" (0 :: Z 2)+ , ringIsoZ "Z' 2 60" (0 :: Z' 2 60)+ , ringIsoZ "Z 3" (0 :: Z 3)+ , ringIsoZ "Z' 3 60" (0 :: Z' 3 60)+ , ringIsoZ "Z 10" (0 :: Z 10)+ , ringIsoZ "Z' 10 60" (0 :: Z' 10 60)+ , ringIsoZ "Z 65535" (0 :: Z 65535)+ , ringIsoZ "Z' 65535 60" (0 :: Z' 65535 60)+ ]++divMulZ ::+ (Show a, Eq a, Integral a, PadicNum a, KnownRadix p, Digit a ~ Mod p)+ => a -> a -> a -> Property+divMulZ t a b = isInvertible b ==> b * (a `div` b) === a++------------------------------------------------------------+pAdicUnitTests :: TestTree+pAdicUnitTests = testGroup "p-adic units."+ [ testCase "13" $ splitUnit (0 :: Z 2) @?= (0, 64)+ , testCase "14" $ splitUnit (1 :: Z 2) @?= (1, 0)+ , testCase "15" $ splitUnit (100 :: Z 2) @?= (25, 2)+ , testCase "16" $ splitUnit (4 `div` 15 :: Z 2) @?= (1 `div` 15, 2)+ ]++------------------------------------------------------------+testSuite :: TestTree+testSuite = testGroup "Integer"+ [+ showTests+ , digitsTests + , equivTest+ , ringIsoZTests+ , pAdicUnitTests+ ]++main = defaultMain testSuite
+ test/Test/Rational.hs view
@@ -0,0 +1,129 @@+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE UndecidableInstances #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-}++module Test.Rational (testSuite) where+++import Math.NumberTheory.Padic+import Test.Base+import GHC.TypeLits hiding (Mod)+import Test.Tasty+import Test.Tasty.HUnit+import Test.Tasty.QuickCheck+import Test.Tasty.ExpectedFailure+import Test.QuickCheck+import Data.Mod+import Data.Ratio+++------------------------------------------------------------+digitsTestQ :: (Show n, Eq n, PadicNum n) => n -> n -> Property+digitsTestQ t n = valuation n == 0 ==> fromDigits (digits n) === n++digitsTests = testGroup "Conversion to and from digits"+ [ testProperty "Q 2" $ digitsTestQ (0 :: Q 2)+ , testProperty "Q 13" $ digitsTestQ (0 :: Q 13)+ , testProperty "Q 257" $ digitsTestQ (0 :: Q 257)+ , testCase "1" $ firstDigit (1 :: Q 3) @?= 1+ , testCase "-1" $ firstDigit (-1 :: Q 3) @?= 2+ , testCase "2" $ firstDigit (0 :: Q 3) @?= 0+ , testCase "3" $ firstDigit (9 :: Q 3) @?= 0+ , testCase "4" $ firstDigit (9 :: Q 10) @?= 9+ ]+ +------------------------------------------------------------+equivTest :: TestTree+equivTest = testGroup "Equivalence tests"+ [ testCase "5" $ (432100000 :: Q' 10 5) @?= 0+ , testCase "6" $ (0 :: Q' 10 5) @/= 543210000+ , testCase "7" $ (1/7 :: Q' 10 5) @?= 57143+ , testCase "8" $ (1/7 :: Q' 10 5) @?= 657143+ , testCase "9" $ (1/7 :: Q' 10 5) @/= 67143+ ]++------------------------------------------------------------+showTests = testGroup "String representation"+ [ testCase "0" $ show (0 :: Q 3) @?= "0.0"+ , testCase "3" $ show (3 :: Q 3) @?= "10.0"+ , testCase "-3" $ show (-3 :: Q 3) @?= "(2)0.0"+ , testCase "123" $ show (123 :: Q 10) @?= "123.0"+ , testCase "123" $ show (123 :: Q 2) @?= "1111011.0"+ , testCase "-123" $ show (-123 :: Q 10) @?= "(9)877.0"+ , testCase "1/2" $ show (1/2 :: Q 2) @?= "0.1"+ , testCase "-1/2" $ show (-1/2 :: Q 2) @?= "(1).1"+ , testCase "1/15" $ show (1/15 :: Q 3) @?= "(1210).2"+ , testCase "1/700" $ show (1/700 :: Q 10) @?= "(428571).43"+ , testCase "100/7" $ show (100/7 :: Q 10) @?= "(285714)300.0"+ , testCase "0.1" $ show (0.1 :: Q 10) @?= "0.1"+ , testCase "0.01" $ show (0.01 :: Q 10) @?= "0.01"+ , testCase "1/23" $ show (1/23 :: Q 10) @?= "…565217391304347826087.0"+ , testCase "1/23" $ show (1/23 :: Q' 17 5) @?= "… 8 14 13 5 3 . 0"+ , testCase "123456" $ show (123456 :: Q' 257 4) @?= "1 223 96 . 0"+ , testCase "123456" $ show (-123456 :: Q' 257 6) @?= "… 256 256 256 255 33 161 . 0"+ ]++ringIsoQ ::+ ( KnownRadix m+ , Fractional n+ , Real n+ , PadicNum n+ , Digit n ~ Mod m+ , Arbitrary n+ , Show n+ )+ => TestName -> n -> TestTree+ringIsoQ s t = testGroup s + [ testProperty "Q <-> Qp" $ homo0 phi psi t+ , testProperty "add Q -> Qp" $ homo1 phi (+) (+) t+ , testProperty "add Qp -> Q" $ homo2 phi psi (+) (+) t+ , testProperty "mul Q -> Qp" $ homo1 phi (*) (*) t+ , testProperty "mul Qp -> Q" $ homo2 phi psi (*) (*) t+ , testProperty "negation Qp" $ invOp phi (+) negate (const True) t+ , testProperty "inversion Qp" $ invOp phi (*) (1 /) isInvertible t+ , ringLaws t+ ]++phi :: (Fractional n, Real n) => SmallRational -> n+phi (SmallRational r) = fromRational r+psi :: (Fractional n, Real n) => n -> SmallRational+psi = SmallRational . toRational ++ringIsoQTests = testGroup "Ring isomorphism"+ [ ringIsoQ "Q 2" (0 :: Q' 2 68)+ , ringIsoQ "Q 3" (0 :: Q' 3 45)+ , ringIsoQ "Q 5" (0 :: Q' 5 29)+ , ringIsoQ "Q 7" (0 :: Q' 7 26)+ , ringIsoQ "Q 13" (0 :: Q 13)+ , ringIsoQ "Q 257" (0 :: Q 257)+ ]++------------------------------------------------------------+ +pAdicUnitTests :: TestTree+pAdicUnitTests = testGroup "p-adic units."+ [ testCase "8" $ splitUnit (0 :: Q' 2 13) @?= (0, 13)+ , testCase "9" $ splitUnit (1 :: Q 2) @?= (1, 0)+ , testCase "10" $ splitUnit (100 :: Q 2) @?= (25, 2)+ , testCase "11" $ splitUnit (1/96 :: Q 2) @?= (1 `div` 3, -5)+ , testCase "12" $ splitUnit (-1/96 :: Q 2) @?= (-1 `div` 3, -5)+ ]+++------------------------------------------------------------+testSuite :: TestTree+testSuite = testGroup "Rational"+ [+ showTests+ , digitsTests + , equivTest+ , ringIsoQTests+ , pAdicUnitTests+ ]++test = defaultMain testSuite