diff --git a/ats-src/numerics-internal.dats b/ats-src/numerics-internal.dats --- a/ats-src/numerics-internal.dats +++ b/ats-src/numerics-internal.dats @@ -28,7 +28,7 @@ $UN.castvwtp0(z) end -// Fast integer exponentiation. Modified from an example in the manual. +// Fast integer exponentiation. fun exp {n:nat} .. (x : int, n : int(n)) : int = case+ x of | 0 => 0 diff --git a/ats-src/numerics.dats b/ats-src/numerics.dats --- a/ats-src/numerics.dats +++ b/ats-src/numerics.dats @@ -2,10 +2,16 @@ #include "share/atspre_staload.hats" #include "ats-src/numerics-internal.dats" +#include "$PATSHOMELOCS/atscntrb-libgmp/mylibies.hats" staload "contrib/atscntrb-hx-intinf/SATS/intinf_vt.sats" +staload "contrib/atscntrb-libgmp/SATS/gmp.sats" extern +fun mpz_free : (&mpz >> mpz?) -> void = + "mac#" + +extern fun is_prime_ats { n : nat | n > 0 } : int(n) -> bool = "mac#" @@ -25,3 +31,6 @@ implement fib_ats (m) = fib_gmp(m) + +implement mpz_free (x) = + mpz_clear(x) diff --git a/bench/Bench.hs b/bench/Bench.hs --- a/bench/Bench.hs +++ b/bench/Bench.hs @@ -28,7 +28,7 @@ , bench "Ext.smallOmega" $ nf (Ext.smallOmega :: Int -> Int) 91 ] , bgroup "factorial" - [ bench "factorial" $ nfIO (factorial 160) + [ bench "factorial" $ nf factorial 160 , bench "Ext.factorial" $ nf Ext.factorial (160 :: Integer) ] , bgroup "σ" @@ -36,23 +36,19 @@ , bench "Ext.sigma" $ nf (Ext.sigma 1) (115 :: Int) ] , bgroup "doubleFactorial" - [ bench "doubleFactorial" $ nfIO (doubleFactorial 79) + [ bench "doubleFactorial" $ nf doubleFactorial 79 , bench "Ext.doubleFactorial" $ nf Ext.doubleFactorial (79 :: Integer) ] , bgroup "choose" - [ bench "choose" $ nfIO (choose 322 16) + [ bench "choose" $ nf (choose 322) 16 , bench "Ext.binomial" $ nf (Ext.binomial 322) (16 :: Int) ] , bgroup "catalan" - [ bench "catalan" $ nfIO (catalan 300) + [ bench "catalan" $ nf catalan 300 , bench "Ext.catalan" $ nf Ext.catalan (300 :: Int) ] - {- , bgroup "jacobi" -} - {- [ bench "jacobi" $ whnf (jacobi 80) 111 -} - {- , bench "Ext.jacobi" $ whnf ((Ext.jacobi :: Int -> Int -> Ext.JacobiSymbol) 80) 111 -} - {- ] -} , bgroup "fibonacci" - [ bench "fibonacci" $ nfIO (fibonacci 200) + [ bench "fibonacci" $ nf fibonacci 200 , bench "hsFibonacci" $ nf (hsFibonacci :: Int -> Integer) 200 ] ] diff --git a/fast-arithmetic.cabal b/fast-arithmetic.cabal --- a/fast-arithmetic.cabal +++ b/fast-arithmetic.cabal @@ -1,5 +1,5 @@ name: fast-arithmetic -version: 0.3.3.0 +version: 0.3.3.1 synopsis: Fast functions on integers. description: Fast functions for number theory and combinatorics with a high level of safety guaranteed by [ATS](http://www.ats-lang.org/). homepage: https://github.com/vmchale/fast-arithmetic#readme diff --git a/src/Numeric/Combinatorics.hs b/src/Numeric/Combinatorics.hs --- a/src/Numeric/Combinatorics.hs +++ b/src/Numeric/Combinatorics.hs @@ -18,6 +18,7 @@ import Foreign.C import Foreign.Ptr import Foreign.Storable +import System.IO.Unsafe (unsafePerformIO) foreign import ccall unsafe double_factorial_ats :: CInt -> Ptr GMPInt foreign import ccall unsafe factorial_ats :: CInt -> Ptr GMPInt @@ -29,16 +30,16 @@ -- -- > λ:> mapM catalan [0..9] -- > [1,1,2,5,14,42,132,429,1430,4862] -catalan :: Int -> IO Integer -catalan = conjugateGMP catalan_ats +catalan :: Int -> Integer +catalan = unsafePerformIO . conjugateGMP catalan_ats -- | See [here](http://mathworld.wolfram.com/BinomialCoefficient.html). -choose :: Int -> Int -> IO Integer -choose = (gmpToInteger <=<) . (peek .* on choose_ats fromIntegral) +choose :: Int -> Int -> Integer +choose = unsafePerformIO .* (gmpToInteger <=<) . (peek .* on choose_ats fromIntegral) -factorial :: Int -> IO Integer -factorial = conjugateGMP factorial_ats +factorial :: Int -> Integer +factorial = unsafePerformIO . conjugateGMP factorial_ats -- | See [here](http://mathworld.wolfram.com/DoubleFactorial.html). -doubleFactorial :: Int -> IO Integer -doubleFactorial = conjugateGMP double_factorial_ats +doubleFactorial :: Int -> Integer +doubleFactorial = unsafePerformIO . conjugateGMP double_factorial_ats diff --git a/src/Numeric/Integer.hs b/src/Numeric/Integer.hs --- a/src/Numeric/Integer.hs +++ b/src/Numeric/Integer.hs @@ -13,13 +13,14 @@ import Foreign.C import Foreign.Ptr import Numeric.Common +import System.IO.Unsafe (unsafePerformIO) foreign import ccall is_prime_ats :: CInt -> CBool foreign import ccall unsafe fib_ats :: CInt -> Ptr GMPInt -- | Indexed starting at @0@. -fibonacci :: Int -> IO Integer -fibonacci = conjugateGMP fib_ats +fibonacci :: Int -> Integer +fibonacci = unsafePerformIO . conjugateGMP fib_ats -- | O(√n) isPrime :: Int -> Bool diff --git a/test/Spec.hs b/test/Spec.hs --- a/test/Spec.hs +++ b/test/Spec.hs @@ -27,10 +27,10 @@ prop "should agree with the pure Haskell function" $ \n -> n < 1 || f n == g n -check :: (Eq a, Show a) => String -> (Int -> IO a) -> (Integer -> a) -> Int -> SpecWith () +check :: (Eq a, Show a) => String -> (Int -> a) -> (Integer -> a) -> Int -> SpecWith () check s f g n = describe s $ it ("should work for n=" ++ show n) $ - f n >>= (`shouldBe` g (fromIntegral n)) + f n `shouldBe` g (fromIntegral n) main :: IO () main = hspec $ parallel $ do @@ -49,8 +49,8 @@ it "should match the arithmoi function" $ toInt (Ext.jacobi (15 :: Int) 19) `shouldBe` toInt (Ext.jacobi (15 :: Int) 19) describe "totient" $ - prop "should be equal to m-1 for m prime" $ - \m -> m < 1 || not (isPrime m) || totient m == m - 1 + prop "should be equal to p-1 for p prime" $ + \p -> p < 1 || not (isPrime p) || totient p == p - 1 describe "derangement" $ prop "should be equal to [n!/e]" $ \n -> n < 1 || n > 18 || (derangement n :: Integer) == floor ((fromIntegral (Ext.factorial (fromIntegral n :: Int) :: Integer) :: Double) / exp 1 + 0.5)