arithmetic 1.3 → 1.4

raw patch · 5 files changed

+97/−3 lines, 5 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

• arithmetic.cabal +5/−3
• src/Arithmetic/Pell.hs +71/−0
• src/Arithmetic/Quadratic.hs +10/−0
• src/Arithmetic/Utility.hs +3/−0
• src/Test.hs +8/−0

@@ -1,5 +1,5 @@ name: arithmetic-version: 1.3+version: 1.4 category: Number Theory synopsis: Natural number arithmetic license: MIT@@ -12,8 +12,9 @@   This package implements a library of natural number arithmetic functions,   including Montgomery multiplication, the Miller-Rabin primality test,   Lucas sequences, the Williams p+1 factorization method, continued fraction-  representations of natural number square roots, the Jacobi symbol and the-  Tonelli-Shanks algorithm for finding square roots modulo a prime.+  representations of natural number square roots, the Jacobi symbol, the+  Tonelli-Shanks algorithm for finding square roots modulo a prime, and+  the Chakravala method for solving Pell's equation.  Library   build-depends:@@ -32,6 +33,7 @@     Arithmetic.Lucas,     Arithmetic.Modular,     Arithmetic.Montgomery,+    Arithmetic.Pell,     Arithmetic.Polynomial,     Arithmetic.Prime,     Arithmetic.Prime.Factor,

@@ -0,0 +1,71 @@+{- |+module: Arithmetic.Pell+description: Pell's equation (a^2 = n*b^2 + 1)+license: MIT++maintainer: Joe Leslie-Hurd <joe@gilith.com>+stability: provisional+portability: portable+-}+module Arithmetic.Pell+where++import OpenTheory.Primitive.Natural++import Arithmetic.Utility+import qualified Arithmetic.Modular as Modular+import qualified Arithmetic.Quadratic as Quadratic++-------------------------------------------------------------------------------+-- Using the Chakravala method to find the fundamental solution of+-- Pell's equation+--+--   a^2 = n*b^2 + 1+--+-- (where n is not a square).+-------------------------------------------------------------------------------++chakravala :: Natural -> [(Natural,Natural,Natural)]+chakravala n =+    if sqrtN * sqrtN == n then []+    else let a = minM 0 1 in reduce a 1+  where+    reduce a b = (a,b,k) : (if k == 1 && a2 > nb2 then [] else reduce a' b')+      where+        a2 = a * a+        nb2 = n * b * b+        k =  distance a2 nb2+        a' = (a * m + n * b) `div` k+        b' = (a + b * m) `div` k+        j = case Modular.divide k (Modular.negate k a) b of+              Just i -> i+              Nothing -> error "pell: couldn't divide"+        m = minM j k++    sqrtN = Quadratic.rootFloor n++    minM j k =+        if n - m_0 * m_0 <= m_1 * m_1 - n then m_0 else m_1+      where+        m_0 = sqrtN - Modular.subtract k sqrtN j+        m_1 = m_0 + k++-------------------------------------------------------------------------------+-- Finding all integer solutions of Pell's equation+-------------------------------------------------------------------------------++solutions :: Natural -> [(Natural,Natural)]+solutions n =+    (a0,b0) : (if null l then [] else go a1 b1)+  where+    a0 = 1+    b0 = 0+    l = chakravala n+    (a1,b1,_) = last l+    go a b = (a,b) : go (a1 * a + n * b1 * b) (a1 * b + b1 * a)++solution :: Natural -> (Natural,Natural)+solution n =+    case solutions n of+      _ : ab : _ -> ab+      _ -> error "Pell's equation a^2 = n*b^2 + 1 has no nontrivial integer solution when n is square"

@@ -12,6 +12,7 @@  import OpenTheory.Primitive.Natural import qualified Data.List as List+import qualified Data.Maybe as Maybe  import Arithmetic.Utility import qualified Arithmetic.ContinuedFraction as ContinuedFraction@@ -33,6 +34,15 @@     if sqrtn * sqrtn == n then sqrtn else sqrtn + 1   where     sqrtn = rootFloor n++destSquare :: Natural -> Maybe Natural+destSquare n =+    if sqrtn * sqrtn == n then Just sqrtn else Nothing+  where+    sqrtn = rootFloor n++isSquare :: Natural -> Bool+isSquare = Maybe.isJust . destSquare  rootContinuedFraction :: Natural -> ContinuedFraction.ContinuedFraction rootContinuedFraction n =

@@ -14,6 +14,9 @@ import OpenTheory.Natural.Divides import qualified OpenTheory.Natural.Bits as Bits +distance :: Natural -> Natural -> Natural+distance m n = if m <= n then n - m else m - n+ functionPower :: (a -> a) -> Natural -> a -> a functionPower f =     loop

@@ -27,6 +27,7 @@ import qualified Arithmetic.Prime.Factor as Factor import qualified Arithmetic.Modular as Modular import qualified Arithmetic.Montgomery as Montgomery+import qualified Arithmetic.Pell as Pell import qualified Arithmetic.Polynomial as Polynomial import qualified Arithmetic.Quadratic as Quadratic import qualified Arithmetic.Ring as Ring@@ -511,6 +512,12 @@     p = Polynomial.fromCoefficients r (map (Ring.fromNatural r) ps)     q = Polynomial.fromCoefficients r (map (Ring.fromNatural r) (qs ++ [1])) +propPellEquation :: Natural -> Bool+propPellEquation n =+    Quadratic.isSquare n || a * a == n * b * b + 1+  where+    (a,b) = Pell.solution n+ {- np = (0 :: Natural) ps = ([] :: [Natural])@@ -585,4 +592,5 @@        check "Polynomial quotient remainder" propPolynomialQuotientRemainder        check "Polynomial quotient remainder monic"          propPolynomialQuotientRemainderMonic+       check "Pell equation solution" propPellEquation        return ()