diff --git a/continued-fractions.cabal b/continued-fractions.cabal --- a/continued-fractions.cabal +++ b/continued-fractions.cabal @@ -1,34 +1,44 @@ name: continued-fractions -version: 0.9.1.1 +version: 0.10.0.0 stability: provisional -cabal-version: >= 1.6 +cabal-version: >= 1.10 build-type: Simple author: James Cook maintainer: James Cook + Alexandre Rodrigues Baldé license: PublicDomain -homepage: /dev/null +homepage: https://github.com/rockbmb/continued-fractions category: Math, Numerical synopsis: Continued fractions. description: A type and some functions for manipulating and evaluating continued fractions. -tested-with: GHC == 6.8.3, - GHC == 6.10.4, - GHC == 6.12.1, GHC == 6.12.3, - GHC == 7.0.4, - GHC == 7.2.1, - GHC == 7.4.1-rc1 +tested-with: GHC==8.0.2, GHC==8.2.2, GHC==8.4.3, GHC==8.6.1 source-repository head type: git - location: https://github.com/mokus0/continued-fractions.git + location: https://github.com/rockbmb/continued-fractions.git -Library +library hs-source-dirs: src exposed-modules: Math.ContinuedFraction - build-depends: base >= 3 && <5 - if impl(ghc >= 7.2) - ghc-options: -trust base + build-depends: base >= 4.9 && <5 + ghc-options: -trust base -Wall + default-language: Haskell2010 + +test-suite continued-fractions-test + type: exitcode-stdio-1.0 + hs-source-dirs: test + ghc-options: -threaded -Wall + main-is: Tests.hs + build-depends: base >= 4.9 && < 5, + containers, + continued-fractions, + test-framework, + test-framework-quickcheck2, + QuickCheck >= 2.10 + other-modules: Math.CFTests + default-language: Haskell2010 diff --git a/src/Math/ContinuedFraction.hs b/src/Math/ContinuedFraction.hs --- a/src/Math/ContinuedFraction.hs +++ b/src/Math/ContinuedFraction.hs @@ -1,8 +1,10 @@ +{-# LANGUAGE CPP #-} {-# LANGUAGE ParallelListComp #-} -{-# LANGUAGE CPP #-} + #if defined(__GLASGOW_HASKELL__) && __GLASGOW_HASKELL__ >= 702 {-# LANGUAGE Safe #-} #endif + module Math.ContinuedFraction ( CF , cf, gcf @@ -27,7 +29,7 @@ ) where import Control.Arrow ((***)) -import Data.List (tails, mapAccumL) +import Data.List (tails) -- * The 'CF' type and basic operations @@ -44,10 +46,10 @@ -- > data CF a -- > = CFInfinity -- eval CFInfinity = ∞ -- > | CFCont a a (CF a) -- eval (CFCont p q x) = p + q / eval x --- +-- -- |A continued fraction. Constructed by 'cf' or 'gcf'. -data CF a +data CF a = CF a [a] -- ^ Not exported. See 'cf', the public constructor. | GCF a [(a,a)] @@ -85,23 +87,23 @@ ) instance Functor CF where - fmap f (CF b0 cf) = CF (f b0) (map f cf) - fmap f (GCF b0 gcf) = GCF (f b0) (map (f *** f) gcf) + fmap f (CF b0 cf') = CF (f b0) (map f cf') + fmap f (GCF b0 gcf') = GCF (f b0) (map (f *** f) gcf') -- |Extract the partial denominators of a 'CF', normalizing it if necessary so -- that all the partial numerators are 1. asCF :: Fractional a => CF a -> (a, [a]) -asCF (CF b0 cf) = (b0, cf) +asCF (CF b0 cf') = (b0, cf') asCF (GCF b0 []) = (b0, []) -asCF (GCF b0 cf) = (b0, zipWith (*) bs cs) +asCF (GCF b0 cf') = (b0, zipWith (*) bs cs) where - (a:as, bs) = unzip cf - cs = recip a : [recip (a*c) | c <- cs | a <- as] + (a:as, bs) = unzip cf' + cs = recip a : [recip (a' * c) | c <- cs | a' <- as] -- |Extract all the partial numerators and partial denominators of a 'CF'. asGCF :: (Num a, Eq a) => CF a -> (a,[(a,a)]) -asGCF (CF b0 cf) = (b0, [(1, b) | b <- cf]) -asGCF (GCF b0 gcf) = (b0, takeWhile ((/=0).fst) gcf) +asGCF (CF b0 cf') = (b0, [(1, b) | b <- cf']) +asGCF (GCF b0 gcf') = (b0, takeWhile ((/=0).fst) gcf') -- |Truncate a 'CF' to the specified number of partial numerators and denominators. truncateCF :: Int -> CF a -> CF a @@ -238,10 +240,10 @@ steed orig = scanl (+) b0 dxs where - (b0, (a1,b1):gcf) = asGCF orig + (b0, (a1,b1):gcf') = asGCF orig - dxs = a1/b1 : [(b * d - 1) * dx | (a,b) <- gcf | d <- ds | dx <- dxs] - ds = 1/b1 : [recip (b + a * d) | (a,b) <- gcf | d <- ds] + dxs = a1/b1 : [(b * d - 1) * dx | (_,b) <- gcf' | d <- ds | dx <- dxs] + ds = 1/b1 : [recip (b + a * d) | (a,b) <- gcf' | d <- ds] -- |Evaluate the convergents of a continued fraction using Lentz's method. -- Only valid if the denominators in the following recurrence never go to @@ -290,7 +292,7 @@ lentzWith :: (Fractional a, Eq a) => (a -> b) -> (b -> b -> b) -> (b -> b) -> CF a -> [b] lentzWith f op inv (CF 0 ( a :rest)) = map inv (lentzWith f op inv (CF a rest)) lentzWith f op inv (GCF 0 ((a,b):rest)) = map (op (f a) . inv) (lentzWith f op inv (GCF b rest)) -lentzWith f op inv c = scanl opF (f b0) (zipWith (*) cs ds) +lentzWith f op _ c = scanl opF (f b0) (zipWith (*) cs ds) where opF x y = op x (f y) (b0, cs, ds) = lentzRecurrence c @@ -327,7 +329,7 @@ modifiedLentzWith :: (Fractional a, Eq a) => (a -> b) -> (b -> b -> b) -> (b -> b) -> a -> CF a -> [[b]] modifiedLentzWith f op inv z (CF 0 ( a :rest)) = map (map inv ) (modifiedLentzWith f op inv z (CF a rest)) modifiedLentzWith f op inv z (GCF 0 ((a,b):rest)) = map (map (op (f a) . inv)) (modifiedLentzWith f op inv z (GCF b rest)) -modifiedLentzWith f op inv z orig = separate (scanl opF (False, f b0) cds) +modifiedLentzWith f op _ z orig = separate (scanl opF (False, f b0) cds) where (b0, cs, ds) = modifiedLentzRecurrence z orig cds = zipWith mult cs ds @@ -339,7 +341,7 @@ -- a new one every time it encounters (True,_). separate [] = [] separate ((_,x):xs) = case break fst xs of - (xs, ys) -> (x:map snd xs) : separate ys + (xs', ys) -> (x:map snd xs') : separate ys -- precondition: b0 /= 0 modifiedLentzRecurrence :: (Fractional a, Eq a) => a -> CF a -> (a,[(Bool, a)],[(Bool, a)]) @@ -368,4 +370,4 @@ -- in the series. sumPartialProducts :: Num a => [a] -> CF a sumPartialProducts [] = cf 0 [] -sumPartialProducts (x:xs) = gcf 0 ((x,1):[(negate x, 1 + x) | x <- xs]) +sumPartialProducts (x:xs) = gcf 0 ((x, 1):[(negate x', 1 + x') | x' <- xs]) diff --git a/test/Math/CFTests.hs b/test/Math/CFTests.hs new file mode 100644 --- /dev/null +++ b/test/Math/CFTests.hs @@ -0,0 +1,255 @@ +{-# LANGUAGE ExtendedDefaultRules #-} + +{-# OPTIONS_GHC -fno-warn-type-defaults #-} +-- Orphan instance are disabled because of @Arbitrary/CoArbitrary@ instances for @CF@. +{-# OPTIONS_GHC -fno-warn-orphans #-} + +module Math.CFTests + ( cfTests + ) where + +import Control.Applicative +import Data.List +import qualified Data.Set as S +import Math.ContinuedFraction +import Test.QuickCheck +import Test.Framework (Test, testGroup) +import Test.Framework.Providers.QuickCheck2 (testProperty) + +default (Integer, Rational, Double) + +eps :: Double +eps = 2.220446049250313e-16 + +instance Arbitrary a => Arbitrary (CF a) where + arbitrary = do + mode <- arbitrary + if mode + then liftA2 cf arbitrary arbitrary + else liftA2 gcf arbitrary arbitrary + +instance (CoArbitrary a, Num a, Eq a) => CoArbitrary (CF a) where + coarbitrary cF = coarbitrary b0 . coarbitrary terms + where + (b0, terms) = asGCF cF + +cfTests :: [Test] +cfTests = + [ testGroup "asCF" asCF_tests + , testGroup "asGCF" asGCF_tests + , testGroup "truncateCF" truncateCF_tests + , testGroup "partitionCF" partitionCF_tests + , testGroup "evenCF" evenCF_tests + , testGroup "oddCF" oddCF_tests + , testGroup "equiv" equiv_tests + , testGroup "setDenominators" setDenominators_tests + , testGroup "setNumerators" setNumerators_tests + , testGroup "convergents" convergents_tests + , testGroup "steed" steed_tests + , testGroup "lentz" lentz_tests + , testGroup "lentzWith" lentzWith_tests + , testGroup "modifiedLentz" modifiedLentz_tests + , testGroup "modifiedLentzWith" modifiedLentzWith_tests + , testGroup "sumPartialProducts" sumPartialProducts_tests + ] + +asCF_tests :: [Test] +asCF_tests = + [ testProperty "preserves convergents" prop_asCF_preserves_convergents + ] + +prop_asCF_preserves_convergents :: (Eq a, Fractional a) => CF a -> Bool +prop_asCF_preserves_convergents orig = convergents orig == convergents new + where + new = uncurry cf (asCF orig) + +asGCF_tests :: [Test] +asGCF_tests = + [ testProperty "preserves convergents" prop_asCF_preserves_convergents + ] + +truncateCF_tests :: [Test] +truncateCF_tests = + [ testProperty "truncates convergents" prop_truncateCF_truncates_convergents + ] + +prop_truncateCF_truncates_convergents + :: (Eq a, Fractional a) => CF a -> NonNegative Int -> Bool +prop_truncateCF_truncates_convergents orig (NonNegative n) = convergents truncated `isPrefixOf` convergents orig + where + truncated = truncateCF n orig + +partitionCF_tests :: [Test] +partitionCF_tests = + [ testProperty "preserves convergents" prop_partitionCF_preserves_convergents + ] + +prop_partitionCF_preserves_convergents :: (Ord a, Fractional a) => CF a -> Property +prop_partitionCF_preserves_convergents cF = + length (snd (asGCF cF)) > 1 + ==> (origConvergents == S.union evenConvergents oddConvergents) + where + origConvergents = S.fromList $ convergents cF + evenConvergents = S.fromList $ convergents evenCf + oddConvergents = S.fromList $ convergents oddC + + (evenCf, oddC) = partitionCF cF + +evenCF_tests :: [Test] +evenCF_tests = + [ testProperty "preserves last convergent" prop_evenCF_preserves_last_convergent + ] + +prop_evenCF_preserves_last_convergent :: (Fractional a, Eq a) => CF a -> Bool +prop_evenCF_preserves_last_convergent orig = origResult == evenResult + where + origResult = last $ convergents orig + evenResult = last $ convergents (evenCF orig) + +oddCF_tests :: [Test] +oddCF_tests = + [ testProperty "preserves last convergent" prop_oddCF_preserves_last_convergent + ] + +prop_oddCF_preserves_last_convergent :: (Fractional a, Eq a) => CF a -> Bool +prop_oddCF_preserves_last_convergent orig = origResult == oddResult + where + origResult = last $ convergents orig + oddResult = last $ convergents (oddCF orig) + +equiv_tests :: [Test] +equiv_tests = + [ testProperty "preserves convergents" prop_equiv_preserves_convergents + ] + +prop_equiv_preserves_convergents :: (Fractional a, Eq a) => [a] -> CF a -> Bool +prop_equiv_preserves_convergents cs orig = convergents orig == convergents new + where + new = equiv cs orig + +setDenominators_tests :: [Test] +setDenominators_tests = + [ testProperty "preserves convergents" prop_setDenominators_preserves_convergents + ] + +prop_setDenominators_preserves_convergents :: (Fractional a, Eq a) => [a] -> CF a -> Bool +prop_setDenominators_preserves_convergents denoms orig = convergents orig == convergents new + where + new = setDenominators denoms orig + +setNumerators_tests :: [Test] +setNumerators_tests = + [ testProperty "preserves convergents" prop_setNumerators_preserves_convergents + ] + +prop_setNumerators_preserves_convergents :: (Fractional a, Eq a) => [a] -> CF a -> Bool +prop_setNumerators_preserves_convergents nums orig = convergents orig == convergents new + where + new = setNumerators nums orig + +convergents_tests :: [Test] +convergents_tests = + [ testProperty "not null" prop_convergents_not_null + ] + +prop_convergents_not_null :: CF Rational -> Bool +prop_convergents_not_null = not . null . convergents + +steed_tests :: [Test] +steed_tests = + [ testProperty "not null" prop_steed_not_null + ] + +prop_steed_not_null :: CF Rational -> Bool +prop_steed_not_null = not . null . steed + +lentz_tests :: [Test] +lentz_tests = + [ testProperty "not null" prop_lentz_not_null + ] + +prop_lentz_not_null :: CF Rational -> Bool +prop_lentz_not_null = not . null . lentz + +lentzWith_tests :: [Test] +lentzWith_tests = + [ testProperty "not null" prop_lentzWith_not_null + , testProperty "lentzWith log agrees with map log . lentz" prop_lentzWith_log_sane + ] + +prop_lentzWith_not_null :: CF Rational -> Bool +prop_lentzWith_not_null = not . null . lentzWith id (*) (1/) + +prop_lentzWith_log_sane :: CF Double -> Bool +prop_lentzWith_log_sane cF = + let signLog x = (signum x, log (abs x)) + addSignLog (xS,xL) (yS,yL) = (xS*yS, xL+yL) + negateSignLog (s,l) = (negate s, l) + + (sX, x) ~= (sY, y) + = sX == sY + && (absErr <= eps * max 1 n + || relErr <= eps * max 1 n + || any isNaN [absErr, relErr]) + where + absErr = abs (x-y) + relErr = absErr / max (abs x) (abs y) + + n = genericLength a + a = map signLog (lentz cF) + b = lentzWith signLog addSignLog negateSignLog cF + in and (zipWith (~=) a b) + +modifiedLentz_tests :: [Test] +modifiedLentz_tests = + [ testProperty "first convergent" prop_modifiedLentz_first_convergent + , testProperty "not null" prop_modifiedLentz_result_not_null + , testProperty "sublists not null" prop_modifiedLentz_sublists_not_null + ] + +prop_modifiedLentz_first_convergent :: (Fractional a, Eq a) => CF a -> Property +prop_modifiedLentz_first_convergent x = + b0 /= 0 ==> (head . head) (modifiedLentz 1e-30 x) == b0 + where b0 = fst (asGCF x) + +prop_modifiedLentz_sublists_not_null :: CF Rational -> Bool +prop_modifiedLentz_sublists_not_null = not . any null . modifiedLentz 1e-30 + +prop_modifiedLentz_result_not_null :: CF Rational -> Bool +prop_modifiedLentz_result_not_null = not . null . modifiedLentz 1e-30 + +modifiedLentzWith_tests :: [Test] +modifiedLentzWith_tests = + [ testProperty "sanity test with log" prop_modifiedLentzWith_log_sane + ] + +prop_modifiedLentzWith_log_sane :: CF Double -> Bool +prop_modifiedLentzWith_log_sane cF = + let signLog x = (signum x, log (abs x)) + addSignLog (xS,xL) (yS,yL) = (xS*yS, xL+yL) + negateSignLog (s,l) = (negate s, l) + + (sX, x) ~= (sY, y) + = sX == sY + && (absErr <= eps * max 1 n * max 1 m + || relErr <= eps * max 1 n * max 1 m + || any isNaN [absErr, relErr]) + where + absErr = abs (x-y) + relErr = absErr / max (abs x) (abs y) + + tiny = 1e-30 + n = genericLength a; m = genericLength (concat a) + a = map (map signLog) (modifiedLentz tiny cF) + b = modifiedLentzWith signLog addSignLog negateSignLog tiny cF + in and (zipWith (~=) (concat a) (concat b)) + +sumPartialProducts_tests :: [Test] +sumPartialProducts_tests = + [ testProperty "preserves partial sums" + prop_sumPartialProducts_preserves_partial_sums + ] + +prop_sumPartialProducts_preserves_partial_sums :: (Eq a, Fractional a) => [a] -> Bool +prop_sumPartialProducts_preserves_partial_sums xs = + scanl (+) 0 (tail $ scanl (*) 1 xs) == convergents (sumPartialProducts xs) diff --git a/test/Tests.hs b/test/Tests.hs new file mode 100644 --- /dev/null +++ b/test/Tests.hs @@ -0,0 +1,9 @@ +#!/usr/bin/env runhaskell +module Main where + +import Test.Framework (defaultMain) + +import Math.CFTests (cfTests) + +main :: IO () +main = defaultMain cfTests